HOME

GRASS 5.7: Vector networking tutorial


Tutorial HOME | Table of contents
README
GRASS 5.7 is currently under development. In case the examples described here do not work properly, you are kindly invited to send us further examples and/or code bugfixes/enhancements.

Example: shortest path vector networking with d.path
[for this exercise we use the Spearfish data set from the GRASS web site]

A module based on the DGLib vector network library is 'd.path' (see screenshot which calculates shortest path on a vector network map. Traveling costs may be either line lengths, or costs saved as attributes in a database table (supported are cost assignments for both arcs and nodes).

We want to use the 'roads' map of the Spearfish demo location (if needed, convert to 5.7 format with 'v.convert'). A simple example is to use the vector length as costs:

d.vect roads
d.path map=roads
More realistic shortest path queries may consider the maximum speed according to the road type. If desired to use different traveling costs for the two directions of a road, forward and backward costs are read from two different table columns. In our example we simply assume the same costs (same table column) for both directions. Node that you have to hit a node with the mouse in the GRASS monitor to select start/end point when using d.path. As the column 'travelcost' is not yet present in the table 'road' we can use StarOffice or OpenOffice to add it (see here). The we can continue and insert values according to the road type (be sure to quote strings when using for where statement):
echo "select * from roads" | db.select
... shows that the travelcost column is filled with 0.0. To assign new travel costs, we run:
echo "update roads set travelcost=5 where cat=1" | db.execute
echo "update roads set travelcost=20 where cat=2" | db.execute
echo "update roads set travelcost=40 where cat=3" | db.execute
echo "update roads set travelcost=60 where cat=4" | db.execute
echo "update roads set travelcost=80 where cat=5" | db.execute
echo "select * from roads" | db.select
Then we can run the 'd.path' with our assigned travel costs:
d.path map=roads afcol=travelcost abcol=travelcost
To see the forward directions of the vector lines, use the 'dir' for the 'display' parameter in 'd.vect'. It will add small arrows indicating the forward directions when plotting the vector map.

Reachability of Schools: v.net.iso
Reachability of Schools (find better name)

We need a street map and a map of school locations. A nice, free vector data set is FRIDA, the vector map of Osnabrück, a city in northern Germany. Additionally you can download the related ASTER DEM at 30m resolution (ARC ASCII GRID raster format, reprojected to Gauss-Krüger, 700kb).
After downloading the dataset (SHAPE files) and extracting the related GRASS location (corrected 1 Oct 2004), we can start. Enter GRASS 5.7 with that osnabrueck LOCATION and user1 mapset.

Translation aid for the map names:

1. Import of SHAPE maps

The SHAPE registration with 'v.external' is one option, but only intended for visualization. Since we want to work with the data, we import them into GRASS with 'v.in.ogr'.

#line/point layers:

 
v.in.ogr dsn=/ssi0/ssi/neteler/grassdata/osnabrueck/maps output=strassen layer=strassen
v.in.ogr dsn=/ssi0/ssi/neteler/grassdata/osnabrueck/maps output=poi layer=poi

#area layer:

 
v.in.ogr dsn=/ssi0/ssi/neteler/grassdata/osnabrueck/maps output=gewaesserflaechen layer=gewaesserflaechen
-> WARNING: 2. centroid found in area 30...etc
-> "WARNING: Area boundaries are not cleaned by this module. Run v.clean on
         imported vector on new map ."
The module 'v.in.ogr' imports non-topological GIS formats into the topological GRASS. This is relevant for areas, where boundaries are duplicated in the SHAPE file. As suggested, we run 'v.clean' to make the imported map fully topological:
 
v.clean input=gewaesserflaechen output=gewaesserflaechen_clean tool=rmdupl,bpol
No more errors appear.

1.b) Add new field to be able to assign two tables later - DIRTY

Later we want to add two tables to the resulting map. For that we have to create a new field (fields are internal ID columns to link attribute tables to GRASS vector maps). By default a vector map has only one field. To add a second field column in the map, we run:
 
v.category poi out=poi_2f field=2 op=add
#verify:
v.category poi_2f field=1,2 op=print  
-> field 1 == field 2

2. Extract school places from POI map (but the map with two fields)

 
v.extract in=poi_2f out=schools type=point where="poiName='Schule' or poiName='Schulen'"
   #TODO: Fix Usulaschule - SQL regex 
