[GRASS-SVN] r45496 - in grass-addons/vector: . v.in.geodesic

svn_grass at osgeo.org svn_grass at osgeo.org
Tue Mar 1 02:47:52 EST 2011


Author: hamish
Date: 2011-02-28 23:47:52 -0800 (Mon, 28 Feb 2011)
New Revision: 45496

Added:
   grass-addons/vector/v.in.geodesic/
   grass-addons/vector/v.in.geodesic/Makefile
   grass-addons/vector/v.in.geodesic/description.html
   grass-addons/vector/v.in.geodesic/v.in.geodesic
Log:
+new module to create a vector line along a great circle path

Copied: grass-addons/vector/v.in.geodesic/Makefile (from rev 43872, grass-addons/vector/v.points.cog/Makefile)
===================================================================
--- grass-addons/vector/v.in.geodesic/Makefile	                        (rev 0)
+++ grass-addons/vector/v.in.geodesic/Makefile	2011-03-01 07:47:52 UTC (rev 45496)
@@ -0,0 +1,7 @@
+MODULE_TOPDIR = ../..
+
+PGM = v.in.geodesic
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script

Copied: grass-addons/vector/v.in.geodesic/description.html (from rev 43872, grass-addons/vector/v.points.cog/description.html)
===================================================================
--- grass-addons/vector/v.in.geodesic/description.html	                        (rev 0)
+++ grass-addons/vector/v.in.geodesic/description.html	2011-03-01 07:47:52 UTC (rev 45496)
@@ -0,0 +1,69 @@
+<h2>DESCRIPTION</h2>
+
+<em>v.in.geodesic</em> creates a new vector map containing a great circle
+line. The user can select to use either a starting and ending point,
+or a starting point with initial azimuth and line length.
+<P>
+The initial azimuth is measured in degrees clockwise from north.
+Line length and segment distance are both measured in meters.
+<!-- or ellipsoid units if that's not in meters (???) -->
+<P>
+If the <b>distance</b> parameter is given, polyline segments will be broken
+at that distance interval. Otherwise the line will be broken into
+<b>nsteps</b> equally spaced segments.
+
+
+<H2>NOTES</H2>
+
+Currently this module only works from a Lat/Long location.
+Use <em>v.proj</em> to reproject the results if needed.
+<P>
+This module requires the <tt>geod</tt> program from
+ <a href="http://proj.osgeo.org">PROJ.4</a> to be installed.
+ (often this is already installed with the GRASS package)
+
+
+<h2>EXAMPLES</h2>
+
+Create a geodesic line between two points:
+<div class="code"><pre>
+start_lat=-49
+start_lon=125
+end_lat=-32
+end_lon=180
+
+v.in.geodesic map=geodesic_line coor="$start_lon,$start_lat,$end_lon,$end_lat"
+</pre></div>
+
+
+Create a geodesic line from a starting point, initial azimuth of ENE
+(67.5 deg), and 500 km long:
+<div class="code"><pre>
+start_lat=-40.5
+start_lon=172.5
+azimuth=67.5
+length=500000
+
+v.in.geodesic map=geodesic_line line="$start_lon,$start_lat,$azimuth,$length"
+</pre></div>
+
+
+<H2>SEE ALSO</H2>
+<em>
+<A HREF="d.geodesic.html">d.geodesic</A><br>
+<A HREF="d.rhumbline.html">d.rhumbline</A><br>
+<A HREF="m.cogo.html">m.cogo</A><br>
+<A HREF="r.profile.html">r.profile</A><br>
+<A HREF="r.transect.html">r.transect</A>
+</em>
+
+<H2>AUTHOR</H2>
+
+Hamish Bowman<BR>
+<i>Dept. Marine Science<br>
+University of Otago<br>
+Dunedin, New Zealand</i>
+<BR>
+
+<p>
+<i>Last changed: $Date$</i>

