[GRASS5] Re: [GRASSLIST:8999] r.statistics vs v.rast.stats vs ??

[Hamish]

> I am trying to calculate statistics such as mean elevation, slope,
> maximum drainage area etc for a very large number of small plots. The
> plots are in a vector file and I have been trying to use v.rast.stats
> on the various rasters that represent elevation, slope, drainage area,
> etc.
> 
> The problem is that the rasters are very high resolution, the plots
> are quite small, and I have more than 600 to process. If you look at
> the code for v.rast.stats, it converts each plot into a raster to use
> as a mask and then runs r.univar on the entire map (now mostly null
> values), once for each statistic. In my situation it can take several
> minutes for each plot. In other words hours for each of the above
> rasters that I am interested in and more than 2 days for the whole
> set.
> 
> Is there a more efficient way of calculating these statistics? I have
> tried to use r.statistics and can't make sense of the inputs and
> outputs. Alternatively, if I could control the region size in the
> critical part of the v.rast.stats script I might be able to improve
> the speed of the script.


for i in $CATSLIST ; do
 echo "Processing category $i ($CURRNUM/$NUMBER)"
 g.remove MASK > /dev/null 2> /dev/null
 r.mapcalc "MASK=if(${VECTOR}_${TMPNAME} == $i, 1, null())" 2> /dev/null
 #n, min, max, range, mean, stddev, variance, coeff_var, sum
 n=`r.univar -g $RASTER | cut -d'=' -f2 | tr '\n' ' ' | cut -d' ' -f1`
 min=`r.univar -g $RASTER | cut -d'=' -f2 | tr '\n' ' ' | cut -d' ' -f2`
 max=`r.univar -g $RASTER | cut -d'=' -f2 | tr '\n' ' ' | cut -d' ' -f3`
 range=`r.univar -g $RASTER | cut -d'=' -f2 | tr '\n' ' ' | cut -d' ' -f4`
 mean=`r.univar -g $RASTER | cut -d'=' -f2 | tr '\n' ' ' | cut -d' ' -f5`
 stddev=`r.univar -g $RASTER | cut -d'=' -f2 | tr '\n' ' ' | cut -d' ' -f6`
 variance=`r.univar -g $RASTER | cut -d'=' -f2 | tr '\n' ' ' | cut -d' ' -f7`
 cf_var=`r.univar -g $RASTER | cut -d'=' -f2 | tr '\n' ' ' | cut -d' ' -f8`
 sum=`r.univar -g $RASTER | cut -d'=' -f2 | tr '\n' ' ' | cut -d' ' -f9`


a couple of things,

1) r.univar is running 9 times when it only needs to be run once:
eval `r.univar -g $RASTER`

Then $cf_var gets renamed $coeff_var, needs to be carried through.


2) if g.region zoom= respects MASK, we can do:
g.region save=${TMPNAME}_rastres

for i in $CATSLIST ; do
 echo "Processing category $i ($CURRNUM/$NUMBER)"
 g.remove MASK > /dev/null 2> /dev/null
 r.mapcalc "MASK=if(${VECTOR}_${TMPNAME} == $i, 1, null())" 2> /dev/null

+g.region zoom=$RASTER

 #n, min, max, range, mean, stddev, variance, coeff_var, sum
 eval `r.univar -g $RASTER`

+g.region ${TMPNAME}_rastres

 ...
done

g.remove region=${TMPNAME}_rastres



but g.region zoom= can take a _long_ time so maybe no advantage/worse?

Would need to check/reset NSRES and EWRES to raster's values or 
is that redundant? (what about g.region -a?) Some checking needs
to be done that we're not screwing up the cells.. at minimum
a before/after check to see that the same answer comes out both
ways.


3) rm twice??
cleanup()
{
 #restore settings:
 g.region region=${TMPNAME}
 g.remove region=${TMPNAME} > /dev/null

twice??

4) isn't  -p > /dev/null  circular?
g.region nsres=$NSRES ewres=$EWRES -ap > /dev/null


5) cleanup should happen before   echo "Done."


Hamish