diff --git a/geomagio/adjusted/AdjustedMatrix.py b/geomagio/adjusted/AdjustedMatrix.py index d4c77ceccf1083442615ba3427102fe6c0316479..b5513877fdd621d4cd33755d78c59d7da5b7d599 100644 --- a/geomagio/adjusted/AdjustedMatrix.py +++ b/geomagio/adjusted/AdjustedMatrix.py @@ -69,30 +69,13 @@ class AdjustedMatrix(BaseModel): predicted = self.matrix @ stacked_ordinates metrics = [] elements = ["X", "Y", "Z", "dF"] - expected = np.vstack( - ( - absolutes, - ChannelConverter.get_computed_f_using_squares(*absolutes), - ) - ) - predicted = np.vstack( - ( - predicted[0:3], - ChannelConverter.get_computed_f_using_squares(*predicted[0:3]), - ) - ) - for i in range(len(elements)): - diff = expected[i] - predicted[i] - metrics.append( - Metric( - element=elements[i], - absmean=abs(np.nanmean(diff)), - stddev=np.std(diff), - ) - ) - metrics.append( - get_metric( - element=elements[i], expected=expected[i], actual=predicted[i] - ) - ) - return metrics + expected = list(absolutes) + [ + ChannelConverter.get_computed_f_using_squares(*absolutes) + ] + predicted = list(predicted[0:3]) + [ + ChannelConverter.get_computed_f_using_squares(*predicted[0:3]) + ] + return [ + get_metric(element=elements[i], expected=expected[i], actual=predicted[i]) + for i in range(len(elements)) + ] diff --git a/geomagio/adjusted/Affine.py b/geomagio/adjusted/Affine.py index 83e5aaa316e7d533a094afae5f27aec97bbfe49a..0a7a979c3a8b9629113d7df6b88e2e3261be90d2 100644 --- a/geomagio/adjusted/Affine.py +++ b/geomagio/adjusted/Affine.py @@ -28,7 +28,7 @@ class Affine(BaseModel): observatory: str = None starttime: UTCDateTime = Field(default_factory=lambda: UTCDateTime() - (86400 * 7)) - endtime: UTCDateTime = UTCDateTime() + endtime: UTCDateTime = Field(default_factory=lambda: UTCDateTime()) update_interval: Optional[int] = 86400 * 7 transforms: List[Transform] = [ RotationTranslationXY(memory=(86400 * 100), acausal=True), @@ -69,10 +69,8 @@ class Affine(BaseModel): or (epoch_end is None or r.time < epoch_end) ] M = self.calculate_matrix(time, readings) - # if readings are trimmed by bad data, mark the valid interval - if M: - M.starttime = epoch_start - M.endtime = epoch_end + M.starttime = epoch_start + M.endtime = epoch_end time += update_interval Ms.append(M) @@ -104,9 +102,9 @@ class Affine(BaseModel): readings=readings, time=time.timestamp, ) - # return None if no valid observations + # raise ValueError if no valid observations if np.sum(weights) == 0: - return AdjustedMatrix(time=time) + raise ValueError(f"No valid observations for: {time}") M = transform.calculate( ordinates=inputs, absolutes=absolutes, weights=weights diff --git a/geomagio/algorithm/AdjustedAlgorithm.py b/geomagio/algorithm/AdjustedAlgorithm.py index 3dbd5f766d55bf879417f4a7265044db56f1d505..9782e9713097912371865078503536ac00811204 100644 --- a/geomagio/algorithm/AdjustedAlgorithm.py +++ b/geomagio/algorithm/AdjustedAlgorithm.py @@ -53,22 +53,19 @@ class AdjustedAlgorithm(Algorithm): raise FileNotFoundError("statefile not found") # Adjusted matrix defaults to identity matrix matrix_size = len([c for c in self.get_input_channels() if c != "F"]) + 1 - if "pier_correction" in data.keys(): - m = AdjustedMatrix(**data) - matrix = m.matrix - pier_correction = m.pier_correction - elif "PC" in data.keys(): + if "pier_correction" in data: + self.matrix = AdjustedMatrix(**data) + elif "PC" in data: matrix = np.eye(matrix_size) # read data from legacy format for row in range(matrix_size): for col in range(matrix_size): matrix[row, col] = np.float64(data[f"M{row+1}{col+1}"]) pier_correction = np.float64(data["PC"]) + self.matrix = AdjustedMatrix(matrix=matrix, pier_correction=pier_correction) else: raise ValueError("pier correction not found in statefile") - self.matrix = AdjustedMatrix(matrix=matrix, pier_correction=pier_correction) - def save_state(self): """Save algorithm state to a file. File name is self.statefile. diff --git a/test/adjusted_test/adjusted_test.py b/test/adjusted_test/adjusted_test.py index 3bab69ac0da268ff1dff81ffb13eab1dcd041c3c..2ac3747660af93332bdfc4cd17c6fe29e106d666 100644 --- a/test/adjusted_test/adjusted_test.py +++ b/test/adjusted_test/adjusted_test.py @@ -3,6 +3,7 @@ import json import numpy as np from numpy.testing import assert_equal, assert_array_almost_equal from obspy.core import UTCDateTime +import pytest from geomagio.adjusted import AdjustedMatrix from geomagio.adjusted.Affine import Affine, get_epochs @@ -118,19 +119,23 @@ def test_BOU201911202001_infinite_weekly(): def test_BOU201911202001_invalid_readings(): - readings = [] starttime = UTCDateTime("2019-11-01T00:00:00Z") - result = Affine( - observatory="BOU", - starttime=starttime, - endtime=UTCDateTime("2020-01-31T23:59:00Z"), - transforms=[ - RotationTranslationXY(memory=np.inf, acausal=True), - TranslateOrigins(memory=np.inf, acausal=True), - ], - update_interval=None, - ).calculate(readings=readings,)[0] - assert result == AdjustedMatrix(time=starttime) + with pytest.raises( + ValueError, match=f"No valid observations for: {starttime}" + ) as error: + readings = [] + result = Affine( + observatory="BOU", + starttime=starttime, + endtime=UTCDateTime("2020-01-31T23:59:00Z"), + transforms=[ + RotationTranslationXY(memory=np.inf, acausal=True), + TranslateOrigins(memory=np.inf, acausal=True), + ], + update_interval=None, + ).calculate( + readings=readings, + ) def test_BOU201911202001_short_acausal():