GEOPHYSICAL RESEARCH LETTERS, VOL. 32, L02307, doi:10.1029/2004GL021897, 2005

Horizontal coseismic deformation of the 2003 Bam (Iran) earthquake measured from SPOT-5 THR satellite imagery R. Binet and L. Bollinger De´partment Analyse et Surveillance de l’Environnement, Commissariat a` l’Energie Atomique, Bruye`res-le-Chaˆtel, France Received 2 November 2004; revised 16 December 2004; accepted 5 January 2005; published 28 January 2005.

[1] On 26 December 2003, a destructive Mw 6.5 earthquake occurred in Kerman province, Iran, killing more than 30,000 people around the ancient city of Bam. We map the fault zone and determine the slip distribution using sub-pixel correlation of images taken by the SPOT-5 THR sensor before and after the main shock. Our results show that the surface displacement occurred on a blind fault oriented N174 and located South of Bam. A 15 km along-strike profile depicts a distribution with a maximum displacement of 1.2 ± 0.15 m near the centre of the fault. The mean slip that we infer (
0.77 ± 0.05 m) is much higher than the ground displacement measured in the field (
0.2 m). We suggest that field-measured displacements are underestimated due to slip accommodation within an 500 m wide shear zone. Furthermore, we demonstrate for the first time the relevancy of 2.5 m SPOT-5 THR images for sub-pixel ground displacement measurements. Citation: Binet, R., and L. Bollinger (2005), Horizontal coseismic deformation of the 2003 Bam (Iran) earthquake measured from SPOT-5 THR satellite imagery, Geophys. Res. Lett., 32, L02307, doi:10.1029/ 2004GL021897.

1. Introduction [2] Bam, an ancient city in Kerman province, Iran, lies at the southeastern end of a large active right-lateral strike-slip fault system, the Gowk fault (Figure 1). Active faulting along this fault system is related to the differential shortening accommodating collision of the Arabian and Eurasian plates. About 15 mm/yr are expected to be accommodated on such right lateral strike slip structures between central Iran and the Kerman plateau, Afghanistan [Vernant et al., 2004]. Prior to the 2003 Bam earthquake, three large earthquakes have recently struck the region, namely Golbaf 1981-Mw 6.6, Sirch 1981-Mw 7.1 and Fandoqa 1998-Mw 6.6. Seismic waveform inversions, surface rupture measurements, determination of coseismic ground displacement by radar interferometry [Berberian et al., 2001] as well as evaluation of aseismic slip on the eastward Shahdad foldand-thrust belt basal detachment [Fielding et al., 2004] lead to a very detailed seismotectonic setting along the fault (Figure 1). However, the lack of instrumental as well as historical seismic event along the southern edge of the Gowk fault suggests the presence of a long standing seismic gap [e.g., Ambraseys and Melville, 1982]. [3] On 26 December 2003, a destructive Mw 6.5 earthquake occurred near Bam killing more than 30,000 people in the city. Macroseismic evidence and re-location of the Copyright 2005 by the American Geophysical Union. 0094-8276/05/2004GL021897$05.00

main shock do not indicate that the rupture occurred on the main Gowk fault segment spatially located between recent earthquakes but on a secondary fault farther East and less distant from the city. The previously mapped Bam fault [e.g., Walker and Jackson, 2002], which follows a morphological high that develops on the southeastern side of the city seemed a reasonable alternative structure whose geometry is consistent with the focal mechanism. Radar interferometric decorrelation effects as well as surface displacement mapped using Envisat radar images in Interferometric Synthetic Aperture Radar (InSAR) processing reveal the potential trace of the seismogenic fault [Talebian et al., 2004]. This identification generated further field investigations allowing the mapping of a right lateral displacement of up to 20 cm following N174 trace of a previously blind structure [Talebian et al., 2004]. This ground rupture is comprised of a series of discontinuous 50 – 100 m en echelon surface breaks. Both traces, (i.e., deduced from the InSAR and the second field mapping) correspond to the Arg-e-Bam segment illuminated by the aftershocks of the main event [Suzuki et al., 2004]. [4] Because the detailed distribution of coseismic fault slip is a key parameter for seismic hazard assessment as well as for rupture modeling, complementary constraints using independent techniques are useful for horizontal surface displacement measurement. For this purpose, we apply for the first time a sub-pixel correlation method to a pair of SPOT-5 THR images acquired before and after the main shock to measure the differential ground displacement. Such measurement is complementary to the InSAR deformation measurements and field-biased horizontal surface displacement measurements since it maps the details of the horizontal ground deformations near the fault. We finally evaluate the horizontal slip on the fault and compare it to the available field measurements as well as with the near field displacement model deduced from InSAR off-fault measurement.

2. Methodology [5] Using a pair of SPOT panchromatic images, [Van Puymbroeck et al., 2000] showed that sub-pixel correlation could provide fault slip measurements with an accuracy of 0.1 pixel (1 meter) [Michel and Avouac, 2002; Dominguez et al., 2003]. We adapted this correlator to ingest SPOT-5 THR images. The increased sampling rate between SPOT-4 (10 m) and SPOT-5 THR (2.5 m) should increase the offset accuracy by a factor of 4. Moreover, SPOT-5 imaging technology is an advancement on that of SPOT-4 due to its higher geometric fidelity. Firstly, the focal plane geometry is provided, and secondly the attitude data accuracy is

L02307

1 of 4

L02307

BINET AND BOLLINGER: HORIZONTAL COSEISMIC DEFORMATION

L02307

Table 1. SPOT-5 THR Parametersa Geometrical Optimisation Acquisition

Date

IA(°)

Roll (rad)

Preseismic 16/10/2003 5.35 5.41 104 Postseismic 28/02/2004 2.06 5.42.104

Pitch (rad)

Yaw (rad)

2.09.104 5.27.105

3.46.104 8.89.106

a

I.A. stands for incidence angle.

Figure 1. Seismotectonic setting along the Gowk fault. Shaded relief topographic map of Bam area derived from SRTM Digital Elevation Model. Focal mechanisms from Harvard CMT catalogue and [Berberian et al., 2001]. Location of the postseismic aftershock cluster on the Arg e Bam fault from [Suzuki et al., 2004] (blue box). Solid and dashed black boxes stand respectively for SPOT-5 images and Figure 2 boundaries. increased by the star tracking sensor of SPOT-5. These geometric enhancements should allow the measurement of low spatial frequency components of the ground deformation, which could give crucial information on earthquake rupture parameters. Finally, the THR specific image acquisition and processing [Latry and Rouge, 2003] provide an optimal sampling. Consequently, two major sources of correlation bias are avoided because on the one hand the images do not experience aliasing and on the other hand the images can be properly interpolated with a sinus cardinal kernel. [6] The approach of [Van Puymbroeck et al., 2000] first consist in resampling the SPOT images into a map projec-

tion so that the remaining image pixel offsets are only due to the earthquake ground deformation. In a second step, the coseismic offset map and its error estimate are computed from the phase shift of the Fourier transform of a sliding window. [7] In order to benefit of SPOT-5 geometric properties, we developed a specific camera model which takes into account the satellite attitude data and the focal plane geometry. Despite the enhancement of the satellite absolute attitude measurement, the pointing accuracy is corrupted by the unknown dilatation of the instrument [Bouillon et al., 2003]. An optimization of the satellite attitude mean angles (roll, pitch, yaw) is therefore required. For this purpose, ground control points are computed automatically by means of correlation of the image with a shaded relief digital elevation model (DEM). [8] Since the phase only filter correlator attributes the same weight to all spatial frequencies and since the high spatial frequencies of the THR oversampled images are null (Figure R1, available as auxiliary material1), the [Van Puymbroeck et al., 2000] correlator has to be adapted to the THR images. The circular top-hat function we apply to weight the spectra is empirically deducted from the circular shape of the image’s power spectrum. The correlator biases are less than 0.02 pixels and its accuracy highly depends on the temporal decorrelation of the images.

3. Disparities and Slip Distribution [9] The correlation window size used in this study is 256 pixels and the window step is 64 pixels in both directions. While no significant offset break is detectable on the EW component (Figure 2b; Table 1), the NS

Figure 2. Disparity map showing detail of the ground displacement induced by the Bam earthquake computed using images 1 and 2 (Table 1). (a) SPOT image; (b) EW component. Notice stereoscopic induced offset biases over Baravat. (c) NS component of the displacement field. 2 of 4

L02307

BINET AND BOLLINGER: HORIZONTAL COSEISMIC DEFORMATION

Figure 3. Disparities analysis. (a) Pitch difference affecting the image lines. (b) NS pixel offsets computed with a 128 pixels window size over the whole images and rotated for visual convenience. Notice the artifacts correlated with the satellite attitude drift (plane component) and the pitch errors (line structured). The high resolution image solution as well as an extract focusing on oscillation artifacts is available in auxiliary data material (auxiliary Figure R2). component reveals a 15 km long N174 right lateral strike slip south of Bam city (Figure 2c). [10] Attempts in reducing the window size lead to a noisy offset map near the fault because of the temporal decorrelation and of the low contrast of the ground. Such offsets computed over the whole image with a window size of 128 pixels are shown in Figure 3b. [11] Whatever window size is chosen, the offsets are corrupted by several biases. We observe a bilinear offset bias whose amplitude is 0.5 pixel over the whole image, which is probably due to a satellite drift bias. Line correlated biases in offset measurements (0.1 pixel) are also visible, well correlated with the differential satellite pitch (Figure 3a), which might be due to undersampled attitude measurements. A periodic pattern whose amplitude is less than 0.1 pixel is also visible in the disparity map (Figure R2, available as auxiliary material), either not correlated with the offsets, nor with the image lines. We suggest that the image processing of the SPOT-5 THR product [Latry and Rouge, 2003] may be the source of such oscillations. DEM errors also induce stereoscopic offset biases only on the EW component because the along track viewing angle of THR images is constant. [12] Despite these biases, the fault slip measurement is not corrupted thanks to the NS orientation of the fault. Cumulated offset profiles across the fault were computed by setting manually the fault geometry (FF0 on Figure 2c), dividing the fault by one kilometre long segments, and computing a weighted mean profile in each segment. These transects enable us to estimate a throw profile depicted on Figure 4b.

4. Discussion [ 13 ] The zone of differential ground displacement matches with the location of the fault segment determined by InSAR as well as with the available field data [Talebian 1 Auxiliary material is available at ftp://ftp.agu.org/apend/gl/ 2004GL021897.

L02307

et al., 2004] and the postseismic cluster location [Suzuki et al., 2004]. Slip profile AA0 which is perpendicular to the rupture trace (Figure 4a) depicts a dramatic increase in slip accommodation through a shear zone with a width of
800 m. Taking into account the convolution effect of the sliding window, the shear zone would be about 500 m (Figure 4a). [14] The 1.2 m integrated displacement measured on that profile leads to a mean displacement gradient around 2  103. We propose that the difference between the field estimations (20 cm) and our measurements (up to 1.2 m) is mainly due to a low frequency deformation accommodation. This might be possible if the displacement is taken up on anastomosing faults, mole tracks, distributed en-echelon tension cracks, as well as on diffuse damping zones induced by a medium with low cohesion. The resulting low deformation gradient might be beyond capacity for field mapping in desert areas and therefore should lead to an underestimation of the integrated ground deformation. [15] The high spatial density of the disparity map allows us to measure the along-strike throw profile. This displacement profile depicts a typical pattern of slip distribution with a maximum displacement Dmax
1.2 ± 0.15 m near the centre of the fault (Figure 4b). The displacement tapers off gradually reaching symmetrical gradients down to the presumed tips of the fault. However, the weak resolution at

Figure 4. Slip distribution measured from the NS pixel offset map. (a) Along the A-A0 fault-perpendicular profile. This profile is the sum of individual offset profiles along a 1 km fault segment weighted by the correlation coefficient of the offset measurement. Dashed line: offset profile of a simulated slip distribution of amplitude 1.2 m located only on the main fault with the same window size. The difference of the two profiles suggests the presence of a shear zone. The correlator yields independent measurements every 250 meters. (b) Along strike throw profile. We estimate a linear weighted least square fit of the offsets for the East and West compartments of the fault. The slip is given by the offset between the two linear fits at the fault location. The slip error is roughly estimated by taking the standard deviation of the misfit.

3 of 4

L02307

BINET AND BOLLINGER: HORIZONTAL COSEISMIC DEFORMATION

the tips does not enable us to decree on the absence of inflections or bell-shape of the throw profile. We finally find an average slip Davg
0.77 ± 0.05 m giving a Davg/ Dmax
0.6 ± 0.1 and a Davg/L
5.1  105 well within the single slip events characteristic range (104 – 105) [e.g., Scholz, 1982]. Assuming a seismogenic crust thickness around 15 km, consistent with the maximum depth of the aftershocks [Suzuki et al., 2004], these estimates are in good agreement with empirical relationships among magnitude, rupture areas and surface displacement [Wells and Coppersmith, 1994], making this event a typical strike-slip event. [16] It should be noted that despite their nearly similar characteristics (depth, Mw, mechanism, geological nature of the ground), adjacent earthquakes along the Gowk fault had highly disparate slip magnitudes (e.g., 1 – 3 m for Fandoqa Mw6.8 vs 0.5 m for Sirch Mw 7.1 and 0.05 m for Golbaf Mw 6.6) [Berberian et al., 2001]. As well as Bam, Sirch and Golbaf surface displacements might have been underestimated in the field due to deformation accommodation. [17] Finally, this first SPOT-5 THR disparity mapping allowed us to measure coseismic ground displacement up to 1.2 m on a 500 m wide shear zone with a 0.3 m accuracy which is
4 times better than with the SPOT-4 case [Van Puymbroeck et al., 2000]. This technique complements the InSAR and GPS measurements since it better constrains the subsurface fault geometry.

References Ambraseys, N. N., and C. P. Melville (1982), A History of Persian Earthquakes, Cambridge Univ. Press, New York. Berberian, M., J. A. Jackson, C. Baker, E. Fielding, B. E. Parsons, K. Priestley, M. Qorashi, M. Talebian, R. Walker, and T. J. Wright (2001), The 14 March 1998 Fandoqa earthquake (Mw 6.6) in Kerman province, S.E. Iran: Re-rupture of the 1981 Sirch earthquake fault, triggering of slip on adjacent thrusts, and the active tectonics of the Gowk fault zone, Geophys. J. Int., 146, 371 – 398. Bouillon, A., E. Breton, F. De Lussy, and R. Gachet (2003), SPOT-5 geometric image quality, in IGARSS 2003: Learning From Earth’s

L02307

Shapes and Sizes: 2003 IEEE International Geoscience and Remote Sensing Symposium: Proceedings, vol. 1, pp. 303 – 305, IEEE Press, Piscataway, N. J. Dominguez, S., J. P. Avouac, and R. Michel (2003), Horizontal coseismic deformation of the 1999 Chi-Chi earthquake measured from SPOT satellite images: Implications for the seismic cycle along the western foothills of central Taiwan, J. Geophys. Res., 108(B2), 2083, doi:10.1029/ 2001JB000951. Fielding, E. J., T. J. Wright, J. Muller, B. E. Parsons, and R. Walker (2004), Aseismic deformation of a fold-and-thrust belt imaged by synthetic aperture radar interferometry near Shahdad, southeast Iran, Geology, 21(7), 577 – 580, doi:10.1130/G204752.1. Latry, C., and B. Rouge (2003), Super resolution: Quincunx sampling and fusion processing, in IGARSS 2003: Learning From Earth’s Shapes and Sizes: 2003 IEEE International Geoscience and Remote Sensing Symposium: Proceedings, vol. 1, pp. 315 – 317, IEEE Press, Piscataway, N. J. Michel, R., and J. P. Avouac (2002), Deformation due to the 17 August Izmit, Turkey, earthquake measured from SPOT images, J. Geophys. Res., 107(B4), 2062, doi:10.1029/2000JB000102. Scholz, C. H. (1982), Scaling laws for large earthquakes: Consequences for physical models, Bull. Seismol. Soc. Am., 72, 1 – 14. Suzuki, S., T. Nakamura, H. Sadeghi, T. Matsushima, Y. Ito, S. Hosseini, A. Jafar Gandomi, M. Maleki, and S. Fatemi Aghda (2004), 3D geometric structure of the blind fault for the 2003 Bam earthquake (southeast Iran) inferred from the aftershock distribution: Existence of the Arg-e-Bam fault proposed, Eos Trans. AGU, 85(47), Fall Meet. Suppl., Abstract S11A-1001. Talebian, M., et al. (2004), The 2003 Bam (Iran) earthquake: Rupture of a blind strike-slip fault, Geophys. Res. Lett., 31, L11611, doi:10.1029/ 2004GL020058. Van Puymbroeck, N., R. Michel, R. Binet, J. P. Avouac, and J. Taboury (2000), Measuring earthquakes from optical satellite images, Appl. Opt., 39(20), 1 – 14. Vernant, P., et al. (2004), Present-day crustal deformation and plate kinematics in the Middle East constrained by GPS measurements in Iran and northern Oman, Geophys. J. Int., 157(1), 381 – 398. Walker, R., and J. A. Jackson (2002), Offset and evolution of the Gowk fault, S.E. Iran: A major intra-continental strike-slip system, J. Struct. Geol., 24, 1677 – 1698. Wells, D. L., and K. J. Coppersmith (1994), New empirical relationships among magnitude, rupture length, rupture width, rupture area, and surface displacement, Seismol. Soc. Am. Bull., 84, 974 – 1002.



R. Binet and L. Bollinger, Laboratoire de De´tection et de Ge´ophysique, CEA, BP12, F-91680 Bruye`res-le-Chaˆtel, France. ([email protected]; [email protected])

4 of 4