CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_6_1_2_SLHC2_patch1/src/TopQuarkAnalysis/TopEventProducers/interface/TtEvtBuilder.h

Go to the documentation of this file.
00001 #ifndef TtEvtBuilder_h
00002 #define TtEvtBuilder_h
00003 
00004 #include <vector>
00005 
00006 #include "FWCore/Framework/interface/Event.h"
00007 #include "FWCore/Framework/interface/EDProducer.h"
00008 #include "FWCore/Framework/interface/Frameworkfwd.h"
00009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00010 #include "AnalysisDataFormats/TopObjects/interface/TtEvent.h"
00011 
00028 template <typename C>
00029 class TtEvtBuilder : public edm::EDProducer {
00030 
00031  public:
00032 
00034   explicit TtEvtBuilder(const edm::ParameterSet&);
00036   ~TtEvtBuilder(){};
00037   
00038  private:
00039 
00042   virtual void produce(edm::Event&, const edm::EventSetup&);
00044   virtual void fillSpecific(C&, const edm::Event&);
00045 
00046  private:
00047 
00049   int verbosity_;
00051   std::vector<edm::InputTag> hyps_;
00053   edm::InputTag genEvt_;
00056   int decayChnTop1_;
00057   int decayChnTop2_;
00058 
00061   edm::ParameterSet kinFit_;
00062   edm::InputTag fitChi2_;
00063   edm::InputTag fitProb_;
00066   edm::ParameterSet hitFit_;
00067   edm::InputTag hitFitChi2_;
00068   edm::InputTag hitFitProb_;
00069   edm::InputTag hitFitMT_;
00070   edm::InputTag hitFitSigMT_;
00073   edm::ParameterSet kinSolution_;
00074   edm::InputTag solWeight_;  
00075   edm::InputTag wrongCharge_;   
00078   edm::ParameterSet genMatch_;
00079   edm::InputTag sumPt_;
00080   edm::InputTag sumDR_;
00083   edm::ParameterSet mvaDisc_;
00084   edm::InputTag meth_;
00085   edm::InputTag disc_;
00086 };
00087 
00088 template <typename C>
00089 TtEvtBuilder<C>::TtEvtBuilder(const edm::ParameterSet& cfg) :
00090   verbosity_   (cfg.getParameter<int>                        ("verbosity"    )),
00091   hyps_        (cfg.getParameter<std::vector<edm::InputTag> >("hypotheses"   )),
00092   genEvt_      (cfg.getParameter<edm::InputTag>              ("genEvent"     )),
00093   decayChnTop1_(cfg.getParameter<int>                        ("decayChannel1")),
00094   decayChnTop2_(cfg.getParameter<int>                        ("decayChannel2"))
00095 {
00096   // parameter subsets for kKinFit
00097   if( cfg.exists("kinFit") ) {
00098     kinFit_  = cfg.getParameter<edm::ParameterSet>("kinFit");
00099     fitChi2_ = kinFit_.getParameter<edm::InputTag>("chi2");
00100     fitProb_ = kinFit_.getParameter<edm::InputTag>("prob");
00101   }
00102   // parameter subsets for kHitFit
00103   if( cfg.exists("hitFit") ) {
00104     hitFit_  = cfg.getParameter<edm::ParameterSet>("hitFit");
00105     hitFitChi2_ = hitFit_.getParameter<edm::InputTag>("chi2");
00106     hitFitProb_ = hitFit_.getParameter<edm::InputTag>("prob");
00107     hitFitMT_ = hitFit_.getParameter<edm::InputTag>("mt");
00108     hitFitSigMT_ = hitFit_.getParameter<edm::InputTag>("sigmt");
00109   }
00110   // parameter subsets for kKinSolution
00111   if( cfg.exists("kinSolution") ) {
00112     kinSolution_  = cfg.getParameter<edm::ParameterSet>("kinSolution");
00113     solWeight_    = kinSolution_.getParameter<edm::InputTag>("solWeight");
00114     wrongCharge_  = kinSolution_.getParameter<edm::InputTag>("wrongCharge");
00115   }  
00116   // parameter subsets for kGenMatch
00117   if( cfg.exists("genMatch") ) {
00118     genMatch_ = cfg.getParameter<edm::ParameterSet>("genMatch");
00119     sumPt_    = genMatch_.getParameter<edm::InputTag>("sumPt");
00120     sumDR_    = genMatch_.getParameter<edm::InputTag>("sumDR");
00121   }
00122   // parameter subsets for kMvaDisc
00123   if( cfg.exists("mvaDisc") ) {
00124     mvaDisc_ = cfg.getParameter<edm::ParameterSet>("mvaDisc");
00125     meth_    = mvaDisc_.getParameter<edm::InputTag>("meth");
00126     disc_    = mvaDisc_.getParameter<edm::InputTag>("disc");
00127   }
00128   // produces a TtEventEvent for:
00129   //  * TtSemiLeptonicEvent 
00130   //  * TtFullLeptonicEvent
00131   //  * TtFullHadronicEvent
00132   // from hypotheses and associated extra information
00133   produces<C>();
00134 }
00135 
00136 template <typename C>
00137 void
00138 TtEvtBuilder<C>::produce(edm::Event& evt, const edm::EventSetup& setup)
00139 {
00140   C ttEvent;
00141 
00142   // set leptonic decay channels
00143   ttEvent.setLepDecays( WDecay::LeptonType(decayChnTop1_), WDecay::LeptonType(decayChnTop2_) );
00144 
00145   // set genEvent (if available)
00146   edm::Handle<TtGenEvent> genEvt;
00147   if( evt.getByLabel(genEvt_, genEvt) )
00148     ttEvent.setGenEvent(genEvt);
00149 
00150   // add event hypotheses for all given 
00151   // hypothesis classes to the TtEvent
00152   typedef std::vector<edm::InputTag>::const_iterator EventHypo;
00153   for(EventHypo h=hyps_.begin(); h!=hyps_.end(); ++h){
00154     edm::Handle<int> key; 
00155     evt.getByLabel(h->label(), "Key", key);
00156 
00157     edm::Handle<std::vector<TtEvent::HypoCombPair> > hypMatchVec; 
00158     evt.getByLabel(*h, hypMatchVec);
00159 
00160     typedef std::vector<TtEvent::HypoCombPair>::const_iterator HypMatch;
00161     for(HypMatch hm=hypMatchVec->begin(); hm != hypMatchVec->end(); ++hm){
00162       ttEvent.addEventHypo((TtEvent::HypoClassKey&)*key, *hm);
00163     }
00164   }
00165 
00166   // set kKinFit extras
00167   if( ttEvent.isHypoAvailable(TtEvent::kKinFit) ) {
00168     edm::Handle<std::vector<double> > fitChi2;
00169     evt.getByLabel(fitChi2_, fitChi2);
00170     ttEvent.setFitChi2( *fitChi2 );
00171     
00172     edm::Handle<std::vector<double> > fitProb;
00173     evt.getByLabel(fitProb_, fitProb);
00174     ttEvent.setFitProb( *fitProb );
00175   }
00176 
00177   // set kHitFit extras
00178   if( ttEvent.isHypoAvailable(TtEvent::kHitFit) ) {
00179     edm::Handle<std::vector<double> > hitFitChi2;
00180     evt.getByLabel(hitFitChi2_, hitFitChi2);
00181     ttEvent.setHitFitChi2( *hitFitChi2 );
00182     
00183     edm::Handle<std::vector<double> > hitFitProb;
00184     evt.getByLabel(hitFitProb_, hitFitProb);
00185     ttEvent.setHitFitProb( *hitFitProb );
00186     
00187     edm::Handle<std::vector<double> > hitFitMT;
00188     evt.getByLabel(hitFitMT_, hitFitMT);
00189     ttEvent.setHitFitMT( *hitFitMT );
00190     
00191     edm::Handle<std::vector<double> > hitFitSigMT;
00192     evt.getByLabel(hitFitSigMT_, hitFitSigMT);
00193     ttEvent.setHitFitSigMT( *hitFitSigMT );
00194   }
00195 
00196   // set kGenMatch extras
00197   if( ttEvent.isHypoAvailable(TtEvent::kGenMatch) ) {
00198     edm::Handle<std::vector<double> > sumPt;
00199     evt.getByLabel(sumPt_, sumPt);
00200     ttEvent.setGenMatchSumPt( *sumPt );
00201 
00202     edm::Handle<std::vector<double> > sumDR;
00203     evt.getByLabel(sumDR_, sumDR);
00204     ttEvent.setGenMatchSumDR( *sumDR );
00205   }
00206 
00207   // set kMvaDisc extras
00208   if( ttEvent.isHypoAvailable(TtEvent::kMVADisc) ) {
00209     edm::Handle<std::string> meth;
00210     evt.getByLabel(meth_, meth);
00211     ttEvent.setMvaMethod( *meth );
00212 
00213     edm::Handle<std::vector<double> > disc;
00214     evt.getByLabel(disc_, disc);
00215     ttEvent.setMvaDiscriminators( *disc );
00216   }
00217 
00218   // fill data members that are decay-channel specific
00219   fillSpecific(ttEvent, evt);
00220 
00221   // print summary via MessageLogger if verbosity_>0
00222   ttEvent.print(verbosity_);
00223 
00224   // write object into the edm::Event
00225   std::auto_ptr<C> pOut(new C);
00226   *pOut=ttEvent;
00227   evt.put(pOut);
00228 }
00229 
00230 template <>
00231 void TtEvtBuilder<TtFullHadronicEvent>::fillSpecific(TtFullHadronicEvent& ttEvent, const edm::Event& evt)
00232 {
00233 }
00234 
00235 template <>
00236 void TtEvtBuilder<TtFullLeptonicEvent>::fillSpecific(TtFullLeptonicEvent& ttEvent, const edm::Event& evt)
00237 {
00238 
00239   // set kKinSolution extras  
00240   if( ttEvent.isHypoAvailable(TtEvent::kKinSolution) ) {
00241     edm::Handle<std::vector<double> > solWeight;
00242     evt.getByLabel(solWeight_, solWeight);
00243     ttEvent.setSolWeight( *solWeight );
00244     
00245     edm::Handle<bool> wrongCharge;
00246     evt.getByLabel(wrongCharge_, wrongCharge);
00247     ttEvent.setWrongCharge( *wrongCharge );   
00248   }
00249 
00250 }
00251 
00252 template <>
00253 void TtEvtBuilder<TtSemiLeptonicEvent>::fillSpecific(TtSemiLeptonicEvent& ttEvent, const edm::Event& evt)
00254 {
00255 
00256   typedef std::vector<edm::InputTag>::const_iterator EventHypo;
00257   for(EventHypo h=hyps_.begin(); h!=hyps_.end(); ++h){
00258     edm::Handle<int> key; 
00259     evt.getByLabel(h->label(), "Key", key);
00260 
00261     // set number of real neutrino solutions for all hypotheses
00262     edm::Handle<int> numberOfRealNeutrinoSolutions;
00263     evt.getByLabel(h->label(), "NumberOfRealNeutrinoSolutions", numberOfRealNeutrinoSolutions);
00264     ttEvent.setNumberOfRealNeutrinoSolutions((TtEvent::HypoClassKey&)*key, *numberOfRealNeutrinoSolutions);
00265 
00266     // set number of considered jets for all hypotheses
00267     edm::Handle<int> numberOfConsideredJets;
00268     evt.getByLabel(h->label(), "NumberOfConsideredJets", numberOfConsideredJets);
00269     ttEvent.setNumberOfConsideredJets((TtEvent::HypoClassKey&)*key, *numberOfConsideredJets);
00270   }
00271 
00272 }
00273 
00274 #endif