Multitaper Misfit

Warning

Please refer to [Laske1996] for a more rigorous mathematical derivation of this misfit function.

The misfit \(\chi_P(\mathbf{m})\) measures frequency-dependent phase differences estimated with multitaper methods. For a given Earth model \(\mathbf{m}\) and a single receiver, \(\chi_P(\mathbf{m})\) is given by

\[\chi_P (\mathbf{m}) = \frac{1}{2} \int_0^W W_P(w) \left| \frac{ \tau^{\mathbf{d}}(w) - \tau^{\mathbf{s}}(w, \mathbf{m})} {\sigma_P(w)} \right|^ 2 dw\]

\(\tau^\mathbf{d}(w)\) is the frequency-dependent phase measurement of the observed data; \(\tau^\mathbf{s}(w, \mathbf{m})\) the frequency-dependent phase measurement of the synthetic data.

The function \(W_P(w)\) denotes frequency-domain taper corresponding to the frequency range over which the measurements are assumed reliable. \(\sigma_P(w)\) is associated with the traveltime uncertainty introduced in making measurements, which can be estimated with cross-correlation method, or Jackknife multitaper approach.

The adjoint source for the same receiver is given by

\[f_P^{\dagger}(t) = \sum_k h_k(t)P_j(t),\]

in which \(h_k(t)\) is one (the \(k\) th) of multi-tapers.

\[\begin{split}P_j(t) = 2\pi W_p(t) * \Delta \tau(t) * p_j(t) \\ P_j(w) = 2\pi W_p(w) \Delta \tau(w) * p_j(w) \\ p_j(w) = \frac{iw s_j}{\sum_k(iw s_k)(iw s_k)^*} \\ \Delta \tau(w) = \tau^{\mathbf{d}}(w) - \tau^{\mathbf{s}}(w, \mathbf{m})\end{split}\]

Usage

adjsrc_type = "multitaper"

The following code snippet illustrates the basic usage of the cross correlation traveltime misfit function. See the corresponding Config object for additional configuration parameters.

import pyadjoint

obs, syn = pyadjoint.get_example_data()
obs = obs.select(component="Z")[0]
syn = syn.select(component="Z")[0]

config = pyadjoint.get_config(adjsrc_type="multitaper", min_period=20.,
                              max_period=100.)

adj_src = pyadjoint.calculate_adjoint_source(config=config,
                                             observed=obs, synthetic=syn,
                                             windows=[(800., 900.)]
                                             )