From abc2c8352a33a51208880aebf2ea508944df53ad Mon Sep 17 00:00:00 2001 From: kienle <kienle@passcal.nmt.edu> Date: Thu, 3 Aug 2023 15:51:01 -0600 Subject: [PATCH] Reorganize RT130 read functions --- .../model/reftek/from_rt2ms/core.py | 2 +- sohstationviewer/model/reftek/reftek.py | 2 +- .../reftek/rt130_experiment/dt_packet.py | 3 +- .../reftek/rt130_experiment/eh_et_packet.py | 3 +- .../model/reftek/rt130_experiment/header.py | 94 +++++++++ .../model/reftek/rt130_experiment/reftek.py | 196 ------------------ .../reftek/rt130_experiment/reftek_helper.py | 128 +++++++++++- .../reftek/rt130_experiment/soh_packets.py | 3 +- 8 files changed, 216 insertions(+), 215 deletions(-) create mode 100644 sohstationviewer/model/reftek/rt130_experiment/header.py delete mode 100644 sohstationviewer/model/reftek/rt130_experiment/reftek.py diff --git a/sohstationviewer/model/reftek/from_rt2ms/core.py b/sohstationviewer/model/reftek/from_rt2ms/core.py index 525522a9d..b94b84536 100644 --- a/sohstationviewer/model/reftek/from_rt2ms/core.py +++ b/sohstationviewer/model/reftek/from_rt2ms/core.py @@ -24,7 +24,7 @@ from obspy.io.reftek.util import _decode_ascii, _parse_long_time from sohstationviewer.model.mseed_data.record_reader_helper import Unpacker from sohstationviewer.model.reftek.from_rt2ms import packet from sohstationviewer.model.reftek.from_rt2ms.soh_packet import Packet -from sohstationviewer.model.reftek.rt130_experiment.reftek import \ +from sohstationviewer.model.reftek.rt130_experiment.reftek_helper import \ ( read_rt130_file, convert_packet_to_obspy_format, ) diff --git a/sohstationviewer/model/reftek/reftek.py b/sohstationviewer/model/reftek/reftek.py index df995d7ec..9d0dc8cac 100755 --- a/sohstationviewer/model/reftek/reftek.py +++ b/sohstationviewer/model/reftek/reftek.py @@ -7,7 +7,7 @@ from typing import Tuple, List, Union import traceback import numpy as np -from sohstationviewer.model.reftek.rt130_experiment.reftek import Reftek130 +from sohstationviewer.model.reftek.rt130_experiment.header import Reftek130 from sohstationviewer.model.reftek.from_rt2ms import ( core, soh_packet, packet) from sohstationviewer.model.reftek.log_info import LogInfo diff --git a/sohstationviewer/model/reftek/rt130_experiment/dt_packet.py b/sohstationviewer/model/reftek/rt130_experiment/dt_packet.py index e4357d1cc..a433dcd9b 100644 --- a/sohstationviewer/model/reftek/rt130_experiment/dt_packet.py +++ b/sohstationviewer/model/reftek/rt130_experiment/dt_packet.py @@ -3,8 +3,7 @@ from dataclasses import dataclass from sohstationviewer.model.mseed_data.record_reader_helper import \ Unpacker -from sohstationviewer.model.reftek.rt130_experiment.reftek_helper import \ - PacketHeader +from sohstationviewer.model.reftek.rt130_experiment.header import PacketHeader def decode_uncompressed(packet: bytes, data_format: str, unpacker: Unpacker): diff --git a/sohstationviewer/model/reftek/rt130_experiment/eh_et_packet.py b/sohstationviewer/model/reftek/rt130_experiment/eh_et_packet.py index a91039767..6e043bfec 100644 --- a/sohstationviewer/model/reftek/rt130_experiment/eh_et_packet.py +++ b/sohstationviewer/model/reftek/rt130_experiment/eh_et_packet.py @@ -3,8 +3,7 @@ import dataclasses from sohstationviewer.model.mseed_data.record_reader_helper import \ Unpacker from sohstationviewer.model.reftek.from_rt2ms.core import eh_et_payload_end_in_packet -from sohstationviewer.model.reftek.rt130_experiment.reftek_helper import \ - PacketHeader +from sohstationviewer.model.reftek.rt130_experiment.header import PacketHeader def read_eh_et_packet(packet: bytes, unpacker: Unpacker): diff --git a/sohstationviewer/model/reftek/rt130_experiment/header.py b/sohstationviewer/model/reftek/rt130_experiment/header.py new file mode 100644 index 000000000..2bd346897 --- /dev/null +++ b/sohstationviewer/model/reftek/rt130_experiment/header.py @@ -0,0 +1,94 @@ +import dataclasses + +from obspy import UTCDateTime + +from sohstationviewer.model.mseed_data.record_reader_helper import \ + Unpacker +from sohstationviewer.model.reftek.rt130_experiment.reftek_helper import ( + RT130ParseError, +) + +from sohstationviewer.model.reftek.from_rt2ms import core + + +@dataclasses.dataclass +class PacketHeader: + packet_type: str + experiment_number: int + unit_id: str + time: UTCDateTime + byte_count: int + packet_sequence: int + + +def parse_rt130_time(year: int, time_bytes: bytes) -> UTCDateTime: + time_string = time_bytes.hex() + # The time string has the format of DDDHHMMSSTTT, where + # D = day of year + # H = hour + # M = minute + # S = second + # T = microsecond + day_of_year, hour, minute, second, millisecond = ( + int(time_string[0:3]), + int(time_string[3:5]), + int(time_string[5:7]), + int(time_string[7:9]), + int(time_string[9:12]) + ) + # RT130 only stores the last two digits of the year. Because the + # documentation for RT130 does not define a way to retrieve the full year, + # we use Obspy's method. Accordingly, we convert 0-49 to 2000-2049 and + # 50-99 to 1950-1999. + if 0 <= year <= 49: + year += 2000 + elif 50 <= year <= 99: + year += 1900 + converted_time = UTCDateTime(year=year, julday=day_of_year, hour=hour, + minute=minute, second=second, + microsecond=millisecond * 1000) + return converted_time + + +def get_rt130_packet_header(rt130_packet: bytes, + unpacker: Unpacker) -> PacketHeader: + try: + # Because RT130 data is always big-endian, it is more convenient to + # use str.decode() than the unpacker. + packet_type = rt130_packet[:2].decode('ASCII') + except UnicodeError: + print(f'Cannot decode packet type.') + print('The given file does not appear to be a valid RT130 file.') + raise RT130ParseError + valid_packet_types = ['AD', 'CD', 'DS', 'DT', 'EH', 'ET', 'OM', 'SH', 'SC', + 'FD'] + if packet_type not in valid_packet_types: + print(f'Invalid packet type found: {packet_type}') + print('The given file does not appear to be a valid RT130 file.') + raise RT130ParseError + + experiment_number = int(rt130_packet[2:3].hex()) + year = int(rt130_packet[3:4].hex()) + # A call to str.upper() is needed because bytes.hex() makes any + # hexadecimal letter (i.e. ABCDEF) lowercase, while we want them to be + # uppercase for display purpose. + unit_id = rt130_packet[4:6].hex().upper() + time_bytes = rt130_packet[6:12] + packet_time = parse_rt130_time(year, time_bytes) + byte_count = int(rt130_packet[12:14].hex()) + packet_sequence = int(rt130_packet[14:16].hex()) + + return PacketHeader(packet_type, experiment_number, unit_id, packet_time, + byte_count, packet_sequence) + + +class Reftek130(core.Reftek130): + pass + + +waveform_file = '/Users/kle/PycharmProjects/sohstationviewer/tests/test_data/RT130-sample/2017149.92EB/2017150/92EB/1/000000015_0036EE80' +soh_file = '/Users/kle/PycharmProjects/sohstationviewer/tests/test_data/RT130-sample/2017149.92EB/2017150/92EB/0/000000000_00000000' +import time +start = time.perf_counter() +(Reftek130.from_file(waveform_file).to_stream()) +# print('Time taken:', time.perf_counter() - start) diff --git a/sohstationviewer/model/reftek/rt130_experiment/reftek.py b/sohstationviewer/model/reftek/rt130_experiment/reftek.py deleted file mode 100644 index 34b3785e5..000000000 --- a/sohstationviewer/model/reftek/rt130_experiment/reftek.py +++ /dev/null @@ -1,196 +0,0 @@ -import os -from typing import Callable, Dict, Union - -import numpy -import numpy as np -from obspy import UTCDateTime - -from sohstationviewer.model.mseed_data.record_reader_helper import \ - Unpacker -from sohstationviewer.model.reftek.rt130_experiment.dt_packet import \ - ( - read_dt_packet, DTPacket, -) -from sohstationviewer.model.reftek.rt130_experiment.eh_et_packet import \ - ( - EHETPacket, read_eh_et_packet, eh_et_payload_end_in_packet, -) -from sohstationviewer.model.reftek.rt130_experiment.reftek_helper import ( - packet_reader_placeholder, RT130ParseError, PacketHeader, -) - -from sohstationviewer.model.reftek.from_rt2ms import core -from sohstationviewer.model.reftek.rt130_experiment.soh_packets import \ - ( - read_soh_packet, SOHPacket, -) - - -def parse_rt130_time(year: int, time_bytes: bytes) -> UTCDateTime: - time_string = time_bytes.hex() - # The time string has the format of DDDHHMMSSTTT, where - # D = day of year - # H = hour - # M = minute - # S = second - # T = microsecond - day_of_year, hour, minute, second, millisecond = ( - int(time_string[0:3]), - int(time_string[3:5]), - int(time_string[5:7]), - int(time_string[7:9]), - int(time_string[9:12]) - ) - # RT130 only stores the last two digits of the year. Because the - # documentation for RT130 does not define a way to retrieve the full year, - # we use Obspy's method. Accordingly, we convert 0-49 to 2000-2049 and - # 50-99 to 1950-1999. - if 0 <= year <= 49: - year += 2000 - elif 50 <= year <= 99: - year += 1900 - converted_time = UTCDateTime(year=year, julday=day_of_year, hour=hour, - minute=minute, second=second, - microsecond=millisecond * 1000) - return converted_time - - -def get_rt130_packet_header(rt130_packet: bytes, - unpacker: Unpacker) -> PacketHeader: - try: - # Because RT130 data is always big-endian, it is more convenient to - # use str.decode() than the unpacker. - packet_type = rt130_packet[:2].decode('ASCII') - except UnicodeError: - print(f'Cannot decode packet type.') - print('The given file does not appear to be a valid RT130 file.') - raise RT130ParseError - valid_packet_types = ['AD', 'CD', 'DS', 'DT', 'EH', 'ET', 'OM', 'SH', 'SC', - 'FD'] - if packet_type not in valid_packet_types: - print(f'Invalid packet type found: {packet_type}') - print('The given file does not appear to be a valid RT130 file.') - raise RT130ParseError - - experiment_number = int(rt130_packet[2:3].hex()) - year = int(rt130_packet[3:4].hex()) - # A call to str.upper() is needed because bytes.hex() makes any - # hexadecimal letter (i.e. ABCDEF) lowercase, while we want them to be - # uppercase for display purpose. - unit_id = rt130_packet[4:6].hex().upper() - time_bytes = rt130_packet[6:12] - packet_time = parse_rt130_time(year, time_bytes) - byte_count = int(rt130_packet[12:14].hex()) - packet_sequence = int(rt130_packet[14:16].hex()) - - return PacketHeader(packet_type, experiment_number, unit_id, packet_time, - byte_count, packet_sequence) - - -def read_rt130_file(file_name: str, unpacker: Unpacker): - # RT130 data looks to be all big-endian (logpeek assumes this, and it has - # been working pretty well), so we don't have to do any endianness check. - - # if os.path.getsize(file_name) % 1024: - # warnings.warn('The size of the data is not a multiple of 1024; it ' - # 'might be truncated.') - - packets = [] - - with open(file_name, 'rb') as rt130_file: - for i in range(os.path.getsize(file_name) // 1024): - packet = rt130_file.read(1024) - packet_header = get_rt130_packet_header(packet, unpacker) - - waveform_handlers: Dict[str, Callable] = { - 'EH': read_eh_et_packet, - 'ET': read_eh_et_packet, - 'DT': read_dt_packet, - } - - soh_handlers: Dict[str, Callable] = dict.fromkeys( - ['AD', 'CD', 'DS', 'FD', 'OM', 'SC', 'SH'], read_soh_packet - ) - - packet_handlers = { - **waveform_handlers, **soh_handlers - } - - packet_handler = packet_handlers.get( - packet_header.packet_type, packet_reader_placeholder - ) - return_val = packet_handler(packet, unpacker) - if packet_header.packet_type == 'DT': - packet_type = DTPacket - elif packet_header.packet_type in ['EH', 'ET']: - packet_type = EHETPacket - else: - packet_type = SOHPacket - - extended_header, data = return_val - current_packet = packet_type(packet_header, extended_header, data) - packets.append(current_packet) - - return packets - - -def convert_packet_to_obspy_format(packet: Union[EHETPacket, DTPacket], - unpacker: Unpacker): - # We want to convert the packet to a tuple. In order to make it easier to - # maintain, we first convert the packet to a dictionary. Then, we grab the - # values of the dictionary as tuple to get the final result. - converted_packet = {} - converted_packet['packet_type'] = packet.header.packet_type - converted_packet['experiment_number'] = packet.header.experiment_number - # Obspy only stores the last two digits of the year. - converted_packet['year'] = packet.header.time.year % 100 - converted_packet['unit_id'] = packet.header.unit_id - converted_packet['time'] = packet.header.time.ns - converted_packet['byte_count'] = packet.header.byte_count - converted_packet['packet_sequence'] = packet.header.packet_sequence - converted_packet['event_number'] = packet.extended_header.event_number - converted_packet[ - 'data_stream_number'] = packet.extended_header.data_stream_number - converted_packet['channel_number'] = packet.extended_header.channel_number - converted_packet[ - 'number_of_samples'] = packet.extended_header.number_of_samples - converted_packet['flags'] = packet.extended_header.flags - converted_packet['data_format'] = packet.extended_header.data_format - - if converted_packet['packet_type'] == 'DT': - # Obspy stores the data as list of 1-byte integers. We store the - # data as an arbitrary length integer, so we need to do some - # conversion. To make converting the resulting tuple to an element - # of a structured array of type PACKET_FINAL_DTYPE easier, we set - # the size of the payload to be 4. This only affect data with format - # 16, and as long as we are careful in self.to_stream, we don't even - # have to make a special case when decoding (note: this is possible - # because of a peculiarity of the 2's complement encoding). - data_size = 4 - format_char = 'B' - converted_packet['payload'] = numpy.empty(1000, np.uint8) - packet_data = list(unpacker.unpack( - f'{data_size}{format_char}', - packet.data.to_bytes(data_size, 'big', signed=True) - )) - converted_packet['payload'][:4] = packet_data - elif converted_packet['packet_type'] in ['EH', 'ET']: - eh_et_payload_size = eh_et_payload_end_in_packet - 24 - converted_packet['payload'] = numpy.empty(1000, np.uint8) - packet_data = numpy.frombuffer(packet.data, np.uint8) - converted_packet['payload'][:eh_et_payload_size] = packet_data - else: - converted_packet['payload'] = numpy.frombuffer(packet.data, np.uint8) - return tuple(converted_packet.values()) - - -class Reftek130(core.Reftek130): - pass - - -waveform_file = '/Users/kle/PycharmProjects/sohstationviewer/tests/test_data/RT130-sample/2017149.92EB/2017150/92EB/1/000000015_0036EE80' -soh_file = '/Users/kle/PycharmProjects/sohstationviewer/tests/test_data/RT130-sample/2017149.92EB/2017150/92EB/0/000000000_00000000' -import time -start = time.perf_counter() -(Reftek130.from_file(waveform_file).to_stream()) -# print('Time taken:', time.perf_counter() - start) diff --git a/sohstationviewer/model/reftek/rt130_experiment/reftek_helper.py b/sohstationviewer/model/reftek/rt130_experiment/reftek_helper.py index c7e7dbff5..fa19960ca 100644 --- a/sohstationviewer/model/reftek/rt130_experiment/reftek_helper.py +++ b/sohstationviewer/model/reftek/rt130_experiment/reftek_helper.py @@ -1,7 +1,26 @@ -import dataclasses -from typing import Any +import os +from typing import Any, Dict, Callable, Union -from obspy import UTCDateTime +import numpy +import numpy as np + +from sohstationviewer.model.mseed_data.record_reader_helper import Unpacker +from sohstationviewer.model.reftek.from_rt2ms.core import \ + eh_et_payload_end_in_packet +from sohstationviewer.model.reftek.rt130_experiment.dt_packet import \ + ( + read_dt_packet, DTPacket, +) +from sohstationviewer.model.reftek.rt130_experiment.eh_et_packet import \ + ( + read_eh_et_packet, EHETPacket, +) +from sohstationviewer.model.reftek.rt130_experiment.header import \ + get_rt130_packet_header +from sohstationviewer.model.reftek.rt130_experiment.soh_packets import \ + ( + read_soh_packet, SOHPacket, +) def packet_reader_placeholder(*args: Any, **kwargs: Any) -> None: @@ -16,11 +35,98 @@ class RT130ParseError(Exception): pass -@dataclasses.dataclass -class PacketHeader: - packet_type: str - experiment_number: int - unit_id: str - time: UTCDateTime - byte_count: int - packet_sequence: int +def read_rt130_file(file_name: str, unpacker: Unpacker): + # RT130 data looks to be all big-endian (logpeek assumes this, and it has + # been working pretty well), so we don't have to do any endianness check. + + # if os.path.getsize(file_name) % 1024: + # warnings.warn('The size of the data is not a multiple of 1024; it ' + # 'might be truncated.') + + packets = [] + + with open(file_name, 'rb') as rt130_file: + for i in range(os.path.getsize(file_name) // 1024): + packet = rt130_file.read(1024) + packet_header = get_rt130_packet_header(packet, unpacker) + + waveform_handlers: Dict[str, Callable] = { + 'EH': read_eh_et_packet, + 'ET': read_eh_et_packet, + 'DT': read_dt_packet, + } + + soh_handlers: Dict[str, Callable] = dict.fromkeys( + ['AD', 'CD', 'DS', 'FD', 'OM', 'SC', 'SH'], read_soh_packet + ) + + packet_handlers = { + **waveform_handlers, **soh_handlers + } + + packet_handler = packet_handlers.get( + packet_header.packet_type, packet_reader_placeholder + ) + return_val = packet_handler(packet, unpacker) + if packet_header.packet_type == 'DT': + packet_type = DTPacket + elif packet_header.packet_type in ['EH', 'ET']: + packet_type = EHETPacket + else: + packet_type = SOHPacket + + extended_header, data = return_val + current_packet = packet_type(packet_header, extended_header, data) + packets.append(current_packet) + + return packets + + +def convert_packet_to_obspy_format(packet: Union[EHETPacket, DTPacket], + unpacker: Unpacker): + # We want to convert the packet to a tuple. In order to make it easier to + # maintain, we first convert the packet to a dictionary. Then, we grab the + # values of the dictionary as tuple to get the final result. + converted_packet = {} + converted_packet['packet_type'] = packet.header.packet_type + converted_packet['experiment_number'] = packet.header.experiment_number + # Obspy only stores the last two digits of the year. + converted_packet['year'] = packet.header.time.year % 100 + converted_packet['unit_id'] = packet.header.unit_id + converted_packet['time'] = packet.header.time.ns + converted_packet['byte_count'] = packet.header.byte_count + converted_packet['packet_sequence'] = packet.header.packet_sequence + converted_packet['event_number'] = packet.extended_header.event_number + converted_packet[ + 'data_stream_number'] = packet.extended_header.data_stream_number + converted_packet['channel_number'] = packet.extended_header.channel_number + converted_packet[ + 'number_of_samples'] = packet.extended_header.number_of_samples + converted_packet['flags'] = packet.extended_header.flags + converted_packet['data_format'] = packet.extended_header.data_format + + if converted_packet['packet_type'] == 'DT': + # Obspy stores the data as list of 1-byte integers. We store the + # data as an arbitrary length integer, so we need to do some + # conversion. To make converting the resulting tuple to an element + # of a structured array of type PACKET_FINAL_DTYPE easier, we set + # the size of the payload to be 4. This only affect data with format + # 16, and as long as we are careful in self.to_stream, we don't even + # have to make a special case when decoding (note: this is possible + # because of a peculiarity of the 2's complement encoding). + data_size = 4 + format_char = 'B' + converted_packet['payload'] = numpy.empty(1000, np.uint8) + packet_data = list(unpacker.unpack( + f'{data_size}{format_char}', + packet.data.to_bytes(data_size, 'big', signed=True) + )) + converted_packet['payload'][:4] = packet_data + elif converted_packet['packet_type'] in ['EH', 'ET']: + eh_et_payload_size = eh_et_payload_end_in_packet - 24 + converted_packet['payload'] = numpy.empty(1000, np.uint8) + packet_data = numpy.frombuffer(packet.data, np.uint8) + converted_packet['payload'][:eh_et_payload_size] = packet_data + else: + converted_packet['payload'] = numpy.frombuffer(packet.data, np.uint8) + return tuple(converted_packet.values()) diff --git a/sohstationviewer/model/reftek/rt130_experiment/soh_packets.py b/sohstationviewer/model/reftek/rt130_experiment/soh_packets.py index b009e14e0..0731ba95b 100644 --- a/sohstationviewer/model/reftek/rt130_experiment/soh_packets.py +++ b/sohstationviewer/model/reftek/rt130_experiment/soh_packets.py @@ -4,8 +4,7 @@ from obspy.io.reftek.util import bcd from sohstationviewer.model.mseed_data.record_reader_helper import \ Unpacker -from sohstationviewer.model.reftek.rt130_experiment.reftek_helper import \ - PacketHeader +from sohstationviewer.model.reftek.rt130_experiment.header import PacketHeader import numpy -- GitLab