src/edu/berkeley/qfat/geom/Plane.java

   1 package edu.berkeley.qfat.geom;
   2 import javax.media.opengl.*;
   3
   4 // FEATURE: represent a plane by the matrix that projects points onto that plane?
   5 public class Plane implements AffineConstraint {
   6
   7     // ax+by+cz=d
   8     public final float a, b, c, d;
   9     public final Point p;
  10
  11     public Vec norm() { return new Vec(a, b, c).norm(); }
  12
  13     public Plane(Point p, Vec norm) {
  14         norm = norm.norm();
  15         this.a = norm.x;
  16         this.b = norm.y;
  17         this.c = norm.z;
  18         this.d = p.x*norm.x+p.y*norm.y+p.z*norm.z;
  19         this.p = p;
  20     }
  21
  22     public Plane(float a, float b, float c, float d) {
  23         // FIXME, what if c==0?
  24         this(a, b, c, d, new Point(0, 0, d/c));
  25     }
  26
  27     /** provided at least one of a,b,c is nonzero, return the Plane representing ax+by+cz=d */
  28     public Plane(float a, float b, float c, float d, Point p) {
  29         this.a = a;
  30         this.b = b;
  31         this.c = c;
  32         this.d = d;
  33         this.p = p;
  34     }
  35
  36     public Point intersect(Plane p1, Plane p2) {
  37         Vec norm = norm();
  38         Plane p3 = this;
  39         float z = p1.norm().dot(p2.norm().cross(p3.norm()));
  40         if (Math.abs(z) < 0.0001) return null;  // planes do not intersect at a point
  41         Vec v1 = p2.norm().cross(p3.norm()).times(-1 * p1.d);
  42         Vec v2 = p3.norm().cross(p1.norm()).times(-1 * p2.d);
  43         Vec v3 = p1.norm().cross(p2.norm()).times(-1 * p3.d);
  44         return Point.ZERO.plus(v1.plus(v2).plus(v3).times(1/z));
  45     }
  46
  47     public AffineConstraint intersect(AffineConstraint c, float epsilon) {
  48         if (c instanceof Plane) {
  49             Plane p = (Plane)c;
  50
  51             // same plane
  52             if (Math.abs(p.a*this.d-this.a*p.d) <= epsilon &&
  53                 Math.abs(p.b*this.d-this.b*p.d) <= epsilon &&
  54                 Math.abs(p.c*this.d-this.c*p.d) <= epsilon &&
  55                 Math.abs(p.d*this.d-this.d*p.d) <= epsilon)
  56                 return this;
  57
  58             // parallel planes
  59             if (Math.abs(p.norm().cross(norm()).mag()) <= epsilon)
  60                 return AffineConstraint.NONE;
  61
  62             Vec u = norm().cross(p.norm());
  63             Point point = null;
  64             if (Math.abs(u.z) >= Math.abs(u.y) && Math.abs(u.z) >= Math.abs(u.y)) {
  65                 point = new Point( (this.b*p.d - p.b*this.d)/(this.a*p.b - p.a*this.b),
  66                                    (this.d*p.a - p.d*this.a)/(this.a*p.b - p.a*this.b),
  67                                    0);
  68             } else if (Math.abs(u.y) >= Math.abs(u.z) && Math.abs(u.y) >= Math.abs(u.x)) {
  69                 point = new Point( (this.c*p.d - p.c*this.d)/(this.a*p.c - p.a*this.c),
  70                                    0,
  71                                    (this.d*p.a - p.d*this.a)/(this.a*p.c - p.a*this.c));
  72             } else {
  73                 point = new Point( 0,
  74                                    (this.c*p.d - p.c*this.d)/(this.b*p.c - p.b*this.c),
  75                                    (this.d*p.b - p.d*this.b)/(this.b*p.c - p.b*this.c));
  76             }
  77
  78             return new Line(point, u);
  79
  80         } else if (c instanceof Line) {
  81             Line l = (Line)c;
  82             if (Math.abs(a + b*l.m + this.c*l.n)     <= epsilon &&
  83                 Math.abs(b * l.c   + this.c*l.d + d) <= epsilon)
  84                 return l;
  85
  86             // FIXME: parallel case
  87
  88             float x = (d - b*l.c - this.c*l.d) / (a + b*l.m + this.c*l.n);
  89             float y = l.m*x+l.c;
  90             float z = l.n*x+l.d;
  91             return new Point(x,y,z);
  92             //throw new RuntimeException("not yet implemented");
  93
  94         } else
  95             return c.intersect(this, epsilon);
  96     }
  97
  98     public AffineConstraint multiply(Matrix m) {
  99         return new Plane( m.times(p), m.times(norm()).norm() );
 100     }
 101
 102     public Point getProjection(Point p) {
 103         Point ret = norm().times(p.minus(this.p).dot(norm())).times(-1).plus(p);
 104         return ret;
 105     }
 106     public String toString() {
 107         return "[plane "+a+"x+"+b+"y+"+c+"z="+d+"]";
 108     }
 109 }