#include <MuonRadiationCorrWeightProducer.h>

Inheritance diagram for MuonRadiationCorrWeightProducer:

Classes
struct	lutEntryType

Public Member Functions
	MuonRadiationCorrWeightProducer (const edm::ParameterSet &)

void	produce (edm::Event &, const edm::EventSetup &)

	~MuonRadiationCorrWeightProducer ()

Public Member Functions inherited from edm::EDProducer
	EDProducer ()

ModuleDescription const &	moduleDescription () const

virtual	~EDProducer ()

Public Member Functions inherited from edm::ProducerBase
void	callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)

	ProducerBase ()

void	registerProducts (ProducerBase , ProductRegistry , ModuleDescription const &)

std::function< void(BranchDescription const &)>	registrationCallback () const
	used by the fwk to register list of products More...

virtual	~ProducerBase ()

Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfo >	consumesInfo () const

	EDConsumerBase ()

ProductHolderIndexAndSkipBit	indexFrom (EDGetToken, BranchType, TypeID const &) const

void	itemsMayGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const

void	itemsToGet (BranchType, std::vector< ProductHolderIndexAndSkipBit > &) const

std::vector < ProductHolderIndexAndSkipBit > const &	itemsToGetFromEvent () const

void	labelsForToken (EDGetToken iToken, Labels &oLabels) const

void	modulesDependentUpon (const std::string &iProcessName, std::vector< const char * > &oModuleLabels) const

void	modulesWhoseProductsAreConsumed (std::vector< ModuleDescription const * > &modules, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const

bool	registeredToConsume (ProductHolderIndex, bool, BranchType) const

bool	registeredToConsumeMany (TypeID const &, BranchType) const

void	updateLookup (BranchType iBranchType, ProductHolderIndexHelper const &)

virtual	~EDConsumerBase ()

Private Types
typedef std::vector < lutEntryType * >	vlutEntryType

Private Attributes
std::map< std::string, std::string >	lutDirectoriesOthers_

std::string	lutDirectoryRef_

std::map< std::string, vlutEntryType >	lutEntriesOthers_

vlutEntryType	lutEntriesRef_

int	maxWarnings_

double	maxWeight_

double	minWeight_

int	numOthers_

int	numWarnings_

edm::InputTag	srcMuonsAfterRad_

edm::InputTag	srcMuonsBeforeRad_

int	verbosity_

Additional Inherited Members
Public Types inherited from edm::EDProducer
typedef EDProducer	ModuleType

Public Types inherited from edm::ProducerBase
typedef ProductRegistryHelper::TypeLabelList	TypeLabelList

Public Types inherited from edm::EDConsumerBase
typedef ProductLabels	Labels

Static Public Member Functions inherited from edm::EDProducer
static const std::string &	baseType ()

static void	fillDescriptions (ConfigurationDescriptions &descriptions)

static void	prevalidate (ConfigurationDescriptions &descriptions)

Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType >	consumes (edm::InputTag const &tag)

EDGetToken	consumes (const TypeToGet &id, edm::InputTag const &tag)

template<BranchType B>
EDGetToken	consumes (TypeToGet const &id, edm::InputTag const &tag)

ConsumesCollector	consumesCollector ()
	Use a ConsumesCollector to gather consumes information from helper functions. More...

template<typename ProductType , BranchType B = InEvent>
void	consumesMany ()

void	consumesMany (const TypeToGet &id)

template<BranchType B>
void	consumesMany (const TypeToGet &id)

template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType >	mayConsume (edm::InputTag const &tag)

EDGetToken	mayConsume (const TypeToGet &id, edm::InputTag const &tag)

template<BranchType B>
EDGetToken	mayConsume (const TypeToGet &id, edm::InputTag const &tag)

Detailed Description

Compute reweighting factor and uncertainty related to correction applied to energy and momentum of embedded taus in order to compensate for muon -> muon + photon radiation of reconstructed muons

Authors: Christian Veelken

Version

Revision:: 1.1

Id:: MuonRadiationCorrWeightProducer.h,v 1.1 2013/01/31 16:15:37 veelken Exp

Definition at line 35 of file MuonRadiationCorrWeightProducer.h.

Member Typedef Documentation

typedef std::vector<lutEntryType*> MuonRadiationCorrWeightProducer::vlutEntryType

private

Definition at line 118 of file MuonRadiationCorrWeightProducer.h.

Constructor & Destructor Documentation

MuonRadiationCorrWeightProducer::MuonRadiationCorrWeightProducer ( const edm::ParameterSet & cfg )

explicit

Definition at line 15 of file MuonRadiationCorrWeightProducer.cc.

References Exception, edm::ParameterSet::exists(), edm::ParameterSet::getParameter(), edm::ParameterSet::getParameterNamesForType(), analyzePatCleaning_cfg::inputFile, writeAntiElectronDiscrMVA_cfi::inputFileName, lutDirectoriesOthers_, lutDirectoryRef_, lutEntriesOthers_, lutEntriesRef_, maxWarnings_, maxWeight_, minWeight_, numOthers_, srcMuonsAfterRad_, srcMuonsBeforeRad_, AlCaHLTBitMon_QueryRunRegistry::string, edm::FileInPath::Unknown, and verbosity_.

   : numOthers_(0),
     numWarnings_(0)
 {
   srcMuonsBeforeRad_ = cfg.getParameter<edm::InputTag>("srcMuonsBeforeRad"); 
   srcMuonsAfterRad_ = cfg.getParameter<edm::InputTag>("srcMuonsAfterRad"); 
 
   edm::FileInPath inputFileName = cfg.getParameter<edm::FileInPath>("inputFileName");
   if ( inputFileName.location() == edm::FileInPath::Unknown)
     throw cms::Exception("MuonRadiationCorrWeightProducer") 
       << " Failed to find File = " << inputFileName << " !!\n";
   std::auto_ptr<TFile> inputFile(new TFile(inputFileName.fullPath().data()));
 
   typedef std::vector<double> vdouble;
   vdouble binningMuonEn = cfg.getParameter<vdouble>("binningMuonEn");
   int numBinsMuonEn = binningMuonEn.size() - 1;
   if ( !(numBinsMuonEn >= 1) ) throw cms::Exception("Configuration")
     << " Invalid Configuration Parameter 'binningMuonEn', must define at least one bin !!\n";
 
   lutDirectoryRef_ = cfg.getParameter<std::string>("lutDirectoryRef");
   for ( int iBinMuPlusEn = 0; iBinMuPlusEn < numBinsMuonEn; ++iBinMuPlusEn ) {
     double minMuPlusEn = binningMuonEn[iBinMuPlusEn];
     double maxMuPlusEn = binningMuonEn[iBinMuPlusEn + 1];
     for ( int iBinMuMinusEn = 0; iBinMuMinusEn < numBinsMuonEn; ++iBinMuMinusEn ) {
       double minMuMinusEn = binningMuonEn[iBinMuMinusEn];
       double maxMuMinusEn = binningMuonEn[iBinMuMinusEn + 1];
       lutEntryType* lutEntry = 
         new lutEntryType(*inputFile, lutDirectoryRef_, minMuPlusEn, maxMuPlusEn, minMuMinusEn, maxMuMinusEn);
       lutEntriesRef_.push_back(lutEntry);
     }
   }
   
   edm::ParameterSet cfgNameOthers = cfg.getParameter<edm::ParameterSet>("lutDirectoryOthers");
   std::vector<std::string> names_others = cfgNameOthers.getParameterNamesForType<std::string>();
   for ( std::vector<std::string>::const_iterator name_others = names_others.begin(); 
         name_others != names_others.end(); ++name_others ) {
     std::string lutDirectoryOther = cfgNameOthers.getParameter<std::string>(*name_others);
     lutDirectoriesOthers_[*name_others] = lutDirectoryOther;
     for ( int iBinMuPlusEn = 0; iBinMuPlusEn < numBinsMuonEn; ++iBinMuPlusEn ) {
       double minMuPlusEn = binningMuonEn[iBinMuPlusEn];
       double maxMuPlusEn = binningMuonEn[iBinMuPlusEn + 1];
       for ( int iBinMuMinusEn = 0; iBinMuMinusEn < numBinsMuonEn; ++iBinMuMinusEn ) {
     double minMuMinusEn = binningMuonEn[iBinMuMinusEn];
     double maxMuMinusEn = binningMuonEn[iBinMuMinusEn + 1];
     lutEntryType* lutEntry = 
           new lutEntryType(*inputFile, lutDirectoryOther, minMuPlusEn, maxMuPlusEn, minMuMinusEn, maxMuMinusEn);
         lutEntriesOthers_[*name_others].push_back(lutEntry);
       }
     }    
     ++numOthers_;
   }
   if ( !(numOthers_ >= 1) )
     throw cms::Exception("MuonRadiationCorrWeightProducer") 
       << " Invalid Configuration Parameter 'lutNameOthers': No alternative models defined !!\n";
   
   minWeight_ = cfg.getParameter<double>("minWeight"); 
   maxWeight_ = cfg.getParameter<double>("maxWeight");
 
   verbosity_ = ( cfg.exists("verbosity") ) ?
     cfg.getParameter<int>("verbosity") : 0;
 
   maxWarnings_ = 10;
 
   produces<double>("weight");
   produces<double>("weightUp");
   produces<double>("weightDown");
 }

MuonRadiationCorrWeightProducer::~MuonRadiationCorrWeightProducer ( )

Definition at line 83 of file MuonRadiationCorrWeightProducer.cc.

References lutEntriesOthers_, and lutEntriesRef_.

 {
   for ( vlutEntryType::iterator it = lutEntriesRef_.begin();
     it != lutEntriesRef_.end(); ++it ) {
     delete (*it);
   }
   
   for ( std::map<std::string, vlutEntryType>::iterator it1 = lutEntriesOthers_.begin();
     it1 != lutEntriesOthers_.end(); ++it1 ) {
     for ( vlutEntryType::iterator it2 = it1->second.begin();
       it2 != it1->second.end(); ++it2 ) {
       delete (*it2);
     }
   }
 }

Member Function Documentation

void MuonRadiationCorrWeightProducer::produce	(	edm::Event &	evt,
		const edm::EventSetup &	es
	)

virtual

Implements edm::EDProducer.

Definition at line 107 of file MuonRadiationCorrWeightProducer.cc.

References gather_cfg::cout, relativeConstraints::error, findMuons(), lutEntriesOthers_, lutEntriesRef_, bookConverter::max, maxWarnings_, maxWeight_, minWeight_, numOthers_, numWarnings_, edm::Event::put(), mathSSE::sqrt(), Clusterizer1DCommons::square(), srcMuonsAfterRad_, srcMuonsBeforeRad_, verbosity_, and puppiForMET_cff::weight.

 {
   if ( verbosity_ ) {
     std::cout << "<MuonRadiationCorrWeightProducer::produce>:" << std::endl;
     std::cout << " srcMuonsBeforeRad = " << srcMuonsBeforeRad_ << std::endl;
     std::cout << " srcMuonsAfterRad = " << srcMuonsAfterRad_ << std::endl;
   }
 
   reco::Candidate::LorentzVector genMuonPlusP4_beforeRad;
   bool genMuonPlus_beforeRad_found = false;
   reco::Candidate::LorentzVector genMuonMinusP4_beforeRad;
   bool genMuonMinus_beforeRad_found = false;
   findMuons(evt, srcMuonsBeforeRad_, genMuonPlusP4_beforeRad, genMuonPlus_beforeRad_found, genMuonMinusP4_beforeRad, genMuonMinus_beforeRad_found);
 
   reco::Candidate::LorentzVector genMuonPlusP4_afterRad;
   bool genMuonPlus_afterRad_found = false;
   reco::Candidate::LorentzVector genMuonMinusP4_afterRad;
   bool genMuonMinus_afterRad_found = false;
   findMuons(evt, srcMuonsAfterRad_, genMuonPlusP4_afterRad, genMuonPlus_afterRad_found, genMuonMinusP4_afterRad, genMuonMinus_afterRad_found);
 
   bool genMuonPlus_found = (genMuonPlus_beforeRad_found && genMuonPlus_afterRad_found);
   bool genMuonMinus_found = (genMuonMinus_beforeRad_found && genMuonMinus_afterRad_found);
   
   double weight = 1.;
   double weightErr = 1.;
   if ( genMuonPlus_found && genMuonMinus_found ) {
     if ( verbosity_ ) {
       std::cout << " muon+: En = " << genMuonPlusP4_afterRad.E() << ", dEn = " << (genMuonPlusP4_beforeRad.E() - genMuonPlusP4_afterRad.E()) << std::endl;
       std::cout << " muon-: En = " << genMuonMinusP4_afterRad.E() << ", dEn = " << (genMuonMinusP4_beforeRad.E() - genMuonMinusP4_afterRad.E()) << std::endl;
     }
     
     int error = 0;
 
     double pRef = 0;
     bool pRef_found = false;
     for ( vlutEntryType::iterator lutEntry = lutEntriesRef_.begin();
       lutEntry != lutEntriesRef_.end(); ++lutEntry ) {
       if ( (*lutEntry)->isWithinBounds(genMuonPlusP4_afterRad.E(), genMuonMinusP4_afterRad.E()) ) {
     pRef = (*lutEntry)->getP(genMuonPlusP4_beforeRad, genMuonPlusP4_afterRad, genMuonMinusP4_beforeRad, genMuonMinusP4_afterRad);
     if ( verbosity_ ) std::cout << "pRef = " << pRef << std::endl;
     pRef_found = true;
       }
     }
     if ( !pRef_found ) {
       if ( numWarnings_ < maxWarnings_ ) {
     edm::LogWarning ("<MuonRadiationCorrWeightProducer>")
       << "Failed to find entry in reference LUT for" 
       << " muon+: En = " << genMuonPlusP4_afterRad.E() << ", dEn = " << (genMuonPlusP4_beforeRad.E() - genMuonPlusP4_afterRad.E()) << ";"
       << " muon-: En = " << genMuonMinusP4_afterRad.E() << ", dEn = " << (genMuonMinusP4_beforeRad.E() - genMuonMinusP4_afterRad.E()) << " !!" << std::endl;
     ++numWarnings_;
       }
       error = 1;
     }
     
     std::map<std::string, double> pOthers; // key = model name
     for ( std::map<std::string, vlutEntryType>::iterator lutEntriesOther = lutEntriesOthers_.begin();
       lutEntriesOther != lutEntriesOthers_.end(); ++lutEntriesOther ) {
       bool pOther_found = false;
       for ( vlutEntryType::iterator lutEntry = lutEntriesOther->second.begin();
         lutEntry != lutEntriesOther->second.end(); ++lutEntry ) {
     if ( (*lutEntry)->isWithinBounds(genMuonPlusP4_afterRad.E(), genMuonMinusP4_afterRad.E()) ) {
       pOthers[lutEntriesOther->first] = (*lutEntry)->getP(genMuonPlusP4_beforeRad, genMuonPlusP4_afterRad, genMuonMinusP4_beforeRad, genMuonMinusP4_afterRad);
       if ( verbosity_ ) std::cout << "pOthers[" << lutEntriesOther->first << "] = " << pOthers[lutEntriesOther->first] << std::endl;
       pOther_found = true;
     }
       }
       if ( !pOther_found ) {
     if ( numWarnings_ < maxWarnings_ ) {
       edm::LogWarning ("<MuonRadiationCorrWeightProducer>")
         << "Failed to find entry in LUT = " << lutEntriesOther->first << " for"
         << " muon+: En = " << genMuonPlusP4_afterRad.E() << ", dEn = " << (genMuonPlusP4_beforeRad.E() - genMuonPlusP4_afterRad.E()) << ";"
         << " muon-: En = " << genMuonMinusP4_afterRad.E() << ", dEn = " << (genMuonMinusP4_beforeRad.E() - genMuonMinusP4_afterRad.E()) << " !!" << std::endl;
       ++numWarnings_;
     }
     error = 1;
       }
     }
 
     if ( !error ) {
       double pMean = pRef;
       for ( std::map<std::string, double>::const_iterator pOther = pOthers.begin();
         pOther != pOthers.end(); ++pOther ) {
     pMean += pOther->second;
       }
       pMean /= (1 + numOthers_);
       if ( verbosity_ ) std::cout << "pMean = " << pMean << std::endl;
       
       double pErr2 = square(pRef - pMean);
       for ( std::map<std::string, double>::const_iterator pOther = pOthers.begin();
         pOther != pOthers.end(); ++pOther ) {
     pErr2 += square(pOther->second - pMean);
       }
       pErr2 /= numOthers_;
       if ( verbosity_ ) std::cout << "pErr = " << sqrt(pErr2) << std::endl;
 
       if ( pRef > 0. ) {
     weight = pMean/pRef;
     weightErr = sqrt(pErr2)/pRef;
       } else {
     weight = maxWeight_;
     weightErr = maxWeight_;
       }
     }
   } else {
     if ( numWarnings_ < maxWarnings_ ) {
       edm::LogWarning ("<MuonRadiationCorrWeightProducer>")
     << "Failed to find generator level taus matching reconstructed muons !!" << std::endl;
       ++numWarnings_;
     }
   }
 
   if ( weight < minWeight_ ) weight = minWeight_;
   if ( weight > maxWeight_ ) weight = maxWeight_;
   double weightUp = weight + weightErr;
   if ( weightUp > (2.*maxWeight_) ) weightUp = 2.*maxWeight_;
   double weightDown = std::max(weight - weightErr, 0.);
   
   if ( verbosity_ ) {
     std::cout << "--> weight = " << weight << " + " << (weightUp - weight) << " - " << (weight - weightDown) << std::endl;
   }
 
   std::auto_ptr<double> weightPtr(new double(weight));
   evt.put(weightPtr, "weight");
   std::auto_ptr<double> weightUpPtr(new double(weightUp));
   evt.put(weightUpPtr, "weightUp");
   std::auto_ptr<double> weightDownPtr(new double(weightDown));
   evt.put(weightDownPtr, "weightDown");
 }