HOME

GRASS 5.7 Usage Examples - basics


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.

Convert from/to GRASS 5.0 vector maps
Convert 5.0 vector maps to 5.7

You can run 5.7 (at the same time) on a 5.0 location without problems.
To convert 5.0 vector data to 5.7, run:

v.convert in=vectormap_from_50 out=vectormap
v.db.connect -p vectormap
g.region -p vect=vectormap
d.vect vectormap
d.m
Vector maps from 5.0 and 5.7 do not interfere (different directories), so you can use the same names.

Convert 5.7 vector map back to GRASS 5.0

v.out.ascii -o newmap.clean out=newmap.clean
Then go to 5.0 and run "v.in.ascii -s".

Geometry storage in various formats
There are various possibilities to store vector geometry in GRASS 5.7:

Read on here about storing vector geometry in GRASS 5.7.

Attribute storage in (external) databases
There are various possibilities to store attributes in GRASS 5.7:

Read on here about storing attributes in GRASS 5.7.

Example: Import 3D DXF/DWG vector data and NVIZ 3D vector visualization
WARNING: v.in.dwg requires OpenDWG toolkit, which is proprietary software! To get this toolkit you must become at least "Associate Member" of OpenDWG Alliance (http://www.opendwg.org), it is free-of-charge (you have to fill one form). If you distribute GRASS binaries with linked OpenDWG, you are violating the license of GRASS. Please help us to create a better DWG library or to convince the OpenDWG Allicance to change the license. [for this exercise we use the Spearfish data set from the GRASS web site]
[The WATRTOWR.DXF can be found here]

Since GRASS 5.7 supports 3D vector data, we can easily import DXF and DWG files and display them with NVIZ:

v.in.dwg in=WATRTOWR.DXF out=watertowerXY
The topology should be build for 'faces' which are (filled) 3D vector polygons. The individual vector colors and DXF/DWG types are imported and stored in the related database table (DBF file per default). To know a bit more about the map, we run:
v.info watertowerXY
To look at the imported XY vector map (the DXF/DWG file is usually unprojected) with NVIZ, we have to generate a raster DEM for this XY region. Later we will transform the vector object to UTM coordinates in the Spearfish region.
As a values for the elevation it is a good idea to select the minimum (bottom) height of the imported 3D vector object(s). This value can be reported with 'v.info' (see above, B: bottom value = zmin):
g.region vect=watertowerXY -p
r.mapcalc "dem_flat=-1874"
nviz el=dem_flat vect=watertowerXY

Transformation to UTM coordinates:

Most DXF/DWG drawings are done within XY coordinates. To transform them to a national grid, we can use 'v.transform' with a 4 point transformation. We select UTM coordinates where to place the water tower and query the elevation (d.what.rast elevation.dem). We decide that the xy-extension of the water tower be 20m and it's height 100m (defined later).
The required transformation points file for 'v.transform' may look like this (L: left, R: right, U: upper, L: lower, N, S, W, E):

# XY UTM
#-------
# UL NW
# UR NE
# LR SW
# LL SE
-584 585  598000 4920770
580  585  598020 4920770
580  -600 598020 4920750
-584 -600 598000 4920750
These values are stored in the ASCII file 'wt.points'. Then the transformation can be performed directly within the Spearfish location. The 'zcale' parameter is useful to rescale the vector object height to map units (e.g. 100m tower height) from XY to real-world length. You can use the output of 'v.info' to calculate a value for 'zcale':
zcale=target_height/(T-B)
The 'zshift' parameter is used to shift the object vertically on top of the elevation model (use the bottom 'B' value shown by 'v.info' and add the local elevation). The '-t' flag automatically shifts the DXF/DWG object(s) to sea level which is useful if the vector objects were composited with partially negative z values:
v.transform -t in=watertowerXY out=watertowerUTM points=wt.points zscale=0.04 zshift=1320
g.region rast=elevation.dem
nviz el=elevation.dem vect=watertowerUTM
You will have to zoom heavily to find the water tower as it appears to be very small compared to the surrounding mountains.

Example: Extract vector data from map to new map with SQL statements

#note: strings must be quoted:
v.extract markveggy.shp output=markveggy.1 where="VEGTYPE = 'PS'"

#selection with OR:
# Take care to use parenthesis on both sides of or ( otherwise the result m
# be bogus). Single clause of one comparison should also be enclosed in
# brackets standing on one of the sides of 'or' like in:
#     select * from tab where (c1 < 5) or (c2 >1)

v.extract markveggy.shp out=markveggy.1 new=1 where="(VEGTYPE = 'Wi') or (VEGTYPE = 'PS') or (PRIME_TYPE='Wi')"
#to use the new map, you should link attributes to it. This can
#be the table from the parent map. Edit the DB file and add a row for
#the new map (here: markveggy.1)
#Then check connection:
v.db.connect -p markveggy.1

#Query new map:
d.what.vect

#Or you can display with labels attached to the vector: First we check
# from which column labels to display:
db.columns tab=markveggy
d.vect markveggy.1 att=vegtype display=attr lcolor=red

Example: Vector map import data using OGR
Import of vector maps supported by OGR (such as SHAPE, UK NTF, SDTS, TIGER, .IHO S-57 (ENC), MapInfo File, DGN, GML, AVCBin, REC, ...) can be done with 'v.in.ogr'. It is important to note that a second cleaning step is required to clean such non-topological formats (SHAPE, MIF etc) to topology for vector areas:

v.in.ogr d=wdb2_unep.mif out=wdb2 lay=wdb2_unep

If the map contains areas, we must clean off collinear boundary vectors as well as multiple centroids with 'v.clean'. This is needed, if you see for example the message
WARNING: 2. centroid found in area 2

v.clean in=wdb2 out=wdb2_clean tool=rmdac 
Duplicate boundaries are corrected with:
v.clean in=wdb2_clean out=wdb2_clean2 tool=rmdupl
d.vect wdb2_clean2

Example: Vector map export data using OGR
#Export of GRASS 5.7 vector map to SHAPE format (generates /tmp/testogr.shp #and related files):

v.out.ogr input=multi typ=line dsn=/tmp layer=testogr

#Export to GML format (generates /tmp/testogr..gml file with layer 'testogr'):

v.out.ogr input=multi typ=line dsn=/tmp/testogr.gml layer=testogr format=GML

Example: SQL queries
Example: we want to read attributes from a connected DBF table (see above how to connect attribute tables):

#Select from table 'markveggy' (in database 'grass57test') all rows where 'VEGTYPE = "IFA"':
#NOTE: string attributes must be quoted:
echo "select VEGCNP_ID from markveggy where VEGTYPE = 'IFA'" | db.select -hc

#NOTE: for areas such selections only work when area centroids are present!
d.erase
d.vect markveggy.vegtype where="VEGTYPE = 'IFA'"
d.what.vect -a
Example: we want to read attributes from a PostgreSQL table connected through ODBC (see above how to connect attribute tables):

echo "select VEGCNP_ID from markveggy where VEGTYPE = 'IFA'" | db.select -hc

#NOTE: for areas such selections only work when area centroids are present!
d.erase
d.vect markveggy.vegtype where="VEGTYPE = 'Wiii'"
d.what.vect -a
You see that it is the same! After connecting a GRASS 5.7 vector map, all db.* modules work in the same way independent from the connected RDBMS.

Hints for GRASS-PostGIS: see e.g. this thread.

Example: Generating a vector map from a point data set
Assume you have a CSV (comma separated values) file and want to generate a map from that. That's rather easy. Convert the CSV file to dBase format with StarOffice or another software. Then store the DBF file into that directory, where you keep your DBF files for the current LOCATION (maybe where your SHAPE files live). The map be "BOTSD90.dbf". Connect the DBMI to that directory:

db.connect driver=shp database=/ssi0/ssi/neteler/grassdata/botswanaLL/shp
db.columns BOTSD90
should display the columns. there must be of course an X and a Y column with coordinates inside. Then try to query the table to check if it's really working:
echo "select X_COORD, Y_COORD, shp_fid from BOTSD90" | db.select -h

Note that "shp_fid" column is generated on the fly. Now we can generate the vector points map (which was the sites file in GRASS 5.0):

echo "select X_COORD, Y_COORD, shp_fid from BOTSD90" | db.select -h |v.in.ascii output=BOTSD90map
There we are!
d.vect BOTSD90map

Finally we have to connect the table to the new vector map. Edit the DB file (if not there, generate it in the MAPSET directory) and add something like:

BOTSD90map 1 BOTSD90 shp_fid /ssi0/ssi/neteler/grassdata/botswanaLL/shp shp
(map field_1 table IDcolumn directory driver)

Query some points:

d.what.vect -a BOTSD90map
For amusement we export the map as SHAPE map:
v.out.ogr BOTSD90map  type=point dsn=./ layer=BOTSD90map
Cheers.

Example: Copying a database table from one DBMS to another within GRASS 5.7
Here we copy a DBF file from current directory into a PostgreSQL database table on machine 'pgserver':

db.copy from_driver=dbf from_database=./ from_table=census1990wet \
 to_driver=pg to_database="host=pgserver,dbname=grass57test" \
 to_table=census1990wetB
Simultaneously we can watch the data arriving in PostgreSQL (launch select command several times):
psql -h pgserver
grass57test=> select count(*) from census1990wetB;
This may take some time depending on network speed and table size.

Distances from one climatic station to others (point to point distances)
Assume you have a map of climatic stations and you want to know all distances from a certain station to all others. Easy to solve:

#show all meteo stations:
d.vect meteostations displ=shape,attr att=name
#highlight station of interest:
d.vect meteostations displ=shape,attr att=name where="name='Dos Gaggio'" col=yellow

#extract station of interest (vector point) into a new map:
v.extract meteostations out=stationdg  where="name='Dos Gaggio'"
d.vect stationdg col=blue lcol=blue

Now we can calculate the distances from all stations to the selected station:

#(take care for the 'max' parameter for maximum distance):
v.distance -p from=meteostations to=stationdg col=dist up=dist

#print with distances sorted:
v.distance -p from=meteostations to=stationdg col=dist up=dist | sort -n -k 2 -t '|'
[...]
from_cat|dist
82|1580.509095
518|5560.323732
54|7458.190196
81|9026.606893
83|9093.747797
647|9820.778075
57|10017.483616
[...]
With 'd.what.vect' or another 'db.select' query could can look up the names for the station IDs.

Clipping vector maps to a particular region
Clipping vector maps to a particular region is easily done:

#set region of interest with: d.zoom/g.region

#create area map with single area of current region:
v.in.region out=clipregion

#cut map according to that area map:
v.overlay ainput=clipregion binput=soils out=soils_clipped operator=and

#look at the clipped vector map:
d.vect soils_clipped
More vector clipping and overlay examples.


Further Links (related software, SQL reference etc).


© 2002-2004 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) $