JFindBadChannels.cc

Go to the documentation of this file.

#include <string>
#include <iostream>
#include <vector>

#include "JCalibrate/JParameters.hh"

#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JModuleRouter.hh"
#include "JDetector/JPMTIdentifier.hh"

#include "JSupport/JSupport.hh"
#include "JSupport/JMeta.hh"

#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"


/**
 * \file
 * Auxiliary program to update the PMT status (bad channel) in the detector file
 * based on the small ToT fraction computation (i.e. the noise to signal ratio).
 * \author acreusot
 */

int main(const int argc, const char * const argv[])
{
  using namespace std;
  using namespace JPP;

  string detectorFile;
  string inputFile;
  double stfCut;
  int debug;
  bool overwriteDetector;

  try {
    JParser<> zap("Auxiliary program to find the bad channels (PMTs) according to the small ToT fraction (i.e. the noise to signal ratio).");
    zap['a'] = make_field(detectorFile, "detector file");
    zap['f'] = make_field(inputFile, "output file of JIntegrateToT");
    zap['t'] = make_field(stfCut, "small ToT fraction cut-off") = 0.5;
    zap['d'] = make_field(debug) = 1;
    zap['A'] = make_field(overwriteDetector, "overwrite detector file with bad channel boolean.");
    zap(argc, argv);
  }
  catch(const exception& error) {
    FATAL(error.what() << endl);
  }

  // Load detector file and update it to latest version
  JDetector detector;
  try {
    load(detectorFile, detector);
  }
  catch(const JException& error) {
    FATAL(error);
  }
  detector.comment.add(JMeta(argc,argv));
  if (detector.setToLatestVersion()) {
    NOTICE("Set detector version to " << detector.getVersion() << endl);
  }
  const JModuleRouter moduleRouter(detector);

  // read PMT ToT parameters from the data file
  ifstream istr(inputFile.c_str());
  vector<JParameters> vParameters;
  for (JParameters curPar; istr >> curPar;) {
    DEBUG("Computing small ToT fraction for DOM" << curPar.domId << " - PMT " << curPar.pmtId << endl);
    vParameters.push_back(curPar);
  }
  istr.close();

  // find and store the bad channels 
  if (overwriteDetector) {
    for (vector<JParameters>::const_iterator par = vParameters.begin(); par != vParameters.end(); ++par) {
      if (moduleRouter.hasModule(par->domId)) {
        // if the histogram is empty (HRV or PMT off) then ratio is set to 1
        const double ratio  = par->signal ? par->noise/par->signal : 1;
        const bool   status = ratio > stfCut ? true : false;

        const JModuleAddress& address = moduleRouter.getAddress(par->domId);

        detector.getModule(address).getPMT(par->pmtId).setStatus(status);

        if (status == true) {
          NOTICE("bad channel: pmt " << par->pmtId << ", dom "
                 << par->domId << " with a STF of " << ratio << endl);
        } 
      }
    }
    NOTICE("Store calibration data on file " << detectorFile << endl);
    store(detectorFile, detector);
  }
}