Geophysical Journal International Royal

Astronomical
Society

Geophys. J. Int. (2024) 237, 949–964 https://doi.org/10.1093/gji/ggae085

Advance Access publication 2024 March 05
GJI Seismology

Crustal structure of Borneo, Makassar Strait and Sulawesi from

ambient noise tomography

N. Heryandoko,1,2 A. D. Nugraha,3,4 Z. Zulfakriza,3,4 S. Rosalia,3 T. Yudistira,3

S. Rohadi,2 D. Daryono,2 P. Supendi,2,5 N. Nurpujiono,2 F. Yusuf,2 F. Fauzi,2
A. Lesmana,1 Y. M. Husni,1 B. S. Prayitno,1,2 R. Triyono,1,2 S. P. Adi,2 D. Karnawati,2
T. Greenfield,5 N. Rawlinson5 and S. Widiyantoro3,6
1 Geophysical Engineering Study Program, Faculty of Mining and Petroleum Engineering, Institut Teknologi Bandung, Bandung 40132, Indonesia. Email:
nova.heryandoko@bmkg.go.id

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

2 Badan Meteorologi Klimatologi dan Geofisika (BMKG), Jakarta 10610, Indonesia .
3 Global Geophysics Research Group, Faculty of Mining and Petroleum Engineering, Institut Teknologi Bandung, Bandung 40132, Indonesia. Email:

andridn104@gmail.com
4 Research Center for Disaster Mitigation (RCDM), Institut Teknologi Bandung, Bandung 40132, Indonesia, .
5 Bullard Laboratories, Department of Earth Sciences, University of Cambridge, Cambridge CB3 0EZ, United Kingdom
6 Faculty of Engineering, Maranatha Christian University, Bandung 40164, Indonesia .

Accepted 2024 March 4. Received 2024 February 15; in original form 2023 September 22

SUMMARY
Borneo and Sulawesi are two large islands separated by the Makassar Strait that lie within
the complex tectonic setting of central Indonesia. The seismic structure beneath this region is
poorly understood due to the limited data availability. In this study, we present Rayleigh wave
tomography results that illuminate the underlying crustal structure. Group velocity is retrieved
from dispersion analysis of Rayleigh waves extracted from the ambient noise field by cross-
correlating long-term recordings from 108 seismic stations over a period of 8 months. We then
produce a 3-D shear wave velocity model via a two-stage process in which group velocity
maps are computed across a range of periods and then sampled over a dense grid of points
to produce pseudo-dispersion curves; these dispersion curves are then separately inverted for
1-D shear wave velocity (Vs), with the resultant models combined and interpolated to form
a 3-D model. In this model, we observed up to ± 1.2 km s−1 lateral Vs heterogeneities as a
function of depth. Our models illuminate a strong low shear wave velocity (Vs) anomaly at
shallow depth (≤ 14 km) and a strong high Vs anomaly at depths of 20–30 km beneath the
North Makassar Strait. We inferred the sediment basement and Moho depth from our 3-D Vs
model based on iso-velocity constrained by the positive vertical gradient of the Vs models.
The broad and deep sedimentary basement at ∼14 ± 2 km depth beneath the North Makassar
Strait is floored by a shallow Moho at ∼22 ± 2 km depth, which is the thinnest crust in the
study area. To the east of this region, our model reveals a Moho depth of ∼45 ± 2 km beneath
Central Sulawesi, the thickest crust in our study area, which suggests crustal thickening since
the late Oligocene. Moreover, the presence of high near-surface Vs anomalies with only slight
changes of velocity with increasing depth in southwest Borneo close to Schwaner Mountain
confirm the existence of a crustal root beneath this region.
Key words: Tomography; Crustal imaging; Moho depth; Seismic interferometry; Seismic
noise.

the Indo-Australian plate with the Eurasian and Philippine plates

1 I N T RO D U C T I O N
has led to a complex tectonic setting in this region, which to date
Borneo and Sulawesi are located in the middle of the Indonesian has only seen limited studies into its deeper structure. Borneo is
Archipelago, where they are separated by the relatively narrow divided among three countries, namely Brunei and Malaysia in the
(∼200–300 km wide) Makassar Strait. The ongoing collision of north and Indonesia in the south. The Indonesian territory is called

C The Author(s) 2024. Published by Oxford University Press on behalf of The Royal Astronomical Society. This is an Open Access
article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0/), which
permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited.
949
950 N. Heryandoko et al.

Kalimantan (Fig. 1), which occupies nearly 75 per cent of the total et al. (2022) and Pranata et al. (2020), but its application remains
land area. Borneo is located at the front of a laterally extruded block limited in Borneo, Makassar Strait and Sulawesi, with the exception
of the Indo-Australian-SE Asian collision zone (Tapponnier et al. of a focused study in northern Borneo (Greenfield et al. 2022). The
1982; Hall et al. 2008) and is relatively stable, with few earthquakes aim of our study is to recover broad scale crustal structure through
and a small number of inactive Quaternary volcanoes (Hall et al. the determination of 3-D shear wave velocity in order to contribute
2008). Sulawesi, located to the east of Borneo, has four arms that to a better understanding of the complex tectonic setting of central
together form its characteristic K-shape. Sulawesi was assembled by Indonesia.
the collision of multiple continental fragments. For example, Cen-
tral and Southeast Sulawesi, Banggai-Sula, and Buton are believed
to have originated in the northern part of the Australian continent 2 D ATA
(Metcalfe 1988, 1990; Audley-Charles 1991; Surono and Bachri We use continuous data collected from a dense seismic network,
2002). Borneo and Sulawesi were part of a single land mass in the with an average station spacing of 58 km, distributed throughout
Late Mesozoic, but were separated during the Cenozoic by rifting in Borneo and Sulawesi (Fig. 2a). This network consists of broadband
the Makassar Strait (Katili 1978; Hamilton 1979; Hall et al 2009). STS-2 instruments from Indonesia’s permanent seismic network
The Makassar Strait has two basins, the North and South Makassar (network code IA) operated by the Meteorological, Climatologi-
Basins, which are separated by the Paternoster Platform (PP). The cal and Geophysical Agency (BMKG), four of Malaysia’s perma-

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

North Basin lies beneath a water depth of ∼2.5 km and is connected nent broad-band STS-2 stations (network code MY) operated by
to the Celebes Sea (CBS) to the north. The South Basin sits beneath the Malaysian Meteorological Service, and two temporary seis-
∼2 km of water and is connected to the shallow shelf area of the mic networks deployed in Borneo and Sulawesi between 2019 and
Java Sea to the south. On the western side of the North Makas- 2020 (Greenfield et al. 2018). Together, these networks provide 108
sar Basin is the Kutai Basin (KB), the largest and deepest basin three-component seismic stations distributed across Borneo and Su-
in Indonesia with a depth to basement of 15 km (Rose & Hartono lawesi. The temporary seismic networks, which feature intermediate
1978); it is also one of the richest hydrocarbon reserves in Indone- broad-band stations (Guralp 6TDs and Nanometric Trillium Com-
sia (Hall et al. 2009). The opening of Makassar Strait resulted from pacts) are especially important for improving the ray coverage of
rifting during the Middle Eocene (∼48 Ma, Guntoro 1999; Moss Borneo, which has fewer permanent seismic stations. For the CCF,
& Chambers 1999; Calvert & Hall 2007), which field observations we used eight months of continuous data recorded on the verti-
and seismic imaging have shown, produced a highly asymmetrical cal component, which allowed us to retrieve good quality Rayleigh
and extremely wide rift in the Late Eocene (Hall et al. 2009). wave Green’s functions.
On the eastern side of North Makassar Basin is the Palu-Koro
Fault, a major sinistral strike-slip fault-oriented north-northwest
(NNW)—south-southeast (SSE) across onshore Central Sulawesi 3 M AT E R I A L S A N D M E T H O D S
(CS) and into the CBS. This fault caused the Mw 7.5 earthquake on
2018 September 28, which unexpectedly produced a large tsunami ANT proceeds by first calculating the ambient noise CCFs between
that inflicted major damage and casualties in the city of Palu and station pairs. These can be considered to be the EGFs corresponding
surrounding areas (Omira et al. 2019). to Earth structure between the stations (Weaver & Lobkis 2001;
During the last few decades, the increased availability of high- Fichtner et al. 2008). To improve the EGF signal-to-noise ratio
quality seismic data has been instrumental for the application of (SNR), we stack daily CCFs and then measure the Rayleigh wave
seismic tomography, for example, Zheng et al. (2021), Greenfield group velocity dispersion by applying wavelet domain analysis (see
et al. (2022) and Chen et al. (2023). Seismic tomography is ideal Section 3.1). Using the dispersion measurements, we then construct
for mapping subsurface velocity changes, which can be used to a group velocity map for each period and extract pseudo-dispersion
make inferences about the structure and tectonic history of a study curves on a 1◦ × 1◦ grid. These pseudo-dispersion curves are then
region (Widiyantoro & Van Der Hilst 1997; Curtis et al. 1998; Laske inverted for 1-D shear wave velocity (Vs) at each gridpoint, with
et al. 2013; Greenfield et al. 2022). In this study, we use ambient a final 3-D shear wave velocity model of the crust obtained by
seismic noise tomography (ANT) to image the lithosphere beneath stitching all the separate 1-D models together via interpolation.
Borneo, Makassar Strait and Sulawesi. The basis of ANT is the
extraction of empirical Green’s Functions (EGF) from long-term
3.1 Dispersion curve measurement
recordings of ambient noise via cross-correlation (CCF, Shapiro
& Campilo 2004). The primary advantage of this method is that We use a modified version of the Python code NoisePy (Jiang et al.
we can obtain a high-resolution image of crustal structure using a 2020) to compute the dispersion curve between station pairs. The
reasonably dense seismic network across a large region, for example, code was designed specifically for large-scale ambient noise seis-
beneath the Australian continent—see Saygin & Kennett (2010, mology and provides routines for downloading data, processing and
2012) and Chen et al. (2023), without relying on earthquakes or dispersion analysis (Jiang et al. 2020). NoisePy implements both
active sources (Suemoto et al. 2020). Body wave information is the wavelet domain scheme as in Fichtner et al. (2008) and time–
typically challenging to extract from ambient noise, so most studies frequency analysis (FTAN; e.g. Bensen et al. 2007) for the group
are restricted to the use of surface waves, which can provide good velocity estimation (Jiang et al. 2020). As discussed in Fichtner
horizontal resolution, but depth resolution is limited. Moreover, the et al. (2008), the wavelet domain approach is equivalent to the
exploitation of group and phase velocity dispersion means that S- FTAN in the case that the width of the sliding window does not
wave velocity tends to be the final product of ANT, with P-wave depend on time or frequency. We applied a standard ambient noise
velocity being too poorly constrained by the data to permit retrieval. processing workflow (Bensen et al. 2007) as follows: the data were
Ambient noise and other surface wave imaging methods have been resampled into daily segments at a sampling rate of 2 Hz, before
widely used in the Indonesian region, for example, Zulfakriza et al. the instrument response, mean and trend were removed from the
(2014, 2020), Martha et al. (2017), Yudistira et al. (2021), Rosalia data; a bandpass filter between 0.02 and 1 Hz was then applied.
 Crustal structure of Borneo, Makassar Strait and Sulawesi 951

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 1. A simplified geological map of central Indonesia (data from Steinshouer et al. 1999); the red square in the insert map indicates the approximate
study area of Borneo and Sulawesi and the red dashed line in the main part of the figure is the horizontal boundary of the shear wave velocity (Vs) model. The
fault lines (thin black lines) are obtained from the Geological Agency of Indonesia (https://geologi.esdm.go.id).

Figure 2. (a) Distribution of seismic stations used in this study; the blue triangles show the location of Indonesia’s permanent stations, the yellow triangles
show the distribution of Malaysia’s permanent stations and the red triangles show the distribution of temporary stations in Borneo and Sulawesi. (b) Eight
months of stacked CCFs in the period band 8–50 s. From this figure, we can see clear Rayleigh wave EGFs emerging at all interstation distances.

This sampling rate was chosen to reduce computation time, since station S and R in the Fourier domain defined by
we only consider the ambient noise signal in the period band 8–50 s.
We used a combination of spectral whitening with a running mean
average (rma) across 10 sample points and waveform CCF between US∗ ( f ) UR ( f )
CCFSR ( f ) = , (1)
{|US ( f )|} {|UR | ( f )}
952 N. Heryandoko et al.

where US (f) is the Fourier transform of the seismogram at virtual CCFs or EGFs prior to dispersion curve analysis. To obtain quan-
source station S, ∗ denotes the complex conjugate, UR (f) is the titative uncertainty estimates of the dispersion measurements, we
Fourier transform of the seismogram at virtual receiver station R, | compute the standard deviation of each period extending from 8
| denotes the absolute value and {} represents a smoothing of the to 50 s with a sampling of 1 s. Fig. 4(c) shows an example of
spectrum using rma (Jiang et al. 2020). The CCFs were computed the dispersion data uncertainties calculated from the station air of
for each station pair with a lag-time of ± 800 s. After being trans- PKKI and TTSI. We can observe that the standard deviation tends
formed back into the time domain, we combined all the daily CCFs to increase with increasing period. Using this approach, we can cal-
from each station pair using a linear stacking method. Fig. 2(b) culate the dispersion data uncertainties for all station pairs. Fig. 4(d)
shows the stacked CCFs extracted from all interstation pairs, which shows the average standard deviations obtained from all interstation
reveals very clear EGFs across the entire distance range. pairs.
The dispersion plot for each interstations pair was obtained using
the wavelet domain analysis tool included in NoisePy. It uses a
Morlet wavelet with default parameters for frequency spacing, scale
of the wavelet and number of scales as defined in the pycwt package
3.3 Group velocity anomalies
(https://pycwt.readthedocs.io, Jiang et al. 2020). To provide high-
quality dispersion data, we first identified and picked the best EGFs We used the Fast-Marching Surface Tomography code (Rawlinson

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

by using the SNR calculated as the ratio of the maximum absolute & Sambridge 2005) to produce group velocity maps as a function
amplitude of the signal window and the root mean square (RMS) of period across the study area. In this step, we only consider the
of a noise window as suggested by Ritzwoller & Feng (2019). The dispersion data ranging from 8 to 45 s with an increment of 1 s to pro-
signal window is taken between the minimum time, tmin , and the duce group velocity maps, since our data have larger uncertainties
maximum time, tmax , where tmin and tmax are defined by 4.5 and at periods above 45 s. We estimated the traveltime between station
1.5 km s−1 , respectively, for a given interstation distance. The noise pairs by picking the corresponding interstation dispersion curve as
window is taken between tmax and the end of the time-series. We a function of period and converting path-average group velocity to
used an SNR threshold of five, as used in a previous ANT study time. To obtain a starting model for each period, we computed the
(Zheng et al. 2021). This threshold of five was chosen to satisfy the uniform velocity field that best satisfied the data set. The resulting
cut-off between the number and the quality of EGF data. Using this starting model exhibits an approximately linear increase in velocity
threshold, the number of EGFs is reduced by about 40 per cent from with period, from 2.7 km s−1 at 8 s period to 3.6 km s−1 at 45 s
the total EGFs obtained from all the station pairs; however, this step period. The inversion was conducted for the target area spanning
leaves only high-quality EGFs, with in this case 1700 subsequently Borneo, the Makassar Strait and Sulawesi. To test the ability of our
used in the dispersion curve analysis. Finally, the dispersion curves data set to recover features of a certain size, evaluate which regions
were extracted from the maximum amplitude of the group velocity are well-recovered and choose optimal grid spacing, we perform
measurements at each period for each station pair. The quality of checkerboard tests at all periods before running the inversion using
dispersion data obtained in these measurements were then assessed, real data. Fig. 5 shows example checkerboard tests at 10, 20, 30 and
as discussed in Section 3.2. 40 s period using the station-pair distribution shown in Fig. 5 (bot-
tom). We test checkerboard anomalies with 0.5◦ (Fig. 3, top) and
1◦ grid spacing (Fig. 5, middle). From the results of these tests, the
1◦ checkerboard anomalies are recovered across most of the study
3.2 Dispersion measurement uncertainties
area, while the 0.5◦ checkerboard anomalies were only recovered
We perform quality checks on the interstation dispersion measure- across Sulawesi. Consequently, we used a 1◦ grid for the inversion
ments prior to group velocity analysis in order to weight the data in of our observed data set, while acknowledging that some finer scale
the inversion step. Here, we examine the spatial and temporal vari- structure in Sulawesi that may be constrained by the data will not
ability of the dispersion curves, for example, Bensen et al. (2007), be recovered. We perform trade-off tests to determine appropriate
as summarized in Figs 3 and 4. Fig. 3 shows dispersion curves damping and smoothing parameters for the real-data inversion. A
extracted from four interstation paths, all with similar geometries, damping value of 15 was found to provide a suitable balance be-
suggestive of similar subsurface properties and standard deviations tween satisfying the data and producing models that do not differ
ranging from 0.01 to 0.1 km s−1 (shown by the blue-shaded area) significantly from the starting model. A smoothing value of 3 was
across all periods. However, this assessment is only valid when chosen to provide smooth models that still exhibit a satisfactory
a cluster of stations subtends a small angle to a relatively distant reduction in data variance (see Fig. S3, Supporting Information).
station, which only occurs for a subset of measurements (Bensen Although there is some variation in path coverage between different
et al. 2007). Thus, we estimate the dispersion uncertainties quanti- periods, the chosen values of damping and smoothing were found
tatively using the temporal variability of dispersion measurements. to be suitable across the period range that was used.
To do so, we create three separate stacks that comprise the first and Following the procedure outlined by Barmin et al. (2001), we
second four months of stacked CCFs and the full eight months of identify and remove the traveltime outliers with an absolute resid-
stacked CCFs (Fig. 4a). Fig. 4(b) shows that the dispersion curve ual between the observed and predicted traveltime greater than two
generated from the second four months of stacked CCFs, which times the standard deviation. Our final Rayleigh wave group ve-
has lower SNR (with an SNR of 4.6, as shown by the red line), locity model for each period was obtained after 12 iterations and
leads to a significant deviation from the full eight month disper- the average RMS traveltime misfit reduction is about 90 per cent.
sion curve. However, the dispersion curve obtained from higher Fig. 6(a) shows the number of paths at each period after removing
SNR data recorded in the first four months (with an SNR of 7.2, as the traveltime outliers. Examples of the resultant group velocity map
shown by the blue line) is observed to be quite similar to the curve at each period are shown in Fig. 7. From these figures, prominent
obtained from the eight months of stacked CCFs. This shows that low-velocity anomalies are clearly revealed at periods between 10
the SNR check is suitable for assessing the quality of the stacked to 20 s. Furthermore, prominent high-velocity anomalies at longer
 Crustal structure of Borneo, Makassar Strait and Sulawesi 953

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 3. Example of the spatial variability of dispersion curve measurements. The inset map shows the cluster of four interstation paths; station PKKI in
Borneo and four stations in Sulawesi (TTSI, PMSI, SRSI and BKSI) used in this analysis. The black lines show the dispersion curves of four interstation
pairs. The blue-shaded area shows the standard deviation of the mean (blue dashed line) as a function of period. Additional examples of the stacked CCFs and
dispersion curve measurements at a number of station pairs can be seen in Fig. S1 in the Supporting Information.

Figure 4. Example of the temporal variability of dispersion curve measurements. (a) The stacked CCFs from all eight months (top), first four months (middle)
and the second four months (bottom) obtained from interstation pair PKKI and TTSI, which are separated by 656 km. (b) The dispersion measurements from
the three time periods with the same colour schema as in (a). (c) The uncertainty of the dispersion measurements calculated using the standard deviation of the
three curves in (b) as a function of period. (d) The average uncertainty of all dispersion measurements obtained from every interstation pair used in this study.
954 N. Heryandoko et al.

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 5. Checkerboard resolution test at 10, 20, 30 and 40 s period. Top: the recovered checkerboard with 0.5◦ × 0.5◦ grid spacing. Middle: the recovered
checkerboard with 1◦ × 1o grid spacing. The black dashed line represents the boundary of the region in which Vs structure was inverted for. Bottom: interstation
paths for all stations that yielded EGFs with an SNR > 5.0 [see Fig. S2, Supporting Information for other checkerboard resolution tests (1◦ × 1◦ ) at different
periods].

periods (periods ≥ 20 s) are clearly shown, that is, in the North Uobs and calculated group velocity, Ucal , is given by:
Makassar Strait (NMS), PP and south Borneo.
In addition, we also performed the inversion using a finer grid of N
∂U
Uobs − Ucal = Vsi , (2)
0.5◦ . The group velocity maps produced using the finer grid (0.5◦ ) ∂ Vsi
i=1
are shown in Fig. S4 (Supporting Information), which reveals some
of the extra details that can be recovered. However, our primary where N is the number of layers, Vsi is updated model from the
interest is to interpret similar-scale features across the entire study prior and the partial derivative of group velocity with respect to
region, and therefore prefer to use the coarser grid in the inversion shear wave velocity are derived by following Rodi et al. (1975):
to produce a single model that is everywhere well resolved.  
∂U U U ∂c U 2 ∂ ∂c
= 2− − 2T , (3)
∂ Vs c c ∂ Vs c ∂ T ∂ Vs
3.4 Shear wave velocity inversion where c is the phase velocity and T is the period. When group
To obtain a 1-D depth-dependent Vs model, we extract pseudo- velocities are computed, the phase velocities are also computed
dispersion curves on a 1◦ grid spacing covering the study area. 116 at the two periods (1 ± h)T. The required phase velocity partial
of 340 gridpoints were selected based on well-recovered areas from derivative can be calculated as follows:
checkerboard tests (see Fig. 5). The pseudo-dispersion curves at ⎛    ⎞
∂c
∂ ∂c ∂ Vs
T + hT − ∂∂cVs T − hT
the selected gridpoints (Fig. 6b) were then used as inputs for shear =T⎝ ⎠. (4)
wave velocity inversions. Strong low-velocity anomalies are clearly ∂ T ∂ Vs 2hT
shown in the period ranges of 10–15 s (Fig. 6b). These anoma-
lies are robust because the group velocity maps (Fig. 7) are well The parameter h is given a value of 0.005. By introducing the
constrained across these periods, as explained in Section 3.3, and regularization parameters of weight, W , and damping, σ , to eq. (2),
also contain distinct low-velocity anomalies. These curves are then the updated model from the prior, Vs , can be retrieved using the
inverted using an iterative least-squares technique contained in the generalized inverse as follows:
package surf96 (Herrmann & Ammon 2002; Herrmann 2013). Ac-
−1
cording to Herrmann & Ammon (2002), the difference in observed, Vs = W −1 V λ2 + σ 2 I λU T U, (5)
 Crustal structure of Borneo, Makassar Strait and Sulawesi 955

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 6. (a) The histogram of interstation path number at each period. (b) The pseudo-dispersion curves at 116 selected gridpoints (red circle in inset map)
obtained from group velocity maps (Fig. 7).

where the U , λ and V are decomposition matrices of group velocity synthetic Vs model (plotted as green line in Fig. 8b). The synthetic
partial derivatives, ∂ U/∂ Vs . Vs model was created using a simple 1-D Vs model obtained from
The prior model in each case was parametrized as a depth de- ak135 (Kennet et al. 1995). We chose an optimum damping value to
pendent seismic structure from the surface to 60 km depth, with obtain the best 1-D Vs model (see Fig. S5, Supporting Information).
homogeneous layers 2 km thick from the surface to 40 km depth We also perturb the synthetic data 1000 times by adding Gaussian
and 5 km thick below 40 km depth, since our data do not have much noise with a standard deviation ranging from 0.14 km s−1 at 8 s
sensitivity below 40 km. We use a monotonic Vs constraint as a period to 0.24 km s−1 at 45 s period, and then invert each synthetic
prior model by setting Vs to linearly increase with increasing depth pseudo-dispersion curve using the same parametrization as before
(see blue lines in Figs 8b and d) as also done in, for example, Feng to obtain the Vs model. As shown in Fig. 8(b), the final 1-D Vs model
& Ritzwoller (2019), Feng (2021a, b) and Feng & Diaz (2023). The (red line) is similar to the average Vs model taken from the ensemble
Vs in the crust increases from 3 km s−1 at 2 km depth to 4.0 km s−1 of 1000 models, indicating that our results are robust. In addition,
at 36 km depth, which we assume to be the Moho. Below the Moho, the grey area shows the standard error of these models, which varies
Vs increases from 4.0 km s−1 at 36 km depth to 4.2 km s−1 at 60 km from 0.03 km s−1 at shallow depths to 0.1 km s−1 at greater depths
depth. We use Brocher (2005) regression to estimate the prior Vp and below the Moho. However, as shown in Fig. 8(b), the velocity
model and the density of each layer. jump across the mid-crust (Conrad discontinuity) and the Moho, as
To estimate the uncertainties and to obtain an appropriate damp- shown in the synthetic model, cannot be resolved by the inversions,
ing value, we perform an inversion of the synthetic group velocity even though the dispersion misfits are relatively small, as shown in
(plotted as black dots with error bars in Fig. 8a) computed from the Fig. 8(a). Previous studies of the Moho using seismic surface waves
956 N. Heryandoko et al.

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 7. Group velocity distribution maps from 10 to 35 s obtained in this study. The prominent low anomaly at periods 10–20 s can be most clearly observed
in the NMS.

reported that Vs inversions suffer from a trade-off between the Moho 15 km), which may indicate that the real Moho could be shallower
depth and the uppermost mantle Vs (Lebedev et al. 2013). Here, we than the CRUST1.0 Moho at this location.
perform an inversion of the synthetic group velocity that explicitly
includes the Moho by applying a small weight in the neighbourhood
of the Moho to permit significant change in the Vs model. Both the 3.5 Crustal interface measurements
depth and velocity contrast at the discontinuity can be subsequently Once the 1-D Vs models have been obtained—as described above—
modified during the inversion. Similar to Figs 8(a)–(c) show the syn- they are combined into a pseudo-3-D model using a smooth inter-
thetic (black dot with error bar) and calculated dispersion curves polator. To delineate the primary interfaces beneath a study area,
(red line) while Fig. 8(d) shows the Vs model (red line) produced Chen et al. (2023) suggested using the gradient maxima to deter-
by the inversion of synthetic group velocity with the reduced regu- mine the sediment–basement interface and the Moho depth. This
larization constraint at 30 km depth. The velocity jump at the given method assumes a ‘typical’ crust that features a distinct sediment–
Moho depth can therefore be recovered by this inversion scheme basement interface and Moho, as seen in Fig. 9 at points A and C.
and the error reduced (see the grey lines in Fig. 8d). Since we use Moho depths from the CRUST1.0 model, the location
For the real data inversions, we use Moho depth information of velocity jumps is variable in depth. However, it is difficult to
from the CRUST1.0 model (Laske et al. 2013) at 116 selected find the gradient maxima for gradients that are not pronounced,
gridpoints covering the study area. Fig. 9 shows an example of the for example, at point B. As a consequence, we use a different
Vs models obtained from the inversion using Moho depth constraints approach based on prescribed velocity contours—in this case, a
of CRUST1.0 and Vs gradients. The gradient is calculated from each value of 3.0 ± 0.1 km s−1 for the sediment–basement interface
Vs solution model. The impulsive gradients at the Moho are clearly and 3.9 ± 0.1 km s−1 for the Moho provided the velocity gradient
shown at points A and C, indicating velocity jumps at this depth. is positive (velocity increases with depth). Mooney et al. (1998)
However, at coordinate point B, a weak gradient occurs at the Moho, produced a P-wave velocity (Vp) classification for sediment lay-
while a rapid velocity change occurs at shallower depths (around ers: 2.0–3.0 km s−1 for soft sedimentary layers on continents and
 Crustal structure of Borneo, Makassar Strait and Sulawesi 957

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 8. The Vs inversion tests using synthetic group velocity dispersion. (a) The black dots with error bars represent the synthetic group velocity at periods
ranging from 8 to 45 s computed from the synthetic velocity model plotted as a green line in (b). The red line is the calculated dispersion curve computed from
the Vs model plotted as a red line in (b). (b) The green line is the velocity model used to generate synthetic group velocity as a function of period. The red line is
the result of the shear wave (Vs) inversion as a function of depth. The grey lines are produced by 1000 separate inversions that include random Gaussian noise
added to the input dispersion curves with a standard deviation at each period equal to that obtained from the group velocity inversion. The dashed blue line is
the starting model for the inversion. (c) and (d) Similar in (a) and (b); here, the Vs model (red line in d) was obtained using an inversion with regularization
reduced in the neighbourhood of the Moho at 30 km depth. The sensitivity kernels for these inversion tests can be seen in Fig. S6 (Supporting Information).

4.0–5.3 km s−1 for hard sedimentary layers. Using Brocher (2005) heterogeneities as a function of depth. The results illuminate sig-
regression, the equivalent for i-wave velocity (Vs) is 0.6–1.4 km s−1 nificant low and high Vs anomalies across the region. We observe a
for soft sediments, and 2.3–3.2 km s−1 for hard sediments. Based strong low near-surface Vs anomaly zone (LSVZ) at shallow depths
on these results, we thus assume that a Vs of 3 km s−1 is the velocity (≤ 14 km) beneath NMS as delineated by the red dashed lines in
of the sediment–basement interface. For the Moho, we assume that Fig. 10. Other LSVZs are also observed beneath Tomini Bay (TB),
a Vs of 3.9 km s−1 is the highest velocity in the crust—based on the Bone Gulf (BG, including the South Arm of Sulawesi, SS), the
ak135 model of Kennet et al. (1995)—and can therefore be used as eastern part of the Southeast Arm of Sulawesi (SES) and Sarawak
a proxy for the Moho. Both the basement interface and Moho can (northwest Borneo). The LSVZ beneath the North Makassar Strait
be clearly recovered in ∼98 per cent of our 1-D Vs models. is strongest at 8 km depth but can be observed down to a depth
of 14 km (Fig. 10). This prominent LSVZ extends from KB, East
Borneo, to the west coast of Sulawesi in an east–west direction
4 R E S U LT S and spans a distance of ∼300 km; it also extends in a north–south
We show horizontal slices through the recovered Vs model between direction from PP to Mangkalihat Peninsula (MP), a distance of
8 and 40 km depth in Fig. 10 and the corresponding errors. From ∼400 km. This is the broadest and deepest LSVZ observed in this
our horizontal slices, we observe up to ± 1.2 km s−1 of lateral Vs region and also has the lowest Vs of ∼2.2 km s−1 observed at 8 km
958 N. Heryandoko et al.

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 9. Examples of 1-D Vs models (red lines) sampled at particular points as shown in the map (top right). The black lines on the right are the Vs gradients
calculated from each Vs solution model. The dashed blue lines are the Moho depth obtained from CRUST1.0.

depth. In contrast, we observe a high near-surface velocity zone the study area. The average depth-to-basement estimated from our
(HSVZ) at 8 km depth (Fig. 10) in Southwest Borneo and CS. Vs model is ∼6 km, with a minimum of 2 km and maximum of
Moreover, we also observe high deep Vs anomaly zones (HDVZs) 14 km, while the average Moho depth is ∼36 km, with a minimum
at depths of 30 km (as indicated by the blue dashed lines in Fig. 10). of 24 km and maximum of 45 km.
These HDVZs tend to be co-located with LSVZs at shallower depth From the vertical Vs profiles, the crustal structure beneath the
(∼8 km), that is, in the North Makassar Strait, TB, SS and the SES, study area can be readily observed. The sediment, the crust and
indicating rapid velocity changes with increasing depth. The errors the lithospheric mantle layer are clearly depicted in our model,
in the model are highest in the North Makassar Strait where strong which illuminates strong LSVZs in the upper crustal layer and in
Vs anomalies are found. the sedimentary layer, that is, in the neighbourhood of the NMS,
Fig. 11 shows six vertical profiles of our Vs model through Bor- including the KB in the west and the MP in the north, as well as TB
neo, the Makassar Strait and Sulawesi, along with their associated and BG. This suggests the likely presence of sedimentary material
error. The black solid lines denote the iso-velocity obtained in this with a Vs ≤ 3.0 km s−1 (Mooney et al. 1998) beneath the NMS, the
study, while the red and black dashed lines reveal the sediment base- SS, the BG and TB.
ment and Moho respectively. The iso-velocities of 3.8–4.0 km s−1
are largely in the vicinity of the CRUST1.0 Moho model (plotted as
black dashed lines); hence, it appears reasonable to assume that the 5 DISCUSSION
iso-velocity contour of 3.9 ± 0.1 km s−1 is appropriate for approx-
The widespread LSVZs beneath Sarawak Basin (SB), NMS, SS,
imately locating the Moho from our Vs models. However, in some
BG, TB and the eastern margin of the SES indicates the presence
regions, our Moho model deviates significantly from CRUST1.0,
of extensive sedimentary basins across the study area. The broad
for example, in the NMS and in CS to the southwest of TB—see
and thick LSVZs observed by this study shows good agreement
Fig. 11 cross-section B-B and E-E . To delineate the geometry of
with previous geology and geophysical studies. For example, the
the crustal interfaces, we compute the 2-D and 3-D depth contour
LSVZ observed beneath Sarawak confirms the presence of SB,
profiles that correspond to the sediment basement and Moho depth
which has been discussed previously by, for example, Hamilton
as shown in Fig. 12. From this figure, we can infer the depth variation
(1979) and Moss (1998). The SB incorporates the Rajang–Crocker
of the sediment basement (Fig. 12a) and the Moho (Fig. 12b) across
Group of Sarawak, which is interpreted as remnant oceanic basin
 Crustal structure of Borneo, Makassar Strait and Sulawesi 959

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 10. Horizontal slices through the final absolute shear wave velocity (Vs) model at selected depths along with the average standard error (bottom right)
across all depths. The red dashed lines at 8 and 14 km depths are iso-velocity contours at 3 km s−1 and the blue dashed lines are iso-velocity contours at
3.9 km s−1 and outline the boundaries of high Vs anomalies. The labelled lines on the 20 km depth slice shows the cross-section lines along which vertical
profiles of absolute Vs structure are plotted in Fig. 11.

fill created by deposition of an extensive Santonian–Palaeocene The most broad and deep sedimentary basement extending to a
(∼86–59 Ma) submarine fan on a passive margin or subduction depth of ∼14 ± 2 km can be observed in the NMS (Fig. 12a). The
accretionary complexes (Hamilton 1979; Moss 1998). The other NMS Basin, which includes the KB, has an asymmetrical shape
example, the LSVZs in SS and BG, confirms the findings of previous and extends laterally by ∼300 km in the east–west direction and
geology studies that reported the BG as a Neogene interarm basin ∼400 km in the north–south direction. This basin is bounded by
(∼23–4 Ma) dominated by extension and underlain by a variety Central Borneo Mountain (CBM) in the west, CS in the east, PP
of pre-Neogene rocks. Proven petroleum plays reported in the BG in the south, and the southwestern margin of the CBS in the north-
include gas fields in the East Sengkang Basin in the SS (Camplin west. The widespread LSVZ (Vs < 3.0 km s−1 ) is underlain by a
& Hall 2014). Our models also show widespread LSVZs in the deep sediment–basement interface at ∼14 ± 2 km depth (Fig. 12a)
eastern margin of the SES. Previous geology studies reported oil and is co-located with a shallowing of the Moho (Fig. 12b). This in-
shows and discoveries in and around the island of Buton at the tip of dicates the presence of thin crust beneath the NMS, likely resulting
the SE Arm (Camplin & Hall 2014). Oil seeps are reported offshore from rifting during the Middle Eocene (∼48 Ma). The thinnest crust
near Kolaka (Thompson et al. 1991). Furthermore, the LSVZ in the beneath the NMS is underlain by a Moho at ∼22 ± 2 km depth,
TB remains enigmatic and less studied; however, based on marine which is slightly thicker than the 20 km elastic thickness estimated
magnetic studies, the TB is underlain by an oceanic-like crust with by Cloke et al. (1999), who assumed that the crust beneath the NMS
suspected thermal extension, thinning and differential subsidence was 48 Myr-old oceanic crust. Previous gravity modelling (Guntoro
in the basement (Kusnida et al. 2009). According to Hall (2011), the 1999) showed that the lowest free-air anomaly of −40 mGal is in
additional extension which contributed to subsidence in this region the middle of the NMS basin and increases to + 70 mGal towards
was caused by subduction rollback at the North Sulawesi Trench. the continental shelves of Borneo and Sulawesi. Guntoro (1999)
960 N. Heryandoko et al.

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 11. Six vertical profiles of absolute Vs structure across eastern Borneo, Makassar Strait and Sulawesi, as shown by labelled lines in Fig. 10 on the 20 km
depth slice. The vertical resolution is 2 km, and the horizontal resolution is 1◦ . The topography and bathymetry profile and the region name are plotted at the
top of each velocity profile, and the seismicity along each profile is plotted as black circles. The red and the black dashed lines show the sediment thickness
and the Moho depth respectively, obtained from the CRUST1.0 model. The black thin lines show iso-velocity contours. The white lines in cross-section B-B
and C-C denote the SLAB2 model (Hayes et al. 2018; Haynie et al. 2020) of North Sulawesi subduction. The red lines that overlap with topography represent
the standard error of the resulting Vs model averaged over all depths, as also shown in Fig. 10 (bottom right).

postulated that there is a change in the thickness of the crust, due which also caused uplift. During the Miocene or earliest Pliocene,
to deformation by extensional thinning. Furthermore, most of Su- the Banggai–Sula continental fragment was accreted onto Eastern
lawesi is characterized by crust rarely thicker than 35 km (Fig. 12b). Sulawesi. A transpressional regime between 5 and 2 Ma then gener-
However, our models show a significant increase in crustal thick- ated further uplift in CS estimated to be 200–700 m Ma−1 (Wilson
ness from ∼22 ± 2 km beneath the NMS basin depocentre to ∼45 & Moss 1999). The vertical extent of LDVZ explains the thickened
± 2 km beneath CS. The crust to the west of the Makassar Strait crust beneath CS and the rapid uplift the region has experienced
is not as thick (30–35 km); this would suggest that the change in since the late Oligocene. The rapid uplift of CS and the tropical lati-
crustal thickness is not related to the extension of the NMS but is tudes have caused very high erosion rates and by extension the thick
instead a pre-existing feature. sedimentary package infilling the NMS and imaged in our velocity
As shown in Fig. 12, the deep Moho beneath CS is caused by model. Due to the limit of vertical resolution at greater depths in
low deeper Vs anomalies (LDVZ) that extend into the lithospheric our model, we cannot say whether the subducted oceanic crust is
mantle. To understand what this region of low velocities could repre- still attached. However, there is no indication in global, for example,
sent, we briefly describe how Sulawesi was assembled. As discussed LLNL-G3Dv3 (Simmons et al. 2012), or regional (Zenonos et al.
earlier, the island of Sulawesi is a juxtaposition of continental frag- 2019; Wehner et al. 2022a,b) tomography models of a subducting
ments from both Laurasia and Gondwana. The process that formed slab in this region. It is likely therefore that the slab has detached,
Sulawesi started in the Oligocene (∼30 Mya) with the continen- further increasing uplift rates in CS (van Hunen & Allen 2011).
tal fragment of Banggai–Sula drifting westwards toward Sulawesi, The uplift and gravitational collapse of Sulawesi is still ongoing
driven by westward subduction beneath the East Sulawesi Ophiolite and driven by the transpressional stress regime and internal defor-
(which comprises ES and SES). At this time, East Sulawesi and mation along major strike-slip faults. The broad westward vergent
West Sulawesi (comprising SS and the western part of CS) were fold and thrust belt in offshore West Sulawesi has been documented
not connected (Wilson & Moss 1999). During the late Oligocene, by Sandwell & Smith (2009) and Puspita et al. (2005). The fold
the East Sulawesi Ophiolite was accreted onto West Sulawesi. This and thrust belt is evidence of shelf failures caused by major uplift
terminated the rifting of the Makassar Stait (Guntoro 1999) and sed- and high exhumation (Hall 2011) and are still active and capable
iments began to accumulate in the subsiding basins. The supply of of producing large magnitude earthquakes (Greenfield et al. 2021;
sediments was increased by the collision of East and West Sulawesi, Meilano et al. 2023). Another example of a vertically extensive
 Crustal structure of Borneo, Makassar Strait and Sulawesi 961

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

Figure 12. The 2-D and 3-D depth contour profile of (a) the sedimentary basement and (b) the Moho model, obtained in this study.

LDVZ, can be found beneath the CBM, which explains the thick- Hall 2012). The SM itself is characterized by Cretaceous granitoid
ened crust in the region. Previous geology studies, for example, plutons intruding the older metamorphic rocks (Hall et al. 2008).
Hall et al. (2008) and Hall (2011), suggested that CBM crust is
mechanically weak and has experienced shortening and uplift in 6 C O N C LU S I O N S
the Early Miocene as a result of compressional forces applied by
adjacent regions of stronger crust. We determined the crustal shear wave velocity structure of Borneo,
The presence of a shallow Moho or thin crust suggests the lack Makassar Strait and Sulawesi using Rayleigh wave group velocity
of a crustal root as discussed by Szanyi et al. (2021) and Simonová ANT. The primary conclusions of our study are as follows:
et al. (2019). A previous ANT study in the Eastern Alps presented (i) We observed lateral Vs heterogeneities of up to ± 1.2 km s−1
by Szanyi et al. (2021) found that higher velocities were observed in the crust and uppermost mantle. The major sedimentary basins in
at shallower depths beneath basin areas, suggesting the presence the study area are well depicted from our Vs models, for example,
of thin crust and hence the absence of a crustal root. Our results in the NMS, TB, BG (including the SS), the eastern part of the SES
beneath NMS, TB and BG show similar crustal characteristics, and in Sarawak (northwest Borneo).
suggesting the lack of a crustal root beneath these locations (see (ii) Our models illuminate a strong low shear wave velocity (Vs)
Fig. 11). In particular, we observe low near-surface Vs anomalies anomaly at shallow depth (≤ 14 km) and a strong high Vs anomaly
(Vs < 3.0 km s−1 ) pointing to the presence of extensive sedimen- at depths of 22–30 km beneath NMS.
tary basins, which are underlain by a strong positive Vs gradient, (iii) We inferred the depth of the sediment basement and Moho
suggesting the presence of a shallow Moho and hence a thin crust. from our 3-D Vs model using iso-velocity contours along with a
Moreover, Szanyi et al. (2021) also suggest that the ∼homogeneous requirement for a positive vertical gradient in Vs. The broad and
structure embedded in the higher velocity mantle material, which deep sedimentary basement at ∼14 ± 2 km depth beneath the NMS
is characterized by slight changes of velocity with depth, indicates is floored by a shallow Moho at ∼22 ± 2 km depth, which is the
the presence of a crustal root. From our results, the presence of thinnest crust in the study area.
HSVZs in southwest Borneo, close to Schwaner Mountain (SM), (iv) Our models reveal a significant change in crustal thickness
confirms the existence of a crustal root, since this region is also to the east, from ∼22 ± 2 km beneath the NMS basin depocentre
characterized by slight changes in velocity with increasing depth to ∼45 ± 2 km beneath CS, which likely reflects crustal thickening
(see Fig. 11). Borneo has a Palaeozoic continental core in its west since the late Oligocene.
to southwest regions, surrounded by ophiolitic, island arc and mi- (v) The presence of high near-surface Vs anomalies and only
crocontinental crust accreted during the Mesozoic between 252 and slight changes in velocity with increasing depth in southwest Borneo
72 Ma (Hamilton 1979), which is interpreted as a block derived close to SM confirms the presence of a crustal root beneath this
from the Gondwana margin of NW Australia (Hall et al. 2009; region.
962 N. Heryandoko et al.

S U P P O RT I N G I N F O R M AT I O N Calvert, S.J. & Hall, R.(2007). Cenozoic Evolution of the Lariang and
Karama regions, North Makassar Basin, western Sulawesi, Indonesia.
Supplementary data are available at GJIRAS online. Pet. Geosci., 13, 353–368. DOI: 10.1144/1354-079306-757.
Camplin, D.J. & Hall, R.(2014). Neogene history of Bone Gulf, Sulawesi, In-
donesia. Mar. Pet. Geol., 57, 88. DOI: 10.1016/j.marpetgeo.2014.04.014.
AC K N OW L E D G M E N T S Chen, Y., Saygin, E., Kennett, B., Qashqai, M.T., Hauser, J., Lumley, D. &
We thank two anonymous-peer-reviewers and the editor Prof Hua- Sandiford, M.(2023). Next-generation seismic model of the Australian
crust from synchronous and asynchronous ambient noise imaging. Nat.
jian Yao, for constructive comments that greatly improved the qual-
Commun., 14, 1192. DOI: 10.1038/s41467-023-36514-z.
ity of this manuscript. We are grateful to Meteorological Clima-
Cloke, I.R., Milsom, J. & Blundell, D.J.B.(1999). Implications of gravity data
tological and Geophysical Agency (BMKG) and Malaysian Me- from east Kalimantan and the Makassar Straits: a solution to the origin of
teorological Service (Ministry of Science, Tech and Innovation, the Makassar Straits? J. Asian Earth Sci., 17, 61–78. DOI: 10.1016/S0743-
MOSTI-MMD-MET) for providing raw seismic data from Borneo 9547(98)00056-7.
and Sulawesi that are available at https://geof .bmkg.go.id/webdc3 Crameri, F.(2018). Scientific colour maps. Zenodo. http://doi.org/10.5281/
and FDSN Web Services at IRIS (http://service.iris.edu). This work zenodo.1243862. Accessed July 18, 2023.
is supported by BMKG and Hibah Riset Institut Teknologi Ban- Curtis, A., Trampert, J., Snieder, R. & Dost, B.(1998). Eurasian fun-
dung (ITB, https://www.itb.ac.id) 2022, awarded to ADN, KE- damental mode surface wave phase velocities and their relation-

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

MENDIKBUDRISTEK (Penelitian Dasar Kemitraan), contract no: ship with tectonic structures. J. geophys. Res., 103, 26919–26947.
DOI: 10.1029/98JB00903.
083/E5/PG.02.00.PT/2022, and matching funds from NFIS awarded
Feng, L.(2021a). Amphibious shear wave structure beneath the Alaska-
to SW and British Council Newton Impact fund G107511. Most of
aleutian subduction zone from ambient noise tomography. Geochem.
the figures were plotted using the Generic Mapping Tools (Wes- Geophys. Geosyst., 22, e2020GC009438. DOI: 10.1029/2020GC009438.
sel et al. 2013) and most of the color figures were plotted using Feng, L.(2021b). High-resolution crustal and uppermost mantle struc-
scientific colour maps from Crameri (2018). ture beneath central Mongolia from Rayleigh waves and receiver
functions. J. geophys. Res.: Solid Earth, 126, e2020JB021161.
DOI: 10.1029/2020JB021161.
D ATA AVA I L A B I L I T Y Feng, L. & Diaz, J.(2023). A high-resolution shear velocity model of the
crust and uppermost mantle beneath westernmost Mediterranean includ-
The data sets used and analysed during the current study are :
ing radial anisotropy. J. geophys. Res.: Solid Earth 128, e2023JB026868.
(i) The waveform data of the Malaysian Seismic Network (ht DOI: 10.1029/2023JB026868.
tps://www.fdsn.org/networks/detail/MY) with station code; KKM, Feng, L. & Ritzwoller, M.H.(2019). A 3-D shear velocity model of
KSM, LDM and SBM are available from the Incorporated Research the crust and uppermost mantle beneath Alaska including apparent
radial anisotropy. J. geophys. Res.: Solid Earth 124, 10468–10497.
Institutions for Seismology (IRIS) Data Services (http://ds.iris.edu/
DOI: 10.1029/2019JB018122.
mda/MY). Fichtner, A., Kennett, B.L., Igel, H. & Bunge, H.-P.(2008), Theoreti-
(ii) The waveform data from Indonesian Seismic Network in cal background for continental-and global-scale full-waveform inver-
Borneo and Sulawesi is restricted access due to the institutional sion in the time–frequency domain. Geophys. J. Int., 175, 665–685.
regulations of BMKG. The data are available at https://geof .bmkg. DOI: 10.1111/j.1365-246X.2008.03923.x.
go.id/webdc3, with permission from BMKG required to access the Greenfield, T., Copley, A.C., Caplan, C., Supendi, P., Widiyantoro,
data. S. & Rawlinson, N.(2021). Crustal deformation and fault strength
(iii) Several IA stations are also available in the GEOFON Pro- of the Sulawesi Subduction Zone. Tectonics, 40, e2020TC006573.
gram (GEOFON), (https://www.fdsn.org/networks/detail/IA) DOI: 10.1029/2020TC006573.
(iv) The waveform data from the temporary network deployed in Greenfield, T., Gilligan, A., Pilia, S., Cornwell, D.G., Tongkul, F., Widiyan-
toro, S. & Rawlinson, N.(2022). Post-subduction tectonics of Sabah,
Borneo (https://doi.org/10.7914/SN/9G 2018) and Sulawesi during
northern Borneo, inferred from surface wave tomography. Geophys. Res.
2019–2020 are restricted; requests for these data should be sent to
Lett., 49, e2021GL096117. DOI: 10.1029/2021GL096117.
TG or NR for Borneo data and to ZZ or SW for Sulawesi data. Greenfield, T., Widiyantoro, S. & Rawlinson, N.(2018). Kalimantan Tempo-
(v) The 3-D shear wave velocity (3-D Vs) model produced in rary Network [Data set]. International Federation of Digital Seismograph
this study and the data underlying this paper are available at https: Networks, DOI: 10.7914/SN/9G 2018. Accessed January 10, 2021.
//doi.org/10.5281/zenodo.10662743. Guntoro, A.(1999). The formation of the Makassar Strait and the separation
between SE Kalimantan and SW Sulawesi. J. Asian Earth Sci., 17, 79–98.
DOI: 10.1016/S0743-9547(98)00037-3.
Hall, R.(2011), Australia–SE Asia collision: plate tectonics and crustal
REFERENCES flow, the SE Asian gateway: history and tectonics of the Australia–Asia
Audley-Charles, M.G.(1991). Tectonics of the New Guinea area. collision. Geol. Soc. Lond. Spec. Publ., 355, 75–109.
Annu. Rev. Earth planet. Sci., 19, 17–41. DOI: 10.1146/an- DOI: 10.1144/SP355.5.
nurev.ea.19.050191.000313. Hall, R.(2012). Late Jurassic-cenozoic reconstructions of the Indone-
Barmin, M.P., Ritzwoller, M.H. & Levshin, A.L.(2001). A fast and reliable sian region and the Indian ocean. Tectonophysics, 570-571, 1–41.
method for surface wave tomography. Pure appl. Geophys., 158, 1351– DOI: 10.1016/j.tecto.2012.04.021.
1375. DOI: 10.1007/PL00001225. Hall, R., Cloke, I.R., Nur’aini, S., Puspita, S.D., Calvert, S.J. & Elders,
Bensen, G.D., Ritzwoller, M.H., Barmin, M.P., Levshin, A.L., Lin, F., C.F.(2009). The North Makassar Strait: what lies beneath? Pet. Geosci.,
Moschetti, M.P., Shapiro, N.M. & Yang, Y., 2007. Processing seismic 15, 147–158. DOI: 10.1144/1354-079309-829.
ambient noise data to obtain reliable broad-band surface wave dispersion Hall, R., Van, H., Marco, W.A. & Spakman, W.(2008). Impact of India–
measurements. Geophys. J. Int., 169, 1239–1260. DOI: 10.1111/j.1365- Asia collision on SE Asia: the record in Borneo. Tectonophysics,
246X.2007.03374.x. Asia out of Tethys: Geochronol. Tectonic Sediment. Records, 451,
Brocher, T.M.(2005). Empirical relations between elastic wavespeeds and 366–389
density in the earth’s crust. Bull. seism. Soc. Am., 95, 2081–2092. Hamilton, W.(1979). Tectonics of the Indonesian region. US Geol. Surv.
DOI: 10.1785/0120050077. Professional Paper, 1078. DOI: 10.3133/pp1078.
 Crustal structure of Borneo, Makassar Strait and Sulawesi 963

Hayes, G.P., Moore, G.L., Portner, D.E., Flamme, H., Furtney, Rawlinson, N. & Sambridge, M.(2005). The fast marching method: an effec-
M. & Smoczyk, G.M.(2018). Slab2, a comprehensive subduction tive tool for tomographic imaging and tracking multiple phases in complex
zone geometry model. Science, 362, 58–61. DOI: 10.1126/sci- layered Media. Explor. Geophys., 36, 341–350. DOI: 10.1071/EG05341.
ence.aat4723. Ritzwoller, M.H. & Feng, L.(2019). Overview of pre-and post-processing
Haynie, K.L., Hayes, G.P., Martinez, E., Fee, J., Las- of ambient noise correlations, in Seismic Ambient Noise, Cambridge Uni-
towka, L. & Guy, M.R.(2020). IN037-05–Cloud-based versity Press, 5, 144–187.
web application for updating USGS 3D subduction Rodi, W.L., Glover, P., Li, T.M.C. & Alexander, S.S.(1975). A fast,
zone geometry models. AGU Abstract, 2020. https://agu.conf ex.c accurate method for computing group-velocity partial derivatives for
om/agu/f m20/meetingapp.cgi/Paper/680941. Rayleigh and Love modes. Bull. seism. Soc. Am., 65, 1105–1114.
Herrmann, R.B.(2013). Computer Programs in seismology: an evolv- DOI: 10.1785/BSSA0650051105.
ing tool for instruction and research. Seismol. Res. Lett., 84, 1081. Rosalia, S., Widiyantoro, S., Cummins, P.R., Yudistira, T., Nugraha, A.D.,
DOI: 10.1785/0220110096. Zulfakriza, Z. & Setiawan, A.(2022). Upper crustal shear-wave velocity
Herrmann, R.B. & Ammon, C.J.(2002). Computer Programs in Seismology: structure beneath Western Java, Indonesia from seismic ambient noise
Surface Waves, Receiver Functions and Crustal Structure, Version 3.30. tomography. Geosci. Lett., 9, 1. DOI: 10.1186/s40562-021-00208-5.
https://www.eas.slu.edu/eqc/eqc cps/CPS/CPS330/cps330c.pdf. Rose, R. & Hartono, P.(1978). Geological evolution of the Tertiary
Jiang, C., Yuan, C. & Denolle, M.(2020). NoisePy: a new high-performance Kutei–Melawi Basin, Kalimantan, Indonesia. In: Proceedings Indonesian
python tool for seismic ambient noise seismology. Seismol. Res. Lett., 91, Petroleum Association, 7th Annual Convention, pp. 225–252

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

1853–1866. DOI: 10.1785/0220190364. Sandwell, D.T. & Smith, W.H.F.(2009). Global marine gravity from retracked
Katili, J.A.(1978). Past and present geotectonic position of Sulawesi, Geosat and ERS-1 altimetry: ridge segmentation v. spreading rate. J.
Indonesia. Tectonophysics, 45, 289–322. DOI: 10.1016/0040- geophys. Res., 114, B01411. DOi: 10.1029/2008JB006008.
1951(78)90166-X. Saygin, E. & Kennett, B.L.N.(2010). Ambient seismic noise to-
Kennett, B.L.N., Engdahl, E.R. & Buland, R.(1995). Constraints on seismic mography of Australian continent. Tectonophysics, 5, 116–125.
velocities in the earth from travel times. Geophys. J. Int, 122, 108–124. DOI: 10.1016/j.tecto.2008.11.013.
DOI: 10.1111/j.1365-246X.1995.tb03540.x. Saygin, E. & Kennett, B.L.N.(2012). Crustal structure of Australia from
Kusnida, D.S. & Nirwana, B. (2009). Basement configuration of the Tomini ambient seismic noise tomography. J. geophys. Res., 117, B01304.
Basin deduced from Marine Magnetic interpretation. J. Geol. Indonesia, DOI: 10.1029/2011JB008403.
4(4), 269–274. Shapiro, N.M. & Campillo, M.(2004). Emergence of broadband Rayleigh
Laske, G., Masters, G., Ma, Z. & Pasyanos, M.(2013). Update on CRUST1.0 waves from correlations of the ambient seismic noise. Geophys. Res. Lett.,
- A 1-degree Global Model of Earth’s Crust. Geophys. Res. Abstracts, 15, 31, L07614. DOI: 10.1029/2004GL019491.
Abstract EGU2013-2658. https://igppweb.ucsd.edu/∼gabi/crust1/laske-e Simmons, N.A., Myers, S.C., Johannesson, G. & Matzel, E.(2012). LLNL-
gu13-crust1.pdf. G3Dv3: global P wave tomography model for improved regional and
Lebedev, S., Adam, J.M.C. & Meier, T.(2013). Mapping the Moho teleseismic travel time prediction. J. geophys. Res., 117, B10302.
with seismic surface waves: a review, resolution analysis, and DOI: 10.1029/2012JB009525. .
recommended inversion strategies. Tectonophysics, 609, 377–394. Simonová, B., Zeyen, H. & Bielik, M., (2019). Continental lithospheric
DOI: 10.1016/j.tecto.2012.12.030. structure from the East European Craton to the Pannonian Basin based
Martha, A.A., Cummins, P., Saygin, E., Widiyantoro, S. & Mastury- on integrated geophysical modelling. Tectonophysics, 750, 289–300. .
ono. (2017). Imaging of upper crustal structure beneath East Java– Steinshouer, D.W., Qiang, J., McCabe, P.J. & Ryder, R.T.(1999). Maps show-
Bali, Indonesia with ambient noise tomography. Geosci. Lett., 4. ing geology, oil and gas fields, and geologic provinces of the Asia Pacific
DOI:10.1186/s40562-017-0080-9 region: U.S. Geological Survey open-file report 97-470-F. 16 p.,
Meilano, I., Salman, R., Susilo, S., Shiddiqi, H.A., Supendi, P., Lythgoe Suemoto, Y., Ikeda, T., Tsuji, T. & Lio, Y.(2020). Identification of a nascent
et al..(2023). The 2021 MW 6.2 Mamuju, West Sulawesi, Indonesia earth- tectonic boundary in the San in area, southwest Japan, using a 3D S wave
quake: partial rupture of the Makassar Strait thrust. Geophys. J. Int., 233, velocity structure obtained by ambient noise surface wave tomography.
1694–1707. DOI: 10.1093/gji/ggac512. Earth Planets Space, 72, 15. DOI: 10.1186/s40623-020-1139-y.
Metcalfe, I.(1988). Origin and assembly of Southeast Asian con- Surono & Bachri, S. (2002). Stratigraphy, sedimentation and palaeogeo-
tinental terranes. Geol. Soc. Lond. Spec. Publ., 37, 101–118. graphic significance of the Triassic Meluhu Formation, southeast arm of
DOI: 10.1144/GSL.SP.1988.037.01.08. Sulawesi, Eastern Indonesia. J. Asian Earth Sci., 20, 177–192
Metcalfe, I.(1990). Allochthonous terrane processes in Southeast Asia. Phi- Szanyi, G., Gráczer, Z., Balázs, B. & Kovács, I.J., AlpArray Working
los. Trans. R. Soc. Lond. A, 331, 625–640. Group(2021). The transition zone between the Eastern Alps and the Pan-
Mooney, W.D., Laske, G. & Masters, T.G.(1998). CRUST5.1: a nonian basin imaged by ambient noise tomography. Tectonophysics, 805,
global crustal model at 5o x 5o . J. geophys. Res., 103, 727–747. 228770. DOI: 10.1016/j.tecto.2021.228770.
DOI: 10.1029/97JB02122. Tapponnier, P., Peltzer, G., Le Dain, A.Y., Armijo, R. & Cobbold, P.(1982).
Moss, S.J.(1998). Embaluh Group turbidites in Kalimantan: evolution Propagating extrusion tectonics in Asia: new insights from simple ex-
of a remnant oceanic basin in Borneo during the late Cretaceous periments with plasticine. Geology, 10, 611–616. DOI: 10.1130/0091-
to Palaeogen. J. geol. Soc. Lond., 155, 509–524. 7613(1982)10611:PETIAN>2.0.CO;2.
DOI: 10.1144/gsjgs.155.3.0509. Thompson, M., Reminton, C., Purnomo, J. & Macgregor, D.(1991). Detec-
Moss, S.J. & Chambers, J.L.C.(1999). Tertiary facies architecture in the tion of liquid hydrocarbon seepage in Indonesian offshore frontier basins
Kutai Basin, Kalimantan, Indonesia. J. Asian Earth Sci., 17, 157–181. using Airborne Laser Fluorosensor (ALF); the results of a Pertamina/BP
DOI: 10.1016/S0743-9547(98)00035-X. joint study. In: Indonesian Petroleum Association, Proceedings 20th An-
Omira, R. et al. 2019. The September 28th, 2018, tsunami In Palu-Sulawesi, nual Convention,pp. 663–689
Indonesia: a post-event field survey. Pure appl. Geophys., 176, 1379– Van Hunen, J. & Allen, M.B.(2011). Continental collision and slab break-off:
1395. DOI: 10.1007/s00024-019-02145-z. a comparison of 3-D numerical models with observations. Earth planet.
Pranata, B. et al. 2020. Shear wave velocity structure beneath Bandung Sci. Lett., 302, 27–37. DOI: 10.1016/j.epsl.2010.11.035.
basin,West Java, Indonesia from ambient noise tomography. Geophys. J. Weaver, R.L. & Lobkis, O.I.(2001). Ultrasonics without a source: ther-
Int., 220, 1045–1054. mal fluctuation correlations at MHz frequencies. Phys. Rev. Lett., 87,
Puspita, S.D., Hall, R. & Elders, C.F.(2005). Structural styles of the off- DOI:10.1103/PhysRevLett.87.134301
shore West Sulawesi fold belt, North Makassar Straits, Indonesia. In: Wehner, D., Blom, N., Rawlinson, N., Daryono, D., Böhm, C., Miller, M.S.,
Proceedings Indonesian Petroleum Association, 30th Annual Convention, Supendi, P. & Widiyantoro, S.(2022a). SASSY21: a 3-D seismic structural
pp. 519–542. model of the lithosphere and underlying mantle beneath Southeast Asia
964 N. Heryandoko et al.

from multi-scale adjoint waveform tomography. J. geophys. Res.: Solid Yudistira, T. et al.(2021). Imaging of a magma system beneath the Merapi
Earth, 127, e2021JB022930. DOI: 10.1029/2021JB022930. Volcano complex, Indonesia, using ambient seismic noise tomography.
Wehner, D., Rawlinson, N., Greenfield, T., Daryono, D., Miller, M.S., Geophys. J. Int., 226, 511–523. DOI: 10.1093/gji/ggab104.
Supendi, P., Lü, C.C. & Widiyantoro, S.(2022b). SASSIER22: Zenonos, A., Siena, L.D., Widiyantoro, S. & Rawlinson, N.(2019). P and
full-waveform tomography of the eastern Indonesian region S wave travel time tomography of the SE Asia-Australia collision zone.
that includes topography, bathymetry, and the fluid ocean. Phys. Earth planet. Inter., 293, 106267. DOI: 10.1016/j.pepi.2019.05.010.
Geochem. Geophys. Geosyst., 23, e2022GC010563. Zheng, L., Fan, X., Zhang, P., Hao, J., Qian, H. & Zheng, T.(2021). Detection
DOI: 10.1029/2022GC010563. of urban hidden faults using group velocity ambient noise tomography
Wessel, P., Smith, W.H.F., Scharroo, R., Luis, J. & Wobbe, F.. (2013). Generic beneath Zhenjiang area, China. Sci. Rep., 11, 987. DOI: 10.1038/s41598-
mapping tools: improved version released. EOS Trans. Am. geophys. Un., 020-80249-6.
94, 409–410. DOI: 10.1002/2013EO450001. Zulfakriza, Z. et al.(2020). Tomographic imaging of the Agung-Batur Vol-
Widiyantoro, S. & Van Der Hilst, R.V. (1997). Mantle structure be- cano Complex, Bali, Indonesia, from the ambient seismic noise. Field.
neath Indonesia inferred from high-resolution tomographic imaging. Front. Earth Sci., 8, 43. DOI: 10.3389/feart.2020.00043.
Geophys. J. Int., 130, 167–182. DOI: 10.1111/j.1365- Zulfakriza, Z., Saygin, E., Cummins, P.R., Widiyantoro, S., Nugraha,
246X.1997.tb00996.x. A.D., Lühr, B.-G. & Bodin, T.. (2014). Upper crustal structure
Wilson, M.E.J. & Moss, S.J.(1999). Cenozoic palaeogeographic evolution of central Java, Indonesia, from transdimensional seismic ambient
of Sulawesi and Borneo. Palaeogeogr. Palaeoclimatol. Palaeoecol., 145, noise tomography. Geophys. J. Int., 197, 630–635.

Downloaded from https://academic.oup.com/gji/article/237/2/949/7623031 by guest on 14 May 2024

303–337. DOI: 10.1016/S0031-0182(98)00127-8. DOI: 10.1093/gji/ggu016.


C The Author(s) 2024. Published by Oxford University Press on behalf of The Royal Astronomical Society. This is an Open Access
article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0/), which
permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited.