Partial answer
Both the spatial and temporal domains are semi-infinite. You can leverage this and use, for instance, Laplace transform in $t$ to find
$$s^2 W(x;s) - s \phi(x)-\psi(x)-c^2 W_{xx}(x;s) = F(x;s) $$
where $W(x;s)$ and $F(x;s)$ are the Laplace transforms of $w(x,t)$ and $f(x,t)$, respectively, defined as
$$W(x;s) = \int^\infty_0 e^{-st}w(x,t)\,\mathrm{d}t $$
Now you can see this equation as a constant-coefficient inhomogenous second order ODE for $W(x;s)$ as a function of $x$
$$ W_{xx} - (s/c)^2 W = - (F+s\phi + \psi)/c^2 \equiv Q(x;s), $$
which is easy to work with, given the initial conditions in $x$. Now we have to work out the solution for $W$. Let $W_1 = e^{sx/c}$ and $W_2 = e^{-sx/c}$ be the solution of the homogenous part of the equation. Then, according to Bender & Orszag (p15), variation of parameters gives the complete solution as
$$ W(x;s) = A_1 W_1 + A_2 W_2 + W_p $$
where the particular solution $W_p$ turns out to be
$$ W_p(x;s) = -W_1 \int^{x}_0 \frac{Q(\xi;s) W_2(\xi;s)}{R(\xi)}\,\mathrm{d}\xi + W_2 \int^{x}_0 \frac{Q(\xi;s) W_1(\xi;s)}{R(\xi)}\,\mathrm{d}\xi $$
where $R = W_1 \, \partial_x W_2 - W_2 \, \partial_x W_1$ is the Wronksian of $(W_1,W_2)$. You now may impose the initial conditions on $x$ and inverse Laplace transform your solution, which, of course depends on the shape of $f$ through the definition of $Q$.
An alternative
View the PDE as
$$ (\partial_t + c \partial_x)(\partial_t - c \partial_x) w = f(x,t) $$
and define $u = (\partial_t - c \partial_x) w$ such that $(\partial_t + c \partial_x) u = f$, the latter being amenable for the use of the method of characteristics, which reads
$$ \mathrm{d}x/c = \mathrm{d}t = \mathrm{d}u/f, $$
which tells us that $x - ct = c_1$ is a characteristic line. On the other hand, one must integrate the relation $\mathrm{d}t = \mathrm{d}u/f(ct + c_1, t)$, or equivalently:
$$ u = c_2 + \int f(c t + c_1, t) \, \mathrm{d}t $$
put $c_2$ as a function of $c_1$ to have $ u = H(x-ct) + \int f(x, t) \, \mathrm{d}t $ where $H$ is an arbitrary function of its argument. Similarly, according to the definition of $u$, one gets the 1st-order PDE for w:
$$ u = H(x-ct) + \int f(x, t) \, \mathrm{d}t = (\partial_t - c \partial_x) w $$ from which the method of characteristics provides
$$ -\mathrm{d}x/c = \mathrm{d}t = \frac{\mathrm{d}w}{H(x-ct) + \int f(x, t) \, \mathrm{d}t} $$ which provides $w = M(x-ct) + N(x+ct) + P$ is as solution, where $M$ and $N$ are arbitrary functions and $P$ is a particular solution depending on $f(x,t)$. Initial conditions should be imposed and convenient limits should be considered for the integrals (but I'm not sure about satisfying the conditions at $ x = 0$).
I hope you find this helpful.