Another option is just to throw everything into the equation for a 2D rotated ellipse and see if the result is less than one.
So a point is inside the ellipse if the following inequality is true
Where (xp,yp) are the point coordinates and (x0, y0) is the center of the ellipse.
I implemented a small Mathematica program demonstrating that this indeed works:
Here it is in action:
And here is the code:
ellipse[x_, y_, a_, b_, [Alpha]_, x0_: 0, y0_: 0] :=
(((x - x0)*Cos[[Alpha]] + (y - y0)*Sin[[Alpha]])/a)^2
+ (((x - x0)*Sin[[Alpha]] - (y - y0)*Cos[[Alpha]])/b)^2;
Manipulate[
RegionPlot[
ellipse[x, y, a, b, [Alpha] [Degree], Sequence @@ pos] < 1, {x, -5, 5}, {y, -5, 5},
PlotStyle -> If[ellipse[Sequence @@ p, a, b, [Alpha] [Degree], Sequence @@ pos] <= 1, Orange, LightBlue],
PlotPoints -> 25]
, {{a, 2}, 1, 5, Appearance -> "Labeled"}
, {{b, 4}, 2, 5, Appearance -> "Labeled"}
, {[Alpha], 0, 180, Appearance -> "Labeled"}
, {{p, {3, 1}}, Automatic, ControlType -> Locator}
, {{pos, {0, 0}}, Automatic, ControlType -> Locator}]