khalil motaghi

This study presents a 3-D model of the shear wave velocity (Vs) for SE Iran, containing the Makran subduction zone and termination of the Zagros collision zone. The model is developed using surface wave tomography and the associated local dispersion curves that were inverted to the Vs velocity maps. Our findings provide insights into the tectonics of the SE Iran region and the geometry of the subducting lithosphere. To obtain Rayleigh waveforms, we extracted data from 458 teleseismic earthquakes recorded by 58 seismic stations in the Makran region between June 2016 and May 2019. Events with magnitudes greater than 5.5, epicentral distances between 30º and 120º, and depths shallower than 50 km were considered for two-plane wave tomography. This resulted in phase velocity maps at nine periods ranging between 25 s and 125 s. Our results suggest the presence of a high-velocity anomaly at the northern coast of the Gulf of Oman extending ~100 km northward up to the Qasr-e Qand fault and the southern edge of the Jaz Murian depression as it was observed at shorter periods of 25-40 s. Further west of Makran and north of the Strait of Hormuz, there is also a high-velocity anomaly at the southwest of the Zendan-Minab-Palami fault at periods of 25-40 s. At periods of 50-125 s, this anomaly is observed at the northern latitudes beneath the Sanandaj-Sirjan zone at the termination of the Zagros collision zone. Our analysis suggests that this anomaly reveals underthrusting of the Arabian lithosphere under Sanandaj-Sirjan zone and Urumieh-Dokhtar magmatic arc. Additionally, we detected a low-velocity anomaly at periods of 25-33 s showing a crustal root potentially generated by underthrusting of the Arabian lithosphere in this region. At the second inversion step, we employed a nonlinear inversion of the local dispersion curves to construct a 3-D Vs model. Our results indicate that there is a flat high-velocity anomaly in the middle of the study region (in the western Makran and under the depression of Jaz Murian), indicating a horizontal oceanic lithosphere, potentially remaining from a truncated oceanic lithosphere between the collision and subduction zones. Furthermore, our findings suggest the presence of a near-horizontal oceanic lithosphere in the eastern part of the study area, extending at a low angle underneath the northern edge of the Jaz Murian depression and subsequently subducting beneath the volcanic arc with a steep angle. Finally, we identified a low Vs anomaly at the crustal-scale in the Sistan suture zone, which stretches towards the north-northwest and is limited by the large N-S directed strike-slip faults.

The potential field data like Gravity data, Magnetic data, Self-potential data, and other natural source data are accessible but hard to interpret and model. Numerous research presents ideas for modeling the potential data; however, they are mostly based on inverse or forward modeling needing a priori constraints and information. In complicated geology and tectonic setting, we do not have convenient access to a priori information to define the constraints. So, we must develop a pure geophysical interpretation method without geological constraints and minimum complexity. In this research, the Normalized Full Gradient’s ability to find the gravity anomaly model was studied. The Normalized Full Gradient method is an effective method for determining anomalous bodies, such as the distribution of oil and gas fields or structural boundaries. The Normalized Full Gradient method depends on the downward analytical continuation of normalized full gradient values of gravity data. Analytical continuation discriminates certain structural anomalies which cannot be distinguished in the observed gravity field. The Normalized Full Gradient of the gravity anomaly is often used rather than the gravity anomaly itself for detecting underground spaces because it is stable and indicates the locations of source bodies. The weakness of the Normalized Full Gradient is the 3D modeling limitation, as we can only calculate the 2D response in practice. On the other hand, the responses can not describe the negative and positive parts of the anomaly. But a unique advantage of the Normalized Full Gradient is, that it does not need the primary information for gravity data modeling. In this research, we used the Normalized Full Gradient for the large-scale Bouguer gravity anomaly interpretation. Bouguer gravity anomaly wavelengths contain information about density distributions of upper mantle and lithosphere structures. A gravity profile is most often a combination of relatively sharp anomalies that must be of shallow origin and very deep and large anomalies with a regional nature.The study profile has a 400 (km) length from SW to NE of North Western Iran, and Sahand, North Tabriz Fault, and Sabalan are the most important structures in the study area. The Normalized Full Gradient synthetic model data study provides the opportunity for the real data recovered model judgment. So, we first showed the Normalized Full Gradient recovered model of the synthetic data test and then based on the resolution of the Normalized Full Gradient used, it provides the lithospheric density of North Western Iran. The result shows the low-density mantle wedge which is probably is beneath the North Tabriz Fault that is responsible for the formation of distinct lithosphere conditions. So, the wedge can explain the complicated tectonic setting of North Western Iran. The mantle wedge has more than 50 (km) wide and more than 40 (km) depth. seemingly, this mantle wedge directly affects the North Tabriz Fault, Sahand, and Sabalan in the shallower depth. The sharpest effect is for the North Tabriz Fault in the shallower part of the mantle wedge and following the shape of the wedge, we can see a sharper effect on the Sabalan in comparison with Sahand.

