Let $u := x + \frac{1}{y}$, $v := y + \frac{1}{x}$ for $x > 0, y > 0$. Then we have
$$f(u+v) = f(u) + f(v) \quad (1)$$
Here $uv = xy + 2 + \frac{1}{xy} \ge 4$. First we show that for any $u, v > 0$ with $uv \ge 4$ can be represented as
$$u = x + \frac{1}{y}, v = y + \frac{1}{x} \quad (2)$$
for some $x > 0, y > 0$. From (2), we deduce
$$u = x + \frac{1}{y} \implies x = u - \frac{1}{y} = \frac{uy - 1}{y} \implies v = y + \frac{1}{x} = \frac{uy^2}{uy - 1} \implies uy^2 - uvy + v = 0.$$

The last quadratic equation has a real solution $y > 0$ because $D = (uv)^2 - 4uv \ge 0$ and $u, v > 0$. (Observe that there's no negative real solution for the equation.) Similarly, there exists real $x > 0$, which satisfies (2) together with the above $y$. Therefore, (1) holds for all $u > 0, v > 0$ $uv \ge 4$.
We verify that (1) is valid for all $u, v > 0$. Indeed, for given $u, v > 0$, we may choose a real number $w > 0$ such that
$$(u + v)w \ge 4, \quad u(v + w) \ge 4, \quad vw \ge 4.$$
Then from
$$\begin{aligned} f(u + v + w) &= f(u + v) + f(w) \\ &= f(u) + f(v + w) \\ &= f(u) + f(v) + f(w). \end{aligned}$$
we may conclude that $f(u + v) = f(u) + f(v)$ for all $u, v \ge 0$.
Now we put $h(x) = f(x) - f(1)x = f(x) - 2008x$ and show that $h(x)$ is a zero function. It is clear that $h$ satisfies (1) and
$$h(x + 1) = h(x) + h(1) = h(x) \quad (3)$$
for all $x \in \mathbb{R}^+$. It follows from the condition $|f(x)| \le x^2 + 1004^2$ that
$$-(x + 1004)^2 \le h(x) \le (x - 1004)^2,$$
which implies that $h(x)$ is bounded on $(0,1]$ and hence $h(x)$ is bounded on $\mathbb{R}^+$ according to (3). Since $h(nx) = nh(x)$ for all positive integer $n$, we must have $h \equiv 0$. Therefore we have $f(x) = 2008x$. $\square$