S23A-1859
Seismic Tomography of Crustal P and S in Eurasia
We present inversion results for Pg and Sg/Lg travel times for Eurasia using data from the LANL Research Knowledge Base. This database consists of integrated local, regional, and global catalogs with arrivals sorted and merged by event. The epicentral ground truth of each event is tested using the Bondar criteria and we weight events by the square root of their ground truth level with a minimum allowed weight of 1. We assign a GT level of 50 km to those events that do not pass any of the Bondar criteria. We use an inversion method analogous to Pn tomography, whereby we assume a straight ray between source and receiver and project travel time anomalies onto this great circle path. Events are restricted to depths less than 33 km and stations to distances greater than .5 degrees. To accommodate depth uncertainty and near-receiver velocity effects we solve for both event and site terms and damp the site term sum to zero. We employ first difference regularization over a one by one degree grid, and solve the set of equations using the LSQR conjugate gradient method. We invert for Pg and Sg/Lg slowness separately over a region from 0 to 80 degrees north latitude and -20 to 195 degrees east longitude. For Pg, the dataset is comprised of 3,709 stations, 407,131 events, and 1,453,318 travel times. The starting RMS error is reduced through inversion by 30% with respect to IASP91. For Sg/Lg, we use 3,084 stations, 266,751 events, and 1,049,125 travel times and achieve a 26% reduction in RMS error. While the Sg/Lg site terms are about twice as large as those for Pg, the Pg and Sg/Lg site terms show similar geographic patterns. Resolution tests with added Gaussian noise of 3 seconds standard deviation are performed for checkerboards having squares of 2, 4, and 10 degrees width. Pg velocities are highest in subduction zones and lowest in eastern Europe and around the Caspian Sea. The ratio of Pg to Lg velocities is high at zones of convergence, including the Pacific Rim and the Tibetan Plateau. Low Vp/Vs ratios are seen in parts of central Europe and southern Russia. Spot checks of our Vp/Vs ratios show consistency with values from constrained receiver functions, active source studies, and tomography studies around Eurasia. Our observed Pg and Sg/Lg velocities are consistent with average upper crustal velocities in apriori models such as Crust2.0 and the DOE Unified Model
S23A-1860
A 3D Tomographic Model of Asia Based on Pn and P Travel Times from GT Events
Increasingly, nuclear explosion monitoring is focusing on detection, location, and identification of small events recorded at regional distances. Because Earth structure is highly variable on regional scales, locating events accurately at these distances requires the use of region-specific models to provide accurate travel times. Improved results have been achieved with composites of 1D models and with approximate 3D models with simplified upper mantle structures, but both approaches introduce non-physical boundaries that are problematic for operational monitoring use. Ultimately, what is needed is a true, seamless 3D model of the Earth. Towards that goal, we have developed a 3D tomographic model of the P velocity of the crust and mantle for the Asian continent. Our model is derived by an iterative least squares travel time inversion of more than one million Pn and teleseismic P picks from some 35,000 events recorded at 4,000+ stations. We invert for P velocities from the top of the crust to the core mantle boundary, along with source and receiver static time terms to account for the effects of event mislocation and unaccounted for fine-scale structure near the receiver. Because large portions of the model are under-constrained, we apply spatially varying damping, which constrains the inversion to update the starting model only where good data coverage is available. Our starting crustal model is taken from the a priori crust and upper mantle model of Asia developed through National Nuclear Security Administration laboratory collaboration, which is based on various global and regional studies, and we substantially increase the damping in the crust to discourage changes from this model. Our starting mantle model is AK135. To simplify the inversion, we fix the depths of the major mantle discontinuities (Moho, 410 km, 660 km). 3D rays are calculated using an implementation of the Um and Thurber ray pseudo-bending approach, with full enforcement of Snell's Law in 3D at the major discontinuities. Due to the highly non-linear nature of our ray tracer, we are forced to substantially damp the inversion in order to converge on a reasonable model. We apply both horizontal and vertical regularization to produce smooth models with velocity feature scale lengths that are consistent with established conventions for mantle velocity structure. To investigate the importance of using true 3D rays for the inversion, as opposed to proxy rays through a reference model, we compare our model and ray paths with the model and ray paths resulting from inverting the same data set using rays traced through a 1D reference model. Finally, we validate the model by performing several inversions with random portions of the data set omitted and then testing the predictive capability of the model against those portions compared with AK135. We test the location performance of the model by relocating the GT events using our model and using AK135. Sandia is a multiprogram laboratory operated by Sandia Corporation, a Lockheed Martin Company, for the United States Department of Energy's National Nuclear Security Administration under Contract DE-AC04- 94AL85000.
S23A-1861
Geometry of plate boundary around Japan-Kuril Trench junction derived from 3D seismic tomography by using LT OBS network
We obtained precise hypocenter distribution and 3-D seismic velocity structure around the junction of Japan and Kuril Trenches by a seismic tomographic inversion technique using data from Long-Term OBS (LTOBS) network. We conducted two observations; the first had a period of eight months using 17 LOBSs from Mar. 2005 and the second used 42 LTOBSs from Nov. 2007 to June 2008. We picked P and S arrivals from the records of 59 LTOBSs and 74 land stations for approximately 2000 events. We constructed an initial velocity structure model from 2-D P-wave velocity structure estimated by the seismic refraction experiments (Nakahigashi et al., 2008). Initial Vp/Vs is set to 1.73. The hypocenters form a landward dipping plane. The planar distribution can be interpreted to be associated with the subduction in the study region. And the subducting oceanic plate was imaged as a low velocity region (7 km/s) in the mantle wedge with P-wave velocity of 8 km/s. Because a large velocity difference between the subducting plate and the overriding landward plate is not detected at depths less than 30 km, we determine positions of the plate boundary using the planar hypocenter distribution. High velocity anomaly is found in and around the rupture area of the 1968 Tokachi-Oki earthquake. There is a possibility that the asperity in this region is characterized by high velocity anomaly and small microseismicity near the plate boundary.
S23A-1862
Solution Methods for 3D Tomographic Inversion Using A Highly Non-Linear Ray Tracer
To develop 3D velocity models to improve nuclear explosion monitoring capability, we have developed a 3D tomographic modeling system that traces rays using an implementation of the Um and Thurber ray pseudo- bending approach, with full enforcement of Snell's Law in 3D at the major discontinuities. Due to the highly non-linear nature of the ray tracer, however, we are forced to substantially damp the inversion in order to converge on a reasonable model. Unfortunately the amount of damping is not known a priori and can significantly extend the number of calls of the computationally expensive ray-tracer and the least squares matrix solver. If the damping term is too small the solution step-size produces either an un-realistic model velocity change or places the solution in or near a local minimum from which extrication is nearly impossible. If the damping term is too large, convergence can be very slow or premature convergence can occur. Standard approaches involve running inversions with a suite of damping parameters to find the best model. A better solution methodology is to take advantage of existing non-linear solution techniques such as Levenberg-Marquardt (LM) or quasi-newton iterative solvers. In particular, the LM algorithm was specifically designed to find the minimum of a multi-variate function that is expressed as the sum of squares of non-linear real-valued functions. It has become a standard technique for solving non-linear least squared problems, and is widely adopted in a broad spectrum of disciplines, including the geosciences. At each iteration, the LM approach dynamically varies the level of damping to optimize convergence. When the current estimate of the solution is far from the ultimate solution LM behaves as a steepest decent method, but transitions to Gauss- Newton behavior, with near quadratic convergence, as the estimate approaches the final solution. We show typical linear solution techniques and how they can lead to local minima if the damping is set too low. We also describe the LM technique and show how it automatically determines the appropriate damping factor as it iteratively converges on the best solution. Sandia is a multiprogram laboratory operated by Sandia Corporation, a Lockheed Martin Company, for the United States Department of Energy's National Nuclear Security Administration under Contract DE-AC04- 94AL85000.
S23A-1863
Multi-resolution global tomography using multiple tessellation tiers
We are currently developing a new 3-D global tomography modeling framework with monitoring purposes at the forefront. The overall monitoring goal is to develop a seamless, 3-D global model capable of accurately locating seismic events from combined teleseismic and regional travel time observations. Such a model will eliminate the need to splice independent regional and/or global models for the purpose of travel time prediction since the singular model will be complete and self-consistent. In order to optimize the number of free parameters, the modeling framework must allow for flexible multiple resolution capabilities (depth and spatial). The model design must also directly account for variable discontinuous structures from a-priori information to maintain travel time prediction accuracy. With these general criteria in mind, we have designed a tessellation procedure whereby triangular sides of a pre-determined polyhedron are recursively subdivided into smaller triangles. Each recursion level (tier) produces 4 triangular regions (daughters) within each larger triangular region (parents) and the vertices are then normalized to the unit sphere. Upon building each higher tier tessellation, parent-daughter triangle information is stored and indexed to allow for rapid model referencing at any desired tessellation tier (resolution level). Therefore sensitivity kernels can be computed for all available tessellation tiers while simultaneously triangulating positions along a given ray path, allowing for easily adjustable multi-resolution model inversion setup. Mixing kernel weights from multiple tessellation tiers also provides for an efficient regularization scheme. Sets of model nodes are defined along the tessellated vertices and discontinuous structures (such as the Moho and subducted slabs) are treated via logical definitions and structured data array design with only a slight hindrance on computational efficiency. A geodetic reference ellipsoid (GRS80) is built into the model and 3-D ray tracing is performed to compute sensitivity kernels and travel times leaving only minimal travel time corrections to be applied. Preliminary modeling results using teleseismic P and a collection of Pn regional times from Eurasia will be demonstrated. This work performed under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC52-07NA27344. LLNL-ABS-406714
S23A-1864
Joint Inversion for 3-dimensional P-, S-velocity structures and radial anisotropy in the mantle along the Tethyan margin
We construct a new 3-D S-velocity model by jointly inverting regional S and Rayleigh waveform fits, teleseismic S and SKS arrival times, fundamental-mode Rayleigh-wave group velocities, and Moho depth estimates from receiver functions, refraction/reflection lines, and gravity surveys for the region which extends from the mid-Atlantic ridge to the Hindu Kush. The joint inversion benefits from both better resolution and wider data coverage than when using only individual data sets, since regional S and Rayleigh waveforms and Rayleigh-wave group velocity data have good resolution for the upper mantle while teleseismic S and SKS arrival time data do for the lower mantle. Because S and SKS waves are also sensitive to some integral measure of upper mantle structure, combining these complementary data sets leads to a better separation of upper and lower mantle structure than possible with teleseismic body waves alone. The independent constraints on Moho depth are incorporated to better separate crustal structure from mantle structure than possible with integral data alone. Resolution tests confirm that the joint inversion yields good resolution ranging from the Moho to 1400 km depth. We convert the 3-D S-velocity model to a 3-D P-velocity model using teleseismic P arrival times. We also estimate radial anisotropy in the uppermost mantle by comparing constraints from Love waveform fits and Love-wave group velocities to the estimated 3-D S-velocity model. We discuss features of our new models, which includes oceanic structure, cratons, subducting slabs, low- velocity mantle plumes, rifts, and basins and characteristics of radial anisotropy in the uppermost mantle.
S23A-1865
Joint-inversion of ZH ratio and Phase Velocity for Crustal S-wave Velocity Structure in LA Basin
Improving our understanding of shallow S-wave velocity structure is important for accurate ground motion prediction. We attempt to improve S-wave crustal structures in the urban area of Los Angeles basin and its surroundings. We use a joint-inversion of two kinds of data sets; the ratio between vertical and horizontal eigenfunctions (ZH ratio) and phase velocities of Rayleigh waves. We found that SCEC Community Velocity Model (CVM-H) is systematically different from S-wave velocity structure derived from ZH ratios and phase velocities. Our inversion results showed that most velocity structures need slower S-wave velocity at shallow depth in order to fit ZH ratio and phase velocity data. This discrepancy cannot be ignored since this is expected to cause differences in ground motions by a factor of few, although this factor changes from frequency to frequency. For ZH ratio data, we take advantage of the fact that stable ZH ratios can be obtained between 0.13-0.37 Hz. This relatively long frequency range enables us to observe different patterns depending on local structure. Within this frequency range, ZH ratios show remarkably sharp sensitivities at upper 5km. As for phase velocity measurement, we use Green's functions obtained from seismic noise by cross- correlation. This enables us to find phase velocity on paths between every possible pair of stations in southern California. We took advantage of having seismic noise with a high S/N between 0.13-0.2 Hz. For this range of frequency, phase velocity is sensitive down to 20 km. Jointly inverting ZH ratios and phase velocities allows us to resolve the entire upper crust while retaining sensitivities at very shallow depth. This method was applied to LA basin area and surroundings with good variance reduction for fitting data. We will show that our inversion results for the LA basin area require slower S-wave structures than CVM-H at shallow depth.
S23A-1866
Crustal Structure Of Southern Africa
In this study, receiver functions and Rayleigh wave group velocities have been jointly inverted to model subsurface structure across southern Africa. The broadband data used for this study come from the Southern African Seismic Experiment (SASE), which consisted of 82 stations deployed between 1997 and 1999, the Kimberly Broadband Array, which was also part of the SASE, GSN/IMS stations, and permanent AfricaArray stations. In all, data from 103 stations were modeled, and from these, good quality results were obtained for 91 stations. The group velocities used for the inversion were taken from the group velocity maps published by Pasyanos (2007), and the inversion method used was taken from Julià et al. (2000, 2003). The shear wave velocity (Vs) results across the geological terrains of southern Africa show several types of anomalous crustal structure (1) Lower than average Vs (< 4.0 km/s) characterizes the lower crust in the central part of the Kimberley terrain. (2) Higher than average Vs (4.2 km/s) characterizes the lower crust in the Namaqua-Natal and Kheiss terrains. (3) Higher than average Vs (Vs > 3.6 km/s) characterizes the upper crust (~ 10 km thick) in the Namaqua-Natal terrain. The low Vs in the lower crust of the Kimberley terrain could be due to thinning of the crust from extension during basin forming events that affected the western part of the Kaapvaal Craton between c. 3.0 Ga and 1.8 Ga. The high Vs in the lower crust of some terrains and the high Vs in the upper crust of the Namaqua-Natal terrain could be attributable to mafic material preserved during suturing events.
S23A-1867
Shear Wave Velocity Structure Of The Bushveld Complex
In this study, upper crustal structure of the Bushveld Complex, which is an ultramafic intrusion located within the northern part of the Kaapvaal craton, is investigated using high frequency Rayleigh wave tomography. Past geophysical studies (e.g. Webb et al., 2004) show that the mafic rocks of the western and eastern limbs of the Bushveld Complex may be continuous and connected at depth. The objective of this study is to derive velocity images (Rayleigh wave group and shear wave velocities) of the shallow crustal structure of the Bushveld Complex (BC) to depths of ~ 15 km and evaluate if the mafic units connect at depth beneath the central part of the complex. Broadband seismic data from mining-induced (mostly from gold mines around the Wits basin) and regional earthquakes (ML ~ 2.5) around the BC have been used to measure event-to-station Rayleigh wave group velocities in the period range of 1- 15 sec. About 46 stations from the SASE (South African Seismic Experiment, 1997 - 1999) and global digital broadband network are used to map the crustal structure of the BC. Preliminary inversions of the group velocity measurements at 150 x 150 km grid are consistent with the presence of high velocity crustal rocks beneath the western and eastern limbs of the BC at depths of ~ 6 - 8 km. Additional semsic data especially from regional earthquakes are being analyzed to improve resolution of our models for the central region of the Bushveld Complex.
S23A-1868
Simultaneous Joint Inversion of Multiple Geophysical Data Sets and 3D Tomography
We present the results of our efforts to model the 3D seismic shear wave velocity structure by simultaneously and jointly inverting multiple geophysical observations. It has long been recognized that derivation of geophysical models for a given observable, such as body wave models or surface wave models, may provide good predictors of behavior for the specific parameter modeled, but often provides poor prediction capabilities for other parameters, even when we have empirical or theoretical knowledge of how the parameters relate to one another. Due to the non-unique properties of inversion methods, often we may find a solution for one data type but we must acknowledge that, although it can predict behavior of that data type, the solution is ambiguous and other models may serve equally well. Thus a joint, simultaneous inversion of all pertinent and available data types for a region may provide us with a far more robust model than obtained by using only a single parameter. Moreover, the extent of our well-resolved model may be significantly expanded due to the variable coverage provided by different measurements, and its features may be better constrained due to the different resolving powers of different data types. We are exploring the 3D seismic structure of a small portion of the Tarim Basin in northwest China using simultaneous joint inversion of surface wave group velocities, teleseismic P-wave receiver functions, S-wave travel times, and Bouguer gravity anomaly observations derived from the Gravity Recovery And Climate Experiment (GRACE) mission. An iterative, conjugate gradient-based least squares inversion is used to jointly model the four different data sets, using shear-velocity variations as the primary model parameters. Improved knowledge of the shear velocity structure of the Asian continent is of fundamental importance for understanding the geodynamic evolution and formation of continents and the processes acting within and on the continental lithosphere.
S23A-1869
An estimation technique of Rayleigh wave phase velocities using arrays with arbitrary geometry
The mictotremor survey method (MSM) is one of the most practical techniques to estimate velocity structure of shear waves in sedimentary layers. In the MSM, the velocity models are determined by inversion analysis of the Rayleigh wave phase velocity dispersion curve observed from microtremors. In most of the cases, the phase velocity dispersion curve is obtained by either the spatial autocorrelation (SPAC) technique or the frequency-wavenumber (F-K) technique applied to array measurements of microtremors. These techniques place significant restrictions on the array geometry and number of stations required, which limits the applicability of MSM, especially in urban areas. We have derived a new technique for estimating phase velocities of Rayleigh waves. This new technique (the direct estimation method: DEM) enables to the use of flexible array configurations and a minimal number of stations. Moreover, the DEM can be applied to records from existing station arrays, such as those in an earthquake monitoring network. In the DEM, microtremors detected by arrays with arbitrary geometry can be represented by complex coherence functions (CCFs: Shiraishi et. al. 2006) of the Rayleigh wave. The CCF is derived from analytic solution of Lamb's problem, and it consists of the Bessel function of the first kind J0(ωr/c) (ω: angular frequency, r: distance between the stations, c: phase velocity), which is well-known function and is used in the SPAC technique to estimate phase velocity. The phase velocities can be estimated by solving the equations with the least squares approach to minimize the residual error between the observed and the theoretical values. A field experiment has been carried out to verify the effectiveness of the DEM, and the phase velocities obtained by the DEM with an array of arbitrary geometry are in excellent agreement with those obtained using the SPAC technique.
S23A-1870
S-Wave Velocity Across Central Mexico Using High Resolution Surface Wave Tomography
The shear wave velocity structure across central Mexico is determined by surface wave dispersion from a dense linear seismic experiment "Mesoamerican Subduction Experiment" (MASE). MASE consisted of 100 portable broadband stations deployed along a line crossing Central Mexico from the Pacific Coast to almost the Gulf of Mexico. Regional records were used to obtain Rayleigh-wave group velocity maps for periods from 5 to 50 s and they show a dramatic variation of velocity (~40%), especially for periods larger of 20 s. Local dispersion curves were reconstructed for each station and inverted to find S-wave velocity by using a simulated annealing algorithm. The results, from inversion, show a significant change, particularly in the lower crust, between the backarc, volcanic arc and forearc regions. The crust in the forearc is thicker and faster than the backarc region. Just below the active Trans Mexican Volcanic Belt (TMVB) (300 km from the coast) is presently a low velocity spot (~3.4 km/s) suggesting presence of anomalous material (probably related to a mantle wedge) as deep as 50 km. The results also show a poorly resolved slab and wedge which correspond to the ones in a model reported recently. The results are supported with consistency checks and resolution tests.
S23A-1871
Crustal and Mantle Structure of the Jalisco Block of western Mexico from Surface Wave Tomography
How a subduction system evolves to a transform system is a fundamental question in plate tectonics and is still not well understood. It is believed that fragmentation of both the subducting and overriding plates is one of the key steps of the evolution. Such a process is occurring in the Rivera subduction zone of the western Mexico. In this region, the Rivera plate detached from the Cocos plate a few million years ago and is presently subducting beneath the Jalisco block, which is separating from the North American plate along the Colima graben and Tepic-Zacoalco Rift to the east and north, respectively. Seismic imaging of the velocity structure and boundaries of the Rivera-Cocos-Jalisco microplate system is thus essential to understanding the regional tectonics as well as the broad question of subduction-to-transform transition. Fifty broadband seismic stations were temporarily deployed between January 2006 and June 2007 to explore the geodynamic processes of this plate boundary region under the MARS (MApping the Rivera Subduction zone) project. Here we present preliminary results on the crustal and upper mantle structure of the Jalisco block based on the analysis of surface waves recorded by the seismic array. Phase velocities at periods from 5 to 35 s were measured by cross correlating continuous records of two-station pairs. High SNR (signal-to-noise ratio) Green's functions were obtained along 474 ray paths and were inverted to generate the phase velocity map of the region. There is a good agreement between the phase velocity maps and the regional tectonic structure. The strongest low velocity anomaly is located beneath the Colima volcano system and is shown in all period bands. The NS trending Colima graben and the NW trending Tepic-Zacoalco Rift are also featured by distinct low velocity lineaments. We also analyzed the fundamental mode Rayleigh waves of 116 teleseismic Mw>6 events in the 20-100s period bands using the two-plane wave method. The resulting phase velocity maps were combined with those derived from ambient noise data to yield dispersion curves covering a period range of 5 to 100 s. These dispersion measurements will be used to estimate the S-wave velocity structure of the crust and upper mantle to ~200 km depth beneath the Jalisco block.
S23A-1872
Spatial Variation in Crustal Thickness in Northern Italy from Surface Waves
The Italian peninsula has long been regarded as a prime example of slab rollback and retreat, extending the overriding Eurasian plate to the southeast beside the retreating Adria microplate. However, the situation in northern Italy appears much more complicated, as the continental Adria microplate lacks the deep intraslab earthquakes typically seen at other subduction zones. Also, Tuscany and the Northern Apennines have little volcanism or history of large subduction thrust earthquakes, as compared with regions further south along the peninsula. Deciphering the tectonics of the region has great implication for the seismic hazard of northern Italy, however, there are competing and very different models for the crustal structure and thickness of the region making it difficult to develop a tectonic model. Here we use surface waves gathered during the RETREAT experiment to examine geographical variations in crustal thickness over the northern part of the Italian peninsula. The dispersion curves are determined using a two-station approach, using only station pairs and events with theoretical ray paths within 30 degrees of the great circle path between the two stations. Individual phase velocity measurements are computed using a multi-wavelet analysis. Dispersion curves are checked for quality control and then stacked based on station pair and the midpoint between the two stations. From the dispersion curves alone it is clear that the Tuscan crust is much thinner than that of the Apennines, which we would expect for the subduction case. This thinner crust is in agreement with recently proposed models from receiver functions generated by other researchers in the RETREAT group. Those geographical bins with sufficient data to generate relatively smooth dispersion curves are then inverted for crustal thickness using Bayesian techniques and compared with other published models.
S23A-1873
Structure of the crust and the uppermost mantle beneath the Kamchatka peninsula and surrounding basins from inversion of surface wave group velocities
We used records of the 2006 Olyutorskoe (Korakia) earthquake (Mw=7.6) and its aftershocks by broadband stations of the Kamchatka seismological network and GSN stations to measure dispersion curves for various paths crossing the Kamchatka peninsula, the Komandor Basin, and the Sea of Okhotsk. Group velocities were measured at periods between 7 and 70 s by applying frequency-time analysis (Levshin et al., 1989). We inverted the obtained measurements for path-average 1D shear-velocity profiles with a Monte-Carlo method. Our results indicate that the crust is about 35 km thick beneath Kamchatka and significantly thinner beneath the Sea of Oknotsk (20-25 km) and the Komandor Basin (15-20 km). We also report a clear difference in seismic velocities in the uppermost mantle between the northern and the southern Kamchatka. The uppermost mantle beneath the east coast south of the Cape Kamchatka is characterized by relatively high seismic velocities due to the presence of the subducted Pacific plate. At the same time, the uppermost mantle is relatively slow beneath the northern Kamchatka, excluding a presence of currently subducting oceanic lithosphere north of the Kamchatka-Aleutian junction.
S23A-1874
Mapping the Moho Beneath Northwest Canada: Crustal Structure Across CANOE
The Canadian Northwest Experiment (CANOE) traverses a wide variety of continental settings, allowing the study of continental assembly over a time span of nearly 4 Ga. CANOE consists of nearly sixty broadband seismometers with branches running approximately east-west from the Slave Craton across a series of Proterozoic orogens to the Canadian Cordillera and from Edmonton, Alberta across the Churchill Province to the Rocky Mountain Front. A third section runs roughly north-south parallel to the Rocky Mountain Front. A pilot group of twelve stations were installed in May 2003, and the remaining stations were operational beginning in May 2004. All stations were removed in September 2005. In this study, we obtain P-s differential times from Moho conversions using teleseismic receiver function analysis. In general, the preliminary results show larger differential times (on the order of 5 seconds) beneath stations on the craton and in the Churchill Province, at the eastern portion of the study area, and smaller differential times (on the order of 3-4 seconds) along the Rocky Mountain Front and beneath the Canadian Cordillera in the west. Crustal thickness is calculated using the differential times from Moho conversions and a local velocity model resulting from ambient-noise cross-correlation tomography beneath the array. The resulting crustal thicknesses range from approximately 25-45 km. In addition to direct P-s Moho conversions, we examine the receiver functions for arrivals from Moho multiples, and mid-crustal, lithospheric, and lithosphere- asthenosphere conversions.
S23A-1875
Mesh Generation for Short-Period Seismic Wave Propagation Based Upon the Spectral- Element Method: Southern California.
The Spectral Element Method (SEM) has been successfully applied to simulate ground motion in Southern California and the Los Angeles Basin for period up to 2 sec. Nowadays, simulations at shorter period are computationally feasible, but they require both a realistic geological model and a detailed unstructured hexahedral mesh. Aiming to include the effect on the seismic propagation due to subsurface geology, topography and low- velocity sedimentary basins, we have generated a 3D unstructured hexahedral mesh of Southern California suitable for shorter period simulations. The grid honors an updated description of the Salton Sea and of the sedimentary basins of Los Angeles, San Fernando and Ventura. We discuss some criteria to determine the geological details that need to be honored, analyzing 2D cross sections of the region for simulations at 10 Hz. We have focused in particular on a profile crossing the Santa Monica Mountains, investigating a possible lensing effect due to a large overthrust. The grid is nearly automatically generated on a massively parallel machine by GEOCUBIT, a Python script collection based upon CUBIT (Sandia Laboratory, www.cubit.sandia.gov). CUBIT is an advanced 3D unstructured hexahedral mesh generator that offers tremendous opportunities in assessing the quality of a mesh both in terms of geometrical complexity and numerical accuracy. For 3D simulations, we have applied SPECFEM3D, which accommodates anisotropy attenuation, free surface topography, fluid-solid boundaries and absorbing boundary conditions. For the 2D simulations we use SPECFEM2D, developed by Roland Martin, Dimitri Komatitsch, Céline Blitz and Nicolas Le Goff (2008).
S23A-1876
Tera3D: A Tera-Scale Full-3D Waveform Tomography (F3DT) in Southern California
We are automating our full-3D waveform tomography (F3DT) and seismic source inversion algorithm based on the scattering-integral (SI) method and applying the automated algorithm to iteratively improve the 3D SCEC Community Velocity Model Version 4.0 (CVM4) in Southern California. In F3DT, the starting model as well as the derived model perturbation is 3D in space and the sensitivity kernels are calculated using the full physics of 3D wave propagation. The SI implementation of F3DT is based on explicitly constructing and storing the sensitivity kernels for individual misfit measurements. Compared with other F3DT implementations, the primary advantages of the SI method are its high computational efficiency and the ease to incorporate 3D Earth structural models into real-time seismic source parameter inversions. In a previous study, we have successfully applied the SI method to improve the 3D seismic velocity structure model (SCEC CVM3) and seismic source models in the Los Angeles region. In this presentation we will report our recent progresses on automating the complete F3DT workflow and extending our analysis to a much larger region in Southern California.
S23A-1877
Using Fréchet kernels to Investigate 3D Wave Excitation and Propagation in Southern California
In this study we use 3D Fréchet kernels to gain insight into wave propagation in Southern California. Three- dimensional structural complexity in Southern California gives rise to waveforms whose properties are unaccounted for by simple 1D models. Fréchet kernels, or sensitivity kernels, are a useful tool which can be used to analyze these phases because kernels can be calculated for any waveform from any part of the seismogram; a priori knowledge of phase type is unnecessary. Based on the 3D SCEC Community Velocity Model version 4.0 (CVM-S4), we compute Fréchet kernels by a scattering-integral method, which involves convolving point-source generated earthquake wavefields with the receiver Green tensors (RGTs), produced by three orthogonal unit point forces acting at the receiver locations. In this presentation we provide examples of sensitivity kernels in both phase and amplitude for waveforms from magnitude 3.5-5.5 earthquakes and discuss the insight they provide into the complexities of 3D wave propagation in Southern California.
S23A-1878
A Study on Efficient Procedure of Non-linear Waveform Inversion --- Applied to Determine The Velocity Structure of The Taipei Basin
The main purpose of this research is to perform a nonlinear waveform inversion to investigate the velocity structures of a basin. We have built up an efficient procedure which incorporates with parallel computing to accomplish the above goal. First, we adopted the nonlinear search algorithm to explore the parameter space and construct a velocity model to calculate the theoretical travel-time for each station. Second, we obtain acceptable models which have minimum travel-time error. Finally, these models are used to waveform simulation by adopting the numerical calculation. By minimizing the misfit between synthetic and observed seismograms, through the nonlinear search algorithm, we are able to obtain more precise velocity structures. By combining all these methods with parallel computer clusters, it allows us to perform the waveform inversion with less time consuming. In order to exam the performance of the procedure described previously, we applied it to investigate the velocity structure of the Taipei basin. Our results indicate that the average travel time error is less than the original data. More importantly, the highly agreement between the observed and the synthetic waveforms suggest the reliability of the velocity structure. The results are not only giving better understanding the complex structure beneath the Taipei basin , but also can allow us to identify waveform characteristics caused by velocity heterogeneity and topographic effect. Furthermore, the geometry of the Taipei basin can be clearly seen from the velocity profiles, such as the depth of the basin increases from east to west, where the deepest part (depth=680m) is near the west edge of the basin. By comparing with the velocity structures obtained from seismic reflection method, we believe that our approach can obtain more precise S-wave velocity structure of the Taipei basin. Thus, the velocity structures obtained from this study can serve as an initial model to invert more detailed structures of the Taipei basin for future research.
S23A-1879
High order non-periodic homogenization for wave propagation in complex 1D and 2D media.
In many cases, in the seismic wave propagation modeling context, scales much smaller than the minimum wavelength are present in the earth model we wish to propagate in. For many numerical methods, like the spectral element method, these small scales need to be meshed in order to be accurate, which is a challenge leading to high numerical cost. The purpose of the work presented here is to understand and to build the effective medium and equations allowing to average the small scales of the original medium without losing the accuracy of the wavefield computation. For this kind of problems, two scale homogenization in periodic and random media is a well-known effective method. Nevertheless, applications of two scale homogenization to non-periodic deterministic media is still a challenge. We first develop a method similar to periodic homogenization for non-periodic media in the layered media case. The order 0 homogenization gives the result that was obtained by Backus in 1962 which implies that order 0 homogenized model is transversely isotropic even though the original model is isotropic. It also appears the order 0 is not enough to obtain surface waves with correct group and phase velocities and that higher order homogenization terms up to 2 are often required. High order homogenization implies to modify the wave equation and the boundary conditions. We then extend this method to media with rapid variation in two directions. Different numerical tests showing the accuracy and convergence of such a method will be presented.
S23A-1880
Nonstandard FDTD for Accurate Modeling of Seismic Wave Propagation in 2D
Finite-difference method in time-domain (FDTD) is one of the most common techniques used for modeling of seismic wave propagation. The algorithm is popular because it is simple and easy to program. In the FDTD, the numerical solutions do not coincide with the theoretical solutions unless the temporal and spatial discretization are sufficiently fine due to the numerical dispersion and grid anisotropy from the FDTD schemes. In this study, we develop a FDTD scheme called nonstandard FDTD (NS-FDTD) for 2D elastic (P- SV) wave computations, which was originally proposed in computational electromagnetics (e.g. Cole, 1997, IEEE Trans. MTT). We implement the NS-FDTD through the following two steps: we first replace the spatial and temporal grid spacings by their frequency optimized counterparts called the correction functions, and we then introduce a finite-difference form of the Laplacian. The nonstandard scheme efficiently reduces numerical dispersion and grid anisotropy to improve the computational accuracy. The high accurate nonstandard versions of the FDTD algorithms are only slightly more complicated than the standard ones, so that existing computer programs could be easily modified to run the nonstandard ones.
S23A-1881
Seismic Modeling Of Reservoir Heterogeneity Scales: An Application To Gas Hydrate Reservoirs
Natural gas hydrates, a type of inclusion compound or clathrate, are composed of gas molecules trapped within a cage of water molecules. The occurrence of gas hydrates in permafrost regions has been confirmed by core samples recovered from the Mallik gas hydrate research wells located within Mackenzie Delta in Northwest Territories of Canada. Strong vertical variations of compressional and shear sonic velocities and weak surface seismic expressions of gas hydrates indicate that lithological heterogeneities control the distribution of hydrates. Seismic scattering studies predict that typical scales and strong physical contrasts due to gas hydrate concentration will generate strong forward scattering, leaving only weak energy captured by surface receivers. In order to understand the distribution of hydrates and the seismic scattering effects, an algorithm was developed to construct heterogeneous petrophysical reservoir models. The algorithm was based on well logs showing power law features and Gaussian or Non-Gaussian probability density distribution, and was designed to honor the whole statistical features of well logs such as the characteristic scales and the correlation among rock parameters. Multi-dimensional and multi-variable heterogeneous models representing the same statistical properties were constructed and applied to the heterogeneity analysis of gas hydrate reservoirs. The petrophysical models provide the platform to estimate rock physics properties as well as to study the impact of seismic scattering, wave mode conversion, and their integration on wave behavior in heterogeneous reservoirs. Using the Biot-Gassmann theory, the statistical parameters obtained from Mallik 5L-38, and the correlation length estimated from acoustic impedance inversion, gas hydrate volume fraction in Mallik area was estimated to be 1.8%, approximately 2x108 m3 natural gas stored in a hydrate bearing interval within 0.25 km2 lateral extension and between 889 m and 1115 m depth. With parallel 3-D viscoelastic Finite Difference (FD) software, we conducted a 3D numerical experiment of near offset Vertical Seismic Profile. The synthetic results implied that the strong attenuation observed in the field data might be caused by the scattering.
S23A-1882
Envelope Synthesis on the Free Surface of a Random Elastic Medium Based on the Markov Approximation
Short-period seismograms mostly consist of the waves scattered by random inhomogeneities of the solid earth. For P-waves, the apparent duration is broadened and the transverse amplitude is excited with travel distance increasing. Usually, we observe the seismic waves on the free surface. For the homogeneous medium case, the vertical incident P-wave amplitude is doubled at the free surface; however, for the inhomogeneous medium case it has not been clear how the free surface affects wave amplitudes since the ray directions are widely distributed because of scattering. Here, we develop the synthesis of vector-wave envelopes on the free surface of a random medium. When the wave length is shorter than the correlation distance of 3-D infinite random media, characterized by a Gaussian autocorrelation function, Sato (2006) derives analytical solutions of vector-wave envelopes on the basis of the Markov approximation, which is a stochastic extension of the split stem method to solve the parabolic wave equation. We extend his method for the synthesis of vector-wave envelopes on the free surface of a random medium. In the Markov approximation, Mean square (MS) envelopes are calculated by using the Fourier transform of two-frequency mutual coherence function (TFMCF) with respect to the angular frequency. The TFMCF is statistically defined on the transverse plane, which is perpendicular to the global ray direction. The Fourier transform of the TFMCF with respect to the transverse coordinates gives the angular spectrum that shows the distribution of ray directions. This angular spectrum has a sharp peak in the global ray direction just after the direct wave arrival and gradually increases its width with the lapse time increasing. We can calculate vector-wave envelopes in the infinite space, simply projecting the angular spectrum to each component and integrating it in the wavenumber space. We calculate the vector-wave envelopes on the free surface by multiplying the amplification factor of the free surface to the angular spectrum instead. We numerically simulate vector-wave envelopes on the free surface of random elastic media with 100km thickness for the vertical incidence of a plane P-wavelet where the MS envelope of the incident wavelet is a delta function. Averaged P-wave and S-wave velocities are 6 km/s and 3.46 km/s, respectively. MS fractional fluctuation of the velocity is 5% and the correlation distance is 5 km. The MS envelopes on the free surface are not simply 4 times larger than those in the corresponding infinite media. The peak amplitude of the horizontal component of the MS envelope on the free surface is 4.8 times larger than that of the MS envelope in the infinite media and in the latter part of the envelope, the amplification rate is less than 4, but the vertical component MS envelope on the free surface is nearly 4 times of that in the infinite media. The peak delay time of the horizontal component on the free surface is about 0.1 s earlier than that in the infinite media.
S23A-1883
Crustal Velocity Structure of the Rwenzori Region, Uganda, From Isotropic and Anisotropic Travel-Time Tomography
The Rwenzori region in Uganda forms part of the western branch of the East-African rift system. During the period from May 2006 to October 2007, a temporary seismological network was operating in the area to constrain the seismic structure of the crust and upper mantle from recordings of local and teleseismic earthquakes. Here we use local earthquake tomography to resolve the 3D velocity structure of the crust down to a depth of about 20 km. Results based on P and S-wave arrival times exhibit a pronounced negative velocity anomaly near the western flank of the Rwenzori Mts. in the upper 5 km. At the surface, this area is characterized by geothermal activity in relation to the Buranga hot springs. As depth increases the negative anomaly gradually widens and shifts to the East. The results are supported by sensitivity tests to constrain the resolving power of the data. We also made an attempt to derive the 3D anisotropic velocity structure in the region. For the anisotropic inversion, we assume that the elastic properties of the crust can be characterized by a simplified form of transverse isotropy, which can be defined by four parameters: a fast and slow velocity and two angles to determine the orientation of the fast velocity axis. Average isotropic velocity variations obtained from the 3D anisotropic model are in good agreement with the results from the purely isotropic inversion. We find that the fast axes in the northern region of the Rwenzori Mts. are dominantly oriented NS, approximately parallel to the strike of mayor faults. Several tests are performed to constrain the uncertainties of the model.
S23A-1884
Seismicity of the Rwenzori Region, Western Uganda
The 5000m high Rwenzori Mountains are situated within the western branch of the East African Rift System (EARS), close to the equator at the border between Uganda and the Democratic Republic of Congo. They represent a basement block within the rift whose origin and relation to the development of the EARS are focus of the RiftLink project (www.riftlink.org). To investigate crustal and upper mantle structure in conjunction with seismic activity on a regional and local scale, a temporary seismic network was deployed over an area of roughly 80 x 140 km and operated from May 2006 to the end of September 2007. The analysis of the registered data revealed high microseismic activity in the region. On average more than 800 events per month were located during the registration period with local magnitudes ranging from -0.5 up to 5.1. Few earthquakes are located within the Rwenzori massive itself. Most of the events occur east and west of the mountains with a pronounced concentration of activity at a depth of about 15 km. Vertical sections across the northern parts of the Rwenzories show, that west of the mountains (towards the rift valley) the focal depths range from 10 to 20 km, whereas the hypocenters go as deep as 30 km on the eastern side. This is in good agreement with Moho depths that were derived from receiver functions and are close to 22 km west and 30 km east of the Rwenzories. There is one exception, however. Approximately 30 km east of the northern mountain ridge, we located a cluster of 7 events exhibiting an anomalous depth of about 60 km that occurred within 20 days in September 2006. These events are unique, up to now we located no other earthquakes at similar depths. P-wave polarities were used to determine fault plane solutions of events that were recorded by an adequate number of stations. Nearly all source mechanisms reveal normal faulting with strike directions more or less parallel to the rift axis and extension forces perpendicular to it. However, there is a group of events whose strike directions seem to be systematically tilted counter-clockwise supporting a numerical model of Koehn et al. (2007) who explain the Rwenzori block as a micro-plate that was captured during the approach of two rift segments and is rotating clockwise. Koehn, D., Aanyu, K., Haines, S., Sachau, T., 2007, Rift nucleation, rift propagation and the creation of basement micro-plates within active rifts, Tectonophysics, in press, doi:10.1016/j.tecto.2007.10.003