The Iranian Plateau is a part of the Alpine–Himalayan orogenic belt located in the western part of Asia. Convergence of the Arabian Plate and Eurasia from the late Cretaceous to the present has generated significant lithospheric deformations such as crustal shortening and thickening in the Plateau and surrounding mountain ranges including the Zagros Fold–Thrust Belt, Alborz and Kopeh Dagh. The convergence across Zagros is accommodated through different mechanisms of diffused shortening and/or thrusting of the Arabian lithosphere beneath Central Iran. This study focuses on the velocity structure of the eastern part of the Iranian Plateau which has not been studied well yet. Our study region also contains Binalud and Kopeh-Dagh deformation domains which were built up by the north-eastern collisional boundary between the Plateau and Eurasia and a small part of the Urumieh-Dokhtar magmatic arc which was the volcanic arc of the past Neotethyan subduction. The lithosphere-asthenosphere system beneath east of Iran is investigated by employing earthquake surface wave tomography. A total of 5862 teleseismic Rayleigh waveforms from 368 events recorded at three permanent networks during a period of three years were used to produce 2-D high-resolution phase velocity maps. We employed a two-plane wave tomography approach to generate phase velocity maps at period ranges of 25–111 s. From a published study of ambient noise tomography, we extracted Rayleigh wave dispersion data at 8–20 s periods to improve resolution in the crust and then inverted them for a 3-D S-wave velocity model. A 3-D velocity model was then constructed by a nonlinear Bayesian Markov chain Monte-Carlo algorithm of local node-wise dispersion data into S-wave velocity models down to a depth of 200 km. The most prominent resolved feature by our 3-D velocity model is a low-velocity asthenospheric channel at 70 and 150 km depths overlaid by a thin lithosphere. We believe that in the lack of an isostatic compensated crustal root in the Iranian Plateau, this feature is supporting high elevation (~1000 m) topography covering the Iranian Plateau. A Moho map for the study region is obtained by mapping the geometry of 4.0 km/s S-wave velocity contour in the 3-D velocity model. It shows that most of the study region is covered by a less deformed crust with a thickness of ~36 km. Two crustal roots are observed, one beneath the Urumieh-Dokhtar magmatic arc and the other beneath the north-eastern part of the study region where an array of the Neotethys suture zones is marked by ophiolite outcrops. Lithospheric scale deformation in a sequence of Neotethys suture zones is high probably responsible for the crustal thickening in NE Iran.

Tabriz city with 1.6 m population is located in the NW Iran. Estimation of empirical attenuation relations for the region is a key for realistic seismic hazard assessments. Recently, Motaghi et al. (2016) have employed microseismicity of the region to estimate an attenuation relationship for the NW Iran. Their regression has left significant amplitude residuals at different stations where neighboring stations have similar average residuals confirming important lateral structural variations ignored in the attenuation estimations. For instance, stations located around the NTF have systematic negative residuals consistent with location of a thick sedimentary basin in the region. Such systematic residual patterns have been reported before in Canada and North America by Atkinson (2004) or in Alborz Mountains, north Iran, by Motaghi and Ghods (2012). These observations have inspired us to conduct 2-D amplitude residual tomography similar to widely used local travel time tomography (e.g., Rawlinson and Sambridge, 2003). We assume that attenuation variations are only affected by anelastic coefficient variations generated by change of rock properties. Thus, we have discretized the study region and considered anelastic attenuation coefficient (inverse of quality factor) as an unknown parameter for the inversion procedure. 

Methodology and Approaches:

We formulized a linear relationship between the amplitude residual (as datum) and anelastic coefficient variation (as unknown parameter), and then, we analyzed seismograms recorded from 943 local earthquakes with magnitudes between 1.6 and 5.2, and azimuthal gap less than 250o . The data were gathered by 35 seismic stations located in the NW Iran to carry out the inversion. We used a weighted damped least squares approach (e.g. Menke, 1989) in which weight matrices for the data and model parameters were used. The weights for the data come from the signal to noise ratio calculated for each signal. The weights for the model parameters come from the number of rays per block. This helps to exclude model parameters associated with blocks not crossed by rays. An optimal damping value of 23 was selected from the trade-off curve between the total residuals and the weighted model variances. The inversion procedure estimated variations of intrinsic attenuation coefficient relative to the average value (0.0012). We converted the obtained values into change of quality factor, ΔQ, which was more usual in the literature.

