import java.util.*;

class Mesh {                      // class for a mesh
  int curfile = -1;
  float orgxsize = 0;
  float orgysize = 0;
  float orgzsize = 0;
  float curscale = 1;
    
  int fileMax=2; 
  String [] filename= new String[fileMax];
   
  int maxArray = 1000;
  TriangleFace[] faces = new TriangleFace[maxArray*3];        //# of possible faces in mesh
  TriangleVertex[] vertices = new TriangleVertex[maxArray];           //# of possible verts in mesh
  TriangleEdge[] edges = new TriangleEdge[maxArray*2];//2 * (3*n*n - 4*n + 1)];    //# of possible edges in mesh
  int triCount = 0;
  int vCount = 0;
  int eCount = 0;
  boolean inflated;
  
  pt ptCenter; 
  Mesh() {
    for(int i=0; i<maxArray; i++) { faces[i] = new TriangleFace(); }
    for(int i=0; i<maxArray; i++) { vertices[i] = new TriangleVertex(); }
    for(int i=0; i<maxArray; i++) { edges[i] = new TriangleEdge(); }
    
    curfile = -1;
    filename[0]="bunny.vts"; filename[1]="horse.vts"; 
    ptCenter = new pt();
  } 
  void setMesh(Mesh m)
  {
    curfile = m.curfile > -1 ? m.curfile-1 : -1;
    curscale = m.curscale;
    println("set mesh : " + curfile + ", " + curscale);
    loadMesh();   
  }
  
  void makeLarger()
  {
    curscale *= 1.25;
    for(int k=0; k<vCount; k++) 
    {
      vertices[k].vert.mul(1.25);
    }
  }
  
  void makeSmaller()
  {
    curscale *= 0.8;
    // rearrage vertex's coordinate
    for(int k=0; k<vCount; k++) 
    {
      vertices[k].vert.mul(0.8);
    }
 }
  
  void AddFace(int a, int b, int c) {
    faces[triCount].set(a, b, c);
    
    pt A = vertices[a].vert;
    pt B = vertices[b].vert;
    pt C = vertices[c].vert;
    
    // Determine if face is valid
    //if (A.z != 0.0 || B.z != 0.0 || C.z != 0.0) { faces[triCount].valid = true; vertices[a].valid = true; vertices[b].valid = true; vertices[c].valid = true; }
    
    // Calculate face normal
    vec AB = new vec(B.x-A.x, B.y-A.y, B.z-A.z);
    vec AC = new vec(C.x-A.x, C.y-A.y, C.z-A.z);
    faces[triCount].faceNormal = AB.cross(AC);
    faces[triCount].faceNormal.div(faces[triCount].faceNormal.norm());  //normalize    
    triCount++;
  }
  
  void AddEdge(int a, int b) {
    int t1=-1, t2=-1;
    
    for ( int i = 0; i < eCount; i++ )
    {
       if(edges[i].contains(a,b))
         return;
    }
    
    for (int i=0; i<triCount; i++) {      
      if (faces[i].contains(a,b)) {
        if (t1==-1) { t1 = i; } else if (t2==-1) { t2 = i; break; }
       }
    }    
    
    if (t1!=-1 && t2!=-1) {
      edges[eCount].set(a,b,t1,t2);
      eCount++; 
    }
    else
    {
       println("error in add edges");
    }
  }
  
  vec edgeTangent(TriangleEdge e) {
    pt A = vertices[e.A].vert;
    pt B = vertices[e.B].vert;
    vec tangent = new vec(B.x-A.x, B.y-A.y, B.z-A.z);
    tangent.makeUnit();
    return tangent;    
  }

