# Bulletin of the Seismological Society of America

- Copyright © 2006 Bulletin of the Seismological Society of America

## Abstract

We have analyzed the aftershocks (*M*_{L} <4.5) following the 1999 Izmit earthquake (*M*_{w} 7.4) to infer the frequency-dependent attenuation characteristics of both *P* and *S* waves, in the frequency range from 1 to 10 Hz and in the distance range from 10 to 140 km. A linear-predictive model is assumed to describe the spectral amplitudes in terms of attenuation and source contributions. The results show that both *P* and *S* waves undergo a strong attenuation along ray paths shorter than 40 km, while the secondary arrivals significantly contribute to the spectral amplitudes over the distance range from 40 to 60 km, as also confirmed by the computation of synthetic seismograms. For longer ray paths, the decrease in attenuation suggests an increase in the propagation efficiency with depth. Finally, the spectral attenuation curves are flattened, or sloped upward at low frequencies in the range from 100 to 140 km, due to the contemporary arrivals of direct waves and postcritical reflections from the Moho. In terms of geometrical spreading and anelastic attenuation, the attenuation in the range from 10 to 40 km is well described by a spreading coefficient *n* = 1 for both *P* and *S* waves, and the quality factors can be approximated by *Q _{S}*(

*f*) = 17

*f*

^{0.80}for 1 ≤

*f*≤ 10 Hz and

*Q*(

_{P}*f*) = 56

*f*

^{0.25}for 2.5 ≤

*f*≤ 10 Hz. For ray paths in the range from 60 to 80 km, the attenuation weakens but the interaction between seismic waves and propagation medium is more complex. The multilapse time window analysis (mltwa) is applied to quantify the amount of scattering loss and intrinsic absorption for

*S*waves. The seismic albedo

*B*

_{0}decreases from 0.5 at 1 Hz to 0.3 at 10 Hz, while the total quality factor

*Q*increases from about 56 to 408. The multiple lapse time-window analysis (mltwa) results provide only an average estimate of the attenuation properties in the range from 10 to 80 km. In fact, by neglecting the variation of attenuation with depth, the mltwa results underestimate attenuation for distances less than 40 km, and do not capture the significant features caused by the integrated energy of the secondary arrivals observed in the range from 40 to 60 km.

_{T}## Introduction

This article presents the results of a crustal-attenuation study performed in northwestern Turkey (Fig. 1). The analyzed region was struck by the Izmit earthquake (*M*_{w} 7.4) on 17 August 1999 and by the Düzce earthquake (*M*_{w} 7.1) on 12 November 1999. We analyzed the aftershocks (*M*_{L} <4.5) of the Izmit earthquake to infer the frequency-dependent attenuation characteristics of both *P* and *S* waves in the frequency range from 1 to 10 Hz. Previous studies were devoted to characterizing the seismic attenuation in Turkey. Evidence of a high-attenuation crust and upper mantle in the region embracing the North Anatolian fault (naf) was presented by Kadinsky-Cade *et al.* (1981). They showed a pattern of inefficient *S _{n}* propagation in northern Turkey and suggested the presence of asthenospheric-type material beneath the crust in this region exhibiting a high attenuation in the uppermost mantle. The development of seismological observatories in northwestern Turkey—for example, the Kandilli Observatory and Earthquake Research Institute (koeri), Sakarya-Bolu (sabo), and German Task Force (gtf) networks—have made also possible the study of attenuation properties at a local scale. Examples (Fig. 1) are the investigations made in the Erzincan region, struck by as

*M*

_{s}6.8 earthquake in March 1992 (Akinci and Eyidogan, 1996; Grosser

*et al.*, 1998; Akinci and Eyidogan, 2000), in the Marmara area (Gündüz

*et al.*, 1998; Kaslilar-Özcan

*et al.*, 2002; Horasan and Boztepe-Güney, 2004), in the Bursa region (Akyol

*et al.*, 2002), and in the western portion of the naf (Akinci, Del Pazzo, and Ibanez, 1995; Akinci, Ibanez,

*et al.*, 1995; Akinci

*et al.*, 2004). In general,

*Q*(

*f*) obtained by considering either the

*P*- or

*S*-wave window or coda confirms the presence of a low

*Q*crust in the region. Considering a power function of the type

*Q*(

*f*) =

*Q*

_{0}

*f*

*, most studies found values of*

^{α}*Q*

