CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_7/src/GeneratorInterface/GenFilters/src/MCLongLivedParticles.cc

Go to the documentation of this file.
00001 
00002 #include "GeneratorInterface/GenFilters/interface/MCLongLivedParticles.h"
00003 
00004 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
00005 #include <iostream>
00006 
00007 using namespace edm;
00008 using namespace std;
00009 
00010 
00011 MCLongLivedParticles::MCLongLivedParticles(const edm::ParameterSet& iConfig) :
00012   hepMCProductTag_(iConfig.getParameter<edm::InputTag>("hepMCProductTag")) {
00013    //here do whatever other initialization is needed
00014   theCut = iConfig.getUntrackedParameter<double>("LengCut",10.);
00015 }
00016 
00017 
00018 MCLongLivedParticles::~MCLongLivedParticles()
00019 {
00020  
00021    // do anything here that needs to be done at desctruction time
00022    // (e.g. close files, deallocate resources etc.)
00023 
00024 }
00025 
00026 
00027 // ------------ method called to skim the data  ------------
00028 bool MCLongLivedParticles::filter(edm::Event& iEvent, const edm::EventSetup& iSetup)
00029 {
00030 
00031   using namespace edm;
00032   
00033   Handle<HepMCProduct> evt;
00034   
00035   iEvent.getByLabel(hepMCProductTag_, evt);
00036   
00037   bool pass = false;
00038   
00039   const HepMC::GenEvent * generated_event = evt->GetEvent();
00040   HepMC::GenEvent::particle_const_iterator p;
00041   
00042   for (p = generated_event->particles_begin(); p != generated_event->particles_end(); p++)
00043     { 
00044 
00045       if((*p)->production_vertex()!=0&&(*p)->end_vertex()!=0)
00046         {
00047           float dist = sqrt((((*p)->production_vertex())->position().x()-((*p)->end_vertex())->position().x())*(((*p)->production_vertex())->position().x()-((*p)->end_vertex())->position().x())+
00048                             (((*p)->production_vertex())->position().y()-((*p)->end_vertex())->position().y())*(((*p)->production_vertex())->position().y()-((*p)->end_vertex())->position().y()));
00049           if(dist>theCut)
00050             pass=true;
00051         }
00052       
00053       if((*p)->production_vertex()==0&&!(*p)->end_vertex()!=0)
00054         {
00055           if(((*p)->end_vertex())->position().perp()>theCut)
00056             pass=true;
00057         }
00058       
00059       if(pass)
00060         return pass;
00061     }
00062   
00063   return pass;
00064 }
00065