| 1 | //
 | 
|---|
| 2 | // surf.h
 | 
|---|
| 3 | //
 | 
|---|
| 4 | // Copyright (C) 1996 Limit Point Systems, Inc.
 | 
|---|
| 5 | //
 | 
|---|
| 6 | // Author: Curtis Janssen <cljanss@limitpt.com>
 | 
|---|
| 7 | // Maintainer: LPS
 | 
|---|
| 8 | //
 | 
|---|
| 9 | // This file is part of the SC Toolkit.
 | 
|---|
| 10 | //
 | 
|---|
| 11 | // The SC Toolkit is free software; you can redistribute it and/or modify
 | 
|---|
| 12 | // it under the terms of the GNU Library General Public License as published by
 | 
|---|
| 13 | // the Free Software Foundation; either version 2, or (at your option)
 | 
|---|
| 14 | // any later version.
 | 
|---|
| 15 | //
 | 
|---|
| 16 | // The SC Toolkit is distributed in the hope that it will be useful,
 | 
|---|
| 17 | // but WITHOUT ANY WARRANTY; without even the implied warranty of
 | 
|---|
| 18 | // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 | 
|---|
| 19 | // GNU Library General Public License for more details.
 | 
|---|
| 20 | //
 | 
|---|
| 21 | // You should have received a copy of the GNU Library General Public License
 | 
|---|
| 22 | // along with the SC Toolkit; see the file COPYING.LIB.  If not, write to
 | 
|---|
| 23 | // the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
 | 
|---|
| 24 | //
 | 
|---|
| 25 | // The U.S. Government is granted a limited license as per AL 91-7.
 | 
|---|
| 26 | //
 | 
|---|
| 27 | 
 | 
|---|
| 28 | #ifndef _math_isosurf_surf_h
 | 
|---|
| 29 | #define _math_isosurf_surf_h
 | 
|---|
| 30 | 
 | 
|---|
| 31 | #ifdef __GNUC__
 | 
|---|
| 32 | #pragma interface
 | 
|---|
| 33 | #endif
 | 
|---|
| 34 | 
 | 
|---|
| 35 | #ifdef HAVE_CONFIG_H
 | 
|---|
| 36 | #include <scconfig.h>
 | 
|---|
| 37 | #endif
 | 
|---|
| 38 | 
 | 
|---|
| 39 | #include <map>
 | 
|---|
| 40 | #include <set>
 | 
|---|
| 41 | #include <vector>
 | 
|---|
| 42 | 
 | 
|---|
| 43 | #include <math/isosurf/triangle.h>
 | 
|---|
| 44 | #include <math/isosurf/volume.h>
 | 
|---|
| 45 | #include <util/render/render.h>
 | 
|---|
| 46 | 
 | 
