Merge branch 'master' into master

AUTHORS.md

@@ -22,3 +22,4 @@ It is not a complete list, in the sense that many others have influenced the dev
 * Timothy Willey
 * Thoralf Dietrich
 * Torsten Dahm
+* Wasja Bloch

CHANGELOG.md

@@ -5,6 +5,9 @@ All notable changes to Pyrocko will be documented in this file.
 The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/).
 
 ## [Unreleased]
+*empty*
+
+## [2021.04.02]
 
 ### Added
 - RectangularSource: added opening_fraction to model tensile dislocations
@@ -18,6 +21,8 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/).
 - Jackseis: Adding `--output-steim` option to control compression. Default
   compression changed to STEIM2.
 - YAML files can now include other YAML files, when loaded through guts.
+- Moment tensor objects can now also be initialized from east-north-up
+  coordinates.
 
 ### Fixed
 - Fix plotting issues in cake.
@@ -36,6 +41,8 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/).
   problems when included in XeLaTeX).
 - QuakeML: Not strictly requiring preferred origin to be set anymore when
   extracting Pyrocko event objects.
+- Snuffler now asks for confirmation when the user attempts to close the
+  window.
 
 ## [2020.10.26]
 

examples/beachball_example03.py

@@ -1,6 +1,7 @@
 from matplotlib import pyplot as plt
 from pyrocko import moment_tensor as pmt
 from pyrocko import plot
+from pyrocko.plot import beachball
 
 
 fig = plt.figure(figsize=(4., 2.))
@@ -11,7 +12,7 @@ axes.set_ylim(0., 2.)
 axes.set_axis_off()
 
 for i, beachball_type in enumerate(['full', 'deviatoric', 'dc']):
-    plot.beachball.plot_beachball_mpl(
+    beachball.plot_beachball_mpl(
             pmt.as_mt((124654616., 370943136., -6965434.0,
                        553316224., -307467264., 84703760.0)),
             axes,

examples/plot_directivity.py

@@ -38,6 +38,8 @@ resp = plot_directivity(
     engine, rect_source, store_id,
     distance=300*km, dazi=5., component='R',
     plot_mt='full', show_phases=True,
-    phase_begin='first{stored:begin}-10%',
-    phase_end='last{stored:end}+20',
+    phases={
+        'First': 'first{stored:begin}-10%',
+        'Last': 'last{stored:end}+20'
+    },
     quantity='displacement', envelope=True)

maintenance/pip/build_wheels.sh

@@ -5,11 +5,13 @@ set -e
 rm -rf wheels_temp wheels
 
 mkdir wheels_temp
+mv pyproject.toml pyproject.toml.orig
+cp maintenance/pip/pyproject-build-pip.toml pyproject.toml
 for x in /opt/python/* ; do
     "$x/bin/python" -c 'import sys ; sys.exit(not (sys.version_info[:2] == (2, 7) or (3, 5, 0) <= sys.version_info < (3, 9, 0)))' || continue
     #"$x/bin/python" -c 'import sys ; sys.exit(not (sys.version_info[:2] == (3, 7)))' || continue
     "$x/bin/pip" install --upgrade pip
-    "$x/bin/pip" install --no-cache-dir -r requirements.txt
+    "$x/bin/pip" install --no-cache-dir -r maintenance/pip/requirements-build-pip.txt
     "$x/bin/pip" wheel -v . -w wheels_temp
 done
 
@@ -18,3 +20,5 @@ for wheel in wheels_temp/pyrocko-*.whl ; do
     auditwheel repair "$wheel" --plat $PLAT -w wheels
     rm "$wheel"
 done
+
+mv pyproject.toml.orig pyproject.toml

maintenance/pip/pyproject-build-pip.toml

@@ -0,0 +1,8 @@
+# Build dependencies as advocated by PEP 518.
+#
+#   https://www.python.org/dev/peps/pep-0518/
+#
+# Supersedes use of `setup_requires` in `setup.py`.
+
+[build-system]
+requires = ["setuptools", "wheel", "numpy>=1.8,<1.17"]

maintenance/pip/requirements-build-pip.txt

@@ -0,0 +1,8 @@
+numpy<1.17
+scipy
+pyyaml
+matplotlib
+progressbar2
+requests
+nose
+PyOpenGL

maintenance/vagrant/centos-7/inside.sh

@@ -22,15 +22,15 @@ rm -f "$outfile_py2"
 cd $HOME
 sudo yum -y install git make gcc mesa-libGL
 
-if [ ! -f "anaconda3.sh" ] ; then
-    curl 'https://repo.anaconda.com/archive/Anaconda3-5.2.0-Linux-x86_64.sh' -o anaconda3.sh
+if [ ! -f "miniconda3.sh" ] ; then
+    curl 'https://repo.anaconda.com/miniconda/Miniconda3-latest-Linux-x86_64.sh' -o miniconda3.sh
 fi
 
-if [ ! -d "anaconda3" ] ; then
-    sh anaconda3.sh -u -b
+if [ ! -d "miniconda3" ] ; then
+    sh miniconda3.sh -u -b
 fi
 
-export PATH="/home/vagrant/anaconda3/bin:$PATH"
+export PATH="/home/vagrant/miniconda3/bin:$PATH"
 
 conda install -y -q -c pyrocko pyrocko
 

setup.py

@@ -26,7 +26,7 @@ except ImportError:
     CustomBDistWheelCommand = None
 
 packname = 'pyrocko'
-version = '2020.10.26'
+version = '2021.04.02'
 
 
 class NotInAGitRepos(Exception):

src/gui/pile_viewer.py

@@ -3751,17 +3751,18 @@ def MakePileViewerMainClass(base):
             self.update()
 
         def myclose(self, return_tag=''):
+            self.return_tag = return_tag
+            self.window().close()
+
+        def cleanup(self):
+            self.about_to_close.emit()
             self.timer.stop()
             if self.follow_timer is not None:
                 self.follow_timer.stop()
-            self.window().close()
 
             for snuffling in list(self.snufflings):
                 self.remove_snuffling(snuffling)
 
-            self.about_to_close.emit()
-            self.return_tag = return_tag
-
         def set_error_message(self, key, value):
             if value is None:
                 if key in self.error_messages:
@@ -4044,6 +4045,9 @@ class PileViewer(qw.QFrame):
         self.viewer.about_to_close.connect(
             self.save_inputline_history)
 
+    def cleanup(self):
+        self.viewer.cleanup()
+
     def get_progressbars(self):
         return self.progressbars
 

src/gui/snuffler.py

@@ -17,7 +17,7 @@ from os.path import join as pjoin
 from optparse import OptionParser
 
 
-from pyrocko import pile
+from pyrocko import pile as pile_mod
 from pyrocko import util
 from pyrocko import model
 from pyrocko import config
@@ -74,12 +74,14 @@ def snuffle(pile=None, **kwargs):
     :param want_markers: bool, whether markers should be returned
     :param launch_hook: callback function called before snuffler window is
         shown
+    :param instant_close: bool, whether to bypass close window confirmation
+        dialog
     '''
     from .snuffler_app import SnufflerWindow, \
         setup_acquisition_sources, PollInjector
 
     if pile is None:
-        pile = pile.make_pile()
+        pile = pile_mod.make_pile()
 
     app = get_snuffler_instance()
 
@@ -298,7 +300,7 @@ def snuffler_from_commandline(args=None):
     if options.gui_toolkit_qt5:
         config.override_gui_toolkit = 'qt5'
 
-    this_pile = pile.Pile()
+    this_pile = pile_mod.Pile()
     stations = []
     for stations_fn in extend_paths(options.station_fns):
         stations.extend(model.station.load_stations(stations_fn))

src/gui/snuffler_app.py

@@ -570,10 +570,12 @@ class SnufflerWindow(qw.QMainWindow):
 
     def __init__(
             self, pile, stations=None, events=None, markers=None, ntracks=12,
-            follow=None, controls=True, opengl=False):
+            follow=None, controls=True, opengl=False, instant_close=False):
 
         qw.QMainWindow.__init__(self)
 
+        self.instant_close = instant_close
+
         self.dockwidget_to_toggler = {}
         self.dockwidgets = []
 
@@ -728,9 +730,23 @@ class SnufflerWindow(qw.QMainWindow):
     def return_tag(self):
         return self.get_view().return_tag
 
+    def confirm_close(self):
+        ret = qw.QMessageBox.question(
+            self,
+            'Snuffler',
+            'Close Snuffler window?',
+            qw.QMessageBox.Cancel | qw.QMessageBox.Ok,
+            qw.QMessageBox.Ok)
+
+        return ret == qw.QMessageBox.Ok
+
     def closeEvent(self, event):
-        event.accept()
-        self.closing = True
+        if self.instant_close or self.confirm_close():
+            self.closing = True
+            self.pile_viewer.cleanup()
+            event.accept()
+        else:
+            event.ignore()
 
     def is_closing(self):
         return self.closing

src/moment_tensor.py

@@ -798,7 +798,7 @@ class MomentTensor(Object):
 
         .. math::
 
-            M0 = \\frac{1}{\\sqrt{2}}\\sqrt{\\sum_{i,j} |M_{ij}|^2}
+            M_0 = \\frac{1}{\\sqrt{2}}\\sqrt{\\sum_{i,j} |M_{ij}|^2}
 
         The scalar moment is calculated based on the Euclidean (Frobenius) norm
         (Silver and Jordan, 1982). The scalar moment returned by this function

src/plot/beachball.py

@@ -9,8 +9,8 @@ from math import pi as PI
 import logging
 import numpy as num
 
-from matplotlib.collections import PathCollection
-from matplotlib.path import Path
+from matplotlib.collections import PatchCollection
+from matplotlib.patches import Polygon
 from matplotlib.transforms import Transform
 from matplotlib.colors import LinearSegmentedColormap
 
@@ -45,7 +45,7 @@ class FixedPointOffsetTransform(Transform):
         self.has_inverse = False
         self.trans = trans
         self.dpi_scale_trans = dpi_scale_trans
-        self.fixed_point = num.asarray(fixed_point, dtype=num.float)
+        self.fixed_point = num.asarray(fixed_point, dtype=num.float64)
 
     def transform_non_affine(self, values):
         fp = self.trans.transform(self.fixed_point)
@@ -251,7 +251,7 @@ def numpy_rtp2xyz(rtp):
     r = rtp[:, 0]
     theta = rtp[:, 1]
     phi = rtp[:, 2]
-    vecs = num.empty(rtp.shape, dtype=num.float)
+    vecs = num.empty(rtp.shape, dtype=num.float64)
     vecs[:, 0] = r*num.sin(theta)*num.cos(phi)
     vecs[:, 1] = r*num.sin(theta)*num.sin(phi)
     vecs[:, 2] = r*num.cos(theta)
@@ -260,7 +260,7 @@ def numpy_rtp2xyz(rtp):
 
 def numpy_xyz2rtp(xyz):
     x, y, z = xyz[:, 0], xyz[:, 1], xyz[:, 2]
-    vecs = num.empty(xyz.shape, dtype=num.float)
+    vecs = num.empty(xyz.shape, dtype=num.float64)
     vecs[:, 0] = num.sqrt(x**2+y**2+z**2)
     vecs[:, 1] = num.arctan2(num.sqrt(x**2+y**2), z)
     vecs[:, 2] = num.arctan2(y, x)
@@ -268,7 +268,7 @@ def numpy_xyz2rtp(xyz):
 
 
 def circle_points(aphi, sign=1.0):
-    vecs = num.empty((aphi.size, 3), dtype=num.float)
+    vecs = num.empty((aphi.size, 3), dtype=num.float64)
     vecs[:, 0] = num.cos(sign*aphi)
     vecs[:, 1] = num.sin(sign*aphi)
     vecs[:, 2] = 0.0
@@ -310,7 +310,7 @@ def eig2gx(eig, arcres=181):
                     continue
 
                 xtheta = num.arctan(num.sqrt(Y))
-                rtp = num.empty(xphi.shape+(3,), dtype=num.float)
+                rtp = num.empty(xphi.shape+(3,), dtype=num.float64)
                 rtp[:, 0] = 1.
                 if sign > 0:
                     rtp[:, 1] = xtheta
@@ -474,7 +474,7 @@ def choose_transform(axes, size_units, position, size):
         raise BeachballError(
             'invalid argument for size_units: %s' % size_units)
 
-    position = num.asarray(position, dtype=num.float)
+    position = num.asarray(position, dtype=num.float64)
 
     return transform, position, size
 
@@ -491,7 +491,7 @@ def mt2beachball(
         projection='lambert',
         view='top'):
 
-    position = num.asarray(position, dtype=num.float)
+    position = num.asarray(position, dtype=num.float64)
     size = size or 1
     mt = deco_part(mt, beachball_type, view)
 
@@ -588,28 +588,29 @@ def plot_beachball_mpl(
             verts = project(poly, projection)[:, ::-1] * size + position[NA, :]
             if alpha == 1.0:
                 data.append(
-                    (Path(verts[::decimation]), color, color, linewidth))
+                    (verts[::decimation], color, color, linewidth))
             else:
                 data.append(
-                    (Path(verts[::decimation]), color, 'none', 0.0))
+                    (verts[::decimation], color, 'none', 0.0))
 
         for poly in lines_upper:
             verts = project(poly, projection)[:, ::-1] * size + position[NA, :]
             data.append(
-                (Path(verts[::decimation]), 'none', edgecolor, linewidth))
-
-    paths, facecolors, edgecolors, linewidths = zip(*data)
-    path_collection = PathCollection(
-        paths,
-        facecolors=facecolors,
-        edgecolors=edgecolors,
-        linewidths=linewidths,
-        alpha=alpha,
-        zorder=zorder,
-        transform=transform)
+                (verts[::decimation], 'none', edgecolor, linewidth))
+
+    patches = []
+    for (path, facecolor, edgecolor, linewidth) in data:
+        patches.append(Polygon(
+            xy=path, facecolor=facecolor,
+            edgecolor=edgecolor,
+            linewidth=linewidth,
+            alpha=alpha))
+
+    collection = PatchCollection(
+            patches, zorder=zorder, transform=transform, match_original=True)
 
-    axes.add_artist(path_collection)
-    return path_collection
+    axes.add_artist(collection)
+    return collection
 
 
 def mts2amps(mts, projection, beachball_type, grid_resolution=200, mask=True,
@@ -622,12 +623,12 @@ def mts2amps(mts, projection, beachball_type, grid_resolution=200, mask=True,
     x = num.linspace(-1., 1., nx)
     y = num.linspace(-1., 1., ny)
 
-    vecs2 = num.zeros((nx * ny, 2), dtype=num.float)
+    vecs2 = num.zeros((nx * ny, 2), dtype=num.float64)
     vecs2[:, 0] = num.tile(x, ny)
     vecs2[:, 1] = num.repeat(y, nx)
 
     ii_ok = vecs2[:, 0]**2 + vecs2[:, 1]**2 <= 1.0
-    amps = num_full(nx * ny, num.nan, dtype=num.float)
+    amps = num_full(nx * ny, num.nan, dtype=num.float64)
 
     amps[ii_ok] = 0.
     for mt in mts:

src/plot/directivity.py

@@ -172,8 +172,8 @@ def hillshade_seismogram_array(
 def plot_directivity(
         engine, source, store_id,
         distance=300*km, azi_begin=0., azi_end=360., dazi=1.,
-        phase_begin='first{stored:any_P}-10%',
-        phase_end='last{stored:any_S}+50',
+        phases={'P': 'first{stored:any_P}-10%',
+                'S': 'last{stored:any_S}+50'},
         quantity='displacement', envelope=False,
         component='R', fmin=0.01, fmax=0.1,
         hillshade=True, cmap=None,
@@ -255,17 +255,15 @@ def plot_directivity(
     targets, azimuths = get_azimuthal_targets(
         store_id, source, distance, azi_begin, azi_end, dazi,
         components='R', quantity=quantity)
+    ref_target = targets[0]
     store = engine.get_store(store_id)
-    mt = source.pyrocko_moment_tensor(store=store, target=targets[0])
+    mt = source.pyrocko_moment_tensor(store=store, target=ref_target)
 
     resp = engine.process(source, targets, nthreads=nthreads)
     data, times = get_seismogram_array(
         resp, fmin, fmax,
         component=component, envelope=envelope)
 
-    timing_begin = Timing(phase_begin)
-    timing_end = Timing(phase_end)
-
     nucl_depth = source.depth
     nucl_distance = distance
 
@@ -287,8 +285,11 @@ def plot_directivity(
         nucl_distance -= anch_x * source.length/2.
         nucl_depth -= anch_y*num.sin(source.dip*d2r) * source.width/2.
 
-    tbegin = store.t(timing_begin, (nucl_depth, nucl_distance))
-    tend = store.t(timing_end, (nucl_depth, nucl_distance))
+    timings = [Timing(p) for p in phases.values()]
+    phase_times = [store.t(t, source, ref_target) for t in timings]
+
+    tbegin = min(phase_times)
+    tend = max(phase_times)
     tsel = num.logical_and(times >= tbegin, times <= tend)
 
     data = data[:, tsel].T
@@ -303,7 +304,7 @@ def plot_directivity(
     if envelope:
         cmw.set_clim(0., vmax)
 
-    ax.set_theta_zero_location("N")
+    ax.set_theta_zero_location('N')
     ax.set_theta_direction(-1)
 
     strike_label = mt.strike1
@@ -369,37 +370,32 @@ def plot_directivity(
         mesh.set_facecolor(color)
 
     if show_phases:
-        _phase_begin = Timing(phase_begin)
-        _phase_end = Timing(phase_end)
-
-        for p in (_phase_begin, _phase_end):
-            p.offset = 0.
-            p.offset_is_slowness = False
-            p.offset_is_percent = False
-
-        tphase_first = store.t(_phase_begin, (nucl_depth, nucl_distance))
-        tphase_last = store.t(_phase_end, (nucl_depth, nucl_distance))
-
+        label_theta = 270.
         theta = num.linspace(0, 2*num.pi, 360)
-        tfirst = num_full_like(theta, tphase_first)
-        tlast = num_full_like(theta, tphase_last)
 
-        ax.plot(theta, tfirst, color='k', alpha=.3, lw=1.)
-        ax.plot(theta, tlast, color='k', alpha=.3, lw=1.)
+        for label, phase_str in phases.items():
+            phase = Timing(phase_str)
 
-        ax.text(
-            num.pi*7/5, tphase_first, '|'.join(_phase_begin.phase_defs),
-            ha='left', color='k', fontsize='small')
+            phase.offset = 0.
+            phase.offset_is_slowness = False
+            phase.offset_is_percent = False
 
-        ax.text(
-            num.pi*6/5, tphase_last, '|'.join(_phase_end.phase_defs),
-            ha='left', color='k', fontsize='small')
+            time = store.t(phase, source, ref_target)
+            times = num_full_like(theta, time)
 
-    description = ('Component {component:s}\n'
-                   'Distance {distance:g} km').format(
-                   component=component, distance=distance / km)
+            ax.plot(theta, times, color='k', alpha=.3, lw=1., ls='--')
+
+            ax.text(
+                label_theta*d2r, time, label,
+                ha='left', color='k', fontsize='small')
+            label_theta += 30.
 
     if show_description:
+        description = (
+            'Component {component:s}\n'
+            'Distance {distance:g} km').format(
+            component=component, distance=distance / km)
+
         if fmin and fmax:
             description += '\nBandpass {fmin:g} - {fmax:g} Hz'.format(
                 fmin=fmin, fmax=fmax)

test/base/test_obspy_compat.py

@@ -23,20 +23,20 @@ class ObsPyCompatTestCase(unittest.TestCase):
         fn = common.test_data_file('test1.mseed')
 
         stream = obspy.read(fn)
-        stream.snuffle(launch_hook=close_win)
+        stream.snuffle(launch_hook=close_win, instant_close=True)
 
         trace = stream[0]
-        trace.snuffle(launch_hook=close_win)
+        trace.snuffle(launch_hook=close_win, instant_close=True)
 
     @common.require_gui
     def test_obspy_fiddle(self):
         fn = common.test_data_file('test1.mseed')
 
         stream = obspy.read(fn)
-        stream2 = stream.fiddle(launch_hook=close_win)  # noqa
+        stream2 = stream.fiddle(launch_hook=close_win, instant_close=True)  # noqa
 
         trace = stream[0]
-        trace2 = trace.fiddle(launch_hook=close_win)  # noqa
+        trace2 = trace.fiddle(launch_hook=close_win, instant_close=True)  # noqa
 
     def test_to_obspy_trace(self):
         traces = io.load(common.test_data_file('test1.mseed'))

test/gui/test_gui.py

@@ -67,7 +67,7 @@ class GUITest(unittest.TestCase):
 
         fpath = common.test_data_file('test2.mseed')
         p = make_pile(fpath, show_progress=False)
-        cls.win = SnufflerWindow(pile=p)
+        cls.win = SnufflerWindow(pile=p, instant_close=True)
         cls.pile_viewer = cls.win.pile_viewer
         cls.viewer = cls.win.pile_viewer.viewer
         pv = cls.pile_viewer