Q:

Weitzenböck identity, and how to define the curved-space Laplacian of a spinor field?

I am following the notes by Freed about the Dirac operator in section 2. I am trying to understand the Weitzenbock equality that for a curved-space spinor's Dirac operator $D\!\!\!/\ $ and the associated 'Laplacian' $\Delta = \nabla^*\nabla $ satisfy
$$D\!\!\!/\ D\!\!\!/\  - \Delta = \frac{R}{4}$$
I'm trying to reproduce this formula using different notation, but I'm having difficulty defining the Laplacian $\Delta$. I know how to define the curved space Dirac operator in terms of an orthonormal frame $e^\mu_a$ and a basis of gamma matrices $\gamma^a$ as:
$$D\!\!\!/\ = \gamma^a e^\mu_a(\partial_\mu + \omega_\mu^I T^I) =: \gamma^a e^\mu_a \nabla_\mu$$
where I'm writing the $SO(n)$ spin connection as $\omega_\mu^I T^I$ for the Lie algebra generators $T^I$ and coefficents $\omega$. However, I'm confused on how the connection Laplacian should be defined. Using the naive definition of 
$$\tilde{\Delta} = (\partial_\mu + \omega_\mu^I T^I)(\partial^\mu + \omega^{\mu I} T^I) = e^a_\mu(\partial_a + \omega_a^I T^I)e^{b\mu}(\partial_b + \omega_b^I T^I) $$
does not reproduce the right identity, and actually gives a differential operator rather than a scalar for $D\!\!\!/\ D\!\!\!/\  - \Delta$.
I have a feeling that the above expression for the Laplacian doesn't make sense for spinor fields, since implicitly I'm raising the spinor index using $g^{\mu\nu}$, whereas everything should in principle be doable via tetrads. What is the correct covariant way to write the operators $\nabla, \nabla^*$ and $\Delta$? 
****EDIT****:
I figured it out, and typed it out below.

A:

First, my expression for the Laplacian above was off, and the definition of $\nabla^*$ and $\nabla^*\nabla$ is a bit more subtle than I thought. $\nabla^*$ is the adjoint of $\nabla$. In particular, some vector $X$, $\nabla^*_X$ is the adjoint of $\nabla_X = X^\mu \nabla_\mu$. This can be computed from trying to integrate $\langle w, \nabla_X v \rangle = \int_M d(vol) w_\nu X^\mu \nabla_\mu v^\nu$ by parts to get 
$$\langle w, \nabla_X v \rangle = \int_M d(vol) w_\nu X^\mu \nabla_\mu v^\nu = -\int_M d(vol) v_\nu \nabla_\mu (X^\mu w^\nu) =: \langle \nabla^*_X w,  v \rangle$$
so that $\nabla_X^* w^\nu = -\nabla_\mu(X^\mu w^\nu)$ in local coordinates.
Note that all of these definitions were covariant derivatives with respect to a fixed orthonormal frame field. The operator $\nabla^*\nabla$ is actually referring to $tr(\nabla^*\nabla) =: \Delta$, where the trace means given some local orthonormal basis $e^\mu_a$, we define the laplacian as 
$$(\Delta w)^\mu = \sum_a(\nabla^*_{e_a}\nabla_{e_a}) = \sum_a -\nabla_\nu(e_a^\nu e_a^\rho \nabla_\rho w^\mu)$$
which we can rewrite as (using Freed's convention that $\{\gamma^a,\gamma^b\}= -2 \delta^{ab}$)
$$(\Delta w)^\mu = \frac{1}{2}\nabla_\nu(e_a^\nu \{\gamma^a,\gamma^b\} e_b^\rho \nabla_\rho w^\mu)$$
Now, we want to compare $D\!\!\!/\ D\!\!\!/\ $ with $\Delta$. If we work in Riemann normal coordinates centered at $x=0$, we note that it's possible to choose $e^\mu_a = \delta^\mu_a + O(x^2)$, so that $\omega_\mu = 0 + O(x)$. This means we can freely rewrite the above expression as 
$$(\Delta w)^\mu = \frac{1}{2} e_a^\nu \{\gamma^a,\gamma^b\} \nabla_\nu(e_b^\rho \nabla_\rho w^\mu) = \frac{1}{2} \{\gamma^\nu,\gamma^\rho\} \nabla_\nu \nabla_\rho w^\mu $$
where we abuse notation and write $\mu = 0,...,d$ to represent the flat space coordinates locally. And we also rewrite $D\!\!\!/\ D\!\!\!/\ $ as
$$(D\!\!\!/\ D\!\!\!/\ w)^\mu = \gamma^\nu \gamma^\rho \nabla_\nu \nabla_\rho w^\mu$$
so that 
$$D\!\!\!/\ D\!\!\!/\ - \Delta = -\frac{1}{2}[\gamma^\nu,\gamma^\rho] \nabla_\nu\nabla_\rho = -\frac{1}{2}\gamma^\nu \gamma^\rho [\nabla_{\nu},\nabla_{\rho}]$$
From here, it becomes straightforward, since $[\nabla_{\nu},\nabla_{\rho}]$ gives the curvature tensor. One tricky point is the eventual need for the identity $\gamma^i\gamma^j\gamma^k\gamma^l R_{ijkl} = 2R$ (which is explained in Freed's notes clearly). This can be shown by noting that the Riemann tensors symmetries require $i \neq j, k \neq l$. The case where $j \neq k \neq l$ are all distinct vanishes by the Bianchi identity $R_{ijkl} + R_{iklj} + R_{iljk} = 0$. The case where $j=k$ gives $\sum_j \gamma^i \gamma^j \gamma^j \gamma^l R_{ijjl} = - \gamma^i \gamma^l R_{il} = R$. the case $j=l$ similarly gives $R$, so the total sum is $2R$. Putting this all together will give the desired result after dealing appropriately with the factors of 2 from the spin connection.

