[GRASS-user] Comparing shaded relief and ESRIs hillshade

[Patrick S.]

A little feedback of current findings for r.shaded.relief in GRASS6.4.2 
in case anybody should have the same problem as we just had. We needed 
to create comparable results between the hillshade algorithm (ArcGIS) 
and the shaded.relief algorithm (GRASS) and compared the output. Values 
of shaded relief differ essentially and are generally smaller.

Based on 
http://help.arcgis.com/en/arcgisdesktop/10.0/help/index.html#//009z000000z2000000.htm 
I could correct the calculation to create the same results. Replacing 
the calculation in shaded relief with the syntax below creates an 
identical range of values. However, this can also been done through 
upscaling the output of r.shaded.relief with a factor of 2.549999108 
(shaded.relief/hillshade=0.392157). Would it make sense to integrate 
this as an option?

The previous reported problem of negative values in the results are 
still existing. "if(cang<0,0,cang)" should remove those- doesn't it? 
Comments on this are very welcome.

patrick


r.mapcalc << EOF
$elev_out = eval( \\
  x=($zmult*$elev[-1,-1] + 2*$zmult*$elev[0,-1] + $zmult*$elev[1,-1] \\
    -$zmult*$elev[-1,1] - 2*$zmult*$elev[0,1] - 
$zmult*$elev[1,1])/(8.*ewres()*$scale) , \\
  y=($zmult*$elev[-1,-1] + 2*$zmult*$elev[-1,0] + $zmult*$elev[-1,1] \\
    -$zmult*$elev[1,-1] - 2*$zmult*$elev[1,0] - 
$zmult*$elev[1,1])/(8.*nsres()*$scale) , \\
  slope=atan(sqrt(x*x + y*y)), \\
  a=round(atan(x,y)), \\
  a=if(isnull(a),1,a), \\
  aspect=if(x!=0||y!=0,if(a,a,360.)), \\
  zenith=90-$alt, \\
  azimuth=360-$az+90, \\
  cang=255*((cos(zenith)*cos(slope)) + (sin(zenith)*sin(slope)* 
cos($az-aspect))), \\
  if(cang<0,0,cang), \\
  if(isnull(cang), null(), cang))
EOF