_{0}less than 100, and values of

*α*close to or even greater than 1. Examples are

*Q*= 35

_{S}*f*

^{0.83}and

*Q*= 29

_{c}*f*

^{1.03}(for a lapse time of 20 sec) for the Erzincan region (Akinci and Eyidogan, 1996) over a distance range from 10 to 40 km and in the frequency range from 1.5 to 24 Hz, and

*Q*= 41

_{c}*f*

^{1.08}and

*Q*= 50

_{S}*f*

^{1.09}for the Marmara region (Gündüz

*et al.*, 1998) for the range from 20 to 110 km and the frequency range from 1.5 to 24 Hz. Studies on

*Q*(

_{c}*f*) that implemented the single-scattering model generally found an increase in

*Q*(

_{c}*f*) with lapse time, suggesting an increase in

*Q*with depth (Akinci and Eyidogan, 1996; Gündüz

*et al.*, 1998). Recently, a comparison between the attenuation characteristics of the western and eastern parts of the naf was performed (Akinci

*et al.*, 2004) using attenuative dispersion measurements of

*P*waves. A slightly lower

*Q*in the western naf with respect to eastern naf was observed.

We first present results from the application of the generalized inversion technique (Castro *et al.*, 1990; Anderson, 1991) to estimate the rate of decay with distance of the spectral amplitudes in northwestern Turkey for both *P* and *S* waves. Second, the attenuation-distance curves are parameterized in terms of geometrical spreading and anelastic attenuation, and the results for *Q _{P}*(

*f*) and

*Q*(

_{S}*f*) are compared over the frequency range 1–10 Hz. Finally, we apply the multilapse time window analysis, mltwa (Hoshiba, 1991; Fehler

*et al.*, 1992) to evaluate the relative contribution of intrinsic absorption and scattering to total

*S*-wave attenuation. We follow the numerical approach introduced by Hoshiba (1991), based on the Monte Carlo technique, to synthesize the theoretical integrals of the energy density for the isotropic scattering and spatial uniformity of attenuation. Moreover, using the method developed by Hoshiba (1997), synthetic seismograms and integrated seismic-wave energy are also computed for a simple conceptual model consisting of two layers over a half-space, with the aim of discussing the discrepancy between the observed and thoretically predicted integrals of seismic energy against distance.

## Data

We analyze 245 aftershocks of the 1999 Izmit earthquake recorded by the sabo and gtf networks (Grosser *et al.*, 1998; Baumbach *et al.*, 2003). The equipment used in this study consists of 1-Hz geophones (Mark L4-3D), a 24- bit digitizer with a sampling rate of 100 samples per second, and Global Positioning System (gps) timing. The analyzed earthquakes have magnitudes between 1 and 4.5, the hypocentral distances range from 10 to 142 km, and most events have a source depth between 5 and 15 km. The statistically estimated average horizontal and vertical errors are 1.8 and 2 km, respectively.

The spectra are collected from 14 stations (Fig. 2 and Table 1). We analyze the spectra of *S*- and *P*-wave windows with widths of 5 sec. If the travel-time difference between an *S* wave (*t _{S}*) and

*P*wave (

*t*) is less than 5 sec, we use the analyzed

_{P}*P*window with a width of (

*t*−

_{S}*t*) sec. Each window is cosine tapered (5%) and Fourier transformed. Instrumental corrections are applied and the spectral amplitudes are smoothed using the Konno–Ohmachi window,

_{P}*b*= 20 (Konno and Ohmachi, 1998). The signal-to-noise ratio (s/n) is greater than about 3 over the frequency range 1 to 10 Hz for

*S*waves, and from 2 to 10 Hz for

*P*waves. The two horizontal (east–west and north–south) and the vertical (

*Z*) spectra are then vectorially summed.

## Method

The attenuation-with-distance curves are estimated by applying the generalized spectral inversion (e.g., Castro *et al.*, 1990). The velocity spectrum *U*(*f,r _{ij}*) at frequency

*f*of event

*j*recorded at station

*i*is modeled as 1 2 where

*r*is the hypocentral distance,

_{ij}*Ŝ*(

_{j}*f*) is a scalar depending upon the spectral amplitude of source

*j*, and

*A*(

*f,r*) is a nonparametric description of attenuation. We discuss here the results for

_{ij}*A*(

*f,r*) obtained solving the first step, whereas the results of the second step, which was performed to split the residual spectra

_{ij}*R*(

_{ij}*f*) =

*U*(

*f,r*)/

_{ij}*A*(

*f,r*) into source

_{ij}*S*(

_{j}*f*) and site

*Z*(

_{i}*f*) spectra, are discussed elsewhere (Parolai

*et al.*, 2004).

The distance range from 10 to 142 km is discretized into 44 bins 3 km wide. The attenuation is constrained to assume a value equal to 1 at 10 km, irrespective of frequency, and the *A*(*f,r*) is constrained to be a smooth function of distance by requiring a small second derivative with respect to *r*. Equation (1) is numerically solved for each frequency separately using the least squares algorithm (lsqr) (Paige and Saunders, 1982). For all events and all stations, equation (1) can be expressed in a matrix form by 3 where **d** is a vector in the data space, **x** is a vector in the model space, and **[ A]** is a matrix relating

**x**to

**d**that only depends upon the source-to-station geometry. The matrix

**[**also includes rows relevant to the constraints. A description of its structure can be found in Castro

*A*]*et al.*(1990). We only recall that the solution vector

**is formed by [**

*x**a*

_{0}, …,

*a*,

_{M}*ŝ*

_{0}, …,

*ŝ*

_{Nev−1}], where

*a*(

_{i}*i*= 0,…,

*M*) are the unknown nonparametric attenuation coefficients, one for each distance knot

*i*, and

*ŝ*(

_{j}*j*= 0,…,

*N*− 1) are related to the unknown spectral amplitude of source

_{ev}*j*.

*M*and

*N*are the number of discrete distance intervals and the number of considered earthquakes, respectively.

_{ev}## Reliability of the Results

We compute both the unit covariance matrix **[C]** and the resolution matrix **[R]** (Menke, 1989) for assessing the effectiveness of the inversion method. These matrices are computed by performing the singular value decomposition, svd (Press *et al.*, 1994) of **[A]**, given by **[A]** = **[U][Λ][V] ^{T}**.

Figures 3 and 4 show the unit covariance and resolution matrices obtained by considering only singular values *λ _{i}* of the diagonal matrix

**[Λ]**greater than 0.01 times the maximum singular value. The resolution matrix is fairly identical, with diagonal values being close to 1. The first 45 columns are those relevant to the attenuation coefficients

*a*There is only a weak smearing of each

_{i}.*ŝ*along the attenuation coefficients, since the off-diagonal values are close to 10

_{j}^{−3}. Therefore, the resolution matrix indicates that each model parameter can be satisfactorily resolved, and the mutual dependence between the term including the source and the term describing the attenuation is weak. The columns of the unit covariance matrix

**[C]**in Figure 3 relevant to

*ŝ*(parameter index > 45) are close to being identical, with standard deviations (for unit data variance) between 0.3 and 1. The bottom frame of Figure 3 shows an example of columns of the unit covariance matrix relevant to

_{j}*a*The standard deviations are smaller than 0.3 (times the data uncertainty). The standard deviations are greater for the longest (>70 km) and shortest (<20 km) distances, depending on the lower sampling of the corresponding bins. In the distance range well constrained by data (20–70 km), the standard deviations of the model parameters are smaller than 0.1, and the spreading of the off-diagonal elements is narrowed.

_{i}.

## Nonparametric Attenuation *A*(*f,r*)

Figure 5 shows *P*- and *S*-wave attenuation-distance curves for frequencies ranging from 1 to 10 Hz, obtained by inverting equation (1). The amplitude is strongly attenuated for *r* less than 40 km in the whole analyzed frequency range, while both *P*- and *S*-wave attenuation-distance curves show a bump in the distance range between 40 and 60 km.

