From d9efa15474524ea2ae9bc42cb8762ac2a232c9c3 Mon Sep 17 00:00:00 2001 From: "Paul M. Sutter" Date: Mon, 21 Apr 2025 16:23:25 -0400 Subject: [PATCH] removed mock particle tracking, since that is no longer used --- c_source/prep/prepObservation.cpp | 18 ++------ c_source/pruning/pruneVoids.cpp | 56 ++++++++++-------------- c_source/pruning/pruneVoids.ggo | 1 - c_source/zobov/jozov.c | 16 +------ c_source/zobov/vozinit.c | 21 +++------ c_source/zobov/voztie.c | 46 +++++--------------- python_source/backend/launchers.py | 61 ++++++++++----------------- python_source/voidUtil/catalogUtil.py | 25 ++--------- 8 files changed, 71 insertions(+), 173 deletions(-) diff --git a/c_source/prep/prepObservation.cpp b/c_source/prep/prepObservation.cpp index 6926270..e641fbb 100644 --- a/c_source/prep/prepObservation.cpp +++ b/c_source/prep/prepObservation.cpp @@ -64,9 +64,6 @@ struct ParticleData vector redshift; vector catalogID; vector uniqueID; - // PMS - int mask_index; - // END PMS int edgeFlag = 0; vector pos; double box[3][2]; @@ -263,12 +260,8 @@ void flagEdgeGalaxies(prepObservation_info& args , size_t nEval; //End Maubert - Needed for comobile distance in mock_sphere - // TODO - REMOVE THIS - output_data.mask_index = output_data.id_gal.size(); - // write a small text file with galaxy position (for diagnostic purposes) // TODO - remove this - FILE *fp; fp = fopen("galaxies.txt", "w"); //for (int i = 0; i < output_data.id_gal.size(); i++) { @@ -414,11 +407,9 @@ void saveForZobov(ParticleData& pdata, const string& fname, fp.putAtt("range_y_max", ncDouble, Lmax/100.); fp.putAtt("range_z_min", ncDouble, -Lmax/100.); fp.putAtt("range_z_max", ncDouble, Lmax/100.); - fp.putAtt("mask_index", ncInt, pdata.mask_index); // PMS - fp.putAtt("is_observation", ncInt, 1); // PMS + fp.putAtt("is_observation", ncInt, 1); - int nOutputPart = pdata.mask_index; - //int nOutputPart = pdata.pos.size(); + int nOutputPart = pdata.pos.size(); NcDim NumPart_dim = fp.addDim("numpart_dim", nOutputPart); NcVar v = fp.addVar("particle_ids", ncInt, NumPart_dim); @@ -497,10 +488,7 @@ int main(int argc, char **argv) saveForZobov(output_data, args_info.output_arg, args_info.params_arg); // PMS - TODO REMOVE THIS - FILE *fp = fopen("mask_index.txt", "w"); - fprintf(fp, "%d", output_data.mask_index); - fclose(fp); - + FILE *fp; fp = fopen("total_particles.txt", "w"); fprintf(fp, "%d", output_data.pos.size()); fclose(fp); diff --git a/c_source/pruning/pruneVoids.cpp b/c_source/pruning/pruneVoids.cpp index a3e4be2..3927608 100644 --- a/c_source/pruning/pruneVoids.cpp +++ b/c_source/pruning/pruneVoids.cpp @@ -116,7 +116,7 @@ double expanFun (double z, void * p) { void openFiles(string outputDir, string sampleName, string prefix, string dataPortion, - int mockIndex, int numKept, + int numPartTot, int numKept, FILE** fpZobov, FILE** fpCenters, FILE** fpCentersNoCut, FILE** fpBarycenter, FILE** fpShapes, @@ -128,7 +128,7 @@ void closeFiles(FILE* fpZobov, FILE* fpCenters, FILE* fpSkyPositions); void outputVoids(string outputDir, string sampleName, string prefix, - string dataPortion, int mockIndex, + string dataPortion, int numPartTot, vector voids, bool isObservation, double *boxLen, bool doTrim, bool doCentralDenCut); @@ -187,7 +187,7 @@ int main(int argc, char **argv) { int i, p, p2, numPartTot, numZonesTot, dummy, iVoid; - int numVoids, mockIndex; + int numVoids; double tolerance; FILE *fp; PART *part, *voidPart, *flaggedPart; @@ -202,7 +202,6 @@ int main(int argc, char **argv) { double nearestEdge, redshift; char line[500], junkStr[10]; string outputDir, sampleName, dataPortion, prefix; - int mask_index; double ranges[3][2], boxLen[3], mul; double volNorm, radius; int clock1, clock2, clock3, clock4; @@ -213,7 +212,6 @@ int main(int argc, char **argv) { gsl_eigen_symmv_workspace *eigw = gsl_eigen_symmv_alloc(3); numVoids = args.numVoids_arg; - mockIndex = args.mockIndex_arg; tolerance = args.tolerance_arg; clock1 = clock(); @@ -309,9 +307,6 @@ int main(int argc, char **argv) { interval = 1.*(clock4 - clock3)/CLOCKS_PER_SEC; printf(" Read %d particles (%.2f sec)...\n", numPartTot, interval); - - if (mockIndex == -1) mockIndex = numPartTot; - // read in desired voids clock3 = clock(); printf(" Loading voids...\n"); @@ -393,12 +388,8 @@ int main(int argc, char **argv) { // and volumes printf(" Loading particle volumes...\n"); fp = fopen(args.partVol_arg, "r"); - fread(&mask_index, 1, 4, fp); - if (mask_index != mockIndex) { - printf("NON-MATCHING MOCK INDICES!? %d %d\n", mask_index, mockIndex); - exit(-1); - } - for (p = 0; p < mask_index; p++) { + fread(&numPartTot, 1, 4, fp); + for (p = 0; p < numPartTot; p++) { fread(&temp[0], 1, 4, fp); part[p].vol = temp[0]; } @@ -408,7 +399,7 @@ int main(int argc, char **argv) { printf(" Loading particle edge flags...\n"); int numFlagged = 0; fp = fopen(args.partEdge_arg, "r"); - for (p = 0; p < mask_index; p++) { + for (p = 0; p < numPartTot; p++) { fscanf(fp, "%d", &part[p].edgeFlag); if (part[p].edgeFlag > 0) numFlagged++; } @@ -417,7 +408,7 @@ int main(int argc, char **argv) { // copy flagged galaxies to separate array for faster processing later flaggedPart = (PART *) malloc(numFlagged * sizeof(PART)); int iFlag = 0; - for (p = 0; p < mask_index; p++) { + for (p = 0; p < numPartTot; p++) { if (part[p].edgeFlag > 0) { flaggedPart[iFlag] = part[p]; iFlag++; @@ -428,20 +419,20 @@ int main(int argc, char **argv) { // and finally finally adjacencies printf(" Loading particle adjacencies...\n"); fp = fopen(args.partAdj_arg, "r"); - fread(&mask_index, 1, 4, fp); - if (mask_index != mockIndex) { - printf("NON-MATCHING MOCK INDICES!? %d %d\n", mask_index, mockIndex); + fread(&numPartTot, 1, 4, fp); + if (numPartTot != numPartTot) { + printf("NON-MATCHING MOCK INDICES!? %d %d\n", numPartTot, numPartTot); exit(-1); } int tempInt; - for (p = 0; p < mask_index; p++) { + for (p = 0; p < numPartTot; p++) { fread(&tempInt, 1, 4, fp); part[p].nadj = tempInt; part[p].ncnt = 0; if (part[p].nadj > 0) part[p].adj = (int *) malloc(part[p].nadj * sizeof(int)); } - for (p = 0; p < mask_index; p++) { + for (p = 0; p < numPartTot; p++) { fread(&tempInt, 1, 4, fp); int nin = tempInt; if (nin > 0) { @@ -452,7 +443,7 @@ int main(int argc, char **argv) { // this bit has been readjusted just in case we are // accidentally still linking to mock particles - //if (tempAdj < mask_index) { + //if (tempAdj < numPartTot) { assert(p < tempAdj); //if (part[p].ncnt == part[p].nadj) { // printf("OVERFLOW %d\n", p); @@ -461,7 +452,7 @@ int main(int argc, char **argv) { //} else { part[p].adj[part[p].ncnt] = tempAdj; part[p].ncnt++; - if (tempAdj < mask_index) { + if (tempAdj < numPartTot) { part[tempAdj].adj[part[tempAdj].ncnt] = p; part[tempAdj].ncnt++; } @@ -542,9 +533,8 @@ int main(int argc, char **argv) { partID = zones2Parts[zoneID].partIDs[p]; // something went haywire - if (partID > mask_index || - (part[partID].vol < 1.e-27 && part[partID].vol > 0.)) { - printf("BAD PART!? %d %d %e", partID, mask_index, part[partID].vol); + if (part[partID].vol < 1.e-27 && part[partID].vol > 0.) { + printf("BAD PART!? %d %d %e", partID, numPartTot, part[partID].vol); exit(-1); } @@ -567,7 +557,7 @@ int main(int argc, char **argv) { //printf("NORMAL!! %d %d %e\n", iVoid, partID, part[partID].vol); } for (int iAdj = 0; iAdj < part[partID].ncnt; iAdj++) { - if (part[partID].adj[iAdj] > mockIndex) { + if (part[partID].adj[iAdj] > numPartTot) { printf("CONTAMINATED!! %d %d %d\n", iVoid, partID, iAdj); } } @@ -949,7 +939,7 @@ int main(int argc, char **argv) { for (int i = 0; i < 2; i++) { dataPortion = dataPortions[i]; outputVoids(outputDir, sampleName, prefix, dataPortion, - mockIndex, + numPartTot, voids, args.isObservation_flag, boxLen, true, true); } @@ -959,7 +949,7 @@ int main(int argc, char **argv) { for (int i = 0; i < 2; i++) { dataPortion = dataPortions[i]; outputVoids(outputDir, sampleName, prefix, dataPortion, - mockIndex, + numPartTot, voids, args.isObservation_flag, boxLen, false, false); } @@ -977,14 +967,14 @@ int main(int argc, char **argv) { // ---------------------------------------------------------------------------- void openFiles(string outputDir, string sampleName, string prefix, string dataPortion, - int mockIndex, int numKept, + int numPartTot, int numKept, FILE** fpZobov, FILE** fpCenters, FILE** fpBarycenter, FILE** fpShapes, FILE** fpExtra, FILE** fpSkyPositions) { *fpZobov = fopen((outputDir+"/"+prefix+"voidDesc_"+dataPortion+"_"+sampleName).c_str(), "w"); - fprintf(*fpZobov, "%d particles, %d voids.\n", mockIndex, numKept); + fprintf(*fpZobov, "%d particles, %d voids.\n", numPartTot, numKept); fprintf(*fpZobov, "Void# FileVoid# CoreParticle CoreDens ZoneVol Zone#Part Void#Zones VoidVol Void#Part VoidDensContrast VoidProb\n"); *fpBarycenter = fopen((outputDir+"/"+prefix+"macrocenters_"+dataPortion+"_"+sampleName).c_str(), "w"); @@ -1022,7 +1012,7 @@ void closeFiles(FILE* fpZobov, FILE* fpCenters, // ---------------------------------------------------------------------------- void outputVoids(string outputDir, string sampleName, string prefix, - string dataPortion, int mockIndex, + string dataPortion, int numPartTot, vector voids, bool isObservation, double *boxLen, bool doTrim, bool doCentralDenCut) { @@ -1034,7 +1024,7 @@ void outputVoids(string outputDir, string sampleName, string prefix, openFiles(outputDir, sampleName, prefix, dataPortion, - mockIndex, voids.size(), + numPartTot, voids.size(), &fpZobov, &fpCenters, &fpBarycenter, &fpShapes, &fpExtra, &fpSkyPositions); diff --git a/c_source/pruning/pruneVoids.ggo b/c_source/pruning/pruneVoids.ggo index 1f02fc3..4c711ba 100644 --- a/c_source/pruning/pruneVoids.ggo +++ b/c_source/pruning/pruneVoids.ggo @@ -19,7 +19,6 @@ option "partAdj" - "Adjacency file from ZOBOV" string required option "partEdge" - "Boundary flag file" string required option "zone2Part" - "Particle file from ZOBOV" string required -option "mockIndex" - "Beginning index of mock particles" int required option "numVoids" - "Number of voids" int required diff --git a/c_source/zobov/jozov.c b/c_source/zobov/jozov.c index da3f3bf..1162249 100644 --- a/c_source/zobov/jozov.c +++ b/c_source/zobov/jozov.c @@ -61,8 +61,6 @@ int main(int argc,char **argv) { double *sorter, e1,maxdenscontrast; int *iord; - int mockIndex; - e1 = exp(1.)-1.; if (argc != 8) { @@ -73,7 +71,6 @@ int main(int argc,char **argv) { printf("arg4: output file containing zones in each void\n"); printf("arg5: output text file\n"); printf("arg6: Density threshold (0 for no threshold)\n"); - printf("arg7: Beginning index of mock galaxies\n\n"); exit(0); } adjfile = argv[1]; @@ -85,10 +82,6 @@ int main(int argc,char **argv) { printf("Bad density threshold.\n"); exit(0); } - if (sscanf(argv[7],"%d",&mockIndex) == 0) { - printf("Bad mock galaxy index.\n"); - exit(0); - } printf("TOLERANCE: %f\n", voltol); if (voltol <= 0.) { printf("Proceeding without a density threshold.\n"); @@ -101,8 +94,6 @@ int main(int argc,char **argv) { exit(0); } fread(&np,1, sizeof(int),adj); - if (mockIndex < 0) - mockIndex = np; printf("adj: %d particles\n", np); FF; @@ -158,7 +149,7 @@ int main(int argc,char **argv) { for (i=0;i 1e30) && i < mockIndex) { - //if ((p[i].dens < 1e-30) || (p[i].dens > 1e30)) { - // END PMS + if ((p[i].dens < 1e-30) || (p[i].dens > 1e30)) { printf("Whacked-out volume found, of particle %d: %f\n",i,p[i].dens); p[i].dens = 1.; } diff --git a/c_source/zobov/vozinit.c b/c_source/zobov/vozinit.c index 21fcac9..7a98ee6 100644 --- a/c_source/zobov/vozinit.c +++ b/c_source/zobov/vozinit.c @@ -28,19 +28,17 @@ int main(int argc, char *argv[]) { int numGuards; int b[3]; int numThreads; - int mockIndex; - if (argc != 10) { + if (argc != 9) { printf("Wrong number of arguments.\n"); printf("arg1: position file\n"); printf("arg2: buffer size (default 0.1)\n"); printf("arg3: box size\n"); - printf("arg6: number of divisions (default 2)\n"); - printf("arg7: suffix describing this run\n"); - printf("arg8: number of parallel threads\n"); - printf("arg9: location of voboz executables\n"); - printf("arg10: output directory\n"); - printf("arg11: index of mock galaxies\n\n"); + printf("arg4: number of divisions (default 2)\n"); + printf("arg5: suffix describing this run\n"); + printf("arg6: number of parallel threads\n"); + printf("arg7: location of voboz executables\n"); + printf("arg8: output directory\n"); exit(0); } posfile = argv[1]; @@ -72,11 +70,6 @@ int main(int argc, char *argv[]) { } vobozPath = argv[7]; outDir = argv[8]; - if (sscanf(argv[9],"%d",&mockIndex) != 1) { - printf("That's no mock galaxy index; try again.\n"); - exit(0); - } - /* Read the position file */ np = posread(posfile,&rfloat,1./boxsize); @@ -187,7 +180,7 @@ int main(int argc, char *argv[]) { } } fprintf(scr,"wait\n"); - fprintf(scr,"%s/voztie %d %s %s %d\n", vobozPath, numdiv,suffix, outDir, mockIndex); + fprintf(scr,"%s/voztie %d %s %s %d\n", vobozPath, numdiv,suffix, outDir); fclose(scr); sprintf(systemstr,"chmod u+x %s",scrfile); diff --git a/c_source/zobov/voztie.c b/c_source/zobov/voztie.c index ae21042..637e953 100644 --- a/c_source/zobov/voztie.c +++ b/c_source/zobov/voztie.c @@ -34,10 +34,6 @@ int main(int argc, char *argv[]) { int numRemoved = 0; - // PMS - int mockIndex; - // END PMS - if (argc != 5) { printf("Wrong number of arguments.\n"); printf("arg1: number of divisions (default 2)\n"); @@ -57,10 +53,6 @@ int main(int argc, char *argv[]) { suffix = argv[2]; outDir = argv[3]; - if (sscanf(argv[4],"%d",&mockIndex) != 1) { - printf("That's no mock galaxy index; try again.\n"); - exit(0); - } np = -1; nvpmax = -1; nvpsum = 0; @@ -92,10 +84,6 @@ int main(int argc, char *argv[]) { printf("We have %d particles to tie together.\n",np); fflush(stdout); printf("The maximum number of particles in a file is %d.\n",nvpmax); - // PMS - if (mockIndex == -1) mockIndex = np; - // END PMS - cnt_adj = (int *)malloc(np*sizeof(int)); if (cnt_adj == NULL) printf("Could not allocate cnt_adj.\n"); @@ -137,7 +125,7 @@ int main(int argc, char *argv[]) { for (p=0;p -1.) - if (fabs(vols[orig[p]]-volstemp)/volstemp > 1.5e-3 && orig[p] < mockIndex) { + if (fabs(vols[orig[p]]-volstemp)/volstemp > 1.5e-3) { printf("Inconsistent volumes for p. %d: (%10g,%10g)!\n", orig[p],vols[orig[p]],volstemp); //exit(0); @@ -183,16 +171,10 @@ int main(int argc, char *argv[]) { // PMS : remove mock galaxies and anything adjacent to a mock galaxy printf("\nRemoving mock galaxies...\n"); - // completely unlink mock particles - for (i = mockIndex; i < np; i++) { - vols[i] = 1.e-29; - adjs[i].nadj = 0; - } - // unlink particles adjacent to mock galaxies - for (i = 0; i < mockIndex; i++) { + for (i = 0; i < np; i++) { for (j = 0; j < adjs[i].nadj; j++) { - if (adjs[i].adj[j] > mockIndex) { + if (adjs[i].adj[j] > np) { //printf("KILLING %d\n", i); vols[i] = 1.e-29; adjs[i].nadj = 0; @@ -203,7 +185,7 @@ int main(int argc, char *argv[]) { } // update all other adjacencies - for (i = 0; i < mockIndex; i++) { + for (i = 0; i < np; i++) { int numAdjSaved = 0; for (j = 0; j < adjs[i].nadj; j++) { @@ -217,10 +199,6 @@ int main(int argc, char *argv[]) { adjs[i].nadj = numAdjSaved; } - printf("Removed %d mock galaxies and %d adjacent galaxies.\n", np-mockIndex, - numRemoved); - printf("There are %d galaxies remaining.\n", mockIndex-numRemoved); - // END PMS */ @@ -254,19 +232,19 @@ int main(int argc, char *argv[]) { printf("Unable to open %s\n",adjfile); exit(0); } - fwrite(&mockIndex,1, sizeof(int),adj); + fwrite(&np,1, sizeof(int),adj); /* Adjacencies: first the numbers of adjacencies, and the number we're actually going to write per particle */ // Recount the number of adjacencies after merge - for(i=0;i i) @@ -287,12 +265,8 @@ int main(int argc, char *argv[]) { printf("Unable to open %s\n",volfile); exit(0); } -// PMS - fwrite(&mockIndex,1, sizeof(int),vol); - fwrite(vols,sizeof(float),mockIndex,vol); - //fwrite(&np,1, sizeof(int),vol); - //fwrite(vols,sizeof(float),np,vol); -// END PMS + fwrite(&np,1, sizeof(int),vol); + fwrite(vols,sizeof(float),np,vol); fclose(vol); diff --git a/python_source/backend/launchers.py b/python_source/backend/launchers.py index 26ccdad..e14fe79 100644 --- a/python_source/backend/launchers.py +++ b/python_source/backend/launchers.py @@ -32,7 +32,6 @@ import shutil import glob import subprocess import sys -from pylab import figure from netCDF4 import Dataset from backend.classes import * from backend.cosmologyTools import * @@ -120,7 +119,7 @@ def launchPrep(sample, binPath, workDir=None, inputDataDir=None, sample.zRange[1], sample.zRange[0], sample.zRange[1], "all", sample.omegaM, useComoving=useComoving) - numTracers = int(open(outputDir+"/mask_index.txt", "r").read()) + numTracers = int(open("total_particles.txt", "r").read()) sample.meanPartSep = (1.*numTracers/boxVol/nbar)**(-1/3.) # flag edge galaxies with python routine for now @@ -158,7 +157,7 @@ def launchPrep(sample, binPath, workDir=None, inputDataDir=None, sample.zRange[1], sample.zRange[0], sample.zRange[1], "all", sample.omegaM, useComoving=useComoving) - numTracers = int(open(outputDir+"/mask_index.txt", "r").read()) + numTracers = int(open("total_particles.txt", "r").read()) sample.meanPartSep = (1.*numTracers/boxVol/nbar)**(-1/3.) if os.access(parmFile, os.F_OK): os.unlink(parmFile) @@ -169,8 +168,7 @@ def launchPrep(sample, binPath, workDir=None, inputDataDir=None, if os.access("comoving_distance.txt", os.F_OK): os.system("mv %s %s" % ("comoving_distance.txt", outputDir)) - if os.access("mask_index.txt", os.F_OK): - os.system("mv %s %s" % ("mask_index.txt", outputDir)) + if os.access("total_particles.txt", os.F_OK): os.system("mv %s %s" % ("total_particles.txt", outputDir)) if os.access("galaxies.txt", os.F_OK): @@ -333,10 +331,8 @@ def launchPrep(sample, binPath, workDir=None, inputDataDir=None, if os.access(parmFile, os.F_OK): os.unlink(parmFile) - if os.access("mask_index.txt", os.F_OK): - os.system("mv %s %s" % ("mask_index.txt", outputDir)) + if os.access("total_particles.txt", os.F_OK): os.system("mv %s %s" % ("total_particles.txt", outputDir)) - #os.system("mv %s %s" % ("sample_info.txt", outputDir)) # save this sample's information with open(outputDir+"/sample_info.dat", mode='wb') as output: @@ -382,15 +378,9 @@ def launchZobov(sample, binPath, outputDir=None, logDir=None, continueRun=None, if os.access(vozScript, os.F_OK): os.unlink(vozScript) - if sample.dataType == "observation": - maskIndex = open(outputDir+"/mask_index.txt", "r").read() - totalPart = open(outputDir+"/total_particles.txt", "r").read() - maxDen = mergingThreshold*float(maskIndex)/float(totalPart) - else: - maskIndex = -1 - maxDen = mergingThreshold - if numZobovDivisions == 1: - print(" WARNING! You are using a single ZOBOV division with a simulation. Periodic boundaries will not be respected!") + maxDen = mergingThreshold + if sample.dataType == "simulation" and numZobovDivisions == 1: + print(" WARNING! You are using a single ZOBOV division with a simulation. Periodic boundaries will not be respected!") if not (continueRun and jobSuccessful(logFile, "Done!\n")): for fileName in glob.glob(outputDir+"/part._"+sampleName+".*"): @@ -404,7 +394,7 @@ def launchZobov(sample, binPath, outputDir=None, logDir=None, continueRun=None, cmd = [binPath+"/vozinit", datafile, "0.1", "1.0", str(numZobovDivisions), \ "_"+sampleName, str(numZobovThreads), \ - binPath, outputDir, str(maskIndex)] + binPath, outputDir] log = open(logFile, 'w') subprocess.call(cmd, stdout=log, stderr=log) log.close() @@ -483,21 +473,24 @@ def launchZobov(sample, binPath, outputDir=None, logDir=None, continueRun=None, else: volFileToUse = outputDir+"/vol_"+sampleName+".dat" - # double-check to make sure galaxies were flagged correctly. - # if any non-edge galaxies have bad tesselations, flag them - # also, re-weight the volumes of any edge galaxies to prevent watershed - # from spilling outside of survey region + # Zobov places guard particles at the domain walls; adjacenies to these + # guard particles are not tracked in the final tesselation. + # Therefore, if a galaxy marked as "interior" is missing an adjancency, that + # means it erroneously connected to a guard particle. + # + # Re-weight the volumes of any edge galaxies to prevent watershed + # from spilling outside of survey region. if sample.dataType == "observation": log = open(logFile, 'a') - log.write("\nDouble-checking all edge flags.\n\n") + log.write("\nDouble-checking all edge flags.\n") # read in the edge flag information edgeFile = outputDir+"/galaxy_edge_flags.txt" edgeFlags = np.loadtxt(edgeFile, dtype=np.int32) # load in and count adjacenies - log.write("Loading adjancies...\n") + log.write(" Loading adjancies...\n") part = [] adjFile = outputDir+"adj_"+sampleName+".dat" with open(adjFile, mode="rb") as File: @@ -516,21 +509,18 @@ def launchZobov(sample, binPath, outputDir=None, logDir=None, continueRun=None, part[p].adjs.append(pAdj) part[pAdj].adjs.append(p) - # if the number of expected adjancies != actual adjacenies, this galaxy - # connected to a guard point on the boundary. This is fine for edge - # galaxies, but occasionally interior galaxies slip through the - # boundary handling, so we'll mark those as an edge for p in range(numPart): nHave = len(part[p].adjs) nExpected = nadjPerPart[p] - #print("Working on %d: %d vs %d" % (p, nHave, nExpected)) - if nHave != nExpected and edgeFlags[p] == 0: - log.write("Found an unflagged galaxy: %d. Correcting.\n" % (p)) + + # an interior galaxy slipped through the boundary + if edgeFlags[p] == 0 and nHave != nExpected: + log.write(" Found an unflagged galaxy: %d. Correcting.\n" % (p)) edgeFlags[p] = 4 # re-write the new edge flag information np.savetxt(edgeFile, edgeFlags, fmt="%d") - log.write("\nDone double-checking. New flags saved.\n") + log.write(" Done double-checking. New flags saved.\n") # read in the appropriate volume file log.write("Re-weighting edge galaxy volumes...\n") @@ -599,15 +589,11 @@ def launchPrune(sample, binPath, open(outputDir+"/voidDesc_"+sampleName+".out")) numVoids -= 2 + maxDen = mergingThreshold if sample.dataType == "observation": - mockIndex = open(outputDir+"/mask_index.txt", "r").read() totalPart = open(outputDir+"/total_particles.txt", "r").read() - maxDen = mergingThreshold*float(mockIndex)/float(totalPart) observationLine = " --isObservation" - #periodicLine = " --periodic=''" else: - mockIndex = -1 - maxDen = mergingThreshold observationLine = "" periodicLine = " --periodic='" + getPeriodic(sample) + "'" @@ -646,7 +632,6 @@ def launchPrune(sample, binPath, cmd += " --partEdge=" + outputDir+"galaxy_edge_flags.txt" cmd += " --extraInfo=" + outputDir+"/zobov_slice_"+str(sampleName)+".par" cmd += " --tolerance=" + str(boundaryTolerance) - cmd += " --mockIndex=" + str(mockIndex) cmd += " --maxCentralDen=" + str(maxDen) cmd += " --zMin=" + str(sample.zRange[0]) cmd += " --zMax=" + str(sample.zRange[1]) diff --git a/python_source/voidUtil/catalogUtil.py b/python_source/voidUtil/catalogUtil.py index 4cbea59..c40e067 100644 --- a/python_source/voidUtil/catalogUtil.py +++ b/python_source/voidUtil/catalogUtil.py @@ -49,7 +49,6 @@ def loadPart(sampleDir): ranges[2][0] = getattr(File, 'range_z_min') ranges[2][1] = getattr(File, 'range_z_max') isObservation = getattr(File, 'is_observation') - maskIndex = getattr(File, 'mask_index') File.close() mul = np.zeros((3)) mul[:] = ranges[:,1] - ranges[:,0] @@ -100,24 +99,8 @@ def loadPart(sampleDir): uniqueID = np.fromfile(File, dtype=np.int64,count=Np) chk = np.fromfile(File, dtype=np.int32,count=1) - ## TEST - outputting of mock particles is temporary for debugging - mockPart = [] - if isObservation == 1: - mockPartx = x[maskIndex+1:Np] - mockParty = y[maskIndex+1:Np] - mockPartz = z[maskIndex+1:Np] - - x = x[0:maskIndex]# * 100/300000 - y = y[0:maskIndex]# * 100/300000 - z = z[0:maskIndex]# * 100/300000 - RA = RA[0:maskIndex] - Dec = Dec[0:maskIndex] - redshift = redshift[0:maskIndex] - uniqueID = uniqueID[0:maskIndex] - partData = np.column_stack((x,y,z)) extraData = np.column_stack((RA,Dec,redshift,uniqueID)) - mockPart = np.column_stack((mockPartx, mockParty, mockPartz)) boxLen = mul @@ -144,7 +127,7 @@ def loadPart(sampleDir): isObservationData = isObservation == 1 - return partData, boxLen, volNorm, isObservationData, ranges, extraData, mockPart + return partData, boxLen, volNorm, isObservationData, ranges, extraData # ----------------------------------------------------------------------------- def getVolNorm(sampleDir): @@ -161,7 +144,6 @@ def getVolNorm(sampleDir): ranges[2][0] = getattr(File, 'range_z_min') ranges[2][1] = getattr(File, 'range_z_max') isObservation = getattr(File, 'is_observation') - maskIndex = getattr(File, 'mask_index') File.close() mul = np.zeros((3)) mul[:] = ranges[:,1] - ranges[:,0] @@ -193,7 +175,7 @@ def getVolNorm(sampleDir): selectionFuncFile=sample.selFunFile, useComoving=sample.useComoving) boxVol = props[0] - volNormObs = maskIndex/boxVol + volNormObs = Np/boxVol return volNorm, volNormObs @@ -457,11 +439,10 @@ def loadVoidCatalog(sampleDir, dataPortion="central", loadParticles=True, if loadParticles: print("Loading all particles...") - partData, boxLen, volNorm, isObservationData, ranges, extraData, mockPart = loadPart(sampleDir) + partData, boxLen, volNorm, isObservationData, ranges, extraData = loadPart(sampleDir) numPartTot = len(partData) catalog.numPartTot = numPartTot catalog.partPos = partData - catalog.mockPart = mockPart catalog.part = [] for i in range(len(partData)): catalog.part.append(Bunch(x = partData[i][0],