$h = u+iv$ is a complex-valued harmonic function if
both real and imaginary part are harmonic:
$$
u_{xx} + u_{yy} = 0 \\
v_{xx} + v_{yy} = 0 \\
$$
In the same manner, $$
f(z) = zh(z) = (x + iy)(u(z) + iv(z)) = (xu(z)-yv(z)) + i (xv(z)+yu(z))
$$
is a complex-valued harmonic function if both real and imaginary
part are harmonic functions:
$$
(xu(z)-yv(z))_{xx} + (xu(z)-yv(z))_{yy} = 0 \\
(xv(z)+yu(z))_{xx} + (xv(z)+yu(z))_{yy} = 0
$$
If you expand the last two equations and simplify the result by
using the fact that $h$ is harmonic (the first two equations),
then you can conclude that $h = u+iv$ satisfies the
Cauchy-Riemann differential equations.
Alternatively (as suggested by Daniel Fischer), using the
Laplace operator $\Delta = \dfrac{\partial^2}{\partial x^2}
+ \dfrac{\partial^2}{\partial y^2}$ and the Wirtinger derivatives:
$$
\Delta (zh(z))= 4\dfrac{\partial^2}{\partial \overline{z}\partial z}(zh(z)) =
4\dfrac{\partial}{\partial \overline{z}} \bigl( h(z) + z \dfrac{\partial}{\partial z} h(z) \bigr) \\
= 4\dfrac{\partial}{\partial \overline{z}} h(z) + 4\dfrac{\partial^2}{\partial \overline{z}\partial z} h(z) =
4\dfrac{\partial}{\partial \overline{z}} h(z) + \Delta h(z)
$$
because $\dfrac{\partial z}{\partial z} = 1$ and $\dfrac{\partial z}{\partial \overline{z}} = 0$.
If $h(z)$ and $zh(z)$ are harmonic then
$\Delta h(z) = \Delta (zh(z)) = 0$ and therefore
$\dfrac{\partial}{\partial \overline{z}} h(z) = 0$,
which is equivalent to $h$ satisfying the Cauchy-Riemann
differential equations.