Commit 48ddbb7e authored by Christian Chapman-Bird's avatar Christian Chapman-Bird
Browse files

updated populations

parent 82db3bc3
Loading
Loading
Loading
Loading
+45 −15
Original line number Diff line number Diff line
@@ -152,7 +152,44 @@ def convert_M(M,z, to_source=False):
    else:
        return M*(1+z)

def draw_population_params(num, M_min=5e4, M_max=5e7, M_lam=2, z_max=20, spin_mu=0.98, spin_std=0.01, e_min=0.01,e_max=0.3):
# def draw_population_params(num, M_min=5e4, M_max=5e7, M_lam=2, z_max=20, spin_mu=0.98, spin_std=0.01, e_min=0.01,e_max=0.3):
#     '''
#     Draw sets of EMRI parameters from population distributions.
#     :param num: Number of events to draw.
#     :param M_min: Minimum mass for mass population.
#     :param M_max: Maximum mass for mass population.
#     :param M_lam: Spectral index for mass population power law.
#     :param z_max: Maximum redshift for redshift population.
#     :param spin_mu: Mean for spin population.
#     :param spin_std: Standard deviation for spin population.
#     :return: Numpy array of shape (num, 5) that consists of the rows (M, a, d_L, z, M_s)
#     '''
#
#     #initialise populations
#
#     mpop = mass_pop(M_min=M_min, M_max=M_max, lam=M_lam)
#     redpop = redshift_pop(z_max=z_max)
#
#     # spin population with scipy truncnorm
#     l1, l2 = (0 - spin_mu) / spin_std, (1 - spin_mu) / spin_std  # rescale the limits (from scipy.stats.truncnorm docs)
#     spinpop = truncnorm(l1, l2, loc=spin_mu,scale=spin_std)
#
#     masses = mpop.rvs(size=num)
#     redshifts = redpop.rvs(size=num)
#     spins = spinpop.rvs(size=num)
#
#     lum_dists = np.zeros(num)
#     for i, z in enumerate(redshifts):
#         lum_dists[i] = dL(z)
#     red_masses = convert_M(masses, redshifts, to_source=False)
#
#     eccentricities = np.random.uniform(e_min, e_max, num)
#
#     out = np.column_stack((red_masses, spins, lum_dists, redshifts, masses, eccentricities))
#     return out


def draw_population_params(num, M_min=5e4, M_max=5e7, M_lam=2, spin_mu=0.98, spin_std=0.01, e_min=0.01,e_max=0.3):
    '''
    Draw sets of EMRI parameters from population distributions.
    :param num: Number of events to draw.
@@ -168,24 +205,16 @@ def draw_population_params(num, M_min=5e4, M_max=5e7, M_lam=2, z_max=20, spin_mu
    #initialise populations

    mpop = mass_pop(M_min=M_min, M_max=M_max, lam=M_lam)
    redpop = redshift_pop(z_max=z_max)

    # spin population with scipy truncnorm
    l1, l2 = (0 - spin_mu) / spin_std, (1 - spin_mu) / spin_std  # rescale the limits (from scipy.stats.truncnorm docs)
    spinpop = truncnorm(l1, l2, loc=spin_mu,scale=spin_std)

    masses = mpop.rvs(size=num)
    redshifts = redpop.rvs(size=num)
    second_ms = 10**np.random.uniform(-0.3,2,num)
    spins = spinpop.rvs(size=num)

    lum_dists = np.zeros(num)
    for i, z in enumerate(redshifts):
        lum_dists[i] = dL(z)
    red_masses = convert_M(masses, redshifts, to_source=False)

    eccentricities = np.random.uniform(e_min, e_max, num)

    out = np.column_stack((red_masses, spins, lum_dists, redshifts, masses, eccentricities))
    out = np.column_stack((masses, second_ms, spins, eccentricities))
    return out


@@ -193,10 +222,9 @@ def draw_other_params(num):
    qS = np.arccos(2*np.random.uniform(0,1,size=num)-1)
    qK = np.arccos(2*np.random.uniform(0,1,size=num)-1)
    phiS = np.random.uniform(0,2*np.pi,size=num)
    Y0 = np.cos(np.arcsin(2*np.random.uniform(0,1,size=num)-1))*0.995  # YLIM = 0.998
    Y0 = np.cos(np.arcsin(np.random.uniform(0,np.sin(np.arccos(0.5/0.99)),size=num)))*0.99  # YLIM = 0.998
    t = np.random.uniform(0,4,size=num)
    out = np.column_stack((Y0,qS,phiS,qK,t))

    return out


@@ -209,5 +237,7 @@ def draw_other_params(num):
#     def _pdf(self,e):
#         return (self.b - e)**(11*(e-self.b))*np.exp(11*(e-self.b))

temp = draw_population_params(100)
print(temp)
 No newline at end of file
# M, mu, a, e = draw_population_params(1)[0]
# Y0, qS, phiS, qK, tPl = draw_other_params(1)[0]
#
# print(M, mu, a, e, Y0, qS, phiS, qK, tPl)
 No newline at end of file