begin time "b" in SAC is not handled correctly

[taotaokai]

Hi,

On line 101 of SAC.jl instead of adding "b" time as milliseconds

ts += round(Int64, b*1.0f3),

it should be added as seconds

ts += round(Int64, b*1.0f6).

[taotaokai]

Also on line 48 of SAC.jl it should be BUF.sac_fv[6] = rem(ts, 1000)*1.0f-6.

[jpjones76]

Really? Is the SAC manual wrong? It clearly states "GMT millisecond", not "GMT microsecond"...

http://ds.iris.edu/files/sac-manual/manual/file_format.html

[taotaokai]

Thank you! I do not know if it worths your time, but below is how I found the problem:

1. Create a sample seismogram "seis.sac" in SAC:
   funcgen seis
   write seis.sac
   quit

2. In Julia read in this sample sac file and write out a new file as "seis2.sac":
   using SeisIO
   s = read_data("sac", "seis.sac")
   writesac(s[1], "seis2.sac")

3. compare the SAC headers using saclst (a utility progam from SAC package)
   $ saclst npts delta b e kztime f seis.sac seis2.sac
   seis.sac 1000 0.01 9.46 19.45 10:38:14.000
   seis2.sac 1000 0.01 0.46 10.46 10:38:14.009

Issues:

1. the seismogram has a duration of 0.01*(1000-1)=9.99 seconds. Since seis.sac is generated by SAC, "b" and "e" should mean seconds. However the "e" time is not set correclty in seis2.sac (should be 10.45);
2. in seis2.sac SeisIO moves 9 seconds from "b", but adds to kztime as 9 milliseconds.

I did the following modifications:

1. SAC.jl line 101: ts is in microsecond, so adding "b" to ts should be ts += round(Int64, b*1.0f6)
2. SAC.jl line 48: for the same reason, the sub-millisecond residual time should be put in "b" as BUF.sac_fv[6] = rem(ts, 1000)*1.0f-6
3. SAC.jl line 49: the "e" time should be calculated as BUF.sac_fv[7] = BUF.sac_fv[6] + BUF.sac_fv[1]*(nx-1)

Then I ran the above test again and got the correct output:
$ saclst npts delta b e kztime f seis.sac seis2.sac
seis.sac 1000 0.01 9.46 19.45 10:38:14.000
seis2.sac 1000 0.01 0.000999 9.991 10:38:23.460

[jpjones76]

Hi again. Thanks for the explanation. I'm currently looking for additional test data to investigate this issue and your proposed fix.

[jpjones76]

Thank you! I do not know if it worths your time, but below is how I found the problem:

Hi again. I'm going to push a fix live tonight, based on this comment, with my thanks.

I think you are correct that my SAC writer incorrectly handled a non-null begin time. The solution should lead to perfectly re-readable data. Writing the SAC test seismogram no longer produces any time shifts, even subsample ones. I verified that the issue was fixed using SeisIO.jl master branch and SAC v101.6a; it's now part of my automated tests.

One reason this took so long to resolve is that I could find no real test data with b != 0.0f0. I'm not 100% confident that the test file generated by SAC funcgen reflects how "b" and "e" are set in real data. If other users report that these changes break their code, I might need to make adjustments.

Thank you very much!

[jpjones76]

Fixed in SeisIO v1.2.0.