From 4c6ec4ba825f84b36d778a7e4ff54bca6432fde1 Mon Sep 17 00:00:00 2001 From: "Paul M. Sutter" Date: Sat, 22 Mar 2025 19:28:14 -0400 Subject: [PATCH] added voidID property to particles in catalog output --- python_source/voidUtil/catalogUtil.py | 69 ++++++++++++++++----------- 1 file changed, 42 insertions(+), 27 deletions(-) diff --git a/python_source/voidUtil/catalogUtil.py b/python_source/voidUtil/catalogUtil.py index ead60a2..9adacf0 100644 --- a/python_source/voidUtil/catalogUtil.py +++ b/python_source/voidUtil/catalogUtil.py @@ -295,6 +295,7 @@ class Catalog: boxLen = np.zeros((3)) ranges = np.zeros((3,2)) part = None + partPos = None zones2Parts = None void2Zones = None voids = None @@ -302,13 +303,14 @@ class Catalog: # ----------------------------------------------------------------------------- def loadVoidCatalog(sampleDir, dataPortion="central", loadParticles=True, - untrimmed=False): + loadAdjacencies = True, untrimmed=False): # loads a void catalog # by default, loads parent-level voids with central densities greater than 0.2*mean # sampleDir: path to VIDE output directory # dataPortion: "central" or "all" # loadParticles: if True, also load particle information +# loadAdjacencies: if True, also load particle adjacency information # untrimmed: if True, catalog contains all voids, regardless of density or hierarchy level sys.stdout.flush() @@ -448,6 +450,7 @@ def loadVoidCatalog(sampleDir, dataPortion="central", loadParticles=True, dec = extraData[i][1], redshift = extraData[i][2], uniqueID = extraData[i][3], + voidID = -1, edgeFlag = 0)) @@ -473,46 +476,51 @@ def loadVoidCatalog(sampleDir, dataPortion="central", loadParticles=True, for ivol in range(len(vols)): catalog.part[ivol].volume = vols[ivol] / volNorm - print(" Loading adjacencies...") - adjFile = sampleDir+"adj_"+sample.fullName+".dat" - with open(adjFile, mode="rb") as File: - numPart = np.fromfile(File, dtype=np.int32,count=1)[0] + if loadAdjacencies: + print(" Loading adjacencies...") + adjFile = sampleDir+"adj_"+sample.fullName+".dat" + with open(adjFile, mode="rb") as File: + numPart = np.fromfile(File, dtype=np.int32,count=1)[0] - # this the total number of adjancies per particle - nadjPerPart = np.fromfile(File, dtype=np.int32,count=numPart) + # this the total number of adjancies per particle + nadjPerPart = np.fromfile(File, dtype=np.int32,count=numPart) - # but the file only stores one half of each pair, so we need to match + # but the file only stores one half of each pair, + # so we need to match + for p in range(numPart): + nin = np.fromfile(File, dtype=np.int32, count=1)[0] + for n in range(nin): + pAdj = np.fromfile(File, dtype=np.int32, count=1)[0] + if (p < pAdj): + catalog.part[p].adjs.append(pAdj) + catalog.part[pAdj].adjs.append(p) + + + print(" Sanity checking adjacenies...") for p in range(numPart): - nin = np.fromfile(File, dtype=np.int32, count=1)[0] - for n in range(nin): - pAdj = np.fromfile(File, dtype=np.int32, count=1)[0] - if (p < pAdj): - catalog.part[p].adjs.append(pAdj) - catalog.part[pAdj].adjs.append(p) + catalog.part[p].nadjs = len(catalog.part[p].adjs) + nHave = len(catalog.part[p].adjs) + nExpected = nadjPerPart[p] + # interior galaxies should not connect to + if (nHave != nExpected and catalog.part[p].edgeFlag == 0): + print(" Error for particle %d: Have %d adj, expected %d (flag: %d)" % (p, nHave, nExpected, catalog.part[p].edgeFlag)) + # end load adjacencies - print(" Sanity checking adjacenies...") - for p in range(numPart): - catalog.part[p].nadjs = len(catalog.part[p].adjs) - nHave = len(catalog.part[p].adjs) - nExpected = nadjPerPart[p] - if (nHave != nExpected and catalog.part[p].edgeFlag == 0): - print(" Error for particle %d: Have %d adj, expected %d (flag: %d)" % (p, nHave, nExpected, catalog.part[p].edgeFlag)) - print(" Loading zone-void membership info...") zoneFile = sampleDir+"/voidZone_"+sample.fullName+".dat" catalog.void2Zones = [] with open(zoneFile, mode="rb") as File: - numZonesTot = np.fromfile(File, dtype=np.int32,count=1)[0] - catalog.numZonesTot = numZonesTot - for iZ in range(numZonesTot): + numVoidsTot = np.fromfile(File, dtype=np.int32,count=1)[0] + catalog.numVoidsTot = numVoidsTot + for iVoid in range(numVoidsTot): numZones = np.fromfile(File, dtype=np.int32,count=1)[0] catalog.void2Zones.append(Bunch(numZones = numZones, zoneIDs = [])) - for p in range(numZones): + for iZ in range(numZones): zoneID = np.fromfile(File, dtype=np.int32,count=1)[0] - catalog.void2Zones[iZ].zoneIDs.append(zoneID) + catalog.void2Zones[iVoid].zoneIDs.append(zoneID) print(" Loading particle-zone membership info...") zonePartFile = sampleDir+"/voidPart_"+sample.fullName+".dat" @@ -529,6 +537,13 @@ def loadVoidCatalog(sampleDir, dataPortion="central", loadParticles=True, partID = np.fromfile(File, dtype=np.int32,count=1)[0] catalog.zones2Parts[iZ].partIDs.append(partID) + print(" Matching particles to voids...") + for void in catalog.voids: + voidID = void.voidID + for zoneID in catalog.void2Zones[voidID].zoneIDs: + for partID in catalog.zones2Parts[zoneID].partIDs: + catalog.part[partID].voidID = voidID + print("Done loading catalog.") return catalog