MePolyPaverToMeshPts.cpp

Go to the documentation of this file.

//------------------------------------------------------------------------------
//------------------------------------------------------------------------------

//----- Included files ---------------------------------------------------------

// 1. Precompiled header

// 2. My own header
#include <xmsmesh/meshing/detail/MePolyPaverToMeshPts.h>

// 3. Standard library headers
#include <cfloat>
#include <list>

// 4. External library headers
#include <boost/unordered_set.hpp>

// 5. Shared code headers
#include <xmscore/points/pt.h>
#include <xmsinterp/geometry/geoms.h>
#include <xmsinterp/geometry/GmPtSearch.h>
#include <xmsmesh/meshing/detail/MeIntersectPolys.h>
#include <xmsmesh/meshing/detail/MePolyCleaner.h>
#include <xmsmesh/meshing/detail/MePolyOffsetter.h>
#include <xmsmesh/meshing/MePolyRedistributePts.h>
#include <xmscore/misc/Progress.h>
#include <xmscore/misc/XmError.h>
#include <xmscore/misc/XmConst.h>

// 6. Non-shared code headers

//----- Forward declarations ---------------------------------------------------

//----- External globals -------------------------------------------------------

//----- Namespace declaration --------------------------------------------------
namespace xms
{
// function to access method on MePolyRedistributePtsImpl that we wanted
// removed from the public interface. This makes it so MePolyOffsetterOutput
// is not visible in the public interface and is only part of the "detail"
// implementation.
void mePolyPaverRedistribute(BSHP<MePolyRedistributePts> a_redist,
                             const MePolyOffsetterOutput& a_input,
                             MePolyOffsetterOutput& a_out,
                             int a_polyOffsetIter);

//----- Constants / Enumerations -----------------------------------------------

//----- Classes / Structs ------------------------------------------------------

//----- Internal functions -----------------------------------------------------

//----- Class / Function definitions -------------------------------------------
namespace
{
class Poly
{
public:
  Poly()
  : m_iter(0)
  {
  }
  std::vector<Pt3d> m_outside;
  std::vector<std::vector<Pt3d>> m_inside;
  int m_iter;
};

class MePolyPaverToMeshPtsImpl : public MePolyPaverToMeshPts
{
public:
  MePolyPaverToMeshPtsImpl()
  : m_meshPts(nullptr)
  , m_xyTol(1e-9)
  , m_bias(1)
  , m_polyEnvelopeArea(0)
  , m_polyOffsetIter(1)
  {
  }

  bool PolyToMeshPts(const std::vector<Pt3d>& a_outPoly,
                     const std::vector<std::vector<Pt3d>>& a_inPolys,
                     double a_bias,
                     double a_xyTol,
                     std::vector<Pt3d>& a_meshPts) override;
  //------------------------------------------------------------------------------
  //------------------------------------------------------------------------------
  void SetRedistributor(BSHP<MePolyRedistributePts> a_) override { m_externalRedist = a_; }
  void Setup();
  void TearDown();
  void ProcessStack();
  void AddPolygonToMeshPoints(const Poly& a_poly, bool a_first);
  void DoPave(const Poly& a_poly);
  void CleanPave(const Poly& a_poly);
  void RedistributePts();
  void ClassifyPolys();
  double AreaFromPolyStack();

  BSHP<VecPt3d> m_meshPts;
  std::list<Poly> m_polyStack;
  BSHP<MePolyOffsetter> m_offsetter;
  BSHP<MePolyCleaner> m_cleaner;
  BSHP<MePolyRedistributePts> m_redist;
  BSHP<MePolyRedistributePts> m_externalRedist;
  double m_xyTol;
  double m_bias;
  double m_polyEnvelopeArea;
  int m_polyOffsetIter;
  boost::unordered_set<std::pair<double, double>> m_ptHash;

  std::vector<MePolyOffsetterOutput> m_offsetOutputs;
};

} // unnamed namespace

//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
BSHP<MePolyPaverToMeshPts> MePolyPaverToMeshPts::New()
{
  BSHP<MePolyPaverToMeshPts> ret(new MePolyPaverToMeshPtsImpl);
  return ret;
} // MePolyPaverToMeshPts::New
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
MePolyPaverToMeshPts::MePolyPaverToMeshPts()
{
}
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
MePolyPaverToMeshPts::~MePolyPaverToMeshPts()
{
}

namespace
{
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
bool MePolyPaverToMeshPtsImpl::PolyToMeshPts(const std::vector<Pt3d>& a_outPoly,
                                             const std::vector<std::vector<Pt3d>>& a_inPolys,
                                             double a_bias,
                                             double a_xyTol,
                                             std::vector<Pt3d>& a_meshPts)
{
  m_bias = a_bias;
  m_xyTol = a_xyTol;
  // Error checks
  if (a_outPoly.empty())
    return false;
  for (size_t i = 0; i < a_inPolys.size(); ++i)
  {
    if (a_inPolys[i].empty())
      return false;
  }

  m_polyStack.push_back(Poly());
  m_polyStack.back().m_outside = a_outPoly;
  m_polyStack.back().m_inside = a_inPolys;
  m_polyStack.back().m_iter = 1;

  Setup();
  ProcessStack();
  // fill the output variable
  a_meshPts.swap(*m_meshPts);
  TearDown();
  return true;
} // MePolyPaverToMeshPtsImpl::PolyToMeshPts
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsImpl::Setup()
{
  // Calculate tolerance for point comparison
  m_meshPts = BSHP<VecPt3d>(new VecPt3d());
  m_offsetter = MePolyOffsetter::New();
  m_cleaner = MePolyCleaner::New();
  m_polyEnvelopeArea = AreaFromPolyStack();
  if (!m_externalRedist)
  {
    m_redist = MePolyRedistributePts::New();
    m_redist->SetSizeFuncFromPoly(m_polyStack.front().m_outside, m_polyStack.front().m_inside,
                                  m_bias);
  }
  else
    m_redist = m_externalRedist;
} // MePolyPaverToMeshPtsImpl::Setup
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsImpl::TearDown()
{
  m_meshPts.reset();
  m_offsetter.reset();
  m_cleaner.reset();
  m_redist.reset();
  m_ptHash.clear();
} // MePolyPaverToMeshPtsImpl::TearDown
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsImpl::ProcessStack()
{
  Progress prog("Paving Polygon");

  double area;
  std::list<Poly>::iterator it(m_polyStack.begin());
  bool first(true);
  while (it != m_polyStack.end())
  {
    Poly& p(*it);

    m_polyOffsetIter = p.m_iter;
    // Pave the polygon
    DoPave(p);
    // clean the results from the pave
    CleanPave(p);
    // redistribute points on the polygons
    RedistributePts();
    // clean again after redistributing the points
    CleanPave(p);
    // classify the newly created polys into polygons defined by the "Poly"
    // class and put them on the stack
    ClassifyPolys();
    // add the points from this polygon to the mesh points
    AddPolygonToMeshPoints(p, first);
    first = false;

    // remove this poly from the stack and process the next one
    m_polyStack.erase(it);
    it = m_polyStack.begin();

    // do progress
    area = AreaFromPolyStack();
    prog.ProgressStatus(1.0 - (area / m_polyEnvelopeArea));
  }
} // MePolyPaverToMeshPtsImpl::ProcessStack
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsImpl::AddPolygonToMeshPoints(const Poly& a_poly, bool // a_first
)
{
  auto itEnd = m_ptHash.end();
  std::pair<double, double> pt;
  for (size_t i = 0; i < a_poly.m_outside.size(); ++i)
  {
    pt.first = a_poly.m_outside[i].x;
    pt.second = a_poly.m_outside[i].y;
    if (m_ptHash.find(pt) == itEnd)
    {
      m_ptHash.insert(pt);
      m_meshPts->push_back(a_poly.m_outside[i]);
    }
  }
  for (size_t i = 0; i < a_poly.m_inside.size(); ++i)
  {
    for (size_t j = 0; j < a_poly.m_inside[i].size(); ++j)
    {
      pt.first = a_poly.m_inside[i][j].x;
      pt.second = a_poly.m_inside[i][j].y;
      if (m_ptHash.find(pt) == itEnd)
      {
        m_ptHash.insert(pt);
        m_meshPts->push_back(a_poly.m_inside[i][j]);
      }
    }
  }

  // if (a_first)
  //{
  //  BSHP<GmPtSearch> ptSearch = GmPtSearch::New(true);
  //  ptSearch->VectorThatGrowsToSearch(m_meshPts);
  //  int idx;
  //  for (size_t i=0; i<a_poly.m_outside.size(); ++i)
  //  {
  //    if (!ptSearch->AddPtToVectorIfUnique(a_poly.m_outside[i], m_xyTol, idx))
  //      idx = -2;
  //  }
  //  for (size_t i=0; i<a_poly.m_inside.size(); ++i)
  //  {
  //    for (size_t j=0; j<a_poly.m_inside[i].size(); ++j)
  //    {
  //      if (!ptSearch->AddPtToVectorIfUnique(a_poly.m_inside[i][j], m_xyTol, idx))
  //        idx = -2;
  //    }
  //  }
  //}
  // else
  //{
  //  m_meshPts->insert(m_meshPts->end(), a_poly.m_outside.begin(),
  //    a_poly.m_outside.end());
  //  for (size_t i=0; i<a_poly.m_inside.size(); ++i)
  //  {
  //    m_meshPts->insert(m_meshPts->end(), a_poly.m_inside[i].begin(),
  //      a_poly.m_inside[i].end());
  //  }
  //}

} // MePolyPaverToMeshPtsImpl::AddPolygonToMeshPoints
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsImpl::DoPave(const Poly& a_poly)
{
  m_offsetOutputs.clear();

  MePolyOffsetterOutput offsetOut;
  MePolyOffsetter::polytype pIn(MePolyOffsetter::INSIDE_POLY), pOut(MePolyOffsetter::OUTSIDE_POLY);
  std::vector<Pt3d> pts;
  // pave in from the outer poly
  m_offsetter->Offset(a_poly.m_outside, pOut, offsetOut, m_xyTol);
  m_offsetOutputs.push_back(offsetOut);

  // pave out from the inner polys
  for (size_t i = 0; i < a_poly.m_inside.size(); ++i)
  {
    m_offsetter->Offset(a_poly.m_inside[i], pIn, offsetOut, m_xyTol);
    m_offsetOutputs.push_back(offsetOut);
  }
} // MePolyPaverToMeshPtsImpl::DoPave
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsImpl::CleanPave(const Poly& a_poly)
{
  m_cleaner->SetOriginalOutsidePolygon(a_poly.m_outside);
  // intersect the new inner polys with the other inner polys
  MePolyOffsetterOutput out, out2;
  m_cleaner->IntersectCleanInPolys(m_offsetOutputs, out, m_xyTol);

  // intersect the new inner polys with the outer polys
  m_cleaner->IntersectCleanInOutPolys(out, out2, m_xyTol);
  m_offsetOutputs.clear();
  m_offsetOutputs.push_back(out2);
} // MePolyPaverToMeshPtsImpl::CleanPave
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsImpl::RedistributePts()
{
  XM_ASSERT(m_redist);
  if (!m_redist)
    return;
  MePolyOffsetterOutput o;
  mePolyPaverRedistribute(m_redist, m_offsetOutputs[0], o, m_polyOffsetIter);
  m_offsetOutputs[0] = o;
} // MePolyPaverToMeshPtsImpl::RedistributePts
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsImpl::ClassifyPolys()
{
  if (m_offsetOutputs[0].m_pts.empty())
    return;

  std::vector<std::vector<size_t>> polys;
  MeIntersectPolys ip;
  ip.ClassifyPolys(m_offsetOutputs[0], polys);

  // put these on the stack
  std::vector<Pt3d>& pts(m_offsetOutputs[0].m_pts);
  Poly p;
  p.m_iter = m_polyOffsetIter + 1;
  for (size_t i = 0; i < polys.size(); ++i)
  {
    p.m_outside.resize(0);
    p.m_inside.resize(0);
    // get the outside loop
    std::vector<size_t>& oLoop(m_offsetOutputs[0].m_loops[polys[i][0]]);
    for (size_t j = 0; j < oLoop.size(); ++j)
    {
      p.m_outside.push_back(pts[oLoop[j]]);
    }
    // do inside loops
    for (size_t j = 1; j < polys[i].size(); ++j)
    {
      std::vector<size_t>& iLoop(m_offsetOutputs[0].m_loops[polys[i][j]]);
      p.m_inside.push_back(std::vector<Pt3d>());
      for (size_t k = 0; k < iLoop.size(); ++k)
      {
        p.m_inside.back().push_back(pts[iLoop[k]]);
      }
    }

    m_polyStack.push_back(p);
  }

} // MePolyPaverToMeshPtsImpl::ClassifyPolys
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
static double iEnvelopeArea(const std::vector<Pt3d>& a_pts)
{
  double area;
  Pt3d pMin(XM_DBL_HIGHEST),
    pMax(XM_DBL_LOWEST); // changed to XM_DBL_LOWEST because it is an arbitrarily low value
  for (size_t i = 0; i < a_pts.size(); ++i)
  {
    gmAddToExtents(a_pts[i], pMin, pMax);
  }
  area = (pMax.x - pMin.x) * (pMax.y - pMin.y);
  return area;
} // iEnvelopeArea
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
double MePolyPaverToMeshPtsImpl::AreaFromPolyStack()
{
  double area(0);
  // out polys are positive and in polys are negative
  std::list<Poly>::iterator it(m_polyStack.begin()), itEnd(m_polyStack.end());
  for (; it != itEnd; ++it)
  {
    Poly& p(*it);
    area += iEnvelopeArea(p.m_outside);
  }
  return area;
} // MePolyPaverToMeshPtsImpl::AreaFromPolyStack

} // unnamed namespace

} // namespace xms

#ifdef CXX_TEST

#include <xmsmesh/meshing/detail/MePolyPaverToMeshPts.t.h>

#include <xmscore/testing/TestTools.h>

// namespace xms {
using namespace xms;

//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsUnitTests::testCreateClass()
{
  BSHP<MePolyPaverToMeshPts> p = MePolyPaverToMeshPts::New();
  TS_ASSERT(p);
} // MeIntersectPolysTest::testCreateClass
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsUnitTests::testCase1()
{
  // x =     0                2                 4
  //
  // y=4     3--------4-------5-------6---------7
  //        |                                  |
  //        |                                  |
  //        2                                  8
  //        |                                  |
  //        |                                  |
  // y=2     1                                  9
  //        |                                  |
  //        |                                  |
  //        0                                  10
  //        |                                  |
  //        |                                  |
  // y=0    15--------14------13------12--------11
  //
  std::vector<Pt3d> outPoly = {{0, 1, 0}, {0, 2, 0}, {0, 3, 0}, {0, 4, 0}, {1, 4, 0}, {2, 4, 0},
                               {3, 4, 0}, {4, 4, 0}, {4, 3, 0}, {4, 2, 0}, {4, 1, 0}, {4, 0, 0},
                               {3, 0, 0}, {2, 0, 0}, {1, 0, 0}, {0, 0, 0}};
  std::vector<std::vector<Pt3d>> inPoly;
  MePolyPaverToMeshPtsImpl paver;
  double bias(1), tol(1e-9);
  std::vector<Pt3d> outPts;
  paver.PolyToMeshPts(outPoly, inPoly, bias, tol, outPts);
  std::vector<Pt3d> basePts = {
    {0.0, 1.0, 0.0},       {0.0, 2.0, 0.0},       {0.0, 3.0, 0.0},       {0.0, 4.0, 0.0},
    {1.0, 4.0, 0.0},       {2.0, 4.0, 0.0},       {3.0, 4.0, 0.0},       {4.0, 4.0, 0.0},
    {4.0, 3.0, 0.0},       {4.0, 2.0, 0.0},       {4.0, 1.0, 0.0},       {4.0, 0.0, 0.0},
    {3.0, 0.0, 0.0},       {2.0, 0.0, 0.0},       {1.0, 0.0, 0.0},       {0.0, 0.0, 0.0},
    {0.8660, 0.8660, 0.0}, {0.8660, 1.8740, 0.0}, {0.8660, 2.8820, 0.0}, {1.6220, 3.1340, 0.0},
    {2.6300, 3.1340, 0.0}, {3.1340, 2.6300, 0.0}, {3.1340, 1.6220, 0.0}, {2.8820, 0.8660, 0.0},
    {1.8740, 0.8660, 0.0}, {2.3148, 1.7390, 0.0}, {1.7390, 2.0535, 0.0}, {2.3323, 2.3802, 0.0}};
  TS_ASSERT_DELTA_VECPT3D(basePts, outPts, 1e-4);
} // MePolyPaverToMeshPtsUnitTests::testCase1
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
void MePolyPaverToMeshPtsUnitTests::testCase2()
{
  // x =     0                2                 4
  //
  // y=4     3--------4-------5-------6---------7
  //        |                                  |
  //        |                                  |
  //        2                                  8
  //        |                18--17            |
  //        |                 |  |             |
  // y=2     1                19--16            9
  //        |                                  |
  //        |                                  |
  //        0                                  10
  //        |                                  |
  //        |                                  |
  // y=0    15--------14------13------12--------11
  //
  std::vector<Pt3d> outPoly = {{0, 1, 0}, {0, 2, 0}, {0, 3, 0}, {0, 4, 0}, {1, 4, 0}, {2, 4, 0},
                               {3, 4, 0}, {4, 4, 0}, {4, 3, 0}, {4, 2, 0}, {4, 1, 0}, {4, 0, 0},
                               {3, 0, 0}, {2, 0, 0}, {1, 0, 0}, {0, 0, 0}};
  std::vector<std::vector<Pt3d>> inPoly;
  std::vector<Pt3d> tmp = {{2.25, 2, 0}, {2.25, 2.25, 0}, {2, 2.25, 0}, {2, 2, 0}};
  inPoly.push_back(tmp);
  MePolyPaverToMeshPtsImpl paver;
  double bias(1), tol(1e-9);
  std::vector<Pt3d> outPts;
  paver.PolyToMeshPts(outPoly, inPoly, bias, tol, outPts);
  std::vector<Pt3d> basePts = {
    {0.0, 1.0, 0.0},       {0.0, 2.0, 0.0},       {0.0, 3.0, 0.0},       {0.0, 4.0, 0.0},
    {1.0, 4.0, 0.0},       {2.0, 4.0, 0.0},       {3.0, 4.0, 0.0},       {4.0, 4.0, 0.0},
    {4.0, 3.0, 0.0},       {4.0, 2.0, 0.0},       {4.0, 1.0, 0.0},       {4.0, 0.0, 0.0},
    {3.0, 0.0, 0.0},       {2.0, 0.0, 0.0},       {1.0, 0.0, 0.0},       {0.0, 0.0, 0.0},
    {2.25, 2.0, 0.0},      {2.25, 2.25, 0.0},     {2.0, 2.25, 0.0},      {2.0, 2.0, 0.0},
    {0.8660, 0.8660, 0.0}, {0.8660, 1.7951, 0.0}, {0.8660, 2.7070, 0.0}, {1.3659, 3.1340, 0.0},
    {2.2505, 3.1340, 0.0}, {3.1340, 3.1224, 0.0}, {3.1340, 2.2316, 0.0}, {3.1340, 1.3448, 0.0},
    {2.6848, 0.8660, 0.0}, {1.7732, 0.8660, 0.0}, {1.7835, 2.1250, 0.0}, {2.0069, 1.7912, 0.0},
    {2.4260, 1.8985, 0.0}, {2.4335, 2.3440, 0.0}, {2.0093, 2.4589, 0.0}};
  TS_ASSERT_DELTA_VECPT3D(basePts, outPts, 1e-4);
} // MePolyPaverToMeshPtsUnitTests::testCase1

//} // namespace xms
#endif