using System;

using System.Drawing;

 

namespace Geometry

{

      /// <summary>

      /// Geometric solutions for circles

      /// </summary>

      public class CircleGeometry

      {

            private CircleGeometry()

            {

            }

           

            /// <summary>

            /// Find circles tangent to two lines and a circle

            /// </summary>

            /// <param name="pt11">First point on first line</param>

            /// <param name="pt12">Second point on first line</param>

            /// <param name="pt21">First point on second line</param>

            /// <param name="pt22">Second point on second line</param>

            /// <param name="ptCentre">Centre of circle</param>

            /// <param name="nRadius">Radius of circle</param>

            /// <param name="Centres">Array to receive list of centre points</param>

            /// <param name="Radii">Array to receive list of radii</param>

            /// <returns>Number of solutions</returns>

            public static int Circle2LC(Point pt11, Point pt12, Point pt21, Point pt22,

                  Point ptCentre, int nRadius, ref Point[] Centres, ref int[] Radii)

                  //*************************************************************************

                  // Find circles tangent to 2 lines and a circle

            {

 

                  // Initialise

                  int nSolutions = 0;

                  Centres = new Point[0];

                  Radii = new int[0];

 

                  // Translate so first circle is on origin

                  Point pt11N = new Point(pt11.X - ptCentre.X, pt11.Y - ptCentre.Y);

                  Point pt12N = new Point(pt12.X - ptCentre.X, pt12.Y - ptCentre.Y);

                  Point pt21N = new Point(pt21.X - ptCentre.X, pt21.Y - ptCentre.Y);

                  Point pt22N = new Point(pt22.X - ptCentre.X, pt22.Y - ptCentre.Y);

 

                  // If first line not vertical

                  if(pt11.X != pt12.X)

                  {

                        nSolutions = NormCircle2LC(pt11N, pt12N, pt21N, pt22N, nRadius, ref Centres, ref Radii);

                  }

 

                        // If second line not vertical

                  else if(pt21.X != pt22.X)

                  {

                        nSolutions = NormCircle2LC(pt21N, pt22N, pt11N, pt12N, nRadius, ref Centres, ref Radii);

                  }

 

                        // If both lines vertical, special case

                  else

                  {

                        // Centres are midway between the lines, and the radius is half the separation

                        Point ptTanCirc = new Point();

                        ptTanCirc.X = (int)Math.Round((pt11N.X + pt21N.X) / 2.0);

                        int radTanCirc = (int)Math.Round(Math.Abs(pt11N.X - pt21N.X) / 2.0);

 

                        // Calculate +/- y-coordinate of centres

                        double ycSquare = ((radTanCirc + nRadius) * (radTanCirc + nRadius)) -

                              (ptTanCirc.X * ptTanCirc.X);

                        int yc = (int)Math.Round(Math.Sqrt(ycSquare));

                        ptTanCirc.Y = ptCentre.Y + yc;

 

                        // Add solutions

                        nSolutions = AddSolution(ptTanCirc, radTanCirc, ref Centres, ref Radii);

                        if(yc != 0)

                        {

                              ptTanCirc.Y = ptCentre.Y - yc;

                              nSolutions = AddSolution(ptTanCirc, radTanCirc, ref Centres, ref Radii);

                        }

                  }

 

                  // Transform results back to original coordinate system

                  for(int iSolution=0; iSolution<nSolutions; ++iSolution)

                  {

                        Centres[iSolution].X += ptCentre.X;

                        Centres[iSolution].Y += ptCentre.Y;

                  }

                 

                  // Return number of solutions

                  return nSolutions;

            }

           

            private static int NormCircle2LC(Point pt11, Point pt12, Point pt21, Point pt22, int nRadius,

                                                            ref Point[] Centres, ref int[] Radii)

                  //*******************************************************************************************

                  // Find circles tangent to 2 lines and a circle centred on (0, 0)

                  //     pt11.X <> pt12.X. (If line is vertical, swap parameters in call)

            {

                  double a1, b1, c1, a2, b2, c2, b24ac, b24acRoot;

                  double t1, t2, r3, u, w, s, A, B, C, fRadius;

 

                  // Initialise

                  int nSolutions = 0;

 

                  // Confirm first line not vertical (so b1 <> 0)

                  if(pt11.X == pt12.X) return 0;

 

                  // Get normalised line equations

                  a1 = 0.0;   b1 = 0.0;   c1 = 0.0;

                  if(!LineGeometry.Equation(pt11, pt12, ref a1, ref b1, ref c1)) return 0;

                  a2 = 0.0;   b2 = 0.0;   c2 = 0.0;

                  if(!LineGeometry.Equation(pt21, pt22, ref a2, ref b2, ref c2)) return 0;

 

                  // Take all combinations of t(n) to get circles on other side of lines

                  // and nRadius to get circles on other side of given circle

                  for(int iCase=0; iCase<8; ++iCase)

                  {

                        t1 = ((iCase & 1) == 0) ? 1 : -1;

                        t2 = ((iCase & 2) == 0) ? 1 : -1;

                        r3 = ((iCase & 4) == 0) ? nRadius : -nRadius;

 

                        // Derive constants

                        u = (t1 * b2) - (t2 * b1);

                        w = (b1 * c2) - (b2 * c1);

                        s = (a1 * b2) - (a2 * b1);

                        A = (u*u) - (2*a1*s*u*t1) + (t1*t1*s*s) - (b1*b1*s*s);

                        B = 2.0 * ((u*w) + (c1*a1*s*u) - (a1*s*t1*w) - (c1*t1*s*s) - (r3*b1*b1*s*s));

                        C = (w*w) + (2*a1*s*c1*w) + (c1*c1*s*s) - (b1*b1*r3*r3*s*s);

 

                        // Calculate radii as roots of a quadratic

                        b24ac = (B * B) - (4 * A * C);

                        if((b24ac < 0.0) || (A == 0.0)) continue;

                        b24acRoot = Math.Sqrt(b24ac);

                        for(int iRoot=0; iRoot<2; ++iRoot)

                        {

                              if((b24ac == 0.0) && (iRoot > 0)) continue;

                              if(iRoot == 0) fRadius = (-B + b24acRoot)/(A + A);

                              else fRadius = (-B - b24acRoot)/(A + A);

                              if((fRadius <= 0.0) || (fRadius > (double)Int32.MaxValue)) continue;

                              int radTanCirc = (int)Math.Round(fRadius);

                             

                              // Derive x-coordinate of centre

                              Point ptCentre = new Point();

                              double xc, yc;

                              if(s != 0.0)

                              {

                                    // Calculate centre where lines not parallel

                                    xc = (fRadius * u + w) / s;

                                    if(Math.Abs(xc) > (double)Int32.MaxValue) continue;

                                    ptCentre.X = (int)Math.Round(xc);

                                    yc = ((-a1 * xc) - c1 + (fRadius * t1)) / b1;

                                    if(Math.Abs(yc) > (double)Int32.MaxValue) continue;

                                    ptCentre.Y = (int)Math.Round(yc);

                                   

                                    // Add solution

                                    nSolutions = AddSolution(ptCentre, radTanCirc, ref Centres, ref Radii);   

                              }

                              else  // If lines are parallel, there are 2 solutions

                              {

                                    double Ac = t1 * t1;

                                    double Bc = 2 * a1 * (c1 - (fRadius * t1));

                                    double Cc = ((fRadius * t1) - c1) * ((fRadius * t1) - c1) - (b1 * b2 * (r3 + fRadius) * ( r3 + fRadius));

                                    b24ac = (Bc * Bc) - (4 * Ac * Cc);

                                    if((b24ac < 0.0) || (Ac == 0.0)) continue;

                                    b24acRoot = Math.Sqrt(b24ac);

                                   

                                    for(int xCase=0; xCase<2; ++xCase)

                                    {

                                          // Calculate x-coordinate

                                          if(xCase == 0) xc = (-Bc + b24acRoot) / (Ac + Ac);

                                          else xc = (-Bc - b24acRoot) / (Ac + Ac);

                                          if(Math.Abs(xc) > (double)Int32.MaxValue) continue;

                                          ptCentre.X = (int)Math.Round(xc);

                                         

                                          // Calculate y-coordinate

                                          yc = (-a1 * xc - c1 + fRadius * t1) / b1;

                                          if(Math.Abs(yc) > (double)Int32.MaxValue) continue;

                                          ptCentre.Y = (int)Math.Round(yc);

                                         

                                          // Add solution

                                          nSolutions = AddSolution(ptCentre, radTanCirc, ref Centres, ref Radii);

                                    }          

                              }

                        }

                  }

 

                  // Return number of solutions

                  return nSolutions;

            }

 

