Discriminating Lithospheric and Asthenospheric Anisotropy Beneath Northern Oman: Sharp Contrast Observed at the Semail Gap Fault Zone

To gain a deeper understanding of the extensive and varied lithospheric deformations beneath northern Oman, we examine seismic anisotropy in this region using splitting analysis of teleseismic shear wave data. Our study utilizes data from a dense network consisting of 13 permanent and 45 temporary seismic stations, which were operational for approximately 2.5 years starting from October 2013. By examining the azimuthal distribution of shear wave splitting (SWS) parameters, we were able to divide the study area into three sub‐regions. The stations located to the west of the Hawasina window exhibit relatively azimuthally invariant SWS parameters suggesting a single anisotropic layer. On the other hand, most of the stations located in the central and eastern regions display variations versus back‐azimuth, indicating the potential presence of depth‐dependent anisotropy. The General NW‐SE trend of the Fast Polarization Directions (FPDs) of the one‐layer anisotropy in the west and FPDs of the upper layers in the east is concordant with the strike of the structures resulting from the collision between the continental and oceanic plates. A clear contrast in SWS parameters is observed in the Semail Gap Fault Zone (SGFZ), suggesting that the SGFZ can be a lithospheric‐scale structure that hampers the intrusion of mafic magma from the southeast. Furthermore, the FPDs of the lower layer in the east exhibit an NE‐SW trend, which may be indicative of the large‐scale mantle flow resulting from the present‐day plate motion.


Introduction
The Semail ophiolite in northern Oman is a prime example of an obducted oceanic lithosphere.The structure of the Semail ophiolite and its tectonic evolution have been subject of geoscientific investigations for decades (e.g., Glennie, 1974;Searle & Malpas, 1980).Geological studies (e.g., Ambrose & Searle, 2019;Duretz et al., 2016;Guilmette et al., 2018;Hacker et al., 1996;Nicolas et al., 1994Nicolas et al., , 2000;;Porkoláb et al., 2021;Rioux et al., 2016;Searle & Cox, 1999) have provided insights into the chronological phases and the geological context underlying the formation of oceanic lithosphere, the dynamics of obduction, and the orogenic development in northern Oman.Previous studies provided a comprehensive understanding of tectonic evolution of the continental crust.However, the significance of pre-obduction tectonic processes in shaping the Arabian continental lithosphere, as well as their role in understanding obduction geometry and dynamics in detail, have been demonstrated only recently (e.g., Weidle et al., 2022Weidle et al., , 2023)).
Two major tectonic events had left their mark before the obduction of the Semail ophiolite, the extensional Neoproterozoic assembly of Gondwana followed by compressional foreland basin as a result of the Permo-Triassic breakdown of Pangea (Ali et al., 2013;Jabir et al., 2023).Both these events are reflected in today's lithospheric architecture in northeast Arabia.Previous surface wave tomography studies provided information about the lithospheric thickness and Moho depth of the study region (e.g., Pilia, Hu, et al., 2020;Pilia, Jackson, et al., 2020, 2021;Priestley et al., 2012;Weidle et al., 2022Weidle et al., , 2023)).Global and regional tomography models (e.g., Al-Lazki et al., 2014;Celli, Lebedev, Schaeffer, & Gaina, 2020;Celli, Lebedev, Schaeffer, Ravenna, & Gaina, 2020;Koulakov et al., 2016) suggest that mantle lithospheric thickness decreases from >200 km beneath the neocratonic core of Arabia to ∼100 km at the northern margins (Weidle et al., 2023).Furthermore, the crust northwest of the Semail Gap Fault Zone (SGFZ) is typically continental with an average Moho depth of ∼40 km, whereas it thins to the east and is likely modified by significant mafic reworking during the mafic breakup of Pangea (Weidle et al., 2023).Given the significant change of crustal architecture across the SGFZ in northern Oman, a main question arises whether similar variations occur in the mantle lithosphere.In the absence of detailed tomographic images of the area, observations from seismic anisotropy could provide more insights into structural properties of the continental lithosphere in the area.
The significance of seismic anisotropy in understanding the deformation and flow pattern in the upper mantle has grown in the last two decades.Seismic anisotropy signatures appear in the mantle and crustal rocks due to prolonged deformation and shearing within the mantle and lithosphere (e.g., Long & Silver, 2009;Savage, 1999).
In general, the main mechanisms that contribute to seismic anisotropy include Lattice-preferred orientation (LPO) of anisotropic minerals, in particular olivine (Hess, 1964;McKenzie, 1979), caused by flow-induced strain (e.g., Long & Silver, 2009;Savage, 1999), and Shape-preferred orientation (SPO) in the form of fluid-filled cracks in the crust or of aligned mineral aggregates in a deformed lithospheric mantle (e.g., Holtzman & Kendall, 2010).In active tectonic regions and beneath oceanic lithospheres, seismic anisotropy can be developing because of deformation and mantle flow.However, "frozen-in" anisotropy can also be present in the lithosphere from prior significant deformations, such as continental collision and shear zones (e.g., Barruol et al., 1998;Silver, 1996).Hence, exploring seismic anisotropy can shed light on the evolutionary processes in the lithospheric mantle.Northern Oman is a continental margin that has experienced an extensive developing process of orogeny, rifting and basin formation (e.g., Ali et al., 2013;Jabir et al., 2023;Searle & Cox, 1999).There is, therefore, a strong likelihood that seismic anisotropy is present due to "fossil" fabrics preserved in the lithosphere, shearing at the base of the lithospheric plate related to present-day plate motion and local sub-lithospheric mantle convective flow, or a combination of these mechanisms.Recently, in the NW region of the study area, Pilia et al. (2021) identified a potential depth-dependent anisotropy through SWS analysis.They attributed the observed anisotropy to orogenic processes, particularly in relation to the obduction of the Semail ophiolite and continental collision in the Zagros mountains.
Shear-wave splitting of core-refracted phases (e.g., SKS, SKKS, PKS, named XKS henceforth) is one of the most common manifestations of seismic anisotropy (Silver & Chan, 1991).Shear-waves by entering an anisotropic medium exhibit birefringence and split into two separate orthogonal quasi shear-waves which travel at different wave speeds.The delay time, δt, between the fast and slow components and the polarization direction of the fast component (Fast Polarization Direction, FPD), φ, represent the strength and orientation of seismic anisotropy, respectively (Silver & Chan, 1991).
In the present study, we analyze data from a temporary seismic network consisting of 45 stations completed by 13 permanent stations, which cover the entire region of northern Oman and its major tectonic structures.Shear-wave splitting analysis and teleseismic core-refracted phases as input data are employed to infer seismic anisotropy.In this research, we examine how lithospheric deformation and mantle flow are linked to seismic anisotropy, providing more insights into the past and current status of the lithosphere beneath northern Oman.

