[GRASS-SVN] r59790 - in grass-addons/grass6: raster/r.surf.nnbathy vector vector/v.surf.nnbathy

svn_grass at osgeo.org svn_grass at osgeo.org
Fri Apr 18 19:24:28 PDT 2014


Author: hamish
Date: 2014-04-18 19:24:28 -0700 (Fri, 18 Apr 2014)
New Revision: 59790

Added:
   grass-addons/grass6/vector/v.surf.nnbathy/
   grass-addons/grass6/vector/v.surf.nnbathy/v.surf.nnbathy
Removed:
   grass-addons/grass6/vector/v.surf.nnbathy/r.surf.nnbathy
Modified:
   grass-addons/grass6/raster/r.surf.nnbathy/description.html
   grass-addons/grass6/vector/v.surf.nnbathy/Makefile
   grass-addons/grass6/vector/v.surf.nnbathy/description.html
Log:
+vector & ascii file variant of r.surf.nnbathy (saves some steps)

Modified: grass-addons/grass6/raster/r.surf.nnbathy/description.html
===================================================================
--- grass-addons/grass6/raster/r.surf.nnbathy/description.html	2014-04-19 02:20:45 UTC (rev 59789)
+++ grass-addons/grass6/raster/r.surf.nnbathy/description.html	2014-04-19 02:24:28 UTC (rev 59790)
@@ -3,7 +3,7 @@
 <em>r.surf.nnbathy</em> is an Bourne Shell and Awk script. It is an
 interface between the external <em>nnbathy</em> utility and <em>GRASS</em>.
 <em>nnbathy</em> is a surface interpolation program provided with
-<a href="https://code.google.com/p/nn-c/"">nn</a> - a natural neighbor
+<a href="https://code.google.com/p/nn-c/">nn</a> - a natural neighbor
 interpolation library, written by Pavel Sakov.
 
 <p>
@@ -67,7 +67,10 @@
 
 <h2>SEE ALSO</h2>
 
-<em><a href="v.to.rast.html">v.to.rast</a></em>
+<em>
+<a href="v.to.rast.html">v.to.rast</a>,
+<a href="v.surf.nnbathy.html">v.surf.nnbathy</a> (AddOn)
+</em>
 
 
 <h2>AUTHOR</h2>

Modified: grass-addons/grass6/vector/v.surf.nnbathy/Makefile
===================================================================
--- grass-addons/grass6/raster/r.surf.nnbathy/Makefile	2013-08-27 15:29:31 UTC (rev 57519)
+++ grass-addons/grass6/vector/v.surf.nnbathy/Makefile	2014-04-19 02:24:28 UTC (rev 59790)
@@ -1,6 +1,6 @@
 MODULE_TOPDIR = ../..
 
-PGM = r.surf.nnbathy
+PGM = v.surf.nnbathy
 
 include $(MODULE_TOPDIR)/include/Make/Script.make
 

Modified: grass-addons/grass6/vector/v.surf.nnbathy/description.html
===================================================================
--- grass-addons/grass6/raster/r.surf.nnbathy/description.html	2013-08-27 15:29:31 UTC (rev 57519)
+++ grass-addons/grass6/vector/v.surf.nnbathy/description.html	2014-04-19 02:24:28 UTC (rev 59790)
@@ -1,13 +1,13 @@
 <h2>DESCRIPTION</h2>
 
-<em>r.surf.nnbathy</em> is an Bourne Shell and Awk script. It is an
+<em>v.surf.nnbathy</em> is an Bourne Shell and Awk script. It is an
 interface between the external <em>nnbathy</em> utility and <em>GRASS</em>.
 <em>nnbathy</em> is a surface interpolation program provided with
-<a href="https://code.google.com/p/nn-c/"">nn</a> - a natural neighbor
+<a href="https://code.google.com/p/nn-c/">nn</a> - a natural neighbor
 interpolation library, written by Pavel Sakov.
 
 <p>
-<em>r.surf.nnbathy</em> provides 3 interpolation algorithms. According to
+<em>v.surf.nnbathy</em> provides 3 interpolation algorithms. According to
 <em>nn</em> library documentation these are: Delaunay interpolation
 (<b>alg=l</b>), Watson's algortithm for Sibson natural neighbor
 interpolation (<b>alg=nn</b>) and Belikov and Semenov's algorithm for
