[long, but hopefully interesting & check out the link to the new surface trials wiki page in the middle]
Frank Broniewski wrote: > I am trying to create a high resolution dem from contour lines Until > now all my tries where not successful. At first I tried r.surf.contour, > but since my interpolation region is not rectangular and the contours > are not evenly distributed (rough terrain), the result was > unfortunately not usable really? how did it go wrong? For me it has always worked beautifully, and I'm working in mountains where the terrain is about as rough as you can get. approx 16000x16000 cell region size without problems. > ( but it took around 7 days to compute, that > alone was already impressive ;) ) yes, the search algorithm is very inefficient for areas with a lot of flat lands where it has to search for a long time to find the nearest contour line. (so in that way it actually likes rough terrain) > My contour map is a combinatin of a national contour line map (5m vert. > resolution) and contours from SRTM with 20m vert. resolution. are you sure that the two are in agreement? do the contour lines cross or does the hole in the srtm avoid that? > I created a "hole" in the srtm contours for the national contour map > and patched both together to avoid large gaps with no height values > (mostly for r.surf.contour) um, why not use srtm directly instead of doing srtm -> contours -> dem? srtm's native form is raster, not contour lines. ... and then patch in the hi-res r.surf.* result in the middle. > My region is 17.000 x 13.000 cells wide (5m horiz. resolution). So my that is getting big, but shouldn't be too bad. except perhaps for v.surf.rst (and r.fillnulls which uses it). hmmm... I wonder if it would speed up r.fillnulls & save lots of memory if we added the no-table and no-topology -b and -z flags to the r.to.vect step in it? > current approach is to use small regions (2000 x 2000) to calculate > small subsets of the dem. Because of the algo used by v.surf.rst to > create the dem the neighboring tiles do have different height values > calculated at the borders. So it was not possible to just create the > tiles and patch them together. I suspect this is making it more complicated than it has to be. Perhaps try running v.generalize before v.surf.rst? > My next approach used an overlapping of 20 cells for each tile and a > moving window average to calculate the mean of the overlapping tiles. > The result was quite good, but the moving window approach resulted in > null values where one tile ended and the other started .. > (similar to the slope and aspect maps, where there is a 1 cell null > border around the map in comparison to the input dem). these are caused because it takes into account all the cells around it, and at the edges some of those are unknown (out of region is equivalent to NULL). > Unfortunately I was not able to remove the null values satisfactorily. > r.fillnulls fails because of the large region, and r.resamp.rst does > the job not very well. The stripes are still visible, though filled > with values. the usual way would be to use one of the r.surf.* or r.resamp* modules, then use r.patch with the old map in the foreground and the interpolated map in the background filling in any gaps. > When calculating a derivate from the dem, like aspect, the errors from > filling null values are quite obvious. I take it the vertical jumps between tile edges are rather large? > So to make my long text short: Is there a technique to combine two or > more raster dem with (or without overlapping) with good > transition/intersection (don't know the correct word) between two > tiles? usually it would be a matter of making the overlaps big enough, and then using r.series or something to patch them together taking the average. I guess it would be better bias away from map edges somehow? Not sure how to do that beyond cropping a set number of cells with an inverse r.grow. but again 13k * 17k cells is not that huge so you shouldn't need to tile. (except for r.fillnulls/v.surf.rst). maybe r.surf.idw(2) or r.surf.nnbathy from wiki addons helps here? > If necessary I can illustrate my efforts by creating a web page or > similar ... I am curious to know what the problem with r.surf.contour was. Usually it does a great job as long as you take a moment to work around the integer and NULL/0 issues. I decided to do some tests using different methods on the 1m contour lines from the NC sample dataset. Results are here: http://grass.osgeo.org/wiki/Contour_lines_to_DEM John Tate wrote: > and then patched (r.patch) them all together. The patching should > average out any differences. r.patch doesn't average, it uses the cell from the first-listed raster map which contains a non-NULL cell. Later-listed rasters fill in the gaps until there are no more gaps or no more maps. r.series will merge overlapping maps using an average though. > I then cropped out each 1km tile (1000x1000). > This was done so that the 1km tiles could be combined for specific > areas by different people (e.g. only a 4kmsq area for academic 'a' or a > 6kmsq area 2km away for academic 'b'). why not just use g.region for that? Pablo Carreira: > I had good results interpolating 5m contours with v.surf.rst, but I > spent some time adjusting the parameters to get a smooth surface. v.surf.rst does not like contours very much. The quadtree method concentrates on areas with high density of points, and on contour lines you have high density along the isoline but large gaps in between. see this wiki page for tips on using v.generalize to even things out a bit. (especially Helena's linked page on it) http://grass.osgeo.org/wiki/RST_Spline_Surfaces#Working_with_contour_lines Nick Cahill: > Just to refer to a previous question - when I have had to make > relatively high resolution DEMs (much smaller than yours, only 5000 x > 3000 cells), I found it most effective to use Arc/Info to create a TIN > from contour lines and points, then rasterize that and import the > raster into GRASS, rather than having GRASS interpolate from contours. both v.surf.rst with generalized contours or r.surf.contour should create much nicer surfaces than a rasterized TIN. Both visually and hydro- logically. Sure, it'll take longer but there is a lot more information in those vertices than a simple linear interpolation between them and their nearest neighbors can expose.. why not fill in the gaps between the points with curves instead of a flat surface? Use whatever works for you, but do examine the results in the above "Contour_lines_to_DEM" wiki page. May I suggest to try r.in.xyz -> r.surf.nnbathy? (or nnbathy directly) Since you are mentioned in the source code thank-yous I'm guessing you are already familiar with it. :) The module is available from GRASS wiki addons. or for massive point cloud data perhaps try -> r.in.xyz -> r.mapcalc *= 100000 -> r.surf.contour -> r.mapcalc /= 100000 ? ? or v.surf.bspline (newly fixed by Markus M) > I was never able to get the parameters right in *.surf.rst, and > processing times were very long. True & true. "auto-tune" has been a wish for a while. it probably would take a while to run though. maybe the tips on the RST Spline Surfaces wiki page along with the ones on the main help page can be useful? The results can be very nice once you get it tuned right though. > Arc/Info does the job very quickly and effectively, and doesn't end up > with overshoots and depressions, which were a problem for me. The Achilles heel of splines is that they do not deal well with hard corners (base of cliffs, etc). You get a ringing effect, http://en.wikipedia.org/wiki/Ringing_artifacts you can play with the tension to minimize it, or avoid using it in areas with sharp changes in topography. (or oversample there) > I wish this were an option in GRASS. I would also like to be able to > work with other vector-based surface models in GRASS. Helena has posted some examples recently, maybe she can say a few words about it. but really, try r.in.xyz (or v.to.rast) + r.surf.nnbathy should do the trick. nnbathy's linear ("alg=l") algorithm creates a rasterized TIN if you really want one. (and it's lightning fast) ... and the winner is ... (drum roll please) ... r.surf.contour. At least IMO. there are some "thatching" artifacts in it, but look closely at the main river valley running up the center. it is the only method that creates it correctly. also it's the only one which is able to recreate the road running along the upper left ridge in the sample dataset. Jed Frechette wrote: (2 weeks ago) > I'm curious why you say TINs are used at the "cost of resolution"? What > if I have a surface sampled by a dense irregularly spaced point data > set, e.g. a lidar point cloud. If I create a TIN using a delaunay > triangulation of those points the resulting model will exactly pass > through all of the original data points with linear interpolations > between. and be limited to a simple linear approximation in between points. > As a result it resolves all of the detail in the original data. well, fundamentally a TIN *is* the original data, and no more. the onus is on what reads/samples/rasterizes the TIN to do the linear interpolation at run-time. > In contrast, if I use a raster map to model the surface it will > inevitably lead to some degree of aggregation as a function of the > interpolation method I choose. as long as the aggregation method is statistically sound, that sort of filtering is often very useful and can help reduce the dataset to a manageable size. (is it really useful to preserve 50 returns per sq m?) > The same is true if I want to incorporate linear features such as > breaklines into the model. this is something I am interested in hearing more about. AFAIK there isn't much in GRASS which can do that. But we can change that. > Of course, the question of which type of model is more useful for a > given application is entirely separate. Absolutely. I'd like to hear what the folks doing hydrological modelling etc think & which artifacts the care about. What aspects important for LIDAR? (yeah, it's just a tool so that's an open question..) Especially in light of the new MDF r.watershed improvements and new r.stream.* addon modules. In theory, for an identical grid of points representing a surface, using triangles instead of grid boxes will match the "true" natural surface much better than a q-bert style collection of boxes (I'm not sure it is best to think about rasters as full-cell blocks, vs. just representing where the model surface intersects the cell centers, however). But in practice a raster surface is going to be vastly more computationally efficient, both in terms of processing time and data storage. So you can get a lot more done with them (at the cost of initial interpolation time) due to the much lower overhead of dealing with them; and by using matrices you open the door to lots of neat mathematical tools. The main usefulness of TINs AFAICS are for setting up triangular-mesh model domains where it is important to have spatial resolution vary by orders of magnitude across the region. (to place computational focus on areas of highest interest) Whether a quad-tree mesh will be more efficient there or not probably depends on the nature of the model. regards, Hamish _______________________________________________ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user