I'll describe the whole scenario to avoid asking for how to implement what i think is the solution rather than asking what actually could be a solution to my problem.
I have to generate a raster from sparse points and constraints doing natural neighbor interpolation (NNI). The input is so large that imposes to use an out-of-core solution. I then chosed to partition the input domain into regular overlapping square tiles (hence not a pure partition due to overlap at tile borders), compute the value at each pixel with NNI, and then assemble the final tiles (this time non overlapping) by applying partition-of-unit (PoU) mixing at borders. To do a NNI on a tile, I take all points and constraints assigned to the tile and to its 8 adjacent neighbors, and construct a constrained delaunay triangulation (CDT) to be used to retrieve NNI coordinates. The problem is that, when a region of the dataset is extremely sparse, taking only the 8 neighbors does not guarantee that adjacent tiles will match at borders (even with PoU strategies) due to extremely different triangulations.
The solution I am using right now is to start at the tile T being processed and make a CDT with its points. Then, I proceed performing a breadth-first flood-fill visit on the tile domain and add (all at once) the points of the tile being visited to the current CDT if at least one of its points is "in conflict" with at least one triangle contained in the bounding rect of T (I "add" this visited tile to the CDT). If a new tile is added to the CDT, all its direct 8 neighbors are scheduled for visit (if not already visited). The visit ends when no more tiles are needed: i thought that this would mean that if I add another point (from a tile at the boundary of the current added tiles region) to the CDT, performing NNI on a point inside the bounding rect of T will continue giving the same result because the new point would not affect the triangulation in T.
Now, this is true (the triangulation inside T is not affected), but this does not mean that NNI is the same because I think I should also check that (at least!) all triangles adjacent to the ones overlapping T remain not in conflict with the point tested to be added.
Is this improvement going to ensure NNI is unchanged? Have you a better solution to the original problem?
I'm sure my description of the problem and the current solution is far from optimal, so please comment and give me any suggestion you may think would improve my explaination.
Thank you very much.
EDIT
Maybe a possible solution could be to take the bounding box B of the circumcircles of all the triangles overlapping T and use B as a conservative boundary, that is, I only add points belonging to all tiles overlapping B.
Or maybe the opposite? As the above is not symmetric, maybe it is better to, in advance, calculate B as above for each tile and then, while processing tile T, include in its CDT all tiles whose "B" overlaps T.