@@ -18,20 +18,21 @@
 
 <p>
 The <b>output</b> raster map is a continous surface interpolated from the
-<b>input</b> raster map.
+<b>input</b> or <b>file</b> data.
 
-<p>
+
+<h2>NOTES</h2>
+
 <em>nnbathy</em>, if built with '-DNN_SERIAL' (default as of nn 1.85), is
 able to create a grid of virtually any size. It interpolates and writes one
 output point at a time only. This eliminates the necessity to hold the whole
-output array in memory. However, even then all the input cells are still
+output array in memory. However, even then all the input points are still
 held in the memory.
 
+<p>
 
-<h2>NOTES</h2>
-
 <ul>
-1. Requires <em>GRASS</em> 6.x and <em>nnbathy</em> 1.76 or greater.<br>
+1. Requires <em>GRASS</em> 6 and <em>nnbathy</em> 1.76 or greater.<br>
 
 2. Build <em>nnbathy</em> according to instructions provided with its source
 code and put it somewhere in your $PATH.<br>
@@ -42,37 +43,35 @@
 4. The output raster map non-NULL area is limited to the convex hull
 encompassing all the non-NULL input cells.<br>
 
-5. The output is double floating point raster map.<br>
+5. The output is double precision floating point raster map (DCELL).<br>
 
 6. Natural neighbor is a an <em>exact</em> interpolation algorithm, so all
-non-NULL input cells have their value exactly preserved in the output.<br>
+non-NULL input points have their value exactly preserved in the output.<br>
 
 7. There is circa 0.2 KB memory overhead per each <em>input</em> cell.
 However, the <em>output</em> grid can be of any size, if <em>nnbathy</em> is
 built with -DNN_SERIAL switch.<br>
 
-8. <em>r.surf.nnbathy</em> creates 3 temporary files: ASCII x,y,z lists of
-the input and output cells, and the output list converted into GRASS ASCII
+8. <em>v.surf.nnbathy</em> creates 3 temporary files: ASCII x,y,z lists of
+the input points and output cells, and the output list converted into GRASS ASCII
 format. Then it makes a GRASS raster map from the latter - and only then it
 removes the 3 temp files, when the script terminates. Thus, at the script
 run time several times more disk space might be required, than the final
 GRASS raster map would actually occupy.<br>
-
-<p>
-I'd like to thank Pavel Sakov for his help, and especially for implementing
-serial input processing.
-
 </ul>
 
 
 <h2>SEE ALSO</h2>
 
-<em><a href="v.to.rast.html">v.to.rast</a></em>
+<em><a href="r.surf.nnbathy.html">r.surf.nnbathy</a></em> (GRASS Add-on script)
 
 
 <h2>AUTHOR</h2>
 
-Maciej Sieczka
+Hamish Bowman, Otago University, New Zealand<br>
+Based on <em>r.surf.nnbathy</em> by Maciej Sieczka<br>
+<em><a href="https://code.google.com/p/nn-c/">nnbathy</a></em>
+ written by Pavel Sakov
 
 <p>
 <i>Last changed: $Date$</i>

Deleted: grass-addons/grass6/vector/v.surf.nnbathy/r.surf.nnbathy
===================================================================
--- grass-addons/grass6/raster/r.surf.nnbathy/r.surf.nnbathy	2013-08-27 15:29:31 UTC (rev 57519)
+++ grass-addons/grass6/vector/v.surf.nnbathy/r.surf.nnbathy	2014-04-19 02:24:28 UTC (rev 59790)
@@ -1,265 +0,0 @@
-#!/bin/sh
-
-############################################################################
-#
-# MODULE:     r.surf.nnbathy
-#
-# AUTHOR(S):  Maciej Sieczka
-#
-# PURPOSE:	  Interpolate raster surface using the "nnbathy" natural
-#             neighbor interpolation program.
-#
-# VERSION:	  1.95, developed over GRASS 6.3 CVS
-#
-# COPYRIGHT:  (c) Maciej Sieczka
-#
-# LICENSE:    This program is free software under the GNU General Public
-#             License (>=v2). Read the file COPYING that comes with GRASS
-#             for details.
-#
-#############################################################################
-
-# NOTES:
-#
-# 1. Requires nnbathy executable v 1.75 or later. Follow the instruction in
-#    html manual page to obtain it.
-#
-# 2. When creating the input raster map, make sure it's extents cover
-#    your whole region of interest, the cells you wish to interplate on are
-#    NULL, and the resolution is correct. Same as most GRASS raster modules
-#    this one is region sensitive too.
-
-# CHANGELOG:
-
-# 1.6, 2006.06.15:
-# - first public release
-
-# 1.9, 2007.01.02:
-# - parse g.region -g with eval, not awk
-# - support all interpolation methods nnbathy provides
-# - create history for output raster with r.support
-# - require nnbathy 1.69 (major bug was fixed)
-# - try to detect if nnbathy failed and exit cleanly
-# - documentation extended
-# - minor cleanups
-# - todo: there is a progress indicator switch (-%) in nnbathy now, but using it
-#   slows down the interpolation 3-4 times on my machine, while works OK on nnbathy
-#   author's; if sorted out, I'll use -%
-
-# 1.95, 2007.11.12:
-# - require nnbathy 1.75 - bugfixes and speed improvemenets for large grids;
-#   refer to file CHANGELOG in nn sorce code for details
-
-# 1.96, 2008.25.03:
-# - handle spaces in pathnames
-# - get rid of Bashims
-# - better run-time error trapping
-# - require nnbathy 1.76 - contains a bugfix; refer to file CHANGELOG in nn
-#   sorce code for details
-# - cosmetics
-# - manual cleaned up
-
-#%Module
-#%  description: Interpolate raster using the nnbathy natural neighbor interpolation program.
-#%END
-
-#%option
-#% key: input
-#% type: string
-#% gisprompt: old,cell,raster
-#% description: The raster map to interpolate on
-#% required : yes
-#%END
-
-#%option
-#% key: output
-#% gisprompt: new,cell,raster
-#% type: string
-#% description: Name of the output raster map
-#% required : yes
-#%END
-
-#%option
-#% key: alg
-#% type: string
-#% options: l,nn,ns
-#% answer: nn
-#% description: Interpolation algorithm for nnbathy to use
-#% required : yes
-#%END
-
-
-# called from GRASS?
-if [ -z "$GISBASE" ]; then
-  echo
-  echo "ERROR: You must be in GRASS GIS to run this program." 1>&2
-  echo
-  exit 1
-fi
-
-if [ "$1" != "@ARGS_PARSED@" ]; then
-  exec g.parser "$0" "$@"
-fi
-
-# check if we have awk
-if [ ! -x "`which awk`" ]; then
-  g.message -e '"awk" executable required but not found.' 1>&2
-  exit 1
-fi
-
-# check if we have nnbathy
-if [ ! -x "`which nnbathy`" ]; then
-  echo
-  g.message -e '"nnbathy" executable required but not found. Follow the instructions in r.surf.nnbathy manual to install it.' 1>&2
-  echo
-  exit 1
-fi
-
-# check nnbathy version
-nnv=`nnbathy -v | sed 's/ /\n/g' | sort -nr | head -n1`
-nnv_ok=`echo $nnv | awk '{ if ($0<1.76) {print 0} else {print 1} }'`
-
-if [ $nnv_ok -eq 0 ]; then
-  g.message -e '"nnbathy" version >= 1.76 is required.'
-  exit 1
-fi
-
-# set up temporary files
-TMP="`g.tempfile pid=$$`"
-if [ $? -ne 0 ] || [ -z "$TMP" ]; then
-  echo
-  g.message -e 'Unable to create temporary files.' 1>&2
-  echo
-  exit 1
-fi
-
-# set environment so that awk works properly in all languages
-unset LC_ALL
-LC_NUMERIC=C
-export LC_NUMERIC
-
-PROG=`basename $0`
-
-# define cleanup procedure
-proc_cleanup() {
-  # Reset traps before normal script termination to avoid bogus ERROR message, as
-  # we put a trap on signal 0.
-  trap - 0 2 3 15 
-  rm -f "$TMP" "$TMP.${PROG}.input_xyz" "$TMP.${PROG}.output_xyz" "$TMP.${PROG}.output_grd"
-}
-
-# define run-time error handling procedure
-proc_runtime_error() {
-  echo
-  g.message -e "There was an error at the script's run time. Please try to debug the problem if you can and let me know by email." 1>&2
-  echo
-  exit 1
-}
-
-# define user-break procedure
-proc_user_break() {
-  echo
-  g.message -w "User break!"
-  echo
-  proc_cleanup
-  exit 1
-}
-
-# set trap for when script terminates
-trap "proc_runtime_error" 0
-
-# trap user break
-trap "proc_user_break" 2 3 15
-
-# assign main variables from user's input
-INPUT="$GIS_OPT_INPUT"
-OUTPUT="$GIS_OPT_OUTPUT"
-ALG=$GIS_OPT_ALG
-
-
-
-### DO IT ###
-
-# Make script terminate (ie. emit signal 0) if any statement returns a
-# non-0 value. Then we trap signal 0, which lets handle such errors.
-# However, the trap on signal 0 must be reset before normal script
-# termination to avoid a bogus ERROR message then - this is done in the
-# cleanup procedure.
-set -e
-
-# grab the current region settings
-eval `g.region -gp`
-
-# spit out non-null (-n) raster coords + values to be interpolated
-r.stats --q -1gn input="${INPUT}" > "$TMP.${PROG}.input_xyz"
-
-# set the working region for nnbathy (it's cell-center oriented)
-nn_n=`echo $n | awk -v res="$nsres" '{printf "%.8f",$1-res/2.0}'`
-nn_s=`echo $s | awk -v res="$nsres" '{printf "%.8f",$1+res/2.0}'`
-nn_w=`echo $w | awk -v res="$ewres" '{printf "%.8f",$1+res/2.0}'`
-nn_e=`echo $e | awk -v res="$ewres" '{printf "%.8f",$1-res/2.0}'`
-
-null=NaN
-type=double
-
-# interpolate
-echo
-g.message -w '"nnbathy" is performing the interpolation now. This may take some time.' 1>&2
-echo
-g.message -w "Once it completes 'All done.' message will be printed." 1>&2
-echo
-
-nnbathy -W 0 -P alg=$ALG -n ${cols}x${rows} -x $nn_w $nn_e -y $nn_n $nn_s -i "$TMP.${PROG}.input_xyz" > "$TMP.${PROG}.output_xyz"
-# Y in "r.stats -1gn" output is in descending order, thus -y must be in MAX MIN order, not MIN MAX, for nnbathy not to produce a grid upside-down
-
-# convert the X,Y,Z nnbathy output into a GRASS ASCII grid, then import with r.in.ascii:
-
-# 1 create header
-cat << EOF > "$TMP.${PROG}.output_grd"
-north: $n
-south: $s
-east: $e
-west: $w
-rows: $rows
-cols: $cols
-null: $null
-type: $type
-EOF
-
-# 2 do the conversion
-echo "Converting nnbathy output to GRASS raster." 1>&2
-echo
-
-awk -v cols="$cols" '
-BEGIN {col_cur=1; ORS=" "}
-{
- if (col_cur==cols) {ORS="\n"; col_cur=0; print $3; ORS=" "}
-		    else {print $3}
- col_cur++
-}' "$TMP.${PROG}.output_xyz" >> "$TMP.${PROG}.output_grd"
-
-# 3 import
-r.in.ascii input="$TMP.${PROG}.output_grd" output="${OUTPUT}" > /dev/null
-
-# store comand history in raster's metadata
-
-r.support map=${OUTPUT} history=""
-r.support map=${OUTPUT} history="script run syntax:"
-r.support map=${OUTPUT} history=""
-r.support map=${OUTPUT} history="r.surf.nnbathy alg=$ALG input=${INPUT} output=${OUTPUT}"
-r.support map=${OUTPUT} history=""
-r.support map=${OUTPUT} history="nnbathy run syntax:"
-r.support map=${OUTPUT} history=""
-r.support map=${OUTPUT} history="nnbathy -W 0 -P alg=$ALG -n ${cols}x${rows} "
-r.support map=${OUTPUT} history="-x $nn_w $nn_e "
-r.support map=${OUTPUT} history="-y $nn_n $nn_s "
-r.support map=${OUTPUT} history="-i tmp_in > tmp_out"
-r.support map=${OUTPUT} history=""
-
-### ALL DONE ###
-
-proc_cleanup
-
-echo
-echo "All done." 1>&2
-echo

Copied: grass-addons/grass6/vector/v.surf.nnbathy/v.surf.nnbathy (from rev 57519, grass-addons/grass6/raster/r.surf.nnbathy/r.surf.nnbathy)
===================================================================
--- grass-addons/grass6/vector/v.surf.nnbathy/v.surf.nnbathy	                        (rev 0)
+++ grass-addons/grass6/vector/v.surf.nnbathy/v.surf.nnbathy	2014-04-19 02:24:28 UTC (rev 59790)
@@ -0,0 +1,316 @@
+#!/bin/sh
+
+############################################################################
+#
+# MODULE:     v.surf.nnbathy
+#
+# AUTHOR(S):  Hamish Bowman, Otago University, New Zealand
+#		Based on r.surf.nnbathy by Maciej Sieczka
+#
+# PURPOSE:    Interpolate raster surface using the "nnbathy" natural
+#             neighbor interpolation program.
+#
+# COPYRIGHT:  (c) 2014 Hamish Bowman, Maciej Sieczka, and the GRASS Development Team
+#
+# LICENSE:    This program is free software under the GNU General Public
+#             License (>=v2). Read the file COPYING that comes with GRASS
+#             for details.
+#
+#############################################################################
+
+# NOTES:
+#
+# 1. Requires nnbathy executable v 1.75 or later. Follow the instruction in
+#    html manual page to obtain it.
+#
+# 2. When creating the input raster map, make sure it's extents cover
+#    your whole region of interest, the cells you wish to interplate on are
+#    NULL, and the resolution is correct. Same as most GRASS raster modules
+#    this one is region sensitive too.
+
+
+#%Module
+#% description: Interpolates a raster map using the nnbathy natural neighbor interpolation program.
+#%End
+#%Option
+#% key: input
+#% type: string
+#% gisprompt: old,vector,vector
+#% description: Name of input vector points map
+#% guisection: Selection
+#% required : no
+#%End
+#%Option
+#% key: output
+#% gisprompt: new,cell,raster
+#% type: string
+#% description: Name for output raster map
+#% required : yes
+#%End
+#%Option
+#% key: zcolumn
+#% type: string
+#% required: no
+#% multiple: no
+#% description: Name of the attribute column with values to be used for approximation (if layer>0)
+#% guisection: Selection
+#%End
+#%Option
+#% key: layer
+#% type: integer
+#% required: no
+#% multiple: no
+#% label: Layer number
+#% description: If set to 0, z coordinates are used. (3D vector only)
+#% answer: 1
+#% gisprompt: old_layer,layer,layer_zero
+#% guisection: Selection
+#%End
+#%Option
+#% key: where
+#% type: string
+#% required: no
+#% multiple: no
+#% key_desc: sql_query
+#% label: WHERE conditions of SQL statement without 'where' keyword
+#% description: Example: income < 1000 and inhab >= 10000
+#% guisection: Selection
+#%End
+#%option
+#% key: alg
+#% type: string
+#% options: l,nn,ns
+#% answer: nn
+#% descriptions: l;Linear;nn;Sibson natural neighbor;ns;Non-Sibsonian natural neighbor
+#% description: Interpolation algorithm
+#%End
+#%Option
+#% key: file
+#% type: string
+#% gisprompt: old_file,file,input
+#% label: Name of input text file
+#% description: Containing x,y,z data as three space separated columns
+#% required : no
+#%End
+
+
+if [ -z "$GISBASE" ]; then
+  echo "ERROR: You must be in GRASS GIS to run this program." 1>&2
+  exit 1
+fi
+
+if [ "$1" != "@ARGS_PARSED@" ]; then
+  exec g.parser "$0" "$@"
+fi
+
+if [ -n "$GIS_OPT_INPUT" ] && [ -n "$GIS_OPT_FILE" ] ; then
+  g.message -e "Please specify either the 'input' or 'file' option, not both."
+  exit 1
+fi
+
+if [ -z "$GIS_OPT_INPUT" ] && [ -z "$GIS_OPT_FILE" ] ; then
+  g.message -e "Please specify either the 'input' or 'file' option."
+  exit 1
+fi
+
+if [ -n "$GIS_OPT_FILE" ] && [ ! -e "$GIS_OPT_FILE" ] ; then
+  g.message -e "File '$GIS_OPT_FILE' does not exist."
+  exit 1
+fi
+
+# check if we have awk
+if [ ! -x "`which awk`" ]; then
+  g.message -e '"awk" executable required but not found.'
+  exit 1
+fi
+
+# check if we have nnbathy
+if [ ! -x "`which nnbathy`" ]; then
+  g.message -e '"nnbathy" executable required but not found. Follow the instructions in v.surf.nnbathy manual to install.'
+  exit 1
+fi
+
+# check nnbathy version
+nnv=`nnbathy -v | sed 's/ /\n/g' | sort -nr | head -n1`
+nnv_ok=`echo $nnv | awk '{ if ($0<1.76) {print 0} else {print 1} }'`
+
+if [ $nnv_ok -eq 0 ]; then
+  g.message -e '"nnbathy" version >= 1.76 is required.'
+  exit 1
+fi
+
+# set up temporary files
+TMP="`g.tempfile pid=$$`"
+if [ $? -ne 0 ] || [ -z "$TMP" ]; then
+  g.message -e 'Unable to create temporary files.' 1>&2
+  exit 1
+fi
+
+# set environment so that awk works properly in all languages
+unset LC_ALL
+LC_NUMERIC=C
+export LC_NUMERIC
+
+PROG=`basename $0`
+
+# define cleanup procedure
+proc_cleanup() {
+  # Reset traps before normal script termination to avoid bogus ERROR message, as
+  # we put a trap on signal 0.
+  trap - 0 2 3 15 
+  rm -f "$TMP" "$TMP.$PROG.input_xyz" "$TMP.$PROG.output_xyz" "$TMP.$PROG.output_grd"
+}
+
+# define run-time error handling procedure
+proc_runtime_error() {
+  g.message -e "There was an error at the script's run time. Please try to debug the problem if you can and let the GRASS developers know."
+  exit 1
+}
+
+# define user-break procedure
+proc_user_break() {
+  g.message -e "User break!"
+  proc_cleanup
+  exit 1
+}
+
+# set trap for when script terminates
+trap "proc_runtime_error" 0
+
+# trap user break
+trap "proc_user_break" 2 3 15
+
+
+# assign main variables from user's input
+INPUT="$GIS_OPT_INPUT"
+OUTPUT="$GIS_OPT_OUTPUT"
+ALG="$GIS_OPT_ALG"
+
+
+### DO IT ###
+
+# Make script terminate (ie. emit signal 0) if any statement returns a
+# non-0 value. Then we trap signal 0, which lets handle such errors.
+# However, the trap on signal 0 must be reset before normal script
+# termination to avoid a bogus ERROR message then - this is done in the
+# cleanup procedure.
+set -e
+
+# grab the current region settings
+eval `g.region -g`
+
+if [ -n "$GIS_OPT_INPUT" ] ; then
+  # spit out coords + values to be interpolated
+  ## FIXME: check if data is 3D or use output column
+  # column=  (slow!!)
+  TMPXYZ="$TMP.$PROG.input_xyz"
+  
+  if [ "$GIS_OPT_LAYER" -eq 0 ] ; then
+    # 3D vector points
+    LAYER=""
+    COLUMN=""
+  else
+    LAYER="layer=$GIS_OPT_LAYER"
+    if [ -n "$GIS_OPT_ZCOLUMN" ] ; then
+      COLUMN="columns=$GIS_OPT_ZCOLUMN"
+    else
+      g.message -e "Name of z column required for 2D vector maps."
+      exit 1
+    fi
+  fi
+
+  if [ -n "$GIS_OPT_WHERE" ] ; then
+    v.out.ascii -r input="$INPUT" output="$TMPXYZ" format=point fs=space \
+       dp=15 where="$GIS_OPT_WHERE" $LAYER $COLUMN
+  else
+    v.out.ascii -r input="$INPUT" output="$TMPXYZ" format=point fs=space \
+       dp=15 $LAYER $COLUMN
+  fi
+
+  if [ "$GIS_OPT_LAYER" -gt 0 ] ; then
+    # crop away the cat column
+    mv "$TMPXYZ" "$TMPXYZ.cat"
+    cut -f1,2,4 -d' ' "$TMPXYZ.cat" > "$TMPXYZ"
+    rm "$TMPXYZ.cat"
+  fi
+
+else
+  TMPXYZ="$GIS_OPT_FILE"
+fi
+
+
+# set the working region for nnbathy (it's cell-center oriented)
+nn_n=`echo $n | awk -v res="$nsres" '{printf "%.8f", $1 - res/2.0}'`
+nn_s=`echo $s | awk -v res="$nsres" '{printf "%.8f", $1 + res/2.0}'`
+nn_w=`echo $w | awk -v res="$ewres" '{printf "%.8f", $1 + res/2.0}'`
+nn_e=`echo $e | awk -v res="$ewres" '{printf "%.8f", $1 - res/2.0}'`
+
+null="NaN"
+type="double"
+
+
+#### interpolate
+g.message '"nnbathy" is performing the interpolation now. This may take some time.'
+g.message "Once it completes an 'All done.' message will be printed."
+
+
+nnbathy -W 0 -P alg="$ALG" -n "${cols}x${rows}" \
+	-x "$nn_w" "$nn_e" -y "$nn_n" "$nn_s" \
+	-i "$TMPXYZ" \
+       > "$TMP.$PROG.output_xyz"
+
+
+# Y in "r.stats -1gn" output is in descending order, thus -y must be in
+# MAX MIN order, not MIN MAX, for nnbathy not to produce a grid upside-down
+
+# convert the X,Y,Z nnbathy output into a GRASS ASCII grid, then import with r.in.ascii:
+
+# 1 create header
+cat << EOF > "$TMP.$PROG.output_grd"
+north: $n
+south: $s
+east: $e
+west: $w
+rows: $rows
+cols: $cols
+null: $null
+type: $type
+EOF
+
+# 2 do the conversion
+g.message "Converting nnbathy output to GRASS raster ..."
+
+awk -v cols="$cols" '
+BEGIN {col_cur=1; ORS=" "}
+{
+ if (col_cur==cols) {ORS="\n"; col_cur=0; print $3; ORS=" "}
+		    else {print $3}
+ col_cur++
+}' "$TMP.$PROG.output_xyz" >> "$TMP.$PROG.output_grd"
+
+# 3 import
+r.in.ascii input="$TMP.$PROG.output_grd" output="$OUTPUT" --quiet
+
+# store comand history in raster's metadata
+
+if [ -n "$GIS_OPT_INPUT" ] ; then
+  r.support map="$OUTPUT" history="v.surf.nnbathy alg=$ALG input=$INPUT output=$OUTPUT"
+  #TODO: layer= where= zcolumn=
+else
+  r.support map="$OUTPUT" history="v.surf.nnbathy alg=$ALG input=$FILE output=$OUTPUT"
+fi
+r.support map="$OUTPUT" history=""
+r.support map="$OUTPUT" history="nnbathy run syntax:"
+r.support map="$OUTPUT" history=""
+r.support map="$OUTPUT" history="nnbathy -W 0 -P alg=$ALG -n ${cols}x$rows "
+r.support map="$OUTPUT" history="-x $nn_w $nn_e "
+r.support map="$OUTPUT" history="-y $nn_n $nn_s "
+r.support map="$OUTPUT" history="-i tmp_in > tmp_out"
+r.support map="$OUTPUT" history=""
+
+### ALL DONE ###
+
+proc_cleanup
+
+g.message "Done. <$GIS_OPT_OUTPUT> created."
+



More information about the grass-commit mailing list