001: /*
002: * GeoTools - OpenSource mapping toolkit
003: * http://geotools.org
004: * (C) 2003-2006, GeoTools Project Managment Committee (PMC)
005: * (C) 2003, Institut de Recherche pour le Développement
006: * (C) 1998, Pêches et Océans Canada
007: *
008: * This library is free software; you can redistribute it and/or
009: * modify it under the terms of the GNU Lesser General Public
010: * License as published by the Free Software Foundation;
011: * version 2.1 of the License.
012: *
013: * This library is distributed in the hope that it will be useful,
014: * but WITHOUT ANY WARRANTY; without even the implied warranty of
015: * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
016: * Lesser General Public License for more details.
017: */
018: package org.geotools.math;
019:
020: // J2SE dependencies
021: import java.awt.geom.Line2D;
022: import java.awt.geom.Point2D;
023: import java.io.Serializable;
024:
025: import javax.vecmath.MismatchedSizeException;
026:
027: import org.opengis.util.Cloneable;
028:
029: /**
030: * Equation of a line in a two dimensional space (<var>x</var>,<var>y</var>).
031: * A line has an equation of the form <var>y</var>=<var>a</var><var>x</var>+<var>b</var>.
032: * At the difference of {@link Line2D} (which are bounded by (<var>x1</var>,<var>y1</var>)
033: * and (<var>x2</var>,<var>y2</var>) points), {@code Line} objects extends toward infinity.
034: *
035: * The equation parameters for a {@code Line} object can bet set at construction
036: * time or using one of the {@code setLine(...)} methods. The <var>y</var> value
037: * can be computed for a given <var>x</var> value using the {@link #y} method. Method
038: * {@link #x} compute the converse and should work even if the line is vertical.
039: *
040: * @source $URL: http://svn.geotools.org/geotools/tags/2.4.1/modules/library/metadata/src/main/java/org/geotools/math/Line.java $
041: * @version $Id: Line.java 22443 2006-10-27 20:47:22Z desruisseaux $
042: * @author Martin Desruisseaux
043: *
044: * @since 2.0
045: *
046: * @see Point2D
047: * @see Line2D
048: * @see Plane
049: */
050: public class Line implements Cloneable, Serializable {
051: /**
052: * Serial number for compatibility with different versions.
053: */
054: private static final long serialVersionUID = 2185952238314399110L;
055:
056: /**
057: * Small value for rounding errors.
058: */
059: private static final double EPS = 1E-12;
060:
061: /**
062: * The slope for this line.
063: */
064: private double slope;
065:
066: /**
067: * The <var>y</var> value at <var>x</var>==0.
068: */
069: private double y0;
070:
071: /**
072: * Value of <var>x</var> at <var>y</var>==0.
073: * This value is used for vertical lines.
074: */
075: private double x0;
076:
077: /**
078: * Construct an initially unitialized line. All methods will returns {@link Double#NaN}.
079: */
080: public Line() {
081: slope = y0 = x0 = Double.NaN;
082: }
083:
084: /**
085: * Construct a line with the specified slope and offset.
086: * The linear equation will be <var>y</var>=<var>slope</var>*<var>x</var>+<var>y0</var>.
087: *
088: * @param slope The slope.
089: * @param y0 The <var>y</var> value at <var>x</var>==0.
090: *
091: * @see #setLine(double, double)
092: */
093: public Line(final double slope, final double y0) {
094: this .slope = slope;
095: this .y0 = y0;
096: this .x0 = -y0 / slope;
097: }
098:
099: /**
100: * Set the slope and offset for this line.
101: * The linear equation will be <var>y</var>=<var>slope</var>*<var>x</var>+<var>y0</var>.
102: *
103: * @param slope The slope.
104: * @param y0 The <var>y</var> value at <var>x</var>==0.
105: *
106: * @see #setLine(Point2D, Point2D)
107: * @see #setLine(Line2D)
108: * @see #setLine(double[], double[])
109: */
110: public void setLine(final double slope, final double y0) {
111: this .slope = slope;
112: this .y0 = y0;
113: this .x0 = -y0 / slope;
114: }
115:
116: /**
117: * Set a line colinear with the specified line segment. The line will continues
118: * toward infinity after the line segment extremities.
119: *
120: * @param line The line segment.
121: *
122: * @see #setLine(Point2D,Point2D)
123: */
124: public void setLine(final Line2D line) {
125: setLine(line.getX1(), line.getY1(), line.getX2(), line.getY2());
126: }
127:
128: /**
129: * Set a line through the specified point. The line will continues
130: * toward infinity after the point.
131: *
132: * @param p1 Coordinate of the first point.
133: * @param p2 Coordinate of the second point.
134: *
135: * @see #setLine(Line2D)
136: */
137: public void setLine(final Point2D p1, final Point2D p2) {
138: setLine(p1.getX(), p1.getY(), p2.getX(), p2.getY());
139: }
140:
141: /**
142: * Set a line through the specified point. The line will continues
143: * toward infinity after the point.
144: *
145: * @param x1 Ordinate <var>x</var> of the first point.
146: * @param y1 Ordinate <var>y</var> of the first point.
147: * @param x2 Ordinate <var>x</var> of the second point.
148: * @param y2 Ordinate <var>y</var> of the second point.
149: *
150: * @see #setLine(Point2D,Point2D)
151: * @see #setLine(Line2D)
152: */
153: private void setLine(final double x1, final double y1,
154: final double x2, final double y2) {
155: this .slope = (y2 - y1) / (x2 - x1);
156: this .x0 = x2 - y2 / slope;
157: this .y0 = y2 - slope * x2;
158: if (Double.isNaN(x0) && slope == 0) {
159: // Occurs for horizontal lines right on the x axis.
160: x0 = Double.POSITIVE_INFINITY;
161: }
162: if (Double.isNaN(y0) && Double.isInfinite(slope)) {
163: // Occurs for vertical lines right on the y axis.
164: y0 = Double.POSITIVE_INFINITY;
165: }
166: }
167:
168: /**
169: * Given a set of data points {@code x[0..ndata-1]}, {@code y[0..ndata-1]},
170: * fit them to a straight line <var>y</var>=<var>b</var>+<var>m</var><var>x</var> in
171: * a least-squares senses. This method assume that the <var>x</var> values are precise
172: * and all uncertainty is in <var>y</var>.
173: *
174: * <p>Reference: <a
175: * href="http://shakti.cc.trincoll.edu/~palladin/courses/ENGR431/statistics/node9.html">Linear
176: * Regression Curve Fitting</a>.
177: *
178: * @param x Vector of <var>x</var> values (independant variable).
179: * @param y Vector of <var>y</var> values (dependant variable).
180: * @return Estimation of the correlation coefficient. The closer
181: * this coefficient is to 1, the better the fit.
182: *
183: * @throws MismatchedSizeException if <var>x</var> and <var>y</var> don't have the same length.
184: */
185: public double setLine(final double[] x, final double[] y)
186: throws MismatchedSizeException {
187: final int N = x.length;
188: if (N != y.length) {
189: throw new MismatchedSizeException();
190: }
191: int count = 0;
192: double mean_x = 0;
193: double mean_y = 0;
194: for (int i = 0; i < N; i++) {
195: final double xi = x[i];
196: final double yi = y[i];
197: if (!Double.isNaN(xi) && !Double.isNaN(yi)) {
198: mean_x += xi;
199: mean_y += yi;
200: count++;
201: }
202: }
203: mean_x /= count;
204: mean_y /= count;
205: /*
206: * We have to solve two equations with two unknows:
207: *
208: * 1) mean(y) = b + m*mean(x)
209: * 2) mean(xy) = b*mean(x) + m*mean(x²)
210: *
211: * Those formulas lead to a quadratic equation. However,
212: * the formulas become very simples if we set 'mean(x)=0'.
213: * We can achieve this result by computing instead of (2):
214: *
215: * 2b) mean(dx y) = m*mean(dx²)
216: *
217: * where dx=x-mean(x). In this case mean(dx)==0.
218: */
219: double mean_x2 = 0;
220: double mean_y2 = 0;
221: double mean_xy = 0;
222: for (int i = 0; i < N; i++) {
223: double xi = x[i];
224: double yi = y[i];
225: if (!Double.isNaN(xi) && !Double.isNaN(yi)) {
226: xi -= mean_x;
227: mean_x2 += xi * xi;
228: mean_y2 += yi * yi;
229: mean_xy += xi * yi;
230: }
231: }
232: mean_x2 /= count;
233: mean_y2 /= count;
234: mean_xy /= count;
235: /*
236: * Assuming that 'mean(x)==0', then the correlation
237: * coefficient can be approximate by:
238: *
239: * R = mean(xy) / sqrt( mean(x²) * (mean(y²) - mean(y)²) )
240: */
241: slope = mean_xy / mean_x2;
242: y0 = mean_y - mean_x * slope;
243: return mean_xy
244: / Math.sqrt(mean_x2 * (mean_y2 - mean_y * mean_y));
245: }
246:
247: /**
248: * Translate the line. The slope stay unchanged.
249: *
250: * @param dx The horizontal translation.
251: * @param dy The vertical translation.
252: */
253: public void translate(final double dx, final double dy) {
254: if (slope == 0 || Double.isInfinite(slope)) {
255: x0 += dx;
256: y0 += dy;
257: } else {
258: x0 += dx - dy / slope;
259: y0 += dy - slope * dx;
260: }
261: }
262:
263: /**
264: * Compute <var>y</var>=<var>f</var>(<var>x</var>).
265: * If the line is vertical, then this method returns an infinite value.
266: * This method is final for performance reason.
267: *
268: * @param x The <var>x</var> value.
269: * @return The <var>y</var> value.
270: *
271: * @see #x(double)
272: */
273: public final double y(final double x) {
274: return slope * x + y0;
275: }
276:
277: /**
278: * Compute <var>x</var>=<var>f</var><sup>-1</sup>(<var>y</var>).
279: * If the line is horizontal, then this method returns an infinite value.
280: * This method is final for performance reason.
281: *
282: * @param y The <var>y</var> value.
283: * @return The <var>x</var> value.
284: *
285: * @see #y(double)
286: */
287: public final double x(final double y) {
288: return y / slope + x0;
289: }
290:
291: /**
292: * Returns the <var>y</var> value for <var>x</var>==0.
293: * Coordinate (0, <var>y0</var>) is the intersection point with the <var>y</var> axis.
294: */
295: public final double getY0() {
296: return y0;
297: }
298:
299: /**
300: * Returns the <var>x</var> value for <var>y</var>==0.
301: * Coordinate (<var>x0</var>,0) is the intersection point with the <var>x</var> axis.
302: */
303: public final double getX0() {
304: return x0;
305: }
306:
307: /**
308: * Returns the slope.
309: */
310: public final double getSlope() {
311: return slope;
312: }
313:
314: /**
315: * Returns the intersection point between this line and the specified one.
316: * If both lines are parallel, then this method returns {@code null}.
317: *
318: * @param line The line to intersect.
319: * @return The intersection point, or {@code null}.
320: */
321: public Point2D intersectionPoint(final Line line) {
322: double x, y;
323: if (Double.isInfinite(slope)) {
324: if (Double.isInfinite(line.slope)) {
325: return null;
326: }
327: x = x0;
328: y = x * line.slope + line.y0;
329: } else {
330: if (!Double.isInfinite(line.slope)) {
331: x = (y0 - line.y0) / (line.slope - slope);
332: if (Double.isInfinite(x)) {
333: return null;
334: }
335: } else {
336: x = line.x0;
337: }
338: y = x * slope + y0;
339: }
340: return new Point2D.Double(x, y);
341: }
342:
343: /**
344: * Returns the intersection point between this line and the specified bounded line.
345: * If both lines are parallel or if the specified {@code line} doesn't reach
346: * this line (since {@link Line2D} do not extends toward infinities), then this
347: * method returns {@code null}.
348: *
349: * @param line The bounded line to intersect.
350: * @return The intersection point, or {@code null}.
351: */
352: public Point2D intersectionPoint(final Line2D line) {
353: final double x1 = line.getX1();
354: final double y1 = line.getY1();
355: final double x2 = line.getX2();
356: final double y2 = line.getY2();
357: double x, y;
358: double m = (y2 - y1) / (x2 - x1);
359: if (Double.isInfinite(slope)) {
360: x = x0;
361: y = x * m + (y2 - m * x2);
362: } else {
363: if (!Double.isInfinite(m)) {
364: x = (y0 - (y2 - m * x2)) / (m - slope);
365: } else {
366: x = 0.5 * (x1 + x2);
367: }
368: y = x * slope + y0;
369: }
370: double eps;
371: /*
372: * Ensure that the intersection is in the range of valid x values.
373: */
374: eps = EPS * Math.abs(x);
375: if (x1 <= x2) {
376: if (!(x >= x1 - eps && x <= x2 + eps)) {
377: return null;
378: }
379: } else {
380: if (!(x <= x1 + eps && x >= x2 - eps)) {
381: return null;
382: }
383: }
384: /*
385: * Ensure that the intersection is in the range of valid y values.
386: */
387: eps = EPS * Math.abs(y);
388: if (y1 <= y2) {
389: if (!(y >= y1 - eps && y <= y2 + eps)) {
390: return null;
391: }
392: } else {
393: if (!(y <= y1 - eps && y >= y2 + eps)) {
394: return null;
395: }
396: }
397: return new Point2D.Double(x, y);
398: }
399:
400: /**
401: * Returns the nearest point on this line from the specified point.
402: *
403: * @param point An arbitrary point.
404: * @return The point on this line which is the nearest of the specified {@code point}.
405: */
406: public Point2D nearestColinearPoint(final Point2D point) {
407: if (!Double.isInfinite(slope)) {
408: final double x = ((point.getY() - y0) * slope + point
409: .getX())
410: / (slope * slope + 1);
411: return new Point2D.Double(x, x * slope + y0);
412: } else {
413: return new Point2D.Double(x0, point.getY());
414: }
415: }
416:
417: /**
418: * Compute the base of a isosceles triangle having the specified summit and side length.
419: * The base will be colinear with this line. In other words, this method compute two
420: * points (<var>x1</var>,<var>y1</var>) and (<var>x2</var>,<var>y2</var>) located in
421: * such a way that:
422: * <ul>
423: * <li>Both points are on this line.</li>
424: * <li>The distance between any of the two points and the specified {@code summit}
425: * is exactly {@code sideLength}.</li>
426: * </ul>
427: *
428: * @param summit The summit of the isosceles triangle.
429: * @param sideLength The length for the two sides of the isosceles triangle.
430: * @return The base of the isoscele triangle, colinear with this line, or {@code null}
431: * if the base can't be computed. If non-null, then the triangle is the figure formed
432: * by joining (<var>x1</var>,<var>y1</var>), (<var>x2</var>,<var>y2</var>) and
433: * {@code summit}.
434: */
435: public Line2D isoscelesTriangleBase(final Point2D summit,
436: double sideLength) {
437: sideLength *= sideLength;
438: if (slope == 0) {
439: final double x = summit.getX();
440: final double dy = y0 - summit.getY();
441: final double dx = Math.sqrt(sideLength - dy * dy);
442: if (Double.isNaN(dx)) {
443: return null;
444: }
445: return new Line2D.Double(x + dx, y0, x - dx, y0);
446: }
447: if (Double.isInfinite(slope)) {
448: final double y = summit.getY();
449: final double dx = x0 - summit.getX();
450: final double dy = Math.sqrt(sideLength - dx * dx);
451: if (Double.isNaN(dy)) {
452: return null;
453: }
454: return new Line2D.Double(x0, y + dy, x0, y - dy);
455: }
456: final double x = summit.getX();
457: final double y = summit.getY();
458: final double dy = y0 - y + slope * x;
459: final double B = -slope * dy;
460: final double A = slope * slope + 1;
461: final double C = Math.sqrt(B * B + A * (sideLength - dy * dy));
462: if (Double.isNaN(C)) {
463: return null;
464: }
465: final double x1 = (B + C) / A + x;
466: final double x2 = (B - C) / A + x;
467: return new Line2D.Double(x1, slope * x1 + y0, x2, slope * x2
468: + y0);
469: }
470:
471: /**
472: * Returns a string representation of this line. This method returns
473: * the linear equation in the form <code>"y=m*x+b"</code>.
474: *
475: * @return A string representation of this line.
476: */
477: public String toString() {
478: if (!Double.isInfinite(slope)) {
479: StringBuffer buffer = new StringBuffer("y= ");
480: if (slope != 0) {
481: buffer.append(slope);
482: buffer.append("*x");
483: if (y0 != 0) {
484: buffer.append(" + ");
485: } else {
486: return buffer.toString();
487: }
488: }
489: buffer.append(y0);
490: return buffer.toString();
491: } else {
492: return "x= " + x0;
493: }
494: }
495:
496: /**
497: * Compare this object with the specified one for equality.
498: */
499: public boolean equals(final Object object) {
500: if (object != null && getClass().equals(object.getClass())) {
501: final Line that = (Line) object;
502: return Double.doubleToLongBits(this .slope) == Double
503: .doubleToLongBits(that.slope)
504: && Double.doubleToLongBits(this .y0) == Double
505: .doubleToLongBits(that.y0)
506: && Double.doubleToLongBits(this .x0) == Double
507: .doubleToLongBits(that.x0);
508: } else {
509: return false;
510: }
511: }
512:
513: /**
514: * Returns a hash code value for this line.
515: */
516: public int hashCode() {
517: final long code = Double.doubleToLongBits(slope) + 37
518: * Double.doubleToLongBits(y0);
519: return (int) code ^ (int) (code >>> 32);
520: }
521:
522: /**
523: * Returns a clone of this line.
524: */
525: public Object clone() {
526: try {
527: return super .clone();
528: } catch (CloneNotSupportedException exception) {
529: throw new AssertionError(exception);
530: }
531: }
532: }
|