Coastal & Marine Geology InfoBank

Home FACS Activities Atlas Geology School Related Sites More

USGS InfoBank program -- dfndcint

Skip navigational links
Loading
InfoBank Programs: by Name   by Topic  
Expanded Description
Topic Description
Name
dfndcint
Synopsis
/infobank/programs/share/dfndcint.for
Description
  Subroutine to find the intersection point(s) of two circles given 

  center coordinates (e.g. UTM) and radii (e.g. range).
       
COMMENTS
  This subroutine has been used to determine ship coordinates

  knowing the range and position of two land-based ranging stations.
   
  All calculations are performed with double precision variables.
  Single precision input and output variables are converted within
  the subroutine.
   
    
USAGE
  "a" = argument, "r" = referenced, "s" = set

  
  call dfndcint                         (with the following arguments)
         X1Center                    [variable r*4 r]
         Y1Center                    [variable r*4 r]
         Radius1                     [variable r*4 r]
         X2Center                    [variable r*4 r]
         Y2Center                    [variable r*4 r]
         Radius2                     [variable r*4 variable r*4 s]
         X1Intersection              [variable r*4 s]
         Y1Intersection              [variable r*4 s]
         X2Intersection              [variable r*4 s]
         Y2Intersection              [variable char*4 sr]
    
  Input Variables:
      real * 4 X1Center        ! X-coordinate for center of circle 1
      real * 4 Y1Center        ! Y-coordinate for center of circle 1
      real * 4 Radius1         ! Radius of circle 1
      real * 4 X2Center        ! X-coordinate for center of circle 2
      real * 4 Y2Center        ! Y-coordinate for center of circle 2
      real * 4 Radius2         ! Radius of circle 2
   
  Output Variables:
      real * 4 X1Intersection  ! X-coordinate of intersection point 1
      real * 4 Y1Intersection  ! Y-coordinate of intersection point 1
      real * 4 X2Intersection  ! X-coordinate of intersection point 2
      real * 4 Y2Intersection  ! Y-coordinate of intersection point 2
      character * 4 StatusFlag ! One of the following:
                "TWO " - There are 2 intersection points
                         (i.e., overlapping circles).
                "ONE " - There is 1 intersection point
                         (i.e., tangential circles).
                         Points 1 and 2 will be set to the same value.
                "NONE" - There are no intersection points
                         (i.e., circles are too far apart,
                         one lies within the other, or
                         one has a zero or negative radius).
                         Points 1 and 2 will be set to 0.0
                "MANY" - There are many intersection points
                         (i.e., circles are identical).
                         Points 1 and 2 will be set to 0.0
   
    
FUNCTIONS
  dsqrt               (FORTRAN intrinsic function)  

    
    
DISCLAIMER
  Although this program has been used by the U.S. Geological Survey,

  no warranty, expressed or implied, is made by the Survey as to the
  accuracy and functioning of the program and related program
  material nor shall the fact of distribution constitute any such
  warranty, and no responsibility is assured by the Survey in
  connection therewith.
    
    
AUTHOR
  Carolyn Degnan            .for

  Clint Steele   02/18/88   Changed variable names for clarity
  Carolyn Degnan 08/15/88   Changed from MULTICS double precision
                            references to VMS real * 8.
                            Changed StatusFlag to indicate number
                            of intersection points.
                            Added check and alternate calculations
                            for when YCenter1 = YCenter2.
                            Added internal conversion of input and
                            output variables to double precision.
    
    
