Skip to content
Snippets Groups Projects
Commit d03495b2 authored by Kien Le's avatar Kien Le
Browse files

Reimplement Reftek130.to_stream()

parent 3ebf39b2
No related branches found
No related tags found
No related merge requests found
...@@ -46,7 +46,6 @@ class Reftek130Exception(ObsPyException): ...@@ -46,7 +46,6 @@ class Reftek130Exception(ObsPyException):
class Reftek130(obspy_rt130_core.Reftek130): class Reftek130(obspy_rt130_core.Reftek130):
def to_stream(self, network="", location="", component_codes=None, def to_stream(self, network="", location="", component_codes=None,
include_mp123=False, include_mp456=False, include_mp123=False, include_mp456=False,
headonly=False, verbose=False, headonly=False, verbose=False,
...@@ -56,8 +55,6 @@ class Reftek130(obspy_rt130_core.Reftek130): ...@@ -56,8 +55,6 @@ class Reftek130(obspy_rt130_core.Reftek130):
:param headonly: Determines whether or not to unpack the data or just :param headonly: Determines whether or not to unpack the data or just
read the headers. read the headers.
""" """
if verbose:
print(self)
if not len(self._data): if not len(self._data):
msg = "No packet data in Reftek130 object (file: {})" msg = "No packet data in Reftek130 object (file: {})"
raise Reftek130Exception(msg.format(self._filename)) raise Reftek130Exception(msg.format(self._filename))
...@@ -99,20 +96,6 @@ class Reftek130(obspy_rt130_core.Reftek130): ...@@ -99,20 +96,6 @@ class Reftek130(obspy_rt130_core.Reftek130):
eh = EHPacket(eh_packets[0]) eh = EHPacket(eh_packets[0])
else: else:
eh = EHPacket(et_packets[0]) eh = EHPacket(et_packets[0])
# only C0, C2, 16, 32 encodings supported right now
if eh.data_format == b"C0":
encoding = 'C0'
elif eh.data_format == b"C2":
encoding = 'C2'
elif eh.data_format == b"16":
encoding = '16'
elif eh.data_format == b"32":
encoding = '32'
else:
msg = ("Reftek data encoding '{}' not implemented yet. Please "
"open an issue on GitHub and provide a small (< 50kb) "
"test file.").format(eh.data_format)
raise NotImplementedError(msg)
header = { header = {
"unit_id": self._data['unit_id'][0], "unit_id": self._data['unit_id'][0],
"experiment_number": self._data['experiment_number'][0], "experiment_number": self._data['experiment_number'][0],
...@@ -158,45 +141,12 @@ class Reftek130(obspy_rt130_core.Reftek130): ...@@ -158,45 +141,12 @@ class Reftek130(obspy_rt130_core.Reftek130):
sample_data = np.array([], dtype=np.int32) sample_data = np.array([], dtype=np.int32)
npts = packets_["number_of_samples"].sum() npts = packets_["number_of_samples"].sum()
else: else:
if encoding in ('C0', 'C2'): sample_data = (packets_['payload'][:, :4])
sample_data = _unpack_C0_C2_data(packets_, sample_data = sample_data.view(np.dtype('>i4')).squeeze()
encoding)
elif encoding in ('16', '32'):
# rt130 stores in big endian
dtype = {'16': '>i2', '32': '>i4'}[encoding]
# just fix endianness and use correct dtype
sample_data = np.require(
packets_['payload'],
requirements=['C_CONTIGUOUS'])
# either int16 or int32
sample_data = sample_data.view(dtype)
# account for number of samples, i.e. some packets
# might not use the full payload size but have
# empty parts at the end that need to be cut away
number_of_samples_max = sample_data.shape[1]
sample_data = sample_data.flatten()
# go through packets starting at the back,
# otherwise indices of later packets would change
# while looping
for ind, num_samps in reversed([
(ind, num_samps) for ind, num_samps in
enumerate(packets_["number_of_samples"])
if num_samps != number_of_samples_max]):
# looping backwards we can easily find the
# start of each packet, since the earlier
# packets are still untouched and at maximum
# sample length in our big array with all
# packets
start_of_packet = ind * number_of_samples_max
start_empty_part = start_of_packet + num_samps
end_empty_part = (start_of_packet +
number_of_samples_max)
sample_data = np.delete(
sample_data,
slice(start_empty_part, end_empty_part))
npts = len(sample_data) npts = len(sample_data)
tr = Trace(data=sample_data, header=copy.deepcopy(header)) tr = Trace(data=sample_data, header=copy.deepcopy(header))
tr.stats.npts = packets_['number_of_samples'].sum()
# channel number is not included in the EH/ET packet # channel number is not included in the EH/ET packet
# payload, so add it to stats as well.. # payload, so add it to stats as well..
tr.stats.reftek130['channel_number'] = channel_number tr.stats.reftek130['channel_number'] = channel_number
...@@ -236,22 +186,5 @@ class Reftek130(obspy_rt130_core.Reftek130): ...@@ -236,22 +186,5 @@ class Reftek130(obspy_rt130_core.Reftek130):
continue continue
tr.stats.channel = "MassPos%s" % (channel_number + 1) tr.stats.channel = "MassPos%s" % (channel_number + 1)
# check if endtime of trace is consistent # check if endtime of trace is consistent
t_last = packets_[-1]['time']
npts_last = packets_[-1]['number_of_samples']
try:
if not headonly:
assert npts == len(sample_data)
if npts_last:
assert tr.stats.endtime == UTCDateTime(
ns=t_last) + (npts_last - 1) * delta
if npts:
assert tr.stats.endtime == (
tr.stats.starttime + (npts - 1) * delta)
except AssertionError:
msg = ("Reftek file has a trace with an inconsistent "
"endtime or number of samples. Please open an "
"issue on GitHub and provide your file for"
"testing.")
raise Reftek130Exception(msg)
st += tr st += tr
return st return st
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment