Hello I have a quick question regarding the Dirac gamma matrices or Dirac equation and the Klein-Gordon equation. Recall that the Klein-Gordon is given by the following \begin{equation} \left(\Box + m^2\right)\Psi=0 \end{equation} where $\Box$ is the D'Alembert operator. Then, we have the Dirac equation written in terms of the gamma matrices (or covariant form) \begin{equation} \left( i\gamma^{\mu}\partial_{\mu}-m\right)\Psi=0 \end{equation} I multiplied both sides of the Dirac equation by the $\beta$ matrix, such that I may be able to reduce the gamma matrices, $\gamma^{\mu}$, to the Dirac matrices, $\alpha^{i}$, after which I $\textbf{square}$ both sides of the equation, from which we have the Klein-Gordon equation after imposing certain restrictions on the alpha and beta matrices, i.e. \begin{equation}\left(-i \alpha_{x} \frac{\partial}{\partial x}-i \alpha_{y} \frac{\partial}{\partial y}-i \alpha_{z} \frac{\partial}{\partial z}+\beta m\right)\left(-i \alpha_{x} \frac{\partial}{\partial x}-i \alpha_{y} \frac{\partial}{\partial y}-i \alpha_{z} \frac{\partial}{\partial z}+\beta m\right) \psi=-\frac{\partial^{2} \psi}{\partial t^{2}}\end{equation} \begin{equation}\begin{aligned} \implies -\frac{\partial^{2} \psi}{\partial t^{2}}=&-\alpha_{x}^{2} \frac{\partial^{2} \psi}{\partial x^{2}}-\alpha_{y}^{2} \frac{\partial^{2} \psi}{\partial y^{2}}-\alpha_{z}^{2} \frac{\partial^{2} \psi}{\partial z^{2}}+\beta^{2} m^{2} \psi \\ &-\left(\alpha_{x} \alpha_{y}+\alpha_{y} \alpha_{x}\right) \frac{\partial^{2} \psi}{\partial x \partial y}-\left(\alpha_{y} \alpha_{z}+\alpha_{z} \alpha_{y}\right) \frac{\partial^{2} \psi}{\partial y \partial z}-\left(\alpha_{z} \alpha_{x}+\alpha_{x} \alpha_{z}\right) \frac{\partial^{2} \psi}{\partial z \partial x} \\ &-\left(\alpha_{x} \beta+\beta \alpha_{x}\right) m \frac{\partial \psi}{\partial x}-\left(\alpha_{y} \beta+\beta \alpha_{y}\right) m \frac{\partial \psi}{\partial y}-\left(\alpha_{z} \beta+\beta \alpha_{z}\right) m \frac{\partial \psi}{\partial z} \end{aligned}\end{equation} such that, \begin{equation}\begin{aligned} \alpha_{x}^{2}=\alpha_{y}^{2}=\alpha_{z}^{2}=\beta^{2} &=1 \\ \alpha_{j} \beta+\beta \alpha_{j} &=0 \\ \alpha_{j} \alpha_{k}+\alpha_{k} \alpha_{j} &=0 \quad(j \neq k) \end{aligned}\end{equation}
Now, this works swell for the normal Dirac equation, but what about the Dirac equation which may have additional terms, i.e. the Dirac equation in curved spacetime, in which we have the spin connection, as well as other terms, i.e. the Levi-Civita term. This proves extremely inefficient, as multiplying both sides by the $\beta$ matrix, then squaring both sides gives an enormous amount of terms. Is there any other way in which perhaps we can derive the Klein-Gordon using gamma matrices, or perhaps derive the Klein-Gordon from the Dirac? Thank you for your help