Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Teleseismic imaging of the southeastern Canadian Shield and Cascadia subduction zone Rondenay, Stéphane 2000

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


831-ubc_2001-611663.pdf [ 17.4MB ]
JSON: 831-1.0053120.json
JSON-LD: 831-1.0053120-ld.json
RDF/XML (Pretty): 831-1.0053120-rdf.xml
RDF/JSON: 831-1.0053120-rdf.json
Turtle: 831-1.0053120-turtle.txt
N-Triples: 831-1.0053120-rdf-ntriples.txt
Original Record: 831-1.0053120-source.json
Full Text

Full Text

TELESEISMIC IMAGING OF THE SOUTHEASTERN CANADIAN SHIELD AND CASCADIA SUBDUCTION ZONE by STEPHANE RONDENAY B.Ing., Ecole Polytechnique de Montreal, 1994 M.Sc.A., Ecole Polytechnique de Montreal, 1996 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in THE FACULTY OF GRADUATE STUDIES (Department of Earth and Ocean Sciences) We accept this thesis as conforming to the required standard THE UNIVERSITY OF BRITISH C O L U M B I A December 2000 © Stephane Rondenay, 2000 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for refer-ence and study. I further agree that permission for extensive copying for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of Earth and Ocean Sciences Hie University of British Columbia 2219 Main Mall Vancouver, BC, Canada V6T 1Z4 Date: f W . G[U , 2COO 11 Abstract This thesis is divided into two parts that illustrate different approaches to multichannel analysis of broadband, teleseismic array data. In Part I, a range of established techniques is applied to data from a teleseismic experiment conducted across the southeastern Canadian Shield. Multi-event SiTS-splitting results yield an av-erage delay time of 0.57±0.22 s and a direction of fast polarization of N93° E ± 18° for the region, which is consistent with an earlier interpretation of fossil strain fields in the lithosphere. Profiling of radial receiver functions reveals an abrupt, northward thinning of the crust, some 65 km south-east of the Grenville Front. This structure is interpreted as a subduction suture extending the full length of the Front and punctuating a major, pre-Grenvillian episode of lithospheric assembly. P-and S-velocity models derived from teleseismic travel-time inversion outline a low-velocity, NW-SE striking corridor that crosses the array at latitude 46°N and lies between 50 and 300 km depth. This feature postdates stabilization of the craton and is postulated to represent interaction of the Great Meteor plume with zones of weakness developed during earlier rifting episodes. Part II demonstrates the application of a new technique for formal, multi-parameter, 2-D in-version of scattered waveforms. Upon discussing a range of practical and numerical considera-tions, the method is applied to data from the CASC93 experiment undertaken in central Oregon to investigate the structure of the Cascadia subduction zone. Two major features are imaged in the resulting model. The continental Moho appears at approximately 150 km from the coast as a positive velocity contrast near 40 km depth. The layers of the eastward dipping Juan de Fuca plate continue from 20 km depth below the coast to at least 100 km beneath the High Cascades. A disruption of oceanic crust is imaged near 40 km depth, below which crustal expression is con-siderably diminished. This depth coincides with an increase in plate dip, in addition to enhanced mantle conductivity and seismic reflectivity documented in previous geophysical studies. I inter-pret these observations to represent the effects of metamorphic reactions within the oceanic crust which lead to dehydration and eclogitization. Table of Contents Abstract ii List of Tables viii List of Figures ix Acknowledgements xi Dedication xii 1 Introduction 1 1 TELESEISMIC IMAGING OF THE LITHOSPHERE BENEATH THE ABITIBI-GRENVILLE TRANSECT 4 2 Overview of Part I 5 2.1 Author's contribution 7 3 Geology and geophysical coverage of the study area 8 3.1 Tectonic setting 8 3.2 Geophysical coverage 11 3.2.1 1994 teleseismic project 14 Shear-wave splitting analysis 14 ii i Table of Contents iv Receiver function analysis 18 4 Analysis of the Abitibi 1996 teleseismic data 20 4.1 Data collection 21 4.1.1 Abitibi 1996 experiment 21 4.1.2 Permanent stations 24 4.1.3 Event selection 24 4.2 Travel-time inversion 25 4.2.1 Method 25 4.2.2 Data set 28 4.2.3 Results 28 4.3 Shear wave splitting analysis 35 4.3.1 Method 35 4.3.2 Data set 35 4.3.3 Results 36 4.4 PS conversions 39 4.4.1 Method 39 4.4.2 Data set 40 4.4.3 Results • 40 5 Interpretation of the results 43 5.1 Variations in crustal thickness 43 Table of Contents v 5.1.1 Mono step and Grenville Front 46 5.1.2 Age of subduction 47 5.2 Nature of the low-velocity mantle corridor 49 5.2.1 Proposed mechanisms for the M W N track 50 5.2.2 Physical properties and low velocity 52 5.3 Concluding remarks 56 II MULTI-PARAMETER 2-D INVERSION OF SCATTERED TELESEISMIC BODY-WAVES: IMPLEMENTATION AND APPLICATION TO THE CASCADIA-93 DATA SET 57 6 Overview of Part II 58 6.1 Theory 58 6.2 Numerical examples 66 6.3 Research goal and overview 66 6.4 Author's contribution 67 7 Considerations for practical implementation 68 7.1 2-D travel times and scattered wave amplitudes for oblique incidence 68 7.1.1 Travel times 68 7.1.2 Scattered wave amplitude 70 7.2 Interpolation 70 Table of Contents vi 7.3 Inversion 71 7.3.1 Choice o f inversion parameters 71 7.3.2 Choice of scattering modes 72 7.3.3 Choice of Jacobian 72 7.4 Preprocessing of the data 74 8 Application to CASC93 76 8.1 Tectonic framework 76 8.2 Previous geophysical studies 78 8.3 Experiment and data set 82 8.4 Reference model 88 8.5 Results 88 8.6 Interpretation of the low velocity crustal disruption . 94 8.6.1 Delamination-tectonic layering 94 8.6.2 Aqueous fluids and hydrated minerals 97 8.7 Concluding remarks 99 9 Conclusions 100 References 102 Appendices 113 A Data set for analyses in Chapter 3 113 Table of Contents B Scattering potentials List of Tables 4.1 List of participants for Abitibi 1996 teleseismic experiment 21 4.2 Abitibi 1996 teleseismic stations 23 4.3 SKS splitting results 38 8.1 Events selected for scattered waveform inversion 84 8.2 Reference, 1-D velocity model for CASC93 89 A. 1 Events selected for analyses along the Abitibi-Grenville Transect 114 viii List of Figures 3.1 Tectonic setting and station coverage 9 3.2 Electrical and seismic anisotropies 13 3.3 Average 1-D S-wave velocity model from receiver function analysis 17 4.1 Data sets for teleseismic analyses along the Abitibi-Grenville Transect 25 4.2 Travel-time inversion grid 26 4.3 Map view of the high resolution region for travel-time tomography 27 4.4 Travel-time tomography - resolution test 29 4.5 Travel-time tomography - P-velocity model 31 4.6 Travel-time tomography - S-velocity model 32 4.7 1996 shear-wave splitting results 37 4.8 Profiling of radial (SV) impulse responses 41 5.1 Compilation of the Abitibi-Grenville teleseismic results 44 5.2 Schematic diagram of a Tectonic Wedge 45 6.1 Quasi-linear array of receivers above a 2-D heterogeneous lithosphere 59 6.2 Scattering modes considered in the inversion method 61 6.3 Geometrical quantities considered for scattered teleseismic waveform inversion 62 ix List of Figures x 6.4 Forward and inverse scattering problems 64 8.1 Tectonic map of the Cascadia subduction zone in western Oregon and Washington 77 8.2 Geophysical studies of central Oregon 80 8.3 Distribution of events selected for multichannel inversion 83 8.4 Preprocessing of three-component recording for event #15 85 8.5 Three-component preprocessed data sections for events #15 and #16 87 8.6 Results from single mode inversion of events #15 and #16 91 8.7 Results from single mode inversions of all events 92 8.8 Results from simultaneous inversion of all events 93 8.9 Schematic diagrams representing delamination and dehydration of the oceanic crust 95 xi Acknowledgements First and foremost, I wish to thank my thesis advisor, Michael Bostock, for his guidance, patience, unfailing support, and for always making time to discuss work in progress and catastrophes on the mend... anytime before 16:30. Thanks to the earthquake seismology group, Charly Bank, Andy Frederiksen and Jeff Shragge, for invaluable help with computer code, impromptu tutorials on analysis techniques, and continuously stimulating discussions. Thanks to the members of my supervisory committee, Bob Ellis (co-supervisor), Ron Clowes and Greg Dipple, for providing insightful and constructive input at various stages during the progress of this research. Thanks to Bruce Buffett, Andrew Calvert and Gene Humphreys for their constructive reviews. I am grateful to everyone who participated, directly or indirectly, in the Abitibi 96 project: Don White from the GSC; Tom Hearn, Anca Rosea and Richard Rapine from N M S U ; Bernard Giroux, Hua Wu and Chantal Balthazard from Ecole Polytechnique; Sid Hellman, Paul Friberg, Dave Lentrichia, Carl Ebeling, John Weber and Tom Jackson from PASSCAL; Rick Benson, from IRIS-DMC, for his tremendous help with the Abitibi, CASC93 and LARS data sets; Jean Goutier fromMRNQ. Thanks to everyone in the department for four wonderful years in Vancouver filled with back-country tripping, geophysical air forcing, hiking, jamming, single malt sipping, running, swim-ming, and BOLFing. This research was supported financially by postgraduate scholarships from the Natural Sci-ence and Engineering Research Council of Canada (NSERC) and the Fonds pour la Formation de Chercheurs et l 'Aide a la recherche (Quebec), and NSERC grants to Michael Bostock and Bob El -lis. The Abitibi 96 teleseismic experiment was co-funded by L i T H O P R O B E , the National Science Foundation and the Geological Survey of Canada. xii A ma famille, et a, la plus grande fan de MC*Aubergine. Chapter 1 Introduction Earthquake seismology, the study of earthquakes and the elastic waves they generate, provides a wealth of information on deep-seated, subsurface structure. The earthquake itself results from a sudden readjustment of stresses within Earth materials, producing an elastic disturbance that ra-diates away from the focus in the form of compressional (P) and shear (S) waves. Eventually, this energy reaches the Earth's surface where ground motion (displacement, velocity or acceleration) is recorded as a function of time on seismograms. Seismograms are generally more complex than the time history of stress release near the hypocenter, because the initial P- and S- waves travel at different speeds and may take a variety of paths to reach a given receiver. Along these paths, the waves encounter variations in material properties where they may be reflected, refracted, at-tenuated, diffracted, decelerated or accelerated. Characterization of deep earth structure is thus accomplished by decoding the complexity of the seismic waveforms that reach the surface. Due to the great depth of investigation it affords, earthquake seismology has been instru-mental in defining the anatomy of our planet, from the solid inner core to the thin outer shell forming the crust [see, for example, Lay and Wallace, 1995; their Chapter 1]. At finer scales, seismology is employed to image detailed crust and upper mantle structure related to past and present tectonic processes. To this end, teleseismic techniques for imaging Earth structure in-volve the analysis of waves generated by distant earthquakes (> 30° epicentral distance) and can be divided into two main categories depending on the nature of the phases analysed [cf , Bostock, 1 Chapter 1. Introduction 2 1999]. Teleseismic phases can be classified as either: (1) transmitted, in which case the analysis focusses on the time required for a primary wave to travel from source to receiver; or (2) scat-tered, in which case the analysis focusses on the secondary waves produced through interaction between a primary wave and subsurface heterogeneity. In teleseismic studies, the choice of recording apparatus (i.e., type of seismometers and ar-ray configuration) is governed by both the scale of structure under investigation and the tech-niques/phases used in the analysis. Seismometers can be classified in terms of number of com-ponents (i.e., the number of mutually orthogonal directions in which ground motion is meas-ured), and the frequency range of the signal they register. Single component seismometers are commonly short-period (i.e., period response centered near 1 s) and oriented to measure vertical ground motion. Resulting seismograms are used for P-traveltime (i.e., transmission) analyses. Three-component seismometers, in contrast, record the vector motion at the surface, which is suit-able for both transmission and diffraction/scattering analyses; these seismometers can either be short-period or broadband. Broadband recordings, as the name implies, register a greater propor-tion of the signal's frequency spectrum, thus allowing more comprehensive analysis of the wave-form (e.g., including S, SS, ScS, SKS, and surface waves). Furthermore, the level of geometrical complexity that can be characterized using teleseismic phases is intimately linked to the surface distribution of seismometers. While single stations are sufficient for imaging planar discontinu-ities in the subsurface [e.g., Langston, 1979; Revenaugh and Jordan, 1991; Bostock, 1998], linear and 2-D grid arrays are better suited for identifying more complex lateral heterogeneity [e.g., Van-Decar, 1991]. The work presented in this thesis focuses on the analysis of teleseismic data recorded on three-component broad-band arrays. The thesis is divided into two parts, illustrating two differ-ent aspects of teleseismic imaging of regional crustal and upper mantle structure. Part I provides Chapter 1. Introduction 3 an account of a regional study undertaken along the LiTHOPROBE Abitibi-Grenville Transect, where an array was specifically designed for the application of well established/proven analysis methods. These techniques include travel-time tomography, shear-wave splitting analysis, and receiver function profiling. The resulting images characterize lithospheric structure associated with the accretion of the southeastern Canadian Shield and its subsequent post-stabilization re-working. Part II presents the numerical implementation of a new teleseismic analysis technique, and its application to a field data set. The technique permits the retrieval of elastic properties from complex 2-D structures beneath dense receiver arrays using scattered energy in the teleseis-mic .P-wave coda. It is tested on a data set recorded in 1993-1994 across central Oregon to image the Cascadia subduction zone. This choice of data set was dictated by the explicit geometrical requirements (e.g., dense station spacing, 2-D regional strike, quasi-linear array) on which the method relies. The results reproduce the main features imaged previously beneath the study area, and resolve more subtle structure associated with the subduction process. Parti TELESEISMIC IMAGING OF THE LITHOSPHERE BENEATH THE ABITIBI-GRENVILLE TRANSECT Chapter 2 Overview of Part I Cratons are generally associated with the idea of stability, a property which is due, at least in part, to the deep lithospheric roots that render them resistant to major structural reworking [Jordan, 1978; Hoffman, 1990]. The detailed architecture of the roots is currently a topic of great interest as it holds the clues to two important issues in solid Earth science: (1) the manner of assembly and stabilization of early continental landmasses, and (2) the nature of processes which have served to modify cratonic lithosphere in its poststabilization evolution, for example, those responsible for kimberlite volcanism or the emplacement of giant mafic dike swarms. Over the past 2 decades, considerable progress has been made on these issues using a variety of tools. Near-surface studies (geological mapping, geochemistry, geochronology, and paleomag-netism) have supplied strong evidence that plate tectonics have been active as far back as Archean time [Ludden and Hubert, 1986; Hoffman, 1988; Wit et al., 1992; Kimura et al, 1993]. Re-cent, high-resolution seismic experiments provide supporting documentation in the form of im-ages of mantle reflectors merging with overlying crust in the vicinity of ancient continental sutures [BABEL Working Group, 1990; Calvert et al, 1995; Warner et al, 1996; Cook et al, 1998]. A complex mantle stratigraphy apparently persists throughout the continental lithospheric column as revealed by deeper-probing seismological investigations [Revenaugh and Jordan, 1991; Bo-stock, 1998]. Widespread, i f infrequent, Phanerozoic kimberlite volcanism indicates that continental roots 5 Chapter 2. Overview of Part I 6 are not completely static entities but continue to evolve in the present day. Some authors [Sykes, 1978; Crough et al, 1980] have suggested a genetic relation between kimberlites and hotspot tracks, here defined as linear, anorogenic volcanic chains mapped over oceanic and young conti-nental lithospheres. The termination of these features over thickening continental roots has been cited as evidence supporting the origin of hotspot tracks in the interaction of fixed mantle plumes with an overriding plate [Morgan, 1972; Burke and Wilson, 1976; Crough et al, 1980; Sleep, 1990a; Duncan and Richards, 1991], although compelling arguments may also be made for as-sociation of hotspot tracks with ancient rifts and sub vertical lithospheric discontinuities [Currie, 1976; Sykes, 1978; Anderson, 1998]. Part I of this thesis presents results from a teleseismic experiment conducted across the Su-perior and Grenville Provinces of the Canadian Shield. This study reveals lithospheric structures that are apparently related to both craton assembly during a major Precambrian compressional event and subsequent alteration associated with the Cretaceous Monteregian hotspot track. Chap-ter 3 provides a brief tectonic overview of the study area and a summary of previous geophysical work in the region, with special attention paid to teleseismic analyses performed on a preliminary small-scale passive array deployed in 1994. Chapter 4 focuses on the acquisition of the 1996 tele-seismic dataset, followed by a detailed account of the different analysis techniques applied to the data and their results. Part I concludes with Chapter 5, which discusses the implications of the results for our understanding of the tectonic evolution of the southeastern Canadian Shield. Chapter 2. Overview of Part I 2.1 Author's contribution 7 The research performed in Part I of this thesis was presented in two refereed publications and a series of poster and oral communications at scientific meetings. Chapter 3 contains mainly ex-cerpts from a paper synthesizing teleseismic results obtained prior to 1 9 9 8 along the LITHOPROBE Abitibi-Grenville Transect [Rondenayetal, 2000b]. Among these previous studies, I participated in the analysis and interpretation of the 1 9 9 4 teleseismic data set as part of my M.A.Sc. project at Ecole Polytechnique de Montreal [Rondenay, 1 9 9 6 ] . The synthesis work presented here in-volved an extensive literature review pertaining to ( 1 ) the tectonic history of the study area, and (2) previous geological/geophysical results obtained in the region. This information was used to tie together the teleseismic results, that had been obtained up until that time from completed and ongoing studies, with the tectonic framework of the southeastern Canadian Shield. The paper was written by myself, with important contributions from Donald J. White and Hua Wu to the section on receiver function analysis. Chapters 4 and 5 are derived from a paper presenting the fi-nal results from analyses performed on the Abitibi 1 9 9 6 data set and their tectonic interpretation [Rondenay et al, 2000a]. I conducted every aspect of this research, from the teleseismic data acquisition in the field, to the preprocessing of the data, its analysis using conventional analysis techniques (or variations thereof) and the interpretation of the results. Chapter 3 Geology and geophysical coverage of the study area The main data set considered in Part I of this thesis was collected along a teleseismic array within the southeastern portion of the Canadian Shield (Figure 3.1). This region has been the focus of extensive investigation from most branches of solid earth sciences, as part of the LITHOPROBE Abitibi-Grenville transect [Clowes, 1997]; for a general and more extensive review, see Ludden and Hynes, 2000; and references therein]. In this chapter, I provide a brief tectonic overview of the study area, followed by a summary of pertinent geophysical studies that have been conducted to date in the region. 3.1 Tectonic setting The study area includes a diversity of geological domains ranging in age from late-Archean in the north to Paleozoic in the south, and a number of major tectonic boundaries. The north-ern half of the study area comprises the Superior Province of the Canadian Shield, the largest Archean province in the world [Card, 1990]. The Superior Province consists of a sequence of generally east-west trending volcano-plutonic (granite-greenstone), metasedimentary, plutonic, and high-grade gneiss subprovinces [Card and Ciesielski, 1986]. The Opatica plutonic belt is the northernmost unit of the Superior Province terranes shown in Figure 3.1. It is interpreted as the deeply eroded core of an Archean orogen [Benn et al, 1992; Sawyer and Benn, 1993], onto 8 Chapter 3. Geology and geophysical coverage of the study area 9 * ABI-96 teleseismic stations A MT stations Figure 3.1: Simplified tectonic map of the study area with the location of teleseismic and mag-netotelluric (MT) stations. CNSN and SOSN denote permanent stations of the Canadian Na-tional Seismograph Network and the Southern Ontario Seismic Network, respectively. Upper inset shows the location of the study area with respect to Canada. Lower inset focuses on the area of the 1994 experiment. Chapter 3. Geology and geophysical coverage of the study area 10 which the younger terranes to the south were accreted by a complex process of northward thrust-ing of the upper-crust and north-dipping subduction of the lower crust [Calvert et al., 1995]. To the south, the Abitibi subprovince is composed of alternating east-west trending volcanic (green-stone) and plutonic terranes [Ludden and Hubert, 1986]. This assemblage is believed to have formed from northward accretion (with north-dipping or perhaps northwest oblique subduction) of oceanic plateaus and/or paired arcs and sedimentary prisms [Kimura et al, 1993; Desrochers et al, 1993]. The southernmost Archean terranes traversed by the teleseismic arrays comprise the Pontiac metasedimentary subprovince, interpreted by Kimura et al. [1993] as a late stage syn-orogenic turbidite fan associated with a Late Archean prograding orogen involving the ac-creted Superior continent. To the west of the concentration of stations, the Kapuskasing Structural Zone (KSZ) transects the greenstone belts of the Abitibi and Wawa subprovinces of the Super-ior Province. The northeast trending structure represents a section of high-grade Archean lower crust, thrust southeastward during the Early Proterozoic [1.95-1.85 Ma; Percival and West, 1994]. The southern half of the study area includes the Grenville Province, which comprises the youngest terranes of the Precambrian Canadian Shield. It is composed of late Archean and Proterozoic rocks that were strongly reworked and deformed during the Grenvillian orogenic cycle, between 1160 and 970 Ma [Rivers etal, 1989]. The Grenville and Superior provinces are separated by the northeast-southwest trending Grenville Front (GF), a major crustal discontinuity marking the northern limit of Grenvillian metamorphic and deformational effects. The Grenville is divided into three main belts according to the age ofthe rocks, their origin, and the episodes of deformation they have experienced: the Parautochthonous Belt, the Allochthonous Polycyclic Belt and the Allochthonous Monocyclic Belt [Rivers et al, 1989]. The 1994 and 1996 teleseismic stations located in the Grenville Province cover mainly the Parautochthonous and Allochthonous Chapter 3. Geology and geophysical coverage of the study area 11 Belts. The southernmost stations sit on Paleozoic sedimentary rocks, which were deposited un-comformably on basement rocks of the Grenville Province. Late Jurassic to Early Cretaceous diamond-bearing kimberlites have been documented for at least two locations in the study area: the "Rapide des Quinze" kimberlites (maximum age of 126 Ma [Ji et al, 1996]) and the Kirkland Lake kimberlite field (age ranging between 147 and 158 Ma [Meyer et al, 1994, and references therein]). Kimberlite magmatism in the region has been associated with the same event responsible for the emplacement of the Cretaceous Monteregian-White Mountains-New England Seamounts (MWN) igneous series [Sykes, 1978; Crough et al, 1980; Adams and Basham, 1991 ]). The trend defined by these three igneous series will be referred to here as the M W N hotspot track. These series are composed mainly of alkaline rocks showing a general age progression from west to east, with radiometric ages dating the Monteregian Hills at 124 ± 1 Ma [Folandet al, 1986], the Cretaceous pulse of the White Mountains at 125-100 Ma [Foland andFaul, 1977], and the New England Seamounts at 103-83 Ma [Duncan, 1984]. 3.2 Geophysical coverage In the past decade, the study area has been the locus of intensive geophysical surveys, as part of the LITHOPROBE Abitibi-Grenville Transect. Seismic reflection profiles have shed some light on the tectonic evolution of the crust and the uppermost mantle [Green et al, 1990; Jackson et al, 1990; Bellefleur etal, 1995; Lacroix and Sawyer, 1995; Ludden et al, 1993; Kellett et al, 1994; Calvert et al, 1995]. Seismic refraction profiles shot over the same regions have provided well constrained velocity models [Grandjean etal, 1995; Winardhi and Mereu, 1997]. The aforemen-tioned studies reveal variations in Moho depth along the corridor traversed by the teleseismic ar-ray. A thinning of the crust occurs near the latitude of the Grenville Front, where the thickness is Chapter 3. Geology and geophysical coverage of the study area 12 32-34 km, relative to 39-43 km in the Grenville, and 34-36 km in the Superior Province [ Winardhi andMereu, 1997; Kellett et al, 1994; Grandjean et al, 1995; Martignole and Calvert, 1996]. A number of magnetotelluric (MT) stations were deployed along the Abitibi-Grenville Tran-sect, in order to obtain electrical conductivity models of the crust and the upper mantle [Kellett etal, 1992; Kellett et al, 1994; Mareschal et al, 1995]. One important result of the M T study is the discovery of deep-seated electrical anisotropy beneath the southeastern Canadian Shield. This anisotropy was first detected below the Pontiac subprovince by Kellett et al. [1992], who interpreted it as a deep crustal feature. Later studies helped constrain the nature of electrical an-isotropy, relocating it within the upper mantle, between 50 and 150 km depth [Kellett et al, 1994; Mareschal et al, 1995; Senechal et al, 1996]. The source of the electrical anisotropy is inter-preted by Mareschal et al. [1992] and Mareschal et al. [1995] as due to interconnected films of graphite deposited along grain boundaries. Mareschal et al. [1995] extended the observation of electrical anisotropy to a larger area of the Canadian Shield, including the Grenville and Superior (Pontiac, Abitibi, KSZ, Wawa) provinces. The direction of maximum conductivity is generally E-W, and parallel to the lineations defined by the major deformation zones that have affected the area during Precambrian time. Along the corridor traversed by the teleseismic arrays, the elec-trical anisotropy shows a fairly uniform orientation of N80°E [Kellett et al, 1994; Senechal et al, 1996]. These results are presented in Figure 3.2(a). Gravity and heat flow surveys [Antonuk and Mareschal, 1992; Guillou et al, 1994] com-plete the main geophysical coverage of the study area. The results from these methods help con-strain and improve the tectonic interpretations obtained from other methods, by providing inform-ation on crustal composition and thickness through modelling. Results show a general increase of Bouguer gravity anomaly and crustal heat flux from east (Grenville) to west (Kapuskasing), while the mantle heat flux stays relatively constant across the different provinces and subprovinces of Chapter 3. Geology and geophysical coverage of the study area 13 the southeastern Canadian Shield [10-14 mW m" 2 ; Guillou et al, 1994]. Finally, at the continental scale, 5-velocity models show the Canadian Shield to be charac-terized by a thick lithosphere (i.e., lithospheric root), extending to depths of 200-250 km beneath the Abitibi-Grenville Transect [Grand, 1994; van der Lee and Nolet, 1997]. Figure 3.2: Azimuths and amplitudes of (a) electrical anisotropy, and (b) shear-wave splitting for the 1994 dataset, beneath the study area. For electrical measurements, the solid bar indicates the most conductive direction, and its length is proportional to the phase difference between the most conductive and orthogonal directions. For teleseismic data, the solid bar indicates the polarization direction of the fast shear-wave, and its length is proportional to the delay time between the two split waves. Chapter 3. Geology and geophysical coverage of the study area 3.2.1 1994 teleseismic project 14 Plans for a teleseismic experiment to investigate mantle structure below the L l T H O P R O B E Abitibi-Grenville Transect were originally motivated by the discovery in the area of deep-seated upper mantle electrical anisotropy from M T measurements (see previous section). The teleseis-mic project was designed as a two-phase experiment, performed in 1994 and 1996. The first phase of data acquisition took place from July to November 1994, using 10 short-period seismographs from the French Lithoscope program. Two types of recording instruments were used, namely HADES and TITAN data acquisition systems, which operated at 25 Hz and 31.5 Hz sampling rate, respectively. A l l stations were equipped with portable short-period three-component Lennartz 3D/5s seismometers, with a natural frequency of 0.2 Hz (5 s), permitting recovery of frequencies between 0.2-12.5 Hz. Stations were deployed along a NNW-SSE line crossing the Grenville Front (Figure 3.1), with an average spacing of ~25 km. The array covered parts of the Pontiac subprovince and the northern Grenville province (Figure 3.1). After installation, the teleseismic stations were serviced once every 3 weeks, for collection of the data, synchronization of the HADES recorders internal clocks, and troubleshooting for tech-nical problems. The recovered dataset was processed into a "workable" format at Ecole Poly-technique de Montreal (EPM) and Laboratoire de Geophysique Interne et de Tectonophysique (LGIT). The following sections present the processing techniques applied to the 1994 teleseismic dataset, the results that were obtained, and their interpretation. Shear-wave splitting analysis Shear-wave splitting analysis is a technique that is widely applied to teleseismic data. It provides information on the azimuthal anisotropy of the upper mantle, which is generally interpreted as Chapter 3. Geology and geophysical coverage of the study area 15 due to the lattice preferred orientation (LPO) of olivine [Babuska and Cara, 1991; Vinnik et al, 1986; Silver and Chan, 1988, 1991]. The exact location of the anisotropy within the upper man-tle is poorly constrained. In the case of stable continental regions, like the Canadian Shield, some researchers interpret it to reside primarily within the lithosphere [Silver and Chan, 1988, 1991] whereas others suggest it is located mainly in the flowing asthenosphere [ Vinnik et al, 1992; Vin-nik et al, 1995]. Physically, the anisotropy causes a plane-polarized, vertically propagating shear-wave to split into two quasi S-waves which are orthogonally polarized along the two main axes of anisotropy, and travel at different velocities. The splitting is generally parameterized by two quantities: (i) a time shift, which depends on the thickness of the anisotropic layer and the strength of the anisotropy; (ii) the polarization direction of the fast wave, which characterizes the orient-ation of the anisotropy. Various methods have been developed to recover these parameters using inverse operators (involving rotation and time shifting) in an attempt to reproduce a single polar-ized wave from the two split waves recorded at the surface [e.g., Vinnik et al, 1986; Bowman and Ando, 1987; Silver and Chan, 1991]. These calculations generally assume the existence of a sin-gle horizontal and homogeneous layer of anisotropy. Such an assumption may be simplistic for an area like the Canadian Shield, which has been subjected to numerous episodes of deformation and reworking, but it is often necessary due to limited backazimuthal coverage in the data. Shear-wave splitting analysis was the first method employed on the Abitibi 1994 dataset. The results that were obtained and their interpretation are presented in Senechal et al [1996], Ji et al. [1996] and Rondenay [1996]. Splitting parameters were calculated independently for a series of events and phases (S, SKS, SKKS, ScS). The average values for each station are presented in Figure 3.2(b). For the entire array, the average direction of fast polarization is N101°E ± 10°, and the average splitting is 1.46 ± 0.21 s. Using this time shift and an S-velocity anisotropy of 3.2% inferred from local xenoliths samples, the thickness of the anisotropic layer is estimated Chapter 3. Geology and geophysical coverage of the study area 16 at approximately 200 km [Ji et al., 1996]. Over the study area, the direction of seismic aniso-tropy correlates well (within a consistent bias of +20°) with the calculated direction of electrical anisotropy (N80°E). Given the depth of the electrical anisotropy (50-150 km), a possible link between the two anisotropics would constrain the location of the seismic anisotropy within the thick lithosphere present below the area [Senechal et al, 1996], which is consistent with the in-terpretation of Silver and Chan [1988, 1991]. The two kinds of anisotropy are both interpreted to manifest the same horizontal stress field that originally caused the now fossilized deformation in the subcrustal lithosphere [Senechal et al, 1996]. This deformation resulted in alignment of a grain boundary conductive phase (source of electrical anisotropy) and LPO of olivine (source of seismic anisotropy). The different thickness estimates of the anisotropic layers obtained with electrical (100 km) and teleseismic (200 km) methods can be reconciled i f the stability field of the conductive phase extends to a maximum depth of approximately 150 km [e.g., graphite, Kennedy and Kennedy, 1976]. The seismic anisotropy can then extend to greater depths within the thick cratonic lithosphere present beneath the study area. While Senechal etal. [ 1996] did not consider the consistent ~20° obliquity observed between the two anisotropics (see Figure 3.2(a)-(b)), it was the basis of the interpretation presented in Ji et al. [1996]. Xenolith samples from the "Rapide des Quinze" locality present a similar obliq-uity between the LPO of olivine, which is responsible for the seismic anisotropy, and its shape preferred orientation (SPO), which can be associated with the electrical anisotropy; this obliquity is attributed to finite noncoaxial strain [Ji et al, 1996]. Since M T data and shear-wave splitting analysis provide a comparable characterization of horizontal strain-induced anisotropy within the upper mantle, Ji et al. [1996] speculate that the obliquity between seismic and electrical aniso-tropy can be employed to indicate the movement of transcurrent shear zones in the mantle. Thus, Chapter 3. Geology and geophysical coverage of the study area 17 the 20° separation between the two anisotropics beneath the Abitibi-Grenville may reflect a dex-tral shear sense in the mantle, which is in agreement with the last cycle of regional deformation inferred from the surface geology [Ji et al, 1996, and references therein]. This interpretation fur-ther supports the location of seismic anisotropy within the lithospheric mantle. Average S-wave velocities Velocity (km/s) Figure 3.3: Average 1-D S-wave velocity model obtained from receiver function analysis on the 1994 teleseismic data set. Velocity models from the individual stations have been averaged to emphasize velocity boundaries which are consistently observed across the array. Shown are the receiver function velocity model (solid line), and the average initial velocity model calculated from the refraction P-wave velocity model (dashed line). The dotted lines represent ± 1 standard deviation ofthe average values. Chapter 3. Geology and geophysical coverage of the study area 18 Receiver function analysis Receiver function analysis [Langston, 1979; Owens et al, 1984; Cassidy, 1995] has been applied to the Abitibi 1994 data set to identify 5-wave velocity discontinuities in the lithosphere beneath the array of three-component stations [Wu and White, in preparation]. The analysis comprised three steps: 1) deconvolution of the vertical component of motion from the radial and transverse components to obtain the radial and transverse receiver functions, 2) summation of individual re-ceiver functions at a given station to enhance the signal-to-noise ratio, and 3) modelling of the resultant receiver function waveforms to obtain an 5-wave velocity model. The first step elimi-nates the earthquake source and instrument response effects so that the receiver functions consist primarily of P to S converted arrivals which originate at ^-velocity contrasts beneath the array. Note that the analysis is restricted to the portion of the teleseismic signal between P and PP ar-rivals, in order to avoid contamination by other main phases and surface wave energy. The main characteristics of the resulting velocity model are shown in Figure 3.3. The follow-ing observations can be made: 1) significant velocity contrasts occurring at the top and the base of the lower crust (~25 km and ~34-43 km, respectively) correspond to layer boundaries observed in the refraction model, whereas the deeper boundary at ~85 km depth may provide new evid-ence for deep lithospheric layering in this region; 2) an average crustal value of Poisson's ratio of 0.263 is required to match the third constraint (see above), and the crustal value of Poisson's ratio increases from 0.250 beneath the Grenville Front to a maximum of 0.275 at the south end ofthe array within the Grenville Province; 3) lower crustal Poisson's ratio values of 0.28-0.30 were de-termined, which in combination with lower crustal P-wave velocities of 6.9-7.2 km/s from refrac-tion modelling [Grandjean et al, 1995; Winardhi and Mereu, 1997] are consistent with a lower crust of anorthositic and/or mafic gneiss composition [cf , Table 1-3, Holbrook et al, 1992]. The Chapter 3. Geology and geophysical coverage of the study area 19 southward increase in average crustal Poisson's ratio from the Grenville Front into the Grenville Province is likely associated with the increasing thickness of the high velocity - high Poisson's ratio lower crustal layer. \ Chapter 4 Analyses of the Abitibi 1996 teleseismic data As seen in Chapter 3, a variety of studies performed in recent years has greatly increased our knowledge of lithospheric evolution in the region encompassed by the LITHOPROBE Abitibi-Grenville Transect, an area which has remained tectonically stable for the past 1000 Ma. Among these studies, earthquake seismology is an essential tool in probing deep Earth structure. Coupled with data from xenoliths that are brought to the surface by localized and infrequent kimberlite volcanism, teleseismic signals provide important information on the nature and structure of the lithosphere and the tectonic processes that have contributed to its evolution. This chapter focuses on the second, large-scale teleseismic array which was deployed in the study area. Since the deployment of the array constituted an integral portion of this thesis project, the first section treats several technical aspects of the field experiment which was undertaken in the summer of 1996. The following sections present results obtained from various analysis tech-niques employed mainly on the 1996 data set. For some of the analyses, the data set was supple-mented with recordings from the 1994 experiment and nearby permanent stations. 20 Chapter 4. Analysis of the Abitibi 1996 teleseismic data 21 4.1 Data collection 4.1.1 Abitibi 1996 experiment The second temporary array was deployed between May and November 1996, as part of a collab-orative project involving L I T H O P R O B E , IRIS-PASSCAL, the U.S. National Science Foundation (NSF), and the Geological Survey of Canada (GSC). Table 4.1 list the people who participated in the field experiment. Prior to the actual deployment of the recording instruments, a series of organisational tasks were performed, as listed below: • Scouting of the sites (August 1995); • Training on instruments, computer interfaces and programs (March 1996); • Purchase of materials for construction of the vaults; • Establishment of suitable instrument headquarters in the city of Rouyn-Noranda (QC); • Shipping of the instruments from the US to Canada, including customs clearance. Table 4.1: List of participants for Abitibi 1996 teleseismic experiment. Principal investigators are indicated by an as-terisk. Name Institution S. Rondenay University of British Columbia D. White* Geological Survey of Canada T. Hearn*, A . Rosea, R. Rapine New Mexico State University B. Giroux, H. Wu Ecole Polytechnique de Montreal D. Lentrichia, S. Hellman IRIS - PASSCAL Chapter 4. Analysis of the Abitibi 1996 teleseismic data 22 Three deployment crews were necessary to install the first 25 recording sites. The work was performed over a period of ~2 weeks, from May 10 to 21, 1996. Three additional seismographs were later deployed by myself, for a total of 28 stations forming the Abitibi 1996 teleseismic array. A l l stations consisted of RefTek data acquisition systems and Streckeisen STS-2 portable three-component broadband seismometers, provided by the IRIS-PASSCAL instrument center, then lo-cated at the Lamont-Doherty Earth Observatory. The seismometers were mounted on concrete paving stones cemented to the ground, and encased in insulated, 1 m 3 wooden vaults. Each station recorded two continuous streams of three component (NS, EW, and vertical) seismic data. The first data stream was recorded at a sampling rate of 1 Hz and was mainly used for preprocessing and quality control. The second data stream was recorded at 20 Hz , and was used for the ana-lyses presented in the following sections. The STS-2 seismometers have a flat response to ground velocity ranging from 0.0083 to 50 Hz permitting recovery of frequencies between 0.0083-10 Hz (using a 20 Hz sampling rate). The locations of the broad-band seismic stations are presented in Figure 3.1 and listed in Ta-ble 4.2. Of the total 28 stations that were deployed for the Abitibi 1996 experiment, 26 formed a linear NNW-SSE array which overlapped with the 1994 line (Figure 3.1). A relatively narrow station spacing of ~20 km was used along this line to facilitate work on receiver function profil-ing and local shear-wave splitting analysis. The remaining two stations were deployed west ofthe main array to allow sampling of the lithosphere below the Kapuskasing Uplift and the Kirkland Lake kimberlite field. During the entire duration of the experiment, I was based in the city of Rouyn-Noranda (Quebec), which is situated close to the center of the array. From there, the stations were ser-viced once every 4 weeks to maintain the instruments and collect the data. After the servicing trips, the data were brought back to headquarters, where they were preprocessed and backed up Chapter 4. Analysis of the Abitibi 1996 teleseismic data 23 Table 4.2: Abitibi 1996 teleseismic stations. Up time refers to the percentage of the total deployment period during which stations were operational. Station Latitude Longitude Elevation Deployment period Up time name °N °E km dd/mm-dd/mm(1993) % DEND 49.797 -79.009 0.265 20/05-25/10 100 VSG4 49.627 -79.000 0.275 20/05-25/10 100 KM44 49.435 -79.999 0.268 20/05-25/09 100 N N N N 49.291 -79.157 0.350 20/05-25/10 82 VBOI 49.104 -79.144 0.310 17/05-25/10 100 V S G L 49.011 -79.109 0.290 19/05-19/10 75 C H A Z 48.832 -79.052 0.315 19/05-26/10 100 POUL 48.633 -79.043 0.277 20/05-19/10 100 DEST 48.460 -78.988 0.295 17/05-26/10 100 LDUF 48.320 -78.938 0.305 15/05-14/10 74 B E L L 48.097 -78.942 0.290 16/05-24/10 95 C A R O 47.915 -78.980 0.330 16/05-24/10 95 POP1 47.810 -79.175 0.275 27/05-27/10 95 POP2 47.593 -79.221 0.297 16/05-27/09 82 POP3 47.455 -79.180 0.275 19/05-27/10 100 POP4 47.226 -79.119 0.305 19/05-27/10 96 POP5 47.119 -78.907 0.310 19/05-14/10 91 POP6 46.949 -78.821 0.310 19/05-27/10 100 POP7 46.810 -78.601 0.305 21/05-27/10 100 POP8 46.624 -78.578 0.365 23/05-20/10 96 POP9 46.468 -78.404 0.350 22/05-24/10 100 POP 10 46.290 -78.292 0.290 22/05-27/10 94 BISC 46.027 -78.241 0.420 22/05-30/09 80 A L Q N 45.739 -78.189 0.455 22/05-29/09 57 M A D W 45.530 -78.004 0.375 22/05-25/10 100 PAPL 45.319 -77.809 0.325 22/05-25/10 85 L L A K 48.104 -79.801 0.290 13/07-28/10 100 K A P U 48.128 -82.906 0.410 12/07-26/10 93 on multiple devices. Intensive quality control was performed on the recorded signal in order to de-tect possible internal problems with the instruments in the field. The preprocessed data were then sent to the IRIS - Data Management Center (DMC) in Seattle, for final processing and archival. Chapter 4. Analysis of the Abitibi 1996 teleseismic data 24 This experiment was the first to successfully integrate the PASSCAL Database System in the data acquisition and handling process. The resulting data set consists of 1150 Gb of continuous seismic traces available to the public in SEED format from the IRIS-DMC ( 4.1.2 Permanent stations Additional data from the Canadian National Seismograph Network (CNSN) and the Southern Ontario Seismic Network (SOSN) were considered for some of the analyses presented below. The study area includes two CNSN broadband (~0.03-10 Hz) three-component stations, namely, G A C (Glen Almond, Quebec) and SADO (Sadowa, Ontario). The SOSN consists of six short-period (1 Hz) three-component stations and is operated by the University of Western Ontario (UWO) for Ontario Power Generation. 4.1.3 Event selection Figure 4.1 displays the locations of events that were selected for analysis (see Appendix A for the full list). The map projection is centered on station POP2 of the 1996 array. Travel-time inversion and receiver function analysis use teleseismic events at epicentral distances between 30° and 104°. The lower limit ensures the exclusion of transition-zone-triplicated phases, whereas the upper one defines the onset of the core shadow zone after which only core diffracted waves reach the surface. Shear-wave splitting analysis, in contrast, is commonly performed on core refracted phases (i.e., SKS, SKKS) which are effectively isolated from other waves at epicentral distances greater than 85°. Chapter 4. Analysis of the Abitibi 1996 teleseismic data 25 • S wave travel time inversion A Receiver functions • SKS splitting Figure 4.1: Data sets used for teleseismic analyses along Abitibi-Grenville Transect: (a) P-wave travel-time inversion; (b) S-wave travel-time inversion, receiver functions, and SKS splitting. See text for details. 4.2 Travel-time inversion 4.2.1 Method Inversion of teleseismic body wave travel times was performed to obtain a velocity model of the upper mantle below the study area. The method used was developed by VanDecar [1991] and can be summarized as follows: (1) optimum relative delay times are determined by cross correlation of all pairs of waveforms, for a given event; (2) delay times are inverted for velocity perturba-tions (with respect to the iasp91 one-dimensional (1 -D) radial Earth model [Kennett and Engdahl, 1991]) beneath the teleseismic array; (3) the velocity model is parameterized by splines under Chapter 4. Analysis of the Abitibi 1996 teleseismic data 26 tension constrained at a 3D grid of regular knots as shown in Figure 4.2, with the region of high resolution bounded by the area of the main location map in Figure 4.3 and extending from the surface down to 500 km depth; and (4) linear inversion is performed using conjugate gradients [Hestenes and Stiefel, 1952], simultaneously solving for slowness perturbation, station time cor-rections, and event mislocation. The inversion is an iterative damped least squares problem, with regularization by flattening and smoothing (minimizing L 2 norms of first and second model pa-rameter derivatives), and an iterative downweighting of large residuals [Bostock and VanDecar, 1995]. Further details are given by VanDecar [1991]. 1000 km Figure 4.2: Grid used for travel-time inversion. Highlighted volume represents the region of in-terest. Chapter 4. Analysis of the Abitibi 1996 teleseismic data 2 7 -84° -82° -80° -78° Figure 4.3: Map view of the high resolution region for travel-time tomography, with location of stations that were used in the analysis. See Figure 3.1 for information on tectonic units and tele-seismic stations. Inset shows the location of the study area with respect to central Canada-NE United States; solid line indicates the surface expression of the Monteregian (M)/White Moun-tains (WM)/New England Seamounts (NES) hotspot track; and dashed line represents the ex-pected hotspot track (low-velocity corridor imaged by travel-time inversion) beneath the Cana-dian Shield. Chapter 4. Analysis of the Abitibi 1996 teleseismic data 4.2.2 Data set 28 The data set for P-waves consisted of2796 travel-time picks from 123 events (5.0<mf,<6.6; Fig-ure 4.1a), while for 5-waves, 918 travel-time picks from 40 events (5.0<mb<6.6; Figure 4.1b) were used. Estimates of relative delay times were obtained from these travel-time picks, with av-erage standard deviations of 0.03 s and 0.13 s for P- and S-waves, respectively. These data did not incorporate seismograms from the Abitibi 1994 experiment due to inaccurate timing. The measured delays are mainly from direct P- and S-phases, but a few core-refracted phases were also included. Both data sets provide reasonable epicentral distance and azimuthal coverage with, however, a concentration of events located along the boundaries of the Pacific plate. As a result, there is good ray coverage along a vertical plane coinciding approximately with the teleseismic array. 4.2.3 Results The inversion results are presented in Figures 4.4-4.6 as horizontal depth slices or vertical profiles through the preferred models. The gray/color scale represents the slowness (reciprocal of veloc-ity) anomaly in terms of percent deviation from the iasp91 radial 1-D model. Regions of poor ray coverage (model parameters sampled by less than four rays) are masked in white on black and white plots (Figures 4.4 and 4.6), and in black on colour plots (Figure 4.5). In all cases shown be-low, curves representing the trade off between data misfit and model variance were constructed to select the appropriate regularization parameters. Damping parameters were set at 600 (P) and 450 (5) for flattening and 17,200 (P) and 14,100 (5) for smoothing. A series of resolution tests were performed in order to assess the reliability of the results. These tests involved the construction of synthetic models and calculation of travel-time residuals Chapter 4. Analysis of the Abitibi 1996 teleseismic data 29 synthetic P-wave S-wave •84" *T -8(T -78' -76' •"«' -84' -82" -80" -78" -W -'?4' _ S r -82' -SO' -78" -76" Figure 4.4: Resolution test for P- and 5-wave travel-time inversion. Results are presented as depth slices through the model, with a gray scale giving percent slowness anomaly with respect to iasp91 and regions of poor ray coverage masked in white, (left) Synthetic model, which is composed of alternating ± 5 % spike anomalies, (middle and right) Inversion results for P and S travel times calculated through the synthetic model. Resolution is generally good within the high ray coverage region, except close to the surface, where some of the variations may be included into station static corrections. We also note that the spikes are better recovered by P than by S travel-time inversion, which is probably due to more comprehensive P-wave coverage. Chapter 4. Analysis of the Abitibi 1996 teleseismic data 30 for the same source-receiver combinations as used in the actual inversion. Travel-time inversion was then performed, and the resulting models were compared to the original synthetic ones. Fig-ure 4.4 presents a resolution test performed for synthetic P and S slowness models composed of alternating ( ± 5 % ) single knot spikes. These spikes are positioned at regular intervals within the model volume, and random Gaussian noise (<r=0.03 s for P-waves and cr=0.13 s for S-waves) has been added to the synthetic travel times. We note that the first layer of spikes (50 km depth) is not well resolved in either P or S inversions. This is explained by the proximity of this layer to the surface; the associated travel-time delays are largely absorbed by station static corrections rather than by velocity variations in the model [ VanDecar, 1991]. For the next three layers (183, 317, and 450 km), the reconstructed image is more consistent with the synthetic input spike model. In the P inversion results, the anomalies are properly positioned within most of the volume illumi-nated by the rays, a region which near the surface is restricted to a narrow band centered on the array (e.g., ~220 km wide at a depth of 50 km) and widens substantially with depth (e.g., ~800 km at a depth of 300 km). In the S inversion tests the spikes are adequately recovered only dir-ectly beneath the teleseismic array, and thus resolution of structure is, not unexpectedly, poorer away from the array axis than for P travel-time results. Both P- and S-wave tests produce amp-litudes of recovered spikes that are lower than their original values, and this disparity becomes more pronounced with increasing depth. This observation is explained by (1) the regularization, which favors flatter, smoother solutions over "spiky" ones and (2) the fact that more rays inter-sect the model points close to the surface than at greater depths, thus placing tighter constraints on amplitudes of the near-surface model perturbations. The P-xvave inversion results for the real data are shown in Figure 4.5; the preferred model explains 93% of the root-mean-square (rms) of P-wave delay time residuals (rms reduction from Chapter 4. Analysis of the Abitibi 1996 teleseismic data 31 Figure 4.5: P-wave travel-time inversion results, presented as (a-d) a selection o f horizontal depth slices and (e) a vertical profile through the preferred model. The color scale gives the percent slowness anomaly with respect toiasp91, and regions of poor ray coverage are masked in black. A N W - S E trending low-velocity corridor ( ~ 1 % slowness anomaly; see arrows) is observed, cross-ing the seismic array near latitude 4 6 ° N . The vertical profile shows that the anomaly is nearly vertical and extends between 50 and 300 km depth, with a relatively constant width o f ~ 1 2 0 km. A t a regional scale the anomaly appears to represent the northwestward extrapolation o f the M W N hotspot track. The low-velocity corridor is flanked on either side by two high-velocity anomalies that probably represent the deep lithospheric roots of the Canadian Shield. Chapter 4. Analysis of the Abitibi 1996 teleseismic data 32 a) Depth: 50 km b) Depth: 100 km Figure 4.6: S- wave travel-time inversion results, presented as (a-d) a selection of horizontal depth slices and (e) a vertical profile through the preferred model. Description of the gray scale is as for Figure 4.4. A low-velocity anomaly is observed near latitude 4 6 ° N , extending mainly west of the station array (see arrows). In this case, the anomaly does not show a clear lateral extension as in the F-wave inversion results. This is probably due to the limited data set available for S-waves and the resulting loss in resolution. Chapter 4. Analysis of the Abitibi 1996 teleseismic data 33 13.57 to 1.00 s). The most prominent feature to appear on all depth slices is an elongate, posit-ive 1.0% slowness anomaly (low velocity) that crosses the southern portion of the main line near latitude 46°N. As depth increases and the illuminated area widens, the anomaly takes the shape of a NW-SE low-velocity corridor, with an average width of 120 km (measured at 0% slowness anomaly). Figure 4.5e shows a vertical profile (AA') that cuts the model along the teleseismic ar-ray and thus yields an approximately transverse view of the low-velocity corridor. The anomaly dips slightly to the south, which may be due in part to smearing in that direction, and its width remains fairly constant with depth, from 50 to 300 km beneath the surface. Two higher-than-average velocity zones flank the anomaly to the north and south. Synthetic tests and analysis of station corrections demonstrate that the limited outward extent of these features are artificial, re-flecting an insensitivity of the method to low wavenumber (including DC) components of struc-ture. Thus the recovered models are band-pass-filtered images of true mantle velocity structure, wherein laterally extended anomalies are well resolved only in the vicinity of more rapid, lateral variations. Accordingly, the high-velocity flanks may persist to the north and south and represent ambient, high-velocity levels of the Canadian Shield, as observed in the North American shear velocity models of Grand [ 1994] and van der Lee and Nolet [ 1997]. Consequently, the true devi-ation of the low-velocity zone from ambient mantle may be better represented by a peak to peak anomaly of ~2%. The S-wave inversion results are shown in Figure 4.6, with the preferred model explaining 97% of the rms of the delay time residuals (rms reduction from 63.97 to 1.64 s). This model also shows the presence of a low-velocity anomaly at the same latitude as that observed for P-waves (i.e., near 46°N). However, as expected, the anomaly is less well resolved away from the main axis, especially on the east side of the array. Vertically, the anomaly is better resolved closer to the surface, between 50 and 200 km (see Figure 4.6e), but the best correlation between Chapter 4. Analysis of the Abitibi 1996 teleseismic data 34 the P and S models is obtained around 300 km depth (compare Figure 4.6d and Figure 4.5d). Dir-ectly north of the low-velocity anomaly, the high-velocity zone observed in the P-wave inversion is also present but more extended, both horizontally and vertically. A potential concern regarding the origin of the low-velocity corridor arises from the pres-ence of known Moho topography near this latitude (see section 4.4). More specifically, a 5-km thickening of the crust from 40 to 45 km in the vicinity above the low-velocity anomaly is fol-lowed, to the north, by a shallowing to 35 km depth that coincides with the Grenville Front and the north flanking high-velocity region (see Figures 4.8 and 5.1). Contributions to relative travel-time residuals of ±0.15s (based on the iasp91 1-D velocity model) might be expected from this structure; however, much of this should be absorbed by station static corrections which are solved for as part of the inversion [ VanDecar, 1991 ]. Moreover, velocities above the zone of crustal thin-ning, as determined from refraction profiling [Winardhi and Mereu, 1997], are ~0.3 km/s slower than for the thicker crust to the south, resulting in a vertically integrated net travel-time anomaly near zero across the structures. The interplay of crustal thickness and velocity may be responsi-ble for the pattern of static station corrections revealed in Figure 4.8b, where negative values oc-cur at latitudes corresponding to both thinner (~47.3°N) and thicker crust (~46.3°N). However effective static station corrections may be, three independent tests were run to further investig-ate the influence of strong lateral crustal variations. First, synthetic tests show that Moho relief (without compensating for variation in velocity) does leak some signal into the shallowmost man-tle but cannot reproduce the depth nor lateral extent of the low-velocity corridor observed in the tomographic results. Second, I performed a number of "squeezing" inversions [Lerner-Lam and Jordan, 1987; Saltzer and Humphreys, 1997] using both real and synthetic data (incorporating Moho relief) where different regularizations were imposed to penalize structure at a range of depth Chapter 4. Analysis of the Abitibi 1996 teleseismic data 35 intervals within the reconstructed model. Results show that, for synthetic data, slowness anom-alies associated with crustal structure are concentrated close to the surface (<100 km), whereas the observed data require a majority of the anomaly to reside at depth (100-300 km) within the lithospheric mantle. Third, I performed an inversion from which travel-time data from the five stations overlying the mantle slowness anomaly were excluded. Crustal structure at the latitude of the anomaly will thereby have minimal influence on the reconstructed model at mantle lithos-pheric depths. The results again indicate a well-defined low-velocity corridor trending NW-SE between 100 and 300 km depth, with average width of 120 km. 4.3 Shear wave splitting analysis 4.3.1 Method As discussed in section, shear wave splitting analysis has been extensively employed over the past 15 years to study azimuthal seismic anisotropy in the continental lithosphere. The method of Silver and Chan [1991] is employed here to recover the time shift 8t and polarization direction 4> which characterize anisotropy. In their method an inverse operator is designed that minimizes the transverse energy recorded at the surface (for initial waves of known polarization, e.g., SKS) or the second eigenvalue of a covariance matrix constructed from the two split waves. The method can be extended to treat several events simultaneously at a given station [ Wolfe and Silver, 1998], increasing the robustness of the results with respect to the analysis of single waveforms. 4.3.2 Data set The data set used in this section consists exclusively of teleseismic SKS phases recorded along the Abitibi 1996 array (results from the analysis of the Abitibi 1994 data set are given by Senechal Chapter 4. Analysis of the Abitibi 1996 teleseismic data 36 et al. [1996] and Ji et al. [1996], and discussed in section Of 69 events recorded on the Abitibi 1996 array with m;>>5.5 and minimum epicentral distance of 85° (for which SKS arrives prior to S), only seven events were suitable for shear wave splitting analysis, and all of these oc-curred in the vicinity of the Sea of Japan (Figure 4.1b). Furthermore, only subsets of these seven events were considered at each individual station (see "number of events" in Table 4.3) due to varying levels of noise on the different recordings. Stations such as K A P U and L L A K were af-fected by generally low SKS signal-to-noise ratio, as evidenced by the fact that only one event could be used in either case. 4.3.3 Results Results from multi-event SKS splitting analysis of the Abitibi 1996 data set are presented in Table 4.3 and Figure 4.7. The average value of <f> over the entire array is N93°E± 18°, in agreement with the results from the Abitibi 1994 data set, and corresponds to the general trend (E-W, NE-SW) of the main tectonic features at the surface [Senechal et al, 1996; Rondenay et al, 2000b]. Over the length of the 1996 array, 4> undergoes a systematic rotation from ENE at the northern stations to ESE in the southern portion of the line, which may be attributed to slight variations in orientation of the last episodes of lithospheric-scale deformation to have affected the area [Ronde-nay et al, 2000b]. Thus anisotropy in the northern portion of the array may reflect E-W regional shear zones of Archean age [Ji et al, 1996; Rondenay et al, 2000b], whereas in the south it may result from younger regional deformation associated with the Grenvillian Orogeny. In the latter case, fast directions are almost perpendicular to Grenvillian structures (Figure 4.7), suggesting that the mantle fabric may have recorded the last extensional event of the orogeny [e.g., 1140-1090 Ma; Carr et al, 2000]. Variations in St from 0.25 to 1.05 s over the entire array (average of 0.57±0.22 s) do not correspond well with previous results (i.e., 1.46±0.21 s, see section Chapter 4. Analysis of the Abitibi 1996 teleseismic data 37 Figure 4.7: Multi-event SKS splitting results. Solid bars indicate the polarization direction <f> of the fast shear wave, and their lengths are proportional to the delay time St between the two split waves. The average direction of polarization is E-W and is believed to represent "fossilized" lithospheric strain fields associated with the last major episode of deformation to have affected the region during the Precambrian. Subtle variations in <f> and St occur at the latitude of the tomo-graphic low-velocity corridor (~46°N; thick dashed line indicates surface projection of the cor-ridor), providing evidence for complexity in the underlying structure. Chapter 4. Analysis of the Abitibi 1996 teleseismic data 38 Table 4.3: SKS splitting results, with number of events considered for analysis at each station. Site Number <t>, St, of Events deg s DEND 6 76.00± 8.50 0.90±0.25 VSG4 6 77.00±15.50 0.65±0.23 KM44 2 79.00±25.50 0.35±0.18 N N N N 5 88.00±31.00 0.25±0.20 VBOI 4 134.00±13.00 0.55±0.68 V S G L 2 53.00± 7.00 0.65±1.60 C H A Z 3 120.00±26.50 0.50±0.40 POUL 4 104.00±26.00 0.40±0.20 DEST 7 80.00±19.50 0.55±0.20 LDUF 2 94.00± 8.00 0.55±0.08 K A P U 1 49.00±37.50 1.05±0.43 B E L L 5 110.00ill .00 0.70±0.15 L L A K 1 65.00± 6.00 0.65±0.23 CARO 4 109.00±13.00. 0.75±0.25 POP1 5 109.00±28.50 0.55±0.23 POP2 2 114.00±13.50 0.45±0.28 POP3 5 115.00±29.50 0.45±0.20 POP4 3 91.00±14.50 0.35±0.08 POP5 4 89.00±16.50 0.45±0.12 POP6 5 77.00±27.00 0.40±0.15 POP7 5 86.00±22.00 0.50±0.20 POP8 5 101.00±20.00 0.75±0.28 POP9 5 112.00± 8.50 1.00±0.30 POP 10 3 83.00±21.50 0.40±0.20 BISC 5 • 110.00ill .00 0.70±0.15 A L Q N 3 100.00±19.50 0.45±0.10 M A D W 4 97.00±17.50 0.40±0.10 PAPL 2 • 77.00±16.00 0.55±0.12 The inconsistency in magnitude of average St between the two experiments may be due to the fol-lowing factors: (l)thedifferenceinbackazimuthalcoverage(285°-335° forthe 1996 experiment, Chapter 4. Analysis of the Abitibi 1996 teleseismic data 39 compared to 165°-3 3 9° for the 1994 experiment) coupled with depth-dependent and laterally vari-able anisotropy, which would result in a wider range of delay times for the analysis of the Abitibi 1994 data; and (2) the possibility of contamination by source-side anisotropy introduced through the use of some S and ScS phases in the Abitibi 1994 analysis. Subtle variations in <j> and St occur at the latitude of the low-velocity corridor (~46°N), providing possible evidence for complexity in the underlying structure (i.e., change in mantle anisotropic fabric). 4.4 P S conversions 4.4.1 Method Receiver function analysis is a method commonly used in the study of P to S (Ps) conversions gen-erated at discontinuities in the crust and uppermost mantle [Langston, 1979; Owens et al, 1984; Cassidy, 1995]. It employs the vertical component of the seismic signal as an estimate of the source and deconvolves it from the radial and transverse components in order to isolate Ps con-versions. The approach adopted here contains a series of additional processing steps [Bostock, 1998] and can be summarized as follows: (1) the rotated components of the ground displacement vector recorded at the surface [UR, UT, UZ] are transformed into an upgoing wave vector [P, SV, SH] using the inverse free-surface transfer matrix [Kennett, 1991]; (2) seismograms are binned in latitudinal intervals at several possible depths of conversion [e.g., Dueker andSheehan, 1997], in order to target known discontinuities (Moho, 410, 660), and possible discontinuities within the lithospheric mantle; (3) simultaneous least squares deconvolution of S components by P compo-nents is performed in the frequency domain for all seismograms within individual bins to obtain an impulse response; and (4) move-out corrections [Vinnik, 1977] calculated using the iasp91 model are applied during deconvolution to coherently stack arrivals from greater depths. Chapter 4. Analysis of the Abitibi 1996 teleseismic data 4.4.2 Data set 40 This approach was applied to a joint data set of 568 high-quality seismograms from 49 events, which includes recordings from both Abitibi 1994 and 1996 experiments. The distribution of events is shown in Figure 4.1b. 4.4.3 Results Figure 4.8 shows the radial impulse responses plotted in color as a function of latitude and binned according to Moho (i.e., 40 km) conversion point. A l l data have been projected onto a N-S line, and the configuration of the profile is such that two bins separate each pair of stations, with three extra bins (~ 10 km each) added to either end of the station line. Each vertical column forming the profiles in Figure 4.8 corresponds to the impulse response of an individual bin. The main feature observed on the radial component profile is the Moho, which is represented by a strong positive polarity pulse located at an average depth of 40 km. Moho topography is apparent with a tapered thickening in crust from 40 to 45 km thickness between latitudes 45°N and 46°N fol-lowed by an abrupt thinning of crust to ~35 km at 46.6°N. These observations are in good agree-ment with Moho relief estimated from seismic refraction profiling [Winardhi and Mereu, 1997] (see Figure 4.8b). The refraction Moho in the vicinity of the proposed suture is apparently con-strained primarily by refracted waves (versus wide-angle reflections), leading to a smoothing of topographical variations [ Winardhi and Mereu, 1997]. The scattered wave profile suggests a more discontinuous jump in Moho depth at ~46.6°N, some 65 km south of the Grenville Front. It is important to stress that the image obtained is likely diffuse due to diffraction at the structure. The image might be improved by employing teleseismic migration techniques [e.g., Bostock and Ron-denay, 1999a,b]. However, as illustrated in Part II of this thesis, such techniques have explicit Chapter 4. Analysis of the Abitibi 1996 teleseismic data 41 Normalized pulse amplitude Latitude North Figure 4.8: SV impulse responses mapped to depth as a function o f latitude, binned according to Moho (i.e., 40 km) conversion point: (a) uninterpreted section; (b) interpreted section with depth to Moho inferred from refraction study [Winardhi and Mereu, 1997] indicated by dashed green line, and station static corrections from travel-time inversion plotted above the section. Each vertical column forming the profiles corresponds to the impulse response o f an individual bin. The most prominent feature is the Moho (M), represented by a strong positive polarity pulse, and showing notable variations in topography in the vicinity of the Grenville Front (GF) . Note jump in Moho depth occuring at ~ 4 6 . 6 ° N , separating thicker crust to the south (40-45 km) from thinner crust to the north (~35 km). Chapter 4. Analysis of the Abitibi 1996 teleseismic data 42 geometrical requirements that were not taken into account in the design of the Abitibi 1996 ar-ray. Specifically, a denser station spacing would have been necessary to avoid spatial aliasing of diffracted signals from Moho depths. Further interpretation of the Moho jump and its relation to similar features observed on seismic profiles crossing the Grenville Front in other locations are presented in Chapter 5. A similar analysis (including profiling and stacking) was performed to detect the presence of possible variations in apparent depths to the 410- and 660-km discontinuities. Estimated Ps times from conversions occuring at the 410-km (P410s) and 660-km (P660s) discontinuities were fairly constant along the entire array at 43.5 and 67.0 s, respectively (compared to theoretical values for iasp91 of 44.1 and 68.1 s). These results are in agreement with those of Bostock [1996], suggest-ing a laterally homogeneous transition zone beneath the study area and below average Ps times related to high-velocity lithospheric mantle beneath the Canadian Shield. Chapter 5 Interpretation of the results A compilation of results presented in the previous chapter is shown in Figure 5.1 along with the depth to Moho inferred in the study of Winardhi and Mereu [1997]. Splitting parameters are shown at the position of each station above the profile. In this chapter we discuss the tectonic and geodynamic implications of these interpreted lithospheric structures; specifically, the relation of variations in crustal thickness to the Grenvillian orogeny and the nature of the low-velocity mantle corridor and mechanism responsible for its origin. 5.1 Variations in crustal thickness M y interpretation of Moho topography is shown as a black line in Figure 5.1. The geometry of Moho topography near 46.6°N (i.e., Moho jump, see section 4.4) is reminiscent of structures ob-served on a number of high-resolution reflection seismic profiles traversing other Precambrian and Phanerozoic orogens [Morgan etal, 1994; Calvert et al, 1995; McBride et al, 1995; Warner et al, 1996; Pfiffner et al, 1990; Cook et al, 1998]. In many of these studies a tapered crustal thickening followed by abrupt thinning has been ascribed directly to the suture between collid-ing plates, where the deepening Moho is associated with delaminating continental lower crust 43 Chapter 5. Interpretation of the results 44 Figure 5.1: Compilation of the results obtained in this study, with depth to Moho inferred from refraction study [Winardhi and Mereu, 1997] indicated by dashed green line. Background shows the radial impulse response profile (Figure 4.8a), with interpretation of Moho topography shown in black. Dashed black line indicates the presence of a possible subducted plate stranded within the cratonic lithosphere. A cross section of the tomographic low-velocity corridor as defined by the 0% slowness anomaly is outlined in red. Splitting parameters are shown at the position of each station above the profile. on the incoming plate sutured to thinner crust on the overriding plate. The results of this pro-cess are usually referred to as lithospheric-scale delamination structures or tectonic wedges (Fig-ure 5.2), since the lithosphere of the incoming plate is effectively split into an upper section ex-humed along a crustal ramp and a lower section subducted beneath the foreland [e.g., Cook et ai, 1998]. On the basis of these previous studies a possible interpretation of the Moho jump is Chapter 5. Interpretation of the results 45 a subduction/underthrusting margin within a tectonic wedge whose surface expression is repre-sented by the Grenville Front (i.e., the surficial limit of Grenville deformation), ~65 km to the north. It is further worth noting that a sequence of admittedly weak but suggestive negative po-larity arrivals project northward into the mantle from the Moho structure and thus may represent Ps conversions from eclogitized oceanic or continental crust stranded within the Superior cratonic lithosphere upon cessation of subduction. As discussed below, this interpretation does not nec-essarily imply that the Grenville Front is a Grenvillian age suture. In this section I explore the implications of the tectonic wedging hypothesis for the observed Moho structure by focussing on two main points, namely, (1) the relation between the structure and the present location of the Grenville Front, and (2) the possible age of the structure. crustal ramp Figure 5.2: Schematic diagram showing the geometry of a tectonic wedge (modified from Beau-mont and Quintan [1994]). The incoming plate on the left is split into an upper section exhumed along a crustal ramp and a lower section subducted beneath the foreland at point S. The compres-sive regime produces a conjugate shear zone opposing the crustal ramp in the upper section ofthe incoming plate. Chapter 5. Interpretation of the results 5.1.1 Moho step and Grenville Front 46 If subduction/underthrusting is truly responsible for the observed variations in crustal thickness and the location of the Grenville Front as an associated surface feature, a similar signature should be expected at other points along strike of the Grenvillian orogen. A number of other profiles have been run to the southwest and include a refraction line in the Sudbury (Ontario) area [ Winardhi and Mereu, 1997], and reflection lines in the Great Lakes region (Great Lakes International Mul-tidisciplinary Program on Crustal Evolution [Green et al, 1988]) and Ohio (Consortium for Con-tinental Reflection Profiling [Culotta et al, 1990]). These studies have all involved lines crossing the Grenville Front and all show some evidence for crustal thickening occurring near but on the Grenville side of the Front. Indeed, Culotta et al. [ 1990] interpret the Ohio profile in terms of sub-duction but invoke a more complicated two-sided subduction style to explain the doubly-vergent crustal reflections. The formation of tectonic wedges developed through subduction (continen-tal collision or intracontinental compression) has been extensively modeled by Beaumont and Quintan [1994] for a range of physical conditions. In their geodynamic modeling, they recon-sider both the Ohio and Great Lakes profiles in terms of continental collision models with single-sided, northwest dipping subduction. Along strike to the northeast, there are, unfortunately, no land based reflection profiles which traverse the Grenville Front, but two lines do approach the boundary from the south. Martignole and Calvert [ 1996] identify a ~40-km-wide zone of thinned crust (i.e., thickness of 36 km) northwest of Mont-Laurier (Quebec). This zone is located some 80 km southeast of the Front, but occurs on a portion of the profile oriented largely parallel to orogen strike and thus must represent short-wavelength local structure. Eaton et al. [1995], in contrast, document a 10-km crustal thickening toward the Front, near the Manicouagan reservoir. Furthermore, a narrow, strike-parallel, negative Bouguer anomaly becomes increasingly evident to the northeast and is generally interpreted to represent locally thickened crust [Rivers et al, Chapter 5. Interpretation of the results 47 1989]. The peak anomaly is located systematically to the Grenville side of the Front, suggest-ing that thickening is restricted to Grenville crust and is thus consistent with the model proposed for more southerly portions of the orogen and with a subduction margin extending the full length of a paleocraton plate boundary. The intensification in Bouguer anomaly toward the N E portion of the orogen may therefore reflect a more pronounced delamination (wedging) of Grenville crust, perhaps related to an earlier commencement of collision due to geometry of the incoming plate boundary. 5.1.2 Age of subduction The next task is to consider the likely age of the proposed tectonic wedge and, in particular, its timing with respect to the Grenvillian Orogeny. I consider three possible scenarios, namely, that the structure is post-, syn-, or pre-Grenvillian. The possibility of the structure being post-Grenvillian in age is the least likely option, given the lack of physical evidence for a major episode of tectonic activity postdating the orogeny. In-deed, a subduction/suture extending the length of the orogen could not have occurred after the main orogenic pulses without producing important post-Grenvillian reworking of the terranes lo-cated southeast of the Grenville Front. Many surface studies conducted in the Grenville Province clearly show that its terranes have retained a metamorphic/deformational imprint dating from the Mesoproterozoic, with no trace of subsequent regional-scale tectonism [Carr et al, 2000, and references therein]. It is widely recognized that the Grenville Front is not a Grenvillian age suture [Rivers et al, 1989; Rivers, 1997], or at least not one which was originally activated syntectonically with the main orogenic pulses. Numerous geological and geophysical studies have outlined the extent of Chapter 5. Interpretation of the results 48 reworked Archean rocks as far as 200 km south of the Grenville Front [Martignole and Calvert, 1996; Dickin and McNutt, 1989; Rivers, 1997, and references therein]. Furthermore, there is geochronologic and isotopic evidence for an Andean-style subduction margin along southeast-ern Laurentia, in effect over a period of ~500 Myr between 1.71 and 1.23 Ga [Rivers, 1997]. At that time, the Laurentian margin between the Great Lakes and Labrador Sea was located up to 600 km southeast of the present location of the Grenville Front [Rivers, 1997; Rivers and Corrigan, 2000]. Thus, at the onset of the Grenvillian Orogeny, a major portion of the Laurentian continent extended well into the present-day Grenville Province, implying that the tectonic wedge proposed in this study could not have been the result of a contemporary suture. However, the wedge could be related to a pre-Grenvillian suture reactivated during the last pulse of the northwest propa-gating orogeny. A related model involving step-up shear zone and associated intra-continental subduction/underthrusting has been suggested for the Grenville Front Tectonic Zone in N E On-tario [Haggart et al, 1993], using geochronological results and geodynamic models of Beaumont and Quinlan [1994]. In a more general sense, several authors have tentatively interpreted the Grenville Front Tectonic Zone as a reactivated zone of weakness that was originally associated with the buildup of Laurentia [Martignole and Calvert, 1996; White et al, 2000; Holmden and Dickin, 1995]. A similar form of intra-continental reworking due to the reactivation of long-lived (possibly accretion related) structures has been inferred from seismic profiles across the Arunta Block in central Australia [Goleby et al, 1989; Korsch et al, 1998]. On the basis of the arguments presented above, I propose the existence of an active conti-nental margin that predated the Grenvillian orogen and approximately paralleled its northwest-ernmost boundary. This margin would have been instrumental in defining the location of the Grenville Front, by providing a pre-existing lithospheric ramp for the exhumation of lower crustal elements towards the surface. The Moho step would thus represent a relict subduction margin, Chapter 5. Interpretation of the results 49 possibly reactivated as an underthrust sheet during the Grenvillian orogeny. To the southwest of the Great Lakes, the inferred subduction may have coincided with the Mesoproterozoic act-ive margin of Laurentia, since arc magmatism of that age has been identified in the Proterozoic foreland of the orogen [Rivers and Corrigan, 2000]. The absence of systematic foreland arc mag-matism over the full extent of the Greville Front remains problematic but might be explained by varying subduction geometry along the margin, as is documented on the western margin of South America [Sacks, 1983]. 5.2 Nature of the low-velocity mantle corridor A n outline of the low-velocity mantle corridor, as defined by the 0% slowness anomaly contour and determined from inversion of P-wave travel-time delays, is presented in Figure 5.1. The es-sential geometrical characteristics of this feature that must be considered in assessing its physical significance are as follows: (1) it is quasi-linear in planform and strikes NW-SE; (2) it is con-strained to the mantle at depths <300 km; and (3) it is < 120 km in cross-strike width. The trend of the structure is well resolved (see dashed line in Figure 4.3) as highly oblique to the NE-SW strike of the Grenville Orogen, thus rendering an association between the two features unlikely. Rather, the geometry of the low-velocity corridor is better correlated with an extrapolation of structures associated with more recent tectonic activity, namely, the alkaline volcanism which defines the Monteregian-White Mountains-New England Seamounts (MWN) hotspot track and the Kirkland Lake-Rapide des Quinze kimberlite fields (Figure 4.3). It is this interpretation which I shall pur-sue in more detail by first considering the process responsible for these surficial features and then the possible factors contributing to low mantle velocities. Chapter 5. Interpretation of the results 5.2.1 Proposed mechanisms for the MWN track 50 Two different processes have been proposed to explain magmatism along the hotspot track, namely, (1) a fixed mantle plume and (2) an episode of continental rifting. The plume hypoth-esis is based on five main arguments. First, the geographical configuration of the M W N Creta-ceous formations and a possible relation with the North American plate drifting over a fixed man-tle plume were originally proposed by Morgan [1972], mainly based on the geographical succes-sion of these linear chains [see also Morgan, 1971]. Second, radiometric ages show a clear age progression along the New England Seamounts (83-119 Ma, from east to west [Duncan, 1984]) and between the other groups, which are dated at 100-125 Ma for the Cretaceous White Moun-tains [Foland and Faul, 1977] and at 124 Ma for the Monteregian Hills (single pulse [Foland et al, 1986]). Similar dates from earlier studies were compiled by Crough [1981] and allowed him to trace the plume back to its present location, that of the Great Meteor hotspot (~30°N, 30° W). Third, the spatiotemporal coordinates of the Kirkland Lake (~150 Ma) and Rapide des Quinze (<126Ma) kimberlite fields are in agreement with the continental drift rate of North America proposed for middle to late Mesozoic time [Crough et al, 1980; Crough, 1981; Duncan, 1984]. Fourth, the evidence of hotspot-related uplift along the continental portion of the track was re-ported and analyzed by Crough [1981], and the swell of the marine portion was studied by Sleep [1990b]. Fifth, a similarity of initial Sr and Nd isotopic ratios between the land and sea portions of the track suggests the mantle plume origin, as noted by Foland et al. [1988]. The effect of a mantle plume on the lithosphere has been examined by several authors [Crough, 1981; Sleep, 1990a, 1996, 1997; Duncan and Richards, 1991; Davies, 1994]. Davies [1994] described a process of thermomechanical erosion which involves heating and softening of Chapter 5. Interpretation of the results 51 the lithospheric base by a plume tail, followed by mechanical removal of the material due to con-vective transfer. In the case of thick continental lithosphere, Davies [ 1994] estimated that the time of contact with a plume tail necessary to obtain partial melting is ~25 Myr. For the present study area it would have taken at least 25 Myr for the plume to produce the track of 600 km from its westernmost point (~ 150 Ma emplacement of the Kirkland Lake kimberlites) to the Monteregian Hills (124 Ma). In that case, any point at the base of the lithosphere beneath the study area would have been in contact with the plume tail for no longer than 5 Myr (for a 120-km-diameter plume tail). This time is inadequate for partial melting to occur and would explain the lack of obvious surficial expression on the Canadian Shield. The continental rifting hypothesis, on the other hand, is supported by four main arguments. First, the different Cretaceous formations of eastern North America coincide with the trend of rift zones that may have been activated or reactivated in relation to the opening of the Atlantic Ocem[Kumarapeli, \916;Sykes, 1978; McHone and Butler, 1984; Bedard, 1985; McHone, 1996; Faure et al, 1996]. Indeed, the Monteregian Hills are located near the junction of the Ottawa-Bonnechere Graben, which is believed to be of late Precambrian age but would have been reacti-vated during the Mesozoic, and the St. Lawrence Valley Paleozoic rift system [Kumarapeli, 1976; Sykes, \91%;Bedard, 1985]. The New England Seamounts are considered by some authors as be-ing related to the early Jurassic Kelvin Fracture Zone, again associated with the opening of the North Atlantic [Sykes, 1978; Bedard, 1985]. Second, even i f the White Mountains are not located above a well-defined rift structure, they still show different episodes of magmatism, with three main pulses at 230 Ma, at 200-165 Ma, and in the Cretaceous at 125-100 Ma [Foland and Faul, 1977]. The Early Jurassic middle pulse corresponds to the initiation of rifting between North America and Africa. Furthermore, there are strong geochemical similarities between formations associated with different pulses [McHone, 1996]. Hence McHone [1996] argues that a long-lived Chapter 5. Interpretation of the results 52 lithospheric tectonic zone is better suited to explain the geographical overlap of similar magmas of different ages. Moreover, Bedard [1985] argues that the doming (or uplift) observed in east-ern North America is associated with the first main pulse of magmatism of this continental rift system, rather than to the younger Cretaceous alkaline formations. Third, only the New Eng-land Seamounts show a clear age progression among individual volcanic centers [Bedard, 1985; McHone, 1996]. Fourth, a paleostress study conducted by Fame et al. [1996] outlines two im-portant events with contrasting extension direction during the Cretaceous: (1) an early regional NE-SW directed extension at ~ 140 Ma, which would have reactivated the Ottawa-Bonnechere Graben and extended farther west in the Canadian Shield to induce kimberlite magmatism; and (2) younger, localized N-S directed extension, at 125 Ma, which would have been generated as a response to the earlier event and is responsible for the single pulse of magmatism leading to the Monteregian Hills. The effect of rifting on a thin continental plate is one of extension accompanied by decom-pressional melting at the base of the lithosphere, producing magmatic episodes at the surface [White and WKenzie, 1989; Bedard, 1985]. Within the thicker lithosphere of the Canadian Shield, events of tensional stress like the one discussed by Faure et al. [1996] could produce some extension, generating occasional kimberlite volcanism without major uprise of decompres-sional melt. 5.2.2 Physical properties and low velocity I now consider the nature of the seismic anomaly observed beneath the study area by discussing the possible factors contributing to low mantle velocities, specifically: (1) temperature, (2) com-position, and (3) mineral fabric/anisotropy. It must be recognized, at the outset, that images in Figure 4.5 likely represent minimum estimates of true mantle velocity contrasts as a result ofthe Chapter 5. Interpretation of the results 53 smoothing/damping equations which are applied to regularize the inverse problem. Thus, on the basis of synthetic tests a peak to peak anomaly of 1-2% may manifest a true velocity contrast of 3-4%. Experimental temperature derivatives of seismic velocities show that P velocity anomalies of this order are produced by 400-600 K perturbations [Kumazawa and Anderson, 1969; Sobolev et al, 1996]. This may appear unrealistically large given that such elevated temperatures would inevitably intersect the solidus in some areas, implying large-scale generation of partial melt. By incorporating temperature-induced effects on anharmonicity, mineral reactions, anelasticity, and production of partial melt and their influence on seismic velocities, Sobolev et al. [1996] argue that smaller temperature perturbations of order 100-200 K may suffice to explain a 3% velocity variation. The presence of kimberlites erupted at 158-126 Ma requires that the underlying mantle must have been subject to thermal perturbations in the past; however, a major present-day thermal contribution to the observed anomaly is judged unlikely on several grounds. First, the anomaly width is of order 120 km, and for heat not to have diffused to greater extent over a period of ~ 150 Myr would require the initial thermal disturbance to have occurred along an exceedingly narrow zone (i.e., fracture). Moreover, this disturbance would have to be maintained over an unrealist-ically long time period to produce a measurable anomaly in the present day. In addition, a local heat flow anomaly would be expected over the feature, whereas thus far none has been observed [Mareschal et al, 2000]. For these reasons I dismiss the notion that a present-day major thermal anomaly is responsible for the low-velocity corridor. However, I retain for consideration the pos-sibility of contributions from minor thermal perturbations (i.e., remnant; <100 K) with related anelasticity, anharmonicity, and mineral reactions (e.g., metasomatism). A second viable mecha-nism is the current existence of small pockets of melt within the lithosphere formed through minor degrees of decompression, as contemporary stress fields appear to produce maximum horizontal Chapter 5. Interpretation of the results 54 extension with a component perpendicular to the strike of the Ottawa-Bonnechere Graben [Zo-back, 1992]. Compositional changes can also produce variations in velocity. In this case, the most impor-tant factor is usually considered to be Fe content, with fertile (Fe-rich) mantle peridotite exhib-iting reduced velocity relative to more depleted rock [e.g., Jordan, 1978, 1979]. In the context of the low-velocity anomaly observed here, Fe-rich magmas may have penetrated the lithosphere through either fracturing or thermomechanical erosion and may have solidified along a narrow fracture system and/or a thermomechanically eroded plume channel. However, the velocity vari-ations associated with Fe content (and those ascribed to different rock types, i.e., eclogite versus peridotite) are generally considered not to exceed 1% in magnitude [e.g., Jordan, 1979; Ander-son, 1990; Sobolev et al, 1996]). Thus compositional changes cannot be wholly responsible for the tomographic anomaly. Changes in anisotropic parameters represent a major source of velocity variation within the continental lithospheric mantle [Anderson, 1990] and may be in major part responsible for the low-velocity corridor observed beneath the Canadian Shield. Some evidence for this comes from the variation in splitting parameters over the vicinity of the anomaly, which presumably betrays a change in configuration of mantle fabric (Figure 4.7). However, a potential anisotropic struc-ture of stranded former oceanic lithosphere dipping to the N W near the tomographic anomaly must also be acknowledged. An origin in anisotropy might tend to favor the rifting hypothesis for the M W N track through attendant rheological considerations; however, it is also conceivable that thermomechanical erosion as envisaged by Davies [1994] could promote mineral alignment. The analysis remains somewhat inconclusive as a clearly superior candidate for the physical Chapter 5. Interpretation of the results 55 property of the anomaly fails to emerge. On the basis of simple conduction models it appears un-likely that the anomaly could represent a thermal signature directly associated with the process re-sponsible for M W N magmatism; however, minor degrees of decompressional melting, composi-tion, and anisotropy remain as potential contenders. The narrow width of the anomalous corridor would further appear to favor continental rifting over plume-induced thermomechanical erosion as a viable mechanism of origin. It is, nevertheless, difficult to ignore the arguments summarized earlier which favor the latter explanation, and some form of mantle plume interaction is clearly required to faithfully accommodate all constraints. An interesting hybrid model which involves both elements has been proposed recently by Ebinger and Sleep [1998] to explain continental magmatism and uplift throughout East Africa. Magmatism there occurs along rift zones, with gravity, geochemical data, and small degrees of extension all pointing toward plume-related vol-canism [Ebinger and Sleep, 1998, and references therein]. They show that a single, large plume can explain the timing and location of volcanism in the area, provided that plume material reaches the surface and/or preferentially causes uplift in regions of thinned lithosphere associated with pre-existing continental rifts [see also Sleep, 1996, 1997]. These results are readily adapted to the present work by observing the alignment of the low-velocity corridor with several ancient (Paleozoic) rift systems. Early rifting episodes would likely have weakened and locally thinned the base of the lithosphere beneath the Canadian Shield. Small degrees of partial melt created during subsequent passage of the North American plate over a fixed mantle plume (i.e., that re-sponsible for the Great Meteor hotspot) would naturally migrate into these zones, interacting with continental lithosphere to produce localized kimberlite in thick cratonic areas and more regular, linear igneous series in peripheral regions. Chapter 5. Interpretation of the results 5.3 Concluding remarks 56 The results presented here highlight the importance of two processes fundamental in shaping the character of the cratonic lithosphere, namely, (1) subduction associated with the assembly of con-tinental blocks, and (2) modification of the craton edifice through rifting and plume-lithosphere interaction. More specifically, the Grenville Front is interpreted to represent the surface limit of a pre-Grenvillian zone of tectonic wedging/indentation between two continental blocks formerly separated by N W dipping subduction. Gravity data and seismic profiles crossing other portions of the Grenville Front suggest that the lithospheric wedge may be ubiquitous over the entire length of the exposed orogen. In addition, I present evidence that mechanisms responsible for linear anor-ogenic volcanic chains in thinner continental and oceanic lithospheres are also likely to leave an imprint within thicker cratonic edifices. While previous studies have extrapolated hotspot tracks in areas of thick lithosphere using surface observations of rare kimberlitic intrusions, the results obtained here provide a more direct link between an extensively documented hotspot track and a well-resolved low-velocity corridor underneath the Canadian Shield. Plume-lithosphere interac-tion is implied on the basis of age progression and geochemistry, but a more precise understanding of the process will await targeted 2-D regional studies of these structures. Part II MULTI-PARAMETER 2-D INVERSION OF SCATTERED TELESEISMIC BODY-WAVES: IMPLEMENTATION AND APPLICATION TO THE CASCADIA-93 DATA SET 57 Chapter 6 Overview of Part II In the past few decades, numerous studies have demonstrated that a wealth of information on subsurface structure can be obtained from the analysis of energy in the coda of primary teleseis-mic phases [e.g., Langston, 1979; Owens et al, 1984; Revenaugh and Jordan, 1991; Dueker and Sheehan, 1997; Bostock, 1998]. Moreover, as the number and fidelity of recording instruments has increased, several approaches to multichannel processing of teleseismic waves that are remi-niscent of migration techniques used in controlled source seismology have emerged [Revenaugh, 1995; Bostock and Rondenay, 1999a, b; Ryberg and Weber, 2000; Sheehan et al, 2000]. The research presented here constitutes the third contribution in a multi-authored tripartite series which investigates the inverse scattering/migration problem in earthquake seismology and, in particular, the development of a method that can be practically implemented using currently available technology. In order to provide an appropriate context for the work which is detailed in this second half of the thesis, the next section presents a brief summary of the motivation and theoretical framework underlying the method. 6.1 Theory The theoretical framework for this work is developed in Bostock et al [2000, hereafter referred to as Paper F\, and constitutes an extension of the approach suggested by Bostock and Rondenay 58 Chapter 6. Overview of Part II 59 X V 2-D he te rogeneous ¥ medium receiver array Incident p lane w a v e distant s o u r c e Figure 6.1: Plan view of a quasi-linear array of receivers (inverted triangles) above a 2-D hetero-geneous lithosphere illuminated by an incident teleseismic wavefield. The incident wavefront is assumed to be planar at the scale of the array, with horizontal slowness pl=0?,P2), and conser-vation of the x2 slowness component p2 for all derived waves as required by Snell's law. [1999a, b; hereafter referred to as BR], Both of these studies consider the interaction of an inci-dent teleseismic wavefield with 2-D lithospheric structure and the generation of resultant scat-tered waves. The total wavefield, comprising both incident and scattered waves, is recorded on a dense, quasi-linear array of receivers situated on the Earth's surface (see Figure 6.1). The scat-tered waves, by virtue of their origin, contain important information on subsurface structure which can be recovered using techniques of inverse scattering. Two important aspects of the method outlined in Paper I serve to render it more generally Chapter 6. Overview of Part II 60 applicable to field data sets than its predecessor BR. First, the new technique accommodates in-cident wavefields arriving from arbitrary back azimuths. This is an important consideration in view of the non-uniform and sparse distribution of global seismicity that is the source of the in-cident teleseismic wavefield. Second, the technique allows for the analysis of both forward- and back-scattered waves, the latter being afforded by free-surface reflection/conversion of the inci-dent wave. The forward scattered wavefield comprises primarily transmitted P-to-S conversions which form the basis of the widely used receiver function technique [Langston, 1979; Owens et al, 1984]. Incorporation of the back-scattered wavefield allows one to further exploit a portion of the scattered wavefield which would otherwise have been treated as noise, and which places additional and complementary constraint on subsurface material properties. In Paper I and the present work, the different scattering modes are identified by index q as (see Figure 6.2): q=l, forward scattering of incident P-wave into P; q=2, forward scattering of incident P-wave into S; q~3, free-surface-reflectedP-wave into back-scattered P; q=4, free-surface-reflectedP-wave into back-scattered S; q=5, free-surface-reflected S-wave into back-scattered P; q = 6,7, free-surface-reflected S- wave into back-scattered S. There are 2 modes of back-scattered S-S energy, one for polarization within the plane containing the incident and scattered rays, and another for polariza-tion perpendicular to that plane. The inversion method is based on a high-frequency, single-scattering formulation of the for-ward scattering problem, in which the scattered displacement field for a particular scattering mode interaction q recorded at the surface and resulting from a planar incident wave can be expressed (cf. equation (26) in Paper L) as: A < ( x ' , p° x , t) = K{t) * / dx f*(x, 6)A*(x, x', p0JS[t - T«(x, x', p")], (6.1) where 2-D spatial variables x, x ' represent scatterer and receiver locations, respectively, and lie Chapter 6. Overview of Part II 61 a) b) c) vvvvvvvvvvvvvvvv vvvvvvvvvvvvvvvv vvvvvvvvvvvvvvvv P^q=5 S^q=6,7 Figure 6.2: Schematic diagrams displaying, on vertical sections, the various scattering modes considered in the inversion method: a) forward scattering of incident P-wave (q = 1, 2); b) back scattering of free-surface-reflected P-wave (q = 3, 4); and c) back scattering of free-surface-reflected S-wave (q = 5, 6, 7). within the xi,x3 plane (i.e., 2-D strike parallels a^-coordinate, with x\ orthogonal in horizon-tal plane and x3 vertical, positive downward; see Figure 6.3); and p° =(p°,p2) is the horizon-tal slowness of the incident wave. The convolutional operator K(t) is determined by the pla-nar nature of the incident wavefield and the 2-D geometry, and possesses the Fourier transform K(tx>)=u>2/The scattering potential, includes radiation pattern coefficients W^d), which are dependent on scattering angle 9, and the material property perturbations Amj=(Aa/a, A/3//3, Ap/p), where a, 8 and p are the P-velocity, S-velocity, and density of the reference medium, respectively. The full expressions for / 9 ( x , 6) are given in Appendix B. Quantity T 9 ( x , x', p° ) is the arrival time at receiver x ' of energy scat-tered from x, and »4£ (x, x', p° ) represents the combined geometrical amplitudes of the incident 3 f ( x ^ ) = ^ ^ ) A m i ( 4 (6.2) i=i Chapter 6. Overview of Part II 62 Receiver Array vvvvvvvvuvvvvvvvvyvvvvvvvvvvvvvvvv y ^ 2 \ x \ , X, p° ) T^ O ^ , x, p^) = constant Incident plane wave with horizontal slowness p° Figure 6.3: Geometrical quantities considered for scattered teleseismic waveform inversion. A l l quantities are projected onto the xi-x3 plane, and represented with solid lines when strictly con-fined to this plane, or with dashed lines when they have a non-zero component in the x2 direction. and scattered waves. That is, with q — 1 for example, x', p 0 x) = A°(x3)Ap(x3; 4; p°±)g{z'a, p° ), (6.3) where A°(x3) is the amplitude of the incident wavefield, and Ap(x3; x'3; p°) and xp(x'3, p°) are the amplitude and source polarization of the scattered P-wave Green's function, respectively. The quantity Au£(x', p° , i ) in equation (6.1) is the scattered displacement field for mode interaction q due to all 2-D line scatterers in model space that could contribute energy at the point Chapter 6. Overview of Part II 63 (x',i) in data space (see Figure 6.4a). These line scatterers are located along an isochronal surface for which, by definition, T 9(x, x', p ° ) is constant, as illustrated in Figures 6.3 and 6.4a. As originally recognized by Beylkin [1985], the scattered displacement field in the single-scattering, high-frequency approximation (i.e., equation 6.1) possesses a form which closely re-sembles the Radon transform in that it involves an integration over surfaces (or lines in 2-D). This analogy with the Radon transform (and its formal inverse) provides the basis for the development of a back-projection operator to the scattering problem. The 2-D Radon transform pair is given where the surfaces of integration in (6.4) are planes defined by their normals n and their distances from the origin s, angle ip defines the direction of n, and H{-} denotes the Hilbert transform. The definition of the inverse Radon transform in (6.5) indicates that it is possible to reconstruct a function /(x) at any point x 0 by summing F(n, s) over all planes that pass through this point. Unlike the Radon transform, the isochronal surfaces of integration in (6.1) are curved (see e.g., Figures 6.4a), but can be approximated as planes by evaluating the integral in the vicinity of a point x0. A back-projection operator for the scattering problem can be constructed by identifying Ati^x', p° , t) with F(n, s) in (6.5), with n and s corresponding to V T 9 and t, respectively [see Miller et al., 1987]. This operator, which allows the recovery of material property perturbations by (6.4) (6.5) Chapter 6. Overview of Part II 64 a) FORWARD PROBLEM Hor izonta l d i s tance [km] Hor izonta l d i s tance [km] b) INVERSE PROBLEM Hor izonta l d i s tance [km] Hor izonta l d i s tance [km] Figure 6.4: Schematic representation of forward and inverse scattering problems for q = 2 in terms of model space (i.e., scatterer location) and data space (i.e., station position versus time), a) Forward problem expresses the scattered displacement field measured at a given time t and station location x', due to all point scatterers in model space that could contribute energy to that point in data space, b) Inverse problem involves weighted diffraction stack along the travel-time curves corresponding to a point scatterer at x and recovers material property perturb-ations Ami=(Aa/a,A8 / 8, Ap/p) at that point (shown here for a single source). Chapter 6. Overview of Part II 65 A m ; , is defined through 1 r I V T 9 I 2 where |.4.9|2 = X)n «^n(x> x')-^n(x>X')> a n ^ m e variable rp is the angle between the gradient of total travel time V T 9 and the vertical axis (see Figure 6.3). The time series v£(x', p° , t = Tq) is the result of convolving the scattered displacement field A M 9 ( X ' , p° , i) with an operator whose Fourier transform is — isgn(u;)/V—i^>. This convolutional operator is a frequency filter which, in conjunction with K{UJ) in (6.1), forms the Hilbert transform required by the definition of the inverse Radon transform. Using equation (6.6) and the definition of the scattering potential / 9 ( x , 0) in (6.2), the mate-rial property perturbations at any scatter point in the model Am(x)=[Aa/a, A8/0, Ap/p]T can be obtained by solving an overdetermined linear system, leading to the expression A m = H- X g. (6.7) The elements of g are given by 1 r i i /* r WTq\" < " ( x ) = n / d l p ° I / d * ? w ' ( , ) W ~ ? • 4 S ( X ' X ' K ( x ' . p ° - ( = r ° ) | V T 9 | 2 «(«i ,7) £ ^ W T Z J - £ ^ (x > x ' ) p i . * = T*)>(6-8) o ^ n where the Jacobian (or "Beylkin determinant") has been introduced in equation (6.8) to transform the integration over geometrical quantities ip, 0 (see Figure 6.3) into one over source and receiver variables 7, x[. The elements of the Hessian matrix H are given by Hih{x) = J d |p° x | J d6Y,W?(0, |p°x| , x ) W2(6t |p» I ,x) . (6.9) Chapter 6. Overview of Part II 66 The inverse problem can, accordingly, be viewed simply as a weighted diffraction stack over all sources and receivers where the weights are determined by the analogy with the Radon trans-form. That is, for any given model point x, a scattering potential (i.e., equation 6.8) is computed by summing the energy along traveltime curves corresponding to a scatter point at that location. The material property perturbations are then retrieved through trivial matrix inversion. A schem-atic representation of the inverse problem is displayed in Figure 6.4b for a single incident wave-field. The inverse procedure is repeated for all points in model space to obtain a depth section of velocity and density perturbations. 6.2 Numerical examples In the second paper by Shragge et al. [2000, hereafter referred to as Paper IL\, the potential ofthe method outlined above was tested through a series of applications involving strictly in-plane, 2-D, synthetic data sets. The main results from that study are: 1) forward- and back-scattered modes display contrasting sensitivity to structure; 2) the null space of the inverse problem decreases sub-stantially with improved source coverage; 3) resolution of structure is greatly improved through combined inversion of different scattering modes; and 4) the method is robust to the inclusion of noise and to irregularly sampled data. 6.3 Research goal and overview Part II of this thesis focuses on practical implementation of the method and its application to data from the Cascadia 1993 (CASC93) experiment conducted in central Oregon, across the southern Cascadia subduction zone [Nabelek et al., 1993; Li, 1996]. The goal of this study is to show that inverse scattering theory/diffraction tomography can be successfully applied to teleseismic data Chapter 6. Overview of Part II 67 sets recorded in the field. Practical implementation is discussed in Chapter 7, in terms of: 1) ef-fect of oblique incidence on scattering kinematics and dynamics, 2) selection and weighting of different scattering modes, 3) choice of parameterization, and 4) computational efficiency. The application to CASC93 is treated in Chapter 8, where I first present the tectonic framework and a summary of previous geophysical work in the study area, followed by a brief overview of the experiment and description of the data set. The results of the analysis are then presented, and I conclude by discussing the implications of imaged structures for subduction zone processes in southern Cascadia. 6.4 Author's contribution As mentioned above, the research performed in Part II of this thesis has been submitted as the third installment of a three-paper series which develops, tests and applies this new technique [i.e., Paper I; Paper II; and Rondenay et al., 2000c, hereafter referred to as Paper III]. The material has also been presented in abbreviated form as a series of poster and oral communications at scientific meetings. Chapters 7 and 8 of this thesis are based in major part on Paper III. The research in all of Papers I, II, and III was conducted in close collaboration between myself, Michael Bostock, and Jeffrey Shragge. While the principal author on each of these three papers was responsible for the main contribution, co-authors were actively involved in most aspects of the research. In particular, I assisted in the development of the theory, first in the early work of BR, and then in the more comprehensive approach derived in Paper 1.1 wrote the computer code which implemented the algorithm for oblique incidence, and which accommodates multiple sources and scattering modes. Finally, I conducted most of the processing and interpretation ofthe Cascadia 1993 data set. Chapter 7 Considerations for practical implementation The inversion methodology outlined in the previous chapter requires the evaluation at each image point x of several potential functions <fr(x) which are defined by equation (6.8). These potential functions are simply weighted diffraction stacks of the data v*(x', p° , t) along moveout curves T 9 ( x , x', p ° ) corresponding to scattering from 2-D (line) perturbations in material properties. The weights include, among other factors, the product of incident and scattered wave geometrical amplitudes ^ ( x , x'). In this chapter, I discuss a variety of practical considerations pertaining to the evaluation of these potential functions, namely: (1) the effect of obliquity on T 9 ( x , x', p ° ) and Aq(x, x'); (2) factors relevant to computational efficiency; (3) the choice of inversion para-meters, Jacobian, and scattering modes; and (4) preprocessing of the raw data to obtain quantity <(x ' ,p°,r.). 7.1 2-D travel times and scattered wave amplitudes for oblique incidence 7.1.1 Travel times For convenience in processing, all travel times are normalized to the arrival of the incident P-wave at each station. Expressions for the travel-time curves in the in-plane case are given in Paper II. The generalizations of these expressions to incident wavefields at arbitrary back azimuths are 68 Chapter 7. Considerations for practical implementation 69 given by: r ^ 2 ( x , x ' , P ° ± ) = f3 D Y S -P°1(X'1-X1) r^1-**™)', (7-1) T- 3 ' 4 (x , x', pi) = P d y 3 - fto(a5i _ ^  + J o "Wy/l - k ( y s ) ] 2 ( r f ) a T- 5. 6. 7(X,X', p o ±) = r _ — d*> _^+ J o *(Vz)y/l ~ k(»a)] 2(pf )2 where, for 9 = 1 , 3 , 5 , (7.4) for g = 2,4,6,7. (7.5) The quantity 1^(2/3) can be interpreted as an obliquity-corrected velocity at depth y 3 in the refer-ence medium, and p\ is the component of the scattered wave slowness in the xx direction (p? = pf for q = 1,3,5; and j>\ = pf for q = 2,4,6, 7). Since p 2 is constant across all scattered waves for any individual, incident wavefield, there is no travel-time dependence on x2 in (7.1)-(7.3), and so the calculations are performed entirely with respect to spatial coordinates xx and x3. vq(yz) (32(y*) - (p 2) 2 Chapter 7. Considerations for practical implementation 7.1.2 Scattered wave amplitude 70 The effect of obliquity on the product of amplitudes, A^(x, x'), is contained entirely within the amplitude of the Green's function since the incident wavefield is planar with an amplitude that depends only on the material properties (e.g., equation (21) of Paper I, cf. also equation (6.3)). Thus, the amplitude of the 2-D, P-wave Green's function corresponding to a line source with axial component of forcing is expressed as (see Appendix 1 of Paper 7): AP(x3; x'3; p2) = (7.6) M*s) \\ 7rp(x3)a(x3)p(x'3)[Jp(x3; x>3)}2y/l - PW(x'3) ' where p is the density, and the geometrical spreading for a horizontally layered medium is given by: P cos(j>(x3)cos(j)(x'3) ( dx'x \ 3 {x*> *3) = \\ IdpJ) • (7-7) The quantity <j> is the angle between the ray and the vertical projected into the xx-x3 plane (see Figure 6.3). Obliquity is evident in equation (7.6) through the factor <^j\ — p\a2(x'3) in the de-nominator, and in equation (7.7) through the obliquity-corrected velocity vq(x3). Corresponding expressions for the scattered -S-wave are obtained by replacing a by 0 and superscript P by S in equations (7.6) and (7.7). 7.2 Interpolation In practice, the quantities presented in the previous section (i.e., T 9 ( x , x ' , p ° ) , vq(x3), Jp(x3; x's), pi) will need to be computed at every image point (xi,x3) for multiple sources (as Chapter 7. Considerations for practical implementation 71 identified by the horizontal slowness \p°L | of the incident wave) and multiple receivers (as iden-tified through surface coordinate x'-J. To minimize the computational burden, a set of tables con-taining travel time and geometrical spreading parameters for a regularly sampled range of (xi-x[), xz, p° and p2, are computed in advance using the 1-D reference velocity model. During inversion the required quantities are recovered through table lookup and interpolation, thereby reducing computing time substantially. 7.3 Inversion Several factors may influence the recovery of material property perturbations in the inversion. These factors include the effect of parameter choice on the conditioning of the problem, the se-lection and relative weighting of mode interactions q included in the inversion, and the choice of Jacobian for an irregular, sparse distribution of seismicity. 7.3.1 Choice of inversion parameters As the forward problem in Paper I is developed for purely isotropic heterogeneity, the material property perturbations can be expressed in terms of three independent parameters. Examples of physically meaningful parameter combinations include the Lame constants and density (AA/A, Afi/p, Ap/p), elastic wave velocities and density ( A a / a , A8/8, Ap/p), or impedances and density (AIP/IP, AIs/Is, &p/p)- However, the optimal choice, as dictated by the condition number of the inverse problem, may not correspond to these simple forms and must be deter-mined through eigenvector analysis for general experimental configurations [Forgues and Lam-bare, 1997]. For inversion of the Oregon data set, we shall employ a reduced parameter set com-prising P- and S-velocity perturbations. Density perturbations are not considered as they are, in Chapter 7. Considerations for practical implementation 72 general, poorly resolved and are effectively indeterminate for realistic array geometries and noise conditions. We further note that although appropriate for forward scattering, the choice of veloc-ities is less ideal for back-scattered modes, where impedances are favoured. Velocities are there-fore adopted here largely as a matter of convenience. 7.3.2 Choice of scattering modes An advantage of our approach is the ability to simultaneously consider multiple scattering modes which leads to improved resolution of structure (see Paper IL). A natural relative weighting for different mode interactions enters equation (6.8) through the factor | V T 9 | 2 which measures the sensitivity of a particular mode's travel time to structural location. As a result it can be shown that forward scattering interactions are less sensitive to structure than back-scattering interactions, a situation which is particularly acute for the P-P (<?=1) mode. However, it should also be noted that the derivation of equation (6.8) assumes a complete separation of scattering modes which can be only partially achieved in practice through the approximate orthogonality of P- and S-wave po-larizations. Consequently, different mode interactions cross-contaminate one another to different degrees with relative amplitudes that are controlled in part by the strength of free-surface reflec-tion/conversion coefficients for a particular incident wave. It is desirable, therefore, to further weight contributions from individual scattering modes on the basis of this latter consideration, and we accordingly adopted an empirical scheme which upweights g=2,3,4,6,7 relative to q=l,5. 7.3.3 Choice of Jacobian Another issue requiring consideration involves the choice of Jacobian in equation (6.8), which is employed in transforming the integration over geometrical quantities (ip, 6) to one involving the Chapter 7. Considerations for practical implementation 73 experimental variables (i.e., sources and receivers). This Jacobian can take several forms since there are at least 3 experimental variables. The source parameter in the Jacobian can be repres-ented through either event back azimuth (7), leading to the expression (see equation (72) in Paper I): (7.8) or the magnitude of horizontal slowness (|p° |), which involves exchanging the roles of |p° | and 7 in equation (6.8), such that we employ (see equation (73) in Paper 7): a(3cos<f>p2(p° + pf) (Js)2sm6 a8cos<j>pf 0 ( * i , | P l l ) (Jsy sine (p° + P f ) Ipil + M 2 , ribll Ipil + Ps + p f t f l p i l pfp°s (7.9) In special cases where structural recovery is limited to a single material parameter, angle 6 no longer needs to be considered in the transformation of variables, and the Jacobian reduces to a simple derivative (see equation (74) in Paper I): COS<f) ( J W + pf) ' (pj + P°z)(pjp°r - P°sPf) , 5N | V T 5 | 2 +P1 (7.10) The preferred choice of Jacobian will thus be governed by the distribution of seismicity and the selection of inversion parameters. For the Oregon data set, we employ equation (7.9) for P-P scattering, to avoid downweighting sources which are close to in-plane (p2 ~ 0) and yield a null Jacobian in (7.8). In the case of PS and SS scattering, we employ the single parameter Jacobian because S-velocity perturbation (88/8) is the sole inversion parameter for these scattering modes. Chapter 7. Considerations for practical implementation 7.4 Preprocessing of the data 74 Prior to inversion, the raw data must be preprocessed to effectively remove the incident wavefield, thereby isolating the scattered displacement field (i.e., A«;(x' , p° , t) in equation (6.1)). This pre-processing involves several steps which are presented in detail by BR, and can be summarized as follows: (1) Transform displacement u = [uR,uT,uz]T to a wave- vector space w = [P, SCI, SC2]T using the inverse free-surface transfer matrix [Kennett, 1991] to isolate the incident P-wave: w = F x u, (7.11) where, and, 4 a / 3 2 | p ° | < M / 3 / S 2 / % ( l - 2 | p ° | 2 / 3 2 ) / 5 0 0 0 2 -2ag a (l - 2 \p°L\ 2 02)/6 A03 |p° | 0 = \fot~ •> 0 I2 2 - I P i l . 9/3 = y j 0 -2 - \ p l \ 2 , 6=1 - 4 P i 02 + A p° 84 + A0 (7.12) (7.13) (2) Time-normalize the resulting P-wave section (i.e., align sections with respect to incident P-wave arrival) using multi-channel cross-correlation derived delay times [ VanDecar and Crosson, 1990]. (3) Estimate and separate incident (i) and scattered (SC3) wavefields by principal component (i.e., Chapter 7. Considerations for practical implementation 75 eigenimage) analysis of the P-wave section [Ulrych et al, 1998]. lis effectively obtained by re-constructing the section using only the first two principal components, which represent the com-ponents most strongly correlated from trace to trace, whereas SC3 is assembled using the remain-ing principal components, as an estimate of the scattered field on the section. (4) Reconstitute scattered displacement A u R T Z = [AUR, AUT, AUZ]T from w' = [SC3, SCI, SC2]T using the free-surface transfer matrix: A u R T Z = Fw'. (7.14) (5) Deconvolve source-time function / from reconstituted displacement sections A U R T Z to obtain the scattered impulse response. (6) Rotate horizontal displacement into a reference frame aligned with the inferred strike of the study area (i.e., x2 parallel to 2-D regional strike). (7) Apply the filter F(ui) = — i s g n ( u > ) / t o the scattered displacement field. This final op-eration yields the quantity v*(x', p° , t) which appears in equation (6.8). Chapter 8 Application to CASC93 To date, few passive seismic experiments have employed array configurations suitable for mul-tichannel processing of scattered waves as envisaged here. However, during 1993-1994, re-searchers at Oregon State University undertook the IRIS-PASSCAL Cascadia 1993 (CASC93) experiment with the objective of imaging detailed structure of the Cascadia subduction zone across central Oregon [Nabelek et al, 1993]. This deployment was unprecedented in both num-ber and density of receivers, and thus presents an ideal opportunity to gauge the performance of the inversion algorithm on an actual field data set. Analysis of the CASC93 data set is preceded by an overview of tectonics in the Cascadia subduction zone and previous geophysical work. 8.1 Tectonic framework The Cascadia subduction zone extends along the west coast of North America, from southern Brit-ish Columbia to northern California (Figure 8.1), where the Juan de Fuca plate subducts beneath the North American plate at a rate of ~42 mm/year in a direction approximately N69°E [DeMets et al, 1990]. Along the trench, the plate is relatively young (i.e., 4-10 Myr) and warm [Heaton and Hartzell, 1987; Hyndman and Wang, 1993]. 76 Chapter 8. Application to CASC93 77 Figure 8.1: Simplified tectonic map of the Cascadia subduction zone in western Oregon and Washington [modified from Mooney and Weaver, 1989], with location of the CASC93 teleseis-mic stations indicated by white squares. Black triangles indicate the position of major Cascades strato-volcanoes. B M denotes the Blue Mountains. Inset shows the location of the study area with respect to North America. Arrows point to the location of (1) conductivity model of Wan-namaker et al. [ 1989], (2) velocity model of Trehu et al. [ 1994], and (3) receiver function model of Li [1996], which are presented in Figure 8.2. The region sampled by the CASC93 experiment represents a record of past and present accre-tionary events along the western North American Cordillera, and extends from the coast of cent-ral Oregon well into the back-arc. The surface geology includes, from west to east: 1) the Coast Chapter 8. Application to CASC93 78 Range, which is composed of seamounts and oceanic basalt terranes accreted to North America during Pal eocene to middle Eocene and subsequently covered by marine forearc sediments; these structures are collectively defined as the Siletz terrane by Trehu et al. [1994], and were up-lifted in the west during late Eocene to Oligocene in response to continued underthrusting by the Juan de Fuca plate [Christiansen et al, 1992]; 2) the Willamette Valley, a structural trough within the Siletz terrane, covered with thick sequences of middle Tertiary to Quaternary sediments [Trehu et al, 1994]; 3) the Cascade Range, which is a result of continental arc volcanism related to the sub-duction of the Juan de Fuca plate in the Tertiary (Western Cascades) and Quaternary (High Cas-cades); 4) the Columbia Intermontane, or back-arc, region comprising various structural basins filled with Miocene to Quaternary volcanic and sedimentary strata. These basins form high el-evation plateaus (e.g., Columbia Basin, High Lava Plains), separated in central Oregon by Late Triassic to Early Cretaceous accreted terranes forming the Blue Mountains [Mooney and Weaver, 1989; Christiansen etal, 1992]. 8.2 Previous geophysical studies A variety of geophysical studies have been undertaken to image the Cascadia subduction zone beneath central Oregon. These include a major magnetotelluric survey [EMSLAB Group, 1988; Wannamaker et al, 1989], a series of controlled and passive source seismic studies [e.g., Lang-ston, 1981; Leaver etal, 1984; Michaelson and Weaver, 1986; Rasmussen and Humphreys, 1988; Keachetal, 1989; Harris et al, 1991; VanDecar, 1991; Nabelek et al, 1993; Trehu et al, 1994; Li, 1996; Flueh et al, 1998], and potential field modelling [Couch and Riddihough, 1989; Flem-ing and Trehu, 1999]. Among the first geophysical evidence for the subducting Juan de Fuca plate beneath central Oregon was provided by Langston [1981], using receiver function analysis at WWSSN station COR (Corvallis, OR). This study identified the Moho of the subducting oceanic Chapter 8. Application to CASC93 79 lithosphere as the base of a low velocity zone 45-50 km beneath the station, dipping eastward at 20°. Wannamaker et al. [1989] modelled magnetotelluric data from the E M S L A B experiment to image the subducting crust and/or overlying sediments of the Juan de Fuca plate as a ~5 km thick conductive layer dipping at ~20° and extending to approximately 60 km depth and 90 km east of the coast (Figure 8.2a). They interpreted the high conductivity as due to the presence of fluids within the subducting oceanic crust. Further inland, they identified a second, more conductive, horizon above the dipping slab, starting in the crust/upper mantle beneath the Willamette Val-ley and extending eastward. This conductive body was interpreted as a region of trapped fluids expelled from the oceanic crust as a result of dehydration/devolatilization accompanying high pressure and temperature metamorphic reactions. Controlled-source, seismic surveys also image the subducting Juan de Fuca crust beneath the forearc, both offshore [Flueh et al., 1998] and to greater depth on land using refracted phases and wide-angle reflections [Trehu et al., 1994]. In the latter study, the authors identified the oceanic crust as a 5-8 km thick, low-velocity layer dipping at 13-16° beneath the Coast Range (Figure 8.2b). Moreover, they model a progressive thickening of the upper crust with increasing depth along the subduction zone. This thickening is associated with a series of strong reflectors at 30-40 km depth [originally identified by Reach et al, 1989], which were tentatively interpreted as the incorporation of a block of oceanic crust within the accretionary prism [Trehu et al, 1994]. The general results in Trehu et al. [1994] were used to constrain models of magnetic and gravity anomalies for the offshore portion of the subduction zone in central Oregon [Fleming and Trehu, 1999]. These models indicate that the sedimentary column entering the subduction zone has a Chapter 8. Application to CASC93 8 0 p(Q-m) <10 <30 <1 00 <300 <103 £ 60 Q . Q 80 -123.5 -123 -122.5 -122 -121.5 -121 -120.5 Minimum thickness of Sttetz Maximum ttnefcfess 7.0 vertical exaggeration 2:1 -124 -123.5 -123 -122.5 -122 -121.5 -121 -120.5 c ) ra=6.5 km/s 20L5=0.33 JdF crust E 40-60-D Q 80-100-120-oc=8.1 km/s CT=0.25 v v v v v v w v v a=6.2 km/s o=0.27 ' V V vv ,v vv Continental Moho no vertical exaggeration -124 -123.5 -123 -122.5 -122 -121.5 -121 -120.5 Longitude [°] Figure 8.2: Geophysical studies of central Oregon (see Figure 8.1 for location of profiles), a) Deep-seated features of the E M S L A B resistivity model [modified from Wannamaker et al, 1989]. b) P-velocity model obtained from seismic refraction profile by Trehu et al, 1994, with RZ indicating the location of crustal reflective zones, c) Integrated receiver function model of Li [1996], with average P-wave velocities a and Poisson's ratio a indicated for individual model-ling areas (shaded). The outline of the receiver function model is superimposed on the first two profiles for comparison. JdF denotes the Juan de Fuca plate. Chapter 8. Application to CASC93 81 maximum thickness of 2 km, in close agreement with the estimates of Flueh et al. [1998]. Fur-thermore, demagnetization of subducting oceanic crust is inferred to occur beneath the continen-tal shelf and result, in part, from increased hydrothermal alteration and metamorphic reactions [Fleming and Trehu, 1999]. Using receiver function analysis and constraints from the crustal velocity model of Trehu et al. [1994], Li [1996] provided an interpretation of the CASC93 data set and modern digital data from station COR. He identified direct and multiple P-to-S conversions from the oceanic crust (modelled as a ~7 km wide, low-velocity layer) which were interpreted to indicate an increase in dip of the subducted plate from 10° beneath the Coast Range to 23° beneath the Willamette Valley and the Western Cascades (Figure 8.2c). The model also includes the continental Moho starting ~150 km from the coast and extending eastward at an average depth of 35 km which is well constrained beneath the backarc region. This depth value is considerably smaller than the crustal thickness of 44 km modelled from refraction/wide-angle reflection data along the north-south axis of the Oregon Cascades [Leaver et al, 1984], indicating a possible eastward thinning of the continental crust between the arc and backarc region. Poisson's ratio estimates vary laterally through the model, with values of 0.33 beneath the Coast Range, 0.25 beneath the Willamette Valley and the Cascades, and 0.27 in the backarc [Li, 1996]. Several P-wave, travel-time tomography studies have been performed in the Pacific Northw-est to investigate the deeper structure of the Cascadia subduction zone [Michaelson and Weaver, 1986; Rasmussen and Humphreys, 1988; Harris etal, 1991; VanDecar, 1991]. These studies im-age the subducted Juan de Fuca plate in the upper mantle beneath the High Cascades, as a steeply dipping tabular feature ~100 km wide and characterized by velocities that are 2-4% faster than the surrounding mantle. The anomaly persists to ~300 km depth beneath central Washington and northern Oregon, but only to 150-200 km depth beneath central and southern Oregon. The loss in Chapter 8. Application to CASC93 82 signature below these depths is interpreted as either a change in thermal/compositional properties that reduces the velocity contrast between the slab and the surrounding mantle [Rasmussen and Humphreys, 1988] or a tear in the slab related to changes in subduction rate over the past ~ 1 OMa [VanDecar, 1991]. Finally, it is worth noting that the Cascadia subduction zone in central and southern Oregon is characterized by abnormally low to absent continental forearc and Wadati-Benioff seismicity relative to more northerly latitudes [Weaver andMichaelson, 1985; Fleming and Trehu, 1999]. This difference in seismicity may be related to the northward migration and clockwize rotation of the accretionary complex with respect to pre-Tertiary North America [Wells et al, 1998], and north-south variations in subduction geometry [Fleming and Trehu, 1999]. 8.3 Experiment and data set The CASC93 data set was acquired using 44 portable RefTek data loggers and Guralp-ESP/STS-2 three-component broadband seismometers, deployed in several phases [Li, 1996]. The instru-ments recorded two continuous data streams at 1 Hz and 20 Hz, with the latter resulting in a flat response to ground velocity between 0.0083-10 Hz for STS-2 and 0.03-10 Hz for Guralp-ESP seismometers. Stations were relocated several times during the one-year period, for a total of 69 sites and an average station spacing of ~5 km [Li, 1996]. Between February and April 1994, fewer than 15 stations were operational resulting in limited data coverage for the easternmost por-tion of the array. The selection of data for multichannel inversion was performed qualitatively, with the re-quirement that traces show high signal-to-noise ratio for incident P-waves and converted S- and P-phases following in the P-coda. Figure 8.3 illustrates the final event distribution on a polar Chapter 8. Application to CASC93 83 Figure 8.3: Distribution of events selected for multichannel inversion of the CASC93 data set. Arrows indicate the location of events #15 and 16 (see Table 8.1), which are retained for detailed analysis. projection centered on the array. A l l events (see Table 8.1) fall between 30° and 104° epicentral distance to avoid transition zone triplications and core diffracted P. This selection resulted in an average of 32 station recordings per event, after removal of individual, low-quality traces. Back-azimuthal coverage is governed by the global distribution of plate boundaries, and events evenly span most of the second, third and fourth quadrants (as measured clockwise from north), although there is no representation between 0° and ~120°. Chapter 8. Application to CASC93 84 Table 8.1: Events selected for scattered waveform inversion. NS denotes the number of station recordings per event, and X's mark scattering modes (<j=2,3,4,6,7) retained for simultaneous in-version. # Date Time Lat °N Lon °E Depth km mt, NS q=2 9=3 q=4 9=6,7 1 93:05:15 21:52:25 51.37 -178.67 32 6.2 26 X X 2 93:05:16 21:44:48 -15.29 -173.33 21 6.1 31 X X 3 93:05:17 16:02:53 -5.34 151.99 17 5.7 30 X X 4 93:05:18 10:19:33 19.91 122.45 169 6.4 30 X X 5 93:05:24 23:51:28 -22.67 -66.54 221 6.6 38 X X X 6 93:06:06 13:23:20 15.82 146.60 14 6.0 36 X X 7 93:06:08 13:03:36 51.22 157.83 71 6.4 33 X 8 93:06:08 23:17:41 -31.56 -69.23 113 6.5 37 X X X 9 93:06:18 11:52:51 -29.05 -176.75 16 6.2 38 X X 10 93:06:18 17:57:46 -28.68 -176.89 11 5.9 38 X X 11 93:08:07 17:53:24 -23.87 179.85 523 6.0 38 X X 12 93:08:07 19:42:41 41.99 139.84 14 6.2 38 X X 13 93:08:09 12:42:48 36.38 70.87 215 6.2 38 X X X 14 93:08:11 14:17:37 13.18 145.65 22 5.9 36 X X X 15 93:09:03 12:35:00 14.52 -92.71 27 5.8 41 X X X X 16 93:09:06 03:56:00 -4.64 153.23 49 6.2 39 X X 17 93:09:10 17:28:08 14.43 -92.81 62 5.3 38 X X X X 18 93:09:10 19:12:54 14.72 -92.64 34 6.2 40 X X X X 19 93:09:19 14:10:56 14.36 -93.32 18 5.7 43 X X X X 20 93:09:26 03:31:14 10.00 138.22 10 6.1 40 X X 21 93:09:30 18:27:50 15.42 -94.70 19 5.8 44 X X X X 22 93:10:11 15:54:21 32.02 137.83 351 6.4 41 X X 23 93:10:24 07:52:15 16.75 -98.72 21 6.3 35 X X X X 24 93:11:13 01:18:04 51.93 158.65 34 6.5 40 X X 25 94:01:10 15:53:50 -13.34 -69.45 596 6.4 4 X X 26 94:02:11 21:17:31 -18.77 169.17 206 6.4 14 X X 27 94:02:12 04:16:26 -10.79 -128.80 15 6.3 15 X X X 28 94:02:12 17:58:24 -20.56 169.33 28 6.3 15 X X 29 94:03:09 23:28:06 -18.04 -178.41 563 6.6 13 X X 30 94:03:14 20:51:24 15.99 -92.43 164 5.8 13 X X X 31 94:03:15 03:36:19 11.11 -88.08 15 5.8 13 X X X Chapter 8. Application to CASC93 85 0 10 20 30 40 Time [seconds] Figure 8.4: Preprocessing of three-component recording at station A03 (44.42°N,-123.95°E) for event #15 (Figure 8.3). Seismic traces represent, from top to bottom: radial, transverse and ver-tical components of the raw data (UR, U T : UZ); decomposed P-, SCI- and SC2- components; es-timates of the incident (I) and scattered (SC3) waveflelds determined from principal component analysis; scattered displacement field in the x\, x2 and x3 directions, bandpass filtered between 0.03-0.3 Hz (Aui, Auz,Au3). Chapter 8. Application to CASC93 86 The choice of strike direction (i.e., x2) for the 2-D inversion is facilitated by the observa-tion that the main surface features (i.e., deformation front, Willamette Valley, High Cascades) all run approximately north-south. Thus, xx and x2 are set to coincide with due east and due south, respectively. Figure 8.4 illustrates the main preprocessing steps, as outlined in sec-tion 7.4, on a three-component recording of a Mexican event (#15 in Table 8.1) at station A03 (44.42°N,-123.95°E). The radial, transverse and vertical components of displacement at the sur-face (UR, UT, UZ) are transformed to a wave-vector space which effectively isolates the incident P-wave to a single trace (P), with remaining two components (SC1,SC2) containing only scat-tered energy. P-wave traces from all 41 stations recording the event are then used to estimate the incident (7) and scattered (SC3) wavefields on each trace through principal component ana-lysis. The scattered displacement at the surface is reconstituted from SCI, SC2 and SC3 using the free-surface transfer matrix. The scattered impulse response (Aut, Au2, A%) is obtained using the incident wavefield (/) as a source estimate to deconvolve the reconstituted displacement. A bandpass filter (0.03-0.3 Hz) is applied to the data to optimize signal-to-noise ratio. Figure 8.5 displays a series of three-component preprocessed data sections for events #15 and 16 (New Ireland Region) in Table 8.1, with incident waves illuminating subsurface structure obliquely from the right and left, respectively. Data sections such as these represent the input to the inverse scattering algorithm. Examples of travel-time curves T 9 ( x , x', p° ) for a model point located at -122.55° longitude and 55 km depth are superposed upon the sections. The two events considered here were selected to illustrate the dependence of scattered response on direction of source incidence (see Paper II for discussion). In the case of the Mexico earthquake, the most striking features are a series of dipping conversions observed over the western portion of the ar-ray. The first of these is a forward-scattered wave (q=2) from the subducted plate, which arrives Chapter 8. Application to CASC93 87 Event #16 a) A « ! component -124 -123.5 -123 -122.5 -122 -121.5 -121 -120.5 b) Au2 component Longitude [°] Event #15 d) AMT component -124 -123.5 -123 -122.5 -122 -121.5 -121 -120.5 -124 -123.5 -123 -122.5 -122 -121.5 -121 -120.5 f) Aw3 component -124 -123.5 -123 -122.5 -122 -121.5 -121 -120.5 Longitude [°] Figure 8.5: Three-component preprocessed data sections for events #15 and #16 (Figure 8.3). Diffraction hyperbolas (T 9 (x , x', p°±) in equations 7.1-7.3) for a model point located at -122.55° longitude and 55 km depth are plotted on each section and identified as (from top to bottom) solid lines for q=\,3,5 and dashed lines for 9=2,4,6-7. Arrows point to main conversions from the Juan de Fuca (JdF) plate and continental Moho (M). between 5 and 10 s after direct P. It is followed by a lower frequency back-scattered (q=4) ar-rival between 10 and 30 s. In the eastern portion of the array, between -121.2° and -122.3° (i.e., Chapter 8. Application to CASC93 88 150-275 km east of the coast), a series of laterally coherent arrivals at ~4-6 s represent a forward scattered wave (q=2) from the continental Moho. For the New Ireland event, the direct conversion (q=2) is largely absent in the western portion of the array, while the reverberation (<?=4) is again observed between 10-30 s, as is the near-horizontal Moho signal in the east. The conversions in this case are strongest on the ux profile, as the incident wavefield arrives nearly in-plane, thus u-i corresponds approximately to the radial component. 8.4 Reference model The inversion approach is based on the single-scattering or Born approximation which is valid for small-contrast material property perturbations defined with respect to a 1-D reference model [Paper T]. It is important therefore to select a broadly representative reference model i f structural features are to be accurately imaged. I constructed a 1 -D S-velocity model adapted from both the 1-D (COR) and 2-D (CASC93) modelling results of Li [1996], and which is presented in Ta-ble 8.2. To compute P-velocities, I employ a constant Poisson's ratio of 0.33 for the crust, consis-tent with the value inferred by Langston [ 1979] and Li [ 1996] in the western portion of the model, and a value of 0.3 for the mantle [Holbrook et al, 1992]. Given the selected velocity model and the bandwidth of the preprocessed signal (0.03-0.3 Hz), the minimum scale of discontinuities re-solved by the technique is of the order of 3 km for 5-waves and 5.5 km for P-waves [i.e., ~A/4; see, e.g., Bostock, 1999]. 8.5 Results A l l inversions were performed through a 2-D grid of model points located every 2 km in both xi and xz directions. I begin by presenting inversion results for individual events and scattering Chapter 8. Application to CASC93 89 Table 8.2: Reference, 1-D velocity model for multi-channel inversion of the CASC93 data set. Layer Depth a 0 P number km km/s km/s g/cm3 1 0-5 5.50 2.77 2.6 2 5-8 6.35 3.20 2.6 3 8-11 6.47 3.26 2.6 4 11-16 6.67 3.36 2.9 5 16-21 6.75 3.40 3.1 6 21-25 6.93 3.49 3.1 7 25-40 7.16 3.61 3.1 8 40-300 8.10 4.33 3.5 modes in order to identify their contributions to the final model. Figure 8.6 displays results for the New Ireland and Mexico events (see Figure 8.3 for location), with the grey scale representing neg-ative (slower) to positive (faster) velocity perturbations as white to black shading, respectively. The most prominent feature observed in single mode inversions is a dipping, low-velocity tabu-lar feature in the western portion of the model, associated with the subducting Juan de Fuca plate (see, in particular, Figure 8.6c-h). Results for (#=4) inversions (i.e., Figure 8.6c,g) further image a horizontal discontinuity in the central/eastern portion of the profile, marking a jump from low to high velocity at deeper levels, and inferred to represent the continental Moho. As a result of their relative signal levels, these two features are not equally well resolved within the individual mode inversions. For example, the (q=2) inversion does not recover the Juan de Fuca slab when illuminated from the left (Figure 8.6a), due to near-normal incidence of the direct P-wave. Back-scattered P-P inversion for the same event also displays a poor recovery of structure. In this case, however, the problem may be due to a loss of low-frequency P-P scattered energy in the estim-ated source-time function during preprocessing. Although not shown here, comparable analyses have been conducted for all individual events. This exercise was useful in devising a weighting Chapter 8. Application to CASC93 90 strategy for individual event-mode combinations within the simultaneous inversion. The intent ofthe weighting is to reduce the contributions from scattering modes that display a weakened or null response. Only modes that recovered some signature of the two main structures (i.e., Juan de Fuca crust and continental Moho) were retained for simultaneous inversion (see Table 8.1). In general, structural recovery was consistently poorer for inversions involving q — 2,3 than those for q = 4,6,7; consequently, the latter group were up-weighted by a factor of 10. Simultaneous inversion of all events listed in Table 8.1 was performed using a selection of scattering modes. Individual scattering-mode images are presented in Figure 8.7, and show an ap-preciable improvement in reconstruction over those of the single-event inversions (Figure 8.6). Resulting P- and S- velocity perturbation models for the multimode inversion are displayed in Fig-ure 8.8 with a red to blue colour scale representing negative (slower) to positive (faster) velocity perturbations. The receiver function model of Li [1996] is superposed on the resulting profile in Figure 8.8c-d for comparison. A good correlation is observed between the two models. The dip-ping, low-velocity layer in the western portion of the array has been interpreted by Li [ 1996] as the oceanic crust of the Juan de Fuca plate sandwiched between high-velocity mafic and ultramafic material of the overriding plate and oceanic mantle, respectively. In my resulting model, this layer is well defined from 20 km depth beneath the coast, to 40 km depth beneath the Willamette Val-ley, with a dip of ~12° and an average thickness near 10 km. Below 40 km, the signature of the oceanic crust becomes noticeably weaker. However the oceanic Moho persists to at least 100 km depth, dipping more steeply at ~27°. The continental Moho is imaged at 35-40 km depth from -122.3° (i.e., 150 km east of the coast) to the eastern edge of the model. A n interesting feature in Figure 8.8a-b is the disruption of oceanic crust at approximately 40 km depth below the Willamette Valley, which accompanies the change in dip and signature of the subducting oceanic crust. While this transition in dip is modelled as a simple "elbow" in the Chapter 8. Application to CASC93 91 Event #16 a) P - S Fo rwa rd Sca t te r ing , q=2 (ApVfi) -124 -123.5 - 123 -122.5 121.5 b) P - P B a c k Sca t te r ing , q=3 (Aoc/oc) -124 -123.5 -123 -122.5 -122 -121.5 -121 c) P - S B a c k Sca t te r ing , q=4 (Ap/p) -124 -123.5 -123 -122.5 -122 -121.5 -121 d) S - S B a c k Sca t te r ing , q=6,7 (Ap/p) -124 -123.5 -123 -122.5 -122 Longitude [°] Event #15 e) P - S Fo rwa rd Sca t te r i ng , q=2 (Ap/p) 0 20 E 40 60 0) Q 80 100 120, -124 -123.5 -123 -122.5 121.5 f) P - P B a c k Sca t te r i ng , q=3 jfto/a) 0 20 E 40 j£ £ 60 I 80 100 120. 124 -123.5 -123 -122.5 -122 -121.5 -121 g) P - S B a c k Sca t te r i ng , q=4 (Ap/p) -"124 -123.5 - 123 -122.5 -122 -121.5 -121 h) S - S B a c k Sca t t e r i ng , q=6,7 (Ap/p) 121.5 -121 no vertical "124 -123.5 exaggeration -123 -122.5 -122 Longitude [°] 121.5 Figure 8.6: Results from single mode inversions of events #15 and #16 (see Figure 8.3). Gray scale indicates slow (white) and fast (black) seismic velocities relative to the reference model. Chapter 8. Application to CASC93 92 a) P-S Forward Scattering, q = 2 (ApVfj) os c) P-S Back Scattering, q=4 (\pVp) -123.5 -123 -122.5 -122 Longitude [°] 121.5 b) P-P Back Scattering, q=3 (Aa/oc) 1 2-°124 d) o 20 ? 40 £ 60 I 80 100| 123.5 -123 -122.5 -122 -121.5 -121 S-S Back Scattering, q=6,7 (Apvp) -j^. . — . ^ j . . 120 124 -123.5 no vertical exaggeration -123 -122.5 -122 Longitude [°] -121.5 -121 Figure 8.7: Results from single mode inversions of all events listed in Table 8.1. Gray scale is de-fined as for Figure 8.6. Note the improvement in resolution (e.g., confinement of the low-velocity dipping layer) with increasing sensitivity of total travel time to scatterer location, V T 9 , from for-ward- through back-scattering modes, in a-d. receiver function analysis of [Li, 1996], it appears here to be associated with a local thickening of the effective low-velocity layer and a possible bifurcation of the crust or plate. We further note that the upper boundary of this feature, situated at ~30 km depth beneath the Willamette Valley, coincides with the zone of high reflectivity and possible increased thickness identified by controlled source seismology [see Figure 8.2b; Trehu etal, 1994; Keach etal, 1989]. Synthetic experiments involving layers with differing dips indicate that the reduction in crustal signature with depth is real. Chapter 8. Application to CASC93 93 Figure 8.8: Results from simultaneous inversion of all events and scattering modes listed in Ta-ble 8.1; (a) P- and (b) S-velocity perturbation models, overlain by (c)-(d) contour of receiver func-tion model by Li [1996], and (e)-(f) contour of conductive features (green) modelled by Wanna-maker et al., [1989]. Colour scale indicates slow (red) and fast (blue) seismic velocities relative to the reference model. The western edge of the model coincides with the coastline (C); major physiographic units traversed by the receiver array are the Coast Range (CR), Willamette Valley (WV), Cascade Range (CASC), and Columbia Plateau (CP). Chapter 8. Application to CASC93 8.6 Interpretation of the low velocity crustal disruption 94 In this section, I investigate two possible interpretations for the disruption of oceanic crust and associated change in properties of the subducting slab (i.e., dip increase and weaker signature) apparent in the model of Figure 8.8. These interpretations are (see Figure 8.9): delamina-tion/underplating of the subducting crust resulting in tectonic layering, and the entrapment of aqueous fluids released through metamorphic reactions in the subducting crust. 8.6.1 Delamination-tectonic layering The first interpretation can be further subdivided into two contrasting classes of subduction struc-ture, based on previous geophysical studies of subduction complexes: (1) underplating of sedi-ments and/or mafic material from the oceanic crust to the accretionary wedge along the subduction decollement [Green etal, 1986; Clowes etal, 1987; Hyndman etal, 1990; Calvert and Clowes, 1990; Calvert, 1996; Keach et al, 1989; Trehu et al, 1994]; and (2) continuous, large-scale de-lamination of the subducting crust [Levander et al, 2000; Tsumura et al, 1999]. Underplating of subducted materials to the accretionary prism was first proposed for the northern Cascadia subduction zone, beneath Vancouver Island from LlTHOPROBE studies [e.g., Green et al, 1986; Clowes et al, 1987; Hyndman et al, 1990]. A ubiquitous reflector band la-belled ' E ' was identified that runs subparallel to, and approximately ~5-10 km above the present location of the subduction decollement [see Hyndman et al, 1990, for a review]. The reflector band dips at a slightly lower angle than the subducting plate, and eventually merges seaward into the interplate decollement [Calvert, 1996]. A second layer, labelled ' C , was intermittently im-aged closer to the surface, and separated from ' E ' by a low-reflectivity, high-velocity zone. One interpretation of this layering proposes a succession of underplated marine sediments (reflective Chapter 8. Application to CASC93 95 a) STRUCTURAL LAYERING Longitude [°] b) FLUID MIGRATION AND HYDROUS REACTIONS £ 60 Q. CD Q 80 120 •124 -123.5 •123 -122.5 -122 Longitude [°] •121.5 -121 no vertical exaggeration Continental crust Oceanic crust Continetal lithospheric mantle Oceanic lithospheric mantle Asthenospheric mantle I Mineral hydration front Figure 8.9: Schematic diagrams representing two possible interpretations for the crustal dis-ruption observed in the model of Figure 8.8: a) delamination of the oceanic crust; b) dehydra-tion/devolatilization of the oceanic crust and mineral reactions in the overriding wedge. Chapter 8. Application to CASC93 96 bands ' C and 'E ') and accreted layers of oceanic crust and mantle (high velocity-low reflectivity zone). The underplating of sediments is postulated to operate continuously, whereas the accre-tion of crustal layers is envisaged as an episodic process associated with seaward jumps of the subduction zone. Similar accretionary mechanisms have been proposed for the tectonic regime in central Oregon, based on layering imaged from seismic reflection [Keach et al, 1989] and re-fraction/wide angle reflection [Trehu et al, 1994] studies. In this context, the apparent thickening and disruption of oceanic crust imaged here could represent a zone of continuously accumulated marine sediments or a discretely accreted fragment of oceanic crust [Trehu et al, 1994]. Both mechanisms can account for the observed geometry and velocity contrast; however, it appears im-plausible that either one could be fully responsible for the increase in slab dip and the diminished signature of the more deeply subducted crust. More specifically, the relative volume of removed material is small and unlikely to result in a significant change in buoyancy. Large-scale delamination of the oceanic crust has been inferred from seismic studies of sub-duction zones in southernmost Cascadia [e.g., Levander et al, 2000] and Japan [Tsumura et al, 1999]. In this scenario, a layer of slab material is separated from the underlying plate and directed either towards the surface [e.g., Hidaka Collision Zone; Tsumura et al, 1999] or into low viscos-ity zones such as slab windows [e.g., south-east of the Mendocino Triple junction; Levander et al, 2000]. In northern California [Levander et al, 2000], the Juan de Fuca plate is inferred to enter the subduction zone with a slightly negative net buoyancy with respect to the continental upper mantle. However, upon separation of the crust, the denser mantle portion of the plate sinks at a higher angle [Levander et al, 2000]. This model can explain the increase in slab dip beneath the imaged crustal disruption; however, one might therefore expect the lower velocities to com-pletely disappear beneath 40 km depth, which is not observed. Furthermore, there is no obvious outlet in central Oregon for the delaminated oceanic crust, i.e. no evidence for recently exhumed Chapter 8. Application to CASC93 97 Juan de Fuca oceanic crust in the forearc, and no slab window into which the delaminated crust could feed. 8.6.2 Aqueous fluids and hydrated minerals A growing body of petrological and geophysical evidence suggests that fluids and fluid migra-tion play an important role in shallow subduction zone dynamics [Peacock, 1996; Hacker, 1996; Ernst, 1999; Kurtz etal, 1986; Hyndman, 1988; Wannamaker, 1989; ANCORP Working Group, 1999]. Recent petrological work [Peacock, 1996; Hacker, 1996; Ernst, 1999] indicates that sig-nificant dehydration at shallow subduction zone pressures and temperatures (P-T) results from prograde metamorphic reactions of hydrated minerals such as Na-amphibole and chlorite, con-tained primarily in altered basaltic and gabbroic layers of the crust. Maximum dehydration occurs through the blueschist-to-eclogite facies transition in most steady state subduction regimes [Pea-cock, 1996], but may also occur in a progressive fashion in younger/warmer subduction zones, such as Cascadia, where intermediate facies including greenschist, amphibolite and granulite may be encountered [Hacker, 1996]. Thermal models for the subduction regime in central Oregon [Hyndman and Wang, 1993; Oleskevich et al, 1999] have been used to constrain temperatures along the main subduction decollement. Beneath the coast, the top of the crust is located at 20 km depth (~0.7 GPa) and has an estimated temperature of 450°C, corresponding to P-T con-ditions in the greenschist facies [e.g., Hacker, 1996]. At 40 km depth (~1.3 Gpa) beneath the Willamette Valley, where the crustal disruption and other changes in slab properties are observed, temperatures on the thrust plane are ~550-600°C, yielding P-T conditions that correspond to the amphibolite-to-eclogite facies transition. A magnetotelluric survey undertaken across Vancouver Island by Kurtz et al [1986] was the first to interpret high-conductivity layering subparallel to the plate decollement in terms of Chapter 8. Application to CASC93 98 subduction-entrained fluids. This layering, which coincides with the reflective band ' E ' discussed above, was later interpreted by Hyndman [1988] as a concentration of fluids expelled from the oceanic crust by dehydration. The fluids are inferred to migrate upwards into the overriding lithos-pheric wedge where they encounter lower P-T conditions corresponding to retrograde facies tran-sitions. Hydrous mineral reactions can then occur, forming impermeable boundaries that trap ris-ing fluids and produce high conductivity/reflectivity anomalies in the subduction complex. Hyn-dman [1988] further suggests that progressive silica precipitation due to decreased solubility in rising fluids may also contribute to the formation of such impermeable boundaries. A similar interpretation has been invoked by the E M S L A B group [see Figure 8.2; Wanna-maker et al, 1989] to explain a strong conductivity anomaly beneath central Oregon, in close proximity to the current study. Figure 8.8e-f shows an outline of the resistivity model from Wan-namaker et al. [1989] superimposed on the velocity perturbation model. Recognizing the non-uniqueness of resistivity models, especially in the vertical direction where only the conductance (conductivity-thickness product) of discrete layers can be resolved [Wannamaker et al, 1989], there is reasonably good agreement between the dipping layers imaged by both techniques. We further note that the western limit of enhanced conductivity extending beneath the Willamette Valley and the Cascades, which has been interpreted by Wannamaker et al. [1989] as a zone of trapped fluids, coincides with the location of oceanic crustal thickening and disruption. Consideration of Figure 8.8 in light of this previous work leads to my preferred explanation for the imaged structure in terms of dehydration and eclogitization. The apparent thickening of the oceanic crust and the shallowing of its upper boundary are interpreted to manifest primarily the migration of fluids rather than structural changes within the crust per se. That is, dehydration reactions flux fluids into material above the oceanic crust, lowering its velocity and increasing Chapter 8. Application to CASC93 99 its conductivity. The upper extent of this process is sharply defined by an impermeable bound-ary which merges seaward with the subduction decollement. This boundary may represent a ret-rograde metamorphic front, or it may be structurally controlled. Increased dip of the plate and reduced signature of the oceanic crust at depths immediately below suggest that the dehydration reactions culminate in eclogitization with an accompanying density increase of 10-15% [Hacker, 1996]. The location, beneath the Willamette Valley, where the changes occur is, moreover, con-sistent with the hypothesis of Rogers (1983) that the density increase is also responsible for the development of trench-parallel lowlands (Willamette Valley, Puget Sound, Georgia Strait) in the Cascadia forearc. 8.7 Concluding remarks A new method for imaging detailed lithospheric structures using the teleseismic body-wave coda has been developed as part of a multi-author collaborative effort. In the present work, I have applied the method to a field data set collected over central Oregon and illustrated the contribu-tions made by different scattering modes and source wavefields in the inversion procedure. I have demonstrated that the method can be used to reveal complex structures associated with dynamic processes in the Cascadia subduction zone. In particular, changes in the expression of the subduct-ing Juan de Fuca plate have been related to dehydration and eclogitization of the oceanic crust. It is my hope that the success of this effort will provide motivation for new, ambitious experiments designed to exploit scattered, teleseismic body waves for the resolution of crust and mantle struc-ture at finer scales than previously possible using conventional teleseismic methods. Chapter 9 Conclusions The two parts of this thesis present two different approaches to the analysis of teleseismic data recorded on broadband three-component arrays. The study in Part I involved the design and de-ployment of an array with the objective of understanding lithospheric structure and evolution in the southeastern Canadian Shield, using well established analysis techniques. In Part II, I have discussed a new method for inverting scattered teleseismic waveforms and have demonstrated its successful application to a field data set. Both studies have been successful in realizing the objectives originally set forth. In Part I, analyses of data from the Abitibi-Grenville teleseismic experiment have identified lithospheric structures related to the accretion and post-stabilization reworking of the Canadian Shield. At crustal levels, variations in Moho topography detected near the Grenville front have been correl-ated to comparable structures at other locations along strike of the front, and interpreted as ev-idence for a pre-Grenvillian subduction suture extending over a large portion of the orogen. At deeper levels within the lithospheric mantle, a low-velocity corridor imaged by P-wave travel-time inversion has been interpreted as the signature of processes responsible for the Moteregian-White Mountains-New England Seamounts hotspot track, which have modified the thick lithos-phere of the Canadian Shield. In Part II, the application of an inverse scattering method to the CASC93 data set has success-fully imaged lithospheric and upper mantle structure of the Cascadia subduction zone in central 100 Chapter 9. Conclusions 101 Oregon. In particular, variations in the signature and geometry of the subducting Juan de Fuca plate have been interpreted, along with results from previous geophysical studies, to manifest the occurrence of metamorphic reactions within the oceanic crust. This new treatment of teleseismic array data can thus provide information that is complementary to seismic reflection and refraction studies at a lower overall cost for the experiment. As is the case in most scientific research, these results point to new questions and future re-search opportunities. Three examples of possible future projects are suggested below. First, an investigation of mantle structure below the Canadian Shield in northern Ontario would likely provide additional constraints on the nature of the hotspot track proposed in Part I, in particu-lar, its relation to North American plate motion, ancient zones of structural weakness and in-traplate/kimberlite volcanism. Second, a more detailed examination of Moho structure in the near vicinity of the Grenville Front using seismic reflection and/or dense teleseismic arrays (i.e., for multichannel inversion) is required to fully understand the origin of this feature and its rela-tion to the Grenville Orogeny. Finally, the demonstration that inverse-scattering techniques can be successfully applied to teleseismic data sets will , hopefully, prompt the design and execution of dense-array experiments and allow for more detailed characterization of lithospheric structure than previously possible. 102 References Adams, J., and P. Basham, The seismicity and seismotectonics of eastern Canada, in Neotecton-ics of North America, edited by D. B. Slemmons et al., pp. 261-276, Geol. Soc. of Am., Boulder, Colo., 1991. ANCORP Working Group, Seismic reflection image revealing offset of Andean subduction-zone earthquake locations into oceanic mantle, Nature, 397, 341-344,1999. Anderson, D. L. , Geophysics of the continental mantle: an historical perspective, in Continental Mantle, edited by M . A. Menzies, pp. 1-30, Clarendon, Oxford, England, 1990. Anderson, D. L., The scales of mantle convection, Tectonophysics, 284, 1-17,1998. Antonuk, C. N . , and J. C. Mareschal, Preliminary gravity modelling along the Lithoprobe seismic reflection profiles 28 and 29, northern Abitibi Subprovince, in Proceedings, Lithoprobe Abitibi-Grenville Transect, Rep. 33, pp. 71-75, Lithoprobe, Vancouver, B.C., 1992. B A B E L Working Group, Evidence for early Proterozoic plate tectonics from seismic reflection profiles in the Baltic shield, Nature, 348, 34-38, 1990. Babuska, V , and M . Cara, Seismic Anisotropy in the Earth, 217 pp., Kluwer Acad., Norwell, Mass., 1991. Beaumont, C , and G. Quinlan, A geodynamic framework for interpreting crustal-scale seismic-reflectivity patterns in compressional orogens, Geophys. J. Int., 116, 754-783,1994. Bedard, J. FL, The opening of the Atlantic, the Mesozoic New England Igneous Province, and mechanisms of continental breakup, Tectonophysics, 113, 209-232, 1985. Bellefleur, G., A . Barnes, A. Calvert, C. Hubert, and M . Mareschal, Seismic reflection constraints from Lithoprobe line 29 on the upper crustal structure of the northern Abitibi greenstone belt, Can. J. Earth Sci., 32, 128-134, 1995. Benn, K. , E. W. Sawyer, and J.-L. Bouchez, Orogen parallel and transverse shearing in the Opatica belt, Quebec: Implications for the structure of the Abitibi Subprovince, Can. J. Earth Sci., 29, 2429-2444, 1992. Beylkin, G., Imaging of discontinuities in the inverse scattering problem by inversion of a causal generalized Radon transform,/. Math. Phys., 26, 99-108,1985. Bostock, M . G., A seismic image of the upper mantle beneath the North American craton, Geo-phys. Res. Lett., 23, 1593-1596, 1996. Bostock, M . G., Mantle stratigraphy and evolution of the Slave province, J. Geophys. Res., 103, 21,183-21,200, 1998. Bostock, M . G., Seismic imaging of lithospheric discontinuities and continental evolution, in De-velopments in Geotectonics 24: Composition, deep structure and evolution of continents, edited by R. D. van der Hilst and W. F. McDonough, pp. 1-16, Elsevier, Amsterdam, The Netherlands, 1999. Bostock, M . G., and S. Rondenay, Migration of scattered teleseismic body waves, Geophys. J. Int., 737,732-746, 1999a. REFERENCES 103 Bostock, M . G., and S. Rondenay, Corrigendum to "Migration of scattered teleseismic body waves", Geophys. J. Int., 139, 597, 1999b. Bostock, M . G., S. Rondenay, and J. Shragge, Multi-parameter 2-D inversion of scattered teleseis-mic body-waves-I. Theory for oblique incidence, submitted to J. Geophys. Res., 2000. Bostock, M . G., and J. C. VanDecar, Upper mantle structure of the northern Cascadia subduction zone, Can. J. Earth. Sci., 32, 1-12, 1995. Bowman, J. R., and M . Ando, Shear-wave splitting in the upper-mantle wedge above the Tonga subduction zone, Geophys. J. R. Astron. Soc, 88, 25-41, 1987. Burke, K . C , and J. T. Wilson, Hot spots on the Earth's surface, Sci. Am., 235(2), 46-57, 1976. Calvert, A . J., E. W. Sawyer, W. J. Davis, and J. N . Ludden, Archaean subduction inferred from seismic images of a mantle suture in the Superior Province, Nature, 375, 670-673, 1995. Calvert, A . J., Seismic reflection constraints on imbrication and underplating of the northern Cas-cadia convergent margin, Can. J. Earth Sci., 33, 1294-1307,1996. Calvert, A . J., and R. M . Clowes, Deep, high-amplitude reflections from a major shear zone above the subducting Juan de Fuca plate, Geology, 18, 1091-1094, 1990. Card, K. D., A review of the Superior Province of the Canadian Shield, a product of Archean accretion. Precambrian Res., 48: 99-156, 1990. Card, K. D., and A . Ciesielski, D N A G No. 1: Subdivisions of the Superior Province of the Ca-nadian Shield, Geosci. Can., 13, 5-13, 1986. Carr, S. D., R. M . Easton, R. A . Jamieson, and N . G. Culshaw, Geologic transect across the Grenville Orogen of Ontario and New York, Can. J. Earth Sci.,37, 193-216, 2000. Cassidy, J. F., Review: Receiver function studies in the southern Canadian Cordillera, Can. J. Earth Sci., 32, 1514-1519,1995. Christiansen, R. L., S. R. Yeats, S. A . Graham, W. A. Niem, A . R. Niem, and D. P. Snavely, Post-Laramide geology of the U.S. Cordilleran region, in The Cordilleran Orogen; con-terminous U.S. edited by C. B. Burchfiel, P. W. Lipman and M . L. Zoback, pp. 261-406, Geol. Soc. of Am., Boulder, Colo., 1992. Clowes, R. M . , Lithoprobe Phase IV: Multidisciplinary studies of the evolution of the continent-A progress report, Geosci. Can., 23, 109-123, 1997. Clowes, R. M . , M . T. Brandon, A . G. Green, C. J. Yorath, A . Sutherland Brown, E. R. Kanasewich, and C. Spencer, L iTHOPROBE - sou thern Vancouver Island: Cenozoic sub-duction complex imaged by deep seismic reflections, Can. J. Earth Sci., 24, 31 -51,1987. Cook, F. A. , A . J. van der Velden, K. W. Hall, and B. J. Roberts, Tectonic delamination and sub-crustal imbrication of the Precambrian lithosphere in northwestern Canada mapped by Lithoprobe, Geology, 26, 839-842, 1998. Couch, R. W , and R. P. Riddihough, The crustal structure of the western continental margin of North America, in Geophysical framework of the continental United States, edited by L . C. Pakiser, and W. D. Mooney, pp. 103-128, Geol. Soc. of Am., Boulder, Colo., 1989. REFERENCES 104 Crough, S. T., Mesozoic hotspot epeirogeny in eastern North America, Geology, 9, 2-6, 1981. Crough, S. T., W. J. Morgan, and R. B. Hargraves, Kimberlites: Their relation to mantle hotspots, Earth Planet. Sci. Lett., 50, 260-274, 1980. Culotta, R. C , T. Pratt, and J. Oliver, A tale of two sutures: COCORP's deep seismic surveys of the Grenville province in the eastern U.S. midcontinent, Geology, 18, 646-649,1990. Currie, K . L., The Alkaline Rocks of Canada, 228 pp., GSC Bulletin 239, Ottawa, Ont, 1976. Davies, G. R, Thermomechanical erosion of the lithosphere by mantle plumes, J. Geophys. Res., 99, 15,709-15,722, 1994. DeMets, C , R. G. Gordon, D. R Argus, and S. Stein, Current plate motions, Geophys. J. Int., 101, 423-478, 1990. De Wit, M . , C. Roehring, R. J. Hart, R. A . Armstrong, C. E. J. de Ronde, R. W. E. Green, M . Tredoux, E. Peberdy, and R. A . Hart, Formation of an Archaean continent, Nature, 357, 553-562, 1992. Desrochers, J.-R, C. Hubert, J. N . Ludden, and P. Pilote, Accretion of Archean oceanic plateau fragments in the Abitibi greenstone belt, Canada, Geology 21, 451-454, 1993. Dickin, A . P., and R. H. McNutt, Nd model age mapping of the southeast margin of the Archean foreland in the Grenville province of Ontario, Geology, 17, 299-302, 1989. Dueker, K . G., and A . F. Sheehan, Mantle discontinuity structure from midpoint stacks of conver-ted P to S waves across the Yellowstone hotspot track, J. Geophys. Res., 102, 8313-8327, 1997. Duncan, R. A. , Age progressive volcanism in the New England seamounts and the opening of the central Atlantic Ocean, J. Geophys. Res., 89, 9980-9990,1984. Duncan, R. A. , andM. A . Richards, Hotspots, mantle plumes, flood basalts, and true polar wander, Rev. Geophys., 29, 31-50, 1991. Eaton, D. W , A . Hynes, A. Indares, and T. Rivers, Seismic images of eclogites, crustal-scale ex-tension, and Moho relief in the eastern Grenville Province, Quebec, Geology, 23, 855-858, 1995. Ebinger, C. J., and N . H . Sleep, Cenozoic magmatism throughout east Africa resulting from im-pact of a single plume, Nature, 395, 788-791,1998. E M S L A B Group, The E M S L A B electromagnetic sounding experiment, Eos Trans. AGU, 69, 89, 98-99, 1988. Ernst, W. G., Hornblende, the continent maker-Evolution of H 2 0 during circum-Pacific subduc-tion versus continental collision, Geology, 27, 675-678, 1999. Faure, S., A . Tremblay, and J. Angelier, State of intraplate stress and tectonism of northeastern America since Cretaceous times, with particular emphasis on the New England-Quebec igneous province, Tectonophysics, 255, 111-134, 1996. Fleming, S. W , and A . M . Trehu, Crustal structure beneath the central Oregon convergent margin from potential-field modeling: Evidence for a buried basement ridge in local contact with a seaward dipping backstop, J. Geophys. Res., 104, 20,431-20,447, 1999. REFERENCES 105 Flueh, E. R., M . A . Fisher, J. Bialas, J. R. Childs, D. Klaeschen, N . Kukowski, T. Parsons, D. W. Scholl, U . ten Brink, A . M . Trehu, and N . Vidal, New seismic images of the Cascadia subduction zone from cruise SO108-ORWELL, Tectonophysics, 293, 69-84, 1998. Foland, K . A. , J.-F. Chen, L. A . Gilbert, and A . W. Hofmann, Nd and Sr isotopic signatures of Mesozoic plutons in northeastern North America, Geology, 16, 684-687,1988. Foland, K. A. , and H. Faul, Ages of the White Mountain intrusives-New Hampshire, Vermont, and Maine, USA, ,4m. J. Sci., 277, 888-904, 1977. Foland, K. A. , L. A . Gilbert, C. A . Sebring, and J.-F. Chen, 4 0 A r / 3 9 A r ages for plutons ofthe Monteregian Hills, Quebec: Evidence for a single episode of Cretaceous magmatism, Geol. Soc. Am. Bull, 97, 966-974, 1986. Forgues, E., and G. Lambare, Parameterization study for acoustic and elastic ray+Born inversion, J. Seis. Expl, 6, 253-277, 1997. Goleby, B. R., R. D. Shaw, C. Wright, B. L. N . Kennett, and K. Lambeck, Geophysical evidence for 'thick-skinned' crustal deformation in central Australia, Nature, 337, 325-330,1989. Grand, S. P., Mantle shear structure beneath the Americas and surrounding oceans, J. Geophys. Res., 99, 11,591-11,621, 1994. Grandjean, G., H. Wu, D. J. White, M . Mareschal, and C. Hubert, Crustal velocity models for the Archean Abitibi greenstone belt from seismic refraction data, Can. J. Earth Sci., 32, 149-166, 1995. Green, A . G., R. M . Clowes, C. J. Yorath, C. Spencer, E. R. Kanasewich, M . T. Brandon, and A . Sutherland Brown, Seismic reflection imaging of the subducting Juan de Fuca plate, Nature, 579,210-213, 1986. Green, A . G., B. Milkereit, A . Davidson, C. Spencer, D. R. Hutchinson, W. F. Cannon, M . W. Lee, W. F. Agena, J. C. Behrendt, and W. J. Hinze, Crustal structure of the Grenville Front and adjacent terranes, Geology, 16, 788-792, 1988. Green, A . G., B. Milkereit, L. J. Mayrand, J. N . Ludden, C. Hubert, S. L. Jackson, R. H. Sutcliffe, G. F. West, P. Verpaelst, and A. Simard, Deep structure of an Archean greenstone terrane, Nature, 344, 327-330, 1990. Guillou, L., J.-C. Mareschal, C. Jaupart, C. Gariepy, G. Bienfait, and R. Lapointe, Heat flow, grav-ity and structure of the Abitibi belt, Superior Province, Canada: Implications for mantle heat flow, Earth Planet. Sci. Lett., 122, 103-123, 1994. Hacker, B. R., Eclogite formation and rheology, buoyancy, seismicity, and H 2 O content of oceanic crust, in Subduction: Top to bottom, Geophys. Monogr. 96, edited by E. G. Bebout, D. W. Scholl, S. H. Kirby, and J. P. Piatt, pp. 337-346, A G U , Washington, D . C , 1996. Haggart, M . J., R. A . Jamieson, P. H. Reynolds, T. E. Krogh, C. Beaumont, and N . G. Culshaw, Last gasp of the Grenville Orogeny: Thermochronology of the Grenville Front tectonic zone near Killarney, Ontario, J. Geol, 101, 575-589,1993. Harris, R. A. , H. M . Iyer, and P. B. Dawson, Imaging the Juan de Fuca plate beneath southern Oregon using teleseismic P wave residuals, J. Geophys. Res., 96, 19,879-19889, 1991. REFERENCES 106 Heaton, T. H., and S. H. Hartzell, Earthquake hazards on the Cascadia subduction zone, Science, 236, 162-168, 1987. Hestenes, M . R., and E. Stiefel, Methods of conjugate gradients for solving linear systems, Nat. Bur. Standards J. Res., 49, 409-436, 1952. Hoffman, R R, United plates of America, The birth of a craton: Early Proterozoic assembly and growth of Laurentia, Annu. Rev. Earth Planet. Sci., 16, 543-603, 1988. Hoffman, P. R, Geological constraints on the origin of the mantle root beneath the Canadian Shield, Philos. Trans. R. Soc. London, Ser. A, 331, 523-532,1990. Holbrook, W. S., W. D. Mooney, and N . I. Christensen, The seismic velocity structure of the deep continental crust, in Continental Lower Crust, edited by D. M . Fountain, R. Arculus and R.W. Kay, pp. 1-34, Elsevier, Amsterdam, 1992. Holmden, C , and A. P. Dickin, Paleoproterozoic crustal history of the southwestern Grenville Province: Evidence from Nd isotopic mapping, Can. J. Earth Sci., 32, 472-485, 1995. Hyndman, R. D., Dipping seismic reflectors, electrically conductive zones, and trapped water in the crust over a subducting plate, J. Geophys. Res., 93, 13,391-13,405, 1988. Hyndman, R. D., and K. Wang, Thermal constraints on the zone of major thrust earthquake fail-ure: the Cascadia subduction zone, J. Geophys. Res., 98, 2039-2060, 1993. Hyndman, R. D., C. J. Yorath, R. M . Clowes, and E. E. Davis, The northern Cascadia subduction zone at Vancouver Island: seismic structure and tectonic history, Can. J. Earth Sc., 27, 313-329, 1990. Jackson, S. L., R. H. Sutcliffe, J. N . Ludden, C. Hubert, A . G. Green, B. Milkereit, L. Mayrand, G. F. West, and P. Verpaelst, Southern Abitibi greenstone belt: Archean crustal structure from seismic-reflection profiles, Geology, 18, 1086-1090,1990. Ji, S., S. Rondenay, M . Mareschal, and G. Senechal, Obliquity between seismic and electrical anisotropies as a potential indicator of movement sense for ductile shear zones in the upper mantle, Geology, 24,1033-1036,1996. Jordan, T. H. , Composition and development of the continental tectosphere, Nature, 274, 544-548, 1978. Jordan, T. H. , Mineralogies, densities and seismic velocities of garnet lherzolites and their geo-physical implications, in The Mantle Sample: Inclusions in Kimberlites and Other Vol-canics, Proceedings ofthe Second International Kimberlite Conference, vol. 2, edited by F. R. Boyd and H . O. A. Meyer, pp. 1-14, A G U , Washington, D . C , 1979. Keach, R. W , J. E. Oliver, L. D. Brown, and S. Kaufman, Cenozoic active margin and shallow Cascades structure: COCORP results from western Oregon, Geol. Soc. Am. Bull, 101, 783-794, 1989. Kellett, R. L., A . E. Barnes, and M . Rive, The deep structure of the Grenville Front: A new per-spective from western Quebec, Can. J. Earth Sci., 31, 282-292, 1994. Kellett, R. L., M . Mareschal, and R. D. Kurtz, A model of lower crustal electrical anisotropy for the Pontiac Subprovince of the Canadian Shield, Geophys. J. Int., Ill, 141-150, 1992. REFERENCES 107 Kennedy, C. S., and G. C. Kennedy, The equilibrium boundary between graphite and diamond, J. Geophys. Res., 81, 2467-2470, 1976. Kennett, B. L. N . , The removal of free surface interactions from three-component seismograms, Geophys. J. Int., 104, 153-163, 1991. Kennett, B. L. N . , and E. R. Engdahl, Traveltimes for global earthquake location and phase iden-tification, Geophys. J. Int., 105, 429-465, 1991. Kimura, G., J. N . Ludden, J.-R Desrochers, and R. Hori, A model of ocean-crust accretion for the Superior province, Canada, Lithos, 30, 337-355, 1993. Korsch, R. J., B. R. Goleby, J. H . Leven, and B. J. Drummond, Crustal architecture of central Australia based on deep seismic reflection profiling, Tectonophysics, 288, 57-69, 1998. Kumarapeli, R S., The St. Lawrence rift system, related metallogeny, and plate tectonic models of appalachian evolution, in Metallogeny and Plate Tectonics, Geol. Assoc. Can. Spec. Pap., 14, 301-320, 1976. Kumazawa, M . , and O. L. Anderson, Elastic moduli, pressure derivatives, and temperature de-rivatives of single-crystal olivine and single-crystal forsterite, J. Geophys. Res., 74(25), 5961-5972, 1969. Kurtz, R. M . , J. M . DeLaurier, and J. C. Gupta, A magnetotelluric sounding across Vancouver Island detects the subducting Juan de Fuca plate, Nature, 321, 596-599, 1986. Lacroix, S., andE. W. Sawyer, An Archean fold-thrust belt in the northwestern Abitibi Greenstone Belt: structural and seismic evidence, Can. J. Earth Sci., 32, 97-112, 1995. Langston, C. A. , Structure under Mount Rainier, Washington, inferred from teleseismic body waves, J. Geophys. Res., 84, 4749-4762,1979. Langston, C. A. , Evidence for the subducting lithosphere under southern Vancouver Island and Western Oregon from teleseismic P wave conversions, J. Geophys. Res., 86, 3857-3866, 1981. Lay, T., and T. C. Wallace, Modern global seismology, 521 pp., Academic Press, San Diego, C a l , 1995. Leaver, D. S., W. D. Mooney, and W. M . Kohler, A seismic refraction study of the Oregon Cas-cades, J. Geophys. Res., 89, 3121-3134, 1984. Lerner-Lam, A . L., and T. H. Jordan, How thick are the continents?, J. Geophys. Res., 92,14,007-14,026, 1987. Levander, A. , A . S. Meltzer, T. J. Henstock, and S. P. S. Gullick, Deformation, mass transfer, and crustal recycling at the Mendocino triple junction: what Gorda gives up, North America receives, submitted to Tectonics, 2000. L i , X.-Q. , Deconvolving orbital surface waves for the source duration of large earthquakes and modeling the receiver functions for the earth structure beneath a broadband seismometer array in the Cascadia subduction zone, Ph.D. thesis, 153 pp., Oregon State Univ., Corval-lis, September 1996. REFERENCES 108 Ludden, J., and C. Hubert, Geologic evolution of the Late Archean Abitibi greenstone belt of Canada, Geology, 14, 707-711, 1986. Ludden, J., C. Hubert, A . Barnes, B. Milkereit, and E. Sawyer, A three dimensional perspective of the evolution of Archean crust: L I T H O P R O B E seismic reflection images in the south-western Superior province, Lithos, 30, 357-372, 1993. Ludden, J., and A. Hynes, The Lithoprobe Abitibi-Grenville transect: two billion years of crust formation and recycling in the Precambrian Shield of Canada, Can. J. Earth Sci., 37,459-476, 2000. Mareschal, J. C , C. Jaupart, C. Gariepy, L. Z. Cheng, L. Guillou-Frottier, G. Bienfait, and R. Lapointe, Heat flow and deep thermal structure near the southeastern edge of the Canadian Shield, Can. J. Earth Sci., 37, 399-414, 2000. Mareschal, M . , W. S. Fyfe, J. Percival, and T. Chan, Grain-boundary graphite in Kapuskasing gneisses and implications for lower-crustal conductivity, Nature, 357, 674-676, 1992. Mareschal, M . , R. L. Kellet, R. D. Kurtz, J. N . Ludden, S. Ji, and R. C. Bailey, Archaean cratonic roots, mantle shear zones and deep electrical anisotropy, Nature, 375, 134-137, 1995. Martignole, J., and A. J. Calvert, Crustal-scale shortening and extension across the Grenville province of western Quebec, Tectonics, 15, 376-386, 1996. McBride, J. H. , D. B. Snyder, M . P. Tate, R. W. England, and R. W. Hobbs, Upper mantle reflector structure and origin beneath the Scottish Caledonides, Tectonics, 14, 1351-1367,1995. McHone, J. G., Constraints on the mantle plume model from Mesozoic alkaline intrusions in northeastern North America, Can. Mineral, 34, 325-334, 1996. McHone, J. G., and J. R. Butler, Mesozoic igneous provinces of New England and the opening of the North Atlantic Ocean, Geol. Soc. Am. Bull, 95, 757-765, 1984. Meyer, H. O. A. , M . A. Waldman, and B. L. Garwood, Mantle xenoliths from kimberlite near Kirkland Lake, Ontario, Can. Mineral, 32, 295-306, 1994. Michaelson, C. A. , and C. S. Weaver, Upper mantle structure from teleseismic P wave arrivals in Washington and Northern Oregon, Geophys. Res., 91, 2077-2094,1986. Miller, D., M . Oristaglio, and G. Beylkin, A new slant on seismic imaging: Migration and integral geometry, Geophysics, 52, 943-964, 1987. Mooney, W. D., and C. S. Weaver, Regional crustal structure and tectonics of the Pacific coastal states; California, Oregon, and Washington, in Geophysical framework of the continental United States, edited by L. C. Pakiser, and W. D. Mooney, pp. 129-161, Geol. Soc. of Am., Boulder, Colo., 1989. Morgan, J. V., M . Hadwin, M . R. Warner, P. J. Barton, and R. P. LI. Morgan, The polarity of deep seismic reflections from the lithospheric mantle: Evidence for a relict subduction zone, Tectonophysics, 232, 319-328, 1994. Morgan, W. J., Convection plumes in the lower mantle, Nature, 230, 42-43, 1971. Morgan, W. J., Deep mantle convection plumes and plate motions, Am. Assoc. Pet. Geol. Bull, 55,203-213, 1972. REFERENCES 109 Nabelek, J., X. -Q. L i , S. Azevedo, J. Braunmiller, A . Fabritius, B. Leitner, A . M . Trehu, and G. Zandt, A high-resolution image of the Cascadia subduction zone from teleseismic conver-ted phases recorded by a broadband seismic array, Eos Trans. AGU, 74 (43, fall meeting suppl), 431, 1993. Oleskevich, D. A. , R. D. Hyndman, and K. Wang, The updip and downdip limits to great subduc-tion earthquakes: Thermal and structural models of Cascadia, south Alaska, SW Japan, and Chile, J. Geophys. Res., 104, 14,965-14,991, 1999. Owens, T. J., G. Zandt, and S. R. Taylor, Seismic evidence from an ancient rift beneath the Cum-berland Plateau, Tennessee: A detailed analysis of broadband teleseismic P waveforms, J. Geophys. Res., 89, 7783-7795, 1984. Peacock, S. M . , Thermal and petrologic structure of subduction zones, in Subduction: Top to bot-tom, Geophys. Monogr. 96, edited by E. G. Bebout, D. W. Scholl, S. H. Kirby, and J. P. Piatt, pp. 119-133, A G U , Washington, D . C , 1996. Percival, J. A. , and G. F. West, The Kapuskasing uplift: A geological and geophysical synthesis, Can. J. Earth Sci., 31, 1256-1286, 1994. Pfiffner, O. A. , W. Frei, P. Valasek, M . Stauble, L. Levato, L. DuBois, S. M . Schmid, and S. B. Smithson, Crustal shortening in the Alpine Orogen: Results from deep seismic reflection profiling in the eastern Swiss Alps, line NFP 20-east, Tectonics, 9, 1327-1355,1990. Rasmussen J., and E. Humphreys, Tomographic image of the Juan de Fuca plate beneath Wash-ington and Western Oregon using teleseismic P-wave travel times, Geophys. Res. Lett., 15, 1417-1420, 1988. Revenaugh, J., A scattered-wave image of subduction beneath the Transverse Ranges, Science, 268, 1888-1892, 1995. Revenaugh, J., and T. H. Jordan, Mantle layering from ScS reverberations, 3, The upper mantle, J. Geophys. Res., 96, 19,781-19,810, 1991. Rivers, T., Lithotectonic elements of the Grenville province: Review and tectonic implications, Precambrian Res., 86, 117-154,1997. Rivers, T., and D. Corrigan, Convergent margin on southeastern Laurentia during the Meso-proterozoic: tectonic implications, Can. J. Earth Sci., 37, 359-383, 2000. Rivers, T., J. Martignole, C. F. Gower, and A. Davidson, New tectonic divisions of the Grenville province, southeast Canadian Shield, Tectonics, 8, 63-84,1989. Rogers, G. C , Seismotectonics of British Columbia, Ph.D. thesis, 247 pp., Univ. of British Co-lumbia, Vancouver, Canada, January 1983. Rondenay, S., Anisotropic sismique et electrique dans le manteau superieur, au niveau du Front de Grenville (Canada), M.A.Sc. thesis, 211 pp., Ecole Polytechnique de Montreal, Montreal, April 1996. Rondenay, S., M . G. Bostock, T. M . Hearn, D. J. White, and R. M . Ellis, Lithospheric assembly and modification of the SE Canadian Shield: Abitibi-Grenville experiment, J. Geophys. Res., 105, 13,735-13,754, 2000a. REFERENCES 110 Rondenay, S., M . G. Bostock, T. M . Hearn, D. J. White, H . Wu, G. Senechal, S. Ji, and M . Mareschal, Teleseismic studies of the lithosphere below the Lithoprobe Abitibi-Grenville Transect, Can. J. Earth Sci., 37, 415-426, 2000b. Rondenay, S., M . G. Bostock, and J. Shragge, Multi-parameter 2-D inversion of scattered tele-seismic body-waves-Ill. Application to the Cascadia-93 data set, submitted to J. Geophys. Res., 2000c. Ryberg, T , and M . Weber, Receiver function arrays: a reflection seismic approach, Geophys. J. Int., 141, 1-11,2000. Sacks, I. S., The subduction of young lithosphere, J. Geophys. Res., 88, 3355-3366, 1983. Saltzer, R. L., and E. D. Humphreys, Upper mantle P wave velocity structure of the eastern Snake River Plain and its relationship to geodynamic models of the region, J. Geophys. Res., 102, 11,829-11,841, 1997. Sawyer, E. W., and K. Benn, Structure of the high-grade Opatica belt and adjacent low-grade Abitibi Subprovince and Archaean mountain front, J. Struct. Geol., 15, 1443-1458,1993. Senechal, G., S. Rondenay, M . Mareschal, J. Guilbert, and G. Poupinet, Seismic and electrical anisotropies in the lithosphere across the Grenville Front, Canada, Geophys. Res. Lett., 23, 2255-2258, 1996. Sheehan, A . F., P. M . Shearer, H. J. Gilbert, and K. G. Dueker, Seismic migration processing of P-SV converted phases for mantle discontinuity structure beneath the Snake River Plain, western United States,/. Geophys. Res., 105, 19055-19065,2000. Shragge, J., M . G. Bostock, and S. Rondenay, Multi-parameter 2-D inversion of scattered tele-seismic body-waves-II. Numerical examples, submitted to J. Geophys. Res., 2000. Silver, P. G., and W. W. Chan, Implications for continental structure and evolution from seismic anisotropy, Nature, 335, 34-39,1988. Silver, P. G., and W. W. Chan, Shear wave splitting and subcontinental mantle deformation, J. Geophys. Res., 96, 16,429-16,454, 1991. Sleep, N . H. , Hotspots and mantle plumes: Some phenomenology, J. Geophys. Res., 95, 6715-6736, 1990a. Sleep, N . H. , Monteregian hotspot track: A long-lived mantle plume, J. Geophys. Res., 95, 21,983-21,990, 1990b. Sleep, N . H. , Lateral flow of hot plume material ponded at sublithospheric depths, J. Geophys. Res., 7(97,28,065-28,083, 1996. Sleep, N . H. , Lateral flow and ponding of starting plume material, / . Geophys. Res., 102,10,001-10,012, 1997. Sobolev, S. V., H . Zeyen, G. Stoll, F. Werling, R. Altherr, and K. Fuchs, Upper mantle tempera-tures from teleseismic tomography of French Massif Central including effects of compo-sition, mineral reactions, anharmonicity, anelasticity and partial melt, Earth Planet. Sci. Lett., 139, 147-163, 1996. REFERENCES 111 Sykes, L. R., Intraplate seismicity, reactivation of preexisting zones of weakness, alkaline mag-matism, and other tectonism postdating continental fragmentation, Rev. Geophys., 16(4), 621-688, 1978. Trehu, A . M . , I. Asudeh, T. M . Brocher, J. H . Luetgert, W. D. Mooney, J. L. Nabelek, Y. Na-kamura, Crustal architecture of the Cascadia forearc, Science, 266, 237-243,1994. Tsumura, N . , H . Ikawa, T. Ikawa, M . Shinohara, T. Ito, K. Arita, T. Moriya, G. Kimura, T. Ikawa, Delamination-wedge structure beneath the Hidaka collision zone, central Hokkaido, Ja-pan inferred from seismic reflection profiling, Geophys. Res. Lett., 26, 1057-1060,1999. Ulrych, T. J., M . D. Sacchi, and S. L. M . Freire, Eigenimage processing of seismic sections, in Covariance Analysis of Seismic Signal Processing, edited by R. L. Kirlin, and W. J. Done, SEG Monograph, Tulsa, 1998. van der Lee, S., and G. Nolet, The upper mantle S velocity structure of North America, J. Geophys. Res., 102, 22,815-22,838, 1997. VanDecar, J. C , Upper-mantle structure of the Cascadia subduction zone from non-linear tele-seismic travel-time inversion, Ph.D. thesis, 165 pp., Univ. of Wash., Seattle, June 1991. VanDecar, J. C , and R. S. Crosson, Determination of teleseismic relative phase arrival times using multi-channel cross-correlation and least squares, Bull. Seismol. Soc. Am., 80, 150-159, 1990. Vinnik, L. P., Detection of waves converted from P to SV in the mantle, Phys. Earth Planet. Inter., 15, 39-45, 1977. Vinnik, L. P., R. W. E. Green, and L. O. Nicolaysen, Recent deformations of the deep continental root beneath southern Africa, Nature, 375, 50-52, 1995. Vinnik, L. P., G. L. Kosarev, and L. I. Makeyeva, Azimuthal anisotropy of the Earth's interior from observations of long-period body waves, Izv. Earth Phy, Engl. Transl., 22(11), 955-960, 1986. Vinnik, L. P., L. I. Makeyeva, A . Milev, and A . Y. Usenko, Global patterns of azimuthal anisotropy and deformations in the continental mantle, Geophys. J. Int., Ill, 433-447, 1992. Wannamaker, P. E., J. R. Booker, A. G. Jones, A . D. Chave, J. H . Filloux, H . S. Waff, and L. K . Law, Resistivity cross section through the Juan de Fuca subduction system and its tectonic implications, J. Geophys. Res., 94, 14,127-14,144, 1989. Warner, M . , J. Morgan, P. Barton, P. Morgan, C. Price, and K. Jones, Seismic reflections from the mantle represent relict subduction zones within the continental lithosphere, Geology, 24, 39-42, 1996. Weaver, C. S., and C. A . Michaelson, Seismicity and volcanism in the Pacific Northwest: evi-dence for the segmentation of the Juan de Fuca plate, Geophys. Res. Lett., 12, 215-218, 1985. Wells, R. E., C. S. Weaver, and R. J. Blakely, Fore-arc migration in Cascadia and its neotectonic significance, Geology, 26, 759-762, 1998. REFERENCES 112 White, D. J., D. A . Forsyth, I. Asudeh, S. D. Carr, H. Wu, and R. M . Easton, A seismic-based cross-section of the Grenville Orogen in southern Ontario and western Quebec, Can. J. Earth Sci., 37, 183-192, 2000. White, R., and D. M cKenzie, Magmatism at rift zones: The generation of volcanic continental margins and flood basalts, J. Geophys. Res., 94, 7685-7729,1989. Winardhi, S., and R. F. Mereu, Crustal velocity structure of the Superior and Grenville Provinces of the southeastern Canadian Shield, Can. J. Earth Sci., 34, 1167-1184, 1997. Wolfe, C. J., and P. G. Silver, Seismic anisotropy of oceanic upper mantle: Shear wave splitting methodologies and observations, J. Geophys. Res., 103, 749-771,1998. Zoback, M . L., Stress field constraints on intraplate seismicity in eastern North America, J. Geo-phys. Res., 97, 11,761-11,782, 1992. Appendix A Data set for analyses in Chapter 3 113 Appendix A. Data set for analyses in Chapter 3 114 Table A . 1: Events selected for analyses along the Abitibi-Grenville Transect. Analysis techniques are: TTI-P, P-wave travel time inversion; TTI-S1, S'-wave travel time inversion; RF, receiver func-tions; SWS, shear-wave splitting. Date Time Latitude Longitude Depth m TTI-P TTI-S RF SWS Y Y : M M : D D hh:mm:ss.s °N °E km 94:07:14 18:38:08 55.36 -163.85 156 5.2 X 94:07:21 18:36:30 42.36 132.91 458 6.5 X 94:08:08 07:55:40 -13.84 -68.33 609 5.4 X 94:08:14 01:31:17 44.58 149.98 19 6.2 X 94:08:18 04:42:57 44.63 150.13 14 6.2 X 94:08:19 10:02:51 -26.60 -63.38 558 6.4 X 94:08:20 04:38:50 44.59 149.18 24 6.2 X 94:08:21 15:55:59 56.74 117.91 12 5.8 X 94:08:28 18:37:20 44.77 150.04 18 6.1 X 94:08:30 06:13:35 44.62 150.05 51 6.2 X 94:08:31 09:07:26 43.61 146.01 76 6.0 X 94:09:01 15:15:53 40.43 -125.69 10 6.6 X 94:09:12 06:29:56 -31.10 -71.72 49 5.9 X 94:09:13 10:01:32 7.14 -76.66 13 5.8 X 94:09:27 23:04:51 5.66 -79.18 33 5.3 X 94:10:04 13:22:57 43.65 147.38 28 7.2 X 94:10:07 02:36:08 43.46 147.33 52 6.1 X 94:10:09 07:55:34 43.89 147.95 33 6.5 X 94:10:09 12:24:21 43.76 147.39 39 5.8 X 94:10:10 21:06:54 51.44 -173.88 33 5.6 X 94:10:11 01:37:21 -32.11 -71.45 57 5.7 X 94:10:16 05:09:59 45.74 149.21 102 6.4 X 94:10:25 00:54:34 36.39 70.99 238 5.9 X 94:10:25 21:00:40 -27.49 -63.21 572 5.1 X 94:10:31 22:59:26 9.05 -83.32 56 5.6 X 94:11:04 01:13:22 -9.40 -71.30 617 5.8 X 94:11:05 12:05:29 -9.38 -71.32 605 5.7 X 96:05:18 03:30:19 55.07 168.17 30 5.1 X 96:05:18 07:42:21 -23.95 -68.82 97 5.3 X 96:05:19 21:19:05 28.79 139.49 382 5.0 X 96:05:23 01:57:22 5.90 -77.58 33 5.5 X X X 96:06:02 00:50:37 -9.68 -79.58 33 5.3 X 96:06:02 02:19:30 30.45 -41.97 20 5.4 X 96:06:02 02:48:47 30.53 -41.89 10 5.3 X X 96:06:02 02:52:09 10.79 -42.25 10 6.1 X X X 96:06:03 19:55:31 46.78 153.72 33 5.6 X 96:06:08 02:55:57 41.65 88.69 0 5.9 X 96:06:08 05:26:09 58.35 -31.98 10 5.0 X 96:06:08 23:19:15 51.49 -178.12 33 5.9 X X X Appendix A. Data set for analyses in Chapter 3 Table A. 1: (Continued) 115 Date Time Latitude Longitude Depth rrifc TTI-P TTI-S RF SWP Y Y : M M : D D hh:mm:ss.s °N °W km 96:06:09 01:12:16 17.44 145.45 149 6.0 X 96:06:10 04:03:35 51.56 -177.63 33 6.6 X X X 96:06:10 04:49:31 51.09 -177.82 33 5.4 X 96:06:10 05:00:30 51.15 -177.52 33 5.2 X 96:06:10 05:34:40 51.14 -177.87 33 5.3 X 96:06:10 08:29:15 51.21 -176.52 33 5.1 X 96:06:10 09:41:05 51.21 -176.48 33 5.4 X 96:06:10 15:24:56 51.47 -176.84 26 5.9 X X X 96:06:10 17:44:16 51.57 -177.04 33 5.4 X 96:06:11 10:40:08 51.27 -176.25 33 5.6 X 96:06:11 16:52:12 17.25 -68.28 33 5.5 X 96:06:12 02:16:48 51.42 -178.21 33 5.5 X 96:06:14 09:17:36 44.61 150.27 55 5.4 X 96:06:14 22:45:40 42.47 72.92 33 5.2 X 96:06:21 13:57:10 51.56 159.11 20 6.0 X X X 96:06:22 14:50:07 51.39 159.23 33 5.6 X X X 96:06:22 16:47:12 75.82 134.61 10 5.6 X X X 96:06:23 01:17:56 51.57 159.54 33 5.4 X 96:06:23 12:45:06 51.58 159.50 33 5.3 X X 96:06:23 14:58:30 51.50 159.41 33. 5.1 X 96:06:23 15:35:25 51.54 159.58 33 5.1 X 96:06:26 03:22:03 27.72 139.74 469 5.5 X X X 96:07:03 16:48:27 -23.37 -70.40 33 5.3 X X X 96:07:03 18:59:26 40.68 142.57 49 5.0 X 96:07:04 11:39:39 61.85 -150.83 55 5.6 X 96:07:06 21:36:28 21.96 142.83 241 5.8 X 96:07:07 10:49:59 58.62 157.75 10 5.6 X X X 96:07:10 05:48:19 52.16 -171.14 33 5.3 X 96:07:13 15:10:30 51.35 -177.85 33 5.2 X 96:07:15 16:51:22 18.72 145.62 177 5.9 X 96:07:15 21:23:34 17.60 -100.96 18 5.7 X 96:07:16 03:48:28 56.08 164.99 33 5.8 X X 96:07:16 06:12:28 56.18 164.56 33 5.0 X 96:07:16 10:07:36 1.01 120.25 33 6.0 X 96:07:18 22:55:03 51.49 159.44 33 5.5 X X 96:07:20 00:00:41 36.14 27.10 33 5.7 X X 96:07:21 00:09:37 51.21 -169.03 33 5.3 X 96:07:22 08:30:21 13.07 -88.72 61 5.2 X X X 96:07:24 17:23:56 35.99 68.76 33 5.1 X 96:07:24 20:15:44 41.78 -125.91 10 5.4 X X X 96:07:24 23:53:22 51.70 -177.12 55 5.1 X 96:07:25 03:30:49 51.87 160.01 41 5.1 X X 96:07:26 18:55:50 40.08 20.65 11 5.2 X Appendix A. Data set for analyses in Chapter 3 Table A.l : (Continued) 116 Date Time Latitude Longitude Depth m& TTI-P TTl-S RF SWP Y : M M : D D hh:mm:ss.s °N °W km 96:08:01 21:51:47 -17.81 -70.62 74 5.0 X X 96:08:02 12:55:29 -10.76 161.44 33 6.2 96:08:05 00:25:15 54.58 161.55 33 5.1 X 96:08:05 02:08:58 -15.26 -173.12 41 6.0 X 96:08:05 21:39:16 -1.99 -81.00 33 5.7 X X 96:08:05 22:38:22 -20.69 -178.31 550 6.4 X X 96:08:08 17:10:52 53.06 -167.09 44 5.7 X 96:08:10 18:12:17 38.90 140.53 10 6.0 X 96:08:10 18:54:11 38.93 140.55 10 5.7 X 96:08:10 23:10:45 38.88 140.62 10 5.6 X 96:08:11 01:59:23 38.84 140.50 10 5.2 X 96:08:13 19:33:40 -15.70 -13.20 10 5.6 X X 96:08:14 01:55:02 40.75 35.34 10 5.3 X 96:08:15 06:25:30 -17.34 -71.14 58 5.1 X X 96:08:19 04:19:16 51.45 -178.36 33 5.7 X X 96:08:20 17:19:56 5.36 -82.63 10 5.2 X X 96:08:23 21:56:06 -4.08 -104.37 10 5.0 X X 96:08:23 22:18:57 -3.92 -104.14 10 5.0 X X 96:08:25 14:09:03 -1.08 -78.67 51 5.1 X 96:08:27 14:36:39 -6.93 -12.72 10 5.2 X 96:08:28 17:16:17 9.37 -84.31 33 5.5 X 96:08:29 20:22:15 1.08 -28.18 10 5.1 X 96:08:30 21:13:41 52.31 151.48 580 5.1 X 96:08:31 20:47:21 51.49 -178.21 44 5.6 X 96:09:03 09:57:26 40.20 142.65 33 5.1 X 96:09:04 04:14:03 36.98 2.87 10 5.3 X 96:09:04 19:06:49 9.36 -84.26 33 5.8 X 96:09:05 08:14:14 -22.11 -113.43 10 6.2 X X 96:09:05 09:10:20 -22.05 -113.09 10 5.6 X 96:09:05 09:46:59 -22.05 -113.08 10 5.6 X X 96:09:05 20:44:09 42.80 17.93 10 5.6 X 96:09:08 08:08:13 -15.57 -73.04 98 5.4 X 96:09:09 00:20:39 -31.90 -71.56 39 6.0 X X 96:09:09 03:21:47 -35.40 -104.64 10 5.0 X 96:09:09 04:34:20 30.43 130.72 33 5.5 X 96:09:11 02:37:14 35.53 140.94 55 6.1 X X 96:09:11 06:28:45 4.25 -76.57 109 5.0 X 96:09:17 23:09:24 0.90 -26.37 10 5.3 X 96:09:17 23:47:07 43.64 147.21 33 5.5 X 96:09:18 17:34:20 11.43 -85.47 193 5.3 X 96:09:21 02:53:18 -18.99 -67.53 224 5.0 X 96:09:22 02:21:33 -15.90 -71.72 140 5.0 X 96:09:22 02:29:07 34.04 140.62 55 5.5 X Appendix A. Data set for analyses in Chapter 3 117 Table A. 1: (Continued) Date Time Latitude Longitude Depth rrii TTI-P ms RF SWP Y Y : M M : D D hh:mm:ss.s °N °W km 96:09:24 11:42:18 15.19 -61.44 147 6.0 X X 96:09:28 12:17:59 11.68 -86.47 33 5.2 X 96:09:28 17:07:04 51.50 161.50 33 5.0 X 96:09:28 18:24:12 -18.16 -70.03 56 5.2 X 96:09:29 10:48:18 64.78 -17.56 10 5.3 X X 96:09:30 11:27:51 45.53 151.82 33 5.5 X 96:10:01 00:52:50 51.87 -166.60 10 5.6 X 96:10:01 23:04:12 -12.68 -76.81 62 5.5 X 96:10:02 11:24:48 45.13 151.16 33 6.1 X X X 96:10:02 16:00:20 53.41 160.38 33 5.6 X 96:10:05 07:37:59 16.79 -86.07 10 5.3 X X 96:10:06 15:32:30 13.22 -44.96 10 5.1 X 96:10:06 20:13:09 49.04 -127.88 10 5.8 X X X 96:10:08 07:52:58 52.87 152.53 627 5.2 X 96:10:09 07:12:25 49.73 -129.60 10 5.3 X 96:10:09 13:10:52 34.55 32.12 33 6.4 X X X 96:10:10 01:10:22 34.56 32.21 33 5.4 X 96:10:13 20:56:12 52.62 160.51 33 5.1 X 96:10:15 00:22:38 43.70 147.10 33 5.5 X 96:10:15 09:55:59 44.79 10.78 10 5.3 X 96:10:17 15:38:20 60.11 -152.95 117 5.4 X 96:10:18 09:26:03 27.69 57.5.7 33 5.4 X 96:10:18 13:41:03 -32.43 -70.04 122 5.2 X 96:10:18 16:44:47 33.68 137.40 338 5.4 X 96:10:19 14:44:40 31.88 131.46 22 6.3 X X X 96:10:22 22:15:02 63.34 -145.35 4 5.7 X 96:10:23 12:05:49 44.65 149.46 33 5.5 X 96:10:24 19:31:53 66.98 -173.22 20 6.0 X X X 96:10:24 21:57:37 67.08 -173.30 18 5.0 X 96:10:25 19:59:41 -17.37 -69.98 116 5.5 X Appendix B Scattering potentials This appendix gives the full expressions for the scattering potentials, /*(x,0) = £ ^ ( 0 ) A m , ( x ) , (B.l) 1=1 which include the ^-dependent radiation patterns W~i(9) and material property perturbations Am;(x)=(Aa/a, A/3//3, Ap/p). These expressions are derived in Paper I for different scattering modes: ^^1 + ^9+^(^29-1) .0 L £ . 0° \ a° k /2'4(x, 9) = p° | | 2^sin20 | + = £ ( sin* + ^s in20 /6(x, 61) = / ( 2 cos20) + ^ (cosfl + cos20)) , / 7 ( X >') = A»° ( 2 C O S *) + y ( i + c o s * ) ) • (B.2) (B.3) (B.4) (B.5) (B.6) 118 


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items