TRTCalibrator.cxx File Reference

#include "TRTCalibrator.h"
#include "Calibrator.h"
#include "TRT_CalibData/TrackInfo.h"
#include "TrkTrack/Track.h"
#include "TrkParameters/TrackParameters.h"
#include "TrkEventPrimitives/LocalParameters.h"
#include "EventPrimitives/EventPrimitives.h"
#include "TrkRIO_OnTrack/RIO_OnTrack.h"
#include "InDetRIO_OnTrack/TRT_DriftCircleOnTrack.h"
#include "InDetPrepRawData/TRT_DriftCircle.h"
#include "TrkPrepRawData/PrepRawData.h"
#include "AtlasDetDescr/AtlasDetectorID.h"
#include "InDetIdentifier/TRT_ID.h"
#include "TRT_ConditionsData/RtRelation.h"
#include "TRT_ConditionsData/BasicRtRelation.h"
#include "TRT_ConditionsData/DinesRtRelation.h"
#include "TRT_ConditionsData/BinnedRtRelation.h"
#include "CommissionEvent/ComTime.h"
#include "TRT_ReadoutGeometry/TRT_DetectorManager.h"
#include <TNtuple.h>
#include <TFile.h>
#include <fstream>

Go to the source code of this file.

Functions
bool TRTCalibrator::calibrate	ATLAS_NOT_THREAD_SAFE ()
	Install fatal handler with default options. More...

Function Documentation

ATLAS_NOT_THREAD_SAFE()

bool TRTCalibrator::calibrate ATLAS_NOT_THREAD_SAFE ( )

inline

Install fatal handler with default options.

This is meant to be easy to call from python via ctypes.

Install fatal handler with default options.

getLorentzAngle() Read LorentzAngle from HIST and write out into local DB

getBSErrors() Read BSErrors from Monitoring HIST and write out into local DB

getEfficiency() Read Efficiency from Monitoring HIST and write out into local DB

getRawOccupancy() Read RawOccupancy from Monitoring HIST and write out into local DB

getNoiseOccupancy() Read NoiseOccupancy from HIST and write out into local DB

getNoisyStrip() Find noisy strips from hitmaps and write out into xml/db formats

Definition at line 309 of file TRTCalibrator.cxx.

                                                     {
 
  float sid;
  Identifier ident;
  const InDetDD::TRT_BaseElement* strawelement;
  const InDetDD::TRT_BarrelElement* barrelelement;
  const InDetDD::TRT_EndcapElement* endcapelement;
  int nTRThist=0, ndethist=0, nlayhist=0, nmodhist=0, nbrdhist=0, nchphist=0, nstwhist=0;
  int nTRThistAr=0, ndethistAr=0, nlayhistAr=0, nmodhistAr=0, nbrdhistAr=0, nchphistAr=0, nstwhistAr=0;
  databundle hitdata;
  
  std::string infile=m_hittuple;
  std::string outfile="calibout";
  
  std::string detname="Detector";
  std::string layname="Layer";
  if (!m_SplitBarrel) {
    detname="WholeBarrel";
    layname="WholeBarrelLayer";
  }    
 
  //create Calibrator objects for each sub-level
  Calibrator TRT(0,"TRT",m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Detector(1,detname.data(),m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Layer(2,layname.data(),m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Module(3,"Module",m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Board(4,"Board",m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Chip(5,"Chip",m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Straw(6,"Straw",m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
 
 
  // make a map of calibrators (so that they can be looped over)
  std::map<std::string,Calibrator*> calibrators;
  calibrators["TRT"]=&TRT;
  calibrators["Detector"]=&Detector;
  calibrators["Layer"]=&Layer;
  calibrators["Module"]=&Module;
  calibrators["Board"]=&Board;
  calibrators["Chip"]=&Chip;
  calibrators["Straw"]=&Straw;
 
  //use default configurations unless "user" option is used
  if(m_calsub!="user"){
    m_doRt=m_config[m_calsub].CalibrateRt;
    m_doT0=m_config[m_calsub].CalibrateT0;
    m_doRes=m_config[m_calsub].FitResidual;
    m_beQuiet=m_config[m_calsub].NoHistograms;
    m_doOutPrint=m_config[m_calsub].PrintT0Out;
    m_doRtPrint=m_config[m_calsub].PrintRtOut;
    m_doLogPrint=m_config[m_calsub].PrintLog;    
  }
 
  //configure the calibrators
  for (const std::string& s : m_doRt)        calibrators[s]->dort=true;
  for (const std::string& s : m_doT0)        calibrators[s]->dot0=true;
  for (const std::string& s : m_doRes)       calibrators[s]->dores=true;
  for (const std::string& s : m_beQuiet)     calibrators[s]->bequiet=true;
  for (const std::string& s : m_doOutPrint)  calibrators[s]->printt0=true;
  for (const std::string& s : m_doRtPrint)   calibrators[s]->printrt=true;
  for (const std::string& s : m_doLogPrint)  calibrators[s]->printlog=true;
  for (const std::string& s : m_useBoardRef) calibrators[s]->usebref=true;
  for (std::pair<const std::string, Calibrator*>& p : calibrators) {
    Calibrator* calib = p.second;
    calib->usep0=m_useP0;
    calib->floatp3=m_floatP3;
    calib->useshortstraws=m_DoShortStrawCorrection;
 
    if (m_doRt.size()==0) calib->nort=true;
    if (m_doT0.size()==0) calib->not0=true;
    
    std::set<int> selset;
    GetSubLevels(m_selstring,calib->level+1,&selset); //get the set of selections for level "level+1" from the m_selstring ...
    for (std::set<int>::iterator imod=selset.begin(); imod!=selset.end(); ++imod){
      calib->selection = selset;//... and configure the calibrator with them
    }
  }     
 
  // FOR ARGON:
  Calibrator TRT_Ar     (0,"TRT_Ar"     ,m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Detector_Ar(1,"Detector_Ar",m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Layer_Ar   (2,"Layer_Ar"   ,m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Module_Ar  (3,"Module_Ar"  ,m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Board_Ar   (4,"Board_Ar"   ,m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Chip_Ar    (5,"Chip_Ar"    ,m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
  Calibrator Straw_Ar   (6,"Straw_Ar"   ,m_mint0,m_minrt,m_rtrel,m_rtbinning,m_t0offset);
 
  // make a map of calibrators (so that they can be looped over)
  std::map<std::string,Calibrator*>     calibratorsAr;
  calibratorsAr["TRT"]                  =&TRT_Ar;
  calibratorsAr["Detector"]             =&Detector_Ar;
  calibratorsAr["Layer"]                =&Layer_Ar;
  calibratorsAr["Module"]               =&Module_Ar;
  calibratorsAr["Board"]                =&Board_Ar;
  calibratorsAr["Chip"]                 =&Chip_Ar;
  calibratorsAr["Straw"]                =&Straw_Ar;
 
  //configure the calibrators
  for (const std::string& s : m_doRt)        calibratorsAr[s]->dort=true;
  for (const std::string& s : m_doT0)        calibratorsAr[s]->dot0=true;
  for (const std::string& s : m_doRes)       calibratorsAr[s]->dores=true;
  for (const std::string& s : m_beQuiet)     calibratorsAr[s]->bequiet=true;
  for (const std::string& s : m_doOutPrint)  calibratorsAr[s]->printt0=true;
  for (const std::string& s : m_doRtPrint)   calibratorsAr[s]->printrt=true;
  for (const std::string& s : m_doLogPrint)  calibratorsAr[s]->printlog=true;
  for (const std::string& s : m_useBoardRef) calibratorsAr[s]->usebref=true;
  for (std::pair<const std::string, Calibrator*>& p : calibratorsAr) {
    Calibrator* calib = p.second;
    calib->usep0=m_useP0;
    calib->floatp3=m_floatP3;
    calib->useshortstraws=m_DoShortStrawCorrection;
 
    if (m_doRt.size()==0) calib->nort=true;
    if (m_doT0.size()==0) calib->not0=true;
 
    std::set<int> selset;
    GetSubLevels(m_selstring,calib->level+1,&selset); //get the set of selections for level "level+1" from the m_selstring ... 
    for (std::set<int>::iterator imod=selset.begin(); imod!=selset.end(); ++imod){
      calib->selection = selset;//... and configure the calibrator with them
    }
  }
 
 
 
  m_options = '_' + TRT.GetOptString() + '_' + Detector.GetOptString() + '_' + Layer.GetOptString() + '_' + Module.GetOptString() + '_' + Board.GetOptString() + '_' + Chip.GetOptString() + '_' + Straw.GetOptString();
  if(m_calsub!="user") m_selstring = '_' + TRT.GetSelString() + '_' + Detector.GetSelString() + '_' + Layer.GetSelString() + '_' + Module.GetSelString() + '_' + Board.GetSelString() + '_' + Chip.GetSelString() + '_' + Straw.GetSelString();
  
  //print some info
  ATH_MSG_INFO(" ");
  ATH_MSG_INFO( "INPUT FILE       : " << infile ); 
  ATH_MSG_INFO( "OUTPUT FILE      : " << outfile ); 
  ATH_MSG_INFO( "SELECTION STRING : " << m_selstring ); 
  ATH_MSG_INFO( "OPTION STRING    : " << m_options ); 
  ATH_MSG_INFO( "RT RELATION      : " << m_rtrel ); 
  ATH_MSG_INFO( "MIN STATISTICS   : RT=" << m_minrt << ", T0=" << m_mint0 ); 
  ATH_MSG_INFO( "T0 OFFSET        : " << m_t0offset ); 
  ATH_MSG_INFO(" ");
  ATH_MSG_INFO( TRT.PrintInfo() );
  ATH_MSG_INFO( Detector.PrintInfo() );
  ATH_MSG_INFO( Layer.PrintInfo() );
  ATH_MSG_INFO( Module.PrintInfo() );
  ATH_MSG_INFO( Board.PrintInfo() );
  ATH_MSG_INFO( Chip.PrintInfo() );
  ATH_MSG_INFO( Straw.PrintInfo() );
  ATH_MSG_INFO(" ");
 
  if (m_DoArXenonSep){
  //print some info
    ATH_MSG_INFO(" ");
    ATH_MSG_INFO( "Setup for the ARGON Straws: " );
    ATH_MSG_INFO( "INPUT FILE       : " << infile );
    ATH_MSG_INFO( "OUTPUT FILE      : " << outfile );
    ATH_MSG_INFO( "SELECTION STRING : " << m_selstring );
    ATH_MSG_INFO( "OPTION STRING    : " << m_options );
    ATH_MSG_INFO( "RT RELATION      : " << m_rtrel );
    ATH_MSG_INFO( "MIN STATISTICS   : RT=" << m_minrt << ", T0=" << m_mint0 );
    ATH_MSG_INFO( "T0 OFFSET        : " << m_t0offset );
    ATH_MSG_INFO( " ");
    ATH_MSG_INFO( TRT_Ar.PrintInfo() );
    ATH_MSG_INFO( Detector_Ar.PrintInfo() );
    ATH_MSG_INFO( Layer_Ar.PrintInfo() );
    ATH_MSG_INFO( Module_Ar.PrintInfo() );
    ATH_MSG_INFO( Board_Ar.PrintInfo() );
    ATH_MSG_INFO( Chip_Ar.PrintInfo() );
    ATH_MSG_INFO( Straw_Ar.PrintInfo() );
    ATH_MSG_INFO(" ");
  }
  //read the chip reference t0 values from finedelays.txt
  int rbrd, rchp, rdet, dum;
  float rt0;
  std::map<std::string,float> reft0map;
  std::ifstream t0ref("finedelays.txt",std::ios::in);
  if(t0ref.is_open()) ATH_MSG_INFO(" Opened finedelays.txt ");
  for(int iref=0;iref<208;iref++){
    t0ref >> rbrd >> rchp >> rdet >> dum >> rt0;
    reft0map[std::string(Form("_%i_%i_%i",rdet,rbrd,rchp))]=rt0;
  }
 
  t0ref.close();
 
  
  //set type of rt-relation
  bool isdines = (m_rtrel.value()).find("dines")!=std::string::npos;
  bool isbinned = (m_rtrel.value()).find("binned")!=std::string::npos;
  int rtint;
  if (isdines) rtint=2;
  else if (isbinned) rtint=1;
  else rtint=0;
  ATH_MSG_INFO(" Rt relation is: " << rtint);
  
  int npop,isid;
  std::map<std::string,epdata> ephasemap;
 
  // open the output histogram file
   std::unique_ptr<TFile> histfile(TFile::Open("calibout.root","RECREATE"));
 
  if(!histfile || histfile->IsZombie()) {
      ATH_MSG_ERROR ("Failed to open calibout.root ");
      return false;
  } else {      
      ATH_MSG_INFO ("Opened calibout.root ");
  }
 
    std::ifstream myFile (infile.data(), std::ios::in | std::ios::binary);
    ATH_MSG_INFO( " Opened " << infile << " as binary histogram file " );        
    int ihist=0;
    int ihistAr=0;
 
    while(true){
 
      //read a binary histogram
      //dangerous cast reproduces pre-existing c style cast
      myFile.read (reinterpret_cast<char*>(&npop),sizeof(int)); //number of populated bins
      if (myFile.eof()) break;
      int* chist=new int[2*npop+2]; //the histogram
      //dangerous cast reproduces pre-existing c style cast
      if (npop>0) myFile.read (reinterpret_cast<char*>(chist+2), sizeof(int)*2*npop);
      //dangerous cast reproduces pre-existing c style cast
      myFile.read (reinterpret_cast<char*>(&isid),sizeof(int)); //the straw id 
      sid = (float)isid;
      if(sid<0) continue;
 
      chist[0]=npop;
      chist[1]=isid;
      //get the straw address (barrel/ec, layer, module, ...) based on the straw identifier
      ident=(Identifier)(isid);
      int chip=0,board=-1;
      m_neighbourSvc->getChip(ident,chip); 
      if(abs(m_TRTID->barrel_ec(ident))<2){
        board=m_neighbourSvc->chipToBoardBarrel(chip,m_TRTID->layer_or_wheel(ident));
      }
      else if (chip<12) {
        board=0;
      }
      else {
        chip=chip-20;
        board=1;
      }
      if (m_SplitBarrel) hitdata.det=(int)m_TRTID->barrel_ec(ident); 
      else hitdata.det=abs((int)m_TRTID->barrel_ec(ident));
      if(hitdata.det!=1 && hitdata.det!=-1 && hitdata.det!=2 && hitdata.det!=-2) {
        ATH_MSG_WARNING( " Invalid Identifier read : " << ident << "  detector decoded to be : " << hitdata.det); 
        continue; 
      }
      hitdata.lay=(int)m_TRTID->layer_or_wheel(ident);
      if( (hitdata.det==1 || hitdata.det==-1) && (hitdata.lay<0 || hitdata.lay>2) ) {
        ATH_MSG_WARNING( " Invalid Identifier read : " << ident << "  barrel layer decoded to be : " << hitdata.lay); 
        continue; 
      } else if( (hitdata.det==2 || hitdata.det==-2) && (hitdata.lay<0 || hitdata.lay>13) ) {
        ATH_MSG_WARNING( " Invalid Identifier read : " << ident << "  endcap layer decoded to be : " << hitdata.lay); 
        continue; 
      }
      hitdata.mod=(int)m_TRTID->phi_module(ident); 
      if( hitdata.mod < 0 || hitdata.mod > 31 ) {
        ATH_MSG_WARNING( " Invalid Identifier read : " << ident << "  phi-module decoded to be : " << hitdata.lay); 
        continue; 
      }
 
      hitdata.brd=board; 
      hitdata.chp=chip; 
      hitdata.stl=(int)m_TRTID->straw_layer(ident);
      if( hitdata.stl < 0 || hitdata.stl > m_TRTID->straw_layer_max(ident) ) {
        ATH_MSG_WARNING( " Invalid Identifier read : " << ident << "  straw-layer decoded to be : " << hitdata.stl); 
        continue; 
      }
 
      hitdata.stw=(int)m_TRTID->straw(ident);
      if( hitdata.stw < 0 || hitdata.stw > m_TRTID->straw_max(ident) ) {
        ATH_MSG_WARNING( " Invalid Identifier read : " << ident << "  straw decoded to be : " << hitdata.stw); 
        continue; 
      }
 
      hitdata.sid=isid;
      
      //get the old rt parameters based on straw identifier
      const float* pcal ;
      std::vector<float> rvalues ;
      const float defaultpcal[] = {0,0,0,0} ;
      if (isdines){
        const TRTCond::DinesRtRelation* rtr = dynamic_cast<const TRTCond::DinesRtRelation*>(m_trtcaldbTool->getRtRelation(ident)) ;
        pcal = rtr ? rtr->cal() : defaultpcal ;
      }
 
      else {
        const TRTCond::BasicRtRelation* rtr = dynamic_cast<const TRTCond::BasicRtRelation*>(m_trtcaldbTool->getRtRelation(ident)) ;
        pcal = rtr ? rtr->cal() : defaultpcal ;
      }
      hitdata.rtpar[0]=pcal[0];
      hitdata.rtpar[1]=pcal[1];
      hitdata.rtpar[2]=pcal[2];
      hitdata.rtpar[3]=pcal[3];
      
 
      //build map keys
      MakeBDKeys(hitdata.det, hitdata.lay, hitdata.mod, hitdata.brd, hitdata.chp, hitdata.sid);
      
      //make the level hierachy dictionaries
      m_trt.t[m_Tkey].d[m_Dkey].l[m_Lkey].m[m_Mkey].b[m_Bkey].c[m_Ckey].s[m_Skey].z=0;
      m_trt_acc.t[m_Tkey].d[m_Dkey_acc].l[m_Lkey_acc].m[m_Mkey].b[m_Bkey].c[m_Ckey].s[m_Skey].z=0;
      
      //populate the hit data structure to be added 
 
 
      hitdata.ievt=0; 
      hitdata.tres=0; 
      hitdata.weight=1.0;
      hitdata.res=0; 
      hitdata.t=0; 
      hitdata.r=0;
      hitdata.t0=m_trtcaldbTool->getT0(ident);
 
      strawelement = m_trtmanager->getElement(ident);
      if(hitdata.det==1 || hitdata.det==-1) {
    barrelelement=static_cast<const InDetDD::TRT_BarrelElement*>(strawelement);
        hitdata.x=(barrelelement->center(ident)).x();
        hitdata.y=(barrelelement->center(ident)).y();
        hitdata.z=(barrelelement->center(ident)).z();
      }else{
    endcapelement=static_cast<const InDetDD::TRT_EndcapElement*>(strawelement);
         hitdata.x=(endcapelement->center(ident)).x();
         hitdata.y=(endcapelement->center(ident)).y();
         hitdata.z=(endcapelement->center(ident)).z();
      }
      
      //in the short straws corrections, autodetect if it was applied on the previous step
      if ( m_DoShortStrawCorrection &&  hitdata.lay==0 && hitdata.stl<9){
        //If correction was done in ctes in db (readed), undo the correction: This fixes problem on t0 averaging....
        double t0test1=m_trtcaldbTool->getT0((Identifier)301998432);
        double t0test2=m_trtcaldbTool->getT0((Identifier)302001504);
        if (t0test1 != t0test2) hitdata.t0+=0.75; //short straw compensation
      }
      
      //the reference t0 from finedelays
      hitdata.rt0=reft0map[std::string(Form("_%i_%i_%i",hitdata.det,hitdata.brd,hitdata.chp))];
     
           // Prepare for Xe-Ar mixed conditions:
            int isArgonStraw = 0;
            if (m_TRTStrawSummaryTool->getStatusHT(ident, Gaudi::Hive::currentContext()) != TRTCond::StrawStatus::Good) {
                  isArgonStraw = 1;
            }
 
 
      //add histogram to the Calibrators (A and C side separated)
   if(!m_DoArXenonSep){  
  
      nTRThist += TRT.AddHit(m_Tkey,hitdata,chist,true);
      if (Detector.CheckSelection(hitdata.det)) { //only add the histogram if it is in the selection
        if (m_SplitBarrel) ndethist += Detector.AddHit(m_Dkey,hitdata,chist,true);
        else ndethist += Detector.AddHit(m_Dkey_acc,hitdata,chist,true);
        if (Layer.CheckSelection(hitdata.lay)) {
          if (m_SplitBarrel) nlayhist += Layer.AddHit(m_Lkey,hitdata,chist,true);
          else { nlayhist += Layer.AddHit(m_Lkey_acc,hitdata,chist,true); continue;}
          if (Module.CheckSelection(hitdata.mod)) {
            nmodhist += Module.AddHit(m_Mkey,hitdata,chist,true);
            if (Board.CheckSelection(hitdata.brd)) {
              nbrdhist += Board.AddHit(m_Bkey,hitdata,chist,true);
              if (Chip.CheckSelection(hitdata.chp)) {
                nchphist += Chip.AddHit(m_Ckey,hitdata,chist,true);
                if (Straw.CheckSelection(hitdata.stw))
                  nstwhist += Straw.AddHit(m_Skey,hitdata,chist,true);
              }
            }
          }
        }
      }                         // Here it closes
      ihist++;
   } else  
   if(isArgonStraw==0){                      // Separate Ar and Xe in endcaps. Here Xe
      nTRThist += TRT.AddHit(m_Tkey,hitdata,chist,true);
      if (Detector.CheckSelection(hitdata.det)) { //only add the histogram if it is in the selection
        if (m_SplitBarrel) ndethist += Detector.AddHit(m_Dkey,hitdata,chist,true);
        else ndethist += Detector.AddHit(m_Dkey_acc,hitdata,chist,true);
        if (Layer.CheckSelection(hitdata.lay)) {
          if (m_SplitBarrel) nlayhist += Layer.AddHit(m_Lkey,hitdata,chist,true);
          else { nlayhist += Layer.AddHit(m_Lkey_acc,hitdata,chist,true); continue;}
          if (Module.CheckSelection(hitdata.mod)) {
            nmodhist += Module.AddHit(m_Mkey,hitdata,chist,true);
            if (Board.CheckSelection(hitdata.brd)) {
              nbrdhist += Board.AddHit(m_Bkey,hitdata,chist,true);
              if (Chip.CheckSelection(hitdata.chp)) {
                nchphist += Chip.AddHit(m_Ckey,hitdata,chist,true);
                if (Straw.CheckSelection(hitdata.stw)) {
                  nstwhist += Straw.AddHit(m_Skey,hitdata,chist,true);
                }
              }
            }
          }
        }
      }                         // Here it closes
      ihist++;
   } else {               // ARGON HITS
      nTRThistAr += TRT_Ar.AddHit(m_Tkey,hitdata,chist,true);
      if (Detector_Ar.CheckSelection(hitdata.det)) { //only add the histogram if it is in the selection
        if (m_SplitBarrel) ndethistAr += Detector_Ar.AddHit(m_Dkey,hitdata,chist,true);
        else ndethistAr += Detector_Ar.AddHit(m_Dkey_acc,hitdata,chist,true);
        if (Layer_Ar.CheckSelection(hitdata.lay)) {
          if (m_SplitBarrel) nlayhistAr += Layer_Ar.AddHit(m_Lkey,hitdata,chist,true);
          else { nlayhistAr += Layer_Ar.AddHit(m_Lkey_acc,hitdata,chist,true); continue;}
          if (Module_Ar.CheckSelection(hitdata.mod)) {
            nmodhistAr += Module_Ar.AddHit(m_Mkey,hitdata,chist,true);
            if (Board_Ar.CheckSelection(hitdata.brd)) {
              nbrdhistAr += Board_Ar.AddHit(m_Bkey,hitdata,chist,true);
              if (Chip_Ar.CheckSelection(hitdata.chp)) {
                nchphistAr += Chip_Ar.AddHit(m_Ckey,hitdata,chist,true);
                if (Straw_Ar.CheckSelection(hitdata.stw)) {
                  nstwhistAr += Straw_Ar.AddHit(m_Skey,hitdata,chist,true);
                }
              }
            }
          }
        }
      }                         // Here it closes
    ihistAr++;
    }                   // END ARGON HITS
 
 
    if ((ihist%10000==9999) | (ihist==m_nevents-1)){
      ATH_MSG_INFO( Form("%7i HISTOGRAMS READ, UNITS ADDED: %i %i %2i %3i %3i %4i %6i",ihist+1,nTRThist, ndethist, nlayhist, nmodhist, nbrdhist, nchphist, nstwhist) );
      if(m_DoArXenonSep) ATH_MSG_INFO( Form("%7i Ar HISTOGRAMS READ, UNITS ADDED: %i %i %2i %3i %3i %4i %6i",ihistAr+1,nTRThistAr, ndethistAr, nlayhistAr, nmodhistAr, nbrdhistAr, nchphistAr, nstwhistAr) );
    }
      
      ihist++;
      delete [] chist;
    } // Finish Straw Loop
   
    ATH_MSG_INFO( Form("%7i HISTOGRAMS READ, UNITS ADDED: %i %i %2i %3i %3i %4i %6i",ihist+1,nTRThist, ndethist, nlayhist, nmodhist, nbrdhist, nchphist, nstwhist) );
    if(m_DoArXenonSep) ATH_MSG_INFO( Form("%7i Ar HISTOGRAMS READ, UNITS ADDED: %i %i %2i %3i %3i %4i %6i",ihistAr+1,nTRThistAr, ndethistAr, nlayhistAr, nmodhistAr, nbrdhistAr, nchphistAr, nstwhistAr) );
 
    TRT.UpdateOldConstants();
    Detector.UpdateOldConstants();
    Layer.UpdateOldConstants();
    Module.UpdateOldConstants();
    Board.UpdateOldConstants();
    Chip.UpdateOldConstants();
 
    if (ihistAr>0 ){
        TRT_Ar.UpdateOldConstants();
        Detector_Ar.UpdateOldConstants();
        Layer_Ar.UpdateOldConstants();
        Module_Ar.UpdateOldConstants();
        Board_Ar.UpdateOldConstants();
        Chip_Ar.UpdateOldConstants();
    }
 
    myFile.close();
 
 
  ATH_MSG_INFO(" ");
  ATH_MSG_INFO( TRT.PrintStat() );
  ATH_MSG_INFO( Detector.PrintStat() );
  ATH_MSG_INFO( Layer.PrintStat() );
  ATH_MSG_INFO( Module.PrintStat() );
  ATH_MSG_INFO( Board.PrintStat() );
  ATH_MSG_INFO( Chip.PrintStat() );
  ATH_MSG_INFO( Straw.PrintStat() );
  ATH_MSG_INFO(" ");
 
 
  if(m_DoArXenonSep){
    ATH_MSG_INFO(" ");
    ATH_MSG_INFO( TRT_Ar.PrintStat() );
    ATH_MSG_INFO( Detector_Ar.PrintStat() );
    ATH_MSG_INFO( Layer_Ar.PrintStat() );
    ATH_MSG_INFO( Module_Ar.PrintStat() );
    ATH_MSG_INFO( Board_Ar.PrintStat() );
    ATH_MSG_INFO( Chip_Ar.PrintStat() );
    ATH_MSG_INFO( Straw_Ar.PrintStat() );
    ATH_MSG_INFO(" ");
  }
 
 
 
  //caldata startdata(true,tbins,rbins);
  caldata startdata(true,55,100);
  //caldata startdata(true,64,64);
  startdata.t0=5.0;
  std::map<std::string,TDirectory*> dirmap;
  TDirectory* trtdir=gDirectory;
  TDirectory* detdir=gDirectory;
  TDirectory* laydir=gDirectory;
  TDirectory* moddir=gDirectory;
  TDirectory* brddir=gDirectory;
  TDirectory* chpdir=gDirectory;
 
  std::ofstream rtcalfile(Form("%s_rt.txt",outfile.data()),std::ios::out);
  std::ofstream binrtcalfile(Form("%s_binrt.txt",outfile.data()),std::ios::out);
  std::ofstream t0calfile(Form("%s_t0.txt",outfile.data()),std::ios::out);
 
 
  //if (m_SplitBarrel){
  
  for (std::pair<const std::string, BDdetector>& pt : m_trt.t) {
    
    if (TRT.Skip()) break;
    if (TRT.HasKey(pt.first)) {
      trtdir = TRT.Calibrate(histfile.get(),pt.first,SubLev(m_options,1),&startdata);
      if (TRT.printt0) t0calfile << Form("-3 -1 -1 -1 -1 : %e %e",TRT.data[pt.first].t0,TRT.data[pt.first].t0err) << std::endl;
      if (TRT.printrt) rtcalfile << Form("-3 -1 -1 -1 -1 : %i %e %e %e %e",rtint,TRT.data[pt.first].rtpar[0],TRT.data[pt.first].rtpar[1],TRT.data[pt.first].rtpar[2],TRT.data[pt.first].rtpar[3]) << std::endl;
    }
    for (std::pair<const std::string, BDlayer>& pd : pt.second.d) {
      
      if(Detector.Skip()) break; 
      if(Detector.HasKey(pd.first)){
        detdir = Detector.Calibrate(trtdir,pd.first,SubLev(m_options,2),&TRT.data[pt.first]);
      }
      for (std::pair<const std::string, BDmodule>& pl : pd.second.l) {
 
        if(Layer.Skip()) break;
        if(Layer.HasKey(pl.first)){
          laydir = Layer.Calibrate(detdir,pl.first,SubLev(m_options,3),&Detector.data[pd.first]);
          if (Layer.printt0) t0calfile << Form("%i %i -1 -1 -1 : %e %e",Layer.data[pl.first].det,Layer.data[pl.first].lay,Layer.data[pl.first].t0,Layer.data[pl.first].t0err) << std::endl;
          if (Layer.printrt) rtcalfile    << Form("%i %i -1 -1 -1 : %i %e %e %e %e",Layer.data[pl.first].det,Layer.data[pl.first].lay,rtint,Layer.data[pl.first].rtpar[0],Layer.data[pl.first].rtpar[1],Layer.data[pl.first].rtpar[2],Layer.data[pl.first].rtpar[3]) << std::endl;
          if (!m_SplitBarrel) {
            if (Layer.printt0) t0calfile    << Form("%i %i -1 -1 -1 : %e %e",-Layer.data[pl.first].det,Layer.data[pl.first].lay,Layer.data[pl.first].t0,Layer.data[pl.first].t0err) << std::endl;
            if (Layer.printrt) rtcalfile    << Form("%i %i -1 -1 -1 : %i %e %e %e %e",-Layer.data[pl.first].det,Layer.data[pl.first].lay,rtint,Layer.data[pl.first].rtpar[0],Layer.data[pl.first].rtpar[1],Layer.data[pl.first].rtpar[2],Layer.data[pl.first].rtpar[3]) << std::endl;
          }
        }
        for (std::pair<const std::string, BDboard>& pm : pl.second.m) {
 
          if(Module.Skip()) break; 
          if(Module.HasKey(pm.first)){ 
            moddir = Module.Calibrate(laydir,pm.first,SubLev(m_options,4),&Layer.data[pl.first]);
            if (Module.printt0) t0calfile << Form("%i %i %i -1 -1 : %e %e",Module.data[pm.first].det,Module.data[pm.first].lay,Module.data[pm.first].mod,Module.data[pm.first].t0,Module.data[pm.first].t0err) << std::endl;
            if (Module.printrt) rtcalfile    << Form("%i %i %i -1 -1 : %i %e %e %e %e",Module.data[pm.first].det,Module.data[pm.first].lay,Module.data[pm.first].mod,rtint,Module.data[pm.first].rtpar[0],Module.data[pm.first].rtpar[1],Module.data[pm.first].rtpar[2],Module.data[pm.first].rtpar[3]) << std::endl;
          } 
          for (std::pair<const std::string, BDchip>& pb : pm.second.b) {
 
            if(Board.Skip()) break; 
            if(Board.HasKey(pb.first)){ 
              brddir = Board.Calibrate(moddir,pb.first,SubLev(m_options,5),&Module.data[pm.first]);
            } 
            for (std::pair<const std::string, BDstraw>& pc : pb.second.c) {
 
              if(Chip.Skip()) break; 
              if(Chip.HasKey(pc.first)){ 
                chpdir = Chip.Calibrate(brddir,pc.first,SubLev(m_options,6),&Board.data[pb.first]);
              }
              for (std::pair<const std::string, BDzero>& ps : pc.second.s) {
 
                if(Straw.Skip()) break;
                if(Straw.HasKey(ps.first)){ 
                  Straw.Calibrate(chpdir,ps.first,SubLev(m_options,7),&Chip.data[pc.first]);
                  if (Straw.printt0) t0calfile << Form("%i %i %i %i %i : %e %e",Straw.data[ps.first].det,Straw.data[ps.first].lay,Straw.data[ps.first].mod,Straw.data[ps.first].stl,Straw.data[ps.first].stw,Straw.data[ps.first].t0,Straw.data[ps.first].t0err) << std::endl;
                  if (Straw.printrt) rtcalfile << Form("%i %i %i %i %i : %i %e %e %e %e",Straw.data[ps.first].det,Straw.data[ps.first].lay,Straw.data[ps.first].mod,Straw.data[ps.first].stl,Straw.data[ps.first].stw,rtint,Straw.data[ps.first].rtpar[0],Straw.data[ps.first].rtpar[1],Straw.data[ps.first].rtpar[2],Straw.data[ps.first].rtpar[3]) << std::endl;
                
 
                }
              }
            }
          }
        }
      }
    }
  }
 
  std::map<std::string,TDirectory*> dirArmap;
  TDirectory* trtdirAr=gDirectory;
  TDirectory* detdirAr=gDirectory;
  TDirectory* laydirAr=gDirectory;
  TDirectory* moddirAr=gDirectory;
  TDirectory* brddirAr=gDirectory;
  TDirectory* chpdirAr=gDirectory;
 
  if (m_DoArXenonSep){
 
  for (std::pair<const std::string, BDdetector>& pt : m_trt.t) {
 
    if (TRT_Ar.Skip()) break;
    if (TRT_Ar.HasKey(pt.first)) {
      trtdirAr = TRT_Ar.Calibrate(histfile.get(),pt.first,SubLev(m_options,1),&startdata);
    }
    for (std::pair<const std::string, BDlayer>& pd : pt.second.d) {
 
      if(Detector_Ar.Skip()) break;
      if(Detector_Ar.HasKey(pd.first)){
        detdirAr = Detector_Ar.Calibrate(trtdirAr,pd.first,SubLev(m_options,2),&TRT_Ar.data[pt.first]);
      }
      for (std::pair<const std::string, BDmodule>& pl : pd.second.l) {
 
        if(Layer_Ar.Skip()) break;
        if(Layer_Ar.HasKey(pl.first)){
          laydirAr = Layer_Ar.Calibrate(detdirAr,pl.first,SubLev(m_options,3),&Detector_Ar.data[pd.first]);
       
        }
        for (std::pair<const std::string, BDboard>& pm : pl.second.m) {
 
          if(Module_Ar.Skip()) break;
          if(Module_Ar.HasKey(pm.first)){
            moddirAr = Module_Ar.Calibrate(laydirAr,pm.first,SubLev(m_options,4),&Layer_Ar.data[pl.first]);
            if (Module_Ar.printt0) t0calfile << Form("%i %i %i -1 -1 : %e %e",Module_Ar.data[pm.first].det,Module_Ar.data[pm.first].lay,Module_Ar.data[pm.first].mod,Module_Ar.data[pm.first].t0,Module_Ar.data[pm.first].t0err) << std::endl;
            if (Layer_Ar.printrt) rtcalfile    << Form("%i %i %i -1 -1 : %i %e %e %e %e",Module_Ar.data[pm.first].det,Module_Ar.data[pm.first].lay,Module_Ar.data[pm.first].mod,rtint,Module_Ar.data[pm.first].rtpar[0],Module_Ar.data[pm.first].rtpar[1],Module_Ar.data[pm.first].rtpar[2],Module_Ar.data[pm.first].rtpar[3]) << std::endl;
          }
          for (std::pair<const std::string, BDchip>& pb : pm.second.b) {
            if(Board_Ar.Skip()) break;
            if(Board_Ar.HasKey(pb.first)){
              brddirAr = Board_Ar.Calibrate(moddirAr,pb.first,SubLev(m_options,5),&Module_Ar.data[pm.first]);
            }
            for (std::pair<const std::string, BDstraw>& pc : pb.second.c) {
 
              if(Chip_Ar.Skip()) break;
              if(Chip_Ar.HasKey(pc.first)){
                chpdirAr = Chip_Ar.Calibrate(brddirAr,pc.first,SubLev(m_options,6),&Board_Ar.data[pb.first]);
              }
              for (std::pair<const std::string, BDzero>& ps : pc.second.s) {
 
                if(Straw_Ar.Skip()) break;
                if(Straw_Ar.HasKey(ps.first)){
                  Straw_Ar.Calibrate(chpdirAr,ps.first,SubLev(m_options,7),&Chip_Ar.data[pc.first]);
                  if (Straw_Ar.printt0) t0calfile << Form("%i %i %i %i %i : %e %e",Straw_Ar.data[ps.first].det,Straw_Ar.data[ps.first].lay,Straw_Ar.data[ps.first].mod,Straw_Ar.data[ps.first].stl,Straw_Ar.data[ps.first].stw,Straw_Ar.data[ps.first].t0,Straw_Ar.data[ps.first].t0err) << std::endl;
                  if (Straw_Ar.printrt) rtcalfile << Form("%i %i %i %i %i : %i %e %e %e %e",Straw_Ar.data[ps.first].det,Straw_Ar.data[ps.first].lay,Straw_Ar.data[ps.first].mod,Straw_Ar.data[ps.first].stl,Straw_Ar.data[ps.first].stw,rtint,Straw_Ar.data[ps.first].rtpar[0],Straw_Ar.data[ps.first].rtpar[1],Straw_Ar.data[ps.first].rtpar[2],Straw_Ar.data[ps.first].rtpar[3]) << std::endl;
                }
              }
            }
          }
        }
      }
    }
  }
 
}// Close for ARGON 
 
 
 
  
  //make ntuples for each calibrator 
  if (!TRT.bequiet) TRT.WriteStat(histfile.get());
  if (!Detector.bequiet) Detector.WriteStat(histfile.get());
  if (!Layer.bequiet) Layer.WriteStat(histfile.get());
  if (!Module.bequiet) Module.WriteStat(histfile.get());
  if (!Board.bequiet) Board.WriteStat(histfile.get());
  Chip.WriteStat(histfile.get());
  Straw.WriteStat(histfile.get());
 
  TRT.DumpConstants();
  Detector.DumpConstants();
  Layer.DumpConstants();
  Module.DumpConstants();
  Board.DumpConstants();
  Chip.DumpConstants();
  Straw.DumpConstants();
  
  rtcalfile.close();
  binrtcalfile.close();
  t0calfile.close();
 
        //AR
  //make ntuples for each calibrator 
  if(m_DoArXenonSep){
    if (!TRT_Ar.bequiet)          TRT_Ar.WriteStat(histfile.get());
    if (!Detector_Ar.bequiet)     Detector_Ar.WriteStat(histfile.get());
    if (!Layer_Ar.bequiet)        Layer_Ar.WriteStat(histfile.get());
    if (!Module_Ar.bequiet)       Module_Ar.WriteStat(histfile.get());
    if (!Board_Ar.bequiet)        Board_Ar.WriteStat(histfile.get());
    Chip_Ar.WriteStat(histfile.get());
    Straw_Ar.WriteStat(histfile.get());
 
   TRT_Ar.DumpConstants();
   Detector_Ar.DumpConstants();
   Layer_Ar.DumpConstants();
   Module_Ar.DumpConstants();
   Board_Ar.DumpConstants();
   Chip_Ar.DumpConstants();
   Straw_Ar.DumpConstants();
  }
 
  histfile->ls();  
  ATH_MSG_INFO( "writing out calibout.root");  
  histfile->Write();
 
  return true;
  
}