MeshCleaner.h

Go to the documentation of this file.

#ifndef Cleaner_HEADER
#define Cleaner_HEADER

#ifdef HAVE_CONFIG_H
# include <config.h>
#endif

#include "Mesh.h"
#include <utilities/hash.h>
#include <vector>

#if 0
#define dprintf(...) printf(__VA_ARGS__);
#else
#define dprintf(...)
#endif

template <unsigned ambient, unsigned topological>
class MeshCleaner {
  BOOST_STATIC_ASSERT(ambient >= topological);
  BOOST_STATIC_ASSERT(topological >= 2);

  typedef typename Mesh<ambient, topological>::Delaunay Delaunay;
  typedef typename Delaunay::Simplex Simplex;
  typedef typename Delaunay::Vertex Vertex;
  typedef typename Delaunay::Complex Complex;
  typedef typename Complex::OSimplex OSimplex;
  typedef ::GenericMesh<ambient> GenericMesh;

  struct FindCorners {
    hudson::FastHashSet<Vertex*, hudson::HashPointer<Vertex>, false > corners_;
    std::vector<Simplex> corner_seeds_;
    template <class v_iterator>
    FindCorners(v_iterator begin, v_iterator end) {
      corners_.insert(begin, end);
    }
    static bool canEnter(const Simplex&) { return true; }
    bool choose(const Simplex& s) {
      BOOST_FOREACH(Vertex *v, s) {
        if (corners_.has(v)) {
          corner_seeds_.push_back(s);
          break;
        }
      }
      return false;
    }
  };


  typedef hudson::OrderedHashSet<Simplex, hudson::hash_by_call<Simplex> >
    SimplexSet;

  struct FindAllDead {
    const SimplexSet& live;
    SimplexSet& dead;
    FindAllDead(const SimplexSet& live, SimplexSet& dead) : live(live), dead(dead) { }
    static bool canEnter(const Simplex&) { return true; }
    bool choose(const Simplex& s) {
      if (!live.has(s)) {
        dead.insert(s);
      }
      return false;
    }
  };

  public:
  static void traverseComponent(Delaunay *del, const Simplex& seed,
      SimplexSet& component,
      hudson::SmallList<Simplex>& next_seeds)
  {
    const Complex& c = del->getComplex();
    hudson::FrontStack<Simplex> fringe;
    // std::stack<Simplex> fringe; // inexplicably, this fails.
    fringe.push(seed);
    dprintf("traversing component from seed %s\n", seed.toString().c_str());

    // Run DFS: find all the cells connected to the corners and mark them.
    while (!fringe.empty()) {
      Simplex s = fringe.top();
      fringe.pop();

      if (component.has(s)) { continue; /* already been here */ }
      component.insert(s);
      dprintf("  reachable: %s\n", s.toString().c_str());

      // collect neighbours
      OSimplex os = c.orientSimplex(s);
      for (size_t faceID = 0; faceID <= topological; ++faceID) {
        // Jump into the neighbour..
        OSimplex neigh = c.getNeighbour(faceID, os);

        // If there is no neighbour, stop here.
        if (!neigh.isInternal()) {
          continue;
        }

        // Compute the face.
        OSimplex fi = c.getFace(faceID, os);

        // Enter through this face unless all vertices share a
        // common (d-1)-mesh.  TODO: this was copied from
        // Mesh::isDegenerate; should be factored out.

        // Compute the intersection of the set of (d-1) meshes on which the
        // vertices of this face lie.
        std::vector<const GenericMesh*> meshes;
        std::copy(fi[0]->begin_uppers(topological-1), fi[0]->end_uppers(topological-1),
            std::back_inserter(meshes));

        for(size_t i = 1 ; i <= topological - 1; ++i) {
          // This uses the assumption that the meshes are stored in sorted order.
          std::vector<const GenericMesh*> tmp;
          std::set_intersection(meshes.begin(), meshes.end(),
              fi[i]->begin_uppers(topological - 1), fi[i]->end_uppers(topological - 1),
              std::back_inserter(tmp));
          meshes.swap(tmp);
        }

        // Either all vertices share one common d-1 face, or none.  Not two!
        assert(meshes.size() == 0 || meshes.size() == 1);

        if (meshes.empty()) {
          // If they share no common face, then we can continue on through.
          fringe.push(neigh);
        } else {
          // This is a boundary, stop here and mark it as a seed.
          next_seeds.push_front(neigh);
        }
      }
    }
  }

  template <class v_iterator>
  static bool clean(Delaunay *del, v_iterator begin_corners, v_iterator end_corners,
      hudson::SmallList<Simplex> *optional_output = 0) {
    const Complex& c = del->getComplex();

    dprintf("Performing cleanup on a %u-mesh\n", topological);

    // Search for the corners; seed the DFS there (this takes its own DFS).
    FindCorners fc(begin_corners, end_corners);
    c.dfsBySimplex(fc, c.getHandle());

    dprintf("Computing exterior simplices\n");
    SimplexSet exterior;
    hudson::SmallList<Simplex> interior_seeds;
    BOOST_FOREACH(const Simplex& seed, fc.corner_seeds_) {
      traverseComponent(del, seed, exterior, interior_seeds);
    }

    dprintf("Computing interior simplices\n");
    SimplexSet interior;
    hudson::SmallList<Simplex> ignored_seeds;
    BOOST_FOREACH(const Simplex& seed, interior_seeds) {
      if (exterior.has(seed)) { continue; } // handle walking through a hole
      traverseComponent(del, seed, interior, ignored_seeds);
    }

    if (interior.empty()) {
      // Not watertight -- we killed all the simplices.  Leave the mesh untouched
      // and return false.
      dprintf("Non-watertight!\n");
      return false;
    }

    // Compute the set of dead simplices.
    FindAllDead fd(interior, exterior);
    c.dfsBySimplex(fd, c.getHandle());

    // Get a writeable handle to the complex and slough off the necrosis.
    Complex& c_writable = del->getComplexRW();
    assert(!exterior.has(interior.front()));
    c_writable.setHandle(interior.front());
    BOOST_FOREACH(const Simplex& toslough, exterior) {
      c_writable.checkedRemoveSimplex(toslough);
    }

    if (optional_output) {
      BOOST_FOREACH(const Simplex& s, interior) {
        optional_output->push_front(s);
      }
    }

    // Who owns the vertices?  This seems to imply a memory leak...
    return true;
  }
};

#undef dprintf

#endif

---