/*****
* bezierpatch.cc
* Authors: John C. Bowman and Jesse Frohlich
*
* Render Bezier patches and triangles.
*****/

#include "bezierpatch.h"
#include "predicates.h"

namespace camp {

using ::orient2d;
using ::orient3d;

#ifdef HAVE_LIBGLM

int MaterialIndex;
bool colors;

const double FillFactor=0.1;

void BezierPatch::init(double res)
{
 res2=res*res;

 if(transparent) {
   Epsilon=0.0;
   MaterialIndex=color ? -1-materialIndex : 1+materialIndex;
   pvertex=&vertexBuffer::tvertex;
 } else {
   Epsilon=FillFactor*res;
   MaterialIndex=materialIndex;
   pvertex=&vertexBuffer::vertex;
 }
}

void BezierPatch::render(const triple *p, bool straight, GLfloat *c0)
{
 triple p0=p[0];
 epsilon=0;
 for(unsigned i=1; i < 16; ++i)
   epsilon=max(epsilon,abs2(p[i]-p0));
 epsilon *= DBL_EPSILON;

 triple p3=p[3];
 triple p12=p[12];
 triple p15=p[15];

 triple n0=normal(p3,p[2],p[1],p0,p[4],p[8],p12);
 if(abs2(n0) <= epsilon) {
   n0=normal(p3,p[2],p[1],p0,p[13],p[14],p15);
   if(abs2(n0) <= epsilon) n0=normal(p15,p[11],p[7],p3,p[4],p[8],p12);
 }

 triple n1=normal(p0,p[4],p[8],p12,p[13],p[14],p15);
 if(abs2(n1) <= epsilon) {
   n1=normal(p0,p[4],p[8],p12,p[11],p[7],p3);
   if(abs2(n1) <= epsilon) n1=normal(p3,p[2],p[1],p0,p[13],p[14],p15);
 }

 triple n2=normal(p12,p[13],p[14],p15,p[11],p[7],p3);
 if(abs2(n2) <= epsilon) {
   n2=normal(p12,p[13],p[14],p15,p[2],p[1],p0);
   if(abs2(n2) <= epsilon) n2=normal(p0,p[4],p[8],p12,p[11],p[7],p3);
 }

 triple n3=normal(p15,p[11],p[7],p3,p[2],p[1],p0);
 if(abs2(n3) <= epsilon) {
   n3=normal(p15,p[11],p[7],p3,p[4],p[8],p12);
   if(abs2(n3) <= epsilon) n3=normal(p12,p[13],p[14],p15,p[2],p[1],p0);
 }

 GLuint i0,i1,i2,i3;
 if(color) {
   GLfloat *c1=c0+4;
   GLfloat *c2=c0+8;
   GLfloat *c3=c0+12;

   i0=data.Vertex(p0,n0,c0);
   i1=data.Vertex(p12,n1,c1);
   i2=data.Vertex(p15,n2,c2);
   i3=data.Vertex(p3,n3,c3);

   if(!straight)
     render(p,i0,i1,i2,i3,p0,p12,p15,p3,false,false,false,false,
            c0,c1,c2,c3);
 } else {
   i0=(data.*pvertex)(p0,n0);
   i1=(data.*pvertex)(p12,n1);
   i2=(data.*pvertex)(p15,n2);
   i3=(data.*pvertex)(p3,n3);

   if(!straight)
     render(p,i0,i1,i2,i3,p0,p12,p15,p3,false,false,false,false);
 }

 if(straight) {
   std::vector<GLuint> &q=data.indices;
   triple Pa[]={p0,p12,p15};
   if(!offscreen(3,Pa)) {
     q.push_back(i0);
     q.push_back(i1);
     q.push_back(i2);
   }
   triple Pb[]={p0,p15,p3};
   if(!offscreen(3,Pb)) {
     q.push_back(i0);
     q.push_back(i2);
     q.push_back(i3);
   }
 }
 append();
}

// Use a uniform partition to draw a Bezier patch.
// p is an array of 16 triples representing the control points.
// Pi are the (possibly) adjusted vertices indexed by Ii.
// The 'flati' are flatness flags for each boundary.
void BezierPatch::render(const triple *p,
                        GLuint I0, GLuint I1, GLuint I2, GLuint I3,
                        triple P0, triple P1, triple P2, triple P3,
                        bool flat0, bool flat1, bool flat2, bool flat3,
                        GLfloat *C0, GLfloat *C1, GLfloat *C2, GLfloat *C3)
{
 pair d=Distance(p);
 if(d.getx() < res2 && d.gety() < res2) { // Bezier patch is flat
   triple Pa[]={P0,P1,P2};
   std::vector<GLuint> &q=data.indices;
   if(!offscreen(3,Pa)) {
     q.push_back(I0);
     q.push_back(I1);
     q.push_back(I2);
   }
   triple Pb[]={P0,P2,P3};
   if(!offscreen(3,Pb)) {
     q.push_back(I0);
     q.push_back(I2);
     q.push_back(I3);
   }
 } else { // Patch is not flat
   if(offscreen(16,p)) return;

   /* Control points are indexed as follows:

       Coordinate
      +-----
       Index

       03    13    23    33
      +-----+-----+-----+
      |3    |7    |11   |15
      |     |     |     |
      |02   |12   |22   |32
      +-----+-----+-----+
      |2    |6    |10   |14
      |     |     |     |
      |01   |11   |21   |31
      +-----+-----+-----+
      |1    |5    |9    |13
      |     |     |     |
      |00   |10   |20   |30
      +-----+-----+-----+
       0     4     8     12

   */

   triple p0=p[0];
   triple p3=p[3];
   triple p12=p[12];
   triple p15=p[15];

   if(d.getx() < res2) { // flat in horizontal direction; split vertically
     /*
       P refers to a corner
       m refers to a midpoint
       s refers to a subpatch

       +--------+--------+
       |P3             P2|
       |                 |
       |       s1        |
       |                 |
       |                 |
    m1 +-----------------+ m0
       |                 |
       |                 |
       |       s0        |
       |                 |
       |P0             P1|
       +-----------------+

     */

     Split3 c0(p0,p[1],p[2],p3);
     Split3 c1(p[4],p[5],p[6],p[7]);
     Split3 c2(p[8],p[9],p[10],p[11]);
     Split3 c3(p12,p[13],p[14],p15);

     triple s0[]={p0  ,c0.m0,c0.m3,c0.m5,
                  p[4],c1.m0,c1.m3,c1.m5,
                  p[8],c2.m0,c2.m3,c2.m5,
                  p12 ,c3.m0,c3.m3,c3.m5};

     triple s1[]={c0.m5,c0.m4,c0.m2,p3,
                  c1.m5,c1.m4,c1.m2,p[7],
                  c2.m5,c2.m4,c2.m2,p[11],
                  c3.m5,c3.m4,c3.m2,p15};

     triple n0=normal(s0[12],s0[13],s0[14],s0[15],s0[11],s0[7],s0[3]);
     if(abs2(n0) <= epsilon) {
       n0=normal(s0[12],s0[13],s0[14],s0[15],s0[2],s0[1],s0[0]);
       if(abs2(n0) <= epsilon)
         n0=normal(s0[0],s0[4],s0[8],s0[12],s0[11],s0[7],s0[3]);
     }

     triple n1=normal(s1[3],s1[2],s1[1],s1[0],s1[4],s1[8],s1[12]);
     if(abs2(n1) <= epsilon) {
       n1=normal(s1[3],s1[2],s1[1],s1[0],s1[13],s1[14],s1[15]);
       if(abs2(n1) <= epsilon)
         n1=normal(s1[15],s1[11],s1[7],s1[3],s1[4],s1[8],s1[12]);
     }

     // A kludge to remove subdivision cracks, only applied the first time
     // an edge is found to be flat before the rest of the subpatch is.

     triple m0=0.5*(P1+P2);
     if(!flat1) {
       if((flat1=Straightness(p12,p[13],p[14],p15) < res2)) {
         if(Epsilon)
           m0 -= Epsilon*unit(differential(s1[12],s1[8],s1[4],s1[0]));
       } else m0=s0[15];
     }

     triple m1=0.5*(P3+P0);
     if(!flat3) {
       if((flat3=Straightness(p0,p[1],p[2],p3) < res2)) {
         if(Epsilon)
           m1 -= Epsilon*unit(differential(s0[3],s0[7],s0[11],s0[15]));
       } else m1=s1[0];
     }

     if(color) {
       GLfloat c0[4],c1[4];
       for(size_t i=0; i < 4; ++i) {
         c0[i]=0.5*(C1[i]+C2[i]);
         c1[i]=0.5*(C3[i]+C0[i]);
       }

       GLuint i0=data.Vertex(m0,n0,c0);
       GLuint i1=data.Vertex(m1,n1,c1);

       render(s0,I0,I1,i0,i1,P0,P1,m0,m1,flat0,flat1,false,flat3,C0,C1,c0,c1);
       render(s1,i1,i0,I2,I3,m1,m0,P2,P3,false,flat1,flat2,flat3,c1,c0,C2,C3);
     } else {
       GLuint i0=(data.*pvertex)(m0,n0);
       GLuint i1=(data.*pvertex)(m1,n1);

       render(s0,I0,I1,i0,i1,P0,P1,m0,m1,flat0,flat1,false,flat3);
       render(s1,i1,i0,I2,I3,m1,m0,P2,P3,false,flat1,flat2,flat3);
     }
     return;
   }
   if(d.gety() < res2) { // flat in vertical direction; split horizontally
     /*
       P refers to a corner
       m refers to a midpoint
       s refers to a subpatch

                m1
       +--------+--------+
       |P3      |      P2|
       |        |        |
       |        |        |
       |        |        |
       |        |        |
       |   s0   |   s1   |
       |        |        |
       |        |        |
       |        |        |
       |        |        |
       |P0      |      P1|
       +--------+--------+
                m0
     */

     Split3 c0(p0,p[4],p[8],p12);
     Split3 c1(p[1],p[5],p[9],p[13]);
     Split3 c2(p[2],p[6],p[10],p[14]);
     Split3 c3(p3,p[7],p[11],p15);

     triple s0[]={p0,p[1],p[2],p3,
                  c0.m0,c1.m0,c2.m0,c3.m0,
                  c0.m3,c1.m3,c2.m3,c3.m3,
                  c0.m5,c1.m5,c2.m5,c3.m5};

     triple s1[]={c0.m5,c1.m5,c2.m5,c3.m5,
                  c0.m4,c1.m4,c2.m4,c3.m4,
                  c0.m2,c1.m2,c2.m2,c3.m2,
                  p12,p[13],p[14],p15};

     triple n0=normal(s0[0],s0[4],s0[8],s0[12],s0[13],s0[14],s0[15]);
     if(abs2(n0) <= epsilon) {
       n0=normal(s0[0],s0[4],s0[8],s0[12],s0[11],s0[7],s0[3]);
       if(abs2(n0) <= epsilon)
         n0=normal(s0[3],s0[2],s0[1],s0[0],s0[13],s0[14],s0[15]);
     }

     triple n1=normal(s1[15],s1[11],s1[7],s1[3],s1[2],s1[1],s1[0]);
     if(abs2(n1) <= epsilon) {
       n1=normal(s1[15],s1[11],s1[7],s1[3],s1[4],s1[8],s1[12]);
       if(abs2(n1) <= epsilon)
         n1=normal(s1[12],s1[13],s1[14],s1[15],s1[2],s1[1],s1[0]);
     }

     // A kludge to remove subdivision cracks, only applied the first time
     // an edge is found to be flat before the rest of the subpatch is.

     triple m0=0.5*(P0+P1);
     if(!flat0) {
       if((flat0=Straightness(p0,p[4],p[8],p12) < res2)) {
         if(Epsilon)
           m0 -= Epsilon*unit(differential(s1[0],s1[1],s1[2],s1[3]));
       } else m0=s0[12];
     }

     triple m1=0.5*(P2+P3);
     if(!flat2) {
       if((flat2=Straightness(p15,p[11],p[7],p3) < res2)) {
         if(Epsilon)
           m1 -= Epsilon*unit(differential(s0[15],s0[14],s0[13],s0[12]));
       } else m1=s1[3];
     }

     if(color) {
       GLfloat c0[4],c1[4];
       for(size_t i=0; i < 4; ++i) {
         c0[i]=0.5*(C0[i]+C1[i]);
         c1[i]=0.5*(C2[i]+C3[i]);
       }

       GLuint i0=data.Vertex(m0,n0,c0);
       GLuint i1=data.Vertex(m1,n1,c1);

       render(s0,I0,i0,i1,I3,P0,m0,m1,P3,flat0,false,flat2,flat3,C0,c0,c1,C3);
       render(s1,i0,I1,I2,i1,m0,P1,P2,m1,flat0,flat1,flat2,false,c0,C1,C2,c1);
     } else {
       GLuint i0=(data.*pvertex)(m0,n0);
       GLuint i1=(data.*pvertex)(m1,n1);

       render(s0,I0,i0,i1,I3,P0,m0,m1,P3,flat0,false,flat2,flat3);
       render(s1,i0,I1,I2,i1,m0,P1,P2,m1,flat0,flat1,flat2,false);
     }
     return;
   }
   /*
     Horizontal and vertical subdivision:
     P refers to a corner
     m refers to a midpoint
     s refers to a subpatch

              m2
     +--------+--------+
     |P3      |      P2|
     |        |        |
     |   s3   |   s2   |
     |        |        |
     |        | m4     |
  m3 +--------+--------+ m1
     |        |        |
     |        |        |
     |   s0   |   s1   |
     |        |        |
     |P0      |      P1|
     +--------+--------+
              m0
   */

   // Subdivide patch:
   Split3 c0(p0,p[1],p[2],p3);
   Split3 c1(p[4],p[5],p[6],p[7]);
   Split3 c2(p[8],p[9],p[10],p[11]);
   Split3 c3(p12,p[13],p[14],p15);

   Split3 c4(p0,p[4],p[8],p12);
   Split3 c5(c0.m0,c1.m0,c2.m0,c3.m0);
   Split3 c6(c0.m3,c1.m3,c2.m3,c3.m3);
   Split3 c7(c0.m5,c1.m5,c2.m5,c3.m5);
   Split3 c8(c0.m4,c1.m4,c2.m4,c3.m4);
   Split3 c9(c0.m2,c1.m2,c2.m2,c3.m2);
   Split3 c10(p3,p[7],p[11],p15);

   triple s0[]={p0,c0.m0,c0.m3,c0.m5,c4.m0,c5.m0,c6.m0,c7.m0,
                c4.m3,c5.m3,c6.m3,c7.m3,c4.m5,c5.m5,c6.m5,c7.m5};
   triple s1[]={c4.m5,c5.m5,c6.m5,c7.m5,c4.m4,c5.m4,c6.m4,c7.m4,
                c4.m2,c5.m2,c6.m2,c7.m2,p12,c3.m0,c3.m3,c3.m5};
   triple s2[]={c7.m5,c8.m5,c9.m5,c10.m5,c7.m4,c8.m4,c9.m4,c10.m4,
                c7.m2,c8.m2,c9.m2,c10.m2,c3.m5,c3.m4,c3.m2,p15};
   triple s3[]={c0.m5,c0.m4,c0.m2,p3,c7.m0,c8.m0,c9.m0,c10.m0,
                c7.m3,c8.m3,c9.m3,c10.m3,c7.m5,c8.m5,c9.m5,c10.m5};

   triple m4=s0[15];

   triple n0=normal(s0[0],s0[4],s0[8],s0[12],s0[13],s0[14],s0[15]);
   if(abs2(n0) <= epsilon) {
     n0=normal(s0[0],s0[4],s0[8],s0[12],s0[11],s0[7],s0[3]);
     if(abs2(n0) <= epsilon)
       n0=normal(s0[3],s0[2],s0[1],s0[0],s0[13],s0[14],s0[15]);
   }

   triple n1=normal(s1[12],s1[13],s1[14],s1[15],s1[11],s1[7],s1[3]);
   if(abs2(n1) <= epsilon) {
     n1=normal(s1[12],s1[13],s1[14],s1[15],s1[2],s1[1],s1[0]);
     if(abs2(n1) <= epsilon)
       n1=normal(s1[0],s1[4],s1[8],s1[12],s1[11],s1[7],s1[3]);
   }

   triple n2=normal(s2[15],s2[11],s2[7],s2[3],s2[2],s2[1],s2[0]);
   if(abs2(n2) <= epsilon) {
     n2=normal(s2[15],s2[11],s2[7],s2[3],s2[4],s2[8],s2[12]);
     if(abs2(n2) <= epsilon)
       n2=normal(s2[12],s2[13],s2[14],s2[15],s2[2],s2[1],s2[0]);
   }

   triple n3=normal(s3[3],s3[2],s3[1],s3[0],s3[4],s3[8],s3[12]);
   if(abs2(n3) <= epsilon) {
     n3=normal(s3[3],s3[2],s3[1],s3[0],s3[13],s3[14],s3[15]);
     if(abs2(n3) <= epsilon)
       n3=normal(s3[15],s3[11],s3[7],s3[3],s3[4],s3[8],s3[12]);
   }

   triple n4=normal(s2[3],s2[2],s2[1],m4,s2[4],s2[8],s2[12]);

   // A kludge to remove subdivision cracks, only applied the first time
   // an edge is found to be flat before the rest of the subpatch is.

   triple m0=0.5*(P0+P1);
   if(!flat0) {
     if((flat0=Straightness(p0,p[4],p[8],p12) < res2)) {
       if(Epsilon)
         m0 -= Epsilon*unit(differential(s1[0],s1[1],s1[2],s1[3]));
     } else m0=s0[12];
   }

   triple m1=0.5*(P1+P2);
   if(!flat1) {
     if((flat1=Straightness(p12,p[13],p[14],p15) < res2)) {
       if(Epsilon)
         m1 -= Epsilon*unit(differential(s2[12],s2[8],s2[4],s2[0]));
     } else m1=s1[15];
   }

   triple m2=0.5*(P2+P3);
   if(!flat2) {
     if((flat2=Straightness(p15,p[11],p[7],p3) < res2)) {
       if(Epsilon)
         m2 -= Epsilon*unit(differential(s3[15],s3[14],s3[13],s3[12]));
     } else m2=s2[3];
   }

   triple m3=0.5*(P3+P0);
   if(!flat3) {
     if((flat3=Straightness(p0,p[1],p[2],p3) < res2)) {
       if(Epsilon)
         m3 -= Epsilon*unit(differential(s0[3],s0[7],s0[11],s0[15]));
     } else m3=s3[0];
   }

   if(color) {
     GLfloat c0[4],c1[4],c2[4],c3[4],c4[4];
     for(size_t i=0; i < 4; ++i) {
       c0[i]=0.5*(C0[i]+C1[i]);
       c1[i]=0.5*(C1[i]+C2[i]);
       c2[i]=0.5*(C2[i]+C3[i]);
       c3[i]=0.5*(C3[i]+C0[i]);
       c4[i]=0.5*(c0[i]+c2[i]);
     }

     GLuint i0=data.Vertex(m0,n0,c0);
     GLuint i1=data.Vertex(m1,n1,c1);
     GLuint i2=data.Vertex(m2,n2,c2);
     GLuint i3=data.Vertex(m3,n3,c3);
     GLuint i4=data.Vertex(m4,n4,c4);

     render(s0,I0,i0,i4,i3,P0,m0,m4,m3,flat0,false,false,flat3,C0,c0,c4,c3);
     render(s1,i0,I1,i1,i4,m0,P1,m1,m4,flat0,flat1,false,false,c0,C1,c1,c4);
     render(s2,i4,i1,I2,i2,m4,m1,P2,m2,false,flat1,flat2,false,c4,c1,C2,c2);
     render(s3,i3,i4,i2,I3,m3,m4,m2,P3,false,false,flat2,flat3,c3,c4,c2,C3);
   } else {
     GLuint i0=(data.*pvertex)(m0,n0);
     GLuint i1=(data.*pvertex)(m1,n1);
     GLuint i2=(data.*pvertex)(m2,n2);
     GLuint i3=(data.*pvertex)(m3,n3);
     GLuint i4=(data.*pvertex)(m4,n4);

     render(s0,I0,i0,i4,i3,P0,m0,m4,m3,flat0,false,false,flat3);
     render(s1,i0,I1,i1,i4,m0,P1,m1,m4,flat0,flat1,false,false);
     render(s2,i4,i1,I2,i2,m4,m1,P2,m2,false,flat1,flat2,false);
     render(s3,i3,i4,i2,I3,m3,m4,m2,P3,false,false,flat2,flat3);
   }
 }
}

void BezierTriangle::render(const triple *p, bool straight, GLfloat *c0)
{
 triple p0=p[0];
 epsilon=0;
 for(int i=1; i < 10; ++i)
   epsilon=max(epsilon,abs2(p[i]-p0));

 epsilon *= DBL_EPSILON;

 triple p6=p[6];
 triple p9=p[9];

 triple n0=normal(p9,p[5],p[2],p0,p[1],p[3],p6);
 triple n1=normal(p0,p[1],p[3],p6,p[7],p[8],p9);
 triple n2=normal(p6,p[7],p[8],p9,p[5],p[2],p0);

 GLuint i0,i1,i2;
 if(color) {
   GLfloat *c1=c0+4;
   GLfloat *c2=c0+8;

   i0=data.Vertex(p0,n0,c0);
   i1=data.Vertex(p6,n1,c1);
   i2=data.Vertex(p9,n2,c2);

   if(!straight)
     render(p,i0,i1,i2,p0,p6,p9,false,false,false,c0,c1,c2);
 } else {
   i0=(data.*pvertex)(p0,n0);
   i1=(data.*pvertex)(p6,n1);
   i2=(data.*pvertex)(p9,n2);

   if(!straight)
     render(p,i0,i1,i2,p0,p6,p9,false,false,false);
 }

 if(straight) {
   triple P[]={p0,p6,p9};
   if(!offscreen(3,P)) {
     std::vector<GLuint> &q=data.indices;
     q.push_back(i0);
     q.push_back(i1);
     q.push_back(i2);
   }
 }
 append();
}

// Use a uniform partition to draw a Bezier triangle.
// p is an array of 10 triples representing the control points.
// Pi are the (possibly) adjusted vertices indexed by Ii.
// The 'flati' are flatness flags for each boundary.
void BezierTriangle::render(const triple *p,
                           GLuint I0, GLuint I1, GLuint I2,
                           triple P0, triple P1, triple P2,
                           bool flat0, bool flat1, bool flat2,
                           GLfloat *C0, GLfloat *C1, GLfloat *C2)
{
 if(Distance(p) < res2) { // Bezier triangle is flat
   triple P[]={P0,P1,P2};
   if(!offscreen(3,P)) {
     std::vector<GLuint> &q=data.indices;
     q.push_back(I0);
     q.push_back(I1);
     q.push_back(I2);
   }
 } else { // Triangle is not flat
   if(offscreen(10,p)) return;
   /* Control points are indexed as follows:

      Coordinate
       Index

                                 030
                                  9
                                  /\
                                 /  \
                                /    \
                               /      \
                              /        \
                         021 +          + 120
                          5 /            \ 8
                           /              \
                          /                \
                         /                  \
                        /                    \
                   012 +          +           + 210
                    2 /          111           \ 7
                     /            4             \
                    /                            \
                   /                              \
                  /                                \
                 /__________________________________\
               003         102           201        300
                0           1             3          6


      Subdivision:
                                  P2
                                  030
                                  /\
                                 /  \
                                /    \
                               /      \
                              /        \
                             /    up    \
                            /            \
                           /              \
                       p1 /________________\ p0
                         /\               / \
                        /  \             /   \
                       /    \           /     \
                      /      \  center /       \
                     /        \       /         \
                    /          \     /           \
                   /    left    \   /    right    \
                  /              \ /               \
                 /________________V_________________\
               003               p2                300
               P0                                    P1
   */

   // Subdivide triangle:
   triple l003=p[0];
   triple p102=p[1];
   triple p012=p[2];
   triple p201=p[3];
   triple p111=p[4];
   triple p021=p[5];
   triple r300=p[6];
   triple p210=p[7];
   triple p120=p[8];
   triple u030=p[9];

   triple u021=0.5*(u030+p021);
   triple u120=0.5*(u030+p120);

   triple p033=0.5*(p021+p012);
   triple p231=0.5*(p120+p111);
   triple p330=0.5*(p120+p210);

   triple p123=0.5*(p012+p111);

   triple l012=0.5*(p012+l003);
   triple p312=0.5*(p111+p201);
   triple r210=0.5*(p210+r300);

   triple l102=0.5*(l003+p102);
   triple p303=0.5*(p102+p201);
   triple r201=0.5*(p201+r300);

   triple u012=0.5*(u021+p033);
   triple u210=0.5*(u120+p330);
   triple l021=0.5*(p033+l012);
   triple p4xx=0.5*p231+0.25*(p111+p102);
   triple r120=0.5*(p330+r210);
   triple px4x=0.5*p123+0.25*(p111+p210);
   triple pxx4=0.25*(p021+p111)+0.5*p312;
   triple l201=0.5*(l102+p303);
   triple r102=0.5*(p303+r201);

   triple l210=0.5*(px4x+l201); // =c120
   triple r012=0.5*(px4x+r102); // =c021
   triple l300=0.5*(l201+r102); // =r003=c030

   triple r021=0.5*(pxx4+r120); // =c012
   triple u201=0.5*(u210+pxx4); // =c102
   triple r030=0.5*(u210+r120); // =u300=c003

   triple u102=0.5*(u012+p4xx); // =c201
   triple l120=0.5*(l021+p4xx); // =c210
   triple l030=0.5*(u012+l021); // =u003=c300

   triple l111=0.5*(p123+l102);
   triple r111=0.5*(p312+r210);
   triple u111=0.5*(u021+p231);
   triple c111=0.25*(p033+p330+p303+p111);

   triple l[]={l003,l102,l012,l201,l111,l021,l300,l210,l120,l030}; // left
   triple r[]={l300,r102,r012,r201,r111,r021,r300,r210,r120,r030}; // right
   triple u[]={l030,u102,u012,u201,u111,u021,r030,u210,u120,u030}; // up
   triple c[]={r030,u201,r021,u102,c111,r012,l030,l120,l210,l300}; // center

   triple n0=normal(l300,r012,r021,r030,u201,u102,l030);
   triple n1=normal(r030,u201,u102,l030,l120,l210,l300);
   triple n2=normal(l030,l120,l210,l300,r012,r021,r030);

   // A kludge to remove subdivision cracks, only applied the first time
   // an edge is found to be flat before the rest of the subpatch is.

   triple m0=0.5*(P1+P2);
   if(!flat0) {
     if((flat0=Straightness(r300,p210,p120,u030) < res2)) {
       if(Epsilon)
         m0 -= Epsilon*unit(differential(c[0],c[2],c[5],c[9])+
                            differential(c[0],c[1],c[3],c[6]));
     } else m0=r030;
   }

   triple m1=0.5*(P2+P0);
   if(!flat1) {
     if((flat1=Straightness(l003,p012,p021,u030) < res2)) {
       if(Epsilon)
         m1 -= Epsilon*unit(differential(c[6],c[3],c[1],c[0])+
                            differential(c[6],c[7],c[8],c[9]));
     } else m1=l030;
   }

   triple m2=0.5*(P0+P1);
   if(!flat2) {
     if((flat2=Straightness(l003,p102,p201,r300) < res2)) {
       if(Epsilon)
         m2 -= Epsilon*unit(differential(c[9],c[8],c[7],c[6])+
                            differential(c[9],c[5],c[2],c[0]));
     } else m2=l300;
   }

   if(color) {
     GLfloat c0[4],c1[4],c2[4];
     for(int i=0; i < 4; ++i) {
       c0[i]=0.5*(C1[i]+C2[i]);
       c1[i]=0.5*(C2[i]+C0[i]);
       c2[i]=0.5*(C0[i]+C1[i]);
     }

     GLuint i0=data.Vertex(m0,n0,c0);
     GLuint i1=data.Vertex(m1,n1,c1);
     GLuint i2=data.Vertex(m2,n2,c2);

     render(l,I0,i2,i1,P0,m2,m1,false,flat1,flat2,C0,c2,c1);
     render(r,i2,I1,i0,m2,P1,m0,flat0,false,flat2,c2,C1,c0);
     render(u,i1,i0,I2,m1,m0,P2,flat0,flat1,false,c1,c0,C2);
     render(c,i0,i1,i2,m0,m1,m2,false,false,false,c0,c1,c2);
   } else {
     GLuint i0=(data.*pvertex)(m0,n0);
     GLuint i1=(data.*pvertex)(m1,n1);
     GLuint i2=(data.*pvertex)(m2,n2);

     render(l,I0,i2,i1,P0,m2,m1,false,flat1,flat2);
     render(r,i2,I1,i0,m2,P1,m0,flat0,false,flat2);
     render(u,i1,i0,I2,m1,m0,P2,flat0,flat1,false);
     render(c,i0,i1,i2,m0,m1,m2,false,false,false);
   }
 }
}

std::vector<GLfloat> zbuffer;

void transform(const std::vector<VertexData>& b)
{
 unsigned n=b.size();
 zbuffer.resize(n);

 double Tz0=gl::dView[2];
 double Tz1=gl::dView[6];
 double Tz2=gl::dView[10];
 for(unsigned i=0; i < n; ++i) {
   const GLfloat *v=b[i].position;
   zbuffer[i]=Tz0*v[0]+Tz1*v[1]+Tz2*v[2];
 }
}

// Sort nonintersecting triangles by depth.
int compare(const void *p, const void *P)
{
 unsigned Ia=((GLuint *) p)[0];
 unsigned Ib=((GLuint *) p)[1];
 unsigned Ic=((GLuint *) p)[2];

 unsigned IA=((GLuint *) P)[0];
 unsigned IB=((GLuint *) P)[1];
 unsigned IC=((GLuint *) P)[2];

 return zbuffer[Ia]+zbuffer[Ib]+zbuffer[Ic] <
   zbuffer[IA]+zbuffer[IB]+zbuffer[IC] ? -1 : 1;
}

void sortTriangles()
{
 if(!transparentData.indices.empty()) {
   transform(transparentData.Vertices);
   qsort(&transparentData.indices[0],transparentData.indices.size()/3,
         3*sizeof(GLuint),compare);
 }
}

void Triangles::queue(size_t nP, const triple* P, size_t nN, const triple* N,
                     size_t nC, const prc::RGBAColour* C, size_t nI,
                     const uint32_t (*PP)[3], const uint32_t (*NN)[3],
                     const uint32_t (*CC)[3], bool Transparent)
{
 if(!nN) return;

 data.clear();
 Onscreen=true;
 transparent=Transparent;
 notRendered();

 data.Vertices.resize(nP);

 MaterialIndex=nC ? -1-materialIndex : 1+materialIndex;

 for(size_t i=0; i < nI; ++i) {
   const uint32_t *PI=PP[i];
   uint32_t PI0=PI[0];
   uint32_t PI1=PI[1];
   uint32_t PI2=PI[2];
   triple P0=P[PI0];
   triple P1=P[PI1];
   triple P2=P[PI2];
   const uint32_t *NI=NN[i];
   if(nC) {
     const uint32_t *CI=CC[i];
     prc::RGBAColour C0=C[CI[0]];
     prc::RGBAColour C1=C[CI[1]];
     prc::RGBAColour C2=C[CI[2]];
     GLfloat c0[]={(GLfloat) C0.R,(GLfloat) C0.G,(GLfloat) C0.B,
                   (GLfloat) C0.A};
     GLfloat c1[]={(GLfloat) C1.R,(GLfloat) C1.G,(GLfloat) C1.B,
                   (GLfloat) C1.A};
     GLfloat c2[]={(GLfloat) C2.R,(GLfloat) C2.G,(GLfloat) C2.B,
                   (GLfloat) C2.A};
     transparent |= c0[3]+c1[3]+c2[3] < 3.0;
     data.Vertices[PI0]=VertexData(P0,N[NI[0]],c0);
     data.Vertices[PI1]=VertexData(P1,N[NI[1]],c1);
     data.Vertices[PI2]=VertexData(P2,N[NI[2]],c2);
   } else {
     data.Vertices[PI0]=VertexData(P0,N[NI[0]]);
     data.Vertices[PI1]=VertexData(P1,N[NI[1]]);
     data.Vertices[PI2]=VertexData(P2,N[NI[2]]);
   }
   triple Q[]={P0,P1,P2};
   std::vector<GLuint> &q=data.indices;
   if(!offscreen(3,Q)) {
     q.push_back(PI0);
     q.push_back(PI1);
     q.push_back(PI2);
   }
 }
 append();
}

#endif

} //namespace camp

You are viewing proxied material from ftp.icm.edu.pl. The copyright of proxied material belongs to its original authors. Any comments or complaints in relation to proxied material should be directed to the original authors of the content concerned. Please see the disclaimer for more details.