contrib/brl/bbas/imesh/imesh_detection.cxx
Go to the documentation of this file.
00001 // This is brl/bbas/imesh/imesh_detection.cxx
00002 #include "imesh_detection.h"
00003 //:
00004 // \file
00005 
00006 #include <vcl_cassert.h>
00007 
00008 //: Return the indices of half edges that are on the mesh boundary
00009 // The results are organized into loops
00010 vcl_vector<vcl_vector<unsigned int> >
00011 imesh_detect_boundary_loops(const imesh_half_edge_set& half_edges)
00012 {
00013   vcl_vector<bool> visited(half_edges.size(),false);
00014   vcl_vector<vcl_vector<unsigned int> > loops;
00015   for (unsigned int i=0; i<half_edges.size(); ++i) {
00016     if (visited[i])
00017       continue;
00018     if (half_edges[i].is_boundary()) {
00019       visited[i] = true;
00020       imesh_half_edge_set::f_const_iterator end(i,half_edges);
00021       imesh_half_edge_set::f_const_iterator itr(end);
00022       vcl_vector<unsigned int> loop(1,end->half_edge_index());
00023       ++itr;
00024       for (; itr != end; ++itr) {
00025         visited[itr->half_edge_index()] = true;
00026         loop.push_back(itr->half_edge_index());
00027       }
00028       loops.push_back(loop);
00029     }
00030   }
00031   return loops;
00032 }
00033 
00034 
00035 //: Trace half edges that have been selected into loops
00036 //  \returns true if all half edges form loops
00037 //  The loops are returned in \param loops as vectors of half edge indices
00038 bool
00039 imesh_trace_half_edge_loops(const imesh_half_edge_set& half_edges,
00040                             const vcl_vector<bool>& flags,
00041                             vcl_vector<vcl_vector<unsigned int> >& loops)
00042 {
00043   loops.clear();
00044   vcl_vector<bool> visited(half_edges.size(),false);
00045   for (unsigned int i=0; i<half_edges.size(); ++i) {
00046     if (visited[i] || !flags[i])
00047       continue;
00048 
00049     vcl_vector<unsigned int> loop;
00050     loop.push_back(i);
00051     visited[i] = true;
00052     imesh_half_edge_set::f_const_iterator end(i,half_edges);
00053     imesh_half_edge_set::f_const_iterator itr(end);
00054     ++itr;
00055     for (; itr!=end; ++itr) {
00056       unsigned int j = itr->half_edge_index();
00057       if (!flags[j]) {
00058         imesh_half_edge_set::v_const_iterator vend(itr);
00059         imesh_half_edge_set::v_const_iterator vitr(vend);
00060         ++vitr;
00061         for (; vitr!=vend && !flags[vitr->half_edge_index()]; ++vitr) /*nothing*/;
00062         if (vitr == vend)
00063           return false; // loop reached a dead end
00064 
00065         itr = imesh_half_edge_set::f_const_iterator(vitr);
00066         j = itr->half_edge_index();
00067       }
00068       if (visited[j]) {
00069         if (j == end->half_edge_index())
00070           break; // loop completed
00071         return false; // loop joins a previous loop
00072       }
00073 
00074       visited[j] = true;
00075       loop.push_back(j);
00076     }
00077     loops.push_back(loop);
00078   }
00079   return true;
00080 }
00081 
00082 
00083 //: Return the indices of contour half edges as seen from direction \param dir
00084 //  The results are organized into loops
00085 vcl_vector<vcl_vector<unsigned int> >
00086 imesh_detect_contour_generator(const imesh_mesh& mesh, const vgl_vector_3d<double>& dir)
00087 {
00088   vcl_vector<bool> contours = imesh_detect_contours(mesh,dir);
00089 
00090   vcl_vector<vcl_vector<unsigned int> > loops;
00091   bool valid = imesh_trace_half_edge_loops(mesh.half_edges(),contours,loops);
00092   assert(valid);
00093 
00094   return loops;
00095 }
00096 
00097 
00098 //: Mark contour half edges as seen from direction \param dir
00099 //  For each contour edge the half edge with the face normal opposing dir is marked
00100 vcl_vector<bool>
00101 imesh_detect_contours(const imesh_mesh& mesh, vgl_vector_3d<double> dir)
00102 {
00103   assert(mesh.has_half_edges());
00104   assert(mesh.faces().has_normals());
00105 
00106   normalize(dir);
00107 
00108   const imesh_half_edge_set& half_edges = mesh.half_edges();
00109   const vcl_vector<vgl_vector_3d<double> >& normals = mesh.faces().normals();
00110 
00111   vcl_vector<double> fdotd(normals.size(),0.0);
00112   for (unsigned int i=0; i<normals.size(); ++i)
00113   {
00114     fdotd[i] = dot_product(normalized(normals[i]),dir);
00115   }
00116 
00117   const unsigned int num_edges = half_edges.size()/2;
00118   vcl_vector<bool> is_contour(half_edges.size(),false);
00119   for (unsigned int i=0; i<num_edges; ++i)
00120   {
00121     const imesh_half_edge& he = half_edges[i*2];
00122     const imesh_half_edge& ohe = half_edges[he.pair_index()];
00123     double v1=0.0, v2=0.0;
00124     if (!he.is_boundary())
00125       v1 = fdotd[he.face_index()];
00126     if (!ohe.is_boundary())
00127       v2 = fdotd[ohe.face_index()];
00128     double prod = v1*v2;
00129     if (prod<=0.0) {
00130       if (v1<0.0) is_contour[i*2] = true;
00131       if (v2<0.0) is_contour[i*2+1] = true;
00132     }
00133   }
00134   return is_contour;
00135 }
00136 
00137 
00138 //: Return the indices of contour half edges as seen from center of projection \param pt
00139 //  The results are organized into loops
00140 vcl_vector<vcl_vector<unsigned int> >
00141 imesh_detect_contour_generator(const imesh_mesh& mesh, const vgl_point_3d<double>& pt)
00142 {
00143   vcl_vector<bool> contours = imesh_detect_contours(mesh,pt);
00144 
00145   vcl_vector<vcl_vector<unsigned int> > loops;
00146   bool valid = imesh_trace_half_edge_loops(mesh.half_edges(),contours,loops);
00147   assert(valid);
00148 
00149   return loops;
00150 }
00151 
00152 
00153 //: Mark contour half edges as seen from center of projection \param pt
00154 //  For each contour edge the half edge with the face normal opposing dir is marked
00155 vcl_vector<bool>
00156 imesh_detect_contours(const imesh_mesh& mesh, const vgl_point_3d<double>& pt)
00157 {
00158   assert(mesh.has_half_edges());
00159   assert(mesh.faces().has_normals());
00160 
00161   const imesh_half_edge_set& half_edges = mesh.half_edges();
00162   const vcl_vector<vgl_vector_3d<double> >& normals = mesh.faces().normals();
00163   const imesh_vertex_array<3>& verts = mesh.vertices<3>();
00164 
00165 
00166   const unsigned int num_edges = half_edges.size()/2;
00167   vcl_vector<bool> is_contour(half_edges.size(),false);
00168   for (unsigned int i=0; i<num_edges; ++i)
00169   {
00170     const imesh_half_edge& he = half_edges[i*2];
00171     const imesh_half_edge& ohe = half_edges[he.pair_index()];
00172     vgl_vector_3d<double> dir1 = vgl_point_3d<double>(verts[he.vert_index()]) - pt;
00173     vgl_vector_3d<double> dir2 = vgl_point_3d<double>(verts[ohe.vert_index()]) - pt;
00174     vgl_vector_3d<double> dir = normalized(dir1+dir2);
00175 
00176     double v1=0.0, v2=0.0;
00177     if (!he.is_boundary())
00178       v1 = dot_product(normalized(normals[he.face_index()]),dir);
00179     if (!ohe.is_boundary())
00180       v2 = dot_product(normalized(normals[ohe.face_index()]),dir);
00181 
00182     double prod = v1*v2;
00183     if (prod<=0.0) {
00184       if (v1<0.0) is_contour[i*2] = true;
00185       if (v2<0.0) is_contour[i*2+1] = true;
00186     }
00187   }
00188   return is_contour;
00189 }
00190 
00191 
00192 //: Segment the faces into groups of connected components
00193 vcl_vector<vcl_set<unsigned int> >
00194 imesh_detect_connected_components(const imesh_half_edge_set& he)
00195 {
00196   vcl_vector<vcl_set<unsigned int> > components;
00197   vcl_vector<bool> visited(he.num_faces(),false);
00198   for (unsigned int i=0; i<visited.size(); ++i)
00199   {
00200     if (visited[i]) continue;
00201     components.push_back(imesh_detect_connected_faces(he,i));
00202     for (vcl_set<unsigned int>::const_iterator itr=components.back().begin();
00203          itr != components.back().end(); ++itr)
00204       visited[*itr] = true;
00205   }
00206   return components;
00207 }
00208 
00209 
00210 //: Compute the set of all faces in the same connected component as \a face
00211 vcl_set<unsigned int>
00212 imesh_detect_connected_faces(const imesh_half_edge_set& he, unsigned int face)
00213 {
00214   vcl_set<unsigned int> component;
00215   vcl_vector<unsigned int > stack(1,face);
00216 
00217   while (!stack.empty())
00218   {
00219     unsigned int f = stack.back();
00220     stack.pop_back();
00221     if (component.find(f) != component.end())
00222       continue;
00223 
00224     component.insert(f);
00225     imesh_half_edge_set::f_const_iterator fi = he.face_begin(f);
00226     imesh_half_edge_set::f_const_iterator end = fi;
00227     do {
00228       unsigned int nf = he[fi->pair_index()].face_index();
00229       if (nf != imesh_invalid_idx)
00230         stack.push_back(nf);
00231       ++fi;
00232     } while (fi != end);
00233   }
00234 
00235   return component;
00236 }