JCompareKatoomba.cc

Go to the documentation of this file.

#include <iostream>
#include <iomanip>
 
#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<JEvt> 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<int, TH2D>  H2(new TH2D("string[%]", NULL, 500, -50.0, +50.0, 500, -50.0, +50.0));
 
  JTreeScanner<JEvt, JEvt::JEvaluator> inA(inputFile[0]);
  JTreeScanner<JEvt, JEvt::JEvaluator> 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);
 
            break;
          }
        }
      }
    }
  }
 
 
  TFile out(outputFile.c_str(), "recreate");
 
  out << H2 << *H2;
 
  out.Write();
  out.Close();
}