Commit 8d9e0244 authored by Marcus Haberland's avatar Marcus Haberland
Browse files

progress in combination of measurements: added LISA_x_icegiant.py

parent 55a7b682
# -*- coding: utf-8 -*-
"""
Created on Fri Dec 17 13:24:49 2021
@author: marcu
"""
from scipy.constants import *
import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate, optimize
import LISA
import time
import astropy
import pickle
pi = np.pi
AU = 1495978707e2 # m
pc = 3.0857e16 # 2e5 AU
yr = 24*3600*365.25 # s
M_J = 1.898e27 # kg
M_S = 1.989e30 # kg
r_S = 2*G/c**2 # m
sqth = np.sqrt(3)/2
omega_E = 2*pi/yr
fig = 'Figures/'
arr_yr = np.linspace(0,yr,1000)
from binary import binary
from ig_binary import ig_binary
infile = open('dict_binaries.txt','rb')
binaries = pickle.load(infile)
infile.close()
infile = open('dict_binaries_ig.txt','rb')
binaries_ig = pickle.load(infile)
infile.close()
def isclose(a, b, rel_tol=1e-04, abs_tol=0.0):
return abs(a-b) <= max(rel_tol * max(abs(a), abs(b)), abs_tol)
def comparison(mode='s', weight=1e-4):
bins = binaries[mode]
ig_bins = binaries_ig[mode]
P = np.array([b['pars']['P'] for b in bins])
P_ig = np.array([b['pars']['P'] for b in ig_bins])
P_fin = []
F_fin = []
Tam_fin = []
for i, Pi in enumerate(P):
vals = np.array([bins[i]['pars']['K'],Pi*yr,1.])
for j, Pj in enumerate(P_ig):
if isclose(Pi,Pj):
P_fin.append(Pi)
A = ig_bins[j]['SNR']**2/2
B = bins[i]['SNR']**2/2*bins[i]['binary'].S_n()
F = bins[i]['Fisher']/bins[i]['SNR']**2/(1+weight*A/B)/bins[i]['pars']['M_P']**2 + weight*ig_bins[j]['Fisher']/2/(weight*A+B)/ig_bins[j]['pars']['M_P']**2
F_fin.append(F)
Tam_fin.append(np.sqrt(np.diag(np.linalg.inv(np.array(F_fin[-1]))/np.outer(vals,vals))))
#print('SNR_L = {:.3f}, SNR_IG = {:.3f}, total = {:.3f}\n'.format(bins[i]['SNR'],weight*bins[i]['SNR']+))
Tam_fin = np.array(Tam_fin)
plt.figure(dpi=300)
plt.tight_layout()
plt.title(r'LISA x IceGiant for $S_n^{LISA}(10 mHZ)/S_n^{IG}(10 mHz) = 10^{-4}$')
plt.loglog(P_fin,Tam_fin[:,0],'rx',label=r'$\sigma_K/K$')
plt.loglog(P_fin,Tam_fin[:,1],'bx',label=r'$\sigma_P/P$')
plt.loglog(P_fin,Tam_fin[:,2],'gx',label=r'$\sigma_\varphi$')
plt.xlabel(r'$P$ in yr')
plt.ylabel(r'$\sqrt{(\Gamma^{-1})_{ii}}/\lambda_i\cdot$SNR$\cdot M_P/M_J$')
plt.tight_layout()
plt.grid()
plt.savefig(fig+'rel_uncertainty_added.png')
comparison(weight=1e-4)
plt.show()
\ No newline at end of file
No preview for this file type
No preview for this file type
Supports Markdown
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