In [1]:
import numpy as np
from astropy.io import fits
import os, sys
path = os.path.abspath('../library/')
if path not in sys.path:
    sys.path.append(path)
from cumhpx import do_likelihood
name = 'final'

In [2]:
print('## Local normalisation, compare to Czekaj+14 table 7 ##')
outputDir = '../output/%s_100pc' %(name)
#outputDir = '../output/mockcat_new_SFR_new_IMF_100pc_more_thick'
x = fits.getdata(outputDir + '/' + 'GDR2mock_20.7Gmag.fits')
smaller_sphere = 50
sphere = (4/3)*np.pi*smaller_sphere**3
local_sum = 0
for popid in np.arange(9):
    selection = (x.popid==popid) & (x.logg < 7) & (np.divide(1,x.parallax)<smaller_sphere/1000)
    if (popid!=8):
        print('popid=',popid,' %.1f Msun/pc^3*10^-3' %(sum(x.smass[selection])/sphere*1000))
        local_sum+=sum(x.smass[selection])/sphere*1000
    else:
        print('popid=',popid,' %.1f Msun/pc^3*10^-6' %(sum(x.smass[selection])/sphere*1e6))
print("## total sum of stars: %.1f" %(local_sum))
selection = (x.popid<=6) & (x.logg > 7)& (np.divide(1,x.parallax)<smaller_sphere/1000)
print('thindisc WD',' %.1f Msun/pc^3*10^-3 %.1f mact' %(sum(x.smass[selection])/sphere*1000,sum(x.mact[selection])/sphere*1000))
selection = (x.popid==7) & (x.logg > 7)& (np.divide(1,x.parallax)<smaller_sphere/1000)
print('thick disc WD',' %.1f Msun/pc^3*10^-3 %.1f mact' %(sum(x.smass[selection])/sphere*1000,sum(x.mact[selection])/sphere*1000))
print('####################################')


outputDir = '../output/%s_0.001' %name
fSample = 0.001
do_likelihood(outputDir,fSample= fSample)


## Local normalisation, compare to Czekaj+14 table 7 ##
popid= 0  1.9 Msun/pc^3*10^-3
popid= 1  5.1 Msun/pc^3*10^-3
popid= 2  4.0 Msun/pc^3*10^-3
popid= 3  3.0 Msun/pc^3*10^-3
popid= 4  5.4 Msun/pc^3*10^-3
popid= 5  5.3 Msun/pc^3*10^-3
popid= 6  8.0 Msun/pc^3*10^-3
popid= 7  1.1 Msun/pc^3*10^-3
popid= 8  14.0 Msun/pc^3*10^-6
## total sum of stars: 34.0
thindisc WD  13.7 Msun/pc^3*10^-3 4.3 mact
thick disc WD  2.6 Msun/pc^3*10^-3 0.9 mact
####################################
GDR2mock
nstar before cleaning: 1686509
nstar after g-band cleaning: 1319834
nstar before cleaning: 1319834
nstar after rp-band cleaning: 970233
970233
after smc cut: 969846
after lmc cut: 967966
GDR2
nstar before cleaning: 1253908
nstar after g-band cleaning: 1106420
nstar before cleaning: 1106420
nstar after rp-band cleaning: 930639
930639
after smc cut: 929218
after lmc cut: 919688
hpx 2055 is empty
hpx 2064 is empty
hpx 2065 is empty
hpx 2066 is empty
hpx 2067 is empty
hpx 2068 is empty
hpx 2069 is empty
hpx 2070 is empty
hpx 2071 is empty
hpx 2072 is empty
hpx 2073 is empty
hpx 2076 is empty
hpx 2080 is empty
hpx 2081 is empty
hpx 2082 is empty
hpx 2088 is empty
hpx 2112 is empty
hpx 2336 is empty
hpx 2338 is empty
hpx 2344 is empty
hpx 2346 is empty
cumhpx 3096 is deleted in the model because of no representation in the data, (might affect likelihood)
cumhpx 3108 is deleted in the model because of no representation in the data, (might affect likelihood)
likelihood:  -317508.52247266035
total diff:  124646
too many: 86459, too little: -38187
hpx 2055 is empty
hpx 2064 is empty
hpx 2065 is empty
hpx 2066 is empty
hpx 2067 is empty
hpx 2068 is empty
hpx 2069 is empty
hpx 2070 is empty
hpx 2071 is empty
hpx 2072 is empty
hpx 2073 is empty
hpx 2076 is empty
hpx 2080 is empty
hpx 2081 is empty
hpx 2082 is empty
hpx 2088 is empty
hpx 2112 is empty
hpx 2336 is empty
hpx 2338 is empty
hpx 2344 is empty
hpx 2346 is empty
hpx 2055 is empty
hpx 2064 is empty
hpx 2065 is empty
hpx 2066 is empty
hpx 2067 is empty
hpx 2068 is empty
hpx 2069 is empty
hpx 2070 is empty
hpx 2071 is empty
hpx 2072 is empty
hpx 2073 is empty
hpx 2076 is empty
hpx 2080 is empty
hpx 2081 is empty
hpx 2082 is empty
hpx 2088 is empty
hpx 2112 is empty
hpx 2336 is empty
hpx 2338 is empty
hpx 2344 is empty
hpx 2346 is empty
967966 919688
/home/rybizki/anaconda3/lib/python3.6/site-packages/healpy/pixelfunc.py:304: RuntimeWarning: invalid value encountered in less_equal
  return np.absolute(m - badval) <= atol + rtol * np.absolute(badval)