Results and Conclusions

Our tomogram for the NW Iran showed a lower quality anomaly for Neogene sedimentary basin around the North Tabriz Fault in contact with high quality Cretaceous volcanic and Cambrian metamorphic basement outcrops in the north of the region. Such clear consistency between our tomogram and lithology variations at the surface disappeared in the eastern part of the tomogram where deep (depth = 20–45 km) events occurring in the oceanic-like crust of the South Caspian Basin are observed. Comparing our eastern anomalies with a teleseismic tomography across the same region, we found that low and high quality patches in eastern part of our tomogram are probably generated by thermal effects of Sabalan volcano and oceanic like crust of the South Caspian Basin, respectively. Good consistency between our results and previously reported features confirms that the amplitude residuals are proper data set to detect structural heterogeneities responsible for large amplitude variations in active seismic regions.

The empirical attenuation relationship for spectral amplitudes was calculated to study the attenuation of shear waves of shallow events at close hypocentral distances inside Tarom-Rudbar region, western Alborz. 3122 waveforms (170 shallow events) recorded by 40 seismic passive stations of two local temporary seismic networks were included in this analysis. The selected events have moment magnitudes between 1.8 and 4.2 and epicentral distances of 10 km to 70 km. All events have location accuracy better than 2 km in epicenter and less than 5 km in depth. The good location quality of the events allows us to estimate accurately geometrical spreading of shear waves at close hypocentral distances. By selecting the small events, we can safely treat them as point sources and thus use hypocentral distance as our distance metric. Additionally, for smaller events, we automatically avoid the non-linearity of amplitude of seismic waves with magnitude and its possible trade-off with geometrical spreading. Due to a rather low dependence of geometrical spreading on magnitude (v. NGA-WEST2 models), our approach of using weak-motion data may provide a means to reliable assessment of the geometrical spreading coefficient which can then be used to partially regionalize NGA-WEST2 GMPEs for regions with low rate of seismicity or lack of enough strong motion records.
The shape of the attenuation curve at different frequencies was obtained using non-parametric fit to the data with Robust Lowess algorithm showing a mono-linear curve in the associated distance. Assuming mono-linear attenuation model and using regression, the value of geometrical spreading coefficient in the equation derived from the spectral method was obtained as 1.6 in frequencies higher than 2 Hz. Spectral amplitude attenuation curves show an obvious super-spherical geometrical spreading at close hypocentral distances. We show that the geometrical spreading is strongly super-spherical in close agreement with those used in some of the NGA-WEST2 GMPEs. The observed super-spherical geometrical spreading of seismic waves could drastically change the level of seismic hazard in close hypocentral distances by localizing strong motion to short hypocentral distances. The calculation of geometrical spreading coefficient using data from more frequent small events recorded by dense local networks can be used to partially regionalize the geometric term of GMPEs in regions with small rate of seismicity.
The residuals were averaged on a station-by-station basis to determine station corrections. The calculated station corrections for the study area shows sharp contrast between the northern and southern hills of western Alborz. The stations in the northern hill, mostly in Gilan plain, show higher amplification (positive station corrections) relative to those in the southern hill. The strong amplification of seismic waves has a strong implication for preparation of seismic hazard maps of the densely populated Gilan province

On 20 December 2010, an earthquake with Mw 6.5 occurred in Rigan, a small town in the desert south of Bam city. The earthquake epicenter was in a low population area so, luckily, it caused only few casualties. Five days later 76 aftershocks reported by Iranian Seismological Center (ISC). On 27 January 2011, another earthquake (Mw 6.2) stroke an area at ~ 20 km southwest of the first earthquake. Bam earthquake Mw 6.6 occurred in 2003 with 40,000 victims is one of the deadliest earthquakes in Iran which is located in shear zones at southeast Iran. Considering the active faults distribution of the region and aftershocks of the 2010 Rigan earthquake encouraged us to better investigate and model the post-seismic deformation related to the 2010 earthquake. Post-seismic syudying provides information about rheology of the surrounding region and improves our knowledge about the strain release after the earthquake. In this study, COSMO-SkyMed (from Italian Space Agency, ASI) images spanning the temporal interval between 27 January 2011 and 15 July 2011 are used to investigate the post-seismic deformation following both earthquakes. We applied the Small Baseline Subset (SBAS) algorithm for images to obtain the post-seismic mean velocity map and the relative deformation time series. 109 interferograms, post-seismic mean velocity map and the relative deformation time series obtained Mean velocity map shows that displacements of post-seismic phase are right lateral strike slip same as co-seismic mechanism. Time series analysis reveals a clear post-seismic signal exponentially increasing with time until reaching the rate of more than 8 mm/year which indicates the end of post-seismic phase and following inter-seismic phase, starts with steady stress accumulation. Later, we modeled the post seismic signal considering a dislocation on a finite fault in an elastic and homogeneous half-space that are the assumptions for the Okada (1985) model. Post-seismic results modeled by adopting a two-step approach: (1) a non-linear inversion performed to constrain the fault geometry parameters and considering a uniform slip, then (2) a linear inversion performed to retrieve the slip distribution on the fault plane previously obtained. The fault plane is split into 1×1 km patches along strike and down-dip. Determining fault parameters and slip distribution by Okada model, indicates that the slip is concentrated in downdip of the coseismic depth with 1.2 m slip and also at the edges of the coseismic asperity. This slip distribution indicates that “afterslip” is the mechanism for post-seismic deformation of the Rigan earthquake.

We analyzed the teleseismic data gathered by a broad-band (CHBR) and four short-period (CDK, CNT, KHB, KSM) seismometers, located in western coastal Makran, north of Chabahar, Iran. The data were gathered by the roughly north-south direction quasi-linear profile and used to calculate P (for all stations) and S (only for CHBR) receiver functions utilizing iterative deconvolution technique of Ligorria and Ammon (1999). Because of backazimuth gaps in south and western directions, we used PKiKP and Pdiff phases to calculate receiver functions in a similar processing approach. Calculated P receiver functions are migrated to depth to clarify the geometry of velocity boundaries at the base of sediments and Moho. The result shows that there is a dipping interface lying at a depth of 27 km (beneath CHBR) to 31 km (beneath CDK), which imply a 2.5o dipping Moho boundary beneath the study region. To avoid the trade-off between velocity model and reported depth, we jointly modeled the stacked receiver function, and group velocity dispersion curve for CHBR and the output model was considered for any time to depth migration of receiver functions.
We analyzed the effects of P and S anisotropy on teleseismic converted waves to map the presence, the strike, and the depth of anisotropic structures. High-resolution PRFs are considered for such analysis. The following criteria are considered to select the high-quality receiver function (Schulte-Pelkum and Mahan, 2014): the signal-to-noise ratio of the three components of the seismograms is at least 1.5; the convolution of the PRF with the vertical component of the seismogram reproduces at least 60% of the horizontal component (defined as variance reduction by Ligorria and Ammon, 1999); the PRF shows a positive polarity direct P arrival; the receiver function amplitude does not exceed 1; any arrivals’ pulse length does not exceed 3.5 s. The latter two criteria are employed because very high amplitudes and long oscillatory pulses are typical characteristics of an unstable deconvolution (Schulte-Pelkum and Mahan, 2014). The calculated PRFs were then binned in 5° azimuthal groups with 5° overlap. In CHBR station, we recognized signs of the top (at 1 km depth) and bottom (at 9 km depth) of an anisotropic layer with almost north-south anisotropic symmetry axis. In addition, we recognized a flat interface beneath CHBR station at 27 km depth that is not in consistency with the result of migration to a depth of RFs showing a 2.5o dip Moho at the same place. For this reason, we utilize forward modelling to calculate synthetic PRFs to explain periodic amplitude variation of P to S converted phases with back-azimuths in each station that could be a signature for anisotropic velocity features. The forward modeling indicates that the horizontal interface makes a similar pattern on simulated PRfs as a low angle dipping interface with dip less than 10o.
Migration of S receiver functions reveals a deep velocity discontinuity at depth around 80 to 100 km that might be considered as a shallow lithosphere-asthenosphere boundary beneath the study region.

As seismic energy propagates through the earth medium, its energy (amplitude) decays due to geometrical spreading, intrinsic attenuation and scattering. Owing to anelastic absorption, intrinsic attenuation converts the seismic energy to heat while scattering redistributes the energy at random heterogeneities. Knowledge of the relative contributions of scattering and intrinsic attenuation is important for appropriate subsurface material identification, tectonic interpretations and quantification of the ground motion. Besides, investigating seismic wave attenuation inside lithosphere allows for a more thorough knowledge as to Earth’s deep structures. The attenuation of short-period S waves, expressed as the inverse of the quality factor (Q−1), helps fathom the physical laws related to the propagation of the elastic energy of an earthquake through the lithosphere. Coda wave attenuation is considered as the combination of scattering and anelastic attenuation. In this study, the quality factor of coda wave was estimated in NW Iran making use of single back scattering method of Aki and Chouet (1975). For this purpose, we analyzed 3720 waveforms recorded by 8 short-period stations of Tabriz network from 1996 to 2013. So as to calculate the frequency relationships for Qc, nine frequency bands with central frequencies of 1.5, 2, 3, 4, 6, 8, 12, 16 and 20 Hz were considered and the lateral and depth variations of Q0 (Qc in 1 Hz) were investigated in the research area. In order to study the lateral variations, we chose coda waves recorded in epicentral distances less than 80 km, in a lapse time window of 30 s. The reason for the selection of such short distance (

NW Iran is part of the complex tectonic system caused by the interaction between Arabian plate, Anatolia and Eurasia. The North Tabriz Fault (NTF) is one of the main structural features of the region and is considered to be the eastern termination of the Gailatu-Siah-Chesmeh-Khoy fault (Karakhanian et al., 2004), which merges with the Maku and the Nakhichevan dextral strike-slip faults and continues to move farther east. Part of the northward motion of Arabia is transferred to Anatolia by this complex system of faults and, the oblique orientation of the motion relative to the Zagros mountain range,results in the partitioning of the motion between shortening and thickening in the Caucasus and right-lateral strike-slip motion along the NTF.
In this research, we investigated the laterally two-dimensional velocity structure in the upper crust of NW Iran (mainly around the NTF) using local earthquake P-waves tomography. Several data sets were utilized, including Pg phase pickings of the Tabriz Network permanent stations governed by Institute of Geophysics, University of Tehran (1996 to 2013), temporary seismic stations installed around the North Tabriz Fault by International Institute of Earthquake Engineering and Seismology (IIEES) (April to July 2004) and temporary seismic stations installed by Institute for Advanced Studies in Basic Sciences (IASBS) (2009 to 2011) that merged with our data set so as to improve ray coverage in the eastern parts of the study area. The merged data set, recorded by 72 stations, consisted of more than 20,000 local earthquakes out of which, only 940 earthquakes were good enough to be selected for the local earthquake tomography. The velocity structures were resolved via a simultaneous solution of the coupled hypocenter and velocity model programmed in SIMULPS14. The tomographic images obtained from the linearized inversion are dependent on the initial velocity models and hypocenter locations. We primarily calculated the initial velocity model through the use of the 940 earthquake selected datasets. The time difference between the observed phase arrival time and predicted arrival time was then calculated and called travel time residuals. The residuals were further used as inputs for SIMULPS14 simulator to be converted into velocity model, which would in turn be used to adjust earthquake location parameters. Following four iterations for our inversion process, we obtained a 2D velocity tomogram that clearly showed different velocity structures on the two sides of the NTF. The velocity contrast across the NTF might have been caused by existence of different kinds of rocks on the two side of the fault trace. The North Tabriz Fault is an active and steep strike-slip fault generating strong structural differences around its surficial trace. It is a WNW–ESE trending fault in which the motion is concentrated on the fault at a rate of 7 mm/year. Such a strong rate of sliding explains the clear structure difference on the two sides of the fault. An anomalous low velocity feature can be seen in the central part of NTF. Comparing the velocity tomogram with the geological map of the region, one can observe that there exists a thick sediment basin in the same area. The low velocity anomaly is probably related to the thick, low velocity sediments deposited in that area.

The Zagros mountain belt, situated on the northern margin of the Arabian plate, is one of the youngest continental collision belts. This belt was formed by a collision between the Arabian plate and the Central Iranian micro-continent. In this study, we used data from 38 temporary seismological stations installed on a 400 km long profile from May to November 2003.The trend of the profile is N58°E across northern Zagros and part of the Central Iran. The stations are part of Zagros03 profile (Paul et al., 2010) operated by the International Institute of Earthquake Engineering and Seismology (IIEES) of Iran in collaboration with CNRS - Université Joseph Fourier, France. We examine the structure of the lithosphere, across the profile by analysis of P-wave receiver functions and Rayleigh wave fundamental mode phase velocity dispersion curves. Joint inversion of Rayleigh wave phase velocity dispersion and receiver functions have been used to estimate the velocity structure beneath 28 seismic stations. Receiver functions are time-series computed from the three-component body-wave seismograms and are sensitive to the earth structure near the receiver station. They are composed of P- to S-wave conversions in discontinuities under the stations. These converted waves are isolated by deconvolving the vertical component of a teleseismic P-wave record from its radial component. For each event, a 120s time-window centered at the direct P arrival is selected and used for the calculation of the receiver function. The deconvolution used is the iterative deconvolution method of Ligorria and Ammon (1999). Surface waves arise from the presence (boundary conditions) of the stress-free surface of the Earth, and in the presence of layering, they are dispersed. They provide valuable information on the absolute S-wave velocity, but they are relatively insensitive to sharp velocity contrasts. On the other hand, receiver functions are sensitive to S-wave velocity contrasts, which give rise to converted phases, but allow for a substantial trade-off between the depth and velocity above an impedance change. Combining them in a joint inversion process bridges the resolution gaps associated with each data set. We jointly inverted the stacked receiver function and surface wave dispersion data. We employ the program joint96 which is available in the software package “Computer Program in Seismology” (Herrmann and Ammon, 2003). In this study, we try to calculate the Moho depth and velocity structure in the north Zagros collision zone using the joint inversion of receiver function and surface wave dispersions. Receiver functions are calculated using teleseismic events of magnitude greater than 5.1, located between 30◦ and 95◦ epicentral distances. The fundamental mode Rayleigh-wave group velocities are extracted from thetomographic study conducted by Rahimi et al. (2014). The 1D velocity models resolved by joint inversion are juxtaposed, and a 2D velocity model is obtained. Results obtained from the 2D model reveal that the thickness of the sediments beneath the Zagros is 12 km, the Moho depth beneath this region of Zagros is 43-57 km, which increases towards Sanandaj–Sirjan zone and Urumieh–Dokhtar magmatic arc and reaches an expanse of 62 km and then decreases in the central Iran with a depth of 42 km. The velocity model confirms the presence of a crustal root and a thick high-velocity lithosphere beneath and north of the suture. These evidences imply that the Arabian plate continues to underthrust beneath the Central Iran.
The Zagros mountain belt, situated on the northern margin of the Arabian plate, is one of the youngest continental collision belts. This belt was formed by a collision between the Arabian plate and the Central Iranian micro-continent. In this study, we used data from 38 temporary seismological stations installed on a 400 km long profile from May to November 2003.The trend of the profile is N58°E across northern Zagros and part of the Central Iran. The stations are part of Zagros03 profile (Paul et al., 2010) operated by the International Institute of Earthquake Engineering and Seismology (IIEES) of Iran in collaboration with CNRS - Université Joseph Fourier, France. We examine the structure of the lithosphere, across the profile by analysis of P-wave receiver functions and Rayleigh wave fundamental mode phase velocity dispersion curves. Joint inversion of Rayleigh wave phase velocity dispersion and receiver functions have been used to estimate the velocity structure beneath 28 seismic stations. Receiver functions are time-series computed from the three-component body-wave seismograms and are sensitive to the earth structure near the receiver station. They are composed of P- to S-wave conversions in discontinuities under the stations. These converted waves are isolated by deconvolving the vertical component of a teleseismic P-wave record from its radial component. For each event, a 120s time-window centered at the direct P arrival is selected and used for the calculation of the receiver function. The deconvolution used is the iterative deconvolution method of Ligorria and Ammon (1999). Surface waves arise from the presence (boundary conditions) of the stress-free surface of the Earth, and in the presence of layering, they are dispersed. They provide valuable information on the absolute S-wave velocity, but they are relatively insensitive to sharp velocity contrasts. On the other hand, receiver functions are sensitive to S-wave velocity contrasts, which give rise to converted phases, but allow for a substantial trade-off between the depth and velocity above an impedance change. Combining them in a joint inversion process bridges the resolution gaps associated with each data set. We jointly inverted the stacked receiver function and surface wave dispersion data. We employ the program joint96 which is available in the software package “Computer Program in Seismology” (Herrmann and Ammon, 2003). In this study, we try to calculate the Moho depth and velocity structure in the north Zagros collision zone using the joint inversion of receiver function and surface wave dispersions. Receiver functions are calculated using teleseismic events of magnitude greater than 5.1, located between 30◦ and 95◦ epicentral distances. The fundamental mode Rayleigh-wave group velocities are extracted from thetomographic study conducted by Rahimi et al. (2014). The 1D velocity models resolved by joint inversion are juxtaposed, and a 2D velocity model is obtained. Results obtained from the 2D model reveal that the thickness of the sediments beneath the Zagros is 12 km, the Moho depth beneath this region of Zagros is 43-57 km, which increases towards Sanandaj–Sirjan zone and Urumieh–Dokhtar magmatic arc and reaches an expanse of 62 km and then decreases in the central Iran with a depth of 42 km. The velocity model confirms the presence of a crustal root and a thick high-velocity lithosphere beneath and north of the suture. These evidences imply that the Arabian plate continues to underthrust beneath the Central Iran.

Estimation of seismic wave attenuation due to anelasticity and geometrical spreading has attracted major interests among earthquake engineering community in recent decades. The choice of ground-motion model has a significant impact on hazard estimates in an active seismic zone such as the NW Iran. Estimation of ground motion for a typical frequency range of 0.5–10 Hz is required for the proper design of earthquake resistant structures and facilities and is considered as input for engineering stochastic ground motion relationships. For seismological purposes, appropriate attenuation models make it possible to calculate more accurately the source parameters such as magnitude and seismic moment. The NW Iran has experienced very few large events during the operation of the accelerometer network of the Building and Housing Research Center (BHRC). The BHRC network has been operating since 1973 but has recorded ground acceleration for few events in the study area, because of the low seismicity rate. The availability of the abundant weak-motion waveform data from the short-period local seismograph network of the Institute of Geophysics of the University of Tehran (IGUT) provides an opportunity to derive a new and more reliable ground-motion relationship for small events to complement those of strong-motion results. In this study, we analysed 3514 records of 943 small and moderate events that were recorded by 8 permanent stations of Tabriz network (belonging to the IGUT) and 16 temporary stations of the Institute for Advanced Studies in Basic Sciences (IASBS) to prepare a dataset including week ground-motion spectral amplitudes for different magnitudes and hypocentral distances. We graphically found the distance at which the nature of geometrical spreading attenuation changes significantly using a locally weighted scatter-plot smoothing called robust LOWESS. A bilinear function with a hinge at distance of about 70 km describes the geometric spreading attenuation with distance. Geometrical spreading and intrinsic attenuation coefficients were calculated using nonlinear regression in different frequencies and an average value of was found as geometrical spreading coefficient for distance range of 10–70 km. This value is consistent with geometrical spreading in a layered Earth. The average geometrical spreading coefficient of was found for the frequency range 0.79–5 Hz and the distance range of 70–200 km. This value is smaller than the values reported for other regions in the world (e.g. .09 for Central Alborz: Motaghi and Ghods, 2012; .2 for North Iran: Motazedian, 2006; .2 for SE Canada and the NE United States: Atkinson, 2004; .1 for SE Australia: Allen et al., 2007) and indicates that the velocity contrast in the Moho discontinuity is smaller than that in the other regions. The low-velocity uppermost mantle in NW Iran was manifested by different types of tomographic results obtained for the region. The geometrical spreading coefficient does not change before and after 70 km distance for frequencies ≥ 5 Hz. Thus, the attenuation relationship in this frequency range changed from bilinear to linear function. Using anelastic attenuation coefficients calculated at different frequencies, the shear-wave quality factor, , obtained equal to for frequencies greater than 1.5 Hz. In fact, the values show a U-shaped behavior in all of the frequency ranges and the function that describes it is defined as .
Estimation of seismic wave attenuation due to anelasticity and geometrical spreading has attracted major interests among earthquake engineering community in recent decades. The choice of ground-motion model has a significant impact on hazard estimates in an active seismic zone such as the NW Iran. Estimation of ground motion for a typical frequency range of 0.5–10 Hz is required for the proper design of earthquake resistant structures and facilities and is considered as input for engineering stochastic ground motion relationships. For seismological purposes, appropriate attenuation models make it possible to calculate more accurately the source parameters such as magnitude and seismic moment. The NW Iran has experienced very few large events during the operation of the accelerometer network of the Building and Housing Research Center (BHRC). The BHRC network has been operating since 1973 but has recorded ground acceleration for few events in the study area, because of the low seismicity rate. The availability of the abundant weak-motion waveform data from the short-period local seismograph network of the Institute of Geophysics of the University of Tehran (IGUT) provides an opportunity to derive a new and more reliable ground-motion relationship for small events to complement those of strong-motion results. In this study, we analysed 3514 records of 943 small and moderate events that were recorded by 8 permanent stations of Tabriz network (belonging to the IGUT) and 16 temporary stations of the Institute for Advanced Studies in Basic Sciences (IASBS) to prepare a dataset including week ground-motion spectral amplitudes for different magnitudes and hypocentral distances. We graphically found the distance at which the nature of geometrical spreading attenuation changes significantly using a locally weighted scatter-plot smoothing called robust LOWESS. A bilinear function with a hinge at distance of about 70 km describes the geometric spreading attenuation with distance. Geometrical spreading and intrinsic attenuation coefficients were calculated using nonlinear regression in different frequencies and an average value of was found as geometrical spreading coefficient for distance range of 10–70 km. This value is consistent with geometrical spreading in a layered Earth. The average geometrical spreading coefficient of was found for the frequency range 0.79–5 Hz and the distance range of 70–200 km. This value is smaller than the values reported for other regions in the world (e.g. .09 for Central Alborz: Motaghi and Ghods, 2012; .2 for North Iran: Motazedian, 2006; .2 for SE Canada and the NE United States: Atkinson, 2004; .1 for SE Australia: Allen et al., 2007) and indicates that the velocity contrast in the Moho discontinuity is smaller than that in the other regions. The low-velocity uppermost mantle in NW Iran was manifested by different types of tomographic results obtained for the region. The geometrical spreading coefficient does not change before and after 70 km distance for frequencies ≥ 5 Hz. Thus, the attenuation relationship in this frequency range changed from bilinear to linear function. Using anelastic attenuation coefficients calculated at different frequencies, the shear-wave quality factor, , obtained equal to for frequencies greater than 1.5 Hz. In fact, the values show a U-shaped behavior in all of the frequency ranges and the function that describes it is defined as .

Reliable determination of earthquake moment magnitude is a fundamental problem of seismic hazard assessment. Short period data does not generally extend to sufficientlylow frequencies to allow for the reliable calculation of moment tensor. In the Tehran region, the available seismic data are mostly short-period records. To resolve this difficulty, following Motazedian and Atkinson (2005), the current study uses M1 magnitude, which closely follows moment magnitude for smallto moderate events. In addition, using this special magnitude scale aids in overcoming the difficulty with the different magnitude scales in the catalogue of the Tehran region. The M1 magnitude scale is obtainedfrom the spectral amplitude at 1 Hz. The Eigene-frequency of short period seismometers is always close to 1 Hz, thus M1 can be determined using data from short-period seismometers. Assuming a Brune point-source model, M1 magnitude is equal to the moment magnitude. To calculate M1 magnitude for a given event, the Fourier power spectrum of the S-window of each observed transverse component waveform is calculated. It is then corrected for geometrical and intrinsic attenuation using the attenuation relationship of Mottaghi (2007), developed for the Tehran region, and the applied to a Butterworth filter centered on 1 Hz to calculate spectrum amplitude at 1 Hz,. Next, using a trial and error method, M1 is calculated by fitting the to the corresponding synthetic Brune spectrum amplitudes. The fit is done over a frequency band of 0.4 to 3 Hz while the band is divided to equal logarithmic bins. A shear velocity of 3.73 km/s has been assumed to calculate the synthetic Brune spectrum, along with an average density of g/cm3 and a stress drop of 10 MPa. The fit is strongly sensitive to 0.1 magnitude unit fluctuations. Magnitude of an event is defined by the average of the calculated values of M1 over all stations that recorded the event. Since 1995, the seismicity of the Tehran region has been monitored by Iran Telemetered Seismograph Network (ITSN), a series of small regional subnetworks operated by the Institute of Geophysics, University of Tehran. A review of the operation of the Tehran subnetwork of the ITSN has been provided by Ghods and Sobouti (2005). The data used in this study are1804 records of 179 earthquakes having magnitudes larger than 3 and occurring in the period of 1996-2004. To ensure reasonable location accuracy and lower sensitivity on radiation pattern, the selected events have an azimuthal gap lower than 250 degree. The events were selected to provide relatively homogenous ray coverage inside the study area [Figure 1]. The selected records have hypocentral distribution in the range of 10 to 440 km with rather dense coverage of hypocentral distance of 350 km (Figure 2) The result of this study is a catalogue of 179 events with magnitude M1 greater than 3 which occurred in the Tehran region. Variations of M1 residues versus epicentral distance show insignificant dependency with distance [Figure 5], implying that the attenuation relationship used in this study (Motaghi, 2007) is appropriate for the Tehran region. A comparison of calculated M1 catalog with the corresponding ML catalog (Figure 4) shows a systematic difference between these two scales following a relationshipof. The fit is computed using Least Absolute Residual algorithm and has a RMS of 0.111.