  void loadMesh() {
    // initialize
    triCount = 0;
    vCount = 0;
    eCount = 0;
    //  String [] ss = loadStrings(fn[fni]);
    curfile++;
    String [] ss = loadStrings(filename[curfile%fileMax]);
    
   // println("length = " + ss.length);
    int s=0;   int comma1, comma2;   int a, b, c;
    float[] coordinate = new float[3];
    
    // intialize vertices' information
    vCount = int(ss[s++]);
    print("number of vertices = "+vCount);
    
    // 0 : x, 1 : y, 2 : z
    int[] smallestidx = new int[3];
    int[] largestidx = new int[3];
    float[] smallest = new float[3];
    float[] largest = new float[3];
    for ( int idx = 0; idx < 3; idx++ )
    {
      smallestidx[idx] = -1;
      largestidx[idx] = -1;
      smallest[idx] = 100000;
      largest[idx] = -100000;
    }    
    for(int k=0; k<vCount; k++) 
    {
      int i=k+s; 
      comma1=ss[i].indexOf(',');   
      coordinate[0]=float(ss[i].substring(0, comma1));
      String rest = ss[i].substring(comma1+1, ss[i].length());
      comma2=rest.indexOf(',');    
      coordinate[1]=float(rest.substring(0, comma2)); 
      coordinate[2]=float(rest.substring(comma2+1, rest.length()));
      vertices[k].vert.setTo(new pt(coordinate[0], coordinate[1], coordinate[2]));
      
      for ( int idx = 0; idx < 3 ; idx++ )
      {
        if ( coordinate[idx] < smallest[idx] )
        {
          smallest[idx] = coordinate[idx];
          smallestidx[idx] = k;    
        }
        
        if ( coordinate[idx] > largest[idx] )
        {
          largest[idx] = coordinate[idx];
          largestidx[idx] = k;         
        }
      }
      
   };
   
   // we can calculate the size of this model
   orgxsize = largest[0] - smallest[0];
   orgysize = largest[1] - smallest[1];
   orgzsize = largest[2] - smallest[2];
   float largestsize = orgxsize > orgysize ? orgxsize : orgysize;
   largestsize = largestsize > orgzsize ? largestsize : orgzsize;
   
   float scalevalue = 200 / largestsize; // we need a model which size is 200
   orgxsize = orgxsize * scalevalue;
   orgysize = orgysize * scalevalue;
   orgzsize = orgzsize * scalevalue; 
    // rearrage vertex's coordinate
    // translate all vertex to make a vertex which y value is smallest be origin
    pt ptFirst = new pt();
    ptFirst.setTo(vertices[smallestidx[1]].vert);
    for(int k=0; k<vCount; k++) 
    {
      vertices[k].vert.sub(ptFirst);
      vertices[k].vert.mul(scalevalue*curscale);
    }
   
   // initialize triangles' information
   s=vCount+1;
   int nt = int(ss[s]); 
   println(", number of triangles ="+nt);
   s++;
   for(int k=0; k<nt; k++) 
   {
      int i=k+s;
      comma1=ss[i].indexOf(',');   
      a=int(ss[i].substring(0, comma1));  
      String rest = ss[i].substring(comma1+1, ss[i].length()); comma2=rest.indexOf(',');  
      b=int(rest.substring(0, comma2)); c=int(rest.substring(comma2+1, rest.length()));
      if ( a >= vCount || b >= vCount || c >= vCount )
      {
         println("Error in Triangle : " + a + ", " +  b + ", " + c);        
      }
      AddFace(a,b,c);
      vertices[a].AddAdj(b);
      vertices[a].AddAdj(c);
      vertices[b].AddAdj(c);
      vertices[b].AddAdj(a);
      vertices[c].AddAdj(a);
      vertices[c].AddAdj(b);
   
      //V[3*k]=a;  V[3*k+1]=b;  V[3*k+2]=c;
    }
    
    //Set up edge list
    for (int k=0; k < triCount; k++) {
      TriangleFace f = faces[k];
      AddEdge(f.A, f.B);
      AddEdge(f.B, f.C);
      AddEdge(f.C, f.A);
    }
    
    println("number of edges = " + eCount);
    
    // center of mass
    pt P = CenterOfMass();
    ptCenter.setTo(P);
    //fill(red);
    //noStroke();
    //ptCenter.show(10);
    for(int k=0; k<vCount; k++) 
    {
      vertices[k].vert.transformPtToY(ptCenter);      
    }
    
    //P = CenterOfMass();
   // ptCenter.setTo(P);
    //ptCenter.show(10);
    //ptCenter.write();
    
    
    
      
  }; 
  
  // from http://www.melax.com/volint/
  pt CenterOfMass()
  {
     pt cen = new pt();
     float volume = 0;
     for( int i = 0;  i < triCount; i++ )
     {
       TriangleFace f = faces[i];
       pt A = vertices[f.A].vert;
       pt B = vertices[f.B].vert;
       pt C = vertices[f.C].vert;
       
       float vol = det(A, B, C);
       //println("vol = " + vol);
       pt P = new pt();
       P.add(A);
       P.add(B);
       P.add(C);
       P.mul(vol);
       cen.add(P);
       volume += vol;   
     }
    
    cen.div(volume*4.0f);
    return cen;
  }
  
  
  void draw(color meshcolor)
  {
    fill (meshcolor);   
    //noFill();
    noStroke();
    for (int i=0;i<triCount;i++) 
    {      
        beginShape(TRIANGLES);
        vertices[faces[i].A].vert.vert();
        vertices[faces[i].B].vert.vert();
        vertices[faces[i].C].vert.vert();
        endShape();
    }
    //ptCenter.show(10);
   // ptCenter.write();
  }
};

class MorphResult {
   int morphstep = 100;
   int curstepidx = 0;
   TriangleFace[] faces = new TriangleFace[10000];  //list of triangles -> A,B,C are indexes of the list list stored in resultPair
   boolean direction = true;
   Hashtable  resultTable;

   int triCount;
  
   MorphResult() { 
     triCount = 0; 
     for(int i=0; i<faces.length;i++) faces[i] = new TriangleFace();
    // for(int i=0; i<resultPair.length;i++) { resultPair[i] = new int[2]; resultPair[i][0] = 0; resultPair[i][1] = 0; }
     
     resultTable = new Hashtable(5000);
   }
   
   void reset()
   {
      triCount = 0;
      resultTable.clear();
   }   
   
   int AddFace(int A, int B, int C) { 
     if (triCount >= faces.length) {
       faces = (TriangleFace[])expand(faces);
       for (int i=triCount; i<faces.length; i++) faces[i] = new TriangleFace();
     }
     faces[triCount].A = A; faces[triCount].B = B; faces[triCount].C = C; triCount++; 
     return (triCount - 1);
   }
   int AddVertexPair(Mesh m1, Mesh m2, int A, int B) {
     int keyValue = A * 1000 + B; // our convension
     if ( resultTable.containsKey(new Integer(keyValue)) == false)
     {
         pt A1 = m1.vertices[A].vert;
         pt A2 = m2.vertices[B].vert;
     
         pt[] resultVertex = new pt[morphstep];
         float stepsize = 1.0/morphstep;
         float morph_value = 0;
     
         for ( int j=0; j < morphstep; j++)
         {
           resultVertex[j] = new pt();
           resultVertex[j].setTo(lerp(A1.x,A2.x,morph_value), lerp(A1.y,A2.y,morph_value), lerp(A1.z,A2.z,morph_value));
           morph_value += stepsize;
         }
         resultTable.put(new Integer(keyValue), resultVertex);
     }
     //println("result num = " + resultMap.size());
     return keyValue;    
   }
   
   void draw()
   {
     if ( curstepidx >= morphstep )
     {
       curstepidx = morphstep - 1;
       direction = !direction; 
    }
    else if (curstepidx < 0 )
    {
       curstepidx = 0;
       direction = !direction; 
    }
    
   fill (blue); 
    noStroke();  
    for (int i=0; i<triCount; i++) 
     {
      TriangleFace face = faces[i];
      
      //pt A1 = m1.vertices[r.resultPair[face.A][0]].vert;
      // pt A2 = m2.vertices[r.resultPair[face.A][1]].vert;
      //pt A = new pt(lerp(A1.x,A2.x,morph_t), lerp(A1.y,A2.y,morph_t), lerp(A1.z,A2.z,morph_t));
      
      //pt B1 = m1.vertices[r.resultPair[face.B][0]].vert;
      //pt B2 = m2.vertices[r.resultPair[face.B][1]].vert;
      //pt B = new pt(lerp(B1.x,B2.x,morph_t), lerp(B1.y,B2.y,morph_t), lerp(B1.z,B2.z,morph_t));
      
      //pt C1 = m1.vertices[r.resultPair[face.C][0]].vert;
      //pt C2 = m2.vertices[r.resultPair[face.C][1]].vert;
      //pt C = new pt(lerp(C1.x,C2.x,morph_t), lerp(C1.y,C2.y,morph_t), lerp(C1.z,C2.z,morph_t));
      
      
      pt[] A = (pt[])resultTable.get(new Integer(face.A));
      pt[] B = (pt[])resultTable.get(new Integer(face.B));
      pt[] C = (pt[])resultTable.get(new Integer(face.C));
      
      beginShape(TRIANGLES);    
      A[curstepidx].vert(); B[curstepidx].vert(); C[curstepidx].vert();
      endShape();
    }
    
    if ( direction )
    {
       curstepidx++;
    }
    else   
    {
      curstepidx--;
    }
   }
};

class TriangleVertex {
  pt vert;
  int maxAdj = 50;
  int[] adjacencies = new int[maxAdj];
  int curAdj = 0;
  TriangleVertex() {
    vert = new pt();
    for (int i = 0; i < maxAdj; i++) { adjacencies[i] = -1; }
  } 
  void AddAdj(int adj)
  {
    if ( curAdj >= maxAdj )
    {
      println("error in AddAdj");
    }
    boolean exist = false;
    for ( int i = 0; i < curAdj; i++ )
    {
      if ( adjacencies[i] == adj )
      {
         exist = true;
         break;
      }
    }
    if ( !exist) 
    { 
      adjacencies[curAdj] = adj;
      curAdj++;
    }
  } 
}

class TriangleFace {
  int A, B, C;
  vec faceNormal;
  
  TriangleFace()
  {
     A = 0; B = 0; C = 0;
     faceNormal = new vec(0,0,0);
  }
  
  void set(int a, int b, int c)
  {
    A=a; B=b; C=c;
  }
    
  boolean contains(int p1, int p2) {
    return ((A==p1 || B==p1 || C==p1) && (A==p2 || B==p2 || C==p2));
  }
};

class TriangleEdge {
  int A, B;
  int tri1, tri2;
  //boolean visited;
  
  TriangleEdge()
  {
    A = -1; B = -1; tri1 = -1; tri2 = -1;
  }
  boolean contains(int p1, int p2)
  {
     if ( (A == p1 && B == p2) || ( A == p2 && B == p1) )
     {
        return true;
     }
     
     return false;
  }
  void set(int a, int b, int t1, int t2) {
    A = a; B = b; tri1 = t1; tri2 = t2;
    //println("Edge1 of tri1= " + A + ", " + tri1 + ", Edge2 of tri2=" + B + ", "+tri2);
  }
}

// from http://www.melax.com/volint/
float det(pt A, pt B, pt C)
{
   //print("A=");A.write();
   //print("B=");B.write();
   //print("C=");C.write();
   float v = A.x*(B.y*C.z-B.z*C.y)-B.x*(A.y*C.z-A.z*C.y)+C.x*(A.y*B.z-A.z*B.y);
   //println("det value = " + v);
   return v;
}