In [1]:

    

import matplotlib.pyplot as plt
import numpy as np
import os
from scipy import stats
import glob
from scipy.stats import ks_2samp, kstest
%matplotlib inline


    

In [2]:

    

def load_summary(filename):
    dtype=[('minr', 'f8'),
           ('maxr', 'f8'), 
           ('ca_ratio', 'f8'),
           ('ba_ratio', 'f8'),
           ('a', 'f8'),
           ('center', 'f8'),
           ('width', 'f8'),
           ('mu', 'f8')]
    summary = np.loadtxt(filename, dtype=dtype)    
    return summary


    

In [3]:

    

def load_experiment(input_path="../data/mstar_selected_summary/vmax_sorted/", fixed_number=False, full_data=False):
    files = glob.glob(input_path+"M31_group_*")
    group_id = []
    for f in files:
        i = int(f.split("_")[-5])
        if i not in group_id:
            group_id.append(i)
    print(group_id, len(group_id))

    n_groups = len(group_id)

    if fixed_number:
        n_iter = np.arange(5)
    else:
        n_iter = np.arange(11,16)
    
    fields = ['width','mu', 'a', 'ba_ratio', 'ca_ratio']
    M31_all = {}
    MW_all = {}
    if not full_data:
        for field in fields:
            M31_all[field] = np.ones(n_groups)
            MW_all[field] = np.ones(n_groups)
            M31_all[field+'_sigma'] = np.ones(n_groups)
            MW_all[field+'_sigma'] = np.ones(n_groups)
        
            M31_all[field+'_random'] = np.ones(n_groups)
            MW_all[field+'_random'] = np.ones(n_groups)
            M31_all[field+'_random_sigma'] = np.ones(n_groups)
            MW_all[field+'_random_sigma'] = np.ones(n_groups)
    else:
        for field in fields:
            M31_all[field] = np.empty((0))
            MW_all[field] = np.empty((0))
            M31_all[field+'_random'] = np.empty((0))
            MW_all[field+'_random'] = np.empty((0))
           

    for g in range(n_groups):

        MW_summary = {}
        M31_summary = {}
    
        for i in n_iter:
            if fixed_number:
                filename_MW = os.path.join(input_path,"MW_group_{}_nmax_{}_iter_{}.dat".format(group_id[g], 11, i))
                filename_M31 = os.path.join(input_path,"M31_group_{}_nmax_{}_iter_{}.dat".format(group_id[g],11, i))
            else:
                filename_MW = os.path.join(input_path,"MW_group_{}_nmax_{}_iter_{}.dat".format(group_id[g], i, 0))
                filename_M31 = os.path.join(input_path,"M31_group_{}_nmax_{}_iter_{}.dat".format(group_id[g], i, 0))

            MW_summary[i] = load_summary(filename_MW)
            M31_summary[i] = load_summary(filename_M31)
    
        
        for field in fields:
            a = np.empty((0))
            b = np.empty((0))
            a_random = np.empty((0))
            b_random = np.empty((0))
        
            for i in n_iter:
                data = M31_summary[i]
                a = np.append(a, data[field][0])
                a_random = np.append(a_random, data[field][1:101])
        
                data = MW_summary[i]
                b = np.append(b, data[field][0])
                b_random = np.append(b_random, data[field][1:101])
                #print('a_random {} iter: {} {}'.format(field, i, a_random))
                
            if not full_data:
                M31_all[field][g] = np.average(a)
                MW_all[field][g] = np.average(b)
                M31_all[field+'_sigma'][g] = np.std(a)
                MW_all[field+'_sigma'][g] = np.std(b)
                M31_all[field+'_random'][g] = np.average(a_random)
                MW_all[field+'_random'][g] = np.average(b_random)
                M31_all[field+'_random_sigma'][g] = np.std(a_random)
                MW_all[field+'_random_sigma'][g] = np.std(b_random)
            else:
                M31_all[field] = np.append(M31_all[field], a)
                MW_all[field] = np.append(MW_all[field], b)
                M31_all[field+'_random'] = np.append(M31_all[field+'_random'], a_random)
                MW_all[field+'_random'] = np.append(MW_all[field+'_random'], b_random)
                
    return M31_all, MW_all


    

In [4]:

    

in_path = "../data/obs_summary/vmag_sorted/"
M31_obs_vmag_sorted, MW_obs_vmag_sorted = load_experiment(input_path=in_path, fixed_number=False, full_data=False)

in_path = "../data/illustris1_mstar_selected_summary/vmax_sorted/"
M31_sim_vmax_sorted_illu, MW_sim_vmax_sorted_illu = load_experiment(input_path=in_path, fixed_number=False)

in_path = "../data/illustris1dark_mstar_selected_summary/vmax_sorted/"
M31_sim_vmax_sorted_illudm, MW_sim_vmax_sorted_illudm = load_experiment(input_path=in_path, fixed_number=False)

in_path = "../data/elvis_mstar_selected_summary/vmax_sorted/"
M31_sim_vmax_sorted_elvis, MW_sim_vmax_sorted_elvis = load_experiment(input_path=in_path, fixed_number=False)


    

    




[0] 1
[0, 10, 11, 13, 14, 16, 18, 1, 20, 21, 22, 24, 25, 2, 3, 4, 5, 6, 8, 9] 20
[0, 10, 11, 12, 13, 14, 15, 16, 17, 18, 1, 20, 21, 22, 23, 25, 2, 3, 4, 5, 6, 7, 8, 9] 24
[0, 10, 11, 1, 2, 3, 4, 5, 6, 7, 8, 9] 12

In [5]:

    

print("M31 observations \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, M31_obs_vmag_sorted[field][0], M31_obs_vmag_sorted[field+'_sigma'][0])
    
    
print("\nMW observations \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, MW_obs_vmag_sorted[field][0], MW_obs_vmag_sorted[field+'_sigma'][0])


    

    




M31 observations 

Natural units
 width 59.06998 3.00637879742
Natural units
 ca_ratio 0.4514 0.0409125897494
Natural units
 ba_ratio 0.8226 0.0680370487308

MW observations 

Natural units
 width 21.27636 2.2566692603
Natural units
 ca_ratio 0.2564 0.0501382089828
Natural units
 ba_ratio 0.8064 0.040262141026

In [6]:

    

print("M31 observations \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:    
    normed_mean = (M31_obs_vmag_sorted[field][0] - M31_obs_vmag_sorted[field+'_random'][0])/M31_obs_vmag_sorted[field+'_random_sigma'][0]
    normed_sigma = M31_obs_vmag_sorted[field+'_sigma'][0]/M31_obs_vmag_sorted[field+'_random_sigma'][0]

    print("Normalized units\n", field, normed_mean, normed_sigma)

    
print("\nMW observations \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    
    normed_mean = (MW_obs_vmag_sorted[field][0] - MW_obs_vmag_sorted[field+'_random'][0])/MW_obs_vmag_sorted[field+'_random_sigma'][0]
    normed_sigma = MW_obs_vmag_sorted[field+'_sigma'][0]/MW_obs_vmag_sorted[field+'_random_sigma'][0]

    print("Normalized units\n", field, normed_mean, normed_sigma)


    

    




M31 observations 

Normalized units
 width -0.485012732811 0.243816540616
Normalized units
 ca_ratio -1.03538999178 0.375236832
Normalized units
 ba_ratio -0.0197275054087 0.824447940313

MW observations 

Normalized units
 width -2.48868439056 0.26020821071
Normalized units
 ca_ratio -2.18257393327 0.429985100027
Normalized units
 ba_ratio -0.132487764218 0.472223888737

In [7]:

    

print("M31 observations (spherically randomized)\n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, M31_obs_vmag_sorted[field+'_random'][0], M31_obs_vmag_sorted[field+'_random_sigma'][0])
    
print("\nMW observations (spherically randomized)\n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, MW_obs_vmag_sorted[field+'_random'][0], MW_obs_vmag_sorted[field+'_random_sigma'][0])


    

    




M31 observations (spherically randomized)

Natural units
 width 65.0504274 12.3304956662
Natural units
 ca_ratio 0.56429 0.109031380345
Natural units
 ba_ratio 0.824228 0.0825243722545

MW observations (spherically randomized)

Natural units
 width 42.859605 8.6725520849
Natural units
 ca_ratio 0.510898 0.116604526482
Natural units
 ba_ratio 0.817696 0.085260703633

In [8]:

    

print("M31 illustris simulation \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, 
          np.mean(M31_sim_vmax_sorted_illu[field]), np.std(M31_sim_vmax_sorted_illu[field+'_sigma']))
    
    
print("\nMW illustris simulation \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, 
          np.mean(MW_sim_vmax_sorted_illu[field]), np.std(MW_sim_vmax_sorted_illu[field+'_sigma']))


    

    




M31 illustris simulation 

Natural units
 width 71.403784 4.4144513083
Natural units
 ca_ratio 0.54648 0.0183691787775
Natural units
 ba_ratio 0.79936 0.019404234985

MW illustris simulation 

Natural units
 width 66.462873 1.98695143337
Natural units
 ca_ratio 0.5136 0.0125069647986
Natural units
 ba_ratio 0.81477 0.0199655504824

In [9]:

    

print("M31 illustris simulation DM \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, 
          np.mean(M31_sim_vmax_sorted_illudm[field]), np.std(M31_sim_vmax_sorted_illudm[field+'_sigma']))
    
    
print("\nMW illustris simulation DM \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, 
          np.mean(MW_sim_vmax_sorted_illudm[field]), np.std(MW_sim_vmax_sorted_illudm[field+'_sigma']))


    

    




M31 illustris simulation DM 

Natural units
 width 61.1427566667 3.25311258406
Natural units
 ca_ratio 0.504316666667 0.0231184691878
Natural units
 ba_ratio 0.810625 0.0183390977986

MW illustris simulation DM 

Natural units
 width 62.0739975 3.14903679161
Natural units
 ca_ratio 0.503858333333 0.0172083642044
Natural units
 ba_ratio 0.799941666667 0.0153044751527

In [10]:

    

print("M31 elvis simulation \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, 
          np.mean(M31_sim_vmax_sorted_elvis[field]), np.std(M31_sim_vmax_sorted_elvis[field+'_sigma']))
    
    
print("\nMW elvis simulation \n")
fields = ['width', 'ca_ratio', 'ba_ratio']
for field in fields:
    print("Natural units\n", field, 
          np.mean(MW_sim_vmax_sorted_elvis[field]), np.std(MW_sim_vmax_sorted_elvis[field+'_sigma']))


    

    




M31 elvis simulation 

Natural units
 width 70.61085 1.96390149498
Natural units
 ca_ratio 0.54725 0.017208460608
Natural units
 ba_ratio 0.80925 0.019973160614

MW elvis simulation 

Natural units
 width 67.37044 4.11301024699
Natural units
 ca_ratio 0.494483333333 0.0248056931154
Natural units
 ba_ratio 0.809066666667 0.0158756410076