1# https://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
6from __future__ import absolute_import, print_function, division
8import math
9import numpy as num
11from pyrocko.guts import Bool, Float, StringChoice
12from pyrocko import cake, automap, plot
13from pyrocko.dataset import topo
14from pyrocko.gui.qt_compat import qw, qc
16from pyrocko.gui.vtk_util import TrimeshPipe, faces_to_cells, \
17 cpt_to_vtk_lookuptable
19from .base import Element, ElementState
20from .. import state as vstate
21from pyrocko import geometry
23from .. import common
25guts_prefix = 'sparrow'
28def ticks(vmin, vmax, vstep):
29 vmin = num.floor(vmin / vstep) * vstep
30 vmax = num.ceil(vmax / vstep) * vstep
31 n = int(round((vmax - vmin) / vstep))
32 return vmin + num.arange(n+1) * vstep
35def nice_value_circle(step):
36 step = plot.nice_value(step)
37 if step > 30.:
38 return 30.
40 return step
43class TopoMeshPipe(TrimeshPipe):
45 def __init__(
46 self, tile, cells_cache=None,
47 mask_ocean=False, mask_land=False, **kwargs):
49 vertices, faces = geometry.topo_to_mesh(
50 tile.y(), tile.x(), tile.data,
51 cake.earthradius)
53 self._exaggeration = 1.0
54 self._tile = tile
55 self._raw_vertices = vertices
57 centers = geometry.face_centers(vertices, faces)
59 altitudes = (geometry.vnorm(centers) - 1.0) * cake.earthradius
60 if mask_ocean:
61 mask = num.all(tile.data.flatten()[faces] == 0, axis=1)
62 altitudes[mask] = None
64 if mask_land:
65 mask = num.all(tile.data.flatten()[faces] >= 0, axis=1)
66 altitudes[mask] = None
68 if cells_cache is not None:
69 if id(faces) not in cells_cache:
70 cells_cache[id(faces)] = faces_to_cells(faces)
72 cells = cells_cache[id(faces)]
73 else:
74 cells = faces_to_cells(faces)
76 TrimeshPipe.__init__(
77 self, self._raw_vertices, cells=cells, values=altitudes, **kwargs)
79 def set_exaggeration(self, exaggeration):
80 if self._exaggeration != exaggeration:
81 factors = \
82 (cake.earthradius + self._tile.data.flatten()*exaggeration) \
83 / (cake.earthradius + self._tile.data.flatten())
84 self.set_vertices(self._raw_vertices * factors[:, num.newaxis])
85 self._exaggeration = exaggeration
88class TopoCPTChoice(StringChoice):
89 choices = ['light', 'uniform']
92class TopoState(ElementState):
93 visible = Bool.T(default=True)
94 exaggeration = Float.T(default=1.0)
95 opacity = Float.T(default=1.0)
96 smooth = Bool.T(default=False)
97 cpt = TopoCPTChoice.T(default='light')
98 shading = vstate.ShadingChoice.T(default='phong')
99 resolution_max_factor = Float.T(default=1.0)
100 resolution_min_factor = Float.T(default=1.0)
101 coverage_factor = Float.T(default=1.0)
103 def create(self):
104 element = TopoElement()
105 return element
108class TopoElement(Element):
110 def __init__(self):
111 Element.__init__(self)
112 self._parent = None
113 self.mesh = None
114 self._controls = None
116 # region = (-180., 180, -90, 90)
117 # self._tile = topo.get('ETOPO1_D8', region)
119 self._visible = False
120 self._active_meshes = {}
121 self._meshes = {}
122 self._cells = {}
123 self._cpt_name = None
124 self._lookuptables = {}
126 def get_name(self):
127 return 'Topography'
129 def bind_state(self, state):
130 Element.bind_state(self, state)
131 for var in ['visible', 'exaggeration', 'opacity', 'smooth', 'shading',
132 'cpt', 'resolution_min_factor', 'resolution_max_factor',
133 'coverage_factor']:
134 self.register_state_listener3(self.update, state, var)
136 self._state = state
138 def unbind_state(self):
139 self._listeners.clear()
141 def set_parent(self, parent):
142 self._parent = parent
144 for var in ['distance', 'lat', 'lon']:
145 self.register_state_listener3(
146 self.update, self._parent.state, var)
148 self._parent.add_panel(
149 self.get_name(),
150 self._get_controls(),
151 visible=True,
152 title_controls=[
153 self.get_title_control_remove(),
154 self.get_title_control_visible()])
156 self.update()
158 def unset_parent(self):
159 self.unbind_state()
160 if self._parent:
161 for k in self._active_meshes:
162 self._parent.remove_actor(self._active_meshes[k].actor)
164 self._active_meshes.clear()
165 self._meshes.clear()
166 self._cells.clear()
168 if self._controls:
169 self._parent.remove_panel(self._controls)
170 self._controls = None
172 self._parent.update_view()
173 self._parent = None
175 def select_dems(self, delta, region):
176 if not self._parent:
177 return [], []
179 _, size_y = self._parent.renwin.GetSize()
181 dmin = 2.0 * delta * 1.0 / float(size_y) \
182 / self._state.resolution_max_factor # [deg]
183 dmax = 2.0 * delta * 20.0 \
184 * self._state.resolution_min_factor / float(size_y)
186 result = [
187 topo.select_dem_names(
188 k, dmin, dmax, topo.positive_region(region), mode='highest')
190 for k in ['ocean', 'land']]
192 if not any(result) and delta > 20.:
193 return ['ETOPO1_D8'], ['ETOPO1_D8']
194 else:
195 return result
197 def update_cpt(self, cpt_name):
198 if cpt_name not in self._lookuptables:
199 if cpt_name == 'light':
200 topo_cpt_wet = 'light_sea'
201 topo_cpt_dry = 'light_land'
203 elif cpt_name == 'uniform':
204 topo_cpt_wet = 'light_sea_uniform'
205 topo_cpt_dry = 'light_land_uniform'
207 cpt_wet = automap.read_cpt(topo.cpt(topo_cpt_wet))
208 cpt_dry = automap.read_cpt(topo.cpt(topo_cpt_dry))
209 cpt_combi = automap.cpt_merge_wet_dry(cpt_wet, cpt_dry)
211 lut_combi = cpt_to_vtk_lookuptable(cpt_combi)
212 lut_combi.SetNanColor(0.0, 0.0, 0.0, 0.0)
214 self._lookuptables[cpt_name] = lut_combi
216 def update(self, *args):
218 pstate = self._parent.state
219 delta = pstate.distance * geometry.r2d * 0.5
220 visible = self._state.visible
222 self.update_cpt(self._state.cpt)
224 step = nice_value_circle(
225 max(1./8., min(2**round(math.log(delta) / math.log(2.)), 10.)))
227 lat_min, lat_max, lon_min, lon_max, lon_closed = common.cover_region(
228 pstate.lat, pstate.lon, delta*self._state.coverage_factor, step)
230 if lon_closed:
231 lon_max = 180.
233 region = lon_min, lon_max, lat_min, lat_max
235 dems_ocean, dems_land = self.select_dems(delta, region)
237 lat_majors = ticks(lat_min, lat_max-step, step)
238 lon_majors = ticks(lon_min, lon_max-step, step)
240 wanted = set()
241 if visible:
242 for ilat, lat in enumerate(lat_majors):
243 for ilon, lon in enumerate(lon_majors):
244 lon = ((lon + 180.) % 360.) - 180.
246 region = topo.positive_region(
247 (lon, lon+step, lat, lat+step))
249 for demname in dems_land[:1] + dems_ocean[:1]:
250 mask_ocean = demname.startswith('SRTM') \
251 or demname.startswith('Iceland')
253 mask_land = demname.startswith('ETOPO') \
254 and (dems_land and dems_land[0] != demname)
256 k = (step, demname, region, mask_ocean, mask_land)
257 if k not in self._meshes:
258 tile = topo.get(demname, region)
259 if not tile:
260 continue
262 self._meshes[k] = TopoMeshPipe(
263 tile,
264 cells_cache=self._cells,
265 mask_ocean=mask_ocean,
266 mask_land=mask_land,
267 smooth=self._state.smooth,
268 lut=self._lookuptables[self._state.cpt])
270 wanted.add(k)
272 # prevent adding both, SRTM and ETOPO because
273 # vtk produces artifacts when showing the two masked
274 # meshes (like this, mask_land is never used):
275 break
277 unwanted = set()
278 for k in self._active_meshes:
279 if k not in wanted or not visible:
280 unwanted.add(k)
282 for k in unwanted:
283 self._parent.remove_actor(self._active_meshes[k].actor)
284 del self._active_meshes[k]
286 for k in wanted:
287 if k not in self._active_meshes:
288 m = self._meshes[k]
289 self._active_meshes[k] = m
290 self._parent.add_actor(m.actor)
292 self._active_meshes[k].set_exaggeration(self._state.exaggeration)
293 self._active_meshes[k].set_opacity(self._state.opacity)
294 self._active_meshes[k].set_smooth(self._state.smooth)
295 self._active_meshes[k].set_shading(self._state.shading)
296 self._active_meshes[k].set_lookuptable(
297 self._lookuptables[self._state.cpt])
299 self._parent.update_view()
301 def _get_controls(self):
302 state = self._state
303 if not self._controls:
304 from ..state import state_bind_slider, state_bind_checkbox, \
305 state_bind_combobox
307 frame = qw.QFrame()
308 layout = qw.QGridLayout()
309 frame.setLayout(layout)
311 iy = 0
313 # exaggeration
315 layout.addWidget(qw.QLabel('Exaggeration'), iy, 0)
317 slider = qw.QSlider(qc.Qt.Horizontal)
318 slider.setSizePolicy(
319 qw.QSizePolicy(
320 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
321 slider.setMinimum(0)
322 slider.setMaximum(2000)
323 layout.addWidget(slider, iy, 1)
325 state_bind_slider(self, state, 'exaggeration', slider, factor=0.01)
327 iy += 1
329 # opacity
331 layout.addWidget(qw.QLabel('Opacity'), iy, 0)
333 slider = qw.QSlider(qc.Qt.Horizontal)
334 slider.setSizePolicy(
335 qw.QSizePolicy(
336 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
337 slider.setMinimum(0)
338 slider.setMaximum(1000)
339 layout.addWidget(slider, iy, 1)
341 state_bind_slider(self, state, 'opacity', slider, factor=0.001)
343 iy += 1
345 # high resolution
347 layout.addWidget(qw.QLabel('High-Res Factor'), iy, 0)
349 slider = qw.QSlider(qc.Qt.Horizontal)
350 slider.setSizePolicy(
351 qw.QSizePolicy(
352 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
353 slider.setMinimum(500)
354 slider.setMaximum(4000)
355 layout.addWidget(slider, iy, 1)
357 state_bind_slider(
358 self, state, 'resolution_max_factor', slider, factor=0.001)
360 iy += 1
362 # low resolution
364 layout.addWidget(qw.QLabel('Low-Res Factor'), iy, 0)
366 slider = qw.QSlider(qc.Qt.Horizontal)
367 slider.setSizePolicy(
368 qw.QSizePolicy(
369 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
370 slider.setMinimum(500)
371 slider.setMaximum(4000)
372 layout.addWidget(slider, iy, 1)
374 state_bind_slider(
375 self, state, 'resolution_min_factor', slider, factor=0.001)
377 iy += 1
379 # low resolution
381 layout.addWidget(qw.QLabel('Coverage Factor'), iy, 0)
383 slider = qw.QSlider(qc.Qt.Horizontal)
384 slider.setSizePolicy(
385 qw.QSizePolicy(
386 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
387 slider.setMinimum(500)
388 slider.setMaximum(4000)
389 layout.addWidget(slider, iy, 1)
391 state_bind_slider(
392 self, state, 'coverage_factor', slider, factor=0.001)
394 iy += 1
396 cb = common.string_choices_to_combobox(TopoCPTChoice)
397 layout.addWidget(qw.QLabel('CPT'), iy, 0)
398 layout.addWidget(cb, iy, 1)
399 state_bind_combobox(self, state, 'cpt', cb)
401 iy += 1
403 cb = qw.QCheckBox('Smooth')
404 layout.addWidget(cb, iy, 0)
405 state_bind_checkbox(self, state, 'smooth', cb)
407 cb = common.string_choices_to_combobox(vstate.ShadingChoice)
408 layout.addWidget(cb, iy, 1)
409 state_bind_combobox(self, state, 'shading', cb)
411 iy += 1
413 layout.addWidget(qw.QFrame(), iy, 0, 1, 2)
415 self._controls = frame
417 return self._controls
420__all__ = [
421 'TopoElement',
422 'TopoState',
423]