Added: grass-addons/vector/v.in.geodesic/v.in.geodesic
===================================================================
--- grass-addons/vector/v.in.geodesic/v.in.geodesic	                        (rev 0)
+++ grass-addons/vector/v.in.geodesic/v.in.geodesic	2011-03-01 07:47:52 UTC (rev 45496)
@@ -0,0 +1,228 @@
+#!/bin/sh
+############################################################################
+#
+# MODULE:       v.in.geodesic
+# AUTHOR(S):    M. Hamish Bowman, Dunedin, New Zealand
+# PURPOSE:      Creates a vector map containing a great circle line.
+# COPYRIGHT:    (C) 2011 Hamish Bowman, and the GRASS Development Team
+#
+#  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.
+#
+############################################################################
+#%Module
+#% description: Creates a vector map containing a great circle line.
+#% keywords: vector, geodesic, great circle
+#%End
+#%Option
+#% key: map
+#% type: string
+#% required: yes
+#% multiple: no
+#% key_desc: name
+#% description: Name for output vector map
+#% gisprompt: new,vector,vector
+#%End
+#%Option
+#% key: coor
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: east1,north1,east2,north2
+#% description: Starting and ending coordinates
+#%End
+#%Option
+#% key: line
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: east,north,azimuth,distance
+#% description: Transect definition
+#%End
+
+# fixme: geod uses as radius for +a distance etc. as well as output. (??)
+##%Option
+##% key: units
+##% type: string
+##% required: no
+##% multiple: no
+##% options: meters,kilometers,feet,miles,nautmiles
+##% description: Units of distance used in transect
+##% answer: meters
+##%End
+
+#%Option
+#% key: nsteps
+#% type: integer
+#% required: no
+#% multiple: no
+#% description: Number of vertices to place along the line
+#% answer: 100
+#% options: 1-10000
+#%End
+#%Option
+#% key: distance
+#% type: double
+#% required: no
+#% multiple: no
+#% description: Distance between vertices in map units
+#%End
+
+
+
+if  [ -z "$GISBASE" ] ; then
+    echo "You must be in GRASS GIS to run this program." 1>&2
+    exit 1
+fi
+
+# save command line
+if [ "$1" != "@ARGS_PARSED@" ] ; then
+    CMDLINE="`basename $0`"
+    for arg in "$@" ; do
+        CMDLINE="$CMDLINE \"$arg\""
+    done
+    export CMDLINE
+    exec g.parser "$0" "$@"
+fi
+
+#### check if we have geod
+if [ ! -x "`which geod`" ] ; then
+    g.message -e "\"geod\" from PROJ.4 is required, please install it first"
+    exit 1
+fi
+
+#### check if we have awk
+if [ ! -x "`which awk`" ] ; then
+    g.message -e "\"awk\" is required, please install it first"
+    exit 1
+fi
+
+#### set environment so that awk works properly in all languages
+unset LC_ALL
+LC_NUMERIC=C
+export LC_NUMERIC
+
+
+if [ -z "$GIS_OPT_COOR" ] && [ -z "$GIS_OPT_LINE" ] ; then
+    g.message -e "Either \"coor\" or \"line\" option must be given"
+    exit 1
+fi
+
+if [ -n "$GIS_OPT_COOR" ] && [ -n "$GIS_OPT_LINE" ] ; then
+    g.message -e "Only one of \"coor\" or \"line\" must be given"
+    exit 1
+fi
+
+#nsteps always present currently
+#if [ -n "$GIS_OPT_NSTEPS" ] && [ -n "$GIS_OPT_DISTANCE" ] ; then
+#    g.message -e "Only one of \"nsteps\" or \"distance\" must be given"
+#    exit 1
+#fi
+
+#are we in LatLong location?
+g.proj -p | grep -i name | grep -i Lon  > /dev/null
+if [ $? -eq 1 ] ; then
+   g.message -e "This module only operates in Lat-Long locations"
+   exit 1
+fi
+
+
+#### setup temporary file
+# (if needed)
+#TEMPFILE="`g.tempfile pid=$$`"
+#if [ $? -ne 0 ] || [ -z "$TEMPFILE" ] ; then
+#   g.message -e "Unable to create temporary file" 
+#   exit 1
+#fi
+
+
+# fixme: geod uses as radius for +a distance etc. as well as output. (??)
+#case "$GIS_OPT_UNITS" in
+#   meters)
+#      UNITS=m
+#      ;;
+#   kilometers)
+#      UNITS=km
+#      ;;
+#   feet)
+#      UNITS=ft
+#      ;;
+#   miles)
+#      UNITS=mi
+#      ;;
+#   nautmiles)
+#      UNITS=kmi
+#      ;;
+#   *)
+#      g.message -e "A programming error occurred."
+#      exit 1
+#      ;;
+#esac
+
+# get ellipsoid parameters in +proj4 terms
+PROJ4_ELLPS=`g.proj -j | grep '^+ellps=\|^+a=\|^+b=\|^+rf=\|^+es=' | tr '\n' ' '`
+
+if [ -n "$GIS_OPT_COOR" ] ; then
+   LON1=`echo "$GIS_OPT_COOR" | cut -f1 -d','`
+   LAT1=`echo "$GIS_OPT_COOR" | cut -f2 -d','`
+   LON2=`echo "$GIS_OPT_COOR" | cut -f3 -d','`
+   LAT2=`echo "$GIS_OPT_COOR" | cut -f4 -d','`
+elif [ -n "$GIS_OPT_LINE" ] ; then
+   LON1=`echo "$GIS_OPT_LINE" | cut -f1 -d','`
+   LAT1=`echo "$GIS_OPT_LINE" | cut -f2 -d','`
+   AZIM=`echo "$GIS_OPT_LINE" | cut -f3 -d','`
+   DIST=`echo "$GIS_OPT_LINE" | cut -f4 -d','`
+fi
+
+
+#### TODO:
+# convert D:M:S to Dm's" with sed  (azimuth too) (see m.proj code)
+# convert non-lat/lon to lat/lon with m.proj
+
+
+if [ -n "$GIS_OPT_DISTANCE" ] ; then
+   # default to region resolution? (avg nsres,ewres)
+   STEP="+del_S=$GIS_OPT_DISTANCE"
+else
+   STEP="+n_S=$GIS_OPT_NSTEPS"
+fi
+
+
+if [ -n "$GIS_OPT_COOR" ] ; then
+
+   geod $PROJ4_ELLPS +lon_1="$LON1" +lat_1="$LAT1" \
+       +lon_2="$LON2" +lat_2="$LAT2" "$STEP" -f "%.9f" | \
+     awk '{print $2,$1}' | \
+     v.in.lines input=- output="$GIS_OPT_MAP" --quiet
+
+elif [ -n "$GIS_OPT_LINE" ] ; then
+
+   geod $PROJ4_ELLPS +lon_1="$LON1" +lat_1="$LAT1" \
+       +S="$DIST" +A="$AZIM" "$STEP" -f "%.9f" | \
+     awk '{print $2,$1}' | \
+     v.in.lines input=- output="$GIS_OPT_MAP" --quiet
+   #  +units="$UNITS"
+fi
+
+# convert output back to local projection with m.proj if needed.
+
+
+# assign our new line a cat number
+v.edit -1 map="$GIS_OPT_MAP" type=line tool=catadd id=1 cat=1 --quiet
+
+# save command line history
+v.support map="$GIS_OPT_MAP" cmdhist="$CMDLINE"
+
+LINE_LEN=`v.to.db -p "$GIS_OPT_MAP" option=length \
+   units=meters --quiet | cut -f2 -d'|' `
+#   units="$GIS_OPT_UNITS" --quiet | cut -f2 -d'|' `
+
+g.message message="Line length:  $LINE_LEN m"  # $GIS_OPT_UNITS
+g.message message="`basename $0` complete."


Property changes on: grass-addons/vector/v.in.geodesic/v.in.geodesic
___________________________________________________________________
Added: svn:executable
   + *
Added: svn:mime-type
   + text/x-sh
Added: svn:eol-style
   + native



More information about the grass-commit mailing list