SAC Header Initialization in Obspy 1.0.0

Users Mailing List Archive
1

Hi All,

Our group is in the process of upgrading to v1.0.0. Most of our scripts work as expected (with a few warnings), but there is one operation that is causing some strange behaviour.

What we are doing is reading in a miniseed file, performing some operations on the stream, and then writing it out as a sac file. Prior to writing the sac file, we want to add the station coordinates to the sac header. Previously (ie 0.10.2), we have done this as follows:

for tr in st:
tr.stats.sac=AttribDict()
tr.stats.sac.stla = sta.stla
tr.stats.sac.stlo = sta.stlo
tr.write((fname,format=‘SAC’)

This would work fine, with the resulting sac file containing the correct origin time and B values, as well as the station coordinates.

After upgrading to v1.0.0, the same operation no produces an error when writing the sac file.
“UserWarning: Old header has invalid reftime.”

The write operation will complete successfully, and the SAC file is written to disk. However, the starttime and endtime are no longer correct, with the start time set to a default value of “1970-01-01T00:00:00.000000Z”.

For instance, reading in the sac file that was just written and viewing the stats gives:
trnew=read(fname,format=‘SAC’)
print trnew.stats
network: BP
station: CCRB
location:
channel: BP1
starttime: 1970-01-01T00:00:00.000000Z
endtime: 1970-01-01T23:59:59.950000Z
sampling_rate: 20.0
delta: 0.05
npts: 1728000
calib: 1.0
_format: SAC
sac: AttribDict({u’lpspol’: 0, u’iftype’: 1, u’stla’: 35.957199, u’leven’: 1, u’nvhdr’: 6, u’stlo’: -120.5516, u’depmen’: -12812.165, u’depmin’: -14073.0, u’kstnm’: u’CCRB’, u’depmax’: -10348.0, u’npts’: 1728000, u’lovrok’: 1, u’knetwk’: u’BP’, u’delta’: 0.050000001, u’e’: nan, u’lcalda’: 1, u’kcmpnm’: u’BP1’})

The stla and stlo are in the header as they should be, but the start and end times are incorrect.

On the other hand, if we don’t set the station coordinates, ie, do not initialize the AttribDict() in tr.stats.sac, then the trace writes out correctly with no warnings, and reading the trace in again shows it has the correct start and end times.
trnew2=read(fname,format=‘SAC)
print trnew2.stats
network: BP
station: CCRB
location:
channel: BP1
starttime: 2003-01-01T00:00:00.023800Z
endtime: 2003-01-01T23:59:59.973800Z
sampling_rate: 20.0
delta: 0.05
npts: 1728000
calib: 1.0
_format: SAC
sac: AttribDict({u’knetwk’: u’BP’, u’nzyear’: 2003, u’nzjday’: 1, u’iztype’: 9, u’nzhour’: 0, u’lcalda’: 1, u’iftype’: 1, u’nvhdr’: 6, u’depmin’: -14073.0, u’kcmpnm’: u’BP1’, u’nzsec’: 0, u’depmen’: -12812.165, u’depmax’: -10348.0, u’lovrok’: 1, u’scale’: 1.0, u’delta’: 0.050000001, u’nzmsec’: 23, u’lpspol’: 0, u’b’: 0.00079999998, u’e’: 86399.953, u’leven’: 1, u’kstnm’: u’CCRB’, u’nzmin’: 0, u’npts’: 1728000})

But now we don’t have the station coordinates (or other variables we may wish to populate.

My question is, how do we go about adding metadata to the SAC header fields to a trace that is read in (either using fdsn client or from a local miniseed file), and write the file out as a SAC file?

Thanks very much,
Andrew

2

Hello Andrew,

This is a known issue with the re-written SAC module (see SAC reference time not automatically set for incomplete stats.sac header by jkmacc-LANL · Pull Request #1285 · obspy/obspy · GitHub), and it is being worked on right now. This behavior is due to changes in the way the ObsPy and SAC headers are merged before writing to file, which fixed some things and unfortunately broke others. In the meantime, you can employ the SACTrace class directly to make the manipulations you desire:

from obspy.io.sac import SACTrace

for tr in st:
    sac = SACTrace.from_obspy_trace(tr)
    sac.stla = sta.stla
    sac.stlo = sta.stlo
    sac.write(fname)

Of course, this is just a short-term fix, and you can expect the correct behavior with tr.write soon. Thanks for your patience while the fix is made.

Best,
Jon

3

Hi Jon,

Thanks very much for pointing out the link and suggesting the work around.

I have tried both the workaround you suggest as well as the workaround suggested by Tom.

I had actually tried the workaround you suggested already, but we some strange problems.
Specifically, it works correctly on one of our computers, but produces an error on another.
In the meantime we plan to use obspy_to_sac_header (from obspy.io.sac.util) until the expected behaviour in trace.write() is back.

For completeness sake, below I outline the error we get with the SACTrace.from_obspy_trace method.

Computer #1 is running Windows 10 (x86_64) with python3.4 and obspy v1.0.0. The SACTrace.from_obspy_trace() works as expected.

Computer #2 is running Ubuntu 14.04 (x86_64) with python2.7 and obspy v1.0.0. The SACTrace.from_obspy_trace() does not work.

On both computers, we assign the following values in ipython (both of which are in range) onto a SACTrace trace:
sactrace.stla = 35.45

sactrace.stlo = -100.45

On the Ubuntu machine this results in the following error:
“Latitude of Point 1 out of bounds! (-90 <= lat1 <=90)”

This error prevents a script from continuing, though if running in ipython, we can see that the value is actually assigned despite the error.

Based on tracking the error, it is actually coming from obspy.geodetics.base.gps2_distazimuth.

I’m not sure why though, it produces an error on one machine and not the other.

Thanks again for the help.

Cheers,
Andrew

4

For the sake of completeness, another workaround working directly with the Trace object:

from obspy.io.sac.util import obspy_to_sac_header

for tr in st:
    tr.stats.sac = obspy_to_sac_header(tr.stats)
    tr.stats.sac.stla = sta.stla
    tr.stats.sac.stlo = sta.stlo
    tr.write(fname, format='SAC')

Best, Tom

5

Hmmm. That "out of bounds” error sounds like an issue in its own right. It may be a good idea to file an issue for this. Thanks for mentioning it.

-Jon