2D system of hyperbolic equation (LeVeque) Again I need help. This time it's exercise 18.1 from LeVeque Finite Volume Methods.

Consider the system $q_t+Aq_x+Bq_y=0$ with $A=\begin{pmatrix}3&1\\1&3\end{pmatrix}$ and $B=\begin{pmatrix}0&2\\2&0\end{pmatrix}$.
Show that these matrices are simultaneously diagonalizible and determine the general solution to this system with arbitrary initial data. In particular sketch how the solution evolves in the $x-y$ plane with data
$$\begin{align*}
&q^1(x,y,0)=\left\{\begin{array}{ll}1&\text{ if } x^2+y^2\leq 1\\0&\text{ otherwise }\end{array}\right.\\
&q^2(x,y,0)\equiv 0
\end{align*}$$

So far I've got $A=R\Lambda^xR^{-1}$ and $B=R\Lambda^yR^{-1}$ with
$$\begin{align*}
&R&&=\begin{pmatrix}1&1\\1&-1\end{pmatrix}\\
&R^{-1}&&=\frac{1}{2}\begin{pmatrix}1&1\\1&-1\end{pmatrix}\\
&\Lambda^x&&=\begin{pmatrix}4&0\\0&2\end{pmatrix}\\
&\Lambda^y&&=\begin{pmatrix}2&0\\0&-2\end{pmatrix}
\end{align*}$$
but now I don't know how to move on. Any tips are appreciated :)

Setting $p = R^{-1}q$, we have
$$
p_t + \Lambda_x p_x + \Lambda_y p_y = 0 \, .
$$
Since the $\Lambda^\alpha$ matrices are diagonal, this system is decoupled.
Applying the method of characteristics componentwise yields
$$
p(x,y,t) = R^{-1} q(x,y,t) = \begin{bmatrix}
F_1(\lambda^x_1 x -t, \lambda^y_1 x - \lambda^x_1 y) \\
F_2(\lambda^x_2 x -t, \lambda^y_2 x - \lambda^x_2 y)
\end{bmatrix}
$$
where $F_1$, $F_2$ are two arbitrary functions, and $\lambda_i^\alpha$ denotes the diagonal entries of $\Lambda_i^\alpha$ (see this post, where the notations correspond to $\Phi = F_i$, $a=1$, $b=\lambda_i^x$, and $c=\lambda_i^y$). Now, it remains to apply the intial condition to determine the arbitrary functions:
$$
q(x,y,0) = R\, p(x,y,0) = \begin{bmatrix}
F_1(4 x, 2 x - 4 y) + F_2(2 x, -2 (x + y)) \\
F_1(4 x, 2 x - 4 y) - F_2(2 x, -2 (x + y))
\end{bmatrix} .
$$
From the initial condition $q_1(x,y,0) = \Bbb I_{x^2 + y^2 < 1}$ and $q_2(x,y,0) = 0$ in OP, we deduce that
\begin{aligned}
\tfrac12\Bbb I_{x^2 + y^2 < 1} &= F_2(2x,-2(x+y)) \\
& = F_2(X,Y) = \tfrac12\Bbb I_{X^2 + (Y+X)^2 < 4} \\
\tfrac12\Bbb I_{x^2 + y^2 < 1} &= F_1(4x,2 x - 4 y) \\
& = F_1(\xi,\eta) = \tfrac12\Bbb I_{\xi^2 + (\eta-\xi/2)^2 < 16}
\end{aligned}
where $\Bbb I$ is the indicator function. Finally,
$$
q(x,y,t) = \frac12 \begin{bmatrix}
\Bbb I_{(x-t/4)^2 + (y-t/8)^2 < 1} + \Bbb I_{(x-t/2)^2 + (y+t/2)^2 < 1} \\
\Bbb I_{(x-t/4)^2 + (y-t/8)^2 < 1} - \Bbb I_{(x-t/2)^2 + (y+t/2)^2 < 1}
\end{bmatrix}
$$
with the above expressions.
