Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:26:28

0001 #ifndef PFProducer_PFCandConnector_H_
0002 #define PFProducer_PFCandConnector_H_
0003 
0004 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0005 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0009 
0010 // \author : M. Gouzevitch
0011 // \date : May 2010
0012 
0013 /// Based on a class from : V. Roberfroid, February 2008
0014 
0015 class PFCandConnector {
0016 public:
0017   PFCandConnector() {
0018     bCorrect_ = false;
0019     bCalibPrimary_ = false;
0020 
0021     fConst_.push_back(1), fConst_.push_back(0);
0022     fNorm_.push_back(0), fNorm_.push_back(0);
0023     fExp_.push_back(0);
0024 
0025     dptRel_PrimaryTrack_ = 0.;
0026     dptRel_MergedTrack_ = 0.;
0027     ptErrorSecondary_ = 0.;
0028   }
0029 
0030   void setParameters(const edm::ParameterSet& iCfgCandConnector) {
0031     bool bCorrect, bCalibPrimary;
0032     double dptRel_PrimaryTrack, dptRel_MergedTrack, ptErrorSecondary;
0033     std::vector<double> nuclCalibFactors;
0034 
0035     /// Flag to apply the correction procedure for nuclear interactions
0036     bCorrect = iCfgCandConnector.getParameter<bool>("bCorrect");
0037     /// Flag to calibrate the reconstructed nuclear interactions with primary or merged tracks
0038     bCalibPrimary = iCfgCandConnector.getParameter<bool>("bCalibPrimary");
0039 
0040     if (iCfgCandConnector.exists("dptRel_PrimaryTrack"))
0041       dptRel_PrimaryTrack = iCfgCandConnector.getParameter<double>("dptRel_PrimaryTrack");
0042     else {
0043       edm::LogWarning("PFCandConnector") << "dptRel_PrimaryTrack doesn't exist. Setting a default safe value 0"
0044                                          << std::endl;
0045       dptRel_PrimaryTrack = 0;
0046     }
0047 
0048     if (iCfgCandConnector.exists("dptRel_MergedTrack"))
0049       dptRel_MergedTrack = iCfgCandConnector.getParameter<double>("dptRel_MergedTrack");
0050     else {
0051       edm::LogWarning("PFCandConnector") << "dptRel_MergedTrack doesn't exist. Setting a default safe value 0"
0052                                          << std::endl;
0053       dptRel_MergedTrack = 0;
0054     }
0055 
0056     if (iCfgCandConnector.exists("ptErrorSecondary"))
0057       ptErrorSecondary = iCfgCandConnector.getParameter<double>("ptErrorSecondary");
0058     else {
0059       edm::LogWarning("PFCandConnector") << "ptErrorSecondary doesn't exist. Setting a default safe value 0"
0060                                          << std::endl;
0061       ptErrorSecondary = 0;
0062     }
0063 
0064     if (iCfgCandConnector.exists("nuclCalibFactors"))
0065       nuclCalibFactors = iCfgCandConnector.getParameter<std::vector<double> >("nuclCalibFactors");
0066     else {
0067       edm::LogWarning("PFCandConnector") << "nuclear calib factors doesn't exist the factor would not be applyed"
0068                                          << std::endl;
0069     }
0070 
0071     setParameters(bCorrect, bCalibPrimary, dptRel_PrimaryTrack, dptRel_MergedTrack, ptErrorSecondary, nuclCalibFactors);
0072   }
0073 
0074   void setParameters(bool bCorrect,
0075                      bool bCalibPrimary,
0076                      double dptRel_PrimaryTrack,
0077                      double dptRel_MergedTrack,
0078                      double ptErrorSecondary,
0079                      const std::vector<double>& nuclCalibFactors);
0080 
0081   static void fillPSetDescription(edm::ParameterSetDescription& iDesc);
0082 
0083   reco::PFCandidateCollection connect(reco::PFCandidateCollection& pfCand) const;
0084 
0085 private:
0086   /// Analyse nuclear interactions where a primary or merged track is present
0087   void analyseNuclearWPrim(reco::PFCandidateCollection&, std::vector<bool>&, unsigned int) const;
0088 
0089   /// Analyse nuclear interactions where a secondary track is present
0090   void analyseNuclearWSec(reco::PFCandidateCollection&, std::vector<bool>&, unsigned int) const;
0091 
0092   bool isPrimaryNucl(const reco::PFCandidate& pf) const;
0093 
0094   bool isSecondaryNucl(const reco::PFCandidate& pf) const;
0095 
0096   /// Return a calibration factor for a reconstructed nuclear interaction
0097   double rescaleFactor(const double pt, const double cFrac) const;
0098 
0099   /// Parameters
0100   bool bCorrect_;
0101 
0102   /// Calibration parameters for the reconstructed nuclear interactions
0103   bool bCalibPrimary_;
0104   std::vector<double> fConst_;
0105   std::vector<double> fNorm_;
0106   std::vector<double> fExp_;
0107 
0108   // Maximal accepatble uncertainty on primary tracks to usem them as MC truth for calibration
0109   double dptRel_PrimaryTrack_;
0110   double dptRel_MergedTrack_;
0111   double ptErrorSecondary_;
0112 
0113   /// Useful constants
0114   static const double pion_mass2;
0115   static const reco::PFCandidate::Flags fT_TO_DISP_;
0116   static const reco::PFCandidate::Flags fT_FROM_DISP_;
0117 };
0118 
0119 #endif