diff --git a/gstlal-calibration/gst/lal/gstlal_sensingtdcfs.c b/gstlal-calibration/gst/lal/gstlal_sensingtdcfs.c
index 2c44472679e8f3ac854e469bcc2771d97cfce52a..520a4baed62f782a516c8c3bf2fb0d51488932d4 100644
--- a/gstlal-calibration/gst/lal/gstlal_sensingtdcfs.c
+++ b/gstlal-calibration/gst/lal/gstlal_sensingtdcfs.c
@@ -69,7 +69,7 @@
"audio/x-raw, " \
"format = (string) {"GST_AUDIO_NE(Z64)", "GST_AUDIO_NE(Z128)"}, " \
"rate = " GST_AUDIO_RATE_RANGE ", " \
- "channels = (int) [4, 5], " \
+ "channels = (int) [1, MAX], " \
"layout = (string) interleaved, " \
"channel-mask = (bitmask) 0"
@@ -77,7 +77,7 @@
"audio/x-raw, " \
"format = (string) {"GST_AUDIO_NE(F32)", "GST_AUDIO_NE(F64)"}, " \
"rate = " GST_AUDIO_RATE_RANGE ", " \
- "channels = (int) [4, 6], " \
+ "channels = (int) [2, MAX], " \
"layout = (string) interleaved, " \
"channel-mask = (bitmask) 0"
@@ -144,43 +144,106 @@ static void set_metadata(GSTLALSensingTDCFs *element, GstBuffer *buf, guint64 ou
#define DEFINE_KAPPA_C_0(DTYPE, F_OR_NOT) \
-void kappa_C_0_ ## DTYPE(DTYPE Gres1_real, DTYPE Gres1_imag, DTYPE Gres2_real, DTYPE Gres2_imag, DTYPE Y1_real, DTYPE Y1_imag, DTYPE Y2_real, DTYPE Y2_imag, DTYPE f1, DTYPE f2, DTYPE *kc1, DTYPE *kc2, DTYPE *kc3, DTYPE *kc4) { \
+void kappa_C_0_ ## DTYPE(DTYPE Cinv_tdep1_real, DTYPE Cinv_tdep1_imag, DTYPE Cinv_tdep2_real, DTYPE Cinv_tdep2_imag, DTYPE f1, DTYPE f2, DTYPE *kc1, DTYPE *kc2, DTYPE *kc3, DTYPE *kc4) { \
\
- DTYPE H1, H2, I1, I2, Xi, Zeta, a, b, c, d, e, Delta0, Delta1, p, q; \
- complex DTYPE Q0, S0, complex_kc1, complex_kc2, complex_kc3, complex_kc4; \
- H1 = f1 * f1 * (Gres1_real - Y1_real); \
- H2 = f2 * f2 * (Gres2_real - Y2_real); \
- I1 = f1 * f1 * (Gres1_imag - Y1_imag); \
- I2 = f2 * f2 * (Gres2_imag - Y2_imag); \
+ DTYPE H1, H2, I1, I2, Xi, Zeta; \
+ /*
+ * Due to numerical instabilities, we typecast to long double precision
+ * here to avoid potentially large numerical errors.
+ */ \
+ complex long double Q0, S0, complex_kc1, complex_kc2, complex_kc3, complex_kc4; \
+ long double a, b, c, d, e, Delta0, Delta1, p, q; \
+ H1 = -f1 * f1 * Cinv_tdep1_real; \
+ H2 = -f2 * f2 * Cinv_tdep2_real; \
+ I1 = -f1 * f1 * Cinv_tdep1_imag; \
+ I2 = -f2 * f2 * Cinv_tdep2_imag; \
Xi = (f2 * f2 * H1 - f1 * f1 * H2) / (f1 * f1 - f2 * f2); \
Zeta = I1 / f1 - I2 / f2; \
- a = f2 * f2 * Xi * (H2 + Xi) * Zeta * Zeta; \
- b = f2 * (f2 * f2 - f1 * f1) * (H2 + Xi) * Zeta * I2 + pow ## F_OR_NOT(f2, 4) * (H2 + 2 * Xi) * Zeta * Zeta; \
- c = pow ## F_OR_NOT(f2, 3) * (f2 * f2 - f1 * f1) * Zeta * I2 + pow ## F_OR_NOT(f2 * f2 - f1 * f1, 2) * pow ## F_OR_NOT(H2 + Xi, 2) + pow ## F_OR_NOT(f2, 6) * Zeta * Zeta; \
- d = 2 * f2 * f2 * pow ## F_OR_NOT(f2 * f2 - f1 * f1, 2) * (H2 + Xi); \
- e = pow ## F_OR_NOT(f2, 4) * pow ## F_OR_NOT(f2 * f2 - f1 * f1, 2); \
+ a = (long double) f2 * f2 * Xi * (H2 + Xi) * Zeta * Zeta; \
+ b = (long double) f2 * (f2 * f2 - f1 * f1) * (H2 + Xi) * Zeta * I2 + pow ## F_OR_NOT(f2, 4) * (H2 + 2 * Xi) * Zeta * Zeta; \
+ c = (long double) pow ## F_OR_NOT(f2, 3) * (f2 * f2 - f1 * f1) * Zeta * I2 + pow ## F_OR_NOT(f2 * f2 - f1 * f1, 2) * pow ## F_OR_NOT(H2 + Xi, 2) + pow ## F_OR_NOT(f2, 6) * Zeta * Zeta; \
+ d = (long double) 2 * f2 * f2 * pow ## F_OR_NOT(f2 * f2 - f1 * f1, 2) * (H2 + Xi); \
+ e = (long double) pow ## F_OR_NOT(f2, 4) * pow ## F_OR_NOT(f2 * f2 - f1 * f1, 2); \
Delta0 = c * c - 3 * b * d + 12 * a * e; \
- Delta1 = 2 * pow ## F_OR_NOT(c, 3) - 9 * b * c * d + 27 * b * b * e + 27 * a * d * d - 72 * a * c * e; \
+ Delta1 = 2 * powl(c, 3) - 9 * b * c * d + 27 * b * b * e + 27 * a * d * d - 72 * a * c * e; \
p = (8 * a * c - 3 * b * b) / (8 * a * a); \
- q = (pow ## F_OR_NOT(b, 3) - 4 * a * b * c + 8 * a * a * d) / (8 * pow ## F_OR_NOT(a, 3)); \
- Q0 = cpow ## F_OR_NOT(0.5 * ((complex DTYPE) Delta1 + cpow ## F_OR_NOT((complex DTYPE) (Delta1 * Delta1 - 4 * pow ## F_OR_NOT(Delta0, 3)), 0.5)), 1.0 / 3.0); \
- S0 = 0.5 * cpow ## F_OR_NOT((-2 * p) / 3.0 + (Q0 + Delta0 / Q0) / (3 * a), 0.5); \
+ q = (powl(b, 3) - 4 * a * b * c + 8 * a * a * d) / (8 * powl(a, 3)); \
+ Q0 = cpowl(0.5 * ((complex long double) Delta1 + cpowl((complex long double) (Delta1 * Delta1 - 4 * powl(Delta0, 3)), 0.5)), 1.0 / 3.0); \
+ S0 = 0.5 * cpowl((-2 * p) / 3.0 + (Q0 + Delta0 / Q0) / (3 * a), 0.5); \
+ \
+ complex_kc1 = -b / (4 * a) - S0 + 0.5 * cpowl(-4 * S0 * S0 - 2 * p + q / S0, 0.5); \
+ complex_kc2 = -b / (4 * a) - S0 - 0.5 * cpowl(-4 * S0 * S0 - 2 * p + q / S0, 0.5); \
+ complex_kc3 = -b / (4 * a) + S0 + 0.5 * cpowl(-4 * S0 * S0 - 2 * p - q / S0, 0.5); \
+ complex_kc4 = -b / (4 * a) + S0 - 0.5 * cpowl(-4 * S0 * S0 - 2 * p - q / S0, 0.5); \
+ \
+ /* DEBUGGING */ \
+ /*complex long double fcc1 = (f2 * f2 - f1 * f1) / (complex_kc1 * (f2 * Cinv_tdep2_imag - f1 * Cinv_tdep1_imag)); \
+ complex long double fs_squared1 = complex_kc1 * (Cinv_tdep2_real - Cinv_tdep1_real) / (1.0 / (f2 * f2) - 1.0 / (f1 * f1)); \
+ complex long double fs_over_Q1 = fcc1 * (complex_kc1 * Cinv_tdep1_real - 1.0 - fs_squared1 / (f1 * f1)); \
+ complex long double fcc2 = (f2 * f2 - f1 * f1) / (complex_kc2 * (f2 * Cinv_tdep2_imag - f1 * Cinv_tdep1_imag)); \
+ complex long double fs_squared2 = complex_kc2 * (Cinv_tdep2_real - Cinv_tdep1_real) / (1.0 / (f2 * f2) - 1.0 / (f1 * f1)); \
+ complex long double fs_over_Q2 = fcc2 * (complex_kc2 * Cinv_tdep1_real - 1.0 - fs_squared2 / (f1 * f1)); \
+ complex long double fcc3 = (f2 * f2 - f1 * f1) / (complex_kc3 * (f2 * Cinv_tdep2_imag - f1 * Cinv_tdep1_imag)); \
+ complex long double fs_squared3 = complex_kc3 * (Cinv_tdep2_real - Cinv_tdep1_real) / (1.0 / (f2 * f2) - 1.0 / (f1 * f1)); \
+ complex long double fs_over_Q3 = fcc3 * (complex_kc3 * Cinv_tdep1_real - 1.0 - fs_squared3 / (f1 * f1)); \
+ complex long double fcc4 = (f2 * f2 - f1 * f1) / (complex_kc4 * (f2 * Cinv_tdep2_imag - f1 * Cinv_tdep1_imag)); \
+ complex long double fs_squared4 = complex_kc4 * (Cinv_tdep2_real - Cinv_tdep1_real) / (1.0 / (f2 * f2) - 1.0 / (f1 * f1)); \
+ complex long double fs_over_Q4 = fcc4 * (complex_kc4 * Cinv_tdep1_real - 1.0 - fs_squared4 / (f1 * f1)); \
\
- complex_kc1 = -b / (4 * a) + S0 + 0.5 * cpow ## F_OR_NOT(-4 * S0 * S0 - 2 * p - q / S0, 0.5); \
- complex_kc2 = -b / (4 * a) + S0 - 0.5 * cpow ## F_OR_NOT(-4 * S0 * S0 - 2 * p - q / S0, 0.5); \
- complex_kc3 = -b / (4 * a) - S0 + 0.5 * cpow ## F_OR_NOT(-4 * S0 * S0 - 2 * p + q / S0, 0.5); \
- complex_kc4 = -b / (4 * a) - S0 - 0.5 * cpow ## F_OR_NOT(-4 * S0 * S0 - 2 * p + q / S0, 0.5); \
+ long double Delta = 256 * a * a * a * e * e * e - 192 * a * a * b * d * e * e - 128 * a * a * c * c * e * e + 144 * a * a * c * d * d * e - 27 * a * a * d * d * d * d + 144 * a * b * b * c * e * e - 6 * a * b * b * d * d * e - 80 * a * b * c * c * d * e + 18 * a * b * c * d * d * d + 16 * a * c * c * c * c * e - 4 * a * c * c * c * d * d - 27 * b * b * b * b * e * e + 18 * b * b * b * c * d * e - 4 * b * b * b * d * d * d - 4 * b * b * c * c * c * e + b * b * c * c * d * d; \
+ long double P = 8 * a * c - 3 * b * b; \
+ long double D = 64 * a * a * a * e - 16 * a * a * c * c + 16 * a * b * b * c - 16 * a * a * b * d - 3 * b * b * b * b; \
+ \
+ complex DTYPE condition1 = Cinv_tdep1_real + I * Cinv_tdep1_imag; \
+ complex DTYPE condition2 = Cinv_tdep2_real + I * Cinv_tdep2_imag; \
+ complex long double solution11 = (1 + I * f1 / fcc1) / complex_kc1 * (f1 * f1 + fs_squared1 - I * f1 * fs_over_Q1) / (f1 * f1); \
+ complex long double solution12 = (1 + I * f2 / fcc1) / complex_kc1 * (f2 * f2 + fs_squared1 - I * f2 * fs_over_Q1) / (f2 * f2); \
+ complex long double solution21 = (1 + I * f1 / fcc2) / complex_kc2 * (f1 * f1 + fs_squared2 - I * f1 * fs_over_Q2) / (f1 * f1); \
+ complex long double solution22 = (1 + I * f2 / fcc2) / complex_kc2 * (f2 * f2 + fs_squared2 - I * f2 * fs_over_Q2) / (f2 * f2); \
+ complex long double solution31 = (1 + I * f1 / fcc3) / complex_kc3 * (f1 * f1 + fs_squared3 - I * f1 * fs_over_Q3) / (f1 * f1); \
+ complex long double solution32 = (1 + I * f2 / fcc3) / complex_kc3 * (f2 * f2 + fs_squared3 - I * f2 * fs_over_Q3) / (f2 * f2); \
+ complex long double solution41 = (1 + I * f1 / fcc4) / complex_kc4 * (f1 * f1 + fs_squared4 - I * f1 * fs_over_Q4) / (f1 * f1); \
+ complex long double solution42 = (1 + I * f2 / fcc4) / complex_kc4 * (f2 * f2 + fs_squared4 - I * f2 * fs_over_Q4) / (f2 * f2); \
+ long double diff1 = cabsl(1.0 - condition1 / solution11) + cabsl(1.0 - condition2 / solution12); \
+ long double diff2 = cabsl(1.0 - condition1 / solution21) + cabsl(1.0 - condition2 / solution22); \
+ long double diff3 = cabsl(1.0 - condition1 / solution31) + cabsl(1.0 - condition2 / solution32); \
+ long double diff4 = cabsl(1.0 - condition1 / solution41) + cabsl(1.0 - condition2 / solution42); \
+ complex long double zero1 = a * complex_kc1 * complex_kc1 * complex_kc1 * complex_kc1 + b * complex_kc1 * complex_kc1 * complex_kc1 + c * complex_kc1 * complex_kc1 + d * complex_kc1 + e; \
+ complex long double zero2 = a * complex_kc2 * complex_kc2 * complex_kc2 * complex_kc2 + b * complex_kc2 * complex_kc2 * complex_kc2 + c * complex_kc2 * complex_kc2 + d * complex_kc2 + e; \
+ complex long double zero3 = a * complex_kc3 * complex_kc3 * complex_kc3 * complex_kc3 + b * complex_kc3 * complex_kc3 * complex_kc3 + c * complex_kc3 * complex_kc3 + d * complex_kc3 + e; \
+ complex long double zero4 = a * complex_kc4 * complex_kc4 * complex_kc4 * complex_kc4 + b * complex_kc4 * complex_kc4 * complex_kc4 + c * complex_kc4 * complex_kc4 + d * complex_kc4 + e; \
+ \
+ if(cabsl(Delta) < 1e109) { \
+ g_print("\nDelta=%Le\nP=%Le\nD=%Le\n", Delta, P, D); \
+ g_print("kc1=%Le+1j*%Le\nfcc1=%Le+1j*%Le\nfs_squared1=%Le+1j*%Le\nfs_over_Q1=%Le+1j*%Le\n", creall(complex_kc1), cimagl(complex_kc1), creall(fcc1), cimagl(fcc1), creall(fs_squared1), cimagl(fs_squared1), creall(fs_over_Q1), cimagl(fs_over_Q1)); \
+ g_print("kc2=%Le+1j*%Le\nfcc2=%Le+1j*%Le\nfs_squared2=%Le+1j*%Le\nfs_over_Q2=%Le+1j*%Le\n", creall(complex_kc2), cimagl(complex_kc2), creall(fcc2), cimagl(fcc2), creall(fs_squared2), cimagl(fs_squared2), creall(fs_over_Q2), cimagl(fs_over_Q2)); \
+ g_print("kc3=%Le+1j*%Le\nfcc3=%Le+1j*%Le\nfs_squared3=%Le+1j*%Le\nfs_over_Q3=%Le+1j*%Le\n", creall(complex_kc3), cimagl(complex_kc3), creall(fcc3), cimagl(fcc3), creall(fs_squared3), cimagl(fs_squared3), creall(fs_over_Q3), cimagl(fs_over_Q3)); \
+ g_print("kc4=%Le+1j*%Le\nfcc4=%Le+1j*%Le\nfs_squared4=%Le+1j*%Le\nfs_over_Q4=%Le+1j*%Le\n", creall(complex_kc4), cimagl(complex_kc4), creall(fcc4), cimagl(fcc4), creall(fs_squared4), cimagl(fs_squared4), creall(fs_over_Q4), cimagl(fs_over_Q4)); \
+ if(diff1 > 1e-12) \
+ g_print("Solution 1 failed: diff=%Le\n", diff1); \
+ if(diff2 > 1e-12) \
+ g_print("Solution 2 failed: diff=%Le\n", diff2); \
+ if(diff3 > 1e-12) \
+ g_print("Solution 3 failed: diff=%Le\n", diff3); \
+ if(diff4 > 1e-12) \
+ g_print("Solution 4 failed: diff=%Le\n", diff4); \
+ g_print("a=%Le b=%Le c=%Le d=%Le e=%Le\n", a, b, c, d, e); \
+ g_print("Solution 1: 0 = %Le+I%Le\n", creall(zero1), cimagl(zero1)); \
+ g_print("Solution 2: 0 = %Le+I%Le\n", creall(zero2), cimagl(zero2)); \
+ g_print("Solution 3: 0 = %Le+I%Le\n", creall(zero3), cimagl(zero3)); \
+ g_print("Solution 4: 0 = %Le+I%Le\n\n", creall(zero4), cimagl(zero4)); \
+ } */ \
\
/* Any solution that is clearly complex or negative is not the solution we want */ \
- if(creal ## F_OR_NOT(complex_kc1) > 0.0 && fabs ## F_OR_NOT(cimag ## F_OR_NOT(complex_kc1) / creal ## F_OR_NOT(complex_kc1)) < 1e-3) \
- *kc1 = creal ## F_OR_NOT(complex_kc1); \
+ if(creall(complex_kc1) > 0.0 && fabsl(cimagl(complex_kc1) / creall(complex_kc1)) < 1e-3) \
+ *kc1 = (DTYPE) creall(complex_kc1); \
\
- if(creal ## F_OR_NOT(complex_kc2) > 0.0 && fabs ## F_OR_NOT(cimag ## F_OR_NOT(complex_kc2) / creal ## F_OR_NOT(complex_kc2)) < 1e-3) \
- *kc2 = creal ## F_OR_NOT(complex_kc2); \
- if(creal ## F_OR_NOT(complex_kc3) > 0.0 && fabs ## F_OR_NOT(cimag ## F_OR_NOT(complex_kc3) / creal ## F_OR_NOT(complex_kc3)) < 1e-3) \
- *kc3 = creal ## F_OR_NOT(complex_kc3); \
- if(creal ## F_OR_NOT(complex_kc4) > 0.0 && fabs ## F_OR_NOT(cimag ## F_OR_NOT(complex_kc4) / creal ## F_OR_NOT(complex_kc4)) < 1e-3) \
- *kc4 = creal ## F_OR_NOT(complex_kc4); \
+ if(creall(complex_kc2) > 0.0 && fabsl(cimagl(complex_kc2) / creall(complex_kc2)) < 1e-3) \
+ *kc2 = (DTYPE) creall(complex_kc2); \
+ if(creall(complex_kc3) > 0.0 && fabsl(cimagl(complex_kc3) / creall(complex_kc3)) < 1e-3) \
+ *kc3 = (DTYPE) creall(complex_kc3); \
+ if(creall(complex_kc4) > 0.0 && fabsl(cimagl(complex_kc4) / creall(complex_kc4)) < 1e-3) \
+ *kc4 = (DTYPE) creall(complex_kc4); \
\
return; \
}
@@ -196,9 +259,9 @@ DEFINE_KAPPA_C_0(double, );
#define DEFINE_F_CC_0(DTYPE) \
-DTYPE f_cc_0_ ## DTYPE(DTYPE Gres1_imag, DTYPE Gres2_imag, DTYPE Y1_imag, DTYPE Y2_imag, DTYPE f1, DTYPE f2, DTYPE kappa_C) { \
+DTYPE f_cc_0_ ## DTYPE(DTYPE Cinv_tdep1_imag, DTYPE Cinv_tdep2_imag, DTYPE f1, DTYPE f2, DTYPE kappa_C) { \
\
- return (f2 * f2 - f1 * f1) / (kappa_C * (f1 * (Gres1_imag - Y1_imag) - f2 * (Gres2_imag - Y2_imag))); \
+ return (f2 * f2 - f1 * f1) / (kappa_C * (f2 * Cinv_tdep2_imag - f1 * Cinv_tdep1_imag)); \
}
@@ -212,9 +275,9 @@ DEFINE_F_CC_0(double);
#define DEFINE_F_S_SQUARED_0(DTYPE) \
-DTYPE f_s_squared_0_ ## DTYPE(DTYPE Gres1_real, DTYPE Gres2_real, DTYPE Y1_real, DTYPE Y2_real, DTYPE f1, DTYPE f2, DTYPE kappa_C) { \
+DTYPE f_s_squared_0_ ## DTYPE(DTYPE Cinv_tdep1_real, DTYPE Cinv_tdep2_real, DTYPE f1, DTYPE f2, DTYPE kappa_C) { \
\
- return kappa_C * (Y2_real - Gres2_real - Y1_real + Gres1_real) / (1.0 / (f2 * f2) - 1.0 / (f1 * f1)); \
+ return kappa_C * (Cinv_tdep2_real - Cinv_tdep1_real) / (1.0 / (f2 * f2) - 1.0 / (f1 * f1)); \
}
@@ -228,9 +291,9 @@ DEFINE_F_S_SQUARED_0(double);
#define DEFINE_FS_OVER_Q_0(DTYPE) \
-DTYPE f_s_over_Q_0_ ## DTYPE(DTYPE Gres1_real, DTYPE Y1_real, DTYPE f1, DTYPE kappa_C, DTYPE f_cc, DTYPE f_s_squared) { \
+DTYPE f_s_over_Q_0_ ## DTYPE(DTYPE Cinv_tdep1_real, DTYPE f1, DTYPE kappa_C, DTYPE f_cc, DTYPE f_s_squared) { \
\
- return f_cc * (kappa_C * Y1_real - 1.0 - f_s_squared / (f1 * f1) - kappa_C * Gres1_real); \
+ return f_cc * (kappa_C * Cinv_tdep1_real - 1.0 - f_s_squared / (f1 * f1)); \
}
@@ -239,20 +302,18 @@ DEFINE_FS_OVER_Q_0(double);
#define DEFINE_FIND_BEST_SOLUTION(DTYPE, F_OR_NOT) \
-guint find_best_solution_ ## DTYPE(DTYPE kc1, DTYPE kc2, DTYPE kc3, DTYPE kc4, DTYPE fcc1, DTYPE fcc2, DTYPE fcc3, DTYPE fcc4, DTYPE fs_squared1, DTYPE fs_squared2, DTYPE fs_squared3, DTYPE fs_squared4, DTYPE fs_over_Q1, DTYPE fs_over_Q2, DTYPE fs_over_Q3, DTYPE fs_over_Q4, DTYPE f1, DTYPE f2, DTYPE complex tdep_sensing_at_f1, complex DTYPE tdep_sensing_at_f2) { \
+guint find_best_solution_ ## DTYPE(DTYPE kc1, DTYPE kc2, DTYPE kc3, DTYPE kc4, DTYPE fcc1, DTYPE fcc2, DTYPE fcc3, DTYPE fcc4, DTYPE fs_squared1, DTYPE fs_squared2, DTYPE fs_squared3, DTYPE fs_squared4, DTYPE fs_over_Q1, DTYPE fs_over_Q2, DTYPE fs_over_Q3, DTYPE fs_over_Q4, DTYPE f1, DTYPE f2, DTYPE complex Cinv_tdep1, complex DTYPE tdep_sensing_at_f2) { \
\
/*
- * In general, more than one of the four solutions can be correct. This is because
- * the variable gain of the sensing function is kappa_C / f_s^2, and there are 3
- * time-dependent poles, which depend on f_cc, f_s, and Q. The algorithm does not
- * care which variable name is assigned to which pole or how the gain is
- * distributed, but we need to care because there is a physical meaning to each
- * variable. We require that:
+ * Any of the four solutions can be "correct". The variable gain of the sensing
+ * function is kappa_C / f_s^2, and there are 3 time-dependent poles, which depend
+ * on f_cc, f_s, and Q. The algorithm does not care which variable name is
+ * assigned to which pole or how the gain is distributed, but we need to care
+ * because we associate a physical meaning with each variable. We require that:
* 1. kappa_C > 0
* 2. f_cc is the highest-frequency pole of the sensing function.
* Of the remaining solutions, we choose the one with the smallest RMS error at
- * the two Pcal line frequencies. Errors can be caused by numerical instabilities
- * or by the fact that some solutions may not be physically correct.
+ * the two Pcal line frequencies. Errors can be caused by numerical instabilities.
*/ \
DTYPE error1, error2, error3, error4, lowest_error = 2.0; \
complex DTYPE fs1, fs2, fs3, fs4, Qinv1, Qinv2, Qinv3, Qinv4; \
@@ -263,7 +324,7 @@ guint find_best_solution_ ## DTYPE(DTYPE kc1, DTYPE kc2, DTYPE kc3, DTYPE kc4, D
Qinv1 = fs_over_Q1 / fs1; \
/* Condition 2 */ \
if(fcc1 > cabs ## F_OR_NOT(fs1 / 2 * (Qinv1 + cpow ## F_OR_NOT(cpow ## F_OR_NOT(Qinv1, 2) + 4, 0.5))) && fcc1 > cabs ## F_OR_NOT(fs1 / 2 * (Qinv1 - cpow ## F_OR_NOT(cpow ## F_OR_NOT(Qinv1, 2) + 4, 0.5)))) { \
- error1 = cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f1 / fcc1) / kc1) * ((f1 * f1 + fs_squared1 - I * f1 * fs_over_Q1) / (f1 * f1)) / tdep_sensing_at_f1 - 1), 2); \
+ error1 = cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f1 / fcc1) / kc1) * ((f1 * f1 + fs_squared1 - I * f1 * fs_over_Q1) / (f1 * f1)) / Cinv_tdep1 - 1), 2); \
error1 += cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f2 / fcc1) / kc1) * ((f2 * f2 + fs_squared1 - I * f2 * fs_over_Q1) / (f2 * f2)) / tdep_sensing_at_f2 - 1), 2); \
if(error1 < lowest_error) { \
lowest_error = error1; \
@@ -275,7 +336,7 @@ guint find_best_solution_ ## DTYPE(DTYPE kc1, DTYPE kc2, DTYPE kc3, DTYPE kc4, D
fs2 = cpow ## F_OR_NOT((complex DTYPE) fs_squared2, 0.5); \
Qinv2 = fs_over_Q2 / fs2; \
if(fcc2 > cabs ## F_OR_NOT(fs2 / 2 * (Qinv2 + cpow ## F_OR_NOT(cpow ## F_OR_NOT(Qinv2, 2) + 4, 0.5))) && fcc2 > cabs ## F_OR_NOT(fs2 / 2 * (Qinv2 - cpow ## F_OR_NOT(cpow ## F_OR_NOT(Qinv2, 2) + 4, 0.5)))) { \
- error2 = cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f1 / fcc2) / kc2) * ((f1 * f1 + fs_squared2 - I * f1 * fs_over_Q2) / (f1 * f1)) / tdep_sensing_at_f1 - 1), 2); \
+ error2 = cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f1 / fcc2) / kc2) * ((f1 * f1 + fs_squared2 - I * f1 * fs_over_Q2) / (f1 * f1)) / Cinv_tdep1 - 1), 2); \
error2 += cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f2 / fcc2) / kc2) * ((f2 * f2 + fs_squared2 - I * f2 * fs_over_Q2) / (f2 * f2)) / tdep_sensing_at_f2 - 1), 2); \
if(error2 < lowest_error) { \
lowest_error = error2; \
@@ -287,7 +348,7 @@ guint find_best_solution_ ## DTYPE(DTYPE kc1, DTYPE kc2, DTYPE kc3, DTYPE kc4, D
fs3 = cpow ## F_OR_NOT((complex DTYPE) fs_squared3, 0.5); \
Qinv3 = fs_over_Q3 / fs3; \
if(fcc3 > cabs ## F_OR_NOT(fs3 / 2 * (Qinv3 + cpow ## F_OR_NOT(cpow ## F_OR_NOT(Qinv3, 2) + 4, 0.5))) && fcc3 > cabs ## F_OR_NOT(fs3 / 2 * (Qinv3 - cpow ## F_OR_NOT(cpow ## F_OR_NOT(Qinv3, 2) + 4, 0.5)))) { \
- error3 = cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f1 / fcc3) / kc3) * ((f1 * f1 + fs_squared3 - I * f1 * fs_over_Q3) / (f1 * f1)) / tdep_sensing_at_f1 - 1), 2); \
+ error3 = cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f1 / fcc3) / kc3) * ((f1 * f1 + fs_squared3 - I * f1 * fs_over_Q3) / (f1 * f1)) / Cinv_tdep1 - 1), 2); \
error3 += cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f2 / fcc3) / kc3) * ((f2 * f2 + fs_squared3 - I * f2 * fs_over_Q3) / (f2 * f2)) / tdep_sensing_at_f2 - 1), 2); \
if(error3 < lowest_error) { \
lowest_error = error3; \
@@ -299,7 +360,7 @@ guint find_best_solution_ ## DTYPE(DTYPE kc1, DTYPE kc2, DTYPE kc3, DTYPE kc4, D
fs4 = cpow ## F_OR_NOT((complex DTYPE) fs_squared4, 0.5); \
Qinv4 = fs_over_Q4 / fs4; \
if(fcc4 > cabs ## F_OR_NOT(fs4 / 2 * (Qinv4 + cpow ## F_OR_NOT(cpow ## F_OR_NOT(Qinv4, 2) + 4, 0.5))) && fcc4 > cabs ## F_OR_NOT(fs4 / 2 * (Qinv4 - cpow ## F_OR_NOT(cpow ## F_OR_NOT(Qinv4, 2) + 4, 0.5)))) { \
- error4 = cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f1 / fcc4) / kc4) * ((f1 * f1 + fs_squared4 - I * f1 * fs_over_Q4) / (f1 * f1)) / tdep_sensing_at_f1 - 1), 2); \
+ error4 = cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f1 / fcc4) / kc4) * ((f1 * f1 + fs_squared4 - I * f1 * fs_over_Q4) / (f1 * f1)) / Cinv_tdep1 - 1), 2); \
error4 += cpow ## F_OR_NOT(cabs ## F_OR_NOT(((1 + I * f2 / fcc4) / kc4) * ((f2 * f2 + fs_squared4 - I * f2 * fs_over_Q4) / (f2 * f2)) / tdep_sensing_at_f2 - 1), 2); \
if(error4 < lowest_error) { \
lowest_error = error4; \
@@ -361,9 +422,7 @@ static GstCaps *transform_caps(GstBaseTransform *trans, GstPadDirection directio
/*
* We know the caps on the source pad, and we want to put constraints on
* the sink pad caps. The sink pad caps are complex, while the source pad
- * caps are real. If there are 4 channels on the source pad, there should
- * be 4 on the sink pad as well. If there are 6 on the source pad, there
- * should be 5 on the sink pad. There are no other possibilities.
+ * caps are real with twice as many channels.
*/
for(n = 0; n < gst_caps_get_size(caps); n++) {
GstStructure *str = gst_caps_get_structure(caps, n);
@@ -372,25 +431,14 @@ static GstCaps *transform_caps(GstBaseTransform *trans, GstPadDirection directio
gint channels_out_min, channels_out_max;
channels_out_min = gst_value_get_int_range_min(v);
channels_out_max = gst_value_get_int_range_max(v);
- g_assert_cmpint(channels_out_min, <=, 6);
- g_assert_cmpint(channels_out_max, >=, 4);
- if(channels_out_min >= 5)
- /* Then we know there must be 5 on the sink pad */
- gst_structure_set(str, "channels", G_TYPE_INT, 5, NULL);
- else if(channels_out_max <= 4)
- /* Then we know there must be 4 on the sink pad */
- gst_structure_set(str, "channels", G_TYPE_INT, 4, NULL);
- else
- gst_structure_set(str, "channels", GST_TYPE_INT_RANGE, 4, 5, NULL);
+ g_assert_cmpint(channels_out_max, >=, 2);
+ gst_structure_set(str, "channels", GST_TYPE_INT_RANGE, channels_out_min / 2, channels_out_max / 2, NULL);
} else if(G_VALUE_HOLDS_INT(v)) {
gint channels_out = g_value_get_int(v);
- if(channels_out == 4)
- gst_structure_set(str, "channels", G_TYPE_INT, 4, NULL);
- else if(channels_out == 6)
- gst_structure_set(str, "channels", G_TYPE_INT, 5, NULL);
- else
- GST_ELEMENT_ERROR(trans, CORE, NEGOTIATION, (NULL), ("invalid number of channels in caps"));
+ /* The number of output channels must be even. */
+ g_assert(channels_out % 2 == 0);
+ gst_structure_set(str, "channels", G_TYPE_INT, channels_out / 2, NULL);
} else
GST_ELEMENT_ERROR(trans, CORE, NEGOTIATION, (NULL), ("invalid type for channels in caps"));
@@ -414,9 +462,7 @@ static GstCaps *transform_caps(GstBaseTransform *trans, GstPadDirection directio
/*
* We know the caps on the sink pad, and we want to put constraints on
* the source pad caps. The sink pad caps are complex, while the source pad
- * caps are real. If there are 4 channels on the sink pad, there should
- * be 4 on the source pad as well. If there are 5 on the sink pad, there
- * should be 6 on the source pad. There are no other possibilities.
+ * caps are real with twice as many channels.
*/
for(n = 0; n < gst_caps_get_size(caps); n++) {
GstStructure *str = gst_caps_get_structure(caps, n);
@@ -425,26 +471,12 @@ static GstCaps *transform_caps(GstBaseTransform *trans, GstPadDirection directio
gint channels_in_min, channels_in_max;
channels_in_min = gst_value_get_int_range_min(v);
channels_in_max = gst_value_get_int_range_max(v);
- g_assert_cmpint(channels_in_min, <=, 5);
- g_assert_cmpint(channels_in_max, >=, 4);
- if(channels_in_min == 5)
- /* Then we know there must be 6 channels on the source pad */
- gst_structure_set(str, "channels", G_TYPE_INT, 6, NULL);
- else if(channels_in_max == 4)
- /* Then we know there must be 4 channels on the source pad */
- gst_structure_set(str, "channels", G_TYPE_INT, 4, NULL);
- else
- gst_structure_set(str, "channels", GST_TYPE_INT_RANGE, 4, 6, NULL);
+ gst_structure_set(str, "channels", GST_TYPE_INT_RANGE, 2 * channels_in_min, 2 * channels_in_max, NULL);
} else if(G_VALUE_HOLDS_INT(v)) {
gint channels_in;
channels_in = g_value_get_int(v);
- if(channels_in == 4)
- gst_structure_set(str, "channels", G_TYPE_INT, 4, NULL);
- else if(channels_in == 5)
- gst_structure_set(str, "channels", G_TYPE_INT, 6, NULL);
- else
- GST_ELEMENT_ERROR(trans, CORE, NEGOTIATION, (NULL), ("invalid number of channels in caps"));
+ gst_structure_set(str, "channels", G_TYPE_INT, 2 * channels_in, NULL);
} else
GST_ELEMENT_ERROR(trans, CORE, NEGOTIATION, (NULL), ("invalid type for channels in caps"));
@@ -528,20 +560,7 @@ static gboolean set_caps(GstBaseTransform *trans, GstCaps *incaps, GstCaps *outc
}
/* Requirements for channels */
- if(channels_in == 4) {
- g_assert_cmpint(channels_out, ==, 4);
- if(element->sensing_model != 0 && element->sensing_model != 1) {
- GST_WARNING_OBJECT(element, "When there are 4 input channels, sensing-model must be either 0 or 1. Resetting sensing-model to 0.");
- element->sensing_model = 0;
- }
- } else if(channels_in == 5) {
- g_assert_cmpint(channels_out, ==, 6);
- if(element->sensing_model != 2) {
- GST_WARNING_OBJECT(element, "When there are 5 input channels, sensing-model must be 2. Resetting sensing-model to 2.");
- element->sensing_model = 2;
- }
- } else
- g_assert_not_reached();
+ g_assert(channels_out == 2 * channels_in);
/*
* record stream parameters
@@ -568,60 +587,6 @@ static gboolean set_caps(GstBaseTransform *trans, GstCaps *incaps, GstCaps *outc
}
-/*
- * transform_size{}
- */
-
-
-static gboolean transform_size(GstBaseTransform *trans, GstPadDirection direction, GstCaps *caps, gsize size, GstCaps *othercaps, gsize *othersize) {
-
- GSTLALSensingTDCFs *element = GSTLAL_SENSINGTDCFS(trans);
-
- /*
- * The data types of inputs and outputs are the same, but the number of channels differs.
- * For N output channels, there are N(N+1) input channels.
- */
-
- switch(direction) {
- case GST_PAD_SRC:
- /*
- * We know the size of the output buffer and want to compute the size of the input buffer.
- * The size of the output buffer should be a multiple of unit_size_out.
- */
-
- if(G_UNLIKELY(size % element->unit_size_out)) {
- GST_DEBUG_OBJECT(element, "buffer size %" G_GSIZE_FORMAT " is not a multiple of %" G_GSIZE_FORMAT, size, (gsize) element->unit_size_out);
- return FALSE;
- }
-
- *othersize = size * element->unit_size_in / element->unit_size_out;
-
- break;
-
- case GST_PAD_SINK:
- /*
- * We know the size of the input buffer and want to compute the size of the output buffer.
- * The size of the input buffer should be a multiple of unit_size_in.
- */
-
- if(G_UNLIKELY(size % element->unit_size_in)) {
- GST_ERROR_OBJECT(element, "buffer size %" G_GSIZE_FORMAT " is not a multiple of %" G_GSIZE_FORMAT, size, (gsize) element->unit_size_in);
- return FALSE;
- }
-
- *othersize = size * element->unit_size_out / element->unit_size_in;
-
- break;
-
- case GST_PAD_UNKNOWN:
- GST_ELEMENT_ERROR(trans, CORE, NEGOTIATION, (NULL), ("invalid direction GST_PAD_UNKNOWN"));
- return FALSE;
- }
-
- return TRUE;
-}
-
-
/*
* start()
*/
@@ -695,44 +660,43 @@ static GstFlowReturn transform(GstBaseTransform *trans, GstBuffer *inbuf, GstBuf
/*
* kappa_C is the solution of a quartic equation, so there are 4 solutions
- * for kappa_C, but only one is correct. For each solution, a different
- * value for fcc, fs, and Q will be computed, and the correct solution
- * will be determined at the end.
+ * for kappa_C. For each solution, a different value for fcc, fs, and Q
+ * will be computed, and the best solution will be determined at the end.
*/
kappa_C1 = kappa_C2 = kappa_C3 = kappa_C4 = f_cc1 = f_cc2 = f_cc3 = f_cc4 = 0.0;
switch(element->sensing_model) {
case 0: ;
f_s_squared1 = f_s_squared2 = f_s_squared3 = f_s_squared4 = f_s_over_Q1 = f_s_over_Q2 = f_s_over_Q3 = f_s_over_Q4 = 0.0;
for(i = 0; i < samples; i++) {
- kappa_C_0_float(crealf(indata[4 * i]), cimagf(indata[4 * i]), crealf(indata[4 * i + 1]), cimagf(indata[4 * i + 1]), crealf(indata[4 * i + 2]), cimagf(indata[4 * i + 2]), crealf(indata[4 * i + 3]), cimagf(indata[4 * i + 3]), f1, f2, &kappa_C1, &kappa_C2, &kappa_C3, &kappa_C4);
+ kappa_C_0_float(crealf(indata[2 * i]), cimagf(indata[2 * i]), crealf(indata[2 * i + 1]), cimagf(indata[2 * i + 1]), f1, f2, &kappa_C1, &kappa_C2, &kappa_C3, &kappa_C4);
/* Only compute f_cc, f_s, and Q if we have to. */
if(kappa_C1 != 0.0) {
- f_cc1 = f_cc_0_float(cimagf(indata[4 * i]), cimagf(indata[4 * i + 1]), cimagf(indata[4 * i + 2]), cimagf(indata[4 * i + 3]), f1, f2, kappa_C1);
- f_s_squared1 = f_s_squared_0_float(crealf(indata[4 * i]), crealf(indata[4 * i + 1]), crealf(indata[4 * i + 2]), crealf(indata[4 * i + 3]), f1, f2, kappa_C1);
- f_s_over_Q1 = f_s_over_Q_0_float(crealf(indata[4 * i]), crealf(indata[4 * i + 2]), f1, kappa_C1, f_cc1, f_s_squared1);
+ f_cc1 = f_cc_0_float(cimagf(indata[2 * i]), cimagf(indata[2 * i + 1]), f1, f2, kappa_C1);
+ f_s_squared1 = f_s_squared_0_float(crealf(indata[2 * i]), crealf(indata[2 * i + 1]), f1, f2, kappa_C1);
+ f_s_over_Q1 = f_s_over_Q_0_float(crealf(indata[2 * i]), f1, kappa_C1, f_cc1, f_s_squared1);
}
if(kappa_C2 != 0.0) {
- f_cc2 = f_cc_0_float(cimagf(indata[4 * i]), cimagf(indata[4 * i + 1]), cimagf(indata[4 * i + 2]), cimagf(indata[4 * i + 3]), f1, f2, kappa_C2);
- f_s_squared2 = f_s_squared_0_float(crealf(indata[4 * i]), crealf(indata[4 * i + 1]), crealf(indata[4 * i + 2]), crealf(indata[4 * i + 3]), f1, f2, kappa_C2);
- f_s_over_Q2 = f_s_over_Q_0_float(crealf(indata[4 * i]), crealf(indata[4 * i + 2]), f1, kappa_C2, f_cc2, f_s_squared2);
+ f_cc2 = f_cc_0_float(cimagf(indata[2 * i]), cimagf(indata[2 * i + 1]), f1, f2, kappa_C2);
+ f_s_squared2 = f_s_squared_0_float(crealf(indata[2 * i]), crealf(indata[2 * i + 1]), f1, f2, kappa_C2);
+ f_s_over_Q2 = f_s_over_Q_0_float(crealf(indata[2 * i]), f1, kappa_C2, f_cc2, f_s_squared2);
}
if(kappa_C3 != 0.0) {
- f_cc3 = f_cc_0_float(cimagf(indata[4 * i]), cimagf(indata[4 * i + 1]), cimagf(indata[4 * i + 2]), cimagf(indata[4 * i + 3]), f1, f2, kappa_C3);
- f_s_squared3 = f_s_squared_0_float(crealf(indata[4 * i]), crealf(indata[4 * i + 1]), crealf(indata[4 * i + 2]), crealf(indata[4 * i + 3]), f1, f2, kappa_C3);
- f_s_over_Q3 = f_s_over_Q_0_float(crealf(indata[4 * i]), crealf(indata[4 * i + 2]), f1, kappa_C3, f_cc3, f_s_squared3);
+ f_cc3 = f_cc_0_float(cimagf(indata[2 * i]), cimagf(indata[2 * i + 1]), f1, f2, kappa_C3);
+ f_s_squared3 = f_s_squared_0_float(crealf(indata[2 * i]), crealf(indata[2 * i + 1]), f1, f2, kappa_C3);
+ f_s_over_Q3 = f_s_over_Q_0_float(crealf(indata[2 * i]), f1, kappa_C3, f_cc3, f_s_squared3);
}
if(kappa_C4 != 0.0) {
- f_cc4 = f_cc_0_float(cimagf(indata[4 * i]), cimagf(indata[4 * i + 1]), cimagf(indata[4 * i + 2]), cimagf(indata[4 * i + 3]), f1, f2, kappa_C4);
- f_s_squared4 = f_s_squared_0_float(crealf(indata[4 * i]), crealf(indata[4 * i + 1]), crealf(indata[4 * i + 2]), crealf(indata[4 * i + 3]), f1, f2, kappa_C4);
- f_s_over_Q4 = f_s_over_Q_0_float(crealf(indata[4 * i]), crealf(indata[4 * i + 2]), f1, kappa_C4, f_cc4, f_s_squared4);
+ f_cc4 = f_cc_0_float(cimagf(indata[2 * i]), cimagf(indata[2 * i + 1]), f1, f2, kappa_C4);
+ f_s_squared4 = f_s_squared_0_float(crealf(indata[2 * i]), crealf(indata[2 * i + 1]), f1, f2, kappa_C4);
+ f_s_over_Q4 = f_s_over_Q_0_float(crealf(indata[2 * i]), f1, kappa_C4, f_cc4, f_s_squared4);
}
- /* Determine which solution is correct. */
- best_solution = find_best_solution_float(kappa_C1, kappa_C2, kappa_C3, kappa_C4, f_cc1, f_cc2, f_cc3, f_cc4, f_s_squared1, f_s_squared2, f_s_squared3, f_s_squared4, f_s_over_Q1, f_s_over_Q2, f_s_over_Q3, f_s_over_Q4, f1, f2, indata[4 * i + 2] - indata[4 * i], indata[4 * i + 3] - indata[4 * i + 1]);
+ /* Determine which solution is our favorite. */
+ best_solution = find_best_solution_float(kappa_C1, kappa_C2, kappa_C3, kappa_C4, f_cc1, f_cc2, f_cc3, f_cc4, f_s_squared1, f_s_squared2, f_s_squared3, f_s_squared4, f_s_over_Q1, f_s_over_Q2, f_s_over_Q3, f_s_over_Q4, f1, f2, indata[2 * i], indata[2 * i + 1]);
switch(best_solution) {
case 1:
@@ -789,44 +753,43 @@ static GstFlowReturn transform(GstBaseTransform *trans, GstBuffer *inbuf, GstBuf
/*
* kappa_C is the solution of a quartic equation, so there are 4 solutions
- * for kappa_C, but only one is correct. For each solution, a different
- * value for fcc, fs, and Q will be computed, and the correct solution
- * will be determined at the end.
+ * for kappa_C. For each solution, a different value for fcc, fs, and Q
+ * will be computed, and the best solution will be determined at the end.
*/
kappa_C1 = kappa_C2 = kappa_C3 = kappa_C4 = f_cc1 = f_cc2 = f_cc3 = f_cc4 = 0.0;
switch(element->sensing_model) {
case 0: ;
f_s_squared1 = f_s_squared2 = f_s_squared3 = f_s_squared4 = f_s_over_Q1 = f_s_over_Q2 = f_s_over_Q3 = f_s_over_Q4 = 0.0;
for(i = 0; i < samples; i++) {
- kappa_C_0_double(creal(indata[4 * i]), cimag(indata[4 * i]), creal(indata[4 * i + 1]), cimag(indata[4 * i + 1]), creal(indata[4 * i + 2]), cimag(indata[4 * i + 2]), creal(indata[4 * i + 3]), cimag(indata[4 * i + 3]), f1, f2, &kappa_C1, &kappa_C2, &kappa_C3, &kappa_C4);
+ kappa_C_0_double(crealf(indata[2 * i]), cimagf(indata[2 * i]), crealf(indata[2 * i + 1]), cimagf(indata[2 * i + 1]), f1, f2, &kappa_C1, &kappa_C2, &kappa_C3, &kappa_C4);
/* Only compute f_cc, f_s, and Q if we have to. */
if(kappa_C1 != 0.0) {
- f_cc1 = f_cc_0_double(cimag(indata[4 * i]), cimag(indata[4 * i + 1]), cimag(indata[4 * i + 2]), cimag(indata[4 * i + 3]), f1, f2, kappa_C1);
- f_s_squared1 = f_s_squared_0_double(creal(indata[4 * i]), creal(indata[4 * i + 1]), creal(indata[4 * i + 2]), creal(indata[4 * i + 3]), f1, f2, kappa_C1);
- f_s_over_Q1 = f_s_over_Q_0_double(creal(indata[4 * i]), creal(indata[4 * i + 2]), f1, kappa_C1, f_cc1, f_s_squared1);
+ f_cc1 = f_cc_0_double(cimagf(indata[2 * i]), cimagf(indata[2 * i + 1]), f1, f2, kappa_C1);
+ f_s_squared1 = f_s_squared_0_double(crealf(indata[2 * i]), crealf(indata[2 * i + 1]), f1, f2, kappa_C1);
+ f_s_over_Q1 = f_s_over_Q_0_double(crealf(indata[2 * i]), f1, kappa_C1, f_cc1, f_s_squared1);
}
if(kappa_C2 != 0.0) {
- f_cc2 = f_cc_0_double(cimag(indata[4 * i]), cimag(indata[4 * i + 1]), cimag(indata[4 * i + 2]), cimag(indata[4 * i + 3]), f1, f2, kappa_C2);
- f_s_squared2 = f_s_squared_0_double(creal(indata[4 * i]), creal(indata[4 * i + 1]), creal(indata[4 * i + 2]), creal(indata[4 * i + 3]), f1, f2, kappa_C2);
- f_s_over_Q2 = f_s_over_Q_0_double(creal(indata[4 * i]), creal(indata[4 * i + 2]), f1, kappa_C2, f_cc2, f_s_squared2);
+ f_cc2 = f_cc_0_double(cimagf(indata[2 * i]), cimagf(indata[2 * i + 1]), f1, f2, kappa_C2);
+ f_s_squared2 = f_s_squared_0_double(crealf(indata[2 * i]), crealf(indata[2 * i + 1]), f1, f2, kappa_C2);
+ f_s_over_Q2 = f_s_over_Q_0_double(crealf(indata[2 * i]), f1, kappa_C2, f_cc2, f_s_squared2);
}
if(kappa_C3 != 0.0) {
- f_cc3 = f_cc_0_double(cimag(indata[4 * i]), cimag(indata[4 * i + 1]), cimag(indata[4 * i + 2]), cimag(indata[4 * i + 3]), f1, f2, kappa_C3);
- f_s_squared3 = f_s_squared_0_double(creal(indata[4 * i]), creal(indata[4 * i + 1]), creal(indata[4 * i + 2]), creal(indata[4 * i + 3]), f1, f2, kappa_C3);
- f_s_over_Q3 = f_s_over_Q_0_double(creal(indata[4 * i]), creal(indata[4 * i + 2]), f1, kappa_C3, f_cc3, f_s_squared3);
+ f_cc3 = f_cc_0_double(cimagf(indata[2 * i]), cimagf(indata[2 * i + 1]), f1, f2, kappa_C3);
+ f_s_squared3 = f_s_squared_0_double(crealf(indata[2 * i]), crealf(indata[2 * i + 1]), f1, f2, kappa_C3);
+ f_s_over_Q3 = f_s_over_Q_0_double(crealf(indata[2 * i]), f1, kappa_C3, f_cc3, f_s_squared3);
}
if(kappa_C4 != 0.0) {
- f_cc4 = f_cc_0_double(cimag(indata[4 * i]), cimag(indata[4 * i + 1]), cimag(indata[4 * i + 2]), cimag(indata[4 * i + 3]), f1, f2, kappa_C4);
- f_s_squared4 = f_s_squared_0_double(creal(indata[4 * i]), creal(indata[4 * i + 1]), creal(indata[4 * i + 2]), creal(indata[4 * i + 3]), f1, f2, kappa_C4);
- f_s_over_Q4 = f_s_over_Q_0_double(creal(indata[4 * i]), creal(indata[4 * i + 2]), f1, kappa_C4, f_cc4, f_s_squared4);
+ f_cc4 = f_cc_0_double(cimagf(indata[2 * i]), cimagf(indata[2 * i + 1]), f1, f2, kappa_C4);
+ f_s_squared4 = f_s_squared_0_double(crealf(indata[2 * i]), crealf(indata[2 * i + 1]), f1, f2, kappa_C4);
+ f_s_over_Q4 = f_s_over_Q_0_double(crealf(indata[2 * i]), f1, kappa_C4, f_cc4, f_s_squared4);
}
- /* Determine which solution is correct. */
- best_solution = find_best_solution_double(kappa_C1, kappa_C2, kappa_C3, kappa_C4, f_cc1, f_cc2, f_cc3, f_cc4, f_s_squared1, f_s_squared2, f_s_squared3, f_s_squared4, f_s_over_Q1, f_s_over_Q2, f_s_over_Q3, f_s_over_Q4, f1, f2, indata[4 * i + 2] - indata[4 * i], indata[4 * i + 3] - indata[4 * i + 1]);
+ /* Determine which solution is our favorite. */
+ best_solution = find_best_solution_double(kappa_C1, kappa_C2, kappa_C3, kappa_C4, f_cc1, f_cc2, f_cc3, f_cc4, f_s_squared1, f_s_squared2, f_s_squared3, f_s_squared4, f_s_over_Q1, f_s_over_Q2, f_s_over_Q3, f_s_over_Q4, f1, f2, indata[2 * i], indata[2 * i + 1]);
switch(best_solution) {
case 1:
@@ -1022,8 +985,8 @@ static void gstlal_sensingtdcfs_class_init(GSTLALSensingTDCFsClass *klass)
"Filter/Audio",
"Solves for the time-dependent correction factors of the sensing function using\n\t\t\t "
"the solution described in LIGO DCC document P1900052, Section 5.2.6. It takes\n\t\t\t "
- "the complex inputs Gres^1, Gres^2, Y^1, and Y^2 (and Y^3 if sensing-model is 2),\n\t\t\t "
- "in that order. The outputs are kappa_C, f_cc, f_s, and Q, in that order.\n\t\t\t "
+ "the complex inputs Cinv_tdep^i, for i = 1, 2, ... N where 2N is the number of\n\t\t\t "
+ "TDCFs solved for. The outputs are kappa_C, f_cc, f_s, and Q, in that order.\n\t\t\t "
"Currently, sensing-model=0 is the only sensing model that has been developed.",
"Aaron Viets <aaron.viets@ligo.org>"
);
@@ -1034,7 +997,6 @@ static void gstlal_sensingtdcfs_class_init(GSTLALSensingTDCFsClass *klass)
transform_class->get_unit_size = GST_DEBUG_FUNCPTR(get_unit_size);
transform_class->set_caps = GST_DEBUG_FUNCPTR(set_caps);
transform_class->transform_caps = GST_DEBUG_FUNCPTR(transform_caps);
- transform_class->transform_size = GST_DEBUG_FUNCPTR(transform_size);
transform_class->start = GST_DEBUG_FUNCPTR(start);
transform_class->transform = GST_DEBUG_FUNCPTR(transform);
@@ -1052,19 +1014,19 @@ static void gstlal_sensingtdcfs_class_init(GSTLALSensingTDCFsClass *klass)
" kappa_C\t\t f^2\n\t\t\t"
" -------------- * ------------------------- * C_res(f)\n\t\t\t"
" 1 + i f / f_cc f^2 + f_s^2 - i f f_s / Q\n\n\t\t\t"
- "Complex inputs: Gres1, Gres2, Y1, Y2. See P1900052, Sec. 5.2.6.\n\t\t\t"
+ "Complex inputs: Cinv_tdep1, Cinv_tdep2. See P1900052, Sec. 5.2.6.\n\t\t\t"
"Real outputs: kappa_C, f_cc, f_s^2, f_s / Q\n\n\t\t\t"
"sensing-model=1:\n\t\t\t"
" kappa_C\n\t\t\t"
" -------------- * (\?\?\?)... not yet modeled\n\t\t\t"
" 1 + i f / f_cc\n\n\t\t\t"
- "Complex inputs: Gres1, Gres2, Y1, Y2.\n\t\t\t"
+ "Complex inputs: Cinv_tdep1, Cinv_tdep2.\n\t\t\t"
"Real outputs: kappa_C, f_cc, \?\?, \?\?\n\n\t\t\t"
"sensing-model=2:\n\t\t\t"
" kappa_C\t\t f^2\n\t\t\t"
" -------------- * ------------------------- * (\?\?\?)... not yet modeled.\n\t\t\t"
" 1 + i f / f_cc f^2 + f_s^2 - i f f_s / Q\n\n\t\t\t"
- "Complex inputs: Gres1, Gres2, Y1, Y2, Y3\?\n\t\t\t"
+ "Complex inputs: Cinv_tdep1, Cinv_tdep2, Cinv_tdep3.\?\n\t\t\t"
"Real outputs: kappa_C, f_cc, f_s^2, f_s / Q, \?\?, \?\?",
0, 2, 0,
G_PARAM_READWRITE | G_PARAM_STATIC_STRINGS | G_PARAM_CONSTRUCT
@@ -1104,16 +1066,16 @@ static void gstlal_sensingtdcfs_class_init(GSTLALSensingTDCFsClass *klass)
)
);
g_object_class_install_property(
- gobject_class,
- ARG_CURRENT_KC,
- g_param_spec_double(
- "current-kc",
- "Current kc",
- "Current value of the variable optical gain of the sensing function.",
- -G_MAXDOUBLE, G_MAXDOUBLE, 1.0,
- G_PARAM_READWRITE | G_PARAM_STATIC_STRINGS | G_PARAM_CONSTRUCT
- )
- );
+ gobject_class,
+ ARG_CURRENT_KC,
+ g_param_spec_double(
+ "current-kc",
+ "Current kc",
+ "Current value of the variable optical gain of the sensing function.",
+ -G_MAXDOUBLE, G_MAXDOUBLE, 1.0,
+ G_PARAM_READWRITE | G_PARAM_STATIC_STRINGS | G_PARAM_CONSTRUCT
+ )
+ );
g_object_class_install_property(
gobject_class,
ARG_CURRENT_FCC,
diff --git a/gstlal-calibration/python/calibration_parts.py b/gstlal-calibration/python/calibration_parts.py
index 412abf0b94823e68bbdc39c326f00a91cb584b70..5f0d612ba1c9b8ad77b8ebc267474f2b4a6c3787 100644
--- a/gstlal-calibration/python/calibration_parts.py
+++ b/gstlal-calibration/python/calibration_parts.py
@@ -1118,7 +1118,7 @@ def compute_exact_kappas_from_filters_file(pipeline, X, freqs, EPICS, rate, defa
CAX = []
CAXreal = []
CAXimag = []
- Gres = []
+ Cinv_tdep = []
kappas = []
for i in range(num_lines):
@@ -1199,19 +1199,19 @@ def compute_exact_kappas_from_filters_file(pipeline, X, freqs, EPICS, rate, defa
kappas[i] = mkmultiplier(pipeline, [kappas[i], mkpow(pipeline, kappas[i - len(kappas) // 2], exponent = -1.0)])
kappas[i] = pipeparts.mktee(pipeline, kappas[i])
- # Next, compute kappa_C. This is going to take some work...
+ # Next, compute the sensing function time dependence.
# Start by computing G_res at each frequency, defined in Eq. 5.2.30
for n in range(2):
- Gres_components = []
+ Cinv_tdep_components = [Y[n]]
for j in range(num_stages):
- kappajGresjatn = pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, pipeparts.mkcapsfilter(pipeline, kappas[j], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), matrix = [[EPICS[2 * (n * (1 + num_stages) + 1 + j)], EPICS[1 + 2 * (n * (1 + num_stages) + 1 + j)]]]))
+ minuskappajGresjatn = pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, pipeparts.mkcapsfilter(pipeline, kappas[j], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), matrix = [[-EPICS[2 * (n * (1 + num_stages) + 1 + j)], -EPICS[1 + 2 * (n * (1 + num_stages) + 1 + j)]]]))
i_omega_tau = pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, pipeparts.mkcapsfilter(pipeline, kappas[num_stages + j], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), matrix = [[0, 2.0 * numpy.pi * freqs[n]]]))
i_omega_tau = mkcapsfiltersetter(pipeline, i_omega_tau, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)
phase = pipeparts.mkgeneric(pipeline, i_omega_tau, "cexp")
- Gres_components.append(mkmultiplier(pipeline, list_srcs(pipeline, kappajGresjatn, phase)))
- Gres.append(mkadder(pipeline, Gres_components))
+ Cinv_tdep_components.append(mkmultiplier(pipeline, list_srcs(pipeline, minuskappajGresjatn, phase)))
+ Cinv_tdep.append(mkadder(pipeline, Cinv_tdep_components))
- sensing_inputs = mkinterleave(pipeline, Gres + Y, complex_data = True)
+ sensing_inputs = mkinterleave(pipeline, Cinv_tdep, complex_data = True)
sensing_outputs = pipeparts.mkgeneric(pipeline, sensing_inputs, "lal_sensingtdcfs", sensing_model = 0, freq1 = freqs[0], freq2 = freqs[1], current_fcc = default_fcc, current_fs_squared = default_fs_squared, current_fs_over_Q = default_fs_over_Q)
sensing_outputs = list(mkdeinterleave(pipeline, sensing_outputs, 4))
diff --git a/gstlal-calibration/tests/check_calibration/Makefile b/gstlal-calibration/tests/check_calibration/Makefile
index 49c95b4a52fb7ce73512cf779ada7057e1796180..6018213f60b0f01be9f61000eed69cdee553c180 100644
--- a/gstlal-calibration/tests/check_calibration/Makefile
+++ b/gstlal-calibration/tests/check_calibration/Makefile
@@ -4,16 +4,14 @@
#################################
# which interferometer (H or L)
-IFO = L
+IFO = H
# determines where to look for filters files (e.g., O1, O2, O3, ER10, ER13, ER14, PreER10, PreER13, PreER14)
OBSRUN = O3
-START = $(shell echo 1251014418 | bc)
+START = $(shell echo 1238288418 | bc)
# 1238288418
-# 1269124389
-END = $(shell echo 1251061218 | bc)
-# 1238295618
-# 1269128350
+END = $(shell echo 1238338818 | bc)
+# 1238338818
SHMRUNTIME = 600
# How much time does the calibration need to settle at the start and end?
PLOT_WARMUP_TIME = 256
@@ -29,8 +27,8 @@ FCC_CORR_CONFIGS = Filters/O3/GDSFilters/H1DCS_test_1256655618_v2_FCC_CORR.ini
ALL_CORR_CONFIGS = Filters/O3/GDSFilters/H1DCS_test_1256655618_v2_ALL_CORR.ini
FCCFS_CORR_CONFIGS = Filters/O2/GDSFilters/H1DCS_FccFsCorrections_Cleaning.ini
NO_CORR_CONFIGS = Filters/O3/GDSFilters/H1DCS_test_1256655618_v2_NO_CORR.ini
-DCSEXACTKAPPASCONFIGS = Filters/O3/GDSFilters/L1DCS_test_1249927218_exactKappas.ini
-DCSAPPROXKAPPASCONFIGS = Filters/O3/GDSFilters/L1DCS_test_1249927218_approxKappas.ini
+DCSEXACTKAPPASCONFIGS = Filters/O3/GDSFilters/H1DCS_test_1237831461_exactKappas.ini
+DCSAPPROXKAPPASCONFIGS = Filters/O3/GDSFilters/H1DCS_test_1237831461_approxKappas.ini
# H1DCS_test_1237831461_exactKappas.ini
# H1DCS_C01_1256655618_v2_test.ini
DCSLINESCONFIGS = ../../config_files/O2/H1/tests/H1DCS_AllCorrections_CleaningLines.ini
@@ -44,7 +42,7 @@ GDSOLDCONFIGS = Filters/ER14/GDSFilters/L1GDS_1235491416_old.ini
GDSBETTERCONFIGS = Filters/ER14/GDSFilters/L1GDS_1235491416_better.ini
GDSBESTCONFIGS = Filters/ER14/GDSFilters/L1GDS_1235491416_best.ini
-all: DCSEXACTKAPPAS_pcal2darm_plots DCSEXACTKAPPAS_act2darm_plots exactkappastimeseries
+all: exactkappastimeseries DCSEXACTKAPPAS_act2darm_plots DCSEXACTKAPPAS_pcal2darm_plots
#pcal_DCS_transfer_functions
#TDCFs_pcal2darm_plots
#actuation_timing_plot
@@ -187,7 +185,7 @@ DCSEXACTKAPPAS_pcal2darm_plots: $(IFO)1_easy_raw_frames.cache $(IFO)1_hoft_DCS_A
python3 pcal2darm_timeseries.py --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --ifo $(IFO)1 --raw-frame-cache $(IFO)1_easy_raw_frames.cache --gstlal-frame-cache-list '$(IFO)1_hoft_DCS_APPROXKAPPAS_frames.cache,$(IFO)1_hoft_DCS_EXACTKAPPAS_frames.cache' --config-file-list '$(DCSAPPROXKAPPASCONFIGS),$(DCSEXACTKAPPASCONFIGS)' --pcal-channel-name CAL-PCALY_RX_PD_OUT_DQ --gstlal-channel-list 'DCS-CALIB_STRAINAPPROXKAPPAS,DCS-CALIB_STRAINEXACTKAPPAS' --labels '{\rm Approx},{\rm Exact}' --magnitude-ranges '0.95,1.05;0.95,1.05;0.8,1.2' --phase-ranges '-3.0,3.0;-3.0,3.0;-15.0,15.0' --latex-labels --pcal-time-advance 0.00006103515625 --show-stats
DCSEXACTKAPPAS_act2darm_plots: $(IFO)1_easy_raw_frames.cache $(IFO)1_hoft_DCS_APPROXKAPPAS_frames.cache $(IFO)1_hoft_DCS_EXACTKAPPAS_frames.cache
- python3 act2darm_timeseries.py --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --ifo $(IFO)1 --raw-frame-cache $(IFO)1_easy_raw_frames.cache --gstlal-frame-cache-list '$(IFO)1_hoft_DCS_APPROXKAPPAS_frames.cache,$(IFO)1_hoft_DCS_EXACTKAPPAS_frames.cache' --config-file-list '$(DCSAPPROXKAPPASCONFIGS),$(DCSEXACTKAPPASCONFIGS)' --gstlal-channel-list 'DCS-CALIB_STRAINAPPROXKAPPAS,DCS-CALIB_STRAINEXACTKAPPAS' --act-channel-list "SUS-ETMY_L3_CAL_LINE_OUT_DQ,SUS-ETMX_L2_CAL_LINE_OUT_DQ,SUS-ETMX_L1_CAL_LINE_OUT_DQ" --labels '{\rm Approx},{\rm Exact}' --magnitude-ranges '0.95,1.05;0.95,1.05;0.95,1.05' --phase-ranges '-3.0,3.0;-3.0,3.0;-3.0,3.0' --latex-labels --act-time-advance 0.00006103515625 --show-stats
+ python3 act2darm_timeseries.py --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --ifo $(IFO)1 --raw-frame-cache $(IFO)1_easy_raw_frames.cache --gstlal-frame-cache-list '$(IFO)1_hoft_DCS_APPROXKAPPAS_frames.cache,$(IFO)1_hoft_DCS_EXACTKAPPAS_frames.cache' --config-file-list '$(DCSAPPROXKAPPASCONFIGS),$(DCSEXACTKAPPASCONFIGS)' --gstlal-channel-list 'DCS-CALIB_STRAINAPPROXKAPPAS,DCS-CALIB_STRAINEXACTKAPPAS' --act-channel-list "SUS-ETMX_L3_CAL_LINE_OUT_DQ,SUS-ETMX_L2_CAL_LINE_OUT_DQ,SUS-ETMX_L1_CAL_LINE_OUT_DQ" --labels '{\rm Approx},{\rm Exact}' --magnitude-ranges '0.95,1.05;0.8,1.2;0.8,1.2' --phase-ranges '-3.0,3.0;-15.0,15.0;-15.0,15.0' --latex-labels --act-time-advance 0.00006103515625 --show-stats
TDCFs_pcal2darm_plots: $(IFO)1_easy_raw_frames.cache $(IFO)1_SCALAR_CORR_frames.cache $(IFO)1_FCC_CORR_frames.cache $(IFO)1_FCCFS_CORR_frames.cache $(IFO)1_ALL_CORR_frames.cache
python3 pcal2darm_timeseries.py --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --ifo $(IFO)1 --raw-frame-cache $(IFO)1_easy_raw_frames.cache --gstlal-frame-cache-list '$(IFO)1_SCALAR_CORR_frames.cache,$(IFO)1_FCC_CORR_frames.cache,$(IFO)1_FCCFS_CORR_frames.cache,$(IFO)1_ALL_CORR_frames.cache' --config-file '$(ALL_CORR_CONFIGS)' --pcal-channel-name CAL-PCALY_RX_PD_OUT_DQ --gstlal-channel-list 'DCS-CALIB_STRAIN,DCS-CALIB_STRAIN,DCS-CALIB_STRAIN,DCS-CALIB_STRAIN' --labels '{\rm Scalars},+f_{\rm cc},+f_{\rm s}+Q,+\tau_i' --latex-labels --pcal-time-advance 0.00006103515625 --show-stats
@@ -255,10 +253,12 @@ kappastimeseries_GDS: $(IFO)1_hoft_GDS_frames.cache $(IFO)1_easy_raw_frames.cach
exactkappastimeseries: $(IFO)1_hoft_DCS_APPROXKAPPAS_frames.cache $(IFO)1_hoft_DCS_EXACTKAPPAS_frames.cache
python3 frame_manipulator.py --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --ifo $(IFO)1 --frame-cache $(IFO)1_hoft_DCS_EXACTKAPPAS_frames.cache --output-path txt --channel-list 'DCS-CALIB_KAPPA_TST_REALEXACTKAPPAS,DCS-CALIB_KAPPA_TST_IMAGINARYEXACTKAPPAS,DCS-CALIB_KAPPA_PUM_REALEXACTKAPPAS,DCS-CALIB_KAPPA_PUM_IMAGINARYEXACTKAPPAS,DCS-CALIB_KAPPA_UIM_REALEXACTKAPPAS,DCS-CALIB_KAPPA_UIM_IMAGINARYEXACTKAPPAS,DCS-CALIB_KAPPA_CEXACTKAPPAS,DCS-CALIB_F_CCEXACTKAPPAS,DCS-CALIB_F_S_SQUAREDEXACTKAPPAS,DCS-CALIB_SRC_Q_INVERSEEXACTKAPPAS'
+ python3 compute_fs_over_Q.py --fs-squared-txt $(IFO)1-DCS-CALIB_F_S_SQUAREDEXACTKAPPAS.txt --Qinv-txt $(IFO)1-DCS-CALIB_SRC_Q_INVERSEEXACTKAPPAS.txt --filename $(IFO)1-DCS-CALIB_F_S_OVER_QEXACTKAPPAS.txt
python3 frame_manipulator.py --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --ifo $(IFO)1 --frame-cache $(IFO)1_hoft_DCS_APPROXKAPPAS_frames.cache --output-path txt --channel-list 'DCS-CALIB_KAPPA_TST_REALAPPROXKAPPAS,DCS-CALIB_KAPPA_TST_IMAGINARYAPPROXKAPPAS,DCS-CALIB_KAPPA_PUM_REALAPPROXKAPPAS,DCS-CALIB_KAPPA_PUM_IMAGINARYAPPROXKAPPAS,DCS-CALIB_KAPPA_UIM_REALAPPROXKAPPAS,DCS-CALIB_KAPPA_UIM_IMAGINARYAPPROXKAPPAS,DCS-CALIB_KAPPA_CAPPROXKAPPAS,DCS-CALIB_F_CCAPPROXKAPPAS,DCS-CALIB_F_S_SQUAREDAPPROXKAPPAS,DCS-CALIB_SRC_Q_INVERSEAPPROXKAPPAS'
+ python3 compute_fs_over_Q.py --fs-squared-txt $(IFO)1-DCS-CALIB_F_S_SQUAREDAPPROXKAPPAS.txt --Qinv-txt $(IFO)1-DCS-CALIB_SRC_Q_INVERSEAPPROXKAPPAS.txt --filename $(IFO)1-DCS-CALIB_F_S_OVER_QAPPROXKAPPAS.txt
python3 plot_kappas_from_txt.py --txt-list '$(IFO)1-DCS-CALIB_KAPPA_TST_REALAPPROXKAPPAS.txt,$(IFO)1-DCS-CALIB_KAPPA_PUM_REALAPPROXKAPPAS.txt,$(IFO)1-DCS-CALIB_KAPPA_UIM_REALAPPROXKAPPAS.txt;$(IFO)1-DCS-CALIB_KAPPA_TST_REALEXACTKAPPAS.txt,$(IFO)1-DCS-CALIB_KAPPA_PUM_REALEXACTKAPPAS.txt,$(IFO)1-DCS-CALIB_KAPPA_UIM_REALEXACTKAPPAS.txt:$(IFO)1-DCS-CALIB_KAPPA_TST_IMAGINARYAPPROXKAPPAS.txt,$(IFO)1-DCS-CALIB_KAPPA_PUM_IMAGINARYAPPROXKAPPAS.txt,$(IFO)1-DCS-CALIB_KAPPA_UIM_IMAGINARYAPPROXKAPPAS.txt;$(IFO)1-DCS-CALIB_KAPPA_TST_IMAGINARYEXACTKAPPAS.txt,$(IFO)1-DCS-CALIB_KAPPA_PUM_IMAGINARYEXACTKAPPAS.txt,$(IFO)1-DCS-CALIB_KAPPA_UIM_IMAGINARYEXACTKAPPAS.txt' --labels 'Approx;Exact' --filename actuation_TDCFs.png
python3 plot_kappas_from_txt.py --txt-list '$(IFO)1-DCS-CALIB_KAPPA_CAPPROXKAPPAS.txt;$(IFO)1-DCS-CALIB_KAPPA_CEXACTKAPPAS.txt:$(IFO)1-DCS-CALIB_F_CCAPPROXKAPPAS.txt;$(IFO)1-DCS-CALIB_F_CCEXACTKAPPAS.txt' --labels 'Approx;Exact' --filename sensing_TDCFs.png
- python3 plot_kappas_from_txt.py --txt-list '$(IFO)1-DCS-CALIB_F_S_SQUAREDAPPROXKAPPAS.txt;$(IFO)1-DCS-CALIB_F_S_SQUAREDEXACTKAPPAS.txt:$(IFO)1-DCS-CALIB_SRC_Q_INVERSEAPPROXKAPPAS.txt;$(IFO)1-DCS-CALIB_SRC_Q_INVERSEEXACTKAPPAS.txt' --labels 'Approx;Exact' --filename SRC_TDCFs.png
+ python3 plot_kappas_from_txt.py --txt-list '$(IFO)1-DCS-CALIB_F_S_SQUAREDAPPROXKAPPAS.txt;$(IFO)1-DCS-CALIB_F_S_SQUAREDEXACTKAPPAS.txt:$(IFO)1-DCS-CALIB_F_S_OVER_QAPPROXKAPPAS.txt;$(IFO)1-DCS-CALIB_F_S_OVER_QEXACTKAPPAS.txt' --labels 'Approx;Exact' --filename SRC_TDCFs.png
kappastimeseries: $(IFO)1_hoft_DCS_frames.cache $(IFO)1_C01_frames.cache
python3 frame_manipulator.py --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --ifo $(IFO)1 --frame-cache $(IFO)1_C01_frames.cache --output-path txt --channel-list 'DCS-CALIB_KAPPA_TST_REAL_C01,DCS-CALIB_KAPPA_TST_IMAGINARY_C01,DCS-CALIB_KAPPA_PUM_REAL_C01,DCS-CALIB_KAPPA_PUM_IMAGINARY_C01,DCS-CALIB_KAPPA_UIM_REAL_C01,DCS-CALIB_KAPPA_UIM_IMAGINARY_C01,DCS-CALIB_KAPPA_C_C01,DCS-CALIB_F_CC_C01,DCS-CALIB_F_S_SQUARED_C01,DCS-CALIB_SRC_Q_INVERSE_C01'