The attenuation is significant in the distance range between 60 and 80 km, but it is weaker than those in the first 40 km. Beyond 80 km, the attenuation-distance curves become flattened and then sloped upward after 100 km, in particular at low frequencies. The difference in the amount of attenuation for *P* and *S* waves decreases with frequency, and the *P*-wave amplitude is always less attenuated than the *S*- wave one. Moreover, the logarithm of attenuation, log*A*(*f,r*), decreases with frequency, and the difference between *P*- and *S*-wave attenuation increases with distance, being larger at low frequencies. In Figure 6, the estimated log(attenuation)- distance curves for *P* and *S* waves at three different frequencies (2, 5, and 10 Hz) are compared with the spectral amplitudes of events having 2.5 < *M*_{L} < 3, which is well sampled by the analyzed earthquakes. In this figure, an arbitrary offset was added to amplitude to make the comparison easier. The general features of the inversion results are consistent with the data, in particular the slope of the decay at short distances, the change in slope between 30 and 40 km, and again between 80 and 90 km. Since the attenuation curves for distances larger than 100 km are weakly constrained by data, the upward slope of the curves could mainly be controlled by the smoothing applied. In any case, the consistency of this result with earlier studies (e.g., Boztepe- Güney and Horasan, 2002; Baumbach *et al.*, 2003) suggests that this behavior of the attenuation-distance curves is a characteristic feature of crustal propagation in northwestern Turkey due to critical Moho reflections. Finally, the scatter affecting the data is mainly related to difference in magnitudes among the events, and to site amplification effects.

## Parametrization of *A*(*f,r*)

The attenuation-distance curves *A*(*f,r*) are parameterized in terms of geometrical spreading and anelastic attenuation. For each frequency *f*, the parameterization is as follows: 4 where *r* is the hypocentral distance in km, *c* is a scale factor, *n* is the geometrical spreading coefficient, and *γ _{φ}* is related to the quality factor

*Q*(

_{φ}*f*) (i.e.,

*γ*=

_{φ}*πf*/

*Q*(

_{φ}*f*)

*v*, where

_{φ}*v*is the velocity of the considered phase

_{φ}*φ*). We assume

*v*= 3.5 km/sec (Table 2) and . To reduce the number of free parameters, we performed the fit constraining

_{S}*c*= 10, irrespective of frequency. Since the bump affecting

*A*(

*f,r*) between 40 and 60 km cannot be modeled by equation (4), two regressions have been performed, one over the range 10–38 km and the other over the range 60–80 km.

### Distance 10–38 km

Figure 7 (left frames) shows *n* and *Q*(*f*) versus frequency for both *S* and *P* waves. For both phases, the average value of *n* is close to 1. *Q _{P}*(

*f*) shows a weak frequency dependency, while the frequency dependence of the

*S*-wave quality factor

*Q*(

_{S}*f*) is remarkable. If the power function

*Q*(

*f*) =

*Q*

_{0}

*f*is fitted to

^{α}*Q*(

*f*), the following relations are obtained: 5 6 The obtained

*Q*(

_{S}*f*) confirms that the

*S*waves undergo a strong attenuation over short paths. The value of

*α*being close to 1 is in agreement with that observed in seismically active regions. The value of

*Q*(

_{P}*f*) confirms that

*P*-wave propagation is more efficient than for

*S*waves. In Figure 7 (top-left panel), the estimated

*Q*(

_{S}*f*) and

*Q*(

_{P}*f*), given in equations (5) and (6), are plotted for comparison.

### Distance 60–80 km

The propagation characteristics of *P* and *S* waves are more complex over this range than those observed in the range from 10 to 38 km (Fig. 7, right frames). In particular, the geometrical spreading *n* is frequency dependent, for both *S* and *P* waves. The mean values vary between 0.95 and 1.35 for *S* waves and from 0.75 to 1.28 for *P* waves.

*Q*_{S}/*Q*_{P} versus Distance

Equation (4) can be used to compute the ratio of *Q _{S}*(

*f*) to

*Q*(

_{P}*f*) versus distance: 7 where

*G*(

*f,r*) is the geometrical spreading.

In Figure 8, the ratios of *Q _{S}*(

*f*) to

*Q*(

_{P}*f*) computed from equation (7) for the selected distances (16, 28, and 38 km in the left panel, 61, 70, and 79 km in the right panel) are shown. We used the average geometrical spreading for

*P*and

*S*waves previously obtained (Fig. 7) and assumed a constant ratio . In the same figure, the

*Q*(

_{S}*f*)/

*Q*(

_{P}*f*) ratios computed from equation (7) are compared to the ratio (thick line) between the same quantities but computed considering, for each distance range, the average quality factors that we obtained in the previous section and are shown in Figure 7. For distances between 10 and 38 km,

*Q*(

_{S}*f*)/

*Q*(

_{P}*f*) values computed from equation (7) increase with frequency, reaching values close to one at 10 Hz. Moreover, it varies with distance. In the range 60–80 km, the ratio oscillates from 0.5 to 1.5, and it shows a negligible dependency on distance. In both cases, the ratio between the average

*Q*and

_{S}*Q*is consistent with the distance-dependent ratios computed by equation (7).

_{P}

In Figure 8, the dotted horizontal lines correspond to , and the theoretical anelastic *Q _{P}* versus

*Q*relation given by Anderson and Hart (1978): 8 where .

_{S}For zero bulk loss (i.e., ), we obtain 9 where we used . It is worth noting that the ratio *v _{P}*/

*v*, assumed to be constant, acts as a scaling factor. The

_{S}*Q*/

_{S}*Q*ratio for values of

_{P}*v*/

_{P}*v*being different from can be easily computed by scaling the

_{S}*Q*/

_{S}*Q*ratios shown in Figure 8.

_{P}## mltwa Analysis for *S* Waves

The mltwa (Hoshiba, 1991; Fehler *et al.*, 1992; Sato and Fehler, 1998) is used to estimate the contribution of scattering loss () and intrinsic absorption () to total attenuation (). This method compares the seismic energy integrated in three consecutive time windows and displayed against hypocentral distance, with the integrals predicted by a theoretical model suitable for the multiple-scattering problem. In particular, under the hypothesis of isotropic scattering and uniform distribution of scatterers, Hoshiba (1991) numerically synthesized the integrals of energy for different attenuation characteristics of the medium via a Monte Carlo approach, and his results agree with the analytical solution of Zeng *et al.* (1991). On the contrary, no analytical solutions have been obtained yet for an inhomogeneous medium, such as a vertical layered medium. In this case, only a numerical approach can be used to compute the integrals of energy, as shown by Hoshiba (1997). The attenuation due to intrinsic absorption and scattering is usually characterized by specifying the values of the extinction length *L _{e}* (i.e., the distance over which the primary

*S*-wave energy is decreased by

*e*

^{−1}), and the seismic albedo

*B*

_{0}(i.e., the ratio between the scattering loss and the total energy loss) (Wu, 1985): 10 11 where

*g*=

*ω*(

*Q*)

_{sc}v^{−1}is the scattering coefficient,

*h*=

*ω*(

*Q*)

_{l}v^{−1}is the absorption coefficient,

*ω*= 2

*πf*is the angular frequency, and

*v*is the

*S*-wave velocity. is defined from equations (10) and (11).

The seismic energy is integrated in three consecutive time windows, each 15 sec wide, starting 1 sec before the *S* arrival *t _{S}*. For each frequency, the integrals are evaluated from the square of the amplitude spectrum at the selected frequency. The amplitude spectrum for each time window is smoothed by applying the Konno–Ohmachi window,

*b*= 60 (Konno and Ohmachi, 1998). The correction for the geometrical spreading is applied by multiplying the integrals for 4

*πr*

^{2}, while the effects of a nonspherical radiation pattern are not accounted for in the present study. The energy in each window is then divided by the coda energy integrated in a fixed reference time window (from 45 to 50 sec) in order to normalize all the different records to a common source and site (Aki, 1980). In the following, we refer to

*E*

_{1}(

*r, f*),

*E*

_{2}(

*r, f*), and

*E*

_{3}(

*r, f*) to indicate the results obtained from integrating the energy in the time windows from

*t*to

_{s}*t*+ 15, from

_{s}*t*+ 15 to

_{s}*t*+ 30, and from

_{s}*t*+ 30 to

_{s}*t*+ 45, respectively. Figure 9 shows

_{s}*E*

_{1},

*E*

_{2}, and

*E*

_{3}for the frequencies 1.6, 3.2, and 8 Hz. We computed the average of

*E*

_{1}(

*r*,

*f*),

*E*

_{2}(

*r*,

*f*), and

*E*

_{3}(

*r*,

*f*) over a window 5 km wide (black lines), −1 standard deviation (gray area). The strategy of the mltwa method is to simultaneously minimize the discrepancy between the computed integrals and those predicted by theory. Assuming a spatially uniform medium, we computed the synthetic integrals for a set of

*B*

_{0}and pairs using the Monte Carlo method. A grid search algorithm is then applied to find the best fit to the data, that is, the solution that minimizes the residual in a least-squares sense. The total residual is computed by summing the relative residual over each time window. Figure 10 shows the total residuals normalized to their minimum value, for two frequencies. In the bottom panels of Figure 10, the pairs (

*B*

_{0}, ) corresponding to a residual within 10% from the minimum are also displayed for estimating the uncertainties.

