pruneVoids now cleans out non-leaf voids; adjusted stackVoidsZero launcher to accept custom stacks

This commit is contained in:
P.M. Sutter 2013-03-30 15:24:45 -05:00
parent 765ba3b500
commit 151c605335
2 changed files with 131 additions and 85 deletions

View file

@ -51,6 +51,8 @@
#define CENTRAL_VOID 1 #define CENTRAL_VOID 1
#define EDGE_VOID 2 #define EDGE_VOID 2
using namespace std;
typedef struct partStruct { typedef struct partStruct {
float x, y, z, vol; float x, y, z, vol;
} PART; } PART;
@ -74,12 +76,25 @@ typedef struct voidStruct {
int accepted; int accepted;
int voidType; int voidType;
int parentID, numChildren; int parentID, numChildren;
bool isLeaf, hasHighCentralDen;
gsl_vector *eval; gsl_vector *eval;
gsl_matrix *evec; gsl_matrix *evec;
} VOID; } VOID;
void openFiles(string outputDir, string sampleName, string name,
int mockIndex, int numKept,
FILE** fpZobov, FILE** fpCenters,
FILE** fpCentersNoCut,
FILE** fpBarycenter, FILE** fpDistances, FILE** fpShapes,
FILE** fpSkyPositions);
void closeFiles(FILE* fpZobov, FILE* fpCenters,
FILE* fpCentersNoCut,
FILE* fpBarycenter, FILE* fpDistances, FILE* fpShapes,
FILE* fpSkyPositions);
void outputVoid(int iVoid, VOID outVoid, FILE* fpZobov, FILE* fpCenters, void outputVoid(int iVoid, VOID outVoid, FILE* fpZobov, FILE* fpCenters,
FILE* fpCenterNoCut, FILE* fpCentersNoCut,
FILE* fpSkyPositions, FILE* fpBarycenters, FILE* fpDistances, FILE* fpSkyPositions, FILE* fpBarycenters, FILE* fpDistances,
FILE* fpShapes, bool isObservation, double *boxLen); FILE* fpShapes, bool isObservation, double *boxLen);
@ -112,15 +127,12 @@ int main(int argc, char **argv) {
int i, p, p2, numPartTot, numZonesTot, dummy, iVoid, iZ; int i, p, p2, numPartTot, numZonesTot, dummy, iVoid, iZ;
int numVoids, mockIndex, numKept; int numVoids, mockIndex, numKept;
double tolerance; double tolerance;
FILE *fp, *fpZobovCentral, *fpZobovAll, *fpCentersCentral, *fpCentersAll, FILE *fp, *fpZobov, *fpCenters, *fpCentersNoCut, *fpBarycenter,
*fpCentersNoCutCentral, *fpCentersNoCutAll, *fpBarycenterCentral, *fpDistances, *fpShapes, *fpSkyPositions;
*fpBarycenterAll, *fpDistancesCentral, *fpDistancesAll,
*fpShapesCentral, *fpShapesAll, *fpSkyPositionsCentral,
*fpSkyPositionsAll;
PART *part, *voidPart; PART *part, *voidPart;
ZONE2PART *zones2Parts; ZONE2PART *zones2Parts;
VOID2ZONE *void2Zones; VOID2ZONE *void2Zones;
std::vector<VOID> voids; vector<VOID> voids;
float *temp, junk, voidVol; float *temp, junk, voidVol;
int junkInt, voidID, numPart, numZones, zoneID, partID, maxNumPart; int junkInt, voidID, numPart, numZones, zoneID, partID, maxNumPart;
int coreParticle, zoneNumPart; int coreParticle, zoneNumPart;
@ -128,6 +140,7 @@ int main(int argc, char **argv) {
float centralRad, centralDen; float centralRad, centralDen;
double nearestEdge, redshift; double nearestEdge, redshift;
char line[500], junkStr[10]; char line[500], junkStr[10];
string outputDir, sampleName, name;
int mask_index; int mask_index;
double ranges[3][2], boxLen[3], mul; double ranges[3][2], boxLen[3], mul;
double volNorm, radius; double volNorm, radius;
@ -253,6 +266,9 @@ int main(int argc, char **argv) {
voids[i-1].radius = pow(voidVol/volNorm*3./4./M_PI, 1./3.); voids[i-1].radius = pow(voidVol/volNorm*3./4./M_PI, 1./3.);
voids[i-1].accepted = 1; voids[i-1].accepted = 1;
voids[i-1].isLeaf = true;
voids[i-1].hasHighCentralDen = false;
voids[i-1].eval = gsl_vector_alloc(3); voids[i-1].eval = gsl_vector_alloc(3);
voids[i-1].evec = gsl_matrix_alloc(3,3); voids[i-1].evec = gsl_matrix_alloc(3,3);
@ -640,18 +656,23 @@ int main(int argc, char **argv) {
for (iVoid = 0; iVoid < voids.size(); iVoid++) { for (iVoid = 0; iVoid < voids.size(); iVoid++) {
if (voids[iVoid].numChildren > 0) { if (voids[iVoid].numChildren > 0) {
numAreParents++; numAreParents++;
voids[iVoid].isLeaf = false;
} else { } else {
voids[iGood++] = voids[iVoid]; //voids[iGood++] = voids[iVoid];
voids[iVoid].isLeaf = true;
} }
} }
voids.resize(iGood); //voids.resize(iGood);
printf(" 5th filter: rejected %d that are not leaf nodes\n", numAreParents); //printf(" 5th filter: rejected %d that are not leaf nodes\n", numAreParents);
for (iVoid = 0; iVoid < voids.size(); iVoid++) { for (iVoid = 0; iVoid < voids.size(); iVoid++) {
if (voids[iVoid].centralDen > args.maxCentralDen_arg) { if (voids[iVoid].centralDen > args.maxCentralDen_arg) {
voids[iVoid].accepted = -1; voids[iVoid].accepted = -1;
voids[iVoid].hasHighCentralDen = true;
numHighDen++; numHighDen++;
} else {
voids[iVoid].hasHighCentralDen = false;
} }
} }
@ -669,80 +690,46 @@ int main(int argc, char **argv) {
printf(" We have %d edge voids\n", numEdge); printf(" We have %d edge voids\n", numEdge);
printf(" We have %d central voids\n", numCentral); printf(" We have %d central voids\n", numCentral);
printf(" We have %d too high central density\n", numHighDen); printf(" We have %d too high central density\n", numHighDen);
printf(" We have %d that are not leaf nodes\n", numAreParents);
printf(" Output...\n"); printf(" Output...\n");
fpZobovCentral = fopen((std::string(args.outputDir_arg)+"/voidDesc_central_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fprintf(fpZobovCentral, "%d particles, %d voids.\n", mockIndex, numKept);
fprintf(fpZobovCentral, "Void# FileVoid# CoreParticle CoreDens ZoneVol Zone#Part Void#Zones VoidVol Void#Part VoidDensContrast VoidProb\n");
fpZobovAll = fopen((std::string(args.outputDir_arg)+"/voidDesc_all_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fprintf(fpZobovAll, "%d particles, %d voids.\n", mockIndex, numKept);
fprintf(fpZobovAll, "Void# FileVoid# CoreParticle CoreDens ZoneVol Zone#Part Void#Zones VoidVol Void#Part VoidDensContrast VoidProb\n");
fpBarycenterCentral = fopen((std::string(args.outputDir_arg)+"/barycenters_central_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fpBarycenterAll = fopen((std::string(args.outputDir_arg)+"/barycenters_all_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fpCentersCentral = fopen((std::string(args.outputDir_arg)+"/centers_central_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fprintf(fpCentersCentral, "# center x,y,z (Mpc/h), volume (normalized), radius (Mpc/h), redshift, volume (Mpc/h^3), void ID, density contrast, num part\n");
fpCentersAll = fopen((std::string(args.outputDir_arg)+"/centers_all_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fprintf(fpCentersAll, "# center x,y,z (Mpc/h), volume (normalized), radius (Mpc/h), redshift, volume (Mpc/h^3), void ID, density contrast, num part\n");
fpCentersNoCutCentral = fopen((std::string(args.outputDir_arg)+"/centers_nocut_central_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fprintf(fpCentersNoCutCentral, "# center x,y,z (Mpc/h), volume (normalized), radius (Mpc/h), redshift, volume (Mpc/h^3), void ID, density contrast, num part\n");
fpCentersNoCutAll = fopen((std::string(args.outputDir_arg)+"/centers_nocut_all_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fprintf(fpCentersNoCutAll, "# center x,y,z (Mpc/h), volume (normalized), radius (Mpc/h), redshift, volume (Mpc/h^3), void ID, density contrast, num part\n");
fpDistancesCentral = fopen((std::string(args.outputDir_arg)+"boundaryDistancesCentral_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fpDistancesAll = fopen((std::string(args.outputDir_arg)+"boundaryDistancesAll_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fpSkyPositionsCentral = fopen((std::string(args.outputDir_arg)+"/sky_positions_central_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fprintf(fpSkyPositionsCentral, "# RA, dec, redshift, radius (Mpc/h), void ID\n");
fpSkyPositionsAll = fopen((std::string(args.outputDir_arg)+"/sky_positions_all_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fprintf(fpSkyPositionsAll, "# RA, dec, redshift, radius (Mpc/h), void ID\n");
fpShapesCentral = fopen((std::string(args.outputDir_arg)+"/shapes_central_"+std::string(args.sampleName_arg)+".out").c_str(), "w");
fprintf(fpShapesCentral, "# void ID, eig(1), eig(2), eig(3), eigv(1)-x, eiv(1)-y, eigv(1)-z, eigv(2)-x, eigv(2)-y, eigv(2)-z, eigv(3)-x, eigv(3)-y, eigv(3)-z\n");
fpShapesAll = fopen((std::string(args.outputDir_arg)+"/shapes_all_"+std::string(args.sampleName_arg)+".out").c_str(), "w"); outputDir = string(args.outputDir_arg);
fprintf(fpShapesAll, "# void ID, eig(1), eig(2), eig(3), eigv(1)-x, eiv(1)-y, eigv(1)-z, eigv(2)-x, eigv(2)-y, eigv(2)-z, eigv(3)-x, eigv(3)-y, eigv(3)-z\n"); sampleName = (string(args.sampleName_arg)+".out");
name = "central";
openFiles(outputDir, sampleName, name,
mockIndex, numKept,
&fpZobov, &fpCenters, &fpCentersNoCut, &fpBarycenter,
&fpDistances, &fpShapes, &fpSkyPositions);
for (iVoid = 0; iVoid < voids.size(); iVoid++) { for (iVoid = 0; iVoid < voids.size(); iVoid++) {
if (voids[iVoid].voidType == CENTRAL_VOID) { if (voids[iVoid].voidType == CENTRAL_VOID) {
outputVoid(iVoid, voids[iVoid], fpZobovCentral, fpCentersCentral, outputVoid(iVoid, voids[iVoid], fpZobov, fpCenters,
fpCentersNoCutCentral, fpSkyPositionsCentral, fpCentersNoCut, fpSkyPositions,
fpBarycenterCentral, fpDistancesCentral, fpShapesCentral, fpBarycenter, fpDistances, fpShapes,
args.isObservation_flag, boxLen); args.isObservation_flag, boxLen);
} }
}
closeFiles(fpZobov, fpCenters, fpCentersNoCut, fpBarycenter,
fpDistances, fpShapes, fpSkyPositions);
name = "all";
openFiles(outputDir, sampleName, name,
mockIndex, numKept,
&fpZobov, &fpCenters, &fpCentersNoCut, &fpBarycenter,
&fpDistances, &fpShapes, &fpSkyPositions);
for (iVoid = 0; iVoid < voids.size(); iVoid++) {
if (voids[iVoid].voidType == EDGE_VOID || if (voids[iVoid].voidType == EDGE_VOID ||
voids[iVoid].voidType == CENTRAL_VOID) { voids[iVoid].voidType == CENTRAL_VOID) {
outputVoid(iVoid, voids[iVoid], fpZobovAll, fpCentersAll, outputVoid(iVoid, voids[iVoid], fpZobov, fpCenters,
fpCentersNoCutAll, fpSkyPositionsAll, fpCentersNoCut, fpSkyPositions,
fpBarycenterAll, fpDistancesAll, fpShapesAll, fpBarycenter, fpDistances, fpShapes,
args.isObservation_flag, boxLen); args.isObservation_flag, boxLen);
} }
} }
closeFiles(fpZobov, fpCenters, fpCentersNoCut, fpBarycenter,
fclose(fpZobovCentral); fpDistances, fpShapes, fpSkyPositions);
fclose(fpZobovAll);
fclose(fpCentersCentral);
fclose(fpCentersAll);
fclose(fpCentersNoCutCentral);
fclose(fpCentersNoCutAll);
fclose(fpBarycenterCentral);
fclose(fpBarycenterAll);
fclose(fpDistancesCentral);
fclose(fpDistancesAll);
fclose(fpShapesCentral);
fclose(fpShapesAll);
fclose(fpSkyPositionsCentral);
fclose(fpSkyPositionsAll);
clock2 = clock(); clock2 = clock();
printf(" Time: %f sec (for %d voids)\n", printf(" Time: %f sec (for %d voids)\n",
@ -755,9 +742,57 @@ int main(int argc, char **argv) {
} // end main } // end main
// ----------------------------------------------------------------------------
void openFiles(string outputDir, string sampleName, string name,
int mockIndex, int numKept,
FILE** fpZobov, FILE** fpCenters,
FILE** fpCentersNoCut,
FILE** fpBarycenter, FILE** fpDistances, FILE** fpShapes,
FILE** fpSkyPositions) {
*fpZobov = fopen((outputDir+"/voidDesc_"+name+"_"+sampleName).c_str(), "w");
fprintf(*fpZobov, "%d particles, %d voids.\n", mockIndex, numKept);
fprintf(*fpZobov, "Void# FileVoid# CoreParticle CoreDens ZoneVol Zone#Part Void#Zones VoidVol Void#Part VoidDensContrast VoidProb\n");
*fpBarycenter = fopen((outputDir+"/barycenters_"+name+"_"+sampleName).c_str(), "w");
*fpCenters = fopen((outputDir+"/centers_"+name+"_"+sampleName).c_str(), "w");
fprintf(*fpCenters, "# center x,y,z (Mpc/h), volume (normalized), radius (Mpc/h), redshift, volume (Mpc/h^3), void ID, density contrast, num part, parent ID, number of children\n");
*fpCentersNoCut = fopen((outputDir+"/centers_nocut_"+name+"_"+sampleName).c_str(), "w");
fprintf(*fpCentersNoCut, "# center x,y,z (Mpc/h), volume (normalized), radius (Mpc/h), redshift, volume (Mpc/h^3), void ID, density contrast, num part, parent ID, number of children\n");
*fpDistances = fopen((outputDir+"boundaryDistances_"+name+"_"+sampleName).c_str(), "w");
*fpSkyPositions = fopen((outputDir+"/sky_positions_"+name+"_"+sampleName).c_str(), "w");
fprintf(*fpSkyPositions, "# RA, dec, redshift, radius (Mpc/h), void ID\n");
*fpShapes = fopen((outputDir+"/shapes_"+name+"_"+sampleName).c_str(), "w");
fprintf(*fpShapes, "# void ID, eig(1), eig(2), eig(3), eigv(1)-x, eiv(1)-y, eigv(1)-z, eigv(2)-x, eigv(2)-y, eigv(2)-z, eigv(3)-x, eigv(3)-y, eigv(3)-z\n");
} // end openFiles
// ----------------------------------------------------------------------------
void closeFiles(FILE* fpZobov, FILE* fpCenters,
FILE* fpCentersNoCut,
FILE* fpBarycenter, FILE* fpDistances, FILE* fpShapes,
FILE* fpSkyPositions) {
fclose(fpZobov);
fclose(fpCenters);
fclose(fpCentersNoCut);
fclose(fpBarycenter);
fclose(fpDistances);
fclose(fpShapes);
fclose(fpSkyPositions);
} // end closeFile
// ---------------------------------------------------------------------------- // ----------------------------------------------------------------------------
void outputVoid(int iVoid, VOID outVoid, FILE* fpZobov, FILE* fpCenters, void outputVoid(int iVoid, VOID outVoid, FILE* fpZobov, FILE* fpCenters,
FILE* fpCenterNoCut, FILE* fpSkyPositions, FILE* fpCentersNoCut, FILE* fpSkyPositions,
FILE* fpBarycenters, FILE* fpDistances, FILE* fpShapes, FILE* fpBarycenters, FILE* fpDistances, FILE* fpShapes,
bool isObservation, double *boxLen) { bool isObservation, double *boxLen) {
@ -795,8 +830,8 @@ void outputVoid(int iVoid, VOID outVoid, FILE* fpZobov, FILE* fpCenters,
outCenter[2] = (outVoid.barycenter[2]-boxLen[2]/2.)*100.; outCenter[2] = (outVoid.barycenter[2]-boxLen[2]/2.)*100.;
} }
if (outVoid.accepted == 1) { if (!outVoid.hasHighCentralDen && outVoid.isLeaf) {
fprintf(fpCenters, "%.2f %.2f %.2f %.2f %.2f %.5f %.2f %d %f %d\n", fprintf(fpCenters, "%.2f %.2f %.2f %.2f %.2f %.5f %.2f %d %f %d %d %d\n",
outCenter[0], outCenter[0],
outCenter[1], outCenter[1],
outCenter[2], outCenter[2],
@ -806,10 +841,12 @@ void outputVoid(int iVoid, VOID outVoid, FILE* fpZobov, FILE* fpCenters,
4./3.*M_PI*pow(outVoid.radius, 3), 4./3.*M_PI*pow(outVoid.radius, 3),
outVoid.voidID, outVoid.voidID,
outVoid.densCon, outVoid.densCon,
outVoid.numPart); outVoid.numPart,
outVoid.parentID,
outVoid.numChildren);
} }
fprintf(fpCenterNoCut, "%.2f %.2f %.2f %.2f %.2f %.5f %.2f %d %f %d\n", fprintf(fpCentersNoCut, "%.2f %.2f %.2f %.2f %.2f %.5f %.2f %d %f %d %d %d\n",
outCenter[0], outCenter[0],
outCenter[1], outCenter[1],
outCenter[2], outCenter[2],
@ -819,7 +856,9 @@ void outputVoid(int iVoid, VOID outVoid, FILE* fpZobov, FILE* fpCenters,
4./3.*M_PI*pow(outVoid.radius, 3), 4./3.*M_PI*pow(outVoid.radius, 3),
outVoid.voidID, outVoid.voidID,
outVoid.densCon, outVoid.densCon,
outVoid.numPart); outVoid.numPart,
outVoid.parentID,
outVoid.numChildren);
fprintf(fpSkyPositions, "%.2f %.2f %.5f %.2f %d\n", fprintf(fpSkyPositions, "%.2f %.2f %.5f %.2f %d\n",
atan((outVoid.barycenter[1]-boxLen[1]/2.) / atan((outVoid.barycenter[1]-boxLen[1]/2.) /

View file

@ -420,6 +420,7 @@ def launchVoidOverlap(sample1, sample2, sample1Dir, sample2Dir,
cmd += periodicLine cmd += periodicLine
cmd += " --outfile=" + outputFile cmd += " --outfile=" + outputFile
cmd += " &> " + logFile cmd += " &> " + logFile
open("temp.par",'w').write(cmd)
os.system(cmd) os.system(cmd)
if jobSuccessful(logFile, "Done!\n"): if jobSuccessful(logFile, "Done!\n"):
@ -438,7 +439,7 @@ def launchStack(sample, stack, binPath, thisDataPortion=None, logDir=None,
zobovDir=None, zobovDir=None,
INCOHERENT=False, ranSeed=None, summaryFile=None, INCOHERENT=False, ranSeed=None, summaryFile=None,
continueRun=None, dataType=None, prefixRun="", continueRun=None, dataType=None, prefixRun="",
idList=None): idList=None, rescaleOverride=None):
sampleName = sample.fullName sampleName = sample.fullName
@ -485,13 +486,18 @@ def launchStack(sample, stack, binPath, thisDataPortion=None, logDir=None,
else: else:
nullTestFlag = "" nullTestFlag = ""
if stack.rescaleMode == "rmax": if rescaleOverride == None: rescaleOverride = stack.rescaleMode
if rescaleOverride == "rmax":
rescaleFlag = "rescale" rescaleFlag = "rescale"
else: else:
rescaleFlag = "" rescaleFlag = ""
idListFile = "idlist.temp"
if idList != None: if idList != None:
idListFlag = "idList '" + str(idList).strip('[]') + "'" idListFlag = "idList " + idListFile
idFile = open(idListFile, 'w')
for id in idList: idFile.write(str(id)+"\n")
idFile.close()
rMinToUse = 0. rMinToUse = 0.
rMaxToUse = 100000. rMaxToUse = 100000.
centralRadius = rMaxToUse centralRadius = rMaxToUse
@ -696,8 +702,9 @@ def launchStack(sample, stack, binPath, thisDataPortion=None, logDir=None,
os.system("mv %s %s" % ("normalizations.txt", voidDir+"/")) os.system("mv %s %s" % ("normalizations.txt", voidDir+"/"))
os.system("mv %s %s" % ("boundaryDistances.txt", voidDir+"/")) os.system("mv %s %s" % ("boundaryDistances.txt", voidDir+"/"))
if os.access(parmFile, os.F_OK): if os.access(idListFile, os.F_OK): os.unlink(idListFile)
os.unlink(parmFile)
if os.access(parmFile, os.F_OK): os.unlink(parmFile)
return return