Template write / read merges traces in rare case where traces from same channel align exactly

[flixha]

Describe the bug Consider a template for a 1-component sensor with only a vertical channel, with P- and S-pick on the same trace. If the traces happen to line up exactly, that is, the first sample of the S-trace follows the last sample of the P-trace, then the obspy / mseed library merges the waveforms for the P- and S-part into one trace. The result is a template where some traces suddenly have double the length, and once there is a mix of trace lengths, then of course EQcorrscan will report an error.

To Reproduce see attached template object - trace NS.ARA1.00.BHZ is merged into one trace when reading, even though the template had two traces for NS.ARA1.00.BHZ before writing. Template_sample.gz (.tgz file renamed to .gz due to file type restrictions on Github)

from eqcorrscan.core.match_filter import Template, Tribe
tribe = Tribe().read('Template_sample.gz')
templ = tribe[0]
assert len(list(set([tr.stats.npts for tr in templ.st]))) == 1

Expected behavior The trace should be properly split into the two parts for P and S.

Desktop (please complete the following information):

• Operating System: RHEL 8.3
• Python version: 3.9
• EQcorrscan version: develop

Additional context This problem is not directly caused by EQcorrscan, but rather by the mseed-reading routines. I'm trying to look for some lower-level options that stop how the data are merged into one trace in the miniseed routines, but I don't know the codes so well so this may take some time. I'm not sure yet whether this is an issue of reading or writing the miniseed file. If the mseed library doesn't allow the trace to be split in this case, then we may have to think of a workaround or some extra checks.

Even though I'm continuing to investigate, if you have any ideas, please let me know!

[calum-chamberlain]

Nice spot! That's an ugly bug, thanks for finding it. I think we should probably try to cope with it in here if we can, although I'm not sure how at the moment....

Sorry for the bug!

[flixha]

I see that the record numbers in the mseed-file stillindicate that it was originally two traces:

NS_ARA1_00_BHZ, 000001, D, 512, 56 samples, 20 Hz, 2007,194,09:10:30.027000
NS_ARA1_00_BHZ, 000002, D, 512, 56 samples, 20 Hz, 2007,194,09:10:32.827000
NS_ARA1_00_BHZ, 000003, D, 512, 56 samples, 20 Hz, 2007,194,09:10:35.627000
NS_ARA1_00_BHZ, 000004, D, 512, 56 samples, 20 Hz, 2007,194,09:10:38.427000
...
NS_ARA1_00_BHZ, 000032, D, 512, 56 samples, 20 Hz, 2007,194,09:11:56.827000
NS_ARA1_00_BHZ, 000033, D, 512, 8 samples, 20 Hz, 2007,194,09:11:59.627000
NS_ARA1_00_BHZ, 000001, D, 512, 56 samples, 20 Hz, 2007,194,09:12:00.027000
NS_ARA1_00_BHZ, 000002, D, 512, 56 samples, 20 Hz, 2007,194,09:12:02.827000
NS_ARA1_00_BHZ, 000003, D, 512, 56 samples, 20 Hz, 2007,194,09:12:05.627000
NS_ARA1_00_BHZ, 000004, D, 512, 56 samples, 20 Hz, 2007,194,09:12:08.427000
...
NS_ARA1_00_BHZ, 000033, D, 512, 8 samples, 20 Hz, 2007,194,09:13:29.627000

But I can't see that this information is read into the stream object:

st.select(station='ARA1')[0].stats.mseed
AttribDict({'dataquality': 'D', 'number_of_records': 66, 'encoding': 'FLOAT64',
    'byteorder': '>', 'record_length': 512, 'blkt1001': AttribDict({'timing_quality': 0}),
    'calibration_type': False, 'filesize': 1563648})