Doing a phase-space portrait in matlab

• MATLAB
bagram
So I have this system of equations:

$$\binom{x_{n+1}}{y_{n+1}}=\begin{pmatrix}e^{r} & 0 \\ 0 & e^{-r} \end{pmatrix}\begin{pmatrix}cos(\phi+I_{n}) & -sin(\phi+I_{n}) \\ sin(\phi+I_{n}) & cos(\phi+I_{n}) \end{pmatrix}\begin{pmatrix}x_{n}\\ y_{n} \end{pmatrix}$$

where
$$I_{n}=x_{n}^2+y_{n}^2$$

I have no idea how to plot that in matlab as a phase-space portrait...

Any help would be great