Python PPSD.load_npz方法代码示例

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
    def test_ppsd_w_iris_against_obspy_results(self):
        """
        Test against results obtained after merging of #1108.
        """
        # Read in ANMO data for one day
        st = read(os.path.join(self.path, "IUANMO.seed"))

        # Read in metadata in various different formats
        paz = {
            "gain": 86298.5,
            "zeros": [0, 0],
            "poles": [-59.4313, -22.7121 + 27.1065j, -22.7121 + 27.1065j, -0.0048004, -0.073199],
            "sensitivity": 3.3554 * 10 ** 9,
        }
        resp = os.path.join(self.path, "IUANMO.resp")
        parser = Parser(os.path.join(self.path, "IUANMO.dataless"))
        inv = read_inventory(os.path.join(self.path, "IUANMO.xml"))

        # load expected results, for both only PAZ and full response
        filename_paz = os.path.join(self.path, "IUANMO_ppsd_paz.npz")
        results_paz = PPSD.load_npz(filename_paz, metadata=None)
        filename_full = os.path.join(self.path, "IUANMO_ppsd_fullresponse.npz")
        results_full = PPSD.load_npz(filename_full, metadata=None)

        # Calculate the PPSDs and test against expected results
        # first: only PAZ
        ppsd = PPSD(st[0].stats, paz)
        ppsd.add(st)
        # commented code to generate the test data:
        # ## np.savez(filename_paz,
        # ##          **dict([(k, getattr(ppsd, k))
        # ##                  for k in PPSD.NPZ_STORE_KEYS]))
        for key in PPSD.NPZ_STORE_KEYS_ARRAY_TYPES:
            np.testing.assert_allclose(getattr(ppsd, key), getattr(results_paz, key), rtol=1e-5)
        for key in PPSD.NPZ_STORE_KEYS_LIST_TYPES:
            for got, expected in zip(getattr(ppsd, key), getattr(results_paz, key)):
                np.testing.assert_allclose(got, expected, rtol=1e-5)
        for key in PPSD.NPZ_STORE_KEYS_SIMPLE_TYPES:
            if key in ["obspy_version", "numpy_version", "matplotlib_version"]:
                continue
            self.assertEqual(getattr(ppsd, key), getattr(results_paz, key))
        # second: various methods for full response
        for metadata in [parser, inv, resp]:
            ppsd = PPSD(st[0].stats, metadata)
            ppsd.add(st)
            # commented code to generate the test data:
            # ## np.savez(filename_full,
            # ##          **dict([(k, getattr(ppsd, k))
            # ##                  for k in PPSD.NPZ_STORE_KEYS]))
            for key in PPSD.NPZ_STORE_KEYS_ARRAY_TYPES:
                np.testing.assert_allclose(getattr(ppsd, key), getattr(results_full, key), rtol=1e-5)
            for key in PPSD.NPZ_STORE_KEYS_LIST_TYPES:
                for got, expected in zip(getattr(ppsd, key), getattr(results_full, key)):
                    np.testing.assert_allclose(got, expected, rtol=1e-5)
            for key in PPSD.NPZ_STORE_KEYS_SIMPLE_TYPES:
                if key in ["obspy_version", "numpy_version", "matplotlib_version"]:
                    continue
                self.assertEqual(getattr(ppsd, key), getattr(results_full, key))

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
    def test_ppsd_temporal_plot(self):
        """
        Test plot of several period bins over time
        """
        ppsd = PPSD.load_npz(self.example_ppsd_npz)

        restrictions = {'starttime': UTCDateTime(2011, 2, 6, 1, 1),
                        'endtime': UTCDateTime(2011, 2, 7, 21, 12),
                        'year': [2011],
                        'time_of_weekday': [(-1, 2, 23)]}

        # add some gaps in the middle
        for i in sorted(list(range(30, 40)) + list(range(8, 18)) + [4])[::-1]:
            ppsd._times_processed.pop(i)
            ppsd._binned_psds.pop(i)

        with ImageComparison(self.path_images, 'ppsd_temporal.png',
                             reltol=1.5) as ic:
            fig = ppsd.plot_temporal([0.1, 1, 10], filename=None, show=False,
                                     **restrictions)
            fig.savefig(ic.name)
        with ImageComparison(self.path_images, 'ppsd_temporal.png',
                             reltol=1.5) as ic:
            ppsd.plot_temporal([0.1, 1, 10], filename=ic.name, show=False,
                               **restrictions)

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
    def test_ppsd(self):
        """
        Test PPSD routine with some real data.
        """
        # paths of the expected result data
        file_histogram = os.path.join(
            self.path,
            'BW.KW1._.EHZ.D.2011.090_downsampled__ppsd_hist_stack.npy')
        file_binning = os.path.join(
            self.path, 'BW.KW1._.EHZ.D.2011.090_downsampled__ppsd_mixed.npz')
        file_mode_mean = os.path.join(
            self.path,
            'BW.KW1._.EHZ.D.2011.090_downsampled__ppsd_mode_mean.npz')
        tr, _paz = _get_sample_data()
        st = Stream([tr])
        ppsd = _get_ppsd()
        # read results and compare
        result_hist = np.load(file_histogram)
        self.assertEqual(len(ppsd.times_processed), 4)
        self.assertEqual(ppsd.nfft, 65536)
        self.assertEqual(ppsd.nlap, 49152)
        np.testing.assert_array_equal(ppsd.current_histogram, result_hist)
        # add the same data a second time (which should do nothing at all) and
        # test again - but it will raise UserWarnings, which we omit for now
        with warnings.catch_warnings(record=True):
            warnings.simplefilter('ignore', UserWarning)
            ppsd.add(st)
            np.testing.assert_array_equal(ppsd.current_histogram, result_hist)
        # test the binning arrays
        binning = np.load(file_binning)
        np.testing.assert_array_equal(ppsd.db_bin_edges, binning['spec_bins'])
        np.testing.assert_array_equal(ppsd.period_bin_centers,
                                      binning['period_bins'])

        # test the mode/mean getter functions
        per_mode, mode = ppsd.get_mode()
        per_mean, mean = ppsd.get_mean()
        result_mode_mean = np.load(file_mode_mean)
        np.testing.assert_array_equal(per_mode, result_mode_mean['per_mode'])
        np.testing.assert_array_equal(mode, result_mode_mean['mode'])
        np.testing.assert_array_equal(per_mean, result_mode_mean['per_mean'])
        np.testing.assert_array_equal(mean, result_mode_mean['mean'])

        # test saving and loading of the PPSD (using a temporary file)
        with NamedTemporaryFile(suffix=".npz") as tf:
            filename = tf.name
            # test saving and loading to npz
            ppsd.save_npz(filename)
            ppsd_loaded = PPSD.load_npz(filename)
            ppsd_loaded.calculate_histogram()
            self.assertEqual(len(ppsd_loaded.times_processed), 4)
            self.assertEqual(ppsd_loaded.nfft, 65536)
            self.assertEqual(ppsd_loaded.nlap, 49152)
            np.testing.assert_array_equal(ppsd_loaded.current_histogram,
                                          result_hist)
            np.testing.assert_array_equal(ppsd_loaded.db_bin_edges,
                                          binning['spec_bins'])
            np.testing.assert_array_equal(ppsd_loaded.period_bin_centers,
                                          binning['period_bins'])

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
    def test_ppsd_add_npz(self):
        """
        Test PPSD.add_npz().
        """
        # set up a bogus PPSD, with fixed random psds but with real start times
        # of psd pieces, to facilitate testing the stack selection.
        ppsd = PPSD(stats=Stats(dict(sampling_rate=150)), metadata=None,
                    db_bins=(-200, -50, 20.), period_step_octaves=1.4)
        _times_processed = np.load(
            os.path.join(self.path, "ppsd_times_processed.npy")).tolist()
        # change data to nowadays used nanoseconds POSIX timestamp
        _times_processed = [UTCDateTime(t)._ns for t in _times_processed]
        np.random.seed(1234)
        _binned_psds = [
            arr for arr in np.random.uniform(
                -200, -50,
                (len(_times_processed), len(ppsd.period_bin_centers)))]

        with NamedTemporaryFile(suffix=".npz") as tf1, \
                NamedTemporaryFile(suffix=".npz") as tf2, \
                NamedTemporaryFile(suffix=".npz") as tf3:
            # save data split up over three separate temporary files
            ppsd._times_processed = _times_processed[:200]
            ppsd._binned_psds = _binned_psds[:200]
            ppsd.save_npz(tf1.name)
            ppsd._times_processed = _times_processed[200:400]
            ppsd._binned_psds = _binned_psds[200:400]
            ppsd.save_npz(tf2.name)
            ppsd._times_processed = _times_processed[400:]
            ppsd._binned_psds = _binned_psds[400:]
            ppsd.matplotlib_version = "X.X.X"
            ppsd.save_npz(tf3.name)
            # now load these saved npz files and check if all data is present
            ppsd = PPSD.load_npz(tf1.name, metadata=None)
            ppsd.add_npz(tf2.name)
            # we changed a version number so this should emit a warning
            with warnings.catch_warnings(record=True) as w:
                warnings.simplefilter('always')
                ppsd.add_npz(tf3.name)
                self.assertEqual(len(w), 1)
            np.testing.assert_array_equal(_binned_psds, ppsd._binned_psds)
            np.testing.assert_array_equal(_times_processed,
                                          ppsd._times_processed)
            # adding data already present should also emit a warning and the
            # PPSD should not be changed
            with warnings.catch_warnings(record=True) as w:
                warnings.simplefilter('always')
                ppsd.add_npz(tf2.name)
                self.assertEqual(len(w), 1)
            np.testing.assert_array_equal(_binned_psds, ppsd._binned_psds)
            np.testing.assert_array_equal(_times_processed,
                                          ppsd._times_processed)

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
 def test_exception_reading_newer_npz(self):
     """
     Checks that an exception is properly raised when trying to read a npz
     that was written on a more recent ObsPy version (specifically that has
     a higher 'ppsd_version' number which is used to keep track of changes
     in PPSD and the npz file used for serialization).
     """
     msg = ("Trying to read/add a PPSD npz with 'ppsd_version=100'. This "
            "file was written on a more recent ObsPy version that very "
            "likely has incompatible changes in PPSD internal structure "
            "and npz serialization. It can not safely be read with this "
            "ObsPy version (current 'ppsd_version' is {!s}). Please "
            "consider updating your ObsPy installation.".format(
                PPSD(stats=Stats(), metadata=None).ppsd_version))
     # 1 - loading a npz
     data = np.load(self.example_ppsd_npz)
     # we have to load, modify 'ppsd_version' and save the npz file for the
     # test..
     items = {key: data[key] for key in data.files}
     # deliberately set a higher ppsd_version number
     items['ppsd_version'] = items['ppsd_version'].copy()
     items['ppsd_version'].fill(100)
     with NamedTemporaryFile() as tf:
         filename = tf.name
         with open(filename, 'wb') as fh:
             np.savez(fh, **items)
         with self.assertRaises(ObsPyException) as e:
             PPSD.load_npz(filename)
     self.assertEqual(str(e.exception), msg)
     # 2 - adding a npz
     ppsd = PPSD.load_npz(self.example_ppsd_npz)
     for method in (ppsd.add_npz, ppsd._add_npz):
         with NamedTemporaryFile() as tf:
             filename = tf.name
             with open(filename, 'wb') as fh:
                 np.savez(fh, **items)
             with self.assertRaises(ObsPyException) as e:
                 method(filename)
             self.assertEqual(str(e.exception), msg)

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
    def test_PPSD_save_and_load_npz(self):
        """
        Test PPSD.load_npz() and PPSD.save_npz()
        """
        _, paz = _get_sample_data()
        ppsd = _get_ppsd()

        # save results to npz file
        with NamedTemporaryFile(suffix=".npz") as tf:
            filename = tf.name
            # test saving and loading an uncompressed file
            ppsd.save_npz(filename)
            ppsd_loaded = PPSD.load_npz(filename, metadata=paz)

        for key in NPZ_STORE_KEYS:
            np.testing.assert_equal(getattr(ppsd, key),
                                    getattr(ppsd_loaded, key))

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
    def test_ppsd_spectrogram_plot(self):
        """
        Test spectrogram type plot of PPSD
        """
        ppsd = PPSD.load_npz(self.example_ppsd_npz)

        # add some gaps in the middle
        for i in sorted(list(range(30, 40)) + list(range(8, 18)) + [4])[::-1]:
            ppsd._times_processed.pop(i)
            ppsd._binned_psds.pop(i)

        with ImageComparison(self.path_images, 'ppsd_spectrogram.png',
                             reltol=1.5) as ic:
            fig = ppsd.plot_spectrogram(filename=None, show=False)
            fig.savefig(ic.name)
        with ImageComparison(self.path_images, 'ppsd_spectrogram.png',
                             reltol=1.5) as ic:
            ppsd.plot_spectrogram(filename=ic.name, show=False)

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
    def test_ppsd_save_and_load_npz(self):
        """
        Test PPSD.load_npz() and PPSD.save_npz()
        """
        _, paz = _get_sample_data()
        ppsd = _get_ppsd()

        # save results to npz file
        with NamedTemporaryFile(suffix=".npz") as tf:
            filename = tf.name
            # test saving and loading an uncompressed file
            ppsd.save_npz(filename)
            ppsd_loaded = PPSD.load_npz(filename, metadata=paz)

        for key in PPSD.NPZ_STORE_KEYS:
            if isinstance(getattr(ppsd, key), np.ndarray) or key == "_binned_psds":
                np.testing.assert_equal(getattr(ppsd, key), getattr(ppsd_loaded, key))
            else:
                self.assertEqual(getattr(ppsd, key), getattr(ppsd_loaded, key))

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
    def test_ppsd_spectrogram_plot(self):
        """
        Test spectrogram type plot of PPSD

        Matplotlib version 3 shifts the x-axis labels but everything else looks
        the same. Skipping test for matplotlib >= 3 on 05/12/2018.
        """
        ppsd = PPSD.load_npz(self.example_ppsd_npz)

        # add some gaps in the middle
        for i in sorted(list(range(30, 40)) + list(range(8, 18)) + [4])[::-1]:
            ppsd._times_processed.pop(i)
            ppsd._binned_psds.pop(i)

        with ImageComparison(self.path_images, 'ppsd_spectrogram.png',
                             reltol=1.5) as ic:
            fig = ppsd.plot_spectrogram(filename=None, show=False)
            fig.savefig(ic.name)
        with ImageComparison(self.path_images, 'ppsd_spectrogram.png',
                             reltol=1.5) as ic:
            ppsd.plot_spectrogram(filename=ic.name, show=False)

