Sediment thickness across Australia from passive seismic methods

the

Australia is an old and largely tectonically stable continent.As such, the sedimentary blanket of the continent has remained largely unaffected by tectonic events since the Mesozoic (∼ 250 million years ago).Sediments and regolith now cover 75-80% of the continental surface (Figure 1).A long history of sedimentation stretching as far back as the Archean has produced basin provinces with sediment accumulations of up to 10 km or more (Yeates et al. 1984;Myers et al. 1996).These various basins are rich in natural resources and mineral deposits critical for the Australian economy and the transition to net-zero (Mudd et al. 2018).
Older basins typically found in Western Australia, such as the 400 Ma Canning Basin, are often characterised by multiple phases of sedimentation, resulting in extremely thick sediment layers that can host hydrocarbon resources (e.g.Towner & Gibson 1983).Younger Australian sediment thickness from receiver functions 3 eastern Australia, remain under-explored in some areas, despite being a prolific oil and gas producer (Cotton et al. 2007).
It is important to understand the structure of sedimentary basins for a number of reasons.
From a seismic hazard perspective, sedimentary basins are known to significantly amplify seismic waves and earthquake ground motion (e.g.Molnar et al. 2014;Wirth et al. 2019).
Despite it's general tectonic stability, Australia is still subject to a number of intraplate earthquakes.Although relatively infrequent and sporadic compared to earthquakes at plate boundaries, intraplate earthquakes are known to cause significant damage (Bowman 1992;Bilham 2014;Geoscience Australia 2023).The city of Melbourne, sitting on a sedimentary basin, recently experienced a M w 5.9 earthquake with an aftershock sequence that lasted several months (Mousavi et al. 2023).With many major cities located on sedimentary basins, understanding of the local basin structures and their thicknesses is imperative for seismic hazard assessment.
From a natural resources perspective, sediment basins host mineral deposits, hydrocarbons and groundwater.However, they also obscure the crustal geology below, making estimating the thickness of the sediments a key first step in any exploration.In response to a sharp decline of new resource discoveries in Australia, and with the majority of easily accessible near-surface resources already discovered and exploited, continent-wide initiatives have been undertaken to develop new strategies for discovering economically-viable deposits beneath the sediment and regolith (Collett & McFadden 2014;Kennett 2020).One of the key objectives identified was a need to assess the depth of cover.Accurate measures of the sediment thickness can be obtained through borehole drilling or active seismic experiments.
However, these are expensive and often impractical endeavours in the remoteness of the Australian outback when it is necessary to transport large equipment to locations where sealed roads do not exist.
Motivated by this problem, this work aims to characterise the sediment thickness across Australia using single-station and passive seismic methods.Passive seismic approaches are significantly cheaper, less logistically challenging, and have a smaller environmental footprint than active seismics or boreholes, making them a promising avenue for geophysical exploration in remote Australia.Examples of passive seismic techniques that can be employed to study sedimentary structures include ambient noise tomography (e.g.Yamaya et al. 2021;Chen et al. 2023), wavefield autocorrelations (e.g.Saygin et al. 2017;Romero & Schimmel 2018), receiver function inversion (e.g.Piana Agostinetti et al. 2018 ;Cunningham & Lekic 2020) and spectral methods (e.g.Perron et al. 2018;Schleicher & Pratt 2021).Recently Agrawal et al. (2022) outlined a new simple approach using teleseismic receiver functions calibrated with borehole observations to obtain an empirical relationship between seismic observations and sediment thickness.This approach circumvents the need to know the underlying velocity structure.However, their study area focused only on the state of South Australia, and the applicability of their empirical relation to other regions with different geological histories is not asserted.In this work, we extend their approach to obtain an continent-wide empirical relation.As far as the authors are aware, this is by far the largest study of its kind, using data from ∼ 1, 500 seismic stations and ∼ 84, 000 boreholes with near complete coverage of the continent, thereby allowing for analysis of basins across a wide range of thicknesses and geological ages, and creating a generalisable empirical relation between receiver function data and sediment thickness.

Receiver functions in the presence of sedimentary basins
Teleseismic receiver functions highlight the conversion of seismic energy from P to S, and vice-versa, at interface between layers.The larger the impedance contrast at the layer interface, the stronger the energy conversion will be.This is a local measurement, arising from the discontinuities beneath the recording seismic station with effects from the source and path removed by deconvolution (Burdick & Langston 1977).As such they have been    Piana Agostinetti & Amato 2009 ;Chen et al. 2010;Kennett et al. 2011;Gallacher & Bastow 2012;Li et al. 2014;Chopra et al. 2014;Mandal & Biswas 2016;Spieker et al. 2018;Hao et al. 2023;Mroczek et al. 2023;Kennett et al. 2023).
The strong impedance contrast between the crustal basement rock and overlying lowvelocity sediments creates a strong P s b conversion.This causes the direct P arrival to broaden and decrease in amplitude, resulting in the arrival of the P -to-S conversion at the sediment basement (P s b ) to be the dominant first arrival peak on the radial receiver function (Yeck et al. 2013;Yu et al. 2015).This arrival is delayed relative to the direct P phase (Figure 2).The impedance contrasts at the basement and free surface also trap energy within the sedimentary layer, resulting in a strong reverberatory signal on the radial receiver function (Levander & Hill 1985).This masks the arrival of phase conversions from deeper depths, and thus has been a long-standing issue of studies investigating, for example, the Moho.A number of approaches have been used by such studies to mitigate the effects of the sedimentary layer, including forward modelling (Sheehan et al. 1995), wavefield continuation and decomposition (Langston 2011) and H-κ stacking (Zhu & Kanamori 2000;Yeck et al. 2013).A much simpler approach was proposed by Yu et al. (2015) by applying a resonance removal filter.The construction of the filter is based on the strength of the reverberations and the two-way travel-time of shear waves in the sedimentary layer, t 2 , as measured from the autocorrelation of the radial receiver function (Figure 2).The two-way travel-time (t 2 ) and the delay of the P s b phase relative to the P phase (t P s b ) are both directly related to the sediment thickness by the shear wave velocity within the sediments, and thus can be used to estimate that thickness (Figure 2).
Receiver functions, specifically the P s b conversion delay, have been used to study sedimentary structures previously in South Australia (Agrawal et al. 2023)  Previously in Australia, Clitheroe et al. (2000) inverted receiver functions at 65 stations across the continent for shear wave speed, including a sedimentary layer that was allowed to vary up to 2 km thick.23 of the 65 stations exhibited a sedimentary layer, primarily located in the east.They found the thickest sediments in the Eromanga Basin (Figure 1), with a thickness of 1.71 km.Clitheroe et al. (2000) were also careful to distinguish between "soft" and "compacted" sediment, highlighting that a compacted sediment may have a similar seismic velocity to the underlying crystalline basement and would thus not be identifiable in receiver function data.Their discussion of sediments was thus limited to, for example, uncompacted alluvium, clays and desert sands.We expect a similar distinction to arise  in our results, particularly in the west where many of the old Precambrian basins have undergone extensive metamorphism and compaction since being deposited, which should reduce the impedance contrast between the (meta-)sediments and the crystalline basement.
In our work, discussion of sediment refers to the layers between the surface and the top of the crystalline basement, and thus can include "soft" and "compacted" layers.

Synthetic demonstration
To demonstrate the principles of the effects of sediments on receiver functions, we begin with a simple synthetic experiment with a 1D Earth model, investigating the stability of the t P s b and t 2 measurements for varying sediment thickness.Our Earth model consists of a sediment layer overlying the crust, which has an average crustal velocity (3.6 km/s) and thickness (32 km) for Australia (Salmon et al. 2012(Salmon et al. , 2013)).The thickness of the sedimentary layer is varied from 0 to 3 km thick, with a shear wave velocity determined by a power law compaction model (Chapman & Godin 2001;Cunningham & Lekic 2020) of the form where p = 0.5 implies a linear increase in shear rigidity µ with depth z.The velocity at the surface v 0 and the growth factor a are chosen such that there is a significant velocity difference at the sedimentary basement for the thickest sediments, giving v 0 = 0.1 km s −1 and a = 1 √ km s −1 .The corresponding compressional velocity, v p , is determined using the empirical Brocher relations (Brocher 2005).Figure 3 shows 1D shear velocity profiles, radial receiver functions and their autocorrelations, and the trend of t P s b and t 2 with sediment thickness.We see that the t P s b measurement is 0 s for the case with no sediment, and would correspond to the direct P arrival.The t 2 estimate is unstable until sediment thickness reaches more than 300 m.Such an observation is to be expected.In sediments where the thickness of the layer is small compared to the wavelength of the seismic wave, reverberations are not observed, meaning t 2 cannot be accurately constrained for such small sediment thicknesses, whereas t P s b can obtain a stable measurement.range where the t 2 measurement is stable, it is largely bounded by 2t P s b and 3t P s b , which we can use to eliminate outliers in our real data.Figure S1 highlights the receiver functions and their autocorrelations for sediment thicknesses of 0, 1.5 and 3 km to further demonstrate the change in character of the receiver functions.

Receiver Functions
We obtained waveforms from 40 permanent and temporary seismic networks across Australia, comprising over 1,500 individual seismic stations (see Figure 1, Table S1 in the supplementary information and the Data Availability Statement for full details).Events at epicentral distances of 30 • < ∆ < 95 • and magnitudes M w > 3 were selected (Figure S2).
For temporary networks, all events for the duration of the network's deployment that met the criteria were selected, whereas for permanent networks we selected events across two years of data during 2018 and 2019.Horizontal components were rotated to the RTZ coordinate system, and iterative time domain deconvolution (Ligorría & Ammon 1999) was used to deconvolve the radial from the vertical component to obtain the receiver function.
The receiver functions were then moveout corrected to account for epicentral distance, and filtered to [0.1, 1.0] Hz as in Agrawal et al. (2022).
An automated quality control was applied to the receiver functions.We keep receiver functions with signal-to-noise ratio greater than 2, as calculated on the vertical component receiver function with a noise window up to 5 s before the P arrival and a signal window from [−5, 25] s around the P arrival.We also only keep receiver functions where the largest arrival occurs within 2 s of the direct P arrival.We note that this immediately places a maximum measurable t P s b of 2 s, however this is necessary to avoid automatically picking deeper conversions such as that from the Moho, which typically appears at around 5 s for average continental crust.Finally, we discard stations for which there are fewer than 10 receiver functions remaining following quality control.This is to ensure a sufficient number of receiver functions are stacked to obtain a clear P s b conversion, and also ensure a sufficient azimuthal coverage to mitigate local effects such as dipping layers.A rare but possible case would be a station with receiver functions all coming from a high seismic zone such as the Tonga-Kermadec subduction zone or Japan.2).This is related to sediment thickness D simply by the difference in travel time of S and P in the sedimentary layer, assuming vertical incidence of the wavefront on the sediment-basement interface, an assumption known to introduce negligible errors (Cunningham & Lekic 2020).For shear velocity continuously varying with depth, as one might expect for progressively more compacted sediments, we have where v s S (v s P ) is the shear (compressional) wave velocity in sediments which varies with depth z.
A subset of our stacked receiver functions is shown in Figure 4.There is a clear change in character of the receiver functions as the first peak (P s b phase) is delayed, which is also visible in the results from our synthetic experiment Figures 3 and S1.When the delay t P s b is between 0-0.2 s (i.e.yellow traces), the first arrival is sharp and no reverberations are seen, indicating little to no sediment.In many cases a smaller positive peak can be clearly seen around 5 seconds consistent with the arrival of the Moho converted phase.As t P s b increases, for values between 0.2 to 0.6 seconds (orange to red colours), a reverberatory signal starts to dominate with a single consistent frequency that decays in amplitude over time.On such traces it is no longer possible to clearly pinpoint a particular arrival from the Moho.For large delay times, greater than 0.6 seconds (dark red to black colours) the receiver function becomes increasingly complex as multiple phases and reverberations overlap.In such cases  there is no longer a single consistent frequency present, and thus the reverberations cannot be effectively removed with a resonance filter (e.g.Agrawal et al. 2023).
Figure 5 maps the distribution of t P s b we obtain at seismic stations across Australia.
Explicit values are given in Table S1.Overall the pattern of delay times shows good spatial

Shear wave two-way travel time
The large impedance contrast between sediment and the crystalline crust can trap shear wave energy which reverberates within the sedimentary layer.The reverberations mask energy conversions from deeper discontinuities such as the Moho, and as a result a number of approaches for removing the effect of sedimentary layers from receiver functions have been studied.Yu et al. (2015) suggested a filter approach that depends on the two-way travel time of shear waves in the sedimentary layer, which is directly related to sediment thickness.
We thus use their approach to measure the two-way travel time of shear waves in the sedimentary layer t 2 .The first reverberation of S waves in the sedimentary layer has negative polarity due to the reflection at the free surface and arrives t 2 seconds after the first P s b arrival.Thus t 2 can be measured as the time between the first peak and first trough of the radial receiver function, or, more simply, the time of the first local minimum in the autocorrelation of the radial receiver function (Yilmaz 2001;Yu et al. 2015).This process is demonstrated in Figure 2. Cunningham & Lekic (2019) suggested first fitting a decaying sinusoid to the autocorrelation, and reading t 2 from the fitted function to help determine if reverberations are present in the receiver function and to help automate the process.
Following Agrawal et al. (2023), we instead measured t 2 from the autocorrelation of each individual (not the stacked) receiver function in order to better identify and remove outliers (Figure S3).Considering a shear wave (of velocity v S ) reverberating in a sediment layer of Australian sediment thickness from receiver functions 15 thickness D, t 2 can simply be defined as and thus is related to t P s b as This provides a lower bound for t 2 , but with the upper bound still unconstrained.Agrawal et al. ( 2023) determined the upper bound manually for each station, but this is impractical for the 1,500 stations used in our work.As such, we determine the upper bound based on synthetic tests (see Figure 3) which show that t 2 is well-bounded by 2t P s b ≤ t 2 ≤ 3t P s b .
If t 2 does not fall within this desired range, the autocorrelation likely does not reflect the sediment reverberation process.Thus for each station we discard autocorrelations of receiver functions for which t 2 does not fall within this desired range (Figure S3).We take the median t 2 of the remaining autocorrelations as the final t 2 value for that station.
Following this protocol, Figure 6 maps the distribution of t 2 measurements we obtain at seismic stations across Australia.Explicit values are given in Table S1.As expected, the spatial pattern correlates strongly with Figure 5, although stations with the smallest t P s b delay times (< 0.2 s) tend to not have any valid t 2 observations, reinforcing that the two-way travel time measurement is only meaningful where there is a significant sedimentary layer.
The largest values of t 2 are observed in the Eromanga Basin.Values also tend to be larger at the interior of basins and lower towards the edges.The strong similarity between these two maps, as well as with the borehole-estimated sediment thicknesses (Figure 7) validates the use of receiver functions for estimating sediment layer thickness.Further, t 2 as measured in this study is only accessible in regions with an already significant sediment layer (Figure 3) such that reverberations are visible in the receiver function.depth to each unit.However, the depth to the sediment-basement interface at each borehole is not reported, and thus it is left to the user to interpret the stratigraphic units to find the basement depth.Considering the large number of boreholes in the ABSUC, we need to automate our interpretation of the database to identify the depth to the sedimentbasement interface.To do this, we cross-reference with the Australian Stratigraphic Units Database (ASUD; Geoscience Australia and Australian Stratigraphy Commission 2017) and the Australian Geological Provinces Database (AGPD Raymond et al. 2018).Starting with the deepest stratigraphic unit, the ASUD identifies the geological province the unit belongs to and the AGPD identifies if said province is sedimentary.If so, the depth to the bottom of that unit is taken to be the depth to the basement.If not, the process is repeated moving upwards for all the stratigraphic units until a sedimentary unit is found or the top of the borehole is reached.
Our automated approach is validated against an independent database compiled by the Geological Survey of South Australia (GSSA SARIG 2021), in which basement depth estimates are provided for South Australian boreholes (Figures S5 and S6).Overall there is a good agreement between the two different estimates of sediment thickness from boreholes, with the overall standard deviation of the estimate error at around 50 m, which is much smaller than what we expect to be able to resolve with passive seismic techniques.
In Figure 7 we show the depth estimate we obtain for each borehole.

