Commit 8520e575 authored by Aaron Viets's avatar Aaron Viets

pipeparts: allow additional properties in capsfilter (e.g., name)

calibration_parts.py:  lower resampler quality for lower latency and lower computational cost.  Also edits to exact kappas solution.
parent 21c9850a
Pipeline #73073 failed with stages
in 1 minute and 7 seconds
......@@ -43,7 +43,13 @@ def mkcomplexqueue(pipeline, head, length = 0, min_length = 0):
head = mkqueue(pipeline, head, length = length, min_length = min_length)
head = pipeparts.mktogglecomplex(pipeline, head)
return head
def mkcapsfiltersetter(pipeline, head, caps, **properties):
# Make a capsfilter followed by a capssetter
head = pipeparts.mkcapsfilter(pipeline, head, caps)
#head = pipeparts.mkcapssetter(pipeline, head, caps, replace = True, **properties)
return head
def mkinsertgap(pipeline, head, **properties):
if "bad_data_intervals" in properties:
# Make sure the array property bad-data-intervals is formatted correctly
......@@ -257,7 +263,7 @@ def demodulate(pipeline, head, freq, td, rate, filter_time, filter_latency, pref
freq_update[2].connect("notify::current-average", update_property_simple, head, "current_average", "prefactor_imag", 1)
elif freq_update is not None:
freq_update.connect("notify::current-average", update_property_simple, head, "current_average", "line_frequency", 1)
head = mkresample(pipeline, head, 5, filter_latency == 0.0, rate)
head = mkresample(pipeline, head, 4, filter_latency == 0.0, rate)
if filter_latency != 0:
# Remove the first several seconds of output, which depend on start time
head = pipeparts.mkgeneric(pipeline, head, "lal_insertgap", chop_length = 7000000000)
......@@ -275,7 +281,7 @@ def remove_harmonics(pipeline, signal, f0, num_harmonics, f0_var, filter_latency
mkqueue(pipeline, head).link(elem)
for i in range(1, num_harmonics + 1):
line = pipeparts.mkgeneric(pipeline, head, "lal_demodulate", line_frequency = i * f0)
line = mkresample(pipeline, line, 5, filter_latency == 0, compute_rate)
line = mkresample(pipeline, line, 4, filter_latency == 0, compute_rate)
line_in_witness = lowpass(pipeline, line_in_witness, compute_rate, length = filter_param / (f0_var * i), fcut = 0, filter_latency = filter_latency)
line = mkresample(pipeline, line, 3, filter_latency == 0.0, rate_out)
line = pipeparts.mkgeneric(pipeline, line, "lal_demodulate", line_frequency = -1.0 * i * f0, prefactor_real = -2.0)
......@@ -583,7 +589,7 @@ def linear_phase_filter(pipeline, head, shift_samples, num_samples = 256, gain =
def compute_rms(pipeline, head, rate, average_time, f_min = None, f_max = None, filter_latency = 0.5, rate_out = 16, td = True):
# Find the root mean square amplitude of a signal between two frequencies
# Downsample to save computational cost
head = mkresample(pipeline, head, 5, filter_latency == 0.0, rate)
head = mkresample(pipeline, head, 4, filter_latency == 0.0, rate)
# Remove any frequency content we don't care about
if (f_min is not None) and (f_max is not None):
......@@ -987,13 +993,13 @@ def compute_exact_kappas_from_filters_file(pipeline, X, freqs, EPICS, rate):
if i < 2:
# Then it's a Pcal line
Y = pipeparts.mktee(pipeline, complex_audioamplify(pipeline, X[i], EPICS[2 * (1 + num_stages) * i], EPICS[2 * (1 + num_stages) * i + 1]))
Yreal.append(pipeparts.mktee(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkgeneric(pipeline, Y, "creal"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0" % rate, name = "capsfilter_Yreal_%d" % i)))
Yimag.append(pipeparts.mktee(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkgeneric(pipeline, Y, "cimag"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0" % rate, name = "capsfilter_Yimag_%d" % i)))
Yreal.append(pipeparts.mktee(pipeline, mkcapsfiltersetter(pipeline, pipeparts.mkgeneric(pipeline, Y, "creal"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate, name = "capsfilter_Yreal_%d" % i)))
Yimag.append(pipeparts.mktee(pipeline, mkcapsfiltersetter(pipeline, pipeparts.mkgeneric(pipeline, Y, "cimag"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate, name = "capsfilter_Yimag_%d" % i)))
else:
# It's an actuator line
CAX.append(pipeparts.mktee(pipeline, complex_audioamplify(pipeline, X[i], EPICS[2 * (1 + num_stages) * i], EPICS[2 * (1 + num_stages) * i + 1])))
CAXreal.append(pipeparts.mktee(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkgeneric(pipeline, CAX[-1], "creal"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0" % rate, name = "capsfilter_CAXreal_%d" % i)))
CAXimag.append(pipeparts.mktee(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkgeneric(pipeline, CAX[-1], "cimag"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0" % rate, name = "capsfilter_CAXimag_%d" % i)))
CAXreal.append(pipeparts.mktee(pipeline, mkcapsfiltersetter(pipeline, pipeparts.mkgeneric(pipeline, CAX[-1], "creal"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate, name = "capsfilter_CAXreal_%d" % i)))
CAXimag.append(pipeparts.mktee(pipeline, mkcapsfiltersetter(pipeline, pipeparts.mkgeneric(pipeline, CAX[-1], "cimag"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate, name = "capsfilter_CAXimag_%d" % i)))
# Let's start by computing the V's of Eqs. 5.2.78 and 5.2.79
for j in range(num_stages):
......@@ -1001,30 +1007,32 @@ def compute_exact_kappas_from_filters_file(pipeline, X, freqs, EPICS, rate):
factor2 = pow(freqs[2 + j], -2) - pow(freqs[1], -2)
factor3 = pow(freqs[0], 2) - pow(freqs[2 + j], 2)
factor4 = freqs[0] * (freqs[2 + j] - freqs[1])
Vj = mkadder(pipeline, list_srcs(pipeline, pipeparts.mkaudioamplify(pipeline, Yreal[1], factor1), pipeparts.mkaudioamplify(pipeline, Yreal[0], factor2)))
Vjplus3 = mkadder(pipeline, list_srcs(pipeline, pipeparts.mkaudioamplify(pipeline, Yimag[1], factor3), pipeparts.mkaudioamplify(pipeline, Yimag[0], factor4)))
Vj = mkadder(pipeline, list_srcs(pipeline, pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Yreal[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor1), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Yreal[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor2)))
Vj = pipeparts.mkcapsfilter(pipeline, Vj, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)
Vjplus3 = mkadder(pipeline, list_srcs(pipeline, pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Yimag[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor3), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Yimag[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor4)))
Vjplus3 = pipeparts.mkcapsfilter(pipeline, Vjplus3, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)
MV_matrix[j] = Vj
MV_matrix[num_stages + j] = Vjplus3
# Now let's compute the elements of the matrix M, given by Eqs. 5.2.70 - 5.2.77
# Many of the elements are constant, so make a stream of ones to multiply
if num_stages > 1:
ones = pipeparts.mktee(pipeline, mkpow(pipeline, Yreal[0], exponent = 0.0))
ones = pipeparts.mktee(pipeline, mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, Yreal[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 0.0))
print("exact kappas 10")
for j in range(num_stages):
# Time-dependent matrix elements
factor = pow(freqs[0], -2) - pow(freqs[1], -2)
addend = (pow(freqs[0], -2) - pow(freqs[2 + j], -2)) * EPICS[2 * ((1 + num_stages) + 1 + j)] + (pow(freqs[2 + j], -2) - pow(freqs[1], -2)) * EPICS[2 * (1 + j)] - (pow(freqs[0], -2) - pow(freqs[1], -2)) * EPICS[2 * ((2 + j) * (1 + num_stages) + 1 + j)]
Mjj = pipeparts.mkgeneric(pipeline, pipeparts.mkaudioamplify(pipeline, CAXreal[j], factor), "lal_add_constant", value = addend)
Mjj = pipeparts.mkgeneric(pipeline, pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, CAXreal[j], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor), "lal_add_constant", value = addend)
factor = -2.0 * numpy.pi * freqs[2 + j] * (pow(freqs[0], -2) - pow(freqs[1], -2))
addend = -2.0 * numpy.pi * freqs[1] * (pow(freqs[0], -2) - pow(freqs[2 + j], -2)) * EPICS[1 + 2 * ((1 + num_stages) + 1 + j)] - 2.0 * numpy.pi * freqs[0] * (pow(freqs[2 + j], -2) - pow(freqs[1], -2)) * EPICS[1 + 2 * (1 + j)] + 2.0 * numpy.pi * freqs[2 + j] * (pow(freqs[0], -2) - pow(freqs[1], -2)) * EPICS[1 + 2 * ((2 + j) * (1 + num_stages) + 1 + j)]
Mjjplus3 = pipeparts.mkgeneric(pipeline, pipeparts.mkaudioamplify(pipeline, CAXimag[j], factor), "lal_add_constant", value = addend)
Mjjplus3 = pipeparts.mkgeneric(pipeline, pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, CAXimag[j], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor), "lal_add_constant", value = addend)
factor = -freqs[2 + j] * (pow(freqs[1], 2) - pow(freqs[0], 2))
addend = freqs[1] * (pow(freqs[0], 2) - pow(freqs[2 + j], 2)) * EPICS[1 + 2 * ((1 + num_stages) + 1 + j)] + freqs[0] * (pow(freqs[2 + j], 2) - pow(freqs[1], 2)) * EPICS[1 + 2 * (1 + j)] + freqs[2 + j] * (pow(freqs[1], 2) - pow(freqs[0], 2)) * EPICS[1 + 2 * ((2 + j) * (1 + num_stages) + 1 + j)]
Mjplus3j = pipeparts.mkgeneric(pipeline, pipeparts.mkaudioamplify(pipeline, CAXimag[j], factor), "lal_add_constant", value = addend)
Mjplus3j = pipeparts.mkgeneric(pipeline, pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, CAXimag[j], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor), "lal_add_constant", value = addend)
factor = -2.0 * numpy.pi * pow(freqs[2 + j], 2) * (pow(freqs[1], 2) - pow(freqs[0], 2))
addend = 2.0 * numpy.pi * pow(freqs[1], 2) * (pow(freqs[0], 2) - pow(freqs[2 + j], 2)) * EPICS[2 * ((1 + num_stages) + 1 + j)] + 2.0 * numpy.pi * pow(freqs[0], 2) * (pow(freqs[2 + j], 2) - pow(freqs[1], 2)) * EPICS[2 * (1 + j)] + 2.0 * numpy.pi * pow(freqs[2 + j], 2) * (pow(freqs[1], 2) - pow(freqs[0], 2)) * EPICS[2 * ((2 + j) * (1 + num_stages) + 1 + j)]
Mjplus3jplus3 = pipeparts.mkgeneric(pipeline, pipeparts.mkaudioamplify(pipeline, CAXreal[j], factor), "lal_add_constant", value = addend)
Mjplus3jplus3 = pipeparts.mkgeneric(pipeline, pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, CAXreal[j], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor), "lal_add_constant", value = addend)
# Add these into the matrix
MV_matrix[(1 + j) * 2 * num_stages + j] = Mjj
MV_matrix[(1 + j) * 2 * num_stages + num_stages + j] = Mjjplus3
......@@ -1036,13 +1044,13 @@ def compute_exact_kappas_from_filters_file(pipeline, X, freqs, EPICS, rate):
knotequalj.remove(j)
for k in knotequalj:
factor = (pow(freqs[0], -2) - pow(freqs[2 + j], -2)) * EPICS[2 * ((1 + num_stages) + 1 + k)] + (pow(freqs[2 + j], -2) - pow(freqs[1], -2)) * EPICS[2 * (1 + k)] - (pow(freqs[0], -2) - pow(freqs[1], -2)) * EPICS[2 * ((2 + j) * (1 + num_stages) + 1 + k)]
Mjk = pipeparts.mkaudioamplify(pipeline, ones, factor)
Mjk = pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, ones, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor)
factor = -2.0 * numpy.pi * freqs[1] * (pow(freqs[0], -2) - pow(freqs[2 + j], -2)) * EPICS[1 + 2 * ((1 + num_stages) + 1 + k)] - 2.0 * numpy.pi * freqs[0] * (pow(freqs[2 + j], -2) - pow(freqs[1], -2)) * EPICS[1 + 2 * (1 + k)] + 2.0 * numpy.pi * freqs[2 + j] * (pow(freqs[0], -2) - pow(freqs[1], -2)) * EPICS[1 + 2 * ((2 + j) * (1 + num_stages) + 1 + k)]
Mjkplus3 = pipeparts.mkaudioamplify(pipeline, ones, factor)
Mjkplus3 = pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, ones, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor)
factor = freqs[1] * (pow(freqs[0], 2) - pow(freqs[2 + j], 2)) * EPICS[1 + 2 * ((1 + num_stages) + 1 + k)] + freqs[0] * (pow(freqs[2 + j], 2) - pow(freqs[1], 2)) * EPICS[1 + 2 * (1 + k)] + freqs[2 + j] * (pow(freqs[1], 2) - pow(freqs[0], 2)) * EPICS[1 + 2 * ((2 + j) * (1 + num_stages) + 1 + k)]
Mjplus3k = pipeparts.mkaudioamplify(pipeline, ones, factor)
Mjplus3k = pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, ones, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor)
factor = 2.0 * numpy.pi * pow(freqs[1], 2) * (pow(freqs[0], 2) - pow(freqs[2 + j], 2)) * EPICS[2 * ((1 + num_stages) + 1 + k)] + 2.0 * numpy.pi * pow(freqs[0], 2) * (pow(freqs[2 + j], 2) - pow(freqs[1], 2)) * EPICS[2 * (1 + k)] + 2.0 * numpy.pi * pow(freqs[2 + j], 2) * (pow(freqs[1], 2) - pow(freqs[0], 2)) * EPICS[2 * ((2 + j) * (1 + num_stages) + 1 + k)]
Mjplus3kplus3 = pipeparts.mkaudioamplify(pipeline, ones, factor)
Mjplus3kplus3 = pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, ones, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), factor)
# Add these into the matrix
MV_matrix[(1 + j) * 2 * num_stages + k] = Mjk
MV_matrix[(1 + j) * 2 * num_stages + num_stages + k] = Mjkplus3
......@@ -1052,9 +1060,11 @@ def compute_exact_kappas_from_filters_file(pipeline, X, freqs, EPICS, rate):
print("exact kappas 20")
# Now pass these to the matrix solver to find kappa_T, kappa_P, kappa_U, tau_T, tau_P, and tau_U.
MV_matrix = mkinterleave(pipeline, MV_matrix)
MV_matrix = pipeparts.mkcapsfilter(pipeline, MV_matrix, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=%d" % (rate, (2 * num_stages) * (2 * num_stages + 1)))
kappas = pipeparts.mkgeneric(pipeline, MV_matrix, "lal_matrixsolver")
kappas = list(mkdeinterleave(pipeline, kappas, 2 * num_stages))
for i in range(len(kappas)):
kappas[i] = pipeparts.mkcapsfilter(pipeline, kappas[i], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)
if i >= len(kappas) / 2:
kappas[i] = mkmultiplier(pipeline, [kappas[i], mkpow(pipeline, kappas[i - len(kappas) / 2], exponent = -1.0)])
kappas[i] = pipeparts.mktee(pipeline, kappas[i])
......@@ -1065,84 +1075,84 @@ def compute_exact_kappas_from_filters_file(pipeline, X, freqs, EPICS, rate):
for n in range(2):
Gres_components = []
for j in range(num_stages):
kappajGresjatn = pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, kappas[j], 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, kappas[num_stages + j], matrix = [[0, 2.0 * numpy.pi * freqs[n]]]))
i_omega_tau = pipeparts.mkcapsfilter(pipeline, i_omega_tau, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0" % rate)
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)]]]))
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 = pipeparts.mktee(pipeline, mkadder(pipeline, Gres_components))
Gresreal.append(pipeparts.mktee(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkgeneric(pipeline, Gres, "creal"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0" % rate)))
Gresimag.append(pipeparts.mktee(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkgeneric(pipeline, Gres, "cimag"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0" % rate)))
Gresreal.append(pipeparts.mktee(pipeline, mkcapsfiltersetter(pipeline, pipeparts.mkgeneric(pipeline, Gres, "creal"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)))
Gresimag.append(pipeparts.mktee(pipeline, mkcapsfiltersetter(pipeline, pipeparts.mkgeneric(pipeline, Gres, "cimag"), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)))
print("exact kappas 30")
# Next, let us find the H's, I's, Xi, and zeta, defined by Eqs. 5.2.48, 5.2.49, 5.2.50, and 5.2.58, respectively.
H1 = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, list_srcs(pipeline, Gresreal[0], pipeparts.mkaudioamplify(pipeline, Yreal[0], -1.0))), pow(freqs[0], 2)))
H1 = pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, list_srcs(pipeline, pipeparts.mkcapsfilter(pipeline, Gresreal[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkcapsfilter(pipeline, Yreal[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -1.0))), pow(freqs[0], 2))
print("exact kappas 31")
H2 = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, list_srcs(pipeline, Gresreal[1], pipeparts.mkaudioamplify(pipeline, Yreal[1], -1.0))), pow(freqs[1], 2)))
H2 = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, list_srcs(pipeline, pipeparts.mkcapsfilter(pipeline, Gresreal[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkcapsfilter(pipeline, Yreal[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -1.0))), pow(freqs[1], 2)))
print("exact kappas 32")
I1 = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, list_srcs(pipeline, Gresimag[0], pipeparts.mkaudioamplify(pipeline, Yimag[0], -1.0))), pow(freqs[0], 2)))
I1 = pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, list_srcs(pipeline, pipeparts.mkcapsfilter(pipeline, Gresimag[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Yimag[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -1.0))), pow(freqs[0], 2))
print("exact kappas 33")
I2 = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, list_srcs(pipeline, Gresimag[1], pipeparts.mkaudioamplify(pipeline, Yimag[1], -1.0))), pow(freqs[1], 2)))
I2 = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, list_srcs(pipeline, pipeparts.mkcapsfilter(pipeline, Gresimag[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Yimag[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -1.0))), pow(freqs[1], 2)))
print("exact kappas 34")
Xi = pipeparts.mktee(pipeline, mkadder(pipeline, list_srcs(pipeline, pipeparts.mkaudioamplify(pipeline, H1, pow(freqs[1], 2) / (pow(freqs[0], 2) - pow(freqs[1], 2))), pipeparts.mkaudioamplify(pipeline, H2, -pow(freqs[0], 2) / (pow(freqs[0], 2) - pow(freqs[1], 2))))))
Xi = pipeparts.mktee(pipeline, mkadder(pipeline, list_srcs(pipeline, pipeparts.mkaudioamplify(pipeline, H1, pow(freqs[1], 2) / (pow(freqs[0], 2) - pow(freqs[1], 2))), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, H2, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -pow(freqs[0], 2) / (pow(freqs[0], 2) - pow(freqs[1], 2))))))
print("exact kappas 35")
zeta = pipeparts.mktee(pipeline, mkadder(pipeline, list_srcs(pipeline, pipeparts.mkaudioamplify(pipeline, I1, 1.0 / freqs[0]), pipeparts.mkaudioamplify(pipeline, I2, -1.0 / freqs[1]))))
zeta = pipeparts.mktee(pipeline, mkadder(pipeline, list_srcs(pipeline, pipeparts.mkaudioamplify(pipeline, I1, 1.0 / freqs[0]), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, I2, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -1.0 / freqs[1]))))
print("exact kappas 40")
# Now, we define the coefficients of the quartic equation a kC^4 + b kC^3 + c kC^2 + d kC + e = 0.
a = pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [Xi, mkadder(pipeline, [H2, Xi]), mkpow(pipeline, zeta, exponent = 2.0)]), pow(freqs[1], 2))
a = pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, Xi, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), mkadder(pipeline, [pipeparts.mkcapsfilter(pipeline, H2, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, Xi, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, zeta, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0)]), pow(freqs[1], 2))
print("exact kappas 41")
a = pipeparts.mktee(pipeline, pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, a, matrix = [[1.0, 0.0]])))
print("exact kappas 42")
b = mkadder(pipeline, [pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [mkadder(pipeline, [H2, Xi]), zeta, I2]), freqs[1] * (pow(freqs[1], 2) - pow(freqs[0], 2))), pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [mkadder(pipeline, [H2, pipeparts.mkaudioamplify(pipeline, Xi, 2.0)]), mkpow(pipeline, zeta, exponent = 2.0)]), pow(freqs[1], 4))])
b = mkadder(pipeline, [pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [mkadder(pipeline, [pipeparts.mkcapsfilter(pipeline, H2, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, Xi, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), pipeparts.mkcapsfilter(pipeline, zeta, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, I2, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), freqs[1] * (pow(freqs[1], 2) - pow(freqs[0], 2))), pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [mkadder(pipeline, [pipeparts.mkcapsfilter(pipeline, H2, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Xi, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), 2.0)]), mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, zeta, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0)]), pow(freqs[1], 4))])
print("exact kappas 43")
b = pipeparts.mktee(pipeline, pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, b, matrix = [[1.0, 0.0]])))
print("exact kappas 44")
c = mkadder(pipeline, [pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [zeta, I2]), pow(freqs[1], 3) * (pow(freqs[1], 2) - pow(freqs[0], 2))), pipeparts.mkaudioamplify(pipeline, mkpow(pipeline, mkadder(pipeline, [H2, Xi]), exponent = 2.0), pow(pow(freqs[1], 2) - pow(freqs[0], 2), 2)), pipeparts.mkaudioamplify(pipeline, mkpow(pipeline, zeta, exponent = 2.0), pow(freqs[1], 6))])
c = mkadder(pipeline, [pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, zeta, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, I2, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), pow(freqs[1], 3) * (pow(freqs[1], 2) - pow(freqs[0], 2))), pipeparts.mkaudioamplify(pipeline, mkpow(pipeline, mkadder(pipeline, [pipeparts.mkcapsfilter(pipeline, H2, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, Xi, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), exponent = 2.0), pow(pow(freqs[1], 2) - pow(freqs[0], 2), 2)), pipeparts.mkaudioamplify(pipeline, mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, zeta, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0), pow(freqs[1], 6))])
print("exact kappas 45")
c = pipeparts.mktee(pipeline, pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, c, matrix = [[1.0, 0.0]])))
print("exact kappas 46")
d = pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, [H2, Xi]), pow(freqs[1], 4) * pow(pow(freqs[1], 2) - pow(freqs[0], 2), 2))
d = pipeparts.mkaudioamplify(pipeline, mkadder(pipeline, [pipeparts.mkcapsfilter(pipeline, H2, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, Xi, "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), pow(freqs[1], 4) * pow(pow(freqs[1], 2) - pow(freqs[0], 2), 2))
d = pipeparts.mktee(pipeline, pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, d, matrix = [[1.0, 0.0]])))
e = pow(freqs[1], 4) * pow(pow(freqs[1], 2) - pow(freqs[0], 2), 2)
print("exact kappas 50")
# More parameters we need to define to solve the quartic...
Delta0 = pipeparts.mktee(pipeline, mkadder(pipeline, [mkpow(pipeline, c, exponent = 2.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [b, d]), -3.0, 0.0), complex_audioamplify(pipeline, a, 12 * e, 0.0)]))
Delta0 = pipeparts.mktee(pipeline, mkadder(pipeline, [mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, c, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, b, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, d, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), -3.0, 0.0), complex_audioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), 12 * e, 0.0)]))
print("exact kappas 51")
Delta1 = pipeparts.mktee(pipeline, mkadder(pipeline, [complex_audioamplify(pipeline, mkpow(pipeline, c, exponent = 3.0), 2.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [b, c, d]), -9.0, 0.0), complex_audioamplify(pipeline, mkpow(pipeline, b, exponent = 2.0), 27.0 * e, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [a, mkpow(pipeline, d, exponent = 2.0)]), 27.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [a, c]), -72.0 * e, 0.0)]))
Delta1 = pipeparts.mktee(pipeline, mkadder(pipeline, [complex_audioamplify(pipeline, mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, c, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 3.0), 2.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, b, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, c, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, d, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), -9.0, 0.0), complex_audioamplify(pipeline, mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, b, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0), 27.0 * e, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, d, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0)]), 27.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, c, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), -72.0 * e, 0.0)]))
print("exact kappas 52")
p = pipeparts.mktee(pipeline, mkmultiplier(pipeline, [mkadder(pipeline, [complex_audioamplify(pipeline, mkmultiplier(pipeline, [a, c]), 8.0, 0.0), complex_audioamplify(pipeline, mkpow(pipeline, b, exponent = 2.0), -3.0, 0.0)]), complex_audioamplify(pipeline, mkpow(pipeline, a, exponent = -2.0), 1.0 / 8.0, 0.0)]))
p = pipeparts.mktee(pipeline, mkmultiplier(pipeline, [mkadder(pipeline, [complex_audioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, c, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), 8.0, 0.0), complex_audioamplify(pipeline, mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, b, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0), -3.0, 0.0)]), complex_audioamplify(pipeline, mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = -2.0), 1.0 / 8.0, 0.0)]))
print("exact kappas 53")
q = mkmultiplier(pipeline, [mkadder(pipeline, [mkpow(pipeline, b, exponent = 3.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [a, b, c]), -4.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [mkpow(pipeline, a, exponent = 2.0), d]), 8.0, 0.0)]), complex_audioamplify(pipeline, mkpow(pipeline, a, exponent = -3.0), 1.0 / 8.0, 0.0)])
q = mkmultiplier(pipeline, [mkadder(pipeline, [mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, b, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 3.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, b, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, c, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), -4.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0), pipeparts.mkcapsfilter(pipeline, d, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), 8.0, 0.0)]), complex_audioamplify(pipeline, mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = -3.0), 1.0 / 8.0, 0.0)])
print("exact kappas 54")
Q0 = pipeparts.mktee(pipeline, mkpow(pipeline, complex_audioamplify(pipeline, mkadder(pipeline, [Delta1, mkpow(pipeline, mkadder(pipeline, [mkpow(pipeline, Delta1, exponent = 2.0), complex_audioamplify(pipeline, mkpow(pipeline, Delta0, exponent = 3.0), -4.0, 0.0)]), exponent = 0.5)]), 0.5, 0.0), exponent = 1.0 / 3.0))
Q0 = pipeparts.mktee(pipeline, mkpow(pipeline, complex_audioamplify(pipeline, mkadder(pipeline, [pipeparts.mkcapsfilter(pipeline, Delta1, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), mkpow(pipeline, mkadder(pipeline, [mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, Delta1, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0), complex_audioamplify(pipeline, mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, Delta0, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 3.0), -4.0, 0.0)]), exponent = 0.5)]), 0.5, 0.0), exponent = 1.0 / 3.0))
print("exact kappas 55")
S0 = pipeparts.mktee(pipeline, mkpow(pipeline, mkadder(pipeline, [complex_audioamplify(pipeline, p, -1.0 / 6.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [mkpow(pipeline, a, exponent = -1.0), mkadder(pipeline, [Q0, complex_division(pipeline, Delta0, Q0)])]), 1.0 / 12.0, 0.0)]), exponent = 0.5))
S0 = pipeparts.mktee(pipeline, mkpow(pipeline, mkadder(pipeline, [complex_audioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, p, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -1.0 / 6.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = -1.0), mkadder(pipeline, [pipeparts.mkcapsfilter(pipeline, Q0, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), complex_division(pipeline, pipeparts.mkcapsfilter(pipeline, Delta0, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, Q0, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate))])]), 1.0 / 12.0, 0.0)]), exponent = 0.5))
print("exact kappas 56")
R0 = pipeparts.mktee(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkgeneric(pipeline, mkadder(pipeline, [mkpow(pipeline, b, exponent = 3.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [d, mkpow(pipeline, a, exponent = 2.0)]), 8.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [a, b, c]), -4.0, 0.0)]), "creal"), "audio/x-raw,format=F64LE,channel-mask=(bitmask)0x0"))
R0 = pipeparts.mktee(pipeline, mkcapsfiltersetter(pipeline, pipeparts.mkgeneric(pipeline, mkadder(pipeline, [mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, b, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 3.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, d, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0)]), 8.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, b, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, c, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), -4.0, 0.0)]), "creal"), "audio/x-raw,format=F64LE,channel-mask=(bitmask)0x0,channels=1"))
print("exact kappas 57")
R0sign = pipeparts.mktee(pipeline, pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, mkmultiplier(pipeline, [R0, mkpow(pipeline, pipeparts.mkgeneric(pipeline, R0, "cabs"), exponent = -1.0)]), matrix = [[1.0, 0.0]])))
R0sign = pipeparts.mktee(pipeline, pipeparts.mktogglecomplex(pipeline, pipeparts.mkmatrixmixer(pipeline, mkmultiplier(pipeline, [pipeparts.mkgeneric(pipeline, pipeparts.mkcapsfilter(pipeline, R0, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), "creal"), mkpow(pipeline, pipeparts.mkgeneric(pipeline, pipeparts.mkcapsfilter(pipeline, R0, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), "cabs"), exponent = -1.0)]), matrix = [[1.0, 0.0]])))
print("exact kappas 60")
# Finally, compute kappa_C
kappa_C = mkadder(pipeline, [complex_audioamplify(pipeline, complex_division(pipeline, b, a), -0.25, 0.0), mkmultiplier(pipeline, [R0sign, S0]), complex_audioamplify(pipeline, mkpow(pipeline, mkadder(pipeline, [complex_audioamplify(pipeline, mkpow(pipeline, S0, exponent = 2.0), -4.0, 0.0), complex_audioamplify(pipeline, p, -2.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [q, R0sign, mkpow(pipeline, S0, exponent = -1.0)]), -1.0, 0.0)]), exponent = 0.5), -0.5, 0.0)])
kappa_C = pipeparts.mktee(pipeline, pipeparts.mkcapsfilter(pipeline, pipeparts.mkgeneric(pipeline, kappa_C, "creal"), "audio/x-raw,format=F64LE,channel-mask=(bitmask)0x0"))
kappa_C = mkadder(pipeline, [complex_audioamplify(pipeline, complex_division(pipeline, pipeparts.mkcapsfilter(pipeline, b, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, a, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)), -0.25, 0.0), mkmultiplier(pipeline, [pipeparts.mkcapsfilter(pipeline, R0sign, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkcapsfilter(pipeline, S0, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), complex_audioamplify(pipeline, mkpow(pipeline, mkadder(pipeline, [complex_audioamplify(pipeline, mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, S0, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = 2.0), -4.0, 0.0), complex_audioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, p, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -2.0, 0.0), complex_audioamplify(pipeline, mkmultiplier(pipeline, [q, pipeparts.mkcapsfilter(pipeline, R0sign, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), mkpow(pipeline, pipeparts.mkcapsfilter(pipeline, S0, "audio/x-raw,format=Z128LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), exponent = -1.0)]), -1.0, 0.0)]), exponent = 0.5), -0.5, 0.0)])
kappa_C = pipeparts.mktee(pipeline, mkcapsfiltersetter(pipeline, pipeparts.mkgeneric(pipeline, kappa_C, "creal"), "audio/x-raw,format=F64LE,channel-mask=(bitmask)0x0,channels=1"))
kappas.append(kappa_C)
print("exact kappas 70")
# Next, compute f_cc using Eq. 5.2.45
f_cc = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkpow(pipeline, mkmultiplier(pipeline, [kappa_C, mkadder(pipeline, [pipeparts.mkaudioamplify(pipeline, Gresimag[0], freqs[0]), pipeparts.mkaudioamplify(pipeline, Yimag[0], -freqs[0]), pipeparts.mkaudioamplify(pipeline, Gresimag[1], -freqs[1]), pipeparts.mkaudioamplify(pipeline, Yimag[1], freqs[1])])]), exponent = -1.0), pow(freqs[1], 2) - pow(freqs[0], 2)))
f_cc = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkpow(pipeline, mkmultiplier(pipeline, [kappa_C, mkadder(pipeline, [pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Gresimag[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), freqs[0]), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Yimag[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -freqs[0]), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Gresimag[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -freqs[1]), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Yimag[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), freqs[1])])]), exponent = -1.0), pow(freqs[1], 2) - pow(freqs[0], 2)))
kappas.append(f_cc)
print("exact kappas 80")
# Compute f_s^2 using Eq. 5.2.41
f_s_squared = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [kappa_C, mkadder(pipeline, [Yreal[1], pipeparts.mkaudioamplify(pipeline, Gresreal[1], -1.0), pipeparts.mkaudioamplify(pipeline, Yreal[0], -1.0), Gresreal[0]])]), 1.0 / (pow(freqs[1], -2) - pow(freqs[0], -2))))
f_s_squared = pipeparts.mktee(pipeline, pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [kappa_C, mkadder(pipeline, [pipeparts.mkcapsfilter(pipeline, Yreal[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Gresreal[1], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -1.0), pipeparts.mkaudioamplify(pipeline, pipeparts.mkcapsfilter(pipeline, Yreal[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate), -1.0), pipeparts.mkcapsfilter(pipeline, Gresreal[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)])]), 1.0 / (pow(freqs[1], -2) - pow(freqs[0], -2))))
kappas.append(f_s_squared)
print("exact kappas 90")
# Compute f_s / Q using Eq. 5.2.33
f_s_over_Q = pipeparts.mktee(pipeline, mkmultiplier(pipeline, [f_cc, pipeparts.mkgeneric(pipeline, mkadder(pipeline, [mkmultiplier(pipeline, [kappa_C, Yreal[0]]), pipeparts.mkaudioamplify(pipeline, f_s_squared, -pow(freqs[0], -2)), pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [kappa_C, Gresreal[0]]), -1.0)]), "lal_add_constant", value = -1.0)]))
f_s_over_Q = pipeparts.mktee(pipeline, mkmultiplier(pipeline, [f_cc, pipeparts.mkgeneric(pipeline, mkadder(pipeline, [mkmultiplier(pipeline, [kappa_C, pipeparts.mkcapsfilter(pipeline, Yreal[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), pipeparts.mkaudioamplify(pipeline, f_s_squared, -pow(freqs[0], -2)), pipeparts.mkaudioamplify(pipeline, mkmultiplier(pipeline, [kappa_C, pipeparts.mkcapsfilter(pipeline, Gresreal[0], "audio/x-raw,format=F64LE,rate=%d,channel-mask=(bitmask)0x0,channels=1" % rate)]), -1.0)]), "lal_add_constant", value = -1.0)]))
kappas.append(f_s_over_Q)
print("exact kappas done")
......@@ -1356,10 +1366,10 @@ def clean_data(pipeline, signal, signal_rate, witnesses, witness_rate, fft_lengt
witnesses = list(witnesses)
witness_tees = []
for i in range(0, len(witnesses)):
witnesses[i] = mkresample(pipeline, witnesses[i], 5, False, witness_rate)
witnesses[i] = mkresample(pipeline, witnesses[i], 4, False, witness_rate)
witness_tees.append(pipeparts.mktee(pipeline, witnesses[i]))
resampled_signal = mkresample(pipeline, signal_tee, 5, False, witness_rate)
resampled_signal = mkresample(pipeline, signal_tee, 4, False, witness_rate)
transfer_functions = mkinterleave(pipeline, numpy.insert(witness_tees, 0, resampled_signal, axis = 0))
if noisesub_gate_bit is not None:
transfer_functions = mkgate(pipeline, transfer_functions, noisesub_gate_bit, 1)
......@@ -1373,7 +1383,7 @@ def clean_data(pipeline, signal, signal_rate, witnesses, witness_rate, fft_lengt
else:
minus_noise = pipeparts.mkgeneric(pipeline, highpass(pipeline, witness_tees[i], witness_rate, fcut = high_pass), "lal_tdwhiten", kernel = numpy.zeros(fir_length), latency = fir_length / 2, taper_length = filter_taper_length)
transfer_functions.connect("notify::fir-filters", update_filters, minus_noise, "fir_filters", "kernel", i)
signal_minus_noise.append(mkresample(pipeline, minus_noise, 5, False, signal_rate))
signal_minus_noise.append(mkresample(pipeline, minus_noise, 4, False, signal_rate))
return mkadder(pipeline, tuple(signal_minus_noise))
......
......@@ -99,7 +99,7 @@ if options.write_ASD_txt:
numpy.savetxt("%s_%s_ASD_%d_%d%s.txt" % (ifo, channel_list[i], start_time, end_time, options.filename_suffix), asd_txt)
# Plot ASDs
colors = ['blue', 'limegreen', "royalblue", "deepskyblue", "red", "yellow", "purple", "pink"]
colors = ['blue', 'limegreen', "red", "deepskyblue", "red", "yellow", "purple", "pink"]
if len(ASD_title):
plot = ASDs[0].plot(color = colors[0], title = "%s" % ASD_title, linewidth = 0.75, label = ASD_labels[0])
else:
......
......@@ -10,14 +10,16 @@ OBSRUN = O3
START = $(shell echo 1240617600 - 256 | bc)
#1238997618
END = $(shell echo 1240617600 + 512 + 256 | bc)
END = $(shell echo 1240617600 + 1024 + 256 | bc)
#1240876818
SHMRUNTIME = 36000
# How much time does the calibration need to settle at the start and end?
PLOT_WARMUP_TIME = 256
PLOT_COOLDOWN_TIME = 256
GDSCONFIGS = Filters/O3/GDSFilters/H1GDS_1239476409_test.ini
GDSCONFIGS = Filters/O3/GDSFilters/H1GDS_1239476409_testAllCorrections.ini
LOWLATENCYCONFIGS = Filters/O3/GDSFilters/H1GDS_lowLatency_1239476409.ini
LOWLATENCYASYMMETRICCONFIGS = Filters/O3/GDSFilters/H1GDS_lowLatency_asymmetric_1239476409.ini
DCSCONFIGS = Filters/O3/GDSFilters/H1DCS_C01_1239472998_test.ini
DCSLINESCONFIGS = ../../config_files/O2/H1/tests/H1DCS_AllCorrections_CleaningLines.ini
#../../config_files/O2/H1/tests/H1DCS_AllCorrections_Cleaning.ini
......@@ -29,7 +31,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: noise_subtraction_ASD_DCS
all: highpass_filter_ASD_GDS filters_tf_GDS
###############################################
### These commands should change less often ###
......@@ -73,6 +75,14 @@ $(IFO)1_hoft_GDS_frames.cache: $(IFO)1_easy_raw_frames.cache filters framesdir
GST_DEBUG=3 gstlal_compute_strain --gps-start-time $(START) --gps-end-time $(END) --frame-cache $(IFO)1_easy_raw_frames.cache --output-path Frames/$(OBSRUN)/$(IFO)1/GDS/ --frame-duration=64 --frames-per-file=1 --wings=0 --config-file $(GDSCONFIGS)
ls Frames/$(OBSRUN)/$(IFO)1/GDS/$(IFO)-$(IFO)1GDS-*.gwf | lalapps_path2cache > $@
$(IFO)1_lowLatency_GDS_frames.cache: $(IFO)1_easy_raw_frames.cache filters framesdir
GST_DEBUG=3 gstlal_compute_strain --gps-start-time $(START) --gps-end-time $(END) --frame-cache $(IFO)1_easy_raw_frames.cache --output-path Frames/$(OBSRUN)/$(IFO)1/GDS/ --frame-duration=64 --frames-per-file=1 --wings=0 --config-file $(LOWLATENCYCONFIGS)
ls Frames/$(OBSRUN)/$(IFO)1/GDS/$(IFO)-$(IFO)1GDS_LOWLATENCY-*.gwf | lalapps_path2cache > $@
$(IFO)1_lowLatency_asymmetric_GDS_frames.cache: $(IFO)1_easy_raw_frames.cache filters framesdir
GST_DEBUG=3 gstlal_compute_strain --gps-start-time $(START) --gps-end-time $(END) --frame-cache $(IFO)1_easy_raw_frames.cache --output-path Frames/$(OBSRUN)/$(IFO)1/GDS/ --frame-duration=64 --frames-per-file=1 --wings=0 --config-file $(LOWLATENCYASYMMETRICCONFIGS)
ls Frames/$(OBSRUN)/$(IFO)1/GDS/$(IFO)-$(IFO)1GDS_LOWLATENCY_ASYMMETRIC-*.gwf | lalapps_path2cache > $@
$(IFO)1_hoft_DCS_frames.cache: $(IFO)1_easy_raw_frames.cache filters framesdir
GST_DEBUG=3 gstlal_compute_strain --gps-start-time $(START) --gps-end-time $(END) --frame-cache $(IFO)1_easy_raw_frames.cache --output-path Frames/$(OBSRUN)/$(IFO)1/DCS/ --frame-duration=64 --frames-per-file=1 --wings=256 --config-file $(DCSCONFIGS)
ls Frames/$(OBSRUN)/$(IFO)1/DCS/$(IFO)-$(IFO)1DCS-*.gwf | lalapps_path2cache > $@
......@@ -119,7 +129,7 @@ pcal_DCS_transfer_functions: $(IFO)1_easy_raw_frames.cache $(IFO)1_hoft_DCS_fram
lines_ratio_DCS: $(IFO)1_hoft_DCS_frames.cache
python demod_ratio_timeseries.py --ifo $(IFO)1 --gps-end-time $(PLOT_END) --gps-start-time $(PLOT_START) --denominator-frame-cache $(IFO)1_hoft_DCS_frames.cache --numerator-frame-cache $(IFO)1_hoft_DCS_frames.cache --denominator-channel-name 'DCS-CALIB_STRAIN' --numerator-channel-name 'DCS-CALIB_STRAIN_CLEAN' --frequencies '35.9,36.7,331.9,1083.7;60,120,180' --magnitude-ranges '0.0,0.1;0.0,1.0' --phase-ranges '-180.0,180.0;-180.0,180.0' --plot-titles '$(IFO)1 Calibration Line Subtraction;$(IFO)1 Power Mains Line Subtraction'
filters_tf_GDS: $(IFO)1_hoft_GDS_frames.cache
filters_tf_GDS: $(IFO)1_lowLatency_GDS_frames.cache
python plot_filters_transfer_function.py --tf-frequency-min 0.5 --tf-frequency-max 8192 --ratio-frequency-min 10 --ratio-frequency-max 8192 --ratio-magnitude-min 0.9 --ratio-magnitude-max 1.1 --tf-phase-min -180 --tf-phase-max 180 --ratio-phase-min -10 --ratio-phase-max 10
filters_tf_DCS: $(IFO)1_hoft_DCS_frames.cache
......@@ -175,9 +185,9 @@ noise_subtraction_ASD_GDS: $(IFO)1_hoft_GDS_frames.cache
noise_subtraction_ASD_C00: $(IFO)1_C00_frames.cache
./ASD_plots --ifo $(IFO)1 --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --frame-cache-list '$(IFO)1_C00_frames.cache,$(IFO)1_C00_frames.cache' --channel-list 'GDS-CALIB_STRAIN,GDS-CALIB_STRAIN_CLEAN'
highpass_filter_ASD_GDS:
#$(IFO)1_hoft_GDS_OLD_frames.cache $(IFO)1_hoft_GDS_BETTER_frames.cache $(IFO)1_hoft_GDS_BEST_frames.cache
./ASD_plots --ifo $(IFO)1 --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --frame-cache-list '$(IFO)1_hoft_GDS_OLD_frames.cache,$(IFO)1_hoft_GDS_BETTER_frames.cache,$(IFO)1_hoft_GDS_BEST_frames.cache' --channel-list 'GDS-CALIB_STRAIN,GDS-CALIB_STRAIN,GDS-CALIB_STRAIN' --ASD-fmin 0.1 --ASD-labels '6s Tukey-windowed actuation filter,3.5s Kaiser-windowed actuation filter,Additional high-pass in residual path'
highpass_filter_ASD_GDS:
#$(IFO)1_C00_frames.cache $(IFO)1_lowLatency_GDS_frames.cache $(IFO)1_lowLatency_asymmetric_GDS_frames.cache
./ASD_plots --ifo $(IFO)1 --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --frame-cache-list '$(IFO)1_C00_frames.cache,$(IFO)1_lowLatency_GDS_frames.cache,$(IFO)1_lowLatency_asymmetric_GDS_frames.cache' --channel-list 'GDS-CALIB_STRAIN,GDS-CALIB_STRAIN,GDS-CALIB_STRAIN' --ASD-fmin 0.5 --ASD-labels 'C00,0.95s Latency Symmetric Filters,0.95s Latency Asymmetric Filters'
noise_subtraction_ASD_DCH_DCS: $(IFO)1_hoft_DCS_frames.cache $(IFO)1_clean_C02_frames.cache
./ASD_comparison_plots --ifo $(IFO)1 --gps-start-time $(PLOT_START) --gps-end-time $(PLOT_END) --raw-frame-cache $(IFO)1_clean_C02_frames.cache --calcs-channel-name DCH-CLEAN_STRAIN_C02 --hoft-frame-cache $(IFO)1_hoft_DCS_frames.cache --hoft-channel-name DCS-CALIB_STRAIN_CLEAN
......
......@@ -59,20 +59,23 @@ def lal_gate_01(pipeline, name):
#
# Make a sine wave
src = test_common.test_src(pipeline, buffer_length = buffer_length, rate = rate, test_duration = test_duration, wave = 0, volume = 1.0, freq = frequency, width = 64)
src = test_common.complex_test_src(pipeline, buffer_length = buffer_length, rate = rate, test_duration = test_duration, wave = 0, freq = frequency, width = 64)
# Add a DC offset
head = pipeparts.mkgeneric(pipeline, src, "lal_add_constant", value = DC_offset)
head = pipeparts.mktee(pipeline, head)
pipeparts.mknxydumpsink(pipeline, head, "%s_in.txt" % name)
control = calibration_parts.mkqueue(pipeline, head, min_length = 8)
control = pipeparts.mkgeneric(pipeline, control, "splitcounter", name = "control")
head = pipeparts.mkgeneric(pipeline, head, "splitcounter", name = "before")
control = test_common.test_src(pipeline, buffer_length = buffer_length, rate = rate, test_duration = test_duration, wave = 0, volume = 1.0, freq = frequency, width = 64, src_suffix = "_control")
control = pipeparts.mkgeneric(pipeline, control, "lal_add_constant", value = DC_offset)
#control = pipeparts.mkgeneric(pipeline, control, "splitcounter", name = "control")
#head = pipeparts.mkgeneric(pipeline, head, "splitcounter", name = "before")
# Gate it
head = calibration_parts.mkgate(pipeline, head, control, threshold, attack_length = int(attack_length * rate), hold_length = int(hold_length * rate))
head = pipeparts.mkgeneric(pipeline, head, "splitcounter", name = "after")
pipeparts.mknxydumpsink(pipeline, head, "%s_out.txt" % name)
#head = pipeparts.mkgeneric(pipeline, head, "splitcounter", name = "after")
real, imag = calibration_parts.split_into_real(pipeline, head)
pipeparts.mknxydumpsink(pipeline, real, "%s_real_out.txt" % name)
pipeparts.mknxydumpsink(pipeline, imag, "%s_imag_out.txt" % name)
#
# done
......
......@@ -405,8 +405,8 @@ def mkndssrc(pipeline, host, instrument, channel_name, channel_type, blocksize =
## Adds a <a href="@gstdoc/gstreamer-plugins-capsfilter.html">capsfilter</a> element to a pipeline with useful default properties
def mkcapsfilter(pipeline, src, caps):
return mkgeneric(pipeline, src, "capsfilter", caps = Gst.Caps.from_string(caps))
def mkcapsfilter(pipeline, src, caps, **properties):
return mkgeneric(pipeline, src, "capsfilter", caps = Gst.Caps.from_string(caps), **properties)
## Adds a <a href="@gstpluginsgooddoc/gst-plugins-good-plugins-capssetter.html">capssetter</a> element to a pipeline with useful default properties
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment