Class Distance

  1 /*
  2  * Copyright (c) 2016 Vivid Solutions.
  3  *
  4  * All rights reserved. This program and the accompanying materials
  5  * are made available under the terms of the Eclipse Public License 2.0
  6  * and Eclipse Distribution License v. 1.0 which accompanies this distribution.
  7  * The Eclipse Public License is available at http://www.eclipse.org/legal/epl-v20.html
  8  * and the Eclipse Distribution License is available at
  9  *
 10  * http://www.eclipse.org/org/documents/edl-v10.php.
 11  */
 12 package org.locationtech.jts.algorithm;
 13  
 14 import org.locationtech.jts.geom.Coordinate;
 15 import org.locationtech.jts.geom.Envelope;
 16 import org.locationtech.jts.math.MathUtil;
 17  
 18 /**
 19  * Functions to compute distance between basic geometric structures.
 20  * 
 21  * @author Martin Davis
 22  *
 23  */
 24 public class Distance {
 25  
 26   /**
 27    * Computes the distance from a line segment AB to a line segment CD
 28    * 
 29    * Note: NON-ROBUST!
 30    * 
 31    * @param A
 32    *          a point of one line
 33    * @param B
 34    *          the second point of (must be different to A)
 35    * @param C
 36    *          one point of the line
 37    * @param D
 38    *          another point of the line (must be different to A)
 39    */
 40   public static double segmentToSegment(Coordinate A, Coordinate B,
 41       Coordinate C, Coordinate D)
 42   {
 43     // check for zero-length segments
 44     if (A.equals(B))
 45       return Distance.pointToSegment(A, C, D);
 46     if (C.equals(D))
 47       return Distance.pointToSegment(D, A, B);
 48   
 49     // AB and CD are line segments
 50     /*
 51      * from comp.graphics.algo
 52      * 
 53      * Solving the above for r and s yields 
 54      * 
 55      *     (Ay-Cy)(Dx-Cx)-(Ax-Cx)(Dy-Cy) 
 56      * r = ----------------------------- (eqn 1) 
 57      *     (Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)
 58      * 
 59      *     (Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay)  
 60      * s = ----------------------------- (eqn 2)
 61      *     (Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx) 
 62      *     
 63      * Let P be the position vector of the
 64      * intersection point, then 
 65      *   P=A+r(B-A) or 
 66      *   Px=Ax+r(Bx-Ax) 
 67      *   Py=Ay+r(By-Ay) 
 68      * By examining the values of r & s, you can also determine some other limiting
 69      * conditions: 
 70      *   If 0<=r<=1 & 0<=s<=1, intersection exists 
 71      *      r<0 or r>1 or s<0 or s>1 line segments do not intersect 
 72      *   If the denominator in eqn 1 is zero, AB & CD are parallel 
 73      *   If the numerator in eqn 1 is also zero, AB & CD are collinear.
 74      */
 75   
 76     boolean noIntersection = false;
 77     if (! Envelope.intersects(A, B, C, D)) {
 78       noIntersection = true;
 79     }
 80     else {
 81       double denom = (B.x - A.x) * (D.y - C.y) - (B.y - A.y) * (D.x - C.x);
 82       
 83       if (denom == 0) {
 84         noIntersection = true;
 85       }
 86       else {
 87         double r_num = (A.y - C.y) * (D.x - C.x) - (A.x - C.x) * (D.y - C.y);
 88         double s_num = (A.y - C.y) * (B.x - A.x) - (A.x - C.x) * (B.y - A.y);
 89         
 90         double s = s_num / denom;
 91         double r = r_num / denom;
 92   
 93         if ((r < 0) || (r > 1) || (s < 0) || (s > 1)) {
 94           noIntersection = true;
 95         }
 96       }
 97     }
 98     if (noIntersection) {
 99       return MathUtil.min(
100             Distance.pointToSegment(A, C, D),
101             Distance.pointToSegment(B, C, D),
102             Distance.pointToSegment(C, A, B),
103             Distance.pointToSegment(D, A, B));
104     }
105     // segments intersect
106     return 0.0
107   }
108  
109   /**
110    * Computes the distance from a point to a sequence of line segments.
111    * 
112    * @param p
113    *          a point
114    * @param line
115    *          a sequence of contiguous line segments defined by their vertices
116    * @return the minimum distance between the point and the line segments
117    */
118   public static double pointToSegmentString(Coordinate p, Coordinate[] line)
119   {
120     if (line.length == 0)
121       throw new IllegalArgumentException(
122           "Line array must contain at least one vertex");
123     // this handles the case of length = 1
124     double minDistance = p.distance(line[0]);
125     for (int i = 0; i < line.length - 1; i++) {
126       double dist = Distance.pointToSegment(p, line[i], line[i + 1]);
127       if (dist < minDistance) {
128         minDistance = dist;
129       }
130     }
131     return minDistance;
132   }
133  
134   /**
135    * Computes the distance from a point p to a line segment AB
136    * 
137    * Note: NON-ROBUST!
138    * 
139    * @param p
140    *          the point to compute the distance for
141    * @param A
142    *          one point of the line
143    * @param B
144    *          another point of the line (must be different to A)
145    * @return the distance from p to line segment AB
146    */
147   public static double pointToSegment(Coordinate p, Coordinate A,
148       Coordinate B)
149   {
150     // if start = end, then just compute distance to one of the endpoints
151     if (A.x == B.x && A.y == B.y)
152       return p.distance(A);
153   
154     // otherwise use comp.graphics.algorithms Frequently Asked Questions method
155     /*
156      * (1) r = AC dot AB 
157      *         --------- 
158      *         ||AB||^2 
159      *         
160      * r has the following meaning: 
161      *   r=0 P = A 
162      *   r=1 P = B 
163      *   r<0 P is on the backward extension of AB 
164      *   r>1 P is on the forward extension of AB 
165      *   0<r<1 P is interior to AB
166      */
167   
168     double len2 = (B.x - A.x) * (B.x - A.x) + (B.y - A.y) * (B.y - A.y);
169     double r = ((p.x - A.x) * (B.x - A.x) + (p.y - A.y) * (B.y - A.y))
170         / len2;
171   
172     if (r <= 0.0)
173       return p.distance(A);
174     if (r >= 1.0)
175       return p.distance(B);
176   
177     /*
178      * (2) s = (Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay) 
179      *         ----------------------------- 
180      *                    L^2
181      * 
182      * Then the distance from C to P = |s|*L.
183      * 
184      * This is the same calculation as {@link #distancePointLinePerpendicular}.
185      * Unrolled here for performance.
186      */
187     double s = ((A.y - p.y) * (B.x - A.x) - (A.x - p.x) * (B.y - A.y))
188         / len2;
189     return Math.abs(s) * Math.sqrt(len2);
190   }
191  
192   /**
193    * Computes the perpendicular distance from a point p to the (infinite) line
194    * containing the points AB
195    * 
196    * @param p
197    *          the point to compute the distance for
198    * @param A
199    *          one point of the line
200    * @param B
201    *          another point of the line (must be different to A)
202    * @return the distance from p to line AB
203    */
204   public static double pointToLinePerpendicular(Coordinate p,
205       Coordinate A, Coordinate B)
206   {
207     // use comp.graphics.algorithms Frequently Asked Questions method
208     /*
209      * (2) s = (Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay) 
210      *         ----------------------------- 
211      *                    L^2
212      * 
213      * Then the distance from C to P = |s|*L.
214      */
215     double len2 = (B.x - A.x) * (B.x - A.x) + (B.y - A.y) * (B.y - A.y);
216     double s = ((A.y - p.y) * (B.x - A.x) - (A.x - p.x) * (B.y - A.y))
217         / len2;
218   
219     return Math.abs(s) * Math.sqrt(len2);
220   }
221  
222 }
223