Let $L=\Delta + c$ in 3 dimensions, where $c$ is a positive constant.

I met this modified mean value property of a solution $u$ of $Lu=0$ as $$u(\xi)=\frac{\sqrt{c}\rho}{sin(\sqrt{c}\rho)}\frac{1}{4\pi \rho^2}\int_{\partial B(\xi,\rho)} u(x)d\sigma(x)$$ where $sin(\sqrt{c}\rho) \ne 0$, and $\sigma$ denotes the surface measure.

I'm trying to prove it. The following is what I did.

First I figured out that $$K(x,\xi)=-\frac{cos(\sqrt{c}r)}{4\pi r}, r=|x-\xi|$$ is a fundamental solution of $L$ with pole $\xi$. Then I tried to apply the second Green identity (which is also true for $L$ for sure) $$\int_{\Omega}(LG) u-\int_{\Omega}G (Lu)=\int_{\partial \Omega} \frac{\partial G}{\partial \nu} u- G \frac{\partial u}{\partial \nu}$$ For fixed $\rho$, let $G=K+\frac{cos(\sqrt{c}\rho)}{4\pi \rho}$ and $\Omega=B_{\rho}(\xi)$, so that $G$ vanishes on the boundary of $\Omega$, and $LG=\delta_{\xi}+c\frac{cos(\sqrt{c}\rho)}{4\pi \rho}$, where $\delta_{\xi}$ denotes the Dirac mass at $\xi$. Then the identity above becomes $$u(\xi)+c\frac{cos(\sqrt{c}\rho)}{4\pi \rho}\int_{B_{\rho}(\xi)}u(x)dx=\frac{\sqrt{c}\rho sin(\sqrt{c}\rho)+cos(\sqrt{c}\rho)}{4 \pi \rho^2} \int_{\partial B(\xi,\rho)} u(x)d\sigma(x)$$

Then no matter how hard I tried, I just could not get the desired form.

I also tried to take the derivative with respect to $\rho$ of the desired form, and I got a complicated form which doesn't seem to be zero. Then I took the second derivative, and things got worse...

Can anyone either give the right reference about the solution of this operator $L$, or help me figure out how to obtain the desired form? Thanks!