CMS 3D CMS Logo

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

Go to the documentation of this file.
00001 // File: DataMixingHcalWorker.cc
00002 // Description:  see DataMixingHcalWorker.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 "DataMixingHcalWorker.h"
00018 
00019 
00020 using namespace std;
00021 
00022 namespace edm
00023 {
00024 
00025   // Virtual constructor
00026 
00027   DataMixingHcalWorker::DataMixingHcalWorker() {sel_=0;} 
00028 
00029   // Constructor 
00030   DataMixingHcalWorker::DataMixingHcalWorker(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
00047 
00048     // Hcal 
00049 
00050     HBHErechitCollectionSig_  = ps.getParameter<edm::InputTag>("HBHEProducerSig");
00051     HOrechitCollectionSig_    = ps.getParameter<edm::InputTag>("HOProducerSig");
00052     HFrechitCollectionSig_    = ps.getParameter<edm::InputTag>("HFProducerSig");
00053     ZDCrechitCollectionSig_   = ps.getParameter<edm::InputTag>("ZDCrechitCollectionSig");
00054 
00055     HBHEPileRecHitInputTag_ = ps.getParameter<edm::InputTag>("HBHEPileRecHitInputTag");
00056     HOPileRecHitInputTag_ = ps.getParameter<edm::InputTag>("HOPileRecHitInputTag");
00057     HFPileRecHitInputTag_ = ps.getParameter<edm::InputTag>("HFPileRecHitInputTag");
00058     ZDCPileRecHitInputTag_ = ps.getParameter<edm::InputTag>("ZDCPileRecHitInputTag");
00059 
00060     HBHERecHitCollectionDM_ = ps.getParameter<std::string>("HBHERecHitCollectionDM");
00061     HORecHitCollectionDM_   = ps.getParameter<std::string>("HORecHitCollectionDM");
00062     HFRecHitCollectionDM_   = ps.getParameter<std::string>("HFRecHitCollectionDM");
00063     ZDCRecHitCollectionDM_  = ps.getParameter<std::string>("ZDCRecHitCollectionDM");
00064 
00065 
00066   }
00067                
00068   // Virtual destructor needed.
00069   DataMixingHcalWorker::~DataMixingHcalWorker() { 
00070     delete sel_;
00071     sel_=0;
00072   }  
00073 
00074   void DataMixingHcalWorker::addHcalSignals(const edm::Event &e) { 
00075     // fill in maps of hits
00076 
00077     LogInfo("DataMixingHcalWorker")<<"===============> adding MC signals for "<<e.id();
00078 
00079     // HBHE first
00080 
00081    Handle< HBHERecHitCollection > pHBHERecHits;
00082 
00083    const HBHERecHitCollection*  HBHERecHits = 0;
00084 
00085    if( e.getByLabel( HBHErechitCollectionSig_, pHBHERecHits) ) {
00086      HBHERecHits = pHBHERecHits.product(); // get a ptr to the product
00087      LogDebug("DataMixingHcalWorker") << "total # HBHE rechits: " << HBHERecHits->size();
00088    } 
00089    
00090  
00091    if (HBHERecHits)
00092      {
00093        // loop over rechits, storing them in a map so we can add pileup later
00094        for(HBHERecHitCollection::const_iterator it  = HBHERecHits->begin();     
00095            it != HBHERecHits->end(); ++it) {
00096 
00097          HBHERecHitStorage_.insert(HBHERecHitMap::value_type( ( it->id() ), *it ));
00098          
00099 #ifdef DEBUG     
00100          LogDebug("DataMixingHcalWorker") << "processed HBHERecHit with rawId: "
00101                                       << it->id() << "\n"
00102                                       << " rechit energy: " << it->energy();
00103 #endif
00104 
00105        }
00106      }
00107 
00108    // HO next
00109 
00110    Handle< HORecHitCollection > pHORecHits;
00111 
00112    const HORecHitCollection*  HORecHits = 0;
00113 
00114    if( e.getByLabel( HOrechitCollectionSig_, pHORecHits) ){
00115      HORecHits = pHORecHits.product(); // get a ptr to the product
00116 #ifdef DEBUG
00117      LogDebug("DataMixingHcalWorker") << "total # HO rechits: " << HORecHits->size();
00118 #endif
00119    } 
00120    
00121  
00122    if (HORecHits)
00123      {
00124        // loop over rechits, storing them in a map so we can add pileup later
00125        for(HORecHitCollection::const_iterator it  = HORecHits->begin(); 
00126            it != HORecHits->end(); ++it) {
00127 
00128          HORecHitStorage_.insert(HORecHitMap::value_type( ( it->id() ), *it ));
00129          
00130 #ifdef DEBUG     
00131          LogDebug("DataMixingHcalWorker") << "processed HORecHit with rawId: "
00132                                       << it->id() << "\n"
00133                                       << " rechit energy: " << it->energy();
00134 #endif
00135 
00136        }
00137      }
00138 
00139    // HF next
00140 
00141    Handle< HFRecHitCollection > pHFRecHits;
00142 
00143    const HFRecHitCollection*  HFRecHits = 0;
00144 
00145    if( e.getByLabel( HFrechitCollectionSig_, pHFRecHits) ) {
00146      HFRecHits = pHFRecHits.product(); // get a ptr to the product
00147 #ifdef DEBUG
00148      LogDebug("DataMixingHcalWorker") << "total # HF rechits: " << HFRecHits->size();
00149 #endif
00150    } 
00151    
00152  
00153    if (HFRecHits)
00154      {
00155        // loop over rechits, storing them in a map so we can add pileup later
00156        for(HFRecHitCollection::const_iterator it  = HFRecHits->begin(); 
00157            it != HFRecHits->end(); ++it) {
00158 
00159          HFRecHitStorage_.insert(HFRecHitMap::value_type( ( it->id() ), *it ));
00160          
00161 #ifdef DEBUG     
00162          LogDebug("DataMixingHcalWorker") << "processed HFRecHit with rawId: "
00163                                       << it->id() << "\n"
00164                                       << " rechit energy: " << it->energy();
00165 #endif
00166 
00167        }
00168      }
00169 
00170    // ZDC next
00171 
00172    Handle< ZDCRecHitCollection > pZDCRecHits;
00173 
00174    const ZDCRecHitCollection*  ZDCRecHits = 0;
00175 
00176    if( e.getByLabel( ZDCrechitCollectionSig_, pZDCRecHits) ) {
00177      ZDCRecHits = pZDCRecHits.product(); // get a ptr to the product
00178 #ifdef DEBUG
00179      LogDebug("DataMixingHcalWorker") << "total # ZDC rechits: " << ZDCRecHits->size();
00180 #endif
00181    } 
00182    
00183  
00184    if (ZDCRecHits)
00185      {
00186        // loop over rechits, storing them in a map so we can add pileup later
00187        for(ZDCRecHitCollection::const_iterator it  = ZDCRecHits->begin();       
00188            it != ZDCRecHits->end(); ++it) {
00189 
00190          ZDCRecHitStorage_.insert(ZDCRecHitMap::value_type( ( it->id() ), *it ));
00191          
00192 #ifdef DEBUG     
00193          LogDebug("DataMixingHcalWorker") << "processed ZDCRecHit with rawId: "
00194                                       << it->id() << "\n"
00195                                       << " rechit energy: " << it->energy();
00196 #endif
00197 
00198        }
00199      }
00200     
00201   } // end of addEMSignals
00202 
00203   void DataMixingHcalWorker::addHcalPileups(const int bcr, const EventPrincipal *ep, unsigned int eventNr) {
00204   
00205     LogDebug("DataMixingHcalWorker") <<"\n===============> adding pileups from event  "<<ep->id()<<" for bunchcrossing "<<bcr;
00206 
00207     // fill in maps of hits; same code as addSignals, except now applied to the pileup events
00208 
00209     // HBHE first
00210 
00211     boost::shared_ptr<Wrapper<HBHERecHitCollection>  const> HBHERecHitsPTR =
00212       getProductByTag<HBHERecHitCollection>(*ep, HBHEPileRecHitInputTag_ );
00213 
00214     if(HBHERecHitsPTR ) {
00215 
00216       const HBHERecHitCollection*  HBHERecHits = const_cast< HBHERecHitCollection * >(HBHERecHitsPTR->product());
00217 
00218       LogDebug("DataMixingEMWorker") << "total # HBHE rechits: " << HBHERecHits->size();
00219 
00220       // loop over digis, adding these to the existing maps                                                     
00221       for(HBHERecHitCollection::const_iterator it  = HBHERecHits->begin();
00222           it != HBHERecHits->end(); ++it) {
00223 
00224         HBHERecHitStorage_.insert(HBHERecHitMap::value_type( (it->id()), *it ));
00225 
00226 #ifdef DEBUG
00227         LogDebug("DataMixingEMWorker") << "processed HBHERecHit with rawId: "
00228                                        << it->id().rawId() << "\n"
00229                                        << " rechit energy: " << it->energy();
00230 #endif
00231       }
00232     }
00233 
00234     // HO Next
00235 
00236     boost::shared_ptr<Wrapper<HORecHitCollection>  const> HORecHitsPTR =
00237       getProductByTag<HORecHitCollection>(*ep, HOPileRecHitInputTag_ );
00238 
00239     if(HORecHitsPTR ) {
00240 
00241       const HORecHitCollection*  HORecHits = const_cast< HORecHitCollection * >(HORecHitsPTR->product());
00242 
00243       LogDebug("DataMixingEMWorker") << "total # HO rechits: " << HORecHits->size();
00244 
00245       // loop over digis, adding these to the existing maps                                                     
00246       for(HORecHitCollection::const_iterator it  = HORecHits->begin();
00247           it != HORecHits->end(); ++it) {
00248 
00249         HORecHitStorage_.insert(HORecHitMap::value_type( (it->id()), *it ));
00250 
00251 #ifdef DEBUG
00252         LogDebug("DataMixingEMWorker") << "processed HORecHit with rawId: "
00253                                        << it->id().rawId() << "\n"
00254                                        << " rechit energy: " << it->energy();
00255 #endif
00256       }
00257     }
00258 
00259     // HF Next
00260 
00261     boost::shared_ptr<Wrapper<HFRecHitCollection>  const> HFRecHitsPTR =
00262       getProductByTag<HFRecHitCollection>(*ep, HFPileRecHitInputTag_ );
00263 
00264     if(HFRecHitsPTR ) {
00265 
00266       const HFRecHitCollection*  HFRecHits = const_cast< HFRecHitCollection * >(HFRecHitsPTR->product());
00267 
00268       LogDebug("DataMixingEMWorker") << "total # HF rechits: " << HFRecHits->size();
00269 
00270       // loop over digis, adding these to the existing maps                                                     
00271       for(HFRecHitCollection::const_iterator it  = HFRecHits->begin();
00272           it != HFRecHits->end(); ++it) {
00273 
00274         HFRecHitStorage_.insert(HFRecHitMap::value_type( (it->id()), *it ));
00275 
00276 #ifdef DEBUG
00277         LogDebug("DataMixingEMWorker") << "processed HFRecHit with rawId: "
00278                                        << it->id().rawId() << "\n"
00279                                        << " rechit energy: " << it->energy();
00280 #endif
00281       }
00282     }
00283 
00284     // ZDC Next
00285 
00286     boost::shared_ptr<Wrapper<ZDCRecHitCollection>  const> ZDCRecHitsPTR =
00287       getProductByTag<ZDCRecHitCollection>(*ep, ZDCPileRecHitInputTag_ );
00288 
00289     if(ZDCRecHitsPTR ) {
00290 
00291       const ZDCRecHitCollection*  ZDCRecHits = const_cast< ZDCRecHitCollection * >(ZDCRecHitsPTR->product());
00292 
00293       LogDebug("DataMixingEMWorker") << "total # ZDC rechits: " << ZDCRecHits->size();
00294 
00295       // loop over digis, adding these to the existing maps                                                     
00296       for(ZDCRecHitCollection::const_iterator it  = ZDCRecHits->begin();
00297           it != ZDCRecHits->end(); ++it) {
00298 
00299         ZDCRecHitStorage_.insert(ZDCRecHitMap::value_type( (it->id()), *it ));
00300 
00301 #ifdef DEBUG
00302         LogDebug("DataMixingEMWorker") << "processed ZDCRecHit with rawId: "
00303                                        << it->id().rawId() << "\n"
00304                                        << " rechit energy: " << it->energy();
00305 #endif
00306       }
00307     }
00308 
00309 
00310   }
00311  
00312   void DataMixingHcalWorker::putHcal(edm::Event &e) {
00313 
00314     // collection of rechits to put in the event
00315     std::auto_ptr< HBHERecHitCollection > HBHErechits( new HBHERecHitCollection );
00316     std::auto_ptr< HORecHitCollection > HOrechits( new HORecHitCollection );
00317     std::auto_ptr< HFRecHitCollection > HFrechits( new HFRecHitCollection );
00318     std::auto_ptr< ZDCRecHitCollection > ZDCrechits( new ZDCRecHitCollection );
00319 
00320     // loop over the maps we have, re-making individual hits or digis if necessary.
00321     DetId formerID = 0;
00322     DetId currentID;
00323     float ESum = 0.;
00324     float HBTime = 0.;
00325 
00326     // HB first...
00327 
00328     HBHERecHitMap::const_iterator iHBchk;
00329 
00330     for(HBHERecHitMap::const_iterator iHB  = HBHERecHitStorage_.begin();
00331         iHB != HBHERecHitStorage_.end(); ++iHB) {
00332 
00333       currentID = iHB->first; 
00334 
00335       if (currentID == formerID) { // we have to add these rechits together
00336 
00337         ESum+=(iHB->second).energy();  
00338 
00339       }
00340       else {
00341         if(formerID>0) {
00342           // cutoff for ESum?                                                                                 
00343           HBHERecHit aHit(formerID, ESum, HBTime);
00344           HBHErechits->push_back( aHit );
00345         }
00346         //save pointers for next iteration                                                                    
00347         formerID = currentID;
00348         ESum = (iHB->second).energy();
00349         HBTime = (iHB->second).time();  // take time of first hit in sequence - is this ok?
00350       }
00351 
00352       iHBchk = iHB;
00353       if((++iHBchk) == HBHERecHitStorage_.end()) {  //make sure not to lose the last one  
00354         HBHERecHit aHit(formerID, ESum, HBTime);
00355         HBHErechits->push_back( aHit );
00356       }
00357     }
00358 
00359     // HO next...
00360 
00361     // loop over the maps we have, re-making individual hits or digis if necessary.
00362     formerID = 0;
00363     ESum = 0.;
00364     float HOTime = 0.;
00365 
00366     HORecHitMap::const_iterator iHOchk;
00367 
00368     for(HORecHitMap::const_iterator iHO  = HORecHitStorage_.begin();
00369         iHO != HORecHitStorage_.end(); ++iHO) {
00370 
00371       currentID = iHO->first; 
00372 
00373       if (currentID == formerID) { // we have to add these rechits together
00374 
00375         ESum+=(iHO->second).energy();  
00376 
00377       }
00378       else {
00379         if(formerID>0) {
00380           // cutoff for ESum?                                                                                 
00381           HORecHit aHit(formerID, ESum, HOTime);
00382           HOrechits->push_back( aHit );
00383         }
00384         //save pointers for next iteration                                                                    
00385         formerID = currentID;
00386         ESum = (iHO->second).energy();
00387         HOTime = (iHO->second).time();  // take time of first hit in sequence - is this ok?
00388       }
00389 
00390       iHOchk = iHO;
00391       if((++iHOchk) == HORecHitStorage_.end()) {  //make sure not to lose the last one  
00392         HORecHit aHit(formerID, ESum, HOTime);
00393         HOrechits->push_back( aHit );
00394       }
00395     }
00396 
00397 
00398     // HF next...
00399 
00400     // loop over the maps we have, re-making individual hits or digis if necessary.
00401     formerID = 0;
00402     ESum = 0.;
00403     float HFTime = 0.;
00404     HFRecHit HFOldHit;
00405 
00406     HFRecHitMap::const_iterator iHFchk;
00407 
00408     for(HFRecHitMap::const_iterator iHF  = HFRecHitStorage_.begin();
00409         iHF != HFRecHitStorage_.end(); ++iHF) {
00410 
00411       currentID = iHF->first; 
00412 
00413       if (currentID == formerID) { // we have to add these rechits together
00414 
00415         ESum+=(iHF->second).energy();  
00416 
00417       }
00418       else {
00419         if(formerID>0) {
00420           // cutoff for ESum?                                                                                 
00421           HFRecHit aHit(formerID, ESum, HFTime);
00422           HFrechits->push_back( aHit );
00423         }
00424         //save pointers for next iteration                                                                    
00425         formerID = currentID;
00426         ESum = (iHF->second).energy();
00427         HFTime = (iHF->second).time();  // take time of first hit in sequence - is this ok?
00428       }
00429 
00430       iHFchk = iHF;
00431       if((++iHFchk) == HFRecHitStorage_.end()) {  //make sure not to lose the last one  
00432         HFRecHit aHit(formerID, ESum, HBTime);
00433         HFrechits->push_back( aHit );
00434       }
00435     }
00436 
00437     // ZDC next...
00438 
00439     // loop over the maps we have, re-making individual hits or digis if necessary.
00440     formerID = 0;
00441     ESum = 0.;
00442     float ZDCTime = 0.;
00443     float lowGainEnergy = 0;
00444     ZDCRecHit ZOldHit;
00445 
00446     ZDCRecHitMap::const_iterator iZDCchk;
00447 
00448     for(ZDCRecHitMap::const_iterator iZDC  = ZDCRecHitStorage_.begin();
00449         iZDC != ZDCRecHitStorage_.end(); ++iZDC) {
00450 
00451       currentID = iZDC->first; 
00452 
00453       if (currentID == formerID) { // we have to add these rechits together
00454 
00455         ESum+=(iZDC->second).energy();  
00456         
00457       }
00458       else {
00459         if(formerID>0) {
00460           // cutoff for ESum?                                                                                 
00461           ZDCRecHit aHit(formerID, ESum, ZDCTime, lowGainEnergy);
00462           ZDCrechits->push_back( aHit );
00463         }
00464         //save pointers for next iteration                                                                    
00465         formerID = currentID;
00466         ESum = (iZDC->second).energy();
00467         lowGainEnergy = (iZDC->second).lowGainEnergy();
00468         ZDCTime = (iZDC->second).time();  // take time of first hit in sequence - is this ok?
00469       }
00470       
00471       iZDCchk = iZDC;
00472       if((++iZDCchk) == ZDCRecHitStorage_.end()) {  //make sure not to lose the last one  
00473         ZDCRecHit aHit(formerID, ESum, HBTime, lowGainEnergy);
00474         ZDCrechits->push_back( aHit );
00475       }
00476     } 
00477   
00478    //done merging
00479 
00480     // put the collection of recunstructed hits in the event   
00481     LogInfo("DataMixingHcalWorker") << "total # HBHE Merged rechits: " << HBHErechits->size() ;
00482     LogInfo("DataMixingHcalWorker") << "total # HO Merged rechits: " << HOrechits->size() ;
00483     LogInfo("DataMixingHcalWorker") << "total # HF Merged rechits: " << HFrechits->size() ;
00484     LogInfo("DataMixingHcalWorker") << "total # ZDC Merged rechits: " << ZDCrechits->size() ;
00485 
00486     e.put( HBHErechits, HBHERecHitCollectionDM_ );
00487     e.put( HOrechits, HORecHitCollectionDM_ );
00488     e.put( HFrechits, HFRecHitCollectionDM_ );
00489     e.put( ZDCrechits, ZDCRecHitCollectionDM_ );
00490 
00491     // clear local storage after this event
00492     HBHERecHitStorage_.clear();
00493     HORecHitStorage_.clear();
00494     HFRecHitStorage_.clear();
00495     ZDCRecHitStorage_.clear();
00496 
00497   }
00498 
00499 } //edm