In general, an ellipse in a general position $[h,k]$ (what I needed) is implicitly given as
$$\frac{(x-h)^2}{a^2} + \frac{(y-k)^2}{b^2} = 1$$
where $a$ and $b$ are semiaxes. Rotating the $x$ and $y$ coords yields (you're missing a $-$ sign in the first line)
$$x' = x\cos(\phi) - y\sin(\phi)$$
$$y' = x\sin(\phi) + y\cos(\phi)$$
Plugging this into ellipse equation, you get
$$\frac{([x\cos(\phi) - y\sin(\phi)]-h)^2}{a^2} + \frac{([x\sin(\phi) + y\cos(\phi)]-k)^2}{b^2} = 1$$
Now, you will need this in a $Ax^2 + 2Bxy + Cy^2 + 2Dx + 2Ey + F = 0$ form, as the ellipse matrix constitutes of the coefficients $A,B,C,D,E,$ and $F$. So you rewrite the previous, getting
$$-(a^2 b^2) + b^2 h^2 + a^2 k^2 - 2 b^2 h x \cos(\phi) - 2 a^2 k y \cos(\phi) + b^2 x^2 \cos(\phi)^2 + a^2 y^2 \cos(\phi)^2 - 2 a^2 k x \sin(\phi) + 2 b^2 h y \sin(\phi) + 2 a^2 x y \cos(\phi) \sin(\phi) - 2 b^2 x y \cos(\phi) \sin(\phi) + a^2 x^2 \sin(\phi)^2 + b^2 y^2 \sin(\phi)^2 = 1$$
The coefficients therefore are
$$\begin{array}
AA &=& b^2 \cos^2(\phi) + a^2 \sin^2(\phi) \\
B &=& a^2 \cos(\phi) \sin(\phi) - b^2 \cos(\phi) \sin(\phi)\\
C &=& a^2 \cos^2(\phi) + b^2 \sin^2(\phi)\\
D &=& -b^2 h \cos(\phi) - a^2 k \sin(\phi)\\
E &=& -a^2 k \cos(\phi) + b^2 h \sin(\phi)\\
F &=& -a^2 b^2 + b^2 h^2 + a^2 k^2\end{array}$$
and the final ellipse matrix
$$M = \begin{bmatrix}A & B & D \\
B & C & E \\
D & E & F
\end{bmatrix}$$
By multiplying $\begin{bmatrix}x&y&1\end{bmatrix}M\begin{bmatrix}x&y&1\end{bmatrix}^T$, you get the original $Ax^2 + 2Bxy + Cy^2 + 2Dx + 2Ey + F = 0$ equation, hence the matrix form.
