Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- #include "TTree.h"
- #include "TFile.h"
- #include "TRandom.h"
- #include <iostream>
- #include <string>
- const int NUM_PIXELS = 256;
- string runNumber = ""; /* four digits; pad with 0 if necessary */
- //string runDir = "/data/macros/treesort"; /* no trailing slash */
- string runDir = "/data/files/atlas1593/gretina"; /* no trailing slash */
- string peakDir = "/data/files/bananas";
- bool verbose(false);
- bool debugBananas(false); /* going to debug with pixel 205 */
- /* interpret user-input run list. Return -1 if there's not at least 1 or if */
- /* non-digit/space characters are encountered. Return number of runs otherwise. */
- /* If the boolean "checkexistence" is set, only runs that actually exist are */
- /* loaded into the list or returned in the count. */
- int ihist(2);
- int bananas [3][256] = {{ }}; /* type (AB=0, AC=1, BC=2), pixel */
- int vertices [3][256][5] = {{{ }}}; /* type, pixel, banana # */
- float points [3][5][256][2][99] = {{{{{ }}}}}; /* type, banana, pixel, X/Y, point */
- long int bananaCounts [3][5] = {{ }}; /* track hits in each banana (up to 2^32 hits) */
- int poly = -1; /* which poly is the point in? */
- bool arraysReady(false); /* init done? */
- int bplot(-1), factor(1.0);
- /* initialise arrays for bananas */
- float pix_angle[NUM_PIXELS] = {
- 37.141, 32.895, 28.445, 23.908, 19.653, 16.391, 15.165, 16.580,
- 40.971, 36.941, 32.748, 28.546, 24.730, 21.943, 20.940, 22.100,
- 45.272, 41.556, 37.721, 33.932, 30.574, 28.196, 27.361, 28.328,
- 50.035, 46.724, 43.333, 40.027, 37.149, 35.153, 34.462, 35.263,
- 55.164, 52.326, 49.446, 46.672, 44.293, 42.667, 42.109, 42.756,
- 60.528, 58.205, 55.871, 53.651, 51.772, 50.503, 50.071, 50.572,
- 65.983, 64.180, 62.390, 60.711, 59.308, 58.371, 58.054, 58.422,
- 71.434, 70.127, 68.848, 67.665, 66.691, 66.047, 65.831, 66.082,
- 37.141, 32.895, 28.445, 23.908, 19.653, 16.391, 15.165, 16.580,
- 40.971, 36.941, 32.748, 28.546, 24.730, 21.943, 20.940, 22.100,
- 45.272, 41.556, 37.721, 33.932, 30.574, 28.196, 27.361, 28.328,
- 50.035, 46.724, 43.333, 40.027, 37.149, 35.153, 34.462, 35.263,
- 55.164, 52.326, 49.446, 46.672, 44.293, 42.667, 42.109, 42.756,
- 60.528, 58.205, 55.871, 53.651, 51.772, 50.503, 50.071, 50.572,
- 65.983, 64.180, 62.390, 60.711, 59.308, 58.371, 58.054, 58.422,
- 71.434, 70.127, 68.848, 67.665, 66.691, 66.047, 65.831, 66.082,
- 37.141, 32.895, 28.445, 23.908, 19.653, 16.391, 15.165, 16.580,
- 40.971, 36.941, 32.748, 28.546, 24.730, 21.943, 20.940, 22.100,
- 45.272, 41.556, 37.721, 33.932, 30.574, 28.196, 27.361, 28.328,
- 50.035, 46.724, 43.333, 40.027, 37.149, 35.153, 34.462, 35.263,
- 55.164, 52.326, 49.446, 46.672, 44.293, 42.667, 42.109, 42.756,
- 60.528, 58.205, 55.871, 53.651, 51.772, 50.503, 50.071, 50.572,
- 65.983, 64.180, 62.390, 60.711, 59.308, 58.371, 58.054, 58.422,
- 71.434, 70.127, 68.848, 67.665, 66.691, 66.047, 65.831, 66.082,
- 37.141, 32.895, 28.445, 23.908, 19.653, 16.391, 15.165, 16.580,
- 40.971, 36.941, 32.748, 28.546, 24.730, 21.943, 20.940, 22.100,
- 45.272, 41.556, 37.721, 33.932, 30.574, 28.196, 27.361, 28.328,
- 50.035, 46.724, 43.333, 40.027, 37.149, 35.153, 34.462, 35.263,
- 55.164, 52.326, 49.446, 46.672, 44.293, 42.667, 42.109, 42.756,
- 60.528, 58.205, 55.871, 53.651, 51.772, 50.503, 50.071, 50.572,
- 65.983, 64.180, 62.390, 60.711, 59.308, 58.371, 58.054, 58.422,
- 71.434, 70.127, 68.848, 67.665, 66.691, 66.047, 65.831, 66.082
- };
- int group[NUM_PIXELS] = {
- 1, 1, 1, 0, 0, 0, 0, 0,
- 2, 1, 1, 1, 0, 0, 0, 0,
- 2, 2, 2, 1, 1, 1, 1, 1,
- 3, 2, 2, 2, 1, 1, 1, 1,
- 3, 3, 3, 2, 2, 2, 2, 2,
- 4, 3, 3, 3, 3, 3, 3, 3,
- 4, 4, 4, 4, 3, 3, 3, 3,
- 4, 4, 4, 4, 4, 4, 4, 4,
- 1, 1, 1, 0, 0, 0, 0, 0,
- 2, 1, 1, 1, 0, 0, 0, 0,
- 2, 2, 2, 1, 1, 1, 1, 1,
- 3, 2, 2, 2, 1, 1, 1, 1,
- 3, 3, 3, 2, 2, 2, 2, 2,
- 4, 3, 3, 3, 3, 3, 3, 3,
- 4, 4, 4, 4, 3, 3, 3, 3,
- 4, 4, 4, 4, 4, 4, 4, 4,
- 1, 1, 1, 0, 0, 0, 0, 0,
- 2, 1, 1, 1, 0, 0, 0, 0,
- 2, 2, 2, 1, 1, 1, 1, 1,
- 3, 2, 2, 2, 1, 1, 1, 1,
- 3, 3, 3, 2, 2, 2, 2, 2,
- 4, 3, 3, 3, 3, 3, 3, 3,
- 4, 4, 4, 4, 3, 3, 3, 3,
- 4, 4, 4, 4, 4, 4, 4, 4,
- 1, 1, 1, 0, 0, 0, 0, 0,
- 2, 1, 1, 1, 0, 0, 0, 0,
- 2, 2, 2, 1, 1, 1, 1, 1,
- 3, 2, 2, 2, 1, 1, 1, 1,
- 3, 3, 3, 2, 2, 2, 2, 2,
- 4, 3, 3, 3, 3, 3, 3, 3,
- 4, 4, 4, 4, 3, 3, 3, 3,
- 4, 4, 4, 4, 4, 4, 4, 4
- };
- float eff[NUM_PIXELS] = {
- 0.305, 0.356, 0.409, 0.460, 0.502, 0.531, 0.541, 0.530,
- 0.347, 0.411, 0.479, 0.546, 0.604, 0.643, 0.656, 0.641,
- 0.386, 0.464, 0.548, 0.633, 0.707, 0.758, 0.776, 0.755,
- 0.419, 0.509, 0.609, 0.710, 0.801, 0.864, 0.886, 0.861,
- 0.442, 0.541, 0.651, 0.766, 0.869, 0.942, 0.967, 0.938,
- 0.450, 0.553, 0.668, 0.787, 0.896, 0.973, 1.000, 0.969,
- 0.444, 0.544, 0.655, 0.771, 0.875, 0.949, 0.975, 0.945,
- 0.422, 0.514, 0.615, 0.718, 0.811, 0.876, 0.898, 0.872,
- 0.305, 0.356, 0.409, 0.460, 0.502, 0.531, 0.541, 0.530,
- 0.347, 0.411, 0.479, 0.546, 0.604, 0.643, 0.656, 0.641,
- 0.386, 0.464, 0.548, 0.633, 0.707, 0.758, 0.776, 0.755,
- 0.419, 0.509, 0.609, 0.710, 0.801, 0.864, 0.886, 0.861,
- 0.442, 0.541, 0.651, 0.766, 0.869, 0.942, 0.967, 0.938,
- 0.450, 0.553, 0.668, 0.787, 0.896, 0.973, 1.000, 0.969,
- 0.444, 0.544, 0.655, 0.771, 0.875, 0.949, 0.975, 0.945,
- 0.422, 0.514, 0.615, 0.718, 0.811, 0.876, 0.898, 0.872,
- 0.305, 0.356, 0.409, 0.460, 0.502, 0.531, 0.541, 0.530,
- 0.347, 0.411, 0.479, 0.546, 0.604, 0.643, 0.656, 0.641,
- 0.386, 0.464, 0.548, 0.633, 0.707, 0.758, 0.776, 0.755,
- 0.419, 0.509, 0.609, 0.710, 0.801, 0.864, 0.886, 0.861,
- 0.442, 0.541, 0.651, 0.766, 0.869, 0.942, 0.967, 0.938,
- 0.450, 0.553, 0.668, 0.787, 0.896, 0.973, 1.000, 0.969,
- 0.444, 0.544, 0.655, 0.771, 0.875, 0.949, 0.975, 0.945,
- 0.422, 0.514, 0.615, 0.718, 0.811, 0.876, 0.898, 0.872,
- 0.305, 0.356, 0.409, 0.460, 0.502, 0.531, 0.541, 0.530,
- 0.347, 0.411, 0.479, 0.546, 0.604, 0.643, 0.656, 0.641,
- 0.386, 0.464, 0.548, 0.633, 0.707, 0.758, 0.776, 0.755,
- 0.419, 0.509, 0.609, 0.710, 0.801, 0.864, 0.886, 0.861,
- 0.442, 0.541, 0.651, 0.766, 0.869, 0.942, 0.967, 0.938,
- 0.450, 0.553, 0.668, 0.787, 0.896, 0.973, 1.000, 0.969,
- 0.444, 0.544, 0.655, 0.771, 0.875, 0.949, 0.975, 0.945,
- 0.422, 0.514, 0.615, 0.718, 0.811, 0.876, 0.898, 0.872
- };
- float cos[NUM_PIXELS] = {
- 0.673, 0.709, 0.742, 0.772, 0.795, 0.810, 0.815, 0.809,
- 0.703, 0.744, 0.783, 0.817, 0.845, 0.863, 0.869, 0.862,
- 0.728, 0.774, 0.818, 0.858, 0.891, 0.912, 0.919, 0.911,
- 0.748, 0.799, 0.847, 0.892, 0.929, 0.953, 0.961, 0.951,
- 0.762, 0.815, 0.867, 0.915, 0.954, 0.980, 0.989, 0.979,
- 0.766, 0.821, 0.874, 0.923, 0.964, 0.991, 1.000, 0.989,
- 0.763, 0.816, 0.868, 0.917, 0.956, 0.983, 0.992, 0.981,
- 0.750, 0.801, 0.850, 0.896, 0.932, 0.957, 0.965, 0.955,
- 0.673, 0.709, 0.742, 0.772, 0.795, 0.810, 0.815, 0.809,
- 0.703, 0.744, 0.783, 0.817, 0.845, 0.863, 0.869, 0.862,
- 0.728, 0.774, 0.818, 0.858, 0.891, 0.912, 0.919, 0.911,
- 0.748, 0.799, 0.847, 0.892, 0.929, 0.953, 0.961, 0.951,
- 0.762, 0.815, 0.867, 0.915, 0.954, 0.980, 0.989, 0.979,
- 0.766, 0.821, 0.874, 0.923, 0.964, 0.991, 1.000, 0.989,
- 0.763, 0.816, 0.868, 0.917, 0.956, 0.983, 0.992, 0.981,
- 0.750, 0.801, 0.850, 0.896, 0.932, 0.957, 0.965, 0.955,
- 0.673, 0.709, 0.742, 0.772, 0.795, 0.810, 0.815, 0.809,
- 0.703, 0.744, 0.783, 0.817, 0.845, 0.863, 0.869, 0.862,
- 0.728, 0.774, 0.818, 0.858, 0.891, 0.912, 0.919, 0.911,
- 0.748, 0.799, 0.847, 0.892, 0.929, 0.953, 0.961, 0.951,
- 0.762, 0.815, 0.867, 0.915, 0.954, 0.980, 0.989, 0.979,
- 0.766, 0.821, 0.874, 0.923, 0.964, 0.991, 1.000, 0.989,
- 0.763, 0.816, 0.868, 0.917, 0.956, 0.983, 0.992, 0.981,
- 0.750, 0.801, 0.850, 0.896, 0.932, 0.957, 0.965, 0.955,
- 0.673, 0.709, 0.742, 0.772, 0.795, 0.810, 0.815, 0.809,
- 0.703, 0.744, 0.783, 0.817, 0.845, 0.863, 0.869, 0.862,
- 0.728, 0.774, 0.818, 0.858, 0.891, 0.912, 0.919, 0.911,
- 0.748, 0.799, 0.847, 0.892, 0.929, 0.953, 0.961, 0.951,
- 0.762, 0.815, 0.867, 0.915, 0.954, 0.980, 0.989, 0.979,
- 0.766, 0.821, 0.874, 0.923, 0.964, 0.991, 1.000, 0.989,
- 0.763, 0.816, 0.868, 0.917, 0.956, 0.983, 0.992, 0.981,
- 0.750, 0.801, 0.850, 0.896, 0.932, 0.957, 0.965, 0.955
- };
- bool haveBanana(int pixel, string histType) {
- string thisPix = to_string((long long int) pixel);
- string strIhist = to_string((long long int) ihist);
- string fileName = peakDir + "/" + runNumber + "/ihist-"
- + strIhist + "/" + histType + "/" + runNumber + "-ihist"
- + strIhist + "-" + histType + "-pixel-" + thisPix + ".poly";
- struct stat st;
- if (stat(fileName.c_str(),&st) == -1) {
- return false;
- } else {
- return true;
- }
- } //end haveBanana
- bool fileExists (string input) {
- ifstream f(input.c_str());
- if (f.good()) {
- f.close();
- return true;
- } else {
- f.close();
- return false;
- }
- }
- bool dirExists (string dir) {
- string currentDir = gSystem->pwd();
- bool existence = gSystem->cd(dir.c_str());
- gSystem->cd(currentDir.c_str());
- return existence;
- } //end ifdirexists ()
- bool isNumber (string input) {
- for (int i=0; i<input.length(); i++) {
- if (isdigit((int) input[i]) == 0) {
- return false;
- }
- }
- return true;
- } // end isNumber
- bool isPercent (string input) {
- for (int i=0; i<(input.length()-1); i++) {
- if (isdigit(input[i]) == 0) {
- return false;
- }
- }
- if (input[input.length()-1] != '%') {
- return false;
- }
- }
- bool isNumberOrPercent (string input) {
- for (int i=0; i<input.length()-1; i++) {
- if (isdigit(input[i]) == 0) {
- return false;
- }
- }
- if (input[input.length()-0] != '%') {return false;}
- return true;
- }
- string removeCommas (string input) {
- string output = "";
- int counter = 0;
- for (int i=0; i<input.length(); i++) {
- if (input[i] != ',') {
- output[counter] = input[i];
- output.length()++;
- counter++;
- }
- }
- /* output[output.length] = '\0'; */
- return output;
- }
- string padRun (string input) {
- string paddedrun(input);
- if (input.length() < 4) {
- for (int iii=0; iii<(4-input.length()); iii++) {
- if (iii == input.length()) {break;} /* stop runaway */
- paddedrun = "0" + paddedrun;
- }
- }
- return paddedrun;
- }
- bool setupPeakArrays (string stringRun, int whatIhist, string ABC) {
- int linecount, banana, polys, type;
- float coordX, coordY;
- char throwaway; //throws out the delimiting character when parsing a line
- string line, linedata, isotope, fileName;
- stringstream convert;
- struct stat st;
- ifstream readFile;
- int maxPix = 256;
- for (int pix=0; pix < maxPix; pix++) {
- convert << pix;
- string thisPix = convert.str(); convert.str("");
- convert << whatIhist;
- string whichIhist = convert.str(); convert.str("");
- fileName = peakDir + "/" + stringRun + "/ihist-" + whichIhist + "/" + ABC
- + "/" + stringRun + "-ihist" + whichIhist + "-" + ABC + "-pixel-" + thisPix + ".poly";
- readFile.open(fileName.c_str());
- if (ABC == "AB") {type = 0;} //must match type settings in inPeak()
- if (ABC == "AC") {type = 1;}
- if (ABC == "BC") {type = 2;}
- //Vertex count - store vertex count in memory. Opens each file once,
- //counts vertices, stores in array vertices, then closes file
- if (!readFile.is_open()) {
- if (verbose) {cout << "No polys for " << ABC << pix << "." << endl;}
- } else {
- banana = -1; polys = 0;
- while (!readFile.eof()) {
- getline(readFile,line); istringstream iss(line);
- if ((line[0] == '#') || (line[0] == '&')) {
- linecount = -1;
- } //reset line count when we hit a header line or a vertex count line
- if (line[0] == '&') {
- banana++;
- int vertexCount;
- iss >> throwaway >> vertexCount;
- vertices[type][pix][banana] = vertexCount;
- } //end parsing vertex count (& line)
- if ((line.length() > 0) && (line[0] != '#') && (line[0] != '&') && (line[0] != '@')) {
- linecount++;
- iss >> coordX >> coordY;
- points [type][banana][pix][0][linecount] = coordX;
- points [type][banana][pix][1][linecount] = coordY;
- } //parse standard line
- if (line[0] == '@') {
- iss >> throwaway >> polys;
- bananas [type][pix] = polys;
- } //end parsing poly count (@ line)
- } //done with readFile
- readFile.close();
- readFile.clear();
- } //done with file
- } //end loading pixels into memory
- if (verbose) {
- cout << "Poly debug:" << endl << endl;
- for (int pixel=0; pixel < 256; pixel++) {
- for (int banana=0; banana < bananas [type][pixel]; banana++) {
- for (int vertex=0; vertex < vertices [type][pixel][banana]; vertex++) {
- cout << "Pixel " << pixel << " Banana: " << banana
- << " X: " << setw(7) << points [type][banana][pixel][0][vertex]
- << " Y: " << setw(7) << points [type][banana][pixel][1][vertex] << endl;
- }
- cout << "Pixel " << ABC << pixel << " bananas: " << bananas [type][pixel] << endl << endl;
- }
- }
- } //end debug output
- arraysReady = true;
- } // end setupPeakArrays
- /* This version of the pnpoly function tests whether the point is inside a given banana. */
- bool inBanana (int whatIhist, string ABC, int whatPixel, float InputX, float InputY, int banana) {
- stringstream convert;
- convert.str("");
- convert << whatPixel;
- string whichPixel = convert.str(); convert.str("");
- bool inPoly = false;
- int i, j, c(0), type;
- if (ABC == "AB") {type = 0;} /* must match type settings in setupPeakArrays() */
- if (ABC == "AC") {type = 1;}
- if (ABC == "BC") {type = 2;}
- if ((vertices[type][whatPixel][banana] != 0)) {
- for (i = 0, j = vertices[type][whatPixel][banana]-1;
- i < vertices[type][whatPixel][banana]; j = i++) {
- if (((points[type][banana][whatPixel][1][i] > InputY)
- != (points[type][banana][whatPixel][1][j] > InputY)) &&
- (InputX < (points[type][banana][whatPixel][0][j]-points[type][banana][whatPixel][0][i])
- * (InputY-points[type][banana][whatPixel][1][i]) /
- (points[type][banana][whatPixel][1][j]-points[type][banana][whatPixel][1][i])
- + points[type][banana][whatPixel][0][i]))
- {c = !c;}
- }
- if (c == 1) {
- inPoly = true;
- }
- }
- return inPoly;
- }
- string intStr (long int input) {
- stringstream convert;
- string output;
- convert << input;
- output = convert.str();
- convert.str("");
- return output;
- } //convert int to string
- string addCommas (long long int number) {
- std::string withCommas = intStr(number);
- int pos = withCommas.length() - 3;
- while (pos > 0) {
- withCommas.insert(pos,",");
- pos = pos - 3;
- }
- return withCommas;
- } //format counter with commas
- int strInt (string input) {
- //convert run number to integer
- for (int i=0; i<input.length(); i++) {
- if (isdigit(input[i]) == 0) {
- return -1;
- }
- }
- char junk[256]; int newInt;
- for (int i=0; i<input.length(); i++) {
- junk[i] = input[i];
- } //convert value to an integer
- junk[input.length()] = '\0';
- newInt = atoi(junk);
- for (int i=0; i<256; i++) {
- junk[i] = '\0';
- }
- return newInt;
- } //convert string to int
- string getUsername() {
- struct passwd *passwd;
- passwd = getpwuid(getuid());
- string userName = (string) passwd->pw_name;
- return userName;
- } //end getUsername()
- string showTime() {
- stringstream convert;
- // current date/time based on current system
- time_t now = time(0);
- // convert now to string form
- convert << ctime(&now);
- string datestamp = convert.str(); convert.str("");
- datestamp = datestamp.erase(datestamp.find_last_not_of( "\r\n\t")+1);
- return datestamp;
- } //end currentTime()
- int digitcount (long long int input) {
- stringstream convert;
- string output;
- convert << input;
- output = convert.str();
- convert.str("");
- return output.length();
- }
- /* START */
- int bananatree () {
- double entries(-999); /* will stop after this many; set to -1 to loop over all */
- double maxEntries(-999);
- long long int totalEntries(0);
- string runlist[9];
- for (int i=0; i<9; i++) {
- runlist[i] = "";
- }
- /* shenanigans; converting to int and back removes leading zeroes */
- //string bananaFolder = runNumber;
- int bananaInt = strInt(runNumber);
- //bananaFolder = intStr(bananaInt);
- bool delorean(true);
- int interval(1000000000);
- string bananaFolder = "654-667";
- /* SETUP */
- bool gotem(false);
- Int_t mult = 0;
- Int_t hitP[256] = {-1};
- Int_t hitA[256] = {-1};
- Int_t hitB[256] = {-1};
- Int_t hitC[256] = {-1};
- Int_t hitT[256] = {-1};
- Int_t angle_group[256] = {-1};
- Float_t angle[256] = {-1};
- Float_t efficiency[256] = {-1};
- Float_t G[256] = {-1};
- Float_t d[256] = {-1};
- Int_t gs_det[256] = {-1};
- Long64_t pw_time = 0;
- Long64_t gs_time[256] = {0};
- Long64_t ts_diff[256] = {0};
- Int_t aux_mult = 0;
- Int_t qdc_chan[256] = {0};
- Int_t qdc_value[256] = {0};
- Int_t tdc_chan[256] = {0};
- Int_t tdc_value[256] = {0};
- Long64_t aux_time = 0;
- Float_t t0[256] = {0.0};
- /* Float_t qdc_chan[256] = {0.0}; */
- /* Float_t qdc_val[256] = {0.0}; */
- /* Float_t tdc_chan[256] = {0.0}; */
- /* Float_t tdc_val[256] = {0.0}; */
- Int_t pw_mult(0), gs_mult(0);
- bool b0, b1, b2, b3, b4, bAny, bNone;
- bool in0[256], in1[256], in2[256], in3[256], in4[256], inNone[256];
- bool have0[256], have1[256], have2[256], have3[256], have4[256];
- bool declared;
- if (!declared) {
- gROOT->ProcessLine(".L /data/macros/treesort/gtree.C");
- gtree ttree; /* load class into memory (not doing this causes crash) */
- declared = true; /* value will persist between runs in this session */
- }
- /* get multiple run numbers */
- string runs(""); int runcount(0);
- cout << "Input runs to analyze or 'c' to continue. 'q' quits." << endl;
- bool haveRun(false);
- while (1) {
- string scratch("");
- haveRun = false;
- cout << "Enter (c)ontinue/(q)uit or (run number): ";
- cin >> scratch;
- if ((scratch.length() == 1) && (scratch[0] == 'c')) {
- if (runcount > 0) {
- break;
- } else {
- cout << "Can't continue with no runs. Select at least one or (q)uit." << endl;
- }
- }
- if ((scratch.length() == 1) && (scratch[0] == 'q')) {
- return;
- }
- /* string testDir = runDir + "/Run" + padRun(scratch); */
- string testFile= runDir + "/Run" + padRun(scratch) + ".root";
- if (fileExists(testFile)) {
- if (runcount > 0) {
- for (int w=0; w<runcount; w++) {
- if (scratch.compare(runlist[w]) == 0) {
- cout << "That run is already on the list. ";
- haveRun = true;
- }
- }
- }
- if (!haveRun) {
- runlist[runcount] = scratch;
- runcount++;
- }
- } else {
- cout << "That run doesn't exist (" << testFile << ")." << endl;
- }
- } /* end while */
- /* check all files to see how many entries they have */
- bool fileGood[15][15] = {{false}};
- string segmentName[15][15];
- long long int entryTotal(0);
- long long int allEntries(0);
- int thisrunEntries(0);
- long long int runTotal[15] = {0};
- long int fileStart[15] = {0};
- long int fileEnd[15] = {0};
- long int fileEntries[15] = {0};
- /* cout << "runcount: " << runcount << endl; */
- for (int r=0; r<runcount; r++) {
- if (runlist[r].length() == 0) {cout << "exiting" << endl; break;}
- /* string dir = runDir + "/Run" + padRun(runlist[r]); */
- string dir(runDir);
- entryTotal = 0;
- for (int iii=0; iii<15; iii++) {
- fileStart[iii] = 0; fileEnd[iii] = 0; fileEntries[iii] = 0;
- string part;
- string fileName = "";
- stringstream convert;
- convert << iii;
- part = convert.str();
- part = "_" + part;
- convert.str("");
- if (iii == 0) {part = "";} /* there's no _0 */
- fileName = dir + "/Run" + padRun(runlist[r]) + part + ".root";
- if (!(fileExists(fileName))) {
- fileGood[r][iii] = false;
- /* If the previous file was last file in a group, store total entry count */
- if (iii > 0) {
- if (fileGood[r][iii-1]) {
- runTotal[r] = entryTotal;
- }
- } else {
- runTotal[r] = entryTotal;
- }
- } else {
- fileGood[r][iii] = true;
- segmentName[r][iii] = fileName;
- TFile *entryCheck = TFile::Open(segmentName[r][iii].c_str(),"READ");
- if ((gDirectory->GetListOfKeys()->Contains("teb")) == 0) {
- fileGood[r][iii] = false;
- break; /* skip invalid file */
- } else {
- fileGood[r][iii] = true;
- TTree *thistree = (TTree*)gROOT->FindObject("teb");
- fileEntries[iii] = thistree->GetEntries();
- delete thistree;
- entryTotal = entryTotal + fileEntries[iii];
- if (iii == 0) {fileStart[iii] = 0;} else {fileStart[iii] = fileEnd[iii-1]+1;}
- fileEnd[iii] = (fileStart[iii] + fileEntries[iii])-1;
- }
- entryCheck->Close();
- entryCheck->Clear();
- delete entryCheck;
- }
- } /* existence check done (iii) */
- /* now add up total across all runs */
- if (r == (runcount-1)) {
- for (int t=0; t<15; t++) {
- allEntries = allEntries + runTotal[t];
- }
- }
- } /* end r loop */
- cout << "Entries across all runs: " << addCommas(allEntries) << endl;
- /* find out how many entries to check */
- string entriesStr = "";
- cout << endl;
- cout << "How many entries (or what percentage) to process out of "
- << addCommas(allEntries) << " entries (0 to quit)? "; cin >> entriesStr;
- while (1) {
- if (entriesStr == '0') {
- cout << endl << "Aw, sorry to see you go! Bye!" << endl;
- cout << endl;
- return;
- }
- if ((entries > allEntries) || (entries == -1)) {
- /* cout << "More or -1." << endl; */
- entries = allEntries;
- break;
- }
- if (entriesStr[entriesStr.length()-1] == '%') {
- entriesStr = entriesStr.substr(0, entriesStr.size()-1);
- if (isNumber(entriesStr)) {
- entries = atoi(entriesStr.c_str());
- }
- if (entries >= 100) {
- entries = allEntries;
- } else {
- entries = TMath::Nint((float)(allEntries * (entries/100))); /* round off answer */
- }
- break;
- } else if (isNumber(entriesStr)) {
- entries = atoi(entriesStr.c_str());
- break;
- } else {
- cout << "I don't understand your answer." << endl;
- break;
- }
- } /* got count */
- cout << "Will process " << fixed << setprecision(2) << ((float)entries/(float)allEntries)*100
- << "% of available events (" << fixed << setprecision(0) << addCommas(entries) << ")." << endl;
- interval = 1000;
- /* if ((entries > 99) && (entries <= 999)) {interval = 100;} */
- /* if ((entries > 9) && (entries <= 99)) {interval = 10;} */
- /* if ((entries > 0) && (entries <= 9)) {interval = 1;} */
- cout << "Will notify of progress every " << addCommas(interval) << " event(s)." << endl;
- cout << "A total of " << addCommas(entries) << " events will be analysed." << endl << endl;
- if (debugBananas) {
- }
- /* loop over files */
- for (int r=0; r<runcount; r++) {
- thisrunEntries = 0;
- for (int iExtra=0; iExtra<15; iExtra++) {
- if (!(fileGood[r][iExtra])) {
- /* file will be skipped if that index didn't exist */
- } else {
- TFile *thisFile = TFile::Open(segmentName[r][iExtra].c_str(),"READ");
- if (!(thisFile->IsOpen())) {
- /* cout << "Could not open file." << endl; */
- } else {
- TKey *tebKey = gDirectory->FindKey("teb");
- if (tebKey == 0) {
- cout << "Can't find tree 'teb'!" << endl;
- return;
- }
- gtree ttree; /* load class into memory (not doing this causes crash) */
- TTree *thistree = (TTree*)gROOT->FindObject("teb");
- /* Add tree-friend (is that an Ent?) */
- string friendFile = segmentName[r][iExtra];
- //friendFile = friendFile.substr(0, (friendFile.length()-5)); /* chop off ".root" */
- friendFile = friendFile.substr(runDir.length(), (friendFile.length()-5)); /* chop off ".root" */
- friendFile = friendFile + ".banana-tree-test-" + entriesStr + "-entries.root";
- friendFile = "/data/files/atlas1593/gretina/bananas" + friendFile; /* no trailing slash */
- cout << "Now saving to " << friendFile << endl;
- TFile *ff = new TFile(friendFile.c_str(),"RECREATE");
- TTree *data = new TTree("data","Particle and gamma/Doppler data");
- //data->AutoSave("500");
- data->AutoSave("SaveSelf");
- thistree->SetBranchStatus("*",0); /* turn off branches we're not interested in right now */
- thistree->SetBranchStatus("pwall.raw.hit",1);
- thistree->SetBranchStatus("pwall.raw.hit.pixel",1);
- thistree->SetBranchStatus("pwall.raw.hit.A",1);
- thistree->SetBranchStatus("pwall.raw.hit.B",1);
- thistree->SetBranchStatus("pwall.raw.hit.C",1);
- thistree->SetBranchStatus("pwall.raw.hit.T",1);
- thistree->SetBranchStatus("pwall.timestamp",1);
- thistree->SetBranchStatus("pwall.aux.qdc.channel",1);
- thistree->SetBranchStatus("pwall.aux.qdc.value",1);
- thistree->SetBranchStatus("pwall.aux.tdc.channel",1);
- thistree->SetBranchStatus("pwall.aux.tdc.value",1);
- thistree->SetBranchStatus("pwall.aux.auxtimestamp",1);
- thistree->SetBranchStatus("xtals.cc",1);
- thistree->SetBranchStatus("xtals.doppler",1);
- thistree->SetBranchStatus("xtals.crystalID",1);
- thistree->SetBranchStatus("xtals.timestamp",1);
- thistree->SetBranchStatus("xtals.t0",1);
- thistree->SetBranchAddress("pwall.raw.hit",&mult);
- thistree->SetBranchAddress("pwall.raw.hit.pixel",&hitP);
- thistree->SetBranchAddress("pwall.raw.hit.A",&hitA);
- thistree->SetBranchAddress("pwall.raw.hit.B",&hitB);
- thistree->SetBranchAddress("pwall.raw.hit.C",&hitC);
- thistree->SetBranchAddress("pwall.raw.hit.T",&hitT);
- thistree->SetBranchAddress("pwall.timestamp",&pw_time);
- thistree->SetBranchAddress("pwall.aux.auxtimestamp",&aux_time);
- thistree->SetBranchAddress("xtals.cc",&G);
- thistree->SetBranchAddress("xtals.doppler",&d);
- thistree->SetBranchAddress("xtals.crystalID",&gs_det);
- thistree->SetBranchAddress("xtals.timestamp",&gs_time);
- thistree->SetBranchAddress("xtals.t0",&t0);
- thistree->SetBranchAddress("pwall.aux.qdc.channel",&qdc_chan);
- thistree->SetBranchAddress("pwall.aux.qdc.value",&qdc_value);
- thistree->SetBranchAddress("pwall.aux.tdc.channel",&tdc_chan);
- thistree->SetBranchAddress("pwall.aux.tdc.value",&tdc_value);
- data->Branch("pw_mult", &mult);
- data->Branch("pixel", hitP, "pixel[pw_mult]/I");
- data->Branch("A", hitA, "A[pw_mult]/I");
- data->Branch("B", hitB, "B[pw_mult]/I");
- data->Branch("C", hitC, "C[pw_mult]/I");
- data->Branch("T", hitT, "T[pw_mult]/I");
- data->Branch("angle_group", angle_group, "angle_group[pw_mult]/I");
- data->Branch("angle", angle, "angle[pw_mult]/F");
- data->Branch("efficiency", efficiency, "efficiency[pw_mult]/F");
- data->Branch("pw_time", &pw_time);
- data->Branch("gs_mult", &gs_mult);
- data->Branch("gs_det", gs_det, "gs_det[gs_mult]/I");
- data->Branch("G", G, "G[gs_mult]/F");
- data->Branch("d", d, "d[gs_mult]/F");
- data->Branch("gs_time", gs_time, "gs_time[gs_mult]/L");
- data->Branch("ts_diff", ts_diff, "ts_diff[gs_mult]/L");
- data->Branch("t0", t0, "t0[gs_mult]/F");
- /* data->Branch("have0", have0, "have0[pw_mult]/O"); */
- /* data->Branch("have1", have1, "have1[pw_mult]/O"); */
- /* data->Branch("have2", have2, "have2[pw_mult]/O"); */
- /* data->Branch("have3", have3, "have3[pw_mult]/O"); */
- /* data->Branch("have4", have4, "have4[pw_mult]/O"); */
- data->Branch("aux_mult", &aux_mult);
- data->Branch("qdc_chan", qdc_chan, "qdc_chan[aux_mult]/I");
- data->Branch("qdc_value", qdc_value, "qdc_value[aux_mult]/I");
- data->Branch("tdc_chan", tdc_chan, "tdc_chan[aux_mult]/I");
- data->Branch("tdc_value", tdc_value, "tdc_value[aux_mult]/I");
- data->Branch("aux_time", &aux_time);
- data->Branch("in0", in0, "in0[pw_mult]/O");
- data->Branch("in1", in1, "in1[pw_mult]/O");
- data->Branch("in2", in2, "in2[pw_mult]/O");
- data->Branch("in3", in3, "in3[pw_mult]/O");
- data->Branch("in4", in4, "in4[pw_mult]/O");
- data->Branch("inNone", inNone, "inNone[pw_mult]/O");
- /* data->Branch("b0", &b0); */
- /* data->Branch("b1", &b1); */
- /* data->Branch("b2", &b2); */
- /* data->Branch("b3", &b3); */
- /* data->Branch("b4", &b4); */
- //data->Branch("bAny", &bAny);
- data->Branch("bNone", &bNone);
- setupPeakArrays(bananaFolder, 2, "AC");
- string bananasDir = peakDir + "/" + bananaFolder + "/ihist-" + intStr(ihist);
- if (debugBananas) {
- TCanvas *DebugCanvas = new TCanvas("DebugCanvas", "Debug Canvas", 2);
- int size = 2000; string title = "Debug canvas";
- TH2C *Debug = new TH2C("Debug",title.c_str(),size,0,size,size,0,size);
- Debug->Draw();
- }
- if ((!gotem) && (iExtra == 0)) {
- cout << "Getting " << addCommas(entries) << " events" << " ";
- }
- /* ------------------------ ANALYSE ENTRY -------------------------- */
- /* loop over entries */
- for (int thisentry=fileStart[iExtra]; thisentry<=fileEnd[iExtra]; thisentry++) {
- for (int iii=0; iii<256; iii++) {
- in0[iii] = false; in1[iii] = false; in2[iii] = false; in3[iii] = false; in4[iii] = false; inNone[iii] = false;
- have0[iii] = false; have1[iii] = false; have2[iii] = false; have3[iii] = false; have4[iii] = false;
- G[iii] = 0.0; d[iii] = 0.0; gs_det[iii] = 0; gs_time[iii] = 0; ts_diff[iii] = 0;
- qdc_chan[iii] = 0; qdc_value[iii] = 0; tdc_chan[iii] = 0; tdc_value[iii] = 0;
- angle_group[iii] = 0; angle[iii] = 0; efficiency[iii] = 0.0;
- }
- b0 = false; b1 = false; b2 = false; b3 = false; b4 = false; bAny = false; bNone = false;
- pw_mult = 0; gs_mult = 0; aux_mult = 0; pw_time = 0; aux_time = 0;
- thistree->GetEntry(thisentry); /* grab entry */
- thisrunEntries++;
- totalEntries++;
- /* figure out mults */
- for (int iii=0; iii<256; iii++) {
- if (G[iii] != 0) {gs_mult = iii+1;}
- if (qdc_chan[iii] != 0) {aux_mult = iii+1;}
- }
- pw_mult = mult;
- /* get PW geometry stats */
- for (int iii=0; iii<pw_mult; iii++) {
- angle[iii] = pix_angle[hitP[iii]];
- angle_group[iii] = group[hitP[iii]];
- efficiency[iii] = 1.0/eff[hitP[iii]];
- }
- char ticker;
- if (thisentry%2 == 0) {ticker = '-';} else {ticker = '|';}
- if (thisentry % interval == 0) {
- cout << " Entries so far: " << addCommas(thisentry) << " \r" << flush;
- }
- if ((thisentry == fileStart[iExtra]) && (thisentry < entries)) {
- int last(fileEnd[iExtra]);
- if (entries < fileEnd[iExtra]) {last = entries;}
- cout << endl << "Starting file " << segmentName[r][iExtra] /* << ": " << iExtra */
- << " (processing entries " << addCommas(fileStart[iExtra])
- << "-" << addCommas(last) << ")." << endl;
- cout << "entries: " << addCommas(entries) << " totalEntries: " << addCommas(totalEntries) << endl;
- cout << endl;
- }
- if (totalEntries > entries) {
- gotem = true;
- //return;
- break;
- } else {
- if (thisrunEntries > entries) {
- gotem = true;
- break;
- }
- }
- for (int ii=0; ii<pw_mult; ii++) {
- if (verbose) {
- cout << "ii: " << ii << " pixel: " << setw(3) << hitP[ii]
- << " A: " << setw(4) << hitA[ii]
- << " C: " << setw(4) << hitC[ii]
- << " E: " << thisentry << endl;
- }
- ts_diff[ii] = pw_time - gs_time[ii];
- if (gs_mult > 0) {
- /* check which AC bananas are defined */
- if (vertices[1][hitP[ii]][0] > 0) {have0[ii] = true;} else {have0[ii] = false;}
- if (vertices[1][hitP[ii]][1] > 0) {have1[ii] = true;} else {have1[ii] = false;}
- if (vertices[1][hitP[ii]][2] > 0) {have2[ii] = true;} else {have2[ii] = false;}
- if (vertices[1][hitP[ii]][3] > 0) {have3[ii] = true;} else {have3[ii] = false;}
- if (vertices[1][hitP[ii]][4] > 0) {have4[ii] = true;} else {have4[ii] = false;}
- if (have0[ii]) {in0[ii] = inBanana(2,"AC",hitP[ii],hitA[ii],hitC[ii],0);}
- if (have1[ii]) {in1[ii] = inBanana(2,"AC",hitP[ii],hitA[ii],hitC[ii],1);}
- if (have2[ii]) {in2[ii] = inBanana(2,"AC",hitP[ii],hitA[ii],hitC[ii],2);}
- if (have3[ii]) {in3[ii] = inBanana(2,"AC",hitP[ii],hitA[ii],hitC[ii],3);}
- if (have4[ii]) {in4[ii] = inBanana(2,"AC",hitP[ii],hitA[ii],hitC[ii],4);}
- /* how many bananas? */
- int howmany(0), inBananas(0);
- /* if (have4[ii]) {howmany = 5;} */
- /* else if (have3[ii]) {howmany = 4;} */
- /* else if (have2[ii]) {howmany = 3;} */
- /* else if (have1[ii]) {howmany = 2;} */
- /* else if (have0[ii]) {howmany = 1;} */
- if (in0[ii] && have0[ii]) {inBananas++;}
- if (in1[ii] && have1[ii]) {inBananas++;}
- if (in2[ii] && have2[ii]) {inBananas++;}
- if (in3[ii] && have3[ii]) {inBananas++;}
- if (in4[ii] && have4[ii]) {inBananas++;}
- if (inBananas == 0) {inNone[ii] = true;} else {inNone[ii] = false;}
- }
- /* if (howmany != inBananas) {inNone[ii] = true;} */
- /* if (in0[ii] && have0[ii]) {b0 = true; bAny = true;} */
- /* if (in1[ii] && have1[ii]) {b1 = true; bAny = true;} */
- /* if (in2[ii] && have2[ii]) {b2 = true; bAny = true;} */
- /* if (in3[ii] && have3[ii]) {b3 = true; bAny = true;} */
- /* if (in4[ii] && have4[ii]) {b4 = true; bAny = true;} */
- /* if ((!in0[ii] && have0[ii]) && (!in1[ii] && have1[ii]) && */
- /* (!in2[ii] && have2[ii]) && (!in3[ii] && have3[ii]) && */
- /* (!in4[ii] && have4[ii])) {inNone[ii] = true;} */
- //if (!in0[ii] && !in1[ii] && !in2[ii] && !in3[ii] && !in4[ii]) {bNone = true;}
- if (!bAny) {bNone = true;}
- if (debugBananas) {
- TMarker *mark1 = new TMarker(hitA[ii],hitC[ii], 2);
- TMarker *mark2 = new TMarker(hitA[ii],hitC[ii], 5);
- mark1->SetMarkerColor(20);
- if (in0[ii]) {mark1->SetMarkerColor(2);}
- if (in1[ii]) {mark1->SetMarkerColor(3);}
- if (in2[ii]) {mark1->SetMarkerColor(6);}
- mark2->SetMarkerColor(20);
- if (b0) {mark2->SetMarkerColor(2);}
- if (b1) {mark2->SetMarkerColor(3);}
- if (b2) {mark2->SetMarkerColor(6);}
- mark1->Draw(); mark2->Draw();
- DebugCanvas->Update();
- }
- if (verbose) {
- cout << "Arguments: "
- << "P: " << hitP[ii] << " "
- << "A: " << hitA[ii] << " "
- << "C: " << hitC[ii] << " "
- << endl;
- cout << "Results:" << " "
- << in0[ii] << "/" << b0 << " "
- << in1[ii] << "/" << b1 << " "
- << in2[ii] << "/" << b2 << " "
- << in3[ii] << "/" << b3 << " "
- << in4[ii] << "/" << b4 << " "
- << bAny << endl;
- cout << "Vertices:" << " "
- << "0: " << vertices[1][hitP[ii]][0] << " "
- << "1: " << vertices[1][hitP[ii]][1] << " "
- << "2: " << vertices[1][hitP[ii]][2] << " "
- << "3: " << vertices[1][hitP[ii]][3] << " "
- << "4: " << vertices[1][hitP[ii]][4] << " "
- << endl;
- cout << endl;
- }
- } /* end ii loop */
- /* if (!have0 && !have1 && !have2 && !have3 && !have4 && !bAny) {bNone = true;} */
- /* if ((have0 && !b0) && (have1 && !b1) && (have2 && !b2) && (have3 && !b3) && (have4 && !b4)) {bNone = true;} */
- /* if ((have0 && !b0) && (have1 && !b1) && (have2 && !b2) && (have3 && !b3) && (have4 && !b4)) { */
- /* bNone = true; */
- /* } else { */
- /* bNone = false; */
- /* } */
- data->Fill();
- } /* we still need more entries (END ENTRY LOOP HERE) */
- } /* entry loop end */
- cout << "Done." << endl;
- data->Write(); /* should save to file ff since it was last file opened */
- ff->Close();
- TFile *_file0 = new TFile(friendFile.c_str(),"READ");
- }
- }
- }
- return 1;
- }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement