Insertion of a grid of points

[LeoPizzo1]

I founded this library very promising. Well done and reasonably stable. I see two problems, to use it effectively in everyday tessellation. In the tessellation algorithms, often is useful to prepare a grid of points and superimpose on it the edge loops. Now all the points are inserted in the same way, regardless if they have a "regular" distribution or they are unorganized. It will be really cool to create an "initVertexGrid" function that initialize the data structure from a grid of vertices (vector<vector<>>?) and apply a predefined pattern to the grid of points in input.

[artem-ogre]

Hi, @LeoPizzo1

Thanks for the feedback.

I am not sure I follow exactly what you mean. Do you mean inserting existing triangulation with points and edges as-is at once?

If yes, that would be possible. However it will not be possible to insert points outside of the already-inserted triangulation as there is no super-triangle to catch them.

[LeoPizzo1]

Usually I have a 2d domain that is subdivided in a regular grid.

On this grid I add the edges defining the trimming loops.

The grid has a fixed structure, so it is possible to initialize the triangulation from this, using its concept of regularity, applying a "pattern" of triangles and avoid the huge set of computations that are required by the insertions.

The result is that I will insert only 48 points instead of 384, in this case.

[LeoPizzo1]

I wrote a simple (raw) algorithm that creates the grid topology and them adds the boundaries points and edges. I created a grid class

template<class T, class pt_type = V2d<T> >
class TriangulationGrid {
public:
  TriangulationGrid( size_t nRows, size_t nCols );

and implemented a function in CDT::Triangulation: void initGrid( const CDT::TriangulationGrid<T>& grid ) The boundary points are added via the already existing function void insertVertices(const std::vector<V2d<T> >& vertices);

The image below shows the trend line with up to ~350.000 grid points and a 700 points (edges) single closed boundary. The elapsed time is measured with boost::timer::auto_cpu_timer. The last point in the list shows the difference of the approaches: 22.765625 versus 0.859375 seconds.

[artem-ogre]

Thanks for the input, really appreciate it! I am away from computer for some days. I will look at it as soon as I have time, but probably not sooner than in a week or two.

[artem-ogre]

@LeoPizzo1

I've modified code so that eraseOuterTrianglesAndHoles works with arbitrary target triangulation: grid in your case instead of super-triangle. It is pushed to a feature-branch: https://github.com/artem-ogre/CDT/tree/feature/19-allow-fast-insertion-of-target

You can now use the following workflow:

1. Insert grid in CDT quickly by modifying vertices and triangles directly to bypass triangulation calculations. Just make sure winding and triangle neighbors are set properly (use noNeighbor four boundary triangles)
2. Then you can insert vertices and constraints with insertVertices and insertEdges
3. Finally, finish with a call to eraseOuterTrianglesAndHoles

Please let me know if this works for you.

[LeoPizzo1]

Thanks. I will try it today and let you know.

[artem-ogre]

I found few issues that need to be resolved before using custom targets instead of super-tri is possible. I am slowly progressing on this.

May I ask why do you prefer to use a grid? Do you hope to get better triangles this way (e.g., for finite element analysis)? Because even with grid you still may get 'bad' triangles.

[LeoPizzo1]

Ciao, sorry but I hadn't yet time enough to test the functionality. In an area of our application (isovalues search and display), we need to tessellate some parametric (UV) surfaces in order to have a regular triangle distribution to provide quality results. This involves a lot of work in our framework to prepare a grid of "equidistant" points on the surface in 3d (usually there is no relation between surfaces parametrization and 3d space). We need, afterwards, to superimpose to this grid the discretized trimming loops, imposing the edges. I am analyzing the possibility to use your class, that seems born for doing this.

I need the grid for 2 reasons:

• I decide by myself the "general shape" of the inner triangles, as you stated.
• It is way faster, since the algorithm does not need to create the connectivity by itself.

BTW, I tested, as I mentioned before, and it is working correctly using the supertriangle: I ended up adding several useless triangles, but that's it.

Another interesting property of the grid is that it will be possible to regularize the X-Y axis so that the grid cell is (approximately) square. Operating such a scaling operation will provide a better behaviour of the predicates and also the results tessellation will have a more "pleasant" distribution.

Given our premises (we are "equidistancing" the points in 3d on the surface), regularizing the grid we will obtain a "quasi 3d" tessellation. But this task has to be accomplished outside the library, IMHO.

Last point, to help the developers to integrate in them own systems, it could be interesting to split the private methods in private and protected. I mean that, at least, the functions: void addSuperTriangle(const Box2d& box); void insertVertex(const V2d& pos); void insertEdge(Edge edge); Should be callable by a developer that is deriving itd class from CDT::Triangulation. At least this is what I am needing in my tests.

[artem-ogre]

Please check out this branch: https://github.com/artem-ogre/CDT/tree/feature/19-initialize-with-grid

I've added a file InitializeWithGrid.h containing method initializeWithGrid and re-factored code to allow custom 'super-geometries' other than super-triangle.

Usage example:

Triangulation cdt(...);
// ...
CDT::initializeWithGrid(xmin, xmax, ymin, ymax, 100, 100, cdt);
cdt.insertVertices(...);
cdt.insertEdges(...)
cdt.eraseOuterTrianglesAndHoles(); // 'eraseOuterTriangles' should work too

Output example:

Could you please try it out and report your feedback?

[LeoPizzo1]

Ciao, interesting, but actually it is not what I need/implemented. You correctly create a regular 2d grid, while in my case I need to adapt the 2d metric to the 3d one, thus the grid "lines" are not equally spaced in XY. Moreover, if the surface has discontinuities (this happens with NURBS) I must insert a line of points at the discontinuity.

Given these request, I think that I have to provide on my side to the task of computing the points. In my first rough implementation I considered to use a specific object to feed the tessellation with a grid. Surely several checks must be done on this (i.e.: the ordering must be respected), and the usage of the vector of vectors replaced with a single vector, but as POC it worked.

template<class T, class pt_type = V2d<T> >
class TriangulationGrid {
public:
  TriangulationGrid( ) = delete;
  TriangulationGrid( size_t nRows, size_t nCols )
  {
    grid.resize( nRows, std::vector< pt_type >( nCols ) );
  }

  size_t rows( ) const { return grid.size( ); }
  size_t cols( ) const { return grid.empty( ) ? 0 : grid[0].size( ); }

  size_t size( ) const { return rows( )*cols( ); }

  inline pt_type& first( ) { return grid[0][0]; }
  inline const pt_type& first( ) const { return grid[0][0]; }

  inline pt_type& last( ) { return grid[grid.size( ) - 1][grid[0].size( ) - 1]; }
  inline const pt_type& last( ) const { return grid[grid.size( ) - 1][grid[0].size( ) - 1]; }

  inline pt_type& operator()( size_t iRow, size_t iCol ) { return grid[iRow][iCol]; }
  inline const pt_type& operator()( size_t iRow, size_t iCol ) const { return grid[iRow][iCol]; }

  inline std::vector<pt_type> &operator[]( size_t iRow ) { return grid[iRow]; }
  inline const std::vector<pt_type> &operator[]( size_t iRow ) const { return grid[iRow]; }

  CDT::Box2d<T> envelope( ) const
  {
    const T max = std::numeric_limits<T>::max( );
    Box2d<T> box = { {max, max}, {-max, -max} };
    for( auto &vertices : grid ) {
      typedef typename std::vector<V2d<T> >::const_iterator Cit;
      for( Cit it = vertices.begin( ); it != vertices.end( ); ++it ) {
        const V2d<T>& v = *it;
        box.min.x = std::min( v.x, box.min.x );
        box.max.x = std::max( v.x, box.max.x );
        box.min.y = std::min( v.y, box.min.y );
        box.max.y = std::max( v.y, box.max.y );
      }
    }
    return box;
  }

  std::vector<pt_type> linearized( ) const
  {
    std::vector<pt_type> res;
    res.reserve( size( ) );
    {
      // create the nodes list
      for( auto &row : grid ) {
        res.insert( res.end( ), row.begin( ), row.end( ) );
      }
    }
    return res;
  }
protected:
  std::vector< std::vector<pt_type> > grid;
};

A furthermore improvement could be to allow the insertion of not uniform grid, to obtain hexagonal shapes (equilateral triangles), but for this a more sophisticated creation is needed, thus probably there is no way to do it in the class itself.

[artem-ogre]

I am not sure I follow. Do you have a working PoC branch somewhere?

You can insert arbitrary geometry (see initializeWithGrid for inspiration). Feel free to provide your own initializeWithXXX.

[LeoPizzo1]

No branches as far as now. If you are not in a hurry, in the next days (weeks maybe) I will be able to provide one.

[LeoPizzo1]

Ciao,

I made some tests and I discovered that the problem (by now) is in "insertEdges", CDT.h row 392: it considers the 3 supertriangle indexes. I'm not really used to git, it is easier (for me) and faster (for you) to provide the code I used to obtain the performance test above (obviously updated with the changes to test the grid without supertriangle) and some instructions. I promise I'll take some time during this year's end time to study git a little bit. Here the code. testcdt.zip

The test is to create a grid and superimpose a circular trimming loop (just one) with increasing number of points: each perform is done both using the initualization with the grid and using the standard insertion. It stores the time measured by boost::auto_cpu_timer and saves all of them in a file. Feel free to remove what is

There is a function (TestCDT::Do) that creates a folder in temp and store the temporary results (ascii stl and pt files) to make some analysis. There is a class (TriangulationGrid) that initialize herself with the grid of points, provided the extremal points (min and max) and the number of row/columns of data. Then I extended the CDT::Triangulation to accept an initialization with the TriangulationGrid created. To accomplish this task I had to set as protected, instead of private, the scope of CDT::Triangulaton variables (to initialize r-tree).

Row 231 define a variable to choose whether to add or not the supertriangle.

In debug mode, the TestCDT::Do() starts creating a 5 row 4 cols grid, my original test case.

FYI, anyway, verifyTopology does not control that the vertices data contains the triangles to which they belongs.

[LeoPizzo1]

Adding a control about the addition of supertriangle now it seems to work correctly. CDT.zip I run all the tests, up to 700x500 points grid "trimmed" with 700 points circle.

[artem-ogre]

Thanks @LeoPizzo1 I will try to look at it after the holidays.

[artem-ogre]

Hi, sorry I don't have time to properly study the code you send me—not enough spare time for that.

If you take a look at the branch you will see that it is possible to initialize triangulation with absolutely any geometry (the only condition is that the points added with insertVertices afterwards will hit triangles of that geometry).

If you look at initializeWithGrid you will see that the general workflow is:

Triangulation cdt(...);
// add absolutely any geometry by using `cdt.triangles` and `cdt.vertices` directly. 
// it can be irregular/triangular/any other grid or geometry
// Adding constraints by manipulating `fixedEdges` should be possible too but I did not test it
// initializeWithGrid makes a regular grid, but it can be anything else.
...
cdt.initializedWithCustomSuperGeometry(); // this will make sure that R-tree (if used) will be initialized with your custom geometry
// now you can insert more vertices and edges
cdt.insertVertices(...);
cdt.insertEdges(...)
cdt.eraseOuterTrianglesAndHoles(); // 'eraseOuterTriangles' will work too

I fail to understand what is missing in this approach in order for it to be applicable in your case.

The plan is to finalize and merge this branch as it definitely has a lot of value: for example in some hard cases I was able to increase performance 50x by adding a grid which avoided calling extended-precision calculations in predicates.

[LeoPizzo1]

Seems to work like a charm, now. Speed is consistent with my original tests. For me you can merge to trunk.

To complete the issue at 100%, in the real life the X-Y steps might not be regular. So, it could be a good idea to provide, at least, a function accepting the X and Y ordered values.

 template <typename T, typename OutputIt, typename ValueIt> 
 void generateGridVertices( OutputIt outFirst, ValueIt xFirst, const ValueIt xLast, ValueIt yFirst, const ValueIt yLast )
 {
 }

And/or, maybe even better, something like this:

template <typename T, typename OutputIt, typename TValueGetX, typename TValueGetY>
void generateGridVertices( OutputIt outFirst, TValueGetX xGetter, TValueGetY yGetter, const std::size_t xres, const std::size_t yres )
{
  for( std::size_t iy = 0; iy <= yres; ++iy, y += ystep ) {
    for( std::size_t ix = 0; ix <= xres; ++ix, x += xstep ) {
      const T x = xGetter( ix, iy );
      const T y = yGetter( ix, iy );
      Vertex<T> v;
      v.pos = V2d<T>::make( x, y );
      ...
    }
  }
}

[artem-ogre]

Added irregular grid functionality. Function signature is:

/**
 * Make a triangulation that uses irregular grid triangles instead of
 * super-triangle. Irregular grid is given by collections of X- and Y-ticks
 *
 * @tparam T type of vertex coordinates (e.g., float, double)
 * @tparam TXCoordIter iterator dereferencing to X coordinate
 * @tparam TYCoordIter iterator dereferencing to Y coordinate
 * @param xfirst beginning of X-ticks range
 * @param xlast end of X-ticks range
 * @param yfirst beginning of Y-ticks range
 * @param ylast end of Y-ticks range
 * @param out triangulation to initialize with grid super-geometry
 */
template <typename T, typename TXCoordIter, typename TYCoordIter>
void initializeWithIrregularGrid(
    const TXCoordIter xfirst,
    const TXCoordIter xlast,
    const TYCoordIter yfirst,
    const TYCoordIter ylast,
    Triangulation<T>& out)

[artem-ogre]

Merged in PR: #30