[GRASS5] [bug #4248] (grass) r.surf.contour: Treats 0 as NULL and other archaic problems

[Maciek Sieczka via RT]

Hamish wrote:

> I really like the output this module creates. I never liked having
> to throw away data (by hand) to stop v.surf.rst from segmenting and aliasing
> itself along the high density contour lines verticies, especially in low
> lying areas & valley floors where is lots of empty space between contours.

Hamish,

I'm also having too much troubles with v.surf.rst to be able to produce a
decent DEM from elevation contours. Seeking for an alternative I found Pavel
Sakov's "nnbathy" - a natural neighbor interpolator, based on Jonathan
Shewchuk's "triangle". It is even better for me than r.surf.contour - it
handles non continous contour lines fine, as well as elevation points between
contour lines, huge areas of no data, rapid elevation gradients, and is darn
fast. The source code is available.

http://www.marine.csiro.au/~sakov/

It is quite complicated to make use of it for Grass. But since I have managed
to do it, everydoby else will too :).

What I was thinking about is, instead of fixing r.surf.contour, to bring
nnbathy to Grass. IMHO it is a good idea because it provides the same
functionality + much more (rapid gradients and clustered input/missing input
tolerant), and is efficient. Of course I'm not able do do it, but wanted to
trigger your interest.

So far, here's how I proceed to produce a DEM with nnbathy from Grass data and
import it back into Grass:

#contour lines into points
v.to.points

#add x,y,z columns into points vector
v.db.addcol 

# copy Z form layer1 to layer2 in the points vector
v.to.db

# populate x,y with points coordinates in layer2
v.to.db

# vector points into ASCII file
v.db.select -c points layer=2 column=x,y,z fs=" " > points.asc

# remove duplicates from points.asc; sort BTW
sort points.asc | uniq > points.asc.sort.uniq

# create the GRID on which nnbathy is supposed to interpolate
g.region as needed
r.mapcalc 'background=null()'
r.stats -g background nv="" > background.asc

# run nnbathy, eg.
./nnbathy -i points.asc.sort.uniq -o background.asc > dem.grd

# replace "nan" in nnbathy output with something Grass understands, we'll
# filter this out later:

cat dem.grd | sed 's/nan/999999.0/g' > dem.grd.nonan

# import it into Grass vector:
v.in.ascii -t -z input=dem.grd.nonan output=dem_vect fs=space x=1 y=2 z=3

# vector into raster:
v.to.rast input=dem_vect output=dem use=z

# change 999999.0 into null
r.null

I wish it wasn't all that complicated. Anybody interested in implementing
nnbathy in Grass :) ? I believe it is worth it.

Maciek


-------------------------------------------- Managed by Request Tracker