[GRASSLIST:1328] correction to shade.rel.sh

[David Finlayson]

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 Finlayson
david_finlayson at yahoo.com
University of Washington 
Box 351310 
Seattle, WA   98195 - 1310

__________________________________________________
Do You Yahoo!?
Yahoo! Photos - Share your holiday photos online!
http://photos.yahoo.com/