v.info schools
-> one dblink
v.category schools field=1,2 op=print
-> 2 fields
d.erase
d.vect strassen
d.vect schools disp=attr attr=poiName bgcolor=white bcolor=black
d.vect schools col=red icon=basic/diamond
Now we have the schools extracted from the poi map.

3. Patch schools into streets map

 
v.patch in=strassen,schools out=streets_schools
v.info streets_schools
-> dblinks = 0
d.erase
d.vect streets_schools
d.vect streets_schools type=point col=red
But the tables are not linked yet (v.patch won't do it for various reasons). So we add the related two database tables to the patched map. GRASS 5.7 allows us to connect two tables to the same map. NOTE: key must be an integer column!
Now we link the patched map to the DBF tables generated during the 'v.in.ogr' import (which are in mapset/dbf/), we do not link to the original tables:
 
db.connect dr=dbf database='$GISDBASE/$LOCATION_NAME/$MAPSET/dbf/'
db.describe -c strassen
ncols:6
Column 1: cat
Column 2: strShapeID
Column 3: strID
Column 4: strTypID
Column 5: strSpuren
Column 6: strEbene
-> the 'cat' column was added by 'v.in.ogr'
 
v.db.connect streets_schools dr=dbf data='$GISDBASE/$LOCATION_NAME/$MAPSET/dbf/' \
  table=strassen field=1 key=cat

v.db.connect streets_schools dr=dbf data='$GISDBASE/$LOCATION_NAME/$MAPSET/dbf/' \
  table=schools field=2 key=cat
v.db.connect -p streets_schools

Display the patched map with the linked two tables:

 
d.erase
d.vect streets_schools
Display attributes of second table (lfield=2):
 
d.vect streets_schools type=point col=red disp=attr attrcol=poiName lfield=2
d.what.vect streets_schools

Connecting schools to street network

Nodes must be on the lines network for vector analysis. To fulfill this request, we have to generate the connection between the schools and the street network by generating a connection map. This connection map (vector lines between the school locations and the nearest street) we will patch into the schools_streets map created above.

First we generate the connection map between schools and closest streets:

v.distance -p from=schools to=strassen output=streets_schools_connect upload=dist column=dist

d.erase
d.vect strassen col=grey
d.vect schools col=red
d.vect streets_schools_connect col=blue

Then we merge these line connections into the schools_streets map (we have to snap nodes and break street lines at intersection of the connecting lines):

v.patch in=streets_schools,streets_schools_connect out=schools_net

#note: next will also remove the few bridges in the streets_schools map...
v.clean in=schools_net out=schools_net_clean tool=snap,break thresh=1
g.remove vect=schools_net

#since g.rename doesn't exist:
g.copy vect=schools_net_clean,schools_net
g.remove vect=schools_net_clean

d.erase
d.vect schools_net

Now the map containing schools and streets is prepared for network analysis.

Vector network analysis: Reachability of selected schools

First we want to see all streets and the categories (IDs) of the schools to select some of them:
d.vect schools_net
#show all schools with their cats:
d.vect schools_net ty=point icon="basic/circle" size=5 col=blue
d.vect schools_net ty=point disp=cat
Now we want to know the reachability of the schools with cats 201, 204, 209. These cats are defined as 'ccats' in 'v.net.iso' (if you want to calculate for all schools, just use ccats=1-10000 or similar):
#next is very slow due to an unsolved bug in graphlib (therefore cache disabled):
#We check for three schools as we give a ccats range:
v.net.iso input=schools_net output=schools_iso ccats=200-210 costs=1000,2500,5000,10000
This calculates the reachability map based on the vector length for 1km, 2.5km, 5km and 10km. Another option were to use vector line attributes from a DBMS table which may contain velocity, street types etc.
Finally we look at the result:
d.erase

#show all streets and schools:
d.vect schools_net

#show reachabilities as calculated above:
d.vect schools_iso col=green  cats=1
d.vect schools_iso col=yellow cats=2
d.vect schools_iso col=orange cats=3
d.vect schools_iso col=red    cats=4

#highlight the selected schools:
d.vect schools disp=attr attr=poiName bcolor=black cat=200-210
d.vect schools_net ty=point cat=200-210 col=red icon="basic/circle" size=5

#nice frame etc:
d.grid -g s=1000
d.barscale -t bc=white tc=black at=60,5

City of Osnabrück - reachability of selected schools (click to enlarge):
City of Osnabrück - reachability of selected schools

Creating subnets: v.net.alloc
Allocates (creates) subnets for nearest centers (direction from center), e.g. zones of responsibility for distributed police stations.

TODO better example (below is close to nonsense):

v.net.alloc input=schools_net output=schools_alloc ccats=200-210

#check for present categories (school IDs)
v.category schools_alloc option=print | sort -nu
0
201
204
209

d.erase

#show all streets and schools:
d.vect schools_net

#show reachabilities as calculated above:
d.vect schools_alloc col=green  cats=201
d.vect schools_alloc col=yellow cats=204
d.vect schools_alloc col=orange cats=209

#highlight the selected schools:
d.vect schools disp=attr attr=poiName bcolor=black cat=200-210
d.vect schools_net ty=point cat=200-210 col=red icon="basic/circle" size=5
 
#nice frame etc:
d.grid -g s=1000lloc col=orange cats=209
d.barscale -t bc=white tc=black at=60,5
City of Osnabrück - patches of selected schools (click to enlarge):
City of Osnabrück - patched of selected schools

v.net.steiner: optimized connection of nodes on a given vector network
Consider a number of hospitals in a city. Medical specialists are usually geographically distributed, and to improve the share of knowledge, applications in telemedicine shall be established. This requires broadband cable connection between the hospitals. As a constraint, to minimize the costs, the cables should follow the existing street network. The GIS task is to find the optimal cable distribution using the smallest amount of wire. This problem is known as "Minimum Steiner Tree" on a graph. In GRASS we can perform such network analysis as described in this section.

We use again the City of Osnabrück vector data as above for this fictitious project. First we extract the hospitals from the points of interest map:

v.extract poi out=hospitals where="poiTypID=2"
d.vect hospitals col=green
d.vect strassen

We can see that the hospitals are not yet connected to the streets network, so we have to generate the connecting lines (using the nearest distance search, of course it were more realistic to digitize with 'v.digit' from a street map):

v.distance -p from=hospitals to=strassen output=hospitals_conn_streets upload=dist column=dist
d.vect hospitals col=yellow
d.zoom
d.vect hospitals_conn_streets col=red
d.vect strassen 

Patch together the hospital points map, the streets map and the map of connecting lines:

v.patch in=hospitals,hospitals_conn_streets,strassen out=h_network
#note that patch doesn't insert nodes where lines are connecting.
#so we clean it (threshold of 1 millimeter):
v.clean h_network out=h_network_clean tool=snap,break thresh=0.001

#display patched map, with hospital points highlighted:
d.vect h_network_clean
d.vect h_network_clean ty=point col=red

Print the categories for the hospital points:

v.category option=print input=h_network_clean type=point
#the range is from 40-215

Calculate the Steiner tree for broadband cable connection network (slow at time to DGLib cache bug):

v.net.steiner in=h_network_clean output=broadband_map afield=1 nfield=1 nsp=-1 tcats=40-300

Let's study the result:

d.erase
d.vect strassen
d.vect hospitals col=yellow
d.vect broadband_map col=red

City of Osnabrück - optimal broadband connection of hospitals for telemedicine applications (click to enlarge):
City of Osnabrück - optimal broadband connection of hospitals

v.net.salesman
Creates a network cycle by connecting given nodes (Traveling salesman problem).

Think that you should introduce medical equipment to all hospitals in the town. You want to minimize the travel route. This is the traveling salesman problem, which can be addressed in GRASS.

TODO (in fact identical preparation up to v.net.steiner from above)

Now we calculate the shortest route to visit all hospitals:

v.net.salesman in=h_network_clean output=salesman nfield=1 ccats=40-300

City of Osnabrück - traveling salesman (click to enlarge):
City of Osnabrück - traveling salesman


Further Links (related software, SQL reference etc).


© 2003 Markus Neteler
Comments about this page | FAQ | Download | Support | Docs | Programming | Back 5.7 Tutorial Home
Last change: $Date: 2008-03-27 21:31:14 +0000 (Thu, 27 Mar 2008) $