diff --git a/Calculation_test.py b/Calculation_test.py index 013e2949dbf03c88770a6c6a217546a49072d4bd..b62c7a3036813356be6dafe77e1f90ddc8add270 100644 --- a/Calculation_test.py +++ b/Calculation_test.py @@ -43,30 +43,3 @@ for filename in os.listdir("etc/residual"): err_msg="Baselines not within 4 decimals", verbose=True, ) - # gather original and resulting diagnostics - o_diagnostics = original.diagnostics - r_diagnostics = original.diagnostics - # test mean mark values - assert_almost_equal( - o_diagnostics.mean_mark, - r_diagnostics.mean_mark, - decimal=2, - err_msg="Baselines not within 4 decimals", - verbose=True, - ) - # test magnetic azimuth values - assert_almost_equal( - o_diagnostics.magnetic_azimuth, - r_diagnostics.magnetic_azimuth, - decimal=2, - err_msg="Baselines not within 4 decimals", - verbose=True, - ) - # test meridian values - assert_almost_equal( - o_diagnostics.meridian, - r_diagnostics.meridian, - decimal=2, - err_msg="Baselines not within 4 decimals", - verbose=True, - ) diff --git a/geomagio/residual/Calculation.py b/geomagio/residual/Calculation.py index 1ec21a885c207603d469bf536927c7885b931d9b..6d58c05401ee6d50655b478fb980321c4496b1ec 100644 --- a/geomagio/residual/Calculation.py +++ b/geomagio/residual/Calculation.py @@ -30,7 +30,7 @@ def calculate(reading: Reading, adjust_reference: bool = True) -> Reading: # reference measurement, used to adjust absolutes reference = reading[mt.WEST_DOWN][0] # calculate inclination - inclination, f, mean = calculate_I( + inclination, f, i_mean = calculate_I( hemisphere=reading.hemisphere, measurements=reading.measurements ) corrected_f = f + reading.pier_correction @@ -38,15 +38,21 @@ def calculate(reading: Reading, adjust_reference: bool = True) -> Reading: absoluteH, absoluteZ = calculate_HZ_absolutes( corrected_f=corrected_f, inclination=inclination, - mean=mean, + mean=i_mean, reference=adjust_reference and reference or None, ) - absoluteD, diagnostics = calculate_D_absolute( + absoluteD, d_mean, m_mean = calculate_D_absolute( azimuth=reading.azimuth, h_baseline=absoluteH.baseline, measurements=reading.measurements, reference=adjust_reference and reference or None, ) + # populate diagnostics object with averaged measurements + diagnostics = Diagnostics( + inclination_measurement=i_mean, + declination_measurement=d_mean, + mark_measurement=m_mean, + ) # calculate scale scale_value = None scale_measurements = reading[mt.NORTH_DOWN_SCALE] @@ -89,14 +95,14 @@ def calculate_D_absolute( mean = average_measurement(measurements, DECLINATION_TYPES) reference = reference or mean # average mark - average_mark = average_measurement(measurements, MARK_TYPES).angle + average_mark = average_measurement(measurements, MARK_TYPES) # adjust based on which is larger mark_up = average_measurement(measurements, [mt.FIRST_MARK_UP]).angle mark_down = average_measurement(measurements, [mt.FIRST_MARK_DOWN]).angle if mark_up < mark_down: - average_mark += 90 + average_mark.angle += 90 else: - average_mark -= 90 + average_mark.angle -= 90 # declination measurements declination_measurements = [ average_measurement(measurements, [t]) for t in DECLINATION_TYPES @@ -118,15 +124,9 @@ def calculate_D_absolute( shift = -180 # add subtract average mark angle from average meridian angle and add # azimuth to get the declination baseline - d_b = (meridian - average_mark) + azimuth + shift + d_b = (meridian - average_mark.angle) + azimuth + shift # calculate absolute d_abs = d_b + np.degrees(np.arctan(reference.e / (reference.h + h_baseline))) - # for diagnostics - magnetic_azimuth = average_mark - mean.angle - if average_mark > 180: - average_mark -= 90 - if meridian > 180: - meridian -= 90 return ( Absolute( @@ -137,11 +137,8 @@ def calculate_D_absolute( starttime=mean.time, endtime=mean.endtime, ), - Diagnostics( - meridian=mean.angle, - mean_mark=average_mark, - magnetic_azimuth=magnetic_azimuth, - ), + mean, + average_mark, ) diff --git a/geomagio/residual/Diagnostics.py b/geomagio/residual/Diagnostics.py index 32d773c1f76931a43eb2a57be8b06b0f2316b71d..7cd1bde862bcf895d8c3243fd0f48a782905583e 100644 --- a/geomagio/residual/Diagnostics.py +++ b/geomagio/residual/Diagnostics.py @@ -1,5 +1,6 @@ from typing import Optional +from .Measurement import Measurement from pydantic import BaseModel @@ -8,11 +9,11 @@ class Diagnostics(BaseModel): Attributes ---------- - meridian: calculated from declination measurements - mean_mark: average mark angles from measurements - magnetic_azimuh: after adjustment + inclination_measurement: Average of inclination measurements + declination_measurement: Average of declination measurements + mark_measurement: Average of mark measurements """ - meridian: float = None - mean_mark: float = None - magnetic_azimuth: float = None + inclination_measurement: Measurement + declination_measurement: Measurement + mark_measurement: Measurement diff --git a/geomagio/residual/SpreadsheetAbsolutesFactory.py b/geomagio/residual/SpreadsheetAbsolutesFactory.py index 63abcb20b9f12baa071d73bfb4d31513259d7702..84d3aefafdd5c74803815e2e4c3b6d4bcf000be9 100644 --- a/geomagio/residual/SpreadsheetAbsolutesFactory.py +++ b/geomagio/residual/SpreadsheetAbsolutesFactory.py @@ -6,7 +6,12 @@ from obspy.core import UTCDateTime import openpyxl from .Absolute import Absolute -from .Calculation import DECLINATION_TYPES, MARK_TYPES, average_measurement +from .Calculation import ( + DECLINATION_TYPES, + MARK_TYPES, + INCLINATION_TYPES, + average_measurement, +) from .Diagnostics import Diagnostics from .Measurement import Measurement from .MeasurementType import MeasurementType as mt @@ -409,18 +414,14 @@ class SpreadsheetAbsolutesFactory(object): """ Gather diagnostics from list of measurements """ - mean_mark = average_measurement(measurements, MARK_TYPES).angle - - meridian = average_measurement(measurements, DECLINATION_TYPES).angle - - magnetic_azimuth = mean_mark - meridian - if meridian > 180: - magnetic_azimuth = mean_mark - (meridian - 90) - if mean_mark > 180: - mean_mark -= 90 - return Diagnostics( - meridian=meridian, mean_mark=mean_mark, magnetic_azimuth=magnetic_azimuth, + inclination_measurement=average_measurement( + measurements, INCLINATION_TYPES + ), + declination_measurement=average_measurement( + measurements, DECLINATION_TYPES + ), + mark_measurement=average_measurement(measurements, MARK_TYPES), )