stf deconvolution

[martinvandriel]

because when the srf reader was removed, there was no stf reader anymore. should be reactivated after implementing the new stf reader,

[sstaehler]

@martinvandriel Is that still an issue?

[martinvandriel]

no idea. Do we have stfs in the kerner now?

[sstaehler]

no

[martinvandriel]

don't we need them?

[auerl]

Optimally, we should convolve the fwd field with the stf (if available), but i wonder if it would result in a huge difference in the kernel? Probably not the highest priority, I would guess.

[sstaehler]

I hope I fixed the STF deconvolution in 3dcf9df2dd5f4ee7e045c1e4a3a85c0ff0a06289 On convolving with the "real" STF: It's trivial, if we have the "real" STF with the correct sampling rate. If not, should we use resample_stf (https://github.com/sstaehler/kerner/blob/master/source.f90#L216)? How does that even work? It's just a linear interpolation, isn't it?

[martinvandriel]

True. So many lines for a linear interpolation, great. If the sampling is far enough from Nyquist, this should be fine. Otherwise, I would consider using the lanczos from intstaseis: https://github.com/krischer/instaseis/blob/master/instaseis/lanczos.py

[sstaehler]

Yes, I just discovered that the core of the lanczos interpolation in instaseis is in Fortran anyway and I'm moving that to the Kerner.

Is there a default file format for STFs in instaseis?

[sstaehler]

Started to work on it 0b69d9bf1a9de7ad111952bd23ca2513039e688f

I would remove the part where multiple sources are possible. Or do we want to have kernels for finite faults?

[sstaehler]

Okay, STF deconvolution and reconvolution is basically done and I can reproduce Instaseis seismograms in the Kerner. For the user, two big changes

• he/she has to supply a file with an event STF all the time
• the filters and the STF have to be chosen more carefully. The STF deconvolution can easily be unstable if neither filter nor the event STF decrease steeply enough towards higher frequencies. For Karin's Gabor filters, I would suggest that their central period has to be at least 2 times the mesh period, since they are extremely broad.

[martinvandriel]

Why not use an additional lowpassfilter to ensure stability? Just the same one for all kernels. Needs to be applied in the measurement as well of course. Butterworth filters should be readily available, so this seems to be mostly an additional input parameter...

[kasra-hosseini]

From the measurement point of view, here are the filtering steps:

1. Bandpass filter between 0.012 - 3 Hz...in fact, the corner frequencies are 0.008, 0.012, 3, 4, so it is one between 0.012 and 3 ---> this will be applied to the real data
2. re-sampling to 10Hz, since we already applied the filter in the first stage, no aliasing is expected.
3. Convolution with the STF
4. Gabor filters on both real and synthetic waveforms, and then doing the measurement.

I know that it is clear, just to summarize what we do in the measurement.

[sstaehler]

@martinvandriel Yes, I actually implemented this additional order 16 lowpass filter in https://github.com/sstaehler/kerner/blob/master/filtering.f90#L298 and most of the times it solves the problem. The combined filter response is written out, so that it can be used to do the measurements accordingly. Problem is that even that is sometimes not enough, since the spectrum of the Gaussian STF in AxiSEM falls with exp(-f^2), while the Butterworth is only f^(-2*n). Anyway, in these cases, the code crashes and the user has to do something.

@kasra-hosseini So, what we would need to add in MCKernel is the highpass filter at 100s, What we would need to change in the measurements is that we replace the Gabor by Gabor+order-16-Butterworth at the Mesh period.

[martinvandriel]

since the spectrum of the Gaussian STF in AxiSEM falls with exp(-f^2)

so do the Gabor filters...

[sstaehler]

Actually with exp( -(ln(f-fc))^2 ) That might be the problem, they are incredibly wide on the high frequency side

[martinvandriel]

What seems to work in instaseis is to deconvolve the gaussian from axisem and reconvolve with the response of a butterworth lowpass (this is how I did the benchmarks in the instaseis paper). So I don't see, why it would be unstable in the kerner case you discribe above, which should be more stable by having an additional bandlimited stf and additional Gabor filter.

[sstaehler]

Okay, I'll check that.