1# https://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
6import math
7import numpy as num
9from pyrocko.guts import Bool, Float, StringChoice
10from pyrocko import cake, automap, plot
11from pyrocko.dataset import topo
12from pyrocko.gui.qt_compat import qw, qc
14from pyrocko.gui.vtk_util import TrimeshPipe, faces_to_cells, \
15 cpt_to_vtk_lookuptable
17from .base import Element, ElementState
18from .. import state as vstate
19from pyrocko import geometry
21from .. import common
23guts_prefix = 'sparrow'
26def ticks(vmin, vmax, vstep):
27 vmin = num.floor(vmin / vstep) * vstep
28 vmax = num.ceil(vmax / vstep) * vstep
29 n = int(round((vmax - vmin) / vstep))
30 return vmin + num.arange(n+1) * vstep
33def nice_value_circle(step):
34 step = plot.nice_value(step)
35 if step > 30.:
36 return 30.
38 return step
41class TopoMeshPipe(TrimeshPipe):
43 def __init__(
44 self, tile, cells_cache=None,
45 mask_ocean=False, mask_land=False, **kwargs):
47 vertices, faces = geometry.topo_to_mesh(
48 tile.y(), tile.x(), tile.data,
49 cake.earthradius)
51 self._exaggeration = 1.0
52 self._tile = tile
53 self._raw_vertices = vertices
55 centers = geometry.face_centers(vertices, faces)
57 altitudes = (geometry.vnorm(centers) - 1.0) * cake.earthradius
58 if mask_ocean:
59 mask = num.all(tile.data.flatten()[faces] == 0, axis=1)
60 altitudes[mask] = None
62 if mask_land:
63 mask = num.all(tile.data.flatten()[faces] >= 0, axis=1)
64 altitudes[mask] = None
66 if cells_cache is not None:
67 if id(faces) not in cells_cache:
68 cells_cache[id(faces)] = faces_to_cells(faces)
70 cells = cells_cache[id(faces)]
71 else:
72 cells = faces_to_cells(faces)
74 TrimeshPipe.__init__(
75 self, self._raw_vertices, cells=cells, values=altitudes, **kwargs)
77 def set_exaggeration(self, exaggeration):
78 if self._exaggeration != exaggeration:
79 factors = \
80 (cake.earthradius + self._tile.data.flatten()*exaggeration) \
81 / (cake.earthradius + self._tile.data.flatten())
82 self.set_vertices(self._raw_vertices * factors[:, num.newaxis])
83 self._exaggeration = exaggeration
86class TopoCPTChoice(StringChoice):
87 choices = ['light', 'uniform']
90class TopoState(ElementState):
91 visible = Bool.T(default=True)
92 exaggeration = Float.T(default=1.0)
93 opacity = Float.T(default=1.0)
94 smooth = Bool.T(default=False)
95 cpt = TopoCPTChoice.T(default='light')
96 shading = vstate.ShadingChoice.T(default='phong')
97 resolution_max_factor = Float.T(default=1.0)
98 resolution_min_factor = Float.T(default=1.0)
99 coverage_factor = Float.T(default=1.0)
101 def create(self):
102 element = TopoElement()
103 return element
106class TopoElement(Element):
108 def __init__(self):
109 Element.__init__(self)
110 self._parent = None
111 self.mesh = None
112 self._controls = None
114 # region = (-180., 180, -90, 90)
115 # self._tile = topo.get('ETOPO1_D8', region)
117 self._visible = False
118 self._active_meshes = {}
119 self._meshes = {}
120 self._cells = {}
121 self._cpt_name = None
122 self._lookuptables = {}
124 def get_name(self):
125 return 'Topography'
127 def bind_state(self, state):
128 Element.bind_state(self, state)
129 self.talkie_connect(
130 state,
131 ['visible', 'exaggeration', 'opacity', 'smooth', 'shading', 'cpt',
132 'resolution_min_factor', 'resolution_max_factor',
133 'coverage_factor'],
134 self.update)
136 self._state = state
138 def set_parent(self, parent):
139 self._parent = parent
141 for var in ['distance', 'lat', 'lon']:
142 self.talkie_connect(
143 self._parent.state, var, self.update)
145 self._parent.add_panel(
146 self.get_title_label(),
147 self._get_controls(),
148 visible=True,
149 title_controls=[
150 self.get_title_control_remove(),
151 self.get_title_control_visible()])
153 self.update()
155 def unset_parent(self):
156 self.unbind_state()
157 if self._parent:
158 for k in self._active_meshes:
159 self._parent.remove_actor(self._active_meshes[k].actor)
161 self._active_meshes.clear()
162 self._meshes.clear()
163 self._cells.clear()
165 if self._controls:
166 self._parent.remove_panel(self._controls)
167 self._controls = None
169 self._parent.update_view()
170 self._parent = None
172 def select_dems(self, delta, region):
173 if not self._parent:
174 return [], []
176 _, size_y = self._parent.renwin.GetSize()
178 dmin = 2.0 * delta * 1.0 / float(size_y) \
179 / self._state.resolution_max_factor # [deg]
180 dmax = 2.0 * delta * 20.0 \
181 * self._state.resolution_min_factor / float(size_y)
183 result = [
184 topo.select_dem_names(
185 k, dmin, dmax, topo.positive_region(region), mode='highest')
187 for k in ['ocean', 'land']]
189 if not any(result) and delta > 20.:
190 return ['ETOPO1_D8'], ['ETOPO1_D8']
191 else:
192 return result
194 def update_cpt(self, cpt_name):
195 if cpt_name not in self._lookuptables:
196 if cpt_name == 'light':
197 topo_cpt_wet = 'light_sea'
198 topo_cpt_dry = 'light_land'
200 elif cpt_name == 'uniform':
201 topo_cpt_wet = 'light_sea_uniform'
202 topo_cpt_dry = 'light_land_uniform'
204 cpt_wet = automap.read_cpt(topo.cpt(topo_cpt_wet))
205 cpt_dry = automap.read_cpt(topo.cpt(topo_cpt_dry))
206 cpt_combi = automap.cpt_merge_wet_dry(cpt_wet, cpt_dry)
208 lut_combi = cpt_to_vtk_lookuptable(cpt_combi)
209 lut_combi.SetNanColor(0.0, 0.0, 0.0, 0.0)
211 self._lookuptables[cpt_name] = lut_combi
213 def update(self, *args):
215 pstate = self._parent.state
216 delta = pstate.distance * geometry.r2d * 0.5
217 visible = self._state.visible
219 self.update_cpt(self._state.cpt)
221 step = nice_value_circle(
222 max(1./8., min(2**round(math.log(delta) / math.log(2.)), 10.)))
224 lat_min, lat_max, lon_min, lon_max, lon_closed = common.cover_region(
225 pstate.lat, pstate.lon, delta*self._state.coverage_factor, step)
227 if lon_closed:
228 lon_max = 180.
230 region = lon_min, lon_max, lat_min, lat_max
232 dems_ocean, dems_land = self.select_dems(delta, region)
234 lat_majors = ticks(lat_min, lat_max-step, step)
235 lon_majors = ticks(lon_min, lon_max-step, step)
237 wanted = set()
238 if visible:
239 for ilat, lat in enumerate(lat_majors):
240 for ilon, lon in enumerate(lon_majors):
241 lon = ((lon + 180.) % 360.) - 180.
243 region = topo.positive_region(
244 (lon, lon+step, lat, lat+step))
246 for demname in dems_land[:1] + dems_ocean[:1]:
247 mask_ocean = demname.startswith('SRTM') \
248 or demname.startswith('Iceland')
250 mask_land = demname.startswith('ETOPO') \
251 and (dems_land and dems_land[0] != demname)
253 k = (step, demname, region, mask_ocean, mask_land)
254 if k not in self._meshes:
255 tile = topo.get(demname, region)
256 if not tile:
257 continue
259 self._meshes[k] = TopoMeshPipe(
260 tile,
261 cells_cache=self._cells,
262 mask_ocean=mask_ocean,
263 mask_land=mask_land,
264 smooth=self._state.smooth,
265 lut=self._lookuptables[self._state.cpt])
267 wanted.add(k)
269 # prevent adding both, SRTM and ETOPO because
270 # vtk produces artifacts when showing the two masked
271 # meshes (like this, mask_land is never used):
272 break
274 unwanted = set()
275 for k in self._active_meshes:
276 if k not in wanted or not visible:
277 unwanted.add(k)
279 for k in unwanted:
280 self._parent.remove_actor(self._active_meshes[k].actor)
281 del self._active_meshes[k]
283 for k in wanted:
284 if k not in self._active_meshes:
285 m = self._meshes[k]
286 self._active_meshes[k] = m
287 self._parent.add_actor(m.actor)
289 self._active_meshes[k].set_exaggeration(self._state.exaggeration)
290 self._active_meshes[k].set_opacity(self._state.opacity)
291 self._active_meshes[k].set_smooth(self._state.smooth)
292 self._active_meshes[k].set_shading(self._state.shading)
293 self._active_meshes[k].set_lookuptable(
294 self._lookuptables[self._state.cpt])
296 self._parent.update_view()
298 def _get_controls(self):
299 state = self._state
300 if not self._controls:
301 from ..state import state_bind_slider, state_bind_checkbox, \
302 state_bind_combobox
304 frame = qw.QFrame()
305 layout = qw.QGridLayout()
306 frame.setLayout(layout)
308 iy = 0
310 # exaggeration
312 layout.addWidget(qw.QLabel('Exaggeration'), iy, 0)
314 slider = qw.QSlider(qc.Qt.Horizontal)
315 slider.setSizePolicy(
316 qw.QSizePolicy(
317 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
318 slider.setMinimum(0)
319 slider.setMaximum(2000)
320 layout.addWidget(slider, iy, 1)
322 state_bind_slider(self, state, 'exaggeration', slider, factor=0.01)
324 iy += 1
326 # opacity
328 layout.addWidget(qw.QLabel('Opacity'), iy, 0)
330 slider = qw.QSlider(qc.Qt.Horizontal)
331 slider.setSizePolicy(
332 qw.QSizePolicy(
333 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
334 slider.setMinimum(0)
335 slider.setMaximum(1000)
336 layout.addWidget(slider, iy, 1)
338 state_bind_slider(self, state, 'opacity', slider, factor=0.001)
340 iy += 1
342 # high resolution
344 layout.addWidget(qw.QLabel('High-Res Factor'), iy, 0)
346 slider = qw.QSlider(qc.Qt.Horizontal)
347 slider.setSizePolicy(
348 qw.QSizePolicy(
349 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
350 slider.setMinimum(500)
351 slider.setMaximum(4000)
352 layout.addWidget(slider, iy, 1)
354 state_bind_slider(
355 self, state, 'resolution_max_factor', slider, factor=0.001)
357 iy += 1
359 # low resolution
361 layout.addWidget(qw.QLabel('Low-Res Factor'), iy, 0)
363 slider = qw.QSlider(qc.Qt.Horizontal)
364 slider.setSizePolicy(
365 qw.QSizePolicy(
366 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
367 slider.setMinimum(500)
368 slider.setMaximum(4000)
369 layout.addWidget(slider, iy, 1)
371 state_bind_slider(
372 self, state, 'resolution_min_factor', slider, factor=0.001)
374 iy += 1
376 # low resolution
378 layout.addWidget(qw.QLabel('Coverage Factor'), iy, 0)
380 slider = qw.QSlider(qc.Qt.Horizontal)
381 slider.setSizePolicy(
382 qw.QSizePolicy(
383 qw.QSizePolicy.Expanding, qw.QSizePolicy.Fixed))
384 slider.setMinimum(500)
385 slider.setMaximum(4000)
386 layout.addWidget(slider, iy, 1)
388 state_bind_slider(
389 self, state, 'coverage_factor', slider, factor=0.001)
391 iy += 1
393 cb = common.string_choices_to_combobox(TopoCPTChoice)
394 layout.addWidget(qw.QLabel('CPT'), iy, 0)
395 layout.addWidget(cb, iy, 1)
396 state_bind_combobox(self, state, 'cpt', cb)
398 iy += 1
400 cb = qw.QCheckBox('Smooth')
401 layout.addWidget(cb, iy, 0)
402 state_bind_checkbox(self, state, 'smooth', cb)
404 cb = common.string_choices_to_combobox(vstate.ShadingChoice)
405 layout.addWidget(cb, iy, 1)
406 state_bind_combobox(self, state, 'shading', cb)
408 iy += 1
410 layout.addWidget(qw.QFrame(), iy, 0, 1, 2)
412 self._controls = frame
414 return self._controls
417__all__ = [
418 'TopoElement',
419 'TopoState',
420]