diff --git a/python_sample/build_dipole_ksz_from_galaxies.py b/python_sample/build_dipole_ksz_from_galaxies.py index 4874136..d5ecb04 100644 --- a/python_sample/build_dipole_ksz_from_galaxies.py +++ b/python_sample/build_dipole_ksz_from_galaxies.py @@ -99,7 +99,7 @@ def build_unit_vectors(N): return ux,uy,uz -def generate_from_catalog(dmin,dmax,Nside,perturb=0.0,y=0.0,do_random=False): +def generate_from_catalog(dmin,dmax,Nside,perturb=0.0,y=0.0,do_random=False,do_hubble=False): import progressbar as pbar cat = np.load("2m++.npy") @@ -114,6 +114,8 @@ def generate_from_catalog(dmin,dmax,Nside,perturb=0.0,y=0.0,do_random=False): ksz_template = np.zeros(Npix, dtype=np.float64) ksz_mask = np.ones(Npix, dtype=np.uint8) + if do_hubble: + ksz_hubble_template = np.zeros(ksz_template.size, dtype=np.float64) for i in pbar.ProgressBar(maxval = cat.size, widgets=[pbar.Bar(), pbar.ETA()])(cat): if do_random: @@ -153,10 +155,17 @@ def generate_from_catalog(dmin,dmax,Nside,perturb=0.0,y=0.0,do_random=False): idx_masked = idx0[idx_masked] ksz_template[idx] += m ksz_mask[idx_masked] = 0 + if do_hubble: + ksz_hubble_template[idx] += m*DA ne.evaluate('ksz_template*ksz_normalization', out=ksz_template) + ne.evaluate('ksz_hubble_template*ksz_normalization', out=ksz_hubble_template) - return ksz_template, ksz_mask + result =ksz_template, ksz_mask + if do_hubble: + return result + ( ksz_hubble_template,) + else: + return result def get_args(): parser=argparse.ArgumentParser(description="Generate Skymaps from CIC maps") @@ -172,6 +181,7 @@ def get_args(): parser.add_argument('--y',type=float,default=0.0) parser.add_argument('--random', type=bool, default=False) parser.add_argument('--perturb', type=float, default=0) + parser.add_argument('--hubble_monopole', type=bool, default=False) return parser.parse_args() def main(): @@ -184,11 +194,19 @@ def main(): print("Generating map...") - proj,mask = generate_from_catalog(args.depth_min,args.depth_max,args.Nside,perturb=args.perturb,y=args.y,do_random=args.random) + r = generate_from_catalog(args.depth_min,args.depth_max,args.Nside,perturb=args.perturb,y=args.y,do_random=args.random,do_hubble=args.hubble_monopole) + hubble_map = None + if args.hubble_monopole: + proj,mask,hubble_map = r + else: + proj,mask = r if args.degrade > 0: proj *= mask proj = hp.ud_grade(proj, nside_out=args.degrade) + if hubble_map is not None: + hubble_map *= mask + hubble_map = hp.ud_grade(hubble_map, nside_out=args.degrade) mask = hp.ud_grade(mask, nside_out=args.degrade) Nside = args.degrade else: @@ -203,6 +221,9 @@ def main(): hp.write_map(args.ksz_map + "_y.fits", proj*y) hp.write_map(args.ksz_map + "_z.fits", proj*z) + if args.hubble_monopole: + hp.write_map(args.ksz_map + "_hubble.fits", hubble_map) + hp.mollview(proj*100*1e6, fig=1, coord='GG', cmap=plt.cm.coolwarm, title='', min=args.minval, max=args.maxval) diff --git a/python_sample/icgen/gen_ic_from_borg.py b/python_sample/icgen/gen_ic_from_borg.py index 8716b2f..09e0bbf 100644 --- a/python_sample/icgen/gen_ic_from_borg.py +++ b/python_sample/icgen/gen_ic_from_borg.py @@ -14,7 +14,7 @@ cosmo['omega_B_0']=0.049 cosmo['SIGMA8']=0.8344 cosmo['ns']=0.9624 -supergen=4 +supergen=1 zstart=99 astart=1/(1.+zstart) halfPixelShift=False diff --git a/src/loadGadget.cpp b/src/loadGadget.cpp index f24f495..cbd458d 100644 --- a/src/loadGadget.cpp +++ b/src/loadGadget.cpp @@ -67,7 +67,7 @@ void loadGadgetHeader(UnformattedRead *f, GadgetHeader& h, SimuData *data, int i data->Hubble = h.HubbleParam = f->readReal64(); f->endCheckpoint(true); - long NumPart = 0, NumPartTotal = 0; + ssize_t NumPart = 0, NumPartTotal = 0; for(int k=0; k<6; k++) { NumPart += h.npart[k]; @@ -112,7 +112,7 @@ SimuData *CosmoTool::loadGadgetMulti(const char *fname, int id, return 0; } - long NumPart = 0, NumPartTotal = 0; + ssize_t NumPart = 0, NumPartTotal = 0; try { diff --git a/src/loadSimu.hpp b/src/loadSimu.hpp index a80886c..b098351 100644 --- a/src/loadSimu.hpp +++ b/src/loadSimu.hpp @@ -74,8 +74,8 @@ namespace CosmoTool float Omega_M; float Omega_Lambda; - long NumPart; - long TotalNumPart; + ssize_t NumPart; + ssize_t TotalNumPart; long *Id; float *Pos[3]; float *Vel[3];