Crust and subduction zone structure of Southwestern Mexico

Crust and Subduction Zone Structure of Southwestern Mexico Sandy Kurniawan Suhardja,1 Stephen P. Grand,1 David Wilson,2 Marco Guzman-Speziale,3 Juan Martin Gomez-Gonzalez,3 Tonatiuh Dominguez-Reyes,4 and James Ni,5 1 The University of Texas at Austin, United States 2 US Geological Survey, United States 3 Universidad Nacional Autonoma de Mexico, Mexico 4 Observatorio Vulcanologico, Universidad de Colima, Mexico 5 New Mexico State University, United States This article has been accepted for publication and undergone full peer review but has not been through the copyediting, typesetting, pagination and proofreading process which may lead to differences between this version and the Version of Record. Please cite this article as doi: 10.1002/2014JB011573 © 2015 American Geophysical Union. All rights reserved. Abstract Southwestern Mexico is a region of complex active tectonics with subduction of the young Rivera and Cocos plates to the south and widespread magmatism and rifting in the continental interior. Here we use receiver function analysis on data recorded by a 50 station temporary deployment of seismometers known as the MARS (MApping the Rivera Subduction zone) array to investigate crustal structure as well as the nature of the subduction interface near the coast. The array was deployed in the Mexican states of Jalisco, Colima, and Michoacan. Crustal thickness varies from 20 km near the coast to 42 km in the continental interior. The Rivera Plate has steeper dip than the Cocos plate and is also deeper along the coast than previous estimates have shown. Inland, there is not a correlation between the thickness of the crust and topography indicating that the high topography in northern Jalisco and Michoacan is likely supported by buoyant mantle. High crustal Vp/Vs ratios (greater than 1.82) are found beneath the trenchward edge of magmatism including below the Central Jalisco Volcanic Lineament and the Michoacan-Guanajuato Volcanic Field implying a new arc is forming closer to the trench than the Trans Mexican Volcanic Belt. Elsewhere in the region, crustal Vp/Vs ratios are normal. The subducting Rivera and Cocos plates are marked by a dipping shear wave low velocity layer. We estimate the thickness of the low velocity layer to be 3 to 4 km with an unusually high Vp/Vs ratio of 2.0 to 2.1 and a drop in S velocity of 25%. We postulate that the low velocity zone is the upper oceanic crust with high pore pressures. The low velocity zone ends from 45 to 50 km depth and likely marks the basalt to eclogite transition. © 2015 American Geophysical Union. All rights reserved. Introduction The tectonic evolution of western North America has been strongly influenced by successive fragmentation events of the ancient Farallon plate, as various segments of the East Pacific rise approached the paleo-trench off shore. When a spreading center encounters a subduction zone, a major change in stress and plate boundaries occurs. The change in stress affects the tectonics of the overriding plate which may include unusual volcanism, fragmentation of the overriding plate, and micro- plate capture (Stock and Lee, 1994). However, the details of how this occurs and the controlling dynamics are still poorly understood. Currently, the East Pacific Rise is close to encountering the Southwestern Mexico subduction zone near the Colima Rift (Figure 1). In this area the Rivera plate detached from the Cocos plate since 10 Myr (Demets and Traylen, 2000) and is subducting beneath Jalisco, Mexico. The exact boundary between the Cocos and Rivera plates is uncertain although it is likely associated with the El Gordo Graben offshore from the Colima Rift (Bandy et al., 1995). The oceanic plates subducting have ages between 10 and 13 Myr except for a small section near the El Gordo Graben that may be as young as 4 Myr (Lawver et al., 2013). Associated with the subduction of the Rivera Plate is the overriding Jalisco Block. The Jalisco Block is bounded by the Colima Rift on the east and the TepicZacoalco Rift to the north. Luhr et al. (1985) proposed that the Jalisco Block was produced by an eastward jump of the East Pacific Rise and that the Colima Rift may evolve to an ocean spreading center in the future with the Jalisco Block captured by the Pacific Plate. This model predicts right lateral strike slip motion along the TepicZacaolco Rift. Rosas-Elguera et al. (1996), however, find evidence of extension along the Tepic-Zacaolco Rift and postulate that the relative motion of the Jalisco Block © 2015 American Geophysical Union. All rights reserved. with respect to the North America plate is due to applied stresses at the plate boundaries as well as the differing motions of the Rivera and Cocos plates respectively. Widespread magmatism is also present in Southwestern Mexico. Within the Tepic-Zacaolco Rift are several stratovolcanoes that make up the Trans-Mexican Volcanic Belt (Figure 1). An exception to this trend is Colima Volcano that lies closer to the trench. Magmatism has been migrating trenchward during the past 10 Ma (Gomez-Tuena at al., 2007; Ferrari et al. 2012) such that the current front of the magmatic arc is 60 to 80 km south of the Tepic-Zacaolco Rift. In the Jalisco Block, the arc front is marked by a line of volcanic fields known as the Central Jalisco Volcanic Lineament (CJVL). Bandy et al. (2001) have dated volcanic rocks along the CJVL and found a northwest trend of younger rocks ranging from 3 Myr near Ayutla to Holocene ages in the Mascota fields (Figure 1). East of the Jalisco Block is the Michoacan-Guanajuato Volcanic Field (MGVF), a recent volcanic field with more than 1000 Quaternary eruptive centers (Gomez-Tuena et al., 2007). In spite of the unique tectonic and magmatic activity within Southwestern Mexico, there has been little detailed work done on the crustal and lithospheric structure of the region. A detailed knowledge of crustal structure provides constraints for the dynamic link between the tectonics of the over-riding plate and the evolution of the down-going Rivera and Cocos slabs. Regional surface wave studies (Gaite et al., 2012) provide estimates of crustal thickness although such studies commonly have strong tradeoffs between Moho depth and velocity. Gravity studies have been the main source of information on crustal thickness in Southwestern Mexico (UrrutiaFucuguachi & Flores-Ruiz; 1996) but also rely on estimates of densities of different rock units to determine crustal thickness. In this study we apply receiver function © 2015 American Geophysical Union. All rights reserved. (RF) analysis to a large data set recorded by two temporary seismic arrays deployed in Southwestern Mexico as well as two stations of the Mexican National Network. We use the H-κ stacking method developed by Zhu and Kanamori (2000) to estimate both crustal thickness (H) as well as the average Vp/Vs ratio (κ) of the crust beneath each station. This approach worked well in the interior of the region where the Moho is relatively flat but did not produce good results near the coast where clear dipping structures are evident. We also use Common Conversion Point (CCP) analysis to improve images near the coast where the dipping oceanic plates produce clear signals in the RF’s. The aim of this study is to determine the crustal structure for both continental and oceanic crust beneath southwestern Mexico Data Receiver functions were computed from data acquired from the Mapping the Rivera Subduction Zone (MARS) array. The MARS array consisted of 50 broad band sensors that were deployed throughout Southwestern Mexico (Figure 1) for 18 months beginning in January, 2006. Stations were spaced from 35 to 50 km apart with a combination of Strekeisen STS-2 and Guralp 3T and ESP sensors. PASSCAL data acquisition systems recorded the data at 40 sps. The MARS experiment involved two U. S. institutions, the University of Texas at Austin and New Mexico State University, in collaboration with two Mexican institutions, Centro de Geociencias, UNAM and the Volcanic Observatory at the Universidad de Colima. We also use data collected by the Colima Volcano Deep Seismic Experiment (CODEX). CODEX consisted of a deployment of 20 short period instruments around Colima Volcano (Figure 1) and overlapped the MARS experiment for five months in time (Gardine et al., 2007). Finally, we also used two broadband seismic stations from the Mexican National Seismic network located within Southwestern Mexico. © 2015 American Geophysical Union. All rights reserved. Method Ideally, receiver functions are time series that show pulses corresponding to converted P to S waves recorded by a particular station. Since the conversion of P to S occurs at sharp boundaries in shear velocity in the subsurface, receiver functions can map out subsurface layering beneath seismic stations. A simple approach to isolate the converted S waves on a seismogram is to deconvolve the vertical from the radial component (Langston, 1977; Ammon, 1991). Data from earthquakes with magnitudes greater than 5.7 and epicentral distances between 30 to 100 degrees were collected from the seismic arrays discussed above. The locations of the hypocenters and the station locations are shown in Figure 2. Most of the good quality data come from Tonga, South America and the Aleutian islands. Seismic traces were cut from 20 seconds before the incoming P wave to 100 seconds after the P wave to ensure all converted phases to a depth of 100 km are included. The data were then rotated to radial and tangential components. A signal to noise ratio check was done by comparing the power in the seismic traces 20 sec before and after the predicted arrival time of the P wave. Only seismograms with signal to noise ratios higher than 2 for both the P and SV components, were used to ensure high quality data before deconvolution. Overall, the signal to noise ratio criteria removed 30-40% of all collected traces. For example, at station MA35, from 120 events collected, 44 were removed due to high noise levels. The deconvolution was done in the frequency domain by using the water-level stabilization method and a low-pass Gaussian filter to remove high frequency noise (Langston, 1977). The receiver function H(ω) is calculated from: © 2015 American Geophysical Union. All rights reserved. H (ω ) = R(ω )Z * (ω ) G (ω ) max Z (ω )Z * (ω ), c* max(Z ( w) Z * ( w)} { (Eq. 1)  −ω2   ……………………………………...(Eq. 2) G (ω ) = exp 2   4α  where, ω is angular frequency, Z(ω) is the fourier transform of the P-component waveform, R(ω) is the transform of the SV component, and Z*(ω) is the complex conjugate of Z(ω). G(ω) is a Gaussian filter that has zero phase distortion and a lack of side-lobes (Langston, 1979). It is added to eliminate high frequency noise in the RF. The frequency content of the Gaussian filter is controlled by the parameter α. Finally, the denominator of equation 1 is either Z(ω)Z*(ω) or a constant c times the maximum of Z(ω)Z*(ω) depending on which value is larger. This is necessary to account for holes in the spectra of Z(ω) and therefore division by a small number which is inherently unstable. The constant c is called the “water level” and serves to stabilize the deconvolution. The values of α and c were chosen by trial and error where we tried to make the receiver function as sharp as possible but also tried to minimize noise. All receiver functions were computed using a water-level parameter c of 0.001 and a Gaussian smoothing parameter α of 3.5 resulting in receiver functions with a dominant period near 3s. A final visual check was also performed. Good receiver functions are identified by having a sharp P wave signal with little energy arriving earlier. Low quality RF’s tend to have anomalously high amplitude signals at later times or very wide side lobes. We eliminated these data before further data processing. An example of receiver functions for station MA18, plotted as a function of epicentral distance, is shown in Figure 3. Most of the receiver functions have similar © 2015 American Geophysical Union. All rights reserved. signals with a peak at 0 seconds (the P wave) followed by negative side lobes. Almost all RFs from MA18 show a strong positive pulse at about 5 second. After 5 second, there is variability in the RFs but we can still observe some continuity of arrivals versus distance. Station MA18 is located far from the coast and is expected to have a relatively flat Moho. At this station, we collected 35 high quality RFs and the clear arrival at about 5 sec is interpreted as the P- to S-wave conversion from the Moho. We also performed a stack of all the MA18 receiver functions in the depth domain. Starting with a 1-D velocity model and horizontal slowness information, each RF can be interpolated into the depth domain to correct for move out and then stacked. On the right side of Figure 3, a stacked RF in the depth domain shows high amplitude at 0 km followed by a strong positive amplitude at 39 km depth interpreted as an arrival from the Moho. Strong signals are also seen at 150 km (16.3 s in the time domain) and 180 km (22 s in the time domain) which are likely PpPs and PsPs/PpSs multiples from the Moho but which could be due to conversions from deeper mantle discontinuities. Receiver Function Imaging Ideally receiver functions show pulses as a function of time that are due to converted P to S waves from interfaces in the subsurface beneath a seismic station. To convert these pulses to depths to interfaces beneath the station requires knowledge of the P and S velocities in the rocks above the interface. Zhu and Kanamori (2000) introduced a method (the H-κ method) that can minimize the ambiguity due to the trade-off between depth and velocity. The time difference between the converted wave and direct P wave depends on the angle of incidence of the incoming P wave, © 2015 American Geophysical Union. All rights reserved. the subsurface P and S velocities and the depth of the interface. The timing of multiple reflected waves (PpPs and PpSs+PsPs) from the conversion depth also depends on incidence angle but differ in their dependence on depth of converter and velocity. Assuming a flat layer over a half space and a uniform velocity model, the estimated arrival time for the Ps wave as well as the multiples can be calculated. The essence of the H-κ method is that for each station, all receiver functions are summed at the times corresponding to the Ps arrival time as well as the arrival times for the multiples for various choices of thickness H of the layer as well as Vp/Vs ratio (Zhu and Kanamori, 20000). For correct values of H and Vp/Vs ratio, the summation should be a maximum. The approach is essentially a grid search but assumes the structure beneath a station is a single layer (the crust) overlying a half space. The H-κ method stacks all receiver function amplitudes at one station including the Ps converted wave and its two multiples. The direct Ps converted wave has the largest amplitude and is the least sensitive to dips in the crust-mantle boundary thus the Ps converted wave is often given heavier weight in the stacking process. The stacking can be written as: S (H , κ ) = ∑ w1 r (t1 ) + w2 r (t 2 ) − w3 r (t 3 ) ………………..(Eq. 3) Where, H is crustal thickness, κ is Vp/Vs ratio, r i (t) is the ith receiver function at times t 1 , t 2 , and t 3 which are the predicted times for T Ps , T PpPs and T PsSs+PsPs . We used a weighting of 0.5, 0.3, and 0.2 for w 1 , w 2 , and w 3, respectively. Zhu and Kanamori (2000) used 0.7, 0.2 and 0.1 for the weighting in their study of Southern California. We tried several weighting schemes and found a slightly better result with the weights © 2015 American Geophysical Union. All rights reserved. given above. However, overall we didn’t see a dramatic change in results using different weighting schemes. Crustal thickness and Vp/Vs ratio of the crust beneath Southwestern Mexico were analyzed using the H-κ stacking method described above. The bounds for the grid search were set to be 20 to 55 km for Moho depth and 1.65 to 2.00 for Vp/Vs ratio. Figure 4, illustrates the results for station MA18 which had 35 high quality receiver functions. We contoured the value of the stacked receiver functions as a function of Moho depth and Vp/Vs ratio with the highest value being the likely true parameters. The contour plot shows a clear maximum with a realistic number for crustal thickness (39.2 km) and crustal Vp/Vs ratio (1.82). The predicted Moho Ps arrival times agree with the receiver function signals showing a strong positive converted wave at 5 s. Predicted times for the multiples are also plotted. The contour plots using the H-κ method show strong peaks and thus good resolution for most inland stations far from the coast but show poorer results for stations closer to the coast. This is not surprising as it is likely there are strongly dipping structures near the coast where the subducting plate is near the surface. For the coastal stations the peaks in the H-κ contour plots are not sharp and just taking the maximum value from the plots results in unrealistic values for Moho depth and Vp/Vs ratio. The H-κ method can be extended to dipping layers (ex. Rossi et al., 2006) but given the complicated nature of the crustal structure near the coast (see below) we did not attempt more sophisticated inversions. For the coastal stations where the H-κ method failed we used a velocity model with a Vp/Vs ratio of 1.78, the global crustal average according to Chevrot and van der Hilst (2000), to stack receiver functions in the depth domain. The peak in the stacks was then measured to give an estimate of crustal thickness. Near the coast the first peak in the RFs is usually a negative pulse indicating a drop in velocity © 2015 American Geophysical Union. All rights reserved. with depth. This has been observed in other subduction zones (Bostock, 2013) and has been interpreted as over pressured oceanic basalt crust. Thus, we picked these negative peaks as the Moho. The RF’s with initial negative pulses are discussed further below. The contour of the H-κ stack amplitude allows one to calculate uncertainties by measuring the flatness at the maximum point as discussed in Zhu and Kanamori (2000). For coastal stations the H-κ method did not give satisfactory results. Therefore, to estimate uncertainties for all data we employed a bootstrapping technique (Efron and Tibshirani, 1991) to estimate the 2σ uncertainty for each stacked receiver function. We sampled 100 random populations of receiver functions for each station, for this purpose, and then processed them as discussed in the text. For stations where the H-κ worked we processed the bootstrapped data with that technique. For stations where the H-κ method did not work we processed data using a Vp/Vs ratio of 1.75 and 1.85 to give uncertainty in Moho depth due to uncertainties in Vp/Vs ratio. We also performed a bootstrapping analysis with a constant Vp/Vs ratio and added the two uncertainties together. The results of the RF analysis are shown in Figures 5 and 6 and listed in Table 1 alo ng with the 2 σ errors. In Figure 5 point measurements of crustal thickness are interpolated and plotted in map view. Figure 6 shows the interpolated Vp/Vs ratios in the region where we were able to make reliable measurements. Crustal thickness varies from 20 km along the coast to a maximum of 42 km inland. There is a difference in crustal structure between the Jalisco block to the northwest and Michoacan to the southeast. Beneath the Jalisco block, which is underthrust by the Rivera plate, the crust is 25-30 km thick along the coast and rapidly reaches 40 km thick within 50 km of the coast. To the southeast of the Colima Rift, however, the © 2015 American Geophysical Union. All rights reserved. crustal thickness along the coast is 20-25 km thick and gradually thickens inland reaching 40 km more than 100 km from the coast. In Figure 7 we plot the crustal thickness estimates at each station on a map of topography. Note that there is not a good correlation between crustal thickness and topography. The thickest crust is roughly in the middle of our study region, parallel to the strike of the subduction zone, but has rather subdued topography. Further inland, the topography increases to near 2 km yet the crust actually thins a little to 37-39 km thick. This indicates the high topography is likely supported by buoyant mantle as has been observed in other backarc regions (Hyndman et al., 2005). This observation also supports a model of slab rollback and perhaps tearing as proposed by Ferrari (2004) and later by Yang et al. (2009). In this model hot asthenospheric mantle fills the space previously occupied by slab to relatively shallow depths providing a source of buoyancy in the mantle beneath the high topography. If this interpretation is correct then the high topography would be relatively recent i.e. following the time of the roll back estimated to be during the last 5 Myr or so (Ferrari, 2004). This model is also consistent with observations of current and recent tectonic activity in the northern Jalisco Block including high rock uplift rates determined from stream incision profiles (Castillo et al., 2014), seismic activity in the northern Jalisco Block (Pacheco et al., 1999; NunezCornu et al., 2002) and Pliocene active deformation seen from structural mapping (Ferrari and Rosas-Elguera, 2000). Another interesting aspect of the crustal thickness map is the lack of significant thinning of crust in the Colima Rift (Figure 5). We do not have the station density to look at the Rift crustal structure in detail but stations close to the Rift such as MA26 and MA24 do not show thinner crust than adjacent stations further from the Rift. We conclude that there has been little extension in the Rift. This agrees with the © 2015 American Geophysical Union. All rights reserved. conclusion of Rosas-Elguera et al. (1996) that the southern Colima Rift is actually a broad slowly extending zone that has been active only since the late Pliocene. We do see a trend to the north of thinning crust from 40-42 km thick to 36-37 km thick. However, given that our deployment did not reach the surface manifestation of the Tepic-Zoacala Rift, it is difficult to place constraints on the overall extension across that Rift. We found the average Vp/Vs ratio of the crust to be 1.81 (Figure 6) which is slightly higher than the global average of 1.78 (Christensen, 1996; Chevrot & Van der Hilst, 2000). There is a large range in values, however, from 1.72 to 1.87. Two regions, in particular, have abnormally high crustal Vp/Vs ratios. The first region is located in the central part of the Jalisco block, close to the volcanoes shown as red triangles in Figure 6. Four stations located at the southwestern edge of the region have Vp/Vs ratios of 1.85 or greater. It is interesting that to the north of these stations the Vp/Vs ratio is closer to normal. The second region, located to the northeast of Colima Volcano, shows a band of high Vp/Vs ratios ranging from 1.85 to 1.87. The eastern part of this band is within the Michoacan-Guanajuato Volcanic Field (MGVF), a wide region of volcanic activity that began about 2.8 Ma and continues today (GomezTuena et al., 2007). The band of high Vp/Vs crust extends further west than the MGVF, towards the northern Colima Rift. Bandy et al. (1995) located the RiveraCocos boundary at depth beneath this region based on the occurrence of thermal springs and crustal seismicity in the region although Yang et al. (2009) place the boundary further west. Average crustal Vp/Vs can be used to interpret the petrology and physical state of the crust. Christensen (1996) showed from laboratory experiments that Vp/Vs ratio does not vary much with changes in temperature and pressure for pressures © 2015 American Geophysical Union. All rights reserved. greater than 100-200 Mpa. The main factor that controls the Vp/Vs ratio in the crust is the presence of melting or fluids and the mineralogy. The relative abundance of quartz and plagioclase feldspar has a dominant effect on Vp/Vs (Christensen, 1996): for felsic quartz-rich rocks such as granite, Vp/Vs is 1.71; intermediate rocks have a Vp/Vs ratio near 1.78 and mafic plagioclase-rich rocks such as gabbro have a Vp/Vs ratio near 1.87. The average composition for continental crust is close to andesite or diorite (Anderson, 1989) and laboratory measurements by Carmichael (1982) confirmed that Vp/Vs for diorite at crustal pressures ranges from 1.75 to 1.79. We found the region just to the east of Colima volcano and the northern part of the Jalisco block have “normal” continental crust Vp/Vs ratios (Figure 6). The high Vp/Vs ratio regions could indicate a very mafic crust or that the crust has high pore pressure fluids or partial melt. The surface geology of the region is largely granitic or siliceous volcanic rock with no indication of a particularly mafic crustal composition (Valencia et al., 2013). Thus we conclude that the high Vp/Vs ratios we observe along two bands in our study area are due to partial melt or high fluid content within the crust. The high Vp/Vs ratios in the Jalisco block are just seaward of a trend of volcanoes (red triangles in Figure 6) known as the Central Jalisco Volcanic Lineament (CJVL). It has been known for some time that there has been an overall trenchward migration of magmatism in Jalisco for the past 10 Ma (Ferrari et al., 2001). The CJVL forms the front of this migrating magmatism. Ferrari (2004) and Yang et al. (2009) proposed that this migration was due to slab rollback of the Rivera plate allowing asthenosphere to warm the lithosphere progressively trenchward. Arc magmas and fluids could also be rising from the deepening slab. The highest Vp/Vs ratios we find in the Jalisco block are just trenchward of the CJVL volcanoes and are not collocated with any volcanoes. This indicates the crust is being broadly heated, with possible © 2015 American Geophysical Union. All rights reserved. partial melting, or has had the addition of fluids trenchward of the recent volcanism. Thermal modeling of Rivera plate subduction shows dehydration of the slab just trenchward of the volcanic front (Ferrari et al., 2012) supporting the latter model. We further note that magnetotelluric results show high crustal conductivity in this region (Corbo-Camargo et al., 2013) consistent with high crustal fluid content The second region with anomalously high crustal Vp/Vs is located in the north east of our study area partly within the Michoacan-Guanajuato Volcanic Field. This region has late Pliocene – Quaternary mafic and intermediate volcanism less than 3 Myr (Gomez-Tuena et al., 2007). It is a 40,000 km2 area with more than 1000 Quaternary eruptive centers and thus a crust with distributed partial melt is not surprising (Gomez-Tuena et al., 2007). We find high crustal Vp/Vs ratios to the west of the MGVF, as well. There is no volcanism in this region so high fluid content is likely the cause which is also consistent with the location of the beginning of slab dehydration in the thermal model derived by Ferrari et al. (2012). Cross sections In order to study the down-dip variation of the Rivera and Cocos subduction zone, we plot RF’s along two linear profiles relative to an estimate of depth to the subducting Rivera and Cocos plates. We divide the data into two groups (Figure 8) corresponding to the Jalisco Block with the Rivera plate subducting beneath it and stations to the east below which is the subducting Cocos plate. We make this division because of the difference in crustal structure discussed above as well as the difference in slab dip reported by Pardo and Suarez (1995) between the two regions. We use estimates of slab interface depth from Pardo and Suarez (1995). For each of our stations, we plotted the stacked receiver function as a function of depth to slab (Figure © 2015 American Geophysical Union. All rights reserved. 9). The stacking was done assuming the Vp/Vs ratio found in our H-κ analysis. For stations that had a null result from the H-κ analysis we assumed a Vp/Vs ratio of 1.78. An estimated depth from H-κ analysis is also placed on the figure. Overall, the stacked RF’s show good agreement with the H-κ analysis results although there are slight differences that may be due to a more complex crustal velocity structure than the simple homogenous layer assumed in H-κ analysis. Both panels in Figure 9 show a northerly dipping structure in the southwest of the array that is visible to near 60 km depth whereas in the northeast there is a relatively simple Moho near 40 km depth with no visible dipping structures. We plot the RFs to 80 km depth because at deeper depths crustal multiples are visible and make interpretation problematic. We interpret the dipping structures as the subducting Rivera and subducting Cocos plates (Figure 9a,b), respectively. Note that for the subducting Cocos plate the top of the descending plate agrees with the results of Pardo and Suarez (1995) but that we find the Rivera plate to be about 10 km systematically deeper than Pardo and Suarez (1995) claim. This discrepancy is not surprising in that seismicity is sparse in the Rivera subduction zone and instrumentation was lacking in the past. The techniques used above, H-κ analysis and stacking of receiver functions, work well for flat lying structures but are problematic for dipping structures. Our analyses to this point have assumed all the arrivals from a given boundary at a single station are produced at the same depth. To account for the different location of conversion points corresponding to different back azimuths to a station we have also processed our receiver functions using the Common Conversion Point (CCP) stacking method. Individual receiver functions were back projected along the raypath corresponding to the receiver function. The back projected amplitudes are stacked in © 2015 American Geophysical Union. All rights reserved. lateral and vertical bins resulting in a 3D image of convertors in the sub surface. This approach has been used in numerous studies (e.g. Dueker and Sheehan (1997), Schulte-Pelkum et al (2005),) but requires a high density of stations with numerous receiver functions such that each bin is sampled by numerous traces from different directions. The MARS array is relatively sparse so we collapsed our images onto 2-D lines so as to better compare the results with Figure 10. We divided the data into the same two groups as before (Figure 8) and back projected the receiver functions onto lines RR’ and CC’ respectively. This was done by projecting each convergence point to its correct location and then projecting that location onto the respective lines along a perpendicular to the lines. Figure 10 shows the results of our CCP stacking along the two lines, RR’ and CC’. We used a bin size radius of 12 km for the CC’ line and 14 km for the RR’ line. The results can be compared to Figure 9 although it should be noted that we use distance from the trench on the horizontal axis here but depth to slab was used in Figure 9 as the horizontal axis. The red color in Figure 10 corresponds to a positive pulse in the receiver functions and thus to a jump in shear velocity with increasing depth. The blue color corresponds to a decrease in shear velocity with depth. The CCP results are similar to the stacking results shown in Figure 9 in that a clear dipping structure is seen in the southwest and a sharp relatively flat continental Moho is seen in the northeast for both profiles. The dipping slab structure is seen to about 50 km depth and then disappears. Slab Structure The cross sections shown in Figures 9 and 10 show several interesting features. First, the dipping structure in the southwest is marked by a negative pulse © 2015 American Geophysical Union. All rights reserved. over a positive pulse (blue over red in Figure 10 –dashed lines in Figure 9) indicating that there is a shear velocity drop with depth associated with the subducting slabs underlain by a sharp jump in velocity. Figure 11 shows the raypaths and a hypothetical velocity model for this situation, as well as an example of a receiver function response. Here we assume the dipping crustal layer has a low velocity layer underlain by a jump in velocity to normal crustal velocity and then finally another jump in velocity representing the oceanic Moho. The Ps converted wave going from low velocity into high velocity has a negative polarity. The converted waves from the deeper layers go from higher velocity at depth to lower velocity above and thus have positive polarities. There are also multiple reflections and conversions within the layers. All these waves arrive close in time resulting in a two sided pulse in the receiver function similar to what we observe in the coastal RF data. One might interpret the slow shear velocity layer to be subducting sediment, however, the gravity study of Manea et al. (2003) suggests that the thickness of the sedimentary column over the Rivera plate probably does not exceed ~20 m, and that it gradually increases eastward along the trench. Moreover, the Deep Sea Drilling Project at site 487, located ~11 km offshore of Guerrero State (lat 15°51.210′ N and long 99°10.518′ W) , found a sedimentary column of ~100 m of Quaternary hemipelagic sediments which overlay ~70 m of late Miocene to Pliocene pelagic. A layer just a few 100 m thick would be too thin to resolve at the periods of the RFs. A slow velocity layer associated with subducting oceanic lithosphere has been observed in several other subduction zones (Bostock; 2013). Bostock et al. (2002) used a scattered wave inversion technique on seismic data recorded over the Juan de Fuca plate in central Oregon. Their images show a layer of slow velocity that is associated with the subducting plate to a depth near 45 km. Bostock et al. (2002) © 2015 American Geophysical Union. All rights reserved. interpret the slow anomaly as the basaltic oceanic crust subducting beneath continental forearc. In Figure 9b it can be seen that the RF peak corresponding to the top of the slow velocity occurs at the depth predicted for the subducting plate by Pardo and Suarez (1995) and is even deeper in Figure 9a. Thus, we interpret the drop in velocity found in the receiver functions to mark the top of the subducting oceanic crust in accord with Bostock et al. (2002). The estimate of thickness of the low velocity layer depends on the Vp/Vs ratio of P and S waves within the slow layer. Using normal values of Vp/Vs ratio the thickness of the layer in Bostock et al. (2002) corresponds to a layer roughly 8 km thick and thus they associated the layer with the entire oceanic basalt crust. However, Audet et al. (2009), using multiple reflections within the slow layer for data collected above Cascadia, found that the Vp/Vs ratio is abnormally high within the low velocity layer, 2.4 to 2.8, and thus the low velocity layer is thinner than previously thought and is on the order of 3 to 5 km thick. Recently, Kim et al. (2010) and Song et al. (2009) examined data from a dense linear profile across central Mexico. They also found a thin 2-4 km thick layer of anomalously slow S velocity with anomalously high Vp/Vs ratio dipping beneath the continent. Audet et al. (2009) as well as Kim et al. (2010) and Song et al. (2009) identify the slow layer as the top half of the subducting oceanic crust. To better constrain the structure producing the RFs, we calculated synthetic receiver functions for various models using a reflectivity code (Levin and Park, 1997). The three parameters we investigated are the drop in shear velocity at the basalt interface, the Vp/Vs ratio within the low velocity layer, and the thickness of the layer. The velocity model we used has 4 layers, a continental crust, a low velocity layer identified with the subducting oceanic crust, a lower oceanic crust, and the mantle. © 2015 American Geophysical Union. All rights reserved. We constrain the continental crust, lower oceanic crust and mantle to have a Vp/Vs ratio of 1.78 and vary the thickness, shear velocity, and Vp/Vs ratio of the low velocity layer to match our data. The velocity of the lower oceanic crust and mantle are taken as normal, 6.8 km/sec Vp and 3.8 km/sec Vs for deep ocean crust and 8.0 km/sec Vp and 4.5 km/sec Vs for the mantle. Our receiver functions along the coast vary considerably with different amounts of noise such that no single model will fit all the data. For the purposes of waveform modeling we stacked all coastal RFs with negative first motions near the Moho arrival time with back azimuths near 3000. The stacked RF is then an average of converted waves from the subducting Rivera and Cocos plates from one direction. Figure 12 shows the results of our simulations with varying thicknesses of the low velocity layer compared to the stacked coastal RF. We varied the thickness of the low velocity zone from one to six km thick. For each thickness, we adjusted the S velocity drop in the low velocity zone as well as the Vp/Vs ratio to best match the stacked data. The S velocity drop controls the amplitude of the first arriving negative pulse and the Vp/Vs ratio controls the timing of the following positive pulses. Our optimal model is a 3 km thick layer with a drop in shear velocity of 25% and a Vp/Vs ratio within the low velocity layer of 2.0-2.1, well above the range for normal rocks at this depth. This is similar to the results of Audet et al. (2009) although we clearly do not have tight constraints on the slow velocity layer’s properties due to variable receiver functions in our profile as well as the fact that multiple variables affect the results. Basalt is usually considered a high velocity component of crust, particularly within continental crust, therefore it is surprising that it is showing up as a slow velocity anomaly in our study as well as in the work of Bostock et al. (2002). Audet et © 2015 American Geophysical Union. All rights reserved. al. (2009) studied this issue and concluded that the basalt would have the observed low velocities if it had pervasive water present that was overpressured. They propose that the top of the oceanic crust formed of pillow basalts and sheeted dikes of gabbro, contains significant water as well as hydrated minerals. As the plate subducts some of the minerals dehydrate creating high pore fluid pressures that cause the dramatic drop in shear velocity observed. This would imply the boundary between the subducting oceanic crust and the overriding continental material is impermeable as well as the deeper crust. Kim et al (2010) and Song et. al. (2009) proposed a similar interpretation for their images in central Mexico although Kim et al. (2010) also suggest hydrous minerals such as talc may be an important contributor to the velocity drop. We observe the disappearance of the slow dipping converted phase at about 45 km (Figure 10). We interpret this to be the depth where basalt begins to transform to eclogite. When this occurs there is a relatively large volume change that may make the subducting crust and its boundary more permeable and thus allow it to lose water to the overlying mantle wedge. If the slow velocity layer is due to high pore pressure fluids then its disappearance may show the depth at which the water can escape into the mantle wedge. Eclogite also has high seismic velocity and thus its contrast with the surrounding mantle will be less than the contrast of the slow basalt with mantle causing the disappearance of the positive converted wave at the oceanic Moho. The slow velocity layer found by Bostock et al. (2002) also disappears at a similar depth. The similarity of subduction zone structure between Cascadia and southwestern Mexico is striking but both regions have young ocean (5Myr – 15Myr) subducting beneath continent. One difference between the regions is the thickness of sedimentary cover in the subduction zone. Off shore southwestern Mexico estimates © 2015 American Geophysical Union. All rights reserved. of the sedimentary thickness range from 20 to 200 m, as mentioned above. However, the subduction zone offshore Cascadia has a thick sedimentary prism and far more flux of sediment into the subduction zone (Rea and Ruff, 1996). The similarity of the two regions shows sediment flux is not a significant factor in the permeability of the oceanic crust interface nor the temperature structure of the slab as the eclogite phase transition seems to occur at similar depths as well. It also supports the model for the slow velocity layer being basalt crust and not a layer of subducting sediment. A second unusual feature of the CCP images of Southwestern Mexico (Figure 10) is the lack of any strong Moho signal in the middle of the images, both for the Rivera plate system and the Cocos plate system. A clear continental Moho is visible to the east and clear dipping convertors are seen to the west but in the middle there are no clear signals corresponding to the Moho. Above the Rivera slab, we measure a 20 km width of weak Moho signal and above the Cocos slab we measure about a 50 km width without a clear Moho. Again, this is similar to what is observed beneath Oregon. Bostock et al. (2002) explained this observation by postulating that the mantle wedge at this location is serpentinized due to the water release from the subducting slab. Serpentine has very low shear velocity and in fact has slower shear velocity than typical lower crustal mineral assemblages (Christensen (1996). Bostock et al. (2002) show that a peridotite with 50% serpentization will have shear velocity similar to lower crustal rocks. With no shear velocity contrast between crust and mantle no converted P to S wave will be created at the Moho. The implication is that this mantle does not participate in the mantle flow associated with the subducting slab and is hydrated by released water from the subducting crust. Our data support this interpretation. It is interesting that the width of the zone with weak Moho is quite different between the Rivera and Cocos subduction zones. The Rivera plate has a © 2015 American Geophysical Union. All rights reserved. steeper subduction angle than the Cocos plate thus creating a thicker mantle wedge as a function of distance inland. A thicker mantle wedge can be entrained in the downward flow easier than a thin wedge. These observations may provide constraints on viscosity of the mantle wedge. Conclusions We have analyzed the lithospheric structure beneath Southwestern Mexico using the receiver function technique. Along the coast, receiver functions show the dipping Rivera and Cocos plates subducting beneath Mexico. We find the Rivera plate along the coast is about 10 km deeper than previously estimated by Pardo and Suarez (1995). The receiver functions also show that the subducting slabs have a 3-4 km thick extremely slow shear velocity layer that is similar to what has been found in several other subduction zones. The thickness of the layer leads us to conclude this is part of the subducting oceanic crust that is under high pore fluid pressure and not a layer of sediment. Further detailed modeling of this feature would require deployment of a dense profile of seismic stations. The contrast of a better imaged profile in Mexico with Cascadia should shed light on the role of thick versus thin accretionary prisms in the subduction process. Inland receiver functions show a clear Moho across much of our study area as well as variations in crustal Vp/Vs ratio. There is no correlation of topography with crustal thickness implying mantle buoyancy causes much of the high relief in Jalisco and Michoacan. High crustal Vp/Vs ratios are seen just trenchward of the arc front or in regions of recent magmatism. This correlation suggests that broad swaths of crust are heated or infiltrated by fluids caused by migrating arc magmatism towards the trench as the Rivera and perhaps Cocos slabs roll back. © 2015 American Geophysical Union. All rights reserved. Acknowledgements We thank all participants in the MARS experiment, especially Alejandro Martinez for leading the field work including siting and maintaining the instruments. Discussions with Luca Ferrari were greatly appreciated and beneficial. Two anonymous reviewers made suggestions that greatly improved the manuscript. Financial support for this study was provided by the National Science Foundation through grant EAR-0335782, the Geology Foundation of the Jackson School of Geosciences at the University of Texas, and grant IN 117205 from Universidad Nacional Autonoma de Mexico. Instrumentation and field support was provided by the IRIS-PASSCAL Instrumentation Center. The Observatorio Vulcanologico of the Universidada de Colima kindly provided space throughout the course of the field work and helped with customs. Finally, we thank Mike West for sharing data from the CODEX experiment. The MARS seismic data are available for download at the IRIS Data Management Center. References Ammon, C. J., The isolation of receiver effects from teleseismic P waveforms, (1991), Bull. Seism. Soc. Am., 81, 2504-2510. Anderson, D. L., (1989), Theory of the Earth.Blackwell Scientific Publications, Boston, p. 366 Audet, P., M. G. Bostock, N. I. Christensen, and S. M. Peacock, (2009), Seismic evidence for overpressured subducted oceanic crust and megathrust fault sealing, Nature, 457, 76–78. Bandy, W., C. Mortera-Gutierrez, J. Urrutia-Fucugauchi, and C. W. T. Hilde. (1995), The subducted Rivera-Cocos plate boundary: Where is it, what is its © 2015 American Geophysical Union. All rights reserved. relationship to the Colima rift?, Geophysical Research Letter, 22, 3075-3078. Bandy, W. L., J. Urrutia-Fucugauchi, F. W. McDowell, and O. Morton-Bermea, (2001), K-Ar ages of four mafic lavas from the Central Jalisco Volcanic Lineament: Supporting evidence for a NW migration of volcanism within the Jalisco block, western Mexico, Geofisica Internacional, v. 40, 259-269. Bostock, M. G., (2013), The Moho in subduction zones, Tectonophys., 609, pp. 547557. Bostock, M. G., Hyndman, D. R., Rondenay, S and Peacock M. S., (2002), An inverted continental Moho and serpentinization of the forearcmantle, Nature, 417, 536–539 Carmichael, R.S., (1982), Handbook of Physical Properties of Rocks. CRC Press, Boca Raton, Fla. Castillo, M., E. Munoz-Salinas, and L. Ferrari, (2014), Response of a landscape to tectonics using channel steepness indices (k sn ) and OSL: A case of study from the Jalisco Block, Western Mexico, Geomorphology, 221, 204-214. Chevrot, S., and R. D. Van der Hilst, (2000), The Poisson's ratio of the Australian crust: geological and geophysical implications, Earth and Planetary Science Letters, v. 183/1-2, p. 121-132. Christensen, N. I., (1996), Poisson's ratio and crustal seismology, Journal of Geophysical Research, 101(B2), 3139–3156. Corbo-Camargo, F., J. A. Arzate-Flores, R. Alvarez-Bejar, J. J. Aranda-Gomez, and V. Yutsis, (2013), Subduction of the Rivera plate beneath the Jalisco block as imaged by magnetotelluric data, Revista Mexicana de Ciencias Geologicas, 30., 268-281. © 2015 American Geophysical Union. All rights reserved. DeMets, C. and S. Traylen, (2000), Motion of the Rivera plate since 10Ma relatives to the Pacific and North American plates and the mantle, Tectophysics, 318,119159. Dueker, K. G. and A. F. Sheehan, (1997), Mantle discontinuity structure from midpoint stacks of converted P to S waves across the Yellowstone hotspot track, Journal of Geophysical Research, 102, 8313-8327. Efron, B. and R. Tibshirani, (1991), Statistical data analysis in the computer age, Science, 253, 390-395. Ferrari, L. and J. Rosas-Elguera, (2000), Late Miocene to Quaternary extension at the northern boundary of the Jalisco block, western Mexico: The Tepic-Zacoalco rift revisited, Geol. Soc. Am. Special Paper, 334. Ferrari, L., C. Petrone, and L. Francalanci (2001), Generation of oceanic-island basalt-type volcanism in the western Trans-Mexican volcanic belt by slab rollback, asthenosphere infiltration, and variable flux melting, Geology, 29, 507-510. Ferrari, L., (2004). Slab detachment control on mafic volcanic pulse and mantle heterogeneity in central Mexico, Geology, 32(1): 77 – 80. Ferrari, L., M. T. Orozco-Esquivel, V. Manea, M. Manea, (2012), The dynamic history of the Trans-Mexican Volcanic Belt and the Mexico subduction zone. Tectonophysics, Invited review paper, doi:10.1016/j.tecto.2011.09.018. Gaite, B., A. Iglesias, A. Villaseñor, M. Herraiz, and J. F. Pacheco., (2012), Crustal structure of Mexico and surrounding regions from seismic ambient noise tomography, Geophysical Journal International, 188(3), 1413-1424, doi: 10.1111/j.1365-246X.2011.05339.x. Gardine, M., T. Dominguez, M. West, S. Grand, S. Suhardja. (2007), The Deep © 2015 American Geophysical Union. All rights reserved. Seismic Structure of Volcan de Colima, Mexico Eos Trans. AGU, 88(23), Fall Meeting Supplementary, Abstract T51A-02. Gómez-Tuena, A., Ma. T. Orozco-Esquivel, and L. Ferrari, ( 2007), Igneous petrogenesis of the Trans-Mexican Volcanic Belt, Geological Society of America Special Paper, 422. Hyndman, R. D., C. A. Currie, and S. P. Mazzotti, (2005), Subduction zone backarcs, mobile belts, and orogenic heat, GSA Today, 15 no. 2, 4-10. Kim, Y., R. W. Clayton, and J. M. Jackson, (2010), Geometry and seismicproperties of the subducting Cocos Plate in central Mexico, Journal of Geophysical Research., 115, B06310. Lawver, L. A., I. W. D. Dalziel, I. O. Norton, L. M. Gahagan, and J. Davis, (2013), The Plates 2013 Atlas of Plate Reconstructions (500 Ma to Present ay), Plates Progress Report No. 359-0413, University of Texas Technical Report No. 199. Langston, C. A., (1977), Corvallis, Oregon, Crustal and upper mantle structure from teleseismic P and S waves, Bull. Seismol. Soc. Am., 67, 713–724. Levin, V. and J. Park, (1997), P-SH conversions in a flat-layered medium with anisotropy of arbitrary orientation. Geophysical Journal International, 131, pp 253-266. Luhr, J. F., F. J. Allan, E. S. I. Carmichael, A. S. Nelson, and T. Hasenaka, (1985), Active rifting in southwestern Mexico: Manifestation of an incipient eastward spreading ridge jump, Geology, 13, 54-57. Manea, M., V. Manea, and V. Kostoglodov,, (2003), Sediment fill in the Middle America Trench inferred from gravity anomalies, Geofísica Internacional, v. 42, no. 4, 603–612. © 2015 American Geophysical Union. All rights reserved. Nunez-Cornu, F. J., R. L. Marta, F. A. Nava, G. Reyes-Davila, and C. SuarezPlacencia, (2002), Charateristics of seismicity in the coast and north of Jalisco Block, Mexico, Phys. Earth Planet. Int., 132, 141-155. Pacheco, J. F., C. A. Morter-Gutierrez, H. Delgado, S. K. Singh, R. W. Valenzuela, N. M. Shapiro, M. A. Santoyo, A. Hurtado, R. Barron, and E. Guttierrez-Moguel, (1999), Tectonic significance of an earthquake sequence in the Zacoalco halfgraben, Jalisco, Mexico , J. South Am. Earth Sci., 12, 557-565. Pardo, M. and G. Suarez, (1995), Shape of the subducted Rivera and Cocos plate in the southern Mexico: Seismic and tectonic implications, Journal of Geophysical Research, 100, 12, 357-12,373. Rea, D.K., and L.J. Ruff, (1996), Composition and mass flux of sediment entering the world's subduction zones: Implications for global sediment budgets, great earthquakes, and volcanism, Earth and Planetary Science Letters, 140, 1-12. Rosas-Elguera, J., L. Ferrari, V. H. Garduno-Monroy and J. Urrutia-Fucugauchi, (1996), Continental boundaries of the Jalisco block and their influence in the Pliocene Quaternary kinematics of western Mexico, Geology, 24, 921-924. Rossi, G., G. A. Abers, S. Rondenay, and D. H. Christensen, (2006), Unusual mantle Poisson’s rato, subduction, and crustal structure in central Alaska, J. Geophys. Res., 111, doi:10.1029/2005JB003956. Stock J. M. and J. Lee, (1994), Do microplates in subduction zones leave a geological record? Tectonics, 13, 6, 1472-1487. Schulte-Pelkum, V., G. Monsalve, A. Sheehan, M. R. Pandey, S. Sapkota, R. Bilham, F. Wu. (2005), Imaging the Indian subcontinent beneath the Himalayan. Nature 435, 1222-1225. Song, T. H., D. V. Helmberger, , M. R. Brudzinski, R. W. Clayton, P. Davis, X. © 2015 American Geophysical Union. All rights reserved. Perez-Campos, and K. S. Singh, (2009), Subducting slab ultra-slow velocity layer coincident with silent earthquakes in southern Mexico, Science, 324, 502-506. Urrutia Fucugauchi, J., J. H. Flores Ruiz, (1996), Bouguer gravity anomalies and regional crustal structure in central Mexico, International Geology Review, 38, 176-194. Valencia, V. A., K. Righter, J. Rosas-Elguera, M. Lopez-Martinez, and M. Grove, (2013), The age and composition of the pre-Cenozoic basement of the Jalisco Block: implications for and relation to the Guerrero composite terrane, Contrib Mineral Petrol, 166, 801-824. Yang, T., S. P. Grand, D. Wilson, M. Guzman-Speziale, M. J. Gomez-Gonzalez, T. Dominguez-Reyes, and J. Ni., (2009), Seismic structure beneath the Rivera subduction zone from finite-frequency seismic tomography, Journal of Geophysical Research, 114, B01302. Zhu, L. and H. Kanamori, (2000), Moho depth variation in southern California from teleseismic receiver functions, Journal of Geophysical Research, 105(B2), 2969–2980. © 2015 American Geophysical Union. All rights reserved. Figure 1. Map of Jalisco Block (JB) and adjacent regions, showing plate boundaries, volcanoes (open circles with number : 1, Volcan Colima; 2, Navado De Colima; 3, Volcan Cantaro; 4, Sierra La Primavera; 5, Volcan Tequila; and 6, Volcan Ceboruco); the Central Jalisco Volcanic Lineament (CJVL)[Bandy et al, 2001](red triangles: A: Ayutla, LV: Los Volcanoes, TA: Talpa de Allende, M: Mascota, and SS: San Sebastian, from southeast to northwest). Black circles are MARS stations. Red circles are CODEX station. Brown circle are Mexican stations. SCR, Southern Colima Rift; NCG Northern Colima Graben; TMVB, Trans-Mexican Volcanic Belt; MGVF, Michoacán-Guanajuato Volcanic Field; MAT, Middle American Trench; RI, Rivera Plate; CO, Cocos Plate; and PA, Pacific Plate. © 2015 American Geophysical Union. All rights reserved. Figure 2. Distribution of earthquakes used for the teleseismic Receiver Function study. Red triangles show station locations. Black circles show earthquake epicenters. © 2015 American Geophysical Union. All rights reserved. Figure 3. Receiver functions for station MA18 with time on the vertical axis and epicentral distance on the horizontal axis. On the right is a linear stack of the receiver functions in the depth domain using a 1-D velocity model to adjust for timing variations due to different incidence angles. © 2015 American Geophysical Union. All rights reserved. Figure 4. (a) is contour of stack amplitudes for station MA18 as a function of crustal thickness H and Vp/Vs ratio. The grid search calculates a stack amplitude of receiver functions for all ranges of crustal thickness (20-45 km) and Vp/Vs ratio (1.65-2.). The final result is taken by choosing the highest amplitude from the contour and uncertainty is calculated by measuring the flatness of the contour peak. (b) is a plot of all receiver functions from all events for this station The green line is the predicted Ps times assuming a crustal thickness of 39.2 km and a Vp/Vs ratio of 1.82. Similarly, the red line is predicted time for PpPs and the blue line is for PpSs. © 2015 American Geophysical Union. All rights reserved. Figure 5. An interpolation of crustal thickness measurements using the H-K method and picked Ps arrivals. The measured crustal thickness beneath individual stations is also given with units of kms. Crustal thickness varies between ~18 km on the coast to 43 km in land. CO is Colima Volcano and MGVF is an active magmatic region known as the Michoacán-Guanajuato Volcanic Field. The red triangles show magmatic centers along the Central Jalisco Volcanic Lineament. © 2015 American Geophysical Union. All rights reserved. Figure 6. An interpolated map of Vp/Vs ratio of the crust as well as the individual measurements. The Colima and Tepic-Zacoalco Rifts are shown as dashed lines. Young volcanoes along the Jalisco Volcanic Lineament are shown as red triangles. CO is Colima Volcano and MGVF is the Michoacán-Guanajuato Volcanic Field. © 2015 American Geophysical Union. All rights reserved. Figure 7. Map of the topography of Southwestern Mexico with crustal thickness measurements at individual stations in km. © 2015 American Geophysical Union. All rights reserved. Figure 8. The groups of stations used for the Rivera and Cocos Line RF plots. R-R’ is the line used for CCP imaging of the Rivera plate and C-C’ is used for CCP imaging of the Cocos plate. © 2015 American Geophysical Union. All rights reserved. Figure 9. (a) is plot of stacked RF’s with respect to depth to the slab for the Rivera Line. (b) is plot of stacked RF’s with respect to depth of the slab for the Cocos Line. The vertical axis is depth and horizontal axis is the depth to the slab taken from Pardo and Suarez (1995). The dashed red line is our interpretation of the dipping oceanic crust. The red stars are the estimates of continental crust depth using the H-K method. © 2015 American Geophysical Union. All rights reserved. Figure 10. (a) is the final CCP image result above the Rivera plate line (line R-R’ on fig 8). The vertical axis is depth in km and the horizontal axis is the distance calculated along the line. Note the starting point is close to the coast. Blue colors show negative receiver function pulses and red positive. (b) is the final CCP image for the Cocos plate (line C-C’). © 2015 American Geophysical Union. All rights reserved. Figure 11. (a) Synthetic seismogram response for the 4 layer velocity model shown in (b). The red line in (a) shows impulsive spikes corresponding to different converted arrivals while the blue line shows the corresponding receiver function at the longer periods used in this study.(c)-(h) show the raypaths for various converted and reflected waves with labels corresponding to the arrivals in (a). © 2015 American Geophysical Union. All rights reserved. Figure 12. The left figures show synthetic receiver function match real stack RF data at different thickness model and the right figures show the optimized velocity model to match the data. © 2015 American Geophysical Union. All rights reserved. Table 1. Crustal thickness and Vp/Vs ratio with 2σ uncertainties. For stations where the H-κ method worked we give both crustal thickness and Vp/Vs ratio (A). For the other stations we have no independent control on Vp/Vs ratio and just give the crustal thickness with associated uncertainties (B). A B station name longitude latitude thickness (km) VpVs ratio station name longitude latitude thickness (km) MA16 -103.254 19.9778 36.6 ± 1.7 1.84 ± 0.04 MA01 -103.9112 18.9014 31 ± 1.9 MA17 -102.033 19.4898 39.1 ± 1.4 1.86 ± 0.04 MA02 -103.6734 18.6273 18 ± 1.4 MA18 -102.3997 19.5489 39.3 ± 0.6 1.82 ± 0.02 MA04 -103.299 18.6885 31 ± 2.2 MA20 -102.616 19.7235 37.8 ± 6.9 1.82 ± 0.07 MA05 -103.125 18.9236 33 ± 3.5 MA21 -102.938 19.854 37.3 ± 0.5 1.87 ± 0.01 MA06 -102.8798 18.1461 22 ± 3.7 MA22 -102.909 19.59 40.3 ± 1.8 1.77 ± 0.04 MA07 -102.7943 18.8933 34 ± 3.0 MA23 -103.1034 19.6671 39.0 ± 0.5 1.77 ± 0.02 MA08 -103.0045 18.5318 29 ± 2.4 MA24 -103.5858 19.8719 37.0 ± 2.2 1.77 ± 0.04 MA09 -102.6567 18.0625 19 ± 0.8 MA27 -103.1435 20.1206 36.9 ± 5.7 1.81 ± 0.09 MA11 -102.344 18.4271 32 ± 3.6 MA28 -102.597 20.1786 40.6 ± 4.0 1.82 ± 0.07 MA12 -102.1906 18.7812 34 ± 5.0 MA29 -102.39 19.9038 41.2 ± 0.5 1.76 ± 0.02 MA14 -103.419 19.2385 36 ± 1.1 MA40 -103.942 20.2198 36.1 ± 6.3 1.73 ± 0.16 MA19 -103.5555 18.9097 32 ± 1.1 MA41 -104.2539 20.1853 36.8 ± 1.1 1.81 ± 0.03 MA25 -104.0746 19.6592 41 ± 2.9 MA42 -104.5175 20.777 38.2 ± 5.9 1.81 ± 0.13 MA26 -103.9354 19.3136 40 ± 3.0 MA43 -104.7815 20.2457 39.6 ± 5.0 1.83 ± 0.10 MA30 -104.2688 19.1269 33 ± 2.1 MA44 -104.2644 19.9627 41.5 ± 3.5 1.85 ± 0.05 MA31 -104.1818 19.467 38 ± 1.7 MA45 -104.2197 19.8121 42.0 ± 6.2 1.87 ± 0.14 MA32 -104.534 19.6162 37 ± 2.0 MA46 -103.968 19.7891 38.6 ± 2.6 1.81 ± 0.06 MA33 -104.5723 19.2432 30 ± 1.2 MA47 -103.8553 19.7118 38.3 ± 3.0 1.86 ± 0.15 MA34 -104.7898 19.6616 40 ± 2.5 MA48 -103.4362 19.5324 41.8 ± 4.9 1.78 ± 0.08 MA35 -104.6349 19.8796 40 ± 4.0 MA49 -103.3109 19.4638 42.0 ± 3.4 1.78 ± 0.03 MA36 -104.8932 19.3594 23 ± 1.1 MA50 -103.1648 19.4062 42.5 ± 0.6 1.77 ± 0.02 MA37 -105.3185 19.899 24 ± 1.5 MA51 -103.018 19.3679 42.5 ± 4.1 1.77 ± 0.05 MA38 -104.9836 19.9624 36 ± 2.1 MA53 -102.7644 19.2382 42.0 ± 6.2 1.67 ± 0.19 CJIG -105.043 19.499 23 ± 0.9 MA54 -103.8605 20.0884 36.1 ± 3.5 1.80 ± 0.08 MMIG -103.3456 18.2885 17 ± 1.0 © 2015 American Geophysical Union. All rights reserved.