# 需要导入模块: from obspy.signal.spectral_estimation import PPSD [as 别名]
# 或者: from obspy.signal.spectral_estimation.PPSD import load_npz [as 别名]
    def test_ppsd_w_iris_against_obspy_results(self):
        """
        Test against results obtained after merging of #1108.
        """
        # Read in ANMO data for one day
        st = read(os.path.join(self.path, 'IUANMO.seed'))

        # Read in metadata in various different formats
        paz = {'gain': 86298.5, 'zeros': [0, 0],
               'poles': [-59.4313, -22.7121 + 27.1065j, -22.7121 + 27.1065j,
                         -0.0048004, -0.073199], 'sensitivity': 3.3554*10**9}
        resp = os.path.join(self.path, 'IUANMO.resp')
        parser = Parser(os.path.join(self.path, 'IUANMO.dataless'))
        inv = read_inventory(os.path.join(self.path, 'IUANMO.xml'))

        # load expected results, for both only PAZ and full response
        filename_paz = os.path.join(self.path, 'IUANMO_ppsd_paz.npz')
        results_paz = PPSD.load_npz(filename_paz, metadata=None)
        filename_full = os.path.join(self.path,
                                     'IUANMO_ppsd_fullresponse.npz')
        results_full = PPSD.load_npz(filename_full, metadata=None)

        # Calculate the PPSDs and test against expected results
        # first: only PAZ
        ppsd = PPSD(st[0].stats, paz)
        ppsd.add(st)
        # commented code to generate the test data:
        # ## np.savez(filename_paz,
        # ##          **dict([(k, getattr(ppsd, k))
        # ##                  for k in PPSD.NPZ_STORE_KEYS]))
        for key in PPSD.NPZ_STORE_KEYS_ARRAY_TYPES:
            np.testing.assert_allclose(
                getattr(ppsd, key), getattr(results_paz, key), rtol=1e-5)
        for key in PPSD.NPZ_STORE_KEYS_LIST_TYPES:
            for got, expected in zip(getattr(ppsd, key),
                                     getattr(results_paz, key)):
                np.testing.assert_allclose(got, expected, rtol=1e-5)
        for key in PPSD.NPZ_STORE_KEYS_SIMPLE_TYPES:
            if key in ["obspy_version", "numpy_version", "matplotlib_version"]:
                continue
            self.assertEqual(getattr(ppsd, key), getattr(results_paz, key))
        # second: various methods for full response
        # (also test various means of initialization, basically testing the
        #  decorator that maps the deprecated keywords)
        for metadata in [parser, inv, resp]:
            with warnings.catch_warnings(record=True) as w:
                warnings.simplefilter('always')
                ppsd = PPSD(st[0].stats, paz=metadata)
            self.assertEqual(len(w), 1)
            self.assertIs(w[0].category, ObsPyDeprecationWarning)

            with warnings.catch_warnings(record=True) as w:
                warnings.simplefilter('always')
                ppsd = PPSD(st[0].stats, parser=metadata)
            self.assertEqual(len(w), 1)
            self.assertIs(w[0].category, ObsPyDeprecationWarning)

            ppsd = PPSD(st[0].stats, metadata)
            ppsd.add(st)
            # commented code to generate the test data:
            # ## np.savez(filename_full,
            # ##          **dict([(k, getattr(ppsd, k))
            # ##                  for k in PPSD.NPZ_STORE_KEYS]))
            for key in PPSD.NPZ_STORE_KEYS_ARRAY_TYPES:
                np.testing.assert_allclose(
                    getattr(ppsd, key), getattr(results_full, key), rtol=1e-5)
            for key in PPSD.NPZ_STORE_KEYS_LIST_TYPES:
                for got, expected in zip(getattr(ppsd, key),
                                         getattr(results_full, key)):
                    np.testing.assert_allclose(got, expected, rtol=1e-5)
            for key in PPSD.NPZ_STORE_KEYS_SIMPLE_TYPES:
                if key in ["obspy_version", "numpy_version",
                           "matplotlib_version"]:
                    continue
                self.assertEqual(getattr(ppsd, key),
                                 getattr(results_full, key))