Psmap flooding example

• List of all ps.map example scripts

About

This is a simple ps.map script demonstrating variable text & queries.

It is based on this mailing list post and considers a hypothetical dynamic flood risk warning map generator.

It builds on the Spearfish catchments generated by the r.watershed help page example. All you need to specify is catchments at risk, the script generates the rest.

Screenshots

screenshot of whole catchment warning: http://bambi.otago.ac.nz/hamish/grass/screenshots/flood_warning_basins.png

---

screenshot of low-lying areas warning: [http://grass.ibiblio.org/grass63/screenshots/images/flood_warning_lowlying.png http://grass.ibiblio.org/grass63/screenshots/images/flood_warning_lowlying_small.jpg]

---

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.

Next we designate the low lying areas. Ideally this would come from a real-time hydrologic model, but as this is an example we just make something up:

# create danger polygons for each catchment
#  (homebrew example of predictor for at-risk low-lying 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 convert 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"?