[GRASS-user] surface interpolation with breaklines

GRASS GRASS user
1

Hello,

trying to generate a DEM from a photogrammetric vector input file
(elevation points, contour lines, and breaklines), I decided to run the
v.triangle addon. The main reason is the local topography has many
ridges and deep gulleys and I want to incorporate a breaklines vector
map in my interpolation process.

Unfortunately v.triangle output is a TIN, and I wish I could turn it
into a raster elevation map.

Would anyone help me figure out how I could perform this conversion,
given that I would rely on a linear interpolation method within
triangles faces (the result would be close to what r.surf.nnbathy does,
except nor this addon does support breaklines).

Exporting a vtk file out of grass I tried to explore capabilities of
other software like Paraview and MEshLab, but without results...

Thank you,
Vincent

2

W dniu 03.02.2014 18:55, Vincent Bain pisze:

Would anyone help me figure out how I could perform this conversion,
given that I would rely on a linear interpolation method within
triangles faces (the result would be close to what r.surf.nnbathy does,
except nor this addon does support breaklines).

Can you rasterize the output triangle 3d sides, and use that as an input for r.surf.nnbathy?

r.surf.nnbathy should create a structure similar to v.triangle product, since both use the same Shewchuk's "triangle" algorithm underneath.

You may try if using only the triangle node vertices alone wouldn't suffice, and extracting those from a 3d vector map should be very easy (v.to.points -n ?).

Maciek

--
Maciej Sieczka
http://www.sieczka.org

3

Hi Maciej,

Thank you for your reply.
I tested your suggestion but the process was /very/ long (v.trianlge
output was roughly 500 000 triangles). Perhaps a wiser solution would be
to v.to.point (i) contour lines, (ii) breaklines, then merge them in a
single raster an run r.surf.nnbathy.

Finally I found on the web a custom solution based on a nice little
python script named tin2raster, you can find it at the address bellow. I
tried to contact Antonio, the guy who wrote it, but did not get a
reply :

http://digilander.libero.it/antonioall/python_tin2raster.html

In my case it was very efficient and fast.
Do you think it could represent a potential new option that could be
applied to v.to.rast in order to handle 3d vector areas ?

Vincent.

Le dimanche 09 février 2014 à 12:53 +0100, Maciej Sieczka a écrit :

W dniu 03.02.2014 18:55, Vincent Bain pisze:

> Would anyone help me figure out how I could perform this conversion,
> given that I would rely on a linear interpolation method within
> triangles faces (the result would be close to what r.surf.nnbathy does,
> except nor this addon does support breaklines).

Can you rasterize the output triangle 3d sides, and use that as an input
for r.surf.nnbathy?

r.surf.nnbathy should create a structure similar to v.triangle product,
since both use the same Shewchuk's "triangle" algorithm underneath.

You may try if using only the triangle node vertices alone wouldn't
suffice, and extracting those from a 3d vector map should be very easy
(v.to.points -n ?).

Maciek

4

W dniu 09.02.2014 13:17, Vincent Bain pisze:

Perhaps a wiser solution would be to v.to.point (i) contour lines,
(ii) breaklines, then merge them in a single raster an run
r.surf.nnbathy.

There is no need to v.to.points isolines before rasterizing them for nnbathy. Unless you want to generalize the contour lines on purpose, e.g. to minimize the staircase artifacts in the output raster map.

As to breaklines, if these are linear features that indicate
surface discontinuity and have no elevation attribute, nnbathy has no
use of them (which you know already). If they are vector isolines, just rasterize them straight away. If these are vector 3d lines, you'd need to rasterize them interpolating their elevation at each cell of the output raster map, which v.to.rast can't do. At least it couldn't a few years ago when I tried it. This could probably be approximated with something like v.to.points -v -i dmax=<raster DEM resolution> + v.to.rast, but that's suboptimal.

Finally I found on the web a custom solution based on a nice little
python script named tin2raster, you can find it at the address
bellow. I tried to contact Antonio, the guy who wrote it, but did not
get a reply :

http://digilander.libero.it/antonioall/python_tin2raster.html

That tool utilizes Vect_tin_get_z
(http://grass.osgeo.org/programming6/tin_8c.html), which claims to be
able to perform such interpolation, not only along a linear 3d feature, but on faces as well. Nice.

I haven't tried the tool, but it seems it does exactly what you are after here.

In my case it was very efficient and fast. Do you think it could
represent a potential new option that could be applied to v.to.rast
in order to handle 3d vector areas ?

I don't know, but I'd like v.to.rast to provide such feature too.

Le dimanche 09 février 2014 à 12:53 +0100, Maciej Sieczka a écrit :

W dniu 03.02.2014 18:55, Vincent Bain pisze:

Would anyone help me figure out how I could perform this
conversion, given that I would rely on a linear interpolation
method within triangles faces (the result would be close to what
r.surf.nnbathy does, except nor this addon does support
breaklines).

Can you rasterize the output triangle 3d sides, and use that as an
input for r.surf.nnbathy?

r.surf.nnbathy should create a structure similar to v.triangle
product, since both use the same Shewchuk's "triangle" algorithm
underneath.

You may try if using only the triangle node vertices alone
wouldn't suffice, and extracting those from a 3d vector map should
be very easy (v.to.points -n ?).

--
Maciej Sieczka
http://www.sieczka.org

5

Le dimanche 09 février 2014 à 15:55 +0100, Maciej Sieczka a écrit :

W dniu 09.02.2014 13:17, Vincent Bain pisze:

> Perhaps a wiser solution would be to v.to.point (i) contour lines,
> (ii) breaklines, then merge them in a single raster an run
> r.surf.nnbathy.

There is no need to v.to.points isolines before rasterizing them for
nnbathy. Unless you want to generalize the contour lines on purpose,
e.g. to minimize the staircase artifacts in the output raster map.

oops, sorry you are right : I forgot to say I first turn isolines to
points (v.to.points, with the vertices option).

As to breaklines, if these are linear features that indicate
surface discontinuity and have no elevation attribute, nnbathy has no
use of them (which you know already). If they are vector isolines, just
rasterize them straight away. If these are vector 3d lines, you'd need
to rasterize them interpolating their elevation at each cell of the
output raster map, which v.to.rast can't do.

Here is my mistake. I have a 3d breaklines vector, and I expected
v.to.rast being able to (i) get vertices z values and (ii) interpolate z
between vertices...

At least it couldn't a few
years ago when I tried it. This could probably be approximated with
something like v.to.points -v -i dmax=<raster DEM resolution> +
v.to.rast, but that's suboptimal.

> Finally I found on the web a custom solution based on a nice little
> python script named tin2raster, you can find it at the address
> bellow. I tried to contact Antonio, the guy who wrote it, but did not
> get a reply :
>
> http://digilander.libero.it/antonioall/python_tin2raster.html

That tool utilizes Vect_tin_get_z
(http://grass.osgeo.org/programming6/tin_8c.html), which claims to be
able to perform such interpolation, not only along a linear 3d feature,
but on faces as well. Nice.

I don't know, but I'd like v.to.rast to provide such feature too.

understanding c libs content is definitely beyond my skills ! but in the
end, having a v.to.rast module able to handle 3d polylines and faces
would be really nice.

Best,
Vincent.

6

Thanks you Antonio,

in my present case, half a cell offset did not matter in any way, but I
notice your modification.
Your script might be well appreciated in the community....
My suggestion was to complete the v.to.rast module so it could handle 3d
geometries (3d polylines and faces), unfortunately I am only good at
scripting in bash (my knowledge of python and C are too poor to pretend
doing that).

Bye,
Vincent.

Le mardi 11 février 2014 à 09:05 +0100, antonioall@libero.it a écrit :
Hello Vincent,

Excuse me for my bad English.
Certainly you can cite me and my job.
Attention: the program considers as reference the angle in low at

left

of the cell and not the central point.
Reading the documentation of Grass found in net, I don't understand

what

reference is used for the cells.
If the reference of the cell is its central point, the lines 115 and

116

of the script become:
             x = xref + j * xres + 0.5 * xres
             y = yref + i * yres + 0.5 * yres
If you find useful information on the matter, please update me.

Good job,
Antonio Alliegro

Le dimanche 09 février 2014 à 16:24 +0100, Vincent Bain a écrit :

Le dimanche 09 février 2014 à 15:55 +0100, Maciej Sieczka a écrit :
> W dniu 09.02.2014 13:17, Vincent Bain pisze:
>
> > Perhaps a wiser solution would be to v.to.point (i) contour lines,
> > (ii) breaklines, then merge them in a single raster an run
> > r.surf.nnbathy.
>
> There is no need to v.to.points isolines before rasterizing them for
> nnbathy. Unless you want to generalize the contour lines on purpose,
> e.g. to minimize the staircase artifacts in the output raster map.

oops, sorry you are right : I forgot to say I first turn isolines to
points (v.to.points, with the vertices option).

>
> As to breaklines, if these are linear features that indicate
> surface discontinuity and have no elevation attribute, nnbathy has no
> use of them (which you know already). If they are vector isolines, just
> rasterize them straight away. If these are vector 3d lines, you'd need
> to rasterize them interpolating their elevation at each cell of the
> output raster map, which v.to.rast can't do.

Here is my mistake. I have a 3d breaklines vector, and I expected
v.to.rast being able to (i) get vertices z values and (ii) interpolate z
between vertices...

> At least it couldn't a few
> years ago when I tried it. This could probably be approximated with
> something like v.to.points -v -i dmax=<raster DEM resolution> +
> v.to.rast, but that's suboptimal.
>
> > Finally I found on the web a custom solution based on a nice little
> > python script named tin2raster, you can find it at the address
> > bellow. I tried to contact Antonio, the guy who wrote it, but did not
> > get a reply :
> >
> > http://digilander.libero.it/antonioall/python_tin2raster.html
>
> That tool utilizes Vect_tin_get_z
> (http://grass.osgeo.org/programming6/tin_8c.html), which claims to be
> able to perform such interpolation, not only along a linear 3d feature,
> but on faces as well. Nice.
>

> I don't know, but I'd like v.to.rast to provide such feature too.

understanding c libs content is definitely beyond my skills ! but in the
end, having a v.to.rast module able to handle 3d polylines and faces
would be really nice.

Best,
Vincent.

_______________________________________________
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user