JInterDomCal_utils.hh

Go to the documentation of this file.

#ifndef __CALIB_UTIL__
#define __CALIB_UTIL__

#include "Detector.hh"
#include "NBRun.hh"
#include "Physics.hh"

/**
 * \author rgruiz
 */



/**
 * It uses the pmts from a source module and from a target module. 
 * For each possible pair source-target, it computes the a factor that needs to be added to the t0s of all the PMTs in the target module. Finally it computes the median  
 *
 * \param det detector with old t0s
 * \param srcs SUPERPMTs of the emitter DOM
 * \param tgts SUPERPMTs of the receiver DOM
 * \return Median correction factor
 */
inline double computeMedianCorrection(JDetector* det , vector < SuperPMT* > srcs , vector < SuperPMT* > tgts){

  double c_w = c_water_genhen() ;
  
  //double c_w = c_water(470 , 3000);

  vector<double> corrections ;
  
  const JPMTRouter pmtRouter(*det);

  for (auto & src : srcs){
    
    JPMT src_pmt = src->getPMT() ;

    NBPulse* src_pulse = src->getNBPulse() ;

    for(auto & tgt : tgts){
    
      JPMT tgt_pmt = tgt->getPMT() ;

      NBPulse* tgt_pulse = tgt->getNBPulse() ;

      double distance = get_PMT_distance(src_pmt , tgt_pmt) ;

      double t_flight = distance / c_w ;

      double delta_traw = tgt_pulse->getPeakTime() - src_pulse->getPeakTime() ;

      double delta_t0 = t_flight - delta_traw ;

      double delta_t0_old = pmtRouter.getPMT(JObjectID(tgt_pmt.getID())).getT0() - pmtRouter.getPMT(JObjectID(src_pmt.getID())).getT0() ; 

      corrections.push_back(delta_t0 - delta_t0_old) ;

    }

  }
  
  sort (corrections.begin() , corrections.end()) ;

  int index = int(0.5*corrections.size()) ;

  return corrections[index] ;

}
 
/**
 * Produces a new detector file by updating the existing t0s of the PMTs 
 *
 * \param detector Old detector
 * \param SuperMods The SuperModules used to calibrate the detector
 * \return calibrated detector
 */
inline JDetector Update_Detector(JDetector detector , const vector< SuperModule* > & SuperMods){
  
  JDetector new_det = detector ;

  for(auto & sm : SuperMods){

    JModuleLocation tgt_loc = sm->getDom().getLocation() ;

    if(sm->get_good_sources().size()>0){
      
      pair < SuperModule* , vector < SuperPMT* > > src = sm->get_closest_good_source();

      vector < SuperPMT* > srcs = src.first->get_ref_pmts() ;

      vector < SuperPMT* > tgts = src.second ;

      double correction = computeMedianCorrection(&new_det , srcs , tgts) ;
       
      int nPMTs = getNumberOfPMTs(detector)/getNumberOfModules(detector) ;

      for(int j=0 ; j < nPMTs ; j++){

        new_det.getModule(tgt_loc).getPMT(j).setT0(new_det.getModule(tgt_loc).getPMT(j).getT0() + correction ) ;
        
      } 
      
    }
    
  }
    
  return new_det; 

}


#endif