From 5244a0297bbeaad819179e669ee445559fee2a0f Mon Sep 17 00:00:00 2001 From: Tyler Sutterley Date: Wed, 23 Jul 2025 11:44:57 -0700 Subject: [PATCH 1/4] fix: use numpy interp for 2015 convention mean pole values --- timescale/eop.py | 38 +++++++++++--------------------------- 1 file changed, 11 insertions(+), 27 deletions(-) diff --git a/timescale/eop.py b/timescale/eop.py index 498a047..2e3ebf6 100644 --- a/timescale/eop.py +++ b/timescale/eop.py @@ -1,7 +1,7 @@ #!/usr/bin/env python u""" eop.py -Written by Tyler Sutterley (04/2024) +Written by Tyler Sutterley (07/2025) Utilities for maintaining and calculating Earth Orientation Parameters (EOP) PYTHON DEPENDENCIES: @@ -14,6 +14,7 @@ utilities.py: download and management utilities for syncing files UPDATE HISTORY: + Updated 07/2025: use numpy interp for 2015 convention mean pole values Updated 04/2024: fixed annotations with multiple types Updated 04/2023: using pathlib to define and expand paths add wrapper function for interpolating daily EOP values @@ -354,14 +355,9 @@ def iers_mean_pole( # reduce to 1971 to end date ii, = np.nonzero(table[:, 0] >= 1971) table = np.copy(table[ii,:]) - # reduce to yearly values + # reduce to yearly values following 2015 conventions jj, = np.nonzero((table[:, 0] % 1) == 0.0) table = np.copy(table[jj,:]) - end_time = table[-1, 0] + 0.2 - # final shape of the table - nrows, *_ = np.shape(table) - else: - end_time = np.inf # allocate for output arrays x = np.full_like(input_epoch, kwargs['fill_value']) y = np.full_like(input_epoch, kwargs['fill_value']) @@ -383,26 +379,14 @@ def iers_mean_pole( y[t] = 0.358891 - 0.0006287*dx flag[t] = True # Conventional mean pole model in IERS Conventions 2015 - # must be below maximum valid date within file (e.g. 2015.2 for 2015) - elif (convention == '2015') and (epoch >= 1975) and (epoch < end_time): - # find epoch within mean pole table - i = 1 - j = nrows+1 - while (j > (i+1)): - k = (i+j)//2 - if (epoch < table[k, 0]): - j = k - else: - i = k - # calculate differential from point in table - dx = epoch - table[i, 0] - if (i == (nrows-1)): - x[t] = table[i, 1] + dx*(table[nrows-1, 1]-table[nrows-2, 1]) - y[t] = table[i, 2] + dx*(table[nrows-1, 1]-table[nrows-2, 2]) - else: - x[t] = table[i, 1] + dx*(table[i+1, 1]-table[i, 1]) - y[t] = table[i, 2] + dx*(table[i+1, 2]-table[i, 2]) - flag[t] = True + # epoch must be within the dates in the mean pole file + elif (convention == '2015') and (epoch >= 1975): + # interpolate using times in table + x[t] = np.interp(epoch, table[:,0], table[:,1], + left=kwargs['fill_value'], right=kwargs['fill_value']) + y[t] = np.interp(epoch, table[:,0], table[:,2], + left=kwargs['fill_value'], right=kwargs['fill_value']) + flag[t] = (x[t] != kwargs['fill_value']) # Secular pole model in IERS Conventions 2018 elif (convention == '2018'): # calculate secular pole From ad15f6edd749d73fca606d7265a8ee4d5989f828 Mon Sep 17 00:00:00 2001 From: Tyler Sutterley Date: Wed, 23 Jul 2025 12:00:06 -0700 Subject: [PATCH 2/4] Update eop.py --- timescale/eop.py | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/timescale/eop.py b/timescale/eop.py index 2e3ebf6..e7e8e62 100644 --- a/timescale/eop.py +++ b/timescale/eop.py @@ -352,12 +352,11 @@ def iers_mean_pole( # read mean pole file input_file = pathlib.Path(kwargs['file']).expanduser().absolute() table = np.loadtxt(input_file) - # reduce to 1971 to end date - ii, = np.nonzero(table[:, 0] >= 1971) + # Reduce table following 2015 conventions: + # 1. trim dates prior to 1971 + # 2. only keep rows falling on exact years + ii, = np.nonzero((table[:, 0] >= 1971) & ((table[:, 0] % 1) == 0.0)) table = np.copy(table[ii,:]) - # reduce to yearly values following 2015 conventions - jj, = np.nonzero((table[:, 0] % 1) == 0.0) - table = np.copy(table[jj,:]) # allocate for output arrays x = np.full_like(input_epoch, kwargs['fill_value']) y = np.full_like(input_epoch, kwargs['fill_value']) From dcd990f4ff02c1dbc8d0a27a115a638fed41953c Mon Sep 17 00:00:00 2001 From: Tyler Sutterley Date: Wed, 23 Jul 2025 16:09:06 -0700 Subject: [PATCH 3/4] Update eop.py --- timescale/eop.py | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/timescale/eop.py b/timescale/eop.py index e7e8e62..145b04a 100644 --- a/timescale/eop.py +++ b/timescale/eop.py @@ -15,6 +15,7 @@ UPDATE HISTORY: Updated 07/2025: use numpy interp for 2015 convention mean pole values + add Desai et al. (2015) secular pole model Updated 04/2024: fixed annotations with multiple types Updated 04/2023: using pathlib to define and expand paths add wrapper function for interpolating daily EOP values @@ -314,7 +315,7 @@ def iers_mean_pole( ): """ Calculates the angular coordinates of the IERS Conventional Mean Pole (CMP) - or IERS Secular Pole (2018 convention) :cite:p:`Petit:2010tp` + or IERS Secular Pole (2018 convention) :cite:p:`Petit:2010tp,Desai:2015jr` Parameters ---------- @@ -327,6 +328,7 @@ def iers_mean_pole( - ``'2003'`` - ``'2010'`` - ``'2015'`` + - ``'Desai'`` - ``'2018'`` input_file: str or pathlib.Path Full path to mean-pole.tab file provided by IERS @@ -386,6 +388,12 @@ def iers_mean_pole( y[t] = np.interp(epoch, table[:,0], table[:,2], left=kwargs['fill_value'], right=kwargs['fill_value']) flag[t] = (x[t] != kwargs['fill_value']) + # Secular pole model in Desai et al. (2015) + elif (convention == 'Desai'): + # calculate secular pole using equation 10a/b of Desai (2015) + x[t] = 0.05097 + 0.00062*(epoch - 2000.0) + y[t] = 0.33449 + 0.00348*(epoch - 2000.0) + flag[t] = True # Secular pole model in IERS Conventions 2018 elif (convention == '2018'): # calculate secular pole From 160ba6add948440d72898cc06186ff4bf2147a23 Mon Sep 17 00:00:00 2001 From: Tyler Sutterley Date: Wed, 23 Jul 2025 17:19:13 -0700 Subject: [PATCH 4/4] Update eop.py --- timescale/eop.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/timescale/eop.py b/timescale/eop.py index 145b04a..16ee8b9 100644 --- a/timescale/eop.py +++ b/timescale/eop.py @@ -306,7 +306,7 @@ def update_finals_file( return # IERS mean or secular pole conventions -_conventions = ('2003', '2010', '2015', '2018') +_conventions = ('2003', '2010', '2015', 'Desai', '2018') # read table of mean pole values, calculate angular coordinates at epoch def iers_mean_pole( input_epoch: np.ndarray,