[GRASS-dev] Fwd: [GRASS GIS] #3246: r.in.srtm: patch to support SRTM water bodies import (SWBD)
Luca Delucchi
lucadeluge at gmail.com
Fri Feb 10 08:29:33 PST 2017
I'm not able to upload the new patch because there is an error [0].
I'm attaching the new patch here
[0] https://trac.osgeo.org/osgeo/ticket/1875
---------- Forwarded message ----------
From: GRASS GIS <trac at osgeo.org>
Date: 10 February 2017 at 17:17
Subject: Re: [GRASS GIS] #3246: r.in.srtm: patch to support SRTM water
bodies import (SWBD)
To:
#3246: r.in.srtm: patch to support SRTM water bodies import (SWBD)
--------------------------+---------------------------------
Reporter: stjo | Owner: grass-dev@…
Type: enhancement | Status: new
Priority: normal | Milestone: 7.2.1
Component: Python | Version: svn-releasebranch72
Resolution: | Keywords: r.in.srtm
CPU: Unspecified | Platform: All
--------------------------+---------------------------------
Comment (by lucadelu):
I updated the patch, but I'm not able to test it because there are no raw
files in https://dds.cr.usgs.gov/srtm/version2_1/SWBD/. In the zip files
there are only shapefile.
--
Ticket URL: <https://trac.osgeo.org/grass/ticket/3246#comment:1>
GRASS GIS <https://grass.osgeo.org>
--
ciao
Luca
www.lucadelu.org
-------------- next part --------------
Index: r.in.srtm.html
===================================================================
--- r.in.srtm.html (revision 70529)
+++ r.in.srtm.html (working copy)
@@ -1,41 +1,73 @@
<h2>DESCRIPTION</h2>
-<em>r.in.srtm</em> imports SRTM hgt files into GRASS.
+<em>r.in.srtm</em> imports SRTM hgt files or SRTM SWBD (SRTM Water Body Data) raw
+files into GRASS GIS.
-SRTM Version 1 and improved Version 2 data sets can be downloaded from
-NASA at this site:<br>
-<a href="http://dds.cr.usgs.gov/srtm/">http://dds.cr.usgs.gov/srtm/</a>
+SRTM Version 1 and improved Version 2 data sets can be downloaded from
+NASA at <a href="http://dds.cr.usgs.gov/srtm/">this site</a>.<br>
<p>
-Gap-filled SRTM Version 3 data can be downloaded from USGS at this site:<br>
-<a href="http://e4ftl01.cr.usgs.gov/SRTM/SRTMGL3.003/2000.02.11/">http://e4ftl01.cr.usgs.gov/SRTM/SRTMGL3.003/2000.02.11/</a>
+Gap-filled SRTM Version 3 (1 arc seconds at 30m) data can be downloaded from USGS
+at <a href="http://e4ftl01.cr.usgs.gov/SRTM/SRTMGL3.003/2000.02.11/">this site</a>.
+<p>
+SRTM SWBD (SRTM Water Body Data; 1 arc seconds at 30m) can be downloaded from USGS at
+<a href="http://e4ftl01.cr.usgs.gov/SRTM/SRTMSWBD.003/2000.02.11/">this site</a>.
+
+<p>
+The default output is the input name, without extension. It is also possible to
+import the zip files directly.
+
<h2>NOTES</h2>
SRTM tiles are of 1 degree by 1 degree size. The SRTM filename contains the
-coordinates which refer to the <b>center</b> of the lower left pixel (e.g., N51E010:
+coordinates which refer to the <b>center</b> of the lower left pixel (e.g., N51E010:
lower left cell center at 10E, 51N). To identify a tile name, a grid can be easily
visualized in the GRASS monitor:
<div class="code"><pre>
+# in latitude-longitude location
+d.mon wx0
d.grid size=1
</pre></div>
-To import TOPEX/SRTM30 PLUS data, use <em>r.in.bin</em>.
+Hint: To import TOPEX/SRTM30 PLUS data, use <em>r.in.bin</em>.
+<h2>EXAMPLES</h2>
+
+<h3>Import of SRTM</h3>
+Import of SRTM V3 (1 arc seconds at 30m):
+<div class="code"><pre>
+r.in.srtm -1 input=N00E006.hgt output=N00E006
+</pre></div>
+
+Import of SRTM V3 zip file (1 arc seconds at 30m):
+<div class="code"><pre>
+r.in.srtm -1 input=N00E006.SRTMGL1.hgt.zip output=N00E006
+</pre></div>
+
+<h3>Import of SRTM water bodies</h3>
+
+Import of SRTM SWBD (1 arc seconds at 30m):
+<div class="code"><pre>
+# using predefined output
+r.in.srtm -1 -w input=N00E006.raw
+</pre></div>
+
<h2>SEE ALSO</h2>
<em>
-<a href="r.in.bin.html">r.in.bin</a>
+<a href="r.import.html">r.import</a>,
+<a href="r.in.bin.html">r.in.bin</a>,
+<a href="r.in.gdal.html">r.in.gdal</a>
</em>
<p>The <a href="http://www2.jpl.nasa.gov/srtm/">Shuttle Radar Topography Mission</a>
homepage at NASA's JPL.
-<br>
-The <a href="http://pub7.bravenet.com/forum/537683448/">SRTM Web Forum</a>
<h2>REFERENCES</h2>
-M. Neteler, 2005. <a href="http://grass.osgeo.org/newsletter/GRASSNews_vol3.pdf">SRTM and VMAP0 data in OGR and GRASS.</a> <i><a href="http://grass.osgeo.org/newsletter/">GRASS Newsletter</a></i>, Vol.3, pp. 2-6, June 2005. ISSN 1614-8746.
+M. Neteler, 2005. <a href="http://grass.osgeo.org/newsletter/GRASSNews_vol3.pdf">SRTM and VMAP0 data in OGR and GRASS.</a>
+<i><a href="http://grass.osgeo.org/newsletter/">GRASS Newsletter</a></i>, Vol.3, pp. 2-6, June 2005. ISSN 1614-8746.
<h2>AUTHORS</h2>
@@ -42,6 +74,7 @@
Markus Neteler<br>
Improved by W. Kyngesburye and H. Bowman<br>
-Update for SRTM V3 by Markus Metz
+Update for SRTM V3 by Markus Metz<br>
+Update for SRTM SWBD by Jonas Strobel
<p><i>Last changed: $Date$</i>
Index: r.in.srtm.py
===================================================================
--- r.in.srtm.py (revision 70529)
+++ r.in.srtm.py (working copy)
@@ -2,13 +2,14 @@
#
############################################################################
#
-# MODULE: r_in_aster.py
+# MODULE: r.in.srtm.py
# AUTHOR(S): Markus Neteler 11/2003 neteler AT itc it
# Hamish Bowman
# Glynn Clements
-# PURPOSE: import of SRTM hgt files into GRASS
+# Jonas Strobel
+# PURPOSE: import of SRTM hgt files or SRTM SWBD raw files into GRASS
#
-# COPYRIGHT: (C) 2004, 2006 by the GRASS Development Team
+# COPYRIGHT: (C) 2004-2017 by the GRASS Development Team
#
# This program is free software under the GNU General Public
# License (>=v2). Read the file COPYING that comes with GRASS
@@ -22,6 +23,7 @@
# April 2006: links updated from ftp://e0dps01u.ecs.nasa.gov/srtm/
# to current links below
# October 2008: Converted to Python by Glynn Clements
+# December 2016: Added flag to read SRTM SWBD (SRTM Water Body Dataset) raw binary files
#########################
# Derived from:
# ftp://e0srp01u.ecs.nasa.gov/srtm/version1/Documentation/Notes_for_ARCInfo_users.txt
@@ -47,14 +49,22 @@
# http://earth-info.nga.mil/GandG/wgs84/gravitymod/egm96/intpt.html
#
#########################
+# Test cases
+# - SRTM HGT 1 arc second
+# - SRTM HGT 3 arc second
+# - SRTM HGT.zip 1 arc second
+# - SRTM HGT.zip 3 arc second
+# - SRTM SWBD 1 arc second
+# - SRTM SWBD.zip 1 arc second
+#########################
#%Module
-#% description: Imports SRTM HGT files into raster map.
+#% description: Imports SRTM HGT or RAW files into raster map.
#% keyword: raster
#% keyword: import
#%End
#%option G_OPT_F_INPUT
-#% description: Name of SRTM input tile (file without .hgt.zip extension)
+#% description: Name of SRTM HGT / SRTM RAW input tile
#%end
#%option G_OPT_R_OUTPUT
#% description: Name for output raster map (default: input tile)
@@ -64,6 +74,10 @@
#% key: 1
#% description: Input is a 1-arcsec tile (default: 3-arcsec)
#%end
+#%flag
+#% key: w
+#% description: Import SRTM SWBD (SRTM Water Body Data)
+#%end
tmpl1sec = """BYTEORDER M
LAYOUT BIL
@@ -117,7 +131,7 @@
def cleanup():
if not in_temp:
return
- for ext in ['.bil', '.hdr', '.prj', '.hgt.zip']:
+ for ext in ['.bil', '.hdr', '.prj', '.hgt.zip', '.raw', '.raw.zip']:
grass.try_remove(tile + ext)
os.chdir('..')
grass.try_rmdir(tmpdir)
@@ -131,6 +145,7 @@
input = options['input']
output = options['output']
one = flags['1']
+ water = flags['w']
# are we in LatLong location?
s = grass.read_command("g.proj", flags='j')
@@ -140,17 +155,23 @@
# use these from now on:
infile = input
- while infile[-4:].lower() in ['.hgt', '.zip']:
+ while infile[-4:].lower() in ['.hgt', '.zip', '.raw']:
infile = infile[:-4]
(fdir, tile) = os.path.split(infile)
if not output:
tileout = tile
+ grass.debug("No output set... using name: " + tileout)
else:
tileout = output
- zipfile = infile + ".hgt.zip"
- hgtfile = os.path.join(fdir, tile[:7] + ".hgt")
+ if not water:
+ zipfile = infile + ".hgt.zip"
+ hgtfile = os.path.join(fdir, tile[:7] + ".hgt")
+ else:
+ zipfile = infile + ".raw.zip"
+ hgtfile = os.path.join(fdir, tile[:7] + ".raw")
+
if os.path.isfile(zipfile):
# check if we have unzip
if not grass.find_program('unzip'):
@@ -162,11 +183,16 @@
tenv['UNZIP'] = '-qq'
if grass.call(['unzip', '-t', zipfile], env=tenv) != 0:
grass.fatal(_("'%s' does not appear to be a valid zip file.") % zipfile)
+ is_zip = True
+ elif not water:
+ os.path.isfile(hgtfile)
+ # try and see if it's already unzipped
+ is_zip = False
- is_zip = True
elif os.path.isfile(hgtfile):
# try and see if it's already unzipped
is_zip = False
+
else:
grass.fatal(_("File '%s' or '%s' not found") % (zipfile, hgtfile))
@@ -176,26 +202,39 @@
os.mkdir(tmpdir)
if is_zip:
- shutil.copyfile(zipfile, os.path.join(tmpdir, tile + ".hgt.zip"))
+ if not water:
+ shutil.copyfile(zipfile, os.path.join(tmpdir, tile + ".hgt.zip"))
+ else:
+ shutil.copyfile(zipfile, os.path.join(tmpdir, tile + ".raw.zip"))
+
else:
- shutil.copyfile(hgtfile, os.path.join(tmpdir, tile + ".hgt"))
+ if not water:
+ shutil.copyfile(hgtfile, os.path.join(tmpdir, tile + ".hgt"))
+ else:
+ shutil.copyfile(hgtfile, os.path.join(tmpdir, tile + ".raw"))
# change to temporary directory
os.chdir(tmpdir)
in_temp = True
- zipfile = tile + ".hgt.zip"
- hgtfile = tile[:7] + ".hgt"
+ if not water:
+ zipfile = tile + ".hgt.zip"
+ hgtfile = tile[:7] + ".hgt"
+ else:
+ zipfile = tile + ".raw.zip"
+ hgtfile = tile[:7] + ".raw"
+
bilfile = tile + ".bil"
if is_zip:
# unzip & rename data file:
- grass.message(_("Extracting '%s'...") % infile)
+ grass.message(_("Extracting '%s'...") % zipfile)
if grass.call(['unzip', zipfile], env=tenv) != 0:
grass.fatal(_("Unable to unzip file."))
- grass.message(_("Converting input file to BIL..."))
- os.rename(hgtfile, bilfile)
+ if not water:
+ grass.message(_("Converting input file to BIL..."))
+ os.rename(hgtfile, bilfile)
north = tile[0]
ll_latitude = int(tile[1:3])
@@ -210,17 +249,27 @@
if east == "W":
ll_longitude *= -1
- # Calculate Upper Left from Lower Left
- ulxmap = "%.1f" % ll_longitude
- # SRTM90 tile size is 1 deg:
- ulymap = "%.1f" % (ll_latitude + 1)
+ if water:
+ # Calculate Upper Left from Lower Left
+ ulxmap = "%.1f" % (ll_longitude + 1)
- if not one:
- tmpl = tmpl3sec
+ # SRTM90 tile size is 1 deg:
+ ulymap = "%.1f" % (ll_latitude + 1)
+
else:
+ # Calculate Upper Left from Lower Left
+ ulxmap = "%.1f" % ll_longitude
+
+ # SRTM90 tile size is 1 deg:
+ ulymap = "%.1f" % (ll_latitude + 1)
+
+ if one or water:
grass.message(_("Attempting to import 1-arcsec data."))
tmpl = tmpl1sec
+ else:
+ tmpl = tmpl3sec
+
header = tmpl % (ulxmap, ulymap)
hdrfile = tile + '.hdr'
outf = file(hdrfile, 'w')
@@ -233,20 +282,38 @@
outf.write(proj)
outf.close()
- try:
- grass.run_command('r.in.gdal', input=bilfile, out=tileout)
- except:
- grass.fatal(_("Unable to import data"))
+ if not water:
+ try:
+ grass.run_command('r.in.gdal', input=bilfile, output=tileout)
+ except:
+ grass.fatal("Unable to import data")
+ else:
+ # If water, these operations are required
+ swbd_res = 0.000277777777777778 # 0:00:01
+
+ n = float(ulymap) + (0.5 * swbd_res)
+ s = float(ll_latitude) - (0.5 * swbd_res)
+ e = float(ulxmap) + (0.5 * swbd_res)
+ w = int(ll_longitude) - (0.5 * swbd_res)
+
+ try:
+ grass.run_command('r.in.bin', input=hgtfile, output=tileout,
+ bytes=1, north=n, south=s, east=e, west=w,
+ rows=3601, cols=3601)
+ except:
+ grass.fatal(_("Unable to import data"))
+
# nice color table
grass.run_command('r.colors', map=tileout, color='srtm')
# write cmd history:
grass.raster_history(tileout)
-
grass.message(_("Done: generated map ") + tileout)
- grass.message(_("(Note: Holes in the data can be closed with 'r.fillnulls' using splines)"))
+ if not water:
+ grass.message(_("(Note: Holes in the data can be closed with 'r.fillnulls' using splines)"))
+
if __name__ == "__main__":
options, flags = grass.parser()
atexit.register(cleanup)
More information about the grass-dev
mailing list