lucene-java-3.5.0/lucene/contrib/spatial/src/java/org/apache/lucene/spatial/geometry/shape/Ellipse.java

   1 /**
   2  * Licensed to the Apache Software Foundation (ASF) under one or more
   3  * contributor license agreements.  See the NOTICE file distributed with
   4  * this work for additional information regarding copyright ownership.
   5  * The ASF licenses this file to You under the Apache License, Version 2.0
   6  * (the "License"); you may not use this file except in compliance with
   7  * the License.  You may obtain a copy of the License at
   8  *
   9  *     http://www.apache.org/licenses/LICENSE-2.0
  10  *
  11  * Unless required by applicable law or agreed to in writing, software
  12  * distributed under the License is distributed on an "AS IS" BASIS,
  13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14  * See the License for the specific language governing permissions and
  15  * limitations under the License.
  16  */
  17
  18 package org.apache.lucene.spatial.geometry.shape;
  19
  20
  21 /**
  22  * Ellipse shape. From C++ gl.
  23  *
  24  * <p><font color="red"><b>NOTE:</b> This API is still in
  25  * flux and might change in incompatible ways in the next
  26  * release.</font>
  27  */
  28 public class Ellipse implements Geometry2D {
  29   private Point2D center;
  30
  31   /**
  32    * Half length of major axis
  33    */
  34   private double a;
  35
  36   /**
  37    * Half length of minor axis
  38    */
  39   private double b;
  40
  41   private double k1, k2, k3;
  42
  43   /**
  44    * sin of rotation angle
  45    */
  46   private double s;
  47
  48   /**
  49    * cos of rotation angle
  50    */
  51   private double c;
  52
  53   public Ellipse() {
  54     center = new Point2D(0, 0);
  55   }
  56
  57   private double SQR(double d) {
  58     return d * d;
  59   }
  60
  61   /**
  62    * Constructor given bounding rectangle and a rotation.
  63    */
  64   public Ellipse(Point2D p1, Point2D p2, double angle) {
  65     center = new Point2D();
  66
  67     // Set the center
  68     center.x((p1.x() + p2.x()) * 0.5f);
  69     center.y((p1.y() + p2.y()) * 0.5f);
  70
  71     // Find sin and cos of the angle
  72     double angleRad = Math.toRadians(angle);
  73     c = Math.cos(angleRad);
  74     s = Math.sin(angleRad);
  75
  76     // Find the half lengths of the semi-major and semi-minor axes
  77     double dx = Math.abs(p2.x() - p1.x()) * 0.5;
  78     double dy = Math.abs(p2.y() - p1.y()) * 0.5;
  79     if (dx >= dy) {
  80       a = dx;
  81       b = dy;
  82     } else {
  83       a = dy;
  84       b = dx;
  85     }
  86
  87     // Find k1, k2, k3 - define when a point x,y is on the ellipse
  88     k1 = SQR(c / a) + SQR(s / b);
  89     k2 = 2 * s * c * ((1 / SQR(a)) - (1 / SQR(b)));
  90     k3 = SQR(s / a) + SQR(c / b);
  91   }
  92
  93   /**
  94    * Determines if a line segment intersects the ellipse and if so finds the
  95    * point(s) of intersection.
  96    *
  97    * @param seg
  98    *            Line segment to test for intersection
  99    * @param pt0
 100    *            OUT - intersection point (if it exists)
 101    * @param pt1
 102    *            OUT - second intersection point (if it exists)
 103    *
 104    * @return Returns the number of intersection points (0, 1, or 2).
 105    */
 106   public int intersect(LineSegment seg, Point2D pt0, Point2D pt1) {
 107     if (pt0 == null)
 108       pt0 = new Point2D();
 109     if (pt1 == null)
 110       pt1 = new Point2D();
 111
 112     // Solution is found by parameterizing the line segment and
 113     // substituting those values into the ellipse equation.
 114     // Results in a quadratic equation.
 115     double x1 = center.x();
 116     double y1 = center.y();
 117     double u1 = seg.A.x();
 118     double v1 = seg.A.y();
 119     double u2 = seg.B.x();
 120     double v2 = seg.B.y();
 121     double dx = u2 - u1;
 122     double dy = v2 - v1;
 123     double q0 = k1 * SQR(u1 - x1) + k2 * (u1 - x1) * (v1 - y1) + k3
 124         * SQR(v1 - y1) - 1;
 125     double q1 = (2 * k1 * dx * (u1 - x1)) + (k2 * dx * (v1 - y1))
 126         + (k2 * dy * (u1 - x1)) + (2 * k3 * dy * (v1 - y1));
 127     double q2 = (k1 * SQR(dx)) + (k2 * dx * dy) + (k3 * SQR(dy));
 128
 129     // Compare q1^2 to 4*q0*q2 to see how quadratic solves
 130     double d = SQR(q1) - (4 * q0 * q2);
 131     if (d < 0) {
 132       // Roots are complex valued. Line containing the segment does
 133       // not intersect the ellipse
 134       return 0;
 135     }
 136
 137     if (d == 0) {
 138       // One real-valued root - line is tangent to the ellipse
 139       double t = -q1 / (2 * q2);
 140       if (0 <= t && t <= 1) {
 141         // Intersection occurs along line segment
 142         pt0.x(u1 + t * dx);
 143         pt0.y(v1 + t * dy);
 144         return 1;
 145       } else
 146         return 0;
 147     } else {
 148       // Two distinct real-valued roots. Solve for the roots and see if
 149       // they fall along the line segment
 150       int n = 0;
 151       double q = Math.sqrt(d);
 152       double t = (-q1 - q) / (2 * q2);
 153       if (0 <= t && t <= 1) {
 154         // Intersection occurs along line segment
 155         pt0.x(u1 + t * dx);
 156         pt0.y(v1 + t * dy);
 157         n++;
 158       }
 159
 160       // 2nd root
 161       t = (-q1 + q) / (2 * q2);
 162       if (0 <= t && t <= 1) {
 163         if (n == 0) {
 164           pt0.x(u1 + t * dx);
 165           pt0.y(v1 + t * dy);
 166           n++;
 167         } else {
 168           pt1.x(u1 + t * dx);
 169           pt1.y(v1 + t * dy);
 170           n++;
 171         }
 172       }
 173       return n;
 174     }
 175   }
 176
 177   public IntersectCase intersect(Rectangle r) {
 178     // Test if all 4 corners of the rectangle are inside the ellipse
 179     Point2D ul = new Point2D(r.MinPt().x(), r.MaxPt().y());
 180     Point2D ur = new Point2D(r.MaxPt().x(), r.MaxPt().y());
 181     Point2D ll = new Point2D(r.MinPt().x(), r.MinPt().y());
 182     Point2D lr = new Point2D(r.MaxPt().x(), r.MinPt().y());
 183     if (contains(ul) && contains(ur) && contains(ll) && contains(lr))
 184       return IntersectCase.CONTAINS;
 185
 186     // Test if any of the rectangle edges intersect
 187     Point2D pt0 = new Point2D(), pt1 = new Point2D();
 188     LineSegment bottom = new LineSegment(ll, lr);
 189     if (intersect(bottom, pt0, pt1) > 0)
 190       return IntersectCase.INTERSECTS;
 191
 192     LineSegment top = new LineSegment(ul, ur);
 193     if (intersect(top, pt0, pt1) > 0)
 194       return IntersectCase.INTERSECTS;
 195
 196     LineSegment left = new LineSegment(ll, ul);
 197     if (intersect(left, pt0, pt1) > 0)
 198       return IntersectCase.INTERSECTS;
 199
 200     LineSegment right = new LineSegment(lr, ur);
 201     if (intersect(right, pt0, pt1) > 0)
 202       return IntersectCase.INTERSECTS;
 203
 204     // Ellipse does not intersect any edge : since the case for the ellipse
 205     // containing the rectangle was considered above then if the center
 206     // is inside the ellipse is fully inside and if center is outside
 207     // the ellipse is fully outside
 208     return (r.contains(center)) ? IntersectCase.WITHIN
 209         : IntersectCase.OUTSIDE;
 210   }
 211
 212   public double area() {
 213     throw new UnsupportedOperationException();
 214   }
 215
 216   public Point2D centroid() {
 217     throw new UnsupportedOperationException();
 218   }
 219
 220   public boolean contains(Point2D pt) {
 221     // Plug in equation for ellipse, If evaluates to <= 0 then the
 222     // point is in or on the ellipse.
 223     double dx = pt.x() - center.x();
 224     double dy = pt.y() - center.y();
 225     double eq=(((k1 * SQR(dx)) + (k2 * dx * dy) + (k3 * SQR(dy)) - 1));
 226
 227     return eq<=0;
 228   }
 229
 230   public void translate(Vector2D v) {
 231     throw new UnsupportedOperationException();
 232   }
 233
 234 }