CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_3/src/Calibration/EcalCalibAlgos/src/ElectronRecalibSuperClusterAssociator.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 //
00003 //
00004 
00005 // user include files
00006 
00007 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
00008 #include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
00009 #include "DataFormats/EgammaReco/interface/BasicClusterFwd.h"
00010 #include "Calibration/EcalCalibAlgos/interface/ElectronRecalibSuperClusterAssociator.h"
00011 
00012 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
00013 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
00014 #include "DataFormats/EgammaCandidates/interface/GsfElectronCoreFwd.h"
00015 #include "DataFormats/EgammaCandidates/interface/GsfElectronCore.h"
00016 
00017 #include <iostream>
00018 
00019 //#define DEBUG 1
00020 
00021 using namespace reco;
00022 using namespace edm;
00023  
00024 ElectronRecalibSuperClusterAssociator::ElectronRecalibSuperClusterAssociator(const edm::ParameterSet& iConfig) 
00025 {
00026 #ifdef DEBUG
00027   std::cout<< "ElectronRecalibSuperClusterAssociator::ElectronRecalibSuperClusterAssociator" << std::endl;
00028 #endif
00029 
00030   //register your products
00031   produces<GsfElectronCollection>();
00032   produces<GsfElectronCoreCollection>() ;
00033   produces<SuperClusterCollection>();
00034   
00035   scProducer_ = iConfig.getParameter<std::string>("scProducer");
00036   scCollection_ = iConfig.getParameter<std::string>("scCollection");
00037 
00038   scIslandProducer_ = iConfig.getParameter<std::string>("scIslandProducer");
00039   scIslandCollection_ = iConfig.getParameter<std::string>("scIslandCollection");
00040 
00041   electronProducer_ = iConfig.getParameter<std::string > ("electronProducer");
00042   electronCollection_ = iConfig.getParameter<std::string > ("electronCollection");
00043 #ifdef DEBUG
00044   std::cout<< "ElectronRecalibSuperClusterAssociator::ElectronRecalibSuperClusterAssociator::end" << std::endl;
00045 #endif
00046 }
00047 
00048 ElectronRecalibSuperClusterAssociator::~ElectronRecalibSuperClusterAssociator()
00049 {
00050 }
00051 
00052 // ------------ method called to produce the data  ------------
00053 void ElectronRecalibSuperClusterAssociator::produce(edm::Event& e, const edm::EventSetup& iSetup) 
00054 {
00055 #ifdef DEBUG
00056   std::cout<< "ElectronRecalibSuperClusterAssociator::produce" << std::endl;
00057 #endif
00058   // Create the output collections   
00059   std::auto_ptr<GsfElectronCollection> pOutEle(new GsfElectronCollection);
00060   std::auto_ptr<GsfElectronCoreCollection> pOutEleCore(new GsfElectronCoreCollection);
00061   std::auto_ptr<SuperClusterCollection> pOutNewEndcapSC(new SuperClusterCollection);
00062 
00063   reco::SuperClusterRefProd rSC = e.getRefBeforePut<SuperClusterCollection>();
00064   edm::Ref<SuperClusterCollection>::key_type idxSC = 0;
00065 
00066   //Get Hybrid SuperClusters
00067   Handle<reco::SuperClusterCollection> pSuperClusters;
00068   e.getByLabel(scProducer_, scCollection_, pSuperClusters);
00069   if (!pSuperClusters.isValid()) {
00070     std::cerr << "Error! can't get the product SuperClusterCollection "<< std::endl;
00071   }
00072   const reco::SuperClusterCollection* scCollection = pSuperClusters.product();
00073   
00074 #ifdef DEBUG
00075   std::cout<<"scCollection->size()"<<scCollection->size()<<std::endl;
00076 #endif
00077   
00078   //Get Island SuperClusters
00079   Handle<reco::SuperClusterCollection> pIslandSuperClusters;
00080   e.getByLabel(scIslandProducer_, scIslandCollection_, pIslandSuperClusters);
00081   if (!pIslandSuperClusters.isValid()) {
00082     std::cerr << "Error! can't get the product IslandSuperClusterCollection "<< std::endl;
00083   }
00084   const reco::SuperClusterCollection* scIslandCollection = pIslandSuperClusters.product();
00085   
00086 #ifdef DEBUG
00087   std::cout<<"scEECollection->size()"<<scIslandCollection->size()<<std::endl;
00088 #endif
00089 
00090   // Get Electrons
00091   Handle<reco::GsfElectronCollection> pElectrons;
00092   e.getByLabel(electronProducer_, electronCollection_, pElectrons);
00093   if (!pElectrons.isValid()) {
00094     std::cerr << "Error! can't get the product ElectronCollection "<< std::endl;
00095   }
00096   const reco::GsfElectronCollection* electronCollection = pElectrons.product();
00097   
00098   for(reco::GsfElectronCollection::const_iterator eleIt = electronCollection->begin(); eleIt != electronCollection->end(); eleIt++)
00099     {
00100       float DeltaRMineleSCbarrel(0.15); 
00101       float DeltaRMineleSCendcap(0.15); 
00102       const reco::SuperCluster* nearestSCbarrel=0;
00103       const reco::SuperCluster* nearestSCendcap=0;
00104       int iscRef=-1;
00105       int iSC=0;
00106       
00107       // first loop is on EB superClusters
00108       for(reco::SuperClusterCollection::const_iterator scIt = scCollection->begin();
00109           scIt != scCollection->end(); scIt++){
00110 #ifdef DEBUG    
00111         std::cout << scIt->energy() << " " << scIt->eta() << " " << scIt->phi() << " " << eleIt->eta() << " " << eleIt->phi() << std::endl;
00112 #endif
00113         
00114         double DeltaReleSC = sqrt ( pow(  eleIt->eta() - scIt->eta(),2) + pow(eleIt->phi() - scIt->phi(),2));
00115         
00116         if(DeltaReleSC<DeltaRMineleSCbarrel)
00117           {
00118             DeltaRMineleSCbarrel = DeltaReleSC;
00119             nearestSCbarrel = &*scIt;
00120             iscRef = iSC;
00121           }
00122         iSC++;
00123       }
00124       iSC = 0;
00125       
00126       // second loop is on EE superClusters
00127       
00128       for(reco::SuperClusterCollection::const_iterator scItEndcap = scIslandCollection->begin();
00129           scItEndcap != scIslandCollection->end(); scItEndcap++){
00130 #ifdef DEBUG    
00131         std::cout << "EE " << scItEndcap->energy() << " " << scItEndcap->eta() << " " << scItEndcap->phi() << " " << eleIt->eta() << " " << eleIt->phi() << std::endl;
00132 #endif
00133         
00134         double DeltaReleSC = sqrt ( pow(  eleIt->eta() - scItEndcap->eta(),2) + pow(eleIt->phi() - scItEndcap->phi(),2));
00135         
00136         if(DeltaReleSC<DeltaRMineleSCendcap)
00137           {
00138             DeltaRMineleSCendcap = DeltaReleSC;
00139             nearestSCendcap = &*scItEndcap;
00140           }
00141         iSC++;
00142       }
00144 
00145       GsfElectronCoreRefProd rEleCore=e.getRefBeforePut<GsfElectronCoreCollection>();
00146       edm::Ref<GsfElectronCoreCollection>::key_type idxEleCore = 0;
00147 
00148       if(nearestSCbarrel && !nearestSCendcap){
00149         reco::GsfElectronCore newEleCore(*(eleIt->core()));
00150         newEleCore.setGsfTrack(eleIt->gsfTrack());
00151         reco::SuperClusterRef scRef(reco::SuperClusterRef(pSuperClusters, iscRef));
00152         newEleCore.setSuperCluster(scRef);
00153         reco::GsfElectronCoreRef newEleCoreRef(reco::GsfElectronCoreRef(rEleCore, idxEleCore ++));
00154         pOutEleCore->push_back(newEleCore);
00155         reco::GsfElectron newEle(*eleIt,newEleCoreRef,CaloClusterPtr(),
00156 //                                TrackRef(),GsfTrackRefVector());
00157                                   TrackRef(),TrackBaseRef(), GsfTrackRefVector());
00158         newEle.setP4(eleIt->p4()*(nearestSCbarrel->energy()/eleIt->ecalEnergy()));
00159 
00160         pOutEle->push_back(newEle);
00161 #ifdef DEBUG
00162         std::cout << "Association is with EB superCluster "<< std::endl;
00163 #endif  
00164       }  
00165 
00166       if(!nearestSCbarrel && nearestSCendcap)
00167         {
00168 #ifdef DEBUG
00169         std::cout << "Starting Association is with EE superCluster "<< std::endl;
00170 #endif  
00171 
00172         float preshowerEnergy=eleIt->superCluster()->preshowerEnergy(); 
00173 #ifdef DEBUG
00174           std::cout << "preshowerEnergy"<< preshowerEnergy << std::endl;
00175 #endif
00176 
00177           CaloClusterPtrVector newBCRef;
00178           for (CaloCluster_iterator bcRefIt=nearestSCendcap->clustersBegin();bcRefIt!=nearestSCendcap->clustersEnd();++bcRefIt){
00179             CaloClusterPtr cPtr(*bcRefIt);
00180             newBCRef.push_back(cPtr);
00181           }
00182          
00183 
00184           reco::SuperCluster newSC(nearestSCendcap->energy() + preshowerEnergy, nearestSCendcap->position() , nearestSCendcap->seed(),newBCRef , preshowerEnergy );
00185           pOutNewEndcapSC->push_back(newSC);
00186           reco::SuperClusterRef scRef(reco::SuperClusterRef(rSC, idxSC ++));
00187 
00188           reco::GsfElectronCore newEleCore(*(eleIt->core()));
00189           newEleCore.setGsfTrack(eleIt->gsfTrack());
00190           newEleCore.setSuperCluster(scRef);
00191           reco::GsfElectronCoreRef newEleCoreRef(reco::GsfElectronCoreRef(rEleCore, idxEleCore ++));
00192           pOutEleCore->push_back(newEleCore);
00193           reco::GsfElectron newEle(*eleIt,newEleCoreRef,CaloClusterPtr(),
00194 //                                TrackRef(),GsfTrackRefVector());
00195                                   TrackRef(),TrackBaseRef(), GsfTrackRefVector());
00196            
00197           newEle.setP4(eleIt->p4()*(newSC.energy()/eleIt->ecalEnergy())) ;
00198           pOutEle->push_back(newEle);
00199 
00200 #ifdef DEBUG
00201         std::cout << "Association is with EE superCluster "<< std::endl;
00202 #endif  
00203       }  
00204     
00205       if(nearestSCbarrel && nearestSCendcap){
00206         reco::GsfElectronCore newEleCore(*(eleIt->core()));
00207         newEleCore.setGsfTrack(eleIt->gsfTrack());
00208 
00209         
00210         if(DeltaRMineleSCendcap>=DeltaRMineleSCbarrel)
00211           {
00212             reco::SuperClusterRef scRef(reco::SuperClusterRef(pSuperClusters, iscRef));
00213             newEleCore.setSuperCluster(scRef);
00214             reco::GsfElectronCoreRef newEleCoreRef(reco::GsfElectronCoreRef(rEleCore, idxEleCore ++));
00215             pOutEleCore->push_back(newEleCore);
00216             reco::GsfElectron newEle(*eleIt,newEleCoreRef,CaloClusterPtr(),
00217 //                                TrackRef(),GsfTrackRefVector());
00218                                   TrackRef(),TrackBaseRef(), GsfTrackRefVector());
00219             newEle.setP4(eleIt->p4()*(nearestSCbarrel->energy()/eleIt->ecalEnergy()));
00220             pOutEle->push_back(newEle);
00221 
00222 
00223 #ifdef DEBUG
00224             std::cout << "Association is with EB superCluster, after quarrel "<< std::endl;
00225 #endif  
00226           }
00227         else if(DeltaRMineleSCendcap<DeltaRMineleSCbarrel)
00228           {
00229             float preshowerEnergy=eleIt->superCluster()->preshowerEnergy(); 
00230             CaloClusterPtrVector newBCRef;
00231             for (CaloCluster_iterator bcRefIt=nearestSCendcap->clustersBegin();bcRefIt!=nearestSCendcap->clustersEnd();++bcRefIt){
00232               CaloClusterPtr cPtr(*bcRefIt);
00233               newBCRef.push_back(*bcRefIt);}
00234             reco::SuperCluster newSC(nearestSCendcap->energy() + preshowerEnergy,  nearestSCendcap->position() , nearestSCendcap->seed(), newBCRef , preshowerEnergy );
00235             pOutNewEndcapSC->push_back(newSC);
00236             reco::SuperClusterRef scRef(reco::SuperClusterRef(rSC, idxSC ++));
00237             newEleCore.setSuperCluster(scRef);
00238             reco::GsfElectronCoreRef newEleCoreRef(reco::GsfElectronCoreRef(rEleCore, idxEleCore ++));
00239             pOutEleCore->push_back(newEleCore);
00240             reco::GsfElectron newEle(*eleIt,newEleCoreRef,CaloClusterPtr(),
00241 //                                TrackRef(),GsfTrackRefVector());
00242                                   TrackRef(),TrackBaseRef(), GsfTrackRefVector());
00243             newEle.setP4(eleIt->p4()*(newSC.energy()/eleIt->ecalEnergy())) ;
00244             pOutEle->push_back(newEle);
00245 #ifdef DEBUG
00246             std::cout << "Association is with EE superCluster, after quarrel "<< std::endl;
00247 #endif  
00248           }     
00249 
00250       }
00251       
00252 
00253     }
00254   
00255   
00256   
00257 #ifdef DEBUG
00258   std::cout << "Filled new electrons  " << pOutEle->size() << std::endl;
00259   std::cout << "Filled new electronsCore  " << pOutEleCore->size() << std::endl;
00260   std::cout << "Filled new endcapSC  " << pOutNewEndcapSC->size() << std::endl;
00261 #endif  
00262   
00263   // put result into the Event
00264 
00265   e.put(pOutEle);
00266   e.put(pOutEleCore);
00267   e.put(pOutNewEndcapSC);
00268   
00269 }
00270 
00271