3-D Imaging of Steeply Dipping Structure Near the
San Andreas Fault, Parkfield, California*

John N. Louie, Robert W. Clayton

Seismological Laboratory, California Institute of Technology, Pasadena, CA 91125
Louie is now at the Seismological Lab, University of Nevada, Reno.

Ronan J. Le Bras

Etudes-Production Schlumberger, 26 Rue de la Canee, Clamant 92142, France
Le Bras is now with SAIC, La Jolla, Calif.

* Presented at 54th SEG Annual Meeting, Atlanta.

This paper presents an excellent example of a seismic reflection dataset that cannot be fully interpreted using the stacking process, and its underlying assumptions. It was recorded across a rather spectacular lateral heterogeneity-- the San Andreas fault zone. It will be seen how the interpretations of this survey done previously could not image the strongest reflections in the dataset. Then a method will be described for imaging these events and the less restrictive assumptions that underlie it. Its effectiveness will be demonstrated through imaging of synthetic data. Finally, the imaging of the field dataset will provoke conclusions on the relationship of the geological and physical setting both of the area near the San Andreas and of seismic reflection targets in general.

In 1977 COCORP recorded 27 km of deep crustal reflection data on a route crossing the San Andreas fault in Monterey County, California, near the town of Parkfield. This section of the fault has long been of interest to seismologists because of the regular occurrence of moderate earthquakes on it. Evaluations of velocities and other seismic characteristics of the region have attended several studies of seismic activity, such as those by Eaton et al. (1970) and Liu (1983). Analysis of reflection profiles just to the north along the San Andreas, in San Benito County, by Feng and McEvilly (1983) shows that the fault zone is marked by ``extreme lateral heterogeneity.'' This is principally expressed as relatively low velocities within a zone surrounding the fault a few kilometers wide.

Fig. 1: Map of the vicinity of the town of Parkfield, showing the route of the COCORP survey, and fault traces taken from Hanna et al. (1972).

Figure 1 is a map of the Parkfield area showing the survey route and major fault traces in the vicinity, as mapped by Hanna et al. (1972). The data from the COCORP survey were originally processed and interpreted by Long (1981). He made interpretations on the history of the crustal blocks juxtaposed by the fault, based on characteristics observed in a stacked section. His line drawing of the major events in that section is given in Figure 2. Long observed differences in the density of events in different parts of the stacked section, relating changes in event density to crustal discontinuities. He interpreted diffractions at shallow levels of the fault zone as the effect of structures truncated by processes of brittle fracture. The deeper, ``transparent'' part of the zone represents a region of ductile flow. Because of the poor quality of the stack, these conclusions could not be made firm.

Fig. 2: Line drawing by Long (1981) of a COCORP preliminary stack of the Parkfield survey. The station numbers increase from southwest to northeast, with the San Andreas fault near station 140.

A better approach was undertaken by McBride and Brown (1986). They present a complete reworking of the dataset facilitated by the, previously unavailable, detailed control of the data processing and reduction. Pre-correlation and pre-stack balancing, filtering, and editing contribute to an overall improvement in the stacked section. The density of stacked events was again used to make associations with regionally known crustal structures.

This kind of interpretation is limited in that it is based on the validity of the stacking process, which, as Feng and McEvilly (1983) showed, is thoroughly violated in this region. Further, it attempts to assign geologic interpretations to a physical phenomenon, stacked event density, where little experimental control exists on the relation of particular geologic units to observable reflections. It is difficult to show that the stacked event density is not an artifact of the survey procedure or the data analysis. This is especially true in an area where the assumptions of this method of analysis, stacking, could be invalid.

In fact, the methods used by Long (1981), and McBride and Brown (1986), were not meant to be capable of imaging the strongest, most interpretable reflection events in the dataset. While these events cannot be analyzed by stacking, they can be reduced through a simple, though time-consuming, procedure. A multi-offset Kirchhoff sum imaging process does succeed in showing where the major physical boundaries of the fault zone lie.

A large number of field records distributed along the entire line were examined. By far the strongest, most coherent events observable in the multi-offset data appeared near the San Andreas fault. These events show horizontal or reverse moveouts, in that their arrival times decrease with increasing offset, and have amplitudes comparable to the direct arrivals (Figure 3). They are found on all of the reasonably clean gathers near the fault. The timing and apparent velocity of the events change rapidly as the orientation of the survey line changes. These factors suggested that the events may have originated as sidewall reflections from steeply dipping structures. Similar events were observed by Robinson in 1945 from reflection surveys on the Gulf Coast. He interpreted them as horizontally propagating refractions reflected off lateral discontinuities in the refracting structures, probably faults. From a number of surveys in different orientations, he was able to locate the faults by relating the arrival time of the events to the propagation time along the refractor.

Fig. 3: Consecutive raw common shot gathers recorded by COCORP in 1977. No amplitude balancing has been applied to these correlated vibrator recordings. Each of the nine panels contains 48 traces from a 4.7 km long off-end receiver spread being pushed by the vibrators towards the northeast with a minimum offset of 0.45 km. The sources are 0.1 km apart. In the center panel, the receiver closest to the San Andreas fault trace is indicated by the vertical arrow. The open arrows indicate the strongest arrival on the records. Its arrival time decreases with increasing offset.

Some compensation for the sinuosity of the COCORP line can be made by sorting out common midpoint (CMP) gathers of the data. All traces whose midpoint fell within 230 m of a node of a two-dimensional grid of points with the same spacing were sorted into a gather for that node, regardless of the orientation of the shot-receiver pair. Two of these gathers are shown in Figure 4. The unusual arrivals at about 2 seconds at the farther offsets can be found on many gathers. As the figure shows, the apparent velocity of the sidewall reflection changes drastically as the survey crosses the San Andreas. On the southwest side, the receivers are between the vibrators and the fault, so the moveout is negative. On the northeast side, the receivers are farther from the fault than the vibrators, so the moveout is normal.

Fig. 4: Two common midpoint gathers of unequalized traces assembled within two-dimensional bins 460 m wide at different locations near the San Andreas. The open arrows indicate the positions of a strong reflection from the steeply dipping fault zone. The timing and apparent velocity of the sidewall reflection change as the geometric relationship of the shots and receivers relative to the fault changes.

Such events, especially where they have negative moveout, can obviously not be stacked using any physically meaningful stacking velocity. The stacking process would destroy their coherency, rendering them invisible in a stacked section. On the other hand, where the line is oriented such that the sidewall reflection has a normal, positive moveout, it may stack coherently, but its location in the section will be completely incorrect. If the reflection point is not in the plane of the survey line, it will not, in fact, be possible to migrate the stacked reflection to its correct location. Yet some process of imaging the reflector producing these events must be found, since they carry most of the energy in the seismic gathers. They therefore represent the most fundamental physical boundaries in the area.

In obtaining this kind of image, it will be necessary to have a starting idea of the velocity structure in the area. The CMP gathers made from the southwest part of the line did show coherent reflections from near-horizontal structures in the upper 3 km of the crust. Interval velocities calculated from velocity semblances of these gathers indicated a strong velocity gradient in this area similar to that found by Liu (1983) from seismicity analysis. The gradient is shown in Figure 5, which also incorporates deeper velocity information derived from refraction surveys by Eaton et al. (1970). These profiles suggest that the strong lateral heterogeneity across the San Andreas is limited to velocities in the upper 5 km of the crust.

Fig. 5: Crustal velocity profiles derived for the blocks on each side of the San Andreas from refraction surveys (Eaton et al., 1970) and as modified by Liu (1983) with arrival times from the local seismicity. The gradient at the top of Liu's profile was adjusted to agree with velocities derived from the moveout of reflections recorded by the COCORP survey to the southwest of the San Andreas.

The velocity gradient in the uppermost crust explains how reflections from a vertically oriented structure could be recorded by a horizontally oriented receiver spread. Horizontal bending of the raypaths with depth, within such a strong velocity gradient, assures that reflections can be located on structures dipping even more than 90 degrees from the horizontal. The ray bending will, unfortunately, also act to limit the range of depths covered by the recorded reflections.

Approximations to the wave equation

With these three approximations the data can be considered to be a linear superposition of rays from individual point scatterers. The tomographic approximation of the inverse of this superposition, as discussed by Le Bras (1985), is simply the superposition of rays from individually recorded reflections. Thus, the scattering potential of the medium can be estimated as the sum of the reflections recorded by each source-receiver pair, positioned according to the travel time of the rays between the surface points and the subsurface reflector.

Since the purpose here is simply to establish the geometry of the scatterers within the medium, we will ignore the amplitude correction factors due to the angle of incidence on the scatterer, and to the length of the travel path. Further, the scatterer will be represented by the sum of reflection wavelets without cross correlation with the source wavelet, since a source wavelet is not available. Certain restrictions will apply. The data put into the inversion should contain only primary P-to-P wave reflections. In addition, the inversion will not be valid for post-critical angle reflections, or refractions.

With these approximations, the Kirchhoff summation method used here to image the geometric distribution of acoustic reflectivity is very similar to that used by McMechan and Fuis (1987), and outlined by Jain and Wren (1980). Figure 6 shows the geometry of a three-dimensional reflection from a steeply dipping fault zone within a vertical velocity gradient. This method is especially versatile in that the reflectivity of any depth point may be inverted from data recorded from sources and receivers at any location. The ordering of the data and of the inversion points are immaterial, since the tomographic sum may be made in any order.

Fig. 6: Geometry of a sidewall reflection from a fault zone within a medium showing an increase in velocity with depth Z. S and g are the source and receiver positions, respectively.

Figure 7 summarizes the imaging procedure. Unsorted seismogram traces are mapped into a depth section by computing the travel time from the source to the depth point and back to the receiver. A velocity model that varies only with depth is used, as derived from the velocity semblances and the model of Liu (1983; Figure 5). The travel time calculation included turning rays, which allows the imaging of structures with greater than 90 degree dips. To allow for propagation through a laterally heterogeneous medium, the travel time calculation could take the form of raytracing through a variable velocity medium. If, however, the bulk of the travel path can be constrained to a part of the medium in which velocity varies principally with depth, then a simple raytracing through a vertically varying medium can be used for data with any orientation of the source-receiver offset. This allows the travel times to be calculated only once, for the range of offset and depth of the experiment, greatly speeding the imaging process. For this reason the migrations performed here will employ mainly sources and receivers on just one side of a major lateral discontinuity such as the San Andreas.

Fig. 7: Schematic representation of of a Kirchhoff sum algorithm used to form a migrated image of multi-offset reflection data. For each recorded source-receiver pair (S, g), over each point (X, Z) in some depth section, a travel time is calculated. The amplitude of the recorded trace at that time is then summed into the point of the depth section.

Once the travel time down to and up from the depth point has been obtained, the value of the seismogram at that time is summed into the section at the depth point. Spurious arcs due to noise bursts and badly gained traces are easily identified with the help of plots of the wavefront shape for the given velocity model. Since small sections in areas of particular interest can be migrated one at a time, and storage of large numbers of traces is not necessary, this method is economical even on a relatively small computer. If the seismic survey has sufficient 3-D coverage, the reflectors can be easily imaged in 3 dimensions by properly locating the depth sections.

To test the method, synthetics of a simulated survey over an idealized model of a steeply dipping fault zone were calculated. A finite difference solution of the two-dimensional acoustic wave equation was used. This solution included all acoustic multiples, post-critical reflections, and refractions. The velocity model was identical to the one identified in Figure 5, except that a 2 km-wide fault zone having a 20% lower velocity at a given depth, and a sinuous profile in cross section, was introduced to test the ability of the method to resolve vertically complex fault geometries (Figure 8).

Fig. 8: Two-dimensional velocity model representing a cross section 25 km wide by 25 km deep across a fault zone having a sinuous profile. The zone, 2 km wide, has a velocity 20% lower than the adjacent regions. Outside the fault zone the model is that of Liu (1983; Figure 5), which varies only with depth. The dashed box is region imaged in Figures 11 and 12.

Synthetics were calculated with shot and receiver spacings meant to simulate the Parkfield survey, but with poorer coverage. The synthetic gathers, given in Figure 9, show strong reflections from both sides of the near-vertical fault zone, which are quite similar to the arrivals in the Parkfield data. The travel times calculated for the velocity model to points at different depths and offsets, and used for the inversions, are given in Figure 10.

Fig. 9: Two-dimensional acoustic finite difference synthetics generated from the velocity model of Figure 8. The geometry of the synthetic spread is similar to the geometry of the COCORP survey. The common shot gathers are presented here the same way as in Figure 3. Note the reflection showing normal hyperbolic moveout at 2.5 s from shots 1 through 5, and the strong sidewall reflections from the fault zone showing negative moveout and intersecting the first arrival.

Fig. 10: Matrix of travel times from a source in the upper left corner to any point in a cross section of depth Z against offset $DELTA bold X$, calculated by ray tracing through the vertically varying velocity model of Liu (Figure 5). Darker areas have the longer travel times, with the contours on equal travel time.

The effect of the Kirchhoff sum migration can be illustrated with a migration of a single trace. This image is shown by Figure 11. Each part of the trace, after muting of the first arrival, has been back projected into the depth section along contours of equal travel time similar to those in Figure 10. At least one point of this projected image is correctly located. As the back projections of more traces are summed into the image, the correctly located point should be reinforced, and the incorrectly located parts of the arcs canceled by destructive interference.

Fig. 11: Depth section located in dashed box of Figure 8 with the Kirchhoff sum migration of a single trace. Note similarity to shape of travel time contours of Figure 10.

Summing in the back projections of all 15 of the 48 trace synthetic gathers produced the image in Figure 12. This image should be a reconstruction of the part of the velocity model in Figure 8 set off by the dashed line. The sinuous geometry of the zone has been quite well reconstructed. In the lower third of the image, however, the reconstruction is not as complete due to the lack of reflection points on the fault zone at those depths. Because of the strong velocity gradient in the first few kilometers of depth, most of the rays turn horizontally or refract at shallow depths. In fact, some of the most strongly reconstructed points lie along refractors that prevail at particular depths.

Fig. 12: Full Kirchhoff migration of 720 traces of 15 synthetic shot gathers within the dashed box of Figure 8. Note the resolution of the sinuous shape of the fault zone, including portions imaged by rays bent more than 90 degrees from vertical.

The reflection data are back projected into four depth sections, shown on Figure 13. The sections were located where there are heavier concentrations of midpoints, with B and D made parallel to test three-dimensional aspects of the image. All back projections were made using travel times calculated from the velocity profile in Figure 5, which is most appropriate for the region to the southwest of the San Andreas. In migrating all of the sections (except for C) only traces having both sources and receivers to the southwest of the Gold Hill fault were used. Section B was migrated both from traces having the first arrivals muted and from traces without any mutes. The two sections showed little difference, so all of the migrations were run on unmuted data. However, ignoring the nearest offset traces, as done by McBride and Brown (1986), did improve the imaging.

Fig. 13: Map of the vicinity of Parkfield showing the locations of the four depth sections migrated from the COCORP shot records. Each section is 5 km long by 5 km deep.

The back projected depth sections are shown in Figure 14. Many parts of these images are artifacts. Where the geometric coverage is poor, because of the layout of the line or the concentration of raypaths along refracting structures, the arching tails of the individual back-projected events may not be canceled. It is most useful to look for strong images that are at least crossed by similar arcs. Among the images in which we have some confidence is, in Section D, a reflector dipping at least 45 degrees to the southwest underneath the San Andreas fault, possibly extending to the surface near the trace of the Gold Hill fault. It can be found to a depth of at least 4 km. Section B contains a similar reflector, although it is less well defined. This reflector also appears on Section C, with a shallower apparent dip. A subhorizontal reflector is observed at a depth of about 3 km between the surface traces of the San Andreas and an unnamed fault to the southwest on section B. Such a reflector is also present in Section D. A strong near-vertical reflector is shown in Section B, extending from the surface to a depth of about 1 km. This reflector can also be discerned in Section A, where it appears to mark the truncation of a strong sub-horizontal reflector at a depth of 2 km, which continues to the southwest.

Fig. 14: Kirchhoff sum migrations of the COCORP shot records in the depth sections shown in Figure 13. The positions on the surface of the unnamed fault to the southwest of the San Andreas, SWF, the San Andreas fault, SAF, and the Gold Hill fault, GHF, are indicated by arrows. Other arrows indicate the reflectors interpreted in the text.

Hanna et al. (1972) have synthesized the surface geology of the Parkfield area with gravity and magnetic data. Figure 15 interprets the location of the imaged reflectors and their association with mapped surface features and suspected geological relations. In making interpretations of the connection of reflectors with geological boundaries, the most obvious, simple boundaries that would provide the greatest velocity contrast should be stressed. The Gold Hill fault incorporates slices of serpentinized ultramafic rocks from the metamorphic Franciscan complex to the northeast and juxtaposes them against crystalline rocks of intermediate to mafic composition, which crop out to the southwest. Such a contrast would create a strong reflector. Sections B and D indicate that the fault may well dip steeply to the southwest and intersect the vertically dipping San Andreas at a depth of 3 to 4 km. Between the San Andreas and the unnamed fault to the southwest, the intermediate-ultramafic contact may be sub-horizontal at a depth of 3 km. It appears to be truncated by the southwestern fault against the possibly granitic basement of the Salinian block. The same fault also may truncate a 2 km thick Tertiary sedimentary section to the southwest against shallow crystalline rocks of intermediate composition in the fault zone, producing the strong vertical reflector. A strong reflection is apparently not observed from the active trace of the San Andreas itself, indicating that it may cut through relatively uniform mafic to intermediate rocks caught in the fault zone above 3 km.

Fig. 15: Perspective block diagram showing mapped features of the Parkfield region in relation to a sectional view of the positions of the major reflectors imaged by the Kirchhoff sum migrations of the COCORP shot records.

The use of data from the southwest side of the San Andreas, combined with the presence of granitic and sedimentary rocks in the fault zone having affinities to the block to the southwest, allowed the laterally homogeneous velocity model to produce fairly accurate reconstructions, at least above the ultramafic rocks. Most of the major reflectors are probably located to within $+-$0.5 km. Any imaging of deeper reflectors will require the use of data far enough away from the fault zone to avoid the complications in the upper several kilometers.

It has been shown that a three-dimensional Kirchhoff sum migration before stack is capable of imaging steeply dipping reflectors, which produced arrivals that cannot be stacked with conventional CMP methods. Such reflections can be observed in the data from a 1977 COCORP survey across the San Andreas fault near Parkfield, California. Previous workers, who used conventional methods, could not interpret the events, which arose from strong lateral heterogeneities. However, applying the Born, WKBJ, and far-field approximations to the wave equation results in the simplification of the inversion of these reflections to a ray equation back-projection process very similar to Kirchhoff sum migration. This process can be easily implemented on unsorted data to back project reflector images into arbitrarily oriented depth sections of limited size. This process was verified by inversions of full wave acoustic synthetics incorporating the same geometric coverage limitations inherent in the COCORP survey. Although limited in accuracy by poor geometric coverage and the laterally inhomogeneous velocity structure, the method imaged reflectors around the San Andreas fault that are consistent with the known structural features of the area. The imaged relationships suggest that the fault zone has juxtaposed mafic to intermediate crystalline rocks against Cenozoic sediments on the southwest and serpentinized ultramafics on the northeast. The modern trace of the San Andreas does not form such a strong reflector that it can be imaged. Such information can prove useful in efforts to reconstruct the complex history of the motion of crustal blocks caught in the transform zone. Further work will be aimed at testing the method on additional structures where data become available, and at taking advantage of the amplitude information available in the multi-offset domain to further constrain the inversion.

Contribution number 4475, Division of Geological and Planetary Sciences, California Institute of Technology.

Eaton, J. P., O'Neill, M. E., and Murdock, J. N., 1970, Aftershocks of the 1966 Parkfield-Cholame earthquake: a detailed study: Bull., Seismol. Soc. Am., 60, 1157-1197.

Feng, R., and McEvilly, T. V., 1983, Interpretation of seismic reflection profiling data for the structure of the San Andreas fault zone: Bull., Seismol. Soc. Am., 73, 1701-1720.

Hanna, W. F., Burch, S. H., and Dibblee, T. W., 1972, Gravity, magnetics, and geology of the San Andreas fault area near Cholame, California: U. S. Geol. Surv. Profess. Paper 646-C.

Jain, S., and Wren, A. E., 1980, Migration before stack-- procedure and significance: Geophysics, 45, 204-212.

Le Bras, R., 1985, Methods of multiparameter inversion of seismic data using the acoustic and elastic Born approximations: Ph.D. thesis, Calif. Institute of Technology.

Liu, H.-L., 1983, Interpretation of near-source ground motion and implications: Ph.D. thesis, Calif. Institute of Technology.

Long, G. H., 1981, A COCORP deep seismic reflection profile along the San Andreas fault, Parkfield, California: M.Sc. thesis, Cornell Univ.

McBride, J. H., and Brown, L. D., 1986, Reanalysis of the COCORP deep seismic reflection profile across the San Andreas fault, Parkfield, California: Bull., Seismol. Soc. Am., 76, 1668-1686.

McMechan, G. A., and Fuis, G. S., 1987, Ray equation migration of wide-angle reflections from southern Alaska: J. Geophys. Res., 92, 407-420.

Robinson, W. B., 1945, Refraction waves reflected from a fault zone: Geophysics, 10, 535-545.

Wu, R., and Aki, K., 1985a, Scattering characteristics of elastic waves by an elastic heterogeneity: Geophysics, 50, 582-595. .in -0.4i

scanned and converted 9/26/96 by J. Louie.