The Parametric Pseudo-Manifold (PPS) Library 1.0
|
00001 00027 #include "ppsfrompnt.h" // PPSfromPNT 00028 00029 #include <cassert> // assert 00030 #include <cmath> // fabs, sqrt 00031 00032 00047 namespace ppsfrompnt { 00048 00049 using pps::PPS ; 00050 00058 PPSfromPNT::PPSfromPNT( Surface* mesh ) : PPS< Surface >( mesh ) 00059 { 00060 // 00061 // Sets the owner of each half-edge attribute set. 00062 // 00063 for ( EdgeIterator eit = edges_begin() ; !is_done( eit ) ; 00064 move_forward( eit ) ) { 00065 Edge* edge = get_edge( eit ) ; 00066 00067 Halfedge* h1 = edge->get_first_halfedge() ; 00068 Halfedge* h2 = edge->get_second_halfedge() ; 00069 00070 h1->get_attributes().set_owner( h1 ) ; 00071 h2->get_attributes().set_owner( h2 ) ; 00072 } 00073 00074 // 00075 // Create the PN triangle surface. 00076 // 00077 build_pnt_surface() ; 00078 } 00079 00080 00095 void 00096 PPSfromPNT::eval_surface( 00097 Face* face , 00098 double u , 00099 double v , 00100 double w , 00101 double& x , 00102 double& y , 00103 double& z 00104 ) 00105 const 00106 { 00107 // 00108 // Make sure the face has a PN triangle associated with it. 00109 // 00110 PNTriangle* patch = face->get_attributes().get_patch() ; 00111 00112 assert( patch != 0 ) ; 00113 00114 // 00115 // Make sure the coordinates are non-negative and add up to 1. 00116 // 00117 assert( 00118 ( u >= 0 ) && ( u <= 1 ) && ( v >= 0 ) && 00119 ( v <= 1 ) && ( w >= 0 ) && ( w <= 1 ) && 00120 ( fabs( 1 - ( u + v + w ) ) <= 1e-15 ) 00121 ) ; 00122 // 00123 // Evaluate the patch. 00124 // 00125 patch->point( u , v , x , y , z ) ; 00126 } 00127 00128 00136 void 00137 PPSfromPNT::build_pnt_surface() 00138 { 00139 // 00140 // For each face of the given triangle mesh, compute a PN triangle. 00141 // 00142 for ( FaceIterator fit = faces_begin() ; !is_done( fit ) ; 00143 move_forward( fit ) ) { 00144 Face* face = get_face( fit ) ; 00145 00146 Halfedge* h1 = face->get_halfedge() ; 00147 Halfedge* h2 = h1->get_next() ; 00148 Halfedge* h3 = h1->get_prev() ; 00149 00150 Vertex* v1 = h1->get_origin() ; 00151 Vertex* v2 = h2->get_origin() ; 00152 Vertex* v3 = h3->get_origin() ; 00153 00154 double p1[ 3 ] ; 00155 double p2[ 3 ] ; 00156 double p3[ 3 ] ; 00157 00158 p1[ 0 ] = v1->x() ; 00159 p1[ 1 ] = v1->y() ; 00160 p1[ 2 ] = v1->z() ; 00161 00162 p2[ 0 ] = v2->x() ; 00163 p2[ 1 ] = v2->y() ; 00164 p2[ 2 ] = v2->z() ; 00165 00166 p3[ 0 ] = v3->x() ; 00167 p3[ 1 ] = v3->y() ; 00168 p3[ 2 ] = v3->z() ; 00169 00170 double n1[ 3 ] ; 00171 double n2[ 3 ] ; 00172 double n3[ 3 ] ; 00173 00174 compute_vertex_normal_vector( v1 , n1[ 0 ] , n1[ 1 ] , n1[ 2 ] ) ; 00175 compute_vertex_normal_vector( v2 , n2[ 0 ] , n2[ 1 ] , n2[ 2 ] ) ; 00176 compute_vertex_normal_vector( v3 , n3[ 0 ] , n3[ 1 ] , n3[ 2 ] ) ; 00177 00178 // 00179 // Assign the PN triangle to the face 00180 // 00181 face->get_attributes().set_patch( 00182 new PNTriangle( 00183 p1 , 00184 p2 , 00185 p3 , 00186 n1 , 00187 n2 , 00188 n3 00189 ) 00190 ) ; 00191 } 00192 00193 return ; 00194 } 00195 00196 00207 void 00208 PPSfromPNT::compute_vertex_normal_vector( 00209 Vertex* vertex , 00210 double& x , 00211 double& y , 00212 double& z 00213 ) 00214 { 00215 unsigned i = 0 ; 00216 00217 x = 0 ; 00218 y = 0 ; 00219 z = 0 ; 00220 00221 Halfedge* h1 = vertex->get_halfedge() ; 00222 Halfedge* h2 = h1 ; 00223 00224 do { 00225 double xaux ; 00226 double yaux ; 00227 double zaux ; 00228 00229 compute_face_normal_vector( h2->get_face() , xaux , yaux , zaux ) ; 00230 00231 x += xaux ; 00232 y += yaux ; 00233 z += zaux ; 00234 00235 ++i ; 00236 00237 h2 = h2->get_prev()->get_mate() ; 00238 } 00239 while ( h2 != h1 ) ; 00240 00241 x /= double( i ) ; 00242 y /= double( i ) ; 00243 z /= double( i ) ; 00244 00245 double length = sqrt( ( x * x ) + ( y * y ) + ( z * z ) ) ; 00246 00247 x /= length ; 00248 y /= length ; 00249 z /= length ; 00250 00251 return ; 00252 } 00253 00254 00265 void 00266 PPSfromPNT::compute_face_normal_vector( 00267 Face* face , 00268 double& x , 00269 double& y , 00270 double& z 00271 ) 00272 { 00273 Halfedge* h1 = face->get_halfedge() ; 00274 Halfedge* h2 = h1->get_next() ; 00275 Halfedge* h3 = h1->get_prev() ; 00276 00277 double x1 = h1->get_origin()->x() ; 00278 double y1 = h1->get_origin()->y() ; 00279 double z1 = h1->get_origin()->z() ; 00280 00281 double x2 = h2->get_origin()->x() ; 00282 double y2 = h2->get_origin()->y() ; 00283 double z2 = h2->get_origin()->z() ; 00284 00285 double x3 = h3->get_origin()->x() ; 00286 double y3 = h3->get_origin()->y() ; 00287 double z3 = h3->get_origin()->z() ; 00288 00289 x2 -= x1 ; 00290 y2 -= y1 ; 00291 z2 -= z1 ; 00292 00293 x3 -= x1 ; 00294 y3 -= y1 ; 00295 z3 -= z1 ; 00296 00297 x = y2 * z3 - z2 * y3 ; 00298 y = z2 * x3 - x2 * z3 ; 00299 z = x2 * y3 - y2 * x3 ; 00300 00301 double length = sqrt( ( x * x ) + ( y * y ) + ( z * z ) ) ; 00302 00303 x /= length ; 00304 y /= length ; 00305 z /= length ; 00306 00307 return ; 00308 } 00309 00310 } 00311 //end of group class.