removing instrument response

Dear users,

I am trying to remove the instrument response from a seed file I obtained from Breq_fast. I can extract the poles and zeros with rdseed, but I would like to do that in Python to not call a different program first and then run my program in Python afterwards. I tried

seismogram_orig = seismogram.copy()
pre_filt = (0, 0.006, 30.0, 35.0)
t1 = UTCDateTime("2010-09-3T16:30:00.000")
date = t1
respf = "RESP.G.TAM.00.VHZ"
seedresp = {'filename': respf,
               'date': date,
               'units': 'VEL'
          }
seismogram.simulate(paz_remove=None, pre_filt=pre_filt, seedresp=seedresp)

but the resulting seismogram is certainly not correct and also the FFT gives weird results. I attach my code file, the seed file and the response files from rdseed. I would be happy if someone could explain me how to handle seed files in ObsPy…

Regards, Boris

instrresp.py (3.08 KB)

RESP.G.TAM.00.VHZ (52.6 KB)

SAC_PZs_G_TAM_VHZ_00_2006.351.23.30.00.0000_2014.032.02.60.60.99999 (1 KB)

Tohoku.seed (224 KB)

Hi Boris,

please note that ObsPy has better support for the more modern StationXML format so you might consider just getting the data via the fdsn web services:

Other than that your basic approach is correct. You might want to have a look at the bode plot of the response and make sure your pre-filter only selects data within the instrument’s valid frequency band. Your data is sampled at a 10 second interval so upper frequency limits of 30-35 Hz in the pre-filter will not do anything.

I also have little experience working with really long period data but you might want to play around with tapering and linear (or higher-order) detrending of the data and just do in manually before calling tr.simulate(). Parameters in it are passed on to� this function here:

Two comments on your code:

  • ObsPy has Stream/Trace.plot() methods - you don’t need to run your own for quick and dirty plots.

  • NumPy has rfft() and rfftfreqs() functions which will not return the negative frequencies for real valued data - then you don’t have to jump through all the hoops to plot the spectra. Or just cut away the negative frequencies.

Hope it helps!

Lion

One Addition to Lion's detailed reply:

Please also note that if you use Trace/Stream.remove_response() with an
Inventory object, you'll have the option to plot the spectra during the
various steps in response removal:

https://docs.obspy.org/packages/autogen/obspy.core.trace.Trace.remove_response.html?highlight=plot#obspy.core.trace.Trace.remove_response

https://docs.obspy.org/tutorial/code_snippets/seismometer_correction_simulation.html?highlight=plot#using-a-stationxml-file-or-in-general-an-inventory-object

T