Let $s$ be a positive integer and $h_{a_1\dots a_s}$ be a traceless and totally symmetric (real) field which is defined modulo gauge transformations of the form $$\delta_{\xi}h_{a_1\dots a_s}=\partial_{(a_1}\xi_{a_2\dots a_s)} \tag{1}$$ where $\xi_{a_1\dots a_{s-1}}$ is also a traceless and totally symmetric (real) field. The indices range from $a=0,1, 2$ (I am working in 3D).
I have seen that using the gauge symmetry it is possible to impose the gauge condition $$\partial^bh_{ba_1\dots a_{s-1}}=0~. \tag{2}$$
I want to understand properly how this can be done, which is what I will attempt below.
In the gauge condition $(2)$ there are $2s-1$ conditions , which happens to be the precise number of functions embedded in the gauge parameter $\xi_{a(s-1)}$ (taking into account total symmetry and tracelessness). Therefore we can choose the components of $\xi_{a(s-1)}$ such that the gauge transformed $h_{a(s)}$ satisfies $(2)$. I think this is usually the reasoning that people employ. The problem I have with this is, how do we know that such a choice is always possible?
To explore this issue, we could investigate how the divergence of $h$ transforms under the gauge transformations, \begin{align} \delta_{\xi}\big(\partial^bh_{ba_1\dots a_{s-1}}\big)&=\partial^b\partial_{(b}\xi_{a_1\dots a_{s-1})}\\ &=\frac{1}{s}\square\xi_{a_1\dots a_{s-1}}+\frac{(s-1)}{s}\partial^b\partial_{(a_1}\xi_{a_2\dots a_{s-1})b} ~. \tag{3} \end{align} Now, if all of the degrees of freedom that were required to impose the gauge $(2)$ were contained entirely in the field $\xi_{a(s-1)}$, then I can see how one could always choose the $\xi_{a(s-1)}$ so that $h$ satisfies $(2)$. By this I mean that if $(3)$ were of the form \begin{align} \delta_{\xi}\big(\partial^bh_{ba_1\dots a_{s-1}}\big)\propto \square\xi_{a_1\dots a_{s-1}} \tag{4} \end{align} or something along those lines. But the second term on the right hand side of $(3)$ seems to complicate things, since all of the degrees of freedom are not contained wholly in $\xi_{a(s-1)}$.
Could it be that in order to impose the gauge $(2)$, we require also that $\xi_{a(s-1)}$ is divergenceless? Or am I misunderstanding how the gauge fixing procedure works?
I am also interested to know how this gauge fixing procedure would change, if at all, when we lift everything to curved space and promote the regular derivatives to covariant ones. Must there be some sort of integrability condition on the Christoffel symbols to allow this choice to always exist?
EDIT
It seems that if I assume that the momentum of the field does not lie on the lightcone (so that the operator $\square$ is invertible), then I can invert the differential equation to solve for $\xi_{a(s-1)}$.
So for example, for $s=2$, in the position representation we have the following: $$\delta h_{ab}=\partial_{(a}\xi_{b)}\implies h'_{ab}=h_{ab}+\partial_{(a}\xi_{b)}.$$ If we want to impose the gauge condition $\partial^ah'_{ab}=0$, then we need to 'solve' the PDE $$0=\partial^bh_{ab}+\frac 12 \square \xi_a+\frac 12 \partial_a\partial^b\xi_b.\tag{5}$$ Taking the divergence of (5) gives $$0=\partial^a\partial^bh_{ab}+\square\partial^a\xi_a \implies \partial^a\xi_a=-\frac{1}{\square}\partial^a\partial^bh_{ab}. \tag{6}$$ By $\frac{1}{\square}$ I mean the inverse operator of $\square$. Inserting the second part of (6) into (5) yields $$\xi_a=-\frac{2}{\square}\partial^bh_{ab}+\frac{1}{\square^2}\partial_a\partial^b\partial^ch_{bc}\tag{7}.$$ With this choice of $\xi_a$ one can check that $h'_{ab}$ is divergenceless.
One can repeat this process for higher-spins (greater $s$), it will just require a couple more steps (an extra step for each increment in $s$).
Comments on the above process are welcome, in particular, if this valid and if so what impact assuming $\square$ is invertible has on the type of functions, $\xi_{a(s-1)}$, we can consider. I should also note that I want to stay off the mass-shell.