| 123456789101112131415161718192021222324252627282930313233343536373839 |
- #!/bin/sh
- # Test script for r.viewshed based on a synthetic DEM
- # RUN THIS IS ANY GRASS LOCATION, e.g. NC or Spearfish
- # create first hemisphere
- g.region n=1000 s=0 w=0 e=1000 -p res=1
- r.mapcalc 'disk.15031=if(sqrt((col() - 500)^2 + (500 - row())^2)<500,sqrt((col() - 500)^2 + (500 - row())^2),null())'
- r.mapcalc 'hemisphere1=500 * sin(acos (disk.15031/500))'
- # create second hemisphere
- g.region n=500 s=0 w=0 e=500 -p res=1
- r.mapcalc 'disk.14947=if(sqrt((col() - 500)^2 + (500 - row())^2)<500,sqrt((col() - 500)^2 + (500 - row())^2),null())'
- r.mapcalc 'hemisphere2=500 * sin(acos (disk.14947/500))'
- g.remove --q rast=disk.14947,disk.15031
- # merge both
- r.mapcalc "hemisphere=hemisphere1 + hemisphere2"
- d.mon x0
- d.rast hemisphere
- # run r.viewshed
- r.viewshed hemisphere out=hemisphere_viewshed coord=250,250 max=1000000 obs=100 mem=2000 --o
- r.shaded.relief hemisphere --o
- d.his h=hemisphere_viewshed i=hemisphere.shade
- # compare to r.los
- r.los hemisphere out=hemisphere_los coord=250,250 max=1000000 obs=100 --o
- d.mon x1
- d.his h=hemisphere_los i=hemisphere.shade
- r.mapcalc "hemisphere_diff = hemisphere_viewshed - hemisphere_los"
- r.colors hemisphere_diff color=differences
- d.mon x2
- d.rast.leg pos=80 map=hemisphere_diff
- nviz hemisphere col=hemisphere_viewshed
|
