Forward Computation
Radial and vertical components of synthetic seismograms are calculated using the Propagator Matrix method of Thomson and Haskell. In the case of offshore setting, we adopt a boundary condition that requires zero traction on the seafloor. Also, we use equations of acoustic wave propagation rather than imitating the acoustic waves by setting Vs ~ 0 on the equations of elastic waves.
We have two options for deconvolution:
Radial component records are deconvolved by the corresponding vertical components using water-level deconvolution. The resulting receiver functions are compared with observed waveforms to measure misfits.
No deconvolution is carried out for this mode. Synthetic radial component seismograms are directly compared with observed waveforms. You may choose this option if you have Green's function estimations instead of receiver functions (e.g., Kumar et al. 2010; Akuhara et al. 2019).
A Gaussian low-pass filter is applied to the deconvolved seismograms:
where a is a tuning parameter.
The synthetic receiver functions are normalized such that the amplitude of the direct P arrival on the vertical component of becomes unity. For receiver function (i.e., deconv_mode = 1), the normalization constant, c, is given by
,
where a is the parameter for Gaussian low-pass filter, N is the number of element in the time series, and 
is the sampling interval.
For Green's functions (i.e., deconv_mode = 0), the amplitude of direct P arrival is measured from the forward-computed synthetic vertical component record. We use this measured amplitude as a normalization factor.
When computing S waves, resultant seismograms are reversed regarding time and multiplied by -1, following the convention for S receiver functions.