On a clean synthetic, all four receiver-function deconvolutions recover the same Ps peak to within a percent. That agreement is a trap. Drop the signal-to-noise ratio to what a real broadband station delivers on a mediocre teleseism and the methods diverge, not by a little, and not in the order you would guess. The fastest, most-cited method is the first to lie.
One impulse response, four estimators
A receiver function is a deconvolution problem. You have the vertical component, which is your best estimate of the incident P wavelet plus source and path, and you have the radial component, which is that same wavelet reverberating through the crust beneath the station. Divide one out of the other and you get the crust's impulse response: the series of Ps conversions whose delays encode the Moho and whose amplitudes encode the velocity contrasts across it. Four methods dominate practice. Spectral water-level division [1, 2], multitaper deconvolution [3], extended-time multitaper [4], and iterative time-domain deconvolution [5]. They all target the same truth. They fail differently, and the failure order is the whole story.
Water-level: fast, cited, and fragile
Water-level division works in the frequency domain. Divide the radial spectrum by the vertical spectrum, and to stop the divide-by-zero at spectral nulls, clamp the denominator so it never drops below a fixed fraction c of its peak value:
That clamp, the water level, is the method's whole personality. At high SNR it barely matters; the nulls are shallow and the clamp rarely engages. But as noise rises, the vertical spectrum develops deep, noisy nulls exactly where there is no signal, and the water level floods them with regularized garbage that inverse-transforms into ringing across the entire trace. The Ps peak is still there, but it now sits on a bed of oscillation that a stack will happily average into a coherent artifact.
The parameter c makes it worse, because c is a bias-variance dial that almost nobody tunes. It is usually inherited verbatim from a script: 0.01, or 0.001, or whatever the last paper used. Set c very small and the Ps pulse comes back sharp, but the noisy nulls divide through almost unclamped and each one smears an oscillation across the whole receiver function. Set c large and the ringing vanishes, but you have over-regularized the spectrum, broadened the pulse, and biased its amplitude low. Receiver-function amplitude is not decorative. It constrains the velocity contrast across the interface, and in H-k stacking it weights the multiples [6], so an amplitude bias propagates directly into the impedance contrast and the crustal composition inference. The sweet spot between those regimes depends on the station's noise level, the source spectra in your catalog, and the frequency band of the target conversion. A single hard-coded c is the right answer for exactly one station on exactly one event, and the wrong answer everywhere else.
Multitaper holds because it averages independent looks
Multitaper deconvolution replaces the single spectral estimate with an average over several orthogonal taper windows, the discrete prolate spheroidal sequences [3]. Each taper gives a statistically independent estimate of the spectrum, and averaging them suppresses the variance that water-level has no defense against. The payoff is a receiver function that degrades gracefully rather than ringing. Extended-time multitaper adds low-frequency stabilization [4], which matters when your target is a deep, low-frequency conversion like the lithosphere-asthenosphere boundary rather than the sharp Moho. Both cost more compute than water-level, and both are worth it once your station log carries a realistic fraction of noisy events.
Iterative deconvolution wins at the bottom of the SNR range
Iterative time-domain deconvolution builds the receiver function one spike at a time [5]. It cross-correlates the current residual with the vertical component, adds the largest-amplitude spike, subtracts its predicted contribution, and repeats until adding spikes stops improving the fit. Because it constructs the answer additively from the data rather than dividing spectra, it never manufactures ringing from a spectral null. There is no null to divide by. That makes it the most robust method as SNR collapses, which is exactly the regime where your marginal events live. The cost is compute and a stopping-criterion choice, but on a noisy catalog it is the method that keeps the Moho pick honest.
Watch the podium reorder
The instrument below makes the divergence literal. Four lanes reconstruct the same synthetic crustal impulse response: direct P, the Moho Ps conversion near 4.3 s, and the first crustal multiple. The amber dashed trace is the truth, and the bar under each lane scores the reconstruction against it. At high SNR the four bars are level and the traces are nearly indistinguishable; a clean-data benchmark genuinely cannot separate these methods. Drag the SNR down. Water-level grows broadband ringing and its bar drops off the podium first. Multitaper broadens and drifts at low frequency, extended-time multitaper broadens with less drift, and iterative stays spiky and close to the truth, paying only a small amplitude loss and, at the very bottom of the range, one spurious late spike.
Two readouts summarize the state: which method currently leads, and the fidelity spread between the best and worst lane. At high SNR the spread is near zero, which is the quantitative version of "any method will do." At low SNR the spread opens to several tenths, which is the difference between a Moho pick and an artifact.
The ranking is a property of your noise, not the method
The uncomfortable conclusion is that there is no best deconvolution in the abstract. Which method you should run is a function of your station's noise distribution, not of which paper you read last. A pipeline that hard-codes one deconvolution for a continent-scale deployment is optimizing for the average station and quietly corrupting the noisy ones, and because the corruption arrives as coherent ringing rather than obvious garbage, the stack hides it instead of flagging it.
The modern move, and the one EarthScan builds toward, is a deconvolution that conditions on each trace's own noise and picks its regularization per event. The water level was always a proxy for a question the data can answer directly: how much do I trust this trace? Answer it per event and you get water-level speed on clean records and iterative robustness on the rest, automatically, with no magic scalar inherited from a decade-old script.
By the numbers
| Method | Mechanism | Failure mode as SNR falls | Relative cost |
|---|---|---|---|
| Water-level division | spectral division with a fixed clamp c | noisy nulls flood into broadband ringing; c biases amplitude | lowest |
| Multitaper | averages independent taper spectra | graceful broadening, mild low-frequency drift | moderate |
| Extended-time multitaper | multitaper plus low-frequency stabilization | graceful broadening, drift suppressed | moderate |
| Iterative time-domain | additive spike construction, no spectral division | small amplitude loss, occasional spurious late spike | highest |
Takeaways
- All four deconvolutions agree on clean synthetics; the divergence only appears at realistic SNR, so a clean-data benchmark tells you nothing about field robustness.
- Water-level's fixed regularization floor is its failure mode: deep noisy spectral nulls become coherent ringing that stacks into fake structure.
- The water-level parameter c is a bias-variance dial with a station-specific sweet spot; a single inherited default is wrong almost everywhere.
- Multitaper degrades gracefully by averaging independent spectral looks, and extended-time multitaper adds the low-frequency stability that deep targets need.
- Iterative time-domain deconvolution is most robust at low SNR because it builds the answer additively and never divides by a null.
- The best method is a property of your noise distribution, so per-event, data-conditioned regularization beats any hard-coded choice.
Limitations
The instrument's fidelity curves and trace artifacts are schematic illustrations of the documented behaviour of each estimator, not the output of a re-run benchmark, and the synthetic impulse response is a three-phase cartoon of a real crustal response. The SNR axis is a single scalar, whereas real noise is colored, time-variable, and partly coherent between components, so the crossover points shown should be read as ordering, not as calibrated thresholds. The comparison also holds the target fixed: for long-period targets or highly reverberant crust the ranking can shift, which is precisely the article's point. Finally, iterative deconvolution's robustness depends on its stopping criterion; run it too long and it will begin fitting noise spikes of its own.
References
[1] Langston, C. A. (1979). Structure under Mount Rainier, Washington, inferred from teleseismic body waves. J. Geophys. Res., 84, 4749-4762. doi:10.1029/JB084iB09p04749
[2] Clayton, R. W., and Wiggins, R. A. (1976). Source shape estimation and deconvolution of teleseismic bodywaves. Geophys. J. R. Astron. Soc., 47, 151-177. doi:10.1111/j.1365-246X.1976.tb01267.x
[3] Park, J., and Levin, V. (2000). Receiver functions from multiple-taper spectral correlation estimates. Bull. Seismol. Soc. Am., 90, 1507-1520. doi:10.1785/0119990122
[4] Helffrich, G. (2006). Extended-time multitaper frequency domain cross-correlation receiver-function estimation. Bull. Seismol. Soc. Am., 96, 344-347. doi:10.1785/0120050098
[5] Ligorria, J. P., and Ammon, C. J. (1999). Iterative deconvolution and receiver-function estimation. Bull. Seismol. Soc. Am., 89, 1395-1400. doi:10.1785/BSSA0890051395
[6] Zhu, L., and Kanamori, H. (2000). Moho depth variation in southern California from teleseismic receiver functions. J. Geophys. Res., 105, 2969-2980. doi:10.1029/1999JB900322
[7] Source: 2018 PhD thesis, receiver-function imaging of the Moho/LAB; comparison of water-level, multitaper, extended-time multitaper, and iterative time-domain deconvolution.
EarthScan builds noise-conditioned receiver-function pipelines that choose their regularization per event instead of per decade. If your network still runs one water level for every station, we should talk.