|---|
| 47 | namespace sc {
 | 
|---|
| 48 | 
 | 
|---|
| 49 | template <class C, class I>
 | 
|---|
| 50 | inline void
 | 
|---|
| 51 | erase_elements_by_value(C &container, I begin, I end)
 | 
|---|
| 52 | {
 | 
|---|
| 53 |   for (I i=begin; i!=end; i++) {
 | 
|---|
| 54 |       container.erase(*i);
 | 
|---|
| 55 |     }
 | 
|---|
| 56 | }
 | 
|---|
| 57 | 
 | 
|---|
| 58 | class TriangulatedSurface: public DescribedClass {
 | 
|---|
| 59 |   protected:
 | 
|---|
| 60 |     int _verbose;
 | 
|---|
| 61 |     int _debug;
 | 
|---|
| 62 | 
 | 
|---|
| 63 |     int _completed_surface;
 | 
|---|
| 64 | 
 | 
|---|
| 65 |     // sets of objects that make up the surface
 | 
|---|
| 66 |     std::set<Ref<Vertex> > _vertices;
 | 
|---|
| 67 |     std::set<Ref<Edge> > _edges;
 | 
|---|
| 68 |     std::set<Ref<Triangle> > _triangles;
 | 
|---|
| 69 | 
 | 
|---|
| 70 |     // map objects to an integer index
 | 
|---|
| 71 |     std::map<Ref<Vertex>,int> _vertex_to_index;
 | 
|---|
| 72 |     std::map<Ref<Edge>,int> _edge_to_index;
 | 
|---|
| 73 |     std::map<Ref<Triangle>,int> _triangle_to_index;
 | 
|---|
| 74 | 
 | 
|---|
| 75 |     // map integer indices to an object
 | 
|---|
| 76 |     std::vector<Ref<Vertex> > _index_to_vertex;
 | 
|---|
| 77 |     std::vector<Ref<Edge> > _index_to_edge;
 | 
|---|
| 78 |     std::vector<Ref<Triangle> > _index_to_triangle;
 | 
|---|
| 79 | 
 | 
|---|
| 80 |     // mappings between array element numbers
 | 
|---|
| 81 |     int** _triangle_vertex;
 | 
|---|
| 82 |     int** _triangle_edge;
 | 
|---|
| 83 |     int** _edge_vertex;
 | 
|---|
| 84 | 
 | 
|---|
| 85 |     // values for each of the vertices
 | 
|---|
| 86 |     int _have_values;
 | 
|---|
| 87 |     std::vector<double> _values;
 | 
|---|
| 88 | 
 | 
|---|
| 89 |     // what to use to integrate over the surface, by default
 | 
|---|
| 90 |     Ref<TriangleIntegrator> _integrator;
 | 
|---|
| 91 |     // other integrators, in terms of time & accuracy:
 | 
|---|
| 92 |     // _fast_integrator <= _integrator <= _accurate_interator
 | 
|---|
| 93 |     Ref<TriangleIntegrator> _fast_integrator;
 | 
|---|
| 94 |     Ref<TriangleIntegrator> _accurate_integrator;
 | 
|---|
| 95 | 
 | 
|---|
| 96 |     void clear_int_arrays();
 | 
|---|
| 97 | 
 | 
|---|
| 98 |     void complete_ref_arrays();
 | 
|---|
| 99 |     void complete_int_arrays();
 | 
|---|
| 100 | 
 | 
|---|
| 101 |     void recompute_index_maps();
 | 
|---|
| 102 | 
 | 
|---|
| 103 |     void add_triangle(const Ref<Triangle>&);
 | 
|---|
| 104 |     void add_vertex(const Ref<Vertex>&);
 | 
|---|
| 105 |     void add_edge(const Ref<Edge>&);
 | 
|---|
| 106 | 
 | 
|---|
| 107 |     // these members must be used to allocate new triangles and edges
 | 
|---|
| 108 |     // since specializations of TriangulatedSurface might need to
 | 
|---|
| 109 |     // override these to produce triangles and edges with interpolation
 | 
|---|
| 110 |     // data.
 | 
|---|
| 111 |     virtual Triangle* newTriangle(const Ref<Edge>&,
 | 
|---|
| 112 |                                   const Ref<Edge>&,
 | 
|---|
| 113 |                                   const Ref<Edge>&,
 | 
|---|
| 114 |                                   int orientation) const;
 | 
|---|
| 115 |     virtual Edge* newEdge(const Ref<Vertex>&,const Ref<Vertex>&) const;
 | 
|---|
| 116 | 
 | 
|---|
| 117 |     // this map of edges to vertices is used to construct the surface
 | 
|---|
| 118 |     std::map<Ref<Vertex>,std::set<Ref<Edge> > > _tmp_edges;
 | 
|---|
| 119 |   public:
 | 
|---|
| 120 |     TriangulatedSurface();
 | 
|---|
| 121 |     TriangulatedSurface(const Ref<KeyVal>&);
 | 
|---|
| 122 |     virtual ~TriangulatedSurface();
 | 
|---|
| 123 | 
 | 
|---|
| 124 |     // control printing
 | 
|---|
| 125 |     int verbose() const { return _verbose; }
 | 
|---|
| 126 |     void verbose(int v) { _verbose = v; }
 | 
|---|
| 127 | 
 | 
|---|
| 128 |     // set up an integrator
 | 
|---|
| 129 |     void set_integrator(const Ref<TriangleIntegrator>&);
 | 
|---|
| 130 |     void set_fast_integrator(const Ref<TriangleIntegrator>&);
 | 
|---|
| 131 |     void set_accurate_integrator(const Ref<TriangleIntegrator>&);
 | 
|---|
| 132 |     virtual Ref<TriangleIntegrator> integrator(int itri);
 | 
|---|
| 133 |     virtual Ref<TriangleIntegrator> fast_integrator(int itri);
 | 
|---|
| 134 |     virtual Ref<TriangleIntegrator> accurate_integrator(int itri);
 | 
|---|
| 135 | 
 | 
|---|
| 136 |     // construct the surface
 | 
|---|
| 137 |     void add_triangle(const Ref<Vertex>&,
 | 
|---|
| 138 |                       const Ref<Vertex>&,
 | 
|---|
| 139 |                       const Ref<Vertex>&);
 | 
|---|
| 140 |     Ref<Edge> find_edge(const Ref<Vertex>&, const Ref<Vertex>&);
 | 
|---|
| 141 |     virtual void complete_surface();
 | 
|---|
| 142 | 
 | 
|---|
| 143 |     // clean up the surface
 | 
|---|
| 144 |     virtual void remove_short_edges(double cutoff_length = 1.0e-6,
 | 
|---|
| 145 |                                     const Ref<Volume> &vol=0, double isoval=0.0);
 | 
|---|
| 146 |     virtual void remove_slender_triangles(
 | 
|---|
| 147 |                                     int remove_slender, double height_cutoff,
 | 
|---|
| 148 |                                     int remove_small, double area_cutoff,
 | 
|---|
| 149 |                                     const Ref<Volume> &vol=0, double isoval=0.0);
 | 
|---|
| 150 |     virtual void fix_orientation();
 | 
|---|
| 151 |     virtual void clear();
 | 
|---|
| 152 | 
 | 
|---|
| 153 |     // get information from the object sets
 | 
|---|
| 154 |     int nvertex() const { return _vertices.size(); };
 | 
|---|
| 155 |     Ref<Vertex> vertex(int i) const { return _index_to_vertex[i]; };
 | 
|---|
| 156 |     int vertex_index(const Ref<Vertex> &o) {
 | 
|---|
| 157 |       std::map<Ref<Vertex>,int>::iterator i = _vertex_to_index.find(o);
 | 
|---|
| 158 |       if (i != _vertex_to_index.end()) return i->second;
 | 
|---|
| 159 |       return -1;
 | 
|---|
| 160 |     }
 | 
|---|
| 161 |     int nedge() const { return _edges.size(); };
 | 
|---|
| 162 |     Ref<Edge> edge(int i) const { return _index_to_edge[i]; };
 | 
|---|
| 163 |     int edge_index(const Ref<Edge> &o) {
 | 
|---|
| 164 |       std::map<Ref<Edge>,int>::iterator i = _edge_to_index.find(o);
 | 
|---|
| 165 |       if (i != _edge_to_index.end()) return i->second;
 | 
|---|
| 166 |       return -1;
 | 
|---|
| 167 |     }
 | 
|---|
| 168 |     int ntriangle() const { return _triangles.size(); };
 | 
|---|
| 169 |     Ref<Triangle> triangle(int i) const { return _index_to_triangle[i]; }
 | 
|---|
| 170 |     int triangle_index(const Ref<Triangle> &o) {
 | 
|---|
| 171 |       std::map<Ref<Triangle>,int>::iterator i = _triangle_to_index.find(o);
 | 
|---|
| 172 |       if (i != _triangle_to_index.end()) return i->second;
 | 
|---|
| 173 |       return -1;
 | 
|---|
| 174 |     }
 | 
|---|
| 175 | 
 | 
|---|
| 176 |     // information from the index mappings
 | 
|---|
| 177 |     int triangle_vertex(int i,int j) const { return _triangle_vertex[i][j]; };
 | 
|---|
| 178 |     int triangle_edge(int i,int j) const { return _triangle_edge[i][j]; };
 | 
|---|
| 179 |     int edge_vertex(int i,int j) const { return _edge_vertex[i][j]; };
 | 
|---|
| 180 | 
 | 
|---|
| 181 |     // associate values with vertices
 | 
|---|
| 182 |     //void compute_colors(Volume&);
 | 
|---|
| 183 |     void compute_values(Ref<Volume>&);
 | 
|---|
| 184 | 
 | 
|---|
| 185 |     // properties of the surface
 | 
|---|
| 186 |     virtual double flat_area(); // use flat triangles
 | 
|---|
| 187 |     virtual double flat_volume(); // use flat triangles
 | 
|---|
| 188 |     virtual double area();
 | 
|---|
| 189 |     virtual double volume();
 | 
|---|
| 190 | 
 | 
|---|
| 191 |     // output of the surface
 | 
|---|
| 192 |     virtual void print(std::ostream&o=ExEnv::out0()) const;
 | 
|---|
| 193 |     virtual void print_vertices_and_triangles(std::ostream&o=ExEnv::out0()) const;
 | 
|---|
| 194 |     virtual void print_geomview_format(std::ostream&o=ExEnv::out0()) const;
 | 
|---|
| 195 |     virtual void render(const Ref<Render> &render);
 | 
|---|
| 196 | 
 | 
|---|
| 197 |     // print information about the topology
 | 
|---|
| 198 |     void topology_info(std::ostream&o=ExEnv::out0());
 | 
|---|
| 199 |     void topology_info(int nvertex, int nedge, int ntri, std::ostream&o=ExEnv::out0());
 | 
|---|
| 200 | };
 | 
|---|
| 201 | 
 | 
|---|
| 202 | 
 | 
|---|
| 203 | class TriangulatedSurfaceIntegrator {
 | 
|---|
| 204 |   private:
 | 
|---|
| 205 |     Ref<TriangulatedSurface> _ts;
 | 
|---|
| 206 |     int _itri;
 | 
|---|
| 207 |     int _irs;
 | 
|---|
| 208 |     double _r;
 | 
|---|
| 209 |     double _s;
 | 
|---|
| 210 |     double _weight;
 | 
|---|
| 211 |     double _surface_element;
 | 
|---|
| 212 |     Ref<Vertex> _current;
 | 
|---|
| 213 |     SCVector3 _dA;
 | 
|---|
| 214 |     Ref<TriangleIntegrator> (TriangulatedSurface::*_integrator)(int itri);
 | 
|---|
| 215 |     Ref<MessageGrp> _grp;
 | 
|---|
| 216 |   public:
 | 
|---|
| 217 |     TriangulatedSurfaceIntegrator();
 | 
|---|
| 218 |     // the surface cannot be changed until this is destroyed
 | 
|---|
| 219 |     TriangulatedSurfaceIntegrator(const Ref<TriangulatedSurface>&);
 | 
|---|
| 220 |     ~TriangulatedSurfaceIntegrator();
 | 
|---|
| 221 |     // Objects initialized by these operators are not automatically
 | 
|---|
| 222 |     // updated.  This must be done with the update member.
 | 
|---|
| 223 |     // The _grp is not copied.
 | 
|---|
| 224 |     void operator = (const TriangulatedSurfaceIntegrator&);
 | 
|---|
| 225 |     TriangulatedSurfaceIntegrator(const TriangulatedSurfaceIntegrator&i) {
 | 
|---|
| 226 |         operator = (i);
 | 
|---|
| 227 |       }
 | 
|---|
| 228 |     // Return the number of integration points.
 | 
|---|
| 229 |     int n();
 | 
|---|
| 230 |     // Assign the surface.  Don't do this while iterating.
 | 
|---|
| 231 |     void set_surface(const Ref<TriangulatedSurface>&);
 | 
|---|
| 232 |     // returns the number of the vertex in the current triangle
 | 
|---|
| 233 |     int vertex_number(int i);
 | 
|---|
| 234 |     inline double r() const { return _r; }
 | 
|---|
| 235 |     inline double s() const { return _s; }
 | 
|---|
| 236 |     inline double w() const { return _weight*_surface_element; }
 | 
|---|
| 237 |     double surface_element() const { return _surface_element; }
 | 
|---|
| 238 |     double weight() const { return _weight; }
 | 
|---|
| 239 |     const SCVector3& dA() const { return _dA; }
 | 
|---|
| 240 |     Ref<Vertex> current();
 | 
|---|
| 241 |     // Tests to see if this point is valid, if it is then
 | 
|---|
| 242 |     // _r, _s, etc are computed and 1 is returned.
 | 
|---|
| 243 |     int update();
 | 
|---|
| 244 |     // This can be used to loop through unique pairs of points.
 | 
|---|
| 245 |     // The argument should be a TriangulatedSurfaceIntegrator for
 | 
|---|
| 246 |     // the same surface as this.
 | 
|---|
| 247 |     int operator < (TriangulatedSurfaceIntegrator&i) {
 | 
|---|
| 248 |         update();
 | 
|---|
| 249 |         return _itri<i._itri?1:(_itri>i._itri?0:(_irs<i._irs?1:0));
 | 
|---|
| 250 |       }
 | 
|---|
| 251 |     // Goes to the next point.  Does not update.
 | 
|---|
| 252 |     void operator++();
 | 
|---|
| 253 |     inline void operator++(int) { operator++(); }
 | 
|---|
| 254 |     // setting TSI = i sets TSI to begin at the triangle i
 | 
|---|
| 255 |     int operator = (int);
 | 
|---|
| 256 |     int itri() const { return _itri; }
 | 
|---|
| 257 |     int irs() const { return _irs; }
 | 
|---|
| 258 |     // the number of points in the current triangle
 | 
|---|
| 259 |     int n_in_tri() const { return (_ts.pointer()->*_integrator)(_itri)->n(); }
 | 
|---|
| 260 |     void distribute(const Ref<MessageGrp> &);
 | 
|---|
| 261 |     void use_fast_integrator();
 | 
|---|
| 262 |     void use_accurate_integrator();
 | 
|---|
| 263 |     void use_default_integrator();
 | 
|---|
| 264 | };
 | 
|---|
| 265 | 
 | 
|---|
| 266 | class TriangulatedImplicitSurface: public TriangulatedSurface {
 | 
|---|
| 267 |   private:
 | 
|---|
| 268 |     // The surface is defined as an isosurface of the volume vol_.
 | 
|---|
| 269 |     Ref<Volume> vol_;
 | 
|---|
| 270 |     double isovalue_;
 | 
|---|
| 271 | 
 | 
|---|
| 272 |     int fix_orientation_;
 | 
|---|
| 273 |     int remove_short_edges_;
 | 
|---|
| 274 |     double short_edge_factor_;
 | 
|---|
| 275 |     int remove_slender_triangles_;
 | 
|---|
| 276 |     double slender_triangle_factor_;
 | 
|---|
| 277 |     int remove_small_triangles_;
 | 
|---|
| 278 |     double small_triangle_factor_;
 | 
|---|
| 279 |     double resolution_;
 | 
|---|
| 280 | 
 | 
|---|
| 281 |     int order_;
 | 
|---|
| 282 | 
 | 
|---|
| 283 |     int inited_;
 | 
|---|
| 284 |   public:
 | 
|---|
| 285 |     TriangulatedImplicitSurface(const Ref<KeyVal>&);
 | 
|---|
| 286 |     ~TriangulatedImplicitSurface();
 | 
|---|
| 287 | 
 | 
|---|
| 288 |     Ref<Volume> volume_object() const { return vol_; }
 | 
|---|
| 289 |     double isovalue() const { return isovalue_; }
 | 
|---|
| 290 | 
 | 
|---|
| 291 |     void init();
 | 
|---|
| 292 |     int inited() const { return inited_; }
 | 
|---|
| 293 | };
 | 
|---|
| 294 | 
 | 
|---|
| 295 | }
 | 
|---|
| 296 | 
 | 
|---|
| 297 | #endif
 | 
|---|
| 298 | 
 | 
|---|
| 299 | // Local Variables:
 | 
|---|
| 300 | // mode: c++
 | 
|---|
| 301 | // c-file-style: "CLJ"
 | 
|---|
| 302 | // End:
 | 
|---|