PeterMakus / SeisMIC

SeisMIC is a python software suite to monitor velocity changes using ambient seismic noise.
https://petermakus.github.io/SeisMIC/
European Union Public License 1.2
49 stars 12 forks source link

A problem regarding the calculation of EGF in the high-frequency range #66

Closed Yuuki-yuuki closed 4 months ago

Yuuki-yuuki commented 4 months ago

Problem

Hello, Peter: I am very interested in your work and appreciate the program you have released. But I have noticed that in certain frequency bands, especially in the high-frequency range, the calculated cross-correlation function shows some issues at Lag Time = -MAX, as shown in the figure below. I am wondering if this is related to the use of the frequency domain method for calculating cross-correlation in the "correlate.py" file's "_pxcorr_matrix" function. Therefore, I am seeking your advice to understand whether this issue might affect the dvv results for high-frequency signals.

1 2

Proposed solution

No response

PeterMakus commented 4 months ago

Hi! Oh that's interesting and certainly not physical. Actually I am wondering whether this is a weird edge effect of the bandpass filter you applied.

On a different note: In which frequency band did you filter the seismic data before correlating? Because if you did not change the parameters in the tutorial then the data were filtered between 2 and 4 Hz, so I would not really expect any content between 6 and 8 Hz and that might be the cause for this strange signal.

PeterMakus commented 4 months ago

PS The filters I am referring to are set in the params.yaml file at this point:

` corr_args : {'TDpreProcessing':[ {'function':'seismic.correlate.preprocessing_td.taper', 'args':{'type':'cosine_taper', 'p': 0.02}}, {'function':'seismic.correlate.preprocessing_td.detrend', 'args':{'type':'linear'}}, {'function':'seismic.correlate.preprocessing_td.TDfilter', 'args':{'type':'bandpass','freqmin':2,'freqmax':4}},

{'function':'seismic.correlate.preprocessing_td.mute',

                                # 'args':{'taper_len':100.,
                                       # 'threshold':1000, absolute threshold
                                #         'std_factor':3,
                                #         'filter':{'type':'bandpass','freqmin':2,'freqmax':4},
                                #         'extend_gaps':True}},
                              #  {'function':'seismic.correlate.preprocessing_td.clip',
                              #   'args':{'std_factor':2.5}},
                               {'function':'seismic.correlate.preprocessing_td.signBitNormalization',
                                'args': {}}
                               ],
              # Standard functions reside in seismic.correlate.preprocessing_fd
             'FDpreProcessing':[
                                # {'function':'seismic.correlate.preprocessing_fd.spectralWhitening',
                                #  'args':{'joint_norm':False}},
                                {'function':'seismic.correlate.preprocessing_fd.FDfilter',
                                 'args':{'flimit':[1.33, 2, 4, 6]}}
                                #  {'function':seismic.correlate.preprocessing_fd.FDsignBitNormalization,
                                # 'args':{}}
                                ],`
Yuuki-yuuki commented 4 months ago

Hi! Oh that's interesting and certainly not physical. Actually I am wondering whether this is a weird edge effect of the bandpass filter you applied.

On a different note: In which frequency band did you filter the seismic data before correlating? Because if you did not change the parameters in the tutorial then the data were filtered between 2 and 4 Hz, so I would not really expect any content between 6 and 8 Hz and that might be the cause for this strange signal.

Thank you for replying, Peter. In fact, besides the example data in the tutorial, I used three datasets from different sites (including broadband seismometers and short-period seismometers), all with a sampling rate of 100Hz. During preprocessing, I modified the contents of 'params.yaml' to:

sampling_rate : 50
preProcessing : [
                {'function':'seismic.correlate.preprocessing_stream.stream_filter',
                  'args':{'ftype':'bandpass',
                          'filter_option':{'freqmin':0.01, #0.01
                                           'freqmax':24.99}}}

corr_args : {'TDpreProcessing':[
                              {'function':'seismic.correlate.preprocessing_td.taper',
                                'args':{'type':'cosine_taper', 'p': 0.02}},
                                {'function':'seismic.correlate.preprocessing_td.detrend',
                                'args':{'type':'linear'}},
                               {'function':'seismic.correlate.preprocessing_td.TDfilter',
                               'args':{'type':'bandpass','freqmin':1,'freqmax':24.99}},
                                #{'function':'seismic.correlate.preprocessing_td.mute',
                                # 'args':{'taper_len':100.,
                                       # 'threshold':1000, absolute threshold
                                #         'std_factor':3,
                                #         'filter':{'type':'bandpass','freqmin':1,'freqmax':12},
                                #         'extend_gaps':True}},
                              #  {'function':'seismic.correlate.preprocessing_td.clip',
                              #   'args':{'std_factor':2.5}},
                               {'function':'seismic.correlate.preprocessing_td.signBitNormalization',
                                'args': {}}
                               ],

However, in all three cases, the phenomena we mentioned earlier still exist: YORYY N}7VB1X()%)ZW6X

Fortunately, it almost does not affect the calculation results of DVV, but only affects the visualization effect of EGF.

Wait, I noticed that I haven't modified the spectral whitening parameters. Could this be the issue? My FDpreProcessing parameters are still set to the default:

              # Standard functions reside in seismic.correlate.preprocessing_fd
             'FDpreProcessing':[
                                # {'function':'seismic.correlate.preprocessing_fd.spectralWhitening',
                                #  'args':{'joint_norm':False}},
                                {'function':'seismic.correlate.preprocessing_fd.FDfilter',
                                 'args':{'flimit':[1.33, 2, 4, 8]}}
                                #  {'function':seismic.correlate.preprocessing_fd.FDsignBitNormalization,
                                # 'args':{}}
                                ],
PeterMakus commented 4 months ago

Hey! So this is the issue here: The seismic data is filtered 3 times. Once, the daily streams are filtered: with ´stream_filter´ Then, the matrix is filtered in the time domain. That is the option you set in ´TDPreProcessing´. Lastly, it is filtered in the frequency domain (´FDpreProcessing´). And here you set a band that will not let any energy through. So basically this part ´´ {'function':'seismic.correlate.preprocessing_fd.FDfilter', 'args':{'flimit':[1.33, 2, 4, 8]}} ´´ either needs to be commented or you change the filter parameters accordingly :). After doing that, I am pretty sure that the edge effect should disappear.

The reason why I filter three times is to avoid any numerical side effects from the preprocessing and the Fourier Transforms.

Just let me know if this works.

CHeers, Peter

Yuuki-yuuki commented 4 months ago

Yes, you are absolutely right. After adjusting the flimit in 'FDpreProcessing', the edge effect completely disappeared. This indicates that the problem was not due to an error in the algorithm. Thank you for your patient response in helping me resolve this issue. Finally, I would like to ask if there are reasonable guidelines for setting the value of flimit. From what I understand, spectral whitening requires two parameters, the lower limit and the upper limit, corresponding to f2 and f3 in FDfilter. Therefore, I'm not sure how the values of f1 and f4 will affect the results. Could you give me some advice?

LG06UFLM7MK04%WLQ`D3T4L

By the way, I asked about this issue in ChatGPT, and the suggestions provided might be useful as a reference: try to choose f1 and f4 outside the frequency range of your signal of interest to ensure that the filter effectively preserves the main components of the signal. Specifically, f1 should be low enough to ensure that important low-frequency information is not filtered out; f4 should be high enough to ensure that important high-frequency information is not filtered out.

Does this mean that we can simply set f1 to a very low frequency and f4 to the Nyquist frequency?

PeterMakus commented 4 months ago

Hey, Great to hear that it works now!

For clarification, the frequency domain filtering that you changed has nothing do with the spectral whitening. The spectral whitening is only 1. applying the Fourier transform 2. deciding the signal in the frequency domain by its amplitude spectrum, so it's not band limited.

For the Frequency domain filter: As this is not a boxcar shaped filter (which would introduce high frequency artefacts), you can define the slope of the filter using f2 and f3 (the frequencies where the filter slope starts). f1 and f4 are the frequencies where the filter reaches zero.

I will close this issue now as it seems to be resolved :)