Finally, let's compute the spectral measure of $X_\theta - 2 \tan(\theta) P$, which we will call $\rho$. This can be done using rank-one perturbation formulas. Alternatively, $X_\theta$ and $P$ are Boolean independent (see https://mast.queensu.ca/~speicher/papers/boolean.ps). The spectral measure of $-2 \tan(\theta) P$ with respect to $\delta_e$ is $\delta_{-2 \tan(\theta)}$. The theory of Boolean independence tells us that $$ F_\rho(z) - z = (F_\nu(z) - z) + (F_{\delta_{-2\tan(\theta)}}(z) - z) = F_\nu(z) - z + 2 \tan(\theta). $$ Hence, $$ F_\rho(z) = \sqrt{(z + 8 \tan \theta)^2 - 4} + 2 \tan(\theta). $$ Thus, for $z$ in the upper half-plane, $$ G_\rho(z) = \frac{1}{\sqrt{(z + 8 \tan \theta)^2 - 4} + 2 \tan(\theta)}, $$ where $G_\rho(z) = \int (z - t)^{-1} \,d\rho(t)$. Using the fact that this is rank-one perturbation that the norm of the new operator $X_\theta - 2 \tan(\theta) P$ is the maximum of the norm of $X_\theta - 2 \tan(\theta) P$ on the cyclic subspace generated by $\delta_e$ and the norm on the orthogonal complement, which will be bounded by $\lVert X_\theta \rVert$. Thus, to show that the two operators have different norms, it suffices to show that the support radius of $\rho$ is strictly larger than that of $\nu$. Looking at $G_\rho$, the support of the measure $\rho$ agrees with the support of $\nu$ except with the addition of atoms at the points where $(2 \tan \theta)^2 = (z + 8 \tan \theta)^2 - 4$. From this you can test whether the norms agree.

Finally, let's compute the spectral measure of $X_\theta - 2 \tan(\theta) P$, which we will call $\rho$. This can be done using rank-one perturbation formulas. Alternatively, $X_\theta$ and $P$ are Boolean independent (see https://mast.queensu.ca/~speicher/papers/boolean.ps). The spectral measure of $-2 \tan(\theta) P$ with respect to $\delta_e$ is $\delta_{-2 \tan(\theta)}$. The theory of Boolean independence tells us that $$ F_\rho(z) - z = (F_\nu(z) - z) + (F_{\delta_{-2\tan(\theta)}}(z) - z) = F_\nu(z) - z + 2 \tan(\theta). $$ Hence, $$ F_\rho(z) = \sqrt{(z + 8 \tan \theta)^2 - 4} + 2 \tan(\theta). $$ Thus, for $z$ in the upper half-plane, $$ G_\rho(z) = \frac{1}{\sqrt{(z + 8 \tan \theta)^2 - 4} + 2 \tan(\theta)}, $$ where $G_\rho(z) = \int (z - t)^{-1} \,d\rho(t)$. Using the fact that this is rank-one perturbation that the norm of the new operator $X_\theta - 2 \tan(\theta) P$ is the maximum of the norm of $X_\theta - 2 \tan(\theta) P$ on the cyclic subspace generated by $\delta_e$ and the norm on the orthogonal complement, which will be bounded by $\lVert X_\theta \rVert$. Thus, to show that the two operators have different norms, it suffices to show that the support radius of $\rho$ is strictly larger than that of $\nu$. 

Looking at $G_\rho$, the support of the measure $\rho$ agrees with the support of $\nu$ except with the addition of atoms at the point where $(2 \tan \theta)^2 = (z + 8 \tan \theta)^2 - 4$ according to the version of the square root described above. Assume $\tan \theta > 0$ since the other case is symmetric (and recall $\theta$ was assumed to be in $(-\pi/2,\pi/2)$). Then we are looking at a point on the negative side where $$(z + 8 \tan \theta)^2 = 4 + 4 \tan^2 \theta = 4 \sec^2 \theta.$$ Hence, $z = -2 \sec \theta - 8 \tan \theta$ is the location of the atom. Since the support of $\nu$ was already farther on the negative side, we have $$ \lVert X_\theta \rVert = 2 + 8 |\tan \theta| < 2 \sec \theta + 8 |\tan \theta| = \lVert X_\theta - (2 \tan \theta) P \rVert. $$

One can use non-commutative probability theory to compute the spectral distributions of $X_\theta$ and $X_\theta - 2 \tan(\theta) P$ with respect to $\delta_e$. Since $G = \mathbb{Z}_2 * \mathbb{Z}_2$, the operators $T_s$ and $T_t$ are freely independent with respect to the vector state given by $\delta_e$. Each of them has the spectral distribution given by the Bernoulli measure $(1/2)(\delta_{-1} + \delta_1)$. Thus, the spectral distribution of $T_s + T_t$ is the free convolution of two Bernoulli distributions. Using the $R$-transform (see Chapter 4 of https://arxiv.org/pdf/1908.08125.pdf), we can compute that the free convolution is the measure $\mu$ satisfying $$ F_\mu(z) = \sqrt{z^2 - 4}, $$ where $F_\mu(z) = 1 / \int (z - t)^{-1}\,d\mu(t)$, also known as the $F$-transform. Here we use the square root which is defined on the upper half-plane and close to $z$ when $z \to \infty$. This means that $\mu$ is an arcsine law $$ d\mu(x) = \frac{1}{\pi \sqrt{4 - x^2}} \chi_{(-2,2)}(x)\,dx $$ This computation is known; see the last full paragraph at the bottom of page 2 here: https://arxiv.org/pdf/1008.5205.pdf.

Next, let $\nu$ be the spectral measure of $T_s + T_t - 8 \tan(\theta)$. This just shifts the measure $\mu$ by $-8 \tan(\theta)$, and we have $F_\nu(z) = F_\mu(z + 8 \tan(\theta))$. Since the support of $\mu$ is $[-2,2]$, the support of $\nu$ is $[-2 - 8 \tan(\theta), 2 - 8 \tan(\theta)]$ and thus the norm of $X_\theta$ is $2 + 8 |\tan \theta|$, at least on the cyclic subspace generated by $\delta_e$. However, because $\delta_e$ is cyclic for the right shift operators by $s$ and $t$, it is cyclic for the commutant of $X_\theta$ and hence the norm of $X_\theta$ on the entire space agrees with the norm on the cyclic subspace generated by $\delta_e$.

Finally, let's compute the spectral measure of $X_\theta - 2 \tan(\theta) P$, which we will call $\rho$. This can be done using rank-one perturbation formulas. Alternatively, $X_\theta$ and $P$ are Boolean independent (see https://mast.queensu.ca/~speicher/papers/boolean.ps). The spectral measure of $-2 \tan(\theta) P$ with respect to $\delta_e$ is $\delta_{-2 \tan(\theta)}$. The theory of Boolean independence tells us that $$ F_\rho(z) - z = (F_\nu(z) - z) + (F_{\delta_{-2\tan(\theta)}}(z) - z) = F_\nu(z) - z + 2 \tan(\theta). $$ Hence, $$ F_\rho(z) = \sqrt{(z + 8 \tan \theta)^2 - 4} + 2 \tan(\theta). $$ Thus, for $z$ in the upper half-plane, $$ G_\rho(z) = \frac{1}{\sqrt{(z + 8 \tan \theta)^2 - 4} + 2 \tan(\theta)}, $$ where $G_\rho(z) = \int (z - t)^{-1} \,d\rho(t)$. Using the fact that this is rank-one perturbation that the norm of the new operator $X_\theta - 2 \tan(\theta) P$ is the maximum of the norm of $X_\theta - 2 \tan(\theta) P$ on the cyclic subspace generated by $\delta_e$ and the norm on the orthogonal complement, which will be bounded by $\lVert X_\theta \rVert$. Thus, to show that the two operators have different norms, it suffices to show that the support radius of $\rho$ is strictly larger than that of $\nu$. Looking at $G_\rho$, the support of the measure $\rho$ agrees with the support of $\nu$ except with the addition of atoms at the points where $(2 \tan \theta)^2 = (z + 8 \tan \theta)^2 - 4$. From this you can test whether the norms agree.

Similar computations can be done for $\mathbb{Z}_2^{*n}$. Let $T_1$, \dots, $T_n$ be the operators corresponding to the generators. Then each $T_1 + \dots + T_n$ as a spectral distribution $\mu_n$ is the free convolution of $n$ Bernoulli laws. This can be done directly with the $R$-transform, or using Proposition 3.1 of this paper: https://arxiv.org/pdf/math/0703295.pdf. This results in $$ F_{\mu_n}(z) = \frac{(n-2)z + n \sqrt{z^2 - 4(n-1)}}{2(n-1)}. $$ Thus, similar computations can be done as in the $n = 2$ case.

One can use non-commutative probability theory to compute the spectral distributions of $X_\theta$ and $X_\theta - 2 \tan(\theta) P$ with respect to $\delta_e$. Since $G = \mathbb{Z}_2 * \mathbb{Z}_2$, the operators $T_s$ and $T_t$ are freely independent with respect to the vector state given by $\delta_e$. Each of them has the spectral distribution given by the Bernoulli measure $(1/2)(\delta_{-1} + \delta_1)$. Thus, the spectral distribution of $T_s + T_t$ is the free convolution of two Bernoulli distributions. Using the $R$-transform (see Chapter 4 of https://arxiv.org/pdf/1908.08125.pdf), we can compute that the free convolution is the measure $\mu$ satisfying $$ F_\mu(z) = \sqrt{z^2 - 4}, $$ where $F_\mu(z) = 1 / \int (z - t)^{-1}\,d\mu(t)$, also known as the $F$-transform. Here we use the square root which is defined on the upper half-plane and close to $z$ when $z \to \infty$. This means that $\mu$ is an arcsine law $$ d\mu(x) = \frac{1}{\pi \sqrt{4 - x^2}} \chi_{(-2,2)}(x)\,dx $$ This computation is known; see the last full paragraph at the bottom of page 2 here: https://arxiv.org/pdf/1008.5205.pdf.

Next, let $\nu$ be the spectral measure of $T_s + T_t - 8 \tan(\theta)$. This just shifts the measure $\mu$ by $-8 \tan(\theta)$, and we have $F_\nu(z) = F_\mu(z + 8 \tan(\theta))$. Since the support of $\mu$ is $[-2,2]$, the support of $\nu$ is $[-2 - 8 \tan(\theta), 2 - 8 \tan(\theta)]$ and thus the norm of $X_\theta$ is $2 + 8 |\tan \theta|$, at least on the cyclic subspace generated by $\delta_e$. However, because $\delta_e$ is cyclic for the right shift operators by $s$ and $t$, it is cyclic for the commutant of $X_\theta$ and hence the norm of $X_\theta$ on the entire space agrees with the norm on the cyclic subspace generated by $\delta_e$.

Finally, let's compute the spectral measure of $X_\theta - 2 \tan(\theta) P$, which we will call $\rho$. This can be done using rank-one perturbation formulas. Alternatively, $X_\theta$ and $P$ are Boolean independent (see https://mast.queensu.ca/~speicher/papers/boolean.ps). The spectral measure of $-2 \tan(\theta) P$ with respect to $\delta_e$ is $\delta_{-2 \tan(\theta)}$. The theory of Boolean independence tells us that $$ F_\rho(z) - z = (F_\nu(z) - z) + (F_{\delta_{-2\tan(\theta)}}(z) - z) = F_\nu(z) - z + 2 \tan(\theta). $$ Hence, $$ F_\rho(z) = \sqrt{(z + 8 \tan \theta)^2 - 4} + 2 \tan(\theta). $$ Thus, for $z$ in the upper half-plane, $$ G_\rho(z) = \frac{1}{\sqrt{(z + 8 \tan \theta)^2 - 4} + 2 \tan(\theta)}, $$ where $G_\rho(z) = \int (z - t)^{-1} \,d\rho(t)$. Using the fact that this is rank-one perturbation that the norm of the new operator $X_\theta - 2 \tan(\theta) P$ is the maximum of the norm of $X_\theta - 2 \tan(\theta) P$ on the cyclic subspace generated by $\delta_e$ and the norm on the orthogonal complement, which will be bounded by $\lVert X_\theta \rVert$. Thus, to show that the two operators have different norms, it suffices to show that the support radius of $\rho$ is strictly larger than that of $\nu$. Looking at $G_\rho$, the support of the measure $\rho$ agrees with the support of $\nu$ except with the addition of atoms at the points where $(2 \tan \theta)^2 = (z + 8 \tan \theta)^2 - 4$. From this you can test whether the norms agree.

Similar computations can be done for $\mathbb{Z}_2^{*n}$. Let $T_1$, \dots, $T_n$ be the operators corresponding to the generators. Then each $T_1 + \dots + T_n$ as a spectral distribution $\mu_n$ is the free convolution of $n$ Bernoulli laws. This can be done directly with the $R$-transform, or using Proposition 3.1 of this paper: https://arxiv.org/pdf/math/0703295.pdf. This results in $$ F_{\mu_n}(z) = \frac{(n-2)z + n \sqrt{z^2 - 4(n-1)}}{2(n-1)}. $$ Thus, similar computations can be done as in the $n = 2$ case.

Unfortunately, for general groups, it is not as easy to compute because the operators might not be self-adjoint. However, the tools of free probability could still be used in theory to compute the spectral distribution of $X^*X$.