SeismicSource / sourcespec

Earthquake source parameters from P- and S-wave displacement spectra
https://sourcespec.seismicsource.org
Other
49 stars 12 forks source link

Problem in _ssp_wave_arrival._wave_arrival #10

Closed krisvanneste closed 2 years ago

krisvanneste commented 2 years ago

If config.vp_tt or config.vs_tt are left at the default value of None, calculation of travel time and takeoff angle using a constant velocity fails in the _ssp_wave_arrival._wave_arrival function. The previous behaviour was to use config.v[p/s]_source if config.v[p/s] is None. This can be reimplemented by replacing: vel = {'P': config.vp_tt, 'S': config.vs_tt} with vel = {'P': config.vp_tt or config.vp_source, 'S': config.vs_tt or config.vs_source} Or perhaps slightly better: vel = {'P': config.vp_tt or (config.vp_source + config.vp_station)/2., 'S': config.vs_tt or (config.vs_source + config.vs_station)/2.}

Alternatively, it may even be possible to eliminate config.v[p/s]_tt entirely, because in principle we know the travel time from the origin time and the phase picks, and we know the hypocentral distance. So, we could obtain these parameters as follows: vp_tt = trace.stats.hypo_dist / (trace.stats.arrivals['P'][1] - trace.stats.hypo.origin_time) vs_tt = trace.stats.hypo_dist / (trace.stats.arrivals['S'][1] - trace.stats.hypo.origin_time) But unfortunately, this does not work, because the arrival times are not yet in trace.stats.arrivals at the time _wave_arrival is called...

claudiodsf commented 2 years ago

I'm a bit confused the issue you're raising: if vp_tt and/or vs_tt are None, then travel time is computed using the "iasp91" velocity model.

If you want vp_tt to be the average of source and station velocity, you should explicitly provide this averaged value in the config file. We might improve this by stating that if vp_tt and/or vs_tt are negative, then they will be computed as the average between source and station velocity. What do you think?

Concerning the idea of computing average velocity from picks, it cannot work because:

  1. Phase picks are not mandatory and the code computes theoretical arrivals when necessary
  2. Phase picks are always double-checked against theoretical picks and rejected if they differ more than p_arrival_tolerance or s_arrival_tolerance
krisvanneste commented 2 years ago

Claudio,

You are right, the current behaviour is in agreement with the documentation. I'm not sure it would be useful to have a special case for negative values, I can easily compute the average of source and station velocity myself and store it in config.

It's a pity that the solution based on phase picks cannot work. If the provided phase picks are reliable (I use rather large tolerances to be sure my phase picks are not rejected), then you can derive the exact vp_tt and vs_tt (which may be different for each station). I need to have a better look at your code to see if this could be done without side effects.

claudiodsf commented 2 years ago

It's a pity that the solution based on phase picks cannot work. If the provided phase picks are reliable (I use rather large tolerances to be sure my phase picks are not rejected), then you can derive the exact vp_tt and vs_tt (which may be different for each station).

But what for? vp_tt and vs_tt are only used to compute theoretical picks. If phase picks are already provided (and if you use a large tolerance), then in practice you will never need to compute theoretical picks.

If you are concerned about computing travel time (for quality factor), then we could effectively compute it from valid picks, when possible, instead of using theoretical travel time. This is not difficult to implement. But I don't think that it will make a huge difference. Anyway, I can do it 😉.

krisvanneste commented 2 years ago

Claudio,

Yes, it would only be for the benefit of a more correct quality factor. But maybe the difference will not be that big.

claudiodsf commented 2 years ago

Ok. I'll do it and report back here 😉

claudiodsf commented 2 years ago

Fixed by ce77e4b.

krisvanneste commented 2 years ago

Thank you, I can confirm that it works!