Psmap flooding example

From GRASS-Wiki
Revision as of 12:14, 18 November 2007 by HamishBowman (talk | contribs) (init)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigation Jump to search

About

This is simple ps.map script demonstrating variable text & queries. It is based on this mailing list post.


In ps.map there is the vlegend instruction which does this automatically.

It builds on the Spearfish catchments generated in the r.watershed help page example. All you need to specify is catchments at risk, it does the rest.

Screenshots

screenshots of output for whole catchment warning and lowlying areas warning: http://bambi.otago.ac.nz/hamish/grass/screenshots/flood_warning_basins.png http://bambi.otago.ac.nz/hamish/grass/screenshots/flood_warning_lowlying.png

screenshot of Spearfish catchments: http://grass.itc.it/screenshots/raster.php

The script

First we need to generate polygons for low lying ares from each catchment created with the r.watershed help page example:

# create danger polygons for each catchment
#  (homebrew example of predictor for at-risk lowlying areas)
r.cost -k in=slope out=slope_cost_from_rivers start_rast=rwater.course
r.colors slope_cost_from_rivers col=bcyr -e

# choose some threshold. this could be based on hydrologic model output.
#  for now just set for lower half of slope-costs in the map
eval `r.univar -g slope_cost_from_rivers`
r.mapcalc "flood.zones = if(slope_cost_from_rivers < $mean, \
     rwater.basin, null() )"

# convert to vector areas
r.to.vect -s in=flood.zones out=flood_zones feature=area
v.db.dropcol map=flood_zones column=label
v.db.renamecol map=flood_zones column=value,catchment

Next we set the dynamic data, either automatically from model output or though interactive selection & highlight of areas (e.g. a script using the HTMLMAP driver, 'd.what.vect -xft', or d.ask)

# Here we say catchment #6 and #12 are at risk.

FLOOD_BASINS=6,12
TIMESTAMP=`date`


The rest is a static template. Here we pass the command through the shell but the input file could be written to a file (taking care to process environment variables) and then passed to ps.map's input option.

# convert catchment list to SQL query
FLOOD_BASIN_SQL=`echo $FLOOD_BASINS | \
  sed -e 's/^/catchment = /' -e 's/,/ OR catchment = /g'`

# generate PostScript image
ps.map out=flood_warning.ps << EOF
raster elevation.dem
# example of dashed lines
vlines roads
  style 000011
  where label ~ 'highway'
  width 0.25
  label Highways
  end
vlines rwater_course
  color grey
  label Water courses
  end
vareas rwater_basins
#change to "flood_zones" vector map for lowlying areas instead of
#  entire catchment
  where $FLOOD_BASIN_SQL
  color red
  fcolor red
  width 2
  pat $GISBASE/etc/paint/patterns/diag_up.eps
  label Catchements subject to flood risk
  end
vlegend
  end
text 100% -3% Warning issued $TIMESTAMP
  font Helvetica-BoldOblique
  fontsize 12
  ref right
  end
end
EOF

Finally we conver the output PostScript image to PDF and PNG formats.

# convert to PDF
ps2pdf flood_warning.ps
# convert to PNG
pstoimg -antialias -aaliastext -out flood_warning.png -scale 1.3 \
  -crop tb flood_warning.ps

The PDF output looks the best.


Rendering issues

  • EPS hatch pattern renders in gv and PNG with black boxes behind lines??
These do not appear in PDF + acrobat reader.
  • EPS hatch patterns are tiled so fill lines have slightly noticeable breaks.
  • We need to add something to allow skipping vectors in ps.map's vlegends
ideas: labelskip [y|N], test for label= or 'none' or 'skip'?