With the added condition $(1)$, one does have $(*)$.

Assume first that $R=1$; the general case will follow by simple rescaling.

By $(1)$ (with $R=1$), without loss of generality
$$\text{$|X|\le 1$ and $|Y|\le 1$}\tag{!}
$$
(see details on this below, at the end of this answer),
and then $(2)$ holds for all real $K$.

Let $F$ and $G$ be the distribution functions (d.f.'s) of $X$ and $Y$, respectively. If the Lévy distance
$$L(F,G):=\inf\{h>0:F(x-h)-h\le G(x)\le F(x+h)+h \text{ for all real }x\}
$$
between $F$ and $G$ is greater than $\delta>0$, then for some real $t$ without loss of generality $F(t-\delta)-G(t)>\delta$, whence
$$\int_{t-\delta}^t[F(x)-G(x)]\,dx\ge\int_{t-\delta}^t[F(t-\delta)-G(t)]\,dx>\delta^2.
$$
But $\int_{t-\delta}^t[F(x)-G(x)]\,dx=\int_{-\infty}^t[F(x)-G(x)]\,dx-\int_{-\infty}^{t-\delta}[F(x)-G(x)]\,dx.
$ So,
$$\Big|\int_{t-\delta}^t[F(x)-G(x)]\,dx\Big|
=\Big|\int_{-\infty}^t[F(x)-G(x)]\,dx-\int_{-\infty}^{t-\delta}[F(x)-G(x)]\,dx\Big|
$$
$$
=\big|\big(E(t-X)_+-E(t-Y)_+\big)-\big(E(t-\delta-X)_+-E(t-\delta-Y)_+\big)\big|
\le2\epsilon
$$
by $(2)$.
Hence, $\delta^2\le2\epsilon$. That is,
$$L(F,G)\le\sqrt{2\epsilon}.
$$
That is, the supremum of the distance in the direction of the vector $(-1,1)$ between the graphs of $F$ and $G$ is $\le2\sqrt\epsilon$. (By the graph of a nondecreasing function $f$ from an open interval $I$ to $\mathbb R$ here I understand the set $\{(x,y)\colon x\in I,f(x-)\le y\le f(x+)\}$, with the corresponding vertical segments added at all points of discontinuity.)

Hence, the area between the graphs of $F$ and $G$ is $\le2\sqrt\epsilon\,\frac3{\sqrt2}=3\sqrt{2\epsilon}$, because the distance between the straight lines through points $(-1,0)$ and $(1,1)$ parallel to the vector $(-1,1)$ is $\frac3{\sqrt2}$.

Without loss of generality, re-define now $X$ and $Y$ by the formulas $X:=F^{-1}(U)$ and $Y:=G^{-1}(U)$, where $U$ is a r.v. uniformly distributed in the interval $(0,1)$ and
$$F^{-1}(u):=\inf\{x\in\mathbb R\colon F(x)\ge u\}
$$
for $u\in(0,1)$.
Then $F$ and $G$ will indeed be the d.f.'s of $X$ and $Y$, respectively. Also, then we will have
$E|X-Y|=\int_0^1|F^{-1}(u)-G^{-1}(u)|\,du$, which is the area between the graphs of $F^{-1}$ and $G^{-1}$, which latter is the same as the area between the graphs of $F$ and $G$, which we showed to be $\le3\sqrt{2\epsilon}$.

So, if $(!)$ holds, then $(*)$ holds with $f(\epsilon):=3\sqrt{2\epsilon}$.

Without assuming $(!)$, but still assuming $(1)$ and $(2)$,
one can replace $X$ and $Y$ by their truncated versions $X_1:=1\wedge((-1)\vee X)$ and $Y_1:=1\wedge((-1)\vee Y)$, so that $|X_R|\le 1$ and $|Y_R|\le 1$. Then $(2)$ will hold for all real $K$ if $\epsilon$, $\mu$, and $\nu$ are replaced there, respectively, by $2\epsilon$, the distribution of $X_1$, and the distribution of $Y_1$.
Also, $E|X-Y|\le E|X_1-Y_1|+E|X-X_1|+E|Y-Y_1|$, $E|X-X_1|=E(-1-X)_++E(X-1)_+=E(-1-X)I\{X<-1\}+E(X-1)I\{X>1\}\le
E(-X)I\{X<-1\}+EXI\{X>1\}=E|X|I\{|X|>1\}\le\epsilon$ (by $(1)$, with $R=1$), and similarly $E|Y-Y_1|\le\epsilon$.
So,
$(*)$ will hold in general for $R=1$ with
$f(\epsilon):=3\sqrt{2\times2\epsilon}+2\epsilon=6\sqrt{\epsilon}+2\epsilon.$

To rescale from $R=1$ to general $R>0$, replace in the above reasoning $X$, $Y$, $\epsilon$ by $\tilde X:=X/R$, $\tilde Y:=Y/R$, $\tilde\epsilon:=\epsilon/R$, respectively, so that the conditions $(1)$ and $(2)$ hold with $1$ in place of $R$, $\tilde\epsilon$ in place of $\epsilon$, and with the distributions of $\tilde X$ and $\tilde Y$ in place of $\mu$ and $\nu$.
Then, by the above, $E|\tilde X-\tilde Y|\le6\sqrt{\tilde\epsilon}+2\tilde\epsilon$, which can be rewritten as
$$E|X-Y|\le6\sqrt{R\epsilon}+2\epsilon.\tag{!!}$$

It is also seen that, instead of $(1)$, the following weaker condition will suffice:
$$E(-R-X)_++E(X-R)_+\le\epsilon,\quad E(-R-Y)_++E(Y-R)_+\le\epsilon.
$$

Condition $(1.5)$ is not needed.

Let us now show that the upper bound in $(!!)$ is best possible in terms of $R$ and $\epsilon$, up to a universal constant factor.
Let $F$ be the d.f.\ of the uniform distribution on $(-R,R)$. For any natural $n$ and all $k=1,\dots,n$, let $G(x):=F(-R+R\,\frac{2k-1}{n})$ for $x\in[-R+R\,\frac{2k-2}{n},-R+R\,\frac{2k}{n})$, with $G=0$ on $(-\infty,-R)$ and $G=1$ on $[R,\infty)$. Then $G$ is the d.f. of a r.v. $Y$, $|E(X-t)_+-E(Y-t)_+|=|\int_t^\infty[F(x)-G(x)]\,dx|\le\frac12\,\frac{2R}{2n}\frac1n=\frac{R}{2n^2}=:\epsilon$ for all real $t$,
$|E(t-X)_+-E(t-Y)_+|=|\int_{-\infty}^t[F(x)-G(x)]\,dx|\le\frac12\,\frac{2R}{2n}\frac1n=\frac{R}{2n^2}=\epsilon$ for all real $t$, $E|X|I\{|X|>R\}=E|Y|I\{|Y|>R\}=0$, whereas
the least possible value here of $E|X-Y|$ is the Wasserstein distance between $F$ and $G$, which equals
$$d(F,G):=\int_0^1|F^{-1}(u)-G^{-1}(u)|\,du
$$
-- see e.g. (2) in **[1]**; hence, this distance also equals $\int_{\mathbb R}|F(x)-G(x)|\,dx=2n\epsilon=2\sqrt{\frac{R}{2\epsilon}}\,\epsilon=\sqrt{2R\epsilon}$, which shows that the upper bound in $(!!)$ is indeed best possible, up to a universal constant factor, for (say) $\epsilon\le R$.