Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:24:22

0001 #include "RecoLocalCalo/EcalRecProducers/plugins/EcalRecHitWorkerRecover.h"
0002 
0003 #include "FWCore/Framework/interface/EventSetup.h"
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "CondFormats/DataRecord/interface/EcalIntercalibConstantsRcd.h"
0006 #include "CondFormats/DataRecord/interface/EcalTimeCalibConstantsRcd.h"
0007 #include "CondFormats/DataRecord/interface/EcalADCToGeVConstantRcd.h"
0008 
0009 #include "DataFormats/EcalDigi/interface/EcalDigiCollections.h"
0010 #include "DataFormats/EcalDetId/interface/EcalScDetId.h"
0011 
0012 #include "CondFormats/EcalObjects/interface/EcalTimeCalibConstants.h"
0013 
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015 #include "FWCore/Framework/interface/EDProducer.h"
0016 
0017 EcalRecHitWorkerRecover::EcalRecHitWorkerRecover(const edm::ParameterSet& ps, edm::ConsumesCollector& c)
0018     : EcalRecHitWorkerBaseClass(ps, c), ecalScaleTokens_(c), tpgscaleTokens_(c) {
0019   rechitMaker_ = std::make_unique<EcalRecHitSimpleAlgo>();
0020   // isolated channel recovery
0021   singleRecoveryMethod_ = ps.getParameter<std::string>("singleChannelRecoveryMethod");
0022   singleRecoveryThreshold_ = ps.getParameter<double>("singleChannelRecoveryThreshold");
0023   sum8RecoveryThreshold_ = ps.getParameter<double>("sum8ChannelRecoveryThreshold");
0024   killDeadChannels_ = ps.getParameter<bool>("killDeadChannels");
0025   recoverEBIsolatedChannels_ = ps.getParameter<bool>("recoverEBIsolatedChannels");
0026   recoverEEIsolatedChannels_ = ps.getParameter<bool>("recoverEEIsolatedChannels");
0027   recoverEBVFE_ = ps.getParameter<bool>("recoverEBVFE");
0028   recoverEEVFE_ = ps.getParameter<bool>("recoverEEVFE");
0029   recoverEBFE_ = ps.getParameter<bool>("recoverEBFE");
0030   recoverEEFE_ = ps.getParameter<bool>("recoverEEFE");
0031   laserToken_ = c.esConsumes<EcalLaserDbService, EcalLaserDbRecord>();
0032   caloTopologyToken_ = c.esConsumes<CaloTopology, CaloTopologyRecord>();
0033   pEcalMappingToken_ = c.esConsumes<EcalElectronicsMapping, EcalMappingRcd>();
0034   pEBGeomToken_ = c.esConsumes<CaloSubdetectorGeometry, EcalBarrelGeometryRecord>(edm::ESInputTag("", "EcalBarrel"));
0035   caloGeometryToken_ = c.esConsumes<CaloGeometry, CaloGeometryRecord>();
0036   chStatusToken_ = c.esConsumes<EcalChannelStatus, EcalChannelStatusRcd>();
0037   ttMapToken_ = c.esConsumes<EcalTrigTowerConstituentsMap, IdealGeometryRecord>();
0038 
0039   dbStatusToBeExcludedEE_ = ps.getParameter<std::vector<int> >("dbStatusToBeExcludedEE");
0040   dbStatusToBeExcludedEB_ = ps.getParameter<std::vector<int> >("dbStatusToBeExcludedEB");
0041 
0042   logWarningEtThreshold_EB_FE_ = ps.getParameter<double>("logWarningEtThreshold_EB_FE");
0043   logWarningEtThreshold_EE_FE_ = ps.getParameter<double>("logWarningEtThreshold_EE_FE");
0044 
0045   tpDigiToken_ =
0046       c.consumes<EcalTrigPrimDigiCollection>(ps.getParameter<edm::InputTag>("triggerPrimitiveDigiCollection"));
0047 
0048   if (recoverEBIsolatedChannels_ && singleRecoveryMethod_ == "BDTG")
0049     ebDeadChannelCorrector.setParameters(ps);
0050 }
0051 
0052 void EcalRecHitWorkerRecover::set(const edm::EventSetup& es) {
0053   laser = es.getHandle(laserToken_);
0054   caloTopology_ = es.getHandle(caloTopologyToken_);
0055   pEcalMapping_ = es.getHandle(pEcalMappingToken_);
0056   ecalMapping_ = pEcalMapping_.product();
0057   // geometry...
0058   pEBGeom_ = es.getHandle(pEBGeomToken_);
0059   caloGeometry_ = es.getHandle(caloGeometryToken_);
0060   chStatus_ = es.getHandle(chStatusToken_);
0061   geo_ = caloGeometry_.product();
0062   ebGeom_ = pEBGeom_.product();
0063   ttMap_ = es.getHandle(ttMapToken_);
0064   recoveredDetIds_EB_.clear();
0065   recoveredDetIds_EE_.clear();
0066   eventSetup_ = &es;
0067 }
0068 
0069 bool EcalRecHitWorkerRecover::run(const edm::Event& evt,
0070                                   const EcalUncalibratedRecHit& uncalibRH,
0071                                   EcalRecHitCollection& result) {
0072   DetId detId = uncalibRH.id();
0073   uint32_t flags = (0xF & uncalibRH.flags());
0074 
0075   // get laser coefficient
0076   //float lasercalib = laser->getLaserCorrection( detId, evt.time());
0077 
0078   // killDeadChannels_ = true, means explicitely kill dead channels even if the recovered energies are computed in the code
0079   // if you don't want to store the recovered energies in the rechit you can produce LogWarnings if logWarningEtThreshold_EB(EE)_FE>0
0080   // logWarningEtThreshold_EB(EE)_FE_<0 will not compute the recovered energies at all (faster)
0081 
0082   if (killDeadChannels_) {
0083     if ((flags == EcalRecHitWorkerRecover::EB_single && !recoverEBIsolatedChannels_) ||
0084         (flags == EcalRecHitWorkerRecover::EE_single && !recoverEEIsolatedChannels_) ||
0085         (flags == EcalRecHitWorkerRecover::EB_VFE && !recoverEBVFE_) ||
0086         (flags == EcalRecHitWorkerRecover::EE_VFE && !recoverEEVFE_)) {
0087       EcalRecHit hit(detId, 0., 0., EcalRecHit::kDead);
0088       hit.setFlag(EcalRecHit::kDead);
0089       insertRecHit(hit, result);  // insert trivial rechit with kDead flag
0090       return true;
0091     }
0092     if (flags == EcalRecHitWorkerRecover::EB_FE && !recoverEBFE_) {
0093       EcalTrigTowerDetId ttDetId(((EBDetId)detId).tower());
0094       std::vector<DetId> vid = ttMap_->constituentsOf(ttDetId);
0095       for (std::vector<DetId>::const_iterator dit = vid.begin(); dit != vid.end(); ++dit) {
0096         EcalRecHit hit((*dit), 0., 0., EcalRecHit::kDead);
0097         hit.setFlag(EcalRecHit::kDead);
0098         insertRecHit(hit, result);  // insert trivial rechit with kDead flag
0099       }
0100       if (logWarningEtThreshold_EB_FE_ < 0)
0101         return true;  // if you don't want log warning just return true
0102     }
0103     if (flags == EcalRecHitWorkerRecover::EE_FE && !recoverEEFE_) {
0104       EEDetId id(detId);
0105       EcalScDetId sc(1 + (id.ix() - 1) / 5, 1 + (id.iy() - 1) / 5, id.zside());
0106       std::vector<DetId> eeC;
0107       for (int dx = 1; dx <= 5; ++dx) {
0108         for (int dy = 1; dy <= 5; ++dy) {
0109           int ix = (sc.ix() - 1) * 5 + dx;
0110           int iy = (sc.iy() - 1) * 5 + dy;
0111           int iz = sc.zside();
0112           if (EEDetId::validDetId(ix, iy, iz)) {
0113             eeC.push_back(EEDetId(ix, iy, iz));
0114           }
0115         }
0116       }
0117       for (size_t i = 0; i < eeC.size(); ++i) {
0118         EcalRecHit hit(eeC[i], 0., 0., EcalRecHit::kDead);
0119         hit.setFlag(EcalRecHit::kDead);
0120         insertRecHit(hit, result);  // insert trivial rechit with kDead flag
0121       }
0122       if (logWarningEtThreshold_EE_FE_ < 0)
0123         return true;  // if you don't want log warning just return true
0124     }
0125   }
0126 
0127   if (flags == EcalRecHitWorkerRecover::EB_single) {
0128     // recover as single dead channel
0129     ebDeadChannelCorrector.setCaloTopology(caloTopology_.product());
0130 
0131     // channel recovery. Accepted new RecHit has the flag AcceptRecHit=TRUE
0132     bool AcceptRecHit = true;
0133     float ebEn = ebDeadChannelCorrector.correct(
0134         detId, result, singleRecoveryMethod_, singleRecoveryThreshold_, sum8RecoveryThreshold_, &AcceptRecHit);
0135     EcalRecHit hit(detId, ebEn, 0., EcalRecHit::kDead);
0136 
0137     if (hit.energy() != 0 and AcceptRecHit == true) {
0138       hit.setFlag(EcalRecHit::kNeighboursRecovered);
0139     } else {
0140       // recovery failed
0141       hit.setFlag(EcalRecHit::kDead);
0142     }
0143     insertRecHit(hit, result);
0144 
0145   } else if (flags == EcalRecHitWorkerRecover::EE_single) {
0146     // recover as single dead channel
0147     eeDeadChannelCorrector.setCaloTopology(caloTopology_.product());
0148 
0149     // channel recovery. Accepted new RecHit has the flag AcceptRecHit=TRUE
0150     bool AcceptRecHit = true;
0151     float eeEn = eeDeadChannelCorrector.correct(
0152         detId, result, singleRecoveryMethod_, singleRecoveryThreshold_, sum8RecoveryThreshold_, &AcceptRecHit);
0153     EcalRecHit hit(detId, eeEn, 0., EcalRecHit::kDead);
0154     if (hit.energy() != 0 and AcceptRecHit == true) {
0155       hit.setFlag(EcalRecHit::kNeighboursRecovered);
0156     } else {
0157       // recovery failed
0158       hit.setFlag(EcalRecHit::kDead);
0159     }
0160     insertRecHit(hit, result);
0161 
0162   } else if (flags == EcalRecHitWorkerRecover::EB_VFE) {
0163     // recover as dead VFE
0164     EcalRecHit hit(detId, 0., 0.);
0165     hit.setFlag(EcalRecHit::kDead);
0166     // recovery not implemented
0167     insertRecHit(hit, result);
0168   } else if (flags == EcalRecHitWorkerRecover::EB_FE) {
0169     // recover as dead TT
0170 
0171     EcalTrigTowerDetId ttDetId(((EBDetId)detId).tower());
0172     edm::Handle<EcalTrigPrimDigiCollection> pTPDigis;
0173     evt.getByToken(tpDigiToken_, pTPDigis);
0174     const EcalTrigPrimDigiCollection* tpDigis = nullptr;
0175     tpDigis = pTPDigis.product();
0176 
0177     EcalTrigPrimDigiCollection::const_iterator tp = tpDigis->find(ttDetId);
0178     // recover the whole trigger tower
0179     if (tp != tpDigis->end()) {
0180       EcalTPGScale ecalScale{ecalScaleTokens_, *eventSetup_};
0181       //std::vector<DetId> vid = ecalMapping_->dccTowerConstituents( ecalMapping_->DCCid( ttDetId ), ecalMapping_->iTT( ttDetId ) );
0182       std::vector<DetId> vid = ttMap_->constituentsOf(ttDetId);
0183       float tpEt = ecalScale.getTPGInGeV(tp->compressedEt(), tp->id());
0184       float tpEtThreshEB = logWarningEtThreshold_EB_FE_;
0185       if (tpEt > tpEtThreshEB) {
0186         edm::LogWarning("EnergyInDeadEB_FE") << "TP energy in the dead TT = " << tpEt << " at " << ttDetId;
0187       }
0188       if (!killDeadChannels_ || recoverEBFE_) {
0189         // democratic energy sharing
0190 
0191         for (std::vector<DetId>::const_iterator dit = vid.begin(); dit != vid.end(); ++dit) {
0192           if (alreadyInserted(*dit))
0193             continue;
0194           float theta = ebGeom_->getGeometry(*dit)->getPosition().theta();
0195           float tpEt = ecalScale.getTPGInGeV(tp->compressedEt(), tp->id());
0196           if (checkChannelStatus(*dit, dbStatusToBeExcludedEB_)) {
0197             EcalRecHit hit(*dit, tpEt / ((float)vid.size()) / sin(theta), 0.);
0198             hit.setFlag(EcalRecHit::kTowerRecovered);
0199             if (tp->compressedEt() == 0xFF)
0200               hit.setFlag(EcalRecHit::kTPSaturated);
0201             if (tp->sFGVB())
0202               hit.setFlag(EcalRecHit::kL1SpikeFlag);
0203             insertRecHit(hit, result);
0204           }
0205         }
0206       } else {
0207         // tp not found => recovery failed
0208         std::vector<DetId> vid = ttMap_->constituentsOf(ttDetId);
0209         for (std::vector<DetId>::const_iterator dit = vid.begin(); dit != vid.end(); ++dit) {
0210           if (alreadyInserted(*dit))
0211             continue;
0212           EcalRecHit hit(*dit, 0., 0.);
0213           hit.setFlag(EcalRecHit::kDead);
0214           insertRecHit(hit, result);
0215         }
0216       }
0217     }
0218   } else if (flags == EcalRecHitWorkerRecover::EE_FE) {
0219     // Structure for recovery:
0220     // ** SC --> EEDetId constituents (eeC) --> associated Trigger Towers (aTT) --> EEDetId constituents (aTTC)
0221     // ** energy for a SC EEDetId = [ sum_aTT(energy) - sum_aTTC(energy) ] / N_eeC
0222     // .. i.e. the total energy of the TTs covering the SC minus
0223     // .. the energy of the recHits in the TTs but not in the SC
0224     //std::vector<DetId> vid = ecalMapping_->dccTowerConstituents( ecalMapping_->DCCid( ttDetId ), ecalMapping_->iTT( ttDetId ) );
0225     // due to lack of implementation of the EcalTrigTowerDetId ix,iy methods in EE we compute Et recovered energies (in EB we compute E)
0226 
0227     EEDetId eeId(detId);
0228     EcalScDetId sc((eeId.ix() - 1) / 5 + 1, (eeId.iy() - 1) / 5 + 1, eeId.zside());
0229     std::set<DetId> eeC;
0230     for (int dx = 1; dx <= 5; ++dx) {
0231       for (int dy = 1; dy <= 5; ++dy) {
0232         int ix = (sc.ix() - 1) * 5 + dx;
0233         int iy = (sc.iy() - 1) * 5 + dy;
0234         int iz = sc.zside();
0235         if (EEDetId::validDetId(ix, iy, iz)) {
0236           EEDetId id(ix, iy, iz);
0237           if (checkChannelStatus(id, dbStatusToBeExcludedEE_)) {
0238             eeC.insert(id);
0239           }  // check status
0240         }
0241       }
0242     }
0243 
0244     edm::Handle<EcalTrigPrimDigiCollection> pTPDigis;
0245     evt.getByToken(tpDigiToken_, pTPDigis);
0246     const EcalTrigPrimDigiCollection* tpDigis = nullptr;
0247     tpDigis = pTPDigis.product();
0248 
0249     // associated trigger towers
0250     std::set<EcalTrigTowerDetId> aTT;
0251     for (std::set<DetId>::const_iterator it = eeC.begin(); it != eeC.end(); ++it) {
0252       aTT.insert(ttMap_->towerOf(*it));
0253     }
0254 
0255     EcalTPGScale tpgscale(tpgscaleTokens_, *eventSetup_);
0256     EcalTPGScale ecalScale(ecalScaleTokens_, *eventSetup_);
0257     // associated trigger towers: total energy
0258     float totE = 0;
0259     // associated trigger towers: EEDetId constituents
0260     std::set<DetId> aTTC;
0261     bool atLeastOneTPSaturated = false;
0262     for (std::set<EcalTrigTowerDetId>::const_iterator it = aTT.begin(); it != aTT.end(); ++it) {
0263       // add the energy of this trigger tower
0264       EcalTrigPrimDigiCollection::const_iterator itTP = tpDigis->find(*it);
0265       if (itTP != tpDigis->end()) {
0266         std::vector<DetId> v = ttMap_->constituentsOf(*it);
0267 
0268         // from the constituents, remove dead channels
0269         std::vector<DetId>::iterator ttcons = v.begin();
0270         while (ttcons != v.end()) {
0271           if (!checkChannelStatus(*ttcons, dbStatusToBeExcludedEE_)) {
0272             ttcons = v.erase(ttcons);
0273           } else {
0274             ++ttcons;
0275           }
0276         }  // while
0277 
0278         if (itTP->compressedEt() == 0xFF) {  // In the case of a saturated trigger tower, a fraction
0279           atLeastOneTPSaturated =
0280               true;  //of the saturated energy is put in: number of xtals in dead region/total xtals in TT *63.75
0281 
0282           //Alternative recovery algorithm that I will now investigate.
0283           //Estimate energy sums the energy in the working channels, then decides how much energy
0284           //to put here depending on that. Duncan 20101203
0285 
0286           totE += estimateEnergy(itTP->id().ietaAbs(), &result, eeC, v, tpgscale);
0287 
0288           /* 
0289                          These commented out lines use
0290                          64GeV*fraction of the TT overlapping the dead FE
0291                         
0292                       int count = 0;
0293                       for (std::vector<DetId>::const_iterator idsit = v.begin(); idsit != v.end(); ++ idsit){
0294                       std::set<DetId>::const_iterator itFind = eeC.find(*idsit);
0295                       if (itFind != eeC.end())
0296                       ++count;
0297                       }
0298                       //std::cout << count << ", " << v.size() << std::endl;
0299                       totE+=((float)count/(float)v.size())* ((it->ietaAbs()>26)?2*ecalScale_.getTPGInGeV( itTP->compressedEt(), itTP->id() ):ecalScale_.getTPGInGeV( itTP->compressedEt(), itTP->id() ));*/
0300         } else {
0301           totE += ((it->ietaAbs() > 26) ? 2 : 1) * ecalScale.getTPGInGeV(itTP->compressedEt(), itTP->id());
0302         }
0303 
0304         // get the trigger tower constituents
0305 
0306         if (itTP->compressedEt() == 0) {  // If there's no energy in TT, the constituents are removed from the recovery.
0307           for (size_t i = 0; i < v.size(); ++i)
0308             eeC.erase(v[i]);
0309         } else if (itTP->compressedEt() != 0xFF) {
0310           //If it's saturated the energy has already been determined, so we do not want to subtract any channels
0311           for (size_t j = 0; j < v.size(); ++j) {
0312             aTTC.insert(v[j]);
0313           }
0314         }
0315       }
0316     }
0317     // remove crystals of dead SC
0318     // (this step is not needed if sure that SC crystals are not
0319     // in the recHit collection)
0320 
0321     for (std::set<DetId>::const_iterator it = eeC.begin(); it != eeC.end(); ++it) {
0322       aTTC.erase(*it);
0323     }
0324     // compute the total energy for the dead SC
0325     const EcalRecHitCollection* hits = &result;
0326     for (std::set<DetId>::const_iterator it = aTTC.begin(); it != aTTC.end(); ++it) {
0327       EcalRecHitCollection::const_iterator jt = hits->find(*it);
0328       if (jt != hits->end()) {
0329         float energy = jt->energy();  // Correct conversion to Et
0330         float eta = geo_->getPosition(jt->id()).eta();
0331         float pf = 1.0 / cosh(eta);
0332         // use Et instead of E, consistent with the Et estimation of the associated TT
0333         totE -= energy * pf;
0334       }
0335     }
0336 
0337     float scEt = totE;
0338     float scEtThreshEE = logWarningEtThreshold_EE_FE_;
0339     if (scEt > scEtThreshEE) {
0340       edm::LogWarning("EnergyInDeadEE_FE") << "TP energy in the dead TT = " << scEt << " at " << sc;
0341     }
0342 
0343     // assign the energy to the SC crystals
0344     if (!killDeadChannels_ || recoverEEFE_) {  // if eeC is empty, i.e. there are no hits
0345                                                // in the tower, nothing is returned. No negative values from noise.
0346       for (std::set<DetId>::const_iterator it = eeC.begin(); it != eeC.end(); ++it) {
0347         float eta = geo_->getPosition(*it).eta();  //Convert back to E from Et for the recovered hits
0348         float pf = 1.0 / cosh(eta);
0349         EcalRecHit hit(*it, totE / ((float)eeC.size() * pf), 0);
0350 
0351         if (atLeastOneTPSaturated)
0352           hit.setFlag(EcalRecHit::kTPSaturated);
0353         hit.setFlag(EcalRecHit::kTowerRecovered);
0354         insertRecHit(hit, result);
0355 
0356       }  // for
0357     }    // if
0358   }
0359   return true;
0360 }
0361 
0362 float EcalRecHitWorkerRecover::estimateEnergy(int ieta,
0363                                               EcalRecHitCollection* hits,
0364                                               const std::set<DetId>& sId,
0365                                               const std::vector<DetId>& vId,
0366                                               const EcalTPGScale& tpgscale) {
0367   float xtalE = 0;
0368   int count = 0;
0369   for (std::vector<DetId>::const_iterator vIdit = vId.begin(); vIdit != vId.end(); ++vIdit) {
0370     std::set<DetId>::const_iterator sIdit = sId.find(*vIdit);
0371     if (sIdit == sId.end()) {
0372       float energy = hits->find(*vIdit)->energy();
0373       float eta = geo_->getPosition(*vIdit).eta();
0374       float pf = 1.0 / cosh(eta);
0375       xtalE += energy * pf;
0376       count++;
0377     }
0378   }
0379 
0380   if (count == 0) {  // If there are no overlapping crystals return saturated value.
0381 
0382     double etsat = tpgscale.getTPGInGeV(0xFF,
0383                                         ttMap_->towerOf(*vId.begin()));  // get saturation value for the first
0384                                                                          // constituent, for the others it's the same
0385 
0386     return etsat / cosh(ieta) * (ieta > 26 ? 2 : 1);  // account for duplicated TT in EE for ieta>26
0387   } else
0388     return xtalE * ((vId.size() / (float)count) - 1) * (ieta > 26 ? 2 : 1);
0389 }
0390 
0391 void EcalRecHitWorkerRecover::insertRecHit(const EcalRecHit& hit, EcalRecHitCollection& collection) {
0392   // skip already inserted DetId's and raise a log warning
0393   if (alreadyInserted(hit.id())) {
0394     edm::LogWarning("EcalRecHitWorkerRecover") << "DetId already recovered! Skipping...";
0395     return;
0396   }
0397   EcalRecHitCollection::iterator it = collection.find(hit.id());
0398   if (it == collection.end()) {
0399     // insert the hit in the collection
0400     collection.push_back(hit);
0401   } else {
0402     // overwrite existing recHit
0403     *it = hit;
0404   }
0405   if (hit.id().subdetId() == EcalBarrel) {
0406     recoveredDetIds_EB_.insert(hit.id());
0407   } else if (hit.id().subdetId() == EcalEndcap) {
0408     recoveredDetIds_EE_.insert(hit.id());
0409   } else {
0410     edm::LogError("EcalRecHitWorkerRecover::InvalidDetId") << "Invalid DetId " << hit.id().rawId();
0411   }
0412 }
0413 
0414 bool EcalRecHitWorkerRecover::alreadyInserted(const DetId& id) {
0415   bool res = false;
0416   if (id.subdetId() == EcalBarrel) {
0417     res = (recoveredDetIds_EB_.find(id) != recoveredDetIds_EB_.end());
0418   } else if (id.subdetId() == EcalEndcap) {
0419     res = (recoveredDetIds_EE_.find(id) != recoveredDetIds_EE_.end());
0420   } else {
0421     edm::LogError("EcalRecHitWorkerRecover::InvalidDetId") << "Invalid DetId " << id.rawId();
0422   }
0423   return res;
0424 }
0425 
0426 // In the future, this will be used to calibrate the TT energy. There is a dependance on
0427 // eta at lower energies that can be corrected for here after more validation.
0428 float EcalRecHitWorkerRecover::recCheckCalib(float eTT, int ieta) { return eTT; }
0429 
0430 // return false is the channel has status  in the list of statusestoexclude
0431 // true otherwise (channel ok)
0432 // Careful: this function works on raw (encoded) channel statuses
0433 bool EcalRecHitWorkerRecover::checkChannelStatus(const DetId& id, const std::vector<int>& statusestoexclude) {
0434   if (!chStatus_.isValid())
0435     edm::LogError("ObjectNotFound") << "Channel Status not set";
0436 
0437   EcalChannelStatus::const_iterator chIt = chStatus_->find(id);
0438   uint16_t dbStatus = 0;
0439   if (chIt != chStatus_->end()) {
0440     dbStatus = chIt->getEncodedStatusCode();
0441   } else {
0442     edm::LogError("ObjectNotFound") << "No channel status found for xtal " << id.rawId()
0443                                     << "! something wrong with EcalChannelStatus in your DB? ";
0444   }
0445 
0446   for (std::vector<int>::const_iterator status = statusestoexclude.begin(); status != statusestoexclude.end();
0447        ++status) {
0448     if (*status == dbStatus)
0449       return false;
0450   }
0451 
0452   return true;
0453 }
0454 
0455 #include "FWCore/Framework/interface/MakerMacros.h"
0456 #include "RecoLocalCalo/EcalRecProducers/interface/EcalRecHitWorkerFactory.h"
0457 DEFINE_EDM_PLUGIN(EcalRecHitWorkerFactory, EcalRecHitWorkerRecover, "EcalRecHitWorkerRecover");