Solution of Klein-Gordon equation Rindler spacetime The Klein Gordon equation in rindler coordinates is $$e^{-2a\epsilon}(-\partial_{n}^2 + \partial_{\epsilon}^2)\phi + m^2\phi = 0.$$
if $m=0$ this would reduce to the KG equation for inertial coordinates, but for the case where $m$ is nonzero, the differential equation is no longer homogenous. How can we solve it?
Rindler coordinates $(n, \epsilon)$ are $x = \frac{1}{a}e^{a\epsilon}\cosh(an)$ and $t = \frac{1}{a}e^{a\epsilon}\sinh(an)$
 A: Separation of variables: $\phi(n,\epsilon) = N(n)E(\epsilon)$.
This leads the equation to,
$\dfrac{1}{E}\dfrac{d^{2}E}{d\epsilon^2} + m^2e^{2a\epsilon} = \dfrac{1}{N}\dfrac{d^{2}N}{dn^2} = \omega^2$ (let)
$\therefore N \propto e^{-\omega n}$ to avoid singularity at $n \rightarrow \infty$.
The remaining part gives the so-called master equation (not to be confused with the master equation of Physical Kinetics!),
$\dfrac{d^{2}E}{d\epsilon^2} + (m^2e^{2a\epsilon}-\omega^2)E = 0$
Now this equation can be solved in numerous way, popular and often most useful two among them are as follows:

*

*Leaver's solution: This gives the exact solution of this master equation, but the final expression turns out to be a horrible one! It gives solution in the form of a continued fraction, which is difficult to handle (both analytically and numerically) in most of the times. But the thorough process is quite difficult to describe here. you may go through Leaver's original paper ($1985$) "An analytic representation for the quasi-normal modes of Kerr black holes, Proc. R. Soc. Lond., $A402285–298$", where the master equation is being solved in Schwarzchild background. you may follow this to do your problem yourself.


*WKB approximation technique: More physically intuitive, but well approximated result is given by WKB method. (Yes, the quantum mechanical one! Fortunately it also works here and is well fitted for such kind of problems!)
Here you should assume $E(\epsilon) \sim exp[\dfrac{1}{\delta}\sum_{k=0}^ {\infty}\delta^kS_k(\epsilon)]$.
If you take first two terms of this expansion, your solution will come in terms of parabolic cylindrical functions and the frequency $\omega$ will then be quantized (Bohr-Sommerfeld quantization rule). Also the entire of this method is difficult to discuss here entirely.
For further reference and a thorough discussion of the solution of this master equation, see the 3rd chapter of Black Hole Perturbation Theory notes by Prof. Emanuele Berti and the references and citations therein. Here the master equation has been discussed for the Schwarzchild case, but you can easily generalize it to your Rindler case.
A: The real eigenfunctions of the equation in @SCh's  answer
$$
\left(-\frac{d^2}{d\xi^2}+m^2 e^{2\xi}\right) \psi_\nu(\xi)=\nu^2 \psi_\nu(\xi)
$$
are
$$
\psi_\nu(\xi)=\left( \frac{2\nu \sinh  \nu \pi}{\pi^2}\right)^{1/2} K_{i\nu}(m e^{\xi}).
$$
where $K_{i\nu}$ is the Macdonald (Bessel K) function with an imaginary index.
$$
{\rm K}_{i\nu}(x)= \int_0^\infty e^{-x\cosh u}\cos \nu u\, du\nonumber\\
=\frac 12 \left(\frac{x}{2}\right)^{i\nu} \int_0^\infty \exp\left(-t-\frac{x^2}{4t}\right) t^{-i\nu -1}\,dt.\nonumber
$$
The $\psi_\nu$  have been normalized to  obey
$$
\int_{-\infty}^{\infty} \psi_\mu(\xi)\psi_\nu(\xi)\, d\xi = \delta(\nu-\mu).
$$
The completeness and orthogonality of eigenfunctions gives us the  gives us the  the Kontorovich-Lebedev  transform pair
$$
\tilde f(\nu)\equiv K[f](\nu)= \int_0^\infty {\rm K}_{i\nu}(x) f(x)\,dx,\nonumber\\
f(x)= \frac{1}{\pi^2 x} \int_0^\infty 2\nu \sinh \nu\pi\, {\rm K}_{i\nu}(x) \tilde f(\nu)\,d\nu.\nonumber
$$
There are more formulae in my online notes
