OSGErrorQuadTree.cpp

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002  *                           OpenSG NURBS Library                            *
00003  *                                                                           *
00004  *                                                                           *
00005  * Copyright (C) 2001-2006 by the University of Bonn, Computer Graphics Group*
00006  *                                                                           *
00007  *                         http://cg.cs.uni-bonn.de/                         *
00008  *                                                                           *
00009  * contact: edhellon@cs.uni-bonn.de, guthe@cs.uni-bonn.de, rk@cs.uni-bonn.de *
00010  *                                                                           *
00011 \*---------------------------------------------------------------------------*/
00012 /*---------------------------------------------------------------------------*\
00013  *                                License                                    *
00014  *                                                                           *
00015  * This library is free software; you can redistribute it and/or modify it   *
00016  * under the terms of the GNU Library General Public License as published    *
00017  * by the Free Software Foundation, version 2.                               *
00018  *                                                                           *
00019  * This library is distributed in the hope that it will be useful, but       *
00020  * WITHOUT ANY WARRANTY; without even the implied warranty of                *
00021  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU         *
00022  * Library General Public License for more details.                          *
00023  *                                                                           *
00024  * You should have received a copy of the GNU Library General Public         *
00025  * License along with this library; if not, write to the Free Software       *
00026  * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.                 *
00027  *                                                                           *
00028 \*---------------------------------------------------------------------------*/
00029 /*---------------------------------------------------------------------------*\
00030  *                                Changes                                    *
00031  *                                                                           *
00032  *                                                                           *
00033  *                                                                           *
00034  *                                                                           *
00035  *                                                                           *
00036  *                                                                           *
00037 \*---------------------------------------------------------------------------*/
00038 #include "OSGErrorQuadTree.h"
00039
00040 #include "OSGDCTPMesh.h"
00041 #include "OSGBezierTensorSurface.h"
00042 #include "OSGBSplineTensorSurface.h"
00043
00044
00045
00046 #include <stdlib.h>
00047 #include <float.h>
00048
00049 OSG_USING_NAMESPACE
00050
00051
00052 //#define min( a, b ) ( ( ( a ) < ( b ) ) ? ( a ) : ( b ) )
00053 //#define max( a, b ) ( ( ( a ) > ( b ) ) ? ( a ) : ( b ) )
00054
00055 //#define OSG_EUCLID_ERRORS
00056
00057 #ifdef _DEBUG
00058  #ifdef WIN32
00059   #undef THIS_FILE
00060 static char THIS_FILE[] = __FILE__;
00061  #endif
00062 #endif
00063 
00064 bool CErrorQuadTree::m_sbNormalApproximation = false;
00065
00066 CErrorQuadTree::CErrorQuadTree()
00067 {
00068     m_fMaxError = -1.0;
00069 #ifdef OSG_USE_NURBS_PATCH
00070     m_ptRoot = NULL;
00071 #else
00072     m_vvptRoot.clear();
00073     m_ptBPRoot = NULL;
00074 #endif
00075     m_fErrorCutoff = DCTP_EPS;
00076 }
00077
00078
00079 CErrorQuadTree::~CErrorQuadTree()
00080 {
00081 #ifdef OSG_USE_NURBS_PATCH
00082     if(m_ptRoot)
00083     {
00084         DeleteNode(m_ptRoot);
00085  #ifdef OSG_ONE_CHILD_PTR
00086         delete m_ptRoot;
00087         m_ptRoot = NULL;
00088  #endif
00089     }
00090 #else
00091     const unsigned int cui_size_u = m_vvptRoot.size();
00092
00093     if(cui_size_u)
00094     {
00095         const unsigned int cui_size_v = m_vvptRoot[0].size();
00096         unsigned int       ui_u;
00097         unsigned int       ui_v;
00098
00099         for(ui_u = 0; ui_u < cui_size_u; ++ui_u)
00100         {
00101             for(ui_v = 0; ui_v < cui_size_v; ++ui_v)
00102             {
00103                 DeleteNode(m_vvptRoot[ui_u][ui_v]);
00104  #ifdef OSG_ONE_CHILD_PTR
00105                 delete m_vvptRoot[ui_u][ui_v];
00106                 m_vvptRoot[ui_u][ui_v] = NULL;
00107  #endif
00108             }
00109         }
00110
00111         m_vvptRoot.clear();
00112     }
00113
00114     if(m_ptBPRoot)
00115     {
00116         DeleteBPNode(m_ptBPRoot);
00117     }
00118 #endif
00119 }
00120
00121
00122 //#ifndef OSG_ARBITRARY_SPLIT
00123 void CErrorQuadTree::BuildMesh(DCTPMesh *pclMesh,
00124 #ifdef OSG_USE_NURBS_PATCH
00125                                BSplineTensorSurface *pclPatch,
00126  #ifdef OSG_ARBITRARY_SPLIT
00127                                const Vec2d cclMinParam, const Vec2d cclMaxParam,
00128  #endif
00129 #else
00130                                std::vector<std::vector<BezierTensorSurface> > *pvvclPatches,
00131                                const std::vector<double> *cpvdIntervalsU,
00132                                const std::vector<double> *cpvdIntervalsV,
00133 #endif
00134                                float fError, float &rfMinError, float &rfMaxError)
00135 {
00136     pclMesh->reinit();
00137 #ifdef OSG_USE_NURBS_PATCH
00138  #ifdef OSG_ARBITRARY_SPLIT
00139     SetInitialCells(pclMesh, fError, pclPatch, cclMinParam, cclMaxParam);
00140
00141     if( (cclMaxParam[0] - cclMinParam[0] < DCTP_EPS) ||
00142         (cclMaxParam[1] - cclMinParam[1] < DCTP_EPS) )
00143     {
00144         return;
00145     }
00146  #else
00147     SetInitialCells(pclMesh, fError, pclPatch);
00148  #endif
00149 #else
00150     SetInitialCells(pclMesh, fError, pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
00151 #endif
00152 
00153     unsigned int   ui_face;
00154     unsigned int   ui_face_cnt = pclMesh->faces.size();
00155     DCTPFace      *pcl_face;
00156     SFaceTreeCell *pt_finfo;
00157     float          f_act_error;
00158     float          f_prev_error;
00159
00160     rfMinError = -1.0;
00161     rfMaxError = FLT_MAX;
00162
00163     if(fError < m_fErrorCutoff)
00164     {
00165         fError = m_fErrorCutoff;
00166     }
00167
00168     if( (fError >= m_fMaxError) && (m_fMaxError >= 0.0) )
00169     {
00170         // ok, we already have a subdivision for that error...
00171         for(ui_face = 0; ui_face < ui_face_cnt; ++ui_face)
00172         {
00173             pcl_face = pclMesh->faces[ui_face];
00174             if(pcl_face->faceinfo)
00175             {
00176                 f_prev_error = FLT_MAX;
00177 #ifdef OSG_USE_NURBS_PATCH
00178                 f_act_error = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fError;
00179                 while(f_act_error > fError)
00180                 {
00181 //                  std::cerr << ",";
00182  #ifdef OSG_USE_KD_TREE
00183   #ifdef OSG_ARBITRARY_SPLIT
00184                     if( (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fSplitValue < 0.0)
00185   #else
00186                     if( (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->bSplitHoriz)
00187   #endif
00188                     {
00189   #ifdef OSG_ARBITRARY_SPLIT
00190                         pclMesh->SubdivideQuadEW(pcl_face, -(static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fSplitValue);
00191   #else
00192                         pclMesh->SubdivideQuadEW(pcl_face);
00193   #endif
00194                     }
00195                     else
00196                     {
00197   #ifdef OSG_ARBITRARY_SPLIT
00198                         pclMesh->SubdivideQuadNS(pcl_face, (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fSplitValue);
00199   #else
00200                         pclMesh->SubdivideQuadNS(pcl_face);
00201   #endif
00202                     }
00203                     SubdivideNode(pclMesh, pcl_face);
00204                     ++ui_face_cnt;      // one new face
00205  #else
00206                     pclMesh->SubdivideQuad(pcl_face);
00207                     SubdivideNode(pclMesh, pcl_face);
00208                     ui_face_cnt += 3;   // three new faces
00209  #endif
00210                     f_prev_error = f_act_error;
00211                     f_act_error  = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fError;
00212                 }
00213 #else
00214                 if( (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell)
00215                 {
00216                     f_act_error = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fError;
00217                     while(f_act_error > fError)
00218                     {
00219                         pclMesh->SubdivideQuad(pcl_face);
00220                         SubdivideNode(pclMesh, pcl_face);
00221                         ui_face_cnt += 3;   // three new faces
00222                         f_prev_error = f_act_error;
00223                         f_act_error  = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fError;
00224                     }
00225                 }
00226                 else
00227                 {
00228                     f_act_error  = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptBPCell->fError;
00229                     f_prev_error = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptBPCell->fPrevError;
00230 //                  std::cerr << pcl_face->orig_quad[ 3 ]->coords << pcl_face->orig_quad[ 1 ]->coords << std::endl;
00231 //                  std::cerr << f_act_error << ", " << f_prev_error << std::endl;
00232                 }
00233 #endif
00234                 if(f_act_error > rfMinError)
00235                 {
00236                     rfMinError = f_act_error;
00237                 }
00238                 if(f_prev_error < rfMaxError)
00239                 {
00240                     rfMaxError = f_prev_error;
00241                 }
00242             }
00243         }
00244
00245         // delete faceinfo
00246         for(ui_face = 0; ui_face < ui_face_cnt; ++ui_face)
00247         {
00248             pcl_face = pclMesh->faces[ui_face];
00249             pt_finfo = static_cast<SFaceTreeCell*>(pcl_face->faceinfo);
00250             if(pt_finfo)
00251             {
00252                 delete pt_finfo;
00253                 pcl_face->faceinfo = NULL;
00254             }
00255 /*          else
00256             {
00257                 pclMesh->removeFace( ui_face );
00258                 --ui_face;
00259                 --ui_face_cnt;
00260             }*/
00261         }
00262
00263 //      std::cerr << "reconstructed " << ui_face_cnt << " faces" << std::endl;
00264 /*      if( pclMesh->faces.size( ) >= 1000 )
00265         {
00266             if( fError > m_fErrorCutoff )
00267             {
00268                 std::cerr << "reduce error cutoff: " << m_fErrorCutoff << " -> " << fError << std::endl;
00269                 m_fErrorCutoff = fError;
00270             }
00271         }
00272         if( fError <= m_fErrorCutoff )
00273         {
00274             rfMinError = 0.0;
00275         }*/
00276         return;
00277     }
00278
00279 #ifdef OSG_ARBITRARY_SPLIT
00280  #ifndef OSG_USE_NURBS_PATCH
00281   #error
00282  #endif
00283     // create first surface
00284     std::vector<BSplineTensorSurface> vcl_surf;
00285     double                            d_split;
00286     Vec2d                             cl_min;
00287     Vec2d                             cl_add;
00288
00289     pt_finfo                       = (static_cast<SFaceTreeCell*>(pclMesh->faces[0]->faceinfo) );
00290     pt_finfo->bOwnSurface          = true;
00291     pt_finfo->pclBSplineSurface    = new BSplineTensorSurface;
00292     *(pt_finfo->pclBSplineSurface) = *pclPatch;
00293     pclPatch->getParameterInterval_U(cl_min[0], cl_add[0]);
00294     pclPatch->getParameterInterval_V(cl_min[1], cl_add[1]);
00295     cl_add -= cl_min;
00296
00297     d_split = (cclMinParam[0] - cl_min[0]) / cl_add[0];
00298     if(d_split > 0.999)
00299         d_split = 0.999;
00300     if(d_split > DCTP_EPS)
00301     {
00302         pt_finfo->pclBSplineSurface->subDivisionU(vcl_surf, d_split);
00303         *(pt_finfo->pclBSplineSurface) = vcl_surf[1];
00304     }
00305
00306     d_split = (cclMaxParam[0] - cclMinParam[0]) / (cl_add[0] - cclMinParam[0]);
00307     if(d_split < 0.001)
00308         d_split = 0.001;
00309     if(1.0 - d_split > DCTP_EPS)
00310     {
00311         pt_finfo->pclBSplineSurface->subDivisionU(vcl_surf, d_split);
00312         *(pt_finfo->pclBSplineSurface) = vcl_surf[0];
00313     }
00314
00315     d_split = (cclMinParam[1] - cl_min[1]) / cl_add[1];
00316     if(d_split > 0.999)
00317         d_split = 0.999;
00318     if(d_split > DCTP_EPS)
00319     {
00320         pt_finfo->pclBSplineSurface->subDivisionV(vcl_surf, d_split);
00321         *(pt_finfo->pclBSplineSurface) = vcl_surf[1];
00322     }
00323
00324     d_split = (cclMaxParam[1] - cclMinParam[1]) / (cl_add[1] - cclMinParam[1]);
00325     if(d_split < 0.001)
00326         d_split = 0.001;
00327     if(1.0 - d_split > DCTP_EPS)
00328     {
00329         pt_finfo->pclBSplineSurface->subDivisionV(vcl_surf, d_split);
00330         *(pt_finfo->pclBSplineSurface) = vcl_surf[0];
00331     }
00332
00333     vcl_surf.clear();
00334 #endif
00335 
00336     // we need to build a new quadtree
00337 /*  if( ( fError > m_fMaxError * 0.7071068115 ) && ( m_fMaxError >= 0.0 ) )
00338     {
00339         // new error must be 1 / sqrt( 2 ) or less than old...
00340         fError = m_fMaxError * 0.7071068115;
00341     }*/
00342
00343     // init errors if computed for the first time
00344 //  if( m_fMaxError < 0.0 )
00345     {
00346         for(ui_face = 0; ui_face < ui_face_cnt; ++ui_face)
00347         {
00348             pcl_face = pclMesh->faces[ui_face];
00349             if(pcl_face->faceinfo)
00350             {
00351                 if( (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell)
00352                 {
00353                     ComputeError(pcl_face);
00354                 }
00355             }
00356         }
00357     }
00358
00359     // subdivide and get max error
00360     m_fMaxError = -1.0;
00361
00362     for(ui_face = 0; ui_face < ui_face_cnt; ++ui_face)
00363     {
00364         pcl_face = pclMesh->faces[ui_face];
00365         if(pcl_face->faceinfo)
00366         {
00367             f_prev_error = FLT_MAX;
00368 #ifdef OSG_USE_NURBS_PATCH
00369             f_act_error = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fError;
00370             while(f_act_error > fError)
00371             {
00372                 const SErrorTreeCell* cpt_errcell = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell;
00373                 bool                  b_split_ok  = true;
00374
00375 //              std::cerr << ".";
00376  #ifdef OSG_USE_KD_TREE
00377   #ifdef OSG_ARBITRARY_SPLIT
00378                 if(cpt_errcell->fSplitValue < 0.0)
00379   #else
00380                 if(cpt_errcell->bSplitHoriz)
00381   #endif
00382                 {
00383   #ifdef OSG_ARBITRARY_SPLIT
00384                     const double cd_ratio = osgMin(-cpt_errcell->fSplitValue, 1.0f + cpt_errcell->fSplitValue);
00385    #ifdef OSG_UNION_TRI_QUAD
00386                     if( (pcl_face->orig_face[1]->coords[0] - pcl_face->orig_face[0]->coords[0]) * cd_ratio > DCTP_EPS * 100.0)
00387    #else
00388                     if( (pcl_face->orig_quad[1]->coords[0] - pcl_face->orig_quad[0]->coords[0]) * cd_ratio > DCTP_EPS * 100.0)
00389    #endif
00390   #else
00391    #ifdef OSG_UNION_TRI_QUAD
00392                     if( (pcl_face->orig_face[1]->coords[0] - pcl_face->orig_face[0]->coords[0]) * 0.5 > DCTP_EPS * 100.0)
00393    #else
00394                     if( (pcl_face->orig_quad[1]->coords[0] - pcl_face->orig_quad[0]->coords[0]) * 0.5 > DCTP_EPS * 100.0)
00395    #endif
00396   #endif
00397                     {
00398   #ifdef OSG_ARBITRARY_SPLIT
00399                         pclMesh->SubdivideQuadEW(pcl_face, -cpt_errcell->fSplitValue);
00400   #else
00401                         pclMesh->SubdivideQuadEW(pcl_face);
00402   #endif
00403                     }
00404                     else
00405                     {
00406                         b_split_ok = false;
00407                     }
00408                 }
00409                 else
00410                 {
00411   #ifdef OSG_ARBITRARY_SPLIT
00412                     const double cd_ratio = osgMin(cpt_errcell->fSplitValue, 1.0f - cpt_errcell->fSplitValue);
00413    #ifdef OSG_UNION_TRI_QUAD
00414                     if( (pcl_face->orig_face[0]->coords[1] - pcl_face->orig_face[3]->coords[1]) * cd_ratio > DCTP_EPS * 100.0)
00415    #else
00416                     if( (pcl_face->orig_quad[0]->coords[1] - pcl_face->orig_quad[3]->coords[1]) * cd_ratio > DCTP_EPS * 100.0)
00417    #endif
00418   #else
00419    #ifdef OSG_UNION_TRI_QUAD
00420                     if( (pcl_face->orig_face[0]->coords[1] - pcl_face->orig_face[3]->coords[1]) * 0.5 > DCTP_EPS * 100.0)
00421    #else
00422                     if( (pcl_face->orig_quad[0]->coords[1] - pcl_face->orig_quad[3]->coords[1]) * 0.5 > DCTP_EPS * 100.0)
00423    #endif
00424   #endif
00425                     {
00426   #ifdef OSG_ARBITRARY_SPLIT
00427                         pclMesh->SubdivideQuadNS(pcl_face, cpt_errcell->fSplitValue);
00428   #else
00429                         pclMesh->SubdivideQuadNS(pcl_face);
00430   #endif
00431                     }
00432                     else
00433                     {
00434                         b_split_ok = false;
00435                     }
00436                 }
00437                 if(b_split_ok)
00438                 {
00439                     SubdivideBuild(pclMesh, pcl_face);
00440                     ++ui_face_cnt;      // one new face
00441                     f_prev_error = f_act_error;
00442                     f_act_error  = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fError;
00443                 }
00444                 else
00445                 {
00446                     f_prev_error = f_act_error;
00447                     f_act_error  = 0.0;
00448                 }
00449  #else
00450                 pclMesh->SubdivideQuad(pcl_face);
00451                 SubdivideBuild(pclMesh, pcl_face);
00452                 ui_face_cnt += 3;   // three new faces
00453                 f_prev_error = f_act_error;
00454                 f_act_error  = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fError;
00455  #endif
00456             }
00457 #else
00458             if( (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell)
00459             {
00460                 f_act_error = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fError;
00461                 while(f_act_error > fError)
00462                 {
00463                     pclMesh->SubdivideQuad(pcl_face);
00464                     SubdivideBuild(pclMesh, pcl_face);
00465                     ui_face_cnt += 3;   // three new faces
00466                     f_prev_error = f_act_error;
00467                     f_act_error  = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptErrorCell->fError;
00468                 }
00469             }
00470             else
00471             {
00472                 f_act_error  = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptBPCell->fError;
00473                 f_prev_error = (static_cast<SFaceTreeCell*>(pcl_face->faceinfo) )->ptBPCell->fPrevError;
00474             }
00475 #endif
00476             if(f_act_error > m_fMaxError)
00477             {
00478                 m_fMaxError = f_act_error;
00479                 rfMinError  = f_act_error;
00480             }
00481             if(f_prev_error < rfMaxError)
00482             {
00483                 rfMaxError = f_prev_error;
00484             }
00485         }
00486     }
00487
00488     // delete faceinfo and constructed bezier patches
00489     for(ui_face = 0; ui_face < ui_face_cnt; ++ui_face)
00490     {
00491         pcl_face = pclMesh->faces[ui_face];
00492         pt_finfo = static_cast<SFaceTreeCell*>(pcl_face->faceinfo);
00493         if(pt_finfo)
00494         {
00495             if(pt_finfo->bOwnSurface)
00496             {
00497 #ifdef OSG_USE_NURBS_PATCH
00498                 delete pt_finfo->pclBSplineSurface;
00499 #else
00500                 delete pt_finfo->pclBezierSurface;
00501 #endif
00502             }
00503             delete pt_finfo;
00504             pcl_face->faceinfo = NULL;
00505         }
00506 /*      else
00507         {
00508             pclMesh->removeFace( ui_face );
00509             --ui_face;
00510             --ui_face_cnt;
00511         }*/
00512     }
00513
00514 //  std::cerr << "build " << ui_face_cnt << " faces" << std::endl;
00515
00516 /*  if( pclMesh->faces.size( ) >= 1000 )
00517     {
00518         if( fError > m_fErrorCutoff )
00519         {
00520             std::cerr << "reduce error cutoff: " << m_fErrorCutoff << " -> " << fError << std::endl;
00521             m_fErrorCutoff = fError;
00522         }
00523     }
00524     if( fError <= m_fErrorCutoff )
00525     {
00526         rfMinError = 0.0;
00527     }*/
00528
00529 //  std::cerr << "done" << std::endl;
00530
00531     return;
00532 }
00533 //#endif
00534
00535
00536 /*#ifdef OSG_ARBITRARY_SPLIT
00537 void CErrorQuadTree::CalculatePoints( std::vector< Vec2d > *pvclInsert, std::vector< Vec2d > *pvclDelete,
00538                                       BSplineTensorSurface *pclPatch,
00539                                       float fError, float &rfMinError, float &rfMaxError )
00540 {
00541     std::vector< SFaceTreeCell >            vt_build_cells;
00542     unsigned int                    ui_cell;
00543     unsigned int                    ui_cell_cnt = 1;
00544     float                           f_act_error;
00545     float                           f_prev_error;
00546     float                           f_old_error;
00547     bool                            b_insert;
00548     std::vector< BSplineTensorSurface > vcl_surfaces;
00549     Vec2dset                        scl_check;
00550 
00551     if( fError < m_fErrorCutoff )
00552     {
00553         fError = m_fErrorCutoff;
00554     }
00555 
00556     pvclInsert->clear( );
00557     pvclDelete->clear( );
00558 
00559     if( m_fMaxError >= 0.0 )
00560     {
00561         f_old_error = rfMinError;
00562         if( fError <= rfMinError )
00563         {
00564             b_insert = true;
00565         }
00566         else if( fError >= rfMaxError )
00567         {
00568             b_insert = false;
00569         }
00570         else
00571         {
00572             return;
00573         }
00574     }
00575     else
00576     {
00577         f_old_error = FLT_MAX;
00578         b_insert = true;
00579     }
00580 
00581     if( m_ptRoot == NULL )
00582     {
00583         m_ptRoot = new SErrorTreeCell;
00584         m_ptRoot->ptChildren = NULL;
00585         m_ptRoot->fError = -1.0;
00586     }
00587 
00588     vt_build_cells.resize( 1 );
00589     vt_build_cells[ 0 ].bOwnSurface = false;
00590     vt_build_cells[ 0 ].pclBSplineSurface = pclPatch;
00591     vt_build_cells[ 0 ].ptErrorCell = m_ptRoot;
00592     pclPatch->getParameterInterval_U( vt_build_cells[ 0 ].clMin[0], vt_build_cells[ 0 ].clMax[0] );
00593     pclPatch->getParameterInterval_V( vt_build_cells[ 0 ].clMin[1], vt_build_cells[ 0 ].clMax[1] );
00594     ComputeError( &( vt_build_cells[ 0 ] ) );
00595 
00596     rfMinError = -1.0;
00597     rfMaxError = FLT_MAX;
00598 
00599     if( ( fError >= m_fMaxError ) && ( m_fMaxError >= 0.0 ) )
00600     {
00601         // ok, we already have a subdivision for that error...
00602         for( ui_cell = 0; ui_cell < ui_cell_cnt; ++ui_cell )
00603         {
00604             f_prev_error = FLT_MAX;
00605             f_act_error = vt_build_cells[ ui_cell ].ptErrorCell->fError;
00606             while( f_act_error > fError )
00607             {
00608                 vt_build_cells.resize( ui_cell_cnt + 1 );
00609                 vt_build_cells[ ui_cell_cnt ].bOwnSurface = false;
00610                 vt_build_cells[ ui_cell_cnt ].pclBSplineSurface = NULL;
00611                 vt_build_cells[ ui_cell_cnt ].clMin = vt_build_cells[ ui_cell ].clMin;
00612                 vt_build_cells[ ui_cell_cnt ].clMax = vt_build_cells[ ui_cell ].clMax;
00613                 vt_build_cells[ ui_cell ].bOwnSurface = false;
00614                 vt_build_cells[ ui_cell ].pclBSplineSurface = NULL;
00615                 if( vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue < 0.0 )
00616                 {
00617                     vt_build_cells[ ui_cell ].clMax[0] =
00618                         vt_build_cells[ ui_cell ].clMin[0]
00619                         - ( vt_build_cells[ ui_cell ].clMax[0] - vt_build_cells[ ui_cell ].clMin[0] )
00620                         * vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue;
00621                     vt_build_cells[ ui_cell_cnt ].clMin[0] = vt_build_cells[ ui_cell ].clMax[0];
00622                 }
00623                 else
00624                 {
00625                     vt_build_cells[ ui_cell ].clMax[1] =
00626                         vt_build_cells[ ui_cell ].clMin[1]
00627                         + ( vt_build_cells[ ui_cell ].clMax[1] - vt_build_cells[ ui_cell ].clMin[1] )
00628                         * vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue;
00629                     vt_build_cells[ ui_cell_cnt ].clMin[1] = vt_build_cells[ ui_cell ].clMax[1];
00630                 }
00631                 // set new error cells
00632                 vt_build_cells[ ui_cell_cnt ].ptErrorCell = &( vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 1 ] );
00633                 vt_build_cells[ ui_cell ].ptErrorCell = &( vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 0 ] );
00634                 if( ( b_insert ) && ( f_act_error >= f_old_error ) )
00635                 {
00636                     // insert new points
00637                     if( scl_check.insert( vt_build_cells[ ui_cell ].clMax ).second )
00638                     {
00639                         pvclInsert->push_back( vt_build_cells[ ui_cell ].clMax );
00640                     }
00641                     if( scl_check.insert( vt_build_cells[ ui_cell_cnt ].clMin ).second )
00642                     {
00643                         pvclInsert->push_back( vt_build_cells[ ui_cell_cnt ].clMin );
00644                     }
00645                 }
00646                 else
00647                 {
00648                     scl_check.insert( vt_build_cells[ ui_cell ].clMax );
00649                     scl_check.insert( vt_build_cells[ ui_cell_cnt ].clMin );
00650                 }
00651                 ++ui_cell_cnt;      // one new face
00652                 f_prev_error = f_act_error;
00653                 f_act_error = vt_build_cells[ ui_cell ].ptErrorCell->fError;
00654             }
00655             if( f_act_error > rfMinError )
00656             {
00657                 rfMinError = f_act_error;
00658             }
00659             if( f_prev_error < rfMaxError )
00660             {
00661                 rfMaxError = f_prev_error;
00662             }
00663             if( !b_insert )
00664             {
00665                 // delete old points
00666                 while( f_act_error > f_old_error )
00667                 {
00668                     vt_build_cells.resize( ui_cell_cnt + 1 );
00669                     vt_build_cells[ ui_cell_cnt ].bOwnSurface = false;
00670                     vt_build_cells[ ui_cell_cnt ].pclBSplineSurface = NULL;
00671                     vt_build_cells[ ui_cell_cnt ].clMin = vt_build_cells[ ui_cell ].clMin;
00672                     vt_build_cells[ ui_cell_cnt ].clMax = vt_build_cells[ ui_cell ].clMax;
00673                     vt_build_cells[ ui_cell ].bOwnSurface = false;
00674                     vt_build_cells[ ui_cell ].pclBSplineSurface = NULL;
00675                     if( vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue < 0.0 )
00676                     {
00677                         vt_build_cells[ ui_cell ].clMax[0] =
00678                             vt_build_cells[ ui_cell ].clMin[0]
00679                             - ( vt_build_cells[ ui_cell ].clMax[0] - vt_build_cells[ ui_cell ].clMin[0] )
00680                             * vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue;
00681                         vt_build_cells[ ui_cell_cnt ].clMin[0] = vt_build_cells[ ui_cell ].clMax[0];
00682                     }
00683                     else
00684                     {
00685                         vt_build_cells[ ui_cell ].clMax[1] =
00686                             vt_build_cells[ ui_cell ].clMin[1]
00687                             + ( vt_build_cells[ ui_cell ].clMax[1] - vt_build_cells[ ui_cell ].clMin[1] )
00688                             * vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue;
00689                         vt_build_cells[ ui_cell_cnt ].clMin[1] = vt_build_cells[ ui_cell ].clMax[1];
00690                     }
00691                     // set new error cells
00692                     vt_build_cells[ ui_cell_cnt ].ptErrorCell = &( vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 1 ] );
00693                     vt_build_cells[ ui_cell ].ptErrorCell = &( vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 0 ] );
00694                     if( scl_check.insert( vt_build_cells[ ui_cell ].clMax ).second )
00695                     {
00696                         pvclDelete->push_back( vt_build_cells[ ui_cell ].clMax );
00697                     }
00698                     if( scl_check.insert( vt_build_cells[ ui_cell_cnt ].clMin ).second )
00699                     {
00700                         pvclDelete->push_back( vt_build_cells[ ui_cell_cnt ].clMin );
00701                     }
00702                     ++ui_cell_cnt;      // one new face
00703                     f_act_error = vt_build_cells[ ui_cell ].ptErrorCell->fError;
00704                 }
00705             }
00706         }
00707         vt_build_cells.clear( );
00708         return;
00709     }
00710 
00711     // we need to build a new Kd-tree
00712     for( ui_cell = 0; ui_cell < ui_cell_cnt; ++ui_cell )
00713     {
00714         f_prev_error = FLT_MAX;
00715         f_act_error = vt_build_cells[ ui_cell ].ptErrorCell->fError;
00716         while( f_act_error > fError )
00717         {
00718             if( vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue < 0.0 )
00719             {
00720                 vt_build_cells[ ui_cell ].pclBSplineSurface->subDivisionU(
00721                     vcl_surfaces, -vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue );
00722             }
00723             else
00724             {
00725                 vt_build_cells[ ui_cell ].pclBSplineSurface->subDivisionV(
00726                     vcl_surfaces, vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue );
00727             }
00728             vt_build_cells.resize( ui_cell_cnt + 1 );
00729             vt_build_cells[ ui_cell_cnt ].bOwnSurface = true;
00730             vt_build_cells[ ui_cell_cnt ].pclBSplineSurface = new BSplineTensorSurface;
00731             *( vt_build_cells[ ui_cell_cnt ].pclBSplineSurface ) = vcl_surfaces[ 1 ];
00732             vt_build_cells[ ui_cell_cnt ].clMin = vt_build_cells[ ui_cell ].clMin;
00733             vt_build_cells[ ui_cell_cnt ].clMax = vt_build_cells[ ui_cell ].clMax;
00734             if( vt_build_cells[ ui_cell ].bOwnSurface )
00735             {
00736                 delete vt_build_cells[ ui_cell ].pclBSplineSurface;
00737             }
00738             vt_build_cells[ ui_cell ].bOwnSurface = true;
00739             vt_build_cells[ ui_cell ].pclBSplineSurface = new BSplineTensorSurface;
00740             *( vt_build_cells[ ui_cell ].pclBSplineSurface ) = vcl_surfaces[ 0 ];
00741             std::cerr << vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue << std::endl;
00742             std::cerr << vt_build_cells[ ui_cell ].clMin << vt_build_cells[ ui_cell ].clMax << std::endl;
00743             if( vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue < 0.0 )
00744             {
00745                 vt_build_cells[ ui_cell ].clMax[0] =
00746                     vt_build_cells[ ui_cell ].clMin[0]
00747                     - ( vt_build_cells[ ui_cell ].clMax[0] - vt_build_cells[ ui_cell ].clMin[0] )
00748                     * vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue;
00749                 vt_build_cells[ ui_cell_cnt ].clMin[0] = vt_build_cells[ ui_cell ].clMax[0];
00750             }
00751             else
00752             {
00753                 vt_build_cells[ ui_cell ].clMax[1] =
00754                     vt_build_cells[ ui_cell ].clMin[1]
00755                     + ( vt_build_cells[ ui_cell ].clMax[1] - vt_build_cells[ ui_cell ].clMin[1] )
00756                     * vt_build_cells[ ui_cell ].ptErrorCell->fSplitValue;
00757                 vt_build_cells[ ui_cell_cnt ].clMin[1] = vt_build_cells[ ui_cell ].clMax[1];
00758             }
00759             std::cerr << vt_build_cells[ ui_cell ].clMin << vt_build_cells[ ui_cell ].clMax << std::endl;
00760             std::cerr << vt_build_cells[ ui_cell_cnt ].clMin << vt_build_cells[ ui_cell_cnt ].clMax << std::endl << std::endl;
00761             // set new error cells
00762             if( vt_build_cells[ ui_cell ].ptErrorCell->ptChildren == NULL )
00763             {
00764                 vt_build_cells[ ui_cell ].ptErrorCell->ptChildren = new SErrorTreeCell[ 2 ];
00765                 vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 0 ].ptChildren = NULL;
00766                 vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 0 ].fError = -1.0;
00767                 vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 1 ].ptChildren = NULL;
00768                 vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 1 ].fError = -1.0;
00769             }
00770             vt_build_cells[ ui_cell_cnt ].ptErrorCell = &( vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 1 ] );
00771             vt_build_cells[ ui_cell ].ptErrorCell = &( vt_build_cells[ ui_cell ].ptErrorCell->ptChildren[ 0 ] );
00772             // compute error
00773             ComputeError( &( vt_build_cells[ ui_cell ] ) );
00774             ComputeError( &( vt_build_cells[ ui_cell_cnt ] ) );
00775             // insert new points
00776             if( scl_check.insert( vt_build_cells[ ui_cell ].clMax ).second )
00777             {
00778                 pvclInsert->push_back( vt_build_cells[ ui_cell ].clMax );
00779             }
00780             if( scl_check.insert( vt_build_cells[ ui_cell_cnt ].clMin ).second )
00781             {
00782                 pvclInsert->push_back( vt_build_cells[ ui_cell_cnt ].clMin );
00783             }
00784             ++ui_cell_cnt;      // one new face
00785             f_prev_error = f_act_error;
00786             f_act_error = vt_build_cells[ ui_cell ].ptErrorCell->fError;
00787         }
00788         if( f_act_error > rfMinError )
00789         {
00790             rfMinError = f_act_error;
00791         }
00792         if( f_prev_error < rfMaxError )
00793         {
00794             rfMaxError = f_prev_error;
00795         }
00796     }
00797 
00798     // delete build cells and constructed bspline patches
00799     for( ui_cell = 0; ui_cell < ui_cell_cnt; ++ui_cell )
00800     {
00801         if( vt_build_cells[ ui_cell ].bOwnSurface )
00802         {
00803             vt_build_cells[ ui_cell ].pclBSplineSurface;
00804         }
00805     }
00806     vt_build_cells.clear( );
00807 }
00808 #endif*/
00809
00810
00811 #ifdef OSG_ONE_CHILD_PTR
00812 void CErrorQuadTree::DeleteNode(SErrorTreeCell *pclNode)
00813 {
00814     if(pclNode->ptChildren)
00815     {
00816         DeleteNode(&(pclNode->ptChildren[0]) );
00817         DeleteNode(&(pclNode->ptChildren[1]) );
00818  #ifndef OSG_USE_KD_TREE
00819         DeleteNode(&(pclNode->ptChildren[2]) );
00820         DeleteNode(&(pclNode->ptChildren[3]) );
00821  #endif
00822         delete[]  pclNode->ptChildren;
00823         pclNode->ptChildren = NULL;
00824     }
00825 }
00826 #else
00827 void CErrorQuadTree::DeleteNode(SErrorTreeCell *&rpclNode)
00828 {
00829     if(rpclNode->aptChildren[0])
00830     {
00831         DeleteNode(rpclNode->aptChildren[0]);
00832         DeleteNode(rpclNode->aptChildren[1]);
00833  #ifndef OSG_USE_KD_TREE
00834         DeleteNode(rpclNode->aptChildren[2]);
00835         DeleteNode(rpclNode->aptChildren[3]);
00836  #endif
00837     }
00838     delete rpclNode;
00839     rpclNode = NULL;
00840 }
00841 #endif
00842 
00843
00844 //#ifndef OSG_ARBITRARY_SPLIT
00845 void CErrorQuadTree::SetInitialCells(DCTPMesh *pclMesh, float fError,
00846 #ifdef OSG_USE_NURBS_PATCH
00847  #ifdef OSG_ARBITRARY_SPLIT
00848                                      BSplineTensorSurface *pclPatch,
00849                                      const Vec2d cclMinParam, const Vec2d cclMaxParam)
00850  #else
00851 BSplineTensorSurface * pclPatch)
00852  #endif
00853 #else
00854 std::vector<std::vector<BezierTensorSurface> > *pvvclPatches,
00855 const std::vector<double> *cpvdIntervalsU,
00856 const std::vector<double> *cpvdIntervalsV)
00857 #endif
00858 {
00859 #ifdef OSG_USE_NURBS_PATCH
00860     double         d_min_x;
00861     double         d_max_x;
00862     double         d_min_y;
00863     double         d_max_y;
00864     Vec3d          acl_edge_points[4];
00865     DCTPFace      *pcl_face;
00866     SFaceTreeCell *pt_finfo;
00867
00868  #ifdef OSG_ARBITRARY_SPLIT
00869     d_min_x = cclMinParam[0];
00870     d_min_y = cclMinParam[1];
00871     d_max_x = cclMaxParam[0];
00872     d_max_y = cclMaxParam[1];
00873  #else
00874     pclPatch->getParameterInterval_U(d_min_x, d_max_x);
00875     pclPatch->getParameterInterval_V(d_min_y, d_max_y);
00876  #endif
00877     acl_edge_points[0][0] = d_min_x;
00878     acl_edge_points[1][0] = d_max_x;
00879     acl_edge_points[2][0] = d_max_x;
00880     acl_edge_points[3][0] = d_min_x;
00881     acl_edge_points[0][1] = d_max_y;
00882     acl_edge_points[1][1] = d_max_y;
00883     acl_edge_points[2][1] = d_min_y;
00884     acl_edge_points[3][1] = d_min_y;
00885     acl_edge_points[0][2] = 0.0;
00886     acl_edge_points[1][2] = 0.0;
00887     acl_edge_points[2][2] = 0.0;
00888     acl_edge_points[3][2] = 0.0;
00889
00890     if(m_ptRoot == NULL)
00891     {
00892         m_ptRoot = new SErrorTreeCell;
00893  #ifdef OSG_ONE_CHILD_PTR
00894         m_ptRoot->ptChildren = NULL;
00895  #else
00896         m_ptRoot->aptChildren[0] = NULL;
00897         m_ptRoot->aptChildren[1] = NULL;
00898         m_ptRoot->aptChildren[2] = NULL;
00899         m_ptRoot->aptChildren[3] = NULL;
00900  #endif
00901         m_ptRoot->fError = -1.0;
00902     }
00903
00904     pcl_face = pclMesh->AddQuad(acl_edge_points[0], acl_edge_points[1],
00905                                 acl_edge_points[2], acl_edge_points[3], 0.0);
00906     pt_finfo              = new SFaceTreeCell;
00907     pt_finfo->bOwnSurface = false;
00908  #ifdef OSG_ARBITRARY_SPLIT
00909     pt_finfo->pclBSplineSurface = NULL;
00910  #else
00911     pt_finfo->pclBSplineSurface = pclPatch;
00912  #endif
00913     pt_finfo->ptErrorCell = m_ptRoot;
00914     pcl_face->faceinfo    = static_cast<void*>(pt_finfo);
00915 #else
00916     const unsigned int cui_u_size = cpvdIntervalsU->size() - 1;
00917     const unsigned int cui_v_size = cpvdIntervalsV->size() - 1;
00918     Vec3d              acl_edge_points[4];
00919     unsigned int       ui_u;
00920     unsigned int       ui_v;
00921     DCTPFace          *pcl_face;
00922     SFaceTreeCell     *pt_finfo;
00923     unsigned int       ui_face;
00924     unsigned int       ui_face_cnt;
00925     float              f_act_error;
00926
00927     acl_edge_points[0][2] = 1.0;
00928     acl_edge_points[1][2] = 1.0;
00929     acl_edge_points[2][2] = 1.0;
00930     acl_edge_points[3][2] = 1.0;
00931
00932     if(m_vvptRoot.size() == 0)
00933     {
00934 //      std::cerr << "new QuadTree" << std::endl;
00935         m_vvptRoot.resize(cui_u_size);
00936
00937         for(ui_u = 0; ui_u < cui_u_size; ++ui_u)
00938         {
00939             m_vvptRoot[ui_u].resize(cui_v_size);
00940
00941             for(ui_v = 0; ui_v < cui_v_size; ++ui_v)
00942             {
00943                 m_vvptRoot[ui_u][ui_v] = new SErrorTreeCell;
00944  #ifdef OSG_ONE_CHILD_PTR
00945                 m_vvptRoot[ui_u][ui_v]->ptChildren = NULL;
00946  #else
00947                 m_vvptRoot[ui_u][ui_v]->aptChildren[0] = NULL;
00948                 m_vvptRoot[ui_u][ui_v]->aptChildren[1] = NULL;
00949                 m_vvptRoot[ui_u][ui_v]->aptChildren[2] = NULL;
00950                 m_vvptRoot[ui_u][ui_v]->aptChildren[3] = NULL;
00951  #endif
00952                 m_vvptRoot[ui_u][ui_v]->fError = -1.0;
00953             }
00954         }
00955
00956         ComputeBPTree(pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
00957     }
00958
00959     acl_edge_points[0][0] = (*cpvdIntervalsU)[0];
00960     acl_edge_points[1][0] = (*cpvdIntervalsU)[cui_u_size];
00961     acl_edge_points[2][0] = (*cpvdIntervalsU)[cui_u_size];
00962     acl_edge_points[3][0] = (*cpvdIntervalsU)[0];
00963     acl_edge_points[0][1] = (*cpvdIntervalsV)[cui_v_size];
00964     acl_edge_points[1][1] = (*cpvdIntervalsV)[cui_v_size];
00965     acl_edge_points[2][1] = (*cpvdIntervalsV)[0];
00966     acl_edge_points[3][1] = (*cpvdIntervalsV)[0];
00967     pcl_face              = pclMesh->AddQuad(acl_edge_points[0], acl_edge_points[1],
00968                                              acl_edge_points[2], acl_edge_points[3], 0.0);
00969     pt_finfo                   = new SFaceTreeCell;
00970     pt_finfo->bOwnSurface      = false;
00971     pt_finfo->pclBezierSurface = NULL;
00972     pt_finfo->ptErrorCell      = NULL;
00973     pt_finfo->ptBPCell         = m_ptBPRoot;
00974     pcl_face->faceinfo         = ( void* ) pt_finfo;
00975
00976     ui_face_cnt = 1;
00977
00978     for(ui_face = 0; ui_face < ui_face_cnt; ++ui_face)
00979     {
00980         pcl_face = pclMesh->faces[ui_face];
00981         if( ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell)
00982         {
00983             f_act_error = ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->fError;
00984             while(f_act_error > fError)
00985             {
00986                 if( ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->aptChildren[1])
00987                 {
00988                     // xx
00989                     // xx
00990                     pclMesh->SubdivideQuad(pcl_face);
00991                     SubdivideNode(pclMesh, pcl_face);
00992                     ui_face_cnt += 3;
00993                 }
00994                 else if( ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->aptChildren[0])
00995                 {
00996                     // x
00997                     // x
00998                     pclMesh->SubdivideQuadNS(pcl_face);
00999                     SubdivideNode(pclMesh, pcl_face);
01000                     ++ui_face_cnt;
01001                 }
01002                 else if( ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->aptChildren[2])
01003                 {
01004                     // xx
01005                     pclMesh->SubdivideQuadEW(pcl_face);
01006                     SubdivideNode(pclMesh, pcl_face);
01007                     ++ui_face_cnt;
01008                 }
01009                 if( ( ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell) &&
01010                     ( ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->aptChildren[3]) )
01011                 {
01012                     f_act_error = ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->fError;
01013                 }
01014                 else
01015                 {
01016                     f_act_error = 0.0;
01017                 }
01018             }
01019         }
01020         if( ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell == NULL)
01021         {
01022             delete ( ( SFaceTreeCell* ) pcl_face->faceinfo);
01023             pcl_face->faceinfo = NULL;
01024         }
01025     }
01026
01027     for(ui_face = 0; ui_face < ui_face_cnt; ++ui_face)
01028     {
01029         pcl_face = pclMesh->faces[ui_face];
01030         if(pcl_face->faceinfo)
01031         {
01032 /*          std::cerr << "quad:";
01033             std::cerr << pcl_face->orig_quad[ 0 ]->coords << " ";
01034             std::cerr << pcl_face->orig_quad[ 1 ]->coords << " ";
01035             std::cerr << pcl_face->orig_quad[ 2 ]->coords << " ";
01036             std::cerr << pcl_face->orig_quad[ 3 ]->coords << std::endl;
01037             std::cerr << ( ( SFaceTreeCell* ) pcl_face->faceinfo )->ptBPCell->uiLeft << " - ";
01038             std::cerr << ( ( SFaceTreeCell* ) pcl_face->faceinfo )->ptBPCell->uiRight << " , ";
01039             std::cerr << ( ( SFaceTreeCell* ) pcl_face->faceinfo )->ptBPCell->uiBottom << " - ";
01040             std::cerr << ( ( SFaceTreeCell* ) pcl_face->faceinfo )->ptBPCell->uiTop << std::endl;*/
01041  #ifdef OSG_UNION_TRI_QUAD
01042             pclMesh->MoveVertex(
01043                 pcl_face->orig_face[0],
01044                 Vec3d( (*cpvdIntervalsU)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiLeft],
01045                        (*cpvdIntervalsV)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiTop + 1],
01046                        0.0) );
01047             pclMesh->MoveVertex(
01048                 pcl_face->orig_face[1],
01049                 Vec3d( (*cpvdIntervalsU)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiRight + 1],
01050                        (*cpvdIntervalsV)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiTop + 1],
01051                        0.0) );
01052             pclMesh->MoveVertex(
01053                 pcl_face->orig_face[2],
01054                 Vec3d( (*cpvdIntervalsU)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiRight + 1],
01055                        (*cpvdIntervalsV)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiBottom],
01056                        0.0) );
01057             pclMesh->MoveVertex(
01058                 pcl_face->orig_face[3],
01059                 Vec3d( (*cpvdIntervalsU)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiLeft],
01060                        (*cpvdIntervalsV)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiBottom],
01061                        0.0) );
01062  #else
01063             pclMesh->MoveVertex(
01064                 pcl_face->orig_quad[0],
01065                 Vec3d( (*cpvdIntervalsU)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiLeft],
01066                        (*cpvdIntervalsV)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiTop + 1],
01067                        0.0) );
01068             pclMesh->MoveVertex(
01069                 pcl_face->orig_quad[1],
01070                 Vec3d( (*cpvdIntervalsU)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiRight + 1],
01071                        (*cpvdIntervalsV)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiTop + 1],
01072                        0.0) );
01073             pclMesh->MoveVertex(
01074                 pcl_face->orig_quad[2],
01075                 Vec3d( (*cpvdIntervalsU)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiRight + 1],
01076                        (*cpvdIntervalsV)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiBottom],
01077                        0.0) );
01078             pclMesh->MoveVertex(
01079                 pcl_face->orig_quad[3],
01080                 Vec3d( (*cpvdIntervalsU)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiLeft],
01081                        (*cpvdIntervalsV)[( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiBottom],
01082                        0.0) );
01083  #endif
01084 /*          std::cerr << "moved to:";
01085             std::cerr << pcl_face->orig_quad[ 0 ]->coords << " ";
01086             std::cerr << pcl_face->orig_quad[ 1 ]->coords << " ";
01087             std::cerr << pcl_face->orig_quad[ 2 ]->coords << " ";
01088             std::cerr << pcl_face->orig_quad[ 3 ]->coords << std::endl;*/
01089             if( ( ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->aptChildren[3] == NULL) &&
01090                 ( ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->fError > fError) )
01091             {
01092                 // single patch
01093                 ui_u = ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiLeft;
01094                 ui_v = ( ( SFaceTreeCell* ) pcl_face->faceinfo)->ptBPCell->uiBottom;
01095 //              std::cerr << ui_u << "," << ui_v << std::endl;
01096                 pt_finfo                   = ( ( SFaceTreeCell* ) pcl_face->faceinfo);
01097                 pt_finfo->bOwnSurface      = false;
01098                 pt_finfo->pclBezierSurface = &( (*pvvclPatches)[ui_u][ui_v]);
01099                 pt_finfo->ptErrorCell      = m_vvptRoot[ui_u][ui_v];
01100                 pt_finfo->ptBPCell         = NULL;
01101 //              pcl_face->faceinfo = ( void* ) pt_finfo;
01102             }
01103         }
01104     }
01105
01106 #endif
01107 }
01108
01109
01110 void CErrorQuadTree::SubdivideNode(DCTPMesh *pclMesh, DCTPFace *pclFace)
01111 {
01112     const unsigned int cui_face_cnt = pclMesh->faces.size();
01113     unsigned int       ui_child;
01114 #ifdef OSG_USE_KD_TREE
01115     DCTPFace *apcl_faces[2];
01116 #else
01117     DCTPFace *apcl_faces[4];
01118 #endif
01119     SFaceTreeCell *pcl_old_node = static_cast<SFaceTreeCell*>(pclFace->faceinfo);
01120     SFaceTreeCell *pcl_new_node;
01121
01122     pclFace->faceinfo = NULL;
01123
01124 #ifdef OSG_USE_NURBS_PATCH
01125     apcl_faces[0] = pclFace;
01126  #ifdef OSG_USE_KD_TREE
01127     apcl_faces[1] = pclMesh->faces[cui_face_cnt - 1];
01128  #else
01129     apcl_faces[1] = pclMesh->faces[cui_face_cnt - 3];
01130     apcl_faces[2] = pclMesh->faces[cui_face_cnt - 2];
01131     apcl_faces[3] = pclMesh->faces[cui_face_cnt - 1];
01132  #endif
01133 
01134  #ifdef OSG_USE_KD_TREE
01135 
01136     for(ui_child = 0; ui_child < 2; ++ui_child)
01137  #else
01138
01139     for(ui_child = 0; ui_child < 4; ++ui_child)
01140  #endif
01141     {
01142         pcl_new_node = new SFaceTreeCell;
01143 //      std::cerr << "subdivide with ErrorCell" << std::endl;
01144  #ifdef OSG_ONE_CHILD_PTR
01145         pcl_new_node->ptErrorCell = &(pcl_old_node->ptErrorCell->ptChildren[ui_child]);
01146  #else
01147         pcl_new_node->ptErrorCell = pcl_old_node->ptErrorCell->aptChildren[ui_child];
01148  #endif
01149         pcl_new_node->pclBSplineSurface = NULL;
01150         pcl_new_node->bOwnSurface       = false;
01151         apcl_faces[ui_child]->faceinfo  = static_cast<void*>(pcl_new_node);
01152     }
01153
01154 #else
01155     if(pcl_old_node->ptErrorCell)
01156     {
01157         apcl_faces[0] = pclFace;
01158         apcl_faces[1] = pclMesh->faces[cui_face_cnt - 3];
01159         apcl_faces[2] = pclMesh->faces[cui_face_cnt - 2];
01160         apcl_faces[3] = pclMesh->faces[cui_face_cnt - 1];
01161
01162         for(ui_child = 0; ui_child < 4; ++ui_child)
01163         {
01164             pcl_new_node = new SFaceTreeCell;
01165 //          std::cerr << "subdivide with ErrorCell" << std::endl;
01166  #ifdef OSG_ONE_CHILD_PTR
01167             pcl_new_node->ptErrorCell = &(pcl_old_node->ptErrorCell->ptChildren[ui_child]);
01168  #else
01169             pcl_new_node->ptErrorCell = pcl_old_node->ptErrorCell->aptChildren[ui_child];
01170  #endif
01171             pcl_new_node->ptBPCell         = NULL;
01172             pcl_new_node->pclBezierSurface = NULL;
01173             pcl_new_node->bOwnSurface      = false;
01174             apcl_faces[ui_child]->faceinfo = static_cast<void*>(pcl_new_node);
01175         }
01176     }
01177     else if(pcl_old_node->ptBPCell->aptChildren[1])
01178     {
01179         // xx
01180         // xx
01181         apcl_faces[0] = pclFace;
01182         apcl_faces[1] = pclMesh->faces[cui_face_cnt - 3];
01183         apcl_faces[2] = pclMesh->faces[cui_face_cnt - 2];
01184         apcl_faces[3] = pclMesh->faces[cui_face_cnt - 1];
01185
01186         for(ui_child = 0; ui_child < 4; ++ui_child)
01187         {
01188             pcl_new_node = new SFaceTreeCell;
01189 //          std::cerr << "subdivide with BPCell" << std::endl;
01190             pcl_new_node->ptErrorCell = NULL;
01191             pcl_new_node->ptBPCell    = pcl_old_node->ptBPCell->aptChildren[ui_child];
01192 /*          std::cerr << pcl_new_node->ptBPCell->uiLeft << " - ";
01193             std::cerr << pcl_new_node->ptBPCell->uiRight << " , ";
01194             std::cerr << pcl_new_node->ptBPCell->uiBottom << " - ";
01195             std::cerr << pcl_new_node->ptBPCell->uiTop << std::endl;
01196             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 0 ]->coords << " ";
01197             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 1 ]->coords << " ";
01198             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 2 ]->coords << " ";
01199             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 3 ]->coords << std::endl;*/
01200             pcl_new_node->pclBezierSurface = NULL;
01201             pcl_new_node->bOwnSurface      = false;
01202             apcl_faces[ui_child]->faceinfo = static_cast<void*>(pcl_new_node);
01203         }
01204     }
01205     else if(pcl_old_node->ptBPCell->aptChildren[0])
01206     {
01207         // x
01208         // x
01209         apcl_faces[0] = pclFace;
01210         apcl_faces[1] = pclMesh->faces[cui_face_cnt - 1];
01211
01212         for(ui_child = 0; ui_child < 2; ++ui_child)
01213         {
01214             pcl_new_node = new SFaceTreeCell;
01215 //          std::cerr << "subdivide with BPCell" << std::endl;
01216             pcl_new_node->ptErrorCell = NULL;
01217             pcl_new_node->ptBPCell    = pcl_old_node->ptBPCell->aptChildren[3 * ui_child];
01218 /*          std::cerr << pcl_new_node->ptBPCell->uiLeft << " - ";
01219             std::cerr << pcl_new_node->ptBPCell->uiRight << " , ";
01220             std::cerr << pcl_new_node->ptBPCell->uiBottom << " - ";
01221             std::cerr << pcl_new_node->ptBPCell->uiTop << std::endl;
01222             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 0 ]->coords << " ";
01223             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 1 ]->coords << " ";
01224             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 2 ]->coords << " ";
01225             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 3 ]->coords << std::endl;*/
01226             pcl_new_node->pclBezierSurface = NULL;
01227             pcl_new_node->bOwnSurface      = false;
01228             apcl_faces[ui_child]->faceinfo = static_cast<void*>(pcl_new_node);
01229         }
01230     }
01231     else
01232     {
01233         // xx
01234         apcl_faces[0] = pclFace;
01235         apcl_faces[1] = pclMesh->faces[cui_face_cnt - 1];
01236
01237         for(ui_child = 0; ui_child < 2; ++ui_child)
01238         {
01239             pcl_new_node = new SFaceTreeCell;
01240 //          std::cerr << "subdivide with BPCell" << std::endl;
01241             pcl_new_node->ptErrorCell = NULL;
01242             pcl_new_node->ptBPCell    = pcl_old_node->ptBPCell->aptChildren[3 - ui_child];
01243 /*          std::cerr << pcl_new_node->ptBPCell->uiLeft << " - ";
01244             std::cerr << pcl_new_node->ptBPCell->uiRight << " , ";
01245             std::cerr << pcl_new_node->ptBPCell->uiBottom << " - ";
01246             std::cerr << pcl_new_node->ptBPCell->uiTop << std::endl;
01247             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 0 ]->coords << " ";
01248             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 1 ]->coords << " ";
01249             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 2 ]->coords << " ";
01250             std::cerr << apcl_faces[ ui_child ]->orig_quad[ 3 ]->coords << std::endl;*/
01251             pcl_new_node->pclBezierSurface = NULL;
01252             pcl_new_node->bOwnSurface      = false;
01253             apcl_faces[ui_child]->faceinfo = static_cast<void*>(pcl_new_node);
01254         }
01255     }
01256 #endif
01257 
01258     delete pcl_old_node;
01259 }
01260
01261
01262 void CErrorQuadTree::SubdivideBuild(DCTPMesh *pclMesh, DCTPFace *pclFace)
01263 {
01264     const unsigned int cui_face_cnt = pclMesh->faces.size();
01265     unsigned int       ui_child;
01266 #ifdef OSG_USE_KD_TREE
01267     DCTPFace *apcl_faces[2];
01268 #else
01269     DCTPFace *apcl_faces[4];
01270 #endif
01271     SFaceTreeCell *pcl_old_node = static_cast<SFaceTreeCell*>(pclFace->faceinfo);
01272     SFaceTreeCell *pcl_new_node;
01273 #ifdef OSG_USE_NURBS_PATCH
01274  #ifdef OSG_USE_KD_TREE
01275     std::vector<BSplineTensorSurface> vcl_surfaces;
01276  #else
01277     std::vector<std::vector<BSplineTensorSurface> > vvcl_surfaces;
01278  #endif
01279 #else
01280  #ifdef OSG_USE_KD_TREE
01281     std::vector<BezierTensorSurface> vcl_surfaces;
01282  #else
01283     std::vector<std::vector<BezierTensorSurface> > vvcl_surfaces;
01284  #endif
01285 #endif
01286 
01287     apcl_faces[0] = pclFace;
01288 #ifdef OSG_USE_KD_TREE
01289     apcl_faces[1] = pclMesh->faces[cui_face_cnt - 1];
01290 #else
01291     apcl_faces[1] = pclMesh->faces[cui_face_cnt - 3];
01292     apcl_faces[2] = pclMesh->faces[cui_face_cnt - 2];
01293     apcl_faces[3] = pclMesh->faces[cui_face_cnt - 1];
01294 #endif
01295 
01296     pclFace->faceinfo = NULL;
01297 #ifdef OSG_USE_NURBS_PATCH
01298  #ifdef OSG_USE_KD_TREE
01299   #ifdef OSG_ARBITRARY_SPLIT
01300     if(pcl_old_node->ptErrorCell->fSplitValue < 0.0)
01301     {
01302         pcl_old_node->pclBSplineSurface->subDivisionU(vcl_surfaces, -pcl_old_node->ptErrorCell->fSplitValue);
01303     }
01304     else
01305     {
01306         pcl_old_node->pclBSplineSurface->subDivisionV(vcl_surfaces, pcl_old_node->ptErrorCell->fSplitValue);
01307     }
01308   #else
01309     if(pcl_old_node->ptErrorCell->bSplitHoriz)
01310     {
01311         pcl_old_node->pclBSplineSurface->midPointSubDivisionU(vcl_surfaces);
01312     }
01313     else
01314     {
01315         pcl_old_node->pclBSplineSurface->midPointSubDivisionV(vcl_surfaces);
01316     }
01317   #endif
01318  #else
01319     pcl_old_node->pclBSplineSurface->midPointSubDivision(vvcl_surfaces);
01320  #endif
01321 #else
01322  #ifdef OSG_USE_KD_TREE
01323   #error
01324  #else
01325     pcl_old_node->pclBezierSurface->midPointSubDivision(vvcl_surfaces);
01326  #endif
01327 #endif
01328 
01329 #ifdef OSG_ONE_CHILD_PTR
01330  #ifdef OSG_USE_KD_TREE
01331 
01332     for(ui_child = 0; ui_child < 2; ++ui_child)
01333  #else
01334
01335     for(ui_child = 0; ui_child < 4; ++ui_child)
01336  #endif
01337     {
01338         pcl_new_node              = new SFaceTreeCell;
01339         pcl_new_node->ptErrorCell = &(pcl_old_node->ptErrorCell->ptChildren[ui_child]);
01340         pcl_new_node->bOwnSurface = true;
01341  #ifdef OSG_USE_NURBS_PATCH
01342         pcl_new_node->pclBSplineSurface = new BSplineTensorSurface;
01343   #ifdef OSG_USE_KD_TREE
01344                                                           (*pcl_new_node->pclBSplineSurface) = vcl_surfaces[ui_child];
01345   #else
01346         (*pcl_new_node->pclBSplineSurface) = vvcl_surfaces[( (ui_child + 1) >> 1) & 1][1 - (ui_child >> 1)];
01347   #endif
01348  #else
01349         pcl_new_node->pclBezierSurface = new BezierTensorSurface;
01350   #ifdef OSG_USE_KD_TREE
01351                                                           (*pcl_new_node->pclBezierSurface) = vvcl_surfaces[ui_child];
01352   #else
01353         (*pcl_new_node->pclBezierSurface) = vvcl_surfaces[( (ui_child + 1) >> 1) & 1][1 - (ui_child >> 1)];
01354   #endif
01355  #endif
01356         apcl_faces[ui_child]->faceinfo = static_cast<void*>(pcl_new_node);
01357     }
01358
01359     if(pcl_old_node->ptErrorCell->ptChildren == NULL)
01360     {
01361  #ifdef OSG_USE_KD_TREE
01362         pcl_old_node->ptErrorCell->ptChildren = new SErrorTreeCell[2];
01363
01364         for(ui_child = 0; ui_child < 2; ++ui_child)
01365  #else
01366         pcl_old_node->ptErrorCell->ptChildren = new SErrorTreeCell[4];
01367
01368         for(ui_child = 0; ui_child < 4; ++ui_child)
01369  #endif
01370         {
01371             pcl_new_node                          = static_cast<SFaceTreeCell*>(apcl_faces[ui_child]->faceinfo);
01372             pcl_new_node->ptErrorCell             = &(pcl_old_node->ptErrorCell->ptChildren[ui_child]);
01373             pcl_new_node->ptErrorCell->ptChildren = NULL;
01374             pcl_new_node->ptErrorCell->fError     = -1.0;
01375             ComputeError(apcl_faces[ui_child]);
01376         }
01377     }
01378 #else
01379 #ifdef OSG_USE_KD_TREE
01380 
01381     for(ui_child = 0; ui_child < 2; ++ui_child)
01382 #else
01383
01384     for(ui_child = 0; ui_child < 4; ++ui_child)
01385 #endif
01386     {
01387         pcl_new_node                      = new SFaceTreeCell;
01388         pcl_new_node->ptErrorCell         = pcl_old_node->ptErrorCell->aptChildren[ui_child];
01389         pcl_new_node->bOwnSurface         = true;
01390         pcl_new_node->pclBezierSurface    = new BezierTensorSurface;
01391         (*pcl_new_node->pclBezierSurface) = vvcl_surfaces[( (ui_child + 1) >> 1) & 1][1 - (ui_child >> 1)];
01392         apcl_faces[ui_child]->faceinfo    = static_cast<void*>(pcl_new_node);
01393         if(pcl_new_node->ptErrorCell == NULL)
01394         {
01395             pcl_old_node->ptErrorCell->aptChildren[ui_child] = new SErrorTreeCell;
01396             pcl_new_node->ptErrorCell                        = pcl_old_node->ptErrorCell->aptChildren[ui_child];
01397             pcl_new_node->ptErrorCell->aptChildren[0]        = NULL;
01398             pcl_new_node->ptErrorCell->aptChildren[1]        = NULL;
01399  #ifndef OSG_USE_KD_TREE
01400             pcl_new_node->ptErrorCell->aptChildren[2] = NULL;
01401             pcl_new_node->ptErrorCell->aptChildren[3] = NULL;
01402  #endif
01403             pcl_new_node->ptErrorCell->fError = -1.0;
01404             ComputeError(apcl_faces[ui_child]);
01405         }
01406     }
01407
01408 #endif
01409 /*  ( SFaceTreeCell* )( apcl_faces[ 0 ]->faceinfo )->pclBezierSurface = vvcl_surfaces[ 0 ][ 1 ];
01410     ( SFaceTreeCell* )( apcl_faces[ 1 ]->faceinfo )->pclBezierSurface = vvcl_surfaces[ 1 ][ 1 ];
01411     ( SFaceTreeCell* )( apcl_faces[ 2 ]->faceinfo )->pclBezierSurface = vvcl_surfaces[ 1 ][ 0 ];
01412     ( SFaceTreeCell* )( apcl_faces[ 3 ]->faceinfo )->pclBezierSurface = vvcl_surfaces[ 0 ][ 0 ];*/
01413
01414     if(pcl_old_node->bOwnSurface)
01415     {
01416 #ifdef OSG_USE_NURBS_PATCH
01417         delete pcl_old_node->pclBSplineSurface;
01418 #else
01419         delete pcl_old_node->pclBezierSurface;
01420 #endif
01421     }
01422     delete pcl_old_node;
01423 }
01424 //#endif
01425
01426
01427 /*#ifdef OSG_ARBITRARY_SPLIT
01428 void CErrorQuadTree::ComputeError( SFaceTreeCell *ptCell )
01429 #else*/
01430 void CErrorQuadTree::ComputeError(DCTPFace *pclFace)
01431 //#endif
01432 {
01433 //#ifndef OSG_ARBITRARY_SPLIT
01434     SFaceTreeCell *ptCell = static_cast<SFaceTreeCell*>(pclFace->faceinfo);
01435 //#endif
01436
01437     if(ptCell->ptErrorCell->fError >= 0.0)
01438     {
01439         return; // error was already calculated
01440     }
01441
01442     if( ( (pclFace->orig_face[1]->coords[0] - pclFace->orig_face[3]->coords[0]) < 10.0) ||
01443         ( (pclFace->orig_face[1]->coords[1] - pclFace->orig_face[3]->coords[1]) < 10.0) )
01444     {
01445         // too small face!
01446         ptCell->ptErrorCell->fError = 0.0;
01447         return;
01448     }
01449
01450 #ifdef OSG_USE_NURBS_PATCH
01451     BSplineTensorSurface      *pcl_surface = ptCell->pclBSplineSurface;
01452     const std::vector<double> &crvd_knot_u = pcl_surface->getKnotVector_U();
01453     const std::vector<double> &crvd_knot_v = pcl_surface->getKnotVector_V();
01454     const unsigned int         cui_dim_u   = pcl_surface->getDimension_U();
01455     const unsigned int         cui_dim_v   = pcl_surface->getDimension_V();
01456     unsigned int               ui_idx;
01457 #else
01458     BezierTensorSurface *pcl_surface = ptCell->pclBezierSurface;
01459 #endif
01460     const std::vector<std::vector <Vec4d> > &crvvcl_cps = pcl_surface->getControlPointMatrix();
01461     const unsigned int                       cui_m      = crvvcl_cps.size() - 1;
01462     const unsigned int                       cui_n      = crvvcl_cps[0].size() - 1;
01463     Vec3d                                    ccl_b00;
01464     ccl_b00[0] = crvvcl_cps[0][0][0] / crvvcl_cps[0][0][3];
01465     ccl_b00[1] = crvvcl_cps[0][0][1] / crvvcl_cps[0][0][3];
01466     ccl_b00[2] = crvvcl_cps[0][0][2] / crvvcl_cps[0][0][3];
01467     Vec3d ccl_b0n;
01468     ccl_b0n[0] = crvvcl_cps[0][cui_n][0] / crvvcl_cps[0][cui_n][3];
01469     ccl_b0n[1] = crvvcl_cps[0][cui_n][1] / crvvcl_cps[0][cui_n][3];
01470     ccl_b0n[2] = crvvcl_cps[0][cui_n][2] / crvvcl_cps[0][cui_n][3];
01471     Vec3d ccl_bm0;
01472     ccl_bm0[0] = crvvcl_cps[cui_m][0][0] / crvvcl_cps[cui_m][0][3];
01473     ccl_bm0[1] = crvvcl_cps[cui_m][0][1] / crvvcl_cps[cui_m][0][3];
01474     ccl_bm0[2] = crvvcl_cps[cui_m][0][2] / crvvcl_cps[cui_m][0][3];
01475     Vec3d ccl_bmn;
01476     ccl_bmn[0] = crvvcl_cps[cui_m][cui_n][0] / crvvcl_cps[cui_m][cui_n][3];
01477     ccl_bmn[1] = crvvcl_cps[cui_m][cui_n][1] / crvvcl_cps[cui_m][cui_n][3];
01478     ccl_bmn[2] = crvvcl_cps[cui_m][cui_n][2] / crvvcl_cps[cui_m][cui_n][3];
01479     Vec3d        cl_cij;
01480     Vec3d        cl_bij;
01481     double       d_quad_size;
01482     unsigned int ui_i;
01483     unsigned int ui_j;
01484 #ifndef OSG_USE_NURBS_PATCH
01485     const double cd_rn = 1.0 / cui_n;
01486     const double cd_rm = 1.0 / cui_m;
01487 #endif
01488     Vec3d cl_ci0;
01489     Vec3d cl_cin;
01490 #ifndef OSG_USE_NURBS_PATCH
01491     const Vec3d ccl_dcx0 = (ccl_bm0 - ccl_b00) * cd_rm;
01492     const Vec3d ccl_dcxn = (ccl_bmn - ccl_b0n) * cd_rm;
01493 #endif
01494     Vec3d cl_norm;
01495 #ifndef OSG_USE_NURBS_PATCH
01496     const Vec3d ccl_dc0x = (ccl_b0n - ccl_b00) * cd_rn;
01497     const Vec3d ccl_dcmx = (ccl_bmn - ccl_bm0) * cd_rn;
01498     Vec3d       cl_c0j;
01499     Vec3d       cl_cmj;
01500 #endif
01501     int   i_err;
01502     Vec2d cl_uv;
01503 #ifdef OSG_USE_NURBS_PATCH
01504     Vec2d cl_min;
01505     Vec2d cl_add;
01506 #endif
01507 #ifdef OSG_ARBITRARY_SPLIT
01508     Vec2d cl_worst;
01509 #endif
01510     // normal stuff
01511     Vec3d cl_nb00;
01512     Vec3d cl_nb0n;
01513     Vec3d cl_nbm0;
01514     Vec3d cl_nbmn;
01515     Vec3d cl_ncij;
01516     Vec3d cl_nbij;
01517     Vec3d cl_nci0;
01518     Vec3d cl_ncin;
01519 #ifndef OSG_USE_NURBS_PATCH
01520     Vec3d cl_ndcx0;
01521     Vec3d cl_ndcxn;
01522     Vec3d cl_ndc0x;
01523     Vec3d cl_ndcmx;
01524     Vec3d cl_ndcix;
01525     Vec3d cl_ndcxj;
01526     Vec3d cl_nc0j;
01527     Vec3d cl_ncmj;
01528 #endif
01529     Vec3f  cl_normal;
01530     Pnt3f  cl_position;
01531     double d_quadcurv = 0.0;
01532
01533 #ifdef OSG_USE_NURBS_PATCH
01534     pcl_surface->getParameterInterval_U(cl_min[0], cl_add[0]);
01535     pcl_surface->getParameterInterval_V(cl_min[1], cl_add[1]);
01536     cl_add -= cl_min;
01537  #ifdef OSG_ARBITRARY_SPLIT
01538     cl_worst = cl_min;
01539  #endif
01540 #endif
01541 
01542     if(m_sbNormalApproximation)
01543     {
01544 //#ifndef OSG_USE_NURBS_PATCH
01545
01546         Vec3d du, dv;
01547         du[0] = crvvcl_cps[1][0][0];
01548         du[1] = crvvcl_cps[1][0][1];
01549         du[2] = crvvcl_cps[1][0][2];
01550
01551         // we take the difference of two homogenious cps, but as the length of the cross product
01552         // is irrevelant (we're calculating normals) we ignore the weights.
01553         // In addition the weight of cp[1][0] might be zero so we avoid dividing by it (multiply is OK).
01554         du -= ccl_b00 * crvvcl_cps[1][0][3];
01555
01556         dv[0]   = crvvcl_cps[0][1][0];
01557         dv[1]   = crvvcl_cps[0][1][1];
01558         dv[2]   = crvvcl_cps[0][1][2];
01559         dv     -= ccl_b00 * crvvcl_cps[0][1][3];
01560         cl_nb00 = du.cross(dv);
01561
01562         du[0]   = crvvcl_cps[1][cui_n][0];
01563         du[1]   = crvvcl_cps[1][cui_n][1];
01564         du[2]   = crvvcl_cps[1][cui_n][2];
01565         du     -= ccl_b0n * crvvcl_cps[1][cui_n][3];
01566         dv[0]   = crvvcl_cps[0][cui_n - 1][0];
01567         dv[1]   = crvvcl_cps[0][cui_n - 1][1];
01568         dv[2]   = crvvcl_cps[0][cui_n - 1][2];
01569         dv     -= ccl_b0n * crvvcl_cps[0][cui_n - 1][3];
01570         cl_nb0n = du.cross(-dv);
01571
01572         du[0]   = crvvcl_cps[cui_m - 1][0][0];
01573         du[1]   = crvvcl_cps[cui_m - 1][0][1];
01574         du[2]   = crvvcl_cps[cui_m - 1][0][2];
01575         du     -= ccl_bm0 * crvvcl_cps[cui_m - 1][0][3];
01576         dv[0]   = crvvcl_cps[cui_m][1][0];
01577         dv[1]   = crvvcl_cps[cui_m][1][1];
01578         dv[2]   = crvvcl_cps[cui_m][1][2];
01579         dv     -= ccl_bm0 * crvvcl_cps[cui_m][1][3];
01580         cl_nbm0 = -du.cross(dv);
01581
01582         du[0]   = crvvcl_cps[cui_m - 1][cui_n][0];
01583         du[1]   = crvvcl_cps[cui_m - 1][cui_n][1];
01584         du[2]   = crvvcl_cps[cui_m - 1][cui_n][2];
01585         du     -= ccl_bmn * crvvcl_cps[cui_m - 1][cui_n][3];
01586         dv[0]   = crvvcl_cps[cui_m][cui_n - 1][0];
01587         dv[1]   = crvvcl_cps[cui_m][cui_n - 1][1];
01588         dv[2]   = crvvcl_cps[cui_m][cui_n - 1][2];
01589         dv     -= ccl_bmn * crvvcl_cps[cui_m][cui_n - 1][3];
01590         cl_nbmn = du.cross(dv);    // -du X -dv  => du X dv
01591
01592 /*
01593         cl_nb00 = ( crvvcl_cps[ 1 ][ 0 ] - crvvcl_cps[ 0 ][ 0 ] ).cross(
01594                     ( crvvcl_cps[ 0 ][ 1 ] - crvvcl_cps[ 0 ][ 0 ] ) );
01595         cl_nb0n = ( crvvcl_cps[ 1 ][ cui_n ] - crvvcl_cps[ 0 ][ cui_n ] ).cross(
01596                     ( crvvcl_cps[ 0 ][ cui_n ] - crvvcl_cps[ 0 ][ cui_n - 1 ] ) );
01597         cl_nbm0 = ( crvvcl_cps[ cui_m ][ 1 ] - crvvcl_cps[ cui_m ][ 0 ] ).cross(
01598                     ( crvvcl_cps[ cui_m - 1 ][ 0 ] - crvvcl_cps[ cui_m ][ 0 ] ) );
01599         cl_nbmn = ( crvvcl_cps[ cui_m ][ cui_n ] - crvvcl_cps[ cui_m - 1 ][ cui_n ] ).cross(
01600                     ( crvvcl_cps[ cui_m ][ cui_n ] - crvvcl_cps[ cui_m ][ cui_n - 1 ] ) );
01601 */
01602 //      std::cerr << cl_nb00 << std::endl;
01603 //      std::cerr << cl_nb0n << std::endl;
01604 //      std::cerr << cl_nbm0 << std::endl;
01605 //      std::cerr << cl_nbmn << std::endl << std::endl;
01606 /*#else
01607         cl_uv = cl_min;
01608         cl_normal = pcl_surface->computeNormal( cl_uv, i_err, cl_position );
01609         cl_nb00[0] = cl_normal[ 0 ];
01610         cl_nb00[1] = cl_normal[ 1 ];
01611         cl_nb00[2] = cl_normal[ 2 ];
01612 
01613         cl_uv[1] += cl_add[1];
01614         cl_normal = pcl_surface->computeNormal( cl_uv, i_err, cl_position );
01615         cl_nb0n[0] = cl_normal[ 0 ];
01616         cl_nb0n[1] = cl_normal[ 1 ];
01617         cl_nb0n[2] = cl_normal[ 2 ];
01618 
01619         cl_uv[0] += cl_add[0];
01620         cl_normal = pcl_surface->computeNormal( cl_uv, i_err, cl_position );
01621         cl_nbmn[0] = cl_normal[ 0 ];
01622         cl_nbmn[1] = cl_normal[ 1 ];
01623         cl_nbmn[2] = cl_normal[ 2 ];
01624 
01625         cl_uv[1] = cl_min[1];
01626         cl_normal = pcl_surface->computeNormal( cl_uv, i_err, cl_position );
01627         cl_nbm0[0] = cl_normal[ 0 ];
01628         cl_nbm0[1] = cl_normal[ 1 ];
01629         cl_nbm0[2] = cl_normal[ 2 ];
01630 //      std::cerr << cl_nb00 << std::endl;
01631 //      std::cerr << cl_nb0n << std::endl;
01632 //      std::cerr << cl_nbm0 << std::endl;
01633 //      std::cerr << cl_nbmn << std::endl << std::endl;
01634 #endif*/
01635
01636         d_quad_size = cl_nb00.squareLength();
01637         if(d_quad_size > DCTP_EPS * DCTP_EPS)
01638         {
01639             cl_nb00 *= 1.0 / sqrt(d_quad_size);
01640         }
01641         d_quad_size = cl_nb0n.squareLength();
01642         if(d_quad_size > DCTP_EPS * DCTP_EPS)
01643         {
01644             cl_nb0n *= 1.0 / sqrt(d_quad_size);
01645         }
01646         d_quad_size = cl_nbm0.squareLength();
01647         if(d_quad_size > DCTP_EPS * DCTP_EPS)
01648         {
01649             cl_nbm0 *= 1.0 / sqrt(d_quad_size);
01650         }
01651         d_quad_size = cl_nbmn.squareLength();
01652         if(d_quad_size > DCTP_EPS * DCTP_EPS)
01653         {
01654             cl_nbmn *= 1.0 / sqrt(d_quad_size);
01655         }
01656         // AKOS !!!!!!!!!!!!!!!!!!!!!!!!!!!!!
01657 //      std::cerr << cl_nb00 << std::endl;
01658 //      std::cerr << cl_nb0n << std::endl;
01659 //      std::cerr << cl_nbm0 << std::endl;
01660 //      std::cerr << cl_nbmn << std::endl << std::endl;
01661
01662 #ifndef OSG_USE_NURBS_PATCH
01663         cl_ndcx0 = (cl_nbm0 - cl_nb00) * cd_rm;
01664         cl_ndcxn = (cl_nbmn - cl_nb0n) * cd_rm;
01665         cl_ndc0x = (cl_nb0n - cl_nb00) * cd_rn;
01666         cl_ndcmx = (cl_nbmn - cl_nbm0) * cd_rn;
01667 #endif
01668 
01669         if( (cl_nb00.squareLength() < DCTP_EPS * DCTP_EPS) ||
01670             (cl_nb0n.squareLength() < DCTP_EPS * DCTP_EPS) ||
01671             (cl_nbm0.squareLength() < DCTP_EPS * DCTP_EPS) ||
01672             (cl_nbmn.squareLength() < DCTP_EPS * DCTP_EPS) ||
01673             ((ccl_b00 - ccl_b0n).squareLength() < DCTP_EPS * DCTP_EPS) ||
01674             ((ccl_b0n - ccl_bmn).squareLength() < DCTP_EPS * DCTP_EPS) ||
01675             ((ccl_bmn - ccl_bm0).squareLength() < DCTP_EPS * DCTP_EPS) ||
01676             ((ccl_bm0 - ccl_b00).squareLength() < DCTP_EPS * DCTP_EPS))
01677         {
01678             d_quadcurv = 1e100;
01679         }
01680         else
01681         {
01682             d_quadcurv = osgMax(osgMax( (cl_nb00 - cl_nb0n).squareLength() / (ccl_b00 - ccl_b0n).squareLength(),
01683                                         (cl_nb0n - cl_nbmn).squareLength() / (ccl_b0n - ccl_bmn).squareLength() ),
01684                                 osgMax( (cl_nbmn - cl_nbm0).squareLength() / (ccl_bmn - ccl_bm0).squareLength(),
01685                                         (cl_nbm0 - cl_nb00).squareLength() / (ccl_bm0 - ccl_b00).squareLength() ) );
01686
01687 //          d_quadcurv *= 4.0;  // 0.5 pixel geometric, 1 pixel shading error
01688         }
01689
01690         if(d_quadcurv < 1e-4)
01691         {
01692 //          std::cerr << d_quadcurv << std::endl;
01693             d_quadcurv = 1e-4;
01694         }
01695 /*      else
01696         {
01697             std::cerr << d_quadcurv << std::endl;
01698         }*/
01699     }
01700
01701 #ifndef OSG_USE_NURBS_PATCH
01702     cl_ci0 = ccl_b00;
01703     cl_cin = ccl_b0n;
01704     if(m_sbNormalApproximation)
01705     {
01706         cl_nci0 = cl_nb00;
01707         cl_ncin = cl_nb0n;
01708     }
01709 #endif
01710 
01711     for(ui_i = 0; ui_i <= cui_m; ++ui_i)
01712     {
01713 #ifdef OSG_USE_NURBS_PATCH
01714 //      cl_uv[0] += ( crvd_knot_u[ cui_dim_u + ui_i ] - crvd_knot_u[ ui_i ] ) / cui_dim_u;
01715         cl_uv[0] = 0.0;
01716
01717         for(ui_idx = 0; ui_idx < cui_dim_u; ++ui_idx)
01718         {
01719             cl_uv[0] += crvd_knot_u[ui_idx + ui_i + 1];
01720         }
01721
01722         cl_uv[0] /= cui_dim_u;
01723
01724         const double cd_xrel = (cl_uv[0] - cl_min[0]) / cl_add[0];
01725
01726         cl_ci0 = (ccl_bm0 - ccl_b00) * cd_xrel + ccl_b00;
01727         cl_cin = (ccl_bmn - ccl_b0n) * cd_xrel + ccl_b0n;
01728         if(m_sbNormalApproximation)
01729         {
01730             cl_nci0 = (cl_nbm0 - cl_nb00) * cd_xrel + cl_nb00;
01731             cl_ncin = (cl_nbmn - cl_nb0n) * cd_xrel + cl_nb0n;
01732         }
01733 #else
01734         cl_cij = cl_ci0;
01735         cl_c0j = ccl_b00;
01736         cl_cmj = ccl_bm0;
01737
01738         const Vec3d ccl_dcix = (cl_cin - cl_ci0) * cd_rn;
01739
01740         if(m_sbNormalApproximation)
01741         {
01742             cl_ncij = cl_nci0;
01743             cl_nc0j = cl_nb00;
01744             cl_ncmj = cl_nbm0;
01745
01746             cl_ndcix = (cl_ncin - cl_nci0) * cd_rn;
01747         }
01748
01749         cl_uv[0] = ui_i * cd_rm;
01750 #endif
01751 
01752         for(ui_j = 0; ui_j <= cui_n; ++ui_j)
01753         {
01754 #ifdef OSG_USE_NURBS_PATCH
01755             cl_uv[1] = 0.0;
01756
01757             for(ui_idx = 0; ui_idx < cui_dim_v; ++ui_idx)
01758             {
01759                 cl_uv[1] += crvd_knot_v[ui_idx + ui_j + 1];
01760             }
01761
01762             cl_uv[1] /= cui_dim_v;
01763
01764             const double cd_yrel = (cl_uv[1] - cl_min[1]) / cl_add[1];
01765
01766             cl_cij = cl_ci0 + (cl_cin - cl_ci0) * cd_yrel;
01767             if(m_sbNormalApproximation)
01768             {
01769                 cl_ncij = cl_nci0 + (cl_ncin - cl_nci0) * cd_yrel;
01770             }
01771 #else
01772             const Vec3d ccl_dcxj = cl_cmj - cl_c0j;
01773
01774             if(m_sbNormalApproximation)
01775             {
01776                 cl_ndcxj = cl_ncmj - cl_nc0j;
01777             }
01778
01779             cl_uv[1] = ui_j * cd_rn;
01780 #endif
01781 #ifdef OSG_DIFF_TO_BILIN
01782  #ifdef OSG_USE_NURBS_PATCH
01783   #ifdef OSG_CONSERVATIVE_ERROR
01784             cl_bij = crvvcl_cps[ui_i][ui_j];
01785   #else
01786             if(m_sbNormalApproximation)
01787             {
01788                 // TODO Vec3f f�r normale, Pnt3f f�r punkt
01789                 cl_nbij     = pcl_surface->computeNormal(cl_uv, i_err, cl_bij);
01790                 d_quad_size = cl_ncij.squareLength();
01791                 if(d_quad_size > DCTP_EPS * DCTP_EPS)
01792                 {
01793                     cl_norm = cl_ncij * (1.0 / sqrt(d_quad_size) );
01794                 }
01795                 else
01796                 {
01797                     cl_norm = cl_nbij;
01798                 }
01799 //              std::cerr << cl_uv << cl_nbij << cl_norm << std::endl;
01800                 cl_nbij -= cl_norm;
01801             }
01802             else
01803             {
01804                 cl_bij = pcl_surface->compute(cl_uv, i_err);
01805             }
01806   #endif
01807  #else
01808   #ifdef OSG_CONSERVATIVE_ERROR
01809             cl_bij = crvvcl_cps[ui_i][ui_j];
01810   #else
01811             cl_bij = pcl_surface->computewdeCasteljau(cl_uv, i_err);
01812   #endif
01813  #endif
01814             if( (m_sbNormalApproximation) && (cl_nbij.squareLength() > DCTP_EPS * DCTP_EPS) )
01815             {
01816                 // calculate shading error
01817                 d_quad_size = cl_nbij.squareLength() / d_quadcurv;
01818
01819  #ifdef OSG_EUCLID_ERRORS
01820                 // for euclidian norm first calculate distance to next correcly shaded pixel
01821                 d_quad_size += (cl_bij - cl_cij).squareLength();
01822  #endif
01823 
01824                 // the corner points are correcly shaded pixels...
01825                 d_quad_size = osgMin(d_quad_size, (cl_bij - ccl_b00).squareLength() );
01826                 d_quad_size = osgMin(d_quad_size, (cl_bij - ccl_b0n).squareLength() );
01827                 d_quad_size = osgMin(d_quad_size, (cl_bij - ccl_bm0).squareLength() );
01828                 d_quad_size = osgMin(d_quad_size, (cl_bij - ccl_bmn).squareLength() );
01829
01830  #ifdef OSG_EUCLID_ERRORS
01831                 d_quad_size = sqrt(d_quad_size);
01832  #else
01833                 // for maximum norm choose the maximum error
01834                 d_quad_size = sqrt(osgMax( (cl_bij - cl_cij).squareLength(), d_quad_size) );
01835  #endif
01836             }
01837             else
01838             {
01839                 d_quad_size = sqrt( (cl_bij - cl_cij).squareLength() );
01840             }
01841  #ifdef OSG_ARBITRARY_SPLIT
01842             if(osgAbs(d_quad_size - ptCell->ptErrorCell->fError) < DCTP_EPS)
01843             {
01844                 double cd_old_x = osgMin(cl_worst[0] - cl_min[0], cl_min[0] + cl_add[0] - cl_worst[0]);
01845                 double cd_old_y = osgMin(cl_worst[1] - cl_min[1], cl_min[1] + cl_add[1] - cl_worst[1]);
01846                 double cd_new_x = osgMin(cl_uv[0] - cl_min[0], cl_min[0] + cl_add[0] - cl_uv[0]);
01847                 double cd_new_y = osgMin(cl_uv[1] - cl_min[1], cl_min[1] + cl_add[1] - cl_uv[1]);
01848 //              std::cerr << cd_old_x * cd_old_x + cd_old_y * cd_old_y << " " << cd_new_x * cd_new_x + cd_new_y * cd_new_y << std::endl;
01849                 if(cd_old_x * cd_old_x + cd_old_y * cd_old_y <= cd_new_x * cd_new_x + cd_new_y * cd_new_y)
01850                 {
01851                     cl_worst = cl_uv;
01852                 }
01853                 if(d_quad_size > ptCell->ptErrorCell->fError)
01854                 {
01855                     ptCell->ptErrorCell->fError = float(d_quad_size);
01856                 }
01857             }
01858             else
01859  #endif
01860             if(d_quad_size > ptCell->ptErrorCell->fError)
01861             {
01862                 ptCell->ptErrorCell->fError = float(d_quad_size);
01863  #ifdef OSG_ARBITRARY_SPLIT
01864                 cl_worst = cl_uv;
01865  #endif
01866             }
01867 #else
01868  #ifdef OSG_CONSERVATIVE_ERROR
01869             cl_bij = crvvcl_cps[ui_i][ui_j];
01870  #else
01871   #ifdef OSG_USE_NURBS_PATCH
01872             cl_bij = pcl_surface->compute(cl_uv, i_err);
01873   #else
01874             cl_bij = pcl_surface->computewdeCasteljau(cl_uv, i_err);
01875   #endif
01876  #endif
01877             d_quad_size = osgMin(computeDistToTriangle(cl_bij, ccl_b00, ccl_b0n, ccl_bm0),
01878                                  computeDistToTriangle(cl_bij, ccl_bmn, ccl_b0n, ccl_bm0) );
01879 //          std::cerr << d_quad_size << std::endl;
01880             if(d_quad_size > ptCell->ptErrorCell->fError)
01881             {
01882 //              std::cerr << d_quad_size << std::endl;
01883                 ptCell->ptErrorCell->fError = ( float ) d_quad_size;
01884             }
01885             d_quad_size = osgMin(computeDistToTriangle(cl_bij, ccl_b0n, ccl_b00, ccl_bmn),
01886                                  computeDistToTriangle(cl_bij, ccl_bm0, ccl_b00, ccl_bmn) );
01887 //          std::cerr << d_quad_size << std::endl;
01888             if(d_quad_size > ptCell->ptErrorCell->fError)
01889             {
01890 //              std::cerr << d_quad_size << std::endl;
01891                 ptCell->ptErrorCell->fError = ( float ) d_quad_size;
01892             }
01893 #endif
01894 #ifndef OSG_USE_NURBS_PATCH
01895             cl_cij += ccl_dcix;
01896             cl_c0j += ccl_dc0x;
01897             cl_cmj += ccl_dcmx;
01898 #endif
01899         }
01900
01901 #ifndef OSG_USE_NURBS_PATCH
01902         cl_ci0 += ccl_dcx0;
01903         cl_cin += ccl_dcxn;
01904 #endif
01905     }
01906
01907     // add twist vector error
01908 #ifdef OSG_DIFF_TO_BILIN
01909  #ifdef OSG_CONSERVATIVE_ERROR
01910     ptCell->ptErrorCell->fError += ( float ) sqrt( (ccl_b00 - ccl_b0n + ccl_bmn - ccl_bm0).squareLength() ) * 0.25;
01911  #endif
01912 #endif
01913 //  ptCell->ptErrorCell->fError *= 0.1;
01914 //  std::cerr << ptCell->ptErrorCell->fError << std::endl;
01915 #ifdef OSG_USE_KD_TREE
01916  #ifdef OSG_ARBITRARY_SPLIT
01917 //  std::cerr << cl_min << cl_add << cl_worst << std::endl;
01918     if( (osgAbs(cl_worst[0] - cl_min[0]) <= osgMax(DCTP_EPS, 1e-4 * cl_add[0]) ) ||
01919         (osgAbs(cl_worst[0] - (cl_min[0] + cl_add[0]) ) <= osgMax(DCTP_EPS, 1e-4 * cl_add[0]) ) )
01920     {
01921 //      std::cerr << "y";
01922         if( (osgAbs(cl_worst[1] - cl_min[1]) <= osgMax(DCTP_EPS, 1e-4 * cl_add[1]) ) ||
01923             (osgAbs(cl_worst[1] - (cl_min[1] + cl_add[1]) ) <= osgMax(DCTP_EPS, 1e-4 * cl_add[1]) ) )
01924         {
01925             if(osgAbs(cl_add[0]) < osgAbs(cl_add[1]) )
01926             {
01927 //              std::cerr << "my";
01928                 ptCell->ptErrorCell->fSplitValue = 0.5;
01929             }
01930             else
01931             {
01932 //              std::cerr << "mx";
01933                 ptCell->ptErrorCell->fSplitValue = -0.5;
01934             }
01935         }
01936         else
01937         {
01938             ptCell->ptErrorCell->fSplitValue = (cl_worst[1] - cl_min[1]) / cl_add[1];
01939         }
01940     }
01941     else if( (osgAbs(cl_worst[1] - cl_min[1]) <= osgMax(DCTP_EPS, 1e-4 * cl_add[1]) ) ||
01942              (osgAbs(cl_worst[1] - (cl_min[1] + cl_add[1]) ) <= osgMax(DCTP_EPS, 1e-4 * cl_add[1]) ) )
01943     {
01944 //      std::cerr << "x";
01945         ptCell->ptErrorCell->fSplitValue = -(cl_worst[0] - cl_min[0]) / cl_add[0];
01946     }
01947     else
01948     {
01949         double d_hdiff;
01950         double d_vdiff;
01951
01952         cl_uv[0] = cl_worst[0];
01953         cl_uv[1] = cl_min[1];
01954         if(m_sbNormalApproximation)
01955         {
01956             cl_nci0   = pcl_surface->computeNormal(cl_uv, i_err, cl_ci0);
01957             cl_uv[1] += cl_add[1];
01958             cl_ncin   = pcl_surface->computeNormal(cl_uv, i_err, cl_cin);
01959             cl_uv[1]  = cl_worst[1];
01960             cl_ncij   = pcl_surface->computeNormal(cl_uv, i_err, cl_cij);
01961
01962             cl_bij      = cl_ci0 + (cl_cin - cl_ci0) * ( (cl_worst[1] - cl_min[1]) / cl_add[1]);
01963             cl_nbij     = cl_nci0 + (cl_ncin - cl_nci0) * ( (cl_worst[1] - cl_min[1]) / cl_add[1]);
01964             d_quad_size = cl_nbij.squareLength();
01965             if(d_quad_size > DCTP_EPS * DCTP_EPS)
01966             {
01967                 cl_nbij *= 1.0 / sqrt(d_quad_size);
01968             }
01969         }
01970         else
01971         {
01972             cl_ci0    = pcl_surface->compute(cl_uv, i_err);
01973             cl_uv[1] += cl_add[1];
01974             cl_cin    = pcl_surface->compute(cl_uv, i_err);
01975             cl_uv[1]  = cl_worst[1];
01976             cl_cij    = pcl_surface->compute(cl_uv, i_err);
01977
01978             cl_bij = cl_ci0 + (cl_cin - cl_ci0) * ( (cl_worst[1] - cl_min[1]) / cl_add[1]);
01979         }
01980
01981 //      const double    cd_hdiff = computeDistToLine( cl_cij, cl_ci0, cl_cin );
01982 //      const double    cd_hdiff = sqrt( ( cl_cij - cl_bij ).squareLength( ) );
01983         if( (m_sbNormalApproximation) && (cl_ncij.squareLength() > DCTP_EPS * DCTP_EPS) )
01984         {
01985             // calculate shading error
01986             d_hdiff = (cl_ncij - cl_nbij).squareLength() / d_quadcurv;
01987
01988    #ifdef OSG_EUCLID_ERRORS
01989             // for euclidian norm first calculate distance to next correcly shaded pixel
01990             d_hdiff += (cl_cij - cl_bij).squareLength();
01991    #endif
01992 
01993             // the new corner points will be correcly shaded pixels...
01994 //          d_hdiff = osgMin( d_hdiff, ( cl_cij - cl_ci0 ).squareLength( ) );
01995 //          d_hdiff = osgMin( d_hdiff, ( cl_cij - cl_cin ).squareLength( ) );
01996
01997    #ifdef OSG_EUCLID_ERRORS
01998             d_hdiff = sqrt(d_hdiff);
01999    #else
02000             // for maximum norm choose the maximum error
02001             d_hdiff = sqrt(osgMax( (cl_cij - cl_bij).squareLength(), d_hdiff) );
02002    #endif
02003 
02004 /*   #ifdef OSG_EUCLID_ERRORS
02005             d_hdiff = sqrt( ( cl_cij - cl_bij ).squareLength( )
02006                       + ( cl_ncij - cl_nbij ).squareLength( ) / d_quadcurv );
02007    #else
02008             d_hdiff = sqrt( osgMax( ( cl_cij - cl_bij ).squareLength( ),
02009                             ( cl_ncij - cl_nbij ).squareLength( ) / d_quadcurv ) );
02010    #endif*/
02011             d_hdiff += sqrt( (cl_ci0 - cl_cin).squareLength() ) * 0.0001;
02012         }
02013         else
02014         {
02015             d_hdiff = sqrt( (cl_cij - cl_bij).squareLength() )
02016                       + sqrt( (cl_ci0 - cl_cin).squareLength() ) * 0.0001;
02017         }
02018 //      const double    cd_hdiff = sqrt( ( cl_ci0 - cl_cin ).squareLength( ) );
02019
02020
02021
02022 /*      std::cerr.precision( 4 );
02023         if( m_sbNormalApproximation )
02024         {
02025             std::cerr << cl_ci0 << cl_cij << cl_cin << cl_nci0 << cl_ncij << cl_ncin << d_hdiff << std::endl;
02026         }
02027         else
02028         {
02029             std::cerr << cl_ci0 << cl_cij << cl_cin << d_hdiff << std::endl;
02030         }*/
02031
02032         if(m_sbNormalApproximation)
02033         {
02034             cl_uv[0]  = cl_min[0];
02035             cl_nci0   = pcl_surface->computeNormal(cl_uv, i_err, cl_ci0);
02036             cl_uv[0] += cl_add[0];
02037             cl_ncin   = pcl_surface->computeNormal(cl_uv, i_err, cl_cin);
02038
02039             cl_bij      = cl_ci0 + (cl_cin - cl_ci0) * ( (cl_worst[0] - cl_min[0]) / cl_add[0]);
02040             cl_nbij     = cl_nci0 + (cl_ncin - cl_nci0) * ( (cl_worst[0] - cl_min[0]) / cl_add[0]);
02041             d_quad_size = cl_nbij.squareLength();
02042             if(d_quad_size > DCTP_EPS * DCTP_EPS)
02043             {
02044                 cl_nbij *= 1.0 / sqrt(d_quad_size);
02045             }
02046         }
02047         else
02048         {
02049             cl_uv[0]  = cl_min[0];
02050             cl_ci0    = pcl_surface->compute(cl_uv, i_err);
02051             cl_uv[0] += cl_add[0];
02052             cl_cin    = pcl_surface->compute(cl_uv, i_err);
02053
02054             cl_bij = cl_ci0 + (cl_cin - cl_ci0) * ( (cl_worst[0] - cl_min[0]) / cl_add[0]);
02055         }
02056
02057 //      const double    cd_vdiff = computeDistToLine( cl_cij, cl_ci0, cl_cin );
02058 //      const double    cd_vdiff = sqrt( ( cl_cij - cl_bij ).squareLength( ) );
02059         if( (m_sbNormalApproximation) && (cl_ncij.squareLength() > DCTP_EPS * DCTP_EPS) )
02060         {
02061             // calculate shading error
02062             d_vdiff = (cl_ncij - cl_nbij).squareLength() / d_quadcurv;
02063
02064    #ifdef OSG_EUCLID_ERRORS
02065             // for euclidian norm first calculate distance to next correcly shaded pixel
02066             d_vdiff += (cl_cij - cl_bij).squareLength();
02067    #endif
02068 
02069             // the new corner points will be correcly shaded pixels...
02070 //          d_vdiff = osgMin( d_vdiff, ( cl_cij - cl_ci0 ).squareLength( ) );
02071 //          d_vdiff = osgMin( d_vdiff, ( cl_cij - cl_cin ).squareLength( ) );
02072
02073    #ifdef OSG_EUCLID_ERRORS
02074             d_vdiff = sqrt(d_vdiff);
02075    #else
02076             // for maximum norm choose the maximum error
02077             d_vdiff = sqrt(osgMax( (cl_cij - cl_bij).squareLength(), d_vdiff) );
02078    #endif
02079 
02080 /*   #ifdef OSG_EUCLID_ERRORS
02081             d_vdiff = sqrt( ( cl_cij - cl_bij ).squareLength( )
02082                             + ( cl_ncij - cl_nbij ).squareLength( ) / d_quadcurv );
02083    #else
02084             d_vdiff = sqrt( osgMax( ( cl_cij - cl_bij ).squareLength( ),
02085                                     ( cl_ncij - cl_nbij ).squareLength( ) / d_quadcurv ) );
02086    #endif*/
02087             d_vdiff += sqrt( (cl_ci0 - cl_cin).squareLength() ) * 0.0001;
02088         }
02089         else
02090         {
02091             d_vdiff = sqrt( (cl_cij - cl_bij).squareLength() )
02092                       + sqrt( (cl_ci0 - cl_cin).squareLength() ) * 0.0001;
02093         }
02094 //      const double    cd_vdiff = sqrt( ( cl_ci0 - cl_cin ).squareLength( ) );
02095
02096 /*      if( m_sbNormalApproximation )
02097         {
02098             std::cerr << cl_ci0 << cl_cij << cl_cin << cl_nci0 << cl_ncij << cl_ncin << d_vdiff << std::endl;
02099         }
02100         else
02101         {
02102             std::cerr << cl_ci0 << cl_cij << cl_cin << d_vdiff << std::endl;
02103         }*/
02104 //      std::cerr << cd_hdiff << " " << cd_vdiff << std::endl;
02105 //      std::cerr << "--------------------" << std::endl;
02106
02107         if(d_hdiff < d_vdiff)
02108         {
02109             ptCell->ptErrorCell->fSplitValue = -(cl_worst[0] - cl_min[0]) / cl_add[0];
02110         }
02111         else
02112         {
02113             ptCell->ptErrorCell->fSplitValue = (cl_worst[1] - cl_min[1]) / cl_add[1];
02114         }
02115     }
02116 //  std::cerr << cl_min << cl_min + cl_add << std::endl;
02117 //  std::cerr << ptCell->ptErrorCell->fSplitValue << std::endl;
02118  #else
02119     cl_uv[0]  = cl_min[0] + cl_add[0] * 0.5;
02120     cl_uv[1]  = cl_min[1];
02121     cl_ci0    = pcl_surface->compute(cl_uv, i_err);
02122     cl_uv[1] += cl_add[1];
02123     cl_cin    = pcl_surface->compute(cl_uv, i_err);
02124     cl_uv[1] -= cl_add[1] * 0.5;
02125     cl_cij    = pcl_surface->compute(cl_uv, i_err);
02126
02127 //  const double    cd_hdiff = computeDistToLine( cl_cij, cl_ci0, cl_cin );
02128 //  const double    cd_hdiff = sqrt( ( cl_cij - ( cl_ci0 + cl_cin ) * 0.5 ).squareLength( ) );
02129 //  const double    cd_hdiff = sqrt( ( cl_cij - cl_ci0 ).squareLength( ) ) + sqrt( ( cl_cij - cl_cin ).squareLength( ) );
02130     const double cd_hdiff = computeDistToLine(cl_cij, cl_ci0, cl_cin)
02131                             + sqrt( (cl_ci0 - cl_cin).squareLength() ) * 0.001;
02132
02133 //  std::cerr.precision( 4 );
02134 //  std::cerr << cl_ci0 << cl_cij << cl_cin << cd_hdiff << std::endl;
02135
02136     cl_uv[0]  = cl_min[0];
02137     cl_ci0    = pcl_surface->compute(cl_uv, i_err);
02138     cl_uv[0] += cl_add[0];
02139     cl_cin    = pcl_surface->compute(cl_uv, i_err);
02140
02141 //  const double    cd_vdiff = computeDistToLine( cl_cij, cl_ci0, cl_cin );
02142 //  const double    cd_vdiff = sqrt( ( cl_cij - ( cl_ci0 + cl_cin ) * 0.5 ).squareLength( ) );
02143 //  const double    cd_vdiff = sqrt( ( cl_cij - cl_ci0 ).squareLength( ) ) + sqrt( ( cl_cij - cl_cin ).squareLength( ) );
02144     const double cd_vdiff = computeDistToLine(cl_cij, cl_ci0, cl_cin)
02145                             + sqrt( (cl_ci0 - cl_cin).squareLength() ) * 0.001;
02146
02147 //  std::cerr << cl_ci0 << cl_cij << cl_cin << cd_vdiff << std::endl;
02148 //  std::cerr << cd_hdiff << " " << cd_vdiff << std::endl;
02149 //  std::cerr << "--------------------" << std::endl;
02150
02151     ptCell->ptErrorCell->bSplitHoriz = (cd_hdiff < cd_vdiff);
02152  #endif
02153 #endif
02154 
02155 /*  const double cd_size = 0.5 * sqrt( osgMax( osgMax( osgMax( ( ccl_b00 - ccl_b0n ).squareLength( ),
02156                                                                ( ccl_b0n - ccl_bmn ).squareLength( ) ),
02157                                                        osgMax( ( ccl_bmn - ccl_bm0 ).squareLength( ),
02158                                                                ( ccl_bm0 - ccl_b00 ).squareLength( ) ) ),
02159                                                osgMax( ( ccl_b00 - ccl_bmn ).squareLength( ),
02160                                                        ( ccl_b0n - ccl_bm0 ).squareLength( ) ) ) );
02161 
02162     if( ( cd_size > DCTP_EPS ) &&
02163         ( ptCell->ptErrorCell->fError > cd_size ) )
02164     {
02165         // error can't be larger than the current face, except for a torus. ;-)
02166         ptCell->ptErrorCell->fError = cd_size;
02167     }*/
02168 }
02169
02170
02171 #ifndef OSG_USE_NURBS_PATCH
02172 void CErrorQuadTree::ComputeBPTree(std::vector<std::vector<BezierTensorSurface> > *pvvclPatches,
02173                                    const std::vector<double>                     * cpvdIntervalsU,
02174                                    const std::vector<double>                     * cpvdIntervalsV)
02175 {
02176     m_ptBPRoot             = new SBPETreeCell;
02177     m_ptBPRoot->uiBottom   = 0;
02178     m_ptBPRoot->uiTop      = (*cpvdIntervalsV).size() - 2;
02179     m_ptBPRoot->uiLeft     = 0;
02180     m_ptBPRoot->uiRight    = (*cpvdIntervalsU).size() - 2;
02181     m_ptBPRoot->fPrevError = 1e+32;
02182     ComputeBPError(m_ptBPRoot, pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02183     SubdivideBPTree(m_ptBPRoot, pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02184 }
02185
02186
02187 void CErrorQuadTree::ComputeBPError(SBPETreeCell                                  * pclBPNode,
02188                                     std::vector<std::vector<BezierTensorSurface> > *pvvclPatches,
02189                                     const std::vector<double>                     * cpvdIntervalsU,
02190                                     const std::vector<double>                     * cpvdIntervalsV)
02191 {
02192     unsigned int         ui_u_surf;
02193     unsigned int         ui_v_surf;
02194     BezierTensorSurface *pcl_surface;
02195     Vec3d                cl_cij;
02196     Vec3d                cl_bij;
02197     double               d_quad_size;
02198     unsigned int         ui_i;
02199     unsigned int         ui_j;
02200     Vec3d                cl_ci0;
02201     Vec3d                cl_cin;
02202     Vec3d                cl_norm;
02203     Vec3d                cl_c0j;
02204     Vec3d                cl_cmj;
02205     int                  i_err;
02206     Vec2d                cl_uv;
02207
02208
02209 //  std::cerr << pclBPNode->uiLeft << " - " << pclBPNode->uiRight << std::endl;
02210 //  std::cerr << pclBPNode->uiBottom << " - " << pclBPNode->uiTop << std::endl;
02211
02212 //  std::cerr << ( *cpvdIntervalsU ).size( ) << " x " << ( *cpvdIntervalsV ).size( ) << std::endl;
02213
02214     pclBPNode->fError = 0.0;
02215
02216 /*  if( ( pclBPNode->uiLeft == pclBPNode->uiRight ) &&
02217         ( pclBPNode->uiBottom == pclBPNode->uiTop ) )
02218     {
02219         return; // handled as single bezier patch
02220     }*/
02221
02222     pcl_surface = &( (*pvvclPatches)[pclBPNode->uiLeft][pclBPNode->uiBottom]);
02223     const Vec3d ccl_a00 = pcl_surface->getControlPointMatrix()[0][0];
02224
02225     pcl_surface = &( (*pvvclPatches)[pclBPNode->uiRight][pclBPNode->uiBottom]);
02226     const Vec3d ccl_s0n = pcl_surface->getControlPointMatrix()[0][pcl_surface->getControlPointMatrix()[0].size() - 1];
02227     const Vec3d ccl_a0n = ccl_s0n - ccl_a00;
02228
02229     pcl_surface = &( (*pvvclPatches)[pclBPNode->uiLeft][pclBPNode->uiTop]);
02230     const Vec3d ccl_sm0 = pcl_surface->getControlPointMatrix()[pcl_surface->getControlPointMatrix().size() - 1][0];
02231     const Vec3d ccl_am0 = ccl_sm0 - ccl_a00;
02232
02233     pcl_surface = &( (*pvvclPatches)[pclBPNode->uiRight][pclBPNode->uiTop]);
02234     const Vec3d ccl_smn = pcl_surface->getControlPointMatrix()[pcl_surface->getControlPointMatrix().size() - 1][pcl_surface->getControlPointMatrix()[0].size() - 1];
02235     const Vec3d ccl_amn = ccl_smn - ccl_am0 - ccl_a00;
02236
02237     for(ui_v_surf = pclBPNode->uiBottom; ui_v_surf <= pclBPNode->uiTop; ++ui_v_surf)
02238     {
02239 //      std::cerr << "v:" << ui_v_surf << std::endl;
02240         const double cd_mul_v  = ( (*cpvdIntervalsV)[ui_v_surf] - (*cpvdIntervalsV)[pclBPNode->uiBottom]) / ( (*cpvdIntervalsV)[pclBPNode->uiTop + 1] - (*cpvdIntervalsV)[pclBPNode->uiBottom]);
02241         const double cd_mul_vn = ( (*cpvdIntervalsV)[ui_v_surf + 1] - (*cpvdIntervalsV)[pclBPNode->uiBottom]) / ( (*cpvdIntervalsV)[pclBPNode->uiTop + 1] - (*cpvdIntervalsV)[pclBPNode->uiBottom]);
02242
02243         for(ui_u_surf = pclBPNode->uiLeft; ui_u_surf <= pclBPNode->uiRight; ++ui_u_surf)
02244         {
02245 //          std::cerr << "u:" << ui_u_surf << std::endl;
02246             pcl_surface = &( (*pvvclPatches)[ui_u_surf][ui_v_surf]);
02247
02248             const std::vector<std::vector <Vec3d> > &crvvcl_cps = pcl_surface->getControlPointMatrix();
02249             const unsigned int                       cui_m      = crvvcl_cps.size() - 1;
02250             const unsigned int                       cui_n      = crvvcl_cps[0].size() - 1;
02251             const double                             cd_mul_u   = ( (*cpvdIntervalsU)[ui_u_surf] - (*cpvdIntervalsU)[pclBPNode->uiLeft]) / ( (*cpvdIntervalsU)[pclBPNode->uiRight + 1] - (*cpvdIntervalsU)[pclBPNode->uiLeft]);
02252             const double                             cd_mul_un  = ( (*cpvdIntervalsU)[ui_u_surf + 1] - (*cpvdIntervalsU)[pclBPNode->uiLeft]) / ( (*cpvdIntervalsU)[pclBPNode->uiRight + 1] - (*cpvdIntervalsU)[pclBPNode->uiLeft]);
02253             const Vec3d                              ccl_b00    = ccl_a00 + ccl_a0n * cd_mul_v + (ccl_am0 + ccl_amn * cd_mul_v) * cd_mul_u;
02254             const Vec3d                              ccl_b0n    = ccl_a00 + ccl_a0n * cd_mul_vn + (ccl_am0 + ccl_amn * cd_mul_vn) * cd_mul_u;
02255             const Vec3d                              ccl_bm0    = ccl_a00 + ccl_a0n * cd_mul_v + (ccl_am0 + ccl_amn * cd_mul_v) * cd_mul_un;
02256             const Vec3d                              ccl_bmn    = ccl_a00 + ccl_a0n * cd_mul_vn + (ccl_am0 + ccl_amn * cd_mul_vn) * cd_mul_un;
02257             const double                             cd_rn      = 1.0 / cui_n;
02258             const double                             cd_rm      = 1.0 / cui_m;
02259             const Vec3d                              ccl_dcx0   = (ccl_bm0 - ccl_b00) * cd_rm;
02260             const Vec3d                              ccl_dcxn   = (ccl_bmn - ccl_b0n) * cd_rm;
02261             const Vec3d                              ccl_dc0x   = (ccl_b0n - ccl_b00) * cd_rn;
02262             const Vec3d                              ccl_dcmx   = (ccl_bmn - ccl_bm0) * cd_rn;
02263
02264 //          std::cerr << cui_n << "," << cui_m << std::endl;
02265 //          std::cerr << ( void* ) pclBPNode << std::endl;
02266
02267             cl_ci0 = ccl_b00;
02268             cl_cin = ccl_b0n;
02269
02270             for(ui_i = 0; ui_i <= cui_m; ++ui_i)
02271             {
02272                 cl_cij = cl_ci0;
02273                 cl_c0j = ccl_b00;
02274                 cl_cmj = ccl_bm0;
02275
02276                 const Vec3d ccl_dcix = (cl_cin - cl_ci0) * cd_rn;
02277
02278                 cl_uv[0] = ui_i * cd_rm;
02279
02280                 for(ui_j = 0; ui_j <= cui_n; ++ui_j)
02281                 {
02282                     const Vec3d ccl_dcxj = cl_cmj - cl_c0j;
02283
02284                     cl_uv[1] = ui_j * cd_rn;
02285 #ifdef OSG_DIFF_TO_BILIN
02286 //                  cl_bij = crvvcl_cps[ ui_i ][ ui_j ] - cl_cij;
02287                     cl_bij = pcl_surface->computewdeCasteljau(cl_uv, i_err) - cl_cij;
02288 //                  std::cerr << ccl_dcix << "x" << ccl_dcxj << std::endl;
02289 /*                  cl_norm = ccl_dcix.cross( ccl_dcxj );
02290                     d_quad_size = cl_norm.squareLength( );
02291                     if( d_quad_size > 0.0 )
02292                     {
02293 //                      std::cerr << "qs:" << d_quad_size << std::endl;
02294                         cl_norm *= 1.0 / sqrt( d_quad_size );
02295                         d_quad_size = osgAbs( cl_bij.dot( cl_norm ) );
02296                         if( d_quad_size > pclBPNode->fError )
02297                         {
02298                             pclBPNode->fError = ( float ) d_quad_size;
02299                         }
02300                     }
02301                     else*/
02302                     {
02303                         d_quad_size = sqrt(cl_bij.squareLength() );
02304 //                      std::cerr << "qs:" << d_quad_size << std::endl;
02305                         if(d_quad_size > pclBPNode->fError)
02306                         {
02307                             pclBPNode->fError = ( float ) d_quad_size;
02308                         }
02309                     }
02310 #else
02311 //                  std::cerr << cl_uv << std::endl;
02312                     cl_bij = pcl_surface->computewdeCasteljau(cl_uv, i_err);
02313 //                  std::cerr << i_err << std::endl;
02314                     if(cl_uv[0] + cl_uv[1] <= 1.0)
02315                     {
02316                         d_quad_size = computeDistToTriangle(cl_bij, ccl_a00, ccl_s0n, ccl_sm0);
02317                     }
02318                     else
02319                     {
02320                         d_quad_size = computeDistToTriangle(cl_bij, ccl_smn, ccl_s0n, ccl_sm0);
02321                     }
02322 //                  d_quad_size = osgMin( computeDistToTriangle( cl_bij, ccl_a00, ccl_s0n, ccl_sm0 ),
02323 //                                     computeDistToTriangle( cl_bij, ccl_smn, ccl_s0n, ccl_sm0 ) );
02324                     if(d_quad_size > pclBPNode->fError)
02325                     {
02326                         pclBPNode->fError = ( float ) d_quad_size;
02327                     }
02328                     if(cl_uv[0] <= cl_uv[1])
02329                     {
02330                         d_quad_size = computeDistToTriangle(cl_bij, ccl_s0n, ccl_a00, ccl_smn);
02331                     }
02332                     else
02333                     {
02334                         d_quad_size = computeDistToTriangle(cl_bij, ccl_sm0, ccl_a00, ccl_smn);
02335                     }
02336 //                  d_quad_size = osgMin( computeDistToTriangle( cl_bij, ccl_s0n, ccl_a00, ccl_smn ),
02337 //                                     computeDistToTriangle( cl_bij, ccl_sm0, ccl_a00, ccl_smn ) );
02338                     if(d_quad_size > pclBPNode->fError)
02339                     {
02340                         pclBPNode->fError = ( float ) d_quad_size;
02341                     }
02342 //                  std::cerr << ccl_a00 << ccl_a0n << ccl_am0 << ccl_amn << std::endl << std::endl;
02343 #endif
02344                     cl_cij += ccl_dcix;
02345                     cl_c0j += ccl_dc0x;
02346                     cl_cmj += ccl_dcmx;
02347                 }
02348
02349                 cl_ci0 += ccl_dcx0;
02350                 cl_cin += ccl_dcxn;
02351             }
02352         }
02353     }
02354
02355     // add twist vector error
02356 #ifdef OSG_DIFF_TO_BILIN
02357 //  pclBPNode->fError += ( float ) sqrt( ( ccl_amn - ccl_a0n ).squareLength( ) ) * 0.25;
02358 #endif
02359 //  pclBPNode->fError *= 0.1;
02360 //  pcl_node->ptErrorCell->fError = 0.0;
02361 }
02362
02363
02364 void CErrorQuadTree::SubdivideBPTree(SBPETreeCell                                  * pclBPNode,
02365                                      std::vector<std::vector<BezierTensorSurface> > *pvvclPatches,
02366                                      const std::vector<double>                     * cpvdIntervalsU,
02367                                      const std::vector<double>                     * cpvdIntervalsV)
02368 {
02369     unsigned int ui_diff;
02370
02371     ui_diff = 1;
02372     while( (ui_diff <= pclBPNode->uiRight - pclBPNode->uiLeft) ||
02373            (ui_diff <= pclBPNode->uiTop - pclBPNode->uiBottom) )
02374     {
02375         ui_diff = (ui_diff << 1);
02376     }
02377     ui_diff = ( (ui_diff - 1) >> 1);
02378
02379 /*  std::cerr << "subdivide:" << std::endl;
02380     std::cerr << pclBPNode->uiLeft << ", " << pclBPNode->uiRight << std::endl;
02381     std::cerr << pclBPNode->uiBottom << ", " << pclBPNode->uiTop << std::endl;
02382     std::cerr << ui_diff << std::endl;*/
02383
02384     if(pclBPNode->uiBottom + ui_diff < pclBPNode->uiTop)
02385     {
02386         if(pclBPNode->uiLeft + ui_diff < pclBPNode->uiRight)
02387         {
02388             // xx
02389             // xx
02390             pclBPNode->aptChildren[3]             = new SBPETreeCell;
02391             pclBPNode->aptChildren[3]->uiBottom   = pclBPNode->uiBottom;
02392             pclBPNode->aptChildren[3]->uiTop      = pclBPNode->uiBottom + ui_diff;
02393             pclBPNode->aptChildren[3]->uiLeft     = pclBPNode->uiLeft;
02394             pclBPNode->aptChildren[3]->uiRight    = pclBPNode->uiLeft + ui_diff;
02395             pclBPNode->aptChildren[3]->fPrevError = osgMin(pclBPNode->fError, pclBPNode->fPrevError);
02396             pclBPNode->aptChildren[2]             = new SBPETreeCell;
02397             pclBPNode->aptChildren[2]->uiBottom   = pclBPNode->uiBottom;
02398             pclBPNode->aptChildren[2]->uiTop      = pclBPNode->aptChildren[3]->uiTop;
02399             pclBPNode->aptChildren[2]->uiLeft     = pclBPNode->aptChildren[3]->uiRight + 1;
02400             pclBPNode->aptChildren[2]->uiRight    = pclBPNode->uiRight;
02401             pclBPNode->aptChildren[2]->fPrevError = osgMin(pclBPNode->fError, pclBPNode->fPrevError);
02402             pclBPNode->aptChildren[0]             = new SBPETreeCell;
02403             pclBPNode->aptChildren[0]->uiBottom   = pclBPNode->aptChildren[3]->uiTop + 1;
02404             pclBPNode->aptChildren[0]->uiTop      = pclBPNode->uiTop;
02405             pclBPNode->aptChildren[0]->uiLeft     = pclBPNode->uiLeft;
02406             pclBPNode->aptChildren[0]->uiRight    = pclBPNode->aptChildren[3]->uiRight;
02407             pclBPNode->aptChildren[0]->fPrevError = osgMin(pclBPNode->fError, pclBPNode->fPrevError);
02408             pclBPNode->aptChildren[1]             = new SBPETreeCell;
02409             pclBPNode->aptChildren[1]->uiBottom   = pclBPNode->aptChildren[0]->uiBottom;
02410             pclBPNode->aptChildren[1]->uiTop      = pclBPNode->uiTop;
02411             pclBPNode->aptChildren[1]->uiLeft     = pclBPNode->aptChildren[2]->uiLeft;
02412             pclBPNode->aptChildren[1]->uiRight    = pclBPNode->uiRight;
02413             pclBPNode->aptChildren[1]->fPrevError = osgMin(pclBPNode->fError, pclBPNode->fPrevError);
02414             ComputeBPError(pclBPNode->aptChildren[3], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02415             ComputeBPError(pclBPNode->aptChildren[2], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02416             ComputeBPError(pclBPNode->aptChildren[0], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02417             ComputeBPError(pclBPNode->aptChildren[1], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02418             SubdivideBPTree(pclBPNode->aptChildren[3], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02419             SubdivideBPTree(pclBPNode->aptChildren[2], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02420             SubdivideBPTree(pclBPNode->aptChildren[0], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02421             SubdivideBPTree(pclBPNode->aptChildren[1], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02422         }
02423         else
02424         {
02425             // x
02426             // x
02427             pclBPNode->aptChildren[3]             = new SBPETreeCell;
02428             pclBPNode->aptChildren[3]->uiBottom   = pclBPNode->uiBottom;
02429             pclBPNode->aptChildren[3]->uiTop      = pclBPNode->uiBottom + ui_diff;
02430             pclBPNode->aptChildren[3]->uiLeft     = pclBPNode->uiLeft;
02431             pclBPNode->aptChildren[3]->uiRight    = pclBPNode->uiRight;
02432             pclBPNode->aptChildren[3]->fPrevError = osgMin(pclBPNode->fError, pclBPNode->fPrevError);
02433             pclBPNode->aptChildren[2]             = NULL;
02434             pclBPNode->aptChildren[0]             = new SBPETreeCell;
02435             pclBPNode->aptChildren[0]->uiBottom   = pclBPNode->aptChildren[3]->uiTop + 1;
02436             pclBPNode->aptChildren[0]->uiTop      = pclBPNode->uiTop;
02437             pclBPNode->aptChildren[0]->uiLeft     = pclBPNode->uiLeft;
02438             pclBPNode->aptChildren[0]->uiRight    = pclBPNode->uiRight;
02439             pclBPNode->aptChildren[0]->fPrevError = osgMin(pclBPNode->fError, pclBPNode->fPrevError);
02440             pclBPNode->aptChildren[1]             = NULL;
02441             ComputeBPError(pclBPNode->aptChildren[3], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02442             ComputeBPError(pclBPNode->aptChildren[0], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02443             SubdivideBPTree(pclBPNode->aptChildren[3], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02444             SubdivideBPTree(pclBPNode->aptChildren[0], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02445         }
02446     }
02447     else if(pclBPNode->uiLeft + ui_diff < pclBPNode->uiRight)
02448     {
02449         // xx
02450         pclBPNode->aptChildren[3]             = new SBPETreeCell;
02451         pclBPNode->aptChildren[3]->uiBottom   = pclBPNode->uiBottom;
02452         pclBPNode->aptChildren[3]->uiTop      = pclBPNode->uiTop;
02453         pclBPNode->aptChildren[3]->uiLeft     = pclBPNode->uiLeft;
02454         pclBPNode->aptChildren[3]->uiRight    = pclBPNode->uiLeft + ui_diff;
02455         pclBPNode->aptChildren[3]->fPrevError = osgMin(pclBPNode->fError, pclBPNode->fPrevError);
02456         pclBPNode->aptChildren[2]             = new SBPETreeCell;
02457         pclBPNode->aptChildren[2]->uiBottom   = pclBPNode->uiBottom;
02458         pclBPNode->aptChildren[2]->uiTop      = pclBPNode->uiTop;
02459         pclBPNode->aptChildren[2]->uiLeft     = pclBPNode->aptChildren[3]->uiRight + 1;
02460         pclBPNode->aptChildren[2]->uiRight    = pclBPNode->uiRight;
02461         pclBPNode->aptChildren[2]->fPrevError = osgMin(pclBPNode->fError, pclBPNode->fPrevError);
02462         pclBPNode->aptChildren[0]             = NULL;
02463         pclBPNode->aptChildren[1]             = NULL;
02464         ComputeBPError(pclBPNode->aptChildren[3], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02465         ComputeBPError(pclBPNode->aptChildren[2], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02466         SubdivideBPTree(pclBPNode->aptChildren[3], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02467         SubdivideBPTree(pclBPNode->aptChildren[2], pvvclPatches, cpvdIntervalsU, cpvdIntervalsV);
02468     }
02469     else
02470     {
02471         pclBPNode->aptChildren[0] = NULL;
02472         pclBPNode->aptChildren[1] = NULL;
02473         pclBPNode->aptChildren[2] = NULL;
02474         pclBPNode->aptChildren[3] = NULL;
02475     }
02476 }
02477
02478
02479 void CErrorQuadTree::DeleteBPNode(SBPETreeCell *&rpclNode)
02480 {
02481     if(rpclNode->aptChildren[0])
02482     {
02483         DeleteBPNode(rpclNode->aptChildren[0]);
02484     }
02485     if(rpclNode->aptChildren[1])
02486     {
02487         DeleteBPNode(rpclNode->aptChildren[1]);
02488     }
02489     if(rpclNode->aptChildren[2])
02490     {
02491         DeleteBPNode(rpclNode->aptChildren[2]);
02492     }
02493     if(rpclNode->aptChildren[3])
02494     {
02495         DeleteBPNode(rpclNode->aptChildren[3]);
02496     }
02497     delete rpclNode;
02498     rpclNode = NULL;
02499 }
02500 #endif
02501 
02502
02503 #ifndef OSG_DIFF_TO_BILIN
02504 double CErrorQuadTree::computeDistToPlane(const Vec3d cclP, const Vec3d cclV1, const Vec3d cclV2, const Vec3d cclV3)
02505 {
02506     Vec3d  cl_norm;
02507     double d_quad_size;
02508     Vec3d  cl_diff;
02509
02510     cl_diff     = cclP - cclV1;
02511     cl_norm     = (cclV2 - cclV1).cross(cclV3 - cclV1);
02512     d_quad_size = cl_norm.squareLength();
02513 //  std::cerr << cclP << cclV1 << cclV2 << cclV3 << std::endl;
02514     if(d_quad_size > 0.0)
02515     {
02516 //      std::cerr << "qs:" << d_quad_size << std::endl;
02517         cl_norm *= 1.0 / sqrt(d_quad_size);
02518 //      std::cerr << cl_diff << cl_norm << std::endl;
02519         return osgAbs(cl_diff.dot(cl_norm) );
02520     }
02521     else
02522     {
02523 //      std::cerr << cl_diff << std::endl;
02524         return sqrt(cl_diff.squareLength() );
02525     }
02526 }
02527
02528
02529 double CErrorQuadTree::computeDistToTriangle(const Vec3d cclP, const Vec3d cclV1, const Vec3d cclV2, const Vec3d cclV3)
02530 {
02531     Vec3d       cl_norm;
02532     double      d_dist;
02533     Vec3d       cl_diff;
02534     Vec3d       cl_proj;
02535     double      d_temp;
02536     double      d_temp2;
02537     const Vec3d ccl_d1 = cclV2 - cclV1;
02538     const Vec3d ccl_d2 = cclV3 - cclV1;
02539
02540     cl_diff = cclP - cclV1;
02541     cl_norm = ccl_d1.cross(ccl_d2);
02542     d_temp  = cl_norm.squareLength();
02543     if(d_temp > 0.0)
02544     {
02545         cl_norm *= 1.0 / sqrt(d_temp);
02546         // project point to triangle plane
02547         d_dist  = cl_diff.dot(cl_norm);
02548         cl_proj = cclP - cl_norm * d_dist;
02549 //      std::cerr << ( cl_proj - cclV1 ).dot( cl_norm ) << std::endl;
02550         // check if point is inside triangle
02551         cl_diff = cl_proj - cclV1;
02552         d_temp  = ccl_d1.dot(cl_diff) / ccl_d1.squareLength();
02553         d_temp2 = ccl_d2.dot(cl_diff) / ccl_d2.squareLength();
02554         if( (d_temp < 0.0) || (d_temp2 < 0.0) || (d_temp + d_temp2 > 1.0) )
02555         {
02556 //          std::cerr << d_temp << " " << d_temp2 << ": " << osgAbs( d_dist );
02557             d_dist = computeDistToLine(cclP, cclV1, cclV2);
02558             d_temp = computeDistToLine(cclP, cclV1, cclV3);
02559             if(d_temp < d_dist)
02560                 d_dist = d_temp;
02561             d_temp = computeDistToLine(cclP, cclV2, cclV3);
02562             if(d_temp < d_dist)
02563                 d_dist = d_temp;
02564 //          std::cerr << " -> " << d_dist << std::endl;
02565         }
02566         return osgAbs(d_dist);
02567     }
02568     else
02569     {
02570         return sqrt(cl_diff.squareLength() );
02571     }
02572 }
02573 #endif
02574 
02575 double CErrorQuadTree::computeDistToLine(const Vec3d cclP, const Vec3d cclV1, const Vec3d cclV2)
02576 {
02577     const Vec3d ccl_d    = cclV2 - cclV1;
02578     const Vec3d ccl_diff = cclP - cclV1;
02579     double      d_temp;
02580     double      d_temp2;
02581
02582     d_temp = ccl_d.dot(ccl_diff);
02583     if(d_temp <= 0.0)
02584     {
02585 //      std::cerr << "a " << sqrt( ccl_diff.squareLength( ) ) << std::endl;
02586         return sqrt(ccl_diff.squareLength() );
02587     }
02588     d_temp2 = ccl_d.squareLength();
02589     if(d_temp2 <= d_temp)
02590     {
02591 //      std::cerr << "b " << sqrt( ( cclP - cclV2 ).squareLength( ) ) << std::endl;
02592         return sqrt( (cclP - cclV2).squareLength() );
02593     }
02594 //  std::cerr << "c " << sqrt( ( ccl_d * ( d_temp / d_temp2 ) - ccl_diff ).squareLength( ) ) << std::endl;
02595     return sqrt( (ccl_d * (d_temp / d_temp2) - ccl_diff).squareLength() );
02596 }
02597
02598 void CErrorQuadTree::WriteTree(std::ostream &rclFile)
02599 {
02600 #ifdef OSG_ARBITRARY_SPLIT
02601     SErrorTreeCell              *pt_act;
02602     std::vector<SErrorTreeCell*> vpt_stack;
02603
02604     // store max error
02605     rclFile << m_fMaxError << std::endl;
02606
02607     // store error cells
02608     vpt_stack.push_back(m_ptRoot);
02609     while(vpt_stack.size() )
02610     {
02611         pt_act = vpt_stack[vpt_stack.size() - 1];
02612         vpt_stack.pop_back();
02613         rclFile << pt_act->fError << " " << pt_act->fSplitValue << " ";
02614         rclFile << ( (pt_act->ptChildren != NULL) ? true : false) << std::endl;
02615         if(pt_act->ptChildren != NULL)
02616         {
02617             vpt_stack.push_back(&(pt_act->ptChildren[1]) );
02618             vpt_stack.push_back(&(pt_act->ptChildren[0]) );
02619         }
02620     }
02621 #endif
02622 }
02623
02624 void CErrorQuadTree::ReadTree(std::istream &rclFile)
02625 {
02626 #ifdef OSG_ARBITRARY_SPLIT
02627     SErrorTreeCell              *pt_act;
02628     std::vector<SErrorTreeCell*> vpt_stack;
02629     bool                         b_has_children;
02630
02631     // load max error
02632     rclFile >> m_fMaxError >> std::ws;
02633
02634     // load error cells
02635     m_ptRoot = new SErrorTreeCell;
02636     vpt_stack.push_back(m_ptRoot);
02637     while(vpt_stack.size() )
02638     {
02639         pt_act = vpt_stack[vpt_stack.size() - 1];
02640         vpt_stack.pop_back();
02641         rclFile >> pt_act->fError >> std::ws >> pt_act->fSplitValue >> std::ws;
02642         rclFile >> b_has_children >> std::ws;
02643         if(b_has_children)
02644         {
02645             pt_act->ptChildren = new SErrorTreeCell[2];
02646             vpt_stack.push_back(&(pt_act->ptChildren[1]) );
02647             vpt_stack.push_back(&(pt_act->ptChildren[0]) );
02648         }
02649         else
02650         {
02651             pt_act->ptChildren = NULL;
02652         }
02653     }
02654 #endif
02655 }