            private static int AddSolution(Point ptCentre, int nRadius,

                  ref Point[] Centres, ref int[] Radii)

                  //*****************************************************

                  // Add a solution to a list if not already present

                  // In this overload, all values are integral

            {

                  Point[] CentreCopy = new Point[0];

                  int[] RadiiCopy = new int[0];

 

                  // Don't accept zero or negative radius solutions

                  int nSolutions = Centres.Length;

                  if(nRadius <= 0) return nSolutions;

 

                  // Check if present

                  bool bPresent = false;

                  for(int iSolution=0; iSolution<nSolutions; ++iSolution)

                  {

                        if((ptCentre == Centres[iSolution]) &&

                              (nRadius == Radii[iSolution]))

                        {

                              bPresent = true;

                              break;

                        }

                  }

 

                  // If found, keep current solutions

                  if(bPresent) return nSolutions;

 

                  // Otherwise, copy existing solutions to temporary arrays

                  if(nSolutions > 0)

                  {

                        CentreCopy = new Point[nSolutions];

                        RadiiCopy = new int[nSolutions];

                        Array.Copy(Centres, 0, CentreCopy, 0, nSolutions);

                        Array.Copy(Radii, 0, RadiiCopy, 0, nSolutions);

                  }

 

                  // Set up new solution arrays and load

                  Centres = new Point[nSolutions + 1];

                  Radii = new int[nSolutions + 1];

                  if(nSolutions > 0)

                  {

                        Array.Copy(CentreCopy, 0, Centres, 0, nSolutions);

                        Array.Copy(RadiiCopy, 0, Radii, 0, nSolutions);

                  }

                  Centres[nSolutions] = ptCentre;

                  Radii[nSolutions] = nRadius;

 

                  // Return updated number of solutions

                  return (nSolutions + 1);

            }

     

     

      }

}

using System;

using System.Drawing;

 

namespace Geometry

{

      /// <summary>

      /// Geometric solutions for straight lines

      /// </summary>

      public class LineGeometry

      {

            private LineGeometry()

            {

 

            }

 

            /// <summary>

            /// Get the normalised equation of a line between two points in the form: Ax + By + C = 0

            /// </summary>

            /// <param name="pt1">First point</param>

            /// <param name="pt2">Second point</param>

            /// <param name="A">X factor</param>

            /// <param name="B">Y factor</param>

            /// <param name="C">Constant</param>

            /// <returns>true on success</returns>

            public static bool Equation(Point pt1, Point pt2, ref double A, ref double B, ref double C)

            //****************************************************************************************

            // Derive the equation of the line between 2 points in the form Ax + By + C = 0

            //  If the points are coincident, return false

            // In this overload, coordinates are integral values

            {

                  // Validation

                  if((pt1 == pt2)) return false;

 

                  // Get constants

                  int xDiff = pt2.X - pt1.X;

                  int yDiff = pt2.Y - pt1.Y;

                  double rSquare = ((double)xDiff * xDiff) + ((double)yDiff * yDiff);

                  double rInv = 1.0 / rSquare;

 

                  // Derive parameters

                  A = -yDiff * rInv;

                  B = xDiff * rInv;

                  C = ((double)pt1.X * pt2.Y - (double)pt2.X * pt1.Y) * rInv;

 

                  // Normalize the equation for convenience

                  double sMult = 1.0 / Math.Sqrt( A * A + B * B);

                  A *= sMult;

                  B *= sMult;

                  C *= sMult;

 

                  // Return success

                  return true;

            }

 

           

      }

}