Geological Setting
The study area is located at the northeastern edge of the Arabian plate (Figure 1a), composed of Neoproterozoic island arcs formed during the Pan-African orogeny (Allen, 2007;Cozzi et al., 2012).As a result of opening of the Neotethys Ocean, the northern Oman (located in the former Gondwana terrane) experienced rifting (Ali et al., 2013;Ali & Watts, 2009;Jabir et al., 2023) and passive margin tectonics (Pillevuit et al., 1997;Ruban et al., 2007) during the Late Carboniferous.This tectonic evolution was followed by the obduction of the Semail Ophiolite in late Cretaceous (e.g., Tippit et al., 1981), and post-obduction processes, including late Eocene extension which led to the uplift of Oman mountains (Ninkabou et al., 2021;Weidle et al., 2023).
The subsurface structure beneath the study area has undergone multiple phases of deformation during various geological periods since the Permian breakup of Pangea (e.g., Ali et al., 2013;Jabir et al., 2023;Ninkabou et al., 2021).Regional-scale tomographic images (e.g., Pilia, Hu, et al., 2020;Pilia, Jackson, et al., 2020;Shad Manaman et al., 2011) illustrate a transition from thick to thin crust across the northern Oman from west to east.Pilia, Hu, et al. (2020) and Pilia, Jackson, et al. (2020) derived higher velocities in the Gulf of Oman and lower velocities toward the west.They attribute this velocity transition zone to the boundary between the Cretaceous Tethyan oceanic lithosphere and the Arabian passive continental margin.Surface wave tomography investigations (Weidle et al., 2022;Wiesenberg et al., 2022) demonstrate that the Arabian foreland crust, west of the Jabal Akhdar zone, has a thickness of 40-44 km and exhibits seismic properties typical of felsic continental crust.The Vp/Vs ratio in their results depicts a transition from felsic to mafic lithology compositions, which occurs within the area spanning from Jabal Akhdar zone to the SGFZ.Moving eastwards, the Moho depths decrease to 30-35 km, and an increased lower crustal Vs value suggests the presence of mafic intrusions.The crustal thickness beneath the peninsula in the northwestern part of the study area was studied using gravity and magnetic modeling (Ali et al., 2020;Geng et al., 2022;Hansen et al., 2007;Ismaiel, Ali, Watts, & Barkat, 2023), teleseismic receiver function (Ismaiel, Ali, Pilia, et al., 2023), and wide-angle seismic measurements (Pilia et al., 2021).These studies reported crustal thickness vayring in the range of 32-38 km beneath the foreland basin, which increases to 40-52 km beneath the Semail ophiolites and decreases again toward the Gulf of Oman.
Large-scale tomography studies (e.g., Celli, Lebedev, Schaeffer, & Gaina, 2020;Celli, Lebedev, Schaeffer, Ravenna, & Gaina, 2020) have depicted variations in shear velocity within the lithospheric mantle of the Arabian plate.The cold lithospheric mantle of the ancient Arabian craton meets the low-velocity region beneath northern Oman.Al-Lazki et al. (2014) associated the very low Pn velocity in northern Oman with a hot upper mantle, likely indicating the presence of partial melt beneath this region.Furthermore, Weidle et al. (2023) discuss how this hot material emerges from depths greater than 200 km below Africa-Arabia and ascends to approximately 80 km on the northern Oman margins, causing partial remelting and erosion of the lithosphere through small-scale convection at the craton edge.

