Vp/Vs ratios in the Parnaíba Basin from joint active-passive seismic analysis – Implications for continental amalgamation and basin formation

The Phanerozoic intracontinental Parnaíba Basin in northeast Brazil lies atop crust composed of Archaean to Mesoproterozoic cratonic blocks and Neoproterozoic mobile belts. Recently, active and passive source geophysical surveys characterised the crustal structure beneath the basin. We use information from published active-source seismic and new, coincident receiver function (RF) data to obtain Vp/Vs ratios for sedimentary and crustal structure and make inferences about crustal compositions and tectonic evolution. In our approach, sedimentary and crustal Vp/Vs ratios are adjusted to match common conversion point (CCP) images of RFs and known Moho and basement geometry. We use a P-wave model from published wide-angle reflection/refraction (WARR) seismics, and structural features from a deep seismic reflection (DSR) profile. CCP images of the primary RF conversions were used to model the crust, whilst conversions of multiples were used for the sedimentbasement interface. The maximum uncertainties in Vp/Vs are estimated to be 0.15 for the basin and 0.03 for the crust. Vp/Vs ratios in the basin were estimated between 1.7 and 2.2. Lower values correlate with the exposure of older units primarily in the east of the basin, whilst higher values coincide with exposed younger units of the Parnaíba Basin. The obtained crustal Vp/Vs ratios between 1.73 and 1.81 support the previously published segmentation of the crust. In particular, we identified three regions of elevated Vp/Vs ratios, which can be related to proposed Neoproterozoic suture zones underlying the Parnaíba Basin, as well as high velocity lower crust beneath. The high Vp/Vs ratios can be explained by mafic compositions, for example metamorphosed or intruded crust, or fluids and sedimentary rocks entrained into highly deformed crust, typical for modifications related to suture zones. These new deep geophysical models provide important and complementary evidence for crustal amalgamation and the formation of the Parnaíba Basin.


Introduction
The intracontinental Parnaíba Basin in northeast Brazil formed in the Phanerozoic as a consequence of lithospheric subsidence related to postorogenic thermal events and regional lithospheric extension after the Late Neoproterozoic/Cambrian Brasiliano-Pan African orogenesis (Góes and Feijó, 1994;Vaz et al., 2007). The formation of the basin was accompanied by magmatic phases in the Mesozoic thought to be related to the opening of the Central and South Atlantic Ocean (Fodor et al., 1990;Merle et al., 2011;Vaz et al., 2007). The regional basement geology surrounding the Parnaíba Basin implies that its underlying lithosphere constitutes a collage of cratonic blocks and Neoproterozoic mobile belts probably consolidated during the Brasiliano-Pan African orogeny (Cordani et al., 2013b;de Wit et al., 2008). While the lithospheric structure of the Parnaíba Basin could previously only be inferred from regional potential field modelling (e.g. Castro et al., 2016Castro et al., , 2014Lima et al., 2017) and surface geology, in recent years, the crustal and sedimentary structure was subject to detailed geophysical surveys. These surveys were concentrated along an approximately 1400 km long E-W oriented profile across the basin, including deep seismic reflection (DSR), wide-angle reflection-refraction (WARR), magnetotelluric (MT) and receiver function (RF) studies, a summary of which was presented in Daly et al. (2018). The new geophysical images allowed a more detailed characterisation of the crust and uppermost mantle beneath the basin providing new insights about the lithospheric evolution in the region.
We present here an analysis of published DSR and WARR images in combination with newly obtained RF data from the Parnaíba Basin to retrieve information on Vp/Vs ratios for the sedimentary package and crystalline crust and then discuss their relation to possible rock compositions and tectono-magmatic evolution. A combined analysis of these three independent datasets reduces uncertainties observed when analysing such datasets separately. In particular, we expect to benefit from the inherently different sensitivities of the employed methods, as WARR modelling dominantly images horizontal P-wave velocity (Vp) interfaces, DSR images vertical impedance contrasts and RFs are dominantly imaging the vertical S-wave velocity (Vs) structure.
To do this, we developed a novel way of constructing common conversion point (CCP) images of the sedimentary basin, using the conversions of RF multiples. While multiples have previously been used to analyse thickness and Vp/Vs ratio of sedimentary basins (Agostinetti et al., 2018;Assumpção et al., 2009;Licciardi and Agostinetti, 2017;Liu et al., 2018;Srinivas et al., 2013;Yeck et al., 2013;Yu et al., 2015;Zheng et al., 2005), this is to our knowledge the first time that multiples are utilised for basin-scale CCP imaging. In addition, we employ "classic" CCP imaging of the crustal structure using the primary conversions. We obtain Vp/Vs ratios for the sedimentary basin and the crystalline crust by matching basement and Moho interfaces in the CCP images with information from the published DSR and WARR data (see Data and Methods section for details).
We demonstrate the feasibility of this methodology by applying this technique to teleseismic data acquired during the Parnaíba Basin Analysis Project (PBAP) by the Lithospheric Studies Lab (LabLitos) from the University of Brasília and published DSR and WARR seismic profiles (Daly et al., 2014;Soares et al., 2018;Lima et al., 2019). The estimated crustal and sedimentary Vp/Vs ratios are used to complement the existing geophysical images and models and help to refine our understanding of the tectonic evolution of the region. In particular, high crustal Vp/Vs ratios can be related to crustal lithologies typically found in collision and suture zones, or crustal intrusions, and variations of sedimentary Vp/Vs ratios can be related to differently compacted sedimentary rocks.

Geological background and seismic character of the Parnaíba Basin
The Parnaíba Basin is a large intracontinental basin of Phanerozoic age situated in northeast Brazil (Fig. 1). The sedimentary basin covers an area >660,000 km 2 with a thickness of c. 3.5 km in its main depocentre (Góes and Feijó, 1994;Vaz et al., 2007). The basin and its underlying crust and upper mantle have been investigated by recent geophysical studies , including DSR (Daly et al., 2014), WARR (Soares et al., 2018), MT (Solon et al., 2018) and sparse RFs (Coelho et al., 2018;Victor et al., 2020). The depositional and subsidence history of the basin can be separated into four major megasequences (Fig. 2a) from the Late Ordovician to the Cretaceous (Vaz et al., 2007), separated by basin-wide disconformities and marked by two magmatic events of Triassic-Early Jurassic and Cretaceous ages; these form the Mosquito and Sardinha formations, respectively (Fodor et al., 1990;Merle et al., 2011). The thickness of igneous intrusions within the sedimentary strata are estimated to be more than 500 m and concentrated in the centre of the basin (Castro et al., 2018;Daly et al., 2014;Miranda et al., 2018). Although the basin developed dominantly as a sag-basin, there are indications for pre-Carboniferous rifting (Brito Neves et al., 1984;Castro et al., 2016Castro et al., , 2014Daly et al., 2019;Góes et al., 1990;Porto et al., 2018). DSR profiling of the Parnaíba Basin has revealed an asymmetric geometry of the basin, with an elevated, more steeply inclined western margin and a gentle dip toward the eastern edge of the basin (Fig. 2a). WARR modelling estimated P-wave velocity range of 3.3-5.5 km/s for the basin with higher velocities in the deeper successions of the basin at the eastern end of the profile (Soares et al., 2018) (Fig. 2c).
The crystalline crust/basement beneath the sedimentary basin consists of a collage of three main crustal domains (Fig. 2b,c): The Fig. 1. a) Physiographic map of the study area with the location of the seismological array (red triangles), and the coincident WARR seismic line (continuous black line) and DSR line (hatched red line) and existing receiver function studies (red and green triangles). The Parnaíba Basin is outlined by the white line and the main neighbouring tectonic provinces are labelled. b) Distribution of events from the teleseismic dataset (orange circles) used for the CCP imaging. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) Amazonian Craton/Araguaia Belt, the Parnaíba block and the Borborema province (Daly et al., 2014). In the following we briefly describe the seismic structure and character of the crust based on Soares et al. (2018) and Daly et al. (2014) (Fig. 2b,c).
The western Amazonian Craton/Araguaia block comprises the Archaean-Palaeoproterozoic Amazonian Craton and ophiolite-bearing Neoproterozoic Araguaia belt, the latter consisting of metasedimentary folded rocks (Cordani et al., 2013a). The crust is up to 40 km thick in the Amazonian Craton and reaches a maximum thickness of approximately 50 km in the Araguaia-Tocantins suture, beneath the Araguaia belt, which includes a more than 5 km thick high velocity lower crustal (HVLC) body. The crystalline crust consists of a typically 22-23 km thick upper crustal layer (Vp=6.2-6.4 km/s) and a 15-25 km thick lower crustal layer (Vp=6.7-6.9 km/s and 7.2-7.4 km/s within the HVLC body). The Amazonian Craton appears to be horizontally layered with reflective upper and lower crust, while the Araguaia belt shows large-scale, east dipping reflections, which abruptly end at its eastern border, which may be related to the Neoproterozoic Tocantins-Araguaia suture.
The central Parnaíba block is entirely buried beneath the Parnaíba Basin and is formed of Precambrian crust thought to be accreted and stabilised during the Neoproterozoic Brasiliano Orogeny (Brito Neves et al., 1984). Recent geophysical results suggest that the Parnaíba block can be subdivided into at least two sub-domains (Castro et al., 2014;Soares et al., 2018;Solon et al., 2018): The crust in the western Grajaú domain is up to 41 km thick overlain by ~3 km of sediments, with a ~20 km thick upper crust (Vp=6.1-6.45 km/s), and a lower crust of approximately 20 km (Vp=6.7-7.2 km/s), which includes an up to 5 km thick HVLC body. The Grajaú domain is remarkably transparent in the DSR, with no significant reflectivity throughout the entire crustal column, but highly reflective in WARR, which is explained by high Vp/Vs ratios caused by dyke intrusions with a characteristic dimension of 50-200 m (Lima et al., 2019). The Teresina domain in the east of the Parnaíba block has an apparent flat lying Moho at 39 km depth, with a 18-21 km thick upper crust (Vp=6.1-6.2 km/s) and a 16-19 km thick lower crust (Vp=6.6-6.9 km/s) (Soares et al., 2018). A prominent, more than 200 km long mid-crustal reflection is observed at 20-25 km depth, likely related to a large sill intrusion, above a moderately reflective lower crust.
The eastern Borborema province is mainly composed of Neoproterozoic Brasiliano supracrustal folded rocks and reworked Archaean/Palaeoproterozoic basement, both intruded by a number of Neoproterozoic-Early Palaeozoic granites, and intensely sheared by transcurrent, NE-SW, and E-W trending late-to post-orogenic shear zones (Almeida et al., 1981;Vauchez et al., 1995). The crust is composed of a 20-23 km thick upper crust (Vp=6.0-6.4 km/s) with a 10-15 km thick lower crust without any signs of HVLC (Vp=6.7-6.85 km/s). The Moho depth decreases from approximately 41 km to approximately 33  (Daly et al., 2014). The seismological stations used in this study are marked as red triangles. b) Post-stack time migrated DSR image of the Parnaíba Basin displaying the crustal structure of the basement with the main crustal blocks and domains indicated at the top, where the black solid lines represent the boundaries between the blocks that from left to right correspond to the Amazonian Craton/Araguaia Belt, Parnaíba block and Borborema province (Daly et al., 2014). c) 2-D velocity model of the crust and upper mantle obtained from the wide-angle reflection/refraction data, with P-wave velocities indicated by the colour bar and interpretations from Soares et al. (2018). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) km from west to east. The entire crustal column of the Borborema province is characterised by strong reflectivity forming coherent, dominantly east-dipping structures, that could be related to Neoproterozoic contractional structures and Mesozoic-Cenozoic extension.
Three suture zones are proposed to lie along the seismic profile buried beneath the Parnaíba Basin, dating back to the Brasiliano-Pan African (or older) amalgamation of the lithosphere: The Araguaia-Tocantins suture between the Amazonian Craton and the Parnaíba block; the Parnaíba-Borborema suture zone between the Parnaíba block and Borborema province; and lastly the Grajaú-Teresina suture subdividing the Parnaíba block into the Grajaú and Teresina domains (Soares et al., 2018). Although regional tectonic reconstructions require suture zones between the Amazonian craton, the Parnaíba block and the Borborema province, the only clear, exposed structural and lithological evidence exists in the Araguaia Belt, while the other suture zones are inferred from the deep geophysical images and models (e.g. Daly et al., 2018;Soares et al., 2018).

Published and new seismological data
The DSR data used in this study (Daly et al., 2014) consist of a poststack time migrated image of the 1430 km long E-W seismic profile that crossed the Parnaíba Basin and neighbouring blocks near latitude 5 • S (Fig. 2b). The original seismic image has a high trace density, with a common-depth-point (CDP) spacing of 12.5 m and a two-way traveltime time window of 20 s (Daly et al., 2014). There is no direct information intrinsic to the DSR data acquisition to estimate depth uncertainty at the basement horizon. However, the regional DSR profile, which was not tuned to optimise imaging at basement depths, is tied to industry seismic lines that were and which, in turn, are tied to exploration wells in the Parnaiba Basin that penetrate to basement or pre-Silurian strata (e.g. Tozer et al., 2017). The DSR basement interpretation is accordingly considered to be very robust in terms of depth. In any case, we also note that the theoretical vertical resolution of a P-wave assuming the WARR velocity model of Soares et al. (2018) and the lowest frequencies of the DSR data of 6 Hz, would be 150-200 m. Therefore, we regard 200 m uncertainty for the basement interface as a maximum estimate.
The WARR model (Soares et al., 2018) is a 2-D velocity model of the crust and uppermost mantle obtained from forward modelling of deep seismic refraction and wide-angle reflection P-wave data acquired with a vertical-component station each 2 km and 22 shotpoints (~1.5 ton. each). The boundaries and velocity distribution that form the refraction model include five distinct velocity layers: a sedimentary layer (sag basin and half-grabens in the Araguaia belt), upper and lower crustal layers, a crustal "underplate" and mantle, as explained above (Fig. 2c). We estimated an uncertainty of 1 km of the WARR model at Moho depth, using the VMONTECARLO algorithm (Loureiro et al., 2016) (Supplementary material S1). VMONTECARLO employs a Monte Carlo algorithm, which estimates the uncertainty of layered wide-angle models by applying random perturbations of the parameters and uses all refracted and reflected phases. This is in good agreement with the theoretical vertical resolution of a P-wave, which is in the range of 0.8-1 km, applying the lower end of dominant frequencies at approximately 2 Hz for a conservative estimate. Furthermore, we observe a mean picking error of 100 ms, which translates to approximately 0.66 km. The three different estimates are quite consistent and in the following, we will use the maximum value of 1 km for the Moho depth uncertainty. The excellent fit of the two active seismic datasets, both in basement depth (and here additionally tied to wells) and Moho depth (Lima et al., 2019), suggests that they have similar uncertainties and that these uncertainties are rather low.
The passive seismological dataset used consists of continuous recordings by 42 three-component short-period stations of the PBAP array ( Fig. 1a, supplementary material S2) equipped with RefTek DAS130 recorders with Sercel L-4A 3D triaxial sensors. The stations were installed at ~30 km intervals along the WARR profile and operated in different acquisition periods with a duration of at least six months (see Fig. 1, supplementary material S3 and S4).

Receiver functions
The RF method provides estimates of the Earth's impulse response by deconvolving the incident P-wavefield of teleseismic earthquakes from the P-to-S (Ps) converted wavefield (e.g. Langston, 1979;Vinnik, 1977). Deconvolution equalises source and path effects, as well as the instrument response and visualises the Ps conversions as isolated pulses in the time domain (e.g. Ammon, 1991;Clayton and Wiggins, 1976;Langston, 1979;Ligorría and Ammon, 1999;Vinnik, 1977). An RF is composed of the primary conversion of every interface (Ps) and conversions from multiple ray paths (reflections between lower layer interface and the surface), one comprising one P-s conversion (so-called PpPs) and the other comprising two P-S conversions (so-called PsPs or PpSs). In cases of multi-layer models, multiple conversions from shallower layers overprint and potentially disturb the primary conversions at deeper discontinuities .
In this study, RFs were initially derived from a total of 412 events of magnitude (mb) greater than 5 in an epicentral distance between 30 • and 90 • (Fig. 1b). The raw data were processed using the SAC software (Goldstein and Snoke, 2005), including detrending and filtering to eliminate high-and low-frequency noise. We used the processed vertical and radial components in the time window from 60 s before to 90 s after the expected P-wave arrival of each event and generated RFs using a waterlevel deconvolution in the frequency domain implemented in the program pwaveqn (Ammon, 1991) with a waterlevel of 0.001 and a Gaussian factor of 4, retaining frequencies of up to approximately 2 Hz.
The selection of individual events for each station varies according to the acquisition period and signal-to-noise ratio of the waveform recording. First, the processed waveforms with a signal-to-noise ratio smaller than 2 on the vertical and radial component were rejected for RF processing. Secondly, RFs with an amplitude smaller than zero at t=0 s (P-arrival), and RFs with a maximum amplitude of >1 were eliminated from the selection. The applied criteria resulted in 1 to 108 selected high-quality RFs for the individual stations with an average of 37 (see supplementary material S3 and S4). All used RF waveforms and the resulting stacks for each station can be found in the supplementary material S5.
RFs are a well-established seismological technique to investigate crustal and upper mantle discontinuities and seismic properties, using a range of different analyses (e.g. Ammon et al., 1990;Cassidy and Ellis, 1993;Kind et al., 1995;Owens et al., 1987;Sandvol et al., 1998). The fast, simple and computationally inexpensive H-κ stacking technique searches for the crustal thickness-Vp/Vs ratio combination that results in the most constructive summation of amplitudes at the theoretical delay times of primary and multiple conversions, while assuming a fixed crustal P-wave velocity (Zhu and Kanamori, 2000). Inversion of stacked RF waveforms provides one-dimensional crustal seismic velocity models beneath the recording stations, frequently jointly inverted with other data (e.g. Julià et al., 2000;Kiselev et al., 2008;Schiffer et al., 2019). Common Conversion Point (CCP) stacking produces three-dimensional stacks of RFs using the predicted ray-paths of each used earthquake event in an a-priori velocity model (e.g. Dueker and Sheehan, 1997;Jones and Phinney, 1998;Kosarev et al., 1999). Thereby, the RF amplitudes at each time-step are placed at their corresponding theoretical conversion point in three dimensions. These three-dimensional models can then be visualised on two-dimensional sections. The CCP technique exploits the two-or three-dimensional geometry of the medium (if known from prior models) and can provide a more complete structural image of the crust and upper mantle, compared to stacking at individual stations or H-κ stacking.
Under specific circumstances, high-frequency (>1 Hz) RFs can be used to estimate the thickness of sedimentary basins (Agostinetti et al., 2018;Assumpção et al., 2009;Licciardi and Agostinetti, 2017;Liu et al., 2018;Srinivas et al., 2013;Yeck et al., 2013;Yu et al., 2015;Zheng et al., 2005). As mentioned above, in the presence of thick low-velocity layers such as sedimentary basins, RFs are often disturbed by strong multiple reverberations caused by the high seismic contrast between sedimentary rocks and crystalline basement (e.g. Mohsen et al., 2005;Zelt and Ellis, 1999). These multiples can hinder the identification of the basement conversions and may mask other lithospheric discontinuities. In the present approach, we make use of the sedimentary reverberations by depth-converting these to their origin using an a-priori velocity model. As visible in Fig. 5 and Supplementary information 5, the Moho conversion is not compromised by sedimentary at most stations. Instead the multiples are located at upper to mid crustal levels.

