Click Here

JOURNAL OF GEOPHYSICAL RESEARCH, VOL. 113, B09403, doi:10.1029/2007JB005317, 2008

for

Full Article

Uplift of the 2004 Sumatra-Andaman earthquake measured from differential hyperspectral imagery of coastal waters Se´bastien Smet,1,2 Re´mi Michel,1 and Laurent Bollinger1 Received 8 August 2007; revised 1 March 2008; accepted 21 May 2008; published 9 September 2008.

[1] We describe a procedure to measure coseismic change of shallow coastal bathymetry

using multispectral imagery. This technique is applied to HYPERION hyperspectral images acquired along the shallow coast of the North Andaman Islands to estimate the uplift induced by the Mw 9.15, 26 December 2004, earthquake. Attenuation coefficient of the particularly clear coastal waters is estimated from two preevent images with a 22 cm tide level difference. Various sources of noise on the estimate of the uplift resulting from atmospheric correction, data registration, sensor noise, and fundamental assumption of stationary optical properties of the scene with time are studied. Average uplift over the shallow bathymetry covered by the imagery is 0.85 ± 0.10 m, increasing from south to north from 0.56 ± 0.10 to 1.12 ± 0.10 m. The uplift amplitude is consistent with local field measurements. These data place constraints on the width of the megathrust rupture in the Andaman area, estimated to about 160 km, and on the amount of coseismic slip there estimated to about 10 m. Citation: Smet, S., R. Michel, and L. Bollinger (2008), Uplift of the 2004 Sumatra-Andaman earthquake measured from differential hyperspectral imagery of coastal waters, J. Geophys. Res., 113, B09403, doi:10.1029/2007JB005317.

1. Introduction [2] The location of the pivot line (or neutral axis) dividing uplifted from subsiding regions as well as some crude estimates of the uplift gradient are key constraints to evaluate the rupture parameters of large subduction earthquakes (Figure 1). Some of these parameters such as rupture width, slip amplitude and azimuth as well as dipping of the subduction interface have been estimated along the 1600 km ruptured by the Sumatra-Andaman mega-earthquake on the basis of seismological and far field geodetical data set [e.g., Banerjee et al., 2005; Vigny et al., 2005; Subarya et al., 2006; Chlieh et al., 2007]. Some of these studies also benefit from local estimates of ground displacement deduced from near-field geodetical data [Gahalaut et al., 2006] and/or coastal subsidence and uplift observations [Chia et al., 2005; Subarya et al., 2006; Meltzner et al., 2006]. In spite of these efforts, the details of the rupture extent and coseismic slip pattern remain poorly constrained. In this paper we show that uplift or subsidence of the seafloor can be estimated from remote sensing, and that this technique could be an efficient way to complement more traditional ways of measuring coseismic ground deformation. [3] Remote sensing imagery is commonly used to measure coseismic deformation from various techniques including SAR interferometry and correlation of optical images 1 Laboratoire de De´tection et de Geophysique, CEA, Bruye`res-le-Chaˆtel, France. 2 Now at Actimar, Brest, France.

Copyright 2008 by the American Geophysical Union. 0148-0227/08/2007JB005317$09.00

acquired before and after the event [e.g., Massonnet et al., 1993; Michel et al., 1999; Van Puymbroeck et al., 2000; Wright et al., 2004; Zebker et al., 1994]. These techniques are rarely used in coastal environments where subduction earthquakes can produce significant vertical deformation due to difficult applications, induced by poor SAR images correlations or very low frequency deformation difficult to constrain while correlating optical images. However, passive optical spectral imagery of coastal waters yields information about the bathymetry, radiative properties of seawater and of seafloor. Some techniques have been developed to recover these parameters from hyperspectral images [e.g., Pozdnyakov and Grassl, 2003; Mobley, 1994; Adler-Golden et al., 2005; Mobley and Sundman, 2000]. Absolute estimate of the water depth by optical spectral imagery is challenging however because of the strong temporal and geographic variability of the optical properties of water and because of trade-off between parameters [Pozdnyakov and Grassl, 2003; Mobley, 1994; Adler-Golden et al., 2005]. In this study we exploit three hyperspectral images of the northwest Andaman Islands coastal area acquired by the HYPERION satellite before and after the Mw 9.15, December 2004, Sumatra-Andaman earthquake to estimate the tectonic uplift (Figure 2 and Table 1). The three images studied have been chosen because of their similar radiance. [4] In the present study we first provide a tectonic setting of the study area. We then recall the basics of radiative transfer useful for that study, quantify the optical properties of the seawater and measure the bathymetric changes from comparing the images and estimating tide levels at the time of the images acquisition. Error analysis and evidence of temporal stationary of the averaged scenes are then de-

B09403

1 of 10

B09403

SMET ET AL.: DIFFERENTIAL HYPERSPECTRAL BATHYMETRY

B09403

2. Background

recurrent major earthquakes, such as the 1941 Ms 7.7 Andaman or the 1881 Mw 7.8 Carnicobar earthquakes [Ortiz and Bilham, 2003]. Although these historical events have recently allowed the unlocking of significant segments along the trench, their rupture traces have been reruptured on 26 December 2004, participating to the giant Mw 9.15 Sumatra earthquake. The large coseismic moment release [Ammon et al., 2005; Stein and Okal, 2005] is associated with a seismic source corresponding to a 1300 km rupture length inferred from the spatial extent and location of the aftershocks, of the T wave radiations along the trench [Guilbert et al., 2005] as well as characterized from tsunami induced altimetric signal monitored by satellite imagery and tide gauges [Titov et al., 2005]. The far field GPS displacements measured further constrain the amount of slip accommodated along the ruptured area [Vigny et al., 2005]. These far field GPS data poorly resolve the slip amplitudes given to range from 0 to 10 m under the North Andaman. Complementary near-field GPS solutions [Subarya et al., 2006; Gahalaut et al., 2006] improved resolution of the rupture characteristics, as well as local measurements of subsidence (or uplift) estimated from coral heads mean distance below (above) their highest level of survival [Briggs et al., 2006; Meltzner et al., 2006] or minimum coastal uplift derived from optical imagery [Meltzner et al., 2006]. All these data sets have been inverted in detailed rupture models [Subarya et al., 2006; Chlieh et al., 2007]. [6] However, it appears that several local estimates of the rupture width, azimuth and amplitude along strike are still debated since the spatial coverage of the observations is inhomogeneous, a pattern largely induced by the poor spatial sampling possibilities available since most of the rupture area is under water. The N012E orientation of the Hyperion hyperspectral images at 10– 15N latitudes, similar to the Andaman coastline orientation (Figure 2), near the northern edge of the rupture, as well as the presence of a suspected uplift gradient along the coastline due to the presence of a pivot line oblique to the island [Meltzner et al., 2006] offer an opportunity to evaluate the potentiality of hyperspectral differential bathymetry as well as further quantify the rupture parameters in that area.

2.1. Tectonics Context [5] The convergence velocity of about 5 cm a
1 between the Indo-Australian plate and the Sunda block is absorbed by underthrusting of the Indo-Australian plate along the Sunda subduction zone [e.g., Bock et al., 2003]. Thrusting is predominantly accommodated by creep on the subduction interface at depth greater than 50 km [e.g., Simoes et al., 2004]. The shallower part of this interface is locked in the time period separating major earthquakes (the ‘‘interseismic period’’). The hangingwall of the thrust system behaves elastically and deforms generating localized uplift along the coast [Simoes et al., 2004]. The energy stored during the interseismic phase at midcrustal depths is released during large devastating earthquakes. These events are well documented along Sumatra coastal area by the dating of corals microatoll growth, monitoring continuously uplift and subsidence of the coastal region [e.g., Sieh et al., 1999; Natawidjaja et al., 2004]. Further north along the Andaman segment of Sunda trench, the convergence rates is lower although not very well constrained (probably around 12 mm a
1 [Paul et al., 2001]). These islands are also prone to

2.2. Hyperspectral Imagery [7] HYPERION provides nadir-viewing radiance images in the range (400 – 2500 nm) with a 10 nm spectral resolution, a pixel size of about 30  30 m in the ground and a noise level of about 0.09 mW cm
2 sr
1 nm
1 per pixel (EO-1 Validation Report, http://eo1.gsfc.nasa.gov/new/validationReport/) (see Figure 3). Useful spectral bandwidth for water analysis is reduced to (400– 800 nm) because of near total absorption of light for wavelengths greater than about 800 nm [Pozdnyakov and Grassl, 2003]. We further reduced it to the range (570– 690 nm) in order to minimize major and poorly constrained contribution of aerosol at short wavelengths and residual miscompensated contribution of atmospheric water vapor [Zhao and Nakajima, 1997]. [8] Atmospheric corrections of sensor radiances are first performed using ad hoc methods in a multiple scattering theory of atmosphere radiative transfer on the complete Hyperion spectral range (400 –2500 nm). A standard atmospheric correction using Fast Line-of-Sight Atmospheric Analysis of Spectral Hypercubes (FLAASH), a MODTRAN

Figure 1. Inferred surface deformation induced by a rupture at a subduction trench: a pivot line on the vertical component divides uplifted from subsiding areas. Locating the pivot line yields information about the width and dip angle of the ruptured zone. Locating the pivot line, often far out at sea, is not straightforward and may be inferred from uplift gradients estimates. In this study, uplift and its lateral gradient is estimated from hyperspectral imagery along the coastline. scribed. We finally use the new data set obtained to constrain the rupture azimuth, width, dip and amplitude of the 26 December 2004 rupture along the northern section of the North Andaman and further discuss implications of the technique for remote sensing of coastal regions.

2 of 10

B09403

SMET ET AL.: DIFFERENTIAL HYPERSPECTRAL BATHYMETRY

Figure 2

3 of 10

B09403

B09403

SMET ET AL.: DIFFERENTIAL HYPERSPECTRAL BATHYMETRY

B09403

Table 1. Hyperion Hyperspectral Images and Sea Level Used in the Present Studya 1 2 3

Image Reference

Date

Site Longitude

Site Latitude

Sea Level (m)

EO11340512004049110PY EO11340512004065110PY EO11340512005042110KZ

18 February 2004 5 March 2004 11 February 2005

92.8842 92.8842 92.8000

13.3422 13.3422 12.5000

1.42 1.64 1.90

a

Images acquired to minimize effects of seasonal variations on radiative transfer. Theoretical sea level from tide gauge at Port Blair, South Andaman Islands (http://www.shom.fr/).

based model [Kneizys et al., 1996] including the ability to compensate for atmospheric adjacency effect, is performed [Adler-Golden et al., 1998]. A standard maritime aerosol is selected with an aerosol-scale height set to 2.0 km, while the CO2 mixing ratio is set to 390 ppm. Second, we further tested the transmittance obtained with that first approach and optimized the atmospheric correction by applying a gas retrieval methodology using the joint reflectance and gas estimator (JRGE) methodology [Marion et al., 2004]. This technique estimates the variations in gas concentrations relatively to the standard atmospheric model used in the first step (applied here on land and vegetation areas, zone C in Figure 3). It allows us to evaluate the biases on the estimates of the surface reflectance R generated by the first approach, evaluating p.e. a pixel-by-pixel water vapor content with a precision better than 10%.

3. Estimating Bathymetric Changes From Hyperspectral Images 3.1. Estimating Bathymetric Changes From Radiance Images [9] Surface reflectance, R, which can be estimated from hyperspectral images, can provide information on water depth provided that the water is clear and depth shallow enough for the reflectance from the sea bottom to be non negligible. [10] Single scattering theory of seawater column transfer function yield, R ¼ Rg þ Rw þ ðRb
Rw Þe
2Kz ;

ð1Þ

where Rg is the sea surface contribution combining Fresnel reflection and foam spectral component, Rw is the deep water subsurface scattering contribution of photons that did not reach the sea bottom, Rb is the reflectance of the seafloor, K is an effective attenuation coefficient that includes contribution of upwelling and downwelling attenuation, z is the depth of the water column. For image number i (see Table 1 and Figure 4) the depth zi is the sum of the bathymetry b and of the tide ti. In addition, the bathymetry of the postearthquake image is reduced by the tectonic uplift U.

[11] Let us now see how tectonic uplift might be retrieved by comparing images. As a preliminary, we make the assumption that the radiative transfer properties of the shallow water does not vary significantly between images. This assumption seems reasonable for the three images considered here (Figure 3). Moreover, differences in water depth between the studied images are typically metric or submetric [Meltzner et al., 2006; Chia et al., 2005] and the considered waters are exceptionally clear (Figure 3b and tides in Table 1). The difference between the single and the multiple scattering approximation increases with the number of scatters. Hence, the radiative contribution of the difference in water depth between the considered images can be modeled by the single scattering approximation of (1). 3.2. Estimate of the Uplift U From Radiance Images [12] After pixel-to-pixel atmospheric correction, in order to reduce the effect of noise, we first work on pixels averaged over three areas describing the deep waters (area A: 8.4 km2), the shallow waters (area B: 8.7 km2), and the inland (area C: 5.7 km2); areas A, B, and C, respectively, defined in Figure 3a. [13] (Rg + Rw)i is estimated for each image i as the average value of the reflectance within the deep waters area A. This is because of negligible contribution of the sea bottom in this region. Rg + Rw can be considered equal for the deep and shallow waters areas A and B because of the lack of large sediment deposit from rivers and because of the kilometric size of the averaged pixels. In the following R
(Rg + Rw) is referred to as the equivalent reflectance. [14] K is then estimated from (1) and from the equivalent reflectances, R1 and R2 of preearthquake images 1 and 2, respectively, as "
 # R2
Rg þ Rw 2 1
 K¼ log 2ðz1
z2 Þ R1
Rg þ Rw 1 "
 # R2
Rg þ Rw 2 1
 : ¼ log 2ðt1
t2 Þ R1
Rg þ Rw 1

ð2Þ

Figure 2. (a) Tectonic setting of the North Andaman Island. Also shown are minimum and maximum estimates of uplift and subsidence derived from satellite imagery from Meltzner et al. [2006] as well as GPS displacement field from Gahalaut et al. [2006]. Envelope of the pivot line (in gray) is derived from uplift estimates of these studies. Location of the pivot line in the northernmost region of the island is far out at sea and therefore poorly constrained (gray area). (b) Location of the hyperspectral Hyperion scenes considered (red box). Black squares Z1 to Z5 identify the five subzones considered in Figure 3. GPS displacements at Aerial Bay (AB) and Diglipur (DGLP) are from Gahalaut et al. [2006] and Jade et al. [2005], respectively. Green outlined white filled triangles are uplift estimates from Kayanne et al. [2007]. Minimum vertical displacement (blue circles) from Meltzner et al. [2006]. The numbers in red (followed by d for days) correspond to the number of postseismic days integrated in the measure. Location of cross section CC0 (Figure 7). 4 of 10

B09403

SMET ET AL.: DIFFERENTIAL HYPERSPECTRAL BATHYMETRY

B09403

[15] The tectonic uplift U is then estimated from images 1 and 3 acquired, respectively, before and after the earthquake as

Figure 4. Water column depth zi at date i depends on bathymetry b, tides ti, and tectonic uplift U (Table 1).

"
 # R3
Rg þ Rw 3 1
 : U ¼ t3
t1 þ log 2K R1
Rg þ Rw 1

ð3Þ

U is estimated for each wavelength and then averaged over the spectral range (570 –690 nm). 3.3. Errors Analysis [16] The procedure proposed here suffers uncertainties due to various sources of noise (Table 2): sensor noise, registration of images, atmospheric correction, tidal value and the assumption of stationary optical properties of the scene (including the water) with time. The uncertainties of the tidal value taken into account, is given by the Service Hydrographique et Oce´anographique de la Marine (SHOM, http:// www.shom.fr/) to be less than 5% of the tide differential (assuming no atmospheric perturbation), corresponding here to a 0.01 m uncertainty on uplift. [17] The HYPERION data noise level is about 0.09 mW cm
2 sr
1 nm
1 per pixel (EO-1 Validation Report, http:// eo1.gsfc.nasa.gov/new/validationReport/). It corresponds to a signal-to-noise ratio (SNR) equal to about 40 db per pixel of the raw image for shallow water radiance of area B for wavelengths in the range (570– 690 nm). Neglecting systematic errors, the SNR is enhanced following a square root n procedure by averaging the n pixels of shallow water area B. The uncertainty Dz on z is estimated from differentiation of (1) as Dz ¼

DR 1 1 1 pffiffiffi ¼ pffiffiffi : 2KR n 2KSNR n

ð4Þ

[18] Overestimate of Dz is computed from the low attenuation coefficient K of pure water [Pozdnyakov and Grassl, 2003], and we found an error due to sensor radiometry equal to 0.04 m for each pixel of 30  30 m

Figure 3. Quicklook of Hyperion images 1, 2, 3 (lR = 570 nm, lG = 620 nm, and lB = 680 nm; see Table 1 for further description). (a) Deep bathymetry area used to estimate average oceanic scattering and sea surface radiative contribution (area A). Averaged uplift is estimated from area B, a region of shallow bathymetry. Atmospheric correction is checked on inland area C. (b) Zoom along the coast on a clear water – shallow bathymetry region. 5 of 10

SMET ET AL.: DIFFERENTIAL HYPERSPECTRAL BATHYMETRY

B09403

B09403

Table 2. Sources of Errors and Space Averaged Amplitudes of the Uncertaintiesa Source of Error

Comments

Uncertainty on Uplift (m)

Sensor noise Atmospheric miscorrection Data misregistration Ocean tide

SNR about 40 maximum error on water vapor: 10% amplitude about 0.2 pixel less than 5% of the tide differential

0.04 per pixel 0.01 0.03 per pixel 0.01

a Variations in optical properties of the sea yield about 0.08 m of uncertainty on estimated uplift. The uncertainty on the ocean tide determination is evaluated assuming no atmospheric effects on the tide level. Global uncertainty on uplift is ±0.10 m. SNR is signal-to-noise ratio.

of the raw images yielding to negligible errors for the 9636 averaged pixels of shallow water area B. [19] Misregistration in the spatial and spectral wavelength domains yields uncertainty on Dz from (2) and (3). We estimated the misregistration to be equal to about 0.2 pixels (rms) representative from the subpixel correlation procedure used [Van Puymbroeck et al., 2000]. In order to estimate the associated uncertainty Dz on z we shifted one image by 0.2 pixel and we computed z between the raw and shifted images. We found Dz equal to 0.03 m per pixel and thus negligible when averaged over area B. [20] Had hoc procedure estimates the atmospheric contribution from the radiance images and yield atmospheric correction with typical precision better than 10% [Marion et al., 2004]. In order to reduce the effect of residual uncertainty on the atmospheric contribution, we limit our bandwidth to (570 – 690 nm). We estimate uncertainty Dz

on z from two equivalent reflectances of image 1 computed for atmospheric water content varying randomly by an amount of 10%. We computed an estimate of Dz equal to 0.01 m. [21] In order to estimate the effect of changes in optical properties of the water column and sea bottom we simulate the reflectance R30 averaged over area B from (1), image 1, tide parameters t1 and t3 and estimate of the uplift U as

 
 R03 ¼ R1
Rg þ Rw 1 e
2K ðt3
t1
U Þ þ Rg þ Rw 3 :

ð5Þ

  DR R3
R03  Dz ¼ ; ¼ 2KR3 2KR

ð6Þ

and

Figure 5. Averaged reflectances. (a) Averaged inland reflectances (region C) vary, between images 1 and 2, and images 1 and 3, by 6% and 13%, respectively; 13% include both atmospheric miscompensation and vegetation changes, compatible with reported 10% uncertainty on atmospheric correction (transmittance dotted line) (see section 3.3). (b) Averaged deep water reflectance over area A (Figure 1) used to compensate shallow reflectances estimated over area B (c) yields equivalent reflectance (d). Equivalent reflectances show slopes resulting from water attenuation and differences in water depth (tide and uplift). Absorption coefficient derived from images 1 and 2; uplift estimated from images 1 and 3. Modeled postearthquake equivalent reflectance (dashed line) fits data (gray line) within 5%. 6 of 10

B09403

SMET ET AL.: DIFFERENTIAL HYPERSPECTRAL BATHYMETRY

Figure 6. Uplift and temporal water variation. (a) Averaged uplift over wavelengths (570– 690 nm) is 0.85 m, and 0.08 m (rms) variations mainly result from temporal changes of water attenuation (b). Attenuation coefficient K (K12) derived from images 1 and 2; K13 derived from images 1 and 3 and estimated uplift. [22] R30 and R3 differ by an averaged 5% (see section 4) over the considered bandwidth (Figure 5d). Associated Dz is estimated to be equal to 0.08 m for area B. As a conclusion, we can estimate the resulting uncertainty Dz on z as the geometrical mean of reported errors. We found Dz equal to 0.10 m for pixels averaged.

4. Results [23] Figure 5a shows the averaged reflectance computed for inland area C from the three images used in this study. Reflectance varies by an average 6% between images 1 and 2 (preearthquake) and by 13% between images 1 and 3 (acquired, respectively, before and after the earthquake). Those differences may originate partly from changes in optical properties of the vegetation. Even if we consider that incorrect compensation of the atmospheric transfer function is responsible for all the observed difference, the result is still compatible with the 10% of errors described in section 3. [24] Figures 5b and 5c show averaged reflectances for areas A (deep water) and B (shallow water), respectively. The reflectances over deep water only result from the terms Rg + Rw defined in equation (1). Those reflectances vary less between image 1 and 2 than between images 1 and 3. The variations do not depend significantly on wavelength and may originate from changes of scattering properties of the oceanic surface and shallow water. Reflectances in shallow water show slopes characteristic of the nearexponential attenuation of the incident photon with water depth (1). Differences between the plots in Figure 5c

B09403

depend on both the variations of the averaged deep water reflectances of Figure 5b and on the differences in the depth of shallow water column. [25] Figure 5d shows the equivalent reflectance in shallow area B. The equivalent reflectance of image 2 is below that of image 1 because the depth of water column is larger for image 2 because of a greater tide (Table 1). The attenuation coefficient K is derived from those plots from (2). The equivalent reflectance of image 3 is above that of image 1 indicating that the water column for image 3 is smaller than for image 1 because of a combination of tide and uplift effects (Figure 4). [26] Figure 6a shows the uplift U estimated for the wavelength in the range (570– 690 nm). The average is 0.85 m and the standard deviation is 0.08 m (consistent with assumption error). Variation of U with wavelength are partially correlated with the water attenuation suggesting that the optical properties of the water changed between images, contributing to the uncertainty on the estimate of U. Figure 6b illustrates the differences in attenuation coefficient K, respectively estimated from images 1 and 2 (K12) using (2) and estimated from images 1 and 3 (K13) from (3), once U has been estimated and averaged over wavelength (Figure 6a). K12 and K13 differ by a maximum of 20% which is typical of low variations of optical properties of seawater [Pozdnyakov and Grassl, 2003]. Those uncertainties on U and K yield errors on models of equivalent reflectances. Modeled reflectance R30 of image 3 from image 1, K and U differs from R3 by 5% (Figure 5d). Thus a posteriori estimates of variation of optical parameters with time are consistent with the assumptions that the radiative transfer properties of the considered shallow waters did not vary significantly between the three studied images.

5. Comparison With Existing Uplift Data Set [27] We describe in a first subsection local estimates of cumulated coseismic and postseismic shoreline uplift and GPS displacement field available in the first months after the earthquake. We then assume that insignificant afterslip signal variations bias the data set, a necessary assumption to compare altogether the sparse uplift measures available. Our uplift measures are finally used to evaluate a model of integrated and homogeneous slip under our region of interest. We describe in a second subsection the evidences for a significant afterslip, moderating the results interpretation. 5.1. Integrated Coseismic and Postseismic Uplift Estimations [28] The 0.85 ± 0.10 m vertical ground displacement determined in this study is consistent with the first coarse estimates at 1 – 2 ± 1 m along the northern Andaman west coast [Bilham, 2005] (M. Searle, personal communication, 2006) as well as extensive estimates of minimum uplift along the coast covered by our HYPERION scenes at 0.36 ± 0.14 m [Meltzner et al., 2006] (Table 3). Furthermore, GPS measurements of the coseismic displacement lead to estimate nearby uplift at Aerial Bay (AB) as well as further south in Diglipur (DGLP), on the east coast. Uplifts at these sites are 0.49 ± 0.05 m [Gahalaut et al., 2006] and 0.59 ± 0.01 m [Jade et al., 2005] (Figures 2b and 7). More

7 of 10

SMET ET AL.: DIFFERENTIAL HYPERSPECTRAL BATHYMETRY

B09403

B09403

Table 3. Uplift Estimates for the Coastal Region Covered by the Images (Area B) as Well as for Five Subzonesa Zone

z1

Longitude Latitude Minimum uplift (m) from Meltzner et al. [2006] Uplift (m) from Chia et al. [2005] Uplift (m) derived from HYPERION Uplift (m) derived from Elastic Dislocation Model

92.826° 13.23° 1.12 ± 0.10 1.01

z2

z3

z4

z5

92.818° 92.808° 92.808° 92.795° 13.20° 13.15° 13.12° 13.07° 0.36 ± 0.14 0.36 ± 0.14 0.36 ± 0.14 0.37 ± 0.14 about 1 m in average in the range [0.5, 1.5 m] 0.89 ± 0.10 0.80 ± 0.10 0.77 ± 0.10 0.56 ± 0.10 0.98 0.91 0.85 0.78

Average 0.36 ± 0.14 0.85 ± 0.10 0.90

a

Subzones; i.e., zones 1 to 5. Minimum uplift estimates from Meltzner et al. [2006]; modeled uplift from this study (Figure 7).

recently, local estimates of the uplift have been obtained in the area from studies involving measurements of coral microatolls and mussel and oyster bed elevations [Kayanne et al., 2007; Rajendran et al., 2007] (Figure 2b and 7). These authors estimated uplifts around 0.6 ± 0.1 m about 10– 25 km east of the image, in a region of presumed positive uplift gradient toward the west (see Figures 1 and 2). An uplift amplitude above the range (0.49 to 0.6 ± 0.1 m) was therefore expected on the west coast. Similarly, an uplift amplitude below the range (1.1 to 1.3 ± 0.1 m) was expected from surveyed coral microatolls on North Reef island, about 5 – 10 km west of our region of interest. [29] Meltzner et al. [2006] as well as Kayanne et al. [2007], interpolating their uplift observations, propose that the neutral axis should be oriented 35°E and around 15 to 30°E in our region of interest. Assuming a constant slip alongside, these results suggest that the coast oriented at 05 to 10°E may have recorded a northward positive uplift gradient, the northernmost sites (zone z1) being significantly farther from the neutral axis in a region of positive east-west uplift gradient.

[30] The division of our area B into five parts (Figure 3), allows us to evaluate uplift variations along the coast (Table 3). The uplift decreases from north to south, by a significant 0.56 ± 0.2 m (Table 3) that is beyond noise estimates of U, as suggested by the geometry. Assuming a constant geometry along strike (i.e., constant dip and azimuth of the ruptured zone) this gradient suggests that the slip varies along strike (in amplitude or azimuth) or that the azimuth of the downdip end of the ruptured zone is slightly different from the azimuth of the coastline. However, the azimuth of the GPS derived coseismic displacement field appears to be constant along strike in the Andaman region [Gahalaut et al., 2006]. For this reason, we do not suspect much change in the azimuth of the slip at depth. Futhermore, assuming that the increase in uplift from south to north is due to a northward increase in slip on the subduction interface is not likely since North Andaman Island is suspected to be near the tip of the ruptured area. Therefore, the uplift gradient might have been generated by the existence of a clockwise angle in between the azimuth of the downdip end of the ruptured zone and the coastline (Figure 2b).

Figure 7. Comparison between the displacements estimated from this study (squares), measured on the field (triangles), deduced from other remote sensing imagery techniques (circles) and a simplistic best fit elastic dislocation model (lines) projected along CC0 (Figure 2b). GPS velocities are from Gahalaut et al. [2006] and Jade et al. [2005], minimum uplift magnitude is from Meltzner et al. [2006], coral head and oyster bed uplifts are from Kayanne et al. [2007], and differential bathymetry estimates are from this study. Numbers (followed by d for days) correspond to the number of postseismic days integrated in the measure. Thick, dashed, and dotted lines correspond to vertical, E-W, and S-N displacements, respectively, predicted by a best fit (see section 5.1) constant slip and constant geometry elastic dislocation model of the surface displacement (35°N dislocation in an isotropic elastic half-space with a 9 m slip, rake 48°, on a 160 km large 21°E dipping plane). Elastic parameters in the model have been a Poisson coefficient n = 0.5 and a crustal ratio of P tied to classical value, l = m = 0.33  1011 pffiffiyielding ffi over S wave seismic velocities Vp/Vs = 3. 8 of 10

B09403

SMET ET AL.: DIFFERENTIAL HYPERSPECTRAL BATHYMETRY

[31] Assuming the pivot line in our region of interest presents an azimuth around N035 and a distance of about 145 km from the subduction trench as given by [Meltzner et al., 2006], we varied these parameters around these values (±20° and ±20 km every 5° and 5 km). We further determine a range of dip angle, rake (i.e., the angle between the fault azimuth and the slip on the fault plane) as well as slip amplitude that fit both the estimated uplift as well as the vertical and horizontal component of the displacement measured at GPS monuments. The best fitting model is determined from minimizing a reduced c2 criterion, measuring the discrepancy between modeled and observed velocities. The best fitting model (c2 = 2.43, with residuals typically less than about 15 cm, average less than 2 cm, showing no systematic pattern in terms of their geographic distribution) leads to an estimate of a 9 m rupture dislocation at N035, dipping at 21°E over a 160 km width encompassing a rake at 48° (Figure 7; see caption for details). This result should be interpreted with care: this purely elastic dislocation model is based on strong assumptions including a constant geometry of the dislocation, only justified at small scale, as well as a constant integrated coseismic and postseismic slip magnitude among others. Although locally satisfactory, its displacement predictions cannot be extrapolated at larger scales suggesting that some initial assumptions are too simplistic for larger-scale models, arguing for lateral and/or temporal variations of the coseismic and/or postseismic slip. 5.2. Afterslip Biases [32] Our measures of the uplift incorporate 47 days of postseismic slip, given the late acquisition of the 3rd Hyperion image. Unfortunately, no near-field continuous GPS station network was continuously monitoring the postseismic deformation, very few days after the event. First-order estimates of the surface displacement field occurring over the first 50 days after the earthquake has been monitored by far field continuous GPS stations, mainly in Thailand, Malaysia and Sumatra [Vigny et al., 2005]. Their study put in evidence large variations of the afterslip magnitude over that period, locally up to 1.25 times the initial coseismic displacement as well as a more important relative after- to coseismic slip in the northern, including Nicobar-Andaman Islands, than in the southern part of the rupture. Chlieh et al. [2007], on the basis of the inversion of far field GPS solutions, suggest that more than 35% of the coseismic displacement is accommodated by postseismic slip under the Andaman Islands 30 days after the event. Although difficult to ascertain, given the sparse far field data used to constrain it, this postseismic slip seems mainly restricted to the downdip end of the ruptured zone [Chlieh et al., 2007]. This result seems consistent with a possible postseismic uplift measured on a mid-January 2005 12 day period in Port Blair (at PB GPS station), coincident to a 1.2 m a
1 N259° velocity of the station [Gahalaut et al., 2006]. However, near-field subsidence observations in the northern part of North Andaman Island [Kayanne et al., 2007; Rajendran et al., 2007] suggest that significant slip is also accommodated at the updip end of the ruptured zone, probably due to the upward propagation of the slip at the trench [Kayanne et al., 2007]. Eyewitness accounts indicate a 90 day postseismic subsidence of the shoreline as large as

B09403

40% of the coseismic uplift [Kayanne et al., 2007] at the forest Office Camp in Interview Island as well as in Mayabunder. Interpretations of second-order differential uplift estimates integrating various time period within 3 months from the earthquake is therefore hazardous. Furthermore, coseismic slip, 5 day cumulated slip, and 30 day cumulated slip do not appear to be proportional [Vigny et al., 2005; Chlieh et al., 2007]. Any decent correction of the uplift observations data set, integrating 0 to 90 days of afterslip, needed to further constrain biases or uncertainties in our uplift determinations, appears therefore illusive.

6. Discussion and Conclusion [33] This study attempting to retrieve uplift from hyperspectral reflectance data using ocean color techniques illustrates the potential of differential multispectral imagery for the monitoring of earthquake deformations within coastal areas, mainly of interest in subduction context. [34] The availability of Hyperion data in the rupture area of the 26/12/2004 Great Sumatra earthquake offered an opportunity to test such a technique. By averaging about 10,000 pixels (in area B corresponding to shallow bathymetry) in 12 bands we estimate the tectonics uplift and its variations along the coast with a precision of ±0.10 m. The mean uplift amplitude determined in this study, around 0.85 ± 0.10 m and its spatial gradient, are consistent with the sparse uplift data set available in the area. This technique allows to measure pluri-decimetrical uplift and subsidence over large areas, with potentially the ability to revisit regularly the surveyed scenes. [35] Unfortunately, the method described here appears to be applicable in a few opportunistic cases only since its application requires clear waters to limit the influence of signal attenuation, large areas of shallow (typically <30 m) bathymetry depending on the sensors’ resolution, very few spatial and temporal variations of coastal turbidity often driven by wind, as well as high Sun angle, and homogeneous tidal response, among other parameters. It further makes necessary the constitution of large imagery databases with similar acquisition parameters to properly evaluate the accuracy of the uplift measurement, depending on local environmental parameters, often seasonal and sometimes transient, influencing the ocean color. Variations in water quality and bottom albedo at most places [Pozdnyakov and Grassl, 2003] are generally much more important than in the present Andaman context so that enhanced techniques including in situ measurement of optical properties of water, sea bottom albedo as well as the join use of bathymetric Lidar with the hyperspectral data should be instrumented [Guenther et al., 2000]. Finally, variations of coefficient K with time could also be constrained from a regularization procedure that would ensure that U does not depend on the spectral wavelength (Figure 6). [36] It should be noted that hyperspectral imagery is not mandatory for such a multispectral differential bathymetry approach. We primarily benefited from an opportunistically rich Hyperion hyperspectral database, with clear water and cloud free preearthquake and postearthquake scenes. Hyperion imagery gives us here the opportunity to discuss atmospheric transfer function biases and make atmospheric corrections nearby gas absorption bands, including water

9 of 10

B09403

SMET ET AL.: DIFFERENTIAL HYPERSPECTRAL BATHYMETRY

vapor, using the potential use of as much as 240 spectral bands. However, only 12 bands, in the 570 – 690 nm range, weakly influenced by water vapor content of the atmosphere, allows us to work on the shallow bottom reflectance variations given the total absorption of light for wavelengths greater than 800 nm. A range further reduced to minimize major unconstrained aerosol and molecular scattering biases. Most multispectral sensors, although possessing lower spectral resolution, allow working at higher spatial resolution with higher signal-to-noise ratio and may have the capacity to properly measure the bathymetry differential. This suggests that both multi and hyperspectral data offer new potentials to monitor tectonic deformation in coastal areas. [37] Acknowledgments. We are most grateful to two anonymous reviewers and Pedro Elosegui, the Associate Editor, for their comments and suggestions, which greatly helped improve this paper.

References Adler-Golden, S. M., A. Berk, L. S. Bernstein, S. C. Richtsmeier, P. K. Acharya, M. W. Matthew, G. P. Anderson, C. L. Allred, L. S. Jeong, and J. H. Chetwynd (1998), FLAASH, a MODTRAN4 atmospheric correction package for hyperspectral data retrievals and simulations, in 1998. Proceedings of the Seventh JPL Airborne Earth Science Workshop, JPL Publ. 98-21, pp. 9 – 14, Jet Propul. Lab., Pasadena, Calif. Adler-Golden, S. M., P. K. Acharya, A. Berk, M. W. Matthew, and D. Gorodetzky (2005), Remote bathymetry of the littoral zone from AVIRIS, LASH, and Quickbird imagery, IEEE Trans. Geosci. Remote Sens., 43(2), 337 – 347, doi:10.1109/TGRS.2004.841246. Ammon, C. J., et al. (2005), Rupture process of the 2004 SumatraAndaman earthquake, Science, 308, 1133 – 1139, doi:10.1126/ science.1112260. Banerjee, P., F. F. Pollitz, and R. Burgmann (2005), The size and duration of the Sumatra-Andaman earthquake from far-field static offsets, Science, 308, 1769 – 1772, doi:10.1126/science.1113746. Bilham, R. (2005), A flying start, then a slow slip, Science, 308, 1126 – 1127, doi:10.1126/science.1113363. Bock, Y., L. Prawirodirdjo, J. F. Genrich, C. W. Stevens, R. McCaffrey, C. Subarya, S. S. O. Puntodewo, and E. Calais (2003), Crustal motion in Indonesia from Global Positioning System measurements, J. Geophys. Res., 108(B8), 2367, doi:10.1029/2001JB000324. Briggs, R., et al. (2006), Deformation and slip along the Sunda megathrust in the great 200 5 Nias-Simeulue earthquake, Science, 311(5769), 1897 – 1901, doi:10.1126/science.1122602. Chia, Aik Song, Soo Chin Liew, A. W. C. Heng, and Leong Keong Kwoh (2005), Satellite observations of coastline changes in the Andaman islands after the 2004 Sumatra earthquake, in Proceedings of the IEEE International Geoscience and Remote Sensing Symposium 2005. IGARSS ’05. Proceedings. 2005 IEEE International, vol. 3, pp. 1838 – 1840, IEEE, Washington, D. C. Chlieh, M., et al. (2007), Coseismic slip and afterslip due to Mw 9.15, 2004 great Sumatra-Andaman earthquake, Bull. Seismol. Soc. Am., 97, doi:10.1785/0120050631. Gahalaut, V. K., B. Nagarajan, J. K. Catherine, and S. Kumar (2006), Constraints on 2004 Sumatra-Andaman earthquake rupture from GPS measurements in Andaman-Nicobar Islands, Earth Planet. Sci. Lett., 242, 365 – 374, doi:10.1016/j.epsl.2005.11.051. Guenther, G. C., M. W. Brooks, and P. E. LaRocque (2000), New capabilities of ‘‘SHOALS’’ airborne Lidar bathymeter, Remote Sens. Environ., 73, 247 – 255, doi:10.1016/S0034-4257(00)00099-7. Guilbert, J., J. Vergoz, E. Schissele, A. Roueff, and Y. Cansi (2005), Use of hydroacoustic and seismic arrays to observe rupture propagation and source extent of the Mw = 9.0 Sumatra earthquake, Geophys. Res. Lett., 32, L15310, doi:10.1029/2005GL022966. Jade, S., M. B. Ananda, P. D. Kumar, and S. Banerjee (2005), Co-seismic and post-seismic displacements in Andaman and Nicobar Islands from GPS measurements, Curr. Sci., 88(12), 1980 – 1984. Kayanne, H., Y. Ikeda, T. Echigo, M. Shishikura, T. Kamataki, K. Satake, J. N. Malik, S. R. Basir, G. K. Chakrabortty, and A. K. G. Roy (2007), Coseismic and postseismic creep in the Andaman Islands associated with the 2004 Sumatra-Andaman earthquake, Geophys. Res. Lett., 34, L01310, doi:10.1029/2006GL028200.

B09403

Kneizys, F. X., et al. (1996), The MODTRAN 2/3 Report and LOWTRAN 7 MODEL, edited by L. W. Abreu and G. P. Anderson, Rep. NrF19628-91C-0132, 261 pp., Ontar Corp. North Andover, Mass. Marion, R., R. Michel, and C. Faye (2004), Measuring trace gases in plumes from hyperspectral remotely sensed data, IEEE Trans. Geosci. Remote Sens., 42, 854 – 864, doi:10.1109/TGRS.2003.820604. Massonnet, D., M. Rossi, C. Carmona, F. Adragna, G. Peltzer, K. Feigl, and T. Rabaute (1993), The displacement field of the Landers earthquake mapped by radar interferometry, Nature, 364, 138 – 142, doi:10.1038/ 364138a0. Meltzner, A. J., K. Sieh, M. Abrams, D. C. Agnew, K. W. Hudnut, J. P. Avouac, and D. H. Natawidjaja (2006), Uplift and subsidence associated with the great Aceh-Andaman earthquake of 2004, J. Geophys. Res., 111, B02407, doi:10.1029/2005JB003891. Michel, R., J. P. Avouac, and J. Taboury (1999), Measuring near field ground displacements from SAR images: Application to the Landers earthquake, Geophys. Res. Lett., 26, 3017 – 3020, doi:10.1029/ 1999GL900524. Mobley, C. D. (1994), Light and Water. Radiative Transfer in Natural Waters, 592 pp., Academic Press, San Diego, Calif. Mobley, C. D., and L. K. Sundman (2000), Hydrolight 4.1 Users’ Guide, 85 pp., Sequoia Sci., Redmond, Wash. Natawidjaja, D., K. Sieh, S. Ward, H. Cheng, R. L. Edwards, J. Galetzka, and B. Suwargadi (2004), Paleogeodetic records of seismic and aseismic subduction from central Sumatran microatolls, Indonesia, J. Geophys. Res., 109, B04306, doi:10.1029/2003JB002398. Ortiz, M., and R. Bilham (2003), Source area and rupture parameters of the 31 December 188 1 Mw = 7.9 Car Nicobar earthquake estimated from tsunamis recorded in the Bay of Bengal, J. Geophys. Res., 108(B4), 2215, doi:10.1029/2002JB001941. Paul, J., et al. (2001), The motion and active deformation of India, Geophys. Res. Lett., 28, 647 – 650, doi:10.1029/2000GL011832. Pozdnyakov, D., and H. Grassl (2003), Colour of Inland and Coastal Waters: A Methodology for Its Interpretation, 170 pp., Praxis, Chichester, UK. Rajendran, C. P., K. Rajendran, and R. Anu (2007), Crustal deformation and seismic history associated with the 2004 Indian Ocean earthquake: A perspective from the Andaman-Nicobar Islands, Bull. Seismol. Soc. Am., 97(1A), S174 – S191, doi:10.1785/0120050630. Sieh, K., S. N. Ward, D. Natawidjaja, and B. W. Suwargadi (1999), Crustal deformation at the Sumatran subduction zone revealed by coral rings, Geophys. Res. Lett., 26(20), 3141 – 3144, doi:10.1029/1999GL005409. Simoes, M., J. P. Avouac, R. Cattin, and P. Henry (2004), The Sumatra subduction zone: A case for a locked fault zone extending into the mantle, J. Geophys. Res., 109, B10402, doi:10.1029/2003JB002958. Stein, S., and E. A. Okal (2005), Speed and size of the Sumatra earthquake, Nature, 434, 581 – 582, doi:10.1038/434581a. Subarya, C., M. Chlieh, L. Prawirodirdjo, J. P. Avouac, Y. Bock, K. Sieh, A. Meltzner, D. Natawidjaja, and R. McCaffrey (2006), Plate-boundary deformation associated with the great Sumatra-Andaman earthquake, Nature, 440, 46 – 51, doi:10.1038/nature04522. Titov, V., A. B. Rabinovich, H. O. Mofjeld, R. E. Thomson, and F. I. Gonzalez (2005), The global reach of the 26 December 2004 Sumatra tsunami, Science, 309, 2045 – 2048, doi:10.1126/science.1114576. Van Puymbroeck, N., R. Michel, R. Binet, J. P. Avouac, and J. Taboury (2000), Measuring earthquakes from optical satellite images, Appl. Opt., 39(20), 3486 – 3494, doi:10.1364/AO.39.003486. Vigny, C., et al. (2005), Insight into the 2004 Sumatra-Andaman earthquake from GPS measurements in southeast Asia, Nature, 436, 201 – 206, doi:10.1038/nature03937. Wright, T. J., B. E. Parsons, and Z. Lu (2004), Toward mapping surface deformation in three dimensions using InSAR, Geophys. Res. Lett., 31, L01607, doi:10.1029/2003GL018827. Zebker, H. A., P. Rosen, R. M. Goldstein, A. Gabriel, and C. L. Werner (1994), On the derivation of coseismic displacement fields using differential radar interferometry: The Landers earthquake, J. Geophys. Res., 99, 19,617 – 19,634, doi:10.1029/94JB01179. Zhao, F., and T. Nakajima (1997), Simultaneous determination of waterleaving reflectance and aerosol optical thickness from Coastal Zone Color Scanner measurements, Appl. Opt., 36(27), 6949 – 6956, doi:10.1364/ AO.36.006949.























L. Bollinger and R. Michel, CEA-DAM-DASE-SLDG, BP 12, F-91680 Bruye`res-le-Chaˆtel CEDEX, France. ([email protected]; remi.michel@ cea.fr) S. Smet, Actimar, 36 quai de la douane, F-29200 Brest, France. ([email protected])

10 of 10