From e5cc6297fa4f190a80ff36eaf543a8e11021ecc7 Mon Sep 17 00:00:00 2001 From: "P.M. Sutter" Date: Wed, 12 Feb 2014 04:52:58 -0600 Subject: [PATCH] consistent handling of comoving coords in sims --- c_tools/stacking/pruneVoids.cpp | 9 ++- c_tools/stacking/pruneVoids.ggo | 2 + .../pipeline_source/prepareCatalogs.in.py | 4 +- .../void_python_tools/backend/classes.py | 9 +-- .../void_python_tools/backend/launchers.py | 55 +++++++++---------- 5 files changed, 43 insertions(+), 36 deletions(-) diff --git a/c_tools/stacking/pruneVoids.cpp b/c_tools/stacking/pruneVoids.cpp index 6e54cc6..6a1bd20 100644 --- a/c_tools/stacking/pruneVoids.cpp +++ b/c_tools/stacking/pruneVoids.cpp @@ -163,7 +163,7 @@ int main(int argc, char **argv) { gsl_function expanF; expanF.function = &expanFun; struct my_expan_params expanParams; - expanParams.Om = 0.27; + expanParams.Om = args.omegaM_arg; expanParams.w0 = -1.0; expanParams.wa = 0.0; expanF.params = &expanParams; @@ -741,7 +741,12 @@ int main(int argc, char **argv) { } else { voids[iVoid].redshiftInMpc = voids[iVoid].barycenter[2]; - voids[iVoid].redshift = voids[iVoid].barycenter[2]/LIGHT_SPEED*100.; + if (args.useComoving_flag) { + voids[iVoid].redshift = gsl_interp_eval(interp, dL, redshifts, + voids[iVoid].redshiftInMpc, acc); + } else { + voids[iVoid].redshift = voids[iVoid].barycenter[2]/LIGHT_SPEED*100.; + } nearestEdge = 1.e99; diff --git a/c_tools/stacking/pruneVoids.ggo b/c_tools/stacking/pruneVoids.ggo index d7b84ee..77f4521 100644 --- a/c_tools/stacking/pruneVoids.ggo +++ b/c_tools/stacking/pruneVoids.ggo @@ -25,6 +25,8 @@ option "isObservation" - "We are working with observational data" flag off option "useComoving" - "Void positions are in comoving coordinates" flag off +option "omegaM" - "Omega_M for redshift convertion" double optional default="0.27" + option "zMin" - "Minimum redshift of sample" double optional default="0.0" option "zMax" - "Maximum redshift of sample" double optional default="10.0" diff --git a/python_tools/pipeline_source/prepareCatalogs.in.py b/python_tools/pipeline_source/prepareCatalogs.in.py index 1a89846..6242572 100644 --- a/python_tools/pipeline_source/prepareCatalogs.in.py +++ b/python_tools/pipeline_source/prepareCatalogs.in.py @@ -230,7 +230,8 @@ newSample = Sample(dataFile = "{dataFile}", profileBinSize = "auto", includeInHubble = True, partOfCombo = False, - {autoStack} + {autoStack}, + numAPSlices = {numAPSlices}, isCombo = False, boxLen = {boxLen}, usePecVel = {usePecVel}, @@ -372,6 +373,7 @@ newSample.addStack({zMin}, {zMax}, 2*{minRadius}+18, 2*{minRadius}+24, True, Fal omegaM=Om, boxLen=lbox, autoStack=autoStack, + numAPSlices=numAPSlices, usePecVel=useVel, minRadius=minRadius, numSubvolumes=numSubvolumes, diff --git a/python_tools/void_python_tools/backend/classes.py b/python_tools/void_python_tools/backend/classes.py index 5d827aa..4d4b21c 100644 --- a/python_tools/void_python_tools/backend/classes.py +++ b/python_tools/void_python_tools/backend/classes.py @@ -76,11 +76,11 @@ class Sample: profileBinSize = 2 # Mpc autoNumInStack = -1 # set to >0 to automatically generate stacks of size N autoPartInStack = -1 # set to >0 to automatically generate stacks with N particles + numAPSlices = 1 volumeLimited = True includeInHubble = True partOfCombo = False isCombo = False - useComoving = False # if True, convert population to comoving coordinates comboList = [] @@ -88,7 +88,6 @@ class Sample: boxLen = 1024 # Mpc/h usePecVel = False subsample = 1.0 - useLightCone = True numSubvolumes = 1 mySubvolume = 1 @@ -98,13 +97,14 @@ class Sample: nickName="", maskFile="", selFunFile="", zBoundary=(), zRange=(), zBoundaryMpc=(), minVoidRadius=0, fakeDensity=0.01, volumeLimited=True, + numAPSlices=1, includeInHubble=True, partOfCombo=False, isCombo=False, comboList=(), profileBinSize=2.0, skyFraction=0.19, boxLen=1024, usePecVel=False, omegaM=0.27, numSubvolumes=1, mySubvolume=1, dataFormat="sdss", - useComoving=False, + useComoving=True, dataType="observation", - subsample=1.0, useLightCone=True, autoNumInStack=-1, + subsample=1.0, useLightCone=False, autoNumInStack=-1, autoPartInStack=-1): self.dataFile = dataFile self.fullName = fullName @@ -121,6 +121,7 @@ class Sample: self.partOfCombo = partOfCombo self.isCombo = isCombo self.comboList = comboList + self.numAPSlices = numAPSlices self.zobovDir = None self.profileBinSize = profileBinSize self.skyFraction = skyFraction diff --git a/python_tools/void_python_tools/backend/launchers.py b/python_tools/void_python_tools/backend/launchers.py index 9d9b393..7648258 100644 --- a/python_tools/void_python_tools/backend/launchers.py +++ b/python_tools/void_python_tools/backend/launchers.py @@ -459,6 +459,7 @@ def launchPrune(sample, binPath, cmd += observationLine cmd += periodicLine cmd += useComovingFlag + cmd += " --omegaM=" + str(sample.omegaM) cmd += " --outputDir=" + zobovDir cmd += " --sampleName=" + str(sampleName) #cmd += " &> " + logFile @@ -594,8 +595,8 @@ def launchStack(sample, stack, binPath, thisDataPortion=None, logDir=None, if useComoving or not useLightCone: zMin = LIGHT_SPEED/100.*vp.angularDiameter(zMin, Om=omegaM) zMax = LIGHT_SPEED/100.*vp.angularDiameter(zMax, Om=omegaM) - print min(sample.zRange[1],stack.zMax)*LIGHT_SPEED/100., zMax - print max(sample.zRange[0],stack.zMin)*LIGHT_SPEED/100., zMin + print min(sample.zRange[1],stack.zMax)*LIGHT_SPEED/100., zMax, sample.zBoundaryMpc[1] + print max(sample.zRange[0],stack.zMin)*LIGHT_SPEED/100., zMin, sample.zBoundaryMpc[0] else: zMin *= LIGHT_SPEED/100. zMax *= LIGHT_SPEED/100. @@ -1185,7 +1186,8 @@ def launchProfile(sample, stack, voidDir=None, logFile=None, continueRun=None, # ----------------------------------------------------------------------------- -def launchFit(sample, stack, logFile=None, voidDir=None, figDir=None, +def launchFit(sample, stack, useComoving=None, + logFile=None, voidDir=None, figDir=None, continueRun=None, thisDataPortion=None): sampleName = sample.fullName @@ -1259,6 +1261,11 @@ def launchFit(sample, stack, logFile=None, voidDir=None, figDir=None, ret,fits,args = vp.compute_inertia(voidDir, stack.rMax, mode="2d", nBootstraps=500, rMaxInertia=0.7) badChain = False + ## rescale to expected stretch if using comoving coords + #if useComoving or not sample.useLightCone: + # exp = args[1] + # expError = args[2] + #np.save(voidDir+"/chain.npy", ret) np.savetxt(voidDir+"/fits.out", fits) @@ -1378,11 +1385,8 @@ def launchHubble(dataPortions=None, dataSampleList=None, logDir=None, # (not not fully accurate) plots for (iZBin, stack) in enumerate(sample.getUniqueZBins()): - if sample.dataType == "observation": - aveDist = vp.aveStretchCone(stack.zMin, stack.zMax, - skyFrac = sample.skyFraction) - else: - aveDist = vp.aveStretch(stack.zMin, stack.zMax, Om=sample.omegaM) + aveDist = vp.aveStretch(sample, stack.zMin, stack.zMax, + Om=sample.omegaM) aveDistList[iZBin, 0] = stack.zMin aveDistList[iZBin, 1] = aveDist @@ -1427,11 +1431,8 @@ def launchHubble(dataPortions=None, dataSampleList=None, logDir=None, # skyFrac = sample.skyFraction, # voidRedshifts=voidRedshifts) - if sample.dataType == "observation": - aveDist = vp.aveStretchCone(zBin.zMin, zBin.zMax, - skyFrac = sample.skyFraction) - else: - aveDist = vp.aveStretch(zBin.zMin, zBin.zMax, Om=sample.omegaM) + aveDist = vp.aveStretch(sample, zBin.zMin, zBin.zMax, + Om=sample.omegaM) expList[0, iR, iZBin, 2] = aveDist @@ -1439,11 +1440,8 @@ def launchHubble(dataPortions=None, dataSampleList=None, logDir=None, for (iRCheck,rBinCheck) in enumerate(allRBins): if zBin.zMin == zBinCheck.zMin and zBin.zMax == zBinCheck.zMax: if rBin.rMin == rBinCheck.rMin and rBin.rMax == rBinCheck.rMax: - if sample.dataType == "observation": - aveDist = vp.aveStretchCone(zBin.zMin, zBin.zMax, - skyFrac = sample.skyFraction) - else: - aveDist = vp.aveStretch(zBin.zMin, zBin.zMax) + aveDist = vp.aveStretch(sample, zBin.zMin, zBin.zMax, + Om=sample.omegaM) allExpList[iSample, iRCheck, iZCheck, 0] = exp allExpList[iSample, iRCheck, iZCheck, 1] = expError @@ -1470,8 +1468,8 @@ def launchHubble(dataPortions=None, dataSampleList=None, logDir=None, #sys.stderr = open(logFile, 'a') plotTitle = "Sample: "+sample.nickName+", "+thisDataPortion+" voids" if doPlot: - #vp.do_all_obs(zbase, expList, workDir+"/avedistortion_", - vp.do_all_obs(zbase, expList, aveDistList, + #vp.do_all_obs(sample, zbase, expList, workDir+"/avedistortion_", + vp.do_all_obs(sample, zbase, expList, aveDistList, rlist, plotTitle=plotTitle, plotAve=True) figure(1).savefig(figDir+"/hubble_"+setName+"_"+sampleName+"_"+\ thisDataPortion+\ @@ -1540,11 +1538,8 @@ def launchHubble(dataPortions=None, dataSampleList=None, logDir=None, if zBin.zMaxPlot > plotZmax: plotZmax = zBin.zMaxPlot if zBin.zMinPlot < plotZmin: plotZmin = zBin.zMinPlot - if sample.dataType == "observation": - aveDist = vp.aveStretchCone(zBin.zMin, zBin.zMax, - skyFrac = sample.skyFraction) - else: - aveDist = vp.aveStretch(zBin.zMin, zBin.zMax) + aveDist = vp.aveStretch(sample, zBin.zMin, zBin.zMax, + Om=sample.omegaM) aveDistList[iZ, 0] = zBin.zMin aveDistList[iZ, 1] = aveDist @@ -1567,10 +1562,11 @@ def launchHubble(dataPortions=None, dataSampleList=None, logDir=None, #plotTitle = "all samples, "+thisDataPortion+" voids" plotTitle = setName #plotTitle = '' - vp.do_all_obs(zbase, allExpList, aveDistList, + vp.do_all_obs(sample, zbase, allExpList, aveDistList, rlist, plotTitle=plotTitle, sampleNames=shortSampleNames, plotAve=True, mulfac = 1.0, biasLine = 1.16, - plotZmin=plotZmin, plotZmax=plotZmax) + plotZmin=plotZmin, plotZmax=plotZmax, + plotOm1=True) figure(1).savefig(figDir+"/hubble_combined_"+setName+"_"+ \ thisDataPortion+\ ".eps",bbox_inches='tight') @@ -1588,10 +1584,11 @@ def launchHubble(dataPortions=None, dataSampleList=None, logDir=None, #plotTitle = "all samples, "+thisDataPortion+\ # " voids (systematics corrected)" plotTitle = setName + "(systematics corrected)" - vp.do_all_obs(zbase, allExpList, aveDistList, + vp.do_all_obs(sample, zbase, allExpList, aveDistList, rlist, plotTitle=plotTitle, sampleNames=shortSampleNames, plotAve=True, mulfac = 1.16, - plotZmin=plotZmin, plotZmax=plotZmax) + plotZmin=plotZmin, plotZmax=plotZmax, + plotOm1=True) figure(1).savefig(figDir+"/hubble_combined_"+setName+"_"+\ thisDataPortion+\ "_debiased.eps",bbox_inches='tight')