Fault Zone Resistivity Structure and Monitoring at the Taiwan Chelungpu Drilling Project ( TCDP )

The Taiwan Chelungpu-fault drilling project (TCDP) has undertaken scientific drilling and directly sampled the sub-surface rupture of the 1999 Chi-Chi earthquake. Audio-magnetotelluric (AMT) measurements were used to investigate electrical resistivity structure at the TCDP site from 2004 2006. These data show a geoelectric strike direction of N15 E to N30 E. Inversion and forward modeling of the AMT data were used to generate a 1-D resistivity model that has a prominent low resistivity zone (< 10 ohm-m) between depths of 1100 and 1500 m. When combined with porosity measurements, the AMT measurements imply that the ground water has a resistivity of 0.55 ohm-m at the depth of the fault zone. Time variations in the measured AMT data were observed from 2004 2005 with maximum changes of 43% in apparent resistivity and 18 in phase. The change in apparent resistivity is greatest in the 1000 100 Hz frequency band. These frequencies are sensitive to the resistivity structure of the upper 500 m of the hanging wall of the Chelungpu Fault. The decrease in resistivity over time appears to be robust and could be caused by an increase in porosity and a re-distribution of the groundwater.


INTRODUCTION
The Chi-Chi earthquake occurred on 21 September 1999 in the Western Foothills of central Taiwan (Kao and Chen 2000;Shin et al. 2000;Teng et al. 2001).This M w = 7.6 earthquake produced a 90 km long surface rupture and caused severe damage across Taiwan (Tsai et al. 2001;Wu et al. 2004).The coseismic displacement on the Chelungpu fault was one of the largest ever observed (Ma et al. 2000;Ma et al. 2006).The Taiwan Chelungpu drilling project (TCDP) began in 2003 and resulted in a 2000 m well that recovered cores from the fault zone at A-hole and finished in 2005 with two boreholes (A-hole and B-hole) being completed.The Chelungpu fault that caused the Chi-Chi earthquake was observed in the core at a depth of 1111 m (FAZ1111).Another fault zone (Sanyi Fault -FAZ1710) was observed at depths of 1500 -1710 m.
Since the electrical resistivity of rocks is sensitive to the presence of fluids, geophysical methods that remotely sense sub-surface resistivity, such as Magnetotellurics (MT), can be a powerful tool in investigating the fluid distribution in the shallow crust (Chen and Chen 1998;Fiordelisi et al. 2000;Chiang et al. 2005;Santos et al. 2006).The effectiveness of MT in imaging fault zones has been demonstrated by studies of the San Andreas Fault zone in California, the US and elsewhere (Unsworth et al. 1999;Chen et al. 2007).In magnetotellurics, the depth of exploration increases as the signal frequency decreases.Thus for imaging shallow fault zone structure at the TCDP site, the higher frequency audio-magnetotelluric (AMT) method is the most suitable.In this paper, AMT data collected at the TCDP site from 2004 to 2006 are presented.Spatial and temporal variations are described and interpreted in terms of the tectonic setting.

AUDIO-MAGNETOTELLURIC DATA COLLECTION AND PROCESSING
The audio-magnetotelluric (AMT) method is a frequency domain geophysical technique that utilizes naturally occurring electromagnetic signals to measure resistivity in the upper 1 -2 km of the Earth (Vozoff 1972).A detailed description of the MT and AMT field techniques and data analysis methods is provided by Simpson and Bahr (2005).AMT data was recorded at two stations (T01 and T02) near A-hole and B-hole from 2004 to 2006 at the locations shown in Fig. 1.Three components of the magnetic field (H x , H y , and H z ) and two components of the electric field (E x and E y ) were recorded in the frequency band 10000 -1 Hz using a V5-2000 system manufactured by Phoenix Geophysics.AMT data was recorded on three separate occasions at T01 and T02 (December 2004, September 2005, and May 2006) close to the TCDP-A site.The time series data were converted to frequency domain using a statistically robust algorithm (Jones et al. 1989) and the resulting AMT sounding curves are shown in Fig. 2. Note that all AMT curves show increasing apparent resistivity from 10000 to 30 Hz and a decreasing apparent resistivity from 30 to 1 Hz.This implies that a high resistivity layer is present at the surface and a lower resistivity layer is present at depth.Quantitative resistivity models derived from these AMT data are presented later in this paper.
Time variations in subsurface resistivity may be expected in fault zones (Mazzella and Morrison 1974;Madden et al. 1993;Hanekop and Simpson 2006).Studies of a major fault zone following a major earthquake give an opportunity to look for such effects.Detecting subsurface resistivity changes with AMT is challenging, since the natural electromagnetic (EM) signals used in this survey also vary with time.Care must be taken to determine if temporal changes in the measured data are due to changes in the AMT signals or changes in the subsurface resistivity structures.Figure 2 also shows changes in the AMT data at station T01, compared to the initial AMT measurement at the same location in December 2004.Both apparent resistivity and phase curves show a small variation in the transverse magnetic (TM) mode that is computed from electric currents flowing across geological strike.The transverse electric (TE) mode is computed from electric currents flowing along strike and appears to be less stable at the lowest frequency analyzed.The phase variations are more stable than the apparent resistivity since they are less influenced by bias and static shifts.These variations will be analyzed in detail in the next section.

