Mark Srednicki's QFT book presents a regularization of the $\delta$ function in calculating the chiral anomaly (see section 77 of the book). This regularization reads \begin{equation} \delta (x-y)=\lim_{M \rightarrow \infty}\int \frac{d^4k}{(2\pi)^4} e^{(i\gamma ^{\mu}D_{\mu})^2/M^2}\circ e^{-ik(x-y)}, \end{equation} where $D_{\mu}=\partial_{\mu}-igA_{\mu}$.
Now I am trying to applying this method to calculate the chiral anomaly of a massless fermion in a gravity field but without gauge field. The action in the gravitational field is \begin{equation} S=\int d^4x \sqrt{g}\bar{\Psi}i\gamma^{\mu}D_{\mu}\Psi, \end{equation} where $D_{\mu}$ now is $\partial_{\mu}+\frac{1}{2}\omega_{\mu}^{ab}\sigma_{ab}$. Under a chiral transformation \begin{equation} \Psi'(x)=e^{-i\alpha(x)\gamma^5}\Psi(x)=\int d^4y ~\delta(x-y)e^{-i\alpha(y)\gamma^5}\Psi(y), \end{equation} one gets the shift of the path integral measure: \begin{equation} \mathcal{D \Psi'}\mathcal{D}{\bar{\Psi}'}=\mathcal{D \Psi}\mathcal{D}{\bar{\Psi}}\exp\bigg\{2i\int d^4x \alpha(x) Tr[\delta (x-x) \gamma^5]\bigg\}. \end{equation} At this step, I still regularize $\delta$ function as \begin{equation} \begin{aligned} \delta (x-y)&=\lim_{M \rightarrow \infty}\int \frac{d^4k}{(2\pi)^4} e^{(i\gamma ^{\mu}D_{\mu})^2/M^2}\circ e^{-ik(x-y)}\\ &=\lim_{M \rightarrow \infty}\int \frac{d^4k}{(2\pi)^4}e^{-ik(x-y)}\circ e^{-(\gamma ^{\mu}D_{\mu}-ik_{\mu}\gamma^{\mu})^2/M^2}. \end{aligned} \end{equation} We can expand the square as \begin{equation} (\gamma ^{\mu}D_{\mu}-ik_{\mu}\gamma^{\mu})^2=\frac{1}{\sqrt{g}}D_{\mu}\sqrt{g}g^{\mu \nu}D_{\nu}-\frac{R}{4}+\{\gamma ^{\mu}D_{\mu},-ik_{\mu}\gamma^{\mu}\}-k^2. \end{equation} Here, I cannot continue. Could someone please correct me? I have already known the result is a quadratic expression in terms of Riemann tensor and its dual. So I guess the square should produce something like $[D_{\mu},D_{\nu}]=R^{ab}_{\mu \nu}\sigma_{ab}/2$, together with the $\gamma^5$ put in the trace, the result can be recovered.