I'll describe a computational method, assuming the rotation has been translated to be around an axis through the origin.
If you know the axis of rotation $A=(a,b,c)$, then you can find a vector orthogonal to this one. For example, if $a \ne 0$, $b \ne 0$, then $V=(b,-a,0)$ is such a vector. Or, if one component is $0$--say, for example $c=0$, then $V=(0,0,1)$ is orthogonal. And that covers all cases.
A third vector in a right-handed triple is found using the vector cross product: $W=A\times V$. Normalize $V$ and $W$ to unit vectors $\hat{V}$, $\hat{W}$.
Now, feed $\hat{V}$ into your matrix and determine the output in terms of $\hat{V}$, $\hat{W}$ using dot products
$$
R\hat{V} = \{(R\hat{V})\cdot\hat{V}\}\hat{V}+\{(R\hat{V})\cdot\hat{W}\}\hat{W}=\alpha\hat{V}+\beta\hat{W},
$$
and compare to the action of a right-handed rotation about $A$ through an angle $\theta$:
$$
\hat{V} \mapsto \cos\theta \hat{V}+\sin\theta \hat{W}.
$$
If $T$ denotes your translation by $(2,0)$, and $R$ your rotation by $90^\circ$ around the origin, then your map is $Q=R\circ T$. You can find $Q$ either analytically or geometrically.
Analytic solution
$T$ and $R$ are defined by
$$
T(x)=x+(2,0)=(x_1+2,x_2), \, R(x)=(-x_2,x_1) \quad \forall x \in \mathbb{R}^2.
$$
Therefore
$$
Q(x)=R(T(x))=R(x_1+2,x_2)=(-x_2,x_1+2) \quad \forall x=(x_1,x_2)\in \mathbb{R}^2
$$
Notice that $Q$ has exactly one fixed point, i.e. $a=(-1,1)$. In fact $Q(x)=x \iff x=(-1,1)$ and therefore $Q$ is a rotation.
For every $x\in \mathbb{R}^2$ we have
$$
(x-a)\cdot(Q(x)-a)=(x_1+1,x_2-1)\cdot(-x_2+1,x_1+1)=0,
$$
i.e. $(x-a)\perp (Q(x)-a)$, and therefore $Q$ is a rotation by $90^\circ$ around the point $a=(-1,1)$.
Geometric solution
If $\Delta\subset \mathbb{R}^2$ is a straight line, we denote by $S_\Delta$ the reflection about $\Delta$. We can then write
$$
T=S_{L_0}\circ S_{L_1}, \quad R=S_{L_2}\circ S_{L_0},
$$
with
\begin{eqnarray}
L_0&=&\{(x_1,x_2)\in \mathbb{R}^2:\, x_1=0 \},\\
L_1&=&\{(x_1,x_2)\in \mathbb{R}^2:\, x_1=1 \},\\
L_2&=&\{(x_1,x_2)\in \mathbb{R}^2:\, x_2=-x_1 \}
\end{eqnarray}
We have
$$
Q=R\circ T=(S_{L_2}\circ S_{L_0})\circ(S_{L_0}\circ S_{L_1})=S_{L_2}\circ (\underbrace{S_{L_0}\circ S_{L_0}}_{\mbox{identity}})\circ S_{L_1}=S_{L_2}\circ S_{L_1}
$$
Since the two lines $L_1$ and $L_2$ intersect at the point $a=(1,-1)$ and the angle $(L_1,L_2)$ between $L_1$ and $L_2$ is $45^\circ$, it follows that $Q$ is a rotation by $2\times45^\circ=90^\circ$ around the point $a=(-1,1)$.
Best Answer
If the covariance is the $2\times2$ identity matrix, then the density is $$e^{−(x^2+y^2)/2}$$ multiplied by a suitable normalizing constant. If $\begin{bmatrix} X \\ Y \end{bmatrix}$ is a random vector with this distribution, then you rotate that random vector by multiplying on the left by a typical $2\times 2$ orthogonal matrix: $$ G \begin{bmatrix} X \\ Y \end{bmatrix} = \begin{bmatrix} \cos\theta & -\sin\theta \\ \sin\theta & \cos\theta \end{bmatrix}\begin{bmatrix} X \\ Y \end{bmatrix}. $$
If the question is how to "rotate" the probability distribution, then asnwer is that it's invariant under rotations about the origin since it depends on $x$ and $y$ only through the distance $\sqrt{x^2+y^2}$ from the origin to $(x,y)$.
If you multiply on the left by a $k\times2$ matrix $G$, you have $$ \mathbb{E}\left(G\begin{bmatrix} X \\ Y \end{bmatrix}\right) = G\mathbb{E}\begin{bmatrix} X \\ Y \end{bmatrix} $$ and $$ \operatorname{var}\left( G \begin{bmatrix} X \\ Y \end{bmatrix} \right) = G\left(\operatorname{var}\begin{bmatrix} X \\ Y \end{bmatrix}\right)G^T, $$ a $k\times k$ matrix. If the variance in the middle is the $2\times2$ identity matrix and $G$ is the $2\times 2$ orthogonal matrix given above, then it's easy to see that the variance is $$ GG^T $$ and that is just the $2\times 2$ identity matrix. The only fact you need after that is that if you multiply a multivariate normal random vector by a matrix, what you get is still multivariate normal. I'll leave the proof of that as an exercise.