|
/*
* (C) 2004 - Geotechnical Software Services
*
* This code is free software; you can redistribute it and/or
* modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation; either
* version 2.1 of the License, or (at your option) any later version.
*
* This code is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU Lesser General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public
* License along with this program; if not, write to the Free
* Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
* MA 02111-1307, USA.
*/
//package no.geosoft.cc.geometry;
/**
* Collection of geometry utility methods. All methods are static.
*
* @author <a href="mailto:jacob.dreyer@geosoft.no">Jacob Dreyer</a>
*/
public final class Geometry
{
/**
* Return true if c is between a and b.
*/
private static boolean isBetween (int a, int b, int c)
{
return b > a ? c >= a && c <= b : c >= b && c <= a;
}
/**
* Return true if c is between a and b.
*/
private static boolean isBetween (double a, double b, double c)
{
return b > a ? c >= a && c <= b : c >= b && c <= a;
}
/**
* Check if two double precision numbers are "equal", i.e. close enough
* to a given limit.
*
* @param a First number to check
* @param b Second number to check
* @param limit The definition of "equal".
* @return True if the twho numbers are "equal", false otherwise
*/
private static boolean equals (double a, double b, double limit)
{
return Math.abs (a - b) < limit;
}
/**
* Check if two double precision numbers are "equal", i.e. close enough
* to a prespecified limit.
*
* @param a First number to check
* @param b Second number to check
* @return True if the twho numbers are "equal", false otherwise
*/
private static boolean equals (double a, double b)
{
return equals (a, b, 1.0e-5);
}
/**
* Return smallest of four numbers.
*
* @param a First number to find smallest among.
* @param b Second number to find smallest among.
* @param c Third number to find smallest among.
* @param d Fourth number to find smallest among.
* @return Smallest of a, b, c and d.
*/
private static double min (double a, double b, double c, double d)
{
return Math.min (Math.min (a, b), Math.min (c, d));
}
/**
* Return largest of four numbers.
*
* @param a First number to find largest among.
* @param b Second number to find largest among.
* @param c Third number to find largest among.
* @param d Fourth number to find largest among.
* @return Largest of a, b, c and d.
*/
private static double max (double a, double b, double c, double d)
{
return Math.max (Math.max (a, b), Math.max (c, d));
}
/**
* Check if a specified point is inside a specified rectangle.
*
* @param x0, y0, x1, y1 Upper left and lower right corner of rectangle
* (inclusive)
* @param x,y Point to check.
* @return True if the point is inside the rectangle,
* false otherwise.
*/
public static boolean isPointInsideRectangle (int x0, int y0, int x1, int y1,
int x, int y)
{
return x >= x0 && x < x1 && y >= y0 && y < y1;
}
/**
* Check if a given point is inside a given (complex) polygon.
*
* @param x, y Polygon.
* @param pointX, pointY Point to check.
* @return True if the given point is inside the polygon, false otherwise.
*/
public static boolean isPointInsidePolygon (double[] x, double[] y,
double pointX, double pointY)
{
boolean isInside = false;
int nPoints = x.length;
int j = 0;
for (int i = 0; i < nPoints; i++) {
j++;
if (j == nPoints) j = 0;
if (y[i] < pointY && y[j] >= pointY || y[j] < pointY && y[i] >= pointY) {
if (x[i] + (pointY - y[i]) / (y[j] - y[i]) * (x[j] - x[i]) < pointX) {
isInside = !isInside;
}
}
}
return isInside;
}
/**
* Check if a given point is inside a given polygon. Integer domain.
*
* @param x, y Polygon.
* @param pointX, pointY Point to check.
* @return True if the given point is inside the polygon, false otherwise.
*/
public static boolean isPointInsidePolygon (int[] x, int[] y,
int pointX, int pointY)
{
boolean isInside = false;
int nPoints = x.length;
int j = 0;
for (int i = 0; i < nPoints; i++) {
j++;
if (j == nPoints) j = 0;
if (y[i] < pointY && y[j] >= pointY || y[j] < pointY && y[i] >= pointY) {
if (x[i] + (double) (pointY - y[i]) / (double) (y[j] - y[i]) *
(x[j] - x[i]) < pointX) {
isInside = !isInside;
}
}
}
return isInside;
}
/**
* Find the point on the line p0,p1 [x,y,z] a given fraction from p0.
* Fraction of 0.0 whould give back p0, 1.0 give back p1, 0.5 returns
* midpoint of line p0,p1 and so on. F
* raction can be >1 and it can be negative to return any point on the
* line specified by p0,p1.
*
* @param p0 First coordinale of line [x,y,z].
* @param p0 Second coordinale of line [x,y,z].
* @param fractionFromP0 Point we are looking for coordinates of
* @param p Coordinate of point we are looking for
*/
public static double[] computePointOnLine (double[] p0, double[] p1,
double fractionFromP0)
{
double[] p = new double[3];
p[0] = p0[0] + fractionFromP0 * (p1[0] - p0[0]);
p[1] = p0[1] + fractionFromP0 * (p1[1] - p0[1]);
p[2] = p0[2] + fractionFromP0 * (p1[2] - p0[2]);
return p;
}
/**
* Find the point on the line defined by x0,y0,x1,y1 a given fraction
* from x0,y0. 2D version of method above..
*
* @param x0, y0 First point defining the line
* @param x1, y1 Second point defining the line
* @param fractionFrom0 Distance from (x0,y0)
* @return x, y Coordinate of point we are looking for
*/
public static double[] computePointOnLine (double x0, double y0,
double x1, double y1,
double fractionFrom0)
{
double[] p0 = {x0, y0, 0.0};
double[] p1 = {x1, y1, 0.0};
double[] p = Geometry.computePointOnLine (p0, p1, fractionFrom0);
double[] r = {p[0], p[1]};
return r;
}
/**
* Extend a given line segment to a specified length.
*
* @param p0, p1 Line segment to extend [x,y,z].
* @param toLength Length of new line segment.
* @param anchor Specifies the fixed point during extension.
* If anchor is 0.0, p0 is fixed and p1 is adjusted.
* If anchor is 1.0, p1 is fixed and p0 is adjusted.
* If anchor is 0.5, the line is adjusted equally in each
* direction and so on.
*/
public static void extendLine (double[] p0, double[] p1, double toLength,
double anchor)
{
double[] p = Geometry.computePointOnLine (p0, p1, anchor);
double length0 = toLength * anchor;
double length1 = toLength * (1.0 - anchor);
Geometry.extendLine (p, p0, length0);
Geometry.extendLine (p, p1, length1);
}
/**
* Extend a given line segment to a given length and holding the first
* point of the line as fixed.
*
* @param p0, p1 Line segment to extend. p0 is fixed during extension
* @param length Length of new line segment.
*/
public static void extendLine (double[] p0, double[] p1, double toLength)
{
double oldLength = Geometry.length (p0, p1);
double lengthFraction = oldLength != 0.0 ? toLength / oldLength : 0.0;
p1[0] = p0[0] + (p1[0] - p0[0]) * lengthFraction;
p1[1] = p0[1] + (p1[1] - p0[1]) * lengthFraction;
p1[2] = p0[2] + (p1[2] - p0[2]) * lengthFraction;
}
/**
* Return the length of a vector.
*
* @param v Vector to compute length of [x,y,z].
* @return Length of vector.
*/
public static double length (double[] v)
{
return Math.sqrt (v[0]*v[0] + v[1]*v[1] + v[2]*v[2]);
}
/**
* Compute distance between two points.
*
* @param p0, p1 Points to compute distance between [x,y,z].
* @return Distance between points.
*/
public static double length (double[] p0, double[] p1)
{
double[] v = Geometry.createVector (p0, p1);
return length (v);
}
/**
* Compute the length of the line from (x0,y0) to (x1,y1)
*
* @param x0, y0 First line end point.
* @param x1, y1 Second line end point.
* @return Length of line from (x0,y0) to (x1,y1).
*/
public static double length (int x0, int y0, int x1, int y1)
{
return Geometry.length ((double) x0, (double) y0,
(double) x1, (double) y1);
}
/**
* Compute the length of the line from (x0,y0) to (x1,y1)
*
* @param x0, y0 First line end point.
* @param x1, y1 Second line end point.
* @return Length of line from (x0,y0) to (x1,y1).
*/
public static double length (double x0, double y0, double x1, double y1)
{
double dx = x1 - x0;
double dy = y1 - y0;
return Math.sqrt (dx*dx + dy*dy);
}
/**
* Compute the length of a polyline.
*
* @param x, y Arrays of x,y coordinates
* @param nPoints Number of elements in the above.
* @param isClosed True if this is a closed polygon, false otherwise
* @return Length of polyline defined by x, y and nPoints.
*/
public static double length (int[] x, int[] y, boolean isClosed)
{
double length = 0.0;
int nPoints = x.length;
for (int i = 0; i < nPoints-1; i++)
length += Geometry.length (x[i], y[i], x[i+1], y[i+1]);
// Add last leg if this is a polygon
if (isClosed && nPoints > 1)
length += Geometry.length (x[nPoints-1], y[nPoints-1], x[0], y[0]);
return length;
}
/**
* Return distance bwetween the line defined by (x0,y0) and (x1,y1)
* and the point (x,y).
* Ref: http://astronomy.swin.edu.au/pbourke/geometry/pointline/
* The 3D case should be similar.
*
* @param x0, y0 First point of line.
* @param x1, y1 Second point of line.
* @param x, y, Point to consider.
* @return Distance from x,y down to the (extended) line defined
* by x0, y0, x1, y1.
*/
public static double distance (int x0, int y0, int x1, int y1,
int x, int y)
{
// If x0,y0,x1,y1 is same point, we return distance to that point
double length = Geometry.length (x0, y0, x1, y1);
if (length == 0.0)
return Geometry.length (x0, y0, x, y);
// If u is [0,1] then (xp,yp) is on the line segment (x0,y0),(x1,y1).
double u = ((x - x0) * (x1 - x0) + (y - y0) * (y1 - y0)) /
(length * length);
// This is the intersection point of the normal.
// TODO: Might consider returning this as well.
double xp = x0 + u * (x1 - x0);
double yp = y0 + u * (y1 - y0);
length = Geometry.length (xp, yp, x, y);
return length;
}
/**
* Find the angle between twree points. P0 is center point
*
* @param p0, p1, p2 Three points finding angle between [x,y,z].
* @return Angle (in radians) between given points.
*/
public static double computeAngle (double[] p0, double[] p1, double[] p2)
{
double[] v0 = Geometry.createVector (p0, p1);
double[] v1 = Geometry.createVector (p0, p2);
double dotProduct = Geometry.computeDotProduct (v0, v1);
double length1 = Geometry.length (v0);
double length2 = Geometry.length (v1);
double denominator = length1 * length2;
double product = denominator != 0.0 ? dotProduct / denominator : 0.0;
double angle = Math.acos (product);
return angle;
}
/**
* Compute the dot product (a scalar) between two vectors.
*
* @param v0, v1 Vectors to compute dot product between [x,y,z].
* @return Dot product of given vectors.
*/
public static double computeDotProduct (double[] v0, double[] v1)
{
return v0[0] * v1[0] + v0[1] * v1[1] + v0[2] * v1[2];
}
/**
* Compute the cross product (a vector) of two vectors.
*
* @param v0, v1 Vectors to compute cross product between [x,y,z].
* @param crossProduct Cross product of specified vectors [x,y,z].
*/
public static double[] computeCrossProduct (double[] v0, double[] v1)
{
double crossProduct[] = new double[3];
crossProduct[0] = v0[1] * v1[2] - v0[2] * v1[1];
crossProduct[1] = v0[2] * v1[0] - v0[0] * v1[2];
crossProduct[2] = v0[0] * v1[1] - v0[1] * v1[0];
return crossProduct;
}
/**
* Construct the vector specified by two points.
*
* @param p0, p1 Points the construct vector between [x,y,z].
* @return v Vector from p0 to p1 [x,y,z].
*/
public static double[] createVector (double[] p0, double[] p1)
{
double v[] = {p1[0] - p0[0], p1[1] - p0[1], p1[2] - p0[2]};
return v;
}
/**
* Check if two points are on the same side of a given line.
* Algorithm from Sedgewick page 350.
*
* @param x0, y0, x1, y1 The line.
* @param px0, py0 First point.
* @param px1, py1 Second point.
* @return <0 if points on opposite sides.
* =0 if one of the points is exactly on the line
* >0 if points on same side.
*/
private static int sameSide (double x0, double y0, double x1, double y1,
double px0, double py0, double px1, double py1)
{
int sameSide = 0;
double dx = x1 - x0;
double dy = y1 - y0;
double dx1 = px0 - x0;
double dy1 = py0 - y0;
double dx2 = px1 - x1;
double dy2 = py1 - y1;
// Cross product of the vector from the endpoint of the line to the point
double c1 = dx * dy1 - dy * dx1;
double c2 = dx * dy2 - dy * dx2;
if (c1 != 0 && c2 != 0)
sameSide = c1 < 0 != c2 < 0 ? -1 : 1;
else if (dx == 0 && dx1 == 0 && dx2 == 0)
sameSide = !isBetween (y0, y1, py0) && !isBetween (y0, y1, py1) ? 1 : 0;
else if (dy == 0 && dy1 == 0 && dy2 == 0)
sameSide = !isBetween (x0, x1, px0) && !isBetween (x0, x1, px1) ? 1 : 0;
return sameSide;
}
/**
* Check if two points are on the same side of a given line. Integer domain.
*
* @param x0, y0, x1, y1 The line.
* @param px0, py0 First point.
* @param px1, py1 Second point.
* @return <0 if points on opposite sides.
* =0 if one of the points is exactly on the line
* >0 if points on same side.
*/
private static int sameSide (int x0, int y0, int x1, int y1,
int px0, int py0, int px1, int py1)
{
return sameSide ((double) x0, (double) y0, (double) x1, (double) y1,
(double) px0, (double) py0, (double) px1, (double) py1);
}
/**
* Check if two line segments intersects. Integer domain.
*
* @param x0, y0, x1, y1 End points of first line to check.
* @param x2, yy, x3, y3 End points of second line to check.
* @return True if the two lines intersects.
*/
public static boolean isLineIntersectingLine (int x0, int y0, int x1, int y1,
int x2, int y2, int x3, int y3)
{
int s1 = Geometry.sameSide (x0, y0, x1, y1, x2, y2, x3, y3);
int s2 = Geometry.sameSide (x2, y2, x3, y3, x0, y0, x1, y1);
return s1 <= 0 && s2 <= 0;
}
/**
* Check if a specified line intersects a specified rectangle.
* Integer domain.
*
* @param lx0, ly0 1st end point of line
* @param ly1, ly1 2nd end point of line
* @param x0, y0, x1, y1 Upper left and lower right corner of rectangle
* (inclusive).
* @return True if the line intersects the rectangle,
* false otherwise.
*/
public static boolean isLineIntersectingRectangle (int lx0, int ly0,
int lx1, int ly1,
int x0, int y0,
int x1, int y1)
{
// Is one of the line endpoints inside the rectangle
if (Geometry.isPointInsideRectangle (x0, y0, x1, y1, lx0, ly0) ||
Geometry.isPointInsideRectangle (x0, y0, x1, y1, lx1, ly1))
return true;
// If it intersects it goes through. Need to check three sides only.
// Check against top rectangle line
if
|