Here is my code:
https://drive.google.com/file/d/137x6BmLilntRax-4w8Depw9653xDaKbh/view?usp=sharing
Here is my file:
https://drive.google.com/file/d/10bp_I_Vz6ngSbqZ8WryC_USGVXUQrDMU/view?usp=sharing
My ppsd.plot_temporal doesn’t work.
Please help me.
Thanks.
Getting an error downloading that notebook
Here is my code.
Please help me check what problem in my code, thank you.
from obspy import read, Stream
from obspy.signal import PPSD
paz = {‘gain’: 1.,
‘poles’: [-21.991+22.435j, -21.991-22.435j],
‘sensitivity’: 257362500,
‘zeros’: [0j, 0j]}
ppsd = None
st = read(‘CL05.2020.199.0717.EHZ.sac’)
ppsd = PPSD(st[0].stats, metadata=paz)
ppsd.add(st)
ppsd.plot()
ppsd.plot_spectrogram()
ppsd.plot_temporal([0.1, 1, 10])
I don’t see any obvious problems there. Can you try ppsd.add(st, verbose=True) and what is the return value of that command, it should be True if data was added and processed properly.
Maybe the PSD values are just outside of the plotted range there? You can have a look at the underlying PSD data like this:
from obspy.signal.tests.test_spectral_estimation import _get_ppsd
ppsd = _get_ppsd()
print(ppsd.psd_values)
Hmm they should be in the plot above, judging from the range. Can you upload the data somewhere else? Above links don’t work for me
No problem, this is my data file. Thanks.
https://drive.google.com/file/d/19aHwc8v65qiVOwYsklzZUeyjH3bcecaQ/view?usp=drive_link
This seems to be a bug, looks like the time step is not properly incremented and then every second time slice isn’t processed, and when plotting the temporal plot it recognises these gaps and the plot probably ends up with single dots plotted that can’t be seen. I’m looking at it..
Somehow the sample times in the stream experience some kind of subsample creep, which eventually makes the PPSD processing think that every second slice of the trace should not get processed because it overlaps with the segment before.
It is pretty easy to fix in PPSD, even though I still need to double check if that fix causes issues elsewhere.
On the other hand, we need to look into how this happens in the first place because that is at the very core of obspy at the Trace level..
@tinyulin0000 see Ppsd fix processing time info by megies · Pull Request #3387 · obspy/obspy · GitHub
For a temporary solution you could make a new conda environment and install obspy from that branch of the pull request. I am still looking at where this problem is actually coming from, since it seems to come from core functionality
So, after chasing this issue for quite some time, turns out, the problem comes from how SAC is storing sampling rate information, it seems. See the linked obspy issue for more details.
@tinyulin0000 what you could do as a quick short term fix is to fix the sampling rate and get rid of the single-precision inaccuracy. Obviously, you need to take care that the sampling rate is what you expect it to be..
st = read("...")
for tr in st:
assert tr.stats.sampling_rate - 25.0 < 1e-5
tr.stats.sampling_rate = 25.0
...
ppsd.add(st)
Sorry for the late reply and thank you very much !
I also thought the problem result from my data file, but I didn’t know the detail problem.
So I really appreciate your help very much!!!