Skip to content
Snippets Groups Projects
Commit 29fe224a authored by Wernle, Alexandra Nicole's avatar Wernle, Alexandra Nicole
Browse files

Deleting BIQParserFactory (renamed to BIQfactory).

parent 906795e3
No related branches found
No related tags found
1 merge request!247make_iaga app and BIQParserFactory
"""Factory that parses BIQ files into arrays and/or ObsPy Stream."""
import numpy as np
import obspy
from typing import BinaryIO
class BIQParser:
"""Geomag BIQ file fields."""
def __init__(self, filename: BinaryIO):
self.filename = filename
# define the data types for a day of data
self.day_block_dtype = np.dtype(
[
("Station ID", "S4"),
("Year & Day", np.int32),
("Co-Latitude", np.int32),
("Longitude", np.int32),
("Elevation", np.int32),
("Orientation", "S4"),
("Origination", "S4"),
("D Conversion", np.int32),
("Quality", "S4"),
("Instrumentation", "S4"),
("Custom Field 1", np.int32),
("Custom Field 2", np.int32),
("Custom Field 3", "S4"),
("Custom Field 4", "S4"),
("Custom Field 5", np.int32),
("Custom Field 6", "S4"),
("X-minute", np.int32, (1440,)),
("Y-minute", np.int32, (1440,)),
("Z-minute", np.int32, (1440,)),
("G-minute", np.int32, (1440,)),
("X-hourly", np.int32, (24,)),
("Y-hourly", np.int32, (24,)),
("Z-hourly", np.int32, (24,)),
("G-hourly", np.int32, (24,)),
("X-daily", np.int32, (1,)),
("Y-daily", np.int32, (1,)),
("Z-daily", np.int32, (1,)),
("G-daily", np.int32, (1,)),
("K-values", np.int32, (8,)),
("Leftovers", np.int32, (4,)),
]
)
def read_data(self):
with open(self.filename, "rb") as f:
data = np.fromfile(f, dtype=self.day_block_dtype)
# call the check_station_id method to check that the Station ID for each day is the same
self.check_station_id(data)
return data
def check_station_id(self, data):
# get the unique station IDs
station_ids = set(day_data["Station ID"].decode() for day_data in data)
# raise an error if there is more than one unique station ID
if len(station_ids) > 1:
raise ValueError(
"One or more station IDs do not match. File may be corrupted or written incorrectly!"
)
def get_stream(self):
data = self.read_data()
# combine the minute data from all the days into an obspy stream
st = obspy.Stream()
tr = obspy.Trace(
header={
"network": "NT",
"location": "Q0", # THIS IS AN ASSUMPTION!
"station": data[0]["Station ID"].decode().strip(),
"delta": 60,
"starttime": obspy.UTCDateTime(str(data[0]["Year & Day"])),
}
)
for day in data:
# add stats for header information
tr.stats.starttime = obspy.UTCDateTime(str(day["Year & Day"]))
x_tr = tr.copy()
y_tr = tr.copy()
z_tr = tr.copy()
g_tr = tr.copy()
x_tr.stats.channel = "X"
y_tr.stats.channel = "Y"
z_tr.stats.channel = "Z"
g_tr.stats.channel = "G"
# convert to floats so that we're not stuck dealing with ints later
x_tr.data = day["X-minute"].astype(np.float64)
y_tr.data = day["Y-minute"].astype(np.float64)
z_tr.data = day["Z-minute"].astype(np.float64)
g_tr.data = day["G-minute"].astype(np.float64)
st += obspy.Stream(traces=[x_tr, y_tr, z_tr, g_tr])
st.merge()
# data was stored as scaled integers, unscale them and translate the 99999's to nan
for tr in st:
tr.data[tr.data == 999999] = np.nan
tr.data = tr.data / 10
st = self.append_f(st)
return st
def append_f(self, st):
# calculate F from XYZ and G ("delta-F")
f_xyz = np.zeros_like(st.select(channel="G")[0].data)
for tr in st.select(channel="[X|Y|Z]"):
f_xyz = f_xyz + np.square(tr.data)
f_xyz = np.sqrt(f_xyz)
f_tr = obspy.Trace(
f_xyz - st.select(channel="G")[0].data,
header=st.select(channel="G")[0].stats,
)
f_tr.stats.channel = "F"
st.append(f_tr)
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