//**************************************
// 3D POINTS AND VECTOR CLASSES
// written by Jarek Rossignac
// June 2006
//*************************************
class p3 { float x,y,z; 
  p3 (float px, float py, float pz) {x = px; y = py; z = pz;};
  p3 make() {return(new p3(x,y,z));};
  void show(float r) { pushMatrix(); translate(x,y,z); sphere(r); popMatrix();}; 
  void showLineTo (p3 P) {line(x,y,z,P.x,P.y,P.z); }; 
  void setToPoint(p3 P) { x = P.x; y = P.y; z = P.z;}; 
  void setTo(p3 P) { x = P.x; y = P.y; z = P.z;}; 
  void setTo (float px, float py, float pz) {x = px; y = py; z = pz;}; 
  void setToMouse() { x = mouseX; y = mouseY; }; 
  void write() {println("("+x+","+y+","+z+")");};
  void addVec(v3 V) {x += V.x; y += V.y; z += V.z;};
  void addScaledVec(float s, v3 V) {x += s*V.x; y += s*V.y; z += s*V.z;};
  void subVec(v3 V) {x -= V.x; y -= V.y; z -= V.z;};
  void vert() {vertex(x,y,z);};
  void vertext(float u, float v) {vertex(x,y,z,u,v);};
  boolean isInWindow() {return(((x<0)||(x>width)||(y<0)||(y>height)));};
  void label(String s, v3 D) {text(s, x+D.x, y+D.y, z+D.z);  };
  v3 vecTo(p3 P) {return(new v3(P.x-x,P.y-y,P.z-z)); };
  float disTo(p3 P) {return(sqrt( sq(P.x-x)+sq(P.y-y)+sq(P.z-z) )); };
  v3 vecToMid(p3 P, p3 Q) {return(new v3((P.x+Q.x)/2.0-x,(P.y+Q.y)/2.0-y,(P.z+Q.z)/2.0-z )); };
  v3 vecToProp (p3 B, p3 D) {
      v3 CB = this.vecTo(B); float LCB = CB.norm();
      v3 CD = this.vecTo(D); float LCD = CD.norm();
      v3 U = CB.make();
      v3 V = CD.make(); V.sub(U); V.mul(LCB/(LCB+LCD));
      U.add(V);
      return(U);  
      };  
  void addPt(p3 P) {x+=P.x; y+=P.y; z+=P.z;};
  void subPt(p3 P) {x-=P.x; y-=P.y; z-=P.z; };
  void mul(float f) {x*=f; y*=f; y*=f;};
  void pers(float d) { y=d*y/(d+z); x=d*x/(d+z); z=d*z/(d+z); };
  void inverserPers(float d) { y=d*y/(d-z); x=d*x/(d-z); z=d*z/(d-z); };
  boolean coplanar (p3 A, p3 B, p3 C) {return(abs(tetVol(this,A,B,C))<0.0001);};
  boolean cw (p3 A, p3 B, p3 C) {return(tetVol(this,A,B,C)>0.0001);};
  } ;
 
class v3 { float x,y,z; 
  v3 (float px, float py, float pz) {x = px; y = py; z = pz;};
  void setTo (float px, float py, float pz) {x = px; y = py; z = pz;}; 
  v3 make() {return(new v3(x,y,z));};
  void setTo(v3 V) { x = V.x; y = V.y; z = V.z;}; 
  void show (p3 P) {line(P.x,P.y, P.z,P.x+x,P.y+y,P.z+z); }; 
  void add(v3 V) {x += V.x; y += V.y; z += V.z;};
  void addScaled(float m, v3 V) {x += m*V.x; y += m*V.y; z += m*V.z;};
  void sub(v3 V) {x -= V.x; y -= V.y; z -= V.z;};
  void mul(float m) {x *= m; y *= m; z *= m;};
  void div(float m) {x /= m; y /= m; z /= m;};
  void write() {println("("+x+","+y+","+z+")");};
  float norm() {return(sqrt(sq(x)+sq(y)+sq(z)));}; 
  void makeUnit() {float n=this.norm(); if (n>0.0001) {this.div(n);};};
  void back() {x= -x; y= -y; z= -z;};
  boolean coplanar (v3 V, v3 W) {return(abs(mixed(this,V,W))<0.0001);};
  boolean cw (v3 U, v3 V, v3 W) {return(mixed(this,V,W)>0.0001);};
  } ;
  
v3 triNormalFromPts(p3 A, p3 B, p3 C) {v3 N = cross(A.vecTo(B),A.vecTo(C));  return(N); };
float tetVol (p3 A, p3 B, p3 C, p3 D) { return(dot(triNormalFromPts(A,B,C),A.vecTo(D))); };
float dot(v3 U, v3 V) {return(U.x*V.x+U.y*V.y+U.z*V.z); };
v3 cross(v3 U, v3 V) {return(new v3( U.y*V.z-U.z*V.y, U.z*V.x-U.x*V.z, U.x*V.y-U.y*V.x )); };
float mixed(v3 U, v3 V, v3 W) {return(dot(cross(U,V),W)); };
p3 average (p3 A, p3 B) {return(new p3((A.x+B.x)/2 , (A.y+B.y)/2, (A.z+B.z)/2 )); };
p3 average (p3 A, p3 B, p3 C) {return(new p3((A.x+B.x+C.x)/3 , (A.y+B.y+C.y)/3, (A.z+B.z+C.z)/3 )); };
p3 average (p3 A, p3 B, p3 C, p3 D) {return(new p3( (A.x+B.x+C.x+D.x)/4 , (A.y+B.y+C.y+D.y)/4, (A.z+B.z+C.z+D.z)/4 ) ); };
p3 between (p3 A, float s, p3 B) {return(new p3((s-1)*A.x+s*B.x , (s-1)*A.y+s*B.y,(s-1)*A.z+s*B.z )); };
v3 between (v3 A, float s, v3 B) {return(new v3((s-1)*A.x+s*B.x , (s-1)*A.y+s*B.y,(s-1)*A.z+s*B.z )); };
v3 dif(p3 A, p3 B) {return(new v3( B.x-A.x , B.y-A.y , B.z-A.z)); };
v3 dif(v3 U, v3 V) {return(new v3(V.x-U.x,V.y-U.y,V.z-U.z)); };
v3 sum(v3 U, v3 V) {return(new v3(V.x+U.x,V.y+U.y,V.z+U.z)); };
v3 average(v3 U, v3 V) {return(new v3((U.x+V.x)/2,(U.y+V.y)/2,(U.z+V.z)/2)); };
v3 average (v3 A, v3 B, v3 C, v3 D) {return(new v3( (A.x+B.x+C.x+D.x)/4 , (A.y+B.y+C.y+D.y)/4, (A.z+B.z+C.z+D.z)/4 ) ); };