TENSOR DECOMPOSTION
Before converting apparent resistivity and phase data into a model of true resistivity as a function of depth, it is important to understand the dimensionality of the AMT data.The tensor decomposition approach of McNeice and Jones (2001) was applied to the data and the resulting geoelectric strike directions are shown in Fig. 3.The dominant strike direction is in the range N15°E to N30°E and essentially paral-  lel to the regional geology as expected.Tensor decomposition also computes the degree to which the regional impedance data has been distorted by galvanic effects arising from near surface resistivity structures.This distortion is characterized by the so-called twist and shear angles.At the TCDP site, the twist angle was in the range 5°to -10°and shear angle was in the range 20°to -10°.These values are relatively small and indicate that galvanic distortion is not strong.The TM mode generally resolves lateral variations in resistivity better than the TE mode (Simpson and Bahr 2005).In addition, the TE and TM curves are quite similar at station T02 (Fig. 2) indicating that the electrical structure is relatively 1D (Jones and Garcia 2003) and justifies the application of a 1D AMT analysis at the TCDP site.

INVERSION OF THE AMT DATA
The most important step in AMT data analysis is to convert the frequency domain data into a resistivity model of the subsurface.In this study, a combination of forward modeling and automated inversions were used.The inversion was carried out using the Occam algorithm developed by Constable et al. (1987).This is a regularized inversion that overcomes the inherent non-uniqueness of the AMT inverse problem by seeking the smoothest model that fits the measured data.For MT data of limited bandwidth, the spatial smoothing imposed during inversion may lead to artifacts in the final model.Since it may not be possible to detect the base of a conductive layer, an inversion algorithm will smear the conductor to infinite depth and give a conductance for the layer that is too large.To solve this problem, the approach of Li et al. (2003) was used, with the resistivity fixed below a given depth.The inversion models are shown in Fig. 4 and the data misfit values are listed in Table 1.Note that the 2006 AMT data were collected at station T02, and not at the same location as in 2004 -2005 (T01) making a quantitative comparison impossible for the 2006 data.
The resistivity models show the presence of a low resistivity zone at depth, as anticipated by examination of the AMT sounding curves shown in Fig. 2. The AMT derived resistivity models can be validated by comparison with the well log data measured at the TCDP site (Fig. 4).The resistivity logs exhibit a low resistivity layer (10 ohm-m) in the depth range of 1000 -1500 m.The low resistivity zone in Fig. 4 at 1100 -1500 m is broadly coincident with the Chinshui Shale.Note that this low resistivity zone is not coincident with the Chelungpu fault (FAZ1111), but is located slightly below it at depth of 1200 -1500 m.This observation may be explained by the presence of saline aqueous fluids below the Chelungpu fault.Saline aqueous fluids are good electrical conductors and laboratory studies show that electrical conductivity of brines depends on the salinity, pressure and temperature (Sourirajan and Kennedy 1962).Below 300°C the conductivity of most saline fluids is relatively independent of pressure, while at higher temperatures the conductivity decreases with pressure (Quist and Marshall 1968).The low resistivity of this layer could be due to an anomaly in: (a) temperature, (b) pressure, (c) porosity, or (d) clay content.
Based on well data, the temperature at a depth of 1000 m is expected to be ~33°C (Suppe 1981); however, this temperature cannot explain the observed low resistivity (< 10 ohm-m).Pressure variations are also unlikely to explain this zone of low resistivity.Variations in porosity are an alternative explanation for the low resistivity in this layer.The electrical resistivity of a porous, saturated rock can be estimated with Archie's Law which states that the bulk resistivity of a rock r 0 is given by: (1) where r w is the pore water resistivity, f is the porosity, S is the saturation, and n is an exponent, assumed to be close to unity (Archie 1942).The exponent m is called the cementation factor and varies according to the degree of interconnection of the pore space.Experimental studies show that this parameter typically varies from 1.1 to 2.2 for sedimentary rocks (Doveton 1986) and a mid-range value of 1.5 was used in this study.Measurements on the core recovered from the A-hole by Tanaka et al. (2006) reported porosity values of 8% and 30% at 1000 and 1300 m depth, respectively.Porosity typically decreases exponentially with depth (Rubey and Hubbert 1959;Giles et al. 1998).
However, in this case, the increase in porosity with depth is likely associated with the damage zone associated with the Chelungpu Fault zone.A pore water resistivity 0.55 ohm-m was calculated from the porosity values.These porosity values predict bulk resistivity values of 44 ohm-m at 1000 m and 3.34 ohm-m at 1300 m; by using m = 1.5 in Archie's Law.The final factor that could cause the low resistivity is the presence of clay minerals, either in the Chinshui shale or in the fault gouge.The influence of clay on resistivity is well known and has been observed elsewhere (Xiao and Unsworth 2006).The results show elevated porosity and in terms of the damage zone around the fault, it is necessary for these pores to be full of aqueous fluids to explain the low resistivity.It should also be noted that agreement between the MT derived model and well log resistivity is not as close as might be expected.Explanation of this effect may be due to the fact that different volumes are sampled in each method of measuring subsurface resistivity.The well log gives a measurement of a small volume around the well, while the AMT sounding gives an average over a volume with horizontal extent comparable to the depth of exploration.

EVIDENCE FOR TIME VARIATIONS IN RESISTIVITY
Spatial changes in fluid distribution may occur as a consequence of stress changes during the earthquake cycle.This change in distribution may cause changes in electrical resistivity over time.These effects might be expected at the TCDP site since the Chi-Chi earthquake ruptured the fault in 1999.Detecting these resistivity changes is challenging for any geophysical method, since the changes often occur at significant depth.The TCDP represents a good opportunity to look for these changes, since the active fault is quite close to the surface at a depth of just over 1 km.To determine if temporal changes in resistivity have occurred at the TCDP site, AMT measurements were repeated from 2004 to 2006 at the site.Significant care was taken with AMT data collection in the field to eliminate the possibility that changes in the data collection method could change the measurement.
The repeat AMT measurements shown in Fig. 2 show subtle changes in the measured apparent resistivity and phase from 2004 -2006.These changes are significant compared to the maximum error bars of ± 0.8 ohm-m in apparent resistivity and ± 7 degrees in phase.Overall the apparent resistivity has decreased and the phase has increased, compared with the first measurement.The changes are generally larger than the error bars and suggest that time variations in subsurface resistivity have occurred between 2004 and 2005.Note that the strongest change in apparent resistivity is at high frequency (1000 -100 Hz), indicating that any changes in subsurface resistivity are occurring close to the surface.
The AMT data for each year were inverted using the approach described above and the three models are shown in Fig. 4. As anticipated from the frequencies at which the changes in AMT data are observed, the resistivity changes are localized at relatively shallow depth and correspond to a decrease in resistivity in the hanging wall of the fault.To further quantify these changes, a 1-D AMT modeling approach was used.This eliminates the effects of smoothing, which can add artifacts into the resistivity model.A resistivity model was constructed for the 2004 AMT measurement and the fit is shown in Fig. 5.This model is inconsistent with the 2005 AMT data with significant phase misfit in the bands 10000 -1000 Hz and 100 -10 Hz (Fig. 5).Note that the change in phase is significantly greater than the 4°e rror floor used.To fit the 2005 AMT data, a decrease in resistivity in the upper 500 m of the model is required.Similar resistivity models were obtained when both the TM mode and the invariant (average of TE and TM, Fig. 2) were used.
As can be seen from the 1-D resistivity models, this corresponds to a depth of 0 -1000 m (Fig. 4), and is above the depth of the Chelungpu Fault.The decrease in resistivity likely reflects a change in the groundwater distribution in the hanging wall of the fault.The change in resistivity in the upper 400 m is from 300 to 100 ohm-m.This change can be explained by using Archie's Law with an exponent of m = 1.5 to compute the porosity.Increasing the effective porosity from 1.4 to 3.1% can account for the observed decrease in resistivity.These porosity values are reasonable and in agreement with studies of near surface materials (Giles et al. 1998).Alternatively, the change in resistivity could be explained by maintaining a constant porosity, but increasing the saturation.For example, if the porosity was a constant 5%, then an increase in saturation from 0.16 to 0.5 is needed to explain the reduction in resistivity from 300 to 100 ohm-m.
It is possible that these changes in subsurface resistivity are associated with the earthquake.A decrease in resistivity in hanging wall could be due to fracturing associated with deformation.However, to convincingly prove a causal link, more repeat measurements recorded closer to the date of an earthquake are needed.It should also be noted that other effects can lead to a change in the measured apparent resistivity and phase.These include a change in coupling of electrodes, instrument calibration etc.

