UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

A geophysical analysis of the Kapuskasing Structural Zone Boland, Andrew V. 1989

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

Item Metadata

Download

Media
831-UBC_1989_A1 B64.pdf [ 23.63MB ]
Metadata
JSON: 831-1.0052973.json
JSON-LD: 831-1.0052973-ld.json
RDF/XML (Pretty): 831-1.0052973-rdf.xml
RDF/JSON: 831-1.0052973-rdf.json
Turtle: 831-1.0052973-turtle.txt
N-Triples: 831-1.0052973-rdf-ntriples.txt
Original Record: 831-1.0052973-source.json
Full Text
831-1.0052973-fulltext.txt
Citation
831-1.0052973.ris

Full Text

A G E O P H Y S I C A L ANALYSIS OF T H E KAPUSKASING STRUCTURAL ZONE  Andrew V . Boland B. A . Mod. (Experimental Physics) Trinity College, Dublin M . Sc. (Applied Geophysics) University College, Galway  A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF T H E REQUIREMENTS FOR T H E D E G R E E OF D O C T O R OF PHILOSOPHY  in T H E FACULTY OF GRADUATE STUDIES GEOPHYSICS AND ASTRONOMY  We accept this thesis as conforming to the required standard  T H E UNIVERSITY OF BRITISH COLUMBIA  September 1989  © Andrew V . Boland, 1989  In  presenting  degree at the  this  thesis  in  University of  partial  fulfilment  of  of  department  this thesis for or  publication of  by  his  or  her  The University of British Columbia Vancouver, Canada  DE-6 (2/88)  It  this thesis for financial gain shall not  Department of  Oct.  representatives.  ^  for  an advanced  Library shall make it  agree that permission for extensive  scholarly purposes may be  permission.  Date  requirements  British Columbia, I agree that the  freely available for reference and study. I further copying  the  is  granted  by the  understood  that  head of copying  my or  be allowed without my written  Abstract  A crustal scale seismic refraction experiment was conducted over the Kapuskasing Structural Zone, Northern Ontario, in 1984. The zone cuts obliquely across the east-west structural grain of the Superior Province in the Canadian shield and has been proposed as a cross-section of Archean crust exposed by thrust faulting along the Ivanhoe Lake Cataclastic Zone during early Proterozoic time. Five seismic refraction lines of 360-450 km in length were shot over the area to obtain a velocity structure for the uplift region. There were 18 profile shots and two fan shots with a recorder spacing of 2 to 5 km. Single trace processing has been developed for enhancement of secondary arrivals. Homomorphic filtering was used to obtain an estimate of the wavelet and then spectral division and cross-correlation with that wavelet were implemented to enhance the later arrivals on the trace. The improved secondary arrivals were used to locate lower crustal refracted phases and wide-angle reflections which were combined with the first arrivals to construct the velocity models. The wide-angle reflections were also correlated with events from coincident Vibroseis reflection data. The travel-times and amplitudes of the data from the profile lines have been modeled using asymptotic ray theory methods.  We have imaged a low velocity zone, ranging  from 4-5 km to 9-12 km depth, under the Abitibi greenstone belt; it is underlain by a highly reflective zone. There is a considerable deepening of the Moho from 40-43 km to 50-53 km under and to the west of the southern end of the Kapukasing Structural Zone. Analysis of wide-angle reflections on the fan shots has corroborated this thickening of the crust under the structure. A high velocity anomaly of 6.5-6.6 km/s has been imaged in the upper crust down to 20 km depth beneath the Kapuskasing structure with a suggested dip of 15° ± 2° to the west. High Poisson's ratio values (0.26), determined from shear  n  wave arrivals, were imaged in the structure, suggesting a high mafic content for that material. Gravity profile modeling has been performed along the five seismic lines. The observed Bouguer anomaly variations have been matched using the structure as delineated by the seismic modeling and velocity-density values that fall within the scatter of the Nafe-Drake set of data points. 2d frequency domain filtering techniques have been applied to the regional gravity and magnetic data. Directional, bandpass, continuation and derivative filtering have been used to locate faults and terrane boundaries. The main thrust fault, the Ivanhoe Lake Cataclastic Zone (ILCZ), has been more clearly delineated by this filtering and agrees well with the surface geological mapping. A n extension of the Lepage Fault has been outlined and indicates that this normal fault may have been important in the post-thrust stabilising period. Seismic estimates of the soling depth for the ILCZ were further constrained by considering the rheological properties of the constituent lithologies. From the results of previous heat flow and heat production work in this region and the structures inferred from the seismic interpretation, geotherms and strength curves were constructed for the present, a period after the stabilisation of the Superior craton, and the time of uplift. Present day 2  heat flows were successfully matched using a mantle heat flow of 20 mW m ~ . The formation pressures and temperatures of the Kapuskasing granulites were achieved using only a 25 % increase in the mantle heat flow. From analysis of the strength curves for the time of uplift, estimates of the depth to the brittle-ductile transition for dry quartz ranged from 17 to 21 km. This depth is corroborated by the results of the seismic refraction and reflection analysis. The results of this thesis support the Percival and Card thrust model but with a soling depth of 17-21 km. This corresponds to a granulite zone thrust up from the mid to lower crust and is indicative of large scale horizontal tectonic processes in late Archean to early Proterozoic times.  iii  Table of Contents  Abstract  ii  List of Figures  xiii  Acknowledgement  xiv  1  2  Introduction  1  1.1  Formation of the Precambrian Crust  4  1.2  Geological and Tectonic Setting of the K S Z  10  1.3  Previous Geophysical Studies  14  1.4  Present Geophysical Studies  15  1.5  The 1984 Seismic Refraction Experiment  17  1.6  Gravity and Magnetic Data Acquisition  36  1.7  Outline of Thesis  38  Seismic Data Processing  41  2.1  Spectral Characteristics and Bandpass Filtering  43  2.2  Enhancement of Secondary Arrivals  43  2.2.1  Convolutional Model  48  2.2.2  Homomorphic System  50  2.2.3  Wavelet Extraction  52  2.2.4  Processing Scheme  55  2.2.5  Synthetic Data Examples  60.  2.2.6  Real Data Examples  67  2.2.7  Conclusions  72  iv  3  4  Seismic Data Analysis  76  3.1  The Interpretation Method  76  3.2  Travel Time and Amplitude Modeling of P-Wave Arrivals  81  3.2.1  Line 1  84  3.2.2  Line 2  90  3.2.3  Line 3  97  3.2.4  Line 4  106  3.2.5  Line 5  114  3.3  Fan Shot Analysis  126  3.4  Travel Time Modeling of Shear Wave Arrivals  129  3.4.1  Interpretation Approach  130  3.4.2  Results  133  Geopotential Field Data Analysis  145  4.1  Gravity Profile Modeling  145  4.2  2d Potential Field Analysis  154  4.2.1  Potential Field Theory  154  4.2.2  Construction of 2d Frequency Domain Filters  157  4.2.3  Results of Analysis .  161  4.3 5  Correlation of Potential Fields  188  Discussion and Interpretation  195  5.1  195  5.2  Seismic Results 5.1.1  P Wave Structure  •  5.1.2  S Wave Structure  199  5.1.3  Inverse vs Forward Modeling for Line 5  200  5.1.4  Comparison of Refraction and Laboratory Measured Velocities . .  209  Potential Field Results  195  212  v  5.3  5.4 6  5.2.1  Relationship of Density to Velocity Structures  212  5.2.2  Gravity and Magnetic Interpretation  216  Palaeogeotherm and Palaeostrength Curves  219  5.3.1  Present Day Thermal Structure and Strength  219  5.3.2  Post-Stabilisation Thermal Structure and Strength  222  5.3.3  Pre-Uplift Thermal Structure and Strength  225  Tectonic Development  229  Conclusions  236  References  243  A Potential Field Theory  257  vi  List of Figures  1.1  Geology map of Kapuskasing Structural Zone  2  1.2  Proposed geological cross-section for Kapuskasing uplift  3  1.3  Down-sagging of volcano-sedimentary basin  5  1.4  Simatic model for greenstone basin  6  1.5  Greenstone and gneiss vertical relationships  7  1.6  Back-arc greenstone belt formation  8  1.7  Micro-continental collisions resulting in alternating greenstone and gneiss terranes  9  1.8  Models for granulite metamorphism  9  1.9  Location map of previous seismic experiments in the region  12  1.10 1984 seismic refraction survey over Kapuskasing region  18  1.11 Shot B of Line 1  20  1.12 Shot H of Line 1  21  1.13 Shot C of Line 2  22  1.14 Shot J of Line 2  23  1.15 Shot D of Line 3  24  1.16 Shot G of Line 3.  25  1.17 Shot K of Line 3  26  1.18 Shot A of Line 4  27  1.19 Shot C of Line 4  -  28  1.20 Shot E of Line 4  29  1.21 Shot H of Line 5  30  1.22 Shot G of Line 5  31  vii  1.23 Shot E of Line 5  32  1.24 Shot J of Line 4  33  1.25 Shot F of Line 1  34  1.26 Location of receivers in refraction experiment  35  1.27 Bouguer gravity anomaly map  39  1.28 Aeromagneitc anomaly map  40  2.1  Two seismograms and periodograms for Line 5, Shot E  44  2.2  Two seismograms and periodograms for Line 4, Shot A  45  2.3  Two seismograms and periodograms for Line 5, Shot H  46  2.4  Raw and bandpassed sections for Line 4, Shot E  47  2.5  Flow diagram for single trace processing  57  2.6  Example illustrating contributions of each convolutional component to time and cepstral domain  62  2.7  Low quefrency region of convolutional components  63  2.8  Amplitude and phase spectra for convolutional components  64  2.9  Extracted wavelets, amplitude spectra and cepstra for convolutional components and processed versions of those components  66  2.10 Close up of low quefrency region of processed convolution components.  .  67  2.11 Examples of effects of increasing added noise on extracted wavelet and reflectivity  68  2.12 Examples of effects of added noise on wavelet extraction and cross-correlation. 69 2.13 A comparison of the extracted reflectivity, cross-correlation with extracted wavelet, and spike deconvolution for several synthetic seismograms.  . . .  70  2.14 Four real data examples showing the wavelet extraction and the crosscorrelation output  71  2.15 Data section showing effects of automatic gain control and power boosting. 73 2.16 Example of full processing scheme applied to real data section viii  74  3.1  Example of synthetic seismogram procedure  79  3.2  Locations of shot points and receiver lines in refraction experiment.  3.3  Examples of different upper mantle refraction arrivals.  3.4  Ray diagram and travel-time comparison for shot B of Line 1  86  3.5  Real and synthetic data sections for shot B of Line 1  87  3.6  Ray diagram and travel-time comparison for shot H of Line 1  88  3.7  Real and synthetic data sections for shot H of Line 1  89  3.8  Resultant velocity model for Line 1  90  3.9  Ray diagram and travel-time comparison for shot C of Line 2  91  . . .  82 .  83  3.10 Real and synthetic data sections for shot C of Line 2  92  3.11 Ray diagram and travel-time comparison for shot J of Line 2  93  3.12 Blow-up of low amplitude first arrivals for shot J of Line 2  94  3.13 Real and synthetic data sections for shot J of Line 2  95  3.14 Resultant velocity model for Line 2  96  3.15 Ray diagram and travel-time comparison for shot D of Line 3  98  3.16 Real and synthetic data sections for shot D of Line 3  99  3.17 Ray diagram and travel-time comparison for shot G of Line 3  100  3.18 Real and synthetic data sections for shot G of Line 3  101  3.19 Ray diagram and travel-time comparison for shot K of Line 3  102  3.20 Real and synthetic data sections for shot K of Line 3  103  3.21 Resultant velocity model for Line 3  104  3.22 Ray diagram and travel-time comparison for shot A of Line 4  107  3.23 Real and synthetic data sections for shot A of Line 4  108  3.24 Ray diagram and travel-time comparison for shot C of Line 4  109  3.25 Real and synthetic data sections for shot C of Line 4  110  3.26 Ray diagram and travel-time comparison for shot E of Line 4  Ill  3.27 Real and synthetic data sections for shot E of Line 4  112  IX  3.28 Resultant velocity model for Line 4  113  3.29 Ray diagram and tra.vel-time comparison for shot H of Line 5  115  3.30 Annotated ray diagram and travel-times of refracted families for shot H of Line 5  116  3.31 Real and synthetic data sections for shot H of Line 5  117  3.32 Blow-up of first arrivals from shot H of Line 5  118  3.33 Ray diagram and travel-time comparison for shot G of Line 5  120  3.34 Real and synthetic data sections for shot G of Line 5  121  3.35 Strong reflectors imaged by shot G on Line 5  122  3.36 Ray diagram and travel-time comparison for shot E of Line 5  123  3.37 Real and synthetic data sections for shot E of Line 5  124  3.38 Resultant velocity model for Line 5  125  3.39 N M O representations of the fan shot from J into Line 4  127  3.40 N M O representations of the fan shot from F into Line 1  128  3.41 Example of relative amplitudes of P- and S-wave arrivals  131  3.42 Variation in travel-times with change in Poisson's ratio  132  3.43 Shear wave sections with superimposed travel-times for Line 1  134  3.44 Shear wave sections with superimposed travel-times for Line 2  135  3.45 Shear wave sections with superimposed travel-times for shots D and G of Line 3  137  3.46 Shear wave section with superimposed travel-times for shot K of Line 3. .  138  3.47 Shear wave section with superimposed travel-times for shot A of Line 4. . 139 3.48 Shear wave sections with superimposed travel-times for shots C and E of Line 4  140  3.49 Shear wave models for Lines 1,2,3, and 4  141  3.50 Shear wave reflection from shot G of Line 5  143  3.51 Poisson's ratio models for Line 5  144  x  4.1  Calculated and observed gravity anomaly curves for Line 1  147  4.2  Calculated and observed gravity anomaly curves for Line 1  149  4.3  Calculated and observed gravity anomaly curves for Line 1  150  4.4  Calculated and observed gravity anomaly curves for Line 1. . .  151  4.5  Calculated and observed gravity anomaly curves for Line 1  152  4.6  Comparison of chosen velocity-density function to bounds of Nafe-Drake data points and velocity-density points for the final models  153  4.7  Geological map of region covered by gravity and magnetic data  162  4.8  Bouguer gravity anomaly map for the Kapuskasing region  163  4.9  Log power spectrum for the Bouguer gravity anomaly map  165  4.10 Radially averaged power spectrum for the Bouguer gravity anomaly map.  166  4.11 Low pass filtered gravity map. .  167  4.12 Directionally filtered gravity map  169  4.13 East-west trends of gravity map  170  4.14 First derivative map of gravity data  171  4.15 Downward continued gravity map  172  4.16 Second derivative map of gravity data  173  4.17 Upward continued gravity map  174  4.18 Upward continued gravity map  175  4.19 Aeromagnetic anomaly map for the Kapuskasing region  177  4.20 Log power spectrum for the aeromagnetic anomaly map  178  4.21 Radially averaged power spectrum for the aeromagnetic anomaly map.  .  179  4.22 Aeromagnetic map reduced to pole  180  4.23 Bandpass filtered aeromagnetic map  181  4.24 Low pass and directionally filtered aeromagnetic map  • • • • 182  4.25 East-west trends of aeromagnetic map  183  4.26 First derivative map of aeromagnetic data  184  xi  4.27 Downward continued aeromagnetic map  185  4.28 Aeromagnetic map upward continued to 10.0 km elevation  186  4.29 Aeromagnetic map upward continued to 50.0 km elevation  187  4.30 Gravity data filtered by amplitude spectrum of aeromagnetic data  190  4.31 Aeromagnetic data filtered by amplitude spectrum of gravity data  191  4.32 Cross-correlation of gravity and aeromagnetic fields  192  4.33 Cross-correlation for east-west trends of gravity and aeromagnetic fields.  193  5.1  Comparison of proposed geological cross-section and velocity model. . . .  198  5.2  Starting model for inversion.  202  5.3  Residual velocity anomalies for three iterations of the inversion  204  5.4  Travel-time comparisons and ray diagram for the final iteration of the inversion  205  5.5  Velocity models for four stages of inversion  207  5.6  Comparison of final results from forward and inverse modeling  208  5.7  Comparison of laboratory measured and refraction model velocities. . . . 210  5.8  Velocity, density variation, and a models for Lines 1,2, and 3  213  5.9  Velocity, density variation, and a models for Lines 4 and 5  215  5.10 Thermal profiles and strength envelope for present day crust  220  5.11 Thermal profiles and strength envelopes for period immediately after stabilisation of Superior craton  223  5.12 Thermal profiles and strength envelopes for period immediately prior to - 2  the Kapuskasing uplift for a reduced mantle heat flow of 25 mW m . . . 226 5.13 Thermal profiles and strength envelopes for period immediately prior to 2  the Kapuskasing uplift for a reduced mantle heat flow of 20 mW m ~ . . .  228  5.14 Development of crust in southern Superior province in 2765 to 2650 M a period  230  5.15 Collision of Churchill and Superior cratons in 1900 to 1800 M a period. xii  . 232  6.1  Contour maps of depth to material with a velocity of 6.6 km/s and thickness of the crust  6.2  237  Final geophysical and geological models for profile 4 across the northern part of the K S Z  6.3  239  Final geophysical and geological models for profile 5 across the southern part of the K S Z  241  xiii  Acknowledgement  I thank Bob Ellis for his latitudinal support and encouragement during the course of my thesis work. His advice and confidence in my work helped me to follow the long and winding path through the many fields of geophysics. Ron Clowes was always there to act as another advisor in matters of science and, more importantly, politics. I also extend my gratitude to Doug Oldenburg and Hugh Greenwood who, as members of my committee, ensured the verity of my work when it trod on their ground. John Percival of the Geological Survey of Canada deserves honourable mention for the many enlightening discussions on Kapuskasing and the wonderful processes of geology. I sincerely thank the members and students, both past and present, of the Department of Geophysics and Astronomy for the many years of lively debate, honest criticism, and humourous insight. The experience has been rich and rewarding, if not always good for my liver. Finally, I would like to thank my parents, Paddy and Sally, and my wife, Lynn, for the long years of support and love. They have helped lighten the load when the legs were creaking. I hope that they realise the fruits of their labour and that we all get to enjoy them together.  xiv  Chapter 1  Introduction  Since the 1960's, many of the features and movements of the Earth's crust that we have observed in our geological and geophysical analyses have been explained using the concepts of plate tectonics. To have a better understanding of the Earth's formational history, it is important to obtain a link between its accretion and these modern day plate motions.  Information on Archean continental growth mechanisms is limited as  the wealth of our knowledge pertains to events in the Phanerozoic. Two fundamental questions are: did the cratons behave as rigid plates; and did a vertical or horizontal style of tectonics dominate the crustal growth process? A n opportunity to study these problems is afforded to us through analyses of the limited large scale, Archean exposures which provide snapshots of crustal deformation episodes. A n example of such a large scale exposure is found in Northern Ontario and enables us to study several metamorphic grades of Archean crust and an associated tectonic feature. The Kapuskasing Structural Zone (KSZ) cuts obliquely across the predominantly east-west structural grain of the Archean Superior Province (Figure 1.1 ). The surface metamorphic assemblages and gravity profile modeling across the zone suggest that the structure is an exposed cross-section of Archean crust (Figure 1.2) [Percival and Card, 1983]. Similar cross-sections of the crust have been observed in other parts of the world: the Pikwitonei belt of northern Manitoba, the Ivrea zone of northern Italy, the Musgrave Block of central Australia, the Kasila group in western Africa, the Limpopo belt in southern Africa, and the southern edge of the Dharwar Craton in southern India [Fountain and Salisbury, 1981; Coward and Fairhead, 1980; Raase et al, 1986]. In July of 1984, as part of the multidisciplinary L I T H O P R O B E program, a regional scale seismic refraction 1  Chapter 1.  Introduction  84°  2  82°  80°  Figure 1.1. Schematic geology map for the Kapuskasing region in the central Superior Province of the Precambrian Shield after Percival and Fountain [1989]. The Kapuskasing Structural Zone (KSZ) and the four subprovinces of Wawa, Abitibi, Quetico, and Opatica are shown. Several important features are indicated: the Ivanhoe Lake Cataclastic Zone (ILCZ), the Lepage fault (LF), the Saganash Lake fault (SLF), the Wakusimi River fault (WRF), the Fraserdale-Mosoonee block (FMB), the Groundhog River block (GRB), the Chapleau block (CB) and the Val Rita block (VRB). A A ' is the location of the geological cross-section shown in Figure 1.2.  Chapter 1.  Introduction  3  Figure 1.2. Proposed geological cross-section for Kapuskasing uplift after Percival and Fountain [1989]. Location is shown in Figure 1.1.  Chapter 1.  Introduction  4  experiment was conducted over the K S Z to obtain a velocity structure for this exposed Archean cross-section. The Kapuskasing structure cuts across the east-west striping of greenstone and gneiss belts in the Superior craton. The tectonic evolution of these Archean greenstone and gneiss terranes is poorly understood. Most greenstone belts are severely tectonized [de Wit and Ashwal, 1986] which suggests an active environment during and after their formation. A knowledge of the deep structure and geometry of these greenstone/gneiss  belts is  fundamental to our understanding of Archean tectonic processes and the formation of continental cratons. Thus, the geophysical structural analysis presented here provides critical information at depth about these Precambrian structures.  1.1  Formation of the Precambrian Crust  Tectonic models for the formation of greenstone and gneiss belts fall into two broad accretionary categories: vertical and horizontal. Implicit in these models is belief in the existence of a dominating tectonic mechanism in Archean and Proterozoic times. A uniformitarian assumption can be made which carries modern day processes back into the Archean with the only modifications being the scale and speed [eg. Sleep and Windley, 1982]. Alternatively, we can assume a different, i.e. vertical, tectonic regime whereby the crust grew by magma injection and basaltic underplating [e.g. Kroner, 1985]. Common to all Precambrian terranes is the presence of large scale greenstone belts with associated gneiss belts. The relationship of these two belts to each other and to the lesser exposed, high grade granulite terranes is and has been the subject of great interest [Windley, 1975; Goodwin, 1978; Sinitsyn, 1979; Kroner, 1981; deWit and Ashwal, 1986].  Models for Greenstone and Gneiss Belts The following models are examples of vertical accretion in both simatic and sialic environments. Anhaeusser et al. [1969] and Anhaeusser [1971] proposed that the greenstone  Chapter 1.  Introduction  5  belts formed from gravitational down-sagging of a thin sialic crust due to volcanic piles (Figure 1.3 ). This downwarp was then further loaded by sediments. Finally, the volcanosedimentary basin was stoped around its margins by granitic invasion. This explanation  (b)  (cl  Figure 1.3. Diagrammatic models showing development of greenstone belt (from Anhaeusser, 1971). See text for explanation. does not account for the presence of the associated gneiss belts. Glikson [1972] suggested that the evolution of the cratons was intertwined with the generation of greenstone belts from simatic or oceanic crust. In this model, the oceanic crust is deformed into large scale folds (Figure 1.4). The resulting zones of partial melting result in plutonism which accentuates the topographical variations of the folds. The erosion of the topographic highs leads to sedimentary infill of the basins with volcanism associated with the downsagging. Further orogenic activity leads to deformation of both the basin and gneiss terranes. However, the geochemistry of greenstone basins does not agree with a purely simatic origin [Windley, 1986]. Further developments of the down-sagging model [Glikson, 1976; Goodwin, 1981] characterise the gneiss and granulite sequences as the coeval roots of greenstone belts that have been exposed by uplift and erosion (Figure 1.5). The initial formation of the greenstones occurred in a failed rift environment on sialic crust  Chapter 1.  Introduction  1  6  OCEANIC STAGE  Figure 1.4. Development of greenstone and gneiss belts from simatic origins (from Glikson, 1972). See text for explanation. and was accompanied by the development of a highly metamorphosed root. However, a more common approach to the greenstone/gneiss formation problem is to make certain uniformitarian assumptions and look for a horizontal accretionary environment to explain the observations. Tarney et al. [1976J suggested a back-arc location for the greenstones as found today in the Rocas Verdes of Chile. A regular subduction environment can result in extension in the back-arc region. This rifting is not always sufficient to open up an ocean, but it does thin the crustal layer and provide passages for volcanic extrusions to the surface. The resultant basin is downwarped by the volcanic piles, continentally derived sediments and lavas and sediments from the volcanic arc. Then, basin closure deforms the volcano-sedimentary belt with possible stoping by granites at a later stage (Figure 1.6). Windley [1986] further developed this idea to explain the formation of stable continental cratons. He considered the model of Tarney et al. to  Chapter 1.  Introduction  7  Figure 1.5. Coeval development of greenstone and gneiss belts in a sialic crustal environment (from Glikson, 1976). See text for explanation. describe the processes on each mini-plate or micro-continent. Then, if one considered the surface of the earth to consist of many mini-plates moving around, their accretion led to cratonisation (Figure 1.7). The deformation at their boundaries gave rise to high grade terranes of interthrusted tonalitic and basement gneisses with adjacent  greenstone/gneiss  belts. These horizontal tectonic interpretions for the cratonisation period (3.0 to 2.5 bya) have been common [eg. Fyfe, 1974; Bickle et al, 1980; Glikson, 1981; Sleep and Windley, 1982; Nisbet and Fowler, 1983]. However, it can be argued that these processes do not adequately account for the significant growth rates of the late Archean and that magma injection and crustal underplating must contribute considerable volumes to the continental crust [Kroner, 1985].  Formation of Granulite Terranes The above model of Windley [1986] predicts that high grade granulite zones develop at depth in the gneiss terranes during metamorphism in the collision of micro-continents. However, it is difficult to link all occurrences of granulites with such a specific setting.  Chapter 1. Introduction  8  Figure 1.6. The development of greenstone and gneiss belts in a back-arc environment (from Tarney et al, 1976). See text for explanation. Newton [1987] presented the main mechanisms for high grade metamorphism from consideration of their penological features.  The main features of note are: accompanying  crustal thickening, thermal perturbations at depth outlasting causative tectonic events, low water activity preventing melting, and mappable isograds analogous to younger erogenic metamorphic belts. There are four main models for the formation of granulites (Figure 1.8). The first model has a hot spot source for the metamorphic event. Mafic magmatic underplating of a stretched and thinned crust is the dominant process. Anatectic plutons, rising from the heated lower crust, pond at mid-crustal levels (Figure 1.8a). The lower crustal addition is predominantly gabbro and should result in a gravity high.  Chapter 1. Introduction St4im«nti and to vat in ulensional back-arc  9 Sialic'basemcnt'of foldtd gnaistet { Amtttaq-type ) with relic  TonaliUl as horiiontal «h«t intruiioni in moin arc  Shtlf-iype  Figure 1.7. A series of micro-continents which develop back-arc basins are accreted horizontally on to each other to form a craton (from Windley, 1986). See text for explanation. A second model [Dewey and Burke, 1973] advocates accordion style thickening of the foreland arising from a continent-continent collision. Magma from the subduction zone provides heat to the .thickened crust and creates the conditions for granulite metamorphism (Figure 1.8b).  This process does not explain the horizontal strain structures  prevalent in granulites. The final two models involve continental subduction at collision  Bosin Sediments, Evoporites Anotthositei 2  MohO  X]  X^HFr?——Contact  firr^' C0 w<^>  Aiwtectit  Gobbto  HOT-SPOT  S7=>  V Underplot*  *C0  ;  A«(eol«  ~ C ] WHOLESALE  •S/S  CONTINENT SU80UCT10N  Collapsed  Or,  S  ^* \ "pi CONTINENT COLLISION •=-1 ACCOROION THICKENING  ^ N. N.  \ Subduction | Zone Gobbros  D J A-SUBDUCTION  Stdimtnts  ^ ^ ^ l  Figure 1.8. Four models for high pressure and temperature granulite metamorphism (from Newton, 1987). See text for explanation.  zones. The subducted continental mass undergoes melting and supplies magma and heat to the base of the overriding plate [Hodges et al, 1982]. The thickened crustal section produces more radioactive heating and thus may elevate lower crustal temperatures to  Chapter 1.  Introduction  10  the granulite field (Figure 1.8c). This model adequately explains the pressures and temperatures required for the granulite metamorphism and can account for the horizontal strain patterns. Also, after erosion, isostatic rebound can uplift the granulites. A variation of this model is the A-type subduction of Kroner [1981]. The overriding plate is split and the lower lithosphere is dragged down into the mantle (Figure 1.8d). The consequences are the same as for the above model with the additional advantage that this type of event can take place in an intra-cratonic setting.  1.2  Geological and Tectonic Setting of the K S Z  The nucleus of the North American continental craton is the Superior Province which is the world's largest contiguous area of exposed Precambrian rock to remain tectonically stable since the early Proterozoic.  The greenstone belts of its metavolcanic-plutonic  subprovinces contain most of the mineral deposits for which the Superior Province is known. The KSZ is a zone of high grade metamorphic assemblages, trending S S W - N N E for approximately 500 km, that obliquely transects the predominant east-west structural grain of the central Superior Province (Figure 1.1). There is a continuous transition in metamorphic grade from low grade greenschist fades of the Michipicoten belt to the west of the map region, through the amphibolite fades of the Wawa domal gneiss terrane to the high grade upper amphibolite and granulite fades gneisses of the Kapuskasing structure. The Ivanhoe Lake Cataclastic Zone marks the abrupt return to lower grade greenschist and amphibolite facies in the Abitibi subprovince to the east [Percival and Card, 1983]. Previous tectonic interpretations  have postulated that the K S Z is a suture [Gibb,  1978], a failed arm of a triple rift junction [Burke and Dewey, 1973], a broad sinistral transcurrent fault zone [Watson, 1980] and an intracratonic basement uplift [Percival and Card, 1983; Ernst and Halls, 1984; Percival and Fountain, 1989].  Wilson [1968],  Chapter 1. Introduction  11  Gibb and Walcott [1971], and Gibb [1978] have all discussed the Kapuskasing structure in connection with the Circum-Superior suture zone. This suggested a link between the KSZ and the early Proterozoic collision of the Superior and Churchill cratons. Wilson [1968] took the Ungava suture along the prominent Kapuskasing structure but Gibb and Walcott [1971] followed the northern periphery of the Superior Province (Figure 1.9). Later, Gibb [1978] related the Kapuskasing gneisses to rifting due to the inferred dextral slip on the Winisk and Kenyon fault structures. Burke and Dewey [1973] postulated that the KSZ had been the failed arm of a triple-rift junction on two occasions. Firstly, it was the failed southerly arm of the Ungava rift circa 1700 M a and, secondly, as a failed arm of the Keweenawan rifting episode at 1100 M a . Both of these examples were considered to be plume generated triple rift junctions resulting in two active rifts and one failed arm or aulacogen, the KSZ. Watson [1980] reexamined the existing geological and geophysical evidence and proposed that the main faulting was a sinistral transcurrent motion that occurred shortly after the formation of the Superior Province, 2700 Ma. She postulated that the fault penetrated to the mantle and had been a zone of weakness over an extended period between the formation of the Superior craton and the Keewenawan rifting episode at 1100 M a . Misalignments of the fault with the predominant strike direction resulted in strips of granulite being uplifted to the surface in a similar manner to the Great Glen fault of Scotland. However, by the mid-1980's, the geological and geophysical evidence presented below, strongly supported the model of an intracratonic basement uplift with an unknown vertical displacement [Percival and Card, 1983]. The ages of the metasedimentary and metavolcanic suites in the Abitibi and Wawa subprovinces have been determined from zircon dates on deformed metamorphosed volcanics and post-metamorphic plutons [Percival and Krogh, 1983]. The supracrustal sequences and early intrusions developed between 2750 and 2700 M a followed by deformation and regional metamorphism between 2700 and 2680 Ma. The mafic gneisses of the Kapuskasing structure give concordant dates of 2650 and 2627 M a from rounded  Chapter 1.  Introduction  12  Figure 1.9. Location map of previous seismic experiments in the region. The thick solid lines are Precambrian province boundaries. The thick dashed line are the inferred locations of these province boundaries under the Phanerozoic cover. The boundary between the Superior and Churchill provinces is also called the Ungava suture zone. The medium dashed line is the U.S./Canada border. The thin dashed lines are the provincial boundaries. The dotted lines are the locations of seismic profiles from previous experiments.  Chapter 1.  Introduction  13  zircons of metamorphic origin. The common geochronological features of the greenstone belts in the Abitibi and Wawa subprovinces, allied with their lithological similarities, suggest that the two are part of a formerly continuous belt that was disrupted by the Kapuskasing uplift. The uplift is inferred to precede 1907 Ma from the age of alkaline plutons that cut post-thrust normal faults [Percival and McGrath, 1986]. From analysis of isotopic systems in the KSZ region, the major thrust is constrained to the early Proterozoic [Percival et al., 1988]: in the 2250-2200 M a range from biotite/argon data or in the 1950-1900 M a age range from Rb/Sr biotite data. Geobarometric analysis techniques, using garnet-orthopyroxene-plagioclase—quartz assemblages [Newton and Perkins, 1982], were utilized by Percival [1983] and Percival and McGrath [1986] to determine the formation pressures of the various penological suites. The above west-east transition is marked by an increase in final mineral equilibrium pressures from 2-3 kb in the greenschists through 4-6 kb in the amphibolites to 7-9 kb in the granulites. These values suggest a formation depth difference of approximately 20 km between the greenschists and the granulites. The amphibolite/grahulite transition 3  corresponds to an increase in average density from 2700 to 2820 kg m ~ [Percival, 1986] and in seismic velocity, at 600 M P a confining pressure, from 6.55 to 6.90 km/s [Fountain and Salisbury, 1986]. This transition may correspond to the Conrad or other similar discontinuities as imaged seismicaliy at depth in the Superior Province [Berry and Fuclis, 1973; Hall and Hajnal, 1973; Green et al., 1979] and combined with the lithological observations and geobarometric results suggests that the KSZ is an uplifted oblique cross-section through the upper two-thirds of the Archean crust [Percival and Card, 1983; Percival and Fountain, 19891. Thus, defining a geophysical model for the structure will increase our understanding of the formation and deformation of the central Superior Province and its relationship to the surrounding Archean plate motions.  Chapter 1.  1.3  Introduction  14  Previous Geophysical Studies  The KSZ was first recognized as a positive gravity anomaly [Garland, 1950; Innes, 1960] and a positive magnetic anomaly [MacLaren and Charbonneau, 1968].  The variable  character of the anomalies suggest considerable structural and lithological changes along the strike of the zone. Percival and McGrath [1986] have divided the Kapuskasing uplift region into 3 main blocks, each with their, own geopotential field and geobarometric characteristics.  Figure 1.1 shows the relationship of these blocks to the KSZ. The  Fraserdale-Moosonee Block ( F M B ) is the granulite facies outcrop to the north of the study area. The main granulite exposure is split by the Wakusimi River fault into the Groundhog River Block (GRB) to the north and the Chapleau Block (CB) to the south. The Fraserdale-Moosonee and Chapleau Blocks have coextensive positive aeromagnetic and gravity anomalies whereas in between, over the Groundhog River Block, the anomalies diverge by up to 45 km. In this central region, the positive gravity anomaly runs over the Val Rita Block ( V R B ) to the west while the aeromagnetic high, strikes along the northern part of the zone over the G R B . These differences have been explained by on-strike variations of post-thrust normal faulting. Percival and McGrath [1986] have found similarities between the K S Z faulting geometries and those of the Wind River thrust system [Lynn and Thompson, 1983]. Previous seismic studies in the area (Figure 1.9) have concentrated on the Lake Superior Basin, 300 km to the south-west, the Ontario-Manitoba border area, 700 km to the west, and the Grenville Front, 300 km to the south-east. Halls [1982] produced a summary of the work in and around the Lake Superior Basin that showed the crust to be of the order of 50 km thick under the center of the basin, but shallowing towards the Kapuskasing region. Green et al. [1979] produced velocity/depth profiles for the southern region of the Ontario-Manitoba border from expanding spread reflection profiles. Their results showed mid-crustal reflective zones at 19 and 22 km and a crustal thickness of 38 km. Hall and Hajnal [1973], from data collected in the southern Ontario-Manitoba  Chapter 1. Introduction  15  border region and north of Lake Winnipeg, imaged reflective zones at 18.2±4.8 km and estimated a crustal thickness of 34.3±2.8 km. Berry and Fuchs [1973] obtained velocity/depth profiles for the eastern Superior Province from Id modeling of refraction data. Their data suggested the possibility of a low velocity zone at 5.0 km. The main results included imaging of the Conrad discontinuity at 13 to 16 km and the Moho at 40 to 45 km depth. Mereu et al. [1986] analyzed data from the Grenville Province and the Grenville Front and showed the province to have crustal thicknesses between 35 and 40 km with thickening under the Front. Parker [1984] suggested the existence of a low velocity zone at 12 to 17 km depth and obtained crustal thicknesses of 33 to 40 km for the southern Abitibi region. At the same time as the K S Z refraction experiment in 1984, a 10 km long seismic reflection line was shot. This pilot line was located over the high grade metamorphic zone of the C B and across the ILCZ. Cook [1985] has interpreted a set of dipping reflectors as being related to the depth extent of the Ivanhoe Lake Cataclastic Zone. These events dip approximately 40° to the west at depths of 6 to 12 km along the short section.  1.4  Present Geophysical Studies  As mentioned above, the seismic refraction experiment was part of a multidisciplinary investigation of the region. Geophysical and geochemical projects, as part of the Canadian geoscience program, L I T H O P R O B E , were conducted over the zone to obtain a comprehensive geoscientific insight into this anomalous structure. Data from these research programs are still being interpreted. A Vibroseis reflection survey was conducted over the Kapuskasing Structural Zone and the southern Abitibi belt between November 1987 and February 1988. The interpretation of these data is still at a preliminary stage. As part of this survey, several high resolution profiles were included to investigate interesting features in the upper crust. The high resolution data over the ILCZ [Geis et al., 1988] images reflections, associated with the  Chapter 1. Introduction  16  ILCZ, from the surface down to 4.0 s two-way travel-time. Also, there are reflections that are probably related to the Shawmere anorthosite complex within the K S Z and to a zone of decollement beneath the structure itself. The regional data set confirms these results and indicates that the ILCZ soles at 4.0 to 5.0 s two-way vertical travel-time. The apparent dips of the lower crustal sections suggest a structural east-west trend similar to that of the southern Superior Province [Wesf et al., 1988]. There is no evidence of significant reflectivity associated with the Moho. Both natural and controlled source electromagnetic experiments have been run over the KSZ. Kurtz et al. [1988] have examined the magnetotelluric ( M T ) and U T E M responses of the C B and the ILCZ. As is typical for crustal profiles, they observed an increase in conductivity at 25 to 35 km depth. In general, they imaged a uniform horizontal electrical structure for the region except for a weakly conductive zone, dipping to the west in the top 3.0 to 3.5 km, at the ILCZ. A controlled source experiment using the U T E M method [Bailey et al., 1989] has also investigated the C B and ILCZ area. They imaged two near surface subhorizontal conductive zones at 2 and 5 km depth associated with the ILCZ. There was also an increase in conductivity at 15 to 20 km depth. These conductive zones do not appear to correlate with the seismic reflections or the lithology.  5  Most of the crust appears to be extremely resistive ( in the 10 fim range ) and  the conductive zones may relate to post-thrust activities in a sub-horizontal stress field. Chouteau et al. [1988] and Mareschal et al. [1988] have examined the M T response of the G R B further to the north. Their preUminary results indicate that there is not a zone of conductivity associated with the ILCZ, but there is one along the Saganash Lake fault. They have observed an increase in conductivity in the 15 to 30 km depth range. Nkwate and Salisbury [1988] have conducted a detailed gravity survey over the G R B . They have confirmed that the G R B does not have an associated anomalous field but that the high grade surface rocks do have high densities. Thus, the granulites of the G R B form a thin layer at the surface. The Val Rita block (VRB) has an anomalous gravity  Chapter 1. Introduction  17  high but the surface densities are not unusually high. This implies a buried causative body for this anomaly. From preliminary modeling of their data, they infer a crustal thickness of 51 km under the G R B and the V R B . A palaeomagnetic analysis of samples from sites within the C B requires westward tilting of the block, in the post 2550 M a period, in order that the palaeopoles fall on existing polar wander curves [Constanzo-Alvarez and Dunlop, 1988]. Similar analysis of samples from the dyke swarms in the Wawa, KSZ, and Abitibi terranes also require crustal tilting within the KSZ to obtain consistent results for dykes of similar ages [Ernst and Halls, 1984; Halls and Palmer, 1988].  1.5  The 1984 Seismic Refraction  Experiment  As part of phase I of L I T H O P R O B E , five reversed refraction lines of lengths 360 to 450 km were shot, three parallel to the strike of the structure and two in a cross-strike direction (Figure 1.10). In addition, Lines 1 and 4 were fan shot from sites F and J , respectively, and Lines 3, 4 and 5 were center shot. Lines 2 and 5 were shot in two parts in order to reduce the recording spacing and thus obtain better coverage on these important profiles. Thus, with 60 instruments deployed for each shot, this provided 3 km recorder spacing on these two lines and 5 km spacing on Lines 1, 3 and 4. Twenty shots, varying in size from 800 to 2000 kg, were used as the sources. To reduce costs, the shots were detonated in six small lakes and four flooded mine pits. Thus, the shot sites were limited in availability and this factor, combined with the line access difficulties, influenced the selection of geometries for the profiles. Each shot was loaded into several 150 liter plastic barrels and detonated typically at 20 m below the water surface. These explosions were recorded by vertical component seismometers on six different designs of instrument.  Each instrument had similar flat responses in the 1 to 14 Hz range and  thus were only scaled for differences in amplitude response in that range. information about the experiment is provided in Northey and West [1985].  Additional  Chapter 1.  Introduction  18  Figure 1.10. Locations of shot points and receiver lines in the 1984 seismic refraction experiment over the KSZ. The dotted lines represent the mid-points of the shot-receiver offsets for the fan shots from F and J into lines 1 and 4, respectively. Inset shows the location of the study area in North America.  Chapter 1.  Introduction  19  Figures 1.11 to 1.25 show the 13 profile and 2 fan shot lines. They are presented in a trace-normalised, reduced time format. No filtering of the data has been performed. The distances on the profile shots are measured as offset from the northern shot points on lines 1, 2 and 3 and from the western shot points on lines 4 and 5. The fan shots are measured azimuthally in a clockwise direction from due south. The character of the refracted and reflected phases is highly variable over all the sections. This is most likely a result of variable shot site conditions (lakes and flooded quarries) and the petrological and structural differences between the terranes. The refracted first arrivals from the upper crust ( P ) range from strong (Figure 1.11) to mod9  erate (Figure 1.12) to very weak (Figure 1.21) amplitudes. This indicates a variation in velocity gradients over the region. However, in the case of shot H of Line 5 (Figure 1.21), the amplitudes were excessively low and could possibly be the result of diffracted energy arriving first at those receivers. The refracted first arrivals from the mid and lower crust (P ) exhibit the same amplitude variations from strong (Figure 1.17) to weak (Figc  ures 1.11 and 1.15). The upper mantle refracted phases also var}' in strength from strong (Figure 1.23) to weak (Figure 1.15) as well as varying in signal to noise ratio depending on the spectral characteristics of the shot. These mantle phases arrive at distances greater than 200 km and were clearest when the shots had a predominantly low frequency character (c.f. Figure 1.21 to Figure 1.23). The intracrustal reflectivity is highly variable in amplitude and reverbatory nature. Some sections show a few distinct events (Figures 1.11 and 1.15) while others show continuous back-scattering from all levels in the crust (Figures 1.13 and 1.14). The reflections from the crust-mantle boundary vary from a strong and clear arrival (Figure 1.18) to a more diffuse arrival (Figure 1.11). The variations in intracrustal and Moho reflectivity indicate that the terranes vary not only in petrological structure, as seen above in the refracted arrivals, but also in the fine-scale layering of their gneissic structures and shear zones.  Figure 1.11. Profile shot for Line 1 from shotpoint B . This section, and each of the following ones, shows a schematic representation of the experiment layout and a fan of the velocities associated with the slopes of the travel-time branchs.  to  o  o  DISTANCE  (KM) to to  Figure 1.13. Profile shot for Line 2 from shotpoint C.  Figure 1.14. Profile shot for Line 2 from shotpoint J .  to  Figure 1.15. Profile shot for Line 3 from shotpoint D.  DISTANCE  (KM)  Figure 1.16. Profile shot for Line 3 from shotpoint G.  to  Figure 1.17. Profile shot for Line 3 from shotpoint K .  o  DISTANCE  (KM)  Figure 1.18. Profile shot for Line 4 from shotpoint A.  to 00  Figure 1.19. Profile shot for Line 4 from shotpoint C.  DISTANCE  (KM)  Figure 1.20. Profile shot for Line 4 from shotpoint E.  DISTANCE  (KM) CO  o  Figure 1.21. Profile shot for Line 5 from shotpoint H.  o  DISTANCE  (KM)  Figure 1.23. Profile shot for Line 5 from shotpoint E.  Figure 1.24. Fan shot for Line 4 from shotpoint J.  to Figure 1.25. Fan shot for Line 1 from shotpoint F.  35  Chapter 1. Introduction  50.5  50.0 A 49.5  C  D o  49.0 JL  : 48.5  3  E  .... :~:o...  .  G 48.0 Shot Point  O 47.5  Recorder .o  •  .'•  • Kapuskasing Structural Zone  J  47.0 100  I  46.5 85.0  zd  K  I I I I I I I I I I I I I I I I I I I I I I I I  84.0  83.0  •  i  i  I  t  82.0  i  i  i  i  i  i  i  i  I  81.0  i  i  i  i  i  •  Kilometres i < l i • • i  •  80.0  LONGITUDE  Figure 1.26. Location of receivers in 1984 refraction experiment.  79.0  Chapter 1.  Introduction  36  The receivers for the fan shots were at distances greater than 230 km and so the phases of interest are the P (e.g. Figure 1.25) first arrivals and the later P P (e.g. Figure 1.24) n  m  arrivals. Both these phases show variability in travel-time and amplitude across the two fans. This is indicative of significant structure on the crust-mantle boundary. Due to the access problems, receivers and shot points were not always colinear (Figure 1.26) and thus the crooked line geometry has led to intrinsic uncertainties in the travel-times and amplitudes due to out-of-plane effects. Along any given profile, the azimuthal variation was approximately ± 1 0 ° . These uncertainties cannot be quantitatively accounted for in the 2d forward modeling of these data and thus can only be kept under consideration when assessing the possible parameter variations of the final models. The different types of shot location resulted in varying seismic spectral characteristics (e.g., compare Figures 1.21 and 1.23) and thus different transmission properties. The character of the data also reflects the lithological and structural complexity of the region. Localized geological features and different footings for the seismometers, ranging from swamp to hard rock, caused some travel-time and amplitude variations that were more dependent on the local effects (e.g., 170 to 190 km in Figure 1.11 and 90 to 110 km in Figure 1.15) than the general crustal structure and thus were not resolvable by the data. Therefore, only the broad trends in the amplitudes were considered in the modeling of the data.  1.6  Gravity and Magnetic Data Acquisition  Regional gravity and magnetic data sets were obtained from the Geological Survey of Canada (GSC) in order to integrate information from these data with results from the analysis of the seismic data. The gravity data consisted of 19,000 measurements from an 8 by 6 degree region around the Kapuskasing Structure with irregular station spacing varying from 100-200 m to 5-7 km. A data grid with a 5 km digitisation interval was constructed from these  Chapter 1.  Introduction  37  points using a 2d Lagrangian and splined interpolation scheme. The Lagrangian component ensured that no false peaks were introduced and the splined component created a smoother surface.  The magnetic data were collected along flight lines 800 m apart  at an altitude of 300 m and were then digitised by the GSC onto a regular grid with a spacing of 812.5 m. These data were regridded onto the same base grid as the gravity data to facilitate direct comparisons. The region covered by the grid was 80°- 85° W longitude and 4 6 ° - 51° N latitude. Figures 1.27 and 1.28 show the Bouguer gravity and aeromagnetic anomaly maps, respectively, for this regridded region. A n overlay is provided in a pocket at the back of the thesis to facilitate relating of potential anomalies to their causative geological terranes. The main features are gravity highs associated with the Chapleau, Val Rita, and Fraserdale-Moosonee blocks and a magnetic high associated with the Groundhog River block. The east-west boundaries between the Wawa/Abitibi and Quetico/Opatica subprovinces are also prevalent on both anomaly maps. During the gridding procedure, the effects of 2d aliasing had to be considered. The problem is more complex than the Id case as it cannot be viewed as simple folding. The discrete data sets obtained are the result of the true field multiplied by a 2d comb filter with spacing Sx. This is equivalent to the convolution of the true spectrum with a 2d comb filter with a spacing of j^.  Thus, our spatial sampling was such that the  overlap of the repeated spectra was minimised [Clement, 1972]. Because the gravity data were irregularly sampled originally, the new gravity grid was optimised by a careful balance of the Lagrangian and splined interpolation schemes. However, to guard against aliasing effects in the magnetic data, the power spectra of the original data sets were studied to ensure that large amounts of high frequency energy were not abased into the desired passband. When the gridding had been completed, the results were contoured 1  and compared carefully to master G S C maps to ensure that no obvious aliasing had occurred. Inspection of the new power spectra revealed that the magnetic data, because Preliminary magnetic anomaly map (residual total field), National Earth Science series [1:1,000,000], Dept. of Energy, Mines, and Resources.  Chapter 1. Introduction  38  of their high frequency nature, may exhibit some small aliasing effects.  1.7  Outline of Thesis  In this thesis, I shall present the methods, analyses and results of seismic, geopotential and rheological investigations of the Kapuskasing Structural Zone. This chapter has outlined the geological background to this anomalous structure, the importance of understanding these features, their role in Precambrian tectonics, and previous geophysical work in the region. Chapter 2 details the theory and practice of the processing applied to the seismic data to enhance the phases of interest and discusses the shortcomings of the various techniques. Chapter 3 contains the seismic interpretation of both the P-wave and S-wave arrivals. The method of modeling and the interpretational approach are discussed, followed by a description of the analysis of each line. Chapter 4 presents the 2d and 3d analyses of the geopotential field data. Gravity profile modeling was performed along the seismic lines using the velocity models obtained in Chapter 3 as constraints. A brief outline of the formulation of the 3d analysis is then presented, followed by the design of the frequency domain filters. The application of these filters to the real data is then discussed with a presentation of the results. Chapter 5 discusses the seismic and potential field results and relates them to each other and to the constraints imposed on the faulting by geothermal and geostrength considerations. Finally, the tectonic implication of these geophysical models is discussed for the K S Z and the cratonisation of the Superior province in Precambrian times.  Figure 1.27. are in mgals.  Bouguer corrected gravity anomaly map for the K S Z region. The units  Figure 1.28. Aeromagnetic anomaly map for the KSZ region. The units are in nanotesla.  Chapter 2  Seismic Data Processing  Regional scale seismic refraction experiments are conducted to obtain a broadscale velocity structure for the crust and upper mantle. These data sets typically have consisted of one to five profiles of 200 to 400 km length, one to three shots per line, and receiver station spacing of 2 to 10 km. This poor spatial sampling renders the data unsuitable for automated inversions or multi-trace attribute analyses (e.g. f-k, tau-p transformations) which use information from the whole data set [Bessanova et al., 1974; McMechan and Fuis, 1987]. So, in general, only the travel-times and amplitudes for the first arrivals (refractions), the Moho reflection and possibly one intracrustal reflection are used to obtain a velocity structure for the earth [e.g. Meltzer et al., 1987; Zeyen et al., 1985]. More recently, experimental shot and receiver spacings have decreased with the result that more sophisticated inversion techniques can be employed [Kanasewich and Chin, 1985; Chapman, 1987]. In the last fifteen years, crustal refraction analysis has been performed using traveltimes and amplitudes in trial and error 2d forward modeling techniques using synthetic seismograms [e.g. Cerveny et al., 1977]. There has been modeling of the general character of secondary arrivals by one dimensional full reflectivity techniques [e.g. Sandmeier and Wenzel, 1986] but little effort has been devoted to examining these later phases individually. In 2d modeling, the suites of wide-angle intracrustal reflections (secondary arrivals) in the data are usually not analyzed even though they can provide additional knowledge about the structure under examination. This lack of analysis is primarily due to poor spatial sampling, trace amplitudes dominated by first arrivals, and reverbatory waveforms; but, more importantly, it is due to lack of processing applied to the traces. 41  Chapter 2. Seismic Data Processing  42  Even when the travel-times for the secondary events are analyzed [ e.g. Klemperer and Luetgert, 1987; Grad and Luosto, 1987; Catchings and Mooney, 1988], very little processing is applied to enhance these arrivals. In recent years the number of coincident seismic reflection/refraction surveys and reflection profiles over regions previously investigated by refraction studies has increased considerably [Mooney and Brocher, 1987]. Thus, the importance of relating the secondary reflected arrivals on refraction surveys to the prominent events on the reflection profiles has become evident. In this chapter, we show examples of filtering to enhance secondary arrivals and help in travel-time picking. Firstly, the variable spectral characteristics of the data and the usefulness of bandpass filtering will be shown. Then, we will present methods for the enhancement of secondary arrivals by adaptation and application of simple existing, reflection data processing techniques as applicable to single traces. In particular, we have adapted homomorphic filtering [Oppenheim, 1965; Ulrych, 1971] to wavelet extraction and subsequent deconvolution with that estimated wavelet. These processes make no assumptions about the phase character of the wavelet or the spectral character of the earth structure and succeed in enhancing crustal refraction data. Single trace processing was chosen because the spatial sampling in this type of data was too poor for multi-trace attribute analysis. Averaging phase or cepstral properties over several traces to improve the wavelet extraction has been employed previously [Clayton and Wiggins, 1976; Otis and Smith, 1977; Lane, 1982]. However, the recorded seismic waveform depends greatly on the receiver's local environment [Cormier and Spudich, 1984] and in typical crustal scale refraction experiments the footing for the seismometers can vary from solid outcrop to marsh, sometimes from one station to the next. Thus, a stacking approach for these data was deemed unsuitable.  Chapter 2. Seismic Data Processing  2.1  43  Spectral Characteristics and Bandpass Filtering  In order to improve the general character of the seismic sections and to help in phase identification, bandpass filtering was applied to each section. Spectral tests were performed on each section to estimate the seismic signal frequency band. This was then used to choose the optimum filter for the whole section. Both non-zero and zero phase 8-pole Butterworth filters were utilised for this procedure [Kanasewich, 1981]. The spectral character of the sections varied considerably due to different shot sizes and locations and due to variations in lithology. The signal frequency band of the sections varied from 2-7 Hz (Figure 2.1b,d), to 3-10 Hz (Figure 2.2b,d), to 2-17 Hz (Figure 2.3b,d).  The medium and high frequency range examples also show the expected  downward shift of the passband from near offsets (40-60 km) to far offsets (240-260 km). However, for a low frequency signal (Figure 2.1), there is little change in the frequency content with offset. In general, the crystalline terranes of the region have a high Q (1000 to 1200) as evidenced by the small shift of 2-4 Hz in the peak of the amplitude spectrum. Examinations of individual traces, such as above, were used to choose suitable bandpass filters for portions of or whole sections. A n example of a filtered section is shown in Figure 2.4 for Shot E of Line 4. A non-zero phase filter with a 1.0 to 9.5 Hz passband was applied to improve the lateral phase coherency and improved the general appearance of the section. The sets of arrivals marked 1 and 2 have improved lateral coherency. This filtering was also used as an aid in the picking of travel times through noisy portions of the data. A n example of this is shown by the sets of arrivals marked 3 and 4.  2.2  Enhancement of Secondary Arrivals  The filtering method of the previous section did not always help and, in particular, the later arrivals (both refracted and reflected) were not enhanced. Thus, some effort was made to emphasize these later phases by using single trace filtering techniques that are  Chapter 2. Seismic Data Processing  44  Figure 2.1. Seismic traces (a and c) and their corresponding periodograms (b and d) for approximate offsets of 60 km and 240 km. Data are from Shot E of Line 5.  Chapter 2. Seismic Data Processing  45  Figure 2.2. Seismic traces (a and c) and their corresponding periodograms (b and d) for approximate offsets of 60 km and 240 km. Data are from Shot A of Line 4.  Chapter 2. Seismic Data Processing  46  Figure 2.3. Seismic traces (a and c) and their corresponding periodograms (b and d) for approximate offsets of 60 km and 240 km. Data are from Shot H of Line 5.  Chapter 2. Seismic Data Processing  47  o  -20  30  80  130  180 DISTANCE  -20  30  80  130  180 DISTANCE  230  280  330  380  230  280  330  380  (KM)  (KM)  Figure 2.4. (a) Raw data for Shot E of Line 4 shown in common maximum amplitude format, (b) Bandpass-filter (1.0 to 9.5 hz) applied to section from (a).  Chapter 2. Seismic Data Processing  48  more common in the processing of reflection data. A brief overview of the convolutional model is presented, followed by a description of the basis to the wavelet extraction technique, and finally the application of the method to synthetic and real data examples.  2.2.1  Convolutional Model  One class of seismic reflection analysis is based on the convolutional model. In this model, the recorded signal s(t) is the convolution of a source wavelet w(t) with the primary impulse responses of a one dimensional plane layer earth r(i). It disregards multiples and the effects of anelastic attenuation and dispersion on the propagating waveform. For our purposes, the impulse responses of the earth structure are both reflected and refracted events. The mathematical basis of the convolutional model is  s(t) =  (2.1)  w(t)*r(t)  but more realistically the recorded data d(t) are represented by  d(t) = w(t) * r(t) + w(t) * ip(t) + n(t)  where w(t) * VKO  15  (2.2)  convolutional noise such as multiples and conversions and n(t) is  added noise from local receiver effects. To look at amplitude and phase spectra we need to Fourier transform, whereby we obtain  D(f) = W(f).R(f)  where D(f),  W{f),  R(f),  + W(f)Mf)  9(f), and N(f)  + N(f)  (2-3)  are the Fourier transforms of the data,  wavelet, reflectivity, noise impulses, and added noise respectively.  Clearly, we must  Chapter 2. Seismic Data Processing  49  make assumptions about the type and magnitude of the noise in order to obtain a good estimate of the desired wavelet or reflectivity. The process of deconvolution is the recovery of the earth structure from the seismic data, d(t), by making assumptions about the convolutional model and the spectral characteristics of its signal, s(t), and noise, if>(t) and n(t), components.  If the source  wavelet is known then the problem is simplified but, in general, all components are unknown. Existing techniques such as Wiener filtering, spiking, predictive, adaptive and zero phase deconvolution have varying degrees of success when certain assumptions are valid [Claerbout, 1976; Yilmaz, 1988]. Often, the autocorrelation of the wavelet is assumed equivalent to that of the data which requires that the reflectivity and noise are white. Also, from experience with reflection data and the modeling of explosive sources, the wavelet is assumed minimum phase. In the case of crustal refraction data, we are generally interested in post-critical reflections and so the reflectivity whiteness [Fokkema and Ziolkowski, 1987] may be a valid assumption. However, visual inspection of the data show that most of the source functions recorded are not minimum phase but are variable mixed phase wavelets. Thus, an alternative to the above deconvolutional techniques was required that allowed for a mixed phase wavelet. These standard deconvolutional processes attempt to improve the temporal resolution of the data by compressing the source signature to a spike. In our case, the resolution of secondary arrivals is poor due to emergent phases, low frequency content and the interference effects of fine scale layering on pulse shapes. So, we seek only to emphasize these previously obscured events. The secondary refracted phases can then be used as additional constraints in the travel-time modeling of the first arrivals and the wide-angle reflections can be compared to corresponding near vertical ones on coincident Vibroseis reflection data.  Chapter 2. Seismic Data Processing  2.2.2  50  Homomorphic System  It is analytically convenient to look at signals that have been added together. A linear filtering system may then be employed to separate the signal components. Such a system, L, is said to obey a principle of superposition  L[ax {t) + x {t)] = aL[xi(i)] + L[x (t)} x  where x (t) x  and x (t) 2  2  2  are the two signal components and a is a constant.  (2.4)  But in our  case the signals have been convolved, so we would like a system matched to the principle of superposition. Let us consider a class of non-linear system that obeys a generalized principle of superposition, H, such that  # M0** (i)] = # M 0 ] * # M 0 ] 2  (2.5a,b) [aOzi(0] = where A ,  aVH[ (t)] Xl  0 and's? are arbitrary operations. This is a homomorphic system [Oppenheim  and Schafer, 1975]. Now, we must define our vector spaces and transformations for the convolutional model. We need a characteristic system to transform our data to a domain where the components are linearly separable and an inverse characteristic S3'stem to bring us back to the time domain. For our purposes the z transform is the key to the characteristic system.  It allows the transformation from a convolutional system to a  multiplicative one, then the log function can be used to produce an additive system  x{t) = x (t) * x {t) x  X(z) =  2  (2.6a, b,c)  X {z).X {z) 1  2  log[X(z)] = loglX^z)} + log[X (z)) 2  Chapter 2. Seismic Data Processing  51  Taking the log calculation we have  X{z) = log[X(z)} (2.7a, 6) =  where X (z) R  and Xi(z)  X (z)+jX (z) R  I  are the real and imaginary components of the log of the z  transform and j is y/—l. If X(z) is also a z transform then it must be analytic in a region that includes the unit circle. Here, Equation 2.7b becomes,  X(e^)  =  X (en+jX (e^) R  I  3U1  The analyticity of X(z) also implies that X(e ) is a continuous function of a;. Since  X(en  = log \X(e")\ +  g[X(e^)]  jaX  JW  then Xfl(e ) is defined provided no zeros exist on the unit circle in the complex z plane J1  and the continuity of X/(e ") depends on the definition of the complex logarithm. In JW  general, the principal value AKG[X(e )]  will be used which is discontinuous and, in that JU  form, does not satisfy the continuity requirement of the z transform. AKG[X(e ')] is a multi-valued function with jumps of 2 T caused by wrap-around in the log function as it crosses UJ — ±7r and remains on the same Riemann surface.  JU  However, Xi(e ')  can  be calculated and is a smoothly varying function that results from the superposition of its continuous components. Thus, the principal value argument function, A R G , must be JW  unwrapped to eliminate the jumps and produce the smooth analytic function arg[A"(e )]. This phase unwrapping is performed before an inverse z transform is implemented to obtain the cepstral representation. The cepstrum is a time delay domain where ideally our convolutional components exist at different delay times or quefrencies; the wavelet and reflectivity have low and high quefrency characteristics, respectively. For the cepstrum  Chapter 2. Seismic Data Processing  52  to be real, X^e ") must be a periodic function of ui with a period of 2 T and 3  3W  and Xi(e )  3  X^e ")  must be even and odd functions of a; respectively [Oppenheim and Schafer,  1975].  2.2.3  Wavelet Extraction  A homomorphic filtering approach was deemed most suitable for the wavelet extraction as no specific assumptions had to be made about the components of the convolutional model. Combining the information from the previous two sections we get the following expressions [Jin and Eisner, 1984] for the characteristic system for the signal, s(t), oo s{z) = £ * ( 0 * - ' t— — OO  S(z) =  log[S(z)]=log\S(z)\+jarg[S{z)] (2.8a, 6, c,d)  S*{T) = — 2nj  T  l  is{z)z ~ dz JC T  + r (r)  — W^(T)  where the dagger superscript (f) indicates the cepstral domain, r is delay time or quefrency, and the circular integration is defined over z = p+ju;, —7r < LJ < IT. But, when we consider noise effects (Equation 2.2) in the input to Equation 2.8a, a more complicated expression for the cepstrum is obtained due to the convolutional and added noise. If, however, we assume the noise components to be significantly smaller than the signal, then schematically equation 2.8d becomes  T  d (r)  = wi(T)  + r\r)+^(T)  (2.9)  where £*(T) is the cepstral representation of a noise term divided by a signal term and may have components that overlap with the signal. We now have a representation for our convolutional model that is a linear superposition of its components. Ideally, we should  Chapter 2. Seismic Data Processing  53  be able to linearly filter out the desired wavelet cepstral response. Then, using an inverse characteristic system, we obtain the original wavelet, w(t), by the following  •W(z)=  W\T): (2.10a, b,c)  W(z) = exp[W{z)} t  1  (t) = J <fw(z)z - dz 2irj Jc This extracted wavelet can be used in deterministic deconvolution procedures on the original data d(t). For practical implementation, the z transform is evaluated on the unit circle and so a digital Fourier transform is used in the characteristic systems. It is important to consider the effects of the various components on the log output of the characteristic system when the Fourier transform is used. If there is a signal present, the real part or log amplitude will be stable. The imaginary part, or phase, exhibits instability in the presence of noise and is the focus of most of the attention in estimation of the wavelet. Thus, to assess the noise effects on the phase, we return to Equation 2.3 and rearranging it yields  D(f) =  W(f)R(f)  N  1  V) W(f)R(f)  +  R(f)  Taking the log of both sides we get  log[D{f)]=log[W(f)R(f)}+log = log[W(f)}+log[R(f)}  , *(/)  i +  R(f)  , +  N(f) W(f)R(f)  (2.11a, 6)  + log[T(/)]  where (2.12)  Chapter 2. Seismic Data Processing  54  is very important in the correct estimation of the phase of our signal. If \?(/) = N(f)  = 0, then log[T(f)] is zero and our amplitude and phase estimates are uncontam-  inated. This ideal case is not encountered in practice, so we would like to estimate the perturbing effects of both types of noise on the spectral estimates. If the additive noise is non-zero then N(f)  log\T(f)] = log 1 + so that for N(f)  <C W(f)R(f)  the phase should be dependable, i.e. close to the true  one. This is a reasonable assumption, as we would only use data with an obvious signal content. If there is only convolutional noise then  log[T(f)} = log 1 + n i l R(f) so that if  <C R(f) or the desired reflectivity power is considerably greater than that  of the multiple and converted impulse responses, then again the phase result would be dependable. But more realistically both noise sources are non-zero and so the criterion for phase dependability is Arm  *(/)  +  TSS « W(f)  (2.13)  R  (f)  and it follows that  log[D(f)) « log[W(f)R(f)}  + R(f)  mil * ( / ) + W(f)  The phase of D(f) is  R(f)  *(/)  </.,(/) + A sin(A^)  N(f) W(f)  .sin(^ (/)-dv(/)) u  (2.14)  Chapter 2. Seismic Data Processing  55  where <£„(/) is the phase of the left hand side of Equation 2.13, and <$> (f) and </>,(/) are T  the phases of the reflectivity and our desired signal phase. We can be selective and discard data traces on inspection where N(f)  is comparatively large. But 9(f)  is more difficult  to estimate as it includes contributions from multiples, conversions and scattering from a generally inhomogeneous medium. If both N(f)  and 9(f)  are white then the wavelet  extraction is less affected than the reflectivity. However, the noise is normally coloured which results in spurious effects at short times in the cepstrum and therefore an unstable wavelet estimate. Our main assumption was that the wavelet is stationary along a given trace.  In  shield areas, a high Q exists throughout the crust [Der and McElfresh, 1977; Singh and Herrmann, 1983]; thus the effects of dispersion and absorption along a given trace can be assumed nominal and the assumption of stationarity is valid. However, we have ignored that along an individual trace there may be various pre- and post-critical phases of the wavelet from the different reflecting horizons. By taking short time windows of data and using damped weighting, it was hoped to limit the effect of this phase variation on each trace. For the wavelet extraction, the wavelet and reflectivity were assumed linearly separable in the cepstral domain. This is generally true as the wavelet dominates the short delays or quefrencies of the cepstrum and the reflectivity predominates at long delays or quefrencies.  2.2.4  Processing Scheme  A processing scheme based on a convolutional model and homomorphic system as described above was designed for single trace processing of crustal scale seismic refraction data. No assumptions were made about the phase of the wavelet or about the colour of the reflectivity spectrum. The wavelet was assumed to be linearly separable in the cepstral domain [Ulrych, 1971; Ulrych et al., 1972; Tribolet, 1979]. Noise was estimated  Chapter 2. Seismic Data Processing  56  to contaminate all parts of the spectrum and cepstrum. The processing scheme is represented in the flow diagram of Figure 2.5. The following sections describe each step of that scheme.  Preprocessing A variable gain control was applied to boost some of the weaker secondary energy. This improved the output of the deterministic deconvolution with the extracted wavelet. The method used was to divide each point, d(t), along a trace by an energy window, eft), so that <t) = ±  mi-  o ( - e{t) + e d { t )  0  where o{t) was the output trace, w was the length of the window, m was a variable power and e was a water level parameter.  No bandpass filtering was applied as zeros  in the amplitude spectrum would have caused instability in the log function. Instead, a form of data adaptive filtering or power boosting was utilized. The time-domain trace was Fourier transformed and divided into its amplitude and phase components. The phase spectrum was left untouched. The amplitude spectrum was normalized and raised to a power greater than one. Thus, amplitudes close to unity were boosted with respect to the smaller components. The new amplitudes and the original phases were recombined and inverse Fourier transformed to obtain the new signal. The effects of this operation are similar to data adaptive bandpass filtering except that zeros are not introduced into the spectrum. Some data were improved greatly by the above procedures and occasionally, at this stage, the processing was deemed satisfactory for travel-time analysis.  Chapter 2. Seismic Data Processing  57  Input Data  Preprocessing  Picking & Weighting  Characteristic System  Linear Filter  Inverse Characteristic System  Shift & Deweight  Reflectivity  Extracted Wavelet  Deterministic Deconvolution  Output Data  Figure 2.5. Flow diagram for single trace processing. There are two outputs to the linear filter: one contains wavelet information and the other pertains to the reflectivity.  Chapter 2. Seismic Data Processing  58  Picking and Weighting This step started the wavelet extraction procedure. The beginning of a data trace was shifted to the onset of the first significant signal within that time window. The running window used in searching for this onset was assigned a length equal to the predominant wavelength of that data trace. Then, the first window with an energy greater than the calculated root mean square energy of the whole trace was chosen. This shifted trace was weighted with an exponential or ramp damping factor. The exponential damping had the effect of radially scaling the zeros and poles of the data in the z plane and creating a time series with a minimum delay character.  Also, if zeros existed close to the unit  circle (the region of convergence), then they were moved away from this region and the instability was eliminated [Ulrych, 1971]. Application of a damping ramp had the same effect but was more moderate. It did not alter the form of the input data as much as the exponential weighting which tended to bias the phase estimate towards that of the first event on the data trace.  Characteristic System The characteristic system was that as described in equations 2.8 of the wavelet extraction section.  The data were padded with zeros to twice their original length to improve  spectral and cepstral resolution.  The phase unwrapping was performed by a simple  routine for eliminating 27r jumps or the more sophisticated phase adaptive approach developed by Tribolet [1977]. If required, the unwrapped phase could be bandpass filtered as suggested by Buttkus [1975] to eliminate the effects of the reflectivity.  Linear Filter The cepstrum was linearly filtered [Ulrych, 1971; Buttkus, 1975; Tribolet, 1979] to separate the wavelet and reflectivity components. The low quefrency part of the cepstrum was multiplied by a boxcar or cosine bell filter. The length of the cepstral filter for the  Chapter 2. Seismic Data Processing  59  wavelet was initially chosen to be 0.5 to 1.0 times the apparent length of the wavelet in the time domain. This length was selected after visual inspection of the time-domain seismogram. Then, the cepstral filter length was varied in a trial and error manner to establish the optimum length for that seismogram or group of seismograms.  Inverse Characteristic System Equation 2.10 describes this inverse system which returned the data to the time domain. Since the exponential function is single valued, there were no problems with analyticity in this calculation.  Shift and Deweight In the phase unwrapping procedure of the characteristic system, the phase function was forced to be odd (as explained in the homomorphic system section). This involved the removal of a linear phase component which transforms to a time domain shift [Brigham, 1974]. The sign of the linear trend dictated the direction of the shift. The data components were also deweighted and repositioned with respect to a reference time after the shifting applied previously at the picking stage.  Extracted Wavelet and Reflectivity The wavelet extraction was the stable output. The wavelet was isolated by automated scanning of the extended output to discard the low level noise background. The reflectivity output was unstable due to noise sensitivity at long quefrencies.  Deterministic Deconvolution The output wavelet was then used in a deterministic deconvolution procedure with the preprocessed data output. A spectral division in the frequency domain [Ziolkowski, 1984]  Chapter 2. Seismic Data Processing  60  was employed to obtain a reflectivity output  7?m_ U  )  D  W  V) *W W(f)W*(f) + \  where W*(f) was the complex conjugate of W(f)  and A was a water level parameter  used to stabilise the inversion. If the extracted wavelet was subject to ringing this technique did not work well. In that case, the wavelet was simply cross-correlated with the preprocessed data to emphasize the signal arrivals. This operation was equivalent to the spectral division with a high water level.  Although this operation did give a  smeared zero phase output, it did succeed in highlighting the signal content better than the spectral division. There was also an option to take the envelope of the output and look merely for energy trends.  2.2.5  Synthetic Data Examples  The above scheme was applied to several synthetic data examples. Both a sparse and a rich reflectivity were used with added bandpass noise' varying from 0% to 25% . A Berlage wavelet [Berlage, 1930] was used to produce a mixed phase (A) and a minimum phase (B) waveform. The analytic form of this wavelet is  n  at  w(t) = AH(t)t e- cos(2ivf t 0  + (f> ) 0  where H(t) is the Heaviside step function given by  H(t) = 0  t <0  = 1 t >0 and A is an amplitude scaling factor, n is a time exponent, a is an exponential decay rate, / is the oscillatory frequency, and </> is the initial phase angle. 0  0  Chapter 2. Seismic Data Processing  61  Time, cepstral and spectral representations of the input wavelets, reflectivities, noise and their convolutions are shown in Figures 2.6, 2.7, and 2.8 to illustrate the contributions of each of the components to our problem. For a noise-free signal, the wavelet and reflectivity were linearly separable in the cepstral domain (Figures 2.6b(l) and (2)). However, both the convolutional and added noise manifested themselves at low quefrencies (Figures 2.7(4) and (6) ) and so contribute to the cepstrum (Figures 2.6b(7) and 2.7(7)) of the data. This was due to the fact that the rich reflectivities and coloured noise unwrapped to produce a low frequency phase curve (Figures 2.8b(4) and (6)). Even though the noise components had relatively little power (Figures 2.8a(4) and (6)), they had a marked influence on the extracted wavelet because of these disruptive phase effects. Thus, it was important to preprocess the data in an optimal manner to eliminate not only additive but also convolutional noise. Figures 2.9 and 2.10 show the effects of the various preprocessing steps on the extracted waveforms, the unwrapped phases and the cepstra. Without any preprocessing of the noisy data, the scheme was unstable and dominated by noise effects (Figure 2.9a(5)). In fact, the unwrapped phase differed significantly from the original wavelet phase (cf. Figure 2.9b(2) to b(5)). The exponential weighting (Figure 2.9a(6)) gave an improved result by shifting the zeros of the time series away from the region of convergence of the Z transform. However, the phase curve (Figure 2.9b(6)) was still not close to the true one but was more minimum phase with the result that the extracted wavelet was less ringy. A combination of this weighting with power boosting resulted in a good estimation of the wavelet (cf. Figure 2.9a(l) to a(7)) and a phase curve and cepstrum close to the true ones (cf. Figure 2.9b(2) to b(7); Figure 2.10(2) to (7)). The power boosting resulted in a sparse reflectivity and so the low quefrency problems associated with convolutional noise were eliminated. Buttkus [1975] suggested high pass filtering of the unwrapped phase (Figure 2.9b(9)) to improve the wavelet extraction (Figure 2.9a(9)).  However,  this operation is equivalent to taking the long times or high quefrencies of the cepstrum  Chapter 2. Seismic Data Processing  62 a  TIME (Sec)  0.0  1.0  2.0  3.0  4.0  5.0  6.0  7.0  8.0  ' 1I 1I 1I I I 1111I 1 I I I I 111111I I I 11I 1I 1I I I 11111 I  w Q  WW-  w^-  3 4  g S  5 6 7 I I 1 1I 111I 1 111I I I I I M I I 11 I 1l'-M 11M I 11 1I 1  QUEFRENCY (Sec)  1  1  •A'  i  « Q O H  a£  2 3 4 5 6 7  Figure 2.6. Example illustrating contributions- of each convolutional component to (a) time and (b) cepstral domains. Each trace is numbered to correspond to: (1) wavelet A of mixed phase; (2) sparse reflectivity; (3) convolved signal; (4) added convolutional reflectivity noise; (5) signal and convolutional noise; (6) added bandpass noise; (7) signal and both noise sources. The cepstra have been normalised.  Chapter 2. Seismic Data Processing  63  QUEFRENCY (Sec) -0.6  0.0  "~"v\_/^  0.6  2 3  W P  3  <  4  I  j«  I5  ^  6 7  Figure 2.7. Close-up of low quefrency region of Figure 2.6b with same annotation as that figure. The cepstra have been normalised.  64  Chapter 2. Seismic Data Processing  a F R E Q U E N C Y (Hz)  b F R E Q U E N C Y (Hz)  Figure 2.8. (a) Amplitude and (b) phase spectra for convolutional components. Numbers refer to inputs as described in Figure 2.6. The amplitude spectra are true relative amplitudes.  Chapter 2. Seismic Data Processing  65  which we know contain the contributions of the reflectivity and noise. Low pass filtering of the phase (Figure 2.9b(8)) is equivalent to short time or low quefrency filtering of the cepstrum but this, too, gave an unstable output (Figure 2.9a(8)). Phase curves are too sensitive to be directly filtered. Since an equivalent operation exists in the cepstral domain then it is advisable to use that more stable approach. Figures 2.11 and 2.12 show the effects of varying amounts of noise on the wavelet extraction and the deconvolved output for the convolution of wavelet B with a sparse reflectivity.  The cepstrally extracted reflectivity was correct with a noise free signal  (Figure 2.11b) but with noise of 10% or greater (Figures 2.11c and d) the output was unstable.  For noise levels of up to 25% in the signal passband, the cepstral wavelet  extraction was satisfactorily stable (Figures 2.11b, c and d). Spectral division using this wavelet and a small water level parameter had problems when the noise was above 15%. For this reason, cross-correlation or spectral division with a high water level was preferred for work with realistic amounts of noise (Figures 2.12a and b). Figure 2.13 shows a comparison of our deconvolution outputs to a spiking deconvolution (predictive deconvolution with a gap of one) for several examples with no noise and 10% noise. The cepstrally extracted reflectivity reproduced the original reflectivity from the noise free traces (Figures 2.13a(3), b(3) and c(3)). In the presence of noise, the cepstral reflectivity extraction was poor (Figure 2.13a(5)) or gave an output with variable time shifts and polarity (Figure 2.13c(5)). These problems were due to the high sensitivity of the reflectivity output to errors in the unwrapped phase estimates. When the source wavelet was mixed phase (Figure 2.13b), the spiking deconvolution had problems because of the violation of its phase assumption but the cepstral wavelet extraction and subsequent cross-correlation were successful. The cross-correlation performed well in all the cases but in the case of a rich reflectivity (Figure 2.13c), closely spaced events were not resolved. A rich reflectivity cannot be distinguished from convolutional noise, in our case, and so the trace was preprocessed to eliminate the undesirable effects of many  o P"" tsi  co  P  00 O  c+-  a era  °  P  P  P co  o P  3  5  H  3 p 3 Q to p o CO a < 3 o a P- P <T> o p P to era o 00 p£2. !T O  - J  ^ c r a 8P , h— P i—*  p  cn AMPLITUDE  1-1  a o -—- o < p P-> c-ter*P e £ 01 o O Pi-» P cd p 3 ° P P- o 3 p p cr. o P_ P L . <<<1—' co to _ n <p «- a^, a- tS) l-l O3 C Ov 13 o . o CO m CO /-—s re er l_i p P 0> CO n 3 P o a O n O O 01 CD C n n C n pb  1s  p  H O  de  y  n P> Co Co  a co O  s  :  O  ^  o  IS tsi  -  „ ^ S' • — ' o N  co  W  co J  CO  era' ;ur  (T> " 'oo era' P £-  0)  v  era - ' o to era' co co P to co p> P P o(D  co  T3  o> o CO P era P 01  P p  0)  p  •O co  t—» Cn  tt>  PO  P  n o  01  p "  to *2.  o CO O  o  i ?  ft  $. cf  P  5^ co  O  P  G ta  3  M  O  CO  2.  nd  p  bo AMPLITUDE  0> Ol ^ U IO  Chapter 2. Seismic Data Processing  -0.6 "1 1  w  67  Q U E F R E N C Y (Sec) 0.0 0.6 1 1 1 1 1 1  —  OH  =  Jly.  Q  !=> H  ^  —^NA^V^^—  —t—  "1  1 1"  Figure 2.10. Close up of low quefrency region of processed convolution components. The traces are numbered according to the caption in Figure 2.9. The output to the combined damped weighting and power boosting is shown by 7 and compares favourably to the true cepstrum in 2. The cepstra are normalised. small events. Thus, power boosting was used which produced a sparser time series and eliminated the phase unwrapping problems associated with convolutional noise.  2.2.6  Real Data Examples  Seismograms and profile sections from the refraction data were used as examples. Figure 2.14 shows the processing applied to four single traces with varying character. The extracted wavelet was cross-correlated with the data and compared to spiking and zero phase deconvolution applied to the same data. Our processing scheme performed well but the parameters for each example had to be carefully chosen. The output enhanced data is slightly ringy, but has successfully highlighted several late arrivals: Figure 2.14a(3) shows highlighted events at 1.2 s and 2.9 s; Figure 2.14b(3) enhances events at 2.3 s and 3.6 s; Figure 2.14c(3) shows prominent events at 2.4 s and 4.0 s. Figure 2.14d(3) shows an example of where the output to the procedure has produced extensive ringing. The other deconvolution techniques (4 and 5 in each panel) did not produce the desired  Figure 2.11. Examples of effects of increasing added noise on extracted wavelet and reflectivity, (a) From top to bottom: wavelet B , sparse reflectivity, and convolved output. For other three frames from top to bottom: input data, extracted wavelet, output from spectral division of input data with extracted wavelet, and homomorphically extracted reflectivity. Ramp weighting was applied to the inputs. These three frames represent: (b) no noise added; (c) 10 % noise added; (d) 10 % noise added and power boosting applied.  Chapter 2. Seismic Data Processing  69 a  w o  (1.  JL. oo  4-5  90 OO  TIME (Sec)  4-5  90  TIME (Sec)  Figure 2.12. Examples of effects of added noise on wavelet extraction and cross-correlation. For each frame from top to bottom: input data-, extracted wavelet, output from cross-correlation of input data with extracted wavelet, and homomorphically extracted reflectivity. Ramp weighting was applied to the inputs. These three frames represent: signal of figure 2.11a with (a) 10 % noise added, and (b) 25 % noise added. results due to the mixed phases of the waveforms. Before analysis of the wide angle reflections preceded, travel-time modeling of refracted first arrivals from that profile was performed. The resultant velocity model had discontinuities which were used to generate the theoretical travel-times of reflections. Processing was then applied to the profile to highlight the secondary events and the calculated times were used to locate possible reflections in the data. The velocity model was adjusted to agree with the observed reflections in the data while remaining consistent with the refraction modeling. Thus, the enhancement of the secondary arrivals was important in this final stage of the forward modeling. If an observed event did not correspond with calculated times from discontinuities in the velocity model, then a shear zone was asssumed to be the cause of the reflection. Figure 2.15 shows the effects of the preprocessing procedures on part of a real data section.  In this example, full processing was not required to produce an acceptable  Chapter 2. Seismic Data Processing  70 TIME (Sec)  0.0  1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 111111111111111111 i 1111111111111111111111111  w Q P> frill B.  a •< •v—'v  A/  W  /  r  V\A^ ~ V^ "— \f\j^^~^s\J\l\f\^s\r~~~  I 1 1 I I I 1 I I I I I 1 I I I I 0.0  I I  'v—\A/VVW\^-  I I 1 I I I 1 I I I I 1 I [ 1 I 1 I I 1 1 I I I 1  1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 I I I I 1 I 1 1 I 1 I I I I I I I I I I 1 I 1 I I I I 1 I I I i 1 1 I I I 1 I I I I I I  -^j\f^ W \ / V — ^ ~ -T  —v-  V W -  I 1 I ! ! 1 I I II I I I I 1 1 I 1 1 1 I I I 1 1 I I 1 1 1 I 1 I 1 I 1 I I 1 I I I I I 0.0  1.0  2.0  3.0  4.0  5.0  6.0  7.0  8.0  Figure 2.13. A comparison of the extracted reflectivity, cross-correlation with extracted wavelet, and spike deconvolution for (a) a sparse reflectivity and wavelet B , (b) a sparse reflectivity and wavelet A , and (c) a rich reflectivity and wavelet B . Each trace is numbered to correspond to: (1) the original reflectivity; (2) the noise-free seismogram; (3) homomorphically extracted reflectivity; (4) seismogram and 10 % added noise; (5) homomorphically extracted reflectivity; (6) cross-correlation of extracted wavelet and seismogram; (7) spike deconvolution of seismogram.  Chapter 2. Seismic Data Processing  71  TIME (Sec) 0.0  0.6  i  I  1.2 i i  I  1.8 i I  I  I  I  I  2.4 I I  I  3.0 I I  I  3.6 I  I  I I I  o.o  4.;  6 . 1.2. . 1.8. , 2.4 I I I I i I I I I I  |  3.^  |  3.6  . 4.2 . 4.8 . 5.4 I I I 1 I I I I I I I  « Q £»  3 a.  0.0  0.8 II  i  I  1.6  i l l  2.4  'I" I  3.2  I I I ! I I I I I I I I i II  I  4.0  • I  4.8  I -I  I  5.6  I I I I I I I I I I I  I  i i i i i  0.0  •  ! I I I I I I I I I I I I I I I I I I I I I I I I I I I  i  i i i i i  i  i i i i i i i  i ii  i i i i i i i  0.8 1.6 2.4 3.2 4.0 4.8 5.6 I I I 1 I 1 I I I I I I I I I I I I I I I I I I I I I I I  —  I I I I I I I I I I I I I I I I M  I I I I l l l l l l  Figure 2.14. Four real data examples (a,b,c, and d) showing a comparison of the cross-correlated output to other deconvolution outputs. Each trace is numbered to correspond to: (1) input data; (2) extracted wavelet; (3) cross-correlation of extracted wavelet and input; (4) spike deconvolution output; (5) zero phase deconvolution output.  Chapter 2. Seismic Data Processing  72  result. The suite of intra-crustal reflections (a to e) and the Moho reflection (f) have been successfully emphasized and were used as constraints in synthetic seismogram modeling of these data. As with all processing, comparisons were made between the processed and raw data to ensure the arrivals were real and not artifacts of the enhancement technique. In Figure 2.15a, branch a is weak between 35 and 100 km offset but has been successfully boosted in Figure 2.15b. Similarly, branch b between 55 and 100 km, branch c between 48 and 74 km, branch d between 55 and 95 km, and branch e between 90 and 110 km have all been emphasized to help in tracing these phases through increasing offset. The traveltime and amplitude variations along a given reflection are indicative of a heterogeneous structure that has been sparsely sampled. Figure 2.16 shows the full processing sequence applied to another data section. The time window looks at a noisy part of the data behind the first arrivals. Spectral division with a high water level was used as the deterministic deconvolution procedure. The peaks of the responses have delineated four possible wideangle reflections. Events a and b are reflections off the bottom of a low velocity zone, c is a lower crustal event and d is a reflection off the Moho. In all cases, the processed data were compared carefully to the raw data to ensure the verity of the results and to help resolve closely spaced events (a and b). When events are already apparent on the raw data (d), the peak of the spectral division response provides an objective estimate of their travel-times. These reflection travel-times have been used to supplement picks made on enlargements of the raw data sections in the forward travel-time modeling of the data.  2.2.7  Conclusions  Processing applied to crustal scale seismic refraction data is useful and can help constrain forward modeling of those data. Results of synthetic data analysis show that for bandpass noise of up to 25%, a stable wavelet extraction can be performed. A n optimum amount of information about secondary events can be obtained by cross-correlation or stabilized  Chapter 2. Seismic Data Processing  73  Figure 2.15. (a) Raw data section for Shot A of Line 4 plotted with true amplitudes scaled by distance; (b) preprocessing of automatic gain control and power boosting applied to (a). The solid curves represent,the travel-times as calculated through the velocity model for that profile. This section shows a suite of intra-crustal reflections (a to e) terminated by a reflection off the Moho (f).  Chapter 2. Seismic Data Processing  74  Figure 2.16. Example of full processing scheme applied to part of profile from Shot K of Line 3. (a) Raw data with true amplitudes scaled by distance; (b) full processing applied to data of (a) with spectral division. Peaks in resultant waveform mark the event travel-times and are joined by dotted curves.  Chapter 2. Seismic Data Processing  75  spectral division with this wavelet on the preprocessed data. Although the output to this procedure is ringy with the associated loss of temporal resolution, it does succeed in highlighting the secondary events which were previously obscured. True amplitude information is lost in this procedure.  However, it has been shown that interference  effects due to fine scale layering in crystalline rocks are the predominant influence on reflection signal amplitudes [Jones and Nur, 1984; Christensen and Szymanski, 1988]. Thus, the true amplitudes of reflections should not be used in ray theoretical modeling of these secondary events. The real data examples supported these results but showed that careful choices of the processing parameters were required for a good result. This factor was a problem in the simultaneous processing of large numbers of traces. Because of the variability of the data, the same parameters did not work well for all traces. Selection of the parameters for each individual trace was the solution. Damped weighting and power boosting were shown to be essential for a successful wavelet extraction. The weighting stabilized the data for the log transformation (Equation 2.8b and Figure 2.9a(6)) and the power boosting reduced the convolutional noise and improved the phase estimate (Equation 2.14 and Figure 2.9b(7)). The linear filter length and its shape were the most variable parameters. It would be desirable to design an automated scheme for estimating a filter length from the data and its phase. However, the phase unwrapping is a problem as there exists no way of correctly performing this operation on an unknown quantity [Tribolet, 1977; Poggiagliolmi et al., 1982]. Thus, the phase estimate may be incorrectly unwrapped and/or contaminated with noise. Unfortunately there is no way of determining when these problems exist in field data and, for this reason, many have discarded the homomorphic approach to deconvolution. However, with careful application and comparison of the processed output to bandpassed and raw data, more information can be obtained from crustal scale refraction data.  Then, comparisons of reflections in refraction and  vertical reflection data are made possible.  Chapter 3  Seismic Data Analysis  Two-dimensional forward modeling of the seismic refraction data has been performed. Travel-time and amplitude analysis using synthetic seismograms has been applied to the P-wave arrivals of the five profile lines to obtain velocity structures for those crustal cross-sections. The two fan shot lines were used to image structure on the crust-mantle boundary from wide-angle reflections. Because the shots were located in small lakes and flooded quarries, the water-rock interface provided a locus for P to S-wave conversion and hence most sections show evidence of SV-wave recordings on the vertical component seismometers. Travel-time analysis was performed on these shear wave arrivals to obtain variations in Poisson's ratio for the crustal profiles.  3.1  The Interpretation Method  The travel-times and amplitudes of the seismic data were modeled using an algorithm [Zelt and Ellis, 1988] based on zero order asymptotic ray theory [Cerveny et al., 1977]. The refracted arrivals were used to construct the velocity models. Then, wide angle reflections were used to locate zones of reflectivity within these velocity models. If the reflections could not be associated with discontinuities in these velocity models consistent with the refracted data information, then their amplitudes were assumed to indicate the presence of fine scale or lamellar banded zones [Fountain et al., 1987] in the structure. Thus, an attempt was made to obtain greater structural information by using the intracrustal reflection phases as well as the refracted phases.  76  Chapter 3. Seismic Data Analysis  77  The picking of the first arrivals was performed from analysis of the raw data, bandpassed data and instantaneous frequency representations of individual traces as appropriate. Sections were plotted at various sizes and with various perspectives to help in the travel-time determinations. Processing was applied to the data to assist in the picking of the secondary arrivals. The profiles were processed using single trace reflection processing techniques as described in the previous chapter. As with all processing of seismic data, comparisons were made between the processed and raw data to ensure the arrivals were real and not artifacts of the processing technique.  The cumulative error on the  travel-time picks was estimated at ± 50 ms for the first arrivals and ± 100 ms for the later arrivals. The modeling approach was chosen so as to maximize objectivity in the construction of the velocity model. Firstly, a tau-p inversion [Bessanova et al., 1974] was performed on the first arrivals of the end shots for the three lines along strike. From this procedure, we obtained extremal bounds on the tau-p functions for each shot. Then, for each profile, a particular velocity-depth function which satisfied the tau-p function to within the observed errors was obtained using linear programming techniques [Garmany et al., 1979; Au, 1981]. These results were used to construct the starting models for input to a Id forward modeling approach.  Reflecting boundaries with small velocity jumps of 0.05  km/s were included in the model at the refracting layer interfaces that did not already have discontinuities. Then, Id travel-time modeling was performed iteratively to reduce discrepancies between the observed data and the calculated arrivals. modeling was used to fully match the travel-times.  2d travel-time  The approach used was to alter  the model minimally in order to satisfy the observed travel-times.  After modeling of  these strike lines, the velocity-depth profiles at the intersections with the cross-strike lines were used as their starting models. When the travel-times for these two lines had been successfully matched within the accepted error margins, the velocity-depth profiles at the intersections were then used to return to the three strike lines and further constrain  Chapter 3. Seismic Data Analysis  78  those velocity models. Thus, the modeling proceeded until the travel-times for all five lines had been matched satisfactorily and the intersection velocity-depth profiles for those lines were in agreement. The amplitude modeling provided a refinement to these resultant models by constraining the velocity gradients and it, too, was performed iteratively until the general trends in the data had been matched. A n example of the important crustal phases and synthetic seismogram modeling procedure is shown in Figure 3.1. This figure shows a simplified ray diagram with a single ray representing each type of arrival. The main refracted phases propagate through the upper crust ( P ) , the middle and lower crust (P ), and the upper mantle (P„). The upper p  c  mantle phases have been modeled using both turning rays in the upper mantle and head waves along the crust-mantle interface. Reflected phases from throughout the crust (PiP) and from the Moho ( P P ) are modeled as additional constraints. The ray model is conm  structed of a number of subhorizontal layers, each of which can be broken up into blocks with vertical boundaries. These trapezoids have fixed velocities at the top and bottom 1  and thus allow for horizontal gradient variations. Reflections were generated from the layer boundaries and have been included to model the later arrivals. When the data required a change in dip of a layer, the boundary was divided so as to obtain a gradual change in slope and thereby reduce the problems of the ray method associated with sharp block edges. This procedure occasionally gives the impression of over-resolution in the model, but is actually associated with the construction of a smoothly varying model (eg. 200-250 km offset in Figure 3.1a). Figure 3.1b shows the travel-time curves generated by ray tracing with a full complement of rays for each family. The first arrivals are refractions through the crust and upper mantle. Behind these arrivals are the intra-crustal reflective phases terminated with a reflection off the Moho. These calculated travel-times were then compared to picks from the observed data, and suitable adjustments were made to the velocity model to improve the agreement between the real and synthetic data. The x  For details see Zelt and Ellis [1988].  Chapter 3. Seismic Data Analysis  0  50  79  100  150 200 0ISTANCE  250 (KM)  300  350  400  Figure 3.1. Example of S3mthetic seismogram modeling procedure with important phases annotated, (a) Simplified ray diagram with one ray per family; (b) calculated travel-time curves for model of (a) with full complement of rays in each family; (c) synthetic seismogram section corresponding to arrivals of (b). The labeled phases in (a) and (b) are: upper crustal refraction(P ); middle and lower crust refraction(P ); upper mantle refraction(P„); intra-crustal reflector(PiP); Moho reflector(P P). p  c  m  Chapter 3. Seismic Data Analysis  80  errors in the picking and the modeling resulted in uncertainties in the velocities of ±0.05 km/s and in the boundary depths of ± 0.5 km. Using a resolution limit of jjjA, where A is the dominant wavelength of the data, the resolving power varies on average from 0.6 km in upper crust to 1.0 km in the mid-crust and 1.5-1.6 km in the lower crust. However, considering the effects of a realistic heterogeneous crust on travel-times and amplitudes [Ojo and Mereu, 1986], the effective uncertainties are more likely to be ± 0.1 km/s and ± 1.0 km in the upper crust, ± 0.15 km/s and ± 1.5 km in the mid-crust and ± 0.2 km/s and ± 2.0 km in the lower crust and upper mantle. From trial and error tests on the velocity models, variations of 0.1 km/s and 0.2 km/s for a 20 km body in the upper and lower crust, respectively, were found to be detectable. The wavelet used in the construction of a synthetic seismic section (Figure 3.1c ) was a stack of several wavelets extracted from the corresponding real data section. Q values of 1000 to 1200 have been obtained by Der and McElfresh [1977] and Singh and Herrmann [1983] for shield rocks in the Northeastern United States. These values were deemed suitable for the attenuation due to anelasticity but probably did not account adequately for attenuation due to scattering in this faulted region and hence a lower bound Q estimate of 1000 was selected for the modeling procedure. The traces of the synthetic sections were plotted at distances consistent with those of the real data. Not all the travel-time and amplitude modeling was completely satisfactory.  The  discrepancies are due not only to the previously mentioned data problems but also to inherent limitations of the shot-receiver geometries and of the modeling procedure. In regions of modest lateral crustal variations, the reverse and center shot spreads are deemed adequate for resolving a velocity structure but, in geologically complex zones, the upper crustal variations away from the shot points are poorly constrained and their effects on the travel-times and amplitudes are unknown. The velocities of these zones are thus constrained by the more vertically travelling energy at different angles of approach from each shot and by their intersection with other lines. In these complex faulted zones, the  Chapter 3. Seismic Data Analysis  81  infinite frequency approximation of asymptotic ray theory is invalid due to wavefxont effects such as diffractions and multiple scattering [Cormier and Spudich, 1984; Weidmann, 1984]. The Kapuskasing structure is highly tectonized and faulted and so the above effects are most probably present in the data. Specifically, the data from shot H of Line 5 (Figure 1.21) and shot J of Line 2 (Figure 1.14) show these scattering characteristics which on occasion may manifest themselves as extremely low amplitude first arrivals. Some reflection amplitudes were not deemed to have been generated by simple step discontinuities. So, after many modeling iterations, if the high amplitudes in the data had not been matched then the boundary was left as a small velocity discontinuity (eg. 0.05 to 0.1 km/s) consistent with the refraction data and the reflector was assumed to be a complex layered zone. It has been shown that in faulted crystalline terranes the high amplitudes are most likely generated by the constructive interference of reflections from the fine scale layering in gneissic and mylonitic textured rocks [Jones and Nur, 1984; Christensen and Szymanski, 1988].  3.2  Travel Time and Amplitude Modeling of P-Wave Arrivals  In this section each line of the experiment (Figure 3.2) is discussed separately. There is a summary of the main features of the profiles followed by a description of the velocity structure obtained from the modeling of those data. The sections are plotted in reduced time format with a velocity of 7.0 km/s and with their true amplitudes scaled by a suitable distance factor as indicated in the figure captions. The synthetic sections are scaled by the same distance factor as their corresponding real data sections. A l l distances, here and in the following descriptions, are model distances for that line unless otherwise stated. The refraction amplitudes are characterized as weak, moderate or strong. Due to the relatively small amplitudes of the upper mantle refracted phases at far offset, they are sometimes difficult to see on a full section. Thus, as examples, four different upper mantle responses from the data are shown in Figure 3.3. They are plotted in common maximum  Chapter 3. Seismic Data Analysis  82  Figure 3.2. Locations of shot points and receiver lines in refraction experiment. The dotted lines represent the mid-points of the shot-receiver offsets for the fan shots from F and J into Lines 1 and 4 , respectively. Inset shows location of study area in North America. A small-scale schematic representation of the seismic array will be shown on all data modeling figures in Chapters 3 and 4.  Figure 3.3. Examples of different upper mantle refraction arrival qualities: (a) strong and clear response from shot C of Line 2; (b) moderate response from shot J of Line 2; (c) weak response from shot D of Line 3; (d) moderate response from shot K of Line 3. (a) and (b) were bandpass filtered between 1.5 and 9.5 Hz. (c) and (d) were filtered between 1.0 and 7.5 Hz.  Chapter 3. Seismic Data Analysis  84  amplitude format and show the variations in clarity of this phase. Figures 3.3a and 3.3b show strong upper mantle arrivals from Line 2 along the strike of the KSZ. Figures 3.3c and 3.3d show moderate arrivals from Line 3 over the Abitibi Greenstone belt. The reflection responses are classified by their amplitudes and reverberatory nature. The responses vary from diffuse, a relatively low amplitude event extended in time, to sharp, a high amplitude event with only one to two cycles.  In general, the sections  show considerable crustal scattering. As mentioned above, ray theory is inadequate for modeling of this type of energy. When the scattered energy showed lateral coherency over distances of 40 km or greater, it was modeled as a wide-angle reflection. These reflectors are represented in the models as dashed line segments. The resultant velocity models are smoothed for presentation purposes. The 2d smoothing operator is 5.0 and 1.0 km in horizontal and vertical extent, respectively, and smears any existing first order discontinuities in the models. Model regions of little or no ray coverage are shown as stippled zones and so those parts constrained by the data are apparent.  3.2.1  Line 1  Line 1 is the western strike line over the Wawa domal gneiss terrane (Figure 3.2). The northern shot (Figures 3.4 and 3.5 ) is characterized by strong upper crust arrivals fading to weak mid-crustal and lower crustal arrivals. There are some strong intra-crustal reflections and a moderate Moho response. The ray coverage is fairly uniform apart from a small shadow zone at 12 to 20 km depth and 50 to 150 km distance caused by lateral velocity variations in the upper part of the model (Figure 3.4b). The lower crust arrivals manifest themselves as first arrivals at 200 to 260 km offset and indicate the presence of high velocities of 7.4 to 7.6 km/s above the Moho (Figure 3.4a). The travel-times are well matched throughout this profile. The first arrivals of the real and synthetic sections compare favourably (Figure 3.5), but the synthetic reflected amplitudes are too weak aside from the prominent reflections from 20-25 km depth at 3.0 to 3.5 s.  The south-  Chapter 3. Seismic Data Analysis  85  Figure 3.4. Travel-time modeling for shot B of Line 1: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  86  Chapter 3. Seismic Data Analysis  DI STANCE (KM)  Figure 3.5. Synthetic seismogram modeling for shot B of Line 1: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 1 0  Chapter 3. Seismic Data Analysis  87  Figure 3.6. Travel-time modeling for shot H of Line 1: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  88  01 STANCE (KM)  Figure 3.7. Synthetic seismogram modeling for shot H of Line 1: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 1-2  Figure 3.8. Smoothed velocity model for Line 1. The solid curves are the isovelocity contours and the dashed lines represent zones of reflectivity modeled by wide-angle reflections. The stippled region has little or no ray coverage. The stars represent the shot locations, the arrows are the locations of intersections with the crossing lines, and the diamonds are the positions of the velocity versus depth profiles located at each end of the model.  Chapter 3. Seismic Data Analysis  90  em shot (Figure 3.6 and Figure 3.7) shows weak mid-crustal arrivals and a moderate response from the lower crust. Again, there is evidence of high velocities in the lower crust at 200 to 260 km offset (Figure 3.6). As in the northern shot, the travel-times and amplitudes of the first arrivals are matched satisfactorily. In both sections, there are moderate to strong intra-crustal wide-angle reflections at mid-crustal depths. The upper mantle refractions at far offset are weak but distinct and the Moho reflective response is moderate. The near-surface velocities increase from 5.9 km/s to 6.1 km/s towards the southern end (Figure 3.8). The upper crust, with gradients of 0.01 to 0.06 km/s/km, follows this velocity trend apart from the pull-down of the 6.2 and 6.4 km/s contours at 5 to 12 km depth between 90 and 160 km south of shot point B . A pull-up of the travel-times and reduction in amplitudes characterize this region which coincides with the outcrop of the LePage Fault at 130 km. The mid-crust shows higher velocities at the southern end of 6.6 to 6.8 km/s at 16 to 24 km, up from 26 to 28 km in the north. The lower crust has high velocities of 7.2 to 7.7 km/s with low gradients of 0.005 to 0.01 km/s/km. There appear to be strong reflectors at 15 to 25 km depth. The Moho shows little topography and is horizontal at 48 km depth with a moderate reflection response. The upper mantle velocities are uniform at 8.1 to 8.2 km/s.  3.2.2  Line 2  Line 2, the centre strike line, runs along the K S Z structure (Figure 3.2). The northern shot (Figure 3.9 and 3.10) shows strong upper crustal first arrivals with a pull-down of the travel-time curves and decrease in amplitudes as the Wakusimi River fault at 80-110 km (Figure 3.9a and 3.10) is crossed by the profile. This delimits the northern boundary of the Chapleau block. The lower crustal refracted phases do not appear as first arrivals and are probably lost in the later scattered energy, but there are clear upper mantle refracted arrivals at far offset (Figures 3.3a). The mid and lower crust yield high  Chapter 3. Seismic Data Analysis  91  Figure 3.9. Travel-time modeling for shot C of Line 2: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  92  DISTANCE o  (KM)  N  -8()  -20  40  100 160 DISTANCE (KM)  220  280  340  Figure 3.10. Synthetic seismogram modeling for shot C of Line 2: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 1 2 5  Chapter 3. Seismic Data Analysis  93  Figure 3.11. Travel-time modeling for shot J of Line 2: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  94  120 140 DISTANCE (KM)  220  Figure 3.12. Blow-up of low amplitude first arrivals for shot J of Line 2 as southern edge of Chapleau block is crossed. The secondary arrivals have been clipped for presentation purposes so as not to overlap the traces and obscure the arrivals of interest. The calculated travel-time curves are superimposed.  95  Chapter 3. Seismic Data Analysis  01 STANCE (KM)  Figure 3.13. Synthetic seismogram modeling for shot J of Line 2: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r ' . 1  25  Figure 3.14. Smoothed velocity model for Line 2. The solid curves are the isovelocity contours and the dashed lines represent zones of reflectivity modeled by wide-angle reflections. The stippled region has little or no ray coverage. The stars represent the shot locations, the arrows are the locations of intersections with the crossing lines, and the diamonds are the positions of the velocity versus depth profiles located at each end of the model.  Chapter 3. Seismic Data Analysis  97  amplitude secondary energy that is only partially coherent and terminates with a diffuse Moho reflection response. The ray coverage is uniform and the travel-times of the first arrivals are well matched (Figure 3.9). The amplitude comparison is favourable for the first arrivals (Figure 3.10) but the high amplitude reverbatory secondary arrivals could not be matched. The southern shot shows the same characteristics as the northern one. There is a travel-time pull-down and amplitude decrease at 170-180 km (Figures 3.11 and 3.13) as the south-western edge of the KSZ is crossed. A blow-up of these arrivals is shown in Figure 3.12 with all other arrivals clipped. Again, the first arrival modeling is satisfactory but the reverbatory nature of the secondary arrivals is not matched (Figure 3.13). Due to their ringy nature, the sections have been scaled so that, at least, the first arrivals are clear. There is a noticeable increase in near surface velocities from 5.9-6.0 km/s at either end to 6.3-6.4 km/s in the center over the Chapleau block (Figure 3.14). This is accompanied by an upwarping of crustal velocites of 6.6 to 7.0 km/s from 18-33 km to 10-29 km under the K S Z itself. The lower crust shows no obvious lateral velocity variations but, from inspection of the data character, probably has a complex and variable fine layered structure which resulted in the observed scattered energy. The Moho deepens significantly from 44-45 km to 52-53 km to the south under the structure. This is observed from delays in the upper mantle refractions at 220 to 260 km offset. The upper mantle velocities range laterally from 8.2 km/s at the ends to 8.35 km/s at around 100 km with the Moho having a diffuse reflection response.  3.2.3  Line 3  Line 3, the eastern strike line, lies over the Abitibi greenstone belt (Figure 3.2). The northern shot (Figure 3.15 and 3.16) shows upper crustal arrivals that are strong out to 50-60 km. There is a delay in the travel-times at 50 to 80 km offset (Figures 3.15a and 3.16a) with the subsequent arrivals characterized by very weak amplitudes. This  Chapter 3. Seismic Data Analysis  98  CD  O _J  !  i_J  1  i_J  1  ! I i i i I  i  I  CD  Figure 3.15. Travel-time modeling for shot D of Line 3: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  o  99  _ l  1  I  i  I  [  i  i  I  0  40  80  120  160  200  240  280  320  N  DISTANCE  (KM)  360  S  o  b  0  40  80  120  160 200 DISTANCE (KM)  240  280  320  360  Figure 3.16. Synthetic seismogram modeling for shot D of Line 3: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r ' . 1  25  Chapter 3. Seismic Data Analysis  100  Figure 3.17. Travel-time modeling for shot G of Line 3: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  0  40  80  101  120  160  DISTANCE  200  (KM)  240  280  320  360  N  S b  0  40  80  120  1 60 200 01 STANCE (KM)  240  280  320  360  Figure 3.18. Synthetic seismogram modeling for shot G of Line 3: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 0 5  Figure 3.19. Travel-time modeling for shot K of Line 3: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  40  80  103  120  160 200 OISTRNCE (KM)  240  280  320  360  Figure 3.20. Synthetic seismogram modeling for shot K of Line 3: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 1 5  Figure 3.21. Smoothed velocity model for Line 3. The solid curves are the isovelocity contours and the dashed lines represent zones of reflectivity modeled by wide-angle reflections. The stippled region has little or no ray coverage. The stars represent the shot locations, the arrows are the locations of intersections with the crossing lines, and the diamonds are the positions of the velocity versus depth profiles located at each end of the model.  Chapter 3. Seismic Data Analysis  105  suggests the presence of one or more low velocity zones. The ray coverage is sparse through the low velocity zone at 4 to 10 km depth and through another smaller inversion at approximately 20 km depth (Figure 3.15b). The refracted energy from the upper mantle, beyond 200 km offset, shows a moderate response (Figure 3.16). The sections (Figure 3.16) are dominated by the strong upper crustal arrivals and the strong reflections from the mid-crust (3.0 s onset time) and from the Moho. The low amplitudes of the lower crustal event were not well matched and, since the travel-time modeling would not allow further change in the gradient, they may indicate a low Q is required for the deep regions of this crustal section. The center shot again shows characteristics of a low velocity zone with a delay in travel-times and decrease in amplitudes at 50-60 km and 260-270 km (Figure 3.17 and 3.18) There are several distinct intra-crustal reflections from the mid to lower crustal region and the Moho has a moderate reflective response. The southern shot (Figures 3.19 and 3.20) shows no significant deviations in travel-time and amplitude with offset, and thus no evidence for a further continuation of the low velocity zone to the southern end of the line. Again, there is evidence for strong reflectivity at mid-crustal depths (Figure 3.20). The lower crustal refracted arrival shows up at 170 to 180 km offset as a first and, further out, as a secondary arrival. The Moho has a fairly diffuse reflection response. On the end shots, the upper mantle refractions have a moderate response (Figure 3.3c and d). The ray diagrams for all the shots show a non-uniform coverage due to the presence of the low velocity zones. The average crustal velocity is lower than that of the other two strike lines. The upper 5 km has a velocity of 5.9 to 6.2 km/s with a high gradient of 0.05 to 0.1 km/s/km (Figure 3.21). Most noticeably, there is a distinct low velocity zone for most of the profile at a depth of 4-5 km to 9-12 km. There are two strongly reflective zones at 24-26 km and 29-31 km. The lower crust shows no significant lateral variations and has a comparatively lower velocity range of 7.0 to 7.4 km/s. The Moho shallows from 46 to 40 km both to  Chapter 3. Seismic Data Analysis  106  the north and south on the line with an upper mantle velocity of 8.2 km/s.  3.2.4  Line 4  This profile traverses the northern end of the K S Z exposure (Figure 3.2). The western shot (Figures 3.22 and 3.23) shows crustal refracted arrivals that are quite strong out to 100 km offset. There is an amplitude fall off as the Lepage fault is crossed at 100 km (Figure 3.23b). Distinct lower crustal first arrivals are observed at 170 to 220 km offset (Figure 3.22a) and there appears to be weak but coherent intra-crustal reflectivity over the Quetico belt (0 to 100 km) and the Val Rita block (100 to 170 km). The Moho gives a very distinct and sharp reflection response with moderately strong upper mantle refracted arrivals. The ray coverage is fairly uniform with depth. The center shot (Figures 3.24 and 3.25) shows stronger upper crust arrivals to the east and weaker ones to the west (Figure 3.25b). The mid-crustal arrivals are strong to the west and very weak to the east. This directional variation indicates significant differences between the Val Rita block to the west and the Abitibi belt to the east. The ray coverage is better to the west where stronger gradients are indicated at mid-crustal depths (Figure 3.24b). The Moho shows a moderate reflection response for this central portion of the line. The eastern shot is of poorer quality (Figures 3.26 and 3.27) with weak mid and lower crustal refractions. There are strong arrivals from the 25 to 35 km depth range that come in at 4.0 s and 200 to 280 km offset. Again, there is a strong reflection from the Moho. The ray coverage is good particularly in the lower crust. One general feature is the dramatic decrease in frequency at around 180 km. This may be related to scattering as the wavefronts pass through the ILCZ. There is an increase in near surface velocities from 5.9 to 6.2 km/s as one goes from west to east (Figure 3.28). On crossing the Lepage fault at 100 to 130 km, there is a general increase in upper and mid-crustal velocities towards the K S Z with quite high gradients. There appear to be only small variations in the lower crust with a velocity  Chapter 3. Seismic Data Analysis  107  Figure 3.22. Travel-time modeling for shot A of Line 4: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  n -20  i 30  i 80  108  i 130  i 180  DISTANCE  i 230  (KM)  I 280  1 330  I 380  W o  b  -20  30  80  130  180 230 DISTANCE (KM)  280  330  380  Figure 3.23. Synthetic seismogram modeling for shot A of Line 4: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 1 0  Figure 3.24. Travel-time modeling for shot C of Line 4: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  110  Chapter 3. Seismic Data Analysis  -20  -20  30  30  80  80  130  130  180  230-  280  330  380  180 230 01 STANCE (KM)  280  330  380  Figure 3.25. Synthetic seismogram modeling for shot C of Line 4: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 0 7 5  Chapter 3. Seismic Data Analysis  111  Figure 3.26. Travel-time modeling for shot E of Line 4: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  112  Chapter 3. Seismic Data Analysis i nil  -20  30  80  130  180 230 DISTANCE (KM)  iiiiiiiM  280  330  380  Figure 3.27. Synthetic seismogram modeling for shot E of Line 4: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 1 5  Figure 3.28. Smoothed velocity model for Line 4. The solid curves are the isovelocity contours and the dashed lines represent zones of reflectivity modeled by wide-angle reflections. The stippled region has little or no ray coverage. The stars represent the shot locations, the arrows are the locations of intersections with the crossing lines, and the diamonds are the positions of the velocity versus depth profiles located at each end of the model.  Chapter 3. Seismic Data Analysis  114  range of 7.4 to 7.6 km/s. The crust thickens from 42 km in the west to 50 km under the northern end of the Val Rita block and then thins again to 42-45 km under the Abitibi belt. The upper mantle velocities range from 8.1 to 8.2 km/s.  3.2.5  Line 5  This line crosses the southern outcrop (the Chapleau block) of the K S Z (Figure 3.2). Shot H , from the west, displays a strong pull-down of travel-times, very low amplitude first arrivals, and a large amount of high frequency scattered energy (Figure 3.29, 3.30, and 3.31). There is a noticeable pull down of the travel-times at 110 to 170 km over the Chapleau block accompanied by a large decrease in amplitudes (Figure 3.32). Prominent reflectors at 15 and 25 km depth are imaged by the high frequency energy. This energy (predominantly 16 Hz as compared to 4 to 6 Hz in most of the other sections) dies away at 160-180 km offset after the Ivanhoe Lake cataclastic zone has been crossed. As a result, the Moho reflection is diffuse and difficult to trace and the upper mantle arrivals are weak and often buried in noise. The travel-time comparison and ray tracing diagram are shown for the refracted and reflected families in Figure 3.29, and for the refracted ray groups only in Figure 3.30. Because of the significant lateral variations in the upper parts of the velocity model, there is poor ray coverage in the 18 to 28 km depth range. Each calculated refraction family is numbered (Figure 3.30) to show what parts of the model are constraining the first arrival travel-time variations. The synthetic section was constructed with two different wavelets in an attempt to match the unusual spectral character of the real data (Figure 3.31).  One wavelet was chosen for the oscillatory,  high frequency upper crustal arrivals and another, lower in frequency, was selected for arrivals at greater offset from deeper in the crust.  The seismic characteristics of this  section are extremely variable and reflect the geological complexity of this crustal crosssection and hence zero order asymptotic ray theory may have been inadequate for the modeling of parts of this profile.  The center shot (Figure 3.33 and 3.34) shows low  Chapter 3. Seismic Data Analysis  115  Figure 3.29. Travel-time modeling for shot H of Line 5: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  116  Figure 3.30. Annotated ray model and travel-times for shot H of Line 5: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled refracted phases. Each refracted raj' family is numbered in both (a) and (b).  Chapter 3. Seismic Data Analysis  117  DISTANCE (KM)  Figure 3.31. Synthetic seismogram modeling for shot H of Line 5: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 1 6  Chapter 3. Seismic Data Analysis  118  F i g u r e 3.32. Blow-up of low amplitude first arrivals that have been pulled down over the Chapleau block. The later arrivals are clipped for presentation purposes. Calculated travel-time curves are superimposed. The pull-down of the arrivals between 85 and 125 km is noticeable, but then the first arrivals weaken considerably beyond that distance. However, there are some arrivals that can be discerned and they have been indicated by arrow heads.  Chapter 3. Seismic Data Analysis  119  amplitude first arrivals to the east with a travel time delay at 140 km offset (340 km model distance) for the mid-crustal arrivals. There is some intra-crustal reflectivity and a moderately sharp reflection response from the Moho. The western half of the center shot is characterized by a travel-time pull-down of low amplitude events over the K S Z at 20 to 70 km offset and two distinct intra-crustal reflections with remarkably high amplitudes (Figures 3.34 and 3.35). These reflection amplitudes could not be reproduced by a simple step discontinuity in the ray model, so we have assumed fine-scale layering associated with a mylonitic fault zone to be the source of these remarkably strong events [Christensen and Szyrnanski, 1988]. There is poor refracted ray coverage below 10 km to the west and below 20 km to the east. The shot from E (Figure 3.36 and 3.37) in the east is much lower in frequency content and secondary energy than shot H. The upper crustal arrivals exhibit strong amplitudes out to 80 km offset giving way to weaker arrivals from the mid and lower crust. The upper mantle arrivals beyond 200 km offset are strong and distinct and the Moho reflection is strong and continuous (Figure 3.37). There is little intra-crustal reflectivity apparent on this section. There is good refracted ray coverage in the mid and lower crust. The resultant ray tracing model was extremely complex in order to account for the highly variable nature of the travel-times and amplitudes. The near surface velocities grade from 5.9 km/s over the Wawa domal gneisses in the west to 6.5 km/s over the Chapleau block and then back to 6.0 km/s in the east over the Abitibi belt (Figure 3.38). Material with a velocity of 6.6 to 6.7 km/s, and a very low gradient of less than 0.005 km/s/km, rises from 20 km under the Wawa gneisses to 4 km beneath the KSZ. The disturbance of the velocity contours decreases with depth suggesting that this velocity anomaly does not reach into the present day lower crust (see discussion in conclusions). Again, a low velocity zone is required by the data in the Abitibi belt. The lower crustal velocities range from 7.2-7.7 km/s in the west to 7.2-7.4 km/s in the east. The Moho deepens considerably from 48 to 53 km under and to the west of the Chapleau block  Chapter 3. Seismic Data Analysis  120  Figure 3.33. Travel-time modeling for shot G of Line 5: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  121  DISTANCE  (KM)  W o  OISTRNCE (KM)  Figure 3.34. Synthetic seismogram modeling for shot G of Line 5: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 0 7 5  Chapter 3. Seismic Data Analysis  122  o  80  90  100  110  120 130 DISTANCE (KM)  140  150  160  170  Figure 3.35. A blow-up of the strong reflectors imaged from shot G of Line 5 under the KSZ. These are reflectors I and II of Figure 3.38. The solid curves superimposed are the travel-times calculated through the ray model.  Chapter 3. Seismic Data Analysis  123  Figure 3.36. Travel-time modeling for shot E of Line 5: (a) comparison of observed travel-times (symbols) to those calculated from ray model (curves) shown at times with a reduction velocity of 7.0 km/s; (b) simplified ray diagram showing modeled phases.  Chapter 3. Seismic Data Analysis  124  01 STANCE (KM)  Figure 3.37. Synthetic seismogram modeling for shot E of Line 5: (a) synthetic seismogram section with traces plotted at distances consistent with the real data; (b) real data section. Both sections are plotted with a reduction velocity of 7.0 km/s and with true amplitudes scaled by r . 0 7 5  Figure 3.38. Smoothed velocity model for Line 5. The solid curves are the isovelocity contours and the dashed lines represent zones of reflectivity modeled by wide-angle reflections. The stippled region has little or no ray coverage. The stars represent the shot locations, the arrows are the locations of intersections with the crossing lines, and the diamonds are the positions of the velocity versus depth profiles located at each end of the model.  £o  Chapter 3. Seismic Data Analysis  126  with upper mantle velocities varying from 8.1 to 8.2 km/s. Under the KSZ there is much high frequency intra-crustal scattering but there were two distinct reflective zones (I and II) imaged by the center shot.  The tops of these zones are placed at 15-19 km and  23-27 km, respectively, under the K S Z and show small dips of 6° to 10° to the west. Using velocities from this refraction analysis, the top event correlates very well with a prominent event on a brute-stack of a coincident high resolution section from the recently acquired Kapuskasing reflection data [Geis et al., 1988]. The Moho is a strong reflector to the east but gives a diffuse response from the deep root in the west.  3.3  Fan Shot Analysis  Fan shots provide a cross-sectional view of the crust. These shots were used to observe wide angle reflections from the crust-mantle boundary. Both fan shots show lower crustal refraction arrivals between the Moho reflection and refracted first arrivals from the upper mantle. The approximate locations of the imaged points at the source-receiver midpoints are shown in Figure 3.2. Thus, after a velocity structure for the region had been obtained from the preceding synthetic seismogram modeling, root-mean-square  velocity ( K ^ )  values were calculated at these midpoint locations. A static normal move-out (NMO) procedure was applied to the fan shots with this laterally varying V-rm,. The static N M O correction was a constant time shift calculated for the reflector of interest (in this case the Moho).  This reduced the fan shot data traces to a two-way vertical travel-time  representation at the offset midpoints (Figures 3.39 and 3.40).  The N M O time shift  proved insensitive to reasonable variations in the Vrm, as the large offset values were the dominant factor in the correction. In Figure 3.39, the arrow head at an azimuth of 216° marks the onset of the Moho reflection for the equivalent trace on shot J of line 2. This provides a starting point for locating the Moho reflection event desired on the fan shot. Variable gain control, data adaptive filtering, and bandpass filtering were applied to the data to help in locating the  Chapter 3. Seismic Data Analysis  127  Figure 3.39. (a) Static N M O corrected fan shot for J into Line 4; (b) same as (a) with A G C control and data adaptive filtering applied. In both frames, the stippled region represents a window around the Moho reflection event and the arrow marks the travel-time of the Moho reflection on the trace from Line 2 that is coincident with that midpoint.  Chapter 3. Seismic Data Analysis  128  Figure 3.40. (a) Static N M O corrected fan shot for F into Line 1; (b) same as (a) with A G C control and data adaptive filtering applied. In both frames, the stippled region represents a window around the Moho reflection event and the arrow marks the travel-time of the picked Moho reflection from the trace in Figure 3.39 at the approximate intersection of the two sets of midpoints.  Chapter 3. Seismic Data Analysis  129  onset of these arrivals. Also, a wavelet was extracted manually from each fan section and used in a deterministic deconvolution procedure. One wavelet (1.0 to 8.0 Hz amplitude spectrum) was deemed adequate as the spectral content of the arrivals was sufficiently consistent over the whole fan shot.  A spectral division [Ziolkowoski, 1984] with that  wavelet was performed in the frequency domain by a stochastic inversion procedure. This provided an objective picking of secondary arrival onset times that was used together with the raw, variably gained, and bandpassed data. Using the Vrm, values, approximate depth variations for the reflectors were obtained and compared to the results of the synthetic seismogram analysis. The example in Figure 3.39 for shot J into Line 4 shows a deepening of the Moho under and to the west of the KSZ. The decrease of the V„„ values to the 4  east serves to accentuate this travel-time trend when true depth is considered.  The  shot from F into Line 1 images a slight thinning of the crust to the north-east and south under the Abitibi belt. The V  m  , crustal values are fairly constant across this  fan shot and thus the depth trend mirrors the travel-time variations. From analysis of the two fan shots, it appears that the general trends and approximate converted depths corroborate the results of the profile modeling and provide an independent confirmation for deepening of the crust-mantle boundary beneath the KSZ. However, for a complex geological situation, two fan shots are inadequate for a satisfactory reflection mapping of the crust-mantle boundary and the results can only serve as a small component of the refraction interpretation.  3.4  Travel Time Modeling of Shear Wave Arrivals  As mentioned earlier, because the shotpoints were located in small lakes and flooded quarries, some of the explosion energy was converted from P- to S-wave motion. The available recordings were made on single component vertical seismometers and thus that component of the SV-wave motion was recorded by the instruments. Figure 3.41 shows an example of the relative amplitudes of the P- and S-waves for shot A of Line 4. On  Chapter 3. Seismic Data Analysis  130  several of the shots, the conversion of energy was significant and so picking of S-arrivals was rendered relatively easy for offsets up to 120-150 km. Lines 1 to 4 were used in modeling of Poisson's ratio, but the signal-to-noise ratio on Line 5 was too small to discern any refracted shear arrivals.  3.4.1  Interpretation Approach  The first assumption was that the conversion occurred at the shot location as the observed travel-times at the closest receiver stations were consistent with this assumption and the water-rock interface provided a strong locus of conversion. Then, the boundaries of the P-wave velocity models were deemed appropriate for the S-wave modeling and an initial Poisson's ratio (cr) of 0.25 was chosen as a suitable starting value for the crustal models. From laboratory rock measurements and field experiments, a a of 0.25 appears to be an average for the crystalline crust [Birch, 1966; Assumpgdo and Bamford, 1978]. The recorded shear wave motion (SV) could have been effected by the anisotropy of the mafic gneisses [Fountain et al, 1989], but this effect was assumed to be averaged out over whole ray paths. The water content of rock pores and cracks has a strong influence on cr [Nur and Simmons, 1969]; an increase in saturation results in a higher cr. The electromagnetic studies over the region [Mareschal et al., 1988; Kurtz et al., 1988; Bailey et al., 1989] image a generally resistive crust and indicate that the pore water connectivity, if any water filled pores exist, is poor. However, this does not eliminate the possibility of the existence of isolated pores and cracks with water content. For the purposes of this study, the presence of water in the shield rocks was considered minimal and so the felsic/mafic nature of the rocks was considered the dominating influence on cr variations. Felsic and mafic characters are associated with low and high <x, respectively, and have been used to infer penological variations in the lower crust [Holbrook et al., 1988]. If non-connected water cracks and pores exist throughout the shield crust, the observed a values will be shifted upwards but the relative differences should still reflect the felsic/mafic variations of  Chapter 3. Seismic Data Analysis  131  o CM  0  15  30  45  DISTANCE  (KM)  60  75  90  Figure 3.41. Example of relative amplitudes of P- and S-wave arrivals for shot A of Line 4. The low frequency surface wave arrivals, behind the S-wave arrivals, are also marked.  Chapter 3. Seismic Data Analysis  132  Figure 3.42. Variation in travel-times of shear-wave arrivals with change in Poisson's ratio. The ray model and the phases are those used in Figure 3.1.  133  Chapter 3. Seismic Data Analysis  the constituent rocks. Variations in cr can also be associated with the degree of hydration of the minerals [Hall and Ali, 1985], but there is no evidence for such variations from the exposed lithologies. The observed variation of Poisson's ratio within the final models was 0.23 to 0.26. A n example of the significance of this variation for an entire model is shown in Figure 3.42. The model structure for this example is the same one as in Figure 3.1. A smaller cr implies a lower P- to S-wave velocity ratio. For offsets of 100 to 150 km, the travel-time difference can be of the order of 1 s for this range of cr values. In general, the refracted shear arrivals are clear out to 100 km offset. Between 100 and 200 km, the arrivals are erratic and at further offset there is no evidence for refracted arrivals from the upper mantle. Previous refraction experiments have had little success in observing the shear wave equivalent of P  n  [Gajewski et al., 1987]. On our profiles where  P P arrivals are observed, their shear wave equivalents are also observed but generally m  with a more diffuse nature. Intra-crustal reflections are usually obscure and buried in continuous signal; however, where lateral coherence is evident a travel-time curve has been generated. The general amplitude characteristics of the P-wave section were used as guidelines when interpreting arrivals on an equivalent S-wave section. The calculated travel-time curves for the final a models are superimposed on the data sections.  3.4.2  Results  Line 1, over the Wawa domal gneiss terrane, showed strong conversions of the shot energy (Figure 3.43). The northern shot has clear first arrivals out to 160 to 180 km offset. There are two distinct intra-crustal reflection events from 25 and 35 km depth at 4.5 s to 6.0 s in Figure 3.43a. The Moho reflection is of moderate energy but diffuse in character. The southern shot has a lower signal-to-noise ratio with the clearest first arrivals observed at 100 to 200 km offset. The intra-crustal reflectivity shows little coherence aside from an event from 15 km depth at 3.0 to 3.5 s and the Moho reflection is difficult to follow. Line 2, on strike along the KSZ, shows similar degrees of secondary backscattered  134  Chapter 3. Seismic Data Analysis  -40  60  160 DISTANCE (KM)  260  360  DISTANCE (KM)  Figure 3.43. Shear wave sections with superimposed travel-times curves for: (a) shot B, and (b) shot H of Line 1. Both sections are plotted in common maximum amplitude format and with a reducing velocity of 4.0 km/s.  Chapter 3. Seismic Data Analysis  135  energy as the equivalent P-wave sections (Figure 3.44). On both shots, the first arrivals are good out to 100 km offset but rather obscure beyond that point aside from a few key arrivals. There is a noticeable intra-crustal reflection from the 17 to 22 km depth range on each section. High amplitude scattered energy from the lower crust and crust-mantle boundary region is a strong characteristic of this profile. Line 3, over the Abitibi greenstone belt, shows moderate amounts of converted shot energy (Figures 3.45 and 3.46).  The northern shot shows good first arrivals out to  110 km offset and displays some clear reflections from the mid to lower crust and the Moho. The center shot has moderate arrivals out to 80 km offset, shows a good midcrustal reflector at 2.5 s and 40 to 120 km model distance on the northern half, and has a o  moderate reflection response from the Moho. The southern shot shows clear arrivals out to 180 km offset with strong secondary arrivals between 190 and 250 km model distance and a discontinuous Moho reflection. Line 4, crossing the northern part of the KSZ, displays the clearest shear wave arrivals (Figures 3.47 and 3.48).  The western shot shows very clear arrivals out to 200 km.  However, even with such good conversion, it is impossible to discern an upper mantle phase.  There is one apparent intra-crustal event from the 20 km depth range at 4.5  s to 5.0 s and 80 to 150 km offset, and the Moho has a diffuse but strong reflection response. The center shot shows clear arrivals on its western arm with one good intracrustal reflection at 4.0 to 5.0 s and 30 to 110 km model distance, and also a diffuse Moho reflection response. The eastern portion shows less distinct arrivals but again a mid-crustal event is apparent (3.0 to 3.5 s and 230 to 290 km) and there is a moderate Moho reflection response. The eastern shot shows moderate first arrivals out to 120 km offset and a strong but ringy Moho reflection. Figure 3.49 shows the resultant Poisson's ratio (a) models from the above shear-wave modeling and the constraints of the P-wave velocity models. In general, the resolution in the upper 20 to 25 km is ±0.005 but, below this level, the resolution is poor as arrivals  Chapter 3. Seismic Data Analysis  136  Figure 3.44. Shear wave sections with superimposed travel-times curves for: (a) shot C, and (b) shot J of Line 2. Both sections are plotted in common maximum amplitude format and with a reducing velocity of 4.0 km/s.  137  Chapter 3. Seismic Data Analysis  0  120  160 200 OISTRNCE (KM) •  240  280  320  360  N  160 200 •ISTRNCE (KM)  360  Figure 3.45. Shear wave sections with superimposed travel-times curves for: (a) shot D, and (b) shot G of Line 3. Both sections are plotted in common maximum amplitude format and with a reducing velocity of 4.0 km/s.  Chapter 3. Seismic Data Analysis  138  Figure 3.46. Shear wave section with superimposed travel-times curves for shot K of Line 3. The section is plotted in common maximum amplitude format and with a reducing velocity of 4.0 km/s.  Chapter 3. Seismic Data Analysis  139  Figure 3.47. Shear wave section with superimposed travel-times curves for shot A of Line 4. The section is plotted in common maximum amplitude format and with a reducing velocity of 4.0 km/s.  140  Chapter 3. Seismic Data Analysis CM  -20  30  80  130  180 230 DISTANCE (KM)  280  330  380  -20  30  80  130  180 230 01STANCE (KM)  280  330  380  Figure 3.48. Shear wave sections with superimposed travel-times curves for: (a) shot C, and (b) shot E of Line 4. Both sections are plotted in common maximum amplitude format and with a reducing velocity of 4.0 km/s.  Chapter 3. Seismic Data Analysis  N  141  DISTANCE (KM)  -40  60  J  0.23  I)j4  160 L_  )  260  S  N  360  -80  l__  s  DISTANCE (KM)  -20  100  40  160  220  280  340  0.J5 I—  0-  LINE 1  b  LINE 2 8-  N 0  S  , y  DISTANCE (KM) 40 I  80 I  120  160  200  I  I  I  240 I 0.24  ° 280 I  320 I  360  W v  -20  DISTANCE (KM)  y  30  80  130  180  230  I  '  I  l_  280  330  380  I"' (—  Q.  ONLINE 3  Figure 3.49. Shear wave models for Lines 1, 2, 3, and 4 shown in (a), (b), (c), and (d), respectively. Contours represent a values of 0.235, 0.245, and 0.255.  Chapter 3. Seismic Data Analysis  142  from deep in the crust were obscured in noise. Thus, the values for the mid to lower crust are gross averages inferred by the occasional intra-crustal reflector and the Moho reflections. The northern end of Line 1 shows low c values for the top 5.0 km out to 80 to 85 km. This portion of the profile lies over the metasedimentary migmatites of the Quetico belt and the low tr reflects the high quartz content of these rocks. The southern end of the line has low cr values for the top 15 km and overlies metavolcanics.  The  termination of this layer at 160 to 170 km coincides with the lower P-wave velocities as seen in Figure 3.8. Line 2 shows a high Poisson's ratio for the top 10 to 15 km of the central portion of the line. This zone coincides with the high P-wave velocity anomaly (Figure 3.14) and may be indicative of the mafic content of those rocks. Line 3 shows a small zone of anomalously low <r close to the northern shot point. The southern half of the line shows a low cr zone varying from 12 km depth at 160 to 280 km to 3 km depth between 280 km and the end of the line. This region corresponds with the fading out of the low velocity zone for this crustal section (Figure 3.21). The western end of Line 4 shows a low <t zone between -20.0 and 100.0 km down to a depth of 22 km. However, the data demand that a high Poisson's ratio be included immediately beside this zone. This second anomaly extends from 105 to 155 km at a depth of 5 to 20 km and coincides with the upwarping of the P-wave velocity contours (Figure 3.28) under the Val Rita block. Although Line 5 did not display any prominent first arrivals, shot G did successfully image reflector I (Figure 3.35) of the P-wave analysis. A n enlargement of that reflection phase is shown in Figure 3.50 with travel-times superimposed for different or models. Figure 3.50a shows the calculated travel-times for reflections off zones I and II of Figure 3.38 for a cr of 0.25. Event II is discontinuous and has not been imaged well, but I is a clear event and was used to infer lateral variations in cr near the shot point and the ILCZ at 170 km. The agreement between the reflection arrivals and the calculated travel-time curves is poor at far offsets (90 to 125 km model distance) with the calculated curves also delayed at nearer offsets (135 to 150 km model distance). By altering the cr  Chapter 3. Seismic Data Analysis  143  in  80  90  100  110  80  90  100  110  120  130  140  150  160  170  120  130  1 40  150  160  170  DISTANCE (KM)  IT)  DISTANCE (KM)  Figure 3.50. Shear wave reflection from western segment of shot G of Line 5 with calculated travel-times superimposed. Both sections are bandpass filtered between 1.5 and 9.0 Hz. (a) I and II are calculated reflection curves from reflectors marked in Figure 3.38 for & =0.25. (b) I indicates curves for three different cr models; see text for discussion. in the vicinity of the shot, better agreement between the observed and calculated arrivals was obtained (Figure 3.50b). The dotted curve represents a model with low cr values extending from 140 to 220 km and down to the reflector at 13 to 16 km. The arrivals at the far offsets are too early and so the possibility of low Poisson's ratio values extending into the KSZ was eliminated. The dashed curve is for the model shown in Figure 3.51b. The travel-time agreement is acceptable except for the 105 to 120 km range. The best agreement is shown by the solid curve in Figure 3.50b for the model in Figure 3.51a. This low u zone terminates at the outcrop of the ILCZ and may differentiate between a  Chapter 3. Seismic Data Analysis  W -30  20 1  70 1 •  144  DISTANCE (KM) 120  170  1  220  \y  270 1  320  370  1  1  30  420  20 1_  DISTANCE (KM) 70 ' •  120 I  170  220  j»| OA  270 1  *320 1  370 1  420  0.23  •XL  t—  CL  LINE 5  a  LINE 5  Figure 3.51. (a) and (b) are two Poisson's ratio models for Line 5. The form of the calculated travel-time curves in Figure 3.50 for these models is shown under each panel. For discussion see text. The inverted triangles represent the recording spread, the open circle represents the shot point, the dashed line is the imaged portion of the reflector, and the contours represent cr values of 0.235 and 0.245. high quartz content region under the shot and a more mafic region west of the ILCZ in the Chapleau block.  Chapter 4  Geopotential Field Data Analysis  The gravity and magnetic data over the Kapuskasing region were analysed to provide further constraints on a physical model for the structure. Seismic velocity and density are related properties of a rock type and so modeling their subsurface distributions can help to discriminate between different lithologies. The magnetisation distribution, which may be related to the density distribution (see Poisson's relation in section 4.2.1), provides additional information about the mineralogy of the causative bodies and helps delineate fault zones. Two dimensional forward gravity modeling [Talwani et al., 1959] was performed along the seismic refraction profiles. The starting density models were obtained from the seismic models using a velocity-density function based,on the Nafe-Drake relationship [Ludwig et al., 1970]. The 2d areal data sets were processed using a series of frequency domain filtering techniques. This analysis was performed in order to separate regional and residual anomalies, to delineate structural trends, and to help relate the causative bodies of the potential fields to one another and to the seismic structure.  4.1  Gravity Profile Modeling  Any model of the Earth must be consistent with all the available data. Thus, the seismic structure obtained from analysis of the refraction data had to satisfy the constraints of the observed gravity trends along those profiles [e.g. Spence et al., 1985; Mueller and Ansorge, 1986; Catchings and Mooney, 1988]. 2d forward gravity modeling was used to analyse the gravity variations. The model  145  Chapter 4. Geopotential Field Data Analysis  146  was represented by a series of n-sided polygons each of which had a residual density value. These residuals represented the differences in value from a prespecified average crustal density. The initial density values were calculated from the velocity models obtained from the refraction analysis. A fourth order polynomial curve was fit to the Nafe-Drake set of velocity-density points to obtain a relationship between these two parameters [Zelt, 1989]. Then, from this curve, densities were obtained for each velocity value in the seismic structure. The density models retained the same boundaries as the velocity models and the density within each trapezoid was the free variable. The gravitational attractions for each polygon were calculated and then added together to obtain the net effect for the model. The observed anomaly curves were constructed from the Bouguer anomaly map (Figure 1.27) over the region. No regional trend was removed from the data as a study of the anomaly and geology maps indicated that the dominant sources were crustal. In the forward modeling procedure, the initial calculated anomaly curve was compared to the observed data for that profile and densities in the models were altered in a trial and error manner to improve the agreement between the two curves. Deviations from the starting densities were accepted as long as they lay within the limits of the scatter of data points on the Nafe-Drake relationship.  Results Line 1, over the Wawa domal gneisses, shows no significant short wavelength features in the observed gravity profile (Figure 4.1). There exists a broadscale low in the central and southern part of the fine over the metavolcanic/metasedimentary sequences and felsic plutons, respectively. The densities had to be increased in the northern part of the model and lowered slightly in the south. There is a small observed gravity anomaly at the same location as the near-surface velocity anomaly at 140 km. This may be associated with an extension of the Lepage Fault. Line 2 has a large positive anomaly at the centre that is associated with the Chapleau  Chapter 4. Geopotential Field Data. Analysis  147  Figure 4.1. Line 1 gravity profile modeling, (a) Initial (dashed curve), final (solid curve), and observed (symbols) gravity anomaly curves for Line 1. The arrow heads mark the terrane boundaries, (b) Model showing deviations of density values from initial structure. The 'plus sign' shading shows where the density has had to be increased, and the 'minus sign' shading shows where it has been decreased. The P-wave velocity contours (in km/s) are also included for reference purposes. The stars are the shot locations from the refraction experiment and the arrow heads are the locations of intersections with other lines. The short dashed.lines represent the reflectors imaged in the P-wave modeling.  Chapter 4. Geopotential Field Data Analysis  148  block (Figure 4.2). Densities were lowered over the northern part of the block while higher densities were required further south than the high velocity region. This suggested a net shift south of the density anomaly from the velocity anomaly. No other significant variations were demanded by the gravity data. Line 3 has two noticeable short wavelength features (Figure 4.3). The first between 40 and 100 km overlies a metavolcanic/metasedimentary sequence of the Abitibi Greenstone belt. The second between 180 and 220 km is narrower in extent and may be related to a nearby felsic pluton. These upper crustal lateral variations in density may indicate that the low velocity zone of the seismic model is discontinuous on this north-south profile as it cuts across the east-west structural grain of the central Superior Province. The profile for Line 4 is fairly flat apart from the high (20 mgal) at 130 km over the Val Rita block (Figure 4.4). This is bounded by the Lepage and the Sagansh Lake faults and is coincident with a step-up in the velocity contours. Significantly, there is no anomaly over the Groundhog River block. Line 5 is marked by a strong anomaly (70 mgal peak to peak) over the Chapleau block at 100 to 170 km (Figure 4.5), a low on the east side of the ILCZ (170 km) where the density values were adjusted lower at depth, and a positive anomaly region to the - 3  east where higher densities values (2800 to 2900 kg m ) had to be used. The previous modeling across this line [Percival and Card, 1983] produced a simpler model as no velocity constraints were available at that time. However, the broad trends of our model agree well in the upper and mid-crust. The lower crust extension of the high density region in the Percival and Card model is not required by the gravity data nor is an associated high velocity anomaly imaged by the refraction data. In the above modeling, only the upper crust densities were varied. Effects due to variations in properties of the lower crust and thickness of the crust are relatively small due to the inverse square fall-off of gravitational fields. Thus, considering the predominant wavelength of the observed anomalies and the desire to deviate minimally from the initial  Chapter 4. Geopotential Field Data Analysis  149  Figure 4.2. Line 2 gravity profile modeling, (a) Initial (dashed curve), final (solid curve), and observed (symbols) gravity anomaly curves for Line 1. The arrow heads mark the terrane boundaries, (b) Model showing deviations of density values from initial structure. The 'plus sign' shading shows where the density has had to be increased, and the 'minus sign' shading shows where it has been decreased. The P-wave velocity contours (in km/s) are also included for reference purposes. The stars are the shot locations from the refraction experiment and the arrow heads are the locations of intersections with other lines. The short dashed lines represent the reflectors imaged in the P-wave modeling.  Chapter 4. Geopotential Field Data Analysis  150  Figure 4.3. Line 3 gravity profile modeling, (a) Initial (dashed curve), final (solid curve), and observed (symbols) gravity anomaly curves for Line 1. The arrow heads mark the terrane boundaries, (b) Model showing deviations of density values from initial structure. The 'plus sign' shading shows where the density has had to be increased, and the 'minus sign' shading shows where it has been decreased. The P-wave velocity contours (in km/s) are also included for reference purposes. The stars are the shot locations from the refraction experiment and the arrow heads are the locations of intersections with other lines. The short dashed lines represent the reflectors imaged in the P-wave modeling.  Chapter 4. Geopotential Field Data Analysis  151  Figure 4.4. Line 4 gravity profile modeling, (a) Initial (dashed curve), final (solid curve), and observed (symbols) gravity anomaly curves for Line 1. The arrow heads mark the terrane boundaries, (b) Model showing deviations of density values from initial structure. The 'plus sign' shading shows where the density has had to be increased, and the 'minus sign' shading shows where it has been decreased. The P-wave velocity contours (in km/s) are also included for reference purposes. The stars are the shot locations from the refraction experiment and the arrow heads are the locations of intersections with other lines. The short dashed lines represent the reflectors imaged in the P-wave modeling.  Chapter 4. Geopotential Field Data Analysis  152  Figure 4.5. Line 5 gravity profile modeling, (a) Initial (dashed curve), final (solid curve), and observed (symbols) gravity anomaly curves for Line 1. The arrow heads mark the terrane boundaries, (b) Model showing deviations of density values from initial structure. The 'plus sign' shading shows where the density has had to be increased, and the 'minus sign' shading shows where it has been decreased. The P-wave velocity contours (in km/s) are also included for reference purposes. The stars are the shot locations from the refraction experiment and the arrow heads are the locations of intersections with other lines. The short dashed lines represent the reflectors imaged in the P-wave modeling.  Chapter 4. Geopotential Field Data Analysis  153  model, varying the upper crust densities seemed most sensible and required the smallest changes. All figures show the calculated gravity profiles for the initial density models. The disagreement between these curves and the observed data seemed substantial but the required density changes were realistic. The velocity-density pairs which differed from the starting values were plotted along with the polynomial fit to the Nafe-Drake data set (Figure 4.6). The upper and lower bounds of the scatter on the Nafe-Drake points "*  W  i—i—i—i—i—|—i—i—i—i—|—i—i—i—i—|—i—i—i—i—]—i—i—i—i—|—i—i—i—i  I y'\ 3  i  i  L_i i 4  i  i  ! 5  i  i  i  i  I  i  i  i  i  6  L_j i 7  i  i  I 8  i  i  i  i— 9  P wave velocity (km/s) Figure 4.6. Comparison of chosen velocity-density function (solid curve) and upper and lower bounds of Nafe-Drake data points (dashed curves) to velocity-density values (points) for the final models. were also plotted to show that our deviations were within previously observed variations except for two points with velocities of 7.45 and 7.6 km/s. There are few Nafe-Drake points in this region; however, recent laboratory work on gneisses and anorthosites from Kapuskasing [Fountain et al., 1989] has shown a wide value of velocity-density points  Chapter 4. Geopotential Field Data Analysis  154  3  (2900 to 3200 kg m ~ and 7.3 to 7.5 km/s) to exist for this region of the relationship. The densities, at depth, for the Wawa domal gneisses and Kapuskasing Zone high grade gneisses and ultra-mafics agree with the measured values [Fountain and Salisbury, 1986] for those rocks at corresponding confining pressures. The main features of the density models mirror the important velocity anomalies of the seismic structures while satisfying the observed gravity profiles.  4.2  2d Potential Field Analysis  In the following sections, the theoretical basis to the filtering will be presented in a framework of linear operations. Then, the construction and application of the filters will be discussed along with the resultant anomaly maps.  4.2.1  Potential Field Theory  The gravity and magnetic fields are presented in a similar style after the formulation of Gunn [1975]. This presentation (see Appendix A) facilitates comparisons and the application of linear transformations. The spectral representations of the potential fields show that the mathematical expressions describing them are the result of a convolution of factors which depend on the geometry of the causative body, the physical properties of the body, and the type of field being observed. Thus, it is possible to remove or alter these factors by the application of linear filters and create the desired output field.  Chapter 4. Geopotential Field Data Analysis  155  Gravity and Magnetic Fields The spectral representations of the gravitational ( F ) and magnetic ( F ) fields as meas  m  sured at a height, h, above the ground surface are  F (u,v,h) g  F (u,v,h)  = 2irG  = 2ir[Lju -f Mjv + N(u  m  2  +C  h  2+v2  f °p {u,v,z)e-^- ^ ^dz v  Jo 2  2  2  + v )*][lju +mjv + n(u + v )*]. +0  f °m (u,v,z)e-^^dz (u + v p Jo 2  v  2  (4.1a, 6) where G is the gravitational constant and p (u,v,z) v  transforms of the bulk density (p (x,y,z)) v  and m (u,v,z) v  are the Fourier  and magnetic (m (x,y,z)) distributions, rev  spectively. Poisson's Relation Assuming the same causative body exists for both potential fields and that the density and magnitude and angle of magnetisation are constant throughout the body, the gravity and magnetic potentials can be related by Poisson's relation [Grant and West, 1965]  m .VU a  =GV  (4.2)  Pv  which gives m^iL + j M + kN).{M +j^-+k^) ox oy  = GV Pv  oz  where m is the magnetisation vector, U and V are the gravity and magnetic potentials, a  and L, M, and N are the direction cosines of the direction of magnetisation. Then  Chapter 4. Geopotential Field Data Analysis  156  which is in the same form as the original magnetic potential equation (see Appendix A) and indicates that if ^p£- is constant, then we can reduce each potential field to the equivalent form of the other field. This is normally performed as the reduction-to-pole operation on the magnetic field. The form of our spectral equations describing the gravity and magnetic fields is amenable to the application of linear transformations. Considering measurements taken at the surface of the earth (h=0), we get  F (u,v, 0) = 2TG /  p (u,v,z)  g  H(u,v,z)  v  J  F (u,v, m  ° 0) = 2TT DAU,V)  dz  , oo I m (u,v,d)  (4.4a,6)  +  D (U,V)  I(U,V)  2  t  H(u,v,z)  dz  Jo where  2TTG  scaling factors 2TT 2  2  Di(u,v)  = [Lju + Mjv + N(u  D (u,v)  = [lju + mjv + n(u -f v )?]  2  H(u,v,z) p (u,v,z) v  2  uT+va)  = e-'( J >  *  +v )*] 2  factor for direction of magnetisation factor for direction of measurement  depth factor  depth distribution factors  m (u,v,z) J v  I(u,v) =  J-  extra factor distinguishing magnetic from gravity field.  From the convolution theorem, components which are multiplied together in the frequency domain (as above) are convolved in the space domain. Thus, since convolution is a linear operation, the order of convolutions does not matter and we can remove the effect of any of the terms in the expressions by convolution with a suitable inverse filter. These operations are more easily implemented as divisions in the frequency domain and  Chapter 4. Geopotential Field Data Analysis  157  this is the method employed in this analysis.  4.2.2  Construction of 2d Frequency Domain Filters  The formulation of the problem in the previous section facilitates the construction and application of frequency domain filters. In the following discussions, all terms representing the filters are given in their frequency domain form.  1  Preprocessing Before a 2d Fast Fourier Transform ( F F T ) is applied to the data, certain procedures are followed to minimise the problems caused by d.c. offsets and sharply truncated edges. First, a least squares plane z = a + ax + ay 0  a  3  is removed from the data where o ,Oi,and o are the plane coefficients and x ,y and 0  3  z are the cartesian coordinate system. Any residual mean after this operation is also removed.  A cosine bell taper is then applied to the edges of the data set to reduce  truncation effects in the spectrum. Finally, the data is padded to over twice its lengths in both directions to increase resolution in the frequency domain and reduce wrap around effects [Brigham, 1974]. Further edge effects can be reduced by rotating a circular window of the data and summing the rotation corrected spectra [Ricard and Blakely, 1988]. This eliminates spectral biases along the axes but is a time consuming and expensive procedure that is not deemed neccessary here.  Power Spectra 2d linear or log spectral representations are obtained to study the power distribution, spectral directionality, and potential aliasing effects. A radially averaged power spectrum is calculated to estimate approximate depths to causative bodies [Spector and Grant, 1  For space domain representations see Grant and West [1965] and Darby and Davies [1967].  Chapter 4. Geopotential Field Data Analysis  158  1973; Hahn et al., 1976 ] and provide power estimates for the white noise in the Wiener optimisation of the derivative and continuation filters [Clarke, 1969; Gupta and Ramani, 1980].  Reduction to Pole Baranov [1957] showed that it is theoretically possible to reduce a magnetic field to its form at the poles by assuming that the magnetisation vector is constant for all magnetised masses in the area under analysis. This reduction to pole procedure has many formulations [Baranov and Naudy, 1964; Bhattacharyya, 1965; Kanasewich and Agarwal, 1970; Silva, 1986] but the representation as a linear transformation [Gunn, 1975] is best suited to the development of section 4.2.1. Thus, the reduction to pole can be viewed as a linear filtering procedure where the effects of the direction of magnetisation and measurement are removed to reduce the field to a vertically acting field. The inverse filter applied is „, , F {u,v)~ r  2  2  (u +v ) D (u,v) D (u,v) y  where Di(u,v)  and D (u,v) 2  2  are as specified in equations 4.4a and b. This filtered output  can be further reduced to a pseudogravity field using Poisson's relation (see equation 4.20) and compared directly to the gravity field.  Bandpass and Directional 2d bandpass niters are constructed in an unfolded form. This allows for a simple 'doughnut' shape to be calculated using the filter function  r i /  x  f 1,  if  [_ 0,  otherwise  2  2  < (i* +t; )^ < ^  2  where v and v are the low cut and high cut frequencies of the passband. Half-cosine x  2  bell tapers are used to softly truncate the edges of these filters with typical roll-offs of 3  Chapter 4. Geopotential Field Data Analysis  db per 0.005 k m  - 1  159  .  Fuller [1967], and more recently Thorarinsson et al. [1988], transferred the concept of velocity or pie slice filtering to potential field data to enhance the directionality of certain anomalous features. This procedure is applied to the unfolded spectra using the following relationship  F (u,v) p  f 1, =1 [ 0,  if a < arctan(^) < /3 otherwise  where a and /3 are the angles denning the 'pie slice'. The edges are softly truncated by half-cosine bell tapers with typical roll-offs of 3 db per 5°.  Derivative First and second vertical derivatives are used to study the slopes and inflection points of the potential surfaces. The form of this filter is obtained from Laplace's equation [Darby and Davies, 1967] and is F (u,v) d  2  2  n  = [(u + v f>]  where n is the order of the derivative. First derivative analysis indicates the steepness and sharpness of an anomalous field. It is used for locating the edges of bodies in gravity and reduced to pole magnetic data and for phase shifting inclined magnetic fields to align peaks over their causative bodies. The vertical derivative is more sensitive to changes in densities of near-surface rocks than to those buried at depth. Second derivative analysis shows the curvature of the potential field surface. The zero values of the second derivative filtered output are commonly used to to delineate faults and body boundaries.  Continuation Continuation of the potential fields to other datum levels facilitates a separation into regional and residual components.  Inspection of equations 4.4a and b show that the  Chapter 4. Geopotential Field Data Analysis  depth factor, H(u,v),  160  accounts for the geometric fall-off of the field above the causative  surface. Thus, a continuation filter can be represented by  d  u  +  2  2  e ^ " ) f o r upward continuation F (u,v) c  2  e  2  +d(u +v )2  £  Qr  downward continuation  The first case upward continues the field above the measurement datum and has the effect of damping the higher frequencies. This is a stable procedure because the filter decays as the height of continuation increases. The downward continuation operator cancels the effect of the depth factor, H(u,v),  and accentuates the higher frequencies. It is a stable  procedure providing that the depth is less than that to the causative body. However, once that depth has been passed, the output becomes unstable as the exponential function is now large, oscillatory, and non-decaying.  Wiener Optimisation The above forms for the responses of the derivative and continuation niters are ideal ones. When filtering real data, noise effects can be accentuated by the application of the filter. The application of Wiener optimum filter theory to the problem [Clarke, 1969; Gupta and Ramani, 1980; Jacobsen, 1987] helps reduce erratic filter output. The Wiener filter seeks to minimise the least squares error between the output and the desired signal and can be represented in the frequency domain by  F (u,v) w  -  P,(u,v) + P (u,v) n  where P (u,v) t  and P (u,v) n  are the power of the signal and noise components respectively.  Obviously the noise power at each frequency cannot be estimated, and so a white noise process is assumed which allows for the calculation of one noise power for the whole  Chapter 4. Geopotential Field Data Analysis  spectrum. This value for P (u,v) n  161  can be obtained from visual inspection of the radially  averaged log power spectrum of the data set.  Cross-Correlation Combined analysis of gravity and magnetic fields in the frequency domain can provide information about the coherence of two data sets assuming that the same causative bodies exist for both fields. Poisson's relation and the reduction to pole procedure allow us to directly compare the two potential fields. By a cross-correlation operation, a quantative estimate of the degree of coherence can be obtained. This operation is given by  Fcro.. (u, v) = F  where * denotes the complex conjugate.  m  (u, v).F*(u, v)  The data sets can be prenltered using band-  pass and directional filters before being input to this cross-correlation operation. Thus, r causative bodies of a given scale and orientation can be quantitatively compared.  4.2.3  Results of Analysis  Results of 2d filtering of the two potential fields are presented separately here.  An  interpretation of their joint geological significance will be presented in Chapter 5. As mentioned earlier, to facilitate in the relating of potential anomalies to their causative bodies, an overlay of the main geological and structural features is provided in a pocket at the back of the thesis. This map is shown in Figure 4.7 where the geological abbreviations are explained. Gravity Data There are several strong positive anomalies (-40 to 40 mgal) that are associated with the K S Z (Figure 4.8). The Chapleau block (CB) is marked by a 125 km by 40 km  Chapter 4. Geopotential Field Data Analysis  162  Precambrian Province Boundary ^ •** Subprovince Boundary «"**  Inferred Subprovince Boundary Edge of Phanerozoic Cover Thrust Fault  Normal Fault  Figure 4.7. Geological map of region covered by gravity and magnetic data. The main features are: Southern and Grenville provinces to the south; Phanerozoic cover of the Moose River Basin to the north; greenstone belts (Wabigoon, Michipicoten, Abitibi); gneiss belts (English River, Quetico, Opatica); Wawa domal gneisses; Remi Lake Zone \RLZ[\ blocks (Chapleau [CB], Groundhog River [GRB], Val Rita [VRB], Fraserdale Moosonee [FMB]); normal fault (Lepage [LF]); and thrust fault (Ivanhoe Lake Cataclastic Zone [ILCZ]).  Chapter 4. Geopotential Field Data Analysis  163  Figure 4.8. Bouguer gravity anomaly map for the Kapuskasing region. The units are in mgals.  Chapter 4. Geopotential Field Data Analysis  164  anomaly that strikes NE-SW between 83.5° and 82.5° longitude and 47.5° and 48.5° latitude. Notably, the Groundhog River block (GRB) does not have an anomalous gravity signature but, instead, the gravity high of the C B is continued over the Val Rita block (VRB) further to the west.  This anomalous feature continues up to the Fraserdale-  Moosonee block (FMB) between 82.0° and 81.0° longitude and 49.5° and 50.5° latitude. The region with a gravity high between the V R B and the F M B has been designated the Remi Lake Zone (RLZ). A broadscale east-west high is observed under the southern margin of the Moose River Basin at 50.5° latitude to the west of the F M B . The boundary between the Superior and Southern Provinces and the Grenville Front, in the SE corner of the map area, are marked by several positive anomalies The log power spectrum (Figure 4.9) shows a concentration of power below 0.04 km  - 1  and spectral directionality with a pointing vector of 145°± 5°. This suggests the  existence of strong features in the data with dimensions of 25 km or greater and a strike of approximately N45°E. The radially averaged log power spectrum (Figure 4.10) displays average powers for the wavenumber range and shows spectral roll-off of 23-24 db per 0.01 1  km" . A low pass filter with a cut-off of 0.03 k m  - 1  and a cosine taper of 0.01 k m  - 1  was ap-  plied to the data (Figure 4.11). The result shows that the dominant features of the data, as described above, have characteristic wavelengths of 25 km or more. The suggested directionality of the power spectrum is exploited by the application of a directional enhancement filter with angles of 120° to 150° and a taper of 10° (Figure 4.12). This shows the delineation of two strong NE-SW trending edges. The Ivanhoe Lake Cataclastic Zone (ILCZ) is highlighted from 84.0° longitude and 47.0° latitude to 81.5° longitude and 50.0° latitude. The other feature reveals an extension of the Lepage Fault (LF) with a strike parallel to the ILCZ. It extends from 83.0° longitude and 49.0° latitude to 84.5° longitude and 48.25° latitude. This result suggests that the Kapuskasing uplift is fault-bounded further to the north-west of the K S Z itself. A directional filter with angles of 90° to 150°  Figure 4.9. Log power spectrum for the Bouguer gravity anomaly map. The wavenumber axes units are k m . - 1  Chapter 4. Geopotential Field Data Analysis  166  Chapter 4. Geopotential Field Data Analysis  Figure 4.11. Gravity map low pass filtered with a cut-off of 0.03 k m taper of 0.01 k m . The units are in mgals. - 1  167  - 1  and a cosine  Chapter 4. Geopotential Field Data Analysis  168  shows the strongest feature to be the RLZ. A filter with angles of 70° to 110° highlights the strong low wavenumber east-west trends of the Superior Province (Figure 4.13). In particular, the belts in the northern part of the map area are accentuated.  At 49.1°  the boundary between the Abitibi and Quetico belts is prominent. Even though further north is covered by the Moose River Basin (MRB), there is evidence for the Wabigoon volcanic (the low at 49.8°) and the English River gneiss (the high at 50.2°) belts that are outcrop west of the map region. Two ridge highs cross the Wawa and Abitibi belts with a northward displacement on the east side of the KSZ. The east-west high along the southern end of the map follows the Penokean belt at the northern margin of the Southern province. The first derivative operation shows a residual field and emphasises the near surface density variations (Figure 4.14). The edges of the C B and V R B are highlighted and clearly show the north-westerly displacement of the KSZ high. There are highs observed at northern end of C B , the west side of the V R B , in the centre of the R L Z , and in the F M B . The boundary with the Southern Province, and the Grenville Province behind it, is also prominent and agrees well with previous first derivative studies of that area [Gupta and Grant, 1982]. Downward continuation of the field yields information on the possible depths to the tops of the main causative bodies (Figure 4.15). A depth of 5.0 km is deemed the maximum depth of burial for anomalous bodies of 25 to 30 km dimension under the KSZ. Below this continuation depth, the output is unstable and oscillatory. The strong anomalies are at the northern end of the F M B , over the RLZ, and over the C B . Further downward continuation gave unstable results but the C B remains a strong feature suggesting that the anomalous densities continue down to the lower crust. The second derivative filter outlines the main segment of the ILCZ and delineates the bounding faults of the F M B (Figure 4.16). Upward continuation of the data emphasises the regional trends and indicates gross crustal variations (Figure 4.17 and 4.18). At a height of 10 km the C B , V R B and R L Z are a continuous strong feature. At 50 km, there  Figure 4.12. Gravity map directionally filtered between 120.0° and 150.0° with a taper width of 10.0°. The units are in mgals.  Figure 4.13. Gravity map directionally filtered between 70.0° and 110.0° with a taper width of 10.0°. The units are in mgals.  Figure 4.14. First derivative map of gravity data with a Wiener optimisation factor of 1.0. The units are in mgals/km.  Figure 4.15. Gravity map downward continued to a depth of 5.0 km with a Wiener optimisation factor of 1.0. The units are in mgals.  Figure 4.16. Second derivative map of gravity data with a Wiener optimisation factor of 1.0. The units are in mgals/km/km.  Figure 4.17. Gravity map upward continued to a height of 10.0 km with a Wiener optimisation factor of 1.0. The units are in mgals.  Figure 4.18. Gravity map upward continued to a height of 50.0 km with a Wiener optimisation factor of 1.0. The units are in mgals.  Chapter 4. Geopotential Field Data Analysis  176  is a suggestion that the high is centred over the R L Z and F M B and that a broadscale low extends over the southern Wawa and Abitibi terranes and across the southern part of the KSZ. This low may be a manifestation of the deep crustal root imaged in the refraction analysis.  Magnetic Data Figure 4.19 shows the aeromagnetic anomaly map for the Kapuskasing region.  The  Groundhog River block between 82.75° and 82.0° longitude and 48.5° and 49.5° latitude and the F M B between 81.25° and 80.75° longitude and 50.0° and 51.0° latitude are the dominant features of interest. The southern margin of the Moose River Basin, trending east-west between 50.0° and 50.5°, and the boundary with the Southern province in the SE corner also have strong associated anomalies. The log power spectrum (Figure 4.20) shows high levels of background noise which is evident in the high frequency components of the field data. The strong spectral directionality indicates a pointing vector of 107°± 5°. This translates to a predominance of bodies striking E17°N and represents the east-west structural trends of the central Superior Province. The radially averaged log power spectrum (Figure 4.21) again illustrates - 1  the high levels of white noise and a low spectral roll-off of 8 to 9 db per 0.01 k m . The reduced to pole (Figure 4.22) and pseudogravity transformations have very similar features to the original field. This result is due to the steep inclination of the ambient field (~78°) and the assumption that the remanent field is parallel to it. A bandpass of 0.015 k m  - 1  to 0.04 k m  - 1  with a cosine taper of 0.01 k m  - 1  width optimally highlights the  features of interest (Figure 4.23). The G R B , the mapped L F and its extension into the R L Z are prominent in this pass band. Using a low pass filter with a cut-off of 0.05 k m  - 1  and a directional filter with angles of 120° and 150°, the northern portion of the KSZ is highlighted (Figure 4.24). The F M B appears to have two extensions to the south-west: the first extends through the G R B along the ILCZ and extends south of the C B ; the  Figure 4.19. Aeromagnetic anomaly map for the Kapuskasing region. The units are in nanotesla.  Chapter 4. Geopotential Field Data Analysis  178  Chapter 4. Geopotential Field Data Analysis  179  Chapter 4. Geopotential Field Data Analysis  180  Figure 4.22. Aeromagnetic map reduced to pole from an induced and remanent field with an inclination of 78° and a declination of 10°. The units are in nanotesla.  Figure 4.23. and 0.04 k m  - 1  Aeromagnetic map bandpass filtered with corner wavenumbers of 0.015 and a cosine taper of 0.01 k m . The units are in nanotesla. - 1  - 1  Figure 4.24. Low pass (0.04 k m ) and directionaUy (120.0 and 150.0°) filtered aeromagnetic map. The units are in nanotesla.  1  Figure 4.25. Aeromagnetic map with a low pass filter (0.04 km ) and a directional filter between 70.0° and 110.0° with a taper width of 10.0°. The units are in nanotesla.  Figure 4.26. First derivative map of aeromagnetic data with a Wiener optimisation factor of 8.0. The units are in nanotesla/km.  Figure 4.27. Aeromagnetic map downward continued to a depth of 0.3 km with a Wiener optimisation factor of 8.0. The units are in nanotesla.  Chapter 4.  Geopotential Field Data Analysis  186  Figure 4.28. Aeromagnetic map upward continued to 10.0 km elevation with a Wiener optimisation factor of 8.0. The units are in nanotesla.  Chapter  4.  Geopotential  Field Data. Analysis  187  Longitude  Figure 4.29. Aeromagnetic map upward continued to 50.0 km elevation with a Wiener optimisation factor of 8.0. The units are in nanotesla.  Chapter 4.  Geopotential Field Data. Analysis  188  second follows the mapped L F and there is a weak suggestion of an extension to the feature at 83.5° longitude and 48.0° latitude.  There is a strong linear feature between  82.5° and 80.0° longitude and 46.0° and 47.0° longitude. The same low pass filter and a directional filter with angles of 70° and 110° shows several strong east-west trends (Figure 4.25). These are the magnetic highs observed along the boundary with the Southern provinces (at 47.0°), over some metavolcanics in the Abitibi (at 48.25°), along the southern edge of the Quetico belt (at 49.25°), and at the southern margin of the Moose River Basin (50.1°).  First derivative analysis indicates that the near-surface magnetisations  are greatest in the G R B , at the southern edges of the Moose River Basin and Quetico Belt, and at the southern end of the C B (Figure 4.26). Downward continuation of the field to a depth of 0.3 km (i.e., to the surface) and application of a low pass filter with a cut-off at 0.04 k m  - 1  emphasises the G R B and F M B anomalies and the southern edge  of the Wabigoon greenstone belt (Figure 4.27). However, further continuation results in oscillatory output suggesting that many of the causative bodies are at or near the surface. Upward continuation of the field to 10.0 km produces a strong positive anomaly between the G R B and RLZ (Figure 4.28) and emphasises the low immediately to the west of the C B . There is a strong linear high across the edge of the Moose River Basin as well as an anomaly associated with the Grenville Front. Continuing up to 25.0 km results in a similar pattern but above this, and up to 50.0 km, there appears to be a high under the Moose River Basin and F M B area accompanied by a low under the Quetico belt, the southern Wawa and Abitibi belts and the C B (Figure 4.29). This is similar to the pattern seen at this height for the upward continued gravity data and again may be linked with the seismic crustal thickness.  4.3  Correlation of Potential Fields  A combined interpretation of the main features highlighted by filtering of the gravity and magnetic data shall be presented in chapter 5. The purpose of this section is to  Chapter 4. Geopotential Field Data Analysis  189  quantitatively estimate correlations between the two data sets. Data adaptive filtering was performed on both data sets (Figure 4.30 and 4.31). This involved taking the wavenumber representation of one data set and multiplying it by the amplitude spectrum of the other data set.  This produced an output that highlighted  common spectral components but with the phase of the original data. Figure 4.30 shows the gravity data filtered by the amplitude spectrum of the aeromagnetic data. Prominent features are the F M B , the RLZ, the V R B , the northern part of the C B , and a broad low across the Wawa and Abitibi regions under the KSZ. Figure 4.31 shows the aeromagnetic data filtered by the amplitude spectrum of the gravity data. The M R B , the RLZ, the south-eastern part of the CB, and the broad low to the west of the KSZ were emphasised. The lows over the Quetico, Wawa, and Abitibi belts are also seen in the upward continued data (Figures 4.18 and 4.29) and indicate a thickened lower crust with little magnetisation under these terranes.  The variable highs along the whole length of the K S Z indicate  strong density and magnetisation variations in the uplift crust and emphasise the different faulting geometries along strike. The F M B and R L Z are both strong anomalies, but no seismic control on their structures is available. Next, a cross-correlation of the data sets was performed.  Figure 4.32 shows the  result with a continuous high extending east-west 75 km to the north and north-west of the origin. The reduction to pole of the magnetic data increased the lateral continuity along this broad high. The output shows that there is a net displacement of the strong anomalies of the data sets with respect to each other. This could be due to several causes: positioning errors; different causative bodies; or the same causative bodies with remanent magnetisation causing a displacement of the anomaly peaks. The displacement distances are too large to seriously consider the first solution. Different causative bodies is indeed the most reasonable solution from consideration of the outcropping features. However, to investigate the last possibility, the data were filtered prior to cross-correlation to enhance features of interest.  A low pass filter with cut-off of 0.03 k m  - 1  and a directional filter  Figure 4.30.  Gravity data filtered by amplitude spectrum of aeromagnetic data. See  text for explanation.  Figure 4.31.  Aeromagnetic data filtered by amplitude spectrum of gravity data. See  text for explanation.  Figure 4.32. field.  Cross-correlation of unfiltered gravity and the reduced to pole magnetic  Figure 4.33. Both frames show cross-correlation of gravity and magnetic fields with a low pass filter below 0.03 k m and a directional filter between 70°.and 110°: (a) magnetic field reduced to pole for induced and remanent fields with an inclination of 78° and a declination of 10°; (b) same as (a) but for remanent field with an inclination of -10° and a declination of 0°. - 1  Chapter 4.  Geopotential Field Data Analysis  194  with angles of 70° and 110° were chosen to study the large scale east-west trends of both data sets. The resultant cross-correlation shows the strong high displaced into the northwest quadrant (Figure 4.33a). Obviously, there is a large coherence (i.e., autocorrelation values nearer 1.0) between the two data sets if the magnetic data are offset south by approximately 75 km. If a remanent field with an inclination of 0° ± 10° and declination of 0° is considered in the reduction to pole operation, the output of the cross-correlation has a broad low running east-west close to the x-axis (Figure 4.33b). Thus, the correct field parameters have been chosen to overlay the gravity and magnetic data and delineate the common causative structures of the central Superior Province east-west trend. These parameters for the remanent field give a pole on the apparent polar wander path for North America [Constanzo-Alvarez and Dunlop, 1988] that gives a date of magnetic stabilisation at 2400 to 2500 million years ago, after the accretion of the southern Superior province (2700-2600 Ma).  However, doubts exist about the stability of remanence over such  lengths of time (2000 to 2500 Ma) and whether such a large volume of the crust could cool down through its Curie temperature while the field maintained a constant polarity [Skive, 1989]. From inspection of geological maps, a more reasonable interpretation is that the magnetic field is predominantly ambient and that the gravity trends follow the belts while the aeromagnetic trends demarcate the belt boundaries.  Chapter 5  Discussion and Interpretation  The primary objectives of the thesis project have been realised. A seismic structure for the Kapuskasing uplift and the surrounding terranes has been obtained from analysis of the refraction data. To complement and corroborate these results, gravity and magnetic field analysis has been performed to outline density and magnetisation variations over the same region. These geophysical models are discussed in the following sections and related to one another, and to possible thermal and rheological conditions that prevailed after the accretion of the Superior craton and at the time of uplift.  5.1  5.1.1  Seismic Results  P Wave Structure  The resultant velocity models from each line are discussed seperately. The resolution and errors within the models have been discussed already in section 3.1 of Chapter 3. Line 1 images the structure of the Wawa domal gneiss terrane (Figure 3.8).  The velocities  range from 6.2 km/s at the surface to 7.65 km/s above the Moho with no apparent large scale velocity inversions. The crust is 48 ± 1 km thick with fairly strong intra-crustal reflectivity and a distinct Moho. The high reflectivity is probably due to the fine scale layering of gneissic textured rocks in shear zones associated with the faulting of this terrane. Line 2, overlying the KSZ, shows surface velocities varying from 6.0 km/s to as high as 6.4 km/s over the Chapleau block and an increase of velocities down through the crust to 7.6 km/s over a diffuse Moho at 50 km depth (Figure 3.14). There is a suggestion that a small velocity inversion at a depth of 25 to 30 km could be included  195  Chapter 5. Discussion and Interpretation  196  but it is not neccessary. The entire crust is characterized by high amplitude reflective scattering that shows some lateral coherence. The upper mantle velocities are in the 8.25 to 8.35 km/s range whereas on the orthogonal intersecting lines the velocities vary from 8.1 to 8.2 km/s. Hence, there is a suggestion of velocity anisotropy below the crust-mantle boundary. Similar anisotropy has been observed in other Archean terranes [Drummond, 1985] and can be explained by syntectonic recrystallisation of olivine in the upper mantle rocks [Ave 'LaJlemant and Carter, 1970]. In this process, the olivine crystals realign themselves in the prevailing stress field with their slow velocity axis parallel to that field. Thus, applying the thrust fault model to the KSZ implies that the main horizontal stress field would have been east-west and one would expect the orientation of anisotropy that is observed. Line 3 overlies the Abitibi greenstone belt which appears to be less faulted than the other terranes of the study area (Figure 3.21). There is a 4 to 6 km surface layer with velocities ranging from 5.8-6.0 km/s at the top to 6.2-6.4 km/s at the base of this layer. A thin upper layer has also been inferred from gravity and geoelectrical studies for the Giyani and Murchison greenstone belts of South Africa [Kleyweyt et al, 1987] and for other greenstone belts in the Canadian shield [Thomas et al., 1986]. Under this layer, there is a low velocity zone which has a thickness of 5 to 10 km with its base marked by strong reflections. There is evidence of a low velocity zone at similar depths in the southeastern Abitibi belt from the 1982 Abitibi-Grenville Front seismic refraction experiment [Parker, 1984]. From preliminary analysis of their magnetotelluric data, Chouteau at al. [1988] have suggested that a conductive zone exists at mid-crustal depths under the belt east of the Groundhog River Block. These results indicate that this low velocity zone may be continuous over regions of the Abitibi greenstone belt and also that its reflective base may be associated with a conductive zone. The crust below this is quite reflective and has a lower average velocity than the other lines. The velocities reach 7.4-7.5 km/s above the Moho at 40 to 45 km depth. The belt appears to be deeper at the center than  Chapter 5. Discussion and Interpretation  197  at either end. Lines 4 and 5 cut across the terranes of Lines 1, 2 and 3 and consequently are structurally more complex but contribute information about the terrane boundaries.  Im-  portant general trends for the whole study region can be discerned from their velocity structures. The surface velocities decrease from 6.0-6.2 km/s over the greenstone belts to 5.8-6.0 km/s over the domal gneisses and then increase to 6.3-6.5 km/s over the high grade granulite zone and a pluton close to shot point G (Figure 3.38). Within the domal gneisses, the tonalitic gneisses comprise 63% of the volume [Fountain et ah, 1989] and it is this rock that has a characteristic low velocity. Figure 5.1 shows a comparison of the proposed geological model of Percival and Fountain [1989] with the corresponding portion of the velocity structure from Line 5. The seismic profile is oblique to the geological cross-section (see Figures 1.1 and 1.10) and so the true dip of the velocity structure is 3°-5° greater than that imaged. A model of uplifted crust is supported by the seismic results but with a soling depth for the thrust at 17-20 km. Mid-crustal velocities of 6.5-6.6 km/s are brought from 15 to 25 km depth under the greenstone belts to 2-5 km under the KSZ. This corresponds to the granulites uplifted from the mid-crust along the ILCZ. In isostatic response to this high density anomaly in the upper crust, the crust has thickened both under and to the west of the Kapuskasing exposure. The geological model in Figure 5.1a shows metavolcanics at the western edge of the section; however, due to the different locations of the profiles, the metavolcanics are not crossed until further to the west in the section of Figure 5.1b; this figure is only used to illustrate the dipping high velocity anomaly and the thickened crust. The reflectors (I and II) marked in the velocity model were imaged by the center shot (Figure 3.35) and are probably associated with shear zones along the plane of decollement for the thrust. The top reflector (I) has also been imaged by the high resolution seismic reflection profile over the ILCZ [Geis et ai, 1988]. The Ivanhoe Lake cataclastic zone, the Wakusimi River fault and the Lepage fault  198  Chapter 5. Discussion and Interpretation  a  |  D I S T A N C E  (KM)  D I S T A N C E  (KM)  ARCHEAN | Plutont, orthoonetss Mclavolconic rocks Gronulite  w 0  fades  42  85  127 I  E 170  CD  •—<C\I  CO S I S -  1  -8.0=0 CD . CD  03 6 . 6 - 6 3 km/s  Figure 5.1. (a) Schematic geological cross-section of KSZ across Chapleau block after Percival and Fountain [1989]. The section is along the profile marked A A' in Figure 1.1. (b) Velocity structure from a portion of Line 5 resolved onto the geological cross-section.  Chapter 5. Discussion and Interpretation  199  all have strong perturbing effects on the travel-times and amplitudes and hence on the resultant velocity structures. Thus, they demarcate important terrane or block boundaries. The Ivanhoe Lake cataclastic zone (160-170 km in Figure 3.38) is the proposed main thrust fault plane for the KSZ. The Wakusimi River fault (90-100 km in Figure 3.14) confirms the division of the KSZ into two blocks, the Groundhog River Block to the north and the Chapleau Block to the south. The southern block (100-190 km) has an anomalously high velocity upper crust, whereas the northern block (20-100 km) is only weakly anomalous. The Lepage fault, which separates the orthogneiss (and probable buried mafic gneiss) in the Val Rita Block of the Wawa terrane from the Quetico metasedimentary migmatites, is imaged as a major feature. The fault may project down to the mid-crust (100-140 km in Figure 3.8 and 110-130 km in Figure 3.28) and also extend further to the south-west (Figure 4.12) than is presently mapped.  5.1.2  S Wave Structure -  Figure 3.49 showed the resultant Poisson's ratio (o ) models from analysis of the shear wave arrivals. In general, the resolution in the upper 20 to 25 km is ± 0 . 0 0 5 but, below this level, there is little resolution as arrivals traversing deeper paths were obscured in noise.  Thus, the values for the mid to lower crust are gross averages inferred by the  occasional intra-crustal reflector and the Moho reflections.  It was possible to alter the  er of the lower crust (>35 km) by 0.02 and still have calculated travel-times within the observational error. Thus, since the data did not demand a change, the starting value of 0.25 was retained. At the northern end of Line 1, there are low cr values for the top 5.0 km out to 80 to 85 km. This portion of the profile lies over the metasedimentary migmatites of the Quetico belt and the low a reflects the high quartz content of these rocks. The southern end of the line has low a values for the top 15 km and coincides with the metavolcanics of the Michipicoten belt. The termination of this layer at distances of 160-170 km coincides  Chapter 5. Discussion and Interpretation  200  with P-wave velocities of 6.8 km/s as seen in Figure 3.8. Line 2 shows a high Poisson's ratio for the top 10 to 15 km of the central portion of the line. This zone coincides with the high P-wave velocity anomaly (Figure 3.14) and may be indicative of a high mafic content for those rocks. Line 3 shows a small zone of anomalously low er close to the northern shot point. The southern half of the line shows a low cr zone varying from 12 km depth at distances of 160 to 280 km, to 3 km depth between 280 km and the end of the line. This region corresponds with the fading out of the low velocity zone for this crustal section (Figure 3.21). The southern end of Line 4 shows a low <r zone between -20.0 and 100.0 km down to a depth of 22 km. However, the data demands that a high Poisson's ratio be included immediately to the east of this zone. This second anomaly extends from 105 to 155 km at a depth of 5 to 20 km and coincides with the upwarping of the P-wave velocity contours (Figure 3.28) under the Val Rita block. Figure 3.51a shows the small portion of the structure for Line 5 that was successfully imaged. This low cr zone terminates at the outcrop of the ILCZ and may differentiate between a high quartz content region (felsic pluton) under the shot and a more mafic region (granulites) west of the ILCZ in the Chapleau block.  5.1.3  Inverse vs Forward M o d e l i n g for Line 5  Line 5 crosses the southern end of the Kapuskasing structural zone and results of the synthetic seismogram modeling show a high velocity anomaly dipping to the west under the Chapleau block (Figure 3.38 and 5.1).  Data analysis using trial and error forward  modeling can be influenced by the preconceptions of the interpreter, and so corroboration of the main features of the model is desirable. The forward modeling technique of iterative ray tracing employed in the analysis of chapter 3 results in models that satisfy the data, but are not neccessarily unique. In order to check the verity of the dipping anomaly in the model for Line 5, a small scale inversion was performed on the first arrivals from shots H and G along the western half of that profile.  Chapter 5. Discussion and Interpretation  201  A tomographic inversion technique, designed for the iterative inversion of refracted first arrivals [White, 1989a, 1989b], was used to obtain an alternative model. In this method, a starting model is introduced and rays are traced to each receiver from each shot.  The residual time differences between the observed and calculated travel-times  are converted to velocity anomalies knowing the length of each ray path. This velocity anomaly is distributed to the pixels along each ray path. The updated model is then used in the next ray tracing iteration. The procedure continues until a satisfactory fit to the data has been achieved.  Because only a coarsely sampled model and a small  number of data were used, a matrix inversion could be utilised. This technique has the advantage that estimates of the resolution and variance in the model parameters can be obtained but has the disadvantage that the matrix size possible for inversion is limited [Worthington, 1984; White, 1989b].  Starting Model A model with a length of 250 km and a depth of 25 km was considered for the inversion procedure. This was sampled every 10 km along the horizontal direction and every 2.5 km in depth. This regular grid of points remained constant throughout the inversion and only the velocities at those points varied. In the trial and error forward modeling 1  approach, both the velocities and block sizes are varied , but in the inversion only the velocities are changed. Thus, odd block shapes used in the forward modeling can result in travel-time artifacts in the ray tracing; this problem is minimised in the inversion by having a regular grid. The shots were located at model distances of 20 km and 214 km at the surface, and the 60 receivers, with an uneven distribution along the length of this profile, had an average spacing of 3 to 5 km. The starting model was that one used for the trial and error forward modeling and is shown in Figure 5.2. 1  See chapter 3 for details  Figure 5.2b shows the rays  Chapter 5. Discussion and Interpretation  202  Figure 5.2. (a) Comparison of calculated (curves) and observed (symbols) travel-times for the starting model, (b) Ray diagram for starting model of inversion. The model consists of a grid of velocity values that form a regular block model. The stars show the locations of the shot points and the arrow heads show the terrane boundaries.  Chapter 5. Discussion and Interpretation  203  traced through the regular grid model and Figure 5.2a shows the observed travel-times (sjmbols) compared to the calculated ones (curves) for that starting model. From a simple visual inspection of the reversed profile, the travel-time differences indicate that higher velocities are needed at greater depths to the west. This qualitatively supports the model of a high velocity anomaly dipping to the west.  Results of Inversion A total of seven iterations was required before a satisfactory result was obtained. At 2  each iteration, the results were damped and weighted in order to create a smooth model for ray tracing and thus eliminate shadow zone problems. Firstly, the magnitudes of the anomalies were damped to avoid spikes in the velocity grid. Secondly, the pixels were weighted inversely proportional to the number of rays passing through them from each source; this prevented the anomalies being concentrated around the shot points. Finally, the output velocity model was smoothed to eliminate sudden jumps in velocity that might cause dramatic changes in ray paths and result in shadow zones at the surface. Figure 5.3 shows the residual velocity anomalies for three of the iterations.  The  starting model had a root mean square (rms) time difference between observed and calculated travel-times of 0.226 s. After the first iteration, this rms value was reduced to 0.133 s. Figure 5.3a shows the resultant residual velocity field with an anomalous high under the Chapleau block that dips to the west. The third iteration reduced the rms value to 0.083 s and shows the same anomaly pattern as the first iteration along with a low at the surface between 35 and 70 km. The final iteration reduced the rms time difference to 0.072s and the result is shown in Figure 5.3c.  The travel-time comparison and ray  diagram for this final model are shown in Figure 5.4. The travel-time discrepancies are fairly evenly distributed over all the receivers, but there is a suggestion that the anomaly high under the Chapleau block could be larger. The ray diagram shows how the final 2  See White [1989] for details  Chapter 5. Discussion and Interpretation  204  Figure 5.3. Residual velocity structures after the (a) first, (b) second, and (c) seventh iterations. The velocity contours show velocity differences from the starting model. Contours are in km/s. The stars show the locations of the shot points and the arrow heads mark the terrane boundaries as named in Figure 5.2.  Chapter 5. Discussion and Interpretation  205  Figure 5.4. Comparison of calculated (symbols) and observed (vertical bars) travel-times for shots H (a) and G (b) for the final iteration. The length of the bars corresponds to the acceptable errors for each travel-time pick, (c) Ray diagram for both shots. The stars show the locations of the shot points and the arrow heads mark the terrane boundaries as named in Figure 5.2.  Chapter 5. Discussion and Interpretation  206  anomalous field has created a blind zone at depth directly under the CB. However, earlier iterations did pass rays through this region so that it is partly constrained. In Figure 5.3c, the high velocity anomalies are concentrated in the western half of the Chapleau block and at 5 to 20 km depth to the west of the structural outcrop. There are low velocity zones at shot point H, another at the surface between 50 and 110 km, and a westward dipping low velocity zone that outcrops around the Ivanhoe Lake Cataclastic Zone. Figure 5.5 shows the resultant velocity models for each of the above iterations. The starting model is shown for reference in Figure 5:5a.  The first and third iterations  (Figures 5.5b and c, respectively) show the main features of a high velocity anomaly under the Chapleau block and at depth under the Wawa and Michipicoten terranes. This anomaly corresponds to the uplifted mid-crustal granulites that outcrop in the KSZ. Figure 5.5d shows the final velocity model with the same high velocity anomaly as above augmented and extending to a greater depth. The low velocity anomaly between 175 and 200 km distance is related to the ILCZ and that part of the Abitibi against which the upthrust crust is butted. There is a low velocity zone at shot point H which 3  is probably associated with near-surface rhyolites . The near surface low in the Wawa 3  terrane between 50 and 100 km is related to a large outcrop of tonalitic gneiss . Finally, the resultant models from the forward and inverse modeling were compared and are shown in Figure 5.6. The broadscale features of a westward dipping high velocity anomaly under the CB, a low velocity zone related to the ILCZ and neighbouring Abitibi belt, and a low velocity region in the Wawa tonalitic gneisses are seen in both models. However, details of the velocity values and boundary locations do differ and indicate that the details of the dipping anomaly may not be resolvable with the current data set. The ray coverage in both examples does not allow for good definition of the westward extent of the upper crustal anomaly (see Figures 3.30b, 3.33b, and 5.4c). However, the forward 3  From measured laboratory velocities of Fountain et al. [1989]  Figure 5.5. Velocity structure with contours in km/s for (a) starting model, (b) first iteration, (c) third iteration, and (d) seventh and final iteration. The stars show the locations of the shot points and the arrow heads mark the terrane boundaries as named in Figure 5.2.  Chapter 5. Discussion and Interpretation  208  Figure 5.6. (a) Velocity structure after seven iterations of inversion scheme, (b) Portion of velocity structure from forward modeling of Line 5. The stars show the locations of the shot points and the arrow heads mark the terrane boundaries as named in Figure 5.2.  Chapter 5. Discussion and Interpretation  209  model was constrained by other data. The profile was modeled using not only the traveltimes of first arrivals, but the travel-times of secondary arrivals and the amplitudes of all these arrivals. In addition, the intersecting lines at the C B and shot point H provided additional constraints to the velocity model.  The inverse model shows the C B high  velocity anomaly centered on the western portion of that block. This part of the profile overlies the Robson Dome which is a mafic gneiss unit.  These rocks have very high  measured velocities of 6.5-6.6 km/s at 20 MP a [Fountain et al., 1989]. The depth extent of the anomaly, down to 15 to 20 km in the west, has velocities that are 0.1 km/s less on average than those of the forward model. These values may be a truer representation of the velocity structure as the ray model for the inversion is much smoother. The low velocity region at 50 to 100 km is a strong feature in both models and corresponds to the tonalitic gneisses of the Wawa terrane.  5.1.4  C o m p a r i s o n of Refraction and Laboratory Measured Velocities  Laboratory measurements of velocity and density for the rock suites in the Michipicoten, Wawa, and Kapuskasing terranes have been made at confining pressures ranging from 10 to 600 MPa [Fountain and Salisbury, 1986; Fountain et al., 1989]. In general, there is good agreement between the lithological suite averages from the laboratory samples and the corresponding units in the seismic refraction models.  Importantly, a value of 6.55  km/s for the Kapuskasing gneisses at 100 MPa (4-6 km depth) is in very good agreement with values in the dipping high velocity anomaly of Figure 5.6. Fountain et al. [1989] have constructed a vertical section of rock types and velocities through the Michipicoten belt to a depth of 20 km. This section is shown in Figure 5.7 where it is compared with a one dimensional velocity profile down through the Michipicoten belt at shot point H . Because the refracted energy travels nearly horizontally in low velocity gradient terranes such as this, the horizontal velocity values for the anisotropic rocks were used for the comparison. The velocities shown for the rock suites should be reduced by 0.05 km/s in  Chapter 5.  Discussion and Interpretation  210  Figure 5.7. Comparison of laboratory measured and refraction model velocities after Fountain et al. [1989]. The schematic vertical section shows the three main rock types and their associated velocities. The dashed lines, with velocities at the right side, mark the velocity values at those depths under the metavolcanics at shot point H of Line 5 (see Figure 3.38).  Chapter 5. Discussion and Interpretation  211  the 5 to 15 km depth range and by 0.1 km/s in the 15 to 25 km range due to thermal effects [Fountain et al, 1989].  The agreement between the refraction and laboratory  measured velocities is good but indicates that the bottom part of the section is less mafic than expected.  The velocity profile value of 6.6 km/s may demarcate a 'Conrad' type  discontinuity as imaged previously in the shield [Berry and Fuchs, 1973; Green et al., 1979] and also corresponds to the feature exposed at the Wawa/Kapuskasing boundary [Percival, 1986]. The presence of high velocity metabasalts in the surface layer and the low velocity tonalites in the underlying rocks, shows the possibility of a velocity inversion existing near the surface as was observed along Line 3. The materials with high measured velocities (mafic gneisses and anorthosites) seem to be probable constituents for the present day lower crust in the shield. Xenolith studies indicate that the petrology of the lower crust may be dominated by garnetiferous minerals and pyroxenites [Rudnick and Taylor, 1987] and laboratory measurements on these rock suites give velocity values in the 7.0 to 7.6 km/s range [Jackson and Arculus, 1984]. Values in this range were imaged at the base of our seismic models. Fountain et al. [1989] have calculated high reflection coeffficients (up to 0.13) betweeen the mafic gneiss/anorthosite  and the paragneiss/diorite/tonalite  units.  These  high values and the average 100 m layering observed in the exposures could give rise to constructive interference of the reflected seismic energy and account for the high amplitude events seen in the reflection data [Geis et al., 1988; West et al., 1988] and the wide angle reflections in the refraction data (e.g., reflector in Figure 3.35). The high frequency (12 to 16 Hz) wide angle reflections observed in the refraction data could be caused by 4  constructive interference from layering with 100 to 140 m thicknesses. These thicknesses are supported by field observations and are probably an upper limit to values found in these terranes. 4  See Hurich and Smithson [1987] for explanation of calculation of interference effects  Chapter 5. Discussion and Interpretation  5.2  212  Potential Field Results  It is important to relate the subsurface distributions of the various physical parameters to known lithologies.  In the first subsection, the density is related to variations in P-  wave velocity and Poisson's ratio (cr) for each of the profiles. All densities are described relative to the fourth order polynomial fit to the Nafe-Drake relationship (Figure 4.6), as discussed in the previous chapter. The second subsection discusses the significance of the gravity and magnetic features in the filtered anomaly maps and how they relate to the seismic structure and the field mapping.  5.2.1  Relationship of Density to Velocity Structures  The north end of Line 1 has higher densities down to a depth of 8 km and a low a zone down to 5 km (Figure 5.8a, b; see original Figures 3.8, 3.49a, and 4.1 for more detail). This zone is coincident with the metasedimentary migmatites of the Quetico belt and velocities that range from 5.9 km/s at the surface to 6.6 km/s at 25 km depth. This unit is bounded at 100 km by the Lepage fault which has a low velocity zone and a small gravity low associated with it. The southern end of the line (160 to 260 km) has lower densities (down to 25 km) and a low cr (down to 20 km) under the plutons and orthogneiss of the Wawa terrane. The surface velocity is 6.1 km/s and the bottom of the anomalous features coincides with P-wave velocities of 6.7 to 6.8 km/s. Further to the south, the gravity values increase over the metavolcanics of the Michipicoten belt. Lower densities were required on Line 2 (Figure 5.8c, d; see original Figures 3.14, 3.49b, and 4.2 for more detail) out to 100 km model distance. This region coincided with the orthogneisses of the Abitibi belt and the GRB. The entire KSZ (40 to 170 km) has a higher a down to 12-13 km. This corresponds to a zone of uplifted mafic material and is marked by an increase in P-wave velocities (from 6.0 km/s to 6.4 km/s at 5 km depth). Just south of the C B , higher densities were required in an area of mixed gneissic composition in the Wawa terrane where velocities had decreased again.  Chapter 5. Discussion and Interpretation  N 0  DISTANCE (KM) 40  80  120  160  200  240  213  DISTANCE (KM)  S 280  320  360.0  40  120 I  160 l  200 l  240 l  280 l  320 I  Figure 5.8. P-wave velocity models with variations from Nafe-Drake velocity/density relationship for Lines (a) 1, (c) 2, and (e) 3. Stars represent shot locations and arrow heads mark intersections with crossing lines. Contours are in km/s and 'plus sign' and 'minus sign' shading represents zones with densities greater than and less than the polynomial fit to the Nafe-Drake set of data points. Variations in Poisson's ratio for Lines (b) 1, (d) 2, and (f) 3.  360  Chapter 5. Discussion and Interpretation  214  Line 3 has a lower density region out to 50 km and down to 14 km depth that is coincident with a 4-5 km deep low a zone (Figure 5.8e, f; see original Figures 3.21, 3.49c, and 4.3 for more detail). This part of the profile covers the plutons and orthogneiss of the northern Abitibi belt. Surface velocities are in the 5.9 to 6.1 km/s range. There is a high density region from 70 to 150 km model distance and down to 10 km depth which coincides with the metavolcanics and underlying low velocity zone. Another high density zone extends from 180 to 270 km and down to 10 km depth. This appears to correspond to a low <7 zone of 10 to 15 km thickness which runs from 160 to 280 km. This also underlies the metavolcanics but may be influenced by a nearby outcrop of orthogneiss. There is a broad scale lower density region at the southern end of this line, possibly indicating a change in the dominant upper crustal composition for the Abitibi belt. The western end of Line 4 has a large low a zone extending to 20 km depth and ranging from -40 km to 100 km at the surface (Figure 5.9a, b; see original Figures 3.28, 3.49d, and 4.4 for more detail). This marks the extent of the metasedimentary migmatites of the Quetico belt. It is also characterised by velocities that range from 5.9 km/s at the surface to 6.4 km/s at 20 km depth. There is an increase in densities at depth under the V R B which is coincident with the step-up in P-wave velocity (6.1-6.2 km/s to 6.3-6.4 km/s) at distances of 100 to 140 km and down to a depth of 10 km. The 5 to 20 km depth range in this region is marked by a high <r indicating that this buried anomalous body may be associated with the uplifted mafic rocks of the C B . Lower densities are required under the GRB and the orthogneisses of the adjoining Abitibi belt between 180 and 280 km. A higher density region coincides with the metavolcanics of the Abitibi at the eastern end of the model. Line 5, at the western end, shows higher densities under the metavolcanics of the Michipicoten belt (Figure 5.9c, d; see original Figures 3.38, 3.51a, and 4.5 for more detail). There are lower densities in the tonalitic gneisses of the Wawa terrane between 20 and 100 km where lower velocities were also imaged. Higher densities are required in  Figure 5.9. P-wave velocity models with variations from Nafe-Drake velocity/density relationship for Lines (a) 4 and (c) 5. Stars represent shot locations and arrow heads mark intersections with crossing lines. Contours are in km/s and 'plus sign' and 'minus sign' shading represents zones with densities greater than and less than the polynomial fit to the Nafe-Drake set of data points. Variations in Poisson's ratio for Lines (b) 4 and (d) 5.  ^  Chapter 5. Discussion and Interpretation  216  the anomalous high velocity zone (6.6-6.7 km/s) under the C B . To the east of the ILCZ at 170 km, lower and higher densities coincide with the orthogneiss and metavolcanics of the Abitibi, respectively.  A low cr zone is tentatively suggested for the east side of the  ILCZ, against which the thrust fault is butted. This zone may relate to the depth extent of the orthogneiss outcrop on the Abitibi side.  5.2.2  G r a v i t y and Magnetic Interpretation  In this section, the potential field characteristics of each block in the KSZ region are discussed. Then the strong trends, lineations, and fault lineaments are discussed for the region. The Chapleau block is characterised by a strong gravity high, especially over the mafic gneisses and anorthosites in the N E part of the block (Figures 4.11 and 4.14). The ILCZ, which demarcates the eastern margin of the C B , has a small broadscale aeromagnetic response but there are larger magnetic anomalies over the carbonatite complexes at 82.9° longitude and 48.3° latitude, and 83.2° longitude and 47.8° latitude (Figure 4.22). -3  The gneisses of the C B have a moderate susceptibility range (0.3 to 18.8 x 10 ) and - 3  high densities (2840 to 3190 kg m ) [Skive and Fountain, 1988; Fountain et al., 1989]; thus, the observed potential anomalies are as expected. The high densities corroborate the existence of the near surface high velocity anomaly (6.5-6.6 km/s) in the seismic interpretation. No susceptibility values are available for the carbonatites. The G R B does not have a gravity anomaly (Figures 4.11 and  4.14) despite the  fact that there are high densities at the surface [Nkwate and Salisbury, 1988RRR From structural mapping [Leclair, 1988], the G R B appears to be a perched thrust tip and so the high density material forms only a thin layer at the surface. There are no significant velocity anomalies associated with this zone, so the missing gravity anomaly is reasonable. However, there is a very strong aeromagnetic response over this block (Figure 4.22) which is bounded by the ILCZ to the east, the Saganash Lake fault to the west, and the Wakisimi  Chapter 5. Discussion and Interpretation  217  River fault to the south. There are no susceptibility measurements from this area, but from analysis of the gneisses in the C B [Shive and Fountain, 1988], magnetite is the probable cause of the magnetic response. The F M B has both a strong gravity and aeromagnetic signature (Figures 4.11 and 4.22) but was not covered by the seismic experiment. The ILCZ and Lepage fault have notable responses (Figure 4.12) as they converge within the Quetico and Opatica gneiss belts. The gravity high extends over a region to the south and south west of the F M B (Figure 4.15) where there are exposed granulites. The Rerni Lake zone has a gravity high that joins those of the VR.B and the F M B (Figure 4.11). There is also a suggestion of a magnetic anomaly at depth under this zone (Figure 4.29). The V R B is characterised by an extension of the C B gravity high in an arcuate form over this region (Figure 4.14).  There are no high density materials apparent at the  surface; thus a buried source is suggested and is in agreement with the higher velocities imaged at depth under the region.  These results, together with the high a, strongly  suggests a mafic body, similar to that of the C B . There are small aeromagnetic signatures associated with the bordering Lepage and Saganash Lake faults. To the south-west of this block, there is a noticeable low in both the gravity and aeromagnetic responses over the Wawa domal gneisses (Figures 4.15 and 4.27).  It has also been imaged as a low  velocity region at the surface along Line 5. This terrane comprises 63% tonalitic gneiss and it is this material that has a low density, a low velocity, and a low susceptibility. There are strong N E - S W trends in both the gravity and magnetic data sets (Figures 4.12 and 4.24). The main feature of note is the continuity of the anomalies along the full length of the KSZ. The ILCZ is outlined and matches the location of the mapped fault. This lineament is also marked by a strong lateral velocity contrast (6.6-6.3 km/s) along the eastern side of the C B . One important result is the delineation of an extension of the Lepage fault down through the Wawa terrane.  This fault also causes a lateral  discontinuity in the seismic structure, especially on Line 4. The Lepage fault may have  Chapter 5. Discussion and Interpretation  218  important implications in understanding the role of the syn- or post-thrust normal faulting. The east-west trends of the potential fields in the Superior province were successfully accentuated by directional filtering (Figures 4.13 and 4.25).  The results confirm the  mapped trends and also help delineate some buried or more subtle features. The southern edge of the Quetico belt (49.0° to 49.2°) is highlighted by both a gravity and magnetic response.  The gravity high continues north into the belt over the metasedimentary  migmatite. This belt is characterised by material with a low to moderate velocity, a low er, and a high density. There is a gravity low between 49.8° and 50.2° that coincides with an extension of the Wabigoon greenstone belt, which outcrops to the west of the map region, under the cover of the Moose River basin (MRB). There is also a magnetic high which may be associated with the edge of the M R B or the underlying greenstone belt. At the northern edge of this belt, there is an increase in the gravity field that is most likely associated with a continuation of the English River gneiss belt under the M R B cover. Within the Michipicoten and Abitibi greenstone belts, there are some potential field variations that can be attributed to the distribution of metavolcanics. On the west side of the KSZ, there are gravity highs related to the metavolcanics at 47.0° and 48.2°. These belts continue on the east side but seem shifted to the north by 50 to 60 km. This may indicate some sinistral transcurrent motion on the ILCZ. These metavolcanic belts also have high velocities and moderate susceptibility. In the south, there is a magnetic high at 47.0° marking the boundary of the Superior and Southern provinces. A gravity high, immediately to the south of this boundary, delineates the Penokean belt of the Southern province.  Chapter 5. Discussion and Interpretation  5.3  219  Palaeogeotherm and Palaeostrength Curves  The seismic refraction modeling provided estimates of the thicknesses of the main lithological units exposed in the Abitibi, Wawa, and Kapuskasing terranes. Rock heat production and thermal conductivity values have been obtained for the exposed rocks [Ashwal et al., 1987] and so there is control over the thermal character of the upper 20 km of the crust. Surface heat flow values [Cermak and Jessop, 1971] for several locations in the region were then used to obtain estimates on the lower crust heat production and the reduced mantle heat flow. These thermal profiles were used to assess the possibility of magnetisation values at depth and, with layer and composition estimates from the seismic study, to estimate rheological conditions at depth around the KSZ. The thermal values and a knowledge of the pre-uplift crustal structure were used to recreate the rheological conditions in the stabilisation period at 2700 to 2600 Ma and immediately prior to the Kapuskasing uplift at approximately 2000 Ma. Thus, estimates of the mineralogies that controlled the faulting and information on possible mid-crustal detachment zones could be obtained.  5.3.1  Present D a y T h e r m a l Structure and Strength  Using a knowledge of the subsurface structure obtained from the seismic results and measurements of heat flow [Cermak and Jessop, 1971] and rock heat production [Ashwal et al., 1987] from a traverse across the KSZ, static earth geotherms [Chapman, 1986] for the present day crust were calculated.  Each layer in these models was assigned a  name that corresponded to the terrane of that dominant composition.  Thus, the top  layer of metavolcanics and metasediments was called the Abitibi layer; the next layer of orthogneiss and plutons was called the Wawa layer; the mafic gneisses and anorthosites were called the Kapuskasing layer; the lower crustal layer corresponded to a highly mafic and garnetiferous unit; and finally the mantle was an ultramafic unit. surface heat flows of 43 and 33 mW m  - 2  The observed  were well matched by the models for the  220  Chapter 5. Discussion and Interpretation  Figure 5.10. (a) Geothermal profiles for the Wawa (A), Abitibi (B), and Kapuskasing (C) terranes. The surface heat flows in mW m" are shown in the enclosed box. (b) Strength envelope (stippled region) for the Abitibi terrane. The straight line curve is Byerlee's law and the other curves represent power law creep for wet quartzite (1), dry quartzite (2), and anorthosite (3). 2  Chapter 5. Discussion and Interpretation  221  Abitibi and Kapuskasing terranes, respectively (Figure 5.10a). The measured value of 50 2  mW m ~ for the Wawa terrane may have been unduly influenced by high heat production granites that occur within that terrane [Ashwal et al., 1987]. The data were satisfied by using a reduced mantle heat flow of 20 mW m  - 2  . This value is in agreement with values  2  (~ 20 mW m~ ) estimated by Cermak and Jessop [1971] and Ashwal et al. [1987] for - 2  the same region but is less than the value (27 mW m ) calculated for shield averages by Morgan [1985]. The discrepancy can be accounted for by his underestimate of mid and lower crust heat production. The thermal profiles are typical for shield regions [c.f., Chapman, 1986; Morgan, 1985] and give Moho temperatures in the 4 0 0 - 4 7 5 ° C range. This value is less than the Curie temperature for magnetite and thus it would be possible to have magnetic sources in the lower crust. Using the thermal profiles and a knowledge of layer thicknesses from the seismic analysis, strength envelopes were constructed for the crust (Figure 5.10b). These envelopes were used to estimate the levels of shear stress sustainable in a given crustal column with assumed thermal conditions and mineralogical composition. The envelopes were limited by a brittle failure curve [Byerlee, 1968; Sibson, 1974] and a ductile flow curve [Ord and Hobbs, 1989] and hence the mode of deformation at a given depth could be inferred. The maximum shear stress considered possible in the crust was 300 M P a [Ord and Hobbs, 1989; C. Beaumont, pers. comm., 1989]. Thus, Byerlee's law was used up to that stress level and below the corresponding depth, the response of the crust was considered elastic (i.e., no change of shear stress with depth). A typical geological strain rate of 1 0  -14  s  _1  [Kusznir and Park, 1984] was used for all the model calculations. The ductile flow curves were constructed using an empirical relationship for power law creep [Ord and Hobbs, 1989; White and Breton, 1985] for each mineral. Each lithological layer was deemed to be Theologically controlled by one dominant mineral [Kusznir and Park, 1986; Ord and Hobbs, 1989; Braun and Beaumont, 1989]. Within that layer, the intersection of the flow curve for that mineralogy and the brittle curve marked the onset of a stress drop,  Chapter 5. Discussion and Interpretation  222  and therefore the top of a potential detachment zone that extended down to the shear stress minimum. The strength envelopes calculated for each terrane show no potential detachment zones or great strength contrasts and are in keeping with the stability of shield regions (Figure 5.10). The top 12 to 13 km are in a brittle field and below this the crust behaves elastically.  Seismicity studies have revealed activity in the northern KSZ region with  events up to a magnitude of 5 [Forsyth et al., 1983]. However, no estimates of focal depths for the events were obtained, and so no comparison could be made with the brittle field estimate. The strength envelope for the Abitibi belt is shown in Figure 5.10b and displays the typical characteristics discussed above. Thus, any earthquakes in this region should occur in the upper 12 to 13 km. Recent preliminary laboratory work on the strength of the Kapuskasing gneisses [Wi/£s and Carter, 1989] suggests that these rocks behave elastically up to shear stresses of 250 M P a at a confining pressure of 1000 MPa. Thus, our maximum stress assumptions for the crustal column are of the correct order.  5.3.2  Post-Stabilisation T h e r m a l Structure and Strength  The Kapuskasing uplift disturbed a formerly continuous greenstone belt (Michipicoten and Abitibi) [Percival and Card, 1985] and so the structure of the present day Abitibi belt, to the east of the faulted terrane, should be a good representation of the pre-uplift crustal structure. An average one dimensional structure derived from the seismic model of Line 3 and a vertical reconstruction of the oblique exposure [Percival and Card, 1985] were used in the thermal and rheological modeling. The one dimensional thermal model for the Abitibi belt, along with two other related models (see Figure 5.11), were used to study palaeogeotherms and palaeostrength curves for the period after the stabilisation of the southern Superior province at 2700 to 2600 Ma. Values for present day and Archean mantle heat flow are highly variable [Morgan, 1985;  Chapter 5. Discussion and Interpretation  223  Figure 5.11. Thermal profiles and strength envelopes for period immediately after stabilisation of Superior craton. (a) Geothermal profiles for the models in frames b (C), c (B), and d (A). The surface heat flows in mW m are shown in the enclosed box. Strength envelopes (stippled region) for the Abitibi (d) and two other models (b and c). The straight line curve is Byerlee's law and the other curves represent power law creep for wet quartzite (1), dry quartzite (2), anorthosite (3), diopsidite (4), wet dunite (5), and dry dunite (6). - 2  Chapter 5. Discussion and Interpretation  224  Condie, 1984; Davies, 1979; Watson, 1978]. For our purposes, a conservative mantle 2  heat flow of 25 mW m ~ was used for the post-stabilisation period.  This represents  an increase of 25% over the present mantle heat flow as calculated above.  The heat  production of the various lithologica! units was estimated at 2.25 times their present production. This value was obtained from consideration of the present abundances of the radiogenic elements and knowledge of their half-lives [Windley, 1986]; thus, an estimate of their former abundance could be predicted for a closed system. For this period, an extra layer (the K crust) was added as a pre-Kenoran uplift and erosion surface unit [Percival et al.,  1988]. It was assigned heat production and conductivity values the  same as the underlying metasedimentary and metavolcanic unit (Abitibi layer). The Abitibi belt model is shown in Figure 5.lid and two extreme variations on this model are shown in Figures 5.11b and c. In these alternate models, the thicknesses of the Wawa and Kapuskasing layers are varied by 5 km to account for maximum possible variations in these units. Figure 5.11a shows the thermal profiles and surface heat flows for each of these models. The surface heat flow values are reasonable for a young active terrane [Kusznir and Park, 1986]. Only curve A satisfies the granulite formation pressures (600 to 800 MPa) and temperatures (700°C to 800° C) as estimated by Percival [1983] for the Kapuskasing gneisses. This curve corresponds to the Abitibi belt model of Figure 5.lid. However, by increasing the mantle heat flow by a further 50% , the granulite formation conditions can be matched for curves B and C. Figures 5.11b, c, and d show the strength envelopes for each of the crustal models. The Abitibi and Wawa lithological layers are Theologically controlled by wet and dry quartz deformation, respectively.  The Kapuskasing layer, although it is classified as a  mafic gneiss unit, is controlled by the high content of tonalitic gneiss [J.A. Percival, pers. comm., 1989] which is dominated by dry quartz.  The lower crustal layer and  the upper mantle are dominated by pyroxene and olivine deformation, respectively. In Figure 5.11b, there is a zone of weakness between the Abitibi and Wawa layers, with all  Chapter 5. Discussion and Interpretation  225  layers below the Wawa behaving in a ductile manner. However, within the Kapuskasing layer, only the bottom half is in a ductile state for pyroxene (curve 4) and this is not supported by the mapped field evidence [Percival and Card, 1985; Moser, 1988; Bursnall, 1988]. The model in Figure 5.11c also shows a zone of weakness between the Abitibi and Wawa layers with all of the thin Kapuskasing layer now in a ductile state for pyroxene. Finally, the Abitibi model in Figure 5.lid shows the same zone of weakness as above. The Wawa layer is now completely in a dry quartz ductile state and the Kapuskasing layer is in a ductile state for pyroxene. This is corroborated by the ductile deformation structures observed in the field for the two terranes [Percival and Card, 1985; Moser, 1988]. Also, in this model, there is very little strength contrast across the Moho whereas the contrast is greater in the previous two models. There is no evidence that the Moho was a decollememt plane during this period or subsequent tectonic events and so the model of Figure 5.lid seems the most appropriate on consideration of this evidence and the above temperature and strength factors.  5.3.3  Pre-Uplift T h e r m a l Structure and Strength  Geothermal regimes and geostrength curves were constructed for a period immediately prior to the Kapuskasing uplift at circa 2000 Ma. The same three models as the above section were used but with the surface layer (K crust) having been eroded away after the Kenoran orogeny (2650 Ma). Mantle heat flows of 25 and 20 mW m  - 2  were used: the  first value assumed that the mantle had not redirected heat to the surrounding oceanic crust during the cratonisation of the Superior province; the second value assumed that the mantle heat flow had settled to the values currently under the shield. The heat production values utilised were 2.0 times those of the present day [Windley, 1986]. These models were constructed to evaluate possible zones of detachment within the crustal column that could have initiated the faulting of the Kapuskasing uplift. If the approximate layer thicknesses and heat productions were the only constraints available for this analysis,  Chapter 5. Discussion and Interpretation  226  Figure 5.12. Thermal profiles and strength envelopes for period immediately prior to the Kapuskasing uplift, for a reduced mantle heat flow of 25 mW m . (a) Geothermal profiles for the models in frames b (C), c (B), and d (A). The surface heat flows in mW m are shown in the enclosed box. Strength envelopes (stippled region) for the Abitibi (d) and two other models (b and c). The straight line curve is Byerlee's law and the other curves represent power law creep for wet quartzite (1), dry quartzite (2), anorthosite (3), diopsidite (4), wet dunite (5), and dry dunite (6). - 2  - 2  Chapter 5. Discussion and Interpretation  then many models could be constructed to satisfy the data.  227  However, an additional  constraint available was that temperatures at the base of the uplifted portion of crust were 300°C to 350°C immediately prior to the uplift. This is a value given by analysis of blocking temperatures for Rb/Sr and K / A r isotope systems in rocks close to the thrust fault [Percival et al., 1988]. Thus, the palaeostrength models were constrained to a fairly small set of plausible answers. Figure 5.12 shows the palaeogeotherm and palaeostrength curves for a mantle heat flow of 25 mW m  - 2  . The models in Figures 5.12b and c show the onset of a zone of  weakness in the 18 to 19 km depth range. However, from examination of their thermal profiles, the ductile deformation starts at 330°C to 340°C and goes to 420°C. Also, the model in Figure 5.12c has too small a stress difference across that zone.  There is a  small zone of weakness at the crust-mantle interface which could have detached but a preliminary analysis of the seismic reflection data shows no prominent shear zones near the base of the crust [W^esi et al., 1988]. The model in Figure 5.12d has a zone of weakness at the top of the Kapuskasing layer and would not have uplifted enough of that unit to account for the surface exposures. The temperatures of these gneisses would also be too high and the large stress drop at the Moho would have resulted in a different faulting pattern. Figure 5.13 shows the palaeogeotherm and palaeostrength curves for a mantle heat - 2  flow of 20 mW m . The model in Figure 5.13b has a zone of weakness at approximately the right depth and temperatures in the 300°C to 350°C range, but the stress difference across the zone is very small. The model of Figure 5.13c shows no zones of weakness and so cannot be considered. Finally, the model in Figure 5.13d (the Abitibi belt model) shows only one zone of weakness starting at 18 km depth and with a temperature in the 330°C to 360°C range. Also, there is a large stress difference of 250 M P a across the zone, making it a likely candidate for faulting in a tectonic stress field. Thus, from the above results, the models in Figures 5.12b and 5.13d are the most  Chapter 5. Discussion and Interpretation  0.0  10.0  20.0 30.0 D e p t h (Km)  40.0  228  50.0  c o 0  S h e o r S t r e s s (MPa) 150 300 ISO 600 750 900 1050 1200  d 0  0  S h e a r S t r e s s (MPo) 150 300 <50 600 750 900 1050 1200  Figure 5.13. Thermal profiles and strength envelopes for period immediately prior to the Kapuskasing uplift for a reduced mantle heat flow of 20 mW m~ . (a) Geothermal profiles for the models in frames b (C), c (B), and d (A). The surface heat flows in mW m ~ are shown in the enclosed box. Strength envelopes (stippled region) for the Abitibi (d) and two other models (b and c). The straight line curve is Byerlee's law and the other curves represent power law creep for wet quartzite (1), dry quartzite (2), anorthosite (3), diopsidite (4), wet dunite (5), and dry dunite (6). 2  2  Chapter 5. Discussion and Interpretation  229  promising structures. Additional calculations were made and it was found that a model with a Wawa layer of 8 to 9 km thickness, a Kapuskasing layer of 11 to 12 km, and mantle heat flow of 22 to 23 mW m  - 2  was optimal. However, if realistic two dimensional  effects and variations in strain rate are considered, the resolution in the models is not adequate to differentiate between these final models.  The important results from this  analysis are that the main thrust faulting was controlled by a dry quartz rheology. This is corroborated by the widespread presence of tonalitic leucosomes within the Kapuskasing gneisses [Percival and Card, 1985]. Thus, the original detachment zone is limited to the 17 to 20 km depth range by the results of the refraction analysis, a preliminary analysis of the reflection data [West et al., 1988], and the above rheological modeling.  5.4  Tectonic Development  The formation and subsequent tectonic activation of the Archeah crust in the central Superior province will be discussed in this section.  Details of these developments can  be found in Goodwin [1981], Percival and Card [1985], and Hoffman [1988], and so only those aspects related to the results of the geophysical modeling will be presented here.  Formation of Crust The major crust forming event at 2750-2650 M a was followed by metamorphism, intracrustal melting, and plutonism [Percival and Card, 1985]. Figure 5.14 shows a schematic representation of the formation of the Archean crust after Percival and Card [1985]. There was a volcanic and sedimentary succession laid prior to 2765 Ma. This was followed by an eruption of volcanics and deposition of sediments in the 2750 to 2696 Ma range. During and after the end of this activity, with a minimum age of 2707 Ma, large volumes of tonalite were injected into the crust. These were either juvenile magmas from the mantle or the result of partial melting of the lower crustal successions.  The tonalites engulfed  and detached fragments from the overlying greenstone succession and the underlying  Chapter 5. Discussion and Interpretation  230  • V V V V V V V V V V W V ' V V V V V V V V V V V V V W V V V V V V V V V V V V V V V ' ,'y/yyy.vV  vvyvyvvvvvvvv vvvvvvv; ,Early volcanics and sediments -Shawmere anorthosite  > 2.765 m.y.  10 km  (a)  Abitibi. Michipicoten "sequences  • 20km  (b)  Figure 5.14. Development of crust in southern Superior province in 2765 to 2650 M a period after Percival and Card [1985]. (a) Early volcanic and sedimentary piles forming thin sialic crust, (b) Later addition of volcanic extrusions and erosional sediments, (c) Intra-crustal intrusion of tonalites. /•• volcanic and sedimentary sequences. In the seismic refraction models and potential field anomalies, these tonalitic gneisses have the distinctive characteristics of low velocity, low density, and low susceptibility compared to both the overlying greenstones and underlying Kapuskasing sequences. The intraplating mechanism of vertical magmatic accretion by which these tonalites were emplaced would have resulted in varying thicknesses depending on local stress fields and thermal structure. Thus, the variability of this layer as imaged seismically in the greenstone belts of the Canadian shield can be explained.  Assembly and Deformation of Superior Province A major deformation episode during the 2696 to 2680 M a period is recorded in the surface volcanic/sedimentary sequences (Abitibi) and the tonalitic layers (Wawa). The deformation episode may have continued into the 2650 to 2627 Ma range in the old  Chapter 5. Discussion and Interpretation  volcanic/sedimentary sequences (Kapuskasing).  231  This would mean lower crustal defor-  mation occurring 25 to 50 Ma after the intrusion of post-tectonic plutons in the Abitibi and Wawa layers. This high grade metamorphism would require a major burial episode after the tectonism and magmatism of the overlying crust. There is no apparent cause for this event but it is possible that there was slower cooling at depth and thus some magmatism associated with a hotter lower crust [Percival and Card, 1985]. However, a more reasonable interpretation is to link this lower crustal metamorphism with a proposed transcurrent motion on the KSZ at that time [Watson, 1980]. There are northeast lineations in the tonalitic gneisses of the Kapuskasing zone that are the result of a late ductile tectonism related to a sinistral transcurrent motion along the KSZ during the emplacement of the Matachewan dyke swarm (2633 Ma). This event could have reset the zircon dates to the 2650 to 2627 Ma range. This sinistral motion is supported by the 50 km northward displacement of the east-west trending potential field anomalies on the east side of the KSZ (Figures 4.13 and 4.25). This shift would be difficult to relate to the.uplift and so is best explained by this sinistral motion. Goodings [1988] has also inferred a sinistral motion on the KSZ from detailed studies of magnetic lineaments in this region and further to the north. The east-west striping, that shows the above transcurrent motion, relates strongly to the greenstone and gneiss terrane patterns of the southern Superior province. The gravity highs over the greenstone belts are accompanied by lows over the gneiss belts and orthogneiss of the greenstone belts. The magnetic lineaments demarcate the boundaries of the greenstone/gneiss belts, even when buried under Phanerozoic cover. This potential field anomaly support of the mapped geological outcrops and buried boundaries is 5  supportive of the microcontinental accretion theories of Windley [1986]. Goodwin [1981] and Hoffman [1988] have also suggested that this type of accretion may have occurred in the southern Superior province. Hoffman [1988] likened the accretion of the Proterozoic 5  See Chapter 1  Chapter 5. Discussion and Interpretation  232  continent, Laurentia, to that of Eurasia in the Phanerozoic period. This ascribes a mobilistic interpretation to the Hudsonian orogen (circa 1800 Ma) and favours horizontal tectonics as a dominating mechanism for protocontinental accretion in the late Archean and early Proterozoic.  Kapuskasing Uplift  Figure 5.15. Collision of Churchill and Superior cratons in 1900 to 1800 M a period after Gibb [1983]. The stippled and white regions are the Churchill and Superior cratons, respectively. The following are the annotated features: K is the Kapuskasing structure; KS is the Kenyon structure; W R F is the Winisk River fault; TS is the Thelon salient; US is the Ungava salient; and LS is the Labrador sea.  It is still unclear how the Kapuskasing uplift relates to activities around the Superior craton in the 2200 to 1800 Ma range [Percival et al., 1988]. Fom analysis of isotopic systems in the KSZ region, the uplift is constrained to the early Proterozoic [Percival et al, 1988]: in the 2250-2200 Ma range from biotite/argon data or in the 1950-1900 Ma range from Rb/Sr data.  Some doubt exists about the first date because of the  possible presence of excess argon in granulite terranes [Poland, 1979]. The second date is easier to relate to other tectonic events of the period.  In particular, Gibb [1983]  related the Kapuskasing suture to the collision of the Superior and Churchill cratons in Hudsonian time (circa 1800 Ma). Figure 5.15 shows part of the collision development,  Chapter 5. Discussion and Interpretation  233  where the Ungava salient (US) advanced further than the Thompson salient (TS) and resulted in the Kapuskasing thrusting event. However, the movements on the Winisk River and Kenyon faults appears to have the wrong sense of motion. The Winisk River fault is considered a dextral slip motion from analysis of outcrops and potential field anomalies [Gibb, 1975]. However, the evidence for dextral motion on the Kenyon structure is equivocal and possibly sinistral motion on that fault accommodated the difference in convergence between the two salients of the Superior province. Hoffman [1988] outlined the various events that were associated with the trans-Hudson orogen.  Deformation  events in the Fox River belt ( 1 8 8 3 ± 2 Ma) and the Thompson belt (1880 and 1770 Ma) on the western edge of the Superior craton may have been contemporaneous with the Kapuskasing uplift. A series of carbonatite complexes pin the post-thrust normal faults and a minimum possible age for the uplift of 1750 M a is given by the youngest of those intrusions.  However, from analysis of several alkalic intrusions, an age of 1907 Ma is  inferred for the uplift [Percival and McGrath, 1986]. For events to be tectonically linked over large distances, the lithosphere must act as a stress guide and stresses must be large enough to initiate a major deformation event. Stress levels of 100-150 MPa have been shown, from thermo-rheological modeling [Kusznir and Park, 1986], to be large enough to initiate crustal-scale deformation. Braun and Beaumont [1989] have used 2d finite element modeling to show how offset features can be created in a tectonically induced stress regime. Present day stress differences of 100150 MPa have been inferred for a large thrust structure in Central Australia [Lambeck, 1983]. Thus, in the past, stress levels must have been equal to or greater than this value to produce the structure. Considering the timing of events at the edges of the craton, this intracratonic thrust may have been the result of the lithosphere acting as a stress guide to forces operating at the craton margins, 1000 km away [Goleby et al., 1989]. This implies that stresses of 100-150 MPa or more can be transmitted over large distances and hence it is plausible to link events at the western edge of the Superior province (1200 km  234  Chapter 5. Discussion and Interpretation  away) to the Kapuskasing uplift. England and Bickle [1984] have shown from thermal modeling and the analysis of metamorphic terranes that stresses in the Archean crust were comparable to today's tectonic stresses, i.e., 80 to 160 MPa. Thus, both the stress magnitudes and the transfer of those stresses were viable for linking circum-Superior events to the Kapuskasing uplift in the early Proterozoic.  Two dimensional thermal  modeling of the KSZ has set a lower limit of 20 my for the duration of this uplift event [Mareschal and LeQuentrec, 1988]. The seismic and geostrength evidence show that the thrust fault soles at depths between 17 and 21 km. The model for the pre-uplift crust suggests a 5 to 7 km mountain formed by the thrusting.  This resulted in a root forming to isostatically  compensate  this extra mass at the surface. However, thermal modeling indicates that rapid erosion could have kept pace with the differential uplift [LeQuentrec et al., 1989]. Nevertheless, with this fast erosion or after erosion of the mountain, a root is still neccessary due to the presence of high density material in the upper crust under and to the west of the KSZ. The seismic character of this root indicates a mafic origin for the material. Thus, a ductile flow of materials may have occurred at or near the base of the crust as part of the isostatic compensation. This root is probably part of some general underplating [Furlong and Fountain, 1986] that occurred under the shield after the uplift event.  The lower  crustal velocities of the refraction models (7.2 to 7.7 km/s) are in good agreement with those estimated for underplated crust [Furlong and Fountain, 1986]. These underplating episodes might be linked to the exposed carbonatite complexes which are indicative of magmatic events at mantle depths. The magmatic events of underplating and carbonatite intrusions may have been initiated by delamination of the lower lithosphere during the upper crust thrusting episode. This would result in A-type subduction [Kroner, 1981] as described in Chapter 1. The geostrength analysis of section 5.3 shows that the normal faulting that occurred subsequent to the thrust may have been resettling and slumping episodes that faulted along zones of weakness already existing within the crust. The  Chapter 5. Discussion and Interpretation  235  Lepage fault may be one example where a zone of weakness had been created between the Abitibi and Wawa layers in the stabilisation period (Figure 5.1 Id) and was activated subsequently at the time of uplift.  2d finite element geodynamic modeling may help  resolve possible deformation patterns for the lower crust resulting from this uplift.  Chapter 6  Conclusions  A geophysical model has been constructed for the Archean crust that shows a Proterozoic uplift feature rooted in the mid crust. Mid-crustal velocities of 6.5 to 6.6 km/s are brought from 15-25 km depth under the surrounding greenstone belts to 2-5 km under the southern block of the KSZ (Figure 6.1a).  This material is mid to lower Archean  crust uplifted along the Ivanhoe Lake Cataclastic Zone during an intraplate compressive episode in the early Proterozoic. The results indicate that the Chapleau Block and the Val Rita block are remanents of this event that have been delimited by subsequent normal faulting. This general trend is in good agreement with the results of the seismic data analysis of Wu and Mereu [1988]. The surface exposure of this granulite feature obliquely intersects the striped greenstone/gneiss accretion pattern of the Superior province. This east-west striping pattern of alternate greenstone and gneiss belts has been corroborated by the potential field analysis and these belts have been shown to have distinct potential field and seismic characteristics.  This striping pattern is supportive of the microconti-  nental accretion theories of Windley [1986] and Hoffman [1988] for cratonisation in the late Archean period. A crustal model for greenstone belts has been constructed and helps to explain many features observed here and elsewhere in the shield. The feature of medium velocity (6.16.3 km/s) surface greenstones underlain by a variable thickness of low velocity (5.9-6.1 km/s) tonalitic gneisses is supported by both laboratory and field measured velocities. The density and susceptibility of this greentone layer is variable but the tonalites have low densities and susceptibilities. These units are underlain by a mixed layer of mafic gneisses, anorthosites and paragneisses that have higher velocities (6.5-6.9 km/s). The high quartz 236  LONGITUDE  LONGITUDE  Figure 6.1. (a) Map of depth to material with a velocity of 6.6 km/s. Depths are in km. (b) Map of crustal thickness in km. The stippled region is the KSZ and the stars are the shot locations.  to to -4  Chapter 6. Conclusions  238  content of this layer resulted in the development of a zone of weakness at 2000 Ma which was exploited by the thrusting episode. Below this layer, there is a higher concentration of pyroxenites and garnetiferous rocks as inferred from the imaged velocities (>7.0 km/s). This has resulted in a lower crust that is both mafic and strong.  The crust thickens  significantly from 40-43 to 50-53 km under and to the west of the Chapleau Block at the southern end of the Kapuskasing structure (Figure 6.1b); a feature that is also imaged by a time-term analysis of the upper mantle refraction arrivals by Northey and West [1986]. The thrusting of one section of crust over the other may have produced the thickening or it may have developed by continental underplating [Furlong and Fountain, 1986] to form an isostatic root and compensate the presence of the anomalous high velocity and high density material in the upper crust. The proposed geological model of an intra-cratonic uplift is supported by the results of the seismic data and gravity analysis (Figure 6.2 and 6.3).  Figure 6.2 shows the  geophysical and derived geological models for Line 4 which traverses the Groundhog River block of the KSZ. To the west, the metasedimentary migmatites of the Quetico gneiss belt (see Figure 1.1) have a low cr and velocity structure varying from 5.9 km/s at the surface to 6.5 km/s at the base. This unit shows no velocity inversions and thus seems to have a different seismic structure than the greenstone belts.  The increase in  velocities at depth under the Val Rita block (130 to 160 km) is accompanied by a high cr zone and a positive gravity anomaly. This indicates that there are buried mafic rocks under the block associated with the uplift. Detailed geological mapping [Leclair, 1988] and gravity analysis [Nhwate and Salisbury, 1988] support the interpretation of Percival and McGrath [1986] that the GRB (162-172 km) is a perched thrust tip isolated by postthrust normal faulting. There is no significant velocity variation as the G R B is crossed. To the east, there are lower and higher densites in the orthogneiss and metavolcanics, respectively, of the northern Abitibi belt. At depth, the granulite/tonalite layer is defined by velocities of 6.6 to 7.0 km/s. This is underlain by a mafic pyroxene lower crust with  Chapter 6. Conclusions  239  w -20  DISTANCE (KM) 30  80  130  180  *«j p > N - D relationship  III  (7 > 0.2S  ->j p < N - D relationship  H§j  <t < 0.25  Groundhog River Vol Rita  Quetico  DISTANCE 130  180  280  230  330  380  330  380  Abitibi  (KM) 230  280  ....It  pyroxene granulite / tonalite  ^ A ^ | metavolcanics | | • j metasedimentary j  migmatite  j orthogneiss / plutons  pyroxene granulite f§3§8g garnet granulite  Figure 6.2. (a) Geophysical model for Line 4. The contours are P-wave velocity in km/s. Superimposed on this are: variations in density (p) away from the polynomial fit to the Nafe-Drake set of points; deviations of Poisson's ratio (a) from the starting value of 0.25. (b) Geological model derived from the above variations, surface exposures and fault mapping. The heavy solid line is the thrust fault that soles into the detachment zone. The dashed lines are normal faults. Vertical exaggeration is 2:1.  Chapter 6. Conclusions  240  velocities in the 7.0 to 7.4 km/s range. For velocities greater than 7.4 km/s, a separate layer (garnet granulite) was included. These high velocity materials are probably the result of underplating [Furlong and Fountain, 1986]. The model shows this layer to be thickest under the anomalous upper crust of the Val Rita block. Regions with velocities greater than 7.9 km/s were designated ultramafic mantle. Figure 6.3 shows the geophysical and derived geological models for Line 5 which crosses the Chapleau block of the KSZ. The surface greenstones at the east and west ends of the models show higher densities. The up thrust feature has not been faulted as much as the G R B , and so appears as a continuous zone from the surface into the midcrust. The Wawa tonalitic gneisses have low velocities and lower densities; these rocks are a probable constituent for the low velocity zone in the east of the model and in other regions of the greenstone belts. The Abitibi orthogneisses, to the east of the ILCZ at 170 km, have a low a and this suggests a high quartz content for these rocks. A high velocity anomaly, dipping to the west at 1 3 ° ± 3 ° , is imaged beneath the southern block of the KSZ. The anomaly extends from 20 km depth in the west to the surface in the KSZ and is coincident with higher densities in the near surface region between 120 and 170 km. This corresponds to an upthrusting of the granulite zone from shallower depths than initially proposed by Percival and Card [1983].  Although the high velocity feature is imaged  only in the upper half of the present day crust, the poorer resolution at depth may have failed to delineate a downward extension of this anomaly.  Alternatively, the ambient  pressure, temperature and fluid conditions at that depth may, over time, have altered the rocks and eliminated any velocity anomalies. These suppositions would project the thrust plane into the lower crust and infer a thicker Archean crust. However, the seismic reflection data indicates that the anomaly is confined to the upper crust. The placement of the sole of the thrust at 17 to 20 km depth agrees well with the geobarometric and gravity modeling results [Percival and McGrath, 1986] and with an initial interpretation of the reflection data [Wes£ et al., 1988]. Again, the granulite/tonalite layer is underlain  Chapter 6. Conclusions  241  DISTANCE  -30  20  70  Michipieoten  Wawa  120  170  Chapleau  Y^L\ metavolcanics  (KM)  220  270  370  420  Abitibi  pyroxene granulite / tonalite  |t j j | metasedimentary migmatite  | : | | | pyroxene granulite  j  | ^ garnet granulite  j orthogneiss / plutons  320  Figure 6.3. (a) Geophysical model for Line 5. The contours are P-wave velocity in km/s. Superimposed on this are: variations in density (p) away from the polynomial fit to the Nafe-Drake set of points; deviations of Poisson's ratio (cr) from the starting value of 0.25. (b) Geological model derived from the above variations, surface exposures and fault mapping. The heavy solid line is the thrust fault that soles into the detachment zone. The dashed lines are normal faults. Vertical exaggeration is 2:1.  Chapter 6. Conclusions  242  by a more mafic and higher velocity lower crust.  The bottom garnetiferous layer is  thicker here, especially under the KSZ where it is greater than 10 km thick. This can be explained by large episodes of underplating associated with A-type subduction or isostatic compensation. The results of this analysis demonstrate the usefulness of integrating information from geological and geophysical surveys. A self-consistent geophysical model for the Kapuskasing uplift has been constructed and successfully linked to the penological studies and the mapped surface exposures.  A need for a more detailed seismic refraction ex-  periment (40 km shot spacing and 1 km receiver spacing) is apparent and would prove more definitive for a complex geological situation such as the KSZ. 2d finite element rheological modeling [e.g., Braun and Beaumont, 1989] would be useful in constraining possible tectonic explanations and the proposed drilling project [Percival, 1988] will give insights into the nature of the reflectors and the relationship of the velocity anomalies to the depth distribution of lithologies.  References  Anhaeusser, C. R., Cyclic volcanicity and sedimentation in the evolutionary development of Archean greenstone belts of shield areas, Geol. Soc. Australia Spec. Publ., 3, 57-70, 1971. Anhaeusser, C . R., R. Mason, M . J . Viljoen, and R. P. Viljoen, A reappraisal of some aspects of Precambrian shield geology, Geol. Soc. Am. Bull, 80, 2175-2200, 1969. Au, C . Y . D., Crustal Structure from an Ocean Bottom Seismometer Survey of the Nootka Fault Zone, Ph. D. Thesis, Univ. of British Columbia, 1981. Ashwal, L . D., P. Morgan, S. A. Kelley, and J . A. Percival, Heat production in an Archean crustal profile and implications for heat flow and mobilisation of heat producing elements, Earth Planet. Sci. Lett, 85, 439-450, 1987. Assumpcao, M . , and D. Bamford, LISPB - V . Studies of crustal shear waves, Geophys. J. R. Astron. Soc, 54, 61-73, 1978. Ave 'Lallemant, H . G . , and N. I. Carter, Syntectonic recrystallisation of olivine and modes of flow in the upper mantle, Geol. Soc. Am. Bull, 81, 2203-2220, 1970. Bailey, R. C , J . A. Craven, J . C.Macnae, and B. D. Polzer, Imaging of deep fluids in Archaean crust, Nature, 340, 136-138, 1989. Baranov, V . , A new method for interpretation of aeromagnetic maps: pseudo-gravimetric anomalies, Geophysics, 12, 359-383, 1957. Baranov, V . , and H . Naudy, Numerical calculation of the formula of reduction to the magnetic pole, Geophysics, 29, 517-531, 1964. Berlage, H . P. Jnr., Seismometer, Auswertung der diagramme, Handbuk der Geophysik, 2, 299-526, 1930. Berry, M . J . , and G . F. West, An interpretation of the first arrival data of the Lake Superior experiment by the time-term method, Bull. Seismol. Soc. Am., 56, 141171, 1966. Berry, M . J . , and K. Fuchs, Crustal structure of the Superior and Grenville provinces of the northeastern Canadian shield, Bull. Seismol. Soc. Am., 63, 1393-1432, 1973. Bessanova, E . N., V . M . Fishman, V . Z. Ryaboyi, and G . A. Sitnikova, The tau method  243  References  244  for inversion of travel-times — I. Deep seismic sounding data, Geophys. J. R. Astron. Soc, 36, 377-398, 1974. Bickle, M . J . , L . F . Bettaney, C. A . Boulter, D. L . Groves, and P. Morant, Horizontal tectonic interaction of an Archean gneiss belt and greenstones, Pilbara block, W . Australia, Geology, 8, 525-529, 1980. Birch, F. J . , Compressibility and elastic constants, in Handbook of Physical Constants, edited by S. P. Clark, 97-174, G . S. A. Memoir 97, New York, 1966. Bhattacharyya, B. K., Two dimensional harmonic analysis as a tool for magnetic interpretation, Geophysics, 30, 829-857, 1965. Bhattacharyya, B. K., Continuous spectrum of the total-magnetic-field anomaly due to a rectangular prismatic body, Geophysics, 31, 97-121, 1966. Braun, J . , and C . Beaumont, Dynamical models of the role of crustal shear zones in asymmetric continental extension, Earth Planet. Sci. Lett., in press, 1989. Brigham, E . O., The Fast Fourier Transform, 252 pp., Prentice Hall, Inc., New Jersey, 1974. Burke, K., and J . F . Dewey, Plume-generated triple junctions: key indicators in applying plate tectonics to old rocks, J. Geology, 81, 406-433, 1973. Bursnall, J . T . , Contact relatioships between the Kapuskasing structural zone and Abitibi subprovince in the vicinity of Ivanhoe Lake (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Nov., 1988. Buttkus, B., Homomorphic filtering-theory and practice, Geophys. Prospect., 23, 712748, 1975. Byerlee, J . D., Brittle-ductile transition in rocks, J. Geophys. Res., 73, 4741-4750, 1968. Catchings, R. D., and W. D. Mooney, Crustal structure of East Central Oregon: relation between Newberry Volcano and regional crustal structure, J. Geophys. Res., 93, 10081-10094, 1988. Cermak, V . , and A. M . Jessop, Heat flow, heat generation and crustal temperature in the Kapuskasing area of the Candian shield, Tectonophysics, 11, 287-303, 1971. Cerveny, V . , I. A. Molotkov, and I. Psencik, Ray Method in Seismology, 214 pp., University of Karlova, Prague, Czechoslovakia, 1977. Chapman, C. H . , The Radon transform and seismic tomography, in Seismic Tomography, edited by G . Nolet, 25-47, Reidel Publ. Co., Dordrecht, 1987.  References  245  Chapman, D. S., Thermal gradients in the continental crust, in The Nature of the Lower Continental Crust, edited by J . B. Dawson, D. A. Carswell, J . Hall, and K. H . Wedepohl, 63-70, Geol. Soc. Sp. Publ. 24, 1986. Chouteau, M . , M . Mareschal, R. Chakridi and K . Bouchard, Preliminary results of a magnetotelluric survey across the Groundhog River block of the Kapuskasing Zone (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Feb., 1988. Christensen, N. I., and D. L . Szymanski, Origin of reflections from the Brevard Fault Zone, J. Geophys. Res., 93, 1087-1102, 1988. Claerbout, J . F. , Fundamentals of Geophysical Data Processing, 274 pp., McGraw Hill Book Co., 1976. Clarke, G . K. C , Optimum second-derivative and downward-continuation filters, Geophysics, 34, 424-437, 1969. Clayton, R. W., and R. A. Wiggins, Source shape estimation and deconvolution of teleseismic bodywaves, Geophys. J. R. Astron. Soc, 4?, 151-177, 1976. Clement, W. G . , Basic principles of two-dimensional digital filtering, Geophys. Prosp., 21, 125-145, 1972. Condie, K. C , Archean geotherms and supracrustal assemblages, Tectonophysics, 105, 29-41, 1984. Constanzo-Alvarez, V . , and D. J . Dunlop, Palaeomagnetic evidence for post-2.55 Ga tectonic tilting and 1.1 Ga reactivation in the southern Kapuskasing zone, Ontario, Canada, J. Geophys. Res., 93, 9126-9136, 1988. Cook, F. A., Geometry of the Kapuskasing structure from a Lithoprobe pilot reflection survey, Geology, 13, 368-371, 1985. Cormier, V . F., and P. Spudich, Amplification of ground motion and waveform complexity in fault zones: examples from the San Andreas and Calaveras faults, Geophys. J. R. Astron. Soc, 79, 135-152, 1984. Coward, M . P., and J. D. Fairhead, Gravity and structural evidence for the deep structure of the Limpopo belt, Southern Africa, Tectonophysics, 68, 31-43, 1980. Darby, E . K., and E . B. Davies, The analysis and design of two-dimensional filters for two-dimensional data, Geophys. Prosp., 15, 383-406, 1967. Davies, G . F., Thickness and thermal history of continantal crust and root zones, Earth Planet. Sci. Lett., 44, 231-238, 1979.  References  246  Der, Z. A. and T . W. McElfresh, The relationship between anelastic attenuation and regional amplitude anomalies of short-period P waves in North America, Bull. Seismol. Soc. Am., 67, 1303-1317, 1977. Dewey, J . F., and K. C. Burke, Tibetan, Variscan, and Precambrian basement reactivation: products of continental collision, J . Geol., 81, 683-692, 1973. deWit, M . J . , and L . D. Ashwal (eds.), Workshop on Tectonic Evolution of Greenstone Belts, 1-4, LPI Tech. Rept., 86-10, Lunar and Planetary Institute, Houston, 1986. Drummond, B. J . , Seismic P-wave anisotropy in the subcrustal lithosphere of north-west Australia, Geophys. J. R. Astron. Soc, 81, 497-519, 1985. England, P., and M . Bickle, Continental thermal and tectonic regimes during the Archean, J. of Geology, 92, 353-367, 1984. Ernst, R. E . , and H. C. Halls, Palaeomagnetism of the Hearst dike swarm and implications for the tectonic history of the Kapuskasing Structural Zone, Northern Ontario, Can. J. Earth Sci., 21, 1499-1506, 1984. Fokkema, J . T . , and A. Ziolkowski, The critical reflection theorem, Geophysics, 52, 965972, 1987. Foland, K . A . , Limited mobility of argon in a metamorphic terrane, Geochim. mochim. Acta, 43, 793-801, 1979.  Cos-  Forsyth, D., P. Morel, H. Hasegawa, R. Wetmiller, J . Adams, A . Goodacre, D. Nagy, R. Coles, J . Harris, and P. Basham, Comparative study of the geophysical and geological information in the Timiskaming-Kapuskasing area, Earth Phys. Br. Tech. Rec. TR238, 1983. Fountain, D. M . , D. T. McDonagh, and J . M . Gorham, Seismic reflection models of continental crust based on metamorphic terrains, Geophys. J. R. Astron. Soc, 89, 61-66, 1987. Fountain, D. M . , and M . H. Salisbury, Exposed cross sections through the continental crust: implications for crustal structure, petrology and evolution, Earth Planet. Sci. Lett, 56, 263-277, 1981. Fountain, D. M . , and M . H . Salisbury, Seismic properties of the Superior Province crust based on seismic velocity measurements on rocks from the Michipicoten-WawaKapuskasing terranes, Ontario (abstract), Geol. Assoc. Can., Min. Assoc. Can., Can. Geophys. U., Prog., 11, p. 69., 1986. Fountain, D. M . , M . H. Salisbury, and J . A. Percival, Seismic structure of the continent crust based on rock velocity measurements from the Kapuskasing uplift, submitted to  References  247  J. Geophys. Res., 1989. Fuller, B. D., Two-dimensional frequency analysis and design of grid operators, in Mining Geophys., Theory, II, 658-708, S.E.G., Tulsa, Oklahoma, 1967. Furlong, K. P., and D. M . Fountain, Continental crustal underplating: thermal considerations and seismic-petrologic consequences, J. Geophys. Res., 91, 8285-8294, 1986. Fyfe, W . S., Archean tectonics, Nature, 249, 338, 1974. Gajewski, D. K. Fuchs, C . Prodehl, and R. Stangl, Anomalous disappearance of S phases near the continental Moho (abstract), International Union of Geodesy and Geophysics, XIX Assembly, Prog., 1, p. 58, 1987. n  Garland, G . D., Interpretations of gravimetric and magnetic anomalies on traverses in the Canadian Shield in Northern Ontario, Publ. Dom. Obs., Ottawa, 16, 5-57, 1950. Garmany, J . , J . A . Orcutt, and R. L . Parker, Travel time inversion: a geometrical approach, J. Geophys. Res., 84, 3615-3622, 1979. Geis, W . T . , A . G . Green, and F . A . Cook, Preliminary results from the high resolution seismic reflection profile: Chapleau Block, Kapuskasing Structural Zone (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Feb., 1988. Gibb, R. A . , Collision tectonics in the Canadian shield, Earth Planet. Sci. Lett., 27, 378-382, 1975. Gibb, R. A., A gravity survey of James Bay and its bearing on the Kapuskasing gneiss belt, Ontario, Tectonophysics, 45, T7-T13, 1978. Gibb, R. A., Model for suturing of superior and Churchill plates: an example of double indentation tectonics, Geology, 11, 413-417, 1983. Gibb, R. A . , and R. I. Walcott, A Precambrian suture in the Canadian Shield, Earth Planet. Sci. Lett., 10, 417-422, 1971. Glikson, A. Y . , Early Precambrian evidence of a primitive ocean crust and island nuclei of sodic granite, Geol. Soc. Am. Bull, 83, 3323-3344, 1972. Glikson, A . Y . , Stratigraphy and evolution of primary and secondary greenstones: significance of data from shields of the southern hemisphere, in The Early History of the Earth, edited by B. F . Windley, 257-277, J . Wiley and Sons, London, 1976. Glikson, A. Y . , Uniformatarian assumptions, plate tectonics and the Precambrian earth, in Precambrian Plate Tectonics, edited by A . Kroner, 91-104, Elsevier Sci. Publ. Co.,  References  248  Amsterdam, 1981. Goleby, B. R., R. D. Shaw, C. Wright, B. L. N. Kennett, and K. Lambeck, Geophysical evidence for 'thick-skinned' crustal deformation in central Australia, Nature, 337, 325330, 1989. Goodings, C. R., The Kapuskasing structure and its relationship to the Proterozoic movements in the Superior province (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Nov., 1988. Goodwin, A. M . , The nature of Archean crust in the Canadian shield, in Evolution of the Earth's Crust, edited by D. H . Tarling, 175-218, Academic Press, New York, 1978. Goodwin, A . M . , Archean plates and greenstone belts, in Precambrian Plate Tectonics, edited by A. Kroner, 105-135, Elsevier Sci. Publ. Co., Amsterdam, 1981. Grad, M . and U. Luosto, Seismic models of the crust of the Baltic shield along the S V E K A profile in Finland, Annates Geophysicae, 5B, 639-650, 1987. Grant, F. S., and G . F. West, Interpretation Theory in Applied Geophysics, McGraw-Hill Book Co., New York, 1965. Green, A. G . , N. L . Anderson, and 0. G . Stephenson, A n expanding spread seismic reflection survey across the Snake Bay-Kakagi greenstone belt, N W Ontario, Can. J. Earth Sci., 16, 1599-1612, 1979. Gunn, P. J . , Linear transformations of gravity and magnetic fields, Geophys. Prosp., 23, 300-312, 1975. Gupta, V . K . , and F. S. Grant, Mineral-exploration aspects of gravity and aeromagnetic surveys in the Sudbury-Cobalt area, Ontario, in The Utility of Regional Gravity and Magnetic Anomaly Maps, edited by W . J . Hinze, 392-411, S.E.G., Tulsa, Oklahoma, 1985. Gupta, V . K . , and N. Ramani, Some aspects of regional-residual separation of gravity anomalies in a Precambrian terrain, Geophysics, 45, 1412-1426, 1980. Hahn, A., E . G . Kind, and D. C. Mishra, Depth estimation of magnetic sources by means of Fourier amplitude spectra, Geophys. Prosp., 24, 287-308, 1976. Hall, D. H . , and Z. Hajnal, Deep crustal studies in Manitoba, Bull. Seismol. Soc. Am., 63, 885-910, 1973. Hall, J . , and M . Ah, Shear waves in a seismic survey of Lewisian basement: control on lithological variation and porosity, J. Geol. Soc. 1985.  an extra  Lond., 142, 677-688,  References  249  Halls, H . C . , Crustal thickness in the Lake Superior region, Geol. Soc. Am. Mem., 156, 239-243, 1982. Halls, H . C , and H. C. Palmer, A palaeomagnetic and petrographic study of dykes in deep crust: implications for the evolution of the Kapuskasing zone and two early Proterozoic dyke swarms (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Nov., 1988. Hodges, K. V . , J . M . Bartley, and B. C. Burchfiel, Structural evolution of an A-type subduction zone, Lofoten-Rombak area, Northern Scandanavian Caledonides, Tectonics, 1, 441-462, 1982. Hoffman, P. F., United plates of America, the birth of a craton: early Proterozoic assembly and growth of Laurentia, Ann. Rev. Earth Planet. Sci., 16, 543-603, 1988. Holbrook, W . S., D. Gajewski, A . Krammer, and C. Prodehl, A n interpretation of wideangle compressional and shear wave data in south-west Germany: Poisson's ratio and penological implications, J. Geophys. Res., 93, 12081-12106, 1988. Hurich, C . A . , and S. B. Smithson, Compositional variation and the origin of deep crustal reflections, Earth Planet. Sci. Lett, 85, 416-426, 1987. Innes, M . J . S., Gravity and isostasy in Northern Ontario and Manitoba, Publ. Dom. Obs., Ottawa, 21, 265-338, 1960. Jackson, I., and R. J . Arculus, Laboratory wave velocity measurements on lower crustal xenoliths from Calcutteroo, South Australia, Tectonophysics, 101, 185-197, 1985. Jacobsen, B. H . , A case for upward continuation as a standard separation filter for potential-field maps, Geophysics, 52, 1138-1148, 1987. Jin, D. J . , and E . Eisner, A review of homomorphic deconvolution, Rev. Geophys. Space Phys., 22, 255-263, 1984. Jones, T . D., and A. Nur, The nature of seismic reflections from deep crustal fault zones, J. Geophys. Res., 89, 3153-3171, 1984. Kanasewich, E . R., Time Sequence Analysis in Geophysics, 237-279, Univ. of Alberta Press, 1981. Kanasewich, E . R., and R. G . Agarwal, Analysis of combined gravity and magnetic fields in wave number domain, J. Geophys. Res., 75, 5702-5712, 1970. Kanasewich, E . R., and S. K. Chiu, Least squares inversion of spatial seismic refraction data, Bull. Seismol. Soc. Am., 75, 865-880, 1985.  References  250  Klemperer, S. L., and J . H. Luetgert, A comparison of reflection and refraction processing and interpretation methods applied to conventional refraction data from coastal Maine, Bull. Seismol. Soc. Am., 77, 614-630, 1987. Kleywegt, R. J . , E . H . Stettler, G . Brandl, R W. Day, J . H. de Beer and A. W . A. Duvenhage, The structure of the Giyani greenstone belt as derived from geophysical studies, S. -Afr. Tydskr. Geol, 90, 282-295, 1987. Kroner, A . , Evolution, growth and stabilisation of the Precambrian lithosphere, Phys. Chem. Earth, 15, 69-106, 1985. Kroner, A . , Precambrian plate tectonics, in Precambrian Plate Tectonics, edited by A. Kroner, 57-90, Elsevier Sci. Publ. Co., Amsterdam, 1981. Kurtz, R. D., E . R. Niblett, J . A. Craven, R. A. Stevens, and J . C. Macnae, Electromagnetic studies over the Kapuskasing structural zone (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Feb., 1988. Kusznir, N. J . , and R. G . Park, Intraplate lithosphere deformation and the strength of the lithosphere, Geophys. J. R. Astron. Soc, 79, 513-538, 1984. Kusznir, N. J . , and R. G . Park, Continental lithosphere strength: the critical role of lower crustal deformation in The Nature of the Lower Continental Crust, edited by J . B. Dawson, D. A. Carswell, J . Hall, and K . H . Wedepohl, 79-93, Geol. Soc. Sp. Publ. 24, 1986. Lambeck, K . , Structure and evolution of the intracratonic basins of central Australia, Geophys. J . R. Astron. Soc, 74, 843-886, 1983. Lane, M . , Multi-Channel Homomorphic Wavelet Extraction, M . Sc. Thesis, Univ. of British Columbia, 1982. Leclair, A . , Results of geological mapping in the Kapuskasing area (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Nov., 1988. LeQuentrec, M . F . , J . C. Mareschal, and 0. Parphenuk, A finite element model of the thermal evolution of the Kapuskasing structure (abstract), Geol. Assoc. Can., Min. Assoc. Can., Prog., U, p. A103., 1989. Ludwig, J . W . , J . E . Nafe, and C. L . Drake, Seismic refraction, in The Sea, edited by A. E . Maxwell, 53-84, J . Wiley and Sons, New York, 1970. Lynn, H., and G . Thompson, Depth migration and interpretation of the C O C O R P Wind River, Wyoming, seismic reflection data, Geology, 11, 462-469, 1983. Maclaren, A. S., and B. W . Charbonneau, Characteristics of magnetic data over major  References  251  subdivisions of the Canadian Shield, Proc. Geol. Assoc. Can., 19, 57-65, 1968. McMechan, G. A., and G. S. Fuis, Ray equation migration of wide-angle reflections from southern Alaska, J. Geophys. Res., 92, 407-420, 1987. Mareschal, J . C , and M . F. LeQuentrec, Two dimensional modelling of the thermal regime of the crust during the Kapukasing uplift (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Nov., 1988. Mareschal, M . , R. Chakridi, and M . Chouteau, A magnetotelluric survey across the Groundhog River block: progress report on the pseudo Id interpretation (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Nov., 1988. Meltzer, A. S., A. R. Levander, and W . D. Mooney, Upper crustal structure, Livermore valley and vicinity, California Coast ranges, Bull. Seismol. Soc. Am., 77, 1655-1673, 1987. Mereu, R. F . , D. Wang, O. Kuhn, D. A. Forsyth, A. G . Green, P. Morel, G . G . R. Buchbinder, D. Crossley, E . Schwarz, R. duBerger, C. Brooks, and R. Clowes, The 1982 C O C R U S T seismic experiment across the Ottawa-Bonnechere graben and Grenville Front in Ontario and Quebec, Geophys. J. R. Astron. Soc, 84, 491-514, 1986. Mooney, W . D., and T . M . Brocher, Coincident seismic reflection/refraction studies of the continental lithosphere: a global review, Rev. Geophys., 25, 723-742, 1987. Morgan, P., Crustal radiogenic heat production and the selective survival of ancient continental crust, J. Geophys. Res., 90, c561-c570, 1985. Moser, D., Structure of the Wawa gneiss terrane near Chapleau, Ontario, in Current Research, 88-lc, 93-99, Geol. Sur. Canada, 1988. Mueller, St., and J . Ansorge, Long range seismic refraction profiles in Europe, in Reflection Seismology: A Global Perspective, edited by M . Barazangi and L . Brown, 167-181, A G U Geodyn. Ser., 14, 1986. Newton, R. C , Petrologic aspects of Precambrian granulite facies terrains bearing in their origins, in Proterozoic Lithospheric Evolution, edited by A. Kroner, 11-26, Geodyn. Ser. 17, A . G . U . , Washington, D.C., 1987. Newton, R. C , and D. Perkins III, Thermodynamic calibration of geobarometers based on the assemblages garnet-plagioclase- orthopyroxene (clinopyroxene)-quartz, Amer. Min., 67, 203-222, 1982. Nisbet, E . G . , and C. M . R. Fowler, Model for Archean tectonics, Geology, 11, 376-379, 1983.  References  252  Nkwate, E . A., and M . H. Salisbury, New gravity mapping of the Kapuskasing structure in the Val Rita and Groundhog River blocks: prehminary results (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Nov., 1988. Northey, D. J . , and G . F. West, The Kapuskasing Structural zone seismic refraction experiment, DSS Contract Rept. OST84-00060 06ST.23235-3-1478, 117 pp., Geophys. Lab., Physics Dept., Univ. of Toronto, 1985. Northey, D. J . , and G. F. West, A crustal scale seismic refraction experiment over the Kapuskasing Structural Zone (abstract), Geol. Assoc. Can., Min. Assoc. Can., Can. Geophys. U., Prog., 11, p. 108, 1986. Nur, A., and G . Simmons, The effect of saturation on velocity in low porosity rocks, Earth Planet. Sci. Lett., 7, 183-193, 1969. O'Brien, P. N. S., Lake Superior structure - A reinterpretation of the 1963 experiment, J. Geophys. Res., 73, 2669-2689, 1968. Ojo, S. B., and R. F. Mereu, The effect of random velocity functions on the travel-times and amplitudes of seismic waves, Geophys. J. R. Astron. Soc, 84, 607-618, 1986. Oppenheim, A. V . , Superposition in a class of nonlinear systems, Tech. Rept. 432, Res. Lab. Elec, MIT, Cambridge, Mass., 1965. Oppenheim, A. V . , and R. W. Schafer, Digital Signal Processing, 585 pp., Prentice Hall, Inc., New Jersey, 1975. Ord, A., and B. E . Hobbs, The strength of the continental crust, detachment zones and the development of plastic instabilities, Tectonophysics, 158, 269-289, 1989. Otis, R. M . , and R. B. Smith, Homomorphic deconvolution by log spectral averaging, Geophysics, 42, 1146-1157, 1977. Parker, C. L . , Crustal structure of the Abitibi greenstone belt determined from refraction seismology, M . Sc. Thesis, McGill University, Montreal, 1984. Percival, J . A., High grade metamorphism in the Chapleau-Foyelet area, Ontario, Am. Mineral., 68, 667-686, 1983. Percival, J . A., A possible exposed Conrad discontinuity in the Kapuskasing uplift, Ontario, in Reflection Seismology: The Continental Crust, edited by M . Barazangi and L. Brown, 135-141, A G U Geodyn. Ser., 14, 1986. Percival, J . A., The basal fault zone of the intracratonic Kapuskasing uplift, in Canadian Continental Drilling Program, Scientific Drilling: Major Faults, edited by M . J . Drury, Can. Cont. Drill. Prog. Rept. 88-3, 1988.  References  253  Percival, J . A., and K. D. Card, The Archean crust as revealed in the Kapuskasing uplift, Superior Province, Canada, Geology, 11, 323-326, 1983. Percival, J . A . , and K. D. Card, Structure and evolution of Archean crust in central Superior province, Canada, in Evolution of Archean Supracrustal Sequences, edited by L. D. Ayres, P. C. Thurston, K. D. Card, and W . Weber, 179-192, Geol. Assoc. Can., Sp. Paper 28, 1985. Percival, J . A., and D. M . Fountain, Metamorphism and melting at an exposed example of the Conrad discontinuity, Kapuskasing uplift, Canada, in Fluid Movements, Element Transport and the Composition of the Crust, edited by D. Bridgwater, D . Reidel, Dordrecht, in press, 1989. Percival, J . A., and T . E . Krogh, U-Pb zircon geochronology of the Kapuskasing structural zone and vicinity in the Chapleau- Foleyet area, Can. J. Earth Sci., 20, 830-843, 1983. Percival, J . A., and P. H . McGrath, Deep crustal structure and tectonic history of the Northern Kapuskasing uplift of Ontario: An integrated petrological-geophysical study, Tectonics, 5, 553-572, 1986. Percival, J . A., R. R. Parrish, T. E . Krogh, and Z. E . Peterman, When did the Kapuskasing zone come up? (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Feb., 1988. Poggiagliolmi, E . , A . J . Berkhout, and M . M . Boone, Phase unwrapping, possibilities and limitations, Geophys. Prospect., 30, 281-291, 1982. Raase, P., M . Raith, D. Ackermand, and R. K. Lai, Progressive metamorphism of mafic rocks from greenschist to granulite facies in the Dharwar craton of South India, J. Geology, 94, 261-282, 1986. Ricard, Y . , and R. J . Blakely, A method to minimise edge effects in two-dimensional discrete Fourier transforms, Geophysics, 53, 1113-1117, 1988. Rudnick, R. L . , and S. R. Taylor, The composition and petrogenesis of the lower crust: a xenolith study, / . Geophys. Res., 92, 13981-14005, 1987. Sandmeier, K . J . , and F . Wenzel, Synthetic seismograms for a complex crustal model, Geophys. Res. Lett., 13, 22-25, 1986. Shive, P. N., Can remanent magnetisation in the deep crust contribute to long wavelength magnetic anomalies?, Geophys. Res. Lett., 16, 89-92, 1989. Shive, P. N., and D. M . Fountain, Magnetic mineralogy in an Archean crustal cross section: implications for crustal magnetisation, J. Geophys. Res., 93, 12177-12186,  References  254  1988. Sibson, R. H., Frictional constraints on thrust, wrench and normal faults, Nature, 249, 542-544, 1974. Silva, J . B. C , Reduction to the pole as an inverse problem and its application to lowlatitude anomalies, Geophysics, 51, 369-382, 1986. Singh, S., and R. B. Herrmann, Regionalization of crustal coda Q in the continental United States, J. Geophys. Res., 88, 527-538, 1983. Sinitsyn, A. V . , Origin of Precambrian greenstone Taelts, Geotectonics, 13, 423-433, 1979. Sleep, N. H., and B. F. Windley, Archean plate tectonics: constraints and inferences, J. of Geology, 90, 363-379, 1982. Spector, A., and F. S. Grant, Statistical models for interpreting aeromagnetic data, Geophysics, 35, 293-302, 1970. Spence, G . D., R. M . Clowes, and R. M . Ellis, Seismic structure across an active subduction zone of western Canada, / . Geophys. Res., 90, 6754-6772, 1985. Talwani, M . , J . L . Worzel, and M . Landisman, Rapid gravity computations for twodimensional bodies with application to the Mendocino submarine fracture zone, / . Geophys. Res., 64, 49-59, 1959. Tarney, J . , I. W . Dalziel, and M . J . deWit, Marginal basin 'Rocas Verdes' complex from S. Chile: a model for Archean greenstone belt formation, in The Early History of the Earth, edited by B. F. Windley, 131-146, J . Wiley and Sons, London, 1976. Thomas, M . D., L . Loiser, P. C. Thurston, V . K . Gupta, R. A. Gibb, and R. A. F. Grieve, Geophysical characteristics and crustal structure of greenstone terranes, Canadian shield, in Tectonic Evolution of Greenstone Belts, edited by M . J . deWit and L. D. Ashwal, 203-206, L.P.I. Tech. Rpt., 86-10, 1986. Thorarinsson, F . , S. G . Magnusson, and A. Bjornsson, Directional spectral analysis and filtering of geophysical maps, Geophysics, 53, 1587-1591, 1988. Tribolet, J . M . , A new phase unwrapping algorithm, IEEE Trans. Acoust. Speech Signal Process., 25, 170-177, 1977. Tribolet, J . M . , Seismic Applications of Homomorphic Signal Processing, 195 pp., Prentice Hall, Inc., New Jersey, 1979. Ulrych, T. J . , Application of homomorphic deconvolution to seismology, Geophysics, 36, 650-660, 1971.  References  255  Ulrych, T . J . , 0. G. Jensen, R. M . Ellis, and P. G . Somerville, Homomorphic deconvolution of some teleseismic events, Bull. Seismol. Soc. Am., 62, 1269-1281, 1972. Watson, J . , Precambrian thermal regimes, Phil. Trans. R. Soc. Lond. A, 288, 431-440, 1978. Watson, J . , The origin and history of the Kapuskasing structural zone, Ontario, Canada, Can. J. Earth Sci., 17, 866-875, 1980. Weidmann, S., Diffraction phenomena in seismic models, Geophys. J. R. Astron. 79, 207-216, 1984.  Soc,  West, G . F . , A. G . Green, F. Cook, B. Milkereit, P. Hurley, and W. Geis, Kapuskasing seismic reflection survey programme (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Nov., 1988. White, D. J . , Two-dimensional refraction tomography, Geophys. J., 97, 223-246, 1989a. White, D. J . , Shallow Crustal Structure beneath the Juan de Fuca Ridge from 2° Seismic Refraction Tomography, Ph. D. Thesis, Univ. of British Columbia, 1989b. White, S. H., and P. G. Bretan, Rheological controls on the geometry of deep faults and the tectonic delamination of the continental crust, Tectonics, 4, 303-309, 1985. Wilks, K. R., and N. D. Carter, Mechanical behaviour of the Adirondack and Pikwitonei granulites, EOS Transactions, 70, p. 477, American Geophysical Union, 1989. Wilson, J . T . , Comparison of the Hudson Bay arc with some other features, in Science, History and Hudson Bay, edited by C. S. Beals and D. A. Shenstone, Dep. Energy Mines Res., Ottawa, 1968. Windley, B. F . , The Early History of the Earth, 619 pp., J . Wiley and Sons, London, 1975. Windley, B. F . , The Evolving Continents, 399 pp., J . Wiley and Sons, Avon, 1986. Worthington, M . H . , An introduction to geophysical tomography, First Break, 20-26, Nov. 1984. Wu, J . , and R. F. Mereu , Crustal models of the Kapuskasing Structural Zone: Results from the 1984 seismic refraction experiment (ext. abstract), Kapuskasing Structural Zone Transect Workshop, Univ. of Toronto, Feb., 1988. Yilmaz, 0., Seismic Data Processing, 526 pp., Investigations in Geophysics, 2, edited by S. M . Doherty, Society of Exploration Geophysicists, Tulsa, 1988. Zelt, C. A . , Seismic Structure of the Crust and Upper Mantle in the Peace River Arch  References  256  region, Ph. D. Thesis, Univ. of British Columbia, 1989. Zelt, C. A . , and R. M . Ellis, Practical and efficient ray tracing in two-dimensional media for rapid travel-time and amplitude modelling, Can. J. Explor. Geophys., 24, 16-31, 1988. Zeyen, H . J . , E . Banda, J . Gallart, and J . Ansorge, A wide angle seismic reconnaissance survey of the crust and upper mantle in Celtiberian Chain of eastern Spain, Earth Planet. Sci. Lett., 75, 393-402, 1985. Ziolkowoski, A . , Deconvolution, 175 pp., Human Res. Devel. Corp., Boston, 1984.  Appendix A  Potential Field T h e o r y  The basic equations for the gravity and magnetic fields [Grant and West, 1965; Bhattacharyya, 1966; Fuller, 1967] are presented for the sake of completeness. They are presented in a similar style after the formulation of Gunn [1975]. This presentation facilitates comparisons and the application of linear transformations.  Gravity Field The force per unit mass on a point P, a distance r from a mass m, is 777.  F(r) = - G - r 7*  (Al)  Such a field is conservative and so is derivable from a scalar function  F(r) = - W ( r )  G  (A2)  m  where U(r) — ~ is the gravitational potential of m. Since potentials in free space are additive, the gravitational potential due to a mass distribution with density p(ro) at an exterior point P is r  £«,) =  - / 3  v  , j p(r )d r = —G Jv |f-f | 0  '  0  0  where |? — ?o| = \Jr +r% — 2rr cos <f>. If we consider the bulk density distribution p(fo) = pv{a,B,z), then at a point h above this distribution the potential is 2  0  U(x,y,h) = -G  +  [ °°  +0  l°  +  [ °° -  2  [( -a) z  257  Pv(<*,P,z) + (y-a)* +  d  {z-hy]i  a  d / 3  d  z  (  A  4  )  Appendix A. Potential Field Theory  258  where a ,B, and z are the coordinates of the elemental mass and z is positive downwards. The inner part of the equation is in the form of a convolution integral f + oo / / -oo • , - « > \( -a) -oo J — oo \(T. -  o (a B z) ' dadj3 = p(x,y,z)**R(x,y,z-h) + (y - B) + (z - h) ]2 -I- (v - BV 4- (z - hW 2  + oo  K  (A.5)  T  2  2  x  2  '  v  '  where ** denotes two-dimensional convolution and  x  R( >y,  z  - h) =  2  2  [x +  2  y -r(z-h) )-i  When two terms are convolved in the space domain, then they are multiplied in the frequency domain. So, if we take the Fourier transform of both sides we get f + oo U(u,v,h) =—G  p (u,v,z)R(u,v,z  — h)dz  v  Jo  where the Fourier transform pair is defined as I  p (u,v,z) v  f+oo  = —  t+oo ux+vy  /  p (x,y,z)e-*  >  v  dx dy  Z7T •/ —oo J—oo  p (x,y,z) v  = — I  /  (A.7a,b) +  m  p (u,v,z)e * ^  du dv  v  ZTT J— OO J — OO  and U(x, y,z — h) and R(x, y,z — h) are transform pairs with U(u, v,z —h) and R(u, v, z — h), respectively. Bhattacharya (1966) has shown that  f+oo  /  -j(ux+vy)  f+oo  2  lre  /  J-oo  J-oo  r  [a;2  + y  2  +  2  2 -{^-h){u +v )h  e  dxdy =  ( _fc)2]5 2  *  ; 2  2  (u +V )2  (A8) '  Using this integral solution, we obtain  ,+oo U(u,v,h) =-2TTG  -(z-h)(u*+v*)l  e  p (u,v,z)  1—dz Jo (u -\-v )i In the field, we measure the vertical component of gravity and so v  2  2  (A9)  Appendix A. Potential Field Theory  259  which gives a  F {u,v,h) = 2-KG Jo  a  p {u,v,z) -(*-'0(« +'' )*  g  v  e  (A.10)  DZ  as the spectral representation of the gravitational field.  Magnetic Field Let us consider a magnetic dipole and its field at a point P at a distance r, which is much greater than the separation of the poles /. The potential at P due to the two poles is m m m m i 2 • T — | cos 6 r + | cos 6 r  R  m (, I = — 1 + — cos 0 + r \ 2r ml cos 8  I - 1 + — cos 0 + 2r  M cos 6 R  2  where m and M are the magnetic moments of the pole and dipole, respectively. If we now consider a bulk magnetisation distribution m (a, 3, z) at a point (x,y,h) then v  m (a,3,z)cos6 v  V(x,y,h)  -  [ ( z - a ) 2 - + (y - 3)* + (z -  h)*]i  for a single element in the magnetisation distribution. This is equivalent to  Vix.y,h)  — —6k  : [(a, _  a  r  ) 2 + (j, _ 3)2  +  (  z  _  fc)  2]i  where £  5  6  £  ^ is the directional derivative in the direction of magnetisation and L, M , and N are the direction cosines of the direction of magnetisation. For the whole distribution it follows that  V  ( „ , , , * ) = L  r r r——  ,  f(MM  (Jtll)  Appendix A. Potential Field Theory  260  This form of equation is the same as that developed for the gravity case except for the rr term. Thus 6k  V{u,v,h) =  m (u,v,z)—-  2tt— /  —r—dz  v  2  (u +  OK JO  (A.12)  2  V p  If / and F are Fourier transform pairs then  -f- = juF ox  ^f- = jvF by  Therefore equation 4.12 becomes  2  2  2TT\LJU + Mjv + N(u + v )2]  V(u,v,h) = —LJ.  J 2  \  [u -f- v  2  r+°°  )s  ,  s  ,  . a , 2  h  LA /  m (u,v,z)e-^- ^ v  a  +" > dz (A13)  •'0  To convert to magnetic field intensity, we have to differentiate in the direction of the component of the field being measured. This gives  F (u, v, h) = m  os  V(u, v, h)  where — - I— 8s Sx  m  —  8y  6  n 8z  j- is in the direction of the field and I, m, and n are direction cosines in the same direction as the component of the field being measured. Thus we have t  F (u,v,h) m  2  2  2  2  = 2T[LJU + Mjv + N(u + v )$][lju + mjv + n(u + v )*].  -^-—T 2  2  2  (u + v )2 Jo  r°°m {u,v, )e-^ ^dz v  Z  (A14) as the spectral representation of the magnetic field.  Appendix A. Potential Field Theory  261  Poisson's Relation The gravity and magnetic fields are now expressed in terms of equivalent layers. These forms of the equations are not neccessary for the 2d filtering of the data, but they serve to illustrate the principles of the convolutions! model as applied to potential fields. A density distribution p (a, 8, z) has an equivalent surface density distribution p (a, 8, d) on a surface at a depth, d, which gives the same gravity anomaly. This is called the equivalent layer. Similarly, a surface magnetisation distribution m (a,8, d) exists that produces the same magnetic anomaly as the bulk distribution m (ct,8,z). The expressions describing the potentials arising from these surface distributions are v  t  a  v  u( , ,h)=X y  v{x  r  G  p  r  -A^m  r d adB  ^ Lr~r~  v  h) = '  6k J-oo  J-oo  ( A 1 5 )  Tdad3  [(x  -  a)  2  + (y  -  0)  2  (A16)  2  + (d -  v  /i) ]5  ;  By the same approach as above, without having to integrate over z, we can obtain F {u,v,h) g  F (u,v,h) m  fc  = 2irG {u,v,d) Pt  2  = 2TT[LJU + Mjv + N(u  e  8  9  -<*- X« +« )*  (A17)  2  2  2  + v )^][lju + mjv + n(u + v )^]. 2  e  2  -(d-h)(u +« )i  m,(u,v,d)2  (u  +  (A18)  2  v )i  as the spectral representations of the gravity and magnetic fields in terms of equivalent layers. Assuming the same causative body exists for both fields and that the density and magnitude, and angle of magnetisation are constant throughout the body, the gravity and magnetic fields can be related by Poisson's relation [Grant and West, 1965]  m .VU a  = Gp V  -  v  ( A 19)  which gives , *r\ /.b~U , 6U. ^ _ m (iL + j M + kN . i — + j — + k — = G V ox by bz T  v  r  Pv  262  Appendix A. Potential Field Theory  where m is the magnetisation vector. Therefore a  which is in the same form as the original magnetic potential equation and indicates that if is constant, then we can calculate the distributions using  m {u,v,d) = —^ Gp t  p (u,v,d)  (A.21)  a  v  where m (u, v,d) and p,(u, v,d) are the same functions except for a scaling factor. Finally, we have s  m (x,y,d) t  = m  v  p (x,y,d) = p a  v  B(x,y,d) B(x,y,d)  where B is the function denning spatial distribution. The equivalent layer equations are true for geological situations where the direction of magnetisation is constant over the area and the magnitude of magnetisation and density is constant for each body. The equivalent layer is the sum of the individual equivalent layers of each isolated anomaly. The form of our spectral equations describing the gravity and magnetic fields is amenable to the application of linear transformations. Considering measurements taken at the surface of the earth, we get  F (u,v, 0) = 2TTG p,(u,v,d)  H(u,v,d)  g  (A.22a,b) F (u,v, m  0) = 2-K Di(u v) }  D (u,v) 2  I(u,v) m (u,v,d) e  H(u,v,d)  Appendix A. Potential Field Theory  263  where  2-KG  scaling factors 2TT 2  Di(u,v) = [Lju -f Mjv + N(u 2  2  2  2  D2(u,v) = [lju + mjv + n(u -f v ) ] H(u,v,d) = -  d i u 2 + v 2 )  e  *  factor for direction of magnetisation  2  + i> ) ]  factor for direction of measurement  depth factor  p,{u,v,d) equivalent layer factors m,(u,v,d) I(u,v)  extra factor distinguishing magnetic from gravity field. 2  (u  2  + V )?  From the convolution theorem, components which are multiplied together in the frequency domain are convolved in the space domain. This gives us  F = 2irG.p (x,y,d) g  s  * H(x,y,d) (A.2Za,b)  F  m  = 2-K.D (x,y) * D (x,y) * I(x,y) * m^x^y^) * H(x,y,d) l  2  Thus, since convolution is a linear operation, the order of the convolutions does not matter and we can remove the effect of any of the terms in the expressions by convolution with a suitable inverse filter. However, these operations are more easily implemented as divisions in the frequency domain and this is the method employed in this analysis.  

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

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"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0052973/manifest

Comment

Related Items