Data Acquisition
Seismic data used in this study were recorded from network 5H (COOL project), consisting of 40 broadband seismometers across the Oman Mountains, which were in operation for about 2 years and 5 months (from October 2013 to February 2016).Some of these seismometers were occasionally displaced, resulting in 45 station locations.Three different types of seismometers are exploited in this network, CMG-3ESP, Trillium120QA, STS-2, which were able to continuously record at 100 sps (Weidle et al., 2013).This network was deployed with the purpose of illuminating the 3D geometry and internal properties of the world's reference obducted oceanic crust, and the underlying continental lithosphere.The network is complemented by 13 permanent stations (network OM).These stations are equipped with an STS-2.5, Trillium 240, CMG-3TB, or Trillium 120 seismometer, operated by the local Earthquake Monitoring Center in Oman (Figure 1b).

Shear Wave Splitting (SWS)
In the past two decades, XKS phases have been extensively used to study seismic anisotropy (e.g., Kind et al., 1985;Savage & Silver, 1993;Silver & Chan, 1988, 1991).These phases have the advantage of being converted from P-to-S waves at the core-mantle boundary (CMB), meaning that splitting observed at the surface is the result of receiver-side anisotropy.In addition, with a near-vertical angle of incidence, these phases provide a relatively good lateral resolution beneath the monitoring stations.However, the signature of anisotropy in these phases is a vertically integrated effect from the CMB to the surface, and thus they can provide limited constraints on the depth of anisotropy.
Seismograms used in this study were gathered from teleseismic events that occurred in epicentral distances of 85°-140°with a magnitude of Mw > 5.8, resulting in a total number of 10,354 seismograms from 336 events.The National Earthquake Information Center (NEIC) provided the origin time and location of these events.The IASP91 reference model was used to calculate the theoretical arrival times of XKS phases and to determine the time windows for the splitting analysis.
We use the SplitRacer code of Reiss andRümpker (2017a, 2017b), an open-source MATLAB plugin, to perform the SWS analysis.Signals are filtered between 0.02 and 0.2 Hz to improve the signal-to-noise ratio (SNR).The SNR threshold is set to be above 2.0.SplitRacer allows categorizing single splitting measurements, assigning a quality rating labeled as good, average, or poor to each measurement (Figure S1 in Supporting Information S1).This is done manually, taking into consideration various factors such as the elliptical shape of the uncorrected particle motion, the percentage of the energy reduction on the transverse component, the size and shape of the 95% confidence contour, the distribution of φ and δt values from different chosen time windows (here 10 different time windows for each phase), and the correlation between the fast and slow split shear waves components (see also Komeazi et al., 2023, Section 2.2).If the uncorrected particle motion is apparently linear or the energy reduction of the T-component is negligible (less than 10%), the measurement is considered as "null."A null measurement could be attributed to either no anisotropy or very complex structure beneath the station or alignment of the initial polarization direction of the XKS phase with the fast/slow direction within the anisotropic medium.
In the first step, for each station we search for splitting parameters (φ and δt) that minimize the energy of Tcomponent (Silver & Chan, 1991) for each individual XKS phase.These single SWS measurements are then categorized based on their quality and only good and average measurements (Figure S1 in Supporting Information S1) are used to confirm the azimuthal variation of the parameters.SplitRacer also allows for the joint splitting analysis of the XKS waveforms.Joint splitting analysis refers to simultaneous minimization of the total energy of T-components of all XKS phases at one station.In the case of a 90-degree periodicity of the splitting parameters versus back azimuth, indicating a depth-dependent anisotropic media (Rümpker & Silver, 1998;Silver & Savage, 1994), we conduct a joint inversion, assuming a two-layer anisotropy.Due to its fast performance, a modified version of SplitRacer provided by Link et al. (2022) was used to conduct the two-layer joint splitting analysis.The correction of the T-component waveforms for splitting parameters is performed using an inverse splitting operator (Rümpker & Silver, 1998).The splitting parameters (for one-and two-layer models) obtained from the joint inversion for all stations are presented in Tables S1 and S2 in Supporting Information S1.