CONCLUSIONS
AMT data was collected at the TCDP from 2004 to 2006.The resistivity models derived from the data show the presence of a low resistivity layer in the depth range (1000 -1500 m).This is broadly in agreement with the depths of low resistivity layers observed in the TCDP well logs, although the layer imaged with AMT is significantly thicker.This low resistivity (< 10 ohm-m) zone at the depth 1100 -1500 m is likely due to a combination of clay within the Chinshui Shale and elevated porosity associated with the fault zone.Subtle, yet significant time variations were observed from 2004 -2005 in the hanging wall of the Chelungpu fault.While external origins for these resistivity changes cannot be completely excluded, these changes may reflect a post-seismic re-distribution of the groundwater regime.Further geophysical measurements of the Chelungpu Fault and related features are needed to further elucidate this effect.

Fig. 1 .
Fig. 1.Geological map showing the location of the AMT stations at the TCDP site (modified from Lin et al. 2007).CHF: Chunghwa fault, SYF: Sanyi fault, CLF: Chelungpu fault, STF: Shuantung fault.T01 was 500 m southeast of the TCDP site; station T02 was close to A-hole.

Fig. 2 .
Fig. 2. (a) Apparent resistivity and (b) phase data from AMT sounding at the TCDP site.The dashed lines are the transverse (TE) mode and solid lines are the transverse magnetic (TM) mode data.(c) Time variations of apparent resistivity plotted as percentage change and (d) phase variations in degrees.The curves labeled INV denote the invariant which is an average of the TE and TM modes.
Fig. 3. (a) Variations of geoelectric strike direction as a function of frequency using tensor decomposition (b) geoelectric strike all frequencies showing a dominant direction of N15°E.Geographic north is denoted by 0 degrees.

Fig. 4
Fig. 4. 1-D resistivity models derived by smooth 1-D inversions of the TM mode data.Also shown are resistivity logs and geological sequence (Lin et al. 2007).FAZ1111 was identified as the Chelungpu fault and FAZ1710 as the Sanyi fault.Note that 2006 data are not recorded in the same location as the 2004 and 2005 data.

Fig. 5 .
Fig. 5. Resistivity models and model responses obtained by forward modeling of AMT data at the TCDP site.(a) resistivity models for 2004 and 2005 (b) 2004 AMT data and the best fitting resistivity model.Dots show data values and thin lines denote upper and lower limits with error floor applied.The residuals are defined as data misfit divided by the standard error for each datum.An ideal fit is within ±1 residual; (c) 2004 resistivity model and 2005 AMT data showing that the model is inconsistent with the 2005 data; (d) 2005 AMT data and 2005 resistivity model.