Mapping borehole sediment thickness to seismic stations
In order to determine an empirical relationship between t P s b and sediment thickness (surface to crystalline basement), we first need to determine how to assign a sediment thickness beneath a given seismic station based on the borehole data (Figure 7).We first select a radius of 0.25 • around a station to search for available boreholes.This radius is chosen considering the likelihood of geological consistency and typical station spacing, while also allowing for a sufficient number of stations to have an estimate from which an empirical relation can be obtained.Figure S7 shows the distribution of seismic stations that have at least one borehole within the 0.25 • radius -about half of all seismic stations studied.These stations have between 1 and over 3,000 boreholes in their vicinity.Figure S7 also shows where those stations with no boreholes within the allowed radius are located.The furthest distance to the nearest borehole for these stations is over 400 km for stations in southern Tasmania.
Having identified which boreholes are located within 0.25 • of each station, the sediment thickness estimate is then taken to be the deepest such measurement within the prescribed radius.We investigated other estimates, including the mean sediment thickness of all boreholes within the radius, the inverse distance weighted mean, and the sediment thickness of the closest borehole (Figure S8), all of which produced consistent depth estimates.Our preference is for the maximum depth as this would be the most useful in practice considering the cost to drill.This provides a more conservative estimate of the sediment thickness, as it is preferable to over-estimate than underestimate the sediment thickness when drilling towards a target.

Constructing the empirical relationship
Following Equations 1 and 2, sediment thickness should approximately vary as the square of t P s b .As such we fit a quadratic to the binned measurements of t P s b and the sediment thickness as measured by nearby boreholes.where D is in metres and t P s b is in seconds.The other three estimates (mean, inverse distance weighted mean and closest borehole) all produce similar fits.Our preferred fit systematically gives greater thicknesses with increasing t P s b time, as expected, but also has a visually similar rate of change as the others over the range of measured seismic data.We note that the minimum of this relation is at (t P s b , D) = (0.2, 382.62).Thus measurements in the range 0 < t P s b < 0.2 would give decreasing estimates of D with increasing t P s b , which is unphysical and inconsistent with the compaction model.Furthermore, the minimum implies our method cannot resolve the difference between the case of no sediments and a sediment layer thinner than ∼ 400 m.This observation is largely consistent with the limit of consistency between t P s b and t 2 that we found in our synthetic experiment (Figure 3).
In Figure 8, all data from across Australia has been included to formulate a generalised equation that incorporates the full spectrum of geological age.However, Agrawal et al. (2023) suggested hints of an age dependence in their empirical relation, in that older Proterozoic regions did not appear to follow the same trend as younger Phanerozoic regions.Indeed, they found that t P s b did not vary much at stations in Precambrian provinces despite the large range of sediment thicknesses expected.They argue that older, more compacted and metamorphosed sedimentary layers have less of an impedance contrast with the underlying crust, resulting in no clear P s b signal.To investigate this further, we separate our seismic data by geological era.The corresponding geological age at each seismic station is taken to be the age of the oldest basin beneath the station (in the case of overlapping basins).We obtain a similar result to Agrawal et al. (2022) for Precambrian/Proterozoic regions (Figure S10) in which only small delay times (< 0.5 s) are retrieved despite borehole sediment thickness estimates of several kilometres.Also, comparing Figures 5-7 we can identify regions such as the young Cenozoic Murray Basin in the south east that has a strong signal in the seismic data, but less so in the borehole data, where as the older Mesozoic Eromanga Basin in the east has a strong signal in both datasets.As such, we also fit different empirical relations for each era of the Phanerozoic (Figure 9).ern Australia, we find sediment thicknesses of generally less than 0.5 km.In these regions, as previously mentioned, the sediment is too compacted or metamorphosed for a reliable t P s b measurement to be obtained, resulting in a likely underestimate of sediment thickness.
OZSEEBASE on the other hand predicts sediment thicknesses ranging between 5-10 km in these regions.In general, the overall pattern of sedimentary regions predicted by our empirical relation is consistent with that of OZSEEBASE.However, in the younger eastern Australia, OZSEEBASE tends to predict sediments 2-3 km thicker than our empirical relation.A point of difference between our empirical predictions and the OZSEEBASE model can be seen in the western part of the Eucla Basin (see also Figure S4).Our approach clearly delineates the western margin of the basin where it abuts with the Albany-Fraser Orogen.
The sediment thickness clearly increases towards the interior of the basin in our results, whereas OZSEEBASE seems to simply interpolate zero thickness across most of the Eucla Basin.With the original geophysical and geological data used to construct OZSEEBASE not publicly available, it is difficult to determine the reliability of its features and interpret the differences between the model and our empirical predictions.
Recently, Chen et al. ( 2023) generated estimates of sediment thickness across Australia using the 3.1-3.3kms −1 contours of their new ambient noise tomography model of the Australian crust.Their results are highly similar to OZSEEBASE, although again OZSEEBASE tends to predict thicker sediments.The advantage of our approach over that of Chen et al. (2023), with both being passive seismic approaches, is that ours is a result of a direct measurement of the sedimentary layer (the P s b phase) at a single station, rather than an interpolated indirect measurement (i.e. the 3.1-3.3kms −1 contour).Our approach is also much more straight-forward, requiring the picking of a single phase, whereas the ambient noise tomography is a two-stage inversion from ambient noise to shear wave speed.Thus our approach is much better suited for obtaining a fast first order estimate of sediment thickness.
The velocity contouring approach of Chen et al. ( 2023) is likely more applicable than our approach in regions of older, more compact sediment where t P s b is difficult to measure, such as the Kimberley or Canning basins.The empirical relations (Equations 5 and 6) provide a good first estimate of the sediment thickness.With a relatively brief deployment and small footprint of a single seismometer (for less than 6 months), an estimate of sediment thickness can quickly be obtained in a manner that is much cheaper than drilling a borehole and less time consuming than performing an inversion or undertaking full H-κ stacking.Further, Equation 5 can be used with little geological knowledge of the study region, the only requirement being that the underlying basin is no older than the Cambrian period.Equation 6 only requires knowledge of a rough age bracket to obtain a potentially more accurate estimate of sediment thickness.Importantly, Equation 5 can be used even when basins are overlapping.However, the uncertainties may be too large depending on the application and the level of accuracy required, in which case more expensive techniques, such as active seismics, will still need to be employed.Care also needs to be taken when a small t P s b is measured.As previously discussed, the range 0 < t P s b 0.2 does not necessarily follow the trend of a sediment compaction model, and Equation 5 cannot constrain the difference between zero sediments and sediments thinner than ∼ 400 m.Finally, the empirical relations are not useful in Proterozoic regions.The t P s b measurement in these regions is typically underestimated (discussed further in the following).Indeed there are regions in Western Australia where sediments are expected to be over 10 km (Yeates et al. 1984;Myers et al. 1996).These would require t P s b > 2.8 s according to our overall empirical relation.

Unquantified uncertainties in seismic measurements
The fundamental assumption of our approach is that the first, largest peak on the radial receiver function is the signature of the P -to-S conversion at the interface between the sediment and the crystalline basement.This implies that the crystalline basement is the interface with the largest impedance contrast, which may not in fact always be the case.For example, it is conceivable that immediately above the crystalline basement there is an old, well-compacted sediment layer with similar impedance as the basement, above which there is loose sediment of lower velocity.In this case the first large peak on the radial receiver function could conceivably be the signature of the interface between the two sediment layers, and our t P s b measurement would underestimate the overall sediment thickness.This potentially explains our observation of low t P s b at stations on Precambrian basins (Figure S10), and consistent with the findings of Agrawal et al. (2022).Further, the picked peak may also include the signature of the first reverberatory phase P P s b in thin sediments (Cunningham & Lekic 2019, 2020), which may either positively or negatively bias the t P s b measurement (Agrawal et al. 2022).The geographical consistency of our t P s b and t 2 measurements (Figures 5 and 6) however suggests that these biases are not significantly affecting our results.

Unquantified uncertainties in sediment thickness estimates
The greatest source of uncertainty in our final empirical relations is from the sediment thicknesses estimated from boreholes (Figure 7) and how this is mapped to a sediment thickness at a seismic station (Figure S8) for later calibration.The sediment thickness estimated at boreholes depends on the completeness of the various databases combined in this study (Section 3).The uncertainties in the databases themselves (e.g. the thickness of each individual stratigraphic unit, the interpretation of units as sediment etc.) are difficult to quantify.The ABSUC itself contains multiple sources of information for the same borehole, from which a preferred interpretation needs to be picked.Mapping the borehole sediment thickness estimates to seismic stations then introduces further uncertainty.We consider a smaller acceptable radius to identify nearby boreholes than previous work (0.25 • compared to 0.5 • in Agrawal et al. (2022)), but there may still be significant geological differences between the borehole site and the seismic station site within this range.The combination of these uncertainties results in a large spread (standard deviation of ∼ 1000 m) of thickness estimates as a function of t P s b (Figures 8 and S9).Finally, none of our borehole sediment thickness estimates reach the 10 km thick sediments expected in some Australian regions (Yeates et al. 1984;Myers et al. 1996).This points to a general limitation of borehole data, where it is difficult to ascertain if the borehole has indeed drilled all the way down to the crystalline basement.

Future Work
New continent-wide and regional deployments of seismic stations in Australia (e.g.Gorbatov et al. 2020;Murdie et al. 2022;Miller et al. 2023) will be excellent opportunities to use the empirical relations to obtain new sediment thickness estimates in some of the most remote areas of Australia, where few boreholes have been drilled.They may also provide opportunities for the empirical relations to be refined by filling in regional gaps.The continental scale nature of our study, covering rocks of a wide range of ages, also offers an exciting opportunity to potentially generalise these empirical relations to other sedimentary regions around the world, however this does remain to be fully explored.A useful future exercise that tests our hypotheses would be to undertake a blind application of the empirical relations for receiver functions obtained across another continental array, such as the US Transportable Array (IRIS Transportable Array 2003), followed by a comparison with sediment thicknesses estimates from other geophysical methods or borehole compilations.Such an experiment would be valuable in order to validate the relations for future global application.

CONCLUSIONS
We have obtained teleseismic receiver functions at ∼ 1, 500 seismic stations across Australia and measured the delay time of the P -to-S conversion at the sediment-crystalline basement interface.Regions of greater delay time illuminate the sedimentary basins, particularly the Cenozoic and Mesozoic basins in eastern Australia.We find that little to no delay is measured in older Proterozoic basins, due to the additional compaction that these sediments have experienced resulting in a smaller impedance contrast with the underlying basement.By calibrating these delay measurements with sediment thickness estimates from boreholes, we have produced an empirical relationship to establish a first order estimate of sediment thickness in a relatively straight-forward and cost-effective manner.This is a Mesozoic and Cenozoic basins to the east such as the Eromanga Basin over central and

Figure 1 .
Figure 1.(left) Geological map of Australia showing sedimentary provinces coloured by geological period (Raymond et al. 2018).Grey regions indicate Proterozoic basins.White regions are nonsedimentary.Names of major basins are indicated.(right) Map of 3-component seismic stations used in this study.

Figure 2 .
Figure 2. (top) Diagram illustrating the arrival of the direct P phase and the P -to-S converted phase from the basement-sediment interface (P s b ) at a seismic station sitting on a sedimentary layer.(middle) Vertical and radial receiver functions, and measurement of t P s b .All the energy on the vertical component is from the direct P phase.The P s b phase is delayed by t P s b seconds relative to the direct P phase.(bottom) Autocorrelation of the radial receiver function above.The time lag of the first local minimum corresponds to the two-way travel time of a shear wave in the sedimentary layer (t 2 ).

OFigure 3 .
Figure 3. Forward modelling of receiver function delay time and two-way travel time for varying sediment thickness.Forward modelling of the receiver functions is performed using the BayHunter package (Dreiling & Tilmann 2019).Sediment shear wave velocities (left) are modelled using a compaction model from Cunningham & Lekic (2020), with compressional wave velocities obtained using the empirical Brocher relations (Brocher 2005).Colours correspond to sediment thickness.(top right) The radial receiver function and it's autocorrelation for the corresponding shear velocity profiles.(bottom right) Delay time (t P s b , orange) and two-way travel time (t 2 , green) as a function of sediment thickness.Also shown in grey are the lines for 2t P s b and 3t P s b which provide bounds on t 2 .

Figure 4 .
Figure 4. Subset of stacked receiver functions, ordered and coloured by delay time t P s b .All receiver functions are normalised to unit amplitude.Labels are in network.stationformat.The subset is chosen as every 10 th stack from the full population of stacks ordered by t P s b , thereby producing a subset that is representative both in terms of the character of the receiver functions and relative frequency of t P s b values.

Figure 5 .O
Figure 5. Delay times of the P-to-S sediment converted phase relative to the direct P arrival (t P s b ).Grey lines outline major geological provinces with initials indicating sedimentary basins (see Figure 1).

Figure 6 .
Figure 6.Two-way travel time of shear waves in the sedimentary layer (t 2 ).

O
thickness from borehole databases The recent Australian Borehole Stratigraphic Units Compilation (ABSUC; Vizy & Rollet 2023, released by Geoscience Australia) provides a framework to consistently interpret borehole datasets spread across Australia.From 47 different data sources, the ABSUC reports the preferred interpretation of stratigraphic units traversed by 171,367 boreholes, giving the Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024

O
Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 Australian sediment thickness from receiver functions 19 Figure 8 shows the fits for the four different borehole-based sediment thickness estimations.As previously mentioned, our favoured fit is with the maximum sediment thickness at nearby boreholes (blue).The best-fitting quadratic relation for this is D = 1369.64t 2 P s b − 545.52tP s b + 382.62,

OFigure 8 .
Figure 8. Empirical relation between t P s b and our four different estimates of sediment thickness based on borehole data.Error bars in sediment thickness correspond to one standard deviation.Error bars in t P s b correspond to bin widths.Our favoured relation is the more conservative maximum thickness estimate (blue).Fits to the unbinned measurements are shown in Figure S9.

Figure 9 .O
Figure 9. Empirical relations as a function of sedimentary basin age, showing Cenozoic (blue), Mesozoic (green) and Paleozoic (orange) eras.Fits to the unbinned measurements are shown in Figure S10.

OFigure 10 .
Figure 10.(left) Sediment thickness at all seismic stations on sedimentary basins across Australia as determined by our overall empirical relation (Equation 5).Empty triangles denote stations not located on a sedimentary basin.(right) OZSEEBASE (Geognotics 2021) depth-to-basement model interpolated from a range of industrial geophysical methods.Note the different range of the colour bars.

O
Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 et al. valuable new tool that will benefit future exploration work across Australia and potentially other parts of the world in sedimented regions.The relationship is quadratic, reflective of a sediment compaction model in which rigidity increases linearly with depth, and is applicable to Phanerozoic basins.The relationship has also been further refined for the three Phanerozoic eras, with sediment thickness varying faster with delay time for the Palaeozoic than for the Mesozoic and Cenozoic.Australian sediment thickness from receiver functions 29 eAustralia/hiperseis) package(Hassan et al. 2020).All processing codes are available at https://github.com/auggiemarignier/AusSedThick.Receiver function data for all networks are available atMarignier (2023) along with the data in TableS1in various computer-readable formats.All figures have been built using matplotlib (Hunter 2007) and PyGMT(Tian et al. 2023).
Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 (Yeck et al. 2013)t thickness from receiver functions 112.2 Ps b basement conversion delayThe delay time of the P s b basement-sediment conversion, t P s b , is measured from the stacked receiver function at each station.The incidence of seismic waves from distant earthquakes is sub vertical due to decreasing velocity approaching the surface and Snell's Law.With the radial component seismogram recording horizontal motions, no P energy is recorded and thus the first main arrival on the radial component is the P s b basement to sediment conversion(Yeck et al. 2013).The delay of this arrival relative to the P arrival on the vertical component receiver function is what we refer to as t P s b (Figure Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 Thicker sediments, over 2000 m are identified in the Eromanga, Otway, Gippsland, Surat and Perth basins, whereas the Murray Basin seems thin (less than 1000 m).The borehole coverage is extremely sparse in much of Western Australia and the Northern Territory.Some regions of rural New South Wales and Victoria, and all of Tasmania, have no boreholes in the database.As a result, some major sedimentary basins such as the Eucla Basin have no borehole depth estimates.This highlights one of the main motivations for this study.On the other hand, there are locations such as the Perth Basin in which thick sediments are observed in the borehole data, but we currently have very few seismic stations in this region.
Downloaded from https://academic.oup.com/gji/advance-article/doi/10.1093/gji/ggae070/7616930 by guest on 09 March 2024 In general, older rocks tend to show a steeper relation than younger rocks due to an increase in shear-wave velocity as expected due to greater compaction or metamorphism.The sediment thickness of younger (Cenozoic) basins is also better constrained, with the uncertainty of each binned t P s b measurement noticeably smaller than for other eras.It is also worth noting, over the range of t P s b measured in each era, the relation is near linear.