Ps-Splitting Analysis
In order to examine the presence of seismic anisotropy in the crust and its potential effect on the XKS splitting measurements (e.g., Latifi et al., 2018), we perform splitting analysis technique developed by Rümpker et al. (2014) on the Moho-converted P-to-S (Ps) phases observed in receiver function data computed at the same stations used for the XKS splitting analysis (Figure 1b).Teleseismic events in the epicentral distances 30°-95°are used for receiver function calculations.Preprocessing of receiver function analysis includes detrending and demeaning of waveforms, bandpass filtering from 0.02 to 2.0 Hz, and visual inspection of the waveforms to make sure about a clear P-wave arrivals.Receiver functions are then computed using the time-domain deconvolution technique (Ligorría & Ammon, 1999).
The finite-frequency fast and slow shear waves overlap in weakly anisotropic media (T >> δt), resulting in a characteristic sinusoidal move-out as a function of backazimuth in the Ps arrival time on radial receiver functions.We follow a two-step grid search approach to find the splitting parameters that cause the sinusoidal move-out on the radial components of the receiver functions and the presence of energy on the transverse components (Rümpker et al., 2014).First, the best pair of parameters (φ and δt) is obtained by stacking the amplitudes on the radial components along the sinusoidal move-out.In the second step, the estimated parameters are verified to ensure the energy reduction in the transverse (T) component.If necessary, the parameters are modified within the range of the initial estimates (±15°for φ and ±0.3 s for δt) to better remove the T-component energy.Finally, if the estimated parameters fail to significantly (more than 40%) reduce the T-component energy, the measurement is considered unreliable and not reported.See Rümpker et al. (2014) for a detailed description of the approach.

Individual XKS Splitting Measurements
After carefully choosing the input data based on epicentral distance, magnitude, and SNR, a total number of 553 well-defined (quality good or average) individual XKS splitting measurements (1-29 measurements per station) were retained.Figure 2a  In the western part, the measurements performed at stations MOH, ASH, COO26, COO28, and COO29 are characterized by consistent φ and δt values.For the majority of the stations in this region, the orogenic belt and FPDs are aligned in the NNW-SSE/NW-SE direction, except for station BAN, located at the far north of western Oman, with some NE-SW-oriented FPDs.The delay times measured in the western part mainly range around 1.0 s.In the central region, the majority of high-quality measured FPDs fall within the range of N00°E to N045°E (NE-SW) and N125°E to N180°E (NW-SE).These observed variations in the FPDs could indicate the presence of multiple anisotropic layers.Additionally, most of the null measurements obtained in this study are from stations located in the central part, including COO12, 17, 20, 24, and 43.
Generally, more individual SWS measurements are obtained at the eastern part relative to the western part of the study area, not only due to a denser station coverage, but also due to more measurements obtained per station.To investigate whether the different patterns of splitting parameters observed in the central part are due to a different distribution of data in this region, we have selectively plotted single measurements for all stations based only on those events that were recorded by the stations in the central area (Figure S3 in Supporting Information S1).Upon examination, we observed no significant change in the pattern of splitting parameters.Most of the FPDs at the stations located in the eastern part are oriented in the WWN-EES, NW-SE and NNW-SSE directions.However, a few stations also include FPDs with a NE-SW trend.The observed variation of the FPDs at stations such as BID, COO33 and COO53 indicate the structural complexity underneath this region, which could be due to multiple sources of seismic anisotropy.A relatively sharp change in the direction of FPDs occurs at the SGFZ, from dominantly NW-SE orientation in the eastern side to dominantly NE-SW trend on the west.This pattern (also seen at the neighboring stations) might reflect short-wavelength variations in the upper mantle structure.
A further consideration is that some stations possess backazimuth-dependent splitting parameters.This might be due to the presence of multiple layers of anisotropy (Rümpker & Silver, 1998;Silver & Savage, 1994).Clear azimuthal variations can be seen at most of the stations in our study, particularly in the eastern and central parts (Figure S4 in Supporting Information S1).At several stations, the values of the fast polarization display systematic azimuthal variations with a period of ∼90°.Such patterns are consistent with a two-layer anisotropic model throughout the eastern and central part of the study area.As an example, Figure 2b shows the azimuthal variation in splitting parameters at station COO17.However, some stations located in the western parts (BAN, MDH, ASH) suffer from insufficient azimuthal coverage to infer the dependency of the splitting parameters on the backazimuth.Some stations exhibit large delay times (>2 s, up to 2.77 s recorded at COO35) in the eastern part of the study area.These (apparent) large delay times could be the effect of depth-dependent anisotropy.
We additionally investigate whether the observed individual splitting measurements originated from a single layer of anisotropy positioned at a specific depth.For this purpose, we follow the spatial coherency approach of Liu and Gao (2011) modified by Gao and Liu (2012).In this method, individual splitting measurements are projected to the ray piercing points, and then by calculating the spatial variation factors, which is a dimensionless weighted average of delay times and FPDs at each point within the Fresnel zone (Liu & Gao, 2011), the procedure determines the depth at which a maximum level of spatial coherency is reached.For a comprehensive understanding of the approach and the detail description of the parameters and mathematical computations please refer to Liu and Gao (2011).The anisotropy at this depth is taken as the main contributor to the observed shear wave splitting.We divided the study area into grids with a size of 0.1°in both longitude and latitude directions.Then for each grid point the average φ and δt for all the splitting parameters within the Fresnel zone (R f = 85 + 0.2*depth (km) after Liu & Gao, 2011) is calculated.Finally, at each assumed depth the variation factor is computed.Our analysis suggests an optimal depth of 290 km (Figure S5a in Supporting Information S1) for the source of anisotropy assuming a single layer of anisotropy.However, even by projecting the splitting parameters at the optimal depth (Figure S5b in Supporting Information S1), there are still a considerable number of non-coherent splitting parameters.These results confirm that a large-scale one-layer model is not a precise approximation for the XKS observations at most of our stations.

