Plot Earth Orientation Parameters

This (notebook) demonstrates plotting the time series of Earth Orientation Parameters (EOPs) as a function of time. Pole tides are driven by variations in Earth’s rotation axis (e.g. Chandler wobble) with respect to the location of the mean (or secular) pole. This polar motion is forced primarily by fluctuations in the angular momentum of the ocean and atmosphere.

Python Dependencies

• numpy: Scientific Computing Tools For Python

• scipy: Scientific Tools for Python

• matplotlib: Python 2D plotting library

• timescale: Python tools for time and astronomical calculations

import datetime
import timescale
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.offsetbox as offsetbox

Create time series from beginning of EOP dataset

start_date = np.array('1973-01-01', dtype=f'datetime64[D]')
end_date = np.array(datetime.datetime.now(), dtype=f'datetime64[D]')
ts = timescale.from_range(start_date, end_date)

Calculate IERS Mean/Secular Pole

# calculate angular coordinates of mean/secular pole at time
mpx, mpy, fl = timescale.eop.iers_mean_pole(ts.year, convention='2018')
# read and interpolate IERS daily polar motion values
px, py = timescale.eop.iers_polar_motion(ts.MJD, k=1, s=0)
# calculate differentials from mean/secular pole positions
# using the latest definition from IERS Conventions (2010)
mx = px - mpx
my = -(py - mpy)

Plot EOPs and deviations from Secular Pole

f1, ax1 = plt.subplots(num=1, ncols=2, figsize=(8, 4.875))
# plot pole position
sc = ax1[0].scatter(px, py, c=ts.year, cmap='plasma_r', s=0.5)
# plot deviation from mean/secular pole
sc = ax1[1].scatter(mx, my, c=ts.year, cmap='plasma_r', s=0.5)
ax1[1].axhline(0, color='0.4', ls='--', lw=0.5)
ax1[1].axvline(0, color='0.4', ls='--', lw=0.5)
# add axis labels
ax1[0].set_xlabel('X Pole [asec]')
ax1[1].set_xlabel('X Pole [asec]')
ax1[0].set_ylabel('Y Pole [asec]')
labels = ['a)', 'b)']
for i, label in enumerate(labels):
    ax1[i].tick_params(which='both', direction='in')
    at = offsetbox.AnchoredText(label,
        loc=2, pad=0.25, borderpad=0.5, frameon=False,
        prop=dict(size=14,weight='bold',color='k'))
    ax1[i].axes.add_artist(at)
# add title
ax1[0].set_title('IERS Polar Motion')
ax1[1].set_title('IERS Deviation from Secular Pole')
# set axis limits
ax1[0].set_xlim([-0.35, 0.35])
ax1[0].set_ylim([-0.05, 0.65])
ax1[1].set_xlim([-0.35, 0.35])
ax1[1].set_ylim([-0.35, 0.35])
ax1[0].set_aspect('equal')
ax1[1].set_aspect('equal')

# Add colorbar with a colorbar axis
# Add an ax1es at position rect [left, bottom, width, height]
cbar_ax = f1.add_axes([0.07, 0.085, 0.925, 0.04])
# extend = add extension triangles to upper and lower bounds
# options: neither, both, min, max1
cbar = f1.colorbar(sc, cax=cbar_ax, extend='neither',
    drawedges=False, orientation='horizontal')
# rasterized colorbar to remove lines
cbar.solids.set_rasterized(True)
# Add label to the colorbar
cbar.ax.set_title('Time [yr]', rotation=0, y=-1.65, va='top')
cbar.ax.xaxis.set_label_coords(1.075, 0.5)
# add minor ticks
m1 = np.ceil(ts[0].to_calendar().year/2)*2
m2 = np.ceil(ts[-1].to_calendar().year/2)*2
minor_ticks = np.arange(m1, m2, 2)
cbar.set_ticks(minor_ticks, minor=True)
# set tick parameters
cbar.ax.tick_params(which='both', width=1, length=5, direction='in')

# adjust subplot within figure
f1.subplots_adjust(left=0.07,right=0.995,bottom=0.15,top=0.995,wspace=0.125)
plt.show()

Plot EOP Time Series

f2, ax2 = plt.subplots(num=2, nrows=3, sharex=True, figsize=(5, 5))
ax2[0].plot(ts.year, px, label='IERS')
ax2[0].plot(ts.year, mpx, label='Secular')
ax2[1].plot(ts.year, py, label='IERS')
ax2[1].plot(ts.year, mpy, label='Secular')
# calculate TT-UT1 and convert to seconds
tt_ut1 = 86400.0*ts.tt_ut1
ax2[2].plot(ts.year, tt_ut1)
# set axis labels
ax2[0].set_ylabel('X Pole [asec]', labelpad=0)
ax2[1].set_ylabel('Y Pole [asec]', labelpad=8)
ax2[2].set_ylabel('TT-UT1 [s]', labelpad=11)
ax2[2].set_xlabel('Time [yr]')
labels = ['a)', 'b)', 'c)']
for i, label in enumerate(labels):
    ax2[i].tick_params(which='both', direction='in')
    at = offsetbox.AnchoredText(label,
        loc=2, pad=0.0, borderpad=0.5, frameon=False,
        prop=dict(size=12,weight='bold',color='k'))
    ax2[i].axes.add_artist(at)
# add legend
lgd = ax2[0].legend(frameon=False, ncols=1, labelspacing=0.1, borderpad=0.1)
for line in lgd.get_lines():
    line.set_linewidth(6)
# set axis limits
ax2[0].set_xlim([1973, m2])
ax2[2].set_xticks(minor_ticks, minor=True)
# adjust subplot within figure
f2.subplots_adjust(left=0.105,right=0.99,bottom=0.075,top=0.99,hspace=0.1)
plt.show()