Common conversion point imaging
We employed a CCP technique based on Schiffer et al. (2014) using 3D ray-tracing through a prior reference P-wave model assuming different Vp/Vs ratios. We produce CCP images of primary and conversions of multiples, in which the RF amplitudes are placed over a Gaussian distribution at the corresponding conversion point along the exact theoretical three-dimensional ray-path beneath every station. The width of the 3D Gaussian distribution corresponds to the vertical and horizontal seismic resolution. Multiple information from different events and stations at every conversion point are averaged and normalised considering the weight of each amplitude over the individual Gaussian distributions. The method requires a three-dimensional reference velocity model. We use the published P-wave model and information on depth of basement and Moho depth along the profile derived from WARR and DSR data, roughly coincident to the acquired teleseismic data (Daly et al., 2014;Soares et al., 2018). The employed model is two-dimensional, but we create a 2.5D model to approximate threedimensional effects and place the RF signals in a 3D model. Corresponding S-wave velocity models are created by assuming different Vp/ Vs ratios for the sedimentary and crustal layer, which is the base of our modelling procedure for Vp/Vs ratios (see below).
For the mid-to-lower crustal and upper mantle structure, we use a "classic" CCP approach, calculating the theoretical conversion points of primary P-to-S (Ps) conversions in space. In order to enhance large structures in the crustal section, we apply a Gaussian filter with a Gaussian factor of 2 (equivalent to approximately 1 Hz) to the original "high-frequency" RFs generated with a Gaussian factor of 4. In the CCP Note that M2 has a reversed polarity and needs to be multiple by − 1 in order to positively interfere with M2. (d)-(g) Examples of CCP imaging using primaries (d) and multiples (e)-(f). In (f) the negative PpSs+PsPs CCP stack is shown, so that the images have comparable polarity. (g) The sum of multiple CCP images (e) and (f). The stippled lines in (d) to (g) mark the basement depth from a WARR seismic study with uncertainties (Soares et al., 2018). The dotted line in (g) marks the basement interpretation based on this particular CCP image. stack, RF information from usually 2 or 3 stations with overlapping raypaths and Fresnel zones of the used events contribute to one CCP.
We developed a CCP stacking approach using conversions of multiples (PpPs, PsPs, PpSs, see Fig. 3) to image and analyse the sedimentary basin. These multiples have delay times approximately 3-4.5 times larger than the primary conversions (depending on Vp/Vs ratios) and are thereby placed at their predicted conversion point at depth. In a Ps CCP stack, these multiples are wrongly placed at depths 3-4.5 times deeper than the actual convertor (clearly visible in Fig. 5 at 10-20 km depth depending on model). Delay times of conversions of multiples have been used previously to estimate thicknesses and seismic properties of sedimentary basins (Li et al., 2007), as well as in modified H-κ techniques for application to sedimentary basins (Yeck et al., 2013;Yu et al., 2015), but CCP stacking of multiples has only been used to image crustal structure (Wilson et al., , 2003. The multiple CCP approach is in principle the same as for the more common CCP stacking of primary Ps conversions, but allocates the RF amplitudes on a Gaussian 3D-distribution along the multiple ray-paths in three dimensions considering all reflections and conversions (see Fig. 3). For the multiple CCP image of the sedimentary basin, we assume the original dominant frequency content of ~2 Hz (i.e. a Gaussian factor of 4). Considering an average station distance of roughly 30 km, the raypaths and typical Fresnel zones do not create an overlap between the individual stations at depths shallower than 5 km, even for the multiple ray-paths, which span a wider area due to the reflections in the basin.
We produced separate CCP images of the first multiple (PpPs) and another for the second multiples (PsPs+PpSs) and finally produce a stack of both multiple phases (PpPs and PsPs+PpPs), while the second multiple is added with reversed sign, to account for its reversed polarity (Fig. 3d-f). We decided not to stack the primary (Ps) and the multiple phases, as Ps has strong energy close to the P-arrival and might make very shallow interfaces undistinguishable. Seismic discontinuities will be enhanced where PpPs and PsPs+PpPs phases positively interfere.
We estimate uncertainties based on the theoretical vertical seismic resolution of the computed RF, using the dominant frequency content, depth and velocity of the target area. Accordingly, the RFs have a vertical resolution of approximately 0.38-0.65 km at basement depth (assuming P-wave velocities of 3-5 km/s and maximum frequencies of 2 Hz) and 1.65-2 km at Moho depth (assuming P-wave velocities of 6.5-8 km/s and maximum frequencies of 1 Hz) and a horizontal resolution (first Fresnel zone) of 1.8-4.5 km for basement (assuming possible depths of 1-5 km) and 20-27 km for the Moho (assuming representative depths of 35-45 km). The estimated resolutions are consistent with the resulting signal widths on the CCP images (Figs. 4, 5).

Vp/Vs analysis
Analogously to the Constant-Velocity Stack (CVS) technique (e.g. Yilmaz, 2001, chapter 3) applied in the velocity analysis of common mid-point (CMP) gathers to normal move-out correction of seismic reflection data, CCP stacking panels were generated using a range of plausible Vp/Vs ratios (e.g. Castagna et al., 1985;Christensen, 1996) in the sedimentary basin [1.5-2.2] and subsequently in the crust [1.64-1.82]. We produced each panel by assuming different constant Fig. 4. Series of CCP panels using PpPs -PsPs+PpPs multiples stacking from basin-basement interface using constant Vp/Vs ratios ranging from 1.5 to 2.2. Grey line is the basement from DSR and WARR data, the width is equivalent to the uncertainty. The arrows indicate interfaces matching in WARR and RF data in particular regions of the model. Vp/Vs ratios along the entire model, which results in an S-wave velocity model based on Vp/Vs ratio and the reference (WARR) P-wave velocity model for the sedimentary basin and crystalline crust. These panels were visually compared with the sediment-basement and Moho interfaces from the WARR model. The Vp/Vs ratios were finally defined by choosing the CCP panel that yields the best visual correlation between the RF conversions and WARR interfaces in a particular region of the model (see Figs. 4 and 5), i.e. where the known interfaces matched with the centre of the RF signals. The larger wavelengths of the RFs result in uncertainties when trying to match the interfaces derived from the different methods, which are discussed below. For both the sedimentary and the crustal model, we allow laterally varying Vp/Vs ratios, although we only consider relatively large-scale variations in Vp/Vs ratio, especially for the crustal structure. This may result in a regional representation of the Vp/Vs ratios in contrast, for example, to station-by-station estimates from H-κ stacking. The remaining misfits between RFs and WARR interfaces will be discussed below.

Results and discussion
Sedimentary and crustal Vp/Vs ratios in the Parnaíba Basin were obtained by combining CCP stacking of new RF data with coincident existing DSR-WARR profiles. The results are presented in Fig. 6, showing the sedimentary CCP (Fig. 6a), the crustal CCP (Fig. 6b), the final sedimentary Vp/Vs model (Fig. 6 c) and the crustal Vp/Vs model (Fig. 6d).

Sedimentary basin structure
CCP imaging using conversions of RF-multiples exhibits a good agreement with the DSR and WARR images (Fig. 6a) and is able to reproduce some of the internal stratigraphic layering of the basin, as well as basement complexity. We observe a "widening" of some phases and/or deeper-seated phases (Fig. 6a) that do not coincide with the sediment-basement interpretation from active source seismic methods. Some of these phases may represent pre-Silurian rift zones from the early Fig. 5. Series of CCP panels using Ps phase stacking for the crustal structure using constant Vp/Vs ratios ranging from 1.64 to 1.82. Grey lines show sedimentbasement interface (top), a mid-crustal discontinuity (middle) and the top high-velocity lower crustal layer and Moho (bottom). The width of the latter is equivalent to the uncertainty. The arrows indicate interfaces matching in WARR and RF data in particular regions of the model. evolution of the Parnaíba Basin (Castro et al., 2016(Castro et al., , 2014Porto et al., 2018), including the Jaibaras rift at approximately 880-930 km (see Fig. 6a-b). These could, however, also represent local, small-scale anomalies in Vp/Vs ratio not recovered by the relatively smooth Vp/ Vs ratio model resulting from our modelling approach.
Our final Vp/Vs ratio model of the sedimentary basin ranges from 1.7 to 2.2 (Fig. 6c). At the rim of the basin, predominantly in the east, we obtain Vp/Vs ratios of 1.7-1.8 matching well with exposures of the older, Silurian to Early Carboniferous sedimentary successions cropping out in the eastern part of the basin. In the eastern part of the Parnaíba Basin, silicified sandstones of the Serra Grande group crop out and Vp/ Vs ratios as low as 1.5 are plausible for such dry, silicified and lowporous sandstones and quartzites (e.g. Castagna et al., 1985;Christensen, 1996;Domenico, 1984). Values of 1.9-2.0 coincide well with the exposure of younger, perhaps less-compacted Cretaceous units of the central Grajaú sub-basin. Lastly, the highest Vp/Vs ratios of 2.0-2.2 coincide with P-wave velocities of up to 5.5 km/s in the basin in the WARR model, interpreted as the effect of intrabasin intrusions (Soares et al., 2018) and with the location of sill intrusions into the crust (Daly et al., 2014). However, the high Vp/Vs ratios are unlikely to be caused by the intrusions, which are typically associated to Vp/Vs values of approximately 1.8-1.84 for diorites, dolerites or basalts (Carmichael, 2017;Christensen, 1996). Vp/Vs ratios as high as 2.2 are possibly indicative of sedimentary rocks such as shales, mudrocks and carbonates (Carmichael, 2017;Christensen, 1996;Castagna et al., 1985). The generally higher Vp compared to most other parts of the basin is not consistent with the presence of fluids, which would have a declining effect on Vp.
Recent joint inversion of high-frequency RFs and surface waves from 11 stations in the central part of the Parnaíba Basin retrieved a  (Victor et al., 2020). The results show a basement interface often coincident with the DSR, WARR and CCP models, but sometimes a few hundreds of metres deeper. This might (partly) coincide with the deeper phases imaged by the CCP stacking, which we interpret as possible pre-Silurian units. When comparing the CCP image of the basin with the velocity models obtained by joint inversion (Victor et al., 2020), the imaged CCP phases coincide well with the velocity discontinuities in the onedimensional models (Fig. 7a). Positive conversions (red) are located at the depths where downward increasing velocity changes were modelled and vice versa, negative phases (blue) are present at depths where downward decreasing velocity changes (low velocity zones) have been modelled. In this case, the negative (blue) phase in our CCP image at approximately 1.5-2.2 km depth may form a coherent structure along almost the entire length of the basin, which interestingly coincides with a strong, basin wide reflector belonging to the Devonian-Silurian boundary .
Our results for the sedimentary successions are generally consistent with other studies using RFs for the analysis of sedimentary Vp/Vs ratios; for example, Yeck et al. (2013) estimated Vp/Vs ratios in the range of 1.66-2.54 (average 2.05) for selected stations of the US Array and Agostinetti et al. (2018) retrieved sedimentary Vp/Vs ratios of ~1.9-2.4 at a study site in the East African Rift. At the edges of the basin, where the depth of the structures (the basement) are shallower than the seismic resolution (below approximately 0.7 km), we may encounter problems of detecting the right phases. The Vp/Vs ratios at the sides of the basin are therefore less well-constrained (red-hatched areas in Fig. 6c).

Crustal segmentation and composition
Our final crustal Vp/Vs model, after incorporating the sedimentary Vp/Vs model presented above, shows some considerable lateral variation (Fig. 6d). Most crustal Vp/Vs ratios are consistent with the regional lithospheric evolution, while some remain challenging to interpret, as discussed as follows.
(i) The estimated Vp/Vs ratios in the Amazonian Craton range from 1.74 at the western edge of the profile (0-50 km profile distance) to 1.81 near the Araguaia-Tocantins suture zone. Although a Vp/ Vs ratio of 1.74 is close to average values for Archaean and Palaeo-Mesoproterozoic crystalline crust (1.734 and 1.739, respectively; Christensen, 1996;Laske et al., 2013), unpublished H-κ stacking results suggest Vp/Vs ratios as high as 1.76-1.80 in the easternmost Amazonian Craton and Araguaia Belt (Trindade, 2019). The fact that our westernmost estimate is based on the comparison of RF and WARR Moho depths at only one single point at the edge of the WARR model, where it has poor raycoverage and larger uncertainties (Soares et al., 2018, Fig . 6), we deem this result as uncertain. The Amazonian Craton to the east, beneath the Araguaia Belt (approximately 50-300 km profile distance) shows Vp/Vs ratios of 1.78-1.81. The HVLC body modelled by WARR is constrained by positive conversions (red) both at its top and bottom (see Fig. 6d). Such high Vp/Vs ratios can be explained by mafic intrusions, metamorphosed or metasomatic rocks, or alternatively highly deformed crust or entrained volatiles in sedimentary rocks (Christensen, 1996), which can all be expected to be present in suture zones and related high velocity lower crustal bodies.  thus, the Vp/Vs ratios increase as an effect of the intrusions. In the Grajaú domain, we modelled the strongest RF conversions on top of the HVLC body. Matching the lower interface would have resulted in an extremely low Vp/Vs ratio of 1.66 (Fig. 5), which are difficult to explain with any plausible rock composition in this tectonic setting. The top of the HVLC body convincingly forms a strong discontinuity in Vs, visible in RFs, while at the Moho as seen in WARR and DSR data, only weak conversions are visible and not along the whole length of the HVLC body (Fig. 6d). This hints toward a complex nature with possibly gradual increase in Vs from the HVLC to upper mantle. More work is required in order to make stronger inferences about the velocity structure of the HVLC. Based on H-κ stacking, Coelho et al. (2018) found generally similar Vp/Vs ratios in the Grajaú domain, with values of 1.77-1.78 in the centre and more moderate values of 1.73-1.74 at the edges (Fig. 7b). However, the Moho depth obtained by the same H-κ stacking study shows very similar values as the WARR-DSR images, while our CCP approach seems to be more complex and shows the strongest conversion at the top. This mismatch could indicate structural and compositional complexity within this body combined with the different sensitivity of surface waves. Similarly, joint inversion of teleseismic data approximately 50-300 km south of the PBAP profile typically resulted in a gradual crust-mantle transition within the Tocantins-Araguaia suture zone and consistently showed a shallower Moho compared to estimates from H-κ stacking (Trindade, 2019). This supports there being a complex distribution of seismic properties within the HVLC.
(iii) Vp/Vs ratios in the Teresina domain (600-900 km), the eastern constituent of the Parnaíba block, are low in the west (1.73) and slightly higher in the east (1.76). The western part of the Teresina domain hosts mid-crustal intrusions, which we initially expected to increase Vp/Vs ratios, especially because the eastern part with no (or little) evidence for intrusions shows higher Vp/Vs ratios. In the same area (west), we observe a 3-5 km shallower and more complex RF Moho architecture compared to the WARR model (feature 5 in Fig. 6b). Matching this particular part of the Moho with the WARR model would require very low Vp/Vs ratios of approximately 1.68 (compare with Fig. 5), which is rather unlikely because of the tectonic setting with probably deformed and intruded crust. The DSR data show moderate lower crustal reflectivity in this area (Daly et al., 2014). One possibility to explain this misfit is that the RFs are sensing a discontinuity on top or within this reflective layer, whereas the WARR data are strongly affected by velocities nearer the bottom of this layer. Additionally, our model suggests the highest Vp/Vs ratios are found in the sedimentary basin above. Although this is beyond the scope of this paper, one reason for these contrasting features could be a lateral emplacement of the intrusions into the domain, resulting in localised sill intrusions in dominantly non-intruded bedrock, which could retain a generally low Vp/Vs ratio. H-κ stacking results from Coelho et al. (2018) also indicate moderate Vp/Vs ratios of 1.7-1.77 beneath stations in the Teresina domain (Fig. 7b). However, the highest value (1.77) was estimated beneath the westernmost station, while our study obtained the highest Vp/Vs ratio in the easternmost Teresina domain (1.76).
The lowest value obtained from H-κ stacking roughly coincides with the area of apparently shallower Moho in our CCP image. Moho depth from H-κ stacking is otherwise very consistent with WARR, DSR and RF CCP stacking. In the eastern Teresina domain, we again relate the higher Vp/Vs ratios to the proximity of the Parnaíba-Borborema suture zone. The WARR model indicates that this part of the Teresina domain has slightly lower P-wave velocities in the lower crust (6.7 km/s) compared to the western Teresina domain (6.8 km/s) (Soares et al., 2018). High Vp/Vs ratios and lower Vp are strong indicators for the presence of fluids. (iv) Finally, the Borborema province (900-1150 km) in the east of the PBAP profile shows Vp/Vs ratios of 1.74, consistent with averages of Proterozoic crust (1.739-1.744) (Laske et al., 2013). Previous results by Luz et al. (2015) are consistent with the Moho depths obtained by WARR, DSR and RFs with 40 km in the Parnaíba-Borborema suture zone, approximately 50 km north of the profile, and 30-37.5 in the western Borborema province, some tens of kilometres east and north/south of the profile, respectively (Fig. 7b).
We have also identified some structural features in the crustal CCP section. Features 1 and 2 in Fig. 6b are east-dipping structures, that might be related to relict suture zones at least dating back to the Neoproterozoic Brasiliano Orogeny. Feature 1 could relate to the Tocantins-Araguaia suture and the eastward underthrusting of the Amazonian Craton beneath the Parnaíba block. Feature 2 could relate to the proposed Grajaú-Teresina suture in the centre of the Parnaíba block (Soares et al., 2018). There may be further indications for crustal scale westvergent structures in the Araguaia Belt and the western Parnaíba Basin, but the RF conversions are too weak to interpret robustly. Weak conversion patterns in the eastern Parnaíba block and the Borborema block coincide with reflectivity patterns in the DSR data (features 3 and 4 in Fig. 6b). Soares et al. (2018) proposed that the HVLC beneath the Grajaú domain and eastern Araguaia Belt may be related to mafic intrusions analogous to the Mosquito and/or Sardinha Formation. Another possibility to explain the HVLC are metamorphosed lower crust or metasomatised peridotites (e.g. Gernigon et al., 2004;Petersen and Schiffer, 2016;Schiffer et al., 2016), which could be related to the Neoproterozoic amalgamation of Western Gondwana (or even older events). This type of HVLC would make particular sense in the Araguaia Belt, as the body seems to be related to the Tocantins-Araguaia suture and a distinct Moho offset.

Resolution and data sensitivity
We assume uncertainties of 1 km and 2 km at Moho depth for WARR and RF data, respectively, as well as 0.2 km and 0.65 km at the basement depth for DSR and RF data, respectively. When assuming that these uncertainties are normally distributed, we can estimate a total uncertainty by convolving each pair of theoretical probability density functions (pdfs), which results in an uncertainty of ~2.2 km for the Moho and ~0.7 km for the sediment-basement interface. The uncertainties of the Vp/Vs ratio model are assessed visually. While the model values are defined by the best match of the active and passive source data "means", i.e. the centre of the "pdfs", the uncertainties are defined when the upper bound of the RF 68% confidence interval matches the lower bound of the DSR/WARR 68% confidence interval and vice versa.
For the tested range of Vp/Vs ratios in the sedimentary basin [1.5-2.2], we observe a change of ~1.4 km in maximum basement depth in the CCP images. The RF-resolution at basement depth (0.38-0.65 km) is visually confirmed in the CCP images by the width of the RF phases (approximately 1 km). The RF-resolution is rather low with respect to the target depth/thickness of the basin (up to 3.5 km), which is due to the inherently larger wavelengths of the teleseismic RFs, but also artificial smearing of interfaces caused by uncertainties of the reference velocity model used when stacking different RF phases (PsPs+PpSs, PpPs). The width of the RF phases and the DSR uncertainty of 200 m translates into an apparent Vp/Vs uncertainty of approximately 0.15 for the sedimentary basin (assuming conservative data uncertainties).
At crustal level, we observe a change in Moho depth over approximately 10 km for the tested range of Vp/Vs ratios [1.64-1.82] (Fig. 5). The estimated vertical RF resolution at Moho depth (1.65-2 km) is also visually confirmed by the width of the Moho conversion in the crustal CCP image (~3-5 km phase width) (Fig. 6b). The width of the Moho conversion translates to a Vp/Vs ratio uncertainty for the crustal column of approximately 0.02. With the additional Moho uncertainty of approximately 1.0 km from WARR data, lateral changes in Vp/Vs ratio below 0.03 are therefore not considered robust (again, assuming conservative/maximum data uncertainties). This is, however, only the case for the transition between the Grajaú domain to the Teresina domain of the Parnaíba block, and the boundary between eastern Teresina domain/Parnaíba-Borborema suture zone and Borborema province. Here, the "boundaries" can be regarded as gradual transitions.
Although it is noted that there may also be inherent mismatches and inconsistencies when comparing the data due to different sensitivities of the employed methods, we consider the major features discussed reliable within the desired maximum uncertainties. The observation that all three datasets do show very good agreements of basin and crustal structure, adds to our confidence. Careful analysis and discussion of these existing mismatches can add further insight into the nature of complex lithologies and structures.
With regards to comparing Vp/Vs ratios using the CCP panel approach from this study with published results from H-κ stacking and RF inversion (Coelho et al., 2018;Trindade, 2019), we need to consider that CCP imaging exploits the two-dimensional variations beneath the stations, while the latter two methods provide punctual results assuming a one-dimensional structure beneath every station. We do expect mismatches between the different methods. However, all published RF results are generally consistent with our outcomes and interpretations.
One outcome of our modelling approach is that the Vp/Vs model is defined in lateral sections with vertical boundaries. In this approach, we can only compare the interfaces from different datasets at a given profile offset and choose the overlying Vp/Vs ratios accordingly. This works reasonably well, because teleseismic earthquakes generally arrive almost vertically beneath a station and we are taking into account threedimensional effects as completely as possible. Hence, it is defensible to model the Vp/Vs ratios vertically beneath the stations. Nevertheless, the vertical nature of our model may seem somewhat incompatible with the complex geology and especially the dip of suture zones. The Vp/Vs ratios in these "blocks" must be examined and interpreted as average values for all overlying strata and its three-dimensional variation.

Summary and conclusions
We present a novel multi-panel RF CCP stacking procedure, which systematically compares the depth of Moho and sediment-basement interface visualised in CCP sections with those obtained by WARR modelling to obtain Vp/Vs ratios in a simple and quick way. Combined analysis of RF CCP images and WARR-DSR data/models is an effective and robust method to obtain Vp/Vs ratios and to draw additional information from the sedimentary and crustal structure, compared to application of one of the methods alone. In particular, the analysis of different datasets to constrain Vp and Vs independently adds to the robustness of this approach, as RFs are primarily sensitive to Vs, whereas wide-angle and reflection seismic estimates are primarily sensitive to Vp and impedance, respectively. The inherently different sensitivities of the employed seismic methods allow additional inferences about the structure and composition of the crust. This is especially interesting for the lower crust, where we sometimes observe slight mismatches between one and the other two datasets, which is likely the result of lower crustal complexity.
Although the data density and resolution of RFs are not comparable to active source seismic methods, there is potential in using CCP stacking of RF multiples for the imaging of sedimentary basins. This is especially interesting as a cost-efficient tool to characterise the structure of sedimentary basins and underlying crust in areas where no prior deep structural information is available or where access is difficult and regular geophysical exploration expensive, for example when preparing and planning active source experiments.
For our joint analysis of active and passive source seismic data to obtain Vp/Vs ratios of the Parnaíba Basin and crystalline crust beneath, we estimated maximum uncertainties of 0.15 and 0.03, respectively. The results reveal a lateral subdivision of the crust into blocks of characteristic Vp/Vs ratios that can be interpreted in terms of segmentation of the crust and lithospheric amalgamation processes.
The intracontinental Parnaíba Basin itself has Vp/Vs ratios of 1.7-2.2, with low values estimated in areas of exposed older units that are more compacted, in particular on the eastern end of the basin where strongly silicified sandstones are exposed, and the higher values in the centre of the basin, where younger, less compacted units are exposed and where magmatic intrusions where mapped. The sedimentary CCP images correlate with previously published sedimentary stratification derived from high-frequency RF inversion (Victor et al., 2020). The results support a model in which the depocentre has shifted from east to west during evolution of the Parnaíba Basin.
The WARR velocity model as well other geophysical and geological information indicate that the crust beneath the Parnaíba Basin comprises several distinct blocks of differing character, namely the Amazonian Craton to the west, the Parnaíba block in the centre, which itself can be sub-divided into the Grajaú and Teresina domains, as well as the Borborema province in the east. The published seismic data also indicate the presence of HVLC beneath the boundary zone between Amazonian Craton and Parnaíba Basin and continues throughout the Grajaú domain. These features are generally consistent with our Vp/Vs analysis.
In addition, the new CCP stacking Vp/Vs analysis reveals zones of elevated Vp/Vs ratios, spatially related to the locations of known and proposed Neoproterozoic (or older) suture zones and HVLC identified by the WARR data. The regionally highest Vp/Vs ratios (1.76-1.81) are found in the easternmost Amazonian Craton (beneath the Tocantins-Araguaia suture) to the boundary between Grajaú and Teresina domain, in the centre of the Parnaíba block, where another suture zone has been proposed (Soares et al., 2018). The entire Grajaú domain and the Tocantins-Araguaia suture zone are underlain by ~3-8 km thick HVLC and Vp/Vs ratios increase with increasing HVLC thickness and vicinity to suture zones. Although not as distinct, slightly higher Vp/Vs ratios (1.76) are also found close to the Parnaíba-Borborema suture, here without any indications of a HVLC body. All other areas have been estimated to have moderate Vp/Vs ratios (1.73-1.74), typical for Archaean-Proterozoic crust. We relate the high Vp/Vs ratios in suture zones to high-pressure rocks, entrained sediments and fluids, and fractured rocks.
Our results support previously proposed structures and evolutionary scenarios of the Parnaíba Basin and its underlying crust: 1) High Vp/Vs ratios, roughly coincident with known or proposed suture zones are a strong additional indicator for their existence, especially because two of these (Grajaú-Teresina and Parnaíba-Borborema) are completely buried beneath the Parnaiba Basin with no surficial expression. 2) The two highest Vp/Vs ratios (1.81 and 1.79) are estimated for the crust including thick HVLC. This is a good indicator for mafic compositions, metamorphism, intrusions or fluids in the crust. Whether the whole HVLC as inferred from the WARR data forms a single body, or whether it is sub-divided into two, is unclear. Several factors could indicate two different origins of the two parts: Both have different Vp (7.3 and 7.1 km/s), different morphology ("S-shaped" and subhorizontal), the HVLC beneath the Grajaú domain is associated to surface volcanism and crustal intrusions (Lima et al., 2019), the western HVLC body is spatially related to a suture zone, and both are connected by an only ~3 km thin section with only poor raycoverage and no reflections at its lower interface (Soares et al., 2018). The possibility of two separate and generically different bodies suggests that the western part (beneath the Tocantins-Araguaia suture) may be related to collision processes (i.e. metamorphosed, metasomatised lower crust/upper mantle) and the eastern part (beneath the Grajaú domain) is related to Cretaceous magmatism.
3) The successively younger exposed sediments from the east to the west of the Parnaíba Basin, Triassic-Jurassic magmatism at the surface (Mosquito Formation) and the presence of at least the eastern part of the HVLC in the Grajaú domain, directly beneath the Cretaceous depocentre and the Triassic-Jurassic magmatic rocks, strongly suggests that the basin depocentre shift was related to intrusions into the Grajaú-lower crust, followed by cooling, loading and sagging. The fact that both RF inversion and our CCP image of the Parnaíba Basin show velocity decrease/negative conversions at similar depths suggests basin-wide high-velocity layers, potentially related to Cretaceous sills.
This study adds valuable, complementary geophysical information constraining the crustal and sedimentary structure and composition of the Parnaíba Basin and its underlying crystalline crust. These constraints contribute to the understanding of amalgamation processes affecting unexposed crystalline basement in NE Brazil as well as the processes of deep lithospheric processes (represented by the HVLC) responsible for the formation and evolution of the Parnaíba Basin itself.

Credit author statement
The method was developed by ML and CS, with critical input by JS and RS. The manuscript was developed by CS with major contributions by ML, JS, RS and RF, and other substantial input by VA, FL, FR and CT. Figures were prepared by CS, with substantial help by ML and JS. The PBAP project was led by JS and RS, and data collected was majorly assisted by ML, VA, FL, FR and CT.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.