#include #include #include #include #include "TROOT.h" #include "TFile.h" #include "TH1D.h" #include "TH2D.h" #include "TGraph.h" #include "JROOT/JManager.hh" #include "JROOT/JGraph.hh" #include "JROOT/JRootToolkit.hh" #include "JDetector/JDetector.hh" #include "JDetector/JDetectorToolkit.hh" #include "JDynamics/JDynamics.hh" #include "JTools/JQuantile.hh" #include "JAcoustics/JEvt.hh" #include "JAcoustics/JSupport.hh" #include "JSupport/JMultipleFileScanner.hh" #include "JSupport/JMeta.hh" #include "JSupport/JSupport.hh" #include "Jeep/JPrint.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" /** * * Example program to compare dynamic position calibrations. * * \author cgatius */ int main(int argc, char **argv) { using namespace std; using namespace JPP; JMultipleFileScanner calibrationFile1; JMultipleFileScanner calibrationFile2; string detectorFile; string outputFile; double Tstep_s; int debug; try { JParser<> zap("Example program to compare dynamic position calibrations."); zap['a'] = make_field(detectorFile); zap['f'] = make_field(calibrationFile1, "output of JKatoomba[.sh]"); zap['F'] = make_field(calibrationFile2, "output of JKatoomba[.sh]"); zap['o'] = make_field(outputFile) = "compare_dynamic_positions.root"; zap['T'] = make_field(Tstep_s) = 600.0; zap['d'] = make_field(debug) = 1; zap(argc, argv); } catch(const exception &error) { FATAL(error.what() << endl); } JDetector detector; try { load(detectorFile, detector); } catch(const JException& error) { FATAL(error); } getMechanics.load(detector.getID()); JDynamics dynamics1(detector, 0.0); STATUS("loading input from file " << calibrationFile1 << "... " << flush); dynamics1.load(calibrationFile1); STATUS("OK" << endl); Double_t xmin1 = numeric_limits::max(); Double_t xmax1 = numeric_limits::lowest(); for (JDynamics::JPosition::const_iterator string = dynamics1.position.begin(); string != dynamics1.position.end(); ++string) { if (!string->second.empty()) { xmin1 = min(xmin1, string->second.getXmin()); xmax1 = max(xmax1, string->second.getXmax()); } } JDynamics dynamics2(detector, 0.0); STATUS("loading input from file " << calibrationFile1 << "... " << flush); dynamics2.load(calibrationFile2); STATUS("OK" << endl); Double_t xmin2 = numeric_limits::max(); Double_t xmax2 = numeric_limits::lowest(); for (JDynamics::JPosition::const_iterator string = dynamics2.position.begin(); string != dynamics2.position.end(); ++string) { if (!string->second.empty()) { xmin2 = min(xmin2, string->second.getXmin()); xmax2 = max(xmax2, string->second.getXmax()); } } Double_t xmin = numeric_limits::max(); Double_t xmax = numeric_limits::lowest(); //latest start time and earliest end time xmin = max(xmin1, xmin2); xmax = min(xmax1, xmax2); const JFormat_t format(4, 0, std::ios_base::fmtflags(), '0'); map HO; map HA; JManager HC(new TH1D ("M[%].comparetx", NULL, 100, 0, 2),'%', format); JManager HD(new TH1D ("N[%].comparety", NULL, 100, 0, 2),'%', format); for (JDynamics::JPosition::const_iterator string1 = dynamics1.position.begin(); string1 != dynamics1.position.end(); ++string1) { for (JDynamics::JPosition::const_iterator string2 = dynamics2.position.begin(); string2 != dynamics2.position.end(); ++string2) { if (string1->first == string2->first) { TH1D* hc = HC[string1->first]; TH1D* hd = HD[string1->first]; for (Double_t t = xmin; t <= xmax; t = t + Tstep_s) { HO[string1->first].put(t, (string1->second(t).tx - string2->second(t).tx) * 1000); HA[string1->first].put(t, (string1->second(t).ty - string2->second(t).ty) * 1000); hc->Fill((string1->second(t).tx - string2->second(t).tx) * 1000); hd->Fill((string1->second(t).ty - string2->second(t).ty) * 1000); } break; } } } TFile out(outputFile.c_str(), "recreate"); out << HC << HD; for (map::const_iterator i = HO.begin(); i != HO.end(); ++i) { out << JGraph(i->second, MAKE_CSTRING("H[" << i->first << "].comparetx")); } for (map::const_iterator i = HA.begin(); i != HA.end(); ++i) { out << JGraph(i->second, MAKE_CSTRING("K[" << i->first << "].comparety")); } }