This data publication contains (i) a slab model of the Cascadia subduction zone, derived from receiver functions, parameterized as depth to the three interfaces: t (top), c (central) and m (Moho), in NetCDF format; (ii) the station measurements of all parameters in the model in tabular and Raysum model file format; (iii) the raw receiver functions in SAC format; and (iv) auxiliary scripts for loading and plotting the data.

A total of 45,601 individual receiver functions recorded at 298 seismic stations distributed across the Cascadia forearc contributed to the slab model. For each station, 100 s recordings symmetric about the P -wave arrival (i.e. 50 s noise and 50 s signal) of earthquakes with magnitudes between 5.5 and 8, in the distance range between 30 and 100 degree, were downloaded from the Incorporated Research Institutions for Seismology (IRIS) data center, the Northern California Earthquake Data Center (NCEDC), and the Natural Resources Canada Data Center (NRCAN). After quality control, radial and transverse receiver functions were computed through frequency-domain simultaneous deconvolution, with an optimal damping factor found through generalized cross validation.

The continental forearc and subducting slab were parameterized as three layers over a mantle half-space, with the subduction stratigraphy bounding interfaces labeled as t (top), c (central) and m (Moho). Synthetic receiver functions were calculated through ray-theoretical modeling of plane-wave scattering at the model interfaces. The thickness, S -wave velocity (VS) and P - to S -wave velocity ratio (VP/VS) of each layer, as well as the common strike and dip of the bottom two layers and the top of the half space (in total 11 parameters) were optimized simultaneously through a simulated annealing global parameter search scheme. The misfit was defined as the anti-correlation (1 minus the cross-correlation coefficient) between the observed and predicted receiver functions, bandpass filtered between 2 and 20 s period duration.

In total, 171, 143 and 137 quality A nodes were determined to constrain the t, c and m interfaces, respectively. At the trench, 105 nodes at 3 km below the local bathymetry were inserted to constrain the t and c interfaces, and at 6.5 km deeper to constrain the m interface, representing typical sediment and igneous crustal thicknesses. A spline surface was fitted to these nodes to yield margin-wide depth models. The spline coefficients were found using singular value decomposition, with the nominal depth uncertainties supplied as weights. The solution was damped by retaining the 116, 117, and 116 largest singular values for the t, c and m interfaces, respectively, based on analysis of L-curves and the Akaike information criterion.

The data set is the supplemental material to Bloch, W., Bostock, M. G., Audet, P. (2023) A Cascadia Slab Model from Receiver Functions. Geochemistry, Geophysics, Geosystems.