set_injection_parameters.py 3.34 KB
Newer Older
Marc Arene's avatar
Marc Arene committed
1 2 3
import sys
# cf https://docs.python.org/3/tutorial/modules.html
sys.path.append('../')
4 5 6
import numpy as np
import bilby.gw.conversion as conv
from Headers.PN_Coefficients import * # Headers.Constants already imported in PN_Coeff
Marc Arene's avatar
Marc Arene committed
7 8
# import json
from configparser import ConfigParser
9
import Library.python_utils as pu
Marc Arene's avatar
Marc Arene committed
10

11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54


def ComputeChirpTime3p5PN(f_low, m1, m2):
    """
    Computes chirp times at 3.5 PN.

    Inputs:
    -------
        f_low: [float] low frequency cutoff of the detector
        sigpar: [array_like] parameters of the system
    Outputs:
    --------

    """
    Mt = m1 + m2
    Mt_sec = Mt * GEOM           # Total mass in s
    etaM = m1 * m2 / Mt**2
    etaM2 = etaM*etaM

    Tau_0 = tf1 * Mt_sec / etaM
    Tau_2 = tf2 + tf3 * etaM
    Tau_3 = tf4 * PI
    Tau_4 = tf5 + tf6 * etaM + tf7 * etaM2
    Tau_5 = PI*(tf9 + tf8*etaM)
    Tau_6 = tf11 + tf12*PI2 + tf10*(EULER + np.log(4.0)) + ( tf13 + tf14*PI2 )*etaM + tf15*etaM2 + tf16*etaM2*etaM
    Tau_6_log = tf10
    Tau_7 = PI*(tf19 + tf18*etaM + tf17*etaM2 )
    v_0 = pow(PI*Mt_sec*f_low,1.0/3.0)
    v_0_2 = v_0*v_0
    v_0_3 = v_0_2*v_0
    v_0_4 = v_0_3*v_0
    v_0_5 = v_0_4*v_0
    v_0_6 = v_0_5*v_0
    v_0_7 = v_0_6*v_0

    Tau_chirp = Tau_0*pow(v_0,-8.0)*(1.0 + Tau_2*v_0_2 + Tau_3*v_0_3 + Tau_4 *v_0_4 +  Tau_5 *v_0_5 +  (Tau_6 + Tau_6_log*np.log(v_0))*v_0_6 +  Tau_7 *v_0_7 )
    return Tau_chirp

def Compute_fmax(Rmin, m1, m2):
    vmax = np.sqrt(1.0 / Rmin)
    Mt = m1 + m2
    fmax = vmax**3 / (PI * GEOM * Mt)      # fmax signal
    return fmax

55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98
def ini_file_to_dict(file_path='../examples/GW170817.ini'):

    config = ConfigParser()
    config.read(file_path)
    injection_parameters = pu.config_parser_to_dict(config)['parameters']


    injection_parameters = conv.generate_mass_parameters(injection_parameters)
    injection_parameters['reduced_mass'] = (injection_parameters['mass_1'] * injection_parameters['mass_2']) / (injection_parameters['mass_1'] + injection_parameters['mass_2'])


    # Set the duration and sampling frequency of the data segment that we're going
    # to inject the signal into. For the
    # TaylorF2 waveform, we cut the signal close to the isco frequency
    minimum_frequency=30.0
    m1 = injection_parameters['mass_1']
    m2 = injection_parameters['mass_2']
    tc_3p5PN = ComputeChirpTime3p5PN(minimum_frequency, m1, m2)
    # We integrate the signal up to the frequency of the "Innermost stable circular orbit (ISCO)"
    R_isco = 6.      # Orbital separation at ISCO, in geometric units. 6M for PN ISCO; 2.8M for EOB
    m1_min = 1
    m2_min = 1
    f_high = Compute_fmax(R_isco, m1, m2)
    f_high_run = Compute_fmax(R_isco, m1_min, m2_min)
    # import IPython; IPython.embed();
    duration = tc_3p5PN + 2
    # sampling_frequency = 4096
    # sampling_frequency = 4397
    sampling_frequency = 2 * f_high_run
    start_time = injection_parameters['geocent_time'] - tc_3p5PN

    meta_params = dict(
        minimum_frequency = minimum_frequency,
        tc_3p5PN = tc_3p5PN,
        f_high = f_high,
        f_high_run = f_high_run,
        duration = duration,
        sampling_frequency = sampling_frequency,
        start_time = start_time
    )

    injection_parameters['meta'] = meta_params

    return injection_parameters
99 100 101 102 103 104 105 106 107




if __name__=='__main__':

    injection_parameters = ini_file_to_dict()

    pu.print_dict(injection_parameters, indent=3, align_keys=True)