Joint Inversion Modeling
The joint-analysis approach described in Section 3.2 allows for the inversion of the XKS waveforms for one or two layers of anisotropy.Final joint splitting parameters from our XKS analysis are shown in Figure 3, along with more detailed parameters in Tables S1 and S2 in Supporting Information S1.As expected from the backazimuthal variations of the individual splitting measurements, we retrieve mostly two-layer models in the eastern and central parts (Table S2 in Supporting Information S1), whereas the western part is dominated by one-layer models (Table S1 in Supporting Information S1).At stations with two-layer models, T-component energy reduction is significantly improved relative to the one-layer model at the same stations (Figure S6 in Supporting Information S1).Insufficient coverage of back azimuths at a few stations (SOH, COO25, COO28, COO27 and COO23) makes it difficult to draw definitive conclusions about the depth-dependent anisotropy beneath these stations (Figure S7 in Supporting Information S1).The joint XKS splitting delay times are on average around 1.1 s for one-layer models, and 0.83 and 0.76 s for the upper and lower layers, respectively, in the case of two-layer models.The slightly larger delay times in the upper layer compared to the lower layer suggest that the study area is mostly characterized by upper mantle and lithospheric anisotropy.Furthermore, in central parts, the delay times are much smaller than the other parts of the study area.In the eastern part, the fast symmetry axes of the upper layer predominantly trend toward the WWN-EES, consistent with the trend of the FPDs of the one-layer model beneath the stations in the western part.The lower-layer FPDs appear to be well-clustered toward the NNE-SSW/N-S direction.A lateral variation in the upper layer FPDs occurs along the strike of the mountain range changing from a WWN-EES trend in the eastern part to a NE-SW trend in central part at the SGFZ that then changes to a dominantly NW-SE direction in the west.
In the example depicted in Figure 2, the measurement errors, derived from 95% confidence levels (Reiss & Rümpker, 2017a, 2017b), appear substantial.This underscores the necessity for conducting additional tests to validate the assumption of systematic backazimuthal variations.Moreover, Latifi et al. (2018) and Rümpker et al. (2023) argued that the results of layered anisotropy modeling could be highly non-unique and emphasized the need for prior information from the region to validate these findings.In our study, we conducted several tests by incorporating prior information from the study area to fix the FPDs in the upper layer (in the direction of mafic intrusions) or lower layer (in the direction of the absolute plate motion (APM)) (Figures S8 and S9 in Supporting Information S1).Fixing a parameter (e.g., the upper/lower layer FPDs in our tests) leads to avoid decreasing the variance due to more degree of freedom.In this scenario, the reduction in transverse component energy becomes more comparable to the energy reductions observed in the 1-layer modeling.These tests illustrate that fixing the FPDs in the upper and lower layers to the NW-SE and NE-SW directions, respectively, results in more reliable two-layer parameters, consistent with the results presented in this study.
In addition, we conducted F-test analysis (Lomax & Hahs-Vaughn, 2013) for some stations in the eastern area (Table S4 in Supporting Information S1).This analysis involves comparing the variances of FPDs and delay times resulting from one-layer and two-layer modeling.After fitting both models to the data and calculating their residual sum of squares (RSS), the F-statistic is computed.Then, we compare the F-statistic to a critical value from the F-distribution (in our study it is determined based on a significance level of 0.05).The results for these stations demonstrate p-values greater than critical value for both FPDs and delay times, indicating a significant improvement in data fit when employing the two-layer modeling in this area.

