Hi, I am revisiting a previous script I made in producing "Multidirectional, oblique-weighted, shaded-relief". ArcGis has this script: http://blogs.esri.com/Support/blogs/mappingcenter/archive/2008/10/07/updated-hillshade-toolbox.aspx
The commads in ARC 6.0.1 GRID by Rober Mark (USGS) are as follows: shade225 = hillshade(hawaii, 225, 30, shade, 5) shade270 = hillshade(hawaii, 270, 30, shade, 5) shade315 = hillshade(hawaii, 315, 30, shade, 5) shade360 = hillshade(hawaii, 360, 30, shade, 5) h00 = resample(hawaii, 1000) h01 = focalmean(h00) h02 = focalmean(h01) h03 = focalmean(h02) asp = aspect(h03) asp1 = con(isnull(asp), 293, asp) w225 = sqr(sin((asp1 - 225) div deg )) w270 = sqr(sin((asp1 - 270) div deg)) w315 = sqr(sin((asp1 - 315) div deg )) w360 = sqr(sin(asp1 div deg )) setcell minof temp = w225 * shade225 + w270 * shade270 + w315 * shade315 + w360 * shade360 shade4 = int(temp div 2) I tried porting in GRASS shell: g.region rast=$GIS_OPT_input -p # Compute hillshade at azimuth 225, 270, 315 and 360 at 30 degrees sun illumination angle r.shaded.relief map=$GIS_OPT_input shadedmap=shade225 altitude=30 azimuth=225 zmult=1 scale=1 r.shaded.relief map=$GIS_OPT_input shadedmap=shade270 altitude=30 azimuth=270 zmult=1 scale=1 r.shaded.relief map=$GIS_OPT_input shadedmap=shade315 altitude=30 azimuth=315 zmult=1 scale=1 r.shaded.relief map=$GIS_OPT_input shadedmap=shade360 altitude=30 azimuth=360 zmult=1 scale=1 # Resample DEM map r.neighbors in=$GIS_OPT_input out=res1 method=average size="$GIS_OPT_window" # compute aspect map r.slope.aspect elevation=res1 format=degrees prec=float aspect=aspect zfactor=1.0 min_slp_allowed=0.0 # compute weights 225, 270, 315 and 360 r.mapcalc "w225 = sin((aspect - 225)*sin(aspect - 225))" r.mapcalc "w270 = sin((aspect - 270)*sin(aspect - 270))" r.mapcalc "w315 = sin((aspect - 315)*sin(aspect - 315))" r.mapcalc "w360 = sin((aspect)*sin(aspect))" #compute weighted r.mapcalc "weightedshade = (((w225 * shade225) + (w270 * shade270) + (w315 * shade315) + (w360 * shade360))/2)" r.colors map=weightedshade color=grey But I get not so good results, here's the output using nc_sp_08's elev_ned_30m (X0=MDOW, X1= default r.shaded.relief module) http://farm5.static.flickr.com/4148/5037622535_113924312b_b.jpg Anything wrong with my commands? -- cheers, maning ------------------------------------------------------ "Freedom is still the most radical idea of all" -N.Branden wiki: http://esambale.wikispaces.com/ blog: http://epsg4253.wordpress.com/ ------------------------------------------------------ _______________________________________________ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
