CI development

This MR will include many useful unit tests in the continuous integration. This long overdue upgrade will save a lot of manual effort in making plots for reviews, reduce the occurrence of bugs in the code, and make bugs easier to catch right away. Here are the tests that are now included:

• Run 3 versions of the calibration pipeline, and check that the right number of calibrated frames are produced
• Check the ASDs of GDS-CALIB_STRAIN and GDS-CALIB_STRAIN_CLEAN, and compare them to a "standard", nominally correct ASD. Plots are included in the CI job artifacts.
• Check the TDCFs, and compare them to nominally correct TDCFs. Plots included.
• Check the calibration state vector, and compare each bit to a nominally correct state vector. Plots included.
• Measure ratios of h(t) / Pcal and h(t) / Actuation injections, and compare to a nominally correct version. Time-series plots included.
• Use FrDiff to compare GDS-CALIB_STRAIN_CLEAN in frames produced by two otherwise identical pipelines with different start times. They are required to be identical.
• Measure the latency in a simulated low-latency-like configuration. Check that the mean latency is not too large, that there are not frequent or large excursions, and that the latency does not increase steadily with time.

Job artifacts from the most recent commit (at the time of writing) can be seen here, including plots and txt files with the data being plotted.

Getting all stages of the CI to work also required me to address some packaging issues. If I did things correctly, this MR Closes #7 (closed), #12 (closed), #14 (closed), #15 (closed), #26 (closed). I was not able to add @duncanmmacleod as a reviewer, but he created several of these issues.

tests/tests_pytest/error.py

-#!/usr/bin/env python3
-# Copyright (C) 2018  Aaron Viets, Alexander Bartoletti
+# Copyright (C) 2023  Alexander Bartoletti
 #
 # This program is free software; you can redistribute it and/or modify it
 # under the terms of the GNU General Public License as published by the
@@ -16,68 +15,133 @@
 # 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
 
 import numpy as np
+import sys
 
-#def rms(hoft_f='hoft_asd.txt', clean_f='clean_asd.txt'): 
 
 def rms(d_type, typ='A'):
-    pcal_freq_list = [102.1, 1083.3, 1083.7, 17.1, 283.9, 33.4, 410.2, 410.3, 56.4, 77.7]
-    act_freq_list = [15.6, 16.4, 17.6]
+    ifo = 'H1'
+    pcal_freq_list = [102.1, 1083.7, 17.1, 283.9, 33.4, 410.3, 77.7]
+    act_freq_list = [17.6, 16.4, 15.6]
+    act_names = ['TST_exc', 'PUM_exc', 'UIM_exc']
+    TDCF_names = ['KAPPA_TST_REAL', 'KAPPA_TST_IMAGINARY', 'KAPPA_PUM_REAL', 'KAPPA_PUM_IMAGINARY', 'KAPPA_UIM_REAL', 'KAPPA_UIM_IMAGINARY', 'KAPPA_C', 'F_CC', 'F_S_SQUARED', 'SRC_Q_INVERSE']
     st_f_list = []
     f_list = []
     if d_type[0].upper() == 'H':
-            name = 'hoft'
-            st_f_list.append('tests/tests_pytest/ASD_data/hoft_asd_standard.txt')
-            f_list.append('tests/tests_pytest/ASD_data/hoft_asd.txt')
+            name = ['hoft (Approx TDCFs)', 'hoft (Exact TDCFs)']
+            standard_n = 32769
+            st_f_list.append('tests/tests_pytest/ASD_data/Approx_hoft_asd_standard.txt')
+            f_list.append('tests/tests_pytest/ASD_data/Approx_hoft_asd.txt')
+            st_f_list.append('tests/tests_pytest/ASD_data/Exact_hoft_asd_standard.txt')
+            f_list.append('tests/tests_pytest/ASD_data/Exact_hoft_asd.txt')
     elif d_type[0].upper() == 'C':
-            name = 'clean'
-            st_f_list.append('tests/tests_pytest/ASD_data/clean_asd_standard.txt')
-            f_list.append('tests/tests_pytest/ASD_data/clean_asd.txt')
+            name = ['clean hoft (Approx TDCFs)', 'clean hoft (Exact TDCFs)']
+            standard_n = 32769
+            st_f_list.append('tests/tests_pytest/ASD_data/Approx_clean_asd_standard.txt')
+            f_list.append('tests/tests_pytest/ASD_data/Approx_clean_asd.txt')
+            st_f_list.append('tests/tests_pytest/ASD_data/Exact_clean_asd_standard.txt')
+            f_list.append('tests/tests_pytest/ASD_data/Exact_clean_asd.txt')
+    elif d_type[0].upper() == 'T':
+            name = []
+            standard_n = 59
+            for TDCF_name in TDCF_names:
+                    name.append("Approx %s" % TDCF_name)
+                    name.append("Exact %s" % TDCF_name)
+                    f_list.append("tests/tests_pytest/TDCF_data/%s_Approx.txt" % TDCF_name)
+                    st_f_list.append("tests/tests_pytest/TDCF_data/%s_Approx_standard.txt" % TDCF_name)
+                    f_list.append("tests/tests_pytest/TDCF_data/%s_Exact.txt" % TDCF_name)
+                    st_f_list.append("tests/tests_pytest/TDCF_data/%s_Exact_standard.txt" % TDCF_name)
     elif d_type[0].upper() == 'P':
-            #typ = input('Approx or Exact: ')
             if typ[0].upper() == 'A':
-                    name = 'Approx pcal'
+                    name = []
+                    standard_n = 57
                     for freq in pcal_freq_list:
-                            f_list.append('tests/tests_pytest/pcal_data/H1_GDS-CALIB_STRAIN_over_Pcal_0_at_%sHz.txt' % freq)
-                            st_f_list.append('tests/tests_pytest/pcal_data/H1_GDS-CALIB_STRAIN_over_Pcal_standard_0_at_%sHz.txt' % freq)
+                            name.append('Deltal / pcal (Approx TDCFs) at %0.1f Hz' % freq)
+                            f_list.append('tests/tests_pytest/pcal_data/%s_GDS-CALIB_STRAIN_over_Pcal_0_at_%0.1fHz.txt' % (ifo, freq))
+                            st_f_list.append('tests/tests_pytest/pcal_data/%s_GDS-CALIB_STRAIN_over_Pcal_standard_0_at_%0.1fHz.txt' % (ifo, freq))
             if typ[0].upper() == 'E':
-                    name = 'Exact pcal'
+                    name = []
+                    standard_n = 57
                     for freq in pcal_freq_list:
-                            f_list.append('tests/tests_pytest/pcal_data/H1_GDS-CALIB_STRAIN_over_Pcal_1_at_%sHz.txt' % freq)
-                            st_f_list.append('tests/tests_pytest/pcal_data/H1_GDS-CALIB_STRAIN_over_Pcal_standard_1_at_%sHz.txt' % freq)
+                            name.append('Deltal / pcal (Exact TDCFs) at %0.1f Hz' % freq)
+                            f_list.append('tests/tests_pytest/pcal_data/%s_GDS-CALIB_STRAIN_over_Pcal_1_at_%0.1fHz.txt' % (ifo, freq))
+                            st_f_list.append('tests/tests_pytest/pcal_data/%s_GDS-CALIB_STRAIN_over_Pcal_standard_1_at_%0.1fHz.txt' % (ifo, freq))
     elif d_type[0].upper() == 'A':
             if typ[0].upper() == 'A':
-                    name = 'Approx act'
-                    for freq in act_freq_list:
-                            f_list.append('tests/tests_pytest/act_data/H1_GDS-CALIB_STRAIN_over_Act_0_at_%sHz.txt' % freq)
-                            st_f_list.append('tests/tests_pytest/act_data/H1_GDS-CALIB_STRAIN_over_Act_standard_0_at_%sHz.txt' % freq)
+                    name = []
+                    standard_n = 57
+                    for i in range(len(act_freq_list)):
+                            name.append('Deltal / %s (Approx TDCFs) at %0.1f Hz' % (act_names[i], act_freq_list[i]))
+                            f_list.append('tests/tests_pytest/act_data/%s_GDS-CALIB_STRAIN_over_%s_0_at_%0.1fHz.txt' % (ifo, act_names[i], act_freq_list[i]))
+                            st_f_list.append('tests/tests_pytest/act_data/%s_GDS-CALIB_STRAIN_over_%s_standard_0_at_%0.1fHz.txt' % (ifo, act_names[i], act_freq_list[i]))
             if typ[0].upper() == 'E':
-                    name = 'Exact act'
-                    for freq in act_freq_list:
-                            f_list.append('tests/tests_pytest/act_data/H1_GDS-CALIB_STRAIN_over_Act_1_at_%sHz.txt' % freq)
-                            st_f_list.append('tests/tests_pytest/act_data/H1_GDS-CALIB_STRAIN_over_Act_standard_1_at_%sHz.txt' % freq)
+                    name = []
+                    standard_n = 57
+                    for i in range(len(act_freq_list)):
+                            name.append('Deltal / %s (Exact TDCFs) at %0.1f Hz' % (act_names[i], act_freq_list[i]))
+                            f_list.append('tests/tests_pytest/act_data/%s_GDS-CALIB_STRAIN_over_%s_1_at_%0.1fHz.txt' % (ifo, act_names[i], act_freq_list[i]))
+                            st_f_list.append('tests/tests_pytest/act_data/%s_GDS-CALIB_STRAIN_over_%s_standard_1_at_%0.1fHz.txt' % (ifo, act_names[i], act_freq_list[i]))
 
 
-    e_file = open('tests/tests_pytest/error_results.txt', 'a')
-    for f in st_f_list:
-            standard = np.loadtxt(f)
-    for f in f_list:
-            data = np.loadtxt(f)
-    standard = np.transpose(standard)[-1]
-    data = np.transpose(data)[-1]
-    sum_sq = 0
-    for i in range(len(standard-3)):
-            if standard[i] != 0: 
-                    sum_sq += abs(1 - data[i]/standard[i])
-    mean_sq = sum_sq/len(data)
-    rms = np.sqrt(mean_sq)
-    #try:
-    #assert rms < 1
-    #except AssertionError:
-    #s1 = 'bad ' + str(name) + ' data' + '\n'
-    #e_file.write(s1)
-    #finally:
-    s2 = str(name) + ' rms = ' + str(rms) + '\n'
-    e_file.write(s2)
-    e_file.close()
-    assert rms < 1
-#error between 1% and 10^-15
+    assert len(st_f_list) == len(f_list), "Mismatch between number of standard files and number of test files: %d != %d" % (len(st_f_list), len(f_list))
+    for i in range(len(f_list)):
+        standard = np.loadtxt(st_f_list[i])
+        data = np.loadtxt(f_list[i])
+        if standard_n == 57:
+            # Revove the last 64 s of data, which are affected by end of stream
+            data = data[:-64]
+        elif standard_n == 59:
+            # Revove the last 40 s of data, which are affected by end of stream
+            data = data[: -40 * 16]
+        sum_sq = 0
+        std_idx = 0
+        # Loop through the data and find values at the same time or frequency as the standard
+        for j in range(len(data)):
+            if data[j][0] == standard[std_idx][0]:
+                # Then check this data point
+                # If there are 3 columns in the file, then columns 2 and 3 are magnitude and phase
+                if(len(data[j]) == 3):
+                    data_value = data[j][1] * np.exp(1j * np.pi * data[j][2] / 180)
+                    standard_value = standard[std_idx][1] * np.exp(1j * np.pi * standard[std_idx][2] / 180)
+                    sum_sq += (abs(1 - data_value/standard_value))**2
+                else:
+                    data_value = data[j][1]
+                    standard_value = standard[std_idx][1]
+                    if 'SRC_Q' in name[i]:
+                        # Deal with asymptotic behavior and zero crossings
+                        if abs(data_value) > 1 or abs(standard_value) > 1:
+                            data_value = 1.0
+                            standard_value = 1.0
+                        data_value += 10
+                        standard_value += 10
+                    elif 'F_S_SQUARED' in name[i]:
+                        # Deal with zero crossings
+                        data_value += 100
+                        standard_value += 100
+                    elif 'IMAGINARY' in name[i]:
+                        # Deal with zero crossings
+                        data_value += 1
+                        standard_value += 1
+                    sum_sq += (abs(1 - data_value/standard_value))**2
+                std_idx += 1
+                if std_idx == standard_n:
+                    break
+
+        # Make sure the data sets overlapped as expected
+        assert j > 0.8 * len(data), "%s: failed to use enough testing data" % name[i]
+        if standard_n == 32769:
+            # Then it is an ASD, and we should use check every frequency
+            assert std_idx == standard_n, "%s: failed to use all of the standard data" % name[i]
+        else:
+            # It is time series data, so just make sure a reasonable minimum amount was checked
+            assert std_idx > 15, "Not enough data for a reliable test" # ~15-minute minimum requirement
+
+        mean_sq = sum_sq / std_idx
+        rms = np.sqrt(mean_sq)
+        e_file = open('tests/tests_pytest/error_results.txt', 'a')
+        s2 = str(name[i]) + ' rms error = ' + str(rms)
+        e_file.write(s2 + '\n')
+        e_file.close()
+        print(s2, file = sys.stderr)
+
+        assert rms < 0.002, "%s: error too large!" % name[i]
+