Coverage for /usr/local/lib/python3.11/dist-packages/pyrocko/dataset/topo/etopo1.py: 40%
48 statements
« prev ^ index » next coverage.py v6.5.0, created at 2023-10-06 06:59 +0000
« prev ^ index » next coverage.py v6.5.0, created at 2023-10-06 06:59 +0000
1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
6'''
7NOAA `ETOPO1 <https://www.ngdc.noaa.gov/mgg/global/>`_ Global Relief Map.
8'''
10import zipfile
11import os.path as op
13import numpy as num
15from . import tile, dataset
17try:
18 range = xrange
19except NameError:
20 pass
22citation = '''
23Amante, C. and B.W. Eakins, 2009. ETOPO1 1 Arc-Minute Global Relief Model:
24Procedures, Data Sources and Analysis. NOAA Technical Memorandum NESDIS
25NGDC-24. National Geophysical Data Center, NOAA. doi:10.7289/V5C8276M
26'''
28__doc__ += '''
30.. rubric:: Citation
31''' + citation
34class ETOPO1(dataset.TiledGlobalDataset):
35 '''
36 ETOPO1 as a tiled global dataset.
37 '''
39 def __init__(
40 self,
41 name='ETOPO1',
42 data_dir=op.join(op.dirname(__file__), 'data', 'ETOPO1'),
43 base_fn='etopo1_ice_g_i2',
44 raw_data_url=('https://mirror.pyrocko.org/www.ngdc.noaa.gov/mgg'
45 '/global/relief/ETOPO1/data/ice_surface/'
46 'grid_registered/binary/%s.zip')):
48 dataset.TiledGlobalDataset.__init__(
49 self,
50 name,
51 21601, 10801, 1351, 1351,
52 num.dtype('<i2'),
53 data_dir=data_dir,
54 citation=citation)
56 self.base_fn = base_fn
57 self.raw_data_url = raw_data_url
59 def download(self):
60 fpath = op.join(self.data_dir, '%s.zip' % self.base_fn)
61 if not op.exists(fpath):
62 self.download_file(
63 self.raw_data_url % self.base_fn, fpath)
65 self.make_tiles()
67 def make_tiles(self):
68 fpath = op.join(self.data_dir, '%s.zip' % self.base_fn)
70 zipf = zipfile.ZipFile(fpath, 'r')
71 rawdata = zipf.read('%s.bin' % self.base_fn)
72 zipf.close()
73 data = num.frombuffer(rawdata, dtype=self.dtype)
74 assert data.size == self.nx * self.ny
75 data = data.reshape((self.ny, self.nx))[::-1]
77 for ity in range(self.ntilesy):
78 for itx in range(self.ntilesx):
79 tiledata = data[ity*(self.nty-1):(ity+1)*(self.nty-1)+1,
80 itx*(self.ntx-1):(itx+1)*(self.ntx-1)+1]
82 fn = '%s.%02i.%02i.bin' % (self.base_fn, ity, itx)
83 fpath = op.join(self.data_dir, fn)
84 with open(fpath, 'w') as f:
85 tiledata.tofile(f)
87 def get_tile(self, itx, ity):
88 assert 0 <= itx < self.ntilesx
89 assert 0 <= ity < self.ntilesy
91 fn = '%s.%02i.%02i.bin' % (self.base_fn, ity, itx)
92 fpath = op.join(self.data_dir, fn)
93 if not op.exists(fpath):
94 self.download()
96 with open(fpath, 'r') as f:
97 data = num.fromfile(f, dtype=self.dtype)
99 assert data.size == self.ntx*self.nty
100 data = data.reshape((self.ntx, self.nty))
102 return tile.Tile(
103 self.xmin + itx*self.stx,
104 self.ymin + ity*self.sty,
105 self.dx, self.dx, data)