Select Git revision
test_rad_pressure.py
Forked from
finesse / pykat
Source project has a limited visibility.
-
Daniel Brown authored
adding in more support for nodes. Can now access any node using kat.nodes.n1. Added in functions for Node object to set the gauss parameters. Only added gauss for w0 and z so far.
Daniel Brown authoredadding in more support for nodes. Can now access any node using kat.nodes.n1. Added in functions for Node object to set the gauss parameters. Only added gauss for w0 and z so far.
plot_data.py 2.69 KiB
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
import os
from tqdm import tqdm
from oct2py import octave
import glob
filenames = glob.glob("CollectedOutput/*.txt")
plt.style.use('paper')
Tspan = 100 * 86400
def Recovery(Tspan, Depth, twoFstar=60, detectors='H1,L1'):
numDetectors = len(detectors.split(','))
cmd = ("DetectionProbabilityStackSlide('Nseg', 1, 'Tdata', {},"
"'misHist', createDeltaHist(0), 'avg2Fth', {}, 'detectors', '{}',"
"'Depth', {})"
).format(numDetectors*Tspan, twoFstar, detectors, Depth)
return octave.eval(cmd, verbose=False)
def binomialConfidenceInterval(N, K, confidence=0.95):
cmd = '[fLow, fUpper] = binomialConfidenceInterval({}, {}, {})'.format(
N, K, confidence)
[l, u] = octave.eval(cmd, verbose=False, return_both=True)[0].split('\n')
return float(l.split('=')[1]), float(u.split('=')[1])
df_list = []
for fn in filenames:
df = pd.read_csv(
fn, sep=' ', names=['depth', 'h0', 'dF0', 'dF1', 'twoF', 'runTime'])
df['CLUSTER_ID'] = fn.split('_')[1]
df_list.append(df)
df = pd.concat(df_list)
twoFstar = 70
depths = np.unique(df.depth.values)
recovery_fraction = []
recovery_fraction_CI = []
for d in depths:
twoFs = df[df.depth == d].twoF.values
N = len(twoFs)
K = np.sum(twoFs > twoFstar)
print d, N, K
recovery_fraction.append(K/float(N))
[fLower, fUpper] = binomialConfidenceInterval(N, K)
recovery_fraction_CI.append([fLower, fUpper])
yerr = np.abs(recovery_fraction - np.array(recovery_fraction_CI).T)
fig, ax = plt.subplots()
ax.errorbar(depths, recovery_fraction, yerr=yerr, fmt='sr', marker='s', ms=2,
capsize=1, capthick=0.5, elinewidth=0.5,
label='Monte-Carlo result', zorder=10)
fname = 'analytic_data_{}.txt'.format(twoFstar)
if os.path.isfile(fname):
depths_smooth, recovery_analytic = np.loadtxt(fname)
else:
depths_smooth = np.linspace(10, 550, 100)
recovery_analytic = []
for d in tqdm(depths_smooth):
recovery_analytic.append(Recovery(Tspan, d, twoFstar))
np.savetxt(fname, np.array([depths_smooth, recovery_analytic]))
depths_smooth = np.concatenate(([0], depths_smooth))
recovery_analytic = np.concatenate(([1], recovery_analytic))
ax.plot(depths_smooth, recovery_analytic, '-k', label='Theoretical maximum')
ax.set_ylim(0, 1.05)
ax.set_xlabel(r'Signal depth', size=10)
ax.set_ylabel(r'Recovered fraction', size=10)
ax.legend(loc=1, frameon=False)
fig.tight_layout()
fig.savefig('allsky_recovery.png')
total_number_steps = 6 * 50.
fig, ax = plt.subplots()
ax.hist(df.runTime/total_number_steps, bins=50)
ax.set_xlabel('run-time per step [s]')
fig.tight_layout()
fig.savefig('runTimeHist.png')