CCE BIN Clock Calibrations¶
Run programs to compare real time clock with BIN system time and with instrument clocks - Mike McCann 12 May 2016
The Python programs used to get data into formats to be read into this Notebook are in MBARI's CVS on moonjelly: http://moonjelly.shore.mbari.org/cgi-bin/cvsweb.cgi/DPforSSDS/py/ (they have also been copied to the CCE_Processed/Software/ssds directory). They were executed thusly:
cd /u/ssdsadmin/dev/DPforSSDS/py
./bin_clocks_parse.py > /mbari/CCE_Processed/BIN/20151013/clocks/bin_sys_rtc_clocks.csv
./ssds_adcp_extract.py --timestamps --stride 10000 --device_id 1825 \
> /mbari/CCE_Processed/BIN/20151013/clocks/adcp_1825_clocks.csv
./ssds_adcp_extract.py --timestamps --stride 10000 --device_id 1827 \
> /mbari/CCE_Processed/BIN/20151013/clocks/adcp_1827_clocks.csv
./ssds_adcp_extract.py --timestamps --stride 10000 --device_id 1828 \
> /mbari/CCE_Processed/BIN/20151013/clocks/adcp_1828_clocks.csvTo execute this Notebook mount the cifs://atlas.shore.mbari.org/cce_processed share on your compter and cd to the BIN/20151013/clocks directory (e.g. on a Mac that would be cd /Volumes/CCE_Processed/BIN/20151013/clocks) and then execute the jupyter notebook command in a terminal window. This assumes that you have Anaconda2 (or equivalent) installed on your computer.
Examine the real time clock and BIN contoller system clock¶
Read the csv file into a Pandas data frame and show the first 10 records
import pandas as pd
df = pd.read_csv('bin_sys_rtc_clocks.csv', parse_dates=True, index_col='realtimeclock_datetime')
df.head()
Plot over time the seconds difference between the real time clock and the system clock
%matplotlib inline
from bokeh.plotting import figure
from bokeh.io import output_notebook, show
from bokeh.resources import INLINE
output_notebook(resources=INLINE, hide_banner=True)
p = figure(title="BIN Controller Clock Time Difference with Real Time Clock",
width=900, height=300,
x_axis_type='datetime',
x_axis_label='RTC Time (GMT)',
y_axis_label='seconds_difference')
p.line(df.index, df.seconds_difference, line_width=1)
_ = show(p)
Examine the ADCP clocks compared to the BIN controller clock¶
Read the csv files created by the ssds_adcp_extract.py script into Pandas data frames and make a common index name
df_1825 = pd.read_csv('adcp_1825_clocks.csv', parse_dates=True, index_col='1825_siam_datetime')
df_1827 = pd.read_csv('adcp_1827_clocks.csv', parse_dates=True, index_col='1827_siam_datetime')
df_1828 = pd.read_csv('adcp_1828_clocks.csv', parse_dates=True, index_col='1828_siam_datetime')
df_1825.index.names = ['siam_datetime']
df_1827.index.names = ['siam_datetime']
df_1828.index.names = ['siam_datetime']
The ~7 hour time difference indicates that the ADCPs were set to local time. Join the data together and plot the clock drifts on the same figure.
import matplotlib.pyplot as plt
plt.rcParams['figure.figsize'] = (16, 4)
df_adcp = df_1825.join(df_1827, how='outer').join(df_1828, how='outer')
dfp = df_adcp[['1825_seconds_difference', '1827_seconds_difference', '1828_seconds_difference']]
dfp.plot(marker='.', grid=True)
The SIAM system clock drifted about 1200 seconds relative to all the ADCP instrument clocks.
Calculate SIAM system to real time clock regression¶
Define function that returns slope and intercept of linear fit through data
import statsmodels.api as sm
def linear_fit(X, y):
'Return slope, intercept, and y_hat predicted values'
model = sm.OLS(y, sm.add_constant(X))
results = model.fit()
intercept, slope = results.params
print(('slope = {}, intercept = {}').format(slope, intercept))
y_hat = results.predict(sm.add_constant(X))
return slope, intercept, y_hat
Examination of the rows of df reveals that 2016-02-18 19:19:47 is the last time before the BIN started rebooting, pull those data out for the regression. Make a copy of that subset and compute a linear least-squares fit of system (SIAM) clock to real time clock for the time period we have ADCP data in SSDS.
df_subset = df.loc[:'2016-02-18 19:19:47'].copy()
X = df_subset.system_epoch_seconds
y = df_subset.realtimeclock_epoch_seconds
siam2rtc_slope, siam2rtc_intercept, siam2rtc_y_hat = linear_fit(X, y)
The equation:
corrected_time = 0.999882964235 * siam_clock + 169079.899744
Applies to siam times where there is a slope in the BIN Controller Clock Time Difference with Real Time Clock plot. Let's lis the records at the boundaries of the period of system rebooting to get the siam clock times between which no correction is needed.
df.loc['2016-02-18 17:00':'2016-02-18 20:40']
df.loc['2016-02-29 18:45':'2016-02-29 23:40']
The equation:
corrected_time = 0.999882964235 * siam_clock + 169079.899744
Applies to siam times approximately before 2016-02-18 20:19:50 (system_epoch_seconds = 1455826790) and after 2016-02-29 19:57:47 (system_epoch_seconds = 1456779467). In between these times the system continually rebooted with less than 2 second difference from the real time clock.
Find the exact time of the system rebooting period¶
Extract a .csv file from SSDS for a few day period around the start of the time when the system started rebooting:
ssds_ctd_extract.py --device_id 1533 --csv --start 20160218T000000 --end 20160220T000000 > ctd.csv
Load that file into a Pandas dataframe and plot the diff of datetime to find the exact time of the first system reboot
df_ctd = pd.read_csv('ctd.csv', parse_dates=True, index_col='datetime_gmt')
p = figure(title="CTD Time Around Start of System Reboots",
width=900, height=300,
x_axis_type='datetime',
x_axis_label='SIAM Time',
y_axis_label='datetime_es.diff (sec)')
p.line(df_ctd.index, df_ctd.datetime_es.diff(), line_width=1)
_ = show(p)
Examine the rows of data at the time of first reboot
pd.set_option('display.float_format', lambda x: '%.1f' % x)
df_ctd['datetime_es_diff'] = df_ctd.datetime_es.diff()
df_ctd.loc['2016-02-18 20:41:45':'2016-02-18 20:43:00']
Any SIAM time greater than 2016-02-18 20:42:11.581200 (datetime_es = 1455828131) and less than 2016-02-29 19:57:47 (system_epoch_seconds = 1456779467) is essentially the same as the real time clock time.
Examine the error¶
Add the predicted time and the difference from observed to the data frame
df_subset['y_hat'] = siam2rtc_y_hat
df_subset['error'] = df_subset.realtimeclock_epoch_seconds - df_subset.y_hat
Plot the error of the predicted value to confirm the functioning of the algorithm and the quality of the fit
p = figure(title="Predicted Real Time Clock Error",
width=900, height=300,
x_axis_type='datetime',
x_axis_label='RTC Time (GMT)',
y_axis_label='error (seconds)')
p.line(df_subset.index, df_subset.error, line_width=1)
_ = show(p)
# Spot check
df_subset.system_epoch_seconds[900] * siam2rtc_slope + siam2rtc_intercept, df_subset.realtimeclock_epoch_seconds[900]
Calculate ADCP instrument to real time clock lookup, via the SIAM system clock¶
Compute slope and intercepts of the regression of ADCP clocks onto the SIAM system clock
adcp1825siam_slope, adcp1825siam_intercept, _ = linear_fit(
df_adcp['1825_instrument_epoch_seconds'].dropna(),
df_adcp['1825_siam_epoch_seconds'].dropna())
adcp1827siam_slope, adcp1827siam_intercept, _ = linear_fit(
df_adcp['1827_instrument_epoch_seconds'].dropna(),
df_adcp['1827_siam_epoch_seconds'].dropna())
adcp1828siam_slope, adcp1828siam_intercept, _ = linear_fit(
df_adcp['1828_instrument_epoch_seconds'].dropna(),
df_adcp['1828_siam_epoch_seconds'].dropna())
Equation to convert ADCP clock to real time clock time is:
rtc = siam2rtc_slope * (adcp2siam_slope * adcp_clock + adcp2siam_intercept) + siam2rtc_intercept
Compute slope and intercept for each ADCP clock conversion to real time clock
adcp1825rtc_slope = siam2rtc_slope * adcp1825siam_slope
adcp1825rtc_intercept = siam2rtc_slope * adcp1825siam_intercept + siam2rtc_intercept
print(('ADCP 1825 - 300 KHz: slope = {:.10f}, intercept = {:.6f}').format(
adcp1825rtc_slope, adcp1825rtc_intercept))
adcp1827rtc_slope = siam2rtc_slope * adcp1827siam_slope
adcp1827rtc_intercept = siam2rtc_slope * adcp1827siam_intercept + siam2rtc_intercept
print(('ADCP 1827 - 1200 KHz: slope = {:.10f}, intercept = {:.6f}').format(
adcp1827rtc_slope, adcp1827rtc_intercept))
adcp1828rtc_slope = siam2rtc_slope * adcp1828siam_slope
adcp1828rtc_intercept = siam2rtc_slope * adcp1828siam_intercept + siam2rtc_intercept
print(('ADCP 1828 - 600 KHz: slope = {:.10f}, intercept = {:.6f}').format(
adcp1828rtc_slope, adcp1828rtc_intercept))