Effects of Crustal Anisotropy
Out of all the stations where we applied the Ps-splitting method, only 15 stations fulfilled the requirement of sufficient azimuthal coverage, which is necessary to obtain reliable results in Ps-splitting method.In Figure S10 (and Table S3 in Supporting Information S1), these computed splitting parameters are illustrated.At stations situated on the western side of the SGFZ, the orogen-parallel splitting parameters are consistent with the upper crust anisotropy from surface waves (Wiesenberg et al., 2022).On the eastern side of the SGFZ, however, FPDs are rather NE-SW and thus more consistent with lower crustal anisotropy observed from surface waves.This suggests a dominant upper and lower crustal anisotropy in the central and eastern parts of the study area, respectively.
Using crustal anisotropy results from the Ps-splitting (an example is illustrated in Figure S11 in Supporting Information S1), we examine the effect of the crustal anisotropy on the XKS splitting observations at several stations with adequate azimuthal coverage.At these locations, we search for two-layer models by following two schemes.(I) two-layer modeling by fixing the fast direction of the upper layer in the direction of the crustal anisotropy.(II) two-layer modeling after removing the crustal anisotropy effect from the waveforms.
In scheme I, the joint inversion of the XKS waveforms is conducted for a two-layer model in which the parameters of the upper layer are constrained to vary around the parameters obtained from the Ps-splitting.Then the joint inversion gives the delay time of the upper layer and the delay time and FPD of the lower layer.In scheme II, in the first step, the XKS waveforms are corrected according to the splitting parameters obtained by the Ps-splitting analysis of receiver functions.The required mathematical equations to correct the waveforms for a certain φ and δt are explained in detail in Rümpker et al. (2023).After correcting the waveforms for the crustal anisotropy parameters, we applied the joint splitting method for a two-layer anisotropy.
The results of this inversion for both schemes are shown in Figure 4 along with the two-layer joint splitting parameters obtained without constraining the crustal anisotropy for the same stations.Only six stations delivered reliable layered models by constraints from crustal anisotropy.Due to a more effective reduction in average Tcomponent energy provided by scheme II, the splitting parameters derived from this scheme are selected as the final two-layer parameters.As illustrated in Figure 4, for stations COO16 and COO20, situated in the western part of the SGFZ, the FPDs given by scheme II are nearly identical to those derived from the inversions without crustal constraint.Nonetheless, the delay times at both stations increased under scheme II.In contrast, for stations COO32 and COO37, the delay times of the upper layer decreased in scheme II compared to the inversions without crustal constraint.Additionally, at stations COO05 and COO51, distinct anisotropic layers are observed in various schemes.

Discussion
A shear wave splitting analysis along northern Oman is provided in this study by using a dense network of stations.We discuss the possible causes of seismic anisotropy beneath northern Oman in the context of the tectonic setting of the region.We examine the effects of extensive processes of orogeny, rifting and basin formation at a continental margin on the lithosphere beneath the northern Oman.The presented results of the XKS Many of the delay times are too large to be attributed exclusively to the lithosphere.Considering the maximum delay times measured in this study, which are of the order of ∼2.5 s, a maximum thickness of ∼125-250 km can be expected for the anisotropic layer, assuming an average anisotropy strength of 2%-4% (e.g., Helffrich, 1995;Savage, 1999) and an average Versus of ∼4 km/s for the upper mantle (Weidle et al., 2022).Weidle et al. (2023) suggest that the lithospheric thickness beneath our study region is in the order of ∼80-120 km.A uniform anisotropic layer with a mean value of 4% of anisotropy can explain the average delay time value of ∼1.0 s observed in the region.However, as explained in previous section, a single layer of anisotropy cannot explain the lateral and azimuthal variations in the observed XKS splitting measurements.Therefore, we need to discuss the observations by considering both the lithospheric and sub-lithospheric sources of anisotropy.

Journal of Geophysical
In Figure 5, we schematically summarize the inferred depth-dependent anisotropic model that is obtained for the entire region by joint inversion of the XKS data and considering the effect of the crustal anisotropy as constrained by the Ps-splitting analysis.

Lithospheric Deformation
The one-layer joint splitting parameters determined in our study (Figure 3a) can be categorized into NW-SE and NNE-SSW/NE-SW trending FPDs.A prominent mountain parallel pattern of the one-layer FPDs is observed in the study area, which is likely due to the compressional forces resulting from the Neo-Tethys closure that could have developed "fossil" anisotropic fabric in the lithospheric mantle.However, a clear distinction regarding the pattern of anisotropy occurs at the SGFZ.On the eastern side of the SGFZ, FPDs are trending NW-SE, whereas the FPDs are mostly oriented NE-SW west of the SGFZ.Finally, moving from west of Hawasina window toward western part of the study area, the FPDs revert to the NW-SE trend.This patterns of the FPDs portraits the central part to be different relative to the neighboring areas.The variation in our observations in the central part could be due to the complexity arising from the localization within the transition zone (see also Section 2).The entire study area underwent extensive lithospheric modification during the Permian breakup of Pangea-first by extension and breakup to the SE (Indian Ocean) and NW (Neo-Tethys) and second by mafic modification, on eastern domain of the SGFZ, which is characterized by a strongly volcanic margin (e.g., Chauvet et al., 2009;Weidle et al., 2023).In the western parts of the study area, the one-layer FPDs (with NW-SE orientation) are very similar, suggesting the overall uniformity of the anisotropic structure at depth.It is worth noting that regarding the long wavelength of the waves in our data set (∼50 km), it is not possible to detect the local small-scale variations beneath this region.In comparison to the previous study (Pilia et al., 2021), their averaged values agree well with our joint splitting parameters for one-layer models in the same region or stations (Figure S12 in Supporting Information S1).At different stations, mean split times range from 0.6 to 1.6 s, with smallest values suggesting weak anisotropy, a relatively thin anisotropic layer where the crust may have played a significant role, or possibly two perpendicular anisotropic layers that cancel out each other's effects.
The backazimuthal variation of the splitting parameters at some stations, particularly in the central and eastern parts (Figure 2b and Figure S4 in Supporting Information S1), prompted us to examine two-layer models.We also employed the method proposed by Liu and Gao (2011) to assess the consistency of our measurements with a single-layer anisotropy.In comparison to the examples presented in the study by Liu and Gao (2011), our results reveal that the variation factor shows negligible depth dependency.Additionally, the observed discrepancy of the splitting parameters at the optimum depth, where variations are minimized, suggests that a single anisotropic layer cannot adequately explain our observations.The two-layer joint splitting parameters are shown in Figure 3b.A NW-SE trend is observed for the upper layer and NE-SW/NNE-SSW for the lower layer, which we attribute to the lithospheric mantle and asthenosphere, respectively (Figure 5).As delay times of the lower layer increase from west to east, we infer a thickening of the anisotropic layer (asthenosphere) under a thinning lithosphere.
The majority of the FPDs of the upper layers, mostly in the eastern side of the SGFZ, as well as the FPDs of the one-layer models in the entire region are subparallel to the main tectonic trends visible at the surface, resulting from compressional tectonics during Arabia-Eurasia convergence.In addition, on the eastern side of the SGFZ, emplacement of the mafic intrusions in the ∼WNW-ESE direction (Weidle et al., 2023) can partially explain the observed NW-SE apparent anisotropy in this area.Furthermore, the SGFZ is discussed in previous studies as a distinguishing tectonic feature.For instance, the abundance of mafic rocks in the Saih Hatat Dome has been attributed to the plate breakup (Chauvet et al., 2009) and since the Jabal Akhdar Dome region was minimally impacted by the mafic rocks, the SGFZ is thought to be acting as a structural barrier at the time, as proposed by Scharf et al. (2019).In the absence of mafic intrusions on the western side of the SGFZ, the upper-and lower-layer delay times are comparable, suggesting a thicker lithospheric mantle beneath the central parts.Therefore, given the NE-SW trend of the one-layer FPDs in the central part, we argue the crustal anisotropy resulting from the NE-SW extension plays a significant role in this area.However, in the west, no significant backazimuthal variation is observed and the predominant NW-SE trend of the one-layer FPDs arise from the thick and old lithospheric mantle, which has undergone the compressional tectonic regime.
The average delay times for the crustal anisotropy as obtained using the Ps-splitting analysis are in good agreement with observations reported from other regions, which is of order of ∼0.3-0.5 s (Barruol & Mainprice, 1993;Silver, 1996).The orogen parallel pattern of the Ps-splitting results for the crustal anisotropy resemble those for the upper-layer joint splitting results as well as the azimuthal anisotropy estimated by surface wave tomography (Weidle et al., 2023;Wiesenberg et al., 2022).The analysis of crustal anisotropy has revealed more insights into the deep structure of the eastern and central regions.At stations COO32 and COO37 in the eastern region, it appears that the crust and lithospheric mantle have undergone coherent deformations, resulting in smaller delay times for the upper layer when correcting the XKS waveforms for crustal anisotropy.In contrast, stations COO51 and COO05 yield completely different models relative to the case with no crustal constraints, indicating the complexity beneath this area.West of the SGFZ, stations COO16 and COO20 show a nearly orthogonal FPD between the crustal anisotropic layer and the lower layer, thereby weakening the effect of the bottom layers.This could explain the relatively small joint splitting delay times for the stations located west of the SGFZ (e.g., COO24 and 13), even though the individual measurements contain much larger delay times.These findings provide insights into the complex dynamics of the crust and mantle across the study area.

Large-Scale Mantle Flow Due to The Plate Motion and Subduction
A basic assumption regarding the shear (basal drag) in the transition layer between the lithosphere and the underlying asthenosphere is that it is influenced by the APM, causing the orientation of the olivine a-axis to be subparallel to the plate motion (Silver, 1996;Zhang & Karato, 1995).The development of the basal drag fabric is believed to be influenced by plate velocity (Debayle & Ricard, 2013), wherein slow-moving plates are unable to generate a sufficiently strong basal drag fabric to produce significant anisotropy.While the plate convergence speed reached 6 cm/yr during the Cretaceous (Agard et al., 2007), GPS measurements by Vernant et al. (2004) and Khorrami et al. (2019) showed a rate of approximately 27 mm/yr for the current Arabian plate motion.However, there is no evidence to suggest that this reduction in velocity would impact the pervasive anisotropy fabric developed by the basal drag in the direction of APM.
The predominant observed FPDs within the lower anisotropic layer (N7°-55°E) in our study exhibit an approximate alignment with the APM direction (∼N47°E) in northern Oman, assuming the no-net-rotation reference frame for plate motion (Kreemer et al., 2014).The slight discrepancy between the observed FPDs in the lower layer and the APM can be attributed to the northward movement of the lithospheric mantle.This north component of the lithospheric mantle could be due to the northward dipping subduction zone within the Neo-Tethys Ocean (Tavani et al., 2020), which resulted in a northward deviation of the FPDs caused by the APM.The same argument applies to the FPDs observed at stations COO08, COO09, and COO10 in the southeastern region, which are predominantly oriented in NNE-SSW direction.The measured splitting parameters at these stations arise from a combination of mantle flow and minor influences of lithospheric anisotropy.

Conclusions
Using a combination of data from both permanent and temporary seismic networks, we studied seismic anisotropy beneath northern Oman.Most of the stations located in the eastern and central parts of the study area exhibit systematic azimuthal variations of the FPDs of splitting analysis.These azimuthal variations suggest the potential presence of at least two layers of anisotropy beneath this region.On the other hand, for the stations located to NW, our observations suggest dominantly one-layer anisotropic models.Overall, most of the upper layer in the twolayer models in the east, coincident with the fast directions of the one-layer models to west, which in turn are subparallel to the major tectonic trends and boundaries visible at the surface.The one-layer anisotropy parameters as well as the parameters of the upper layer in the two-layer cases mark a sharp boundary at the SGFZ.This, along with the significant variation in the FPDs between some neighboring stations, suggests that the lithospheric structural alignments are the main source of anisotropy of the upper layer.The observed orogen-parallel anisotropy supports the idea of frozen-in lithospheric mantle anisotropy resulting from the compressional forces in northern Oman.In the east, the combination of thin lithosphere and mafic intrusions effects are dominant, making the NW-SE trend.Within the central part, to the west of SGFZ, the presence of a thicker lithosphere empowers the role of crustal anisotropy.A relatively weak apparent anisotropy in the central part might be due to complex lateral and depth variation of the anisotropic structure in the lithosphere.In the west, the NW-SE FPDs highlight the significant role of the thick and old lithosphere of the Arabian plate in the resulting anisotropy.The predominantly NE-SW trend of the lower layer anisotropy in the central and southeastern part of the study area is parallel to the APM.Therefore, we argue that this anisotropic fabric is created by the large-scale mantle flow beneath northern Oman due to the plate motion.

Figure 1 .
Figure 1.(a) An overview map of the tectonic setting in northern Oman (after Scharf et al., 2019).SHD, Saih Hatat dome.JAD, Jabal Akhdar dome.HW, Hawasina Window.SGFZ, Semail Gap Fault Zone, (b) Map showing the temporary (black triangles, network 5H) and permanent (red triangles, network OM) seismic stations used in this study.The inset maps show the event distribution within the epicentral distance of 85°-140°used to extract the XKS waveform data together with a rose diagram illustrating the backazimuthal coverage of the events.The elevation data is taken from Tozer et al. (2019).
depicts the individual SWS measurements at each station.The backazimuths of the events giving null measurements are shown in Figure S2 in Supporting Information S1.Based on the individual FPDs, the study region can be divided into three subregions: western, central, and eastern parts.The SGFZ separates the central part from the eastern part, and the western part is separated from the central part at around longitude 56.5°E, west of Hawasina Window (HW).

Figure 2 .
Figure 2. (a) Individual splitting parameters of XKS waveforms are shown by red bars.They are oriented in the fast polarization direction with length proportional to the delay time.(b) An example of backazimuthal variation at station COO17.The red and blue circles show non-null and null measurements.The green dash-lines are theoretical two-layer modeling parameters derived from joint inversion of XKS waveforms.δt and φ for the upper and lower layer are respectively ( 43.23°, 0.56 s) and (21.01°, 0.49 s).

Figure 3 .
Figure 3. (a) One-layer splitting parameters obtained from joint-splitting analysis of XKS waveforms, shown with red bars, (b) two-layer splitting parameters obtained from joint-splitting analysis of XKS waveforms, shown with yellow (upper layer) and blue (lower layer) bars.On both maps, the bars are oriented in the fast polarization direction with length proportional to the delay time.Black arrows illustrate the direction of the absolute plate motion.

Figure 4 .
Figure 4. (a) Joint splitting parameters for two-layer modeling without any constrictions.(b) XKS waveforms are first corrected by Ps-splitting results, and then twolayer joint splitting analysis is performed on corrected XKS waveforms (scheme II).Ps-splitting parameters indicated by red bars.The upper-and lower-layer joint splitting parameters are shown with yellow and blue bars.The bars are oriented in the obtained fast directions, and their length is proportional to the delay time.(c) The fast polarization directions of the upper layer are constrained to vary around the parameters obtained from the Ps-splitting (scheme I).

Figure 5 .
Figure 5. Schematic visualization summarizing representatively the anisotropy parameters observed beneath the study area using the Ps-splitting and Joint splitting methods.Bars in different colors as described in the legend show the orientation of the fast polarization directions and their length are proportional to the delay times.The schematic geophysical features depicted in the depth cross-section are derived from Weidle et al. (2023, Figure 6).