"""
Plotting of NTS WR-link data
AW2026-06-18

datafiles: http://mikes-monitor.ad.vtt.fi/ftp/qdata/NTS_data/
"""

import numpy as np
import matplotlib.pyplot as plt
import h5py
import datetime

# Read TIC data
# UTC(MIKE)-WR(RISE) is D20-D00
with h5py.File('2026_NTS_data.hdf5', 'r') as hfile:
    t00 = hfile['D0/t'][()]  # timestamp
    dt00 = [datetime.datetime.fromtimestamp(x) for x in t00]  # datetime
    x00 = 1e9*hfile['D0/x'][()]  # time interval in ns

    t20 = hfile['D20/t'][()]  # timestamp
    dt20 = [datetime.datetime.fromtimestamp(x) for x in t20]  # datetime
    x20 = 1e9*hfile['D20/x'][()]  # time interval in ns

    offset = 192.22 # ns
    wr = x20 - np.interp(t20, t00, x00) + offset

    # first few datapoints ar junk, throw them away
    print(len(t20), len(wr))
    dt20 = dt20[10:]
    wr = wr[10:]

# Read WRS-monitoring data
with h5py.File('2026_NTS_wrs_monitoring.hdf5', 'r') as hfile:
    #print(hfile.keys())  # ['109.105.121.138', '109.105.121.139', '109.105.121.146']
    #print(hfile['109.105.121.138'].keys())  # ['t', 'wrsPtpRTT']

    t1 = hfile['/109.105.121.146/t'][()]  # timestamp
    dt1 = [datetime.datetime.fromtimestamp(x) for x in t1]  # datetime
    rtt1 = hfile['/109.105.121.146/wrsPtpRTT'][()]/1e3  # RTT in ns

    t2 = hfile['/109.105.121.138/t'][()]
    dt2 = [datetime.datetime.fromtimestamp(x) for x in t2]
    rtt2 = hfile['/109.105.121.138/wrsPtpRTT'][()]/1e3

    t3 = hfile['/109.105.121.139/t'][()]
    dt3 = [datetime.datetime.fromtimestamp(x) for x in t3]
    rtt3 = hfile['/109.105.121.139/wrsPtpRTT'][()]/1e3


#%%

plt.figure()
plt.subplot(4,1,1)
plt.plot(dt1, rtt1, label='rtt1 %.3f ns'%(np.mean(rtt1)))
plt.grid()
plt.legend()
plt.title('RTT(RISE-SUNET)')

plt.subplot(4,1,2)
plt.plot(dt2, rtt2, label='rtt2 %.3f ns'%(np.mean(rtt2)))
plt.grid()
plt.legend()
plt.ylabel('RTT / ns')
plt.title('RTT(SUNET-CSC)')

plt.subplot(4,1,3)
plt.plot(dt3, rtt3, label='rtt3 %.3f ns'%(np.mean(rtt3)))
plt.grid()
plt.legend()
plt.title('RTT(CSC-MIKES)')

plt.subplot(4,1,4)
plt.plot(dt20,wr,'r.', label='UTC(MIKE)-WR(RISE) %.3f ns'%(np.mean(wr)))
plt.grid()
plt.legend()
plt.ylabel('Phase / ns')
plt.title('UTC(MIKE)-WR(RISE)')