redshift boundary tolerance is now user-selectable

This commit is contained in:
Paul M. Sutter 2024-06-07 11:58:36 +02:00
parent c0a69c4536
commit 7e5a51d931
8 changed files with 26 additions and 20 deletions

View file

@ -851,8 +851,8 @@ int main(int argc, char **argv) {
// toss out voids near non-periodic box edges // toss out voids near non-periodic box edges
iGood = 0; iGood = 0;
for (iVoid = 0; iVoid < voids.size(); iVoid++) { for (iVoid = 0; iVoid < voids.size(); iVoid++) {
if (tolerance*voids[iVoid].maxRadius > voids[iVoid].nearestEdge || if (voids[iVoid].maxRadius > voids[iVoid].nearestEdge ||
tolerance*voids[iVoid].radius > voids[iVoid].nearestEdge) { voids[iVoid].radius > voids[iVoid].nearestEdge) {
numNearZ++; numNearZ++;
} else { } else {
voids[iGood++] = voids[iVoid]; voids[iGood++] = voids[iVoid];

View file

@ -99,6 +99,10 @@ newSample = Sample(
# max and min redshifts where you want to find voids # max and min redshifts where you want to find voids
zRange = (0.1, 0.15), zRange = (0.1, 0.15),
# width of redshift boundaries to flag edge galaxies
# (interpreted as boundaryWidth*(zMax-zMin))
boundaryWidth = 0.1,
# leave this at -1 for mean particle separation, # leave this at -1 for mean particle separation,
# or specify your own in Mpc/h # or specify your own in Mpc/h
minVoidRadius = -1, minVoidRadius = -1,

View file

@ -73,6 +73,7 @@ class Sample:
selFunFile = "czselfunc.all.dr72dim.dat" selFunFile = "czselfunc.all.dr72dim.dat"
zBoundary = (0.0, 0.1) zBoundary = (0.0, 0.1)
zBoundaryMpc = (0., 300) zBoundaryMpc = (0., 300)
boundaryWidth = 0.02
shiftSimZ = False shiftSimZ = False
zRange = (0.0, 0.1) zRange = (0.0, 0.1)
omegaM = 0.27 omegaM = 0.27
@ -100,7 +101,8 @@ class Sample:
def __init__(self, dataFile="", fullName="", dataUnit=1, def __init__(self, dataFile="", fullName="", dataUnit=1,
nickName="", maskFile="", nsideForMask=128, selFunFile="", nickName="", maskFile="", nsideForMask=128, selFunFile="",
zBoundary=(), zRange=(), zBoundaryMpc=(), shiftSimZ=False, zBoundary=(), zRange=(), zBoundaryMpc=(), boundaryWidth=0.1,
shiftSimZ=False,
minVoidRadius=-1, fakeDensity=0.01, volumeLimited=True, minVoidRadius=-1, fakeDensity=0.01, volumeLimited=True,
numAPSlices=1, numAPSlices=1,
includeInHubble=True, partOfCombo=False, isCombo=False, includeInHubble=True, partOfCombo=False, isCombo=False,
@ -119,6 +121,7 @@ class Sample:
self.selFunFile = selFunFile self.selFunFile = selFunFile
self.zBoundary = zBoundary self.zBoundary = zBoundary
self.zBoundaryMpc = zBoundaryMpc self.zBoundaryMpc = zBoundaryMpc
self.boundaryWidth = boundaryWidth
self.shiftSimZ = shiftSimZ self.shiftSimZ = shiftSimZ
self.zRange = zRange self.zRange = zRange
self.minVoidRadius = minVoidRadius self.minVoidRadius = minVoidRadius

View file

@ -74,7 +74,7 @@ def launchPrep(sample, binPath, workDir=None, inputDataDir=None,
edgeMaskFile = outputDir + "/mask_edge_map.fits" edgeMaskFile = outputDir + "/mask_edge_map.fits"
findEdgeGalaxies(datafile, sample.maskFile, edgeGalFile, edgeMaskFile, findEdgeGalaxies(datafile, sample.maskFile, edgeGalFile, edgeMaskFile,
sample.zBoundary[0], sample.zBoundary[1], sample.omegaM, sample.zBoundary[0], sample.zBoundary[1], sample.omegaM,
useComoving) useComoving, sample.boundaryWidth)
if useComoving: if useComoving:
useComovingFlag = "useComoving" useComovingFlag = "useComoving"

View file

@ -130,11 +130,13 @@ def figureOutMask(galFile, nside, outMaskFile):
# figures out which galaxies live on a mask edge, and also writes the edge # figures out which galaxies live on a mask edge, and also writes the edge
# map to an auxillary file # map to an auxillary file
def findEdgeGalaxies(galFile, maskFile, edgeGalFile, edgeMaskFile, def findEdgeGalaxies(galFile, maskFile, edgeGalFile, edgeMaskFile,
zmin, zmax, omegaM, useComoving=False): zmin, zmax, omegaM, useComoving, boundaryWidth):
#if useComoving: if useComoving:
# zmin = LIGHT_SPEED/100.*comovingDistance(zmin, Om=omegaM) zmin = comovingDistance(zmin, Om=omegaM)
# zmax = LIGHT_SPEED/100.*comovingDistance(zmax, Om=omegaM) zmax = comovingDistance(zmax, Om=omegaM)
#zmin = LIGHT_SPEED/100.*comovingDistance(zmin, Om=omegaM)
#zmax = LIGHT_SPEED/100.*comovingDistance(zmax, Om=omegaM)
#else: #else:
# zmin *= LIGHT_SPEED/100. # zmin *= LIGHT_SPEED/100.
# zmax *= LIGHT_SPEED/100. # zmax *= LIGHT_SPEED/100.
@ -166,8 +168,7 @@ def findEdgeGalaxies(galFile, maskFile, edgeGalFile, edgeMaskFile,
isOnMaskEdge = any(mask[p] == 0 for p in neighbors) isOnMaskEdge = any(mask[p] == 0 for p in neighbors)
# check the redshift boundaries # check the redshift boundaries
tol = 0.05 # TODO: made this user-adjustable zbuffer = (zmax-zmin)*boundaryWidth
zbuffer = (zmax-zmin)*tol
isOnHighZEdge = (z >= zmax-zbuffer) isOnHighZEdge = (z >= zmax-zbuffer)
isOnLowZEdge = (z <= zmin+zbuffer) isOnLowZEdge = (z <= zmin+zbuffer)
@ -182,6 +183,7 @@ def findEdgeGalaxies(galFile, maskFile, edgeGalFile, edgeMaskFile,
edgeFile.write("0\n") edgeFile.write("0\n")
edgeFile.close() edgeFile.close()
healpy.write_map(edgeMaskFile, edgeMask, overwrite=True, dtype=np.dtype('float64')) healpy.write_map(edgeMaskFile, edgeMask, overwrite=True,
dtype=np.dtype('float64'))
return return

View file

@ -121,8 +121,7 @@ for sample in dataSampleList:
launchPrune(sample, PRUNE_PATH, launchPrune(sample, PRUNE_PATH,
logFile=logFile, outputDir=outputDir, logFile=logFile, outputDir=outputDir,
useComoving=sample.useComoving, continueRun=continueRun, useComoving=sample.useComoving, continueRun=continueRun,
mergingThreshold=mergingThreshold, mergingThreshold=mergingThreshold)
boundaryTolerance=boundaryTolerance)
# ------------------------------------------------------------------------- # -------------------------------------------------------------------------
if (startCatalogStage <= 4) and (endCatalogStage >= 4): if (startCatalogStage <= 4) and (endCatalogStage >= 4):

View file

@ -143,10 +143,6 @@ numZobovThreads = 2
# 1e-9 (or smaller != 0) -> Do not merge anything # 1e-9 (or smaller != 0) -> Do not merge anything
mergingThreshold = 1.e-9 mergingThreshold = 1.e-9
# when trimming away voids near the bounaries, what multiple of the radius to
# use for safety
boundaryTolerance = 1.0
# simulation information # simulation information
numPart = 512*512*512 numPart = 512*512*512
lbox = 999.983 # Mpc/h lbox = 999.983 # Mpc/h

View file

@ -25,6 +25,7 @@ import numpy as np
import os import os
import pylab as plt import pylab as plt
import backend.cosmologyTools as vp import backend.cosmologyTools as vp
import backend.surveyTools as sp
from voidUtil import getArray, shiftPart, getVoidPart from voidUtil import getArray, shiftPart, getVoidPart
def fill_between(x, y1, y2=0, ax=None, **kwargs): def fill_between(x, y1, y2=0, ax=None, **kwargs):
@ -77,17 +78,18 @@ def plotNumberFunction(catalogList,
if sample.dataType == "observation": if sample.dataType == "observation":
maskFile = sample.maskFile maskFile = sample.maskFile
boxVol = vp.getSurveyProps(maskFile, boxVol = sp.getSurveyProps(maskFile,
sample.zBoundary[0], sample.zBoundary[1], sample.zBoundary[0], sample.zBoundary[1],
sample.zRange[0], sample.zRange[1], "all", sample.zRange[0], sample.zRange[1], "all",
selectionFuncFile=sample.selFunFile)[0] selectionFuncFile=sample.selFunFile,
omegaM=sample.omegaM)[0]
else: else:
boxVol = sample.boxLen*sample.boxLen*(sample.zBoundaryMpc[1] - boxVol = sample.boxLen*sample.boxLen*(sample.zBoundaryMpc[1] -
sample.zBoundaryMpc[0]) sample.zBoundaryMpc[0])
boxVol *= 1.e-9 # Mpc->Gpc boxVol *= 1.e-9 # Mpc->Gpc
bins = 100./binWidth bins = int(100./binWidth)
hist, binEdges = np.histogram(data, bins=bins, range=(0., 100.)) hist, binEdges = np.histogram(data, bins=bins, range=(0., 100.))
binCenters = 0.5*(binEdges[1:] + binEdges[:-1]) binCenters = 0.5*(binEdges[1:] + binEdges[:-1])