In Figure 9, the synthetic curves corresponding to the minimum residual are compared to the integrated energy against distance. In the bottom-right panel, the , , and corresponding to the solutions of minimum residual are plotted against frequency. At 1 Hz, the albedo *B*_{0} is 0.5 and the scattering loss and the intrinsic absorption equally contribute to total attenuation. The value of *Q _{T}* at 1 Hz is 56. The albedo decreases with frequency, and it becomes 0.3 at 10 Hz. It follows that the influence of scattering loss on total attenuation decreases with frequency with respect to the intrinsic absorption. At 10 Hz,

*Q*becomes 408. At 1 Hz,

_{T}*Q*is intermediate between the values of

_{T}*Q*obtained in the previous sections for distances in the range 10–38 km and 60–80 km. The comparison at 1.6 Hz between synthetic and actual integrals for the first time window

_{S}*E*shows that the theoretical curve only approximates the trend of the observed data over the distance range considered. In particular, at distances less than 40 km,

_{1}*E*computed from observed data decreases with distance more rapidly than

_{1}*E*for the synthetics. Since in the mltwa method the decay of

_{1}*E*with distance constrains the total attenuation, the obtained value of underestimates the actual attenuation at short distances. We attribute this underestimation to an inadequacy in the assumed uniformity of the propagation medium, that is, the decreases in the rate of decay for distances over 40 km may reflect increases of propagation efficiency with depth. Moreover, for distances between 40 and 60 km, the integrals over the first window show a bump, similar to the behavior observed for the spectral amplitude versus distance in Figure 5. The bump could be a consequence of reflected arrivals whose energy adds to the energy of the direct arrival. The presence of secondary arrivals in the range 40–60 km that cannot be accounted for using a uniform model introduces a further underestimation of the total attenuation at distances less than 40 km.

_{1}Figure 11 shows the synthetic *SH* velocity seismograms computed by applying the method of Wang (1999) to the 1D model shown in Table 2 and derived from Karahan *et al.* (2001). The source depth is 10 km, typical for the considered earthquakes. For distances larger than 30 km, the results for this simplified model show that secondary arrivals can contribute significantly to the total energy in the first time window. For example, the maximum amplitude of synthetic seismograms in Figure 11 is almost constant for distances between 60 and 80 km, and the peak value corresponds to secondary arrivals for distances larger than 70 km. This result suggests that vertical stratification is necessary to obtain an accurate description of attenuation in northwestern Turkey. The Monte Carlo method can also be used to synthesize the energy integrals for a stratified medium (Hoshiba, 1997). However, the computational burden strongly increases with increasing number of layers and limits the grid search procedure to only a coarse grid of *L _{e}* and

*B*

_{0}values for each layer. In Figure 9 (top-left panel, dotted lines) we also show the best-fit synthetic curves obtained by applying the Hoshiba (1997) method to the model given in Table 2. The synthetic integrals are computed only for 1.6 Hz. The values of

*L*and

_{e}*B*

_{0}are fixed for the half-space and the second layer, and the grid search is applied only to the first layer. In particular, we set

*g*= 0.001 and

*h*= 0.003 for the half-space, and

*g*= 0.002 and

*h*= 0.006 for the second layer. The attenuation coefficients of the uppermost layer for the best-fit solution are

*g*= 0.02 and

*h*= 0.10 (i.e.,

*B*

_{0}= 0.15 and . The value of corresponds to

*Q*= 45. The results for the layered model do not increase the quality of the fit with respect to the results for the uniform model. This is not surprising since a better description of the data would require a less simplified vertical model, and to set up a multigrid analysis for investigating the attenuation and velocity properties of all the layers is beyond the aims of the present work. However, the results for the stratified model can be used to stress the importance of considering vertical layering. The theoretical

_{T}*E*

_{1}integral decreases till 50 km, from which it flattens and then slopes upward after 60 km. Even if the upward trend is delayed toward larger distances, this resembles the behavior of the observed

*E*

_{1}. Moreover, the slope of

*E*

_{1}computed for the layered model at short distances better describes the observed energy decay, even if the values of the integrals are slightly overestimated. The overestimation, considered together with the observed underestimation of

*E*

_{2}, suggests that the seismic albedo has been underestimated.

## Discussions and Conclusions

The results of the previous sections depict a complex interaction between the crustal properties of northwestern Turkey and the seismic-wave propagation. We identified four characteristic distance intervals.

For 10 ≤ *r* ≤ 38 km, the waves propagate through a low *Q*-volume. We obtained fairly constant values for *n* being close to 1, suggesting that the selected windows are mainly dominated by direct *S* or *P* waves. For *f* < 10 Hz, *Q _{S}*(

*f*) = 17

*f*

^{0.80}has a behavior similar to that observed in previous studies along the naf, but with a smaller

*Q*-value at 1 Hz. Gündüz

*et al.*(1998) found

*Q*(

_{S}*f*) = 50

*f*

^{1.09}in the Marmara region, Akyol

*et al.*(2002) found

*Q*(

_{S}*f*) = 46.59

*f*

^{0.67}in the Bursa region, and Akinci and Eyidogan (1996) found

*Q*(

_{S}*f*) = 35

*f*

^{0.83}in the Erzincan region. Recently, Horasan and Boztepe-Güney (2004) estimated the

*Q*values ranging from 13

_{S}*f*

^{1.22}to 94

*f*

^{0.83}in the Sea of Marmara using the earthquakes over the distance range 15–70 km and in the frequency range 1.5–12 Hz. Our

*Q*(

_{P}*f*) shows a weak frequency dependence in the range 2–10 Hz, and is in good agreement with the mean value of 40 found by Akinci

*et al.*(2004) at a distance of 30 km. The

*Q*(

_{S}*f*)/

*Q*(

_{P}*f*) ranges from 0.5 to 1.3 (Fig. 8). Assuming a frequency-independent

*Q*, Liu

*et al.*(1991) suggested that the

*Q*/

_{S}*Q*in the New Madrid seismic zone lies between 0.53 and 1.79, in agreement with Chen

_{P}*et al.*(1994). In the Central Rio Grande rift, Carpenter and Sanford (1985) observed a ratio of between 0.72 and 2.94, and assumed it to be frequency independent, whereas Clouser and Langstone (1991) found a ratio between 0.38 and 1.87 for a basin 1.4 km thick in the Gazli region (Uzbekistan). Exploiting the results achieved in laboratory experiments and extrapolating them to

*in situ*analysis, the

*Q*/

_{S}*Q*in sedimentary rocks is known to be affected by several factors, including pore fluids, pore-fluid saturation levels, confining pressure, crack distribution, and rock type. In particular, many laboratory studies (Spencer, 1979; Toksöz

_{P}*et al.*, 1979; Johnston and Toksöz, 1980; Winkler and Nur, 1982; etc.) suggest that a

*Q*greater than

_{P}*Q*is expected for fully saturated or dry rocks, and less than

_{S}*Q*for partially saturated rocks. Whatever the causes controlling the

_{S}*Q*/

_{S}*Q*ratio are, a distance dependency of this ratio is clear in the range from 10 to 38 km. This confirms that there is a strong heterogeneous structure of upper crust in northwestern Turkey.

_{P}For 38 < *r* ≤ 60 km, two factors can be considered to describe the behavior of spectral attenuation with distance: (a) the dependence of *Q* with depth, which may cause a change in the slope of the attenuation curve, and (b) the arrival of reflected waves. Several studies show that the boundary between the fall-off of the direct waves and the emergence of lower crustal or Moho reflections can lead to fairly constant amplitudes over distance ranges that depend on several factors, such focal depth, crustal thickness, crustal-velocity gradient, and so on. For example, Burger *et al.* (1987) found high amplitudes in the range between 60 and 150 km using a model for the central United States. Somerville and Yoshimura (1990) observed high amplitudes in the range 40–100 km using the Loma Prieta earthquake recorded in the San Francisco and Oakland areas, while Hartzell (1992) found high amplitudes in the range 30–40 km while analyzing the 1989 Loma Prieta earthquake to estimate site response along the San Francisco Peninsula. Within the framework of magnitude calibration studies, Bakun and Joyner (1984) found large amplitudes at distance ranges of 75–125 km for central California, while Hutton and Boore (1987) found similar results at around 60 km for southern California. Mori and Helmberger (1996) noted that these enhanced amplitudes might be attributed to *SmS.* Atkinson and Mereu (1992) showed that the contemporary arrivals of direct waves and postcritical reflections from the Moho produce spectral amplitudes that are approximately constant between 70 and 130 km in southeastern Canada. In Turkey, Boztepe-Güney and Horasan (2002) observed large- amplitude *SmS* phases around 100 km for the data recorded toward the west in the Sea of Marmara.

For 60 < *r* ≤ 80 km, the waves propagate through a volume having an average higher *Q* than in the first 38 km, indicating an increasing propagation efficiency with depth. Studies that applied the coda-wave method to estimate the coda-quality factor *Q _{c}* along the naf found an increase of

*Q*with lapse time. This evidence has been interpreted as the presence of a more attenuating shallow crust with respect to the deeper one. The ratio

_{c}*Q*/

_{S}*Q*does not depend on distance, suggesting a higher homogeneity of the properties that control this ratio. The geometrical spreading for

_{P}*S*waves is between 0.9 and 1.2 for frequencies lower than 4 Hz. For higher frequencies, it ranges between 1.2 and 1.4. The geometrical spreading for

*P*waves is about 0.9 for frequencies less than 7 Hz, and increases up to about 1.2 at 10 Hz. The values of geometrical spreading greater than one were found in vertically layered medium by Frankel

*et al.*(1990). For

*r*> 80 km, first the slope of attenuation decay decreases, then fairly constant values are obtained, until finally, the curves slope upward for distances larger than 100 km.

The time integrals of energy against distance also confirm the depth dependence of attenuation properties in the analyzed area and the presence of secondary arrivals that significantly contribute to spectral amplitude in the range 40–60 km. Since the application of the mltwa to a medium that is not spatially uniform requires both a large computational effort and the availability of a realistic model for the area, we interpreted the observations by generating synthetic integrals for a homogeneous medium, being aware that the results can provide only an average description of the attenuation properties in the range from 10 to 80 km. The results of mltwa analysis indicate that the albedo for the area is fairly low, assuming a value of 0.5 at 1 Hz and decreasing to 0.3 at 10 Hz. Then, the intrinsic absorption contributes equally to scattering loss in determining the total attenuation at low frequencies, but it becomes the dominant mechanism at high frequencies, suggesting that the coda could be mainly generated by a back-scattering mechanism (Menke and Chen, 1984). The decrease in albedo with frequency has been previously observed in other seismogenic regions, such as Japan (Hoshiba, 1993), northern Chile (Hoshiba *et al.*, 2001), southern central Alaska (Dutta *et al.*, 2004), and southern Spain (Akinci, Del Pezzo, and Ibanez 1995). The total attenuation over the range 10–80 km varies from about 0.018 at 1 Hz to 0.0024 at 10 Hz. For frequencies less than 2 Hz, is intermediate between the values of obtained for the ranges 10–40 km and 60–80 km, while for higher frequencies, is close to obtained for the range 60–80 km. The result of this comparison confirms that the total attenuation obtained via the mltwa analysis applied to a homogeneous medium provides an average estimate underestimating the attenuation at shallow depths. A first attempt to compare the observations with the results obtained by applying the mltwa method to a layered model have also been performed. Even if the results suffer from the limits imposed by the computational effort required, they confirm that some peculiarities are present in the observations that can be ascribed to the variation of attenuation and velocity with depth. Since several variables play a role in determining the results of the mltwa analysis, such as the number of layers and their thickness, the velocity and attenuation of each layers, and so on, producing a detailed regional model for the investigated area is a matter that deserves future work, and the features appearing in the attenuation-with-distance curves could be helpful in constraining such a model for the area. In conclusion, the obtained spectral model for attenuation can be exploited to calibrate synthetic attenuation relationships regressing ground-motion parameters predicted by applying, for example, stochastic techniques.

## Acknowledgments

The authors thank R. Wang and M. Hoshiba for making available their computer codes for the calculation of the synthetic seismograms and the Monte Carlo integrals of energy for the multiple scattering problem, respectively. The authors express their gratitude to the Hannover Rückversicherung AG for their significant financial support of the field mission. K. Fleming kindly improved our English. Comments from two anonymous reviewers also improved the manusript. The figures have been drawn using the Generic Mapping Tools software (Wessel and Smith, 2000). Part of this work was conducted during the visits of D. Bindi at the GeoForschungsZentrum Potsdam (gfz-Potsdam) that were partially funded by the gfz-Potsdam.

- Manuscript received 4 March 2005.