/Users/mourrain/Devel/mmx/shape/src/ssi_qsegment.cpp

Go to the documentation of this file.

#include <shape/ssi_qsegment.hpp>
//#include <shape/ssi/chrono.h>

#define ParametricSurface shape::surface_parametric<double>

namespace mmx {
namespace  shape_ssi {

  std::ostream & operator<<( std::ostream& o, const lrnode& l )
  {
    o << l.lmin << " " << l.lmax << " " << l.rmin << " " << l.rmax << std::endl;
    return o;
  };

  lrnode * lrnode::alloc( unsigned size )
    {
      unsigned i, acc;
      i = size, acc =0;
      do { acc += i, i = i / 2 + i % 2; } while ( i > 1 );
      return new lrnode[acc+1];
    };
  
  lrnode * lrnode::make( lrnode * data, unsigned size )
  {
#define Mmin(a,b) std::min(a,b)
#define Mmax(a,b) std::max(a,b)
    
    unsigned bsize = size;
    bool     right = true;
    lrnode * built = data;
    lrnode * curr  = data + size;
    unsigned i,j;
    do
      {
        if ( right )
          { 
            for ( i = 0, j = 0; i < bsize-(bsize%2); i += 2, j++ ) 
              {
                curr[j].umin = built[i].umin;
                curr[j].umax = built[i+1].umax;
                curr[j].lmin = Mmin(built[i].lmin,
                                    built[i+1].lmin);
                curr[j].lmax = Mmax(built[i].lmax,
                                    built[i+1].lmax);
                curr[j].rmin = Mmin(built[i].rmin,
                                    built[i+1].rmin);
                curr[j].rmax = Mmax(built[i].rmax,
                                    built[i+1].rmax);
                curr[j].ls   = built + i;
                curr[j].rs   = built + i + 1;
              };
            
            if ( bsize%2 ) { curr[j++] = built[bsize-1];  };
            built = curr;
            curr += j;
          }
        else
          {
            unsigned next_size = bsize/2 + bsize % 2;
            built += bsize-1;
            curr  += next_size-1;
            for ( i = 0, j = 0; i < bsize-(bsize%2); i += 2, j++ ) 
              {
                curr[-j].umin = built[-i-1].umin;
                curr[-j].umax = built[-i].umax;
                curr[-j].lmin = Mmin(built[-i].lmin,
                                     built[-i-1].lmin);
                curr[-j].lmax = Mmax(built[-i].lmax,
                                           built[-i-1].lmax);
                curr[-j].rmin = Mmin(built[-i].rmin,
                                               built[-i-1].rmin);
                curr[-j].rmax = Mmax(built[-i].rmax,
                                     built[-i-1].rmax);
                curr[-j].ls   = built - i - 1;
                curr[-j].rs   = built - i;
              };
            
            if ( bsize%2 ) { curr[-j] = built[-bsize+1];  };
            built = curr - next_size+1;
            curr++;
          };
        right = !right;
        bsize = bsize/2 + bsize%2;
      }
    while( bsize > 1 );
    return --curr;
  };

  
  void shiftm( vector3 * v, unsigned sz, const aabb3 & box )
  {
    //    double m,M;
    for ( vector3 * src = v; src != v+sz; src ++ )
      for ( int i = 0; i < 3; (*src)[i] -= box[i].m, i ++ );
  };

    
  void scale( vector3 * p, double s ){ 
    for ( int i = 0; i < 4; i ++ )
      for ( int d = 0; d < 3; d ++ )
      p[i][d] *= s;
  };

  void qsegment::scale_conflict( vector3 * a, vector3 * b, qtree qa, qtree qb )
  {

    aabb3 h;
    hull(h,qa->box,qb->box);
    double s = 1.0/lmax(h);
    qa->fill(a,this);
    qb->fill(b,this);
    shiftm(a,4,h);
    shiftm(b,4,h);
    scale(a,s);
    scale(b,s);
  };
  
  void print( std::ostream& o, qtree q )
  {
    o << *q << std::endl;
    if ( q->r) print( o, q->r );
    if ( q->l ) print( o, q->l );
  };
  qnode * qsegment::make( rid_t id, coord_t vmin, coord_t vmax, lrnode * lr )
  { 
    qnode * curr, * l, * r;
    curr = l = r = 0;
    
    if (  vmax  <=  lr->lmin || vmin >=  lr->rmax ) return 0; 
    if (  vmin  <   lr->lmax || vmax >  lr->rmin  )
      {
        if ( lr->umax-lr->umin > vmax-vmin)
          {
            l = make( id, vmin, vmax, lr->ls );
            r = make( id, vmin, vmax, lr->rs );
          }
        else 
          {
            unsigned m = (vmin + vmax)/2;
            l = make( id, vmin, m, lr );
            r = make( id, m, vmax, lr );
          };
        
        if ( l && r ) 
          { 
            curr = new qnode();
            curr->umin = lr->umin;
            curr->umax = lr->umax;
            curr->vmin = vmin;
            curr->vmax = vmax;
            hull(curr->box,l->box,r->box);
            curr->l = l;
            curr->r = r;
            l->father = curr;
            r->father = curr;
            return curr;
          };
        
        if ( l ) { return l; };
        if ( r ) { return r; };
      }
    else 
      {
        curr = new qnode();
        curr->l    = 0;
        curr->r    = 0;
        curr->umin = lr->umin;
        curr->umax = lr->umax;
        curr->vmin = vmin;
        curr->vmax = vmax;
        curr->mbox(this);

        return curr;
      };
    return curr;
  };
  
  qsegment::qsegment( const ParametricSurface * s, unsigned w, unsigned h  ) : lsegment( s, w, h )
  {
    //    _qt = time();
    lrnode * lr_chunk   = lrnode::alloc(h);
    for ( unsigned i = 0; i < this->regions.size(); i++ )
      {
        lr_chunk -= this->regions[i].umin();
        for ( coord_t u = this->regions[i].umin(); u < this->regions[i].umax(); u++  )
          {
            lr_chunk[u].umin =  u;
            lr_chunk[u].umax =  u+1;
            lr_chunk[u].lmax = (lr_chunk[u].lmin = this->regions[i].min(u));
            lr_chunk[u].rmax = (lr_chunk[u].rmin = this->regions[i].max(u));
            lr_chunk[u].ls   =  lr_chunk[u].rs   = 0; 
          };

        lr_chunk += this->regions[i].umin();
        
        lrnode * lrt = lrnode::make(lr_chunk,
                                    this->regions[i].umax()-this->regions[i].umin());
        qtree tmp = make( i, lrt->lmin, lrt->rmax, lrt );
        tmp->father = 0;
        trees.push_back( tmp  );
        
      };
    delete[] lr_chunk;
    //    std::cout << "q = " << (time()-_qt) << std::endl;
  };
  
  qsegment::~qsegment() 
  { 
    for ( unsigned i = 0; i < size(); delete trees[i], i++ ); 
  };
  
}

} //namespace mmx    

# undef ParametricSurface