main/src/georegression/geometry/UtilEllipse_F64.java
public static boolean tangentLines(Point2D_F64 pt , EllipseRotated_F64 ellipse ,
									   Point2D_F64 tangentA , Point2D_F64 tangentB )
	{
		// Derivation:
		// Compute the tangent at only point along the ellipse by computing dy/dx
		//    x*b^2/(y*a^2) or - x*b^2/(y*a^2)  are the possible solutions for the tangent
		// The slope of the line and the gradient are the same, so this is true:
		//   y - y'     -x*b^2
		//  -------  =  -------
		//   x - x'      y*a^2
		//
		//  (x,y) is point on ellipse, (x',y') is pt that lines pass through
		//
		//  that becomes
		//  y^2*a^2 + x^2*b^2 = x'*x*b^2 + y'*y*a^2
		//  use the equation for the ellipse (centered and aligned at origin)
		//  a^2*b^2 =  x'*x*b^2 + y'*y*a^2
		//
		// solve for y
		// plug into ellipse equation
		// solve for x, which is a quadratic equation

		// translate and rotate into ellipse reference frame
		double cphi = Math.cos(ellipse.phi);
		double sphi = Math.sin(ellipse.phi);

		double tmpx = pt.x - ellipse.center.x;
		double tmpy = pt.y - ellipse.center.y;

		double xt =  tmpx*cphi + tmpy*sphi;
		double yt = -tmpx*sphi + tmpy*cphi;

		// solve
		double a2 = ellipse.a*ellipse.a;
		double b2 = ellipse.b*ellipse.b;

		// quadratic equation for the two variants.
		// solving for x
		double aa0 = yt*yt/b2 + xt*xt/a2;
		double bb0 = -2.0*xt;
		double cc0 = a2*(1.0-yt*yt/b2);

		double descriminant0 = bb0*bb0 - 4.0*aa0*cc0;

		// solving for y
		double aa1 = xt*xt/a2 + yt*yt/b2;
		double bb1 = -2.0*yt;
		double cc1 = b2*(1.0-xt*xt/a2);

		double descriminant1 = bb1*bb1 - 4.0*aa1*cc1;

		double x0,y0, x1,y1;
		if( descriminant0 < 0 && descriminant1 < 0 ) {
			return false;
		} else if( descriminant0 > descriminant1 ) {
			if( yt == 0 )
				return false;

			double right = Math.sqrt(descriminant0);

			x0 = (-bb0 + right)/(2.0*aa0);
			x1 = (-bb0 - right)/(2.0*aa0);

			y0 = b2/yt - xt*x0*b2/(yt*a2);
			y1 = b2/yt - xt*x1*b2/(yt*a2);

		} else {
			if( xt == 0 )
				return false;

			double right = Math.sqrt(descriminant1);

			y0 = (-bb1 + right)/(2.0*aa1);
			y1 = (-bb1 - right)/(2.0*aa1);

			x0 = a2/xt - yt*y0*a2/(xt*b2);
			x1 = a2/xt - yt*y1*a2/(xt*b2);
		}

		// convert the lines back into world space
		tangentA.x = x0*cphi - y0*sphi + ellipse.center.x;
		tangentA.y = x0*sphi + y0*cphi + ellipse.center.y;

		tangentB.x = x1*cphi - y1*sphi + ellipse.center.x;
		tangentB.y = x1*sphi + y1*cphi + ellipse.center.y;

		return true;
	}
public static boolean tangentLines(Point2D_F64 pt , EllipseRotated_F64 ellipse ,
									   Point2D_F64 tangentA , Point2D_F64 tangentB )
	{
		// Derivation:
		// Compute the tangent at only point along the ellipse by computing dy/dx
		//    x*b^2/(y*a^2) or - x*b^2/(y*a^2)  are the possible solutions for the tangent
		// The slope of the line and the gradient are the same, so this is true:
		//   y - y'     -x*b^2
		//  -------  =  -------
		//   x - x'      y*a^2
		//
		//  (x,y) is point on ellipse, (x',y') is pt that lines pass through
		//
		//  that becomes
		//  y^2*a^2 + x^2*b^2 = x'*x*b^2 + y'*y*a^2
		//  use the equation for the ellipse (centered and aligned at origin)
		//  a^2*b^2 =  x'*x*b^2 + y'*y*a^2
		//
		// solve for y
		// plug into ellipse equation
		// solve for x, which is a quadratic equation

		// translate and rotate into ellipse reference frame
		double cphi = Math.cos(ellipse.phi);
		double sphi = Math.sin(ellipse.phi);

		double tmpx = pt.x - ellipse.center.x;
		double tmpy = pt.y - ellipse.center.y;

		double xt =  tmpx*cphi + tmpy*sphi;
		double yt = -tmpx*sphi + tmpy*cphi;

		// solve
		double a2 = ellipse.a*ellipse.a;
		double b2 = ellipse.b*ellipse.b;

		// quadratic equation for the two variants.
		// solving for x
		double aa0 = yt*yt/b2 + xt*xt/a2;
		double bb0 = -2.0*xt;
		double cc0 = a2*(1.0-yt*yt/b2);

		double descriminant0 = bb0*bb0 - 4.0*aa0*cc0;

		// solving for y
		double aa1 = xt*xt/a2 + yt*yt/b2;
		double bb1 = -2.0*yt;
		double cc1 = b2*(1.0-xt*xt/a2);

		double descriminant1 = bb1*bb1 - 4.0*aa1*cc1;

		double x0,y0, x1,y1;
		if( descriminant0 < 0 && descriminant1 < 0 ) {
			return false;
		} else if( descriminant0 > descriminant1 ) {
			if( yt == 0 )
				return false;

			double right = Math.sqrt(descriminant0);

			x0 = (-bb0 + right)/(2.0*aa0);
			x1 = (-bb0 - right)/(2.0*aa0);

			y0 = b2/yt - xt*x0*b2/(yt*a2);
			y1 = b2/yt - xt*x1*b2/(yt*a2);

		} else {
			if( xt == 0 )
				return false;

			double right = Math.sqrt(descriminant1);

			y0 = (-bb1 + right)/(2.0*aa1);
			y1 = (-bb1 - right)/(2.0*aa1);

			x0 = a2/xt - yt*y0*a2/(xt*b2);
			x1 = a2/xt - yt*y1*a2/(xt*b2);
		}

		// convert the lines back into world space
		tangentA.x = x0*cphi - y0*sphi + ellipse.center.x;
		tangentA.y = x0*sphi + y0*cphi + ellipse.center.y;

		tangentB.x = x1*cphi - y1*sphi + ellipse.center.x;
		tangentB.y = x1*sphi + y1*cphi + ellipse.center.y;

		return true;
	}
<p>Finds two points on the ellipse that in combination with point 'pt' each define
a line that is tangent to the ellipse.</p>

Notes:<br>
Point 'pt' is assumed to be outside of the ellipse.

@param pt Point which the lines will pass though
@param ellipse The ellipse which the lines will be tangent to
@param tangentA (output) Point on the ellipse where tangent line A hits it
@param tangentB (output) Point on the ellipse where tangent line B hits it