gstlal-calibration: plotting script update: use state vector

tests/check_calibration/approxKappasErrorPlot.py

@@ -327,6 +327,11 @@ if oldfs_squared < 0:
 	oldQ = -1j * oldQ
 oldfs_over_Q = np.real(np.sqrt(complex(oldfs_squared)) / oldQ)
 
+if os.path.exists("%s_approxKappasTimes_%d-%d.txt" % (ifo, int(tstart), int(dursec))):
+	os.remove("%s_approxKappasTimes_%d-%d.txt" % (ifo, int(tstart), int(dursec)))
+if os.path.exists("%s_approxKappasTDCFs_%d-%d.txt" % (ifo, int(tstart), int(dursec))):
+	os.remove("%s_approxKappasTDCFs_%d-%d.txt" % (ifo, int(tstart), int(dursec)))
+
 config_index = 0
 num_tries = 0
 for i in range(len(lines)):
@@ -353,12 +358,18 @@ for i in range(len(lines)):
 	# Try this first
 	frdump0 = os.popen('FrDump -d 4 -f %d -l %d -i %s' % (int(t), int(t), lines[i].split('localhost')[1])).read()
 	if 'Vector:%s:%sCALIB_KAPPA_TST_REAL%s' % (ifo, chan_prefixes[config_index], chan_suffixes[config_index]) in frdump0:
-		frdump = frdump0
+		dqvecmin = int(frdump0.split('Vector:%s:%sCALIB_STATE_VECTOR%s' % (ifo, chan_prefixes[config_index], chan_suffixes[config_index]))[1].split('min ')[1].split(' ')[0])
+		dqvecmax = int(frdump0.split('Vector:%s:%sCALIB_STATE_VECTOR%s' % (ifo, chan_prefixes[config_index], chan_suffixes[config_index]))[1].split('max ')[1].split(' ')[0])
+		if (dqvecmin % 2 == 1 and dqvecmax % 2 == 1) or i == 0:
+			frdump = frdump0
 	else:
 		# Then try reading all the frames in the file
 		frdump0 = os.popen('FrDump -d 4 -i %s' % lines[i].split('localhost')[1]).read()
 		if 'Vector:%s:%sCALIB_KAPPA_TST_REAL%s' % (ifo, chan_prefixes[config_index], chan_suffixes[config_index]) in frdump0:
-			frdump = frdump0
+			dqvecmin = int(frdump0.split('Vector:%s:%sCALIB_STATE_VECTOR%s' % (ifo, chan_prefixes[config_index], chan_suffixes[config_index]))[1].split('min ')[1].split(' ')[0])
+			dqvecmax = int(frdump0.split('Vector:%s:%sCALIB_STATE_VECTOR%s' % (ifo, chan_prefixes[config_index], chan_suffixes[config_index]))[1].split('max ')[1].split(' ')[0])
+			if (dqvecmin % 2 == 1 and dqvecmax % 2 == 1) or i == 0:
+				frdump = frdump0
 		# Otherwise frdump remains what it was before
 
 	ktstr = float(frdump.split('Vector:%s:%sCALIB_KAPPA_TST_REAL%s' % (ifo, chan_prefixes[config_index], chan_suffixes[config_index]))[1].split('mean ')[1].split(' ')[0])
@@ -647,6 +658,23 @@ for i in range(len(lines)):
 	#print("fs_over_Q0=%f" % fs_over_Q0)
 	fs = np.sqrt(complex(fs_squared))
 
+	TimesFile = open("%s_approxKappasTimes_%d-%d.txt" % (ifo, int(tstart), int(dursec)), "a")
+	TimesFile.write('%d\n' % int(t))
+	TDCF_file = open("%s_approxKappasTDCFs_%d-%d.txt" % (ifo, int(tstart), int(dursec)), "a")
+	TDCF_file.write('%d\tApproximate\tExact\n' % int(t))
+	TDCF_file.write('ktst\t\t%f\t%f\n' % (ktst, ktst1))
+	TDCF_file.write('tautst\t\t%e\t%e\n' % (tautst, tautst1))
+	TDCF_file.write('kpum\t\t%f\t%f\n' % (kpum, kpum1))
+	TDCF_file.write('taupum\t\t%e\t%e\n' % (taupum, taupum1))
+	TDCF_file.write('kuim\t\t%f\t%f\n' % (kuim, kuim1))
+	TDCF_file.write('tauuim\t\t%e\t%e\n' % (tauuim, tauuim1))
+	TDCF_file.write('kc\t\t%f\t%f\n' % (kc, kc0))
+	TDCF_file.write('fcc\t\t%f\t%f\n' % (fcc, fcc0))
+	TDCF_file.write('fs_squared\t%f\t%f\n' % (fs_squared, fs_squared0))
+	TDCF_file.write('fs_over_Q\t%f\t%f\n\n' % (np.real(fs * Qinv), fs_over_Q0))
+	TimesFile.close()
+	TDCF_file.close()
+
 	oldkc = kc0
 	oldfcc = fcc0
 	oldfs_squared = fs_squared0
@@ -660,6 +688,11 @@ for i in range(len(lines)):
 	MagArray[i] = 100 * (abs(Error) - 1)
 	PhaseArray[i] = np.angle(Error) * 180 / np.pi
 
+# Save data
+np.savetxt("%s_approxKappasErrorMagArray_%d-%d.txt" % (ifo, int(tstart), int(dursec)), MagArray)
+np.savetxt("%s_approxKappasErrorPhaseArray_%d-%d.txt" % (ifo, int(tstart), int(dursec)), PhaseArray)
+np.savetxt("%s_approxKappasErrorFrequency_%d-%d.txt" % (ifo, int(tstart), int(dursec)), fvec)
+
 # Make plots
 realtstart = datetime.datetime.utcnow().replace(tzinfo=datetime.timezone.utc).timestamp()
 fig, ax = plt.subplots(2 * nrows, 1)