/home/rybizki/anaconda3/lib/python3.6/site-packages/healpy/projaxes.py:1053: RuntimeWarning: invalid value encountered in less
  result.data[result.data<0]=0.0
/home/rybizki/anaconda3/lib/python3.6/site-packages/healpy/projaxes.py:1054: RuntimeWarning: invalid value encountered in greater
  result.data[result.data>1]=1.0
/home/rybizki/anaconda3/lib/python3.6/site-packages/matplotlib/colors.py:1031: RuntimeWarning: invalid value encountered in less_equal
  mask |= resdat <= 0

In [3]:
name_comp = 'id5'
compDir = '../output/%s_0.001/GDR2mock_20.7Gmag.fits' %name_comp
fSample = 0.0025
do_likelihood(outputDir,fSample= fSample, comp_file = compDir, savefile=False)


GDR2mock
nstar before cleaning: 1686509
nstar after g-band cleaning: 1319834
nstar before cleaning: 1319834
nstar after rp-band cleaning: 970233
970233
after smc cut: 969846
after lmc cut: 967966
GDR2
nstar before cleaning: 1589492
nstar after g-band cleaning: 1237075
nstar before cleaning: 1237075
nstar after rp-band cleaning: 907793
907793
after smc cut: 907458
after lmc cut: 905764
hpx 2055 is empty
hpx 2064 is empty
hpx 2065 is empty
hpx 2066 is empty
hpx 2067 is empty
hpx 2068 is empty
hpx 2069 is empty
hpx 2070 is empty
hpx 2071 is empty
hpx 2072 is empty
hpx 2073 is empty
hpx 2080 is empty
hpx 2081 is empty
hpx 2082 is empty
hpx 2083 is empty
hpx 2088 is empty
hpx 2112 is empty
hpx 2336 is empty
hpx 2338 is empty
hpx 2339 is empty
hpx 2344 is empty
hpx 2346 is empty
cumhpx 3103 is deleted in the model because of no representation in the data, (might affect likelihood)
cumhpx 3108 is deleted in the model because of no representation in the data, (might affect likelihood)
likelihood:  -4543.584325346532
total diff:  62277
too many: 62237, too little: -40
hpx 2055 is empty
hpx 2064 is empty
hpx 2065 is empty
hpx 2066 is empty
hpx 2067 is empty
hpx 2068 is empty
hpx 2069 is empty
hpx 2070 is empty
hpx 2071 is empty
hpx 2072 is empty
hpx 2073 is empty
hpx 2080 is empty
hpx 2081 is empty
hpx 2082 is empty
hpx 2083 is empty
hpx 2088 is empty
hpx 2112 is empty
hpx 2336 is empty
hpx 2338 is empty
hpx 2339 is empty
hpx 2344 is empty
hpx 2346 is empty
hpx 2055 is empty
hpx 2064 is empty
hpx 2065 is empty
hpx 2066 is empty
hpx 2067 is empty
hpx 2068 is empty
hpx 2069 is empty
hpx 2070 is empty
hpx 2071 is empty
hpx 2072 is empty
hpx 2073 is empty
hpx 2080 is empty
hpx 2081 is empty
hpx 2082 is empty
hpx 2083 is empty
hpx 2088 is empty
hpx 2112 is empty
hpx 2336 is empty
hpx 2338 is empty
hpx 2339 is empty
hpx 2344 is empty
hpx 2346 is empty
967966 905764
/home/rybizki/anaconda3/lib/python3.6/site-packages/healpy/pixelfunc.py:304: RuntimeWarning: invalid value encountered in less_equal
  return np.absolute(m - badval) <= atol + rtol * np.absolute(badval)
/home/rybizki/anaconda3/lib/python3.6/site-packages/healpy/projaxes.py:1053: RuntimeWarning: invalid value encountered in less
  result.data[result.data<0]=0.0
/home/rybizki/anaconda3/lib/python3.6/site-packages/healpy/projaxes.py:1054: RuntimeWarning: invalid value encountered in greater
  result.data[result.data>1]=1.0
/home/rybizki/anaconda3/lib/python3.6/site-packages/matplotlib/colors.py:1031: RuntimeWarning: invalid value encountered in less_equal
  mask |= resdat <= 0

In [ ]: