#include #include #include #include #include #include #include #include #include "TROOT.h" #include "TFile.h" #include "TH1D.h" #include "TKey.h" #include "TString.h" #include "TRegexp.h" #include "JAcoustics/JTransmission_t.hh" #include "JSupport/JMeta.hh" #include "JSystem/JGlob.hh" #include "JSystem/JStat.hh" #include "Jeep/JContainer.hh" #include "Jeep/JParser.hh" #include "Jeep/JMessage.hh" /** * \author mdejong * Auxiliary program to set disable status of transmission based on time-of-arrival histograms. */ int main(int argc, char **argv) { using namespace std; using namespace JPP; typedef JContainer< set > disable_container; typedef array array_type; vector inputFile; string outputFile; int Nmin; double Tmax_us; array_type Q = { 0.1, 0.9 }; string disableFile; int debug; try { JParser<> zap("Auxiliary program to set disable status of transmission based on time-of-arrival histograms."); zap['f'] = make_field(inputFile, "input file (output from JCanberra)."); zap['o'] = make_field(outputFile) = "disable.root"; zap['N'] = make_field(Nmin, "minimum number of entries") = 1; zap['T'] = make_field(Tmax_us, "maximal time [us]"); zap['Q'] = make_field(Q, "quantiles"); zap['!'] = make_field(disableFile, "disable transmission file") = ""; zap['d'] = make_field(debug, "debug.") = 1; zap(argc, argv); } catch(const exception &error) { FATAL(error.what() << endl); } inputFile = getFilenames(inputFile); disable_container disable; if (disableFile != "" && getFileStatus(disableFile.c_str())) { disable.load(disableFile.c_str()); } disable.comment.add(JMeta(argc, argv)); TFile out(outputFile.c_str(), "recreate"); TH1D ha("ha", NULL, 1000, 0.0, 1.0e4); TH1D hb("hb", NULL, 1000, 0.0, 1.0e3); for (vector::const_iterator i = inputFile.begin(); i != inputFile.end(); ++i) { NOTICE("Processing " << *i << endl) ; TFile in(i->c_str(), "read"); TIter iter(in.GetListOfKeys()); for (TKey* key; (key = (TKey*) iter.Next()) != NULL; ) { if (TString(key->GetName()).EndsWith(".toa")) { TH1* h1 = dynamic_cast(key->ReadObj()); if (h1 != NULL) { TString buffer(h1->GetName()); const char* regexp = "[^0-9][0-9]* [0-9]*[^0-9]"; buffer = buffer(TRegexp(regexp)); buffer = buffer(1, buffer.Length() - 2); if (buffer.Length() > 0 && buffer.IsDigit()) { JTransmission_t id; istringstream(buffer.Data()) >> id; const int N = h1->GetEntries(); if (N > 0 && N < Nmin) { NOTICE("disable: " << setw(2) << id.tx << ' ' << setw(8) << id.rx << ' ' << setw(6) << N << " < " << setw(6) << Nmin << endl); disable.insert(id); } ha.Fill(N); if (N > 0) { vector R(Q.size()); h1->GetQuantiles(Q.size(), R.data(), Q.data()); const double T_us = (*R.rbegin() - *R.begin()) * 1.0e6; hb.Fill(T_us); DEBUG("transmission: " << setw(2) << id.tx << ' ' << setw(8) << id.rx << ' ' << setw(6) << N << ' ' << FIXED(6,1) << T_us << " [us]" << endl); if (T_us > Tmax_us) { NOTICE("disable: " << setw(2) << id.tx << ' ' << setw(8) << id.rx << ' ' << FIXED(6,1) << T_us << " > " << FIXED(6,1) << Tmax_us << " [us]" << endl); disable.insert(id); } } } else { ERROR("Histogram name " << h1->GetName() << " not compatible with regular expression " << regexp << "." << endl); } } } } in.Close(); } if (disableFile != "") { disable.store(disableFile.c_str()); } out.Write(); out.Close(); }