JTuneHV.cc

Go to the documentation of this file.

#include <string>
#include <map>
#include <fstream>
#include <sstream>

#include "JSystem/JDate.hh"

#include "JLang/JUUID.hh"

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

#include "JTools/JRange.hh"

#include "JSupport/JMeta.hh"

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

#include "JDB/JUPI.hh"
#include "JDB/JUPI_t.hh"

#include "JCalibrate/JFitToT.hh"
#include "JCalibrate/JHVInterpolator.hh"

#include "JROOT/JManager.hh"

#include "JSon/JSon.hh"

#include "TFile.h"

namespace {

  /**
   * Auxiliary class for creating a PMT-UPI look-up table from ASCII DB file.
   */
  struct JUPITable
  {
    /**
     * Initialise
     *
     * \param  fileName               file name
     */
    JUPITable(const char* const fileName)
    {
      using namespace std;
      using namespace JPP;
      
      ifstream in(fileName);
      
      for (string buffer; getline(in, buffer); ) {

        istringstream is(buffer);
      
        JUPI entry; 

        if (is >> entry.PBS >> entry.VARIANT >> entry.VERSION >> entry.SERIALNUMBER) {

          if (entry.PBS.is_valid()  &&
              entry.PBS == PBS::PMT &&
              entry.SERIALNUMBER >= 0) {
          
            number2upi[entry.PBS][entry.SERIALNUMBER] = JUPI_t(entry.PBS,
                                                               entry.VARIANT,
                                                               entry.VERSION,
                                                               entry.SERIALNUMBER);
          }
        }
      }
      
      in.close();
    }


    /**
     * Get %UPI.
     *
     * \param  pbs                   %PBS
     * \param  number                serial number
     * \return                       %UPI
     */
    JDATABASE::JUPI_t operator()(const JDATABASE::JPBS_t& pbs, const int number)
    {
      using namespace std;
      using namespace JPP;
      
      map<int, JUPI_t>::const_iterator p = number2upi[pbs].find(number);

      if (p == number2upi[pbs].end()) {

        p = number2upi[pbs].find(number);

        if (p == number2upi[pbs].end()) {
          THROW(JValueOutOfRange, "Invalid PBS / serial number " << pbs << " / " << number);
        }
      }
      
      return p->second;
    }

  protected:
    std::map<JDATABASE::JPBS_t, std::map<int, JDATABASE::JUPI_t> > number2upi;
  };
}


/**
 * \file
 *
 * Program for inter-/extrapolation of nominal high-voltage settings.
 * \author acreusot, mdejong and bjung
 */
int main(int argc, char **argv)
{
  using namespace std;
  using namespace JPP;

  string           inputFile;
  string           detectorFile;
  string           outputFile;
  string           UPIFile;

  string           login;
  string           locationID;
  string           testType;
  int              elapsedTime = 0;
  
  double           HVstepMin   = 2 * 3.14;              //!< Minimal HV step  [V]
  JRange<double>   HVrange  (-1500           + 1e-2,                
                               -80           - 1e-2);   //!< Allowed HV range [V]
  JRange<double>   gainRange(FITTOT_GAIN_MIN + 1e-2,
                             FITTOT_GAIN_MAX - 1e-2);   //!< Allowed gain fit-range
  double           gainTarget  = NOMINAL_GAIN;          //!< Target nominal gain value

  int              debug;
 

  try {

    JProperties properties;

    properties.insert(gmake_property(login));
    properties.insert(gmake_property(locationID));
    properties.insert(gmake_property(testType));
    properties.insert(gmake_property(elapsedTime));
    
    JProperties settings;
    
    settings.insert(gmake_property(HVstepMin));
    settings.insert(gmake_property(HVrange));
    settings.insert(gmake_property(gainRange));
    settings.insert(gmake_property(gainTarget));
    
    JParser<> zap("Example program to find optimal HV-settings.");

    zap['f'] = make_field(inputFile,    "input    file (ROOT format)");
    zap['o'] = make_field(outputFile,   "output   file (json format)");
    zap['a'] = make_field(detectorFile, "detector file");
    zap['b'] = make_field(UPIFile,      "PMT UPI ASCII table");
    zap['#'] = make_field(properties,   "database information")                       = JPARSER::initialised();
    zap['@'] = make_field(settings,     "interpolation settings")                     = JPARSER::initialised();
    zap['d'] = make_field(debug,        "debug")                                      = 1;

    zap(argc, argv);
  }
  catch(const exception &error) {
    FATAL(error.what() << endl);
  }

  if (login.empty() || locationID.empty()) {
    
    FATAL("Missing user information (please specify via -#login and -#locationID");
    
  } else if (testType.empty()) {

    FATAL("Missing database test type (please specify via -#testType)");
  }

  
  const JUUID& UUID = JUUID::rndm();

  JDateAndTimeFormat   format = ISO8601;  
  JDateAndTime         timer;
  
  timer->tm_sec -= elapsedTime;
  
  JDetector            detector;

  try {
    load(detectorFile, detector);
  }
  catch (const exception& error) {
    FATAL(error.what() << endl);
  }

  JModuleRouter router  (detector);
  
  JUPITable     fetchUPI(UPIFile.c_str());
  

  // Read input
  
  TFile input(inputFile.c_str(), "READ");
  
  JManager<string, TMultiGraph> manager(new TMultiGraph(), "%.HVxG");
  
  input >> manager;
  
  input.Close();
  
  
  // Find optimal HV corresponding to a nominal gain of 1.0
  // via linear interpolation in log-log scale
  
  JHVCalibration HVcalibrations;

  JHVInterpolator::setMinHVDistance(HVstepMin);
  JHVInterpolator::setHVRange      (HVrange);
  JHVInterpolator::setGainRange    (gainRange);
  
  NOTICE("Searching optimal high-voltage settings..." << endl <<
         LEFT(30) << "UPI" << CENTER(35) << "(identifer / location)" << RIGHT(35) << "Inter-/Extrapolated high-voltage [-]" << endl);

  for (JManager<string, TMultiGraph>::iterator it = manager.begin(); it != manager.end(); ++it) {

    const string&   name = it->first;
    JHVInterpolator data(*(it->second));
    
    // Retrieve upi, location and serialID
    
    const int      seppos   = name.find('.');
    const int      moduleID = stoi(name.substr(0, seppos));
    const int      tdc      = stoi(name.substr(seppos + 1, seppos + 3));

    const JModule& module   = router.getModule(moduleID);
    const JPMT&    pmt      = router.getPMT(JPMTIdentifier(moduleID, tdc));

    const JUPI_t&  pmtUPI   = fetchUPI(PBS::PMT, pmt.getID());
    const string   location = MAKE_STRING(right << FILL(4,'0') << module.getLocation().getString() << '.' <<
                                          right << FILL(2,'0') << module.getLocation().getFloor()  << '.' <<
                                          right << FILL(2,'0') << tdc);
    
    NOTICE(LEFT(30) << pmtUPI << "(a.k.a. " << name << " / " << location << "): ");
    
    if (data.interpolate(gainTarget)) {
      
      const double result = -data.getTargetHV     (gainTarget);
      const double error  =  data.getTargetHVError(gainTarget);
      
      NOTICE(right << FIXED(20,2) << result << " +/- " << FIXED(4,2) << error << endl);
      
      HVcalibrations.push_back(JHVCalibration_t(pmtUPI, OK_t, result));
      
    } else {

      WARNING("Could not find high-voltage corresponding to target gain; Setting zero." << endl);
      
      HVcalibrations.push_back(JHVCalibration_t(pmtUPI, Failed_t, 0.0));
    }
  }

  // Update graphical output

  TFile output(inputFile.c_str(), "RECREATE");

  output << manager;

  output.Close();
  
  // Write json output
  
  json js;

  js[UUID_t]                     = MAKE_STRING(UUID);
  js[User_t]                     = login;
  js[Location_t]                 = locationID;
  js[Test_t + Type_t]            = testType;
  js[StartTime_t]                = timer.toString(format);
  js[EndTime_t]                  = timer().toString(format);
  js[Tests_t]                    = json(HVcalibrations);

  ofstream ofs(outputFile.c_str());

  ofs << setw(2) << setprecision(8);
  ofs << js;

  ofs.close();
}