uafgeotools / mtuq

moment tensor uncertainty quantification
BSD 2-Clause "Simplified" License
68 stars 22 forks source link

Error in apply_statics #216

Closed Rinku-274 closed 1 year ago

Rinku-274 commented 1 year ago

Hi, I am trying to apply statics in MTUQ by modifying 'apply_statics=True' in mtuq/process_data.py as well as the weight file. The code is running well when I apply 'download_greens_tensors' (downloads Green's tensors from syngine). But the problem lies with 'get_greens_tensors' (extracts green's tensors from FK database). I am attaching the snapshot of the error with this message. Hope you will resolve the issue.

Thanks Rinku static_eror

Rinku-274 commented 1 year ago

Kindly look at this issue.

SeismoFelix commented 1 year ago

Hi @Rinku-274 sorry for the late reply. It seems an issue with the length of the data arrays created before performing the operation that is not working. Before the operation that is failing, you have this:

ZSS = self.select(channel="ZSS")[0].data ZDS = self.select(channel="ZDS")[0].data ZDD = self.select(channel="ZDD")[0].data ZEP = self.select(channel="ZEP")[0].data

Therefore, it seems that at least one of the vertical component Green Functions (GFs) is inconsistent with the others. This means that at least one GF is shorter than the others or may have a smaller sampling rate. Would you mind checking that the GFs you are streaming are consistent? This is, they have the same length, sampling rate, etc. Can you try with 2 stations and see what happens (with one station MTUQ is going to fail, this is an issue we are pending to fix). Good luck!

Rinku-274 commented 1 year ago

Hi @SeismoFelix I checked out the consistency of the Green's functions that I am using. All the GFs are in same length and sampling rate. My problem lies with only with process_bw in GridSearch.DoubleCouple.py process_bw = ProcessData( filter_type='Bandpass', freq_min= 0.5, freq_max= 1.0, pick_type='user_supplied', taup_model=model, window_type='body_wave', window_length=12., capuaf_file=path_weights, apply_statics=True ) Same code does not show any error in process_sw (i.e. if I apply statics for surface wave). I think apply_statics has some issue with process_bw.

Thanks

SeismoFelix commented 1 year ago

I see ...

Still, if the error you are having after calling ProcessData for the body waves is the same that the one in the image you are sharing, the issue has to do with the GFs' array size at the time of performing the operation: array[_i, _j+0, :] = -ZSS/2. * np.cos(2*phi) - ZDD/6. + ZEP/3. and probably the others as well.

You said that for processing the surface waves you do not have any issue. So my guess is when it comes to running ProcessData in the case of the body waves, there must be a point when the code trims (cuts) the GFs according to the parameters you passed. So, for any reason, the length of the arrays for all the elementary sources (ZSS, ZDS, etc.) probably are not the same although your raw GFs are. This probably has to do with the length of the body wave window you defined along with manual picking. The next test I would kindly suggest is to make a security copy of FK.py and add a print(len()) for all the ZSS, ZDS, etc. arrays just before performing the operation that is failing . If you run the test for the surface wave processing, you must see that all the lengths are the same, but not so for the body waves one. Once we can confirm that the bottom problem is the length of the arrays then we can move forward and check how the trim is working in the case of the body wave pre-processing with the parameters you are defining.

Good luck!

Félix

thurinj commented 1 year ago

Hi @Rinku-274 , @SeismoFelix ,

The issue is that the code doesn't support body-wave time shifts, as it was implemented to follow CAP implementation and support its weight files format. Here's the part of the code that parses the static time shift values from the weight files: https://github.com/uafgeotools/mtuq/blob/260b827fe8f4934986efd6172b9cc45eecc34318/mtuq/util/cap.py#L110

The alternative for body waves is to use the user-supplied pick option instead of relying on the computed Tau-p from the model if it is too far off what you are observing in the data.

My two cents on the user-supplied option: It is too cumbersome to fill in all the picked values manually, especially if you have only a couple of stations that are off. I've implemented a "mixed" method that computes the arrival time with Tau-p unless a user-provided arrival time, different from 0, is written in the CAP file. I think it would make sense to put it in the main code repo as it would be an obvious QoL improvement.

Rinku-274 commented 1 year ago

Thanks @thurinj @SeismoFelix for your suggestions and kind replies. @thurinj, as you suggested I printed out the length of all ZSS, ZDS etc and all of them are showing a length of 800. @SeismoFelix, I am using only user-supplied pick option. I think the problem lies with FK green function database (while using 'get_greens_tensors'). Since the same is not showing any error in 'download_greens_tensors' option.

Thanks.

rmodrak commented 1 year ago

If the above remains unresolved, feel free to create a new issue with a more specific description of the problem or a more specific feature request