search

jpjones76 / SeisIO.jl

Julia language support for geophysical time series data
http://seisio.readthedocs.org
Other
47 stars 20 forks source link

Parsing gaps correctly? #36

Closed interseismic closed 4 years ago

interseismic commented 4 years ago

I'm working with gappy MSEED data and am unsure if something I'm seeing is expected behavior or not:

julia> st[1].t

21×2 Array{Int64,2}:
        1  1420070400000000
  4208601          13000000
  4557601           6000000
  4815601          19000000
  5036201          16000000
  5719401          12000000
  5812401          34000000
  5941001           4000000
  6007601           2000000
  8597801          23000000
  8598801          48000000
  8635601          48000000
  8717201          17000000
  8721201          15000000
  8732601           4000000
  9000801          13000000
  9254201          49000000
  9412001           4000000
  9555601          24000000
  9556201          58000000
 17198200                 0

Assuming I'm understanding the structure of this correctly, the first gap is 13 sec long. If I read the same file in Python:

st.print_gaps()
Source            Last Sample                 Next Sample                 Delta           Samples 
WR.PEC..EHZ       2015-01-01T05:50:41.995000Z 2015-01-01T05:50:56.000000Z 14.000000       2800    
WR.PEC..EHZ       2015-01-01T06:20:00.995000Z 2015-01-01T06:20:07.000000Z 6.000000        1200    
WR.PEC..EHZ       2015-01-01T06:41:36.995000Z 2015-01-01T06:41:56.000000Z 19.000000       3800    
WR.PEC..EHZ       2015-01-01T07:00:18.995000Z 2015-01-01T07:00:35.000000Z 16.000000       3200    
WR.PEC..EHZ       2015-01-01T07:57:30.995000Z 2015-01-01T07:57:43.000000Z 12.000000       2400    
WR.PEC..EHZ       2015-01-01T08:05:27.995000Z 2015-01-01T08:06:02.000000Z 34.000000       6800    
WR.PEC..EHZ       2015-01-01T08:16:44.995000Z 2015-01-01T08:16:49.000000Z 4.000000        800     
WR.PEC..EHZ       2015-01-01T08:22:21.995000Z 2015-01-01T08:22:24.000000Z 2.000000        400     
WR.PEC..EHZ       2015-01-01T11:58:14.995000Z 2015-01-01T11:58:38.000000Z 23.000000       4600    
WR.PEC..EHZ       2015-01-01T11:58:42.995000Z 2015-01-01T11:59:31.000000Z 48.000000       9600    
WR.PEC..EHZ       2015-01-01T12:02:34.995000Z 2015-01-01T12:03:23.000000Z 48.000000       9600    
WR.PEC..EHZ       2015-01-01T12:10:10.995000Z 2015-01-01T12:10:28.000000Z 17.000000       3400    
WR.PEC..EHZ       2015-01-01T12:10:47.995000Z 2015-01-01T12:11:03.000000Z 15.000000       3000    
WR.PEC..EHZ       2015-01-01T12:11:59.995000Z 2015-01-01T12:12:04.000000Z 4.000000        800     
WR.PEC..EHZ       2015-01-01T12:34:24.995000Z 2015-01-01T12:34:38.000000Z 13.000000       2600    
WR.PEC..EHZ       2015-01-01T12:55:44.995000Z 2015-01-01T12:56:34.000000Z 49.000000       9800    
WR.PEC..EHZ       2015-01-01T13:09:42.995000Z 2015-01-01T13:09:46.000000Z 3.000000        600     
WR.PEC..EHZ       2015-01-01T13:21:44.995000Z 2015-01-01T13:22:09.000000Z 24.000000       4800    
WR.PEC..EHZ       2015-01-01T13:22:11.995000Z 2015-01-01T13:23:10.000000Z 58.000000       11600

The first gap in ObsPy is instead listed as 14 sec, but the rest match the SeisIO durations. Am I missing something here? WR.PEC.EHZ.mseed.zip

jpjones76 commented 4 years ago

What's the sampling frequency?

interseismic commented 4 years ago

200 Hz

jpjones76 commented 4 years ago

I get the same gaps that are shown in ObsPy, so I'm not sure why we obtain such differences.

Below, t_obspy is a String created by copy-pasting your ObsPy output. 

using DelimitedFiles, SeisIO
S = read_data("mseed", "/home/josh/Downloads/WR.PEC.EHZ.mseed")
t1 = S.t[1]
buf = IOBuffer(t_obspy)
t2_raw = readdlm(buf)
close(buf)
t2 = round.(Int64, t2_raw[:,4].*1000000)

The result:

julia> [t1[2:end-1,2] .- t2]
1-element Array{Array{Int64,1},1}:
 [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]

If I load in your time matrix in the original comment as t, there are significant differences.

julia> [(t-t1)[2:end-1,2]] # enclosing in arrays for more compact display
1-element Array{Array{Int64,1},1}:
 [-1000000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1000000, 0, 0]

I can't tell this is an error in an old SeisIO version, or something potentially very serious, like machine-dependent behavior. Just to cover all bases, so to speak, could you please tell me the following additional information?

interseismic commented 4 years ago

Ok my bad. Accidentally was looking at EHN in SeisIO, but EHZ in ObsPy. Everything matches up now -- sorry about wasting your time on this.