Skip to content
Snippets Groups Projects
Commit 86a7d20c authored by Wilbur, Spencer Franklin's avatar Wilbur, Spencer Franklin
Browse files

Removed previous lines of code that checked metadata for Azimuth and Dip and...

Removed previous lines of code that checked metadata for Azimuth and Dip and replaced it with a new function that rotates the Stream object if channels X, Y, or Z, are requested and returns the appropriate channel after rotation.
parent 134defe3
No related branches found
No related tags found
1 merge request!335FDSN 3D Rotation
......@@ -244,10 +244,7 @@ class FDSNFactory(TimeseriesFactory):
)
except FDSNNoDataException:
data = Stream()
print("The data count is", data.count)
# make sure data is 32bit int
for trace in data:
trace.data = trace.data.astype("i4")
data.merge()
if data.count() == 0 and add_empty_channels:
data += self._get_empty_trace(
......@@ -264,29 +261,12 @@ class FDSNFactory(TimeseriesFactory):
TimeseriesUtility.pad_and_trim_trace(
trace=data[0], starttime=starttime, endtime=endtime
)
# Beneath is necessary code to check the reported azimuth
# for the LF1 (E or Y) and LF2 (H or X) channels and determine if intstrument axes have been reversed.
azi, dip = FDSNFactory._get_orientations(self, data[0], starttime, sncl)
# Checks Azimuth for LF2 channel
if 270 > azi and azi > 90.0 and data[0].stats.channel[-2:] == "F2":
data[0].data *= -1
# Checks Azimuth for LF1 channel
elif azi > 180 and azi < 360 and data[0].stats.channel[-2:] == "F1":
data[0].data *= -1
# Checks Dip for LFZ channel
elif dip > 0 and dip < 180 and data[0].stats.channel[-2:] == "FZ":
data[0].data *= -1
# Remove channel Response:
# The function "remove_response" appears to be doing what we want;
# i.e. convert from counts to NT, but this may be a placeholder
# at least until we see how this function behaves if
# a station has a frequency response.
if data.count != 0:
data.remove_response(output="DEF", zero_mean=False, taper=False)
# Rotate the trace into a right handed coordinate frame.
# This will worrk assuming the metadata is reported correctly.
data = self._rotate_and_return_requested_channel(
sncl, starttime, endtime, channel
)
self._set_metadata(data, observatory, channel, type, interval)
......@@ -358,30 +338,63 @@ class FDSNFactory(TimeseriesFactory):
trace.stats, observatory, channel, type, interval
)
def _get_orientations(
self, trace: Trace, starttime: UTCDateTime, sncl
) -> Tuple[float, float]:
# Retrieve station orientation information using FDSN for each trace
FDSN = FDSNClient("IRIS")
inv = FDSN.get_stations(
network=sncl.network,
station=sncl.station,
channel=sncl.channel,
level="channel",
)
# Construct the channel code using the current trace's information
channel_code = f"{sncl.network}.{sncl.station}.{sncl.location}.{sncl.channel}"
def _rotate_and_return_requested_channel(
self, sncl, starttime: UTCDateTime, endtime: UTCDateTime, requested_channel: str
) -> Stream:
"""
Retrieves the necessary LF channels, rotates them to ZNE, and returns the requested channel.
Channels are mapped as follows: 'X' -> 'LF2', 'Y' -> 'LF1', 'Z' -> 'LFZ'.
"""
# Mapping X, Y, Z to LF channels
channel_map = {"X": "LF2", "Y": "LF1", "Z": "LFZ"}
# Initialize FDSN client and get the necessary data
FDSN = self.Client
# Determine if the requested channel is X, Y, or Z
if requested_channel in channel_map:
# Pull the LF* channels needed for rotation
lf_channels = ["LF1", "LF2", "LFZ"]
inv = FDSN.get_stations(
network=sncl.network,
station=sncl.station,
location=sncl.location, # Use location if necessary
channel=",".join(lf_channels), # Request LF1, LF2, and LFZ
starttime=starttime,
endtime=endtime,
level="response",
)
# Get orientation for the constructed channel code and time
orient = inv.get_orientation(channel_code, starttime)
# # Get the waveform data for LF* channels
st = FDSN.get_waveforms(
network=sncl.network,
station=sncl.station,
location=sncl.location,
channel="LF*",
starttime=starttime,
endtime=endtime,
)
azimuth = orient["azimuth"]
dip = orient["dip"]
# Rotate the stream to ZNE
st.rotate(method="->ZNE", inventory=inv)
print(st)
# Set azimuth in trace.stats.azimuth
trace.stats.azimuth = azimuth
# Now return only the requested channel (mapped to Z, N, or E)
if requested_channel == "X":
return st.select(channel="LFN") # N after rotation
elif requested_channel == "Y":
return st.select(channel="LFE") # E after rotation
elif requested_channel == "Z":
return st.select(channel="LFZ") # Z remains Z
print(st)
# Return both azimuth and dip
return azimuth, dip
else:
# If the requested channel is not X, Y, or Z, just retrieve and return that specific channel
st = FDSN.get_waveforms(
network=sncl.network,
station=sncl.station,
location=sncl.location,
channel=sncl.channel,
starttime=starttime,
endtime=endtime,
)
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