added more flexible catalog loading options. can now automatically differentiate between old and new catalogs. first step in eliminating cumbersome void outputs

python_source/voidUtil/catalogUtil.py

@@ -288,19 +288,53 @@ class Catalog:
   sampleInfo = None
 
 # -----------------------------------------------------------------------------
-def loadVoidCatalog(sampleDir, dataPortion="central", loadParticles=True,
-                    loadAdjacencies = True, untrimmed=False):
+def loadVoidCatalog(sampleDir,
+                    loadParticles = False,
+                    loadAdjacencies = False, 
+                    clearEdges = False, 
+                    clearTree = False,
+                    clearNearBoundaries = False,
+                    maxCentralDen = -1,
+                    replicateOldCentralVoids = 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
+#   clearEdges: if True, remove voids with any edge contamination
+#   clearTree: if True, remove all non-leaf voids
+#   clearNearBoundaries: remove voids where the maximum extent is
+#                        greater than the distance to nearest edge
+#   maxHighCentralDen: if != -1, filters based on based on central density
+
+# NOTE: we are moving away from the cumbersome void catalog outputs of
+#       older versions, and will eventually just output a single catalog.
+#       To replicate the old "central" catalog, choose:
+#         clearEdges = True
+#         clearTree = True
+#         clearNearBoundaries = True
+#         maxCentralDen = 1.e-9
+#
+#       ~or~ set replicateOldCentralVoids = True
 
   sys.stdout.flush()
 
+  print("Loading catalog from ", sampleDir)
+
+  isOldCatalog = os.path.exists(sampleDir+"/mask_index.txt")
+
+  if isOldCatalog and clearNearBoundaries:
+    print("WARNING: Old catalog. Unable to clear near boundaries.")
+  if isOldCatalog and maxCentralDen != -1:
+    print("WARNING: Old catalog. Central density cuts already applied.")
+
+  if replicateOldCentralVoids:
+    clearEdges = True
+    clearTree = True
+    clearNearBoundaries = True
+    maxCentralDen = 1.e-9
+
   catalog = Catalog()
   with open(sampleDir+"/sample_info.dat", 'rb') as input:
     sample = pickle.load(input)
@@ -326,13 +360,17 @@ def loadVoidCatalog(sampleDir, dataPortion="central", loadParticles=True,
   catalog.volNormZobov = volNormZobov
   catalog.volNormObs = volNormObs
 
-  if untrimmed:
-    prefix = "untrimmed_"
-  else:
-    prefix = ""
+  # for new catalogs, we will load by default the whole shebang, then
+  # apply filters later. for old catalogs, we need to pick the right file
+  prefix = "untrimmed_"
+  if isOldCatalog and clearTree: prefix = ""
+
+  dataPortion = "all"
+  if isOldCatalog and clearEdges: dataPortion = "central"
 
   print("Loading voids...")
   fileName = sampleDir+"/"+prefix+"voidDesc_"+dataPortion+"_"+sample.fullName+".out"
+
   catData = np.loadtxt(fileName, comments="#", skiprows=2)
   catalog.voids = []
   for line in catData:
@@ -433,7 +471,17 @@ def loadVoidCatalog(sampleDir, dataPortion="central", loadParticles=True,
      iLine += 1
 
   else:
-    print("  Extra info file not found!")
+    print(" Old catalog: extra info file not found")
+
+  # apply filters to new catalogs
+  if not isOldCatalog:
+    print("Filtering catalog...")
+    if clearEdges: catalog = filterOnType(catalog, CENTRAL_VOID)
+    if clearTree: catalog = filterOnTreeLevel(catalog, level=-1)
+    if clearNearBoundaries: catalog = filterOnNearestEdge(catalog)
+    if maxCentralDen != -1: catalog = filterOnCentralDen(catalog, maxCentralDen)
+
+    print("  After filtering there are %d voids remaining" % catalog.numVoids)
 
   if loadParticles:
     print("Loading all particles...")
@@ -587,39 +635,40 @@ def getVoidPart(catalog, voidID):
 # -----------------------------------------------------------------------------
 # various handy catalog filtering routines
 
-def filterVoidsOnSize(catalog, rMin):
+def filterOnSize(catalog, rMin):
   catalog.voids = [v for v in catalog.voids if v.radius >= rMin]
   catalog.numVoids = len(catalog.voids)
   return catalog
 
-def filterVoidsOnTreeLevel(catalog, level):
-  catalog.voids = [v for v in catalog.voids if v.treeLevel == level]
+def filterOnTreeLevel(catalog, level):
   if level == -1:
     catalog.voids = [v for v in catalog.voids if v.numChildren == 0]
+  else:
+    catalog.voids = [v for v in catalog.voids if v.treeLevel == level]
   catalog.numVoids = len(catalog.voids)
   return catalog
 
-def filterVoidsOnCentralDen(catalog, maxCentralDen):
+def filterOnCentralDen(catalog, maxCentralDen):
   catalog.voids = [v for v in catalog.voids if v.centralDen <= maxCentralDen]
   catalog.numVoids = len(catalog.voids)
   return catalog
 
-def filterVoidsOnCoreDen(catalog, maxCoreDen):
+def filterOnCoreDen(catalog, maxCoreDen):
   catalog.voids = [v for v in catalog.voids if v.coreDens <= maxCoreDen]
   catalog.numVoids = len(catalog.voids)
   return catalog
 
-def filterVoidsOnDenCon(catalog, minDenCon):
+def filterOnDensCon(catalog, minDenCon):
   catalog.voids = [v for v in catalog.voids if v.densCon >= minDenCon]
   catalog.numVoids = len(catalog.voids)
   return catalog
 
-def filterVoidsOnType(catalog, voidType):
+def filterOnType(catalog, voidType):
   catalog.voids = [v for v in catalog.voids if v.voidType == voidType]
   catalog.numVoids = len(catalog.voids)
   return catalog
 
-def filterVoidsOnNearestEdge(catalog, factor=1.0):
+def filterOnNearestEdge(catalog, factor=1.0):
   catalog.voids = [v for v in catalog.voids if \
                    factor*v.maxRadius <= v.nearestEdge]
   catalog.numVoids = len(catalog.voids)