1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
6import zipfile
7import os.path as op
8import os
9import re
10# from pyrocko.util import urlopen
12import numpy as num
14from pyrocko import util, config
15from . import tile, dataset
17citation = '''
18Farr, T. G., and M. Kobrick, 2000, Shuttle Radar Topography Mission produces a
19wealth of data. Eos Trans. AGU, 81:583-585.
21Farr, T. G. et al., 2007, The Shuttle Radar Topography Mission, Rev. Geophys.,
2245, RG2004, doi:10.1029/2005RG000183. (Also available online at
23http://www2.jpl.nasa.gov/srtm/SRTM_paper.pdf)
25Kobrick, M., 2006, On the toes of giants-How SRTM was born, Photogramm. Eng.
26Remote Sens., 72:206-210.
28Rosen, P. A. et al., 2000, Synthetic aperture radar interferometry, Proc. IEEE,
2988:333-382.
30'''
33class AuthenticationRequired(Exception):
34 pass
37class SRTMGL3(dataset.TiledGlobalDataset):
39 def __init__(
40 self,
41 name='SRTMGL3',
42 data_dir=op.join(op.dirname(__file__), 'data', 'SRTMGL3'),
43 raw_data_url='https://mirror.pyrocko.org/e4ftl01.cr.usgs.gov/'
44 'MEASURES/SRTMGL3.003/2000.02.11'):
46 dataset.TiledGlobalDataset.__init__(
47 self,
48 name,
49 432001, 216001, 1201, 1201,
50 num.dtype('>i2'),
51 data_dir=data_dir,
52 citation=citation,
53 region=(-180., 180., -60., 60.))
55 self.raw_data_url = raw_data_url
56 self._available_tilenames = None
57 self.config = config.config()
59 def tilename(self, itx, ity):
60 itx -= 180
61 ity -= 90
62 if ity >= 0:
63 s = 'N%02i' % ity
64 else:
65 s = 'S%02i' % -ity
67 if itx >= 0:
68 s += 'E%03i' % itx
69 else:
70 s += 'W%03i' % -itx
72 return s
74 def available_tilenames(self):
75 if self._available_tilenames is None:
76 fpath = op.join(self.data_dir, 'available.list')
77 if not op.exists(fpath) or os.stat(fpath).st_size == 0:
78 util.ensuredirs(fpath)
79 # remote structure changed, we would have to clawl through
80 # many pages. Now keeping tile index here:
81 self.download_file(
82 'https://mirror.pyrocko.org/e4ftl01.cr.usgs.gov/'
83 'MEASURES/SRTMGL3.003/2000.02.11/available.list', fpath)
85 # url = self.raw_data_url + '/'
86 # f = urlopen(url)
87 # data = f.read().decode()
88 # available = re.findall(
89 # r'([NS]\d\d[EW]\d\d\d)\.SRTMGL3\.hgt', data)
90 #
91 # f.close()
92 #
93 # with open(fpath, 'w') as f:
94 # f.writelines('%s\n' % s for s in available)
96 with open(fpath, 'r') as f:
97 available = [
98 s.strip() for s in f.readlines()
99 if re.match(r'^[NS]\d\d[EW]\d\d\d$', s.strip())]
101 self._available_tilenames = set(available)
103 return self._available_tilenames
105 def tilefilename(self, tilename):
106 return tilename + '.SRTMGL3.hgt.zip'
108 def tilepath(self, tilename):
109 fn = self.tilefilename(tilename)
110 return op.join(self.data_dir, fn)
112 def download_tile(self, tilename):
113 fpath = self.tilepath(tilename)
114 fn = self.tilefilename(tilename)
115 url = self.raw_data_url + '/' + fn
116 try:
117 # we have to follow the oauth redirect here...
118 self.download_file(url, fpath)
119 except Exception as e:
120 raise e
122 def download(self):
123 for tn in self.available_tilenames():
124 fpath = self.tilepath(tn)
125 if not op.exists(fpath):
126 self.download_tile(tn)
128 def get_tile(self, itx, ity):
129 tn = self.tilename(itx, ity)
130 if tn not in self.available_tilenames():
131 return None
132 else:
133 fpath = self.tilepath(tn)
134 if not op.exists(fpath):
135 self.download_tile(tn)
137 zipf = zipfile.ZipFile(fpath, 'r')
138 rawdata = zipf.read(tn + '.hgt')
139 zipf.close()
140 data = num.frombuffer(rawdata, dtype=self.dtype)
141 assert data.size == self.ntx * self.nty
142 data = data.reshape(self.nty, self.ntx)[::-1, ::]
143 return tile.Tile(
144 self.xmin + itx*self.stx,
145 self.ymin + ity*self.sty,
146 self.dx, self.dx, data)
149if __name__ == '__main__':
150 import sys
152 util.setup_logging('pyrocko.topo.srtmgl3', 'info')
154 if len(sys.argv) != 2:
155 sys.exit('usage: python -m pyrocko.topo.srtmgl3 download')
157 if sys.argv[1] == 'download':
158 srtmgl3 = SRTMGL3()
159 srtmgl3.download()