[GRASS-SVN] r54755 - in grass-addons/grass6/raster: . r.surf.nnbathy
svn_grass at osgeo.org
svn_grass at osgeo.org
Tue Jan 22 11:02:27 PST 2013
Author: msieczka
Date: 2013-01-22 11:02:27 -0800 (Tue, 22 Jan 2013)
New Revision: 54755
Added:
grass-addons/grass6/raster/r.surf.nnbathy/
grass-addons/grass6/raster/r.surf.nnbathy/Makefile
grass-addons/grass6/raster/r.surf.nnbathy/description.html
grass-addons/grass6/raster/r.surf.nnbathy/r.surf.nnbathy
Log:
Submitting r.surf.nnbathy addon shell script.
Added: grass-addons/grass6/raster/r.surf.nnbathy/Makefile
===================================================================
--- grass-addons/grass6/raster/r.surf.nnbathy/Makefile (rev 0)
+++ grass-addons/grass6/raster/r.surf.nnbathy/Makefile 2013-01-22 19:02:27 UTC (rev 54755)
@@ -0,0 +1,7 @@
+MODULE_TOPDIR = ../..
+
+PGM = r.surf.nnbathy
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script
Property changes on: grass-addons/grass6/raster/r.surf.nnbathy/Makefile
___________________________________________________________________
Added: svn:mime-type
+ text/x-makefile
Added: svn:eol-style
+ native
Added: grass-addons/grass6/raster/r.surf.nnbathy/description.html
===================================================================
--- grass-addons/grass6/raster/r.surf.nnbathy/description.html (rev 0)
+++ grass-addons/grass6/raster/r.surf.nnbathy/description.html 2013-01-22 19:02:27 UTC (rev 54755)
@@ -0,0 +1,40 @@
+<h2>DESCRIPTION</h2>
+
+<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 interpolation library, written by Pavel Sakov.
+
+<p>
+<em>r.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 non-Sibsonian natural neighbor interpolation (<b>alg=ns</b>). For performing the underlaying Delaunay triangulation in all cases <em>nnbathy</em> uses <em>triangle</em> software by <a href="http://www.cs.berkeley.edu/~jrs/">Jonathan Richard Shewchuk</a>.
+
+<p>
+The <b>output</b> raster map is a continous surface interpolated from the <b>input</b> raster map.
+
+<p>
+<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 held in the memory.
+
+<h2>NOTES</h2>
+
+<ul>
+1. Requires <em>GRASS</em> 6.x 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>
+3. The output raster map extent and resolution match the region settings at which the script was started.<br>
+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>
+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>
+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 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>
+
+<h2>AUTHOR</h2>
+
+Maciej Sieczka
+
+<p>
+<i>Last changed: $Date$</i>
Property changes on: grass-addons/grass6/raster/r.surf.nnbathy/description.html
___________________________________________________________________
Added: svn:mime-type
+ text/html
Added: svn:keywords
+ Author Date Id
Added: svn:eol-style
+ native
Added: grass-addons/grass6/raster/r.surf.nnbathy/r.surf.nnbathy
===================================================================
--- grass-addons/grass6/raster/r.surf.nnbathy/r.surf.nnbathy (rev 0)
+++ grass-addons/grass6/raster/r.surf.nnbathy/r.surf.nnbathy 2013-01-22 19:02:27 UTC (rev 54755)
@@ -0,0 +1,265 @@
+#!/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
Property changes on: grass-addons/grass6/raster/r.surf.nnbathy/r.surf.nnbathy
___________________________________________________________________
Added: svn:executable
+ *
Added: svn:mime-type
+ text/x-sh
Added: svn:eol-style
+ native
More information about the grass-commit
mailing list