ALGORITHM
  If a ship's location at a given time is to be determined from the 

  ranges of 2 ranging stations of known positions, one must set 
  up a problem of finding the coordinates of the intersection points 
  of 2 overlapping circles.
     
  The center of each circle is at the ranging station's location.
  The radius of each circle is that stations's range.
      
  The general equation of a circle is:
     (x-XCenter)**2+(y-YCenter)**2 = RadiusCircle**2
  where (XCenter,YCenter) is the origin of the circle
  and RadiusCircle is the radius of the circle.
     
  Expanding this equation gives:
     x**2-(2*XCenter*x)+XCenter**2+
     y**2-(2*YCenter*y)+YCenter**2-RadiusCircle**2 = 0
     
  The equations for the 2 overlapping circles are:
      x**2-(2*XCenter1*x)+XCenter1**2+
      y**2-(2*YCenter1*y)+YCenter1**2-RadiusCircle1**2 = 0
      x**2-(2*XCenter2*x)+XCenter2**2+
      y**2-(2*YCenter2*y)+YCenter2**2-RadiusCircle2**2 = 0
  where...
     (XCenter1,YCenter1) are the coordinates of center of circle 1,
     RadiusCircle1 is the radius of circle 1,
     (XCenter2,YCenter2) are the coordinates of center of circle 2, 
     RadiusCircle2 is the radius of circle 2.
      
  Since both equal zero, set the 2 equations equal to each other.
     x**2-(2*XCenter1*x)+XCenter1**2+
     y**2-(2*YCenter1*y)+YCenter1**2-RadiusCircle1**2  = 
     x**2-(2*XCenter2*x)+XCenter2**2+
     y**2-(2*YCenter2*y)+YCenter2**2-RadiusCircle2**2
     
  Cancel out mutual terms.
     -(2*XCenter1*x)+XCenter1**2
     -(2*YCenter1*y)+YCenter1**2-RadiusCircle1**2  = 
     -(2*XCenter2*x)+XCenter2**2
     -(2*YCenter2*y)+YCenter2**2-RadiusCircle2**2
     
  Combine terms and rearrange to get the form: cx+dy+e = 0.
   (-2*XCenter1+2*XCenter2)*x+(-2*YCenter1+2*YCenter2)*y+
   (XCenter1**2+YCenter1**2-RadiusCircle1**2-XCenter2**2-YCenter2**2+
   RadiusCircle2**2)  =   0
  where:
    XDifference = 2*(XCenter2-XCenter1)
    YDifference = 2*(YCenter2-YCenter1)
    Value1 = (XCenter1**2+YCenter1**2-RadiusCircle1**2)-
        (XCenter2**2-YCenter2**2+RadiusCircle2**2)
      
  Rearrange this equation to get an expression for y in terms of x: 
      y = (-Value1-XDifference*x)/YDifference
      (NOTE:  If YCenter1 = YCenter2, you must proceed
      with an expression for x in terms of y instead.)
      
  Substitute this expression of y back into one of the circle 
  equations,
  x**2-2*XCenter1*x+XCenter1**2+
  ((-Value1-XDifference*x)/YDifference)**2-
  2*YCenter1*((-Value1-XDifference*x)/YDifference)+
  YCenter1**2-RadiusCircle1**2  =  0
     
  This gives:
  x**2-(2*XCenter1)x+XCenter1**2+
  (Value1**2+(2*Value1*XDifference)*x+
  (XDifference**2)*x**2)/YDifference**2-
  2*YCenter1*-Value1-XDifference*x))YDifference+YCenter1**2-
  RadiusCircle1**2 = 0
     
  Combine like terms:
  (1+XDifference**2/YDifference**2)*x**2+
  (-2*XCenter1+(2*Value1*XDifference)/YDifference**2+
  (2*YCenter1*XDifference)/YDifference)*x+
  (XCenter1**2+Value1**2/YDifference**2+
  (2+YCenter1*Value1)/YDifference+YCenter1**2-RadiusCircle1**2) = 0
     
  This gives the form, Value2*x**2+Value3*x+Value4 = 0
  where
  Value2 = 1+(XDifference**2)/(YDifference**2)
  Value3 = -2*XCenter1+(2*Value1*XDifference)/YDifference**2+
           (2*YCenter1*XDifference)/YDifference
  Value4 = XCenter1**2+Value1**2/YDifference**2+
           (2*YCenter1*Value1)/YDifference+YCenter1**2-
           RadiusCircle1**2
    
  Use the equation, x = (-Value3(+ or -)sqrt(Value3**2-
     (4*Value2*Value4)))/(2*Value2), 
  to find possible x values.
     
  To simplify, let Value5 = Value3**2-4*Value2*Value4 and 
                   Value6 = 2*Value2.

    
Reads
Writes
Opens
Calls

Skip footer navigational links

Coastal and Marine Science Centers:  Pacific   St. Petersburg   Woods Hole  
InfoBank   Coastal and Marine Geology Program   Geologic Information   Ask-A-Geologist   USGS Disclaimer  


Accessibility FOIA Privacy Policies and Notices

Take Pride in America logo USA.gov logo U.S. Department of the Interior | U.S. Geological Survey
URL: http://walrus.wr.usgs.gov/infobank/programs/share/dfndcint.doc.html
Page Contact Information: InfoBank staff
Page Last Modified: Mon Sep 16 03:36:50 PDT 2013  (chd)