CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_5/src/RecoHI/HiJetAlgos/plugins/HiGenCleaner.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 //
00003 // Package:    HiGenCleaner
00004 // Class:      HiGenCleaner
00005 // 
00013 //
00014 // Original Author:  Yetkin Yilmaz
00015 //         Created:  Tue Jul 21 04:26:01 EDT 2009
00016 // $Id: HiGenCleaner.cc,v 1.3 2010/06/01 14:44:22 yilmaz Exp $
00017 //
00018 //
00019 
00020 // system include files
00021 #include <memory>
00022 #include <vector>
00023 
00024 // user include files
00025 #include "FWCore/Framework/interface/Frameworkfwd.h"
00026 #include "FWCore/Framework/interface/EDProducer.h"
00027 
00028 #include "FWCore/Framework/interface/Event.h"
00029 #include "FWCore/Framework/interface/MakerMacros.h"
00030 
00031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00032 #include "FWCore/Utilities/interface/InputTag.h"
00033 
00034 #include "DataFormats/Common/interface/View.h"
00035 #include "DataFormats/JetReco/interface/GenJetCollection.h"
00036 #include "DataFormats/GeometryVector/interface/VectorUtil.h"
00037 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
00038 #include "DataFormats/HepMCCandidate/interface/GenParticleFwd.h"
00039 
00040 #include "DataFormats/Math/interface/Point3D.h"
00041 #include "DataFormats/Math/interface/LorentzVector.h"
00042 
00043 using namespace std;
00044 using namespace edm;
00045 
00046 
00047 //
00048 // class decleration
00049 //
00050 
00051 template <class T2>
00052 class HiGenCleaner : public edm::EDProducer {
00053 public:
00054   typedef std::vector<T2> T2Collection;
00055       explicit HiGenCleaner(const edm::ParameterSet&);
00056       ~HiGenCleaner();
00057 
00058    private:
00059       virtual void produce(edm::Event&, const edm::EventSetup&);
00060       // ----------member data ---------------------------
00061 
00062    InputTag jetSrc_;
00063    double deltaR_;
00064    double ptCut_;
00065    bool makeNew_;
00066    bool fillDummy_;
00067 
00068 };
00069 
00070 //
00071 // constants, enums and typedefs
00072 //
00073 
00074 
00075 //
00076 // static data member definitions
00077 //
00078 
00079 //
00080 // constructors and destructor
00081 //
00082 
00083 template <class T2>
00084 HiGenCleaner<T2>::HiGenCleaner(const edm::ParameterSet& iConfig) :
00085   jetSrc_(iConfig.getParameter<InputTag>( "src")),
00086   deltaR_(iConfig.getParameter<double>("deltaR")),
00087   ptCut_(iConfig.getParameter<double>("ptCut")),
00088   makeNew_(iConfig.getUntrackedParameter<bool>("createNewCollection",true)),
00089   fillDummy_(iConfig.getUntrackedParameter<bool>("fillDummyEntries",true))
00090 {
00091    std::string alias = jetSrc_.label();
00092    produces<T2Collection>().setBranchAlias (alias);
00093 }
00094 
00095 template <class T2>
00096 HiGenCleaner<T2>::~HiGenCleaner()
00097 {
00098    // do anything here that needs to be done at desctruction time
00099    // (e.g. close files, deallocate resources etc.)
00100 }
00101 
00102 
00103 //
00104 // member functions
00105 //
00106 
00107 // ------------ method called to produce the data  ------------
00108 template <class T2>
00109 void
00110 HiGenCleaner<T2>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup)
00111 {
00112    using namespace edm;
00113    using namespace reco;
00114 
00115    auto_ptr<T2Collection> jets;
00116    jets = auto_ptr<T2Collection>(new T2Collection);
00117    
00118    edm::Handle<edm::View<T2> > genjets;
00119    iEvent.getByLabel(jetSrc_,genjets);
00120 
00121    int jetsize = genjets->size();
00122 
00123    vector<int> selection;
00124    for(int ijet = 0; ijet < jetsize; ++ijet){
00125       selection.push_back(-1);
00126    }
00127 
00128    vector<int> selectedIndices;
00129    vector<int> removedIndices;
00130 
00131    for(int ijet = 0; ijet < jetsize; ++ijet){
00132 
00133      const T2* jet1 = &((*genjets)[ijet]);
00134      
00135      if(selection[ijet] == -1){
00136        selection[ijet] = 1;
00137        for(int ijet2 = 0; ijet2 < jetsize; ++ijet2){
00138 
00139          if(ijet2 == ijet) continue;
00140          
00141          const T2* jet2 = &((*genjets)[ijet2]);
00142             
00143          if(Geom::deltaR(jet1->momentum(),jet2->momentum()) < deltaR_){
00144            if(jet1->et() < jet2->et()){
00145              selection[ijet] = 0;
00146              removedIndices.push_back(ijet);
00147              break;
00148            }else{
00149              selection[ijet2] = 0;
00150              removedIndices.push_back(ijet2);
00151            }
00152          }
00153          }
00154      }
00155      
00156      double etjet = ((*genjets)[ijet]).et();
00157       
00158       if(selection[ijet] == 1 && etjet > ptCut_){ 
00159          selectedIndices.push_back(ijet);
00160          jets->push_back(*jet1);
00161       }
00162    }
00163    iEvent.put(jets);
00164 
00165 }
00166 
00167 typedef HiGenCleaner<reco::GenParticle> HiPartonCleaner;
00168 typedef HiGenCleaner<reco::GenJet> HiGenJetCleaner;
00169 
00170 DEFINE_FWK_MODULE(HiPartonCleaner);
00171 DEFINE_FWK_MODULE(HiGenJetCleaner);
00172