StEvtSolutionMaker.cc

Go to the documentation of this file.

//
//

#include <memory>

#include "TopQuarkAnalysis/TopEventProducers/interface/StEvtSolutionMaker.h"

StEvtSolutionMaker::StEvtSolutionMaker(const edm::ParameterSet& iConfig)
{
  // configurables
  electronSrcToken_    = mayConsume<std::vector<pat::Electron> >(iConfig.getParameter<edm::InputTag>("electronSource"));
  muonSrcToken_        = mayConsume<std::vector<pat::Muon> >(iConfig.getParameter<edm::InputTag>("muonSource"));
  metSrcToken_         = consumes<std::vector<pat::MET> >(iConfig.getParameter<edm::InputTag>("metSource"));
  jetSrcToken_         = consumes<std::vector<pat::Jet> >(iConfig.getParameter<edm::InputTag>("jetSource"));
  genEvtSrcToken_      = mayConsume<StGenEvent>(edm::InputTag("genEvt"));
  leptonFlavour_  = iConfig.getParameter< std::string >("leptonFlavour");
  jetCorrScheme_  = iConfig.getParameter<int>          ("jetCorrectionScheme");
  //jetInput_        = iConfig.getParameter< std::string >    ("jetInput");
  doKinFit_       = iConfig.getParameter< bool >       ("doKinFit");
  addLRJetComb_   = iConfig.getParameter< bool >       ("addLRJetComb");
  maxNrIter_      = iConfig.getParameter< int >        ("maxNrIter");
  maxDeltaS_      = iConfig.getParameter< double >     ("maxDeltaS");
  maxF_           = iConfig.getParameter< double >     ("maxF");
  jetParam_       = iConfig.getParameter<int>          ("jetParametrisation");
  lepParam_       = iConfig.getParameter<int>          ("lepParametrisation");
  metParam_       = iConfig.getParameter<int>          ("metParametrisation");
  constraints_    = iConfig.getParameter< std::vector<int> > ("constraints");
  matchToGenEvt_  = iConfig.getParameter< bool >    ("matchToGenEvt");

  // define kinfitter
  if(doKinFit_){
    myKinFitter = new StKinFitter(jetParam_, lepParam_, metParam_, maxNrIter_, maxDeltaS_, maxF_, constraints_);
  }
  // define what will be produced
  produces<std::vector<StEvtSolution> >();
}

StEvtSolutionMaker::~StEvtSolutionMaker()
{
  if (doKinFit_) delete myKinFitter;
}

void StEvtSolutionMaker::produce(edm::Event& iEvent, const edm::EventSetup& iSetup)
{
  //
  //  TopObject Selection
  //

  // select lepton (the TtLepton vectors are, for the moment, sorted on pT)
  bool leptonFound = false;
  edm::Handle<std::vector<pat::Muon> > muons;
  if(leptonFlavour_ == "muon"){
    iEvent.getByToken(muonSrcToken_, muons);
    if( !muons->empty() ) leptonFound = true;
  }
  edm::Handle<std::vector<pat::Electron> > electrons;
  if(leptonFlavour_ == "electron"){
    iEvent.getByToken(electronSrcToken_, electrons);
    if( !electrons->empty() ) leptonFound = true;
  }

  // select MET (TopMET vector is sorted on ET)
  bool metFound = false;
  edm::Handle<std::vector<pat::MET> > mets;
  iEvent.getByToken(metSrcToken_, mets);
  if( !mets->empty() ) metFound = true;

  // select Jets
  bool jetsFound = false;
  edm::Handle<std::vector<pat::Jet> > jets;
  iEvent.getByToken(jetSrcToken_, jets);
  unsigned int maxJets=2;//this has to become a custom-defined parameter (we may want 2 or 3 jets)
  if (jets->size() >= 2) jetsFound = true;

  std::vector<StEvtSolution> *evtsols = new std::vector<StEvtSolution>();
  if(leptonFound && metFound && jetsFound){
    std::cout<<"constructing solutions"<<std::endl;
    for (unsigned int b=0; b<maxJets; b++) {
      for (unsigned int l=0; l<maxJets; l++) {
    if(b!=l){  // to avoid double counting
      StEvtSolution asol;
      asol.setJetCorrectionScheme(jetCorrScheme_);
      if(leptonFlavour_ == "muon")     asol.setMuon(muons, 0);
      if(leptonFlavour_ == "electron") asol.setElectron(electrons, 0);
      asol.setNeutrino(mets, 0);
      asol.setBottom(jets, b);
      asol.setLight(jets, l);

      if(doKinFit_) asol = myKinFitter->addKinFitInfo(&asol);

      /* to be adapted to ST (Andrea)
         if(addLRJetComb_){
         asol.setPtrueCombExist(jetCombProbs[m].getPTrueCombExist(&afitsol));
         asol.setPtrueBJetSel(jetCombProbs[m].getPTrueBJetSel(&afitsol));
         asol.setPtrueBhadrSel(jetCombProbs[m].getPTrueBhadrSel(&afitsol));
         asol.setPtrueJetComb(afitsol.getPtrueCombExist()*afitsol.getPtrueBJetSel()*afitsol.getPtrueBhadrSel());
         }
      */
      evtsols->push_back(asol);
    }
      }
    }

    // if asked for, match the event solutions to the gen Event
    if(matchToGenEvt_){
      /*
    edm::Handle<StGenEvent> genEvt;
    iEvent.getByToken (genEvtSrcToken_,genEvt);
    double bestSolDR = 9999.;
    int bestSol = 0;
    for(size_t s=0; s<evtsols->size(); s++) {
    (*evtsols)[s].setGenEvt(genEvt->particles());
    vector<double> bm = BestMatch((*evtsols)[s], false); //false to use DR, true SpaceAngles
    (*evtsols)[s].setSumDeltaRjp(bm[0]); // dRBB + dRLL
    (*evtsols)[s].setChangeBL((int) bm[1]); // has swapped or not
    (*evtsols)[s].setDeltaRB(bm[2]);
    (*evtsols)[s].setDeltaRL(bm[3]);
    if(bm[0]<bestSolDR) { bestSolDR =  bm[0]; bestSol = s; }
    }
    (*evtsols)[bestSol].setBestSol(true);
      */
    }

    //store the vector of solutions to the event
    std::unique_ptr<std::vector<StEvtSolution> > pOut(evtsols);
    iEvent.put(std::move(pOut));
  }
  else
    {

      std::cout<<"@@@ No calibrated solutions built, because:  " << std::endl;;
      if(jets->size()<maxJets)                    std::cout<<"@ nr jets = " << jets->size() << " < " << maxJets <<std::endl;
      if(leptonFlavour_ == "muon" && !leptonFound)                std::cout<<"@ no good muon candidate"<<std::endl;
      if(leptonFlavour_ == "electron" && !leptonFound)             std::cout<<"@ no good electron candidate"<<std::endl;
      if(mets->empty())                       std::cout<<"@ no MET reconstruction"<<std::endl;

      StEvtSolution asol;
      evtsols->push_back(asol);
      std::unique_ptr<std::vector<StEvtSolution> > pOut(evtsols);
      iEvent.put(std::move(pOut));
    }
}