Maybe we would like to see the body waves more clearly and filter out
the surface waves. And instead of the raw velocity data, lets
apply the full instrument response to correct to displacement.
First, lets switch from the LH? channels, recorded
at one sample per second, to the HH? channels, recorded at 100 sps
to give us a wider frequency band to play with. Since we will need
the response information as well, we will need to change from
using queryChannels()
to queryResponses()
in the stationQuery
. This
will return stationXML with the full response filled in for each
channel.
let stationQuery = new seisplotjs.fdsnstation.StationQuery()
.networkCode('CO')
.stationCode('HODGE')
.locationCode('00')
.channelCode('HH?')
.timeWindow(queryTimeWindow);
Promise.all( [ eventQuery.query(), stationQuery.queryResponses() ] )
Of course 100 sps means that we will have 100 times as many samples, so lets reduce the time window to just the region where the P arrival is, say 300 seconds instead of 1800.
let timeWindow = new seisplotjs.util.StartEndDuration(startTime, null, 300);
Now we insert the filtering code after the seismograms have arrived
but before we plot them. We will insert another
then()
call even though
we could just add to the existing one. We create a butterworth
filter using the sampling rate of the seismogram, with a passband
of .5 to 10 Hz. Removing the mean is usually a good idea, then
we apply the filter. Tapering is important before a deconvolution
and then we transfer
the instrument response for the channel.
}).then( ( [ quakeList, networks, ttimes, seismogramDataList ] ) => {
seismogramDataList.forEach(sdd => {
let butterworth = seisplotjs.filter.createButterworth(
2, // poles
seisplotjs.filter.BAND_PASS,
.5, // low corner
10, // high corner
1/sdd.seismogram.sampleRate // delta (period)
);
let rmeanSeis = seisplotjs.filter.rMean(sdd.seismogram);
let filteredSeis = seisplotjs.filter.applyFilter(butterworth, rmeanSeis);
let taperSeis = seisplotjs.taper.taper(filteredSeis);
let correctedSeis = seisplotjs.transfer.transfer(taperSeis,
sdd.channel.response, .001, .02, 250, 500);
sdd.seismogram = correctedSeis;
});
return Promise.all( [ quakeList, networks, ttimes, seismogramDataList ] );
Then, just to make sure we don't correct the data twice, we disable the gain correction.
seisConfig.doGain = false;
Perhaps just for fun we can also plot the amplitude spectra for the three seismograms. We need to add an additional div to hold them.
<div id="fftplot">
</div>
And an additional style to size it.
div#fftplot {
height: 600px;
}
Then we calculate the fft and plot it. Because the resulting plot is just an svg element created with d3, we can add a plot label by modifying it with normal d3 techniques.
}).then( ( [ quakeList, networks, ttimes, seismogramDataList ] ) => {
let div = seisplotjs.d3.select('div#fftplot');
let fftList = seismogramDataList.map(sdd => seisplotjs.fft.fftForward(sdd.seismogram));
let svg = seisplotjs.fftplot.createOverlayFFTPlot('div#fftplot', fftList, true);
svg.append("g").classed("title", true)
.attr("transform", "translate(700, 100)")
.append("text").classed("title label", true)
.selectAll("tspan")
.data(seismogramDataList)
.enter()
.append("tspan")
.text(sdd => " "+sdd.seismogram.channelCode+" ");
return Promise.all( [ quakeList, networks, ttimes, seismogramDataList ] );
Previous: Predicted phase arrival times
Next: Helicorder