Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Crustal velocity structure in the southern coast belt, British Columbia, from a seismic refraction survey O'Leary, Deirdre 1992

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

Item Metadata


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

Full Text

We accept this thesis as conformingto the required standardCRUSTAL VELOCITY STRUCTURE IN THE SOUTHERN COAST BELT,BRITISH COLUMBIA, FROM A SEISMIC REFRACTION SURVEYbyDeirdre O'LearyB.Sc. (Geophysics), McGill University, 1990A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTER OF SCIENCEinTHE FACULTY OF GRADUATE STUDIESDEPARTMENT OF GEOPHYSICS AND ASTRONOMYTHE UNIVERSITY OF BRITISH COLUMBIADecember 1992©Deirdre O'Leary, 1992In presenting this thesis in partial fulfilment of the requirements for an advanceddegree at the University of British Columbia, I agree that the Library shall make itfreely available for reference and study. I further agree that permission for extensivecopying of this thesis for scholarly purposes may be granted by the head of mydepartment or by his or her representatives. It is understood that copying orpublication of this thesis for financial gain shall not be allowed without my writtenpermission.(SignaturDepartment of  G el:v(1y g 1%cs an et A c÷ron 0 r-t..The University of British ColumbiaVancouver, CanadaDate  Dec. 24 ) Iqq2DE-6 (2/88)ABSTRACTThe Coast Belt is one of five morphogeological belts of the Canadian Cordillera. It wascreated through the complex tectonic processes that accommodated the Mesozoic accretionof the allochthonous Insular superterrane to the Intermontane superterrane, then the leadingedge of western North America, and the subsequent overprinting of the suture zone bygranitic intrusions associated with east-dipping subduction. As part of the LITHOPROBESouthern Cordillera transect, seismic refraction data were recorded along a 350 km longstrike profile in the Coast Belt. An iterative combination of two-dimensional travel timeinversion and amplitude forward modelling was used to interpret crust and upper mantleP-wave velocity structure. The structural model features a thin (0.5 to 3.0 kin) near-surfacelayer with an average velocity of 4.40 km/s, a significant vertical velocity gradient and lateralvariations. This uppermost stratum overlies a crustal velocity structure with three layers ofapproximately 10.0 km each. The upper and middle crust have average velocities of 6.20 and6.30 km/s, respectively. The lower crust consists of an upper and lower unit with averagevelocities of 6.50 and 6.65 km/s, respectively. Beneath the lower crust lies Moho, the crust-mantle boundary, which is modelled as an approximately 2 km thick transitional layer withan average depth to its upper boundary of 34.5 km and a maximum depth of 35.5 km inthe southeast. The transitional layer exhibits lateral velocity variations (7.40 to 7.80 km/s)and overlies an upper mantle which is poorly defined in terms of velocity, although the dataindicate relatively high values (> 8.05 km/s).Interpretations of LITHOPROBE reflection data and other refraction lines which cross,or nearly cross, the Coast Belt profile suggest that this profile is within the region of thecollision zone between the Insular and Intermontane superterranes. These interpretationsindicate two models for the collision zone: (i) a crustal delamination model in which theInsular superterrane was displaced along east-vergent faults over the terranes below and (ii)a crustal imbrication model in which imbrication of Insular and Intermontane rocks occursiithroughout the crust. The latter model involves the presence of thick layers of Insularmaterial beneath the Coast Belt refraction profile. However, the velocity model indicatespredominantly Intermontane material, thereby favoring the crustal delamination model andthe eastward extent of the Insular Belt being west of the refraction profile. From comparisonsof the refraction velocity model with the reflection data, the top of the Moho transition zonecorresponds to the top of a prominent band of reflectivity which constitutes the reflectionMoho. When combined with the seismic reflection as well as other geophysical studies, threelikely sources for wide-angle reflections from the upper, middle, and lower crust, as observedin the refraction data, are structural features (e.g. fault zones), Ethological contrasts (e.g.batholithic granites to accreted volcanics) and the transition zones at the top and bottom ofa region of layered porosity in the crust.Table of ContentsABSTRACT^  iiList of Tables  viiList of Figures ^  viiiAcknowledgements  ixChapter I^INTRODUCTION ^ 11.1 Tectonic History — Terrane Accretion  11.2 Geological and Geophysical Background ^ 61.2.1 Geology of Southern Coast Belt  61.2.2 Previous Geophysical Studies ^  81.3 Experimental Objectives  10Chapter II^DATA ACQUISITION AND PROCESSING ^ 112.1 The Refraction Experiment^  112.1.1 LITHOPROBE Southern Cordillera Transect ^  112.1.2 SCoRE '89 & SCoRE '90 ^  112.1.3 Line 10 — Acquisition Geometry of the Field Data^ 122.1.4 Shot and Receiver Site Positioning ^  132.2 Instrumentation ^  132.3 Processing  142.3.1 Initial Data Reduction ^  142.3.2 Filtering ^  142.3.3 Amplitude Corrections ^  15ivChapter III^DATA ANALYSIS AND INTERPRETATION PROCEDURES . . 173.1 Initial Observations and Interpretations ^  173.1.1 Data Characteristics and Quality  173.1.2 Travel Time Picking ^  363.2 Interpretation Methods  373.2.1 The Modelling Algorithms ^  373.2.2 Inversion mathematics  403.2.3 The Modelling Techniques ^  443.2.4 Modelling the Refraction Data  453.2.4a The Upper Crust ^  453.2.4b The Middle Crust  533.2.4c The Lower Crust and Upper Mantle ^  55Chapter IV^INTERPRETATION^  614.1 The Final Model^  614.1.1 Primary Features of the Velocity Model^  614.1.2 Non-uniqueness of the Model^  634.2 Spatial Resolution and Absolute Parameter Uncertainty ^ 63Chapter V^DISCUSSION AND CONCLUSIONS  665.1 Comparison with LITHOPROBE Refraction Interpretations ^ 665.2 Comparison with Other Geophysical Studies — Reflection, Heat Flow, Gravity,and Geomagnetic ^  715.2.1 Comparison with LITHOPROBE Reflection Data and Interpretations . . . ^ 715.2.2 Heat Flow Data and Interpretations ^  835.2.3 Gravity Data and Interpretations  845.2.4 Geomagnetic Data and Interpretations ^  875.3 Conclusions ^  88Appendix A^Terrane Descriptions^  98viList of Tables1. Noisy traces deleted from true amplitude data sections^  162. Number of observations (with corresponding uncertainties) for each shotpoint . . ^ 313. Source-receiver offsets of identified phases ^  344. Source-receiver offsets and reduced arrival times of critical points of PmP ^ 365. Q structure used in amplitude modelling ^  396. Inversion results for the final velocity model  497. Estimated lateral resolution and absolute parameter uncertainty^ 658. Comparison of interpreted velocity structure for Line 10 and Line 2 of SCoRE'89 ^ 679. Comparison of interpreted velocity structure for Line 10 and Line 3 of SCoRE'89 ^ 6810. Comparison of interpreted velocity structure for Line 10, Line 1 of ScoRE'89 and^Line I of VISP '80    70viiList of Figures1. Map of the Canadian Cordillera showing accreted terranes ^ 22. Map of study area showing Line 10 and other regional seismic profiles ^ 43. Collisional model of the creation of the Coast Belt ^ 54. Location map of Line 10 showing the geological features of the surrounding area ^ 75. Comparisons of observed and calculated data for SP-46^  186. Comparisons of observed and calculated data for SP-47  217. Comparisons of observed and calculated data for SP-48^  238. Comparisons of observed and calculated data for SP-49  259. Comparisons of observed and calculated data for SP-50^  2710. Comparisons of observed and calculated data for SP-51  2911. Final velocity model and schematic interpretation^  4612. Location of nodes used in inversion process  4713. Ray tracing diagrams ^  4814. Refraction comparisons — 1–D velocity versus depth profiles ^ 6915. Reflection data from LITHOPROBE profiles 88-13,88-14, and 88-17 ^ 7316. Comparison of refraction model with results of 88-13 interpretation ^ 7617. Comparison of refraction model with results of 88-14/88-17 interpretation ^ 7718. Reflection interpretations ^  8119. Bouguer gravity anomaly map  85vu 'AcknowledgementsI wish to thank my supervisor Dr. Ron Clowes for his deft "coaching" abilities. I wouldalso like to thank my co-supervisor Dr. Robert M. Ellis and Dr. Garry Clarke for criticallyreviewing this thesis.To my fellow linemates Barry Clint Zelt (#3), Mike Perz (#18), and Brad Isbell (#X):you made the GPRF room a place of final solutions, insane humour, and impeccable safety.My thanks also go to manager John Amor, for managing both the system and my nerves.The other departmental team members have made this a great place to study — I am gratefulto all of them for their friendship.And, mesdames et messieurs, for their constant encouragement, the MOLSON/GRPFthree star selection:1. Les premiere etoiles, the first stars, my MOM and DAD!2. La deuxieme etoile, the second star, DAVE HADLEY!3. La troisieme etoile, the third star, CLAIRE DAT!SCoRE '90 was carried out by about 40 participants under the co-leadership of Dr. E.R. Kanasewich (U. of Alberta) and Dr. R. M. Ellis (U.B.C.). I thank all of them for theirefforts in the field. Dr. Colin Zelt provided the ray-tracing inversion code used in this study.SCoRE '90 data were acquired primarily with funds from the NSERC Collaborative SpecialProject and Program grant in support of LITHOPROBE, with additional funds from theGeological Survey of Canada and NSERC Research grants of E. R. Kanasewich, R. M. Ellisand R. M. Clowes. Financial support for my student stipend and analysis costs derived fromNSERC Research and Infrastructure grants to R. M. Clowes. Additional financial supportwas provided by Petro-Canada Inc.ix1I. INTRODUCTION1.1 Tectonic History — Terrane AccretionThe Cordillera of North America is one of Earth's great orogenic belts. Until the 1960's,it was believed that Cordilleran orogenesis was the result of geosynclinal deposition (Monger1992). The advent of the plate tectonic hypothesis coupled with the completion of regionalstratigraphic analysis of supracrustal rocks led to a new theory of Cordilleran evolution basedupon terrane accretion. This theory asserts that the Cordillera has grown 500 km westwardfrom the ancestral North American margin, which was formed by extension, rifting, andsea floor spreading in Early to Late Proterozoic time (1600-570 Ma). The scope of theprocesses involved in this tectonic growth include (i) the deposition of substantial thicknesses(>10 km) of sedimentary and volcanic rocks in a passive margin setting, beginning in theMiddle Proterozoic, (ii) the amalgamation and accretion of volcanic island arc and oceanicterranes, and associated capacious magmatism through Mesozoic and Cenozoic time, and(iii) very large transcurrent displacements during Cretaceous and Cenozoic time (Gabrielseand Yorath 1989).The term "terrane" refers to a block of the Earth's crust which holds a geological recorddistinct from those of surrounding terranes; the term has no significance pertaining to theorigin of a crustal block, or to its past or present position. By definition, the boundariesbetween adjacent terranes are faults. The description of terranes as "accreted" defines themas those which became a part of ancestral North America at an advanced point of their tectonicdevelopment. "Amalgamated" terranes, or "superterranes", consist of at least two terraneswhich were joined prior to their accretion to North America (Gabrielse and Yorath 1989).The map shown in Figure 1 delineates the terranes which make up the Canadian portionof the Cordillera. This present structure is considered to be the result of the accretion of twoallochthonous superterranes — the Insular superterrane and the Intermontane superterrane —FIGURE 1: Map of the Canadian Cordillera showing the accreted terranes in differentpatterns and labelled with circled letters (from Clowes 1989). The list of terranes is shownat the bottom left, with stars indicating terranes that were formed near the west coast of theancient craton; their relationships to North America are not certain. The other terranes arebelieved to have been formed on oceanic crust and then accreted to the ancient margin. Thelarger rectangle (solid line) shows the study area of the LITHOPROBE Southern CordilleraTransect; the inner rectangle (dotted line) shows the map area of Figure 2. Descriptions ofthe terranes discussed in the text are provided in Appendix A.3to the leading edge of the westward-travelling North American Plate prior to mid-Cretaceoustime (e.g. Monger and Price 1979, Coney et al. 1980). The Insular superterrane, whichwas positioned adjacent to ancestral North America, comprises the Alexander and Wrangelliaterranes, while the Intermontane superterrane to the east consists of the Stikine, Bridge River,Cache Creek, and Quesnel terranes. These terranes are briefly described in Appendix A.The Coast Belt is one of the five morphogeological belts of the Canadian Cordillera(inset, Figure 2); their boundaries are largely coincident with terrane boundaries (Gabrielseand Yorath, 1989). Numerous detailed explanations of the tectonic processes involved increating the Coast Belt through the accretion of geologic terranes have been provided (e.g.Davis et al. 1978, Monger et al. 1982, Gehrels and Saleeby 1985, and van der Heyden 1992).In a currently widely-accepted model, Monger et al. (1982) described the metamorphic andstructural character of the Coast Belt as the result of the mid-Cretaceous collision betweenthe Insular superterrane and the Intermontane superterrane, which was by then attached toNorth America. According to this hypothesis, the Insular-Intermontane collision caused theclosure of an intervening ocean basin, and the development of a collisional suture in itsplace. This suture was in turn overprinted by the granitoid intrusions of the Coast Belt,a magmatic arc linked to the east-dipping subduction beneath the new west coast of theNorth American continent. Most recently, van der Heyden (1992) has proposed a variationof this collisional model, whereby the Coast Belt was formed as an Andean-type magmaticarc superimposed upon a Cordilleran structure resulting from the mid-Jurassic accretion ofa superterrane composed of Wrangellia, Alexander and Stikinia to ancestral North America.Although the collision model of the Coast Belt orogeny (which is illustrated in Figure 3) isgenerally accepted, the nature of the superterrane collision zone has been obscured by themagmatic intrusions of the Coast Belt, as well as the other regional geological complexitiesmentioned below.4FIGURE 2: Location map showing shotpoints (triangles) and receiver locations (smallsquares) for Line 10 of SCoRE '90 and Lines 1, 2, and 3 of SCoRE '89. Other seismicsurveys in the region are shown by clashed lines; these include the refraction Line I of VISP'80, and the LITHOPROBE reflection profiles 88-13, 88-14, 88-17 and 88-18. The maparea is shown by the rectangle outlined with dotted lines in Figure 1. The inset shows themorphogeological belts of the Canadian Cordillera.5FIGURE 3: Schematic illustration of the collision model for the evolution of the Coast Belt:Cretaceous magmatic arc, with associated subduction west of the Insular superterrane (from vander Heyden 1992).61.2 Geological and Geophysical BackgroundThe Coast Belt runs for 1700 km between latitudes 49° N and 62° N, is 100-200 kmwide, and reaches elevations above 4 km (although most peaks are in the range of 2-3km) (Monger and Journeay 1992). It consists mostly of Late Jurassic to Early Tertiarygranites (80-85% by area according to Roddick 1983), and its dominant structural featuresare Mid-Cretaceous to Early Tertiary in age.Due to the vast amount of plutonism, along with the high-grade metamorphism andcomplex deformation that characterize the Coast Belt (CB), the orogenesis and early tectonicpast of this morphogeologic belt have proven difficult to discern.1.2.1 Geology of Southern Coast BeltThe area of interest for this study (see Line 10, Figures 2 and 4) lies in the southernCB between latitudes 49° N and 52° N, more than 200 km east of the Cascadia subductionzone which marks the east-clipping subduction currently taking place beneath the westerncontinental margin. Based upon terrane affiliation and the distribution of dated plutonicrocks (Friedman and Armstrong 1992), the southern CB can be divided into eastern andwestern parts, both of which were traversed by Line 10 (Figure 4).Along the eastern margin of the Insular superterrane, there lies a Late-Cretaceous, west-vergent contractional fault belt which exists in both the northern and southern CB (Mongerand Journeay 1992). The Coast Belt Thrust System (CBTS) of southwestern British Columbiahas been identified as a part of this contractional system (Journeay and Friedman, in press);Line 10 passed through the CBTS. On Figure 4 four major faults of the CBTS in the regionof the present study have been identified.The southwestern CB (WCB, Figure 4) is characterized by low grade metamorphism anda mainly plutonic nature, and therefore probably acted as a rigid crustal block during LateCretaceous west-directed thrust faulting (Journeay and Friedman, in press). To the east, the50' - 50'-121'-124' -123' -122'-125'49'-121'-123'-124' -122'-125'7ASCoRE '90 shot point- 52'orphic AssemblagesAssemblageAccretiones52' 011^Plu o•Intrusions=Undivided MetaGambier OverlaOther Post-TerraneOverlap AssemblaCadwallader T: rraneStikine TerraneBridge River Terr neShuksan Terranearrison TerraneRocks ost-Terrane AccretionNTERBONTANET5V - 51'\100 ilometers49';:Z5fig!FIGURE 4: Location map showing shotpoints (triangles) and receiver locations (small squares) forLine 10. Heavy lines indicate the divisions between the Insular, Coast, and Intermontane belts, andlines show the positions of the major faults of the CBTS (Coast Belt Thrust System) in the region(CCBD denotes the Central Coast Belt Detachment). The Harrison Fault divides the Eastern CoastBelt (ECB) from the Western Coast Belt (WCB) (R. Friedman, pers. comm.). The regional geologyis shown in the area surrounded by the dashed line (Wheeler et al. 1991).8WCB is thought to have been underplated beneath the southeastern CB (ECB) along west-vergent thrust faults of the CBTS (Monger and Journeay 1992) during tectonic adjustmentslinked with the impingement of the Insular superterrane in Late Cretaceous time (Journeay1990). The ECB (Figure 4) is characterized by the presence of oceanic Bridge River terranerocks. The rocks range widely in age from at least 160 to 350 million years, and havetherefore been identified as a former large ocean basin (Monger and Journeay 1992).The Line 10 profile traversed mainly granitic plutons created during post-terrane accretionintrusion, as well as slivers of terrane materials and some rocks excluded from terraneclassification, as shown in the region in Figure 4 enclosed with a dashed line. In thenorthwest, the Stikine and Cadwallader terranes were crossed, as were the Gambier overlapassemblage (comprising arc volcanics and local rift volcanics) and another undifferentiatedoverlap assemblage (characterized by foredeep and southwesterly-derived clastic wedgematerial, foredeep marine shales, and arc volcanics). In the center and southeast, whereit runs along the Harrison fault, Line 10 passed over pods of the Cadwallader terrane, theGambier assemblage, and undivided metamorphic assemblages (Figure 4). The terranes ofthe CB are described briefly in Appendix A.1.2.2 Previous Geophysical StudiesIn comparison with the other belts of the Canadian Cordillera, the Coast Belt (CB) hasuntil recently been a frontier area in terms of geophysical studies. The pioneer seismicresearch effort in the CB was the interpretation of Berry and Forsyth (1975) of refractiondata recorded prior to 1971. However, as the data were acquired with large shot and receiverspacings, the resolution of both lateral velocity variations and changes in velocity with depthwere poor. The velocity model in the vicinity of Line 10 that resulted from this study has acrustal velocity of 6.0 km/s to a depth of 20 km, a lower crust with velocities of 6.1 km/sat the top and 6.9 km/s at the bottom, an upper mantle velocity of 7.8 km/s, and a crust-mantle boundary at 33 km depth. A subsequent surface-wave study by Wickens (1977) in9the southern Canadian Cordillera that was generally consistent with that of Berry and Forsyth(1975) lacked data in the center of the southern CB.Although the Berry and Forsyth (1975) experiment was the only seismic survey coin-cident with the present area, previous refraction measurements for the north (Johnson et al.1972) as well as the east (Cumming et al. 1979) and west (McMechan and Spence 1983,Spence et al. 1985, and Drew and Clowes 1990) have also been interpreted. The two for-mer studies were not important to the current one as the Johnson et al. (1972) effort wasrestricted by the limitations described for the Berry and Forsyth (1975) case, and profilesrecorded similarly to Line 10 have yielded updated results for the region to the east (Zeltet al. 1992a, 1992b).In the context of linking known features of the subsurface region to the west with thestructure determined from Line 10, it is important to consider the Spence et al. (1985)interpretation of the refraction data from a refraction line (VISP '80 Line I, see Figure 2)which ran from Vancouver Island to the western CB. Principal features of the model ofSpence et al. (1985) include upper crustal velocities of 6.4 km/s to 6.7 km/s (to a depth of18 km), and a Moho depth of 38 km; however, the model is considered well-constrainedonly to 18 km depth.Other types of geophysical studies performed in the region of interest include geo-magnetic (Caner 1970), magnetic (Coles and Currie 1977), and gravity (Riddihough 1979)surveys. The resolution of the structural information provided by these potential field datais limited in comparison with the detail which is of interest for the present study.The SCoRE '89 and '90 (Southern Cordillera Refraction Experiment) projects whichencompassed Line 10 were run as part of the LITHOPROBE Southern Cordillera Transect(Clowes 1989, Clowes et al. 1992). LITHOPROBE is a multidisciplinary program designedto supplement the existing geological and geophysical data, in order to elucidate the processesinvolved in the westward growth of the North American continent. The profiles whichiocoincide with the present study are reflection profiles 88-13, 88-14, and 88-17, and refractionLines 1, 2, and 3 of SCoRE '89 (see Figure 2). Heat flow (Lewis et al. 1988, 1992)and magnetotelluric (Jones et al. 1992) studies which formed part of the LITHOPROBESouthern Cordillera Transect activities were also undertaken. The elements of these andthe aforementioned surveys which are of particular relevance will be discussed within theframework of the results determined in this thesis.1.3 Experimental ObjectivesSCoRE '89 and '90 were designed to complement LITHOPROBE seismic reflection dataacquired in 1985 and 1988 (Cook et al. 1988, 1991, 1992; Varsek et al. 1992) as they offereda means of extending the reflection results beyond 2—D interpretations. This would fulfillthe general objective of providing structural models for the southern Canadian Cordillera totie the models which exist for the convergent margin (e.g. McMechan and Spence 1983,Spence et al. 1985, Drew and Clowes 1990) with the interior of the continent.More specifically, the ultimate objectives put forth for the SCoRE '89 and '90 surveyswere to study the velocity field, crust and upper mantle structure, spatial range of terranes,fault geometry and westernmost extent of cratonic basement in the southern CanadianCordillera. This was to be achieved by the integrated interpretation with the seismic reflectiondata and other geophysical and geological studies.The specific objective for the present study is to analyse the refraction data of the SCoRE'90 Line 10 profile in order to map crust-mantle transition and intracrustal boundaries, and tointegrate this interpretation with coincident seismic reflection and refraction data, gravity andheat flow profiles, and magnetotelluric information in order to unravel some of the mysteriessurrounding the formation of the Coast Belt.11II. DATA ACQUISITION AND PROCESSING2.1 The Refraction Experiment2.1.1 LITHOPROBE Southern Cordillera TransectAs part of LITHOPROBE's multidisciplinary scientific studies in the Southern CordilleraTransect, the extensive seismic refraction program referred to as SCoRE was can -ied outduring the summers of 1989 and 1990. Line 10, the focus of this study, was recorded asone part of SCoRE ' SCoRE '89 & SCoRE '90The first phase of the refraction experiment, SCoRE '89 (Figure 2), focussed on theIntermontane Belt and the Coast Belt, with an incursion into the eastern Insular Belt to providea closer tie with the VISP '80 experiment (Zelt et al. 1992a). Two shots also were detonatedwest of Vancouver Island to extend the VISP refraction Line I across the CB. SCoRE '89 wasan international collaboration with participants from Canadian universities, the GeologicalSurvey of Canada, the United States Geological Survey, and Japan. It represented the firstthorough application of the spatial seismic refraction recording (S 2R2) method (Kanasewichand Chiu 1985). By recording shots from each apex of the triangular array, the center ofeach arm, and the center of the array, a set of inline and broadside data was obtained. Inlineprofiles are interpreted as 2—D crustal sections, whereas the suite of broadside data provide thebasis for a tomographic study of the interior of the triangle, in this case focussed on the FraserRiver fault system, the boundary zone between the Coast Belt and the Intermontane Belt.The second phase of the experiment, SCoRE '90, completed studies in the CB includingLine 10, and extended coverage eastward. This also was a cooperative, international effort,involving the University of Alberta, the University of British Columbia, the GeologicalSurvey of Canada, the U.S. Geological Survey , the University of Victoria, and the University12of Saskatchewan. It spanned a major portion of the Canadian Cordillera south of 52° Nlatitude. Line 10, which was approximately 350 km long and ran entirely within the CB,was the westernmost profile. The remainder of the project included east-west Lines 7 (400km long) and 9 (360 km long) which spanned the Intermontane, Omineca Crystalline, andForeland Belts of the Cordillera; north-south Line 8 (340 km long) which formed an along-strike profile in the Omineca Belt; and Line 6 (230 km long) located just north of the USborder and deployed to obtain broadside recording of shots along Line 8 (Burianyk et al.1992).2.1.3 Line 10 — Acquisition Geometry of the Field DataLine 10 extended northwest from Harrison Lake in the western CB to north of ChilcoLake in the eastern CB (Figure 2). Based on the experimental plans of shot points at -,50km intervals and the limitations of access imposed by the available road system, Line 10included 6 shot points (SP) (Figures 2 and 4). Large shots were detonated at the endpoints:2500 kg of explosives were used at SP-46, and 1800 kg at SP-51. An 800 kg shot wasfired at SP-49 near the center of the line, while smaller sources were used at the interveningshotpoint locations (400 kg at SP-50, 200 kg at SP-48 and SP-47) in order to providedetailed information on shallow crustal structure.Three hundred portable seismographs were placed at intervals varying from 1 to 3 kmalong the profile. The receiver locations were determined by marking off the intervals on astraight line drawn between the shot points on a map, and then projecting these sites ontothe available road system for ready access (Burianyk et al. 1992). Line 10 representedthe only possible strike line within the CB along which such a survey could be carried outat reasonable cost. Nevertheless, many stations northwest of SP-49 and SP-50 had to bedeployed by helicopter due to a lack of roads, and those along Chilco Lake were deployedby boat.132.1.4 Shot and Receiver Site PositioningAll of the shot points and 94% of the receiver locations were located using twotypes of GPS (Global Positioning System) navigational instruments: the Trimble AdvancedNavigation Sensor GPS, and the Magellan GPS NAV 100 PRO. Two-dimensional fixes (i.e.3 satellites used in the positioning solution) were used to position most coordinates, as this isthe most accurate approach where altitude is known. Altitudes were read from topographicmaps with a scale of 1:50 000, and the resulting measurement errors were approximately100 m for latitude and longitude, and 20 m for altitude.When poor satellite visibility made GPS positioning difficult, the alternative was to readthe locations from the 1:50 000 topographic maps, which resulted in errors of 100-150for latitude and longitude (B. Isbell, personal communication). Due to a limited supply ofinstruments, Loran C positioning was used for a number of helicopter deployments at thenortheastern portion of the line; these were then corrected using redundant GPS measurementsso that accuracy levels were comparable to that of GPS (Burianyk et al. 1992).2.2 InstrumentationThe recording instruments consisted of 180 GSC/LITHOPROBE EDA model PRS-1 digital refraction seismographs and 120 USGS analog cassette recorders, which wereconnected to Mark Products L-4A 2 Hz vertical seismometers. The sampling rate for thePRS-1 s was 120 samples/s, while the analog USGS data were digitized at 200 samples/s(Burianyk et al. 1992). For convenience of discussion, the data sets recorded on the PRS-lsand the analog systems will be referred to as the GSC and USGS data, respectively. Therelative velocity responses of these instruments have been given by Zelt et al. (1990).Each shot point had one to three drill holes holding charges of Nitropel®, a nitroglycerine-based pellet explosive, at depths of 30-42 m. The detonation time of the charges was con-trolled by a clock device which uses the same satellite-based timing system (GOES) as the14recording instruments. The precision involved is such that instrument timing errors are con-sidered negligible (I. Asudeh, personal communication); in other words, the beginning ofeach recorded trace is assumed to be the exact moment of shot detonation.2.3 Processing2.3.1 Initial Data ReductionDigital field data were initially available on 9—track tape in SEGY-LDS format (Spenceret al. 1989), with the GSC and the USGS data forming separate data sets. The first step inprocessing was conversion of the data recorded by the USGS instruments to the standardizedformat of the GSC data. The data were redigitized to 120 Hz from 200 Hz, gains werecorrected to units of nanometers/s, and reformatting from SEGY format 3 (i.e. 2 byte fixedpoint) to SEGY format 1 (4 byte floating point) was carried out. Finally, 20 seconds ofdata (sample values of zero) were added to make the USGS traces equal in length to theGSC traces.At this stage, the headers of both the GSC and USGS traces were edited to update shotand receiver site locations. Corrections were therefore necessary to dependent header valuessuch as source-receiver offset.This initial data reduction was done by the University of Alberta Seismology Laboratory(Burianyk et al. 1992). Additional processing was designed to extract the most informationpossible from the data.2.3.2 FilteringSpectral analyses were performed on representative signal and noise recordings from therefraction data. This testing revealed that the power spectrum of the desired signal was withinthe 2-12 Hz range, whereas that of the noise was mainly outside of this range. An eight-pole Butterworth bandpass filter (Kanasewich 1981) was applied to the data in an attempt15to improve the signal-to-noise ratio. These filtered sections were used concurrently with therespective unfiltered sections to make travel time picks.2.3.3 Amplitude CorrectionsAmplitude corrections were applied to the data to account for spherical spreading Eachtrace was scaled by a factor proportional to d15 , where d is the shot-receiver distance. Thisprocess enhanced the amplitude of the far-offset traces with respect to those near the source.There was one exception to this procedure: since the amount of energy released at SP-46 wasso large, it was necessary to scale the corresponding traces proportionally to d 1 •. Otherwise,the amplitudes of the near-offset traces would be so large that separate phases could notbe distinguished. There was no amplitude variation observed between traces recorded bythe two types of seismographs, and since the velocity responses of each were nominally thesame, no further amplitude corrections were necessary.The final step in preparing the data for interpretation was singling out traces which wereso noisy that they obscured the signal visible on adjacent traces when the data were displayedwith the amplitude corrections described above. For the purpose of these displays (Figures5a-10a), these traces were deleted (see Table 1); this was not necessary when the data wereplotted in common maximum amplitude format (Figures 5d-10d).16Shot Point Trace Numbers of Deleted TracesSP-46 147-151,161SP-47 55,58,69,75,78,97,147-151,161,193SP-48 55,58,69,91,97,147-151,161SP-49 198,247SP-50 147-151,161,178SP-51 135,147-151,161,178,189Table 1. The trace numbers of those traces which had to be deleted (due to their excessivenoisiness) in Figures 5a — 10a are specified for each shot point. Note that the traces arenumbered 1-290 from left to right for each shot point in the data sections shown in the figures.17III. DATA ANALYSIS AND INTERPRETATION PROCEDURES3.1 Initial Observations and Interpretations3.1.1 Data Characteristics and QualityAcquisition of the field data resulted in a data set of 1734 traces. The six record sectionsare displayed in Figures 5a-10a. Of these 1734 traces, approximately 1100 contain primaryand secondary arrivals which could be identified and timed. Distances in Figures 5-10 referto model distances, with zero distance coinciding with SP-51. References to shot pointsaccompanied by a "NW" or "SE" refer to data recorded to the northwest or southeast ofthe shot point. In the text, offset distances will be mentioned; these refer to distances fromthe shots to the receivers.The branches of the phases that were identified and used in the modelling procedureare labelled in Figures 5a-10a; in addition, Figure 5e shows the two phases which are notobservable on a relative amplitude section unless the traces of SP-46 are scaled proportion-ately to d1.5 . In total, eight phases were identified: Ps, Pg, and Pn, which represent energywhich refracted back to the surface through the near-surface, upper crust, and upper mantle,respectively; and R1, R2, R3, PmP, and R, which denote energy which reflected back tothe surface. The R1 reflector is at the base of the upper crust, R2 is at the base of themid-crust, R3 is within the lower crust, and PmP corresponds to energy reflected from thecrust-mantle boundary (see Figures 5a-10a). R is a sub-Moho reflection phase visible only onthe far-offset traces of the SP-46 section when the appropriate amplitude scaling is applied(see Figures 5d and 5e).Due to low signal-to-noise ratio, inadequate shot-receiver offset, or nonexistence of aparticular phase, not all phases are observed on all record sections. The quality of the datais usually good, with the following exceptions: SP-47 NW (50-110 km offset, Figure 6a),and SP-50 SE (120-180 km and 215-255 km offset, Figure 9a). Table 2 gives a synopsis of18 SE0 50 100DISTANCE (kg/1111■1/1111M, AI^:11.111,150^200 250 300 350it.;,!.a Lr.-:..2. „:I 1 .^g■^I ,PrilLIIIIIMI 11011131111 ilia..mairimmatimmolun minI1E'1 pirgiliiti,1•011til  1.1116111411111114101111111M11111113IGNIGANNAIMIRIMIC15:111111MINIffeffloimiliminummtorast121 0(f) 00 50 100 150 200 250 300 35021FIGURE 6:(a) Observed record section for SP-47 (see Figure 5 caption for plotting parameters andgeneral description). (b) Synthetic section for SP-47 (see Figure 5 caption for scaling and wavelet).0^50^100^150^200^250^300^350DISTANCE (km)FIGURE 6:(c) Comparison of observed and calculated travel times (see Figure 5 caption for plottingstyle). (d) Observed record sections scaled to common maximum amplitude.'110^50^100^150^200^250^300^350DISTANCE (km)FIGURE 7:(a) Observed record section for SP-48 (see Figure 5 caption for plotting parameters andgeneral description). (b) Synthetic section for SP-48 (see Figure 5 caption for scaling and wavelet).FIGURE 7:(c) Comparison of observed and calculated travel times (see Figure 5 caption for plottingstyle). (d) Observed record sections scaled to common maximum amplitude.L/I \J I rlisIN-11-- ‘1\111/FIGURE 8:(a) Observed record section for SP-49 (see Figure 5 caption for plotting parameters andgeneral description). (b) Synthetic section for SP-49 (see Figure 5 caption for scaling and wavelet).26.^.FIGURE 8:(c) Comparison of observed and calculated travel times (see Figure 5 caption for plottingstyle). (d) Observed record sections scaled to common maximum amplitude.27FIGURE 9:(a) Observed record section for SP-50 (see Figure 5 caption for plotting parameters andgeneral description). (b) Synthetic section for SP-50 (see Figure 5 caption for scaling and wavelet).28.^.FIGURE 9: (c) Comparison of observed and calculated travel times (see Figure 5 caption for plottingstyle). (d) Observed record sections scaled to common maximum amplitude.29.^.FIGURE 10:(a)Observed record section for SP-51 (see Figure 5 caption for plotting parameters andgeneral description). (b) Synthetic section for SP-51 (see Figure 5 caption for scaling and wavelet).30FIGURE 10:(c) Comparison of observed and calculated travel times (see Figure 5 caption for plottingstyle). (d) Observed record sections scaled to common maximum amplitude.31Shot point Ps & R1 R2 R3 PmP Pn R TotalPgSP46 SE 1 1(50) (50)SP46 NW 101 12 21 62 33 14 243(63) (144) (150) (143) (92) (93) (100)SP47 SE 48 48(52) (52)SP47 NW 42 14 12 35 105(72) (111) (108) (115) (96)SP48 SE 64 40 6 18 128(68) (100) (153) (120) (89)SP48 NW 60 20 23 71 174(93) (95) (104) (109) (101)SP49 SE 73 14 16 38 141(81) (129) (78) (129) (98)SP49 NW 129 10 41 55 235(68) (127) (120) (100) (87)SP50 SE 65 30 23 28 146(61) (88) (114) (96) (82)SP50 NW 54 7 61(56) (62) (57)SP51 SE 108 56 15 43 222(50) (92) (80) (97) (78)SP51 NW 2 2(50) (50)Total 749 145 108 107 350 33 14 1506(67) (96) (124) (103) (114) (92) (93) (88)Table 2. The number of observations (travel time picks) for each phase and each shot point.Corresponding average uncertainties in ms are given in brackets. The right column givesthe total number of observations (and corresponding average uncertainties) for each shotpoint; the bottom row lists the total number of observations (and corresponding averageuncertainties) for each phase.32Ps is visible on all record sections, except on SP-50 NW (Figure 9a). This near-offset portion of SP-50 NW was the least successfully modelled feature of the data, mainlybecause it was difficult to concurrently model the relatively prominent Ps phase on SP-50SE. From section to section, the Ps phase shows large variations of approximately 2 to 5km/s in apparent velocity (i.e. the inverse of the slope of the travel time-distance curve).The amplitude of Ps is invariably high, which indicates a large velocity gradient in thenear-surface layer.The Pg phase appears on all sections, after the Pg-Ps cross-over distance, as the firstarrival to offsets of 110 — 170 km; its average apparent velocity is 6.2 km/s. The phase isessentially linear with no distinct cross-over points, which indicates that there are no distinctlayers of different velocities and no strong vertical velocity gradients in the upper crust. Aswell, there are no substantial travel time offsets to indicate the presence of thick low velocityzones. The fluctuations which appear within the branch were considered to be caused bynear-surface features such as topography and relatively rapid changes in velocity throughoutthe near-surface layer.With the exception of SP-46, the amplitude pattern of Pg varies minimally from section tosection, which suggests that there are no major lateral changes in the velocity gradient withinthe upper crust. As shown in Figure 5a, however, the amplitude of the phase is exceptionallylarge for SP-46, indicating a strong gradient in this region. Within each coherent Pg branch,the amplitude decreases quite abruptly at certain offsets (ranging from ,30 km on SP-50SE to ,-100 km on SP-46 NW; see Figures 9a and 5a). If near-surface effects can be ruledout, such fluctuations indicate a sudden decrease in velocity gradient at some depth withinthe upper crust.The Pn phase appears as a first arrival on only one section — SP-46, where the largestshot was fired. The Pg-Pn crossover point occurs at an offset of about 180 km and the apparentvelocity of Pn is ,8.0 km/s. This Pn phase is also essentially linear, with fluctuations being33attributed to near-surface variations as in the Pg case.The secondary arrivals from the upper crust, identified as R1, are apparent on 5 recordsections, at the source-receiver offsets shown in Table 3. The R1 arrivals on SP-51 SE,SP-50 NW and SE, and SP-48 SE are virtually parallel to the Pg first arrivals, with anaverage delay of 0.30 s; the R1 arrivals on SP-46 are visible at offsets of 140-200 km,where the Pg arrivals are obscured by noise. The R1 phase is generally prominent — itsamplitude is greater than that of the Pg phase in all cases. The amplitudes of R1 were usedto constrain the velocities at the top of the middle crust.The phase labelled R2 was identified as a wide-angle reflection from the base of themiddle crust and is visible on 7 sections; see Table 3 for the source-receiver offset rangesat which R2 appears. The asymptotic apparent velocity of the R2 phase is ,6.40 km. Theamplitude of the R2 arrivals is invariably smaller that of the R1 arrivals. The R2 phase wasthe principal constraint for the velocities of the middle crust and the amplitudes of the R2arrivals were used to determine the velocity contrast across the middle crust — lower crustinterface. However, it was not possible to obtain any information from the R2 amplitudesof SP-47 NW as this phase is visible on sections plotted in common maximum amplitudeformat (Figure 6d), but cannot be accurately discerned on true amplitude plots (Figure 6a).The secondary arrivals from the lower crust consist of two phases: R3, identified as awide-angle reflection from a horizon within the lower crust, and PmP, identified as a reflectionfrom the crust-mantle boundary. The R3 reflections appear on the seismic sections from 5shot points. The source-receiver offsets at which the respective R3 branches are observedare shown in Table 3. The R3 phase typically appears 0.25-0.40 s before the PmP arrivals(except for SP-49 NW, Figure 8a, where the R3 branch is ,0.8 s before PmP at offsetsbetween 50 and 100 km). The asymptotic velocity of the R3 phase is between ,6.40 km/s34Phase^Shot Point^Source-receiver OffsetR1 SP-46 NW 140 -200 kmSP-48 SE^60 -110 kmSP-50 NW 60 - 110 kmSP-50 SE^60 -30 kmSP-51 SE 80 - 150 kmR2^SP-46 NW^185 - 220 kmSP-47 NW 100 -140 kmSP-48 NW^100 - 115 kmSP-48 SE 70 - 100 kmSP-49 NW^150 - 175 kmSP-49 SE 75 - 100 kmSP-50 SE^80 - 105 km/180 - 210 kmR3^SP-47 NW 140 - 160 kmSP-48 NW^120 - 170 kmSP-49 NW 150 - 175 kmSP-49 SE^50 - 170 kmSP-51 SE 110 -135 kmPmP^SP-46 NW^155 - 275 kmSP-47 NW 120 - 180 kmSP-48 NW^110 - 230 kmSP-48 SE 80 - 115 kmSP-49 NW^90 - 170 kmSP-49 SE 100 - 175 kmSP-50 SE^85 - 100 km/180 210 kmSP-51 SE 70 - 170 kmTable 3. The source-receiver offsets at which each phase was observed is specified for therelevant shot points.35and --,6.50 km/s in each case.The amplitudes of the R3 phase vary: for SP-51, it is larger than R1 and quite prominentover the small distance that it is apparent; for SP-49, the R3 phase is less prominent thanits SP-51 counterpart, but still larger in amplitude than the earlier R2 reflection; and theR3 phase visible on the seismic section from SP-48 is comparable in amplitude to the R2phase, as was assumed to be the case for R3 of SP-47. The overall trend is decreasing R3amplitudes from NW-SE. The R3 phase was used to constrain the velocities in the upperportion of the lower crust, and its amplitudes were used to constrain the velocities belowthe R3 reflecting horizon.The very prominent PmP arrivals are visible on 8 record sections at the offsets shownin Table 3. The asymptotic apparent velocities average 6.7 km/s. The source-receiver offsetand reduced arrival time near each critical point — values which are sensitive to the Mohodepth and the overlying velocity structure — vary significantly, as shown in Table 4. Thesevariations indicate that lateral velocity changes in the base of the lower crust, topographicvariations along the crust-mantle boundary, and lateral changes in the velocity contrast acrossthe Moho are likely to be present. In terms of amplitudes, SP-48 NW is anomalous as it isnot relatively large — this may be a further indication of significant lateral complexity in theMoho. The PmP arrivals were used to constrain the velocity at the base of the lower crustand the topography of the Moho, and their amplitudes were incorporated to give informationabout the velocity contrast across the crust-mantle boundary.The wide-angle reflection R is present on the seismic section from SP 46, starting at315 km shot-receiver offset and 8.0 s reduced travel time (see Figures 5d and 5e). Theasymptotic apparent velocity of this feature is > 8.0 km/s, and its amplitude is bigger thanPn and smaller than the other reflections. These features of the phase led it to be identifiedas a sub-Moho reflection.36Shot Point^Offset of Critical Point (km)^Reduced Arrival Time ofCritical Point (s)SP-51 SE^80^ 6.9SP-50 SE 100 6.9SP-49 NW^100 7.2SP-49 SE 100^ 6.5SP-48 NW^110 7.0SP-48 SE 165 7.0SP-47 NW^80^ 6.8SP-46 NW 105 7.6Table 4. The source-receiver offsets and reduced arrival times of the critical point for PmPis given for each appropriate shot point.3.1.2 Travel Time PickingIdentifying and picking travel times of the phases involved looking at plots (of variousscales) as well as computer screen presentations using interactive software. The procedurefollowed for travel time picking was dependent upon the nature of the phase, as outlinedbelow; the preferred trace format was wiggles in all cases.First, the Ps and Pg first arrivals were picked from large-scale presentations of unfilteredtraces using the interactive software package. This was a relatively straightforward task dueto the high quality of the data and the uniform nature of the waveforms. The Pn first arrivalswere picked in the same way. Each travel time pick was assigned an uncertainty dependentupon the difficulty of pinpointing the first motions. Ranging between 25 and 150 ms, theseuncertainty values were used in the interpretation procedure (see below). The average errorsfor Ps, Pg, and Pn, as well as all other phases, are shown in Table 2.The picking of most of the R1 travel times was accomplished in the same manner. Inregions where the onsets of the arrivals were difficult to discern in comparison with Pg, itproved useful to make initial picks on large-scale displays of the data plotted at a reducing37velocity of 6 km/s, and then make fine adjustments using the interactive software. The traveltimes for R2 were also picked in this manner for sections from all shotpoints except SP-47, forwhich it was necessary to use a filtered section due to high noise levels. For the R3 picking,the method used was the same as that of R1, except that the plots were displayed using an8 km/s reducing velocity, in order to increase the coherency of this higher-velocity phase.For the PmP arrivals, the picking procedure was complicated by variations in the waveform.Plotting the sections in common maximum amplitude format (Figures 5d-10d) was helpfulfor the picking of the PmP phase, as this sometimes enhanced the phase coherency. Commonmaximum amplitude plots were also useful in the picking of the upper mantle arrival R whichappeared on SP-46 NW, as seen by comparing Figure 5d with Figures 5a and 5e.3.2 Interpretation Methods3.2.1 The Modelling AlgorithmsThe modelling procedure applied to the refraction data combined two approaches:1. the iterative modelling algorithm of Zelt and Smith (1991) was used to create a modelof interfaces and velocities via the inversion of travel time data;2. the forward amplitude modelling algorithm of Zelt and Ellis (1988) was then used toadjust the model based upon the additional constraints of amplitude information.Both the travel time algorithm and the forward amplitude modelling algorithm are 2-Dmethods based on the numerical solution of the ray-tracing equations of zero-order asymptoticray theory (terven5r et al. 1977), with first-order asymptotic ray theory incorporated for thecalculation of head-wave amplitudes. The numerical evaluation of the travel time along a raypath uses the Runge-Kutta method (see, for example, Sheriff and Geldart 1982) with errorcontrol as suggested by terveq et al. (1977). The application of Snell's law at each pointwhere a ray intersects a boundary was the final aspect of the ray tracing procedure.38An inversion method for the travel time modelling was chosen over a forward travel timemodelling procedure primarily because it is a less time-consuming, more objective approach.In addition, the Zelt and Smith (1991) routine provides useful estimates for the resolutionof model parameters (i.e. boundary depths and layer velocities). Furthermore, the choseninversion program was written for crustal seismic data with a shot-receiver density too lowto fulfill the requirements of conventional tomographic or full wavefield techniques. In thissense, the algorithm is ideally suited to the Line 10 data.The Zelt and Ellis (1988) algorithm was chosen because it efficiently allows rapid forwardmodelling of amplitudes; although more accurate methods (e.g. finite-difference algorithms)exist, they were not deemed computationally practical for this study. An important aspectof this algorithm is that in order to calculate synthetic seismograms, a value of the physicalconstant Q (the "quality factor") must be assumed in order to account for attenuation, thereduction in amplitude or energy caused by the characteristics of the transmitting media. Thevalues of Q, which is defined as 27/(fraction of energy lost per cycle) (Sheriff and Geldart1982), used to calculate the amplitudes of the synthetic seismograms were not determinedthrough the modelling procedure, but were adapted from the values of Hough and Anderson(1988) determined for the Anza, California region. That study area was very similar to theCB site, as it was situated within the Southern California Batholith. The dominant geology(granitic plutons), age (Jura-Cretaceous), and tectonic origin (subduction off the west coastof North America) of the Southern California Batholith make it a suitable analog for thiscomparison with the CB. The Q structure used in modelling is shown in Table 5.One convenience of using the travel time inversion and the forward modelling algorithmsconcurrently was that the model parameterizations for the two methods were identical. Theparameterization procedure is based upon the positioning of two types of nodes: boundaryand velocity nodes. The boundary nodes are essentially depth nodes, as their horizontalposition is user-specified, with the inversion program accounting for the adjustment of the39Depth From Surface^ Q (dimensionless)1.5 km 1008.0 km 30022.0 km^ 100072.0 km 700Table 5. The Q structure used in modelling is shown in terms of depth from the surface.Notes: each Q represents the value assigned to a homogeneous Q layer with its bottom atthe corresponding depth and its top at the previously-mentioned depth. Values were adaptedfrom Hough and Anderson (1988).depth. The boundary nodes are used to construct the model in terms of layers; a layerboundary is represented by an arbitrary number and spacing of boundary nodes which areconnected by linear interpolation. These horizons must cross the model from left to rightwithout crossing another boundary.The velocity structure within each layer is represented by a series of velocity nodes atthe top and at the bottom of each layer. Again the number and spacing of nodes is arbitrary.The velocity between the upper velocity nodes in a layer is determined through horizontallinear interpolation, as is the velocity between the lower velocity nodes. The final pictureof the layer velocity is achieved by linearly interpolating vertically between the upper andlower velocities in each layer (see Zelt and Smith 1991).The resulting parameterized model is made up of layers composed of different-sizedtrapezoids which are divided by vertical segments that are included wherever there is aboundary or velocity node. These vertical segments within each layer do not representvelocity discontinuities — the velocity is laterally continuous within each layer. The quasi-horizontal boundaries between layers can represent velocity discontinuities if required; theycan also represent second-order discontinuities (i.e. a change in velocity gradient), or pseudo-boundaries used for ease of modelling. A smoothing procedure is applied to reduce the effects40(i.e. scattering and focusing of ray paths and geometrical shadows) of this blocky modelparameterization (Zelt and Smith 1991).Since the structure of the final model is strongly dependent upon the nature of theinversion algorithm, the aspects of the inversion procedure that are important in this contextare discussed below. This derivation of the algorithm's essential equations is an alternativeapproach to the one taken by Zelt and Smith (1991), who solved equation (1) using adamped least-squares technique presented by Aid and Richards (1980), following the methoddescribed by Lutter et al. (1990).3.2.2 Inversion mathematicsIf a velocity field is considered in terms of its discrete form instead of its form as acontinuous velocity field v (x, y), the travel time t along a raypath between a source anda receiver ist= E li/viwhere vi is the velocity of the i th ray segment, li is the length of the ith ray segment, andn is the number of segments comprising the ray path (i.e. the number of cells the raysamples in passing through the velocity field). This travel time expression is linear in termsof slowness (reciprocal velocity), but non-linear in terms of raypath since the path of the raydepends upon the velocity structure it encounters. The following linearized expression fortravel time residual (the difference between observed and calculated travel time) is obtainedby expanding the velocity field about a starting velocity model in a Taylor's series expansionand neglecting the higher order terms:At = AAm^ (1)where At is the travel time residual vector, Am is the model parameter adjustment vector,and A is a matrix of partial derivatives. A contains the elements ati/ami, where ti is the41i th observed travel time, and mi is the j th model parameter selected for inversion (i.e. eithera velocity value or the z-coordinate of a boundary node).To solve equation (1), rays are traced through the current model during a particulariteration and At and A are computed, with partial derivatives calculated analytically usingthe method described in Zelt and Smith (1991). This determines the model update vectorAm, which is added to the current model for further iterations or for stopping if the stoppingcriteria are met.Crustal refraction experiments often give rise to a system (1) which is overdetermined,underconstrained, and ill-conditioned. For such systems, the damped least squares methodof solution is a common choice (e.g. Braile 1973, Kanasewich and Chiu 1985, Chiu andStewart 1987, Phadke and Kanasewich 1990). In solving the linearized data equations (1),Zelt and Smith (1991) considered the overdetermined case by choosing to solve for fewermodel parameters than travel time observations, and they sought a model update vector tosatisfy the system of equations in the least-squares sense via a damped least-squares solution.The relevant objective function needs to be of the form:0 (Am) = II AAm — Atir -I- 62 IIAmirwhere 62 is a scalar trade-off parameter that is specified in order to set the relevant importanceof the terms, and 11.II means the 12 norm. Extremizing 0 with respect to Am yields the linearalgebraic equationsAm= (ATA -Fe2I) ATAt^ (2)Following the method suggested by Lutter et al. (1990), if the trade-off parameter isdefined as E2 = Do-? /0-2 where= the standard deviation used to estimate the uncertainty of the i th travel timemeasurement,0- • = the standard deviation used as an a priori estimate of the uncertainty of the i thmodel parameter,42D = an overall damping parameter usually equal to one,then this becomes the Zelt and Smith (1991) damped least-squares solution to (1):Am= (A TC 1 A-1-DC;„;1 ) -1 A T C l At^ (3)whereC t = the estimated data covariance matrix defined by diag (al ), andCm = the estimated model covariance matrix defined by diag fq ).The values of the standard deviation Cri were specified during travel time picking by assigningan error estimate to each observation, as mentioned in the previous section. The choice ofthe values used for the standard deviation a • will be included in the forthcoming discussionof the modelling procedure.It is important to observe that by removing some of the underdeterminancy of the problemby underparameterizing, the inverse solution becomes more stable, but also more dependenton the chosen model parameterization. It is therefore of utmost importance that a modelparameterization appropriate to the subsurface geology is chosen (i.e. within the capacityof blocky trapezoids to represent the subsurface). The optimal parameterization uses themaximum number of nodes that the data are able to constrain, for if these nodes are prudentlypositioned they promise enhanced knowledge about the subsurface.Zelt and Smith (1991) suggested that the optimal node placement may be determined byrunning a series of tests based on the following criteria:1. the travel time residual must be sufficiently small;2. parameter resolution (as determined by the resolution matrix R) must be adequatelyhigh; and3. rays must be traced to all observations.43Using these three criteria, a preliminary series of inversions was performed to attempt tooptimize the number and placement of the velocity and boundary nodes used in inversion.In practice, this was possible for the boundary nodes; for the velocity nodes, however, whilethe number of nodes chosen could be optimized, in certain cases their placement could not.In these instances, which are discussed later, the choice of velocity node positioning asgoverned by the above guidelines was somewhat arbitrary.The resolution matrix R defined by Zelt and Smith (1991),R. (A TCi-1A+DC;,I) -1 A TCT 1 A(derived in Lutter et al. 1990), has as its diagonal elements numbers that range betweenvalues of 0 to 1. Non-zero diagonal elements indicate the degree of averaging or the lineardependence of the true model space as represented by the inverted model. Model parametersassociated with diagonal elements of the resolution matrix that are ?_ 0.5 are consideredmeaningful and well-resolved (see Zelt and Smith 1991).The criteria determining optimal parameterization lead directly to the stopping criteria:1. the final model chosen is that which gives the preferred trade-off between travel timeresidual and parameter resolution; and2. the final model must permit rays to be traced to all observation points.The crux of (1) is that adding more nodes to a model typically allows the travel time residualto be decreased, but at the same time decreases the overall parameter resolution. It is essentialto realize that the R matrix gives the resolution values for the parameters of the model updatevector Am; this does not provide knowledge about the resolution of the parameters of thefinal model.In terms of tracing rays to all observation points, it is often not possible using a modelobtained from an inversion. Zelt and Smith (1991) cited the cause of this as using moreparameters than the data can adequately constrain. This may lead to shadow zones because44of the potential for unrealistic oscillations in layer velocities or boundary shapes, which causerays to focus and scatter. A final model must therefore be rejected if it is not possible to tracerays to most observation points. Note that the ray take-off angles that connect sources andobserved receiver locations are determined so that ray tracing is performed in the "two-pointray tracing" sense, whereby rays are traced only between sources and observation points.3.2.3 The Modelling TechniquesThe approach used for modelling the data with the two algorithms was one of "layer strip-ping", by which all the observed travel times associated with a given layer are simultaneouslyinverted while the values of the parameters of the layers above are fixed. The model thusdetermined with the Zelt and Smith (1991) inversion algorithm was adjusted, if necessary,for the added constraints of amplitudes using the Zelt and Ellis (1988) forward modellingprocedure. This process of alternating travel time inversion and amplitude forward modellingwas continued until an acceptable fit of observed and calculated traveltimes and amplitudeswas achieved for successively deeper layers. In this way, the final model was constructed.This method clearly simplifies and expedites the inversion process. Testing has shownthat the method is valid, since, for example, PmP and Pn arrivals provide little constraint onupper and middle crustal structure (Zelt and Smith 1991). In practice, a trade-off betweenreducing the travel time residual and improving the computed amplitudes is inevitable. Moreemphasis was placed upon achieving a satisfactory travel time fit, because of the limitationsinvolved in calculating amplitudes using ray theory. Due to the large amplitude variationswithin each phase (which is a common trait of data from refraction surveys on land), thecomparison between computed and observed amplitudes was qualitative.In applying the modelling techniques in this manner, two assumptions were made: (1)the inversion algorithm is capable of resolving trade-offs between velocity and boundarydepth; and (2) the effects of using a 2—D (straight-line) representation of the actual shot andreceiver geometry are negligible. Zelt and Smith (1991) asserted that testing has revealed45assumption (1) to be valid. As well, assumption (2) seems valid in the case of Line 10, withthe possible exception of the vicinity to the southeast of SP-48 (Figures 2 and 4).3.2.4 Modelling the Refraction DataThis section describes the application of the modelling procedure to the Line 10 data setin order to construct the final model (Figure 11a). Figures 5c-10c show the match betweenobserved and calculated travel times, and the synthetic sections shown in Figures 5b-10ballow calculated and observed amplitudes to be compared. The location of the boundary andvelocity nodes used in the inversion are shown in Figure 12, which also shows the numberswhich were assigned to the boundaries and layers in the text. Figure 13 shows the raycoverage within (a) the upper crust and (b) the remainder of the crust and the upper mantle.All boundaries — including the surface, which represents the topography — in aninversion model are smoothed by uniform sampling at 250 points and the application ofa three-point averaging filter (Zelt and Ellis 1988). The sub-horizontal boundaries shown bydashed lines in Figure 1 la represent velocity discontinuities (those which delineate changes invertical velocity gradient have been omitted in the presentation). Since the velocity structureof the model suggested a natural division of the crust into thirds, the upper, middle, andlower crust were defined for the discussion. Table 6 provides a summary of the results ofthe inversion procedure.3.2.4a The Upper Crust The uppermost layer in Figure 1 la represents a thin (P-2.5 km)near-surface stratum. Velocities within the layer were determined by inverting the Ps arrivals.The starting model had 6 velocity nodes on the surface at the shot point locations, whichpreliminary modelling had shown to be the optimal node locations; indeed these were theonly zones of ray coverage. Velocity nodes representing values at the bottom of the surface46FIGURE 11: (a) Color contour plot of the velocities beneath Line 10 (vertical exaggerationapproximately 4:1). Numbers represent average P-wave velocities in km/s between thedashed lines (not nodal values). The distance axis has SP-51 as the origin ("model distance"),and the depth axis has sea level as the origin. The shotpoint locations are shown by solid tri-angles, and positions of refraction Line 2 of SCORE '89 and LITHOPROBE reflection profiles88-13, 88-14, and 88-17 are indicated. Thick solid lines represent areas which are constrain-ed by reflected energy. The transitional Moho layer is indicated by "M". The shaded regionaround the model is unconstrained by the travel time data (The upper mantle reflector hasbeen omitted here.)FIGURE 11: (b) A simplified schematic interpretation of (a) based partly on other studies(see text). Vertical exaggeration is approximately 2:1. Note that the interpreted extent of theCB batholiths does not mle out their presence below this depth.60-0—13---v-------a---556NW0SE350DISTANCE (km)50^1001^1501^200 250(3001^•0,,-----------:7 1---0-- 1 1 °.^- i•C\2o--00-78^el^n --- v10N 11FIGURE 12: Location of depth (squares) and velocity (circles) nodes used in the inversion process. The numbers in the left columnrepresent layer identification numbers used in the text, and those in the right column denote boundary numbers. There are no nodes shownon boundaries 3-5 because they were defined through forward modelling of amplitudes (see text). Boundary 1 represents the topography ofLine 10, layer 1 is the near-surface layer, layers 2-5 constitute the upper crust, layer 6 is the middle crust, layers 7 and 8 make up the lowercrust, layer 9 is the transitional Moho layer, and layers 10 and 11 are upper mantle layers separated by the upper mantle reflector shown asboundary 11. The shotpoint locations are shown by solid triangles. Vertical exaggeration is -3:1.0-DISTANCE (km)150^2000 25050 100 300 35010000200•TV00048(a)(b)DISTANCE (km)2150^ 00NFIGURE 13: Two-point ray tracing diagrams (all shots to all receivers) showing ray coverage in (a)the upper crust, and (b) the remainder of the crust and the upper mantle. The vertical exaggerationin (a) is approximately 9:1, and in (b) approximately 3:1. The shotpoint locations are shown bysolid triangles.49PhaseInvertedRMS TravelTime ResidualNormalizedChi-squaredMisfitNumber ofVelocityNodesNumber ofDepthNodesPs .01 7.7 11 0Pg , layer 2 .01 4.5 12 0Pg, layer 3 .12 1.9 0 0Pg, layer 4 .12 2.6 0 0R1 .06 0.7 0 3R2 .09 1.0 6 4R3 .11 1.3 6 3PmP .11 1.4 6 3Pn .09 1.1 1 0R .03 0.2 0 1Table 6: Inversion results for the final velocity model. The RMS travel time residuals, chi-squared misfits, and number of velocity and boundary nodes are displayed for each phase ofthe inversion procedure. Note that the RMS travel time residual is given in seconds. SeeFigure 12 for the placement of the boundary and velocity nodes.layer were initially placed 1 km beneath the shot point locations, a depth based upon theestimated thickness of the surface layer.The damping factor for the inversion was set at 1. A single value of a i (a,), 0.1 km/swas used for all velocities and a single value of ai (az), 0.1 km was used for all boundarynodes. These parameters were unchanged for each step of the layer-stripping procedure,except for the cases which are described below. The relative values of ay and az determinethe trade-off between the size of velocity and boundary adjustments in the inversion.Since the Ps arrivals were limited in both number and spatial extent, extracting velocitygradient information was not feasible. The velocity gradient in this surface layer was thereforefixed at approximately 0.5 s -1 for the inversion procedure, based on a priori estimates ofthe velocities at the top of and at the base of the layer.50The resulting model had an acceptable travel time fit (i.e. the RMS travel time residualof 0.01 s was of the same order as the uncertainty of the travel time picks), and rays weretraced to all observations. However, the resolution of the velocity nodes — 0.57, 0.84, 0.47,0.42, 0.43, and 0.12 — was relatively low, a reflection of the fact that the shot-receivergeometry was inappropriate for obtaining extensive near-surface information.To determine the base of the Ps layer, and to constrain velocities in the upper crust, the702 Pg arrivals were inverted. The modelling of the Pg arrivals was initiated by placing 6boundary nodes beneath each shot point at the Ps-Pg interface (boundary 2 on Figure 12).Six velocity nodes were also placed at these locations, and six additional velocity nodes wereadded at the necessary intervals to invert for the velocity at the top of the Pg layer. Thisapproach was adopted after initial modelling revealed that an acceptable travel time fit couldnot be obtained with fewer velocity nodes.Through forward modelling it was determined that the deepest turning point of the Pgfirst arrivals was 8 km, so lower velocity nodes were placed at this depth beneath the 12 uppernodes; this estimated depth of penetration was based solely upon the extent to which the Pgphase was visible as a first arrival on the field data. The boundary which appears at 8 kmdepth in Figure 12 (Boundary 5) therefore does not represent a structural feature; it indicatesthe depth to which velocity calculations in the crust were controlled by refracting rays.Preliminary inversion results indicated that inverting for the base of the Ps layer was notgoing to be successful — the boundary node positions oscillated with successive iterationsinstead of converging to a solution, and there were wild undulations in the boundary whichled to shadow zones in ray coverage. As a result, the Ps-Pg boundary was determined throughforward modelling of the Pg arrivals. This procedure required that 6 extra boundary nodesbe placed along the Ps-Pg interface. Note that these boundary nodes (12 in total) do notappear in Figure 12, as they were determined through forward modelling alone.Velocities to 8 km depth were then determined by carrying out the inversion procedure.51During the course of the inversion procedure, some of the idiosyncrasies of the inversionprogram mentioned by Zelt and Smith (1991) became evident. The velocity gradient beneathSP-49 became sufficiently high that the Pg arrivals observed at large offsets for shot point46 were not getting rays traced to them. This problem was solved when the velocity gradientbelow SP-49 at the bottom of the layer was held fixed at a more reasonable value for aportion of the inversion procedure. Also, the velocity below SP 51 was unreasonably large,so it was necessary to hold this value fixed during some iterations. For the final iterations,the damping factor was increased to as high as 25 in order to come up with the model thathad the best travel time fit and had rays traced to the maximum number of observation points.While the inversion for the velocity in the Pg layer was taking place, it was necessaryto continue the modelling of the Ps layer in order to account for "kinks" in the essentiallylinear Pg branch. It was apparent that these small undulations were caused by geologicalfeatures, as their positions did not change on the six different record sections; as such, theywere modelled as zones of higher or lower velocity material in the near surface layer. Theresulting Ps layer was defined on its top and bottom by 11 velocity nodes; the bottom nodesdo not appear in Figure 12 because the gradient in the overlying material was fixed, so thevalues attached to these nodes bottom depend upon values determined for the upper nodesthrough inversion.One of the "kinks" was relatively large — it appears at 50-60 km model distance on thePg branches of SP's 49-51 (Figures 8a-10a). This feature corresponds with eight receiverstations in the northwest which traversed the Cadwallader terrane (Figure 4), a Coast Beltterrane of island-arc clastics and volcanics which may be linked with Stikinia (Wheeler et al.1991); it was successfully modelled as a localized high velocity (5.0-5.6 km/s) zone. Thislocation represents one of two times that the refraction line crossed a sliver of Cadwalladerterrane; however, a similar kink in the data does not occur at the second crossing. Althoughthere is no obvious explanation for this feature of the data, the fact that the second crossing52occurs where Line 10 deviated from a straight-line configuration may be significant.The velocity model determined down to 8 km in this manner was then used as thestarting model for the amplitude forward modelling procedure. The initial model wasquite problematic in terms of its amplitudes — there was a general trend on the syntheticseismograms for the amplitudes to increase with offset which was not seen on the seismicsections plotted from the field data. This indicated that the rays which penetrated more deeplyshould be passing through material with a lower gradient than the more shallow-diving rays.In order to accommodate this reduction in gradient with depth, two new boundaries wereadded to the Pg layer, at approximately 2 km and 4 km depth (see boundaries 3 and 4, Figure12). These boundaries, across which there is no velocity discontinuity, delineate a changein velocity gradient. The gradients were then adjusted until the values gave a satisfactoryfit to the amplitude data. Each time a velocity node value was changed in order to changethe velocity gradient, forward travel time modelling was done to ascertain that there was nosignificant deterioration of the travel time fit.The amplitude fit was acceptable when the second, third, and fourth layers had verticalvelocity gradients of about 0.027 s -1 , 0.009 s -1 , and 0.003 s -1 , respectively. The traveltime data were then inverted again to solve for the velocities at the top of layer 2 with thevelocity gradients fixed at these values. The resolution of the 12 velocity nodes in layer 2is good (> 0.71) in all cases.Keeping the velocities in layers 2 — 4 fixed at their final values, the Ps and Pg arrivalswere inverted again for all 11 of the upper velocity nodes in layer 1. The resulting RMStravel time residual is 0.01 s (Table 6). The resolution values are generally acceptable atthe shot point locations (values were between 0.48 and 0.84), while resolution values atthe intervening nodes range from 0.18 to 0.52. The low resolution values of certain nodesindicates that the model depends strongly on the parameterization.The modelling of the upper crust was completed by inverting for the R1 reflection. Tests53showed that the velocity in layer 5 (Figure 12) is not constrained by the data, so the velocityin this layer was modelled by fixing the vertical velocity gradient at 0.003 s -1 as in layer4, and having no velocity discontinuity between layers 4 and 5. Boundary 5 in Figure 12 istherefore a pseudo-boundary, as the velocities in the layer below are a continuation of thosecalculated for the layer above. Three boundary nodes were placed at a depth of 10 km, andwere inverted to determine the depth to the reflecting horizon; the number and position ofthe boundary nodes were determined to be optimal through preliminary modelling.The inversion algorithm gave a model with the R1 reflector having an average depthof 10.3 km; the reflector is flat throughout most of the model, with an upward dip at thenorthwestern end. The three boundary nodes used in inversion have excellent resolution(0.94, 0.96, and 0.94), and the RMS travel time residual value is 0.06 s for the R1 arrivals(Table 6).In the case of R1, the boundary was extrapolated between the portions which areconstrained by the reflections (see Figure 11a). This procedure was adopted for all boundariesdetermined through the inversion of reflected phases.The final upper crustal model has 23 velocity nodes and 3 depth nodes. The travel timefit was 0.096 when rays were traced for all the phases. Rays could not be traced to 10 ofthe 702 Pg data points for the reasons described earlier, and the corresponding travel timeswere therefore not used in the inversion.3.2.4b The Middle Crust The mid-crustal structure of the model was determined byinverting the R2 reflection. Preliminary modelling showed that the velocity gradient withinlayer 6 is not constrained by the data, so the velocity gradient was again fixed at 0.003 s -1as in layers 4 and 5. This procedure was deemed reasonable because there was no evidencein the data for strong gradients in this layer, and it was consistent with similar studies in theSierra Nevada (eg. Bolt and Gutdeutsch 1982, Pakiser and Brune 1980).54To invert this mid-crustal reflection, 4 boundary nodes were placed at a depth of 20 km.Velocity nodes were placed at the top of layer 6 (boundary 6 in Figure 12) beneath four ofthe shot points, to the exclusion of the two on the ends of the profile. Preliminary modellinghad indicated that this was the optimal number of nodes, as inverting with less nodes led to amodel which did not accurately fit the data, and using more nodes significantly decreased theresolution. The location of the nodes was more arbitrary, since the positioning criteria couldbe met equally by a number of placements, as long as the nodes were roughly evenly-spacedacross the region constrained by the data.The inversion produced an initial model which had the R2 reflector at an average depth of22.2 km, except an upwelling which reached , 19 km depth beneath SP 48, where the resultis constrained by data from both SP-46 and SP-47. For this initial model, the velocitieswere 6.30 km/s beneath SP's 50, 49 and 47, and 6.55 km/s beneath SP 48; the portions oflayer 6 below SP's 51 and 46 are not constrained by the data, as indicated above.With the modelling of layer 6 (Figure 12) completed through the inversion of R2, thevelocity contrast across boundary 6 was specified, so that the forward modelling of the R1reflections was possible. This starting model was problematic — the R1 reflections from thenorthwestern part of the line were too low in amplitude, and those from the southeastern parthad traces in which the amplitudes were too high. The velocities in layer 6 were adjusted untilthe R1 amplitudes on the synthetic seismograms matched those on the original data; holdingthese velocities fixed, the R2 arrivals were inverted to correctly position the boundary nodesat the base of layer 6. During the course of the R1 amplitude modelling, it was necessaryto add a velocity node to the northwest of the one beneath SP-50.The resulting R2 reflector has an average depth of 22.8 km, except for the upwellingbeneath SP-48, which reaches a depth of 19.7 km (Figure 11a). The velocities determinedfor the layer are ,--6.35 km/s beneath SP-50 and SP-49, 6.45 km/s beneath SP-48, and 6.40km/s beneath SP-47. The velocities beneath SP-51 and SP-46 are not constrained due to55inadequate data. The RMS travel time residual is 0.09 s for the R2 arrivals (Table 6), andthe resolution of the boundary nodes is very good: 0.97,0.97,0.94, and 0.70. The resolutionof the 4 velocity nodes used in the inversion is also good (0.71 or greater in all cases).3.2.4c The Lower Crust and Upper Mantle The modelling of the lower crustal velocitystructure was initiated by inverting the R3 arrivals. The vertical velocity gradient was fixedas for the R2 phase. Three boundary nodes were placed at a depth of 30 km, and velocitynodes were placed beneath four of the shot point locations at the top of layer 7 (boundary 7in Figure 12). Preliminary modelling was used to ascertain that this was the optimal numberof nodes in terms of the criteria described above; again, the placement of the nodes wassomewhat arbitrary, as described for the boundary 6 case.The initial model which resulted from the inversion procedure had velocities in layer 7of 6.35, 6.40, 6.50, and 6.45 km/s beneath SP's 50, 49, 48 and 47, respectively, and had aflat R3 reflecting boundary with an average depth of 29.9 km. The portions of the modelbeneath SP's 51 and 46 are not constrained by the data.When the model based upon travel time inversion was used in the forward amplitudemodelling procedure, the amplitudes on the synthetic sections were seen to be generally toosmall for the R2 reflections in comparison with those on the true data. The exceptions werethe R2 reflections from SP's 46 and 47, which compared well in terms of amplitudes. Thevelocities in layer 7 were adjusted until the amplitudes of the R2 reflections on the syntheticseismic sections matched those on the real data sections; with these velocities held fixed, theR3 arrivals were inverted to model the base of the layer.The inversion model which produced an acceptable travel time fit (see Table 6) hasvelocities of ,6.45 km/s beneath SP-50, 6.50 km/s beneath SP-49, and SP-48, and 6.45km/s beneath SP-47 (see Figure 11a). The velocities at the nodes beneath SP's 51 and 46were arbitrarily fixed at the values of the adjacent velocity nodes. The R3 reflecting boundaryis essentially flat at a depth of 30.0 km. The resolution values of the boundary nodes are56excellent (0.98,0.98,0.97), as are those of 3 of the 4 velocity nodes used in inversion (> 0.83).For the velocity node beneath SP-47, the resolution is poor (0.32) because the ray coveragewas not high in the vicinity. However, added confidence may be ascribed to this valuebecause the R2 reflection on the southwestern portion of SP-48 boundary 7 is constrainedat the location of its node, and the amplitude of the reflection therefore provides additionalconstraint.Since both the PmP and Pn arrivals constrain the structure of the Moho, these arrivalswere inverted simultaneously to model boundary 9 on Figure 12. According to the criteriadescribed above for the R2 and R3 modelling procedures, the velocity gradient was fixed, 3boundary nodes were located at the base of layer 8, and velocity nodes were placed beneathfour of the shot point locations at the top of layer 8 (Figure 12).Since the Pn data were so limited (the upper mantle was sampled only in the upperkilometer, and the lateral extent of this was confined to within 150 and 300 km model distance,as shown in Figure 13b), it was necessary to fix the vertical velocity gradient in the uppermantle at a reasonable estimated value before performing the inversion. Amplitude modellingtests led to the selection of 0.003 s -1 — a gradient of 0.004 s-1 produced amplitudes of thePn arrivals on synthetic seismic sections which were much larger than those on the observedsections, and the choice of 0.002 s-1 caused shadow zones in ray coverage.The velocity at the top of the upper mantle is represented by a single velocity node (i.e.this was modelled as a laterally homogeneous layer, since the information available was solimited). Preliminary modelling showed that inversion for the crust-mantle boundary resultedin so much relief that rays could not be traced to all sites for the Pn arrivals; the dampingfactor was increased to 10 to alleviate this problem.Inverting the PmP and Pn arrivals in this manner led to a preliminary model in whichthe velocity values above the Moho (layer 8 in Figure 12) ranged from 6.90 km/s beneathSP-50 to 6.45 km/s beneath SP-47. The crust-mantle boundary had an average depth of5734.5 km, except at the northwestern end of the line, where it sloped upward.Forward modelling of R3 amplitudes led to the following changes for layer 8: the velocitynode below SP-50 was moved roughly 25 km closer to SP-49, and the velocity at this nodewas decreased to 6.70 km/s, and the velocity at the nodes beneath SP's 49 and 47 wereincreased to 6.65 km/s and 6.60 km/s, respectively. The velocity of 6.70 km/s at the nodebeneath SP-48 was not changed. The velocities of the nodes centered under SP's 51 and 46were arbitrarily fixed at the values of the adjacent velocity nodes. This overall decrease invelocity contrast from NW to SE means that the observed trend of decreasing R3 amplitudefrom NW to SE is replicated in the calculated data (compare Figures 5b-10b with Figures5a-10a). With these velocities held fixed, the PmP and Pn arrivals were again simultaneouslyinverted for the boundary nodes at the base of layer 8 and the upper mantle velocity node.The initial model of the crust-mantle boundary which resulted from the inversion processwas essentially flat with an average depth of 34.6, except for an upward-sloping section atthe northwestern portion of the line, which reached a depth of 32.4 km at the northwesternextent of the portion of the boundary constrained by the data. The velocity at the top of theupper mantle was 8.05 km/s. The RMS travel time residual for the PmP and Pn phases was0.110 with this model, the resolution values of the velocity nodes in layer 8 were acceptable(0.49-0.87), and those of the crust-mantle boundary nodes were excellent (> 0.97).When this model was used in the amplitude forward modelling program, the amplitudesof the PmP arrivals on the synthetic seismic sections seemed quite large in comparisonwith those on the true data sections. This situation is most readily explained by one of thefollowing:1. the velocity of the upper mantle is too high,2. the Moho cannot be modelled simply as a single-order velocity discontinuity.In order to see if the first explanation was feasible, tests were performed in which thevelocity of the upper mantle was fixed at every reasonable value lower than 8.05 km/s, i.e.587.6 km/s, 7.65 km/s, etc... These tests showed that this lowering of upper mantle velocitydid not sufficiently decrease the amplitude of PmP.The second explanation is consistent with amplitude modelling in previous studies (e.g.Benz et al. 1990, Zelt and Smith 1991, and Zelt et al. 1992a), which showed that thesimilarity between observed and calculated data increased when the Moho was modelled asa transition zone rather than a velocity discontinuity. The method used in these studies wasto incorporate a transitional Moho layer with velocities ranging from 7.4-7.8 km/s; PmP andPn were then modelled respectively as reflections from the top of, and refractions below,this transitional Moho, which ranged in thickness from 0.5-3.5 km (Zelt et al. 1992a) to2.0-5.0 km (Benz et al. 1990).There were two other aspects of the data which suggested that the modelling shouldbe done in this manner. First, the calculated critical points were ,10-20 km closer to thesource than the observed critical points. As was the case in amplitude modelling, decreasingthe upper mantle velocities did not alleviate the problem. Also, the PmP phase of SP-48NW had a significantly lower amplitude than the others. This amplitude anomaly could bemodelled with lateral velocity variation in either upper mantle if the Moho was modelledas a single-order velocity discontinuity, or within the transitional Moho itself; it was moredifficult to justify the former using this one feature of the data as a basis.The velocities in the transitional Moho (layer 9 in Figure 12) were determined throughthe forward modelling of PmP amplitudes. Velocity nodes where placed at the top of theMoho layer beneath each of the shot point locations for the forward modelling procedure,and the velocity node below SP-50 was assigned a value of 7.8 km/s, below SP's 49 and48 a value of 7.4 km/s, and below SP 47 a value of 7.5 km/s. The nodes beneath SP's51 and 46 where the data do not constrain the amplitudes were given the velocity valuesof the adjacent nodes.Since the Pn data were so limited (33 observed arrivals were picked), it was not possible59to model the base of the transition zone. It was therefore arbitrarily decided that the base ofthe layer (boundary 10 in Figure 12) would parallel the top boundary of the layer in shape,for a transitional Moho of uniform thickness (chosen to be 2 km in Figures 11 and 12).With layer 9 defined in this manner, the Pn arrivals were inverted to determine thevelocity in the upper mantle; the resulting velocity value was 8.15 km/s. Obviously, thisvalue depends upon the arbitrarily chosen thickness of the transition zone. The extent of thisdependence was shown by a simple test in which the thickness of layer 9 was decreased to1 km: inverting the Pn arrivals led to an upper mantle velocity of 8.10 km/s.This modelling of the crust—mantle boundary as a transitional layer improves thecalculated amplitudes of all PmP arrivals, puts the critical point of all the arrivals into theappropriate position, and provides a reasonable means of incorporating the PmP amplitudeanomaly of SP-48 NW into the model. The transitional Moho is therefore included in thefinal model, although its thickness is not constrained. This may provide a more realisticpicture of the true nature of the crust-mantle boundary, which has been interpreted by some(see Mooney and Brocher 1987) as a zone of inter-layered materials of higher (ultramafic)and lower velocities (mafic).The upper mantle velocity structure was completed by inverting the arrivals from thesub-Moho reflection R visible on the seismic section from SP-46. One boundary node wasused in the inversion procedure, which resulted in a sub-Moho reflector at 70.3 km depth inthe final model. The resolution of the one boundary node used is very good (0.88). The depthof this reflection depends upon the velocity chosen for the upper mantle. If, for example, a1 km thick transition zone is chosen, the upper mantle velocity is 8.10 km/s, and the depthto the R reflector is 68.3 km.When amplitude modelling was performed for the R reflection, it was shown that a smallvelocity increase (-Ai km/s) resulted in amplitudes on the synthetic seismograms which arecomparable to those on the true data sections. Similarly, a negative velocity discontinuity60(,0.2 km/s) provides a good fit to the observed data (see Figure 5e cf. Figure 5f). Sincethe relative phase (polarity) of the arrival could not be distinguished, the data do not enablea preferred velocity change to be inferred.61IV. INTERPRETATION4.1 The Final Model4.1.1 Primary Features of the Velocity ModelThe principal features of the velocity model are shown in Figure 1 la. The uppermoststratum is a thin (< 2.5 km) near-surface layer with velocities averaging 3.9 km/s at the topand 4.9 km/s at the bottom. Within this average velocity structure, there are zones of higherand lower velocity responsible for fluctuations in the essentially linear Pg branch. The mostsignificant of these is a high velocity (5.0-5.6 km/s) zone which is located below SP 50, andhas been linked with a feature of the surface geology (the Cadwallader Terrane).The velocities at the top of the upper crust average 6.15 km/s, and those at the bottom6.25 km/s, with higher velocities occurring in the southeast. The depth to the lower boundaryof the upper crust (^40 km) is constrained for the northwestern third and for an ,30 kmsegment at the southeastern part of the model by the R1 phase. In the northwest, the R1reflection is from a relatively large velocity contrast of —6.15-6.20 to 6.35 km/s. The baseof the upper crust in the center of the model where the velocity contrast is intermediate(-0.10 km/s) is not constrained by observations of reflected energy, while the R1 reflectionin the southeast occurs where the velocity contrast between the upper and the middle crusthas decreased significantly to 40.05 km/s (Figure 11a). This suggests that another source ofreflectivity (e.g. fine-scaled layering which cannot be resolved in this experiment) other thansimply a relatively small single-order velocity discontinuity may have a role in generating R1in the southeast. Alternatively, the situation may be an artifact of the modelling procedure(i.e. the blocky model parameterization).The velocities in the middle crust average 6.35-6.45 km/s, with slightly higher valuesbelow SP 48, for an overall velocity structure of little lateral variation. The depth to the R262reflecting horizon defines the base of the middle crust at ,22 km; this boundary is well-constrained for the center portion of the model (i.e. 100-250 km model distance). As in theupper crust, the two northwestern reflected phases are generated by a velocity discontinuityof at least 0.15 km/s, and the southeastern reflection by one of ,0.05 km/s. Since the baseof the middle crust is not similarly constrained by observed R2 reflections in areas with anintermediate velocity contrast (,0.10 km/s), it is again suggested that perhaps other sourcesof reflectivity may be in effect.The lower crust, which has the Moho as its base, is divided into an upper and a lowerlayer by the R3 reflector. The lower crust is also homogeneous in its velocity structure withinthe resolution of the survey, with values averaging 6.45-6.50 km/s in the upper portion and6.60-6.70 km/s in the lower layer.The depth to the top of the transitional Moho averages 34.5 kin, dips to 35.5 km at thesoutheastern end of the model, and shows minimal topography. The velocities within thetransitional layer vary from ,7.80 km/s in the northwest to 7.40 and 7.50 km/s in the centerand southeastern portions of the line, respectively. The top of the crust-mantle boundary iswell-constrained by PmP reflections from 50-300 km (model distance), but its thickness is notconstrained by the data. Modelling the Moho as a transitional layer increases the similaritybetween the observed and the calculated data; as well, if the Moho transition is explained as acomplex series of interlayered high and low velocity materials, this presentation is consistentwith recent studies of the character of the crust-mantle boundary (e.g. Hale and Thompson1982, Braile and Chiang 1986).While the sparsity of the Pn arrivals and the unknown thickness of the Moho rule outconstraining the upper mantle velocity, the available data suggest that velocities are relativelyhigh (e.g. ,8.15 km/s for a 2 km thickness, r-8.10 km/s for a 1 km thickness of the Mohotransitional layer) in this region. An upper mantle reflector which is not shown in Figure1 la presentation was imaged at a depth of ,70 km, a value which is also dependent upon63the Moho thickness.4.1.2 Non-uniqueness of the ModelSince two different models may both produce an acceptable fit to the data, be well-resolved, and appear physically acceptable, non-uniqueness is inherent to the final model.The scope of the non-uniqueness issue in this problem is outlined below.The non-unique nature of the solution is evident from the subjective choice of nodeplacement. Furthermore, for a certain choice of node placement, different a priori estimatesof parameter uncertainties, different damping parameters, and different starting models areassociated with a different final model. Even within a certain parameterization, there is asense of non-uniqueness associated with any nodal value output from the RAYINVR inversion— each calculation has an associated standard error, indicating a range of feasible values.Thus, it must be stressed that the final model presented above is a model which fits the data,not a unique solution of the inverse problem.4.2 Spatial Resolution and Absolute Parameter UncertaintyZelt and Smith (1991) described tests which provide estimates of (i) the spatial resolutionof a model about a specific velocity or boundary node, and (ii) the absolute parameteruncertainty of a particular node (note that the term "spatial resolution" is used here todistinguish this quantity from the numerical resolution discussed in a previous chapter).The first step in estimating spatial resolution is to choose a model parameter and perturbit by a measure on the order of its estimated uncertainty o — the perturbation should belarge enough to change travel time, but small enough that ray path pattern is not disturbedsignificantly. The next step is to trace rays through the perturbed model to calculate thetravel times corresponding to the receiver locations and assign the observed pick uncertaintyat each location. Finally, the parameter is reset to its non-perturbed value, and the calculated64data are inverted to determine all the model parameters that were computed at the same timeas the selected parameter during the inversion for the final model.The spatial resolution about the selected parameter is indicated by the amount that thevalues of adjacent parameters differ from the corresponding values in the non-perturbedmodel. If the model is well-resolved about the selected parameter, then all other parametervalues will resemble the non-perturbed values; poor resolution about selected parameter willresult in the perturbation being "smeared" into the adjacent parameters.To estimate the absolute uncertainty of a model parameter, Zelt and Smith (1991)suggested another test, the first step of which again involves choosing a model parameterto perturb. The second step is to invert all other parameters that were determined at thesame time as the selected parameter during the inversion for the final model, while holdingthe selected parameter fixed at its perturbed position. Finally, the size of the perturbation iscontinually increased until the data cannot be fit as well as they were in the non-perturbedmodel according to 1) the ability to trace rays to all observations, and 2) the result of anF—test comparing the chi-squared values of the two final models. The statistical F—testinvolves a ratio calculated such that values either >>1 or <<1 indicate very significantdifferences between two samples; in this application, it reveals in a statistical sense whetheror not the fit of a perturbed model has deteriorated significantly from that of the final model.The size of the largest perturbation which leads to an equivalent fit to the data givesan indication of the absolute uncertainty of the selected parameter. This test for absoluteparameter uncertainty, however, does not account for the added constraints of forwardamplitude modelling. As acknowledged by Zelt and Smith (1991), this type of analysison every node would take longer than the inversion procedure; they suggest that testing anumber of representative velocity and boundary nodes may be sufficient.The estimated spatial resolution and absolute uncertainty values thus obtained for thefinal model are shown in Table 7. Generally, the resolutions decrease and the uncertaintiesVelocities in near-surface layer (layer 1)Velocities at top of upper crust (layer 2)Depth nodes at base of upper crust (boundary 6)Velocities at top of middle crust (layer 6)Depth nodes at base of middle crust (boundary 7)Velocities at top of lower crust (layer 7)Depth nodes of R3 in lower crust (boundary 8)Velocities within the lower crust (top of layer 8)Depth nodes at crust-mantle boundary (boundary 9)Velocities in the upper mantle0.5 km/s0.1 km/s1.5 km0.1-0.2 km/s1.5 km0.2-0.3 km/s1.5 km0.2-0.3 km/s1.5 km0.1 km/s**20-7010-1535403550503050*65increase with depth in the model. As well, this testing procedure revealed a trade-offbetween resolution and uncertainty based upon node spacing: in comparison with widely-spaced nodes, those which were more closely spaced were shown to have both a decreasein resolution and an improvement in absolute uncertainty.Estimated Estimatedlateral^absoluteresolution^uncertaintyNodes^ (km)* There were not sufficient Pn data to determine this quantity.**based upon an assumed value of Moho thickness.Table 7. Estimated lateral resolution of the final velocity model about the given nodes, andabsolute uncertainty of velocity and depth nodes.66V. DISCUSSION AND CONCLUSIONS5.1 Comparison with LITHOPROBE Refraction InterpretationsThe refraction lines of SCoRE '89 — Lines 1, 2 and 3 in Figure 2 — were all withinthe proximity of Line 10. As such, the interpretations of these profiles which exist at presentmay be used for comparison with the results obtained from the Line 10 study. The data fromLine 1 and Line 3 have been interpreted by Zelt et al. 1992a and 1992b, respectively. Apreliminary interpretation of Line 2 has been made by McLean (in preparation).The Line 2 data modelled by McLean (in preparation) were provided by a profile whichran from the eastern Insular Belt across the Coast Belt and into the western IntermontaneBelt (see Figure 2). The proposed model features upper crustal velocities of - ,6.25 km/s(to -10 km depth), mid-crustal velocities of 6.35 —6.50 km/s to 17 km depth overlyingmaterials of 6.55 km/s average velocity which extend to ,24 kin depth, and lower crustalvelocities of 6.8 km/s. The Moho has an average depth of approximately 34 km, and theunderlying upper mantle has velocity of 7.9 km/s.These preliminary results of the line 2 interpretation are important to this study sinceLine 2 intersected Line 10 (Figure 2; reflection profile 88-13 also intersected Line 10 at thislocation). It is important to ascertain that the interpretations of the two refraction data setsagree at their intersection point, or to provide a plausible explanation for any discrepancy.The Line 2 data interpreted by McLean (in preparation) exhibit the same distinct features—i.e. Pg phase confined to upper crustal region with underlying crust giving rise to threeprominent reflection phases — as the Line 10 data, and the parameterization chosen in themodelling procedure was very similar. This suggests a comparison as shown in Table 8,where the approximate values of the main crustal features of the models at their intersectionpoint (shown in Figure 11a) are displayed. Clearly, the two models agree closely, with the67Model Feature^ Line 2^Line 10average velocity above R1depth to R1 reflectoraverage velocity between R1 & R2depth to R2 reflectoraverage velocity between R2 & R3depth to R3 reflectoraverage velocity between R3 & Mohodepth to Moho6.25 km/s10.0 km6.40 km/s20.0 Ian6.55 km/s24.5 km6.80 km/s35.0 km6.30 km/s10.0 km6.40 km/s22.0 km6.50 km/s29.0 km6.70 km/s34.0 kmTable 8. Velocity values in the crust and depth to reflectors for Line 2 of SCoRE '89 andLine 10 of SCoRE '90 at the intersection point of Line 10 and Line 2 (see Figure 2).exception of the depth to the R3 reflecting horizon. Explanations for this discrepancy includethe 3—D nature of Line 10 at this location being ignored in the 2—D approach, the combinationof the absolute errors resulting from the modelling procedures, the possibility that a differentcycle of the phase was picked (testing has shown that a difference of one cycle in time wouldbe equivalent to an "s1 km shift in depth), and the possibility of a structural variation inthe reflector (although the near-coincidence of the region of R3 constraint on each of thetwo models would necessitate a rather abrupt inflection). Alternatively, the most reasonableexplanation for this 4.5 km difference may be that the source of the McLean (in preparation)R3 reflection is an entirely different horizon than the one defined by the Line 10 data.Line 3 of SCoRE '89 extended from the Insular Belt through the southernmost CoastBelt to the western border of the Intermontane Belt (Figure 2), passing about 20 km southof the southeastern end of Line 10. An interpretation of the in-line data has been presentedby Zelt et al. (1992b). The proposed model shows average velocities in the upper crust (toa depth of , 11 km) of 6.3 km/s in the west and 6.2 km/s in the east; in the middle crust (toa depth of ,22 km) of 6.7 km/s in the west and 6.35 km/s in the east; in the lower crustof 7.0 km/s in the west and 6.8km/s in the east; and in the upper mantle of 7.65 km/s. The68Moho has an average depth of 37 km throughout most of the model, except in the west,where a depth of ,32 km is reached at the Insular Belt-CB border. An important featureof the model is an abrupt east-west variation interpreted in the Harrison Fault region; Zeltet al. (1992b) have identified this as the location of the collision zone between the Insularsuperterrane (which exists as Wrangellia in this region) and the Intermontane superterrane.The southward projection of Line 10 intersects Line 3 in the region of Harrison Fault, the Insular superterrane-Intermontane superterrane boundary. Table 9 shows the velocityvalues in the upper crust (to 10 km depth), middle crust (10-22 km depth), and lower crust(22 km — Moho depth) for Line 10 and Line 3 at the intersection location. The Moho depthwas interpreted for Line 3 to be ,-34.5, as was the case for Line 10. "Line 3W" and "Line3E" refer to the portions of the Line 3 model directly to the west and east of the HarrisonFault, respectively. Table 9 shows that the velocities determined for the southeastern portionof Line 10 are much closer to Zelt et al.'s (1992b) values for the eastern side of the collisionzone, i.e. those of the Intermontane superterrane.Model Featurevelocity of upper crust (km/s)velocity of middle crust (km/s)velocity of lower crust (km/s)Line 3W Line 10 Line 3E6.30 6.20 6.206.70 6.40 6.357.00 6.55 6.8Table 9. Velocity values in the crust for Line 10 of SCoRE '90 and Line 3 of SCoRE '89at the intersection point of Line 10's projection onto Line 3 (see Figure 11a). The uppercrust 0-10 km depth, the middle crust 10-22 km depth, the lower crust ti 22-34.5 kmdepth. The "W" and "E" refer to the portions of Line 3 to the west and east, respectively,of the intersection point. See Figure 14a for a pictorial representation.The Intermontane Belt to the east of the CB was the site of Line 1 of the 1989 SouthernCordillera Refraction Experiment (SCoRE '89, Figure 2). Zelt et al. (1992a) have presenteda 2-D model for the along-strike line which shows (1) upper crustal and mid-crustal velocitiesvarying laterally from 6.0 to 6.4 km/s; (2) lower crustal velocities averaging 6.5km/s and 6.769FIGURE 14: (a) and (b) show comparisons of 1—D average velocity (in km/s) versus depth profilesfor the refraction interpretations indicated (see Tables 9 and 10). The shaded region in the Line I(VISP '80) presentation represents an area that is poorly constrained. The "W" and "E" representthe portions of Line 3 to the west and east, respectively, of the intersection point. The Moho depthsare indicated with a thick horizontal line.70km/s at the top and the base; (3) a transitional Moho of 1-3 km thickness at an average depthof 33 km overlying an upper mantle with velocities of 7.9 km/s in the south and 7.7 km/s inthe north; and (4) an upper mantle reflector at ,--16 km below the Moho. thought to representthe base of the lithosphere, which is considered to be thin throughout the Intermontane Belt.For the calculation of upper mantle velocities, Zelt et al. (1992a) had many observations ofPn to use in the inversion algorithm, while the upper mantle velocity determined using theLine 10 data is not well-constrained, as will be discussed later.The velocities beneath Line 1 of the Intermontane Belt as interpreted by Zelt et al.(1992a) bear a strong resemblance to the Line 10 velocities, as shown in Table 10. Further-more, the interpretation of Spence et al. (1985) of VISP 80 Line I has shown Wrangellianvelocities in the upper 18 km which are much higher (at 6.3-6.7 km/s, see Table 10) thanthe velocities determined through the Line 10 analysis.Crustal Velocity (km/s) Line 10upper crust^6.20 (to 10 km)middle crust 6.40 (to 22 km)lower crust^6.55 (to 34.5 km)Line 1 (SCoRE '89)6.20 (to 12.5 km)6.20-6.40 (to 23 km)6.60 (to 33 km)Line I (VISP '80)6.30 (to 9 km)6.70 (to 18 km)* The Line I model of Spence et al. (1985) was not well-constrained below 18 km depth.Table 10. Velocity values in the crust (to depths approximated in parentheses) for Line 10 ofSCoRE '90, Line 1 of SCoRE '89, and the easternmost part of VISP '80 Line I. See Figure14b for a pictorial representation.Clearly, this overview of both concurrent and previous regional refraction studies showsthat the velocities beneath Line 10 are similar to those derived in studies of the Intermontanesuperterrane, and significantly slower than those interpreted for elements of the Insular supert-errane (i.e. Wrangellia); this information proved to be critical to the tectonic interpretationof results from this study.71The Zelt et al. (1992a) interpretation of the Line 1 upper mantle reflector being the baseof the lithosphere is also significant. The review of geophysical, geological, and petrologicalstudies from the Cordillera by Gough (1986) indicated that the lithosphere is thin throughoutthe Intermontane Belt, and tapers to deeper levels beneath the Coast Belt. This suggests thatthe upper mantle reflector below Line 10 could, at a depth of ",35 km beneath the Moho,represent a western counterpart to Zelt et al.'s reflection from the top of the asthenosphere,an idea supported by the fact that a negative velocity contrast across the R reflecting horizonwas shown to provide as good a fit to the observed data as a positive velocity contrast.Another horizon that was considered as a possible source of the upper mantle reflection Rwas the top of the subducting oceanic crust. However, the depth of the subducting platehas been inferred to be , 100 km at a point -,50 km west of Line 10 based on hypocenters(Rogers et al. 1990) and the presence of the Garibaldi volcanic chain (Dickinson 1970).This precluded the top of the subducted plate as a contingent source of the upper mantlereflection R, which was imaged at a depth of ,70 km. Further, the multidisciplinary reviewof Gough (1986) has suggested that the subducted plate loses its identity before it reachesthe region beneath Line 10.5.2 Comparison with Other Geophysical Studies —Reflection, Heat Flow, Gravity, and Geomagnetic5.2.1 Comparison with LITHOPROBE Reflection Data and InterpretationsIn 1988, nearly 350 km of deep seismic reflection data were acquired in the southern CBthrough a LITHOPROBE seismic reflection transect (profiles 11 to 18; Figure 1 of Varseket al. 1992). Lines 88-13, 88-14, 88-17, and, to a lesser degree, 88-18 were in the vicinityof the present study area (Figure 2). The seismic sections generally show east-dipping uppercrustal (i.e. above P-3.5 s) reflectors, which are truncated by the subhorizontal to west-dippingstructures that characterize the middle and lower crust (Varsek et al. 1992). Monger and72Journeay (1992) have provided a preliminary interpretation of the data from these profiles,and Varsek et al. (1992) have made a more detailed interpretation.Figure 15 shows the reflection data for (a) profile 88-13, and (b) profiles 88-14 and88-17; the location of Line 10 is indicated. As shown in Figure 2, the westernmost limit of88-13 abutted against Line 10 at the southeast side of SP-48, and 88-14 and 88-17 provideda crossing of Line 10, with the eastern end of 88-14 coinciding with the refraction line andthe western limit of 88-17 falling ,11 km to the east of it; the locations of the reflectionlines with respect to Line 10 are shown on Figure 1 la. The proximity of these profiles toLine 10 makes it important to consider the nature of the reflection data in both areas.The section of Figure 15a shows that the western part of profile 13 is characterized byeast-dipping structures (indicated by "A") in the upper crust. The middle crust shows zoneswhich dip slightly to the west ("B") and truncate the east-dipping upper crustal reflectors.Nearly-linear features are apparent in the lower crust, as are the east-dipping reflections("C") which are visible at ,-6.0 and 9.5 s underneath the Line 10 location. The pronouncedpackage of reflections from the Moho (denoted with an "M") begins below this at -10.9 sand has a travel time thickness of ,0.3 s; these values vary laterally along the profile, withthe travel time thickness reaching up to "0.8 s.The overall crustal fabric for profile 88-17 (Figure 15b) is again one of upper crustaleast-dipping reflectors ("A") truncated by west-dipping mid-crustal reflections ("B") whichoverlie flat reflectors in the lower crust. Prominent crustal reflections occur at the Line 10site on the 88-14 reflection section in zones of enhanced reflectivity between approximately1.5 and 3.5 s and 6.0 and 7.0 s, and at roughly 10.0 s (each of these is indicated with an"*"); all seem continuous with reflections in the 88-17 data. Reflections from the Moho arenot visible for 88-14, but for 88-17 they again appear as a package of prominent reflections,with a travel time of , 10.8 s at the top beneath Line 10 (denoted with an "M"). This packageof Moho reflections begins at -10.8 to 11.0 s across the profile, and is of variable thickness.Line 10FIGURE 15:(a) LITHOPROBE reflection profile 88-13, migrated and coherency filtered (Varsek et al. 1992), with the Line 10 locationindicated at the surface. "A" indicates upper crustal east-dipping reflectors, "B" mid-crustal west-dipping reflectors, "C" east-dippingreflections embedded in the generally flat reflections of the lower crust, and "M" the Moho reflections.''''^-''.*:-...2",': . - ,.^., '.. ,^■-7- 7,--'-..--':5...:.-.„ ----•-•%--,-:.--i-..-----i---2,----,--7_..r-..m.:_....._•-•■•4•'...4'.,‘, •n-,—,__.. :.-..7,:k....o..-."--,-.-.-..--7,"'7:2;F:::-. --..%--^-_',' • •\, . . - .^.:--- '. '---`....,,..4.,--F-2.t. , ...-.,--- ..,--.' • ,.._-''. ,,---,"---,•*- -.-'-- - 7^.1:-.^-, -' "--,---. ,_^- ..., --_,-.''''..,---- .%-•..--,.-'"._,-.-.-^.-r-_,-k'-,•,,^-. ..'".:--"----";:,-ir-F:.;_"--^:•••-•-•----._--..-..-,---.,--,...--,----,--•:-" -..,..--"..7,..,---7",.....--------., ----••:,_^..,.• - -^_•2468101214MLine 10BKFS50 km(b)FIGURE 15:(b) LITHOPROBE reflection profiles 88-14 and 88-17, migrated and coherency filtered (Varsek et al. 1992), with the Line 10location indicated at the surface. BKFS indicates the Bralome-Kwoiek Creek fault system, and "*" zones of enhanced reflectivity discussedin the text; "A", "B", and "M" are as defined in (a).75Many of the upper crustal structures visible in the reflection sections undoubtedly arosefrom the Mesozoic crustal shortening adjustments that accommodated terrane accretion. Thesimilarity of the data interpretations (i.e. Monger and Journeay 1992, Varsek et al. 1992)for the upper crust have shown that observed seismic reflectors may indeed be successfullycorrelated with the surface expressions of terrane boundaries and associated terrane properties.However, the extensions of the upper crustal reflection interpretations to lower crustal depthshave not been certain; this is a point to which this discussion returns.Figure 16a shows (i) the 1—D velocity profile of the Line 10 interpretation at the locationof profile 88-13 plotted versus two-way travel time, and (ii) a simplified schematic of the88-13 reflection interpretation of Monger and Journeay (1992) at the location of Line 10.Figure 16b is the counterpart of Figure 16a, the difference being that part (ii) represents theinterpretation of Varsek et al. (1992). Figure 17a shows (i) a simplified schematic of the88-14 reflection interpretation of Monger and Journeay (1992) at the location of Line 10,(ii) the 1—D velocity profile plotted versus two-way travel time of the Line 10 interpretationat the location of the 88-14/88-17 crossing of Line 10, and (iii) a simplified schematic ofthe 88-17 reflection interpretation of Monger and Journeay (1992) at the location of Line10. Figure 17b is the counterpart of Figure 17a, with parts (i) and (ii) representing the workof Varsek et al. (1992). Monger and Journeay (1992) did not define a terrane affiliation forall units; these have been added to Figures 16a and 17a where it was deemed appropriate.The conversion of the final model from depth to two-way travel time for Line 10 velocitieswas accomplished by the integration of reciprocal velocities from the surface downward.Since the refraction model changes very slightly between its intersection with 88-14 and itsintersection with the projection of 88-17 (see Figure 11a), and structural continuity betweenthe southwest end of profile 17 and profile 14 is clear from the reflection data, Figures 17aand 17b show a 2—D representation of the crossing.The first important implication of the Figure 16 and Figure 17 presentations is thatFIGURE 16: Figure 16a shows (i) the 1—D velocity profile plotted versus two-way travel time of theLine 10 interpretation at the location of profile 88-13, and (ii) a simplified schematic of the 88-13reflection interpretation of Monger and Joumeay (1992) at the location of Line 10. Figure 16b is thesame as Figure 16a, except that part (ii) represents the interpretation of Varsek et al. (1992). "Wr"indicates that the Insular superterrane exists as Wrangellia in this region.CadwalladerTerrane1-Gambier(OverlapAssemblage upperGambierOverlapAssemblageInsularSuperterrane idc le-------- R2IntermontaneSuperterraneu16EFlntermontaneSuperterrane101214UpperMantle) (iii)88-17 from Monger andJoumeay (1992)Line 10(a) 88-14(iii)88-17 from Varsek et al.(1992)Cadwallader TerraneLine 10Gambier OverlapAssemblageVelocity (knits)3 4 5 6 7 8 9Line 10Gambier OverlapAssemblageInsularSuperterrane— — — — — —6E8101214(b) 88-14Shuksan TerraneShuksan TerraneUpper Mantle77Velocity (km/s)3 4 5 6 7 8 91^I^I^IJ Line 10F13lower a ist78the conversions of the final model from depth to two-way travel time allow the refraction-imaged reflections to be compared with the reflection data shown in Figure 15. Although thereflection lines intersected Line 10 in the cases of 88-13 and 88-14, a precise correspondenceis not expected due to the following points: (i) the possibility of crustal anisotropy; (ii) thepossibility of lateral averaging in the refraction velocity model; (iii) a 3-D effect is introducedto the 88-13 case due to its intersection point being located on the only portion of Line 10which veers significantly from the straight line between its endpoints; and (iv) the quality ofthe reflection data suffers at the extreme ends due to the loss of fold.At the intersection of Line 10 and 88-13, R1 is not constrained (see Figure 11a); thisis consistent with the fact that there is no definite individual reflection imaged at this level(,4.0 s) in the profile 13 data (Figure 15a). The R2 and R3 reflections were modelled at-6.7 and 9.7 s, respectively (see Figure 16a (i) or Figure 16b (i)), and correspond well withstructures visible in the reflection data (Figure 15a) and described above. The transitionalrefraction Moho was modelled at -10.8 s which roughly coincides with the reflection Moho;further discussion of this comparison appears later in this section.At the location of the 88-14/88-17 crossing of Line 10, R1 is again not constrained (seeFigure 1 la); however, an event in the reflection data occurs at N3.5 s (Figure 15b), the time ofRi. Varsek et al. (1992) have interpreted this reflecting horizon as the subsurface continuationof the northern extension of the Thomas Lake fault (Figure 4); the rocks juxtaposed atthis boundary are of theological properties that could not be resolved with the refractionexperiment. The R2 and R3 reflectors were imaged at r-7.2 and 9.8 s, respectively, by theLine 10 data (see Figure 17a (ii) or Figure 17b (ii)). The reflection data (Figure 15b) showthat the base of a high-reflectivity zone appears at '7.0 s on both the 88-14 and 88-17 sidesof the crossing, and that a strong event occurs similarly at 9.8 s.The second significance of Figures 16 and 17 is that the differences of the Mongerand Journeay (1992) and the Varsek et al. (1992) models show the tenuous nature of the79interpretations of deeper (i.e. below .-3.5 s) reflections, so that while the upper crustalinterpretations are alike, below this the solutions diverge. For example, the source of thereflection corresponding R2 at the westernmost end of 88-13 is identified in Varsek et al.(1992) (Figure 16b (ii)) as the top of the "Insular ramp", a proposed unit containing asuture separating the Insular and the Intermontane superterranes; a very different wedge-likestructure appears in the Monger and Journeay (1992) version (Figure 16a (ii)). Also, whilethe Monger and Journeay (1992) interpretation does not put forward a structural origin forthe R3 reflector, the Varsek et al. (1992) model shows R3 to be generated at the bottomof the Insular ramp.Further examples appear in the comparison of Monger and Journeay's (1992) interpre-tation of the Line 10 crossing of 88-14 and 88-17 (Figure 17a (i) and (iii)) with that ofVarsek et al. (1992) (Figure 17b (i) and (iii)). This comparison shows that the former hasthe R2 reflecting horizon embedded in an Intermontane unit, while the latter has it withinthe Insular superterrane. And while Varsek et al. (1992) did not interpret a structural sourcefor the R3 reflection, Monger and Journeay (1992) have shown this as the top of a thrustramp within the Intermontane superterrane.The crustal structures shown in Figures 16 and 17 do show an important similarity: thejuxtaposition of Insular and InterMontane materials. As well, the presence of rocks classifiedwith the Shuksan terrane, which represents the vestiges of an inter-terrane ocean basin,appears in the Varsek et al. (1992) interpretation. This shows that both of the reflectionstudies have concluded that the Line 10 region was the site of the collision zone betweenthe Insular and Intermontane superterranes, as was the case with the Zelt et al. (1992b)refraction study of SCoRE '89 Line 3.It is not surprising that the interpretation of seismic reflection data across such astructurally complex features has lead to ambiguities in linking deep reflections with knownstructural features; in fact, the anticipation of such ambiguities provided part of the motivation80for a study designed in the manner of SCoRE '90. In the region of this study, such ambiguitieshave led to fundamentally different views of the nature of the superterrane collision whichproduced the CB. Varsek et al. (1992) have constructed a collisional model based uponcrustal imbrication, as they identified the source of reflectivity as the lithological layeringresulting from the interfingering of terranes. The Monger and Journeay (1992) model isone of crustal delamination, whereby the sub-horizontal reflections characterizing the middleand lower crust were interpreted as deep counterparts to east-vergent faults identified in theshallower crust to the east (i.e. the Bralorne-Kwoiek Creek fault system, Figure 15b), withWrangellia being dispaced along this fault system over the top of the terrane material below.Figure 18a shows the interpretation of Monger and Journeay (1992) and Figure 18b that ofVarsek et al. (1992) over the entire width of the reflection profiles.A fundamental difference in these two viewpoints is the extent to which Wrangellia,whose easternmost surface expressions are shown as x 's in Figure 4, has penetrated theCordillera. The Monger and Journeay (1992) proposal has Wrangellia being pushed upwardover the Intermontane superterrane, so that all that remains of it below Line 10 is a wedgeoverlying the rest of the middle and upper crust (Figures 16a and 17a; Figure 18a). Theimbrication model of Varsek et al. (1992) results in Wrangellia existing much farther to theeast in the middle and lower crust, so that thick fingers of it are proposed beneath the Line10 location (Figures 16b and 17b; Figure 18b). The comparison of refraction results in theprevious section has shown that if Wrangellia exists under Line 10, it must be limited tozones whose extent could not be resolved; thus, the crustal delamination model of Mongerand Journeay (1992) is the one which is consistent with results derived in the present study.Third and finally, the summary provided in Figures 16 and 17 allows the correspondenceof the refraction-imaged transitional Moho (layer 9 in Figure 12) with the crust-mantleboundaries in the reflection data to be verified. The Moho has been defined by Steinhart(1967) as "that level in the Earth where compressional wave velocity increased rapidly orLine 10'^-Mc......./.'.-• '. -"41..C .1"4 ...-•4 4--:.*.0.1.,y_"•-•'----':'-114.- .^-.•^•^•....' • • • • ..: '^-.... - :45.....,.:•ii, •-7 -fIN7=-- Insular k- -_----,--7.7-,:-C-=----. .(a) -^e- 6,•••-' -v. ;- c:C- 8• •.:" 7-V7-i-••^— • —^- from Monger and Journeay (1992)—^•2• '^•—^- — •-••^—•^s‘s•-••••-•. •^ ••EI•cmt--.^^-84cC1214from Varsek et al. (1992)8110121450 kmFIGURE 18: Interpretations of the reflection data from the 88-14/88-17 crossing of Line 10 by (a)Monger and Joumeay (1992) and (b) Varsek et al. (1992) overlie the reflection data (see Figure 15).Black lines indicate reflections which have been interpreted in the studies, the dotted area representsInsular superterrane, and the line-shaded area represents Intennontane superterrane (see Figures 16and 17 for more detail in the region of Line 10). Where appropriate, some identifications missingfrom the original studies have been added.82discontinuously to a value between 7.6 and 8.6 km/s"; thus, the Moho must by definition bedetermined through refraction studies, which yield velocity information. For reflection data,the Moho has been identified by Klemperer et al. (1986) as the "deepest, high-amplitude,laterally extensive reflection or group of reflections". The question of how the reflectionMoho relates to the conventional definition remains one of current debate (Jarchow andThompson, 1989).The duration of the zones of Moho reflectivity on the 88-13 (Figure 15a) and 88-17(Figure 15b) sections are consistent with the interpretation of the crust-mantle boundary asan alternating series. of high and low velocity material. As shown in the seismic section of88-13 (Figure 15a), the pronounced zone of Moho reflectivity at the site of Line 10 is between, 10.9 and 11.2s. The top of this zone coincides well with the top of the transitional refractionMoho, which was modelled at -40.8 s. The travel time thickness of s observed forthese well-defined reflections is equivalent to -4.9 km, which may be an indication of thethickness of the transitional Moho. The package of Moho reflections visible on the profilefrom 88-17 (Figure 15b) also begins at a travel time (-41.0 s) which corresponds to the topof the refraction-imaged Moho layer (-41.0 s); the thickness of the reflectivity package in thiscase is more difficult to discern. Monger and Journeay (1992) and Varsek et al. (1992) placedthe Moho at the position of the deepest prominent reflection; this interpretation coincideswith the base of the refraction-imaged transitional Moho, if its thickness is assumed to beabout 2.0 km as inferred from the well-defined 88-13 Moho reflections (Figure 15a).The correlation of the top of the refraction-imaged transitional Moho with the top of thezone of prominent Moho reflectivity observed in the 88-13 and 88-17 sections is suggestedby the data; similarly, the base of this transitional layer may correspond to the deepest extentof this Moho reflectivity, although it is important to note that these observations are subjectto the refraction modelling error (Table 7).835.2.2 Heat Flow Data and InterpretationsHeat flow in the broad region of subduction zones characteristically follows a pattern:a band of low heat flow reaches from the trench to the volcanic arc about 200 km inland,and a much higher than normal heat flow exists for a great distance inland (Keen andHyndman 1979). The heat flow measured across southwestern British Columbia in an areaencompassing the LITHOPROBE Southern Cordillera Transect site exhibits this characteristicpattern (Lewis et al. 1988, 1992). The data show high heat flow values for most of the CB(-75-95 mW/m2), and values almost as high for the Intermontane Belt; in comparison, heatflow values are low (<50 mW/m 2) throughout the Insular Belt.In Phanerozoic regions, including the southern Canadian Cordillera, the lower crusttypically exhibits bands of subhorizontal reflectors. Klemperer (1987), Weyer et al. (1987)and others have indicated that a correlation exists between the depths to the tops of thesecharacteristic bands and regional heat flow. Heat flow studies such as these have commonlyemphasized the depths to the 450 and 730 °C isotherms, which are identified respectivelyas the approximate temperatures of the brittle-ductile transition zone in the crust, and thetransition from a mineralogy in equilibrium with coexisting free water to a dry mineralogy. Ina number of Cordilleran reflection sections, the 450 °C isotherm has corresponded to the topof the characteristic, sub-horizontal reflective bands; the 730 °C isotherm has correspondedto the bottom of these characteristic bands (Lewis et al. 1992). The origin of the reflectionshas been the subject of much discussion (e.g. Warner 1990) and support has been providedfor the model of layered porosity, whereby the region between the transition zones describedabove is identified as the portion of the crust where thin pores can stay open, and fluidtrapped within the pores generates high reflectivity coefficients if the porosity is layered. Inthe LITHOPROBE reflection data from the Line 10 region, however, the top and bottom ofsuch reflective bands are not distinguishable, as a strong reflectivity is present throughoutthe entire crust (Figures 15a and 15b).84Lewis et al. (1992) have shown that the level of heat flow in the Line 10 region was-94 mW/m2 ; corresponding crustal temperatures suggested a depth to the 450 °C isothermof ,10.8 km and to the 730 °C of ,20.0 km; these depths coincide well with the depthsof the R1 and R2 reflectors imaged in the refraction data (see Figure 11a). In the case ofthe R1 reflecting horizon, although its depth (-40.0 km) correlates well with that of the 450°C isotherm, its sporadic nature (Figure 11a) does not closely resemble continuous reflectiondescribed by Lewis et al. (1992) as coming from the top of a zone of layered porosity in thecrust. Also, the fact that the regions where R1 is constrained are equivalent to the portions ofLine 10 which crossed highly-faulted regions (Wheeler and McFeely 1991) suggests faultingas a more likely source of the reflectivity. In terms of the R2 horizon, its depth (-22.0 km),continuous, nearly linear nature, and the relatively small velocity contrast across it (0.05-0.15km/s) support the idea of its correspondence with the transition zone linked with the 730 °Cisotherm. However, the reflective nature of the crust below this depth, which is visible onboth the refraction data (i.e. R3 in Figures 5-10) and the reflection data (Figures 15a and15b), means that the abrupt decrease in reflectivity associated with this transition zone cannotbe uniquely identified, so it is not appropriate to attempt a definitive correlation with the R2reflecting horizon. Nevertheless, the corresponding depths of the R2 horizon and the 730 °Cisotherm do suggest that the interpreted structural boundaries „should be regarded with caution.5.2.3 Gravity Data and InterpretationsThe Bouguer gravity anomaly map for the study area is shown in Figure 19 (fromRiddihough and Seeman 1982); the map shows that this is a region of low gravity values.The gravity field across nearly all continental margins where subduction is taking place alsoshows a characteristic pattern (Keen and Hyndman 1979). This is typified by a band ofvery low gravity over the subduction trench, and a parallel band of high gravity extendinginland for approximately 100 km. Much farther inland, there is usually a large region of lowBouguer gravity; this corresponds to the region of the present study. Riddihough (1979)85FIGURE 19: The Bouguer gravity anomaly map for the study area (from Riddihough and Seeman1982), with contours of 10 mgal overlying the topography. The location of Line 10 is indicated, andthe high gravity anomaly at the southeastern part of the profile is denoted with an "H".86pointed out that this low gravity could be attributed to low mantle density, possibly causedby high temperatures related to the high heat flow in the area.The gravity data coincident with Line 10 show that the profile is one of high negativeanomalies, with lower negative values on the ends (-124 mgal in the northwest, and —87 mgalin the southeast) than in the middle, where values reach -200 mgal. Using the Barton (1986)velocity-density relationship, the final velocity model was converted to a density model forthe purpose of gravity modelling. [Barton (1986) gave the mean, and the maximum andminimum bounds, of the measurements used for the conventional Nafe and Drake (1963)relationship; as well, Barton suggests slightly higher densities for seismic velocities greaterthan 6.5 km/s.]The algorithm available for gravity modelling was the interactive Fortran program SAKI(Webring 1985), which calculates and compares theoretical gravity responses of a structuralmodel with profiles of observed data. The input model was constructed as a set of polygonalprisms of finite lateral extent, with the resulting format being 2.5—D. Gravity modelling in2.5—D is suitable for the extrapolation of values to either side of a line; it is not suitablefor modelling off-line changes. The application of the 2.5—D SAKI modelling algorithmto the data set from Line 10 was therefore limited, as the profile ran along the strike of amajor geological feature, with quite different structures flanking it. Testing indicated that thevelocity model could be converted into a density model which generally provided an accuratefit to the gravity data. The exception occurred in the modelling of the southeastern end ofthe profile, where out-of-plane effects impossible to account for with the 2.5—D modellingprogram rendered the algorithm inappropriate (Figure 19).Although the gravity modelling was not applicable, the increase in gravity at thesoutheastern portion of the line (Figure 19) must at least have a qualitative interpretation.Since the crustal thickness along Line 10 is relatively constant (Figure 11a), a possibleexplanation for the observed increase in values is that the densities of the rocks in this area87increase. The geological map of Wheeler and McFeely (1991) shows that the nature ofthe plutonic rocks which dominate the surroundings of Line 10 (Figure 4) changes in thissoutheastern region: the surface rocks are predominantly Mid-Cretaceous quartz diorites, asopposed to the Early Tertiary granodiorites/quartz diorites and Jura-Cretaceous quartz dioritepackages which predominate at the northwest and central portions of the line, respectively.The gravity study of Dehler (1991) included the measured densities of rock samples fromthe CB. The wide range of density values observed per rock unit (e.g. 2600-2900 kg/m 3 forquartz diorite) lends plausibility to this possible interpretation.5.2.4 Geomagnetic Data and InterpretationsNatural geomagnetic induction studies (e.g. geomagnetic deep soundings, Caner 1970;magnetotellurics, Gough and Majorowicz 1991, Jones et al. 1992) provide information aboutthe electrical conductivity of the crust and upper mantle. The early work of Caner (1970)showed that the entire region of the Cordillera in southern British Columbia to be one of highconductivity, or low resistivity. This was explained in the magnetotelluric (MT) investigationof Gough and Majorowicz (1991) as being caused by the wetness of Cordilleran rocks, whichis linked in part with the ongoing subduction in the region.The Gough and Majorowicz (1991) work comprises one of two recent MT modellingstudies in the general region of Line 10, the other being that of Jones et al. (1992). Bothstudies have interpreted subsurface structure in terms of electrical resistivity versus depth,and have shown an upper crust characterized by zones of high resistivity (i.e. > 300 S2•m)overlying a conducting middle and lower crust; the resistivity variations have been explainedas resulting from variations in fracture densities. The resistive zones in the upper crust havebeen identified as granitoid plutons, with low fracture densities.The depth extent of the batholiths interpreted in the two MT studies ranges between10 and 22 km. The vertical homogeneity of the Line 10 velocity model at this depth levelsuggests that the extent of the CB batholiths has not been determined using the refraction88data. This does not eliminate the possibility of the change in rock properties at the limitof a batholith generating reflectivity, since reflections identified as the R1 and the R2 phasecome from depths ("40.0 and 22.0 km, respectively) that could correspond with the base ofMT-imaged batholiths. A good example of this is provided by the intersection point of one ofthe Gough and Majorowicz (1991) MT profiles with Line 10, which occurs at approximatelythe site of SP-48 (the "E" site on their profile). Gough and Majorowicz image a resistivestructure of which they say there is "little doubt that the ... mass ... is the main bulk of theCoast Plutonic granodiorites"; the structure reaches a depth of ,22.0 km, which correspondswith the depth of the R2 reflecting horizon at this location ( ,21.0 km). Such a scenario isless likely for the R1 reflections, however, since the northwestern region where most of theconstraint on the R1 phase was located (Figure 11a) corresponds to the region of Line 10where the surface rocks were terranes, as opposed to large granitic plutons (Figure 4).5.3 ConclusionsThrough the application of successive iterations of inverse travel time modelling andforward amplitude modelling to the refraction data from Line 10 (Figures 2 and 4), a velocitystructure was developed for the central Coast Belt (Figure 11a). The model features a thin(,2 km) near-surface layer with a velocity of 3.9 km/s at the top and a high (0.5 s -1 )velocity gradient. Zones of high and low velocity overprint this average velocity structure,with the most pronounced high velocity zone possibly corresponding to a geological feature(the Cadwallader terrane). The surface layer overlies three major crustal strata: (i) the uppercrust, which has an average depth extent of 10 km, velocity of 6.2 km/s, and velocity gradientof 0.02 s-1 ; (ii) the middle crust, which extends to 22 km depth and has a velocity of 6.35km/s; and (iii) the lower crust, which consists of an 8 km thick upper unit with a velocityof 6.5 km/s overlying a unit of 4 km thickness and 6.7 km/s velocity. The crustal packagedips slightly downward to the southeast, with the average depth to the Moho being 34.5 km.In order that synthetic seismograms provided a close match to the observed data, the Moho89was modelled as a thin transition layer, with velocities varying laterally from -7.40 to 7.80km/s. The upper mantle velocity ranged from 8.05 to 8.15 km/s in accordance with a Mohothicknesses of between 0 to 2 km (0 km thickness represents a first-order discontinuity.) Thisupper mantle velocity estimate must be regarded with caution, however, as the observed dataused for the inversion were sparse (only 33 observations), and the refraction data do notconstrain the thickness of the transitional Moho. Furthermore, this relatively high velocityestimate is not consistent with the high heat flow, low gravity, and low electrical resistivityobserved in the area. An upper mantle reflection was imaged at a depth ranging from -66to 70 km, in keeping with a Moho thickness of 0 to 2 km.The reflection Moho observed on each of the 88-13 and 88-17 sections appears asa band of enhanced reflectivity with two-way travel time thickness varying from r-0.2 to0.8s (or equivalently ,0.7 to 3.0 km), with an upper extent corresponding with the topof the refraction-imaged Moho. This observation is consistent with the interpretation ofthe refraction Moho as a band of alternating high and low velocity materials (Mooney andBrocher 1987), and suggests that its thickness, which is not constrained by the refractiondata, could be constrained by the thickness of the bands. The Moho thickness of 2 kmshown in the final model (Figure 11a) is consistent with the thickness of the Moho at theLine 10/88-13 intersection point, where the thickness of the pronounced band of reflectivityis clearly-defined.In general, the lateral velocity variations determined along the refraction line are minimal— this may be related to the approximately along-strike orientation of the profile (Figure 2,inset). Despite this lateral velocity homogeneity, the limited study of regional gravity hasindicated that the density of the materials below Line 10 may increase in the southeast. Thismay be consistent with a change in surface geology in this region (Wheeler and McFeely1991), in which the age of the plutonic rocks which dominate the region of Line 10 becomesMid-Cretaceous, as opposed to the predominantly Jura-Cretaceous age of plutons to the north.90The strongest wide-angle reflections present in the SCoRE '90 Line 10 data are generatedfrom the Moho; the origin of the three weaker crustal reflections, however, cannot beascertained through the refraction data alone. The comparisons with other geophysical studiesmade in the preceding section have indicated that three possible types of reflecting horizonsmay exist at the depth of those imaged in the crust: (i) structural features (e.g. fault zones),(ii) the transition zones at the top and bottom of a region of layered porosity within the crust,and (iii) lithological contrasts (e.g. the boundary of a batholith and the surrounding terranematerial). It is probable that the source of the R1 reflections was fault zones; the sourceof the R2 reflections could not be uniquely determined, and the only possible origin of theR3 reflections uncovered in this study are deep thrust faults which have been linked to thecollision of the Insular and Intermontane superterranes (Monger and Journeay 1992). Theupper mantle reflection R may represent the base of the lithosphere.Interpretations of seismic profiles which crossed, or nearly crossed, Line 10 have beenprovided by the SCoRE'89 Line 3 refraction study of Zelt et al. (1992b), and the Mongerand Journeay (1992) and Varsek et al. (1992) studies of LITHOPROBE reflection profiles88-13, 88-14, and 88-17 (Figure 2). These studies have determined that the location ofLine 10 coincides with the collision zone between the Insular and the Intermontane terranes,a geological suture whose surface expression has been obscured by the granitic intrusionsof the Coast Belt.Recent refraction interpretations for the southern Cordillera regions surrounding Line10 include the Spence et al. (1985) study of VISP '80 Line I of the Insular Belt, and theZelt et al. (1992a) modelling of SCoRE '89 data from the Intermontane Belt, as well asthe study of Zelt et al. (1992b) (Figure 2). Spence et al. (1985) and Zelt et al. (1992b)have determined velocities for the Insular superterrane — which exists as Wrangellia in thisregion — that are significantly higher than those interpreted for the Line 10 profile (Tables 9and 10, respectively); in contrast, the velocities interpreted by Zelt et al. (1992a) and Zelt et91al. (1992b) for the Intermontane superterrane are very similar to the Line 10 values (Tables9 and 10). Although the presence of Wrangellia is to be expected in the proposed Insular-Intermontane collision zone under Line 10, the refraction velocity model (Figure 11a) showsthat it must be limited to zones whose extent cannot be resolved using the refraction data.Thus, the results of this study suggest that Line 10 represents approximately the easternmostextent reached by Wrangellia in its complex emplacement against the ancient margin of NorthAmerica, an idea which is expressed in the interpretation shown in Figure 1 lb.The collisional model of Monger and Journeay (1992) involves the upward displacementof Wrangellia over the Intermontane rocks of the ancient North American margin, whilethat of Varsek et al. (1992) consists of the imbrication of Insular and Intermontane material(Figures 18a and 18b). Since the extension of Wrangellia to the east of Line 10 in the Varseket al. (1992) model entails the local presence of thick Wrangellian layers, the Monger andJourneay (1992) interpretation which implies its presence as a wedge confined to the uppermiddle crust is favored by the results of this study.92BIBLIOGRAPHYAM, K. and RICHARDS, P.G. 1980. Quantitative Seismology, vol. 2, W.H. Freeman, NewYork.BARTON, P.J. 1986. The relationship between seismic velocity and density in the continentalcrust — a useful constraint? Geophysical Journal of the Royal Astronomical Society, 87:195-208.BENZ, H.M., SMITH, R.B., and MOONEY, W.D. 1990. Crustal structure of the northwestern Basin and Range Province from the 1986 Program for Array Seismic Studies ofthe Continental Lithospheric seismic experiment. Journal of Geophysical Research, 95:21823-21842.BERRY, M.J. and FORSYTH, D.A. 1975. Structure of the Canadian Cordillera from seismicrefraction and other data. Canadian Journal of Earth Sciences, 12: 182-208.BOLT, B.A., and GUTDEUTSCH, R. 1982. Reinterpretation by ray tracing of a transverseseismic profile through the California Sierra Nevada, Part I. Bulletin of the SeismologicalSociety of America, 72: 889-900.BRAILE, L.W. 1973. Inversion of crustal seismic refraction and reflection data. Journal ofGeophysical Research, 78: 7738-7744.BRAILE, L.W. and CHIANG, C.S. 1986. The continental Mohorovi6i6 discontinuity: Re-sults from near-vertical and wide-angle seismic reflection studies. In Reflection Seismol-ogy: A Global Perspective, Geodynamics Series, edited by M. Barazangi and L. Brown,Washington, D.C., vol.13, pp. 257-272.BURIANYK, M.J.A., KANASEWICH, E.R., ELLIS, R.M., and CLOWES, R.M. 1992.SCoRE '90: Southern Cordillera Refraction Experiment Field Report and Data SetDescription. LITHOPROBE Report No. 21, LITHOPROBE Secretariat, University ofBritish Columbia, Vancouver, 30 pp.CANER, B., 1970. Electrical conductivity in western Canada and petrological implications.Journal of Geomagnetism and Geoelectricity, 22: 113-129tERVENt V., MOLOTKOV, I., and PkNefK, I. 1977. Ray method in seismology.Charles University Press, Prague.CHIU, S.K.L., and STEWART, R.R. 1987. Tomographic determination of three-dimensionalseismic velocity structure using well logs, vertical seismic profiles, and surface seismicdata. Geophysics, 52, 1085-1098.CLOWES, R.M. (editor) 1989. LITHOPROBE Phase III Proposal — The Evolution of aContinent, LITHOPROBE Secretariat, The University of British Columbia, Vancouver.213 pp.CLOWES, R.M., COOK, F.A., GREEN, A.G., KEEN, C.E., LUDDEN, J.N., PERCIVAL,J.A., QUINLAN, G.A., and WEST, G.F. 1992. LITHOPROBE — New perspectives oncrustal evolution. Canadian Journal of Earth Sciences. (In press.)93COLES, R.L., and CURRIE, R.G. 1977. Magnetic anomalies and rock magnetizations in thesouthern Coast Mountains, British Columbia: possible relation to subduction. CanadianJournal of Earth Sciences, 14: 1753-1770.CONEY, P.J. JONES, D.L., and MONGER, J.W.H. 1980. Cordilleran suspect terranes.Nature, 288: 329-333.COOK, F.A., GREEN. A.G., SIMONY, P.S., PRICE, R.A., PARRISH, R.M., MILKEREIT,B., GORDY, P.L., BROWN, R.L., COFLIN, K.C., and PATENAUDE, C., 1988. LITHO-PROBE seismic reflection structure of the southeastern Canadian Cordillera: initial results.Tectonics, 7:157-180.COOK, F.A., VARSEK, J.L., CLOWES, R.M., KANASEWICH, E.R., SPENCER, S.C.,PARRISH, R.R., BROWN, R.L., CARR, S.D., JOHNSON, B.J. and PRICE, R.A. 1991.LITHOPROBE crustal reflection cross section of the Southern Canadian Cordillera: Fore-land thrust and fold belt to Fraser River fault. Tectonics. (In press.)COOK, F.A., VARSEK, J.L., CLOWES, R.M., KANASEWICH, E.R., SPENCER, C.S.,PARRISH, R.R., BROWN, R.L., CARR, S.D., JOHNSON, B.J., and PRICE, R.A. 1992.LITHOPROBE crustal reflection cross section of the southern Canadian Cordillera I:foreland thrust and fold belt to Fraser River fault. Tectonics, 11: 12-35.CUMMING, W.B., CLOWES, R.M., and ELLIS, R.M. 1979. Crustal structure from aseismic refraction profile across southern British Columbia. Canadian Journal of EarthSciences, 16: 1024-1040.DICKINSON, W.R. 1970. Relations of andesites, granites, and derivative sandstones toarc-trench tectonics. Review of Geophysics and Space Physics, 8: 813-860.DAVIS, G.A., J.W.H. MONGER, and BURCHFIEL, B.C. 1978. Mesozoic construction of theCordilleran "collage", central British Columbia to California. In Mesozoic paleographyof the Western United states, Pacific Coast Paleography Symposium 2, edited by D.G.Howell and K.A. McDougall, Society of Economic Paleontologists and Mineralogists,Pacific Section, Los Angeles, California, pp.1-32.DEHLER, S.A. 1991. Intergrated geophysical modelling of the northern Cascadia subductionzone. Ph.D. thesis, University of British Columbia, 151 pp.DREW, J.J., and CLOWES, R.M. 1990. A re-interpretation of the seismic structure across theactive subduction zone of western Canada. In Studies of laterally heterogeneous structuresusing seismic refraction and reflection data, edited by A.G. Green, Geological Survey ofCanada, Paper 89-13, pp.115-132.FRIEDMAN, R.M. and ARMSTRONG, R.L. 1992. Jurassic and Cretaceous geochronologyof the southern Coast Belt, B.C., 49°-51°N. Geological Society of America Memoir.(Submitted.)GABRIELSE, H., and YORATH, C.J. 1989. DNAG #4. The Cordilleran Orogen in Canada.Geoscience Canada, 16: 67-83.94GEHRELS, G.E., and SALEEBY, J.B. 1985. Constraints and 'speculations on the displace-ment and accretionary history of the Alexander-Wrangellia-Peninsular superterrane {ab-stract}. Geological Society of America, Abstracts with Programs, 17: 356.GOUGH, D.I. 1986. Mantle upflow mechanics in the Canadian Cordillera. Journal ofGeophysical Research, 91: 1909-1919.GOUGH, D.I. and MAJOROWICZ 1992. Magnetotelluric soundings, structure, and fluidsin the southern Canadian Cordillera. Canadian Journal of Earth Sciences, 29: 609-620.HALE, L.D. and THOMPSON, G.A. 1982. The seismic reflection character of the continentalMohorovi6id discontinuity. Journal of Geophysical Research, 87: 4625-4635.JARCHOW, C.M., and THOMPSON, G.A. 1989. The nature of the Mohorovitid disconti-nuity. Annual Review of Earth and Planetary Sciences, 17: 475-506.JOHNSON, S.H., COUCH, R.W., GEMPERLE, M., and BANKS, E.R. 1972. Seismicrefraction measurements in southeast Alaska and western British Columbia. CanadianJournal of Earth Sciences, 9: 1756-1765.JONES, A.G., GOUGH, D.I., KURTZ, R.D., DeLAURIER, J.M., BOERNER, D.E.,CRAVEN, J.A., ELLIS, R.G. and McNIECE, G.W. 1992. Electromagnetic images ofregional structure in the southern Canadian Cordillera. Geophysical Research Letters. (Inpress.)JOURNEAY, J.M. 1990. Structural and tectonic framework of the southern Coast Belt,British Columbia. In Current research, part E. Geological Survey of Canada, Paper90-1E, pp. 183-197.JOURNEAY, J.M., and FRIEDMAN, R.M. 1992. Structural setting and geochronometry ofthe Lillooet River Fault system, southwestern British Columbia. Tectonics. (In press.)HOUGH, S.E., and ANDERSON, J.G. 1988. High-frequency spectra observed at Anza,California: implications for Q structure. Bulletin of the Seismological Society of America,78: 692-707.KANASEWICH, E.R. 1981. Time Sequence Analysis in Geophysics, third edition: Univer-sity of Alberta Press, Edmonton.KANASEWICH, E.R., and CHIU, S.K.L. 1985. Least-squares inversion of spatial seismicrefraction data. Bulletin of the Seismological Society of America, 75: 865-880.KEEN, C.E. and HYNDMAN, R.D. 1979. Geophysical review of the continental marginsof eastern and western Canada. Canadian Journal of Earth Sciences, 16: 712-749.KLEMPERER, S.L., HAUGE, T.A., HAUSER, E.C., OLIVER, J.E., and POTTER, C.J.1986. The Moho in the northern Basin and Range province, Nevada, along the COCORP40°N seismic reflection transect. Geological Society of America Bulletin, 97: 603-618.KLEMPERER, S.L. 1987. A relation between continental heat flow and the seismicreflectivity of the lower crust. Journal of Geophysics, 61: 1-11.LEWIS, T.J., BENTKOWSKI, W.H., DAVIS, E.E., HYNDMAN, R.D., SOUTHER, J.G.and WRIGHT, J.A. 1988. Subduction of the Juan de Fuca Plate: thermal consequences.Journal of Geophysical Research, 93: 15207-15225.95LEWIS, T.J., BENTKOWSKI, W.H., and HYNDMAN, R.D. 1992. Crustal temperaturesnear the Lithoprobe Southern Cordillera Transect. Canadian Journal of Earth Sciences,29: 1197-1214.McLEAN, N.A. 1993. Velocity structure from seismic refraction across the Coast Plutonicand Intermontane Belts. M.Sc. Thesis, University of Victoria. (In preparation.)LUTTER, W.J., NOWACK R.L., and BRAILE, L.W. 1990. Seismic imaging of uppercrustal structure using travel times from the PASSCAL Ouachita experiment. Journal ofGeophysical Research, 95: 4621-4631.McMECHAN, G.A., and SPENCE, G.D. 1983. P-wave velocity structure of the Earth's crustbeneath Vancouver Island. Canadian Journal of Earth Sciences, 20: 742-752.MONGER, J.W.H. 1992. Canadian Cordilleran Tectonics: from Geosynclines to CrustalCollage. Canadian Journal of Earth Sciences. (In press.)MONGER, J.W.H., and PRICE, R.A. 1979. Geodynamic evolution of the CanadianCordillera — progress and problems. Canadian Journal of Earth Sciences, 16: 770-791.MONGER, J.W.H., PRICE, R.A., and TEMPELMAN-KLUIT, D.J. 1982. Tectonic accretionand the origin of the two major metamorphic and plutonic welts in the Canadian Cordillera.Geology, 10: 70-75.MONGER, J.W.H. and JOURNEAY, J.M. 1992. A field guide to accompany the PenroseConference on the Tectonic evolution of the Coast Mountains orogen, May 1992. 97 pp.MOONEY, W.D. and BROCHER, T.M. 1987. Coincident seismic reflection/refractionstudies of the continental lithosphere: a global review. Reviews of Geophysics, 25:723-742.NAFE, J., and DRAKE, C. 1963. Physical properties of marine sediments. In The Sea.Edited by M.N. Hill, pp. 794-828.PAKISER, L.C., and BRUNE, J.N. 1980. Seismic models of the root of the Sierra Nevada.Science, 210: 1088-1094.PHADKE, S., and KANASEWICH, E.R., 1990. Seismic tomography to obtain velocitygradients and three-dimensional structure and its application to reflection data on VancouverIsland. Canadian Journal of Earth Sciences, 27: 104-116.RIDDIHOUGH, R.P. 1979. Gravity and structure of an active margin — British Columbiaand Washington. Canadian Journal of Earth Sciences, 16: 350-363.RIDDIHOUGH, R.P. and SEEMAN, D.A. 1982. Juan de Fuca Plate Map: JFP-8 GravityAnomaly, Pacific Geoscience Center, scale 1:2 000 000.RODDICK, J.A. 1983. Geophysical review and composition of the Coast Plutonic Complex,south of latitude 55°N. In Circum-Pacific Plutonic Terranes. Edited by J.A. Roddick.Geological Society of America, Memoir 159, pp. 195-211.ROGERS, G.C., SPINDLER, C. and HYNDMAN, R.D. 1990. Seismicity along the Van-couver Island Lithoprobe Corridor, in Proceedings of the Project Lithoprobe: SouthernCanadian Cordillera Transect Workshop, Calgary, 3-4 March 1990,166-169.96SHERIFF, R.E. and GELDART, L.P. 1982. Exploration Seismology, v.1, CambridgeUniversity Press, New York.SPENCE, G.D., CLOWES, R.M., and ELLIS, R.M. 1985. Seismic structure across the activesubduction zone of western Canada. Journal of Geophysical Research, 90: 6754-6772.SPENCER, C., ASUDEH, I., and COTE, T., 1989. SEGY-LDS Version 2.0 Format ReferenceDocument, Manual Version 1.00, Geological Survey of Canada, Ottawa. 17 pp.STEINHART, J.S. 1967. Mohorovitie discontinuity. In International Dictionary of Geo-physics, edited by S.K. Runcorn, Oxford, vol.2, pp. 991-994.VAN DER HEYDEN, P. 1992. A middle Jurassic to early Tertiary Andean-Sierran arc modelfor the Coast Belt of British Columbia. Tectonics 11: 82-97.VARSEK, J.L., COOK, F.A., CLOWES, R.M., JOURNEAY, J.M., MONGER, J.W.H.,PARRISH, R.R., KANASEWICH, E.R., and SPENCER, C.S. 1992. LITHOPROBE crustalreflection structure of the southern Canadian Cordillera II: Coast Mountains transect.Tectonics. (In press.)WARNER, M. 1990. Basalts, water, or shear zones in the lower continental crust? Tectono-physics, 173: 163-174.WEBRING, M. 1985. SAKI: A Fortran program for generalized linear inversion of gravityand magnetic profiles. US Geological Survey, Open File Report 85-122.WEVER, T., TRAPPE, H., and MEISSNER, R. 1987. Possible relations between crustalreflectivity, crustal age, heat flow, and viscosity of the continents. Annales Geophysicae,5B: 255-266.WHEELER, J.O. and McFEELY, P. (comp.) 1991. Tectonic Assemblage Map of theCanadian Cordillera and adjacent parts of the United States of America, Geological Surveyof Canada, Map 1712A, scale 1:2 000 000.WHEELER, J.O., BROOKFIELD, A.J., GABRIELSE, H., MONGER, J.W.H., TIPPER,H.W. and WOODWORTH, G.J. (comp.) 1991. Terrane Map of the Canadian Cordillera,Geological Survey of Canada, Map 1713A, scale 1:2 000 000.WICKENS, A.J. 1977. The upper mantle of southern British Columbia. Canadian Journalof Earth Sciences, 14: 1100-1115.ZELT, B.C., ELLIS, R.M., and CLOWES, R.M. 1990. SCoRE '89: The Southern CordilleraRefraction Experiment Description of Dataset. LITHOPROBE Report No. 20, LITHO-PROBE Secretariat, University of British Columbia, Vancouver, 38 pp.ZELT, B.C., ELLIS, R.M., CLOWES, R.M., KANASEWICH, E.R., ASUDEH, I., LUET-GERT, J.H., HAJNAL, Z., IKAMI, A., SPENCE, G.D., and HYNDMAN, R.D. 1992a.Crust and upper mantle velocity structure of the Intermontane belt, southern CanadianCordillera. Canadian Journal of Earth Sciences, 29: 1530-1548.ZELT, B.C., ELLIS, R.M., and CLOWES, R.M. 1992b. Crust and upper mantle velocitystructure in the eastern Insular and southernmost Coast belts, British Columbia, Canada.Canadian Journal of Earth Sciences. (Submitted.)97LUX, C.A., and ELLIS, R.M., 1988. Practical and efficient ray tracing in 2—D mediafor rapid traveltime and amplitude forward modelling. Canadian Journal of ExplorationGeophysics, 24: 16-31.ZELT, C.A., and SMITH, R.B., 1991. Seismic traveltime inversion for 2—D crustal velocitystructure. Geophysical Journal International, 108: 16-34.Appendix A Terrane DescriptionsTerranes of the Insular SuperterraneThe Alexander terrane is composed of Upper Proterozoic to Triassic volcanic and sed-imentary rocks in various depositional settings (including ocean arc, back arc, platform, rift,trough, and offshelf) and comagmatic intrusions (Monger and Journeay 1992). Wrangellia isa submarine arc terrane consisting of Middle to Late Triassic tholeiitic basalts and associatedmafic intrusions, calcareous sedimentary rocks of similar age, and Early to Middle Jurassicterrigenous clastics and volcanics (Coney et al. 1980).Terranes of the Intermontane SuperterraneThe following terranes of the Intermontane Superterrane are of the island arc type. TheStikine terrane includes interbedded arc volcanics, with ages ranging from Lower Devonianto Permian, overlain unconformably by Upper Triassic arc volcanics and granitic rocks, andunconformably overlying Lower Jurassic andesite arc volcanics and intercalated sedimentaryrocks (Monger and Journeay 1992). The Bridge River terrane, where metamorphism isnot predominant, is characterized by extremely disruptive ribbon chert, argillite and basalt,subordinate rocks including siltstones and greywackes, as well as mafic intrusives. The ageof the terrane ranges from Mississippian to late Middle Jurassic, and metamorphic gradeis mostly sub-greenschist to greenschist (Monger and Journeay 1992). The Cache Creekterrane is partly disrupted and partly coherent, and is composed of Permian to Lower andMiddle Jurassic chert and pelite, with melanges of basalt, ultramafic rocks, and carbonatesof Middle Pennsylvanian and Permian age. The Quesnel terrane is comprised of UpperTriassic to Early Jurassic arc volcanics, granitic and alkaline intrusions, and clastic rocks;overlying these unconformably are the clastic rocks (ca. 190-160 Ma) of the AshcroftFormation (Monger and Journeay 1992).98Terranes of the Coast BeltThe Cadwallader terrane is an island arc terrane consisting of Permian-aged pioneergreenstone, diorite, trondhjemite, gabbro, and alpine-type ultramafic rocks, overlain byUpper Triassic arc volcanics, carbonates, and distinctive clast conglomerates (of the HurleyFormation). These rocks are again overlain, by Jura-Cretaceous elastic rocks (Monger andJourneay 1992). The Shuksan terrane is an oceanic terrane which includes Upper Triassicand Lower Jurassic oceanic crust and sediments metamorphosed to greenschist and blueschistfacies and Jurassic near-arc marginal basin crust and sediments. The Harrison terranefeatures Middle Triassic cherry argillites and mafic volcanics overlain by thick Lower Jurassicandesitic and dacitic volcanics. Strata of sedimentary and andesitic volcanic rocks of Middleand Late Jurassic age are also present within this island arc terrane (Monger and Journeay1992).99


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