further enhancements to void merger trees

c_tools/analysis/voidOverlap.cpp

@@ -245,6 +245,7 @@ void loadCatalog(const char *partFile, const char *volFile,
   float ranges[3][2];
 
   printf("Loading info...\n");
+printf("HELLO %s\n", infoFile);
   NcFile f_info(infoFile);
   ranges[0][0] = f_info.get_att("range_x_min")->as_double(0);
   ranges[0][1] = f_info.get_att("range_x_max")->as_double(0);

crossCompare/analysis/datasetsToAnalyze.py

@@ -3,16 +3,15 @@
 
 workDir = "/home/psutter2/workspace/Voids/"
 figDir  = "./figs"
+dataDir = "./data"
+ 
+
+CTOOLS_PATH = "../../c_tools/"
 
 sampleDirList = [ 
-                  "multidark/md_ss0.1_pv/sample_md_ss0.1_pv_z0.56_d00/",
+                  "multidark/md_ss0.000175/sample_md_ss0.000175_z0.56_d00/",
-                  "multidark/md_halos_min1.393e12_pv/sample_md_halos_min1.393e12_pv_z0.56_d00/",
+                  "multidark/md_ss0.0004/sample_md_ss0.0004_z0.56_d00/",
-                  "random/ran_ss0.000175/sample_ran_ss0.000175_z0.56_d00/",
+                ]
-                  "random/ran_ss0.1/sample_ran_ss0.1_z0.56_d00/",
-                  "multidark/md_hod_dr9mid_pv/sample_md_hod_dr9mid_pv_z0.56_d00/",
-                  "multidark/md_ss0.000175_pv/sample_md_ss0.000175_pv_z0.56_d00/",
-                  "sdss_dr9/sample_lss.dr9cmassmid.dat/",
-                  "lanl/masked/masked_lanl_hod_dr9mid_pv/sample_masked_lanl_hod_dr9mid_pv_z0.0/" ]
 
 dataPortion = "central"
 

crossCompare/analysis/mergerTree.py

@@ -4,7 +4,6 @@
 
 from void_python_tools.backend import *
 from void_python_tools.plotting import *
-import void_python_tools.apTools as vp
 import imp
 import pickle
 import os
@@ -14,16 +13,10 @@ import argparse
 
 # ------------------------------------------------------------------------------
 
-plotNameBase = "compdist"
+dataNameBase = "mergerTree"
 
-obsFudgeFactor = 1.0 # what fraction of the volume are we *reall* capturing?
+parser = argparse.ArgumentParser(description='Analyze.')
-#obsFudgeFactor = .66 # what fraction of the volume are we *reall* capturing?
+parser.add_argument('--parmFile', dest='parmFile', default='datasetsToAnalyze.py',
-
-parser = argparse.ArgumentParser(description='Plot.')
-parser.add_argument('--show', dest='showPlot', action='store_const',
-                   const=True, default=False,
-                   help='display the plot (default: just write eps)')
-parser.add_argument('--parmFile', dest='parmFile', default='datasetsToPlot.py',
                     help='path to parameter file')
 args = parser.parse_args()
 
@@ -40,82 +33,45 @@ globals().update(vars(parms))
 if not os.access(figDir, os.F_OK):
   os.makedirs(figDir)
 
-dataSampleList = []
+outFileName = dataDir + "/" + dataNameBase + ".dat"
+
+for (iSample, sampleDir) in enumerate(sampleDirList):
+  if iSample == 0: continue
 
-for sampleDir in sampleDirList:
   with open(workDir+sampleDir+"/sample_info.dat", 'rb') as input:
-    dataSampleList.append(pickle.load(input))
+    sample = pickle.load(input)
 
-plt.clf()
+  with open(workDir+sampleDirList[0]+"/sample_info.dat", 'rb') as input:
-plt.xlabel("Void Radius [Mpc/h]")
+    baseSample = pickle.load(input)
-plt.ylabel(r"N > R [$h^3$ Gpc$^{-3}$]")
-plt.yscale('log')
-plt.xlim(xmax=120.)
 
-plotName = plotNameBase
+  print " Working with", sample.fullName, "..."
-allData = []
-
-for (iSample,sample) in enumerate(dataSampleList):
 
   sampleName = sample.fullName
-  lineTitle = sampleName
 
-  if sample.dataType == "observation":
+  binPath = CTOOLS_PATH+"/analysis/voidOverlap"
-    boxVol = vp.getSurveyProps(sample.maskFile, 
+  logFile = os.getcwd()+"/mergerTree.log"
-                               sample.zBoundary[0], sample.zBoundary[1], 
+  stepOutputFileName = os.getcwd()+"/thisStep.out"
-                               sample.zRange[0], sample.zRange[1], "all",
+
-                               selectionFuncFile=None)[0]
+  launchVoidOverlap(baseSample, sample, workDir+sampleDirList[0], 
-                               #selectionFuncFile=sample.selFunFile)[0]
+                    workDir+sampleDir, binPath, 
-    boxVol *= obsFudgeFactor
+                    thisDataPortion="central", logFile=logFile,
+                    continueRun=False, workDir=workDir,
+                    outputFile=stepOutputFileName)
+
+  # attach columns to summary file
+  if iSample == 1:
+    os.system("cp %s %s" % (stepOutputFileName, outFileName))
   else:
-    boxVol = sample.boxLen*sample.boxLen*(sample.zBoundaryMpc[1] - 
+    outFile = fp.open("temp.out", "w")
-                                          sample.zBoundaryMpc[0])
+    stepOutFile1 = fp.open(outFileName, "r")
+    stepOutFile2 = fp.open(stepOutputFileName, "r")
 
-  boxVol *= 1.e-9 # Mpc->Gpc
+    for line1 in stepOutFile1:
-
+      line2 = stepOutFile2.readline()
-  filename = workDir+"/"+sampleDirList[iSample]+"/centers_"+dataPortion+"_"+\
+      outFile.write(line1+" "+line2)
-             sampleName+".out"
-  if not os.access(filename, os.F_OK):
-    print "File not found: ", filename
-    continue
-
-  data = np.loadtxt(filename, comments="#")
-  if data.ndim == 1:
-    print " Too few!"
-    continue
-  data = data[:,4]
-  indices = np.arange(0, len(data), 1)
-  sorted = np.sort(data)
-
-  plt.plot(sorted, indices[::-1]/boxVol, '-',
-           label=lineTitle, color=colorList[iSample],
-           linewidth=linewidth)
-
-  hist, bin_edges = np.histogram(data, bins=100, range=(0,100)) 
-  allData.append(hist)
-
-plt.legend(title = "Samples", loc = "upper right", prop={'size':8})
-#plt.title(plotTitle)
-
-plt.savefig(figDir+"/fig_"+plotName+".pdf", bbox_inches="tight")
-plt.savefig(figDir+"/fig_"+plotName+".eps", bbox_inches="tight")
-plt.savefig(figDir+"/fig_"+plotName+".png", bbox_inches="tight")
-
-dataFile = figDir+"/data_"+plotName+".dat"
-fp = open(dataFile, 'w')
-fp.write("# R [Mpc/h], N [h^3 Gpc^-3]\n")
-fp.write("# ")
-for sample in dataSampleList:
-  fp.write(sample.fullName+" ")
-fp.write("\n")
-for i in xrange(100):
-  fp.write(str(bin_edges[i]) + " ")
-  for iSample in xrange(len(dataSampleList)):
-    fp.write(str(allData[iSample][i])+" ")
-  fp.write("\n")
-fp.close()
-
-if args.showPlot:
-  os.system("display %s" % figDir+"/fig_"+plotName+".png")
 
+    outFile.close()
+    stepOutFile1.close()
+    stepOutFile2.close()
 
+if os.access("mergerTree.log", os.F_OK): os.unlink("mergerTree.log")

python_tools/void_python_tools/backend/launchers.py

@@ -309,6 +309,70 @@ def launchPrune(sample, binPath, thisDataPortion=None,
     print "already done!"
 
  
+# -----------------------------------------------------------------------------
+def launchVoidOverlap(sample1, sample2, sample1Dir, sample2Dir, 
+                      binPath, thisDataPortion=None, 
+                      logFile=None, workDir=None,
+                      continueRun=None, outputFile=None):
+
+  sampleName1 = sample1.fullName
+  sampleName2 = sample2.fullName
+
+  periodicLine = " --periodic='"
+  if sample1.dataType != "observation":
+    if sample1.numSubvolumes == 1: periodicLine += "xy"
+    if sample1.zBoundaryMpc[0] == 0 and \
+       sample1.zBoundaryMpc[1] == sample1.boxLen : periodicLine += "z"
+    periodicLine += "' "
+
+  if not (continueRun and jobSuccessful(logFile, "Done!\n")):
+    cmd = binPath
+    cmd += " --partFile1=" + sample1Dir+"/zobov_slice_" + \
+           str(sampleName1)
+    cmd += " --volFile1=" + sample1Dir+"/vol_" + \
+           str(sampleName1)+".dat"
+    cmd += " --voidFile1=" + sample1Dir+"/voidDesc_" + \
+           thisDataPortion+"_"+str(sampleName1)+".out"
+    cmd += " --infoFile1=" + sample1Dir+"/zobov_slice_" + \
+           str(sampleName1)+".par"
+    cmd += " --barycenterFile1=" + sample1Dir + \
+           "/barycenters_"+thisDataPortion+"_"+str(sampleName1)+".out"
+    cmd += " --zoneFile1=" + sample1Dir+"/voidZone_" + \
+           str(sampleName1)+".dat"
+    cmd += " --zonePartFile1=" + sample1Dir+"/voidPart_" + \
+           str(sampleName1)+".dat"
+
+    cmd += " --partFile2=" + sample2Dir+"/zobov_slice_" + \
+           str(sampleName2)
+    cmd += " --volFile2=" + sample2Dir+"/vol_" + \
+           str(sampleName2)+".dat"
+    cmd += " --voidFile2=" + sample2Dir+"/voidDesc_" + \
+           thisDataPortion+"_"+str(sampleName2)+".out"
+    cmd += " --infoFile2=" + sample2Dir+"/zobov_slice_" + \
+           str(sampleName2)+".par"
+    cmd += " --barycenterFile2=" + sample2Dir + \
+           "/barycenters_"+thisDataPortion+"_"+str(sampleName2)+".out"
+    cmd += " --zoneFile2=" + sample2Dir+"/voidZone_" + \
+           str(sampleName2)+".dat"
+    cmd += " --zonePartFile2=" + sample2Dir+"/voidPart_" + \
+           str(sampleName2)+".dat"
+
+    #cmd += " --useID"
+    cmd += periodicLine
+    cmd += " --outfile=" + outputFile
+    cmd += " &> " + logFile
+    os.system(cmd)
+
+    if jobSuccessful(logFile, "Done!\n"):
+      print "done"
+    else:
+      print "FAILED!"
+      exit(-1)
+
+  else:
+    print "already done!"
+
+
 # -----------------------------------------------------------------------------
 def launchStack(sample, stack, binPath, thisDataPortion=None, logDir=None,
                 voidDir=None, freshStack=True, runSuffix=None,