CC BY
21
RUSSIAN JOURNAL OF EARTH SCIENCES, VOL. 8, ES1003, doi:10.2205/2006ES000194, 2006
Crustal velocity structure under the RUKSA seismic array (Karelia, Russia)
I. M. Aleshin, G. L. Kosarev, and O. Yu. Riznichenko
Institute of Physics of the Earth, Russian Academy of Sciences, Moscow, Russia I. A. Sanina
Institute for Dynamics of Geospheres, Russian Academy of Sciences, Moscow, Russia Received 15 February 2006; revised 1 March 2006; accepted 9 March 2006; published 14 March 2006.
[1] Based on data of three three-component seismographs belonging to the temporary small-aperture Russian Karelia Seismic Array (RUKSA) in the Petrozavodsk region (Karelia), a 1-D velocity model of the crust is constructed by the method of the receiver function.
Waveforms of distant earthquakes recorded by short-period instruments with improved characteristics are used. The data were inverted by the simulated annealing method. The inversion was stabilized by using phase velocities of Rayleigh waves and traveltimes of converted Ps waves from the 410-km boundary determined from broadband records of the SVEKALAPKO seismic array. Anomalously low seismic velocities are discovered in the upper part of the cross section beneath the RUKSA array. Index terms: 7205 Seismology:
Continental crust; 7255 Seismology: Surface waves and free oscillations; 7294 Seismology: Seismic instruments and networks; KEYWORDS: Earth’s crust, small aperture array, method of the receiver function.
Citation: Aleshin, I. M., G. L. Kosarev, O. Yu. Riznichenko, and I. A. Sanina (2006), Crustal velocity structure under the RUKSA seismic array (Karelia, Russia), Russ. J. Earth. Sci., 8, ES1003, doi:10.2205/2006ES000194.
Introduction
[2] Observations with the small aperture RUKSA array (Figure 1) were carried out in 1999 in the Petrozavodsk region (Karelia, Russia) within the framework of the SVEKALAPKO international experiment [Bock, 2001]. Small aperture seismic arrays (sometimes referred to as seismic antennas) have been widely used since the 1980s as an effective tool for the location of seismic events of various origins, primarily, for the monitoring of the Non-Proliferation Treaty. Data of seismic antennas are also effectively used for the study of seismic noise properties and seismic wave scattering, as well as for the recording and identification of local seismic events. One of the immediate tasks involved in the analysis of array records is the determination of the structure of the medium in the area of the temporary array site and, in particular, the construction of a 1-D velocity model of the crust and the upper mantle. The method of the receiver function [Langston, 1979; Vinnik, 1977; Vinnik and Kosarev, 1981] developed for the analysis of three-component seismic records is effective for the detection of main interfaces or zones of higher velocity gradients but
Copyright 2006 by the Russian Journal of Earth Sciences.
ISSN: 1681-1208 (online)
fails to reliably determine absolute values of velocities in model layers. The inverse problem of the velocity structure determination from a known receiver function can be solved by invoking additional information on the parameters to be determined [Ammon, 1990; Kosarev et al., 1987]. To reduce the ambiguity of the inversion, Julia et al. [2000] performed joint inversion of the receiver function and group velocities of Rayleigh waves in the Arabia region. Using joint inversion of P and S wave receiver functions, Vinnik et al. [2004] constructed a 3-D model of the crust and upper mantle (down to 150 km) under the Tien Shan. Residuals of P and S traveltimes were included in the inversion in (L. P. Vinnik et al., in press, 2006). In the present paper, Ps traveltimes from the 410- and 660-km boundaries and Rayleigh phase velocities are used as additional optimization constraints.
Observation Conditions and Initial Data
[3] The temporary small aperture RUKSA array (Institute of Physics of the Earth, Russian Academy of Sciences (IPE RAS)) was located in the Hautavaara settlement area (Karelia) 100 km west of the town of Petrozavodsk, between lakes Ladoga and Onega. The array observations (Figure 1) were conducted by a group of IPE RAS reserchears headed
Figure 1. Map showing the position of the SVEKALAPKO and RUKSA stations. The configuration of the RUKSA array is shown in the inset (the triangles are three-component stations, and the circles are stations recording the vertical component alone). The diameter of the outer ring (the stations C1-C5) is 640 m.
by M. V. Nevsky and S. G. Volosov within the framework of the large-scale SVEKALAPKO seismic tomography experiment for the study of the deep structure and evolution of the Baltic Shield. The array consisted of З three-component and 6 one-component (Z) instruments (SM-3KV seismometers) integrated with the digital recording system of the Ekspress
station. The seismometers were positioned symmetrically on concentric circles of diameters of 640 m (outer ring C) and 300 m (inner ring B). The dynamic range of the recording instrumentation is 120 dB, and the recording sampling rate is 100 Hz. Time control was performed via GPS. The RUKSA recording channel after extending its frequency response to-
Figure 2. The observed receiver functions. Curves 1-12 are receiver functions summarized over the stations C1, C3, and C5 for earthquakes with the corresponding number (see the Table 1). Curves C1, C3, C5, and C1+C3+C5 are receiver functions summarized over all earthquakes for the corresponding stations.
ward lower frequencies [Bashilov et al., 1985] ensures reliable recording in the range 0.5-20.0 Hz.
[4] The location of the array is characterized by a rather low level of microseismic noise. The noise spectrum includes stationary components of anthropogenic origin at frequencies of 12.5 Hz and 25 Hz and an intense episodic component with a maximum at 5 Hz. Upon examining the level and spa-tiotemporal variations in seismic noise at the RUKSA array, a frequency range of 0.5-3.0 Hz was chosen for the subsequent analysis.
[5] Over 30 days of observations, the array recorded 100 seismic events of various origins including the largest earthquakes of 1999 in Turkey (Izmit), Taiwan (Chi-Chi), and Greece. This work uses records of the 12 strongest tele-seismic events of magnitudes mb of 5.0 to 7.7 at epicentral distances of 20° to 70° (Table 1).
Analysis of the Receiver Function Obtained from RUKSA Records
[6] The procedure used for the construction and analysis of the receiver functions from P wave records has been repeatedly described in literature. In this work, we follow mainly the scheme described in [Kind et al., 1995]. The processing of P wave seismograms was performed with the use of the Seismic Handler software package [Stammler, 1993] and included the following operations.
[7] (1) The observed vertical (Z) and two horizontal (NS and EW) components of a P wave were rotated to the ray coordinate system (L, Q, and T). The L axis lies in the ray plane and is directed from the source along the main motion in the P wave. The Q axis lies in the same plane and is perpendicular to the L axis. The horizontal projection of the Q axis is positive in the direction from the source. Together with the L and Q axes, the T axis forms a right-hand orthogonal triple. Displacements along the Q axis depend very weakly on the main motion in the incident P wave, determined mostly by the earthquake source, and are mainly composed of multiple and converted waves arising beneath the station.
[8] (2) For applying the stacking procedure to a signal in order to increase signal-to-noise ratio, records of P waves from different earthquakes must be reduced to a standard source of the impulsive type. To do this, the L component of the signal, considered in this case as a source, is transformed into the standard form with the help of a deconvolution filter. The same filter is then applied to the components Q and T. It is the component Q that is called the receiver function. After this transformation, the noise can be suppressed by stacking Q components of all (or some) events, which yields the averaged receiver function Qobs (t).
[9] Figure 2 presents Q components of all processed waveforms from 12 events (see the Table 1). The microseismic noise level can be estimated in the time interval from 5 s
Table 1. Parameters of earthquakes
No. Date Time Latitude Longitude Distance Azimuth Depth (km) Region
1 7 Sep. 1999 1156:49.3 38.1 23.6 24.6 197 10 Greece
2 19 Sep. 1999 1647:00.1 43.2 46.9 20.8 149 43 Caucasus
3 20 Sep. 1999 0932:42.7 46.3 153.5 62.1 42 33 Kurile Islands
4 20 Sep. 1999 1747:18.4 23.8 121.0 68.6 79 33 Taiwan
5 20 Sep. 1999 1803:44.2 23.6 121.3 69.0 79 33 Taiwan
6 20 Sep. 1999 1811:53.6 23.8 121.2 68.8 79 33 Taiwan
7 20 Sep. 1999 2146:41.8 23.4 121.0 69.0 80 33 Taiwan
8 22 Sep. 1999 0014:39.1 23.7 121.2 68.8 79 26 Taiwan
9 22 Sep. 1999 0049:42.7 23.6 121.1 68.8 79 33 Taiwan
10 24 Sep. 1999 1917:14.5 28.6 51.3 35.7 151 33 Iran
11 25 Sep. 1999 2352:48.6 23.7 121.2 68.7 79 17 Taiwan
12 28 Sep. 1999 0500:42.9 54.6 168.3 58.6 28 33 Komandorski Islands
to -1s before a first arrival (the time t = 0). Differences between records of the stations C1, C3, and C5 are small, except for systematically higher amplitudes at the station C5 as compared to C1 and C3 in the interval from 0 s to 1 s. Waves that experienced conversion, refraction, and scattering in the upper few kilometers under the station are observed precisely in this time interval. The higher C5 amplitudes are well seen on the average C5 trace in Figure 2. The signal-to-noise ratio on average and even some individual traces is large enough to enable further analysis and formulate the inversion problem for the velocity structure determination.
[10] The inferred receiver functions display two main peaks at delay times of 0.27 s and 5.1 s. These peaks were initially identified as P to S conversions at an upper crust boundary and at the Moho. Particle motion analysis showed that oscillations at delay times longer than 10 s are associated with arrivals of scattered surface waves having a characteristic elliptic polarization rather than with multiple reflections of body waves. For this reason, we chose the time interval (-2 s, 8 s) for the receiver function used below for reconstructing the velocity structure of the medium.
Surface Wave Data
[11] As mentioned above, additional observational data are required to stabilize the receiver function inversion results. Phase velocities of surface waves have recently been utilized for this purpose [Julia et al., 2000]. Such data are unavailable from RUKSA records due to the overly small aperture. However, the phase velocities were determined from SVEKALAPKO (without RUKSA) broadband records of many earthquakes [Bruneton et al., 2004]. The method consists in the measurement of traveltimes of the Rayleigh wave fundamental mode between their curvilinear fronts within the SVEKALAPKO array and the inversion of the measured data into the phase velocity field. These data in the form of dispersion curves of the Rayleigh phase velocity at grid nodes covering the SVEKALAPKO array were af-
forded by M. Bruneton, who used them for the construction of the SVEKALAPKO tomography model [Bruneton et al., 2004]. A phase velocity curve in the range of periods from 10.5 s to 190 s was constructed at each grid point. We used the dispersion curve that was obtained by extrapolating the set of these curves (Figure 3) to the point with the RUKSA coordinates by the kriging method [Isaacs and Srivastava, 1989]. In our opinion, this procedure is reliable because the spread of the dispersion curves in the region about 500 km across is small (Figure 3) and the distance from the RUKSA array to the nearest SVEKALAPKO stations does not exceed 100 km. The inversion interval of periods was 10.5110 s.
Traveltime of the Converted Wave From the 410-km Boundary
[12] The SVEKALAPKO average traveltime of the converted wave from the boundary at a depth of 410 km tps = Ps410 — P = 42.3 ± 0.1 s was another observational parameter that was used as a constraint on the velocity structure beneath RUKSA. The measurements were made by the technique proposed for the first time by Vinnik [1977]. The SVEKALAPKO traveltimes Ps660 — P = 66.2 s and Ps410 — P = 42.3 s are smaller than the corresponding values for the standard Earth by 1.8 s. The traveltime Ps660 — Ps410 = 23.9 s coincides with the value of the IASP91 standard Earth. Hence, in the absence of Vp/Vs anomalies, seismic velocities in the SVEKALAPKO area coincide with standard values at depths of 410 km to 660 km and exceed them at depths of 0 km to 410 km. Most earthquake epicenters lie east of SVEKALAPKO, implying that the average Ps410 path from the boundary to the Earth’s surface is shifted toward RUKSA by about 200 km. If the size of the Fresnel zone is additionally taken into account, it becomes clear that the region affecting the propagation of the converted waves includes the RUKSA zone, but the exact trace corresponding to the 42.3-s traveltime cannot be determined. Apparently, this time can be considered as an
Figure З. Dispersion curves (solid lines) of Rayleigh phase velocities from SVEKALAPKO data [Bruneton et al., 2004]. The circles are the velocities extrapolated to the RUKSA point.
average for the SVEKALAPKO+RUKSA territory, i.e. for the SE Baltic Shield.
Construction of a Velocity Structure by the Simulated Annealing Method
[13] We modeled the medium under the array by a set of horizontal layers underlain by a homogeneous half-space. The model parameters m included P and S wave velocities, density, and layer thicknesses. The model was constructed through the minimization of the functional
E(m) = w1
+ W2
t2 — t1
t2
^^(Qobs(t) Qsyn(m,t)) dt
1/2
^2 — ^1
Ш2
^ (Robs(^) ^^syn (Ї-И, t))
1/2
+ W3
A synthetic receiver function Qsyn(m,t) was calculated by the formula [Kind et al., 1995]
Qsyn(m,t) = T1 [ Hq(m,^) Lobs(^)exp(iwt)dw .
2ж J
The response components of the layered medium Hq(m,w) and HL(m,w) were calculated by the Thomson-Haskell
method [Haskell, 1962]. (R(w) is the surface wave dispersion curve, and wa are empirically chosen weights.)
[14] The minimization problem formulated above is strongly nonlinear and ill-conditioned. The traditional method used for solving problems of this type is based on the regularization method [Tikhonov and Arsenin, 1979]. The minimization reduces to the solution of a system of linear equations that can be found very rapidly. However, the model obtained from an implementation of this method [Kosarev et al., 1987] is essentially dependent on the initial approximation whose choice is generally complicated.
[15] Statistical (nongradient) methods (primarily, the Monte Carlo method) “scan” the entire space of parameters in the process of solution and are free from the necessity of specifying the initial approximation. However, the direct application of statistical methods to the search for the solution in a multidimensional space involves time consuming computations, which makes the practical use of these methods unrealistic. A possible way of solving this problem is to seek the solution in three stages: identification of one or more promising regions in the initial space of potential models, careful examination of each of these regions, and “fine adjustment” of the best solution inferred at the preceding stage. Such an approach is realized in the modern modification of the simulated annealing (SA) method [Ingber, 1989]. The software implementation of this method reduces the time of the search for the solution to an acceptable value (no more than one hour on a standard PC for the present case).
[16] To invert data, we used a model consisting of 11 layers on a half-space. Velocities of S waves in the layers and the half-space and thicknesses of some layers were varied, whereas the ratio Vp/Vs was fixed (1.73 in the crust and 1.8 in the mantle) and the density was determined by the Birch formula. The chosen time interval of the receiver function
1
1
obs syn
t410 — t410
Figure 4. Results of the joint inversion of receiver functions, surface wave data, and traveltimes of converted waves from boundaries of the phase transition zone in the mantle. The top diagrams are the Vp and Vs velocity structures beneath the RUKSA array. The thin red lines bound the range of velocity variations in the inversion process. The middle curves are the observed (red) and synthetic (blue) receiver functions. The observed receiver function corresponds to curve C1+C3+C5 in Figure 2. The bottom plots show the observed (red circles) and calculated (blue line) phase velocities of Rayleigh waves. The traveltime of the converted Ps wave from the 410-km boundary is tps410 = 42.3 s (the ray parameter is 6.4 s/°.
minimization (—2 s to 8 s) bounds the overall depth of the model by about 70-80 km. Deeper layers have no effect on the receiver function. However, in order to use traveltimes of waves converted at the 410-km boundary, the velocity structure should be specified to a depth of 410 km. For this purpose, our current 11-layer model (more specifically, the line defining the velocity in the half-space) was continued downward until its intersection with the IASP91 standard model of the Earth. IASP91 velocities were used for the calculation of tps410 in the depth interval from the intersection point to 410 km.
[17] The SA method does not require the specification of an initial approximation. However, it is necessary to de-
fine intervals within which the variable parameters can be changed. The resulting model must lie within these intervals. The variation limits of P and S wave velocities and the inversion results are shown in Figure 4. The P velocity variation limits were chosen in such a way that they included the whole range of velocity values published in the recent literature for the Baltic Shield. For example, the P velocity in the layer at depths of 10-15 km (the first layer of the middle crust) was varied from 6.0 km s-1 to 6.4 km s-1 in agreement with the seismic velocity structures obtained for this region [Goncharov et al., 1991]. The upper part of the model consists of six layers. The thickness of the two uppermost layers was varied from 0.1 km to 1 km, whereas
the thickness of the four remaining layers was fixed (2 km). Such a detailed specification of the upper part of the model was due to the necessity of fitting the intense phase of the first second. The remaining part of the crust includes five layers of variable thicknesses that could change from 5 km to 10 km. the mantle is described by a half-space. During optimization, 6000 to 20,000 models were generated.
[18] To estimate the effect of each type of the observed data on the resulting model, we performed the following experiment (Figure 5). The receiver function alone was first inverted. As a result, we obtained a model shown by the red line in Figure 5. Then the receiver function was supplemented by either surface wave dispersion data (blue line) or the traveltime of converted waves from the 410-km boundary (green line). Finally, the fourth curve (black line) represents a model obtained by using the whole data set. Comparison of models reveals that the data of the receiver function alone underestimate the velocities in both the crust and the mantle. The receiver function data supplemented by either the data of surface waves alone or by the tps410 traveltime data alone yield virtually the same and, apparently, correct velocities in the middle crust at depths of 8 km to 20 km. The surface wave data underestimate velocities in the lower crust and overestimate the mantle velocities. On the contrary, the tps410 data overestimate the lower crust velocities but yield “correct” velocities in the mantle. By correct velocities, we tentatively mean here the velocities obtained from the set of all data (the black line).
Discussion and Conclusion
[19] One of the most noticeable features of the inferred structures is anomalously low velocities in the upper crust. Thus, Vp ~ 4.1 km s-1 and Vs ~ 2.4 km s-1 in the uppermost layer about 1 km thick. Figure 4 shows the results of the average receiver function inversion from data of 12 earthquakes recorded at the stations C1, C3, and C5. If the receiver function is inverted from the data of the station C5 alone (Figure 2), the inverted velocities in the upper layer are even lower: Vp = 3.5 km s-1 and Vs =2 km s-1. Apparently, these are the lowest velocities that have ever been observed in the upper crust of the Karelia region. In any case, such low velocities are absent in the velocity structures presented in [Sharov, 2004]. The velocity Vp ~ 5.5 km s-1 is attained at a depth of about 3 km, and Vp ~ 6.2 km s-1, only at a depth of about 10 km. Low velocities in the upper hundred meters are not surprising. All RUKSA stations were mounted on moraine deposits. However, anomalously low velocities are sometimes also observed in bedrocks of the upper crust [Vinnik et al., 2002]. Grad and Luosto [1992] reported relatively low velocities in the upper crust from the SVEKA profile in Finland. The velocities were determined from observations of short-period (0.5-2.5 Hz) Rayleigh waves. The lowest surface wave velocity is 2.7 km s-1, which yields an S wave velocity of 2.9 km s-1. Based on the analysis of O ’Connel and Budiansky, 1977], Grad and Luosto attributed the relatively low velocities to the fracturing of upper crust rocks provided that the cracks are isolated and filled with water. Our observations yield much lower velocities that
Figure 5. Comparison of velocity structures obtained by data inversion of the receiver function alone (red line), receiver function + Rayleigh waves (blue line), receiver function + tps410 (green line), and receiver function + Rayleigh waves + tps410 (black line).
can also be accounted for by the presence of cracks filled with water except that all or some part of the cracks are interconnected. According to the calculations of O’Connel and Budiansky, an S wave velocity decrease in the case of isolated cracks reaches 25% compared to a monolithic rock. The concentration of cracks was assumed to be 0.5. If the cracks are not isolated, the velocity decrease is 60%. In our case, the observed S velocities amount to ~ 2.4 km s-1, whereas a value of 3.6 km s-1 is expected for rocks of the granite type. The velocity decrease amounts to 40% and, as distinct from the SVEKA profile observations, can accounted for by partially interconnected cracks. This inference could be supported (or disproved) by estimating the crack density from measurements of the number and sizes of fractures in the RUKSA area. Vinnik et al. [2002] noted a 30% decrease in S velocities at some points of the upper crystalline crust in central Tien Shan and interpreted this decrease in terms of the fracturing in the shallow layer. Our data are insufficient to assess to what extent the velocity decrease is affected by moraine deposits, the isolation or interconnectedness of cracks, and the possible presence of particularly low velocity rocks in the volcanogenic sedimentary series. All three factors appear to be significant.
[20] In the lower and upper crust, at depths of 10 km to 40 km, the P wave velocity slowly rises to about 6.6 km s-1. Its average gradient at these depths is an order of magnitude smaller than in the upper crust at depths of 0 km to 10 km. The Moho depth is about 40 km. The P velocity immediately under the Moho is 8.4 km s-1. This value
is much higher than in the IASP91 model but is consistent with modern models of ancient shields. Thus, the velocity structure obtained in this work consists of an essentially heterogeneous upper part where the velocity (converted into P wave velocity values) rises from 3.5 km s-1 at the surface to 6.2 km s-1 at depths of 10-12 km and nearly homogeneous middle and lower parts with minimal velocity contrasts in the corresponding layers. Velocity inversion zones are virtually absent in the final variant of the velocity structure.
[21] Previously we presented and discussed the velocity structure beneath RUKSA obtained from a receiver function alone, without invoking data on surface waves and travel-times of converted waves from the 410-km boundary [Sharov, 2004]. The main features of the structure, the Moho depth and average velocities in the crust and the mantle, are close to those inferred in our study. Distinctions are mainly observed in the upper part of the section and can be attributed to the instability of the inversion using data of the receiver function alone. The extent of misfit (the rms deviation) between the synthetic and observed waveforms of the published structure is 1.5 times greater than in our present work.
[22] Reliable features of the presented models (Figure 4) are the presence of a sharp boundary at a depth of about 40 km (the Moho), anomalously low velocities in the depth interval 0-3 km that are untypical of ancient shields, a rapid rise in the velocity with depth in the upper crust, and weak differentiation in velocities of the middle and lower crust. Significant heterogeneities in the upper part of the crystalline crust to depths of 12-15 km were previously noted in studies of the fine structure of seismic wavefields from small aperture array data [Nevsky and Riznichenko, 1980].
[23] Thus, using data of short-period instrumentation with improved characteristics, a detailed velocity structure is obtained beneath the RUKSA small aperture array. These data are necessary for remote location of seismic events. They can also be used for the construction of a 3-D lithosphere model of the European part of the Russian Federation. This approach is promising if mobile small-aperture arrays and natural sources (remote earthquakes) are used for local mapping of crustal boundaries from seismic data of relatively short period instrumentation.
[24] Acknowledgments. We are grateful to M. Bruneton, who kindly provided us with data on phase velocities of surface waves. Data on tps410 traveltimes from SVEKALAPKO records were obtained by I. M. Aleshin and G. L. Kosarev during their visit to Oulu University at the invitation of E. Kozlovskaya (research grant of the Finnish Academy of Sciences, 2004). Experimental observations with the RUKSA array were supported by the INTAS grant 97-0936 and the Zurich ETH Institute (Switzerland). This work was supported by the Russian Foundation for Basic Research, projects nos. 03-05-64654, 04-05-64634, and 04-07-90362-B.
References
Ammon, C. J. (1990), On the nonuniqueness of receiver function inversions, J. Geophys. Res., 95, 2504.
Bashilov, I. P., T. N. Ershova, M. M. Krekov, M. V. Nevskiy, and A. V. Nikolaev (1985), Experimental digital seismic
array, in Development and Testing of Seismic Instrumentation, Seismic Instruments, issue 17 (in Russian), p. 139, Nauka, Moscow.
Bock, G. (2001), Seismic probing of Fennoscandian lithosphere, EOS Trans. AGU, 82, 50.
Bruneton, M. (2004), Complex lithospheric structure under the central Baltic Shield from surface wave tomography, J. Geophys. Res., 109, B10303, doi:10.1029/2003JB002947.
Goncharov, A. G., K. A. Kalnin, S. B. Lobach-Zuchenko, M. D. Lizinsky, L. N. Platonenkova, and V. P. Chekulaev (1991), Seismogeological characteristic of the Karelia crust, in Problems of Joint Interpretation of Geological and Geophysical Data, edited by V. A. Glebovitsky and N. N. Sharov (in Russian), p. 53, Nauka, Leningrad.
Grad, M., and U. Luosto (1992), Fracturing of the crystalline uppermost crust beneath the SVEKA profile in central Finland, Geophysica, 28(1—2), 53.
Haskell, N. A. (1962), Crustal reflection of plane P and SV waves, J. Geophys. Res., 67, 4751.
Ingber, L. (1989), Very fast simulated annealing, Math. Comput. Model., 12, 967.
Isaacs, E. H., and R. M. Srivastava (1989), Applied Geostatistics, 80 pp., Oxford Univ., New York.
Julia, J., C. J. Ammon, R. B. Herrmann, and Am. M. Correig (2000), Joint inversion of receiver function and surface-wave dispersion observations, Geophys. J. Int., 143, 99.
Kind, R., G. L. Kosarev, and N. V. Petersen (1995), Receiver functions at the stations of the German regional seismic network (GRSN), Geophys J. Int., 121, 191.
Kosarev, G. L., L. I. Makeeva, and L. P. Vinnik (1987), Inversion of teleseismic P-wave particle motions for crustal structure in Fennoscandia, Phys. Earth Planet. Inter., 47, 11.
Langston, C. A. (1979), Structure under Mount Rainier, Washington, inferred from teleseismic body waves, J. Geophys. Res., 84, 4749.
Nevsky, M. V., and O. Yu. Riznichenko (1980), Scattering of seismic waves in the crust from data of seismic arrays, Izv. Akad. Nauk SSSR, Fiz. Zemli (in Russian), 16(6), 26.
O’Connel, R. J., and B. Budiansky (1977), Viscoelastic properties of fluid-saturated cracked solids, J. Geophys. Res., 79, 5412.
Sharov, N. V., (Ed.) (2004), Deep Structure and Seismicity of the Karelia Region and Adjacent Areas (in Russian), 97 pp., KNTs RAN, Petrozavodsk.
Stammler, K. (1993), Seismic handler: Programmable multichannel data handler for interactive and automatic processing of seismic analyses, Comput. Geosci., 19, 135.
Tikhonov, A. N., and V. Yu. Arsenin (1979), Methods of Solving Ill-Conditioned Problems (in Russian), 288 pp., Nauka, Moscow.
Vinnik, L. P. (1977), Detection of waves converted from P to SV in the mantle, Phys. Earth Planet. Inter., 15, 39.
Vinnik, L. P., and G. L. Kosarev (1981), Determination of crustal parameters from observations of teleseismic body waves, Dokl. Akad. Nauk SSSR (in Russian), 261(5), 1091.
Vinnik, L. P., S. Roecker, G. L. Kosarev, S. I. Oreshin, and I. Yu. Kulakov (2002), Crustal structure and dynamics of the Tien Shan, Geophys. Res. Lett., 29(22), 94.
Vinnik, L. P., Ch. Reigber, I. M. Aleshin, G. L. Kosarev, M. K. Kaban, S. I. Oreshin, and S. W. Roecker (2004), Receiver function tomography of the central Tien Shan, Earth Phys. Sci. Lett., 225, 131.
I. M. Aleshin, G. L. Kosarev, O. Yu. Riznichenko, Institute of Physics of the Earth, Russian Academy of Sciences, 10 Bol’shaya Gruzinskaya ul., Moscow, 123995 Russia (ima@ifz.ru)
I. A. Sanina, Institute for Dynamics of Geospheres, Russian Academy of Sciences, 38/6 Leninskii prosp., Moscow, 117334 Russia
