fixed: auto_ptr -> unique_ptr
[opensg.git] / Source / System / NodeCores / Drawables / Nurbs / Internal / OSGQuadTreeCreator.cpp
bloba663259f73416b8b107861575d13d218ffea4333
1 /*---------------------------------------------------------------------------*\
2 * OpenSG NURBS Library *
3 * *
4 * *
5 * Copyright (C) 2001-2006 by the University of Bonn, Computer Graphics Group*
6 * *
7 * http://cg.cs.uni-bonn.de/ *
8 * *
9 * contact: edhellon@cs.uni-bonn.de, guthe@cs.uni-bonn.de, rk@cs.uni-bonn.de *
10 * *
11 \*---------------------------------------------------------------------------*/
12 /*---------------------------------------------------------------------------*\
13 * License *
14 * *
15 * This library is free software; you can redistribute it and/or modify it *
16 * under the terms of the GNU Library General Public License as published *
17 * by the Free Software Foundation, version 2. *
18 * *
19 * This library is distributed in the hope that it will be useful, but *
20 * WITHOUT ANY WARRANTY; without even the implied warranty of *
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU *
22 * Library General Public License for more details. *
23 * *
24 * You should have received a copy of the GNU Library General Public *
25 * License along with this library; if not, write to the Free Software *
26 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. *
27 * *
28 \*---------------------------------------------------------------------------*/
29 /*---------------------------------------------------------------------------*\
30 * Changes *
31 * *
32 * *
33 * *
34 * *
35 * *
36 * *
37 \*---------------------------------------------------------------------------*/
38 #include "OSGQuadTreeCreator.h"
39 #include <stdio.h>
41 OSG_USING_NAMESPACE
43 #ifdef WIN32
44 #pragma warning (disable : 985)
45 #endif
47 #ifdef _DEBUG
48 #ifdef WIN32
49 #undef THIS_FILE
50 static char THIS_FILE[] = __FILE__;
51 #endif
52 #endif
54 double QuadTreeCreator::
55 computeApproximationError(DCTPFace *f)
57 //the approximate error comes from to error,
58 //I: interpolation with a bilinear surface
59 //II: interpolating this bilinear surface with two triangles
61 //first gotta get needed controlpoints from the index-ed quadtreeleaf
62 BezierTensorSurface *btsp =
63 static_cast<BezierTensorSurface*>(f->faceinfo);
65 const DCTPVec4dmatrix & cps = btsp->getControlPointMatrix();
66 const DCTPVec4dmatrix::size_type m = cps.size() - 1;
67 //size of control point matrix in x direction
68 const DCTPVec4dvector::size_type n = cps[0].size() - 1;
69 //size of cpm. in y direction
70 Vec3d b00;
71 Vec3d b0n;
72 Vec3d bm0;
73 Vec3d bmn;
74 b00[0] = cps[0][0][0] / cps[0][0][3];
75 b00[1] = cps[0][0][1] / cps[0][0][3];
76 b00[2] = cps[0][0][2] / cps[0][0][3];
77 b0n[0] = (cps[0][n][0] / cps[0][n][3]);
78 b0n[1] = (cps[0][n][1] / cps[0][n][3]);
79 b0n[2] = (cps[0][n][2] / cps[0][n][3]);
80 bm0[0] = (cps[m][0][0] / cps[m][0][3]);
81 bm0[1] = (cps[m][0][1] / cps[m][0][3]);
82 bm0[2] = (cps[m][0][2] / cps[m][0][3]);
83 bmn[0] = (cps[m][n][0] / cps[m][n][3]);
84 bmn[1] = (cps[m][n][1] / cps[m][n][3]);
85 bmn[2] = (cps[m][n][2] / cps[m][n][3]);
86 //now computin' I. error:
87 //we hafta decide a maxmimum value:
88 double max_quad = -1.0;
89 Vec3d cij, bij;
90 double quad_size;
91 DCTPVec3dmatrix::size_type i;
92 DCTPVec3dvector::size_type j;
93 const double rn = 1.0 / n;
94 const double rm = 1.0 / m;
96 Vec3d ci0 = b00;
97 Vec3d cin = b0n;
98 const Vec3d dc0 = (bm0 - b00) * rm;
99 const Vec3d dcn = (bmn - b0n) * rm;
100 // Vec2d cl_uv;
101 // int i_err = 0;
103 // if( m_bForTrimming )
105 // cl_uv[0] = 0.0;
106 for(i = 0; i <= m; ++i)
108 // cl_uv[1] = 0.0;
109 cij = ci0;
110 const Vec3d dci = (cin - ci0) * rn;
112 for(j = 0; j <= n; ++j)
114 if(cps[i][j][3] < DCTP_EPS)
116 // the weight is zero so we return a practically
117 // infinite error enforcing a subdivision of the patch
118 quad_size = 1.0e300;
120 else
122 bij[0] = (cps[i][j][0] / cps[i][j][3]) - cij[0];
123 bij[1] = (cps[i][j][1] / cps[i][j][3]) - cij[1];
124 bij[2] = (cps[i][j][2] / cps[i][j][3]) - cij[2];
125 // bij = btsp->computewdeCasteljau( cl_uv, i_err ) - cij;
126 quad_size = bij.squareLength();
128 if(quad_size > max_quad)
129 max_quad = quad_size;
130 cij += dci;
131 // cl_uv[1] += rn;
132 } //end for
134 ci0 += dc0;
135 cin += dcn;
136 // cl_uv[0] += rm;
139 const double error_I = sqrt(max_quad);
141 //now get error2
142 const Vec3d linerr_vect = (b00 - b0n + bmn - bm0);
143 const double error_II = sqrt(linerr_vect.squareLength() ) * 0.25;
144 //now the error iz the sum o' I and II:
145 return error_I + error_II;
147 /* else
149 Vec3d cl_norm;
150 const Vec3d cdc0 = ( b0n - b00 ) * rn;
151 const Vec3d cdcm = ( bmn - bm0 ) * rn;
152 Vec3d c0j;
153 Vec3d cmj;
155 for( i = 0; i <= m; ++i )
157 cij = ci0;
158 c0j = b00;
159 cmj = bm0;
161 const Vec3d dci = ( cin - ci0 ) * rn;
163 // std::cerr << dci << std::endl;
165 for( j = 0; j <= n; ++j )
167 const Vec3d dcj = cmj - c0j;
169 bij = cps[ i ][ j ] - cij;
170 cl_norm = dci.cross( dcj );
171 // std::cerr << dci << " x " << dcj << std::endl;
172 // std::cerr << cl_norm << "," << bij << std::endl;
173 quad_size = cl_norm.squareLength( );
174 if( quad_size > 0.0 )
176 cl_norm *= 1.0 / sqrt( quad_size );
177 quad_size = osgabs( bij.dot( cl_norm ) );
178 // quad_size = bij.squareLength();
179 if ( quad_size > max_quad ) max_quad = quad_size;
181 cij += dci;
182 c0j += cdc0;
183 cmj += cdcm;
184 } //end for
185 ci0 += dc0;
186 cin += dcn;
188 // std::cerr << max_quad << std::endl;
189 return max_quad;
194 int QuadTreeCreator::setInitialLeaves(
195 const beziersurfacematrix& patches,
196 const DCTPdvector& intervals_u,
197 const DCTPdvector& intervals_v)
199 //sets up the basic leaves array of the quadtree
200 //FIXME: no size & sanity checks made to input, they may be wrong
201 beziersurfacematrix::size_type u_size = patches.size();
202 beziersurfacevector::size_type v_size = patches[0].size();
203 Vec3d a,b,c,d;
205 for(beziersurfacematrix::size_type u = 0; u < u_size; ++u)
207 for(beziersurfacevector::size_type v = 0; v < v_size; ++v)
209 a[0] = intervals_u[u]; a[1] = intervals_v[v + 1]; a[2] = 0.0;
210 b[0] = intervals_u[u + 1]; b[1] = intervals_v[v + 1]; b[2] = 0.0;
211 c[0] = intervals_u[u + 1]; c[1] = intervals_v[v]; c[2] = 0.0;
212 d[0] = intervals_u[u]; d[1] = intervals_v[v]; d[2] = 0.0;
214 // std::cerr << a << b << c << d << std::endl;
216 DCTPFace *face = qtm->AddQuad(a, b, c, d, 0.0);
217 face->faceinfo = new BezierTensorSurface;
218 #ifdef OSG_ADAPTIVE_QUAD_TREE
219 face->norm = -1.0;
220 #endif
221 *(static_cast<BezierTensorSurface*>(face->faceinfo)) = patches[u][v];
222 BezierTensorSurface temp_surface = patches[u][v];
223 } //end inner for
227 return 0;
230 int QuadTreeCreator::finishSubdivisions(DCTPFace *f)
232 /* BezierTensorSurface *btsp = (BezierTensorSurface*)f->faceinfo;
233 beziersurfacematrix created_surfaces( 0 );
234 if ( btsp->midPointSubDivision( created_surfaces ) ) {
235 std::cerr << "QuadTreeCreator::finisSubdivisions:" <<
236 "Cannot midpointsubdivide BezierSurf'!" << std::endl;
237 return -1;
239 dctpfacevector::size_type idx = qtm->faces.size() - 3;
240 (qtm->faces[ idx ])->faceinfo = new BezierTensorSurface;
241 *((BezierTensorSurface*)(qtm->faces[ idx++ ])->faceinfo) =
242 created_surfaces[ 1 ][ 1 ];
243 (qtm->faces[ idx ])->faceinfo = new BezierTensorSurface;
244 *((BezierTensorSurface*)(qtm->faces[ idx++ ])->faceinfo) =
245 created_surfaces[ 1 ][ 0 ];
246 (qtm->faces[ idx ])->faceinfo = new BezierTensorSurface;
247 *((BezierTensorSurface*)(qtm->faces[ idx++ ])->faceinfo) =
248 created_surfaces[ 0 ][ 0 ];
250 delete btsp;
251 f->faceinfo = new BezierTensorSurface;
252 *((BezierTensorSurface*)f->faceinfo) = created_surfaces[ 0 ][ 1 ];
253 return 0;*/
255 BezierTensorSurface *btsp = static_cast<BezierTensorSurface*>(f->faceinfo);
256 beziersurfacevector created_surfaces;
257 if(btsp->midPointSubDivision(created_surfaces) )
259 std::cerr << "QuadTreeCreator::finisSubdivisions:" <<
260 "Cannot midpointsubdivide BezierSurf'!" << std::endl;
261 return -1;
263 dctpfacevector::size_type idx = qtm->faces.size() - 3;
265 (qtm->faces[idx])->faceinfo = new BezierTensorSurface;
266 *(static_cast<BezierTensorSurface*>((qtm->faces[idx++])->faceinfo)) =
267 created_surfaces[2];
269 /* std::cerr << std::endl;
270 std::cerr << ( qtm->faces[ idx ] )->orig_face[ 0 ]->coords << std::endl;
271 std::cerr << ( qtm->faces[ idx ] )->orig_face[ 1 ]->coords << std::endl;
272 std::cerr << ( qtm->faces[ idx ] )->orig_face[ 2 ]->coords << std::endl;
273 std::cerr << ( qtm->faces[ idx ] )->orig_face[ 3 ]->coords << std::endl;
274 std::cerr << created_surfaces[ 1 ].getControlPointMatrix( )[ 0 ][ 0 ] << std::endl;*/
276 (qtm->faces[idx])->faceinfo = new BezierTensorSurface;
277 *(static_cast<BezierTensorSurface*>((qtm->faces[idx++])->faceinfo)) =
278 created_surfaces[1];
280 (qtm->faces[idx])->faceinfo = btsp;
282 f->faceinfo = new BezierTensorSurface;
283 *(static_cast<BezierTensorSurface*>(f->faceinfo)) = created_surfaces[0];
284 return 0;
288 int QuadTreeCreator::createQuadTree(void)
290 //now this is the essentials of this whole class
291 //all other are 'eye-candys'
293 for(dctpfacevector::size_type f = 0;
294 f < qtm->faces.size(); ++f)
296 DCTPFace *face = qtm->faces[f];
297 while(computeApproximationError(face) > error_epsilon)
299 qtm->SubdivideQuad(face);
300 if(finishSubdivisions(face) )
301 return -1;
303 face->norm = error_epsilon / computeBilinearNorm(face);
306 return 0; //all OK, whew...
309 double QuadTreeCreator::computeBilinearNorm(DCTPFace *face)
311 BezierTensorSurface *btsp =
312 static_cast<BezierTensorSurface*>(face->faceinfo);
313 DCTPVec4dmatrix& cps = btsp->getControlPointMatrix();
314 DCTPVec4dmatrix::size_type m = cps.size() - 1;
315 DCTPVec4dvector::size_type n = cps[0].size() - 1;
316 Vec3d b00;
317 Vec3d b0n;
318 Vec3d bm0;
319 Vec3d bmn;
320 b00[0] = cps[0][0][0] / cps[0][0][3];
321 b00[1] = cps[0][0][1] / cps[0][0][3];
322 b00[2] = cps[0][0][2] / cps[0][0][3];
323 b0n[0] = (cps[0][n][0] / cps[0][n][3]) - b00[0];
324 b0n[1] = (cps[0][n][1] / cps[0][n][3]) - b00[1];
325 b0n[2] = (cps[0][n][2] / cps[0][n][3]) - b00[2];
326 bm0[0] = (cps[m][0][0] / cps[m][0][3]) - b00[0];
327 bm0[1] = (cps[m][0][1] / cps[m][0][3]) - b00[1];
328 bm0[2] = (cps[m][0][2] / cps[m][0][3]) - b00[2];
329 bmn[0] = (cps[m][n][0] / cps[m][n][3]) - b00[0];
330 bmn[1] = (cps[m][n][1] / cps[m][n][3]) - b00[1];
331 bmn[2] = (cps[m][n][2] / cps[m][n][3]) - b00[2];
332 double tmp1,tmp2;
333 tmp1 = sqrt(b0n.squareLength() );
334 tmp2 = sqrt(bm0.squareLength() );
335 if(tmp1 < tmp2)
336 tmp1 = tmp2;
337 tmp2 = sqrt(bmn.squareLength() ) / sqrt(2.f); //1.414213; //divided by length
338 if(tmp1 < tmp2)
339 return tmp2;
340 else
341 return tmp1;
344 void QuadTreeCreator::resetMesh(void)
346 dctpfacevector::iterator fe = qtm->faces.end();
348 for(dctpfacevector::iterator f = qtm->faces.begin(); f != fe; ++f)
349 delete static_cast<BezierTensorSurface*>((*f)->faceinfo);
351 qtm->reinit();