CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_0/src/SimGeneral/DataMixingModule/plugins/DataMixingEMWorker.cc

Go to the documentation of this file.
00001 // File: DataMixingEMWorker.cc
00002 // Description:  see DataMixingEMWorker.h
00003 // Author:  Mike Hildreth, University of Notre Dame
00004 //
00005 //--------------------------------------------
00006 
00007 #include <map>
00008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00009 #include "FWCore/Utilities/interface/EDMException.h"
00010 #include "FWCore/Framework/interface/ConstProductRegistry.h"
00011 #include "FWCore/ServiceRegistry/interface/Service.h"
00012 #include "DataFormats/Common/interface/Handle.h"
00013 #include "DataFormats/Provenance/interface/Provenance.h"
00014 #include "DataFormats/Provenance/interface/BranchDescription.h"
00015 //
00016 //
00017 #include "DataMixingEMWorker.h"
00018 
00019 
00020 using namespace std;
00021 
00022 namespace edm
00023 {
00024 
00025   // Virtual constructor
00026 
00027   DataMixingEMWorker::DataMixingEMWorker() { sel_=0;}
00028 
00029   // Constructor 
00030   DataMixingEMWorker::DataMixingEMWorker(const edm::ParameterSet& ps) : 
00031                                                             label_(ps.getParameter<std::string>("Label"))
00032 
00033   {                                                         
00034 
00035     // get the subdetector names
00036     //    this->getSubdetectorNames();  //something like this may be useful to check what we are supposed to do...
00037 
00038     // create input selector
00039     if (label_.size()>0){
00040       sel_=new Selector( ModuleLabelSelector(label_));
00041     }
00042     else {
00043       sel_=new Selector( MatchAllSelector());
00044     }
00045 
00046     // declare the products to produce, retrieve
00047 
00048     EBProducerSig_ = ps.getParameter<edm::InputTag>("EBProducerSig");
00049     EEProducerSig_ = ps.getParameter<edm::InputTag>("EEProducerSig");
00050     ESProducerSig_ = ps.getParameter<edm::InputTag>("ESProducerSig");
00051 
00052 
00053     EBrechitCollectionSig_ = ps.getParameter<edm::InputTag>("EBrechitCollectionSig");
00054     EErechitCollectionSig_ = ps.getParameter<edm::InputTag>("EErechitCollectionSig");
00055     ESrechitCollectionSig_ = ps.getParameter<edm::InputTag>("ESrechitCollectionSig");
00056 
00057     EBPileRecHitInputTag_ = ps.getParameter<edm::InputTag>("EBPileRecHitInputTag");
00058     EEPileRecHitInputTag_ = ps.getParameter<edm::InputTag>("EEPileRecHitInputTag");
00059     ESPileRecHitInputTag_ = ps.getParameter<edm::InputTag>("ESPileRecHitInputTag");
00060 
00061 
00062     EBRecHitCollectionDM_        = ps.getParameter<std::string>("EBRecHitCollectionDM");
00063     EERecHitCollectionDM_        = ps.getParameter<std::string>("EERecHitCollectionDM");
00064     ESRecHitCollectionDM_        = ps.getParameter<std::string>("ESRecHitCollectionDM");
00065    //   nMaxPrintout_            = ps.getUntrackedParameter<int>("nMaxPrintout",10);
00066 
00067    //EBalgo_ = new EcalRecHitSimpleAlgo();
00068    //EEalgo_ = new EcalRecHitSimpleAlgo();
00069 
00070    // don't think I can "produce" in a sub-class...
00071 
00072    //produces< EBRecHitCollection >(EBRecHitCollectionDM_);
00073    //produces< EERecHitCollection >(EERecHitCollectionDM_);
00074 
00075   }
00076                
00077 
00078   // Virtual destructor needed.
00079   DataMixingEMWorker::~DataMixingEMWorker() { 
00080     delete sel_;
00081     sel_=0;
00082   }  
00083 
00084   void DataMixingEMWorker::addEMSignals(const edm::Event &e) { 
00085     // fill in maps of hits
00086 
00087     LogInfo("DataMixingEMWorker")<<"===============> adding MC signals for "<<e.id();
00088 
00089     // EB first
00090 
00091    Handle< EBRecHitCollection > pEBRecHits;
00092 
00093    const EBRecHitCollection*  EBRecHits = 0;
00094 
00095    if(e.getByLabel(EBProducerSig_, pEBRecHits) ){
00096      EBRecHits = pEBRecHits.product(); // get a ptr to the product
00097      LogDebug("DataMixingEMWorker") << "total # EB rechits: " << EBRecHits->size();
00098    }
00099    
00100  
00101    if (EBRecHits)
00102      {
00103        // loop over rechits, storing them in a map so we can add pileup later
00104        for(EBRecHitCollection::const_iterator it  = EBRecHits->begin(); 
00105            it != EBRecHits->end(); ++it) {
00106 
00107          EBRecHitStorage_.insert(EBRecHitMap::value_type( ( it->id() ), *it ));
00108          
00109          LogDebug("DataMixingEMWorker") << "processed EBRecHit with rawId: "
00110                                       << it->id().rawId() << "\n"
00111                                       << " rechit energy: " << it->energy();
00112 
00113        }
00114      }
00115 
00116    // EE next
00117 
00118    Handle< EERecHitCollection > pEERecHits;
00119 
00120    const EERecHitCollection*  EERecHits = 0;
00121 
00122    
00123    if(e.getByLabel(EEProducerSig_, pEERecHits) ){
00124      EERecHits = pEERecHits.product(); // get a ptr to the product
00125      LogDebug("DataMixingEMWorker") << "total # EE rechits: " << EERecHits->size();
00126    } 
00127    
00128  
00129    if (EERecHits)
00130      {
00131        // loop over rechits, storing them in a map so we can add pileup later
00132        for(EERecHitCollection::const_iterator it  = EERecHits->begin(); 
00133            it != EERecHits->end(); ++it) {
00134 
00135          EERecHitStorage_.insert(EERecHitMap::value_type( ( it->id() ), *it ));
00136 #ifdef DEBUG     
00137          LogDebug("DataMixingEMWorker") << "processed EERecHit with rawId: "
00138                                       << it->id().rawId() << "\n"
00139                                       << " rechit energy: " << it->energy();
00140 #endif
00141 
00142        }
00143      }
00144    // ES next
00145 
00146    Handle< ESRecHitCollection > pESRecHits;
00147 
00148    const ESRecHitCollection*  ESRecHits = 0;
00149 
00150    
00151    if(e.getByLabel( ESProducerSig_, pESRecHits) ){
00152      ESRecHits = pESRecHits.product(); // get a ptr to the product
00153 #ifdef DEBUG
00154      LogDebug("DataMixingEMWorker") << "total # ES rechits: " << ESRecHits->size();
00155 #endif
00156    } 
00157    
00158  
00159    if (ESRecHits)
00160      {
00161        // loop over rechits, storing them in a map so we can add pileup later
00162        for(ESRecHitCollection::const_iterator it  = ESRecHits->begin(); 
00163            it != ESRecHits->end(); ++it) {
00164 
00165          ESRecHitStorage_.insert(ESRecHitMap::value_type( ( it->id() ), *it ));
00166          
00167 #ifdef DEBUG     
00168          LogDebug("DataMixingEMWorker") << "processed ESRecHit with rawId: "
00169                                       << it->id().rawId() << "\n"
00170                                       << " rechit energy: " << it->energy();
00171 #endif
00172 
00173        }
00174      }
00175     
00176   } // end of addEMSignals
00177 
00178   void DataMixingEMWorker::addEMPileups(const int bcr, const EventPrincipal *ep, unsigned int eventNr) {
00179   
00180     LogInfo("DataMixingEMWorker") <<"\n===============> adding pileups from event  "<<ep->id()<<" for bunchcrossing "<<bcr;
00181 
00182     // fill in maps of hits; same code as addSignals, except now applied to the pileup events
00183 
00184     // EB first
00185 
00186     boost::shared_ptr<Wrapper<EBRecHitCollection>  const> EBRecHitsPTR =
00187       getProductByTag<EBRecHitCollection>(*ep, EBPileRecHitInputTag_ );
00188 
00189     if(EBRecHitsPTR ) {
00190 
00191       const EBRecHitCollection*  EBRecHits = const_cast< EBRecHitCollection * >(EBRecHitsPTR->product());
00192 
00193       LogDebug("DataMixingEMWorker") << "total # EB rechits: " << EBRecHits->size();
00194 
00195       // loop over digis, adding these to the existing maps                                         
00196       for(EBRecHitCollection::const_iterator it  = EBRecHits->begin();
00197           it != EBRecHits->end(); ++it) {
00198 
00199         EBRecHitStorage_.insert(EBRecHitMap::value_type( (it->id()), *it ));
00200 
00201 #ifdef DEBUG
00202         LogDebug("DataMixingEMWorker") << "processed EBRecHit with rawId: "
00203                                            << it->id().rawId() << "\n"
00204                                            << " rechit energy: " << it->energy();
00205 #endif
00206       }
00207     }
00208 
00209     // EE Next
00210 
00211     boost::shared_ptr<Wrapper<EERecHitCollection>  const> EERecHitsPTR =
00212       getProductByTag<EERecHitCollection>(*ep, EEPileRecHitInputTag_ );
00213 
00214     if(EERecHitsPTR ) {
00215 
00216       const EERecHitCollection*  EERecHits = const_cast< EERecHitCollection * >(EERecHitsPTR->product());
00217 
00218       LogDebug("DataMixingEMWorker") << "total # EE rechits: " << EERecHits->size();
00219 
00220       // loop over digis, adding these to the existing maps                                         
00221       for(EERecHitCollection::const_iterator it  = EERecHits->begin();
00222           it != EERecHits->end(); ++it) {
00223 
00224         EERecHitStorage_.insert(EERecHitMap::value_type( (it->id()), *it ));
00225 
00226 #ifdef DEBUG
00227         LogDebug("DataMixingEMWorker") << "processed EERecHit with rawId: "
00228                                        << it->id().rawId() << "\n"
00229                                        << " rechit energy: " << it->energy();
00230 #endif
00231       }
00232     }
00233 
00234     // ES Next
00235 
00236     boost::shared_ptr<Wrapper<ESRecHitCollection>  const> ESRecHitsPTR =
00237       getProductByTag<ESRecHitCollection>(*ep, ESPileRecHitInputTag_ );
00238 
00239     if(ESRecHitsPTR ) {
00240 
00241       const ESRecHitCollection*  ESRecHits = const_cast< ESRecHitCollection * >(ESRecHitsPTR->product());
00242 
00243       LogDebug("DataMixingEMWorker") << "total # ES rechits: " << ESRecHits->size();
00244 
00245       // loop over digis, adding these to the existing maps                                         
00246       for(ESRecHitCollection::const_iterator it  = ESRecHits->begin();
00247           it != ESRecHits->end(); ++it) {
00248 
00249         ESRecHitStorage_.insert(ESRecHitMap::value_type( (it->id()), *it ));
00250 
00251 #ifdef DEBUG
00252         LogDebug("DataMixingEMWorker") << "processed ESRecHit with rawId: "
00253                                        << it->id().rawId() << "\n"
00254                                        << " rechit energy: " << it->energy();
00255 #endif
00256       }
00257     }
00258 
00259 
00260   }
00261  
00262   void DataMixingEMWorker::putEM(edm::Event &e) {
00263 
00264     // collection of rechits to put in the event
00265     std::auto_ptr< EBRecHitCollection > EBrechits( new EBRecHitCollection );
00266     std::auto_ptr< EERecHitCollection > EErechits( new EERecHitCollection );
00267     std::auto_ptr< ESRecHitCollection > ESrechits( new ESRecHitCollection );
00268 
00269     // loop over the maps we have, re-making individual hits or digis if necessary.
00270     DetId formerID = 0;
00271     DetId currentID;
00272     float ESum = 0.;
00273     float EBTime = 0.;
00274 
00275     // EB first...
00276 
00277     EBRecHitMap::const_iterator iEBchk;
00278 
00279     for(EBRecHitMap::const_iterator iEB  = EBRecHitStorage_.begin();
00280         iEB != EBRecHitStorage_.end(); ++iEB) {
00281 
00282       currentID = iEB->first; 
00283 
00284       if (currentID == formerID) { // we have to add these rechits together
00285 
00286         ESum+=(iEB->second).energy(); 
00287       }
00288       else {
00289           if(formerID>0) {
00290             // cutoff for ESum?
00291             EcalRecHit aHit(formerID, ESum, EBTime);
00292             EBrechits->push_back( aHit );
00293           }
00294           //save pointers for next iteration
00295           formerID = currentID;
00296           ESum = (iEB->second).energy();
00297           EBTime = (iEB->second).time();  // take time of first hit in sequence - is this ok?
00298       }
00299 
00300       iEBchk = iEB;
00301       if((++iEBchk) == EBRecHitStorage_.end()) {  //make sure not to lose the last one
00302         EcalRecHit aHit(formerID, ESum, EBTime);
00303         EBrechits->push_back( aHit );     
00304       }
00305     }
00306 
00307     // EE next...
00308 
00309     // loop over the maps we have, re-making individual hits or digis if necessary.
00310     formerID = 0;
00311     ESum = 0.;
00312     float EETime = 0.;
00313     
00314     EERecHitMap::const_iterator iEEchk;
00315 
00316     for(EERecHitMap::const_iterator iEE  = EERecHitStorage_.begin();
00317         iEE != EERecHitStorage_.end(); ++iEE) {
00318 
00319       currentID = iEE->first; 
00320 
00321       if (currentID == formerID) { // we have to add these rechits together
00322 
00323         ESum+=(iEE->second).energy(); 
00324       }
00325       else {
00326           if(formerID>0) {
00327             // cutoff for ESum?
00328             EcalRecHit aHit(formerID, ESum, EETime);
00329             EErechits->push_back( aHit );
00330           }
00331           //save pointers for next iteration
00332           formerID = currentID;
00333           ESum = (iEE->second).energy();
00334           EETime = (iEE->second).time();  // take time of first hit in sequence - is this ok?
00335       }
00336 
00337       iEEchk = iEE;
00338       if((++iEEchk) == EERecHitStorage_.end()) {  //make sure not to lose the last one
00339         EcalRecHit aHit(formerID, ESum, EETime);
00340         EErechits->push_back( aHit );     
00341       }
00342     }
00343 
00344     // ES next...
00345 
00346     // loop over the maps we have, re-making individual hits or digis if necessary.
00347     formerID = 0;
00348     ESum = 0.;
00349     float ESTime = 0.;
00350 
00351     ESRecHitMap::const_iterator iESchk;
00352 
00353     for(ESRecHitMap::const_iterator iES  = ESRecHitStorage_.begin();
00354         iES != ESRecHitStorage_.end(); ++iES) {
00355 
00356       currentID = iES->first; 
00357 
00358       if (currentID == formerID) { // we have to add these rechits together
00359 
00360         ESum+=(iES->second).energy(); 
00361       }
00362       else {
00363           if(formerID>0) {
00364             // cutoff for ESum?
00365             EcalRecHit aHit(formerID, ESum, ESTime);
00366             ESrechits->push_back( aHit );
00367           }
00368           //save pointers for next iteration
00369           formerID = currentID;
00370           ESum = (iES->second).energy();
00371           ESTime = (iES->second).time();  // take time of first hit in sequence - is this ok?
00372       }
00373 
00374       iESchk = iES;
00375       if((++iESchk) == ESRecHitStorage_.end()) {  //make sure not to lose the last one
00376         EcalRecHit aHit(formerID, ESum, ESTime);
00377         ESrechits->push_back( aHit );     
00378       }
00379     }
00380 
00381     // done merging
00382 
00383     // put the collection of reconstructed hits in the event   
00384     LogInfo("DataMixingEMWorker") << "total # EB Merged rechits: " << EBrechits->size() ;
00385     LogInfo("DataMixingEMWorker") << "total # EE Merged rechits: " << EErechits->size() ;
00386     LogInfo("DataMixingEMWorker") << "total # ES Merged rechits: " << ESrechits->size() ;
00387 
00388     e.put( EBrechits, EBRecHitCollectionDM_ );
00389     e.put( EErechits, EERecHitCollectionDM_ );
00390     e.put( ESrechits, ESRecHitCollectionDM_ );
00391     
00392     // clear local storage after this event
00393 
00394     EBRecHitStorage_.clear();
00395     EERecHitStorage_.clear();
00396     ESRecHitStorage_.clear();
00397 
00398   }
00399 
00400 } //edm