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 = 'bulge2'

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.1 Msun/pc^3*10^-3
popid= 3  3.0 Msun/pc^3*10^-3
popid= 4  4.7 Msun/pc^3*10^-3
popid= 5  4.3 Msun/pc^3*10^-3
popid= 6  6.8 Msun/pc^3*10^-3
popid= 7  1.4 Msun/pc^3*10^-3
popid= 8  27.6 Msun/pc^3*10^-6
## total sum of stars: 31.3
thindisc WD  12.3 Msun/pc^3*10^-3 3.8 mact
thick disc WD  3.3 Msun/pc^3*10^-3 1.1 mact
####################################
GDR2mock
nstar before cleaning: 1659283
nstar after g-band cleaning: 1283758
nstar before cleaning: 1283758
nstar after rp-band cleaning: 934427
934427
after smc cut: 934042
after lmc cut: 932193
GDR2
nstar before cleaning: 1253908
nstar after g-band cleaning: 1106123
nstar before cleaning: 1106123
nstar after rp-band cleaning: 930765
930765
after smc cut: 929364
after lmc cut: 919945
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)
likelihood:  -317486.79124657647
total diff:  108809
too many: 60527, too little: -48282
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
932193 919945
/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 = 'id4'
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: 1659283
nstar after g-band cleaning: 1283758
nstar before cleaning: 1283758
nstar after rp-band cleaning: 934427
934427
after smc cut: 934042
after lmc cut: 932193
GDR2
nstar before cleaning: 1666125
nstar after g-band cleaning: 1291838
nstar before cleaning: 1291838
nstar after rp-band cleaning: 939496
939496
after smc cut: 939123
after lmc cut: 937234
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 2344 is empty
hpx 2346 is empty
cumhpx 3103 is deleted in the model because of no representation in the data, (might affect likelihood)
likelihood:  -4819.276983472037
total diff:  15792
too many: 5374, too little: -10418
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 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 2344 is empty
hpx 2346 is empty
932193 937234
/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 [ ]: