#include #include #include "TROOT.h" #include "TFile.h" #include "TH1D.h" #include "TH2D.h" #include "TProfile.h" #include "JSupport/JMultipleFileScanner.hh" #include "JSupport/JTreeScanner.hh" #include "JROOT/JRootToolkit.hh" #include "JROOT/JGraph.hh" #include "JROOT/JManager.hh" #include "JAcoustics/JEvt.hh" #include "JAcoustics/JSupport.hh" #include "Jeep/JPrint.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" /** * \file * * Example program to compare acoustic fit results. * \author mdejong */ int main(int argc, char **argv) { using namespace std; using namespace JPP; JMultipleFileScanner inputFile; JLimit_t& numberOfEvents = inputFile.getLimit(); string outputFile; int debug; try { JParser<> zap("Example program to compare acoustic fit results."); zap['f'] = make_field(inputFile, "input file (output of JKatoomba[.sh])"); zap['n'] = make_field(numberOfEvents) = JLimit::max(); zap['o'] = make_field(outputFile) = "comparison.root"; zap['d'] = make_field(debug) = 2; zap(argc, argv); } catch(const exception &error) { FATAL(error.what() << endl); } if (inputFile.size() != 2u) { FATAL("Invalid number of input files; need 2 files for comparison." << endl); } JManager H2(new TH2D("string[%]", NULL, 500, -50.0, +50.0, 500, -50.0, +50.0)); JTreeScanner inA(inputFile[0]); JTreeScanner inB(inputFile[1]); while (inA.hasNext() && inB.hasNext()) { STATUS("event: " << setw(10) << inA.getCounter() << '\r'); DEBUG(endl); JEvt* pA = inA.next(); JEvt* pB = inB.next(); // find the same event based on the start and stop times of the events while (pA->UNIXTimeStop < pB->UNIXTimeStart && inA.hasNext()) { pA = inA.next(); } while (pB->UNIXTimeStop < pA->UNIXTimeStart && inB.hasNext()) { pB = inB.next(); } if (pA->UNIXTimeStart < pB->UNIXTimeStop && pB->UNIXTimeStart < pA->UNIXTimeStop) { for (JEvt::const_iterator iA = pA->begin(); iA != pA->end(); ++iA) { for (JEvt::const_iterator iB = pB->begin(); iB != pB->end(); ++iB) { if (iA->id == iB->id) { const double tx = (iA->tx - iB->tx) * 1.0e3; // [mrad] const double ty = (iA->ty - iB->ty) * 1.0e3; // [mrad] H2 ->Fill(tx, ty); H2[iA->id]->Fill(tx, ty); } } } } } TFile out(outputFile.c_str(), "recreate"); out << H2 << *H2; out.Write(); out.Close(); }