Physical Acoustics Lab

In this post we have shown that the adaptive stacking tcas can reflect pickings done manually. Here we would like to present how tcas may help us in our seismic tomograpy study.

We use a catalog of three events from similar hypocentres travelling to an array of receivers.

The tcas result shows the following improvement in trace alignment. Each row shows the progression from raw data to initial alignment to the final alignment for one event.

We then plot the residual patterns at every stations that came out from tcas calculation. We superimpose pattern from all events. Red color shows positive residual which means larger than average value thus a slow field and vice versa for blue. The preliminary indication is that we have a slow region to the east of the AVF and fast to the westward.

Inversion uses a model with velocity grid cells spread roughly every 0.1 degrees (both NS and EW) and 10 km in depth. The picking error is estimated to be 0.1 s. Regularization is used with damping and smoothing terms of 1 and 2. And we have epsilon equals 1 and eta equals 10. The result shows misfit and chi^2 that are decreasing with iterations. With chi^s converge to value around 1 after two iterations.

The resulting model presents a perturbation of velocity field around the average with similar pattern to the presumption.

more events

We move forward this time by adding more earthquakes as our sources. As a qualification, we consider teleseismic events further than 27 degrees epicentrally from Auckland. Also we only consider events with good quality data by plotting and judging the trace output. We end up with 13 events.

There are reversal in polarity at station HIZ for events prior to 2014, and at OUZ before 2004. This has been addressed prior to applying tcas.

Final tcas alignment for all events shows good results. The output of its residuals/final time-shifts too are well contained between +/- 2 s. So we think tcas has performed as we expect.

The map of the residuals across the array shows more consistency. We again see slow velocity to the east of AVF.

Using the same model, inversion is done with standard deviation set to 0.1 s. If we decide to make all regularization terms equal to 0, chi^2 seems to converge around 7 after 4 iterations. This is equivalent to a variance reduction of 65%.

With regularization parameters for damping and smoothing set to 1 and 1 with epsilon of 4 and eta of 8, velocity anomaly start to show a pattern resembling the dt residual map. Chi^2 converge to 8 after 2 iterations.

station TLZ and TOZ

It seems like the conflicting information at station TLZ and TOZ prevent the fit to go below chi^2 equals to 7. Once we removed these two stations from the inversion, chi^2 converge around 5 after 4 iterations with regularization terms set to zero. The variance is reduced from 0.20 to 0.05 s^2 or equal to 75% reduction.

Similarly with the same regularization applied, chi^2 converge to 7 after 2 iterations but with more coherent resulting model.

inversion with catalog of 100+ events

We keep adding events for consideration in the inversion to get a better image of the subsurface presented by the resulting model. We would like to limit our teleseismic event to be not closer than 27 degrees epicentrally from our most Northwesterly station to triplication of arrival. The resulting catalog consist of 166 events:

We would also like for the catalog to contain only data with reasonable signal quality after lowpass filtering of 2.0 Hz has been done. On top of data availability, this cut the usable source data further but still we have 131 earthquakes which we have better confidence to work on.

After performing tcas to these events, we plot the residual dt map using the output from tcas. The map shows that the trend of slow to the east of AVF is retained. Residuals at station TLZ and TOZ are more consistent too. There are some stations near the centre where the prevalent residuals can be negative and positive. But these are fairly close to the orientation of the presumed split from fast to slow, so we think this is an understandable outcome.

The output of the inversion with zero regularization shows the following.

While by using the same regularization parameters as the previously,

We would like to compare the result we obtained in this post to the same set of events but by using the adaptive stacking TCAS approach.

We plot the raw data, initial alignment by removing the predicted arrival travel-time, and the final alignment after the cross-correlation. We allow the final time-shift to be within +/- 2 seconds.

The residuals contained in the output are consistent considering the same event path. Furthermore they are in relative agreement with the residuals obtained by manual picking with the mean removed from the observed minus predicted values.

Firstly the result shows that tcas is robust as it gives a consistent output for equal path. Secondly it shows that it reflects our manual work by showing similar perturbation pattern in the stations. Therefore there is the feasibility to use tcas to automate pickings using more sources.

events

We have three events in the Celebes Sea south of Philippine that have virtually the same epicenter. So presumably these would yield equal or at least similar paths to our stations.

In details, events locations are:

The ray-paths are projected to be indeed of uniform trajectory. Where we have basically a single path to each station from all three events.

arrivals and residuals

We would like to compare the residuals out of these three events i.e. the differential time of predicted ak135 model and observed arrival.

We will pick the arrival and find the residual by computing the time between that picked time (blue line) and predicted time (red line). We would expect them to be fairly close to one another for each respective station.

The results are combined in a spreadsheet:

Of course since the depth of the event is not exactly the same, this yields to a difference of two seconds for the amount of time it took to reach the stations. We think this is a reasonable difference in consideration.

The observed arrivals are all late arrivals in respect to the predicted arrivals. Value wise, the differences in residuals between events are fairly consistent.

inversion and misfit

We proceed with the inversion using this data. Regularization is turned on. It seems that the most satisfying result is obtained by making the damping in the regularization to equal zero. While we do so, smoothing need only to have a factor of 1 to have this progression in misfit.

The model correspond to this inversion steps is the following.

The negative anomaly is completely expected considering all of our arrival are late with respect to the 1D-Earth model.

We perform the inversion using teleseismic events only for comparison. The stations and picks are identical to the data for teleseismic earthquake in this attempt.

stations

Eighteen stations are used, divided into 11 in central Auckland, and 7 surrounding the former array. Low resolution channels are removed.

events

There are 121 events considered, though 13 are discarded due to bad data. We also have issue with the time window while we were picking.

inversion

The inversion shows negative velocity anomaly in most area in the Northland. The negative velocity in the North is a feature of our joint inversion.

There isn’t any change in travel-time residual distribution between initial and final model. The distribution is centred more in the positive direction which explains the mainly negative velocity anomaly.

Relative travel-time tomographic inversion is performed to the Auckland Volcanic Field (AVF), extending to include the northern parts of the North Island New Zealand. First P-wave arrivals are picked manually from earthquake traces that can be ranged in terms of quality. The travel-time residual between the picked arrival and the background ak135 1D-Earth model is the basis of this tomographic inversion.

stations

Eighteen stations are used in the inversion. The stations can be differentiated into the 11 that are centred to the city of Auckland, and the remaining 7 that surround the former array. For all of the stations, we are not using the low resolution channels.

events

The events that we use here can be separated into local and teleseismic events. The local events make up the sources that are originated from the South of the stations array. Meanwhile the teleseismic events, mostly from the pacifics, cover the illumination from the North.

There are 116 events sourced locally that we looked at. Most of these events have good trace data, but there are 4 events in which picking is not applicable due to the noise. From the remaining 112, some stations don’t show earthquake occurrence in the traces. This happened to the stations on the outside periphery, and the most times with station OUZ. We think that this may be related to the location of the stations. Station OUZ is the furthest station away from the sources, and thus the time window in which the traces are shown could be out of range here. All data that exhibit these situations are not considered for the inversion. Other traces that are regarded too noisy are also not included.

local events

Similarly with the teleseismic events, there are 121 events taken into account, 13 of which are discarded due to bad data. The same situation with the local data is also observed here. Station HIZ and OUZ at the two ends of the array are the most often disregarded from the inversion because of this. Again on top of all of this, if a trace data is too noisy it would also be discarded.

teleseismic events

iterations

To assess the inversion performance, we plot the evolution of standard deviation (in seconds) and chi^2 throughout the iterative steps. We see that after 3 iterations, we won’t get any more significant improvement to the model.The graphs show that the overall improvement from the initial model is less than 10%. Relaxing the damping parameter does reduce the std deviation and chi^2 but only slightly. However doing so might introduce redundant features.

standard deviation and chi^2

inversion result

The following result is obtained after performing 3 iteration steps. We see that the Auckland Volcanic Field is at the edge of a predominantly negative velocity anomaly (in km/s). We are still unsure on what might have caused this observation. But the negative anomaly in Auckland is a characteristic for areas of active volcanic activity. For example in the southern extent of the model, the negative anomalies reflect the locations of Mt. Taranaki, and the Taupo Volcanic Zone.

The red boxes signify the extent of our model that we’re most confident on. However it doesn’t mean that features outside of this box are all artefacts, they are simply less accurate due to the configuration of our sources and stations.

We add here the travel-time residual distribution for the initial and final model after 3 iterations. There is little change as explained by the standard-deviation analysis. The histograms are centred near zero, although we do have travel-times that are very late or very early.

In this travel-time tomographic inversion we plan to use more events to get a significant result in terms of AVF subsurface velocity field.

stations

Eighteen stations are used in the inversion. These are divided into 11 stations located around the central area, and 7 outer surrounding stations. We removed the low resolution channels associated with the station.

events and picks

We have 116 earthquakes around the region that we plan on using.

Almost all events have good trace data, which means that we have a good quality source input that we can work with. Our next step is to determine the first arrival pick with impeccable accuracy. The automatic trigger does a good job with picking the arrivals to a couple of seconds or so, but to consistently get an accuracy to the second or even less, picking by hand is still the go to method.

While we are doing the picking, some number of traces don’t show the event. We presume that this is mostly due to the location of the related station. The relative distances can make the time window in which the traces are shown becoming too short. For example, the northernmost station OUZ is the one that mostly show this condition. For the following results, the corresponding picks are excluded from the inversion, but we can revisit these picks and include them later on.

Even though the number of traces that don’t show an earthquake is small, it is still small things that we can improve on.

inversion

The split in the positive and negative anomaly (in km/s) which we also saw in the inversion with less events can be seen more distinctly this time around. Positive anomaly is prominent to the North and similarly, negative in the South. Another highlight is the negative anomaly right below Auckland (36.6-37.2 S, 174.6-175.2 E) at 50+ km.

We can even see a dip in the anomaly although we cannot be sure how much influence does the smearing effect has here. It is consistent with the primary paths direction from the South-East. However unlike any smearing that we have seen before in previous results, that are usually constricted to the outer and deep part of the model, here it is rather central. For example, the deep lineaments in the NS slice are likely to be smearing, and are in agreement with the three deepest events. Furthermore, our previous analysis from the checker-board test shows that at least for the first 70 km or so, the smearing is not too dominant. Therefore, we may be seeing a jump or a continuation of structure to depth commencing somewhere in the south Auckland.

Here are the plots showing how the root-mean-square travel-time residuals, standard deviation, and chi^2 change through each iteration. The first iteration is when the misfit drops drastically, afterward the plots quickly becoming stable. These plots justify the number of iterations we did as after three or four iterations we already have a convergence.

We have here the histogram of the travel-time residuals for initial and final model. There is no change between the two, but that is consistent with the standard-deviation plot. Most residuals are within a couple of seconds, but we do have some that can be more than 20 s late.

The movie below is a 3D visualization of estimates of the p-wave speed under New Zealand, from supplementary material of this paper. The three orthogonal planes cross approximately 50 km under Auckland.

These graphics are made with the python package mayavi.

Below are short descriptions of research projects available for postgraduate students in our group at the University of Auckland. If you are interested in any of these, please contact us for technical questions. For questions around enrollment at the University of Auckland, please go here. Scholarships may be available for the best applications through the University of Auckland, or the Dodd Walls Centre.

the Auckland volcanic field

We have a Marsden Fund Grant to study the Auckland Volcanic Field, with funding for a PhD and an MSc student in seismology. There are opportunities for field work on land and at sea.

Background:

Auckland, New Zealand’s largest city, sits on an active, monogenetic, and intraplate volcanic field. An international team of geoscientists from New Zealand, the USA, and China have joined forces to investigate why this volcanic field exists. An on- and offshore temporary seismic network will provide the density and the aperture of sensors to create the first high-resolution three-dimensional image of the Auckland Volcanic Field (AVF) from the surface to deep in the Earth’s mantle, beyond the depths of the magma that feeds the AVF. The resulting image will then constrain geodynamic modelling of mantle flow to explain what controls the AVF, and with it, this class of intra-plate volcanic fields.

ResonanT Ultrasound Spectroscopy

Resonant ultrasonic spectroscopy fills an important gap between our ultrasonic and seismic research. Together with the PORO lab, we have rock samples to complement laser ultrasound and a host of other petrophysical data with new RUS results.

Laser ultrasonic rock physics under high pressure and temperature

The PAL and PORO lab join forces by combining our respective strengths in laser ultrasound and rock physics to improve data quality and quantity. In this project, we are building the expertise to do laser ultrasound in a pressure vessel with optical windows. Source/receiver locations are varied under computer control with an arduino-controlled servo rotational stage.

quality control of timber and fruit products with laser ultrasound

Laser ultrasound can be applied to products of interest to a wide community. Current methods of testing the quality of fruit and timber, for example, can often be described by one or more of the following terms: sparse, contacting, expensive, and often destructive. In this project, we aim to explore the opportunities for laser ultrasound in estimating the quality of  fruit and timber in a non-contacting and non-destructive matter.

Numerical modelling and experiments in photoacoustic (medical) imaging

The Auckland Volcanic Field (AVF) is intraplate and active. GEONET seismic stations are used to monitor the field.

Stations distribution. Distinguished between those enclosing the AVF and those in the outer ring.

We propose to perform P- (and possibly S-) wave tomography to infer the structure beneath the AVF, with Nick Rawlinson’s FMTT. Our first test is to establish the resolution that we can expect from this exercise. Below, we plot the hypocentres of the earthquakes used in this test:

We compute synthetic seismic P-wave arrival times for the AVF seismic stations, as if the earth structure under the AVF has a checker-board pattern (see images on the left side of the gallery below). For this test, this will be our “true” model of the Earth. And so the objective of FMTT is to invert these synthetic travel times to recover the “true” model checker-board pattern. Areas where the inverted results show a checker-board pattern suggest a good resolution.

Some of the sources need to be removed (17 out of 649 events) in order to run S wave tomography as the code could not determine the traveltimes of these 17 events from their respective sources to the base of the 3D model. Even still, S wave tomography shows a more prominent recovery overall. A major difference can be seen in the EW slice, where a band of positive anomaly cannot be recovered in the P wave equivalent.

Lateral resolution is better in NS rather than in EW orientation in particular at shallower depth. This is expected because the stations are more widely spread in this direction. In both though, resolution of the recovered checker-board is smeared, especially prominent in the lateral extent of the plot.

Directly below the AVF, we can see a pair of checker-board pattern being reflected. The EW slice at -37.0 shows this the best with minimal smearing. The pattern translates horizontally as well as shown in the depth slice. Also, some checker-board pattern persist even to the depth of >200 km.

Summary

• 50-km features seem to be resolved for the very shallow (20km?) and very deep (200km+) parts under the seismic array.
• S wave tomography shows a slighlty higher resolution, especially deeper  and at greater lateral extent.

If we consider closer events we can obtain better resolution at shallower depth. In order to do this we use Nick Rawlinson’s FMTOMO, which results can be seen in this post.