// MuonReconstructor.cxx
 //
 // Implementation of the Muon Smearer class
 //
 // Only the reconstruct() method is overridden
 //
 //
 // Authors: H.T. Phillips, P.Clarke, E.Richter-Was, P.Sherwood, R.Steward
 //
 // Smearer - Reconstructor Conversion: Alexander Richards (UCL)
 //
 
 #include "AtlfastAlgs/MuonReconstructor.h"
 #include <cmath>
 
 #include "CLHEP/Vector/LorentzVector.h"
 #include "CLHEP/Units/SystemOfUnits.h"  
 #include "CLHEP/Random/JamesRandom.h"
 #include "CLHEP/Random/RandGauss.h"
 #include <iostream>
 
 #include "PathResolver/PathResolver.h"
 
 //-------------------------------------------------------------
   /** @brief Provides parameterised muon smearing.
    *
    * MuonReconstructor honours the IReconstructor interface and
    * inherets implementation from DefaultReconstructor.
    */
 
 namespace Atlfast {
   using std::abs;
   using std::pair;
   
 
   // MuonReconstructor methods
   
   /** Constructor called from Athena code - instantiates TrackSmearer */
   MuonReconstructor::MuonReconstructor(const int aseed, const int lumi, const int keymuo,
                            std::string &resolutionFile, MsgStream& log): 
     IDefaultReconstructor(), DefaultReconstructor(aseed), m_lumi(lumi), m_keymuo(keymuo)
   {    
     makeResolutionCalculator(resolutionFile);
    
     // Make a (default) TrackSmearer
     int randSeed = 12345;
     m_tracksmearer = new TrackSmearer(randSeed,log);
   }
   
   /** Constructor that can be called from outside Athena - no TrackSmearer */
   MuonReconstructor::MuonReconstructor(const int aseed, const int lumi, const int keymuo,
                            std::string &resolutionFile): 
     IDefaultReconstructor(), DefaultReconstructor(aseed), m_lumi(lumi), m_keymuo(keymuo)
   {    
     makeResolutionCalculator(resolutionFile);
   }
 
   MuonReconstructor::~MuonReconstructor()
   {
     delete m_muonResCalculator;
     delete m_tracksmearer;
   }
   
   void MuonReconstructor::makeResolutionCalculator(std::string &resolutionFile){
     
     std::string filename =
       PathResolver::find_file(resolutionFile, "DATAPATH");
     std::cout << "MuonReconstructor getting file " << filename << std::endl;
     
     if (filename!="")
       m_muonResCalculator = new MuonResolutionCalculator(filename.c_str());
     else
       std::cout << "Could not open muon resolution file!!" << std::endl;
 
     return;
     
   }
 
   /** reconstruct method for ReconstructedParticle */ 
   ReconstructedParticle MuonReconstructor::reconstruct( const ReconstructedParticle &particle )const{
 
     if(!particle.truth() && m_keymuo>1){
         std::cout << "No GenParticle in the ReconstructedParticle." << m_keymuo 
                   << "ID smearing requested but no GenParticle - unsmeared particle returned." << std::endl;
         return particle;
     }
 
     HepLorentzVector avec(particle.px(),particle.py(),particle.pz(),particle.e());
     
     double sigmaMS;
     double sigmamuon;
     double sigmaID;
     double sigmatrack;
     double sigma=0.;
     HepLorentzVector bvec(0,0,0,0);
     pair<double, double> sigmapr;
 
     if(m_smearParamSchema==1){
       switch (m_keymuo) {
       case 1:
         while( ((sigmapr=resolMS(avec)).first) <-1.) {}
         sigma = sigmapr.first;
         break;
       case 2:
         while( ((sigmapr=resolID(*particle.truth())).first)<-1.) {}
         sigma = sigmapr.first;
         break;
       case 3:
         while( ((sigmapr=resolMS(avec)).first) <-1.) {}
         sigmaMS    = sigmapr.first;
         sigmamuon  = sigmapr.second;
         while( ((sigmapr=resolID(*particle.truth())).first)<-1.) {}
         sigmaID    = sigmapr.first;
         sigmatrack = sigmapr.second;
         
         sigma=resol(sigmaMS, sigmamuon, sigmaID, sigmatrack);
         break;
       default:
         std::cout << "Invalid value for MuonSmearKey: " << m_keymuo 
                   << ", returning unsmearer muon particle" << std::endl;
         return particle;
         break;      
       }
       
       bvec = avec/(1.0+sigma);
       
       //std::cout << "Unsmeared muon: " << avec << ", smeared muon: " << bvec << std::endl;
       
     } ////end of if (m_smearParamSchema==1)
     
     ReconstructedParticle return_particle(particle);
     return_particle.set_momentum(bvec);
     return return_particle; // bvec is zeroes if m_smearParamSchema!=1
   }
   
   pair<double, double> MuonReconstructor::resolMS(const HepLorentzVector& avec)const{
     
     if (!m_muonResCalculator)
       std::cout << "No MuonResolutionCalculator, resolution set to 0.0" << std::endl;
     
     double sigmamuon = m_muonResCalculator ?
       m_muonResCalculator->calculateResolution(avec) : 0.;
     if (!sigmamuon){
       std::cout << "Zero spectrometer resolution!!!" << std::endl;
       pair<double, double>  sigmapr(0,0);
       return sigmapr;
     }
     
     double aa=m_randGauss->fire();
     double sigmams   = aa*sigmamuon;
     pair<double, double> sigmapr(sigmams, sigmamuon);
     
     return sigmapr;
   }
   
   pair<double, double> MuonReconstructor::resolID(const HepMC::GenParticle& particle)const{
 
     if (!m_tracksmearer){
       std::cout << "No TrackSmearer in MuonReconstructor, resolution set to 0.0" << std::endl;
       std::cout << "To get inner detector resolutions, you must instantiate"
                 << " MuonReconstructor with a MsgStream" << std::endl;
     }
 
     double sigmatrack = m_tracksmearer ? 
       m_tracksmearer->getPtResolution(particle) : 0.;
 
     if (!sigmatrack){
       std::cout << "Zero inner detector resolution!!!" << std::endl;
       pair<double, double>  sigmapr(0,0);
       return sigmapr;
     }
 
     double aa=m_randGauss->fire();
     double sigmatr = aa*sigmatrack;
     pair<double, double> sigmapr(sigmatr, sigmatrack);
       
     return sigmapr;
 
   }
 
   double MuonReconstructor::resol ( double sigmams, double sigmamuon, 
                                     double sigmaid, double sigmatrack)const{
     double wmuon  = 1./(sigmamuon*sigmamuon);
     double wtrack = 1./(sigmatrack*sigmatrack);
     double wtot   = wmuon+wtrack;
     double corr   = (wmuon*(1.0+sigmams)+wtrack*(1.0+sigmaid))/wtot;
     double sigma   = corr-1.0;
     
     return sigma;
   }
   
   //cct: implement the setSmearParameters method
   int MuonReconstructor::setSmearParameters (const std::vector<double>& smearValues){
     m_smearParams=smearValues;
     return 0;
   }
   //cct: implement the setSmearParamSchema method
   int MuonReconstructor::setSmearParamSchema ( const int smearSchema){
     m_smearParamSchema=smearSchema;
     return 0;
   }
 
   
   
 } // end of namespace bracket
 

Due to the LXR bug, the updates fail sometimes to remove references to deleted files. The Saturday's full rebuilds fix these problems
This page was automatically generated by the LXR engine.