#2359: r.stream.distance with a large map
-------------------------------+--------------------------------------------
Reporter: hcho | Owner: grass-dev@…
Type: defect | Status: new
Priority: normal | Milestone: 7.1.0
Component: Raster | Version: svn-trunk
Keywords: r.stream.distance | Platform: All
Cpu: All |
-------------------------------+--------------------------------------------
r.stream.distance creates FCELL output maps, but sometimes with a large
map, the precision of FCELL is not enough. For example, to create a
longest flow path map:
{{{
r.stream.distance -o stream_rast=outlet direction=drain method=upstream
distance=flus
r.stream.distance -o stream_rast=outlet direction=drain method=downstream
distance=flds
r.mapcalc expression="flusds=flus+flds"
eval `r.info -r flusds`
r.mapcalc expression="lfp=if(flusds>=$max-0.1, flusds, null())"
}}}
lfp should be a single connected stream representing the longest flow
path, but because of rather big floating-point errors caused by FCELL, it
can have more than one segment of streams.
In line 434 in distance_calc.c:
{{{
cur_dist = tmp_dist + G_distance(...);
}}}
As tmp_dist becomes large, I get the following result:
{{{
cur_dist = tmp_dist + G_distance(...);
double double_dist = (double)tmp_dist + G_distance(...);
fprintf(stderr, "float: %f = %f + %f\n", cur_dist, tmp_dist,
G_distance(...));
fprintf(stderr, "double: %f = %f + %f\n", double_dist, tmp_dist,
G_distance(...));
Output:
float: 117685.031250 = 117642.601562 + 42.426407
double: 117685.027969 = 117642.601562 + 42.426407
}}}
These errors get accumulated and can cause unexpected problems. I suggest
to change the default type of output maps to DCELL and add a flag (-f) for
FCELL outputs, if needed.
--
Ticket URL: <http://trac.osgeo.org/grass/ticket/2359>
GRASS GIS <http://grass.osgeo.org>