Intersect (vertical) line segment with triangulated surface (topo)

What I want to achieve: I have a topo surface (triangulated, “2.5D” type) and a point, and I want to find the relative position of the point to that topo:

  • is it outside the boundary of the topo or inside?
  • if it is inside: what is the vertical distance between point and topo?

I would say that is some kind of “ray tracing problem”, i.e. I have the point and want to send a “ray of infinite length” vertically through the topo. If I find zero intersections, my point is “outside”, and if I find one, I can easily calculate the distance (above or below).

I am pretty sure that this is a trivial thing with the means of VTK, but I simply do not find it among the many classes: maybe it is simply a question of wording - what to search for?!

Thanks for any helpful hint!

Take a look at this example: https://lorensen.github.io/VTKExamples/site/Cxx/PolyData/InterpolateTerrain/

Phantastic - thanks a lot!

I just tried it, and it looks like the second part of your example does exactly what I was asking for!

And even more: With some customer I am about to possibly try a similar thing with GEOTIFF files: with these I could then directly use the first part of your example, i.e. using the vtkProbeFilter.

So I was right indeed: knowing the right terms is the key to finding the right thing with search machines, but if you don’t know the right search words you have to rely on friendly and helpful people!

Regards, Cornelis

In addition to the above, I asked myself what would happen if the topo is not totally “clean” in the sense that it may have some local overlap: In real life topos are hardly ever 100% perfect in that sense, and with remote sensing technology it is not getting “cleaner”…

So I wanted to know if there is a way to “collect” all the intersection points along my line from the bottom to the top.

I did a test with a “topo” that consists of only two triangles on top of each other, not even connected - because I wanted to know what the vtkCellLocator is doing if there are actually 2 or more intersections along the line:

        double btmPt[3] = {0., 0., -500.},
           topPt[3] = {0., 0., 10000.},
           t,
           xyz[3],
           pcoord[3];
    int subId;
    while(cloc->IntersectWithLine(btmPt, topPt, 0.0001, t, xyz, pcoord, subId))
    {
        std::cout << "intersection found" << std::endl
                  << std::fixed << std::setprecision(2) << xyz[0]
                  << " " << std::fixed << std::setprecision(2) << xyz[1]
                  << " " << std::fixed << std::setprecision(2) << xyz[2] << std::endl
                  << "t " << std::fixed << std::setprecision(2) << t << std::endl
                  << std::fixed << std::setprecision(2) << pcoord[0]
                  << " " << std::fixed << std::setprecision(2) << pcoord[1]
                  << " " << std::fixed << std::setprecision(2) << pcoord[2] << std::endl
                  << "id " << subId << std::endl;
        btmPt[2] = xyz[2] + 0.0001;
    }

The effect was that I was “hitting” the two triangles indeed from bottom to top, and if I switched the order of the two triangles in the polydata it remained the same.

Now I don’t understand this cell locator code good enough to see if this is just “good luck” or by design: the fact that I can “iterate” through the intersections in an order from bottom to top!?

The alternative would be to do some kind of binary search: after every “hit” do another search “down” and one “upwards” until they all return false. Possible but a little bit more complicated.