[GRASSLIST:1355] Re: correction to shade.rel.sh

[Markus Neteler]

On Tue, Jan 02, 2001 at 04:54:59PM -0800, David Finlayson wrote:
> It is possible in perfectly flat areas to get a
> division by zero error in the r.mapcalc portion of the
> shade.rel.sh script.  For some reason, the test for
> this condition fails on my copy of grass5beta10 but I
> can't follow the logic to see just why.  I added a
> test for the null condition and this seems to work. 
> Maybe someone has a better idea?  Here is a clip with
> my addition:
> 
> r.mapcalc << EOF
> shade = eval( \\
>  x=($elev[-1,-1] + 2*$elev[0,-1] + $elev[1,-1] \\
>    -$elev[-1,1] - 2*$elev[0,1] -
> $elev[1,1])/(8.*ewres()) , \\
>  y=($elev[-1,-1] + 2*$elev[-1,0] + $elev[-1,1] \\
>    -$elev[1,-1] - 2*$elev[1,0] -
> $elev[1,1])/(8.*nsres()) , \\
>  slope=90.-atan(sqrt(x*x + y*y)), \\
>  a=round(atan(x,y)), \\
>  a=if(isnull(a),1,a), \\  <-- line I added
>  aspect=if(x||y,if(a,a,360.)), \\  <-- test that fails
>  cang = sin($alt)*sin(slope) + cos($alt)*cos(slope) *
> cos($az-aspect), \\
>  if(cang < 0.,0.,100.*cang), \\
>  if(isnull(cang), 22, 100.*cang))
> EOF
> 

David,

many thanks for your fix. I have updated the CVS.

Markus Neteler