[GRASS-SVN] r40323 - grass-addons/vector/v.transect.kia
svn_grass at osgeo.org
svn_grass at osgeo.org
Fri Jan 8 11:26:40 EST 2010
Author: prea
Date: 2010-01-08 11:26:39 -0500 (Fri, 08 Jan 2010)
New Revision: 40323
Added:
grass-addons/vector/v.transect.kia/v.transect.kia
grass-addons/vector/v.transect.kia/v.transect.kia.html
Log:
importing v.transect.kia addon
Added: grass-addons/vector/v.transect.kia/v.transect.kia
===================================================================
--- grass-addons/vector/v.transect.kia/v.transect.kia (rev 0)
+++ grass-addons/vector/v.transect.kia/v.transect.kia 2010-01-08 16:26:39 UTC (rev 40323)
@@ -0,0 +1,861 @@
+#!/bin/bash
+
+# start stopwatch
+start_time="$(date +%s)"
+
+################################################################################
+#
+# MODULE: v.transect.kia
+# VERSION: 0.11rc3
+# AUTHOR(S): Damiano G. Preatoni <prea at uninsubria.it>,
+# Clara Tattoni <clara.tattoni at gmail.com>
+# PURPOSE: prepare datasets for line transect analysis
+# COPYRIGHT: (C) 2007-2010 Damiano G. Preatoni & Clara Tattoni
+#
+# This program is free software; you can redistribute it and/or modify
+# it under the terms of the GNU General Public License as published by
+# the Free Software Foundation; either version 2 of the License, or
+# (at your option) any later version.
+#
+# This program is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+# GNU General Public License for more details.
+#
+################################################################################
+# Version History
+# 0.4 - dgp first release
+# 0.5 - dgp changed call to r.reclass using piping, to be compatible with GRASS 6.2
+# 0.6 - dgp DEM contrours are now clipped just to lines extent (10 u buffer)
+# added code to show progress
+# silenced out all stdout redirecting to /dev/null (1>/dev/null)
+# 0.7 - ct added a speed routine that uses a MASK based on buffered features.
+# ct moved raster processing inside the loop
+# 0.8 - dgp borrowed v.segment.elev v. 0.7 codebase, adapted to do data import
+# and preprocessing
+# 0.9 - dgp rewritten starting from v.transect.segment.sh v. 0.7, renamed to v.transect.KIA
+# rationalised parameter names, rewritten GUI/CLI integration
+# number of mandatory fields in inout shapefile has been reduced to the very bone
+# moved constants declarations at script top, as well as function declarations
+# added full sqlite support to have full SQL calculation capacities
+# rewritten all text output using g.message
+# cleaned up 3D draping routine
+# added waypoints weighting routine
+# enhanced timer functionalities
+# removed -w flag in v.in.ogr calls to preserve names, mandatory fields are searched for
+# and turned to lowercase
+# tidied up v.0.7 ct "speedup hack" in segmantation routine:
+# now saves current region and restores at end
+# added ctrl-c trapping and cleanup routine
+# added output optioni
+# 0.10 - added test for cut program
+# added sqlite detection logici
+# rewritten IKA calculation logic SQL code
+# enhanced elapsed time display at end
+# 0.11 - changed name to all lowercase,
+# polished out according to GRASS-addons coding standards
+# reviewed g.parser metacomments
+# TODO: check and fix temporary maps cleanup at exit (or when user interrupted)
+# TODO: check and fix any use of "echo" instead of g.message
+#
+################################################################################
+# paths shapefile attribute table must contain:
+# TRANS_ID C,255
+#
+# waypoint shapefile attribute table must contain:
+# TRANS_ID C,255
+# TYPE C,4
+# N N,4,0
+#
+################################################################################
+#%Module
+#% label: v.transect.kia
+#% description: Calculate Kilometric Abundance Indexes on line transect surveys of presence signs.
+#% keywords: vector, elevation, geometry, line transect, abundance index, KIA
+#%End
+#%Option
+#% key: paths
+#% description: a line vector map containing line transect paths
+#% type: string
+#% required: yes
+#% multiple: no
+#% key_desc: <paths map>
+#% gisprompt: old_file,file,input
+#%End
+#%Option
+#% key: waypoints
+#% description: a point vector containing objects recorded along transects
+#% type: string
+#% required: yes
+#% multiple: no
+#% key_desc: <waypoints map>
+#% gisprompt: old_file,file,input
+#%End
+#%Option
+#% key: idfield
+#% answer: trans_id
+#% description: field (present both in paths and waypoints) with transect IDs
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: <transect ID field>
+#% gisprompt: old_dbcolumn,dbcolumn,dbcolumn
+#%End
+#%Option
+#% key: typefield
+#% answer: type
+#% description: field in waypoints map with point type information
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: <waypoint type field>
+#% gisprompt: old_dbcolumn,dbcolumn,dbcolumn
+#%End
+#%Option
+#% key: nfield
+#% answer: n
+#% description: field in waypoints map containing number of items (for clustered waypoints)
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: <waypoint cluster size field>
+#% gisprompt: old_dbcolumn,dbcolumn,dbcolumn
+#%End
+#%Option
+#% key: output
+#% description: output line vector map
+#% type: string
+#% required: yes
+#% multiple: no
+#% key_desc: <output map>
+#% gisprompt: new_file,file,output
+#%End
+#%Option
+#% key: elev
+#% description: digital elevation model to correct path lengths for elevation by draping
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: <digital elevation model>
+#% gisprompt: old,cell,raster
+#%End
+#%Option
+#% key: weights
+#% description: text file with waypoint weights by type
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: <waypoint weight file>
+#% gisprompt: old_file,file,input
+#%End
+#%Option
+#% key: groups
+#% description: a vector (polygon) map to segment transect paths with and calculate partial KIA
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: <paths segmentation map>
+#% gisprompt: old,vector,vector
+#%End
+#%Option
+#% key: class
+#% answer: class
+#% description: field in groups map containing class values
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: <class values field>
+#% gisprompt: old_dbcolumn,dbcolumn,dbcolumn
+#%End
+#%flag
+#%key: s
+#%description: assumes both paths and waypoints are in ESRI Shapefile format, does implicit conversion (using v.in.ogr -oe)
+#%end
+#%flag
+#% key: o
+#% description: overwrite output
+#%end
+
+# handle user breaks
+exitprocedure () {
+ echo "Interrupted by user!"
+ cleanup
+ exit 1
+}
+# shell check for user break (signal list: trap -l)
+trap "exitprocedure" 2 3 15
+
+################################################################################
+# constants
+################################################################################
+# lists of mandatory fields in both paths and waypoints attribute tables
+if [ -n "$GIS_OPT_IDFIELD" ] ; then
+ FIELDS_PATHS="$GIS_OPT_IDFIELD"
+else
+ FIELDS_PATHS=trans_id
+fi
+if [ -n "$GIS_OPT_TYPEFIELD" ] ; then
+ FIELDS_TYPES=type
+else
+ FIELDS_TYPES="$GIS_OPT_TYPEFIELD"
+fi
+if [ -n "$GIS_OPT_NFIELD" ] ; then
+ FIELDS_N=n
+else
+ FIELDS_N="$GIS_OPT_NFIELD"
+fi
+FIELDS_WYPTS="`echo ${FIELDS_PATHS[@]}` $FIELDS_TYPES $FIELDS_N"
+
+# settings for segmentation routine
+MSK_RES=2 # speedup mask resolution in map units, keep it very small (2m)
+MSK_BUFSIZE=1000 # speedup mask buffer width in map units
+MSK_BUFUNITS=meters # speedup mask map units specification
+# temporary file/raster/vector names
+TMP_REGION=region.$$
+SAVED_REGION=region.$$.saved
+SAVED_MASK=MASK.saved
+TMP_RAST=TMP_RAST # raster conversion of current vector path
+TMP_BUF=TMP_BUF # raster buffer, MSK_BUFSIZE wide, around TMP_RAST
+TMP_SINGLETON=TMP_SINGLETON
+TMP_STAGING0=TMP_STAGING0
+TMP_STAGING1=TMP_STAGING1
+TMP_SEGMENT=TMP_SEGMENT
+SEGMENT_PART=SEGMENT_PART
+PATHS_REGION=PATHS_REGION
+SEGMENT_REGION=SEGMENT_REGION
+CLIPPER_REGION=CLIPPER_REGION
+SCRATCH_DB=TMP_sqlite.db
+
+# drop unnecessary fields
+# $1: name of vector map whose fields are to be dropped
+# $2: array with names of fields to keep, if present
+drop_fields () {
+ VECTOR=$1
+ local FIELDS_TO_KEEP
+ FIELDS_TO_KEEP="`echo "$2"`"
+ g.message -i message="dropping unneeded fields in $VECTOR... "
+ FIELDS_TO_CHECK=`v.info -c $VECTOR --quiet | awk -F'|' '{print $2}'`
+ for ftc in $FIELDS_TO_CHECK; do
+ g.message -v message="about to drop $ftc from $VECTOR..."
+ local DROP_IT=1
+ for ftk in ${FIELDS_TO_KEEP[@]}; do
+ if [ "$ftc" == "$ftk" ] ; then # ok, will keep field
+ DROP_IT=0
+ fi
+ ftcl=`echo $ftc | awk '{$1=tolower($1);print}'`
+ if [ "$ftcl" == "$ftk" ] ; then # ok, will keep field, too
+ DROP_IT=0
+ fi
+ if [ ${DROP_IT} == 0 ] ; then
+ g.message -v message="keeping $ftc field in $VECTOR"
+ break
+ fi
+ done
+ # if drop flag is still 1, then drop
+ if [ ${DROP_IT} == 1 ] ; then
+ g.message -w message="dropping $ftc field from $VECTOR"
+ v.db.dropcol map=$VECTOR column=$ftc --quiet
+ fi
+ done
+}
+
+# check (case insensitively) whether all needed fields are present
+# if a field is found, its name is automatically converted into lowercase
+# $1: name of vector map whose fields are to be checked
+# $2: array with names of fields that must be present
+check_fields () {
+ VECTOR=$1
+ local FIELDS_TO_KEEP
+ FIELDS_TO_KEEP="`echo "$2"`"
+ g.message -i message="checking necessary fields in $VECTOR... "
+ FIELDS_TO_CHECK=`v.info -c $VECTOR --quiet | awk -F'|' '{print $2}'`
+ ALL_PRESENT=1
+ for ftk in ${FIELDS_TO_KEEP[@]}; do
+ g.message -v message="looking for field $ftk in $VECTOR... "
+ local FOUND=0
+ for ftc in $FIELDS_TO_CHECK; do
+ g.message -v message="checking $ftc against $ftk"
+ if [ "$ftc" == "$ftk" ] ; then
+ FOUND=1
+ fi
+ ftcl=`echo $ftc | awk '{$1=tolower($1);print}'`
+ g.message -v message="(checking $ftcl against $ftk)"
+ if [ "$ftcl" == "$ftk" ] ; then
+ # turn field name to lowercase, to allow final v.db.join to work
+ # sqlite API seems rather case-insensitive, but GRASS is!
+ v.db.renamecol map=$VECTOR column=$ftc,X$ftc --quiet
+ v.db.renamecol map=$VECTOR column=X$ftc,$ftk --quiet
+ g.message -w message="Column name changed: '$ftc' -> '$ftk'"
+ FOUND=1
+ fi
+ g.message -v message="[$FOUND]"
+ if [ ${FOUND} == 1 ] ; then
+ g.message -v message="checking $VECTOR for required field: $ftk found"
+ break
+ fi
+ done
+ if [ ${FOUND} == 0 ] ; then
+ g.message -e message="required field $ftk is not present in input map $VECTOR."
+ ALL_PRESENT=0
+ fi
+ done
+ if [ ${ALL_PRESENT} == 0 ]
+ then
+ g.message -e message="some required field are missing in input map $VECTOR. Aborting."
+ exit 1
+ fi
+}
+
+# extract keys and values from waypoints weight file
+# file format: K1 {K2 {Kn}} = [<integer> | Q | N]
+# last element is always scalar and is the value
+# other elements are keys for that value
+# WEIGHT_KEYS and WEIGHT_VALUE hold the results
+getkeys () {
+ WEIGHT_KEYS=$(echo $1 | cut -d '=' -f 1 )
+}
+
+getvalue () {
+ WEIGHT_VALUE=$(echo $1 | cut -d '=' -f 2 )
+}
+
+# searches for elements in current mapset
+# $1 is 'element' sensu g.findfile, i.e. windows, vector, cell
+# $2 is 'file' sensu g.findfile
+# TODO check if passing $MAPSET as $3 is actually needed
+findelement () {
+ findelement=0
+ eval `g.findfile element=$1 file=$2 mapset=$MAPSET`
+ if [ -n "$file" ] ; then # found
+ findelement=1
+ else
+ findelement=0
+ fi
+}
+
+# housekeeping function
+cleanup() {
+ # restore any pre-existing MASK
+ findelement cell $SAVED_MASK
+ if [ $findelement == 1 ] ; then
+ g.rename rast=$SAVED_MASK,MASK --quiet > /dev/null
+ r.mask -o input=$SAVED_MASK --quiet > /dev/null
+ fi
+ # reset original region
+ findelement windows $SAVED_REGION
+ if [ $findelement == 1 ] ; then
+ g.region region=$SAVED_REGION --quiet > /dev/null
+ g.remove region=$SAVED_REGION --quiet
+ fi
+ # remove stale temporary maps
+ findelement vector $PATHS_IN
+ if [ $findelement == 1 ] ; then
+ g.remove vect=$PATHS_IN --quiet > /dev/null
+ fi
+ findelement vector $WYPTS_IN
+ if [ $findelement == 1 ] ; then
+ g.remove vect=$WYPTS_IN --quiet > /dev/null
+ fi
+ findelement vector $TMP_SEGMENT
+ if [ $findelement == 1 ] ; then
+ g.remove vect=$TMP_SEGMENT --quiet > /dev/null
+ fi
+ findelement vector $SEGMENT_PART
+ if [ $findelement == 1 ] ; then
+ g.remove vect=$SEGMENT_PART --quiet > /dev/null
+ fi
+ findelement vector $TMP_SINGLETON
+ if [ $findelement == 1 ] ; then
+ g.remove vect=$TMP_SINGLETON --quiet > /dev/null
+ fi
+ # restore previous DB connection, if different from sqlite
+ if [ ${DB_PARAMS[0]} != "sqlite" ] ; then
+ db.connect driver="${DB_PARAMS[0]}" database="${DB_PARAMS[1]}"
+ g.message -v message="restored db connection to ${DB_PARAMS[1]}, using driver ${DB_PARAMS[0]}"
+ rm -rf "$GISDBASE/$LOCATION_NAME/$MAPSET/$SCRATCH_DB"
+ fi
+}
+
+
+################################################################################
+# program starts here!
+if [ -z "$GISBASE" ] ; then
+ # for obvious reasons, here g.message is not used...
+ echo "You must be in GRASS GIS to run this program." 1>&2
+ exit 1
+fi
+
+# check if we have awk
+if [ ! -x "`which awk`" ] ; then
+ g.message -e message="awk required, please install awk or gawk first"
+ exit 1
+fi
+
+# check if we have cut
+if [ ! -x "`which cut`" ] ; then
+ g.message -e message="cut required, please install cut first"
+ exit 1
+fi
+
+if [ "$1" != "@ARGS_PARSED@" ] ; then
+ exec g.parser "$0" "$@"
+fi
+
+# check immediately if output exists
+if [ -f "$GIS_OPT_OUTPUT" ] ; then
+ if [ "$GIS_FLAG_O" == 0 ] ; then
+ g.message -e message="output file exists, and -o option not specified. Exiting."
+ exit 1
+ fi
+fi
+
+# be verbose, if asked to
+if [ "$GIS_FLAG_S" == 1 ] ; then
+ g.message -v message="waypoints, paths and output will be in ESRI shapefile format"
+fi
+g.message -v message="input paths are in $GIS_OPT_PATHS"
+g.message -v message="input waypoints are in $GIS_OPT_WAYPOINTS"
+if [ -n "$GIS_OPT_ELEV" ] ; then
+ g.message -v message="digital elevation model is $GIS_OPT_ELEV"
+else
+ g.message -v message="no digital elevation model has been defined, paths won't be 3D corrected"
+fi
+if [ -n "$GIS_OPT_WEIGHTS" ] ; then
+ g.message -v message="waypoint weight file is $GIS_OPT_WEIGHTS"
+else
+ g.message -v message="no waypoint weight file has been defined, all waypoints will count 1 for KIA calculations"
+fi
+if [ -n "$GIS_OPT_GROUPS" ] ; then
+ g.message -v message="grouping map is $GIS_OPT_GROUPS, partial KIA will be calculated on path segments"
+else
+ g.message -v message="no grouping map has been defined, will calculate a global KIA for each path"
+fi
+
+# get GISDBASE, LOCATION, MAPSET
+GISDBASE=`g.gisenv get=GISDBASE`
+LOCATION_NAME=`g.gisenv get=LOCATION_NAME`
+MAPSET=`g.gisenv get=MAPSET`
+
+# store current db connection
+OIFS="$IFS"
+IFS=$'\n'
+DB_PARAMS=(`db.connect -p | awk -F: '{print $2}'`) # bashism!
+IFS="$OIFS"
+
+g.message -v message="testing for an existing sqlite connection"
+# check if already running with sqlite. if not, do a db.connect
+# $DB_PARAMS[0] -> driver, $DB_PARAMS[1] -> database
+if [ "${DB_PARAMS[0]}" != "sqlite" ] ; then
+ g.message -v message="switching to a sqlite database connectiom"
+ # set up a sqlite connection to work with: all tables need be in sqlite
+ db.connect driver=sqlite database="$GISDBASE/$LOCATION_NAME/$MAPSET/$SCRATCH_DB"
+else
+ g.message -v message="found a sqlite db in use at ${DB_PARAMS[1]}"
+fi
+
+# save current MASK if any
+findelement cell MASK
+if [ "$findelement" == 1 ] ; then
+ g.copy rast=MASK,$SAVED_MASK --quiet > /dev/null
+ r.mask -r input=$SAVED_MASK --quiet > /dev/null
+fi
+
+# save current REGION if any
+g.region save=$SAVED_REGION
+
+# initialise some internal variables
+if [ "$GIS_FLAG_S" == 1 ] ; then
+ PATHS=`basename $GIS_OPT_PATHS .shp`
+ WYPTS=`basename $GIS_OPT_WAYPOINTS .shp`
+else
+ PATHS=$GIS_OPT_PATHS
+ WYPTS=$GIS_OPT_WAYPOINTS
+fi
+PATHS_IN=$PATHS"_in"
+PATHS_3D=$PATHS"_3d"
+WYPTS_IN=$WYPTS"_in"
+WYPTS_3D=$WYPTS"_3d"
+if [ -n "$GIS_OPT_ELEV" ] ; then
+ DEM=$GIS_OPT_ELEV
+fi
+if [ -f "$GIS_OPT_WEIGHTS" ] ; then # if a weights file is defined, preprocess it
+ WEIGHTSFILE=$GIS_OPT_WEIGHTS
+ # now, glob weights file
+ declare -a WEIGHTS_LINES
+ OIFS="$IFS" # save old field separator
+ IFS=$'\n' # set field separator to newline
+ set -f # no globbing, see help set
+ WEIGHTS_LINES=( $(< $WEIGHTSFILE ) ) # read file
+ set +f # reset globbing
+ IFS="$OIFS" # reset field separator
+fi
+
+if [ -n "$GIS_OPT_GROUPS" ] ; then
+ SEGMENT=$GIS_OPT_GROUPS
+fi
+
+# convert shapefiles into GRASS vector, if needed
+if [ "$GIS_FLAG_S" == 1 ] ; then
+ g.message -i message="converting input paths in `dirname $GIS_OPT_PATHS`/$PATHS into $PATHS_IN"
+ v.in.ogr -oe dsn=`dirname $GIS_OPT_PATHS` layer=$PATHS output=$PATHS_IN --overwrite --quiet
+ g.message -i message="converting input waypoints in `dirname $GIS_OPT_WAYPOINTS`/$WYPTS into $WYPTS_IN"
+ v.in.ogr -oe dsn=`dirname $GIS_OPT_WAYPOINTS` layer=$WYPTS output=$WYPTS_IN --overwrite --quiet
+else
+ g.copy vect=$GIS_OPT_PATHS,$PATHS_IN --overwrite --quiet
+ g.copy vect=$GIS_OPT_WAYPOINTS,$WYPTS_IN --overwrite --quiet
+fi
+
+# check if all needed fields are present, if not, complain and exit
+# echo FIELDS_PATHS in another variable to properly pass an array as an argument
+CHECK_PATHS=`echo ${FIELDS_PATHS[@]}`
+check_fields $PATHS_IN "$CHECK_PATHS"
+CHECK_WYPTS=`echo ${FIELDS_WYPTS[@]}`
+check_fields $WYPTS_IN "$CHECK_WYPTS"
+# throw away superfluous fields
+KEEP_PATHS=`echo cat ${FIELDS_PATHS[@]}`
+drop_fields $PATHS_IN "$KEEP_PATHS"
+KEEP_WYPTS=`echo cat ${FIELDS_WYPTS[@]}`
+drop_fields $WYPTS_IN "$KEEP_WYPTS"
+
+# no need to save preexisting REGION or MASK, this is done at program start once and for all
+
+# segment path by $GROUPS map
+if [ -n "$SEGMENT" ] ; then
+ if [ -n "$GIS_OPT_CLASS" ] ; then
+ CLASSFIELD=$GIS_OPT_CLASS
+ else
+ CLASSFIELD=class
+ fi
+ g.message -v message="splitting paths in $PATHS_IN using $SEGMENT map, classified by field $CLASSFIELD"
+ # change region to an extent as wide as $PATHS_IN
+ g.region vect=$PATHS_IN
+ # turn it into a $SEGMENT_REGION vector
+ v.in.region output=$PATHS_REGION --overwrite --quiet
+ # select features in $SEGMENT that overlap $SEGMENT_REGION
+ v.select ainput=$SEGMENT binput=$PATHS_REGION output=$SEGMENT_REGION --overwrite
+ # now, use $SEGMENT_REGION instead of bigger and slower $SEGMENT
+ NUM=0 # processed primitives counter
+ FIRST=true # true if first primitive is being processed
+ # list features in $GIS_OPT_INPUT
+ CATS=`v.category input=$PATHS_IN option=print | grep -v "/"`
+ # count features in $GIS_OPT_INPUT (added in 0.6)
+ NUMCATS=`v.category input=$PATHS_IN option=report --quiet | grep all | awk -F" " '{print $2}'`
+ # save current region
+ # g.region save=$TMP_REGION --overwrite
+ # process features
+ for C in $CATS; do
+ NUM=$((NUM+1)) # changed in 0.6
+ g.message -v message="Processing primitive $C ($NUM of $NUMCATS) ..."
+ # extract a single feature from $GIS_OPT_INPUT
+ v.extract input=$PATHS_IN output=$TMP_SINGLETON type=line new=-1 list=$C --overwrite --quiet > /dev/null
+ # added in 0.7 ct
+ # To have a dramatic speed gain, do calculations only in a small neighborhood of
+ # the current feature (as small as possible!)). Use a raster buffer of the current
+ # feature as a mask to set the calculation region.
+ # Also take care when doing v.overlay, since is the slowest part, not to intersect
+ # all the stuff in it with a minuscule singleton path...
+ # This speeds up a 15 h calculation with v. 0.6 cycle to 20 min
+ # There is a difference in lengths as calculated with 0.6 slow algorithm, but is < 1%.
+ # (kudos to ct!)
+
+ # see http://article.gmane.org/gmane.comp.gis.grass.devel/27329
+
+ # Remove existing mask if any
+ g.message -v message="masking ... "
+ findelement cell MASK
+ if [ $findelement == 1 ] ; then
+ r.mask -r input=MASK --quiet > /dev/null
+ fi
+
+ # old code to remove MASK, superseded by findelement function
+# CHKMASK=`g.findfile element=cell file=MASK | grep '^name' | awk -F "=" '{print $2}'`
+# if [ "$CHKMASK" ] ; then
+# r.mask -r $CHKMASK --quiet
+# g.remove MASK --quiet
+# fi
+
+ # Convert input vect to rast, raster buffer calculation is _ways_ faster
+ g.message -v message="setting region to $TMP_SINGLETON at $MSK_RES resolution..."
+ g.region vect=$TMP_SINGLETON res=$MSK_RES
+ g.message -v message="buffering..."
+ v.to.rast input=$TMP_SINGLETON output=$TMP_RAST use=val layer=1 value=1 --overwrite --quiet # rows=4096 is the default
+ r.buffer input=$TMP_RAST output=$TMP_BUF distances=$MSK_BUFSIZE units=$MSK_BUFUNITS --overwrite --quiet
+ g.remove -f rast=$TMP_RAST --quiet > /dev/null
+ # Set calculated buffer as mask
+ r.mask input=$TMP_BUF maskcats=* -o --quiet > /dev/null
+ g.message -v message="setting region to $TMP_BUF..."
+ g.region rast=$TMP_BUF
+ v.in.region output=$CLIPPER_REGION
+ # limit $SEGMENT to the geometries iverlapped by $TMP_SINGLETON
+ g.message -v message="doing v.select... "
+ v.select ainput=$SEGMENT_REGION binput=$TMP_SINGLETON output=$SEGMENT_PART --overwrite
+ # do v.overlay with that feature
+ g.message -v message="intersecting... "
+ v.overlay ainput=$TMP_SINGLETON atype=line binput=$SEGMENT_PART btype=area output=$TMP_SEGMENT operator=and olayer=1,0,0 --overwrite --verbose #--quiet > /dev/null
+ #g.remove vect=$TMP_SINGLETON --quiet 1>/dev/null
+ g.region vect=$PATHS_IN
+ g.message -v message="fixing fields..."
+ #v.db.dropcol map=$TMP_SEGMENT column=id --quiet
+ v.db.dropcol map=$TMP_SEGMENT column=a_cat --quiet
+ v.db.dropcol map=$TMP_SEGMENT column=b_cat --quiet
+ v.db.addcol map=$TMP_SEGMENT columns="trans_id varchar(100), $CLASSFIELD varchar(100)"
+ echo "update $TMP_SEGMENT set trans_id=a_trans_id, $CLASSFIELD=b_$CLASSFIELD" | db.execute
+ v.db.dropcol map=$TMP_SEGMENT column=a_trans_id --quiet
+ v.db.dropcol map=$TMP_SEGMENT column=b_$CLASSFIELD --quiet
+ # create staging vector if first loop
+ g.message -v message="adding ... "
+ if [ "$FIRST" = "true" ] ; then # if first feature processed, convert it to the staging vector to have a coherent attribute table
+ g.rename vect=$TMP_SEGMENT,$TMP_STAGING0 --overwrite --quiet 1>/dev/null
+ FIRST=false
+ else # add the result to a staging vector (using v.patch)
+ v.patch -e input=$TMP_SEGMENT,$TMP_STAGING0 output=$TMP_STAGING1 --overwrite --quiet 1>/dev/null
+ #g.remove vect=$TMP_SEGMENT --quiet 1>/dev/null
+ #g.remove vect=$TMP_STAGING0 --quiet 1>/dev/null
+ g.rename vect=$TMP_STAGING1,$TMP_STAGING0 --overwrite --quiet 1>/dev/null
+ fi
+ # NUM=$((NUM+1)) # changed in 0.6
+ #echo "done."
+ done
+ g.remove MASK --quiet
+ g.remove rast=$TMP_BUF --quiet
+ g.region vect=$PATHS_IN
+ #substitute segmented paths
+ g.remove vect=$PATHS_IN --quiet
+ g.rename vect=$TMP_STAGING0,$PATHS_IN --quiet
+ # clean up unneded fields in the new $PATHS_IN: must keep just PATHS base fields, plus $CLASSFIELD
+ KEEP_FIELDS=`echo cat $CLASSFIELD ${FIELDS_PATHS[@]}`
+ drop_fields $PATHS_IN "$KEEP_FIELDS"
+ g.region region=$TMP_REGION
+ # 'segment' waypoints
+ TYPE=`v.info -c $SEGMENT |grep $CLASSFIELD | awk -F'|' '{print $2,$1}'`
+ v.db.addcol map=$WYPTS_IN column="$TYPE" --quiet
+ v.what.vect vector=$WYPTS_IN column=$CLASSFIELD qvector=$SEGMENT qcolumn=$CLASSFIELD --quiet
+ g.remove vect=$TMP_SEGMENT --quiet
+ g.remove vect=$TMP_SINGLETON --quiet
+fi # end of segmentation routine
+
+# do 3D correction, if a DEM is available
+if [ -n "$DEM" ] ; then
+ g.message -i message="doing 3D draping..."
+ # set region to $DEM, see @FIXME note below
+ g.region rast=$DEM zoom=$DEM
+ # process paths
+ g.message -v message="draping $PATHS_IN on $DEM..."
+ ##@FIXME sometines setting region here causes v.drape to misbehave. setting region to the whole $DEM works instead.
+ ##g.region vect=$PATHS_IN align=$DEM
+ v.drape input=$PATHS_IN type=line rast=$DEM method=cubic output=$PATHS_3D --overwrite --quiet
+ # add true length to paths attribute table
+ v.db.addcol map=$PATHS_3D columns="length double" --quiet 1>/dev/null
+ v.to.db map=$PATHS_3D option=length column=length --quiet 1>/dev/null
+ # process waypoints
+ g.message -v message="draping $WYPTS_IN on $DEM..."
+ ##g.region vect=$WYPTS_IN align=$DEM
+ v.drape input=$WYPTS_IN type=point rast=$DEM method=cubic output=$WYPTS_3D --overwrite --quiet
+ # clean up points_in and paths_in
+ g.remove -f vect=$PATHS_IN,$WYPTS_IN --quiet
+ g.rename vect=$PATHS_3D,$PATHS_IN --quiet
+ g.rename vect=$WYPTS_3D,$WYPTS_IN --quiet
+ LENGTH_FIELD=length
+else # no 3D correction, anyway we need path lengths
+ v.db.addcol map=$PATHS_IN columns="length double" --quiet 1>/dev/null
+ v.to.db map=$PATHS_IN option=length column=length --quiet 1>/dev/null
+ LENGTH_FIELD=length
+fi
+g.message -v message="path length calculated in field $LENGTH_FIELD"
+
+# waypoints weighting
+# add to all waypoints a weight equal to 1
+g.message -v message="weighing points..."
+v.db.addcol map=$WYPTS_IN columns="weight integer" --quiet
+v.db.update map=$WYPTS_IN column=weight value=1 --quiet
+
+# if a weight table has been specified, give waypoints their weight
+if [ -n "$WEIGHTS_LINES" ] ; then
+ g.message -v message="calculating weights..."
+ # iterate thru WEIGHT_LINES
+ for ((i=0; i < "${#WEIGHTS_LINES[@]}"; i++)); do
+ #echo "line $i: ${WEIGHTS_LINES[${i}]}"
+ # TODO add something to skip lines beginning with #'s
+ getkeys "${WEIGHTS_LINES[${i}]}"
+ getvalue "${WEIGHTS_LINES[${i}]}"
+ for K in ${WEIGHT_KEYS[@]}; do
+ WEIGHT_VALUE=${WEIGHT_VALUE/ /} # drop any space, so the case can go smooth
+ case "$WEIGHT_VALUE" in
+ "N") # use N field values
+ g.message -v message="will assign value from N field to points of type $K"
+ v.db.update map=$WYPTS_IN column=weight qcolumn=n where="type='$K'";
+ ;;
+ "Q") # use quartiles
+ # calculate statistics on field N, for that class of waypoint only
+ stats=($( v.univar -ge map=$WYPTS_IN type=point column=n where="type='$K'" | cut -d '=' -f 2 - ))
+ # last 4 elements out uf 17 are Q1, median, Q3 and 90th percentile
+ Q1=${stats[14]}
+ Q2=${stats[15]}
+ Q3=${stats[16]}
+ g.message -i message="will assign value from quartiles of N field to points of type $K:"
+ g.message -i message=" for 0 <= N < $Q1 -> weight = 1"
+ g.message -i message=" for $Q1 <= N < $Q2 -> weight = 2"
+ g.message -i message=" for $Q2 <= N < $Q3 -> weight = 3"
+ g.message -i message=" for N > $Q3 -> weight = 4"
+ v.db.update map=$WYPTS_IN column=weight value=1 where="type='$K' and n < $Q1"
+ v.db.update map=$WYPTS_IN column=weight value=2 where="type='$K' and (n >= $Q1 AND n < $Q2)"
+ v.db.update map=$WYPTS_IN column=weight value=3 where="type='$K' and (n >= $Q2 AND n < $Q3)"
+ v.db.update map=$WYPTS_IN column=weight value=4 where="type='$K' and n >= $Q3"
+ ;;
+ *)
+ g.message -i message="will assign value $WEIGHT_VALUE to points of type $K"
+ v.db.update map=$WYPTS_IN column=weight value=$((WEIGHT_VALUE)) where="type='$K'"
+ ;;
+ esac
+ done
+ done
+fi
+
+# now, both tables are set up properly, do a query to get KIA values, and store them
+g.message -v message="dropping stale KIA table, if any..."
+echo "drop table if exists kia;" | db.execute --quiet
+# although sqlite offers the capability to create a table on the fly,
+# this causes KIA fields to be of text type, so create the table first
+SQL="create table kia (id text"
+if [ -n "$SEGMENT" ] ; then
+ SQL="$SQL, $CLASSFIELD text"
+fi
+SQL="$SQL, KIA real"
+if [ -n "$WEIGHTS_LINES" ] ; then
+ SQL="$SQL, WKIA real"
+fi
+SQL="$SQL )"
+echo "creating table: $SQL"
+echo $SQL | db.execute --quiet
+
+## assembling SQL statement
+g.message -i message="calculating KIA..."
+
+## no-segmentation SQL (KIA & WKIA)
+if [ -z "$SEGMENT" ] ; then
+ SQL="insert into kia(id, KIA"
+ if [ -n "$WEIGHTS_LINES" ] ; then
+ SQL="$SQL, WKIA"
+ fi
+ SQL="$SQL) select w.[trans_id] as id, "
+ SQL="$SQL cast(count(w.cat)/p.length*1000 as numeric(10,6)) as KIA "
+ if [ -n "$WEIGHTS_LINES" ] ; then
+ SQL="$SQL, cast(sum(w.weight)/p.length*1000 as double numeric(10,6)) as WKIA"
+ fi
+ SQL="$SQL from $WYPTS_IN w "
+ SQL="$SQL left join $PATHS_IN p on p.[trans_id] = w.[trans_id] group by w.[trans_id];"
+ # once the SQL has been bulit, execute it
+ g.message -v message="calculating KIA, SQL: $SQL"
+ echo "$SQL" | db.execute --quiet
+fi
+
+## segmentation SQL
+## example query
+##"insert into kia(id, value, KIA, WKIA) select w.[trans_id] as id, p.value, cast(count(w.cat)/p.length*1000 as numeric(10,6)) as KIA , cast(sum(w.weight)/p.length*1000 as double numeric(10,6)) as WKIA from waypoints_in w left join paths_in p on p.[trans_id] = w.[trans_id] and p.value = w.value group by w.[trans_id],w.value"
+##--insert into kia(id, value, KIA, WKIA) select w.[trans_id] as id, p.value, cast(count(w.cat)/p.length*1000 as numeric(10,6)) as KIA , cast(sum(w.weight)/p.length*1000 as double numeric(10,6)) as WKIA from waypoints_in w left join paths_in p on p.[trans_id] = w.[trans_id] and p.value = w.value group by w.[trans_id],w.value
+##-- FROM wypt
+##--select w.[trans_id], w.value, count(w.cat) as CAT from waypoints_in w group by w.[trans_id],w.value
+##-- FROM paths
+##--select p.[trans_id], p.value, sum(p.length) as length from paths_in p group by p.[trans_id],p.value
+##-- FROM both
+## select p.[trans_id], p.value, sum(p.length) as length count(w.cat) as n from paths_in left join waypoints_in w on p.[trans_id] = w.[trans_id] and p.value = w.value p group by p.[trans_id],p.value
+## BUG: selecting in one shot causes paths being repeated.
+## FIX: three-step process: i) select paths; ii) select wypts; iii) join & calculate
+
+
+
+if [ -n "$SEGMENT" ] ; then
+ ## create and populate a paths summary intermediate table (use sqlite temporary views)
+ echo "drop view if exists paths_sum" | db.execute --quiet
+ SQL="create view paths_sum as select p.[trans_id], p.$CLASSFIELD, sum(p.length) as length from $PATHS_IN p group by p.[trans_id],p.$CLASSFIELD"
+ g.message -v message="creating temporary paths table: $SQL"
+ echo $SQL | db.execute --quiet
+ ## create and populate a waypoints summary intermediate table (use sqlite temporary views)
+ echo "drop view if exists wypts_sum" | db.execute --quiet
+ SQL="create view wypts_sum as select w.[trans_id], w.$CLASSFIELD, count(w.cat) as n"
+ if [ -n "$WEIGHTS_LINES" ] ; then
+ SQL="$SQL ,sum(w.weight) as weight"
+ fi
+ SQL="$SQL from $WYPTS_IN w group by w.[trans_id],w.$CLASSFIELD"
+ g.message -v message="creating temporary paths table: $SQL"
+ echo $SQL | db.execute --quiet
+ ## now, calculate KIA (and optionally WKIA) joining the temporary views
+ SQL="insert into kia(id,"
+ SQL="$SQL $CLASSFIELD, "
+ SQL="$SQL KIA"
+ if [ -n "$WEIGHTS_LINES" ] ; then
+ SQL="$SQL, WKIA"
+ fi
+ SQL="$SQL) select w.[trans_id] as id, "
+ SQL="$SQL p.$CLASSFIELD, "
+ SQL="$SQL cast(w.n/p.length*1000 as numeric(10,6)) as KIA "
+ if [ -n "$WEIGHTS_LINES" ] ; then
+ SQL="$SQL, cast(w.weight/p.length*1000 as double numeric(10,6)) as WKIA"
+ fi
+ SQL="$SQL from paths_sum p "
+ SQL="$SQL left join wypts_sum w on p.[trans_id] = w.[trans_id] and p.$CLASSFIELD = w.$CLASSFIELD"
+ # once the SQL has been bulit, execute it
+ g.message -v message="calculating KIA by groups, SQL: $SQL"
+ echo "$SQL" | db.execute --quiet
+ ## remember to drop the intermediate tables paths_sum
+ echo "drop view if exists paths_sum" | db.execute --quiet
+ echo "drop view if exists wypts_sum" | db.execute --quiet
+fi
+
+# print out results, if --verbose
+if [ $((GRASS_VERBOSE)) -ge 3 ] ; then
+ db.select table=kia
+fi
+
+# add results to geometries
+if [ -n "$SEGMENT" ] ; then
+ echo "alter table $PATHS_IN add column join_fld varchar(255)" | db.execute
+ echo "update $PATHS_IN set join_fld=trans_id || $CLASSFIELD" | db.execute
+ echo "alter table kia add column join_f varchar(255)" | db.execute
+ echo "update kia set join_f=id || $CLASSFIELD" | db.execute
+ db.dropcol -f table=kia column=$CLASSFIELD --quiet
+ v.db.join map=$PATHS_IN column=join_fld otable=kia ocolumn=join_f --quiet
+ v.db.dropcol map=$PATHS_IN column=join_fld
+else
+ v.db.join map=$PATHS_IN column=trans_id otable=kia ocolumn=id --verbose
+ #v.db.dropcol map=$PATHS_IN column=id
+fi
+
+# finally, spit out a result vector
+if [ "$GIS_FLAG_S" == 1 ] ; then
+ if [ -f $GIS_OPT_OUTPUT ] ; then
+ g.message -w message="Shapefile <$GIS_OPT_OUTPUT> already exists and will be overwritten"
+ rm -rf `basename $GIS_OPT_OUTPUT .shp`.*
+ fi
+ if [ -n $DEM ] ; then
+ v.out.ogr -e input=$PATHS_IN type=line dsn=`dirname $GIS_OPT_OUTPUT` olayer=`basename $GIS_OPT_OUTPUT .shp` format=ESRI_Shapefile lco="SHPT=ARCZ" --quiet --overwrite
+ else
+ v.out.ogr -e input=$PATHS_IN type=line dsn=`dirname $GIS_OPT_OUTPUT` olayer=`basename $GIS_OPT_OUTPUT .shp` format=ESRI_Shapefile lco="SHPT=ARCZ" --quiet --overwrite
+ fi
+else
+ findelement vector $GIS_OPT_OUTPUT
+ if [ $findelement == 1 ] ; then
+ g.remove vect=$GIS_OPT_OUTPUT --quiet > /dev/null
+ fi
+ g.copy vect=$PATHS_IN,$GIS_OPT_OUTPUT --overwrite --quiet
+fi
+
+# exit nicely
+cleanup
+
+# use stopwatch
+end_time="$(date +%s)"
+elapsed_secs="$(expr $end_time - $start_time)"
+finished=`date`
+d_secs=$(($elapsed_secs%60))
+d_mins=$((($elapsed_secs/60)%60))
+d_hours=$(($elapsed_secs/3600))
+elapsed_time=`date -d $d_hours:$d_mins:$d_secs +%H:%M:%S`
+g.message -i message="Job finished at $finished. Elapsed time: $elapsed_time."
+
+exit 0
+### End of File ##
+# kate: encoding utf-8; syntax bash; space-indent on; indent-width 2;
+# kate: word-wrap-column 80; word-wrap-marker on; word-wrap-marker-color green;
+# kate: auto-brackets on; indent-mode cstyle;
+
Property changes on: grass-addons/vector/v.transect.kia/v.transect.kia
___________________________________________________________________
Added: svn:executable
+
Added: grass-addons/vector/v.transect.kia/v.transect.kia.html
===================================================================
--- grass-addons/vector/v.transect.kia/v.transect.kia.html (rev 0)
+++ grass-addons/vector/v.transect.kia/v.transect.kia.html 2010-01-08 16:26:39 UTC (rev 40323)
@@ -0,0 +1,144 @@
+<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
+<html>
+<head>
+<title>GRASS GIS manual: v.transect.kia</title>
+<meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1">
+<link rel="stylesheet" href="grassdocs.css" type="text/css">
+</head>
+<body bgcolor="white">
+
+<img src="grass_logo.png" alt="GRASS logo"><hr align=center size=6 noshade>
+
+<h2>NAME</h2>
+<em><b>v.transect.kia</b></em> - v.transect.kia<BR>
+Calculate Kilometric Abundance Indexes on line transect surveys of presence signs.
+<h2>KEYWORDS</h2>
+vector, elevation, geometry, line transect, abundance index, KIA
+<h2>SYNOPSIS</h2>
+<b>v.transect.kia</b><br>
+<b>v.transect.kia help</b><br>
+<b>v.transect.kia</b> [-<b>so</b>] <b>paths</b>=<em><paths map></em> <b>waypoints</b>=<em><waypoints map></em> [<b>idfield</b>=<em><transect ID field></em>] [<b>typefield</b>=<em><waypoint type field></em>] [<b>nfield</b>=<em><waypoint cluster size field></em>] <b>output</b>=<em><output map></em> [<b>elev</b>=<em><digital elevation model></em>] [<b>weights</b>=<em><waypoint weight file></em>] [<b>groups</b>=<em><paths segmentation map></em>] [<b>class</b>=<em><class values field></em>] [--<b>verbose</b>] [--<b>quiet</b>]
+
+<h3>Flags:</h3>
+<DL>
+<DT><b>-s</b></DT>
+<DD>assumes both paths and waypoints are in ESRI Shapefile format, does implicit conversion (using v.in.ogr -oe)</DD>
+
+<DT><b>-o</b></DT>
+<DD>overwrite output</DD>
+
+<DT><b>--verbose</b></DT>
+<DD>Verbose module output</DD>
+<DT><b>--quiet</b></DT>
+<DD>Quiet module output</DD>
+</DL>
+
+<h3>Parameters:</h3>
+<DL>
+<DT><b>paths</b>=<em><paths map></em></DT>
+<DD>a line vector map containing line transect paths</DD>
+
+<DT><b>waypoints</b>=<em><waypoints map></em></DT>
+<DD>a point vector containing objects recorded along transects</DD>
+
+<DT><b>idfield</b>=<em><transect ID field></em></DT>
+<DD>field (present both in paths and waypoints) with transect IDs</DD>
+<DD>Default: <em>trans_id</em></DD>
+
+<DT><b>typefield</b>=<em><waypoint type field></em></DT>
+<DD>field in waypoints map with point type information</DD>
+<DD>Default: <em>type</em></DD>
+
+<DT><b>nfield</b>=<em><waypoint cluster size field></em></DT>
+<DD>field in waypoints map containing number of items (for clustered waypoints)</DD>
+<DD>Default: <em>n</em></DD>
+
+<DT><b>output</b>=<em><output map></em></DT>
+<DD>output line vector map</DD>
+
+<DT><b>elev</b>=<em><digital elevation model></em></DT>
+<DD>digital elevation model to correct path lengths for elevation by draping</DD>
+
+<DT><b>weights</b>=<em><waypoint weight file></em></DT>
+<DD>text file with waypoint weights by type</DD>
+
+<DT><b>groups</b>=<em><paths segmentation map></em></DT>
+<DD>a vector (polygon) map to segment transect paths with and calculate partial KIA</DD>
+
+<DT><b>class</b>=<em><class values field></em></DT>
+<DD>field in groups map containing class values</DD>
+<DD>Default: <em>class</em></DD>
+
+</DL>
+
+
+<h2>DESCRIPTION</h2>
+
+<em>v.transect.kia.sh</em> calculates Kilometric Abundance Indexes (KIA) from
+a point vector representing sightings along a survey transect and a line vector representing transect paths. Both vectors are usually recorded in the field using a GPS.<br>
+Basically, a KIA is evaluated as<br><br>
+<center>KIA = sum(objects recorded along transect) / transect length</center>
+
+<h2>NOTES</h2>
+<em>v.transect.kia.sh</em> needs <tt>sqlite</tt> as a database backend. Make sure that your GRASS installation has been compiled with <tt>sqlite</tt> support. If the current database connection already uses <tt>sqlite</tt>, the existing connection will be used. On the contrary, a temporary <tt>sqlite</tt> connection will be created and used. Previous non-sqlite database connection will be restored at script end.
+
+<p>
+Both paths and waypoints maps must contain in their attribute table a field with the same name and type (<b>idfield</b>), containing transect unique identifiers. This field will be used to assign each waypoint to the transect along which it has been recorded. Automatic assignment e.g. based on the nearest transect will be implemented in the future,
+
+<h3>Weighted KIA</h3>
+In case of different kinds of signs, e.g. accounting for different number of individuals and/or different importance in defining a species' presence, waypoints can be diversified by adding a <em>type</em> field and using it with a <em>weight table</em>. In this case a <em>weighted KIA</em> will be calculated as:<br><br>
+<center>KIA = sum(object recorded along transect * object type weight) / transect length</center>
+<p>
+A <em>weight table</em> is a plain ASCII text file containing weight definitions, one per row, formatted as <em>type code( type code...) = weight</em>.<br>
+Two special cases are possible, i.e. use <em>N</em> or </em>Q</em> (uppercase!) as weights: if <em>weight = N</em>, values contained in the <em>cluster size field</em> <b>nfield</b> will be used as weights (and a field name present in the waypoints map must be specified as value for the <b>nfield</b> parameter). If <em>weight = Q</em> will be used, <b>nfield</b> must be supplied as well, and four classes, weighting respectively 1, 2, 3, and 4 will be automatically constructed, grouping <b>nfield</b> values by quartiles.<br>
+A sample <em>weight table</em> could be as follows:
+<div class="code"><pre>
+A CO I P PE S U = 1
+C = Q
+</pre></div>
+I.e., waypoints of type <em>A</em>, <em>CO</em>, <em>I</em>, <em>P</em>, <em>PE</em>, <em>S</em> and <em>U</em> will count "1", whereas <em>C</em> type waypoint will be count as 1, 2, 3,or 4 after a reclassification of <b>nfield</b> values into quartiles.
+
+<h3>Stratified KIA</h3>
+It is often needed to calculate "partial" KIAs, splitting each path by habitat types or elevation ranges (e.g. to check whether any differences in abundance exist on different habitat types or at different elevations). This is achieved supplying a polygon vector that classifies the area covered by paths. The program assumes that the GRASS vector map specified as <b>groups</b> parameter contains a field (either numeric or textual) named <em>class</em> with class codes for each polygon. If the field has a different name, use the <b>class</b> parameter.
+<p>
+If stratified KIAs are calculated, be warned that the process is significantly longer than other KIA calculations, since <tt>v.overlay</tt> has to be called repeatedly to break each path. It is known that <tt>v.overlay</tt> is slow, and even if called appropriately in couple with <tt>v.select</tt> or setting a region as small as possible, calculation time could increase appreciably.
+
+<h3>3D correction</h3>
+In the case of transect in montane areas, the effective path length could be underestimated and should be corrected by draping paths on a digital elevation model. This is achieved using the <b>elev</b> parameter and supplying a DEM at a siutable resolution. The DEM must be prepared in advance.
+
+
+<h2>EXAMPLES</h2>
+
+<h3>Calculate plain KIA</h3>
+A shapefile will be procuced, containing KIA values, directly from waypoints and tracks shapefile downloaded from a GPS unit.
+
+<div class="code"><pre>
+v.transect.kia -s paths=tracks.shp waypoints=wypts.shp output=kia.shp
+</pre></div>
+
+<h3>KIA by elevation ranges, 3D corrected</h3>
+Paths will be split based on a vector derived from a DEM (with <em>r.reclass</em> plus <em>r.to.vect</em>), which classifies the study area into 200 m wide elevation ranges.
+
+<div class="code"><pre>
+v.transect.kia paths=tracks waypoints=sightings output=kia-200m elev=srtm90 groups=elev200m class=height
+</pre></div>
+
+
+<h2>SEE ALSO</h2>
+
+<em>
+<a href="db.connect.html">db.connect</a>,
+<a href="r.reclass.html">r.reclass</a>,
+<a href="r.to.vect.html">r.to.vect</a>,
+<a href="v.in.ogr.html">v.in.ogr</a>,
+<a href="v.overlay.html">v.overlay</a>,
+<a href="v.select.html">v.select</a>
+</em>
+
+
+<h2>AUTHOR</h2>
+Clara Tattoni, Laboratorio di Ecologia, Dipartimento di Ingegneria Civile ed Ambientale, Università degli studi di Trento, Italy<br>
+Damiano G. Preatoni, Dipartimento Ambiente-Salute-Sicurezza, Università degli Studi dell'Insubria, Varese, Italy
+
+</body>
+</html>
More information about the grass-commit
mailing list