Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:32:25

0001 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0002 #include "FWCore/Framework/interface/Event.h"
0003 #include "FWCore/Framework/interface/EventSetup.h"
0004 #include "DataFormats/Common/interface/Handle.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "FWCore/Utilities/interface/InputTag.h"
0007 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0008 #include "SimDataFormats/GeneratorProducts/interface/LHEEventProduct.h"
0009 #include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h"
0010 #include "SimGeneral/HepPDTRecord/interface/ParticleDataTable.h"
0011 #include "PhysicsTools/HepMCCandAlgos/interface/PDFWeightsHelper.h"
0012 #include "FWCore/Utilities/interface/EDMException.h"
0013 #include "FWCore/ServiceRegistry/interface/Service.h"
0014 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0015 #include "TTree.h"
0016 #include <iostream>
0017 using namespace std;
0018 using namespace edm;
0019 using namespace HepMC;
0020 
0021 class PDFWeightsTest : public edm::one::EDAnalyzer<> {
0022 private:
0023   PDFWeightsHelper pdfweightshelper_;
0024   EDGetTokenT<LHEEventProduct> srcToken_;
0025   EDGetTokenT<LHEEventProduct> srcTokenAlt_;
0026   EDGetTokenT<GenEventInfoProduct> srcTokenGen_;
0027 
0028   unsigned int pdfWeightOffset_;
0029   unsigned int nPdfWeights_;
0030   unsigned int nPdfEigWeights_;
0031 
0032   TTree* tree_;
0033   std::vector<float> pdfweights_;
0034   std::vector<float> pdfeigweights_;
0035   float weight_;
0036 
0037 public:
0038   explicit PDFWeightsTest(const ParameterSet& cfg)
0039       : srcToken_(consumes<LHEEventProduct>(edm::InputTag("externalLHEProducer"))),
0040         srcTokenAlt_(consumes<LHEEventProduct>(edm::InputTag("source"))),
0041         srcTokenGen_(consumes<GenEventInfoProduct>(edm::InputTag("generator"))),
0042         pdfWeightOffset_(cfg.getParameter<unsigned int>("pdfWeightOffset")),
0043         nPdfWeights_(cfg.getParameter<unsigned int>("nPdfWeights")),
0044         nPdfEigWeights_(cfg.getParameter<unsigned int>("nPdfEigWeights")),
0045         pdfweights_(nPdfWeights_),
0046         pdfeigweights_(nPdfEigWeights_) {
0047     edm::Service<TFileService> fs;
0048     tree_ = fs->make<TTree>("tree", "");
0049 
0050     tree_->Branch("pdfrep", &pdfweights_);
0051     tree_->Branch("pdfeig", &pdfeigweights_);
0052     tree_->Branch("weight", &weight_);
0053 
0054     edm::FileInPath mc2hessianCSV = cfg.getParameter<edm::FileInPath>("mc2hessianCSV");
0055     pdfweightshelper_.Init(nPdfWeights_, nPdfEigWeights_, mc2hessianCSV);
0056   }
0057 
0058 private:
0059   void analyze(const Event& evt, const EventSetup& es) override {
0060     Handle<LHEEventProduct> lheInfo;
0061     evt.getByToken(srcToken_, lheInfo);
0062 
0063     if (!lheInfo.isValid()) {
0064       evt.getByToken(srcTokenAlt_, lheInfo);
0065     }
0066 
0067     double nomlheweight = lheInfo->weights()[0].wgt;
0068 
0069     Handle<GenEventInfoProduct> genInfo;
0070     evt.getByToken(srcTokenGen_, genInfo);
0071 
0072     weight_ = genInfo->weight();
0073 
0074     //get the original mc replica weights
0075     std::vector<double> inpdfweights(nPdfWeights_);
0076     for (unsigned int ipdf = 0; ipdf < nPdfWeights_; ++ipdf) {
0077       unsigned int iwgt = ipdf + pdfWeightOffset_;
0078 
0079       //this is the weight to be used for evaluating uncertainties with mc replica weights
0080       pdfweights_[ipdf] = lheInfo->weights()[iwgt].wgt * weight_ / nomlheweight;
0081 
0082       //this is the raw weight to be fed to the mc2hessian convertor
0083       inpdfweights[ipdf] = lheInfo->weights()[iwgt].wgt;
0084     }
0085 
0086     std::vector<double> outpdfweights(nPdfEigWeights_);
0087     //do the actual conversion, where the nominal lhe weight is needed as the reference point for the linearization
0088     pdfweightshelper_.DoMC2Hessian(nomlheweight, inpdfweights.data(), outpdfweights.data());
0089 
0090     for (unsigned int iwgt = 0; iwgt < nPdfEigWeights_; ++iwgt) {
0091       double wgtval = outpdfweights[iwgt];
0092 
0093       //the is the weight to be used for evaluating uncertainties with hessian weights
0094       pdfeigweights_[iwgt] = wgtval * weight_ / nomlheweight;
0095     }
0096 
0097     tree_->Fill();
0098   }
0099 };
0100 
0101 #include "FWCore/Framework/interface/MakerMacros.h"
0102 
0103 DEFINE_FWK_MODULE(PDFWeightsTest);