00001 #include "RecoEcal/EgammaClusterProducers/interface/InterestingDetIdCollectionProducer.h"
00002
00003 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00004 #include "FWCore/Framework/interface/EventSetup.h"
00005 #include "FWCore/Framework/interface/ESHandle.h"
00006
00007 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
00008
00009 #include "DataFormats/EcalDetId/interface/EBDetId.h"
00010 #include "DataFormats/EcalDetId/interface/EEDetId.h"
00011
00012 #include "DataFormats/DetId/interface/DetIdCollection.h"
00013 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
00014 #include "DataFormats/EgammaReco/interface/BasicClusterFwd.h"
00015
00016 #include "Geometry/CaloEventSetup/interface/CaloTopologyRecord.h"
00017 #include "Geometry/CaloTopology/interface/CaloTopology.h"
00018 #include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
00019
00020
00021
00022 InterestingDetIdCollectionProducer::InterestingDetIdCollectionProducer(const edm::ParameterSet& iConfig)
00023 {
00024
00025 recHitsLabel_ = iConfig.getParameter< edm::InputTag > ("recHitsLabel");
00026 basicClusters_ = iConfig.getParameter< edm::InputTag > ("basicClustersLabel");
00027
00028 interestingDetIdCollection_ = iConfig.getParameter<std::string>("interestingDetIdCollection");
00029
00030 minimalEtaSize_ = iConfig.getParameter<int> ("etaSize");
00031 minimalPhiSize_ = iConfig.getParameter<int> ("phiSize");
00032 if ( minimalPhiSize_ % 2 == 0 || minimalEtaSize_ % 2 == 0)
00033 edm::LogError("InterestingDetIdCollectionProducerError") << "Size of eta/phi should be odd numbers";
00034
00035
00036 produces< DetIdCollection > (interestingDetIdCollection_) ;
00037
00038 }
00039
00040
00041 InterestingDetIdCollectionProducer::~InterestingDetIdCollectionProducer()
00042 {}
00043
00044 void InterestingDetIdCollectionProducer::beginJob (const edm::EventSetup& iSetup)
00045 {
00046 edm::ESHandle<CaloTopology> theCaloTopology;
00047 iSetup.get<CaloTopologyRecord>().get(theCaloTopology);
00048 caloTopology_ = &(*theCaloTopology);
00049 }
00050
00051
00052 void
00053 InterestingDetIdCollectionProducer::produce (edm::Event& iEvent,
00054 const edm::EventSetup& iSetup)
00055 {
00056 using namespace edm;
00057 using namespace std;
00058
00059
00060 Handle<reco::BasicClusterCollection> pClusters;
00061 iEvent.getByLabel(basicClusters_, pClusters);
00062
00063
00064 Handle<EcalRecHitCollection> recHitsHandle;
00065 iEvent.getByLabel(recHitsLabel_,recHitsHandle);
00066
00067
00068 std::auto_ptr< DetIdCollection > detIdCollection (new DetIdCollection() ) ;
00069
00070 reco::BasicClusterCollection::const_iterator clusIt;
00071
00072 for (clusIt=pClusters->begin(); clusIt!=pClusters->end(); clusIt++) {
00073
00074
00075 float eMax=0.;
00076 DetId eMaxId(0);
00077
00078 std::vector<DetId> clusterDetIds = (*clusIt).getHitsByDetId();
00079 std::vector<DetId>::iterator posCurrent;
00080
00081 EcalRecHit testEcalRecHit;
00082
00083 for(posCurrent = clusterDetIds.begin(); posCurrent != clusterDetIds.end(); posCurrent++)
00084 {
00085 EcalRecHitCollection::const_iterator itt = recHitsHandle->find(*posCurrent);
00086 if ((!((*posCurrent).null())) && (itt != recHitsHandle->end()) && ((*itt).energy() > eMax) )
00087 {
00088 eMax = (*itt).energy();
00089 eMaxId = (*itt).id();
00090 }
00091 }
00092
00093 if (eMaxId.null())
00094 continue;
00095
00096 const CaloSubdetectorTopology* topology = caloTopology_->getSubdetectorTopology(eMaxId.det(),eMaxId.subdetId());
00097 std::vector<DetId> xtalsToStore=topology->getWindow(eMaxId,minimalEtaSize_,minimalPhiSize_);
00098 std::vector<DetId> xtalsInClus=(*clusIt).getHitsByDetId();
00099
00100 for (unsigned int ii=0;ii<xtalsInClus.size();ii++)
00101 {
00102 if (std::find(xtalsToStore.begin(),xtalsToStore.end(),xtalsInClus[ii]) == xtalsToStore.end())
00103 xtalsToStore.push_back(xtalsInClus[ii]);
00104 }
00105
00106 for (unsigned int iCry=0;iCry<xtalsToStore.size();iCry++)
00107 {
00108 if (
00109 std::find(detIdCollection->begin(),detIdCollection->end(),xtalsToStore[iCry]) == detIdCollection->end()
00110 )
00111 detIdCollection->push_back(xtalsToStore[iCry]);
00112 }
00113 }
00114
00115
00116 iEvent.put( detIdCollection, interestingDetIdCollection_ );
00117
00118 }