Logo Search packages:      
Sourcecode: yafray version File versions

triangletools.cc

#include<limits>
#include "triangletools.h"
#include "vector2d.h"

using namespace std;

__BEGIN_YAFRAY

struct square_t
{
      bool isInside(const point2d_t &p)const
      {
            return
                  p.x >= x1 && p.x <= x2 &&
                  p.y >= y1 && p.y <= y2 ;
      };
      bool isInside(const point3d_t &p)const
      {
            return
                  p.x >= x1 && p.x <= x2 &&
                  p.y >= y1 && p.y <= y2 ;
      };

      PFLOAT x1,x2,y1,y2;
};

class planeEquation_t
{
      public:
            planeEquation_t(const point3d_t &a,const vector3d_t &N)
            {
                  A=-N.x;
                  B=-N.y;
                  C=toVector(a)*N;
                  PFLOAT temp=1.0/N.z;
                  A*=temp;
                  B*=temp;
                  C*=temp;
            };
            PFLOAT getZ(PFLOAT x,PFLOAT y)const { return A*x+B*y+C;};
      protected:
            PFLOAT A,B,C;
};

static bool intersectY(const point2d_t &a,const point2d_t &b,PFLOAT y,PFLOAT x1,PFLOAT x2,point2d_t &res)
{
      res.y=y;
      PFLOAT diffX=b.x-a.x;
      if(diffX==0.0)
      {
            res.x=a.x;
            if(!((res.x>=x1) && (res.x<=x2))) return false;
            if(a.y<b.y) return ((res.y>=a.y) && (res.y<=b.y));
            else return ((res.y>=b.y) && (res.y<=a.y));
      }
      else
      {
            PFLOAT diffY=b.y-a.y;
            if(diffY==0.0) return false;
            PFLOAT lambda=(y-a.y)/diffY;
            if((lambda<0) || (lambda>1.0)) return false;
            res.x=lambda*diffX+a.x;
            return ((res.x>=x1) && (res.x<=x2));
      }
}

static bool intersectX(const point2d_t &a,const point2d_t &b,PFLOAT x,PFLOAT y1,PFLOAT y2,point2d_t &res)
{
      res.x=x;
      PFLOAT diffY=b.y-a.y;
      if(diffY==0.0)
      {
            res.y=a.y;
            if(!((res.y>=y1) && (res.y<=y2))) return false;
            if(a.x<b.x) return ((res.x>=a.x) && (res.x<=b.x));
            else return ((res.x>=b.x) && (res.x<=a.x));
      }
      else
      {
            PFLOAT diffX=b.x-a.x;
            if(diffX==0.0) return false;
            PFLOAT lambda=(x-a.x)/diffX;
            if((lambda<0) || (lambda>1.0)) return false;
            res.y=lambda*diffY+a.y;
            return ((res.y>=y1) && (res.y<=y2));
      }
}


template<class F>
bool applyVectorIntersect(const point2d_t &a,const point2d_t &b,
            const square_t &q,F &func)
{
      point2d_t pr;
      int c=0;
      if(intersectX(a, b, q.x1, q.y1, q.y2, pr)) {if(!func(pr)) return false;c++;}
      if(intersectX(a, b, q.x2, q.y1, q.y2, pr)) {if(!func(pr)) return false;c++;}
      if(c>1) return true;
      if(intersectY(a, b, q.y1, q.x1, q.x2, pr)) {if(!func(pr)) return false;c++;}
      if(c>1) return true;
      if(intersectY(a, b, q.y2, q.x1, q.x2, pr)) if(!func(pr)) return false;
      return true;
}

static bool isInside(const point2d_t &p,
                      const point2d_t &a,
                      const point2d_t &b,
                      const point2d_t &c)
{
             vector2d_t w1(b.y-a.y,-(b.x-a.x));
             vector2d_t pa=p-a;
             vector2d_t ca=c-a;
             if((w1*pa) * (w1*ca) < 0) return false;

             vector2d_t w2(c.y-b.y,-(c.x-b.x));
             vector2d_t pb=p-b;
             vector2d_t ab=a-b;
             if((w2*pb) * (w2*ab) < 0) return false;
             
             vector2d_t w3(a.y-c.y,-(a.x-c.x));
             vector2d_t pc=p-c;
             vector2d_t bc=b-c;
             return ((w3*pc) * (w3*bc) >= 0); 
}


template<class F>
bool intersectApply(const point2d_t &a,const point2d_t &b,
                                                            const point2d_t &c,const square_t &q,F &func)
{
      if(!applyVectorIntersect(a, b, q, func)) return false;
      if(!applyVectorIntersect(b, c, q, func)) return false;
      if(!applyVectorIntersect(c, a, q, func)) return false;
      if (q.isInside(a)) if(!func(a)) return false;
      if (q.isInside(b)) if(!func(b)) return false;
      if (q.isInside(c)) if(!func(c)) return false;

  point2d_t qp ( q.x1, q.y1 );
      if (isInside(qp, a,b,c)) if(!func(qp)) return false;
  qp.set( q.x2, q.y1 );
      if (isInside(qp, a,b,c)) if(!func(qp)) return false;
  qp.set( q.x2, q.y2 );
      if (isInside(qp, a,b,c)) if(!func(qp)) return false;
  qp.set( q.x1, q.y2 );
      if (isInside(qp, a,b,c)) if(!func(qp)) return false;
      return true;
}

struct checkPosition_f
{
      checkPosition_f(const point3d_t &a,const vector3d_t &N,PFLOAT z):
            plane(a,N),Z(z),decision(trianglePosition_t::NONE) {};

      bool operator () (const point2d_t & p)
      {
            PFLOAT z=plane.getZ(p.x,p.y);
            if(z==Z) decision=trianglePosition_t::INTERSECT;
            else if(decision==trianglePosition_t::NONE)
            {
                  if(z<Z) decision=trianglePosition_t::LOWER;
                  else decision=trianglePosition_t::HIGHER;
            }
            else
                  if( ((z<Z) && (decision==trianglePosition_t::HIGHER)) || 
                              ((z>Z) && (decision==trianglePosition_t::LOWER)) ) 
                        decision=trianglePosition_t::INTERSECT;
            
            if(decision==trianglePosition_t::INTERSECT) return false;
            else return true;
      };
      planeEquation_t plane;
      PFLOAT Z;
      int decision;
};

static PFLOAT exitDist(const point3d_t &from,const vector3d_t &ray,PFLOAT dist,const square_t &q)
{
      PFLOAT xmin=-numeric_limits<PFLOAT>::infinity(),xmax=numeric_limits<PFLOAT>::infinity();
      PFLOAT ymin=-numeric_limits<PFLOAT>::infinity(),ymax=numeric_limits<PFLOAT>::infinity();
      
      if(ray.x!=0)
      {
            xmin=(q.x1-from.x)/ray.x;
            xmax=(q.x2-from.x)/ray.x;
            if(xmin>xmax) swap(xmax,xmin);
      }
      else if( (from.y<q.y1) || (from.y>q.y2)) return -1;
      
      if(ray.y!=0)
      {
            ymin=(q.y1-from.y)/ray.y;
            ymax=(q.y2-from.y)/ray.y;
            if(ymin>ymax) swap(ymax,ymin);
      }
      else if( (from.x<q.x1) || (from.x>q.x2)) return -1;

      if(ymax<xmax) xmax=ymax;
      if(ymin>xmin) xmin=ymin;
      if(xmax>dist) return -1;
      if(xmin>xmax) return -1; else return xmax;
}

int perpendicularPosition(const point3d_t &a,const point3d_t &b,const point3d_t &c,
                                                                              const square_t &q,PFLOAT Z)
{
      const point3d_t *min=&a,*mid=&b,*max=&c;
      if(mid->z<min->z) swap(min,mid);
      if(mid->z>max->z) swap(max,mid);
      if(mid->z<min->z) swap(min,mid);
      bool maxinside=q.isInside(*max);
      bool midinside=q.isInside(*mid);
      //bool mininside=q.isInside(*min);
      vector3d_t vmax=*max-*mid;
      vector3d_t vmin=*mid-*min;
      vector3d_t vmm=*max-*min;
      PFLOAT dmax=vmax.normLen();
      PFLOAT dmin=vmin.normLen();
      PFLOAT dmm=vmm.normLen();
      PFLOAT zmax=-numeric_limits<PFLOAT>::infinity();
      PFLOAT zmin=numeric_limits<PFLOAT>::infinity();
      if(maxinside) zmax=max->z;
      else
      {
            PFLOAT exit=exitDist(*mid,vmax,dmax,q);
            if(exit>=0) zmax=mid->z+vmax.z*exit;
            exit=exitDist(*min,vmm,dmm,q);
            if(exit>=0)
            {
                  exit=min->z+vmm.z*exit;
                  if(exit>zmax) zmax=exit;
            }
            if(midinside)
            {
                  exit=exitDist(*min,vmin,dmin,q);
                  if(exit>=0)
                  {
                        exit=min->z+vmin.z*exit;
                        if(exit>zmax) zmax=exit;
                  }
            }
      }
      if(q.isInside(*min)) zmin=min->z;
      else
      {
            vmax=-vmax;
            vmin=-vmin;
            vmm=-vmm;
            PFLOAT exit=exitDist(*mid,vmin,dmin,q);
            if(exit>=0) zmin=mid->z+vmin.z*exit;
            exit=exitDist(*max,vmm,dmm,q);
            if(exit>=0)
            {
                  exit=max->z+vmm.z*exit;
                  if(exit<zmin) zmin=exit;
            }
            if(midinside)
            {
                  exit=exitDist(*max,vmax,dmax,q);
                  if(exit>=0)
                  {
                        exit=max->z+vmax.z*exit;
                        if(exit<zmin) zmin=exit;
                  }
            }
      }
      /*
      if(dmm>1.0)
      {
      cout<<"quad "<<q.x1<<"-"<<q.x2<<" "<<q.y1<<"-"<<q.y2<<endl;
      cout<<"t "<<*min<<" "<<*mid<<" "<<*max<<endl;
      cout<<"min "<<zmin<<" max "<<zmax<<endl;
      getchar();
      }
      */
      if((zmin<Z) && (zmax<Z)) return trianglePosition_t::LOWER;
      if((zmin>Z) && (zmax>Z)) return trianglePosition_t::HIGHER;
      return trianglePosition_t::INTERSECT;
}

int cheapPosition(const triangle_t &tri,const bound_t &bound,PFLOAT Z,int axis)
{
      PFLOAT az=0,bz=0,cz=0;
      int decision=trianglePosition_t::NONE;
      point3d_t min,max;
      bound.get(min,max);
      bool inside=false;
      
      switch(axis)
      {
            case AXISX:
                  az=tri.a->x;bz=tri.b->x;cz=tri.c->x;
                  inside  =((tri.a->y>=min.y) && (tri.a->y<=max.y) && (tri.a->z>=min.z) && (tri.a->z<=max.z));
                  inside&=((tri.b->y>=min.y) && (tri.b->y<=max.y) && (tri.b->z>=min.z) && (tri.b->z<=max.z));
                  inside&=((tri.c->y>=min.y) && (tri.c->y<=max.y) && (tri.c->z>=min.z) && (tri.c->z<=max.z));
                  break;
            case AXISY:
                  az=tri.a->y;bz=tri.b->y;cz=tri.c->y;
                  inside  =((tri.a->x>=min.x) && (tri.a->x<=max.x) && (tri.a->z>=min.z) && (tri.a->z<=max.z));
                  inside&=((tri.b->x>=min.x) && (tri.b->x<=max.x) && (tri.b->z>=min.z) && (tri.b->z<=max.z));
                  inside&=((tri.c->x>=min.x) && (tri.c->x<=max.x) && (tri.c->z>=min.z) && (tri.c->z<=max.z));
                  break;
            case AXISZ:
                  az=tri.a->z;bz=tri.b->z;cz=tri.c->z;
                  inside  =((tri.a->x>=min.x) && (tri.a->x<=max.x) && (tri.a->y>=min.y) && (tri.a->y<=max.y));
                  inside&=((tri.b->x>=min.x) && (tri.b->x<=max.x) && (tri.b->y>=min.y) && (tri.b->y<=max.y));
                  inside&=((tri.c->x>=min.x) && (tri.c->x<=max.x) && (tri.c->y>=min.y) && (tri.c->y<=max.y));
                  break;
      }
      
      int intersection;
      if(inside)
            intersection=trianglePosition_t::INTERSECT;
      else
            intersection=trianglePosition_t::NONE; // possible intersection, we don't know

      if(az==Z) return intersection;
      if(az<Z) decision=trianglePosition_t::LOWER; else decision=trianglePosition_t::HIGHER;
      if(bz==Z) return intersection;
      if((bz>Z) && (decision==trianglePosition_t::LOWER)) return intersection;
      if((bz<Z) && (decision==trianglePosition_t::HIGHER)) return intersection;
      if(cz==Z) return intersection;
      if((cz>Z) && (decision==trianglePosition_t::LOWER)) return intersection;
      if((cz<Z) && (decision==trianglePosition_t::HIGHER)) return intersection;
      
      return decision;
}

int expensivePosition(const triangle_t &tri,const bound_t &bound,PFLOAT Z,int axis)
{
      const point3d_t &a3=*(tri.a),&b3=*(tri.b),&c3=*(tri.c);
      point3d_t bmin,bmax;
      const vector3d_t &n=tri.N();
      bound.get(bmin,bmax);

      point2d_t a,b,c;
      vector3d_t N;
      point3d_t planepoint;
      square_t q;
      
      switch(axis)
      {
            case AXISX:
                  q.x1=bmin.z;
                  q.x2=bmax.z;
                  q.y1=bmin.y;
                  q.y2=bmax.y;
                  N.x=n.z;N.y=n.y;N.z=n.x;
                  if(fabs(N.x)<0.0001) 
                        return perpendicularPosition(point3d_t(a3.z,a3.y,a3.x),
                                    point3d_t(b3.z,b3.y,b3.x),point3d_t(c3.z,c3.y,c3.x),q,Z);
                  a.x=a3.z;a.y=a3.y;
                  b.x=b3.z;b.y=b3.y;
                  c.x=c3.z;c.y=c3.y;
                  planepoint.x=a3.z;planepoint.y=a3.y;planepoint.z=a3.x;
                  break;
            case AXISY:
                  q.x1=bmin.x;
                  q.x2=bmax.x;
                  q.y1=bmin.z;
                  q.y2=bmax.z;
                  N.x=n.x;N.y=n.z;N.z=n.y;
                  if(fabs(N.y)<0.0001) 
                        return perpendicularPosition(point3d_t(a3.x,a3.z,a3.y),
                                    point3d_t(b3.x,b3.z,b3.y),point3d_t(c3.x,c3.z,c3.y),q,Z);
                  a.x=a3.x;a.y=a3.z;
                  b.x=b3.x;b.y=b3.z;
                  c.x=c3.x;c.y=c3.z;
                  planepoint.x=a3.x;planepoint.y=a3.z;planepoint.z=a3.y;
                  break;
            case AXISZ:
                  q.x1=bmin.x;
                  q.x2=bmax.x;
                  q.y1=bmin.y;
                  q.y2=bmax.y;
                  N.x=n.x;N.y=n.y;N.z=n.z;
                  if(fabs(N.z)<0.0001) 
                        return perpendicularPosition(a3,b3,c3,q,Z);
                  a.x=a3.x;a.y=a3.y;
                  b.x=b3.x;b.y=b3.y;
                  c.x=c3.x;c.y=c3.y;
                  planepoint.x=a3.x;planepoint.y=a3.y;planepoint.z=a3.z;
                  break;
      }

      checkPosition_f func(planepoint,N,Z);
      intersectApply(a,b,c,q,func);
      return func.decision;
}

__END_YAFRAY

Generated by  Doxygen 1.6.0   Back to index