further improvements to void overlap stuff

This commit is contained in:
P.M. Sutter 2012-12-18 21:41:06 -06:00
parent 3539eaacb4
commit 3dd1281cc9
2 changed files with 65 additions and 61 deletions

View file

@ -159,7 +159,6 @@ int main(int argc, char **argv) {
} }
*/ */
//for (iVoid1 = 0; iVoid1 < 1; iVoid1++) {
for (iVoid1 = 0; iVoid1 < catalog1.numVoids; iVoid1++) { for (iVoid1 = 0; iVoid1 < catalog1.numVoids; iVoid1++) {
printf(" Working on void %d of %d...\n", iVoid1, catalog1.numVoids); printf(" Working on void %d of %d...\n", iVoid1, catalog1.numVoids);
voidID1 = catalog1.voids[iVoid1].voidID; voidID1 = catalog1.voids[iVoid1].voidID;
@ -169,15 +168,12 @@ int main(int argc, char **argv) {
match = false; match = false;
for (iZ1 = 0; iZ1 < catalog1.void2Zones[voidID1].numZones; iZ1++) { for (iZ1 = 0; iZ1 < catalog1.void2Zones[voidID1].numZones; iZ1++) {
// if (match) break;
zoneID1 = catalog1.void2Zones[voidID1].zoneIDs[iZ1]; zoneID1 = catalog1.void2Zones[voidID1].zoneIDs[iZ1];
for (p1 = 0; p1 < catalog1.zones2Parts[zoneID1].numPart; p1++) { for (p1 = 0; p1 < catalog1.zones2Parts[zoneID1].numPart; p1++) {
// if (match) break;
partID1 = catalog1.zones2Parts[zoneID1].partIDs[p1]; partID1 = catalog1.zones2Parts[zoneID1].partIDs[p1];
for (iZ2 = 0; iZ2 < catalog2.void2Zones[voidID2].numZones; iZ2++) { for (iZ2 = 0; iZ2 < catalog2.void2Zones[voidID2].numZones; iZ2++) {
// if (match) break;
zoneID2 = catalog2.void2Zones[voidID2].zoneIDs[iZ2]; zoneID2 = catalog2.void2Zones[voidID2].zoneIDs[iZ2];
for (p2 = 0; p2 < catalog2.zones2Parts[zoneID2].numPart; p2++) { for (p2 = 0; p2 < catalog2.zones2Parts[zoneID2].numPart; p2++) {
@ -234,10 +230,7 @@ int main(int argc, char **argv) {
newMatch.dist = getDist(catalog1, catalog2, iVoid1, iVoid2, newMatch.dist = getDist(catalog1, catalog2, iVoid1, iVoid2,
periodicX, periodicY, periodicZ); periodicX, periodicY, periodicZ);
catalog1.voids[iVoid1].matches.push_back(newMatch); catalog1.voids[iVoid1].matches.push_back(newMatch);
//printf("MATCHED TO %d\n", iVoid2);
//catalog1.voids[iVoid1].matches.push_back(iVoid2);
} // end if match } // end if match
// break;
} }
} // end p2 } // end p2
@ -245,59 +238,41 @@ int main(int argc, char **argv) {
} // end p1 } // end p1
} // end iZ1 } // end iZ1
} // end iVoid2 } // end iVoid2
} // end iVoid1 } // end match finding
for (iVoid1 = 0; iVoid1 < catalog1.numVoids; iVoid1++) { for (iVoid1 = 0; iVoid1 < catalog1.numVoids; iVoid1++) {
sortMatches(catalog1.voids[iVoid1].matches); sortMatches(catalog1.voids[iVoid1].matches);
} }
// output properties of matches // count up significant matches
fp = fopen(args.outfile_arg, "w");
fprintf(fp, "# void ID, radius, radius ratio, common volume ratio, relative dist, num matches, num significant matches\n");
//for (iVoid1 = 0; iVoid1 < 1; iVoid1++) {
for (iVoid1 = 0; iVoid1 < catalog1.numVoids; iVoid1++) { for (iVoid1 = 0; iVoid1 < catalog1.numVoids; iVoid1++) {
// find closest match
closestMatchDist = 0.; closestMatchDist = 0.;
for (iMatch = 0; iMatch < catalog1.voids[iVoid1].matches.size(); iMatch++) { for (iMatch = 0; iMatch < catalog1.voids[iVoid1].matches.size(); iMatch++) {
rdist = catalog1.voids[iVoid1].matches[iMatch].dist;
//iVoid2 = catalog1.voids[iVoid1].matches[iVoid];
//rdist = getDist(catalog1, catalog2, iVoid1, iVoid2,
// periodicX, periodicY, periodicZ);
commonVolRatio = catalog1.voids[iVoid1].matches[iMatch].commonVol / commonVolRatio = catalog1.voids[iVoid1].matches[iMatch].commonVol /
// catalog1.voids[iVoid1].numPart; // catalog1.voids[iVoid1].numPart;
catalog1.voids[iVoid1].vol; catalog1.voids[iVoid1].vol;
if (commonVolRatio > 0.1) catalog1.voids[iVoid1].numBigMatches++; if (commonVolRatio > 0.1) catalog1.voids[iVoid1].numBigMatches++;
if (commonVolRatio > closestMatchDist) {
closestMatchID = iMatch;
closestMatchDist = commonVolRatio;
}
//if (rdist < closestMatchDist) {
// closestMatchID = iMatch;
// closestMatchDist = rdist;
//}
//printf("CENTER %d %d %e %e %e\n", iVoid1, iVoid2, rdist, rdist/catalog1.voids[iVoid1].radius, rdist/catalog2.voids[iVoid2].radius);
} }
catalog1.voids[iVoid1].numMatches = catalog1.voids[iVoid1].matches.size(); catalog1.voids[iVoid1].numMatches = catalog1.voids[iVoid1].matches.size();
}
// actually print out // output summary
std::string filename;
filename = string(args.outfile_arg);
filename = filename.append("_summary.out");
fp = fopen(filename.c_str(), "w");
//fp = fopen(args.outfile_arg, "w");
fprintf(fp, "# void ID, radius, radius ratio, common volume ratio, relative dist, num matches, num significant matches\n");
for (iVoid1 = 0; iVoid1 < catalog1.numVoids; iVoid1++) {
int voidID = catalog1.voids[iVoid1].voidID; int voidID = catalog1.voids[iVoid1].voidID;
if (catalog1.voids[iVoid1].numMatches > 0) { if (catalog1.voids[iVoid1].numMatches > 0) {
iVoid2 = catalog1.voids[iVoid1].matches[closestMatchID].matchID; iVoid2 = catalog1.voids[iVoid1].matches[0].matchID;
//iVoid2 = catalog1.voids[iVoid1].biggestMatchID;
float rRatio = catalog2.voids[iVoid2].radius / float rRatio = catalog2.voids[iVoid2].radius /
catalog1.voids[iVoid1].radius; catalog1.voids[iVoid1].radius;
float commonVolRatio = commonVolRatio = catalog1.voids[iVoid1].matches[0].commonVol /
catalog1.voids[iVoid1].matches[closestMatchID].commonVol / //catalog1.voids[iVoid1].numPart;
//catalog1.voids[iVoid1].numPart; catalog1.voids[iVoid1].vol;
catalog1.voids[iVoid1].vol; rdist = catalog1.voids[iVoid1].matches[0].dist;
rdist = getDist(catalog1, catalog2, iVoid1, iVoid2,
periodicX, periodicY, periodicZ);
rdist /= catalog1.voids[iVoid1].radius; rdist /= catalog1.voids[iVoid1].radius;
fprintf(fp, "%d %.2f %.2f %.2f %.2f %d %d\n", voidID, fprintf(fp, "%d %.2f %.2f %.2f %.2f %d %d\n", voidID,
@ -312,9 +287,37 @@ int main(int argc, char **argv) {
fprintf(fp, "%d %f 0.0 0.0 0.0 0 0\n", voidID, fprintf(fp, "%d %f 0.0 0.0 0.0 0 0\n", voidID,
catalog1.voids[iVoid1].radius); catalog1.voids[iVoid1].radius);
} }
} } // end printing
fclose(fp); fclose(fp);
// output detail
filename = string(args.outfile_arg);
filename = filename.append("_detail.out");
fp = fopen(filename.c_str(), "w");
int MAX_OUT = 10;
fprintf(fp, "# void ID, match common vol\n");
for (iVoid1 = 0; iVoid1 < catalog1.numVoids; iVoid1++) {
int voidID = catalog1.voids[iVoid1].voidID;
fprintf(fp,"%d: ", voidID);
for (iMatch = 0; iMatch < MAX_OUT; iMatch++) {
if (iMatch < catalog1.voids[iVoid].matches.size()) {
commonVolRatio = catalog1.voids[iVoid1].matches[iMatch].commonVol /
//catalog1.voids[iVoid1].numPart;
catalog1.voids[iVoid1].vol;
fprintf(fp, "%.2f ", commonVolRatio);
} else {
fprintf(fp, "0.00 ");
}
}
fprintf(fp, "\n");
} // end printing detail
fclose(fp);
printf("\nDone!\n"); printf("\nDone!\n");
return 0; return 0;
} // end main } // end main

View file

@ -33,7 +33,7 @@ globals().update(vars(parms))
if not os.access(figDir, os.F_OK): if not os.access(figDir, os.F_OK):
os.makedirs(figDir) os.makedirs(figDir)
outFileName = dataDir + "/" + dataNameBase + ".dat" outFileName = dataDir + "/" + dataNameBase #+ ".dat"
for (iSample, sampleDir) in enumerate(sampleDirList): for (iSample, sampleDir) in enumerate(sampleDirList):
if iSample == 0: continue if iSample == 0: continue
@ -50,7 +50,8 @@ for (iSample, sampleDir) in enumerate(sampleDirList):
binPath = CTOOLS_PATH+"/analysis/voidOverlap" binPath = CTOOLS_PATH+"/analysis/voidOverlap"
logFile = os.getcwd()+"/mergerTree.log" logFile = os.getcwd()+"/mergerTree.log"
stepOutputFileName = os.getcwd()+"/thisStep.out" stepOutputFileName = outFileName + "_" + sampleName + "_"
#stepOutputFileName = os.getcwd()+"/data/overlap_"
launchVoidOverlap(baseSample, sample, workDir+sampleDirList[0], launchVoidOverlap(baseSample, sample, workDir+sampleDirList[0],
workDir+sampleDir, binPath, workDir+sampleDir, binPath,
@ -59,23 +60,23 @@ for (iSample, sampleDir) in enumerate(sampleDirList):
outputFile=stepOutputFileName) outputFile=stepOutputFileName)
# attach columns to summary file # attach columns to summary file
if iSample == 1: #if iSample == 1:
os.system("cp %s %s" % (stepOutputFileName, outFileName)) # os.system("cp %s %s" % (stepOutputFileName, outFileName))
else: #else:
outFile = open("temp.out", "w") # outFile = open("temp.out", "w")
stepOutFile1 = open(outFileName, "r") # stepOutFile1 = open(outFileName, "r")
stepOutFile2 = open(stepOutputFileName, "r") # stepOutFile2 = open(stepOutputFileName, "r")
#
for line1 in stepOutFile1: # for line1 in stepOutFile1:
line1 = line1.rstrip() # line1 = line1.rstrip()
line2 = stepOutFile2.readline() # line2 = stepOutFile2.readline()
outFile.write(line1+" "+line2) # outFile.write(line1+" "+line2)
#
os.system("cp %s %s" % ("temp.out", outFileName)) # os.system("cp %s %s" % ("temp.out", outFileName))
outFile.close() # outFile.close()
stepOutFile1.close() # stepOutFile1.close()
stepOutFile2.close() # stepOutFile2.close()
#if os.access("mergerTree.log", os.F_OK): os.unlink("mergerTree.log") #if os.access("mergerTree.log", os.F_OK): os.unlink("mergerTree.log")
if os.access("temp.out", os.F_OK): os.unlink("temp.out") #if os.access("temp.out", os.F_OK): os.unlink("temp.out")
#if os.access("thisStep.out", os.F_OK): os.unlink("thisStep.out") #if os.access("thisStep.out", os.F_OK): os.unlink("thisStep.out")