82 virtual void endJob()
override;
125 int getChannelStatusMaps();
137 int setEvtTPstatus(
const double &tpCntCut,
const int &chnStatus);
175 int setEvtRecHitstatus(
const double &tpValCut,
const int &chnStatus,
const int &towerTest);
183 : taggingMode_ (iConfig.getParameter<bool>(
"taggingMode") )
184 , debug_ (iConfig.getParameter<bool>(
"debug") )
185 , verbose_ (iConfig.getParameter<int>(
"verbose") )
186 , doEEfilter_ (iConfig.getUntrackedParameter<bool>(
"doEEfilter") )
187 , ebReducedRecHitCollection_ (iConfig.getParameter<edm::
InputTag>(
"ebReducedRecHitCollection") )
188 , eeReducedRecHitCollection_ (iConfig.getParameter<edm::
InputTag>(
"eeReducedRecHitCollection") )
191 , maskedEcalChannelStatusThreshold_ (iConfig.getParameter<int>(
"maskedEcalChannelStatusThreshold") )
192 , etValToBeFlagged_ (iConfig.getParameter<double>(
"etValToBeFlagged") )
193 , tpDigiCollection_ (iConfig.getParameter<edm::
InputTag>(
"tpDigiCollection") )
195 , makeProfileRoot_ (iConfig.getUntrackedParameter<bool>(
"makeProfileRoot") )
196 , profileRootName_ (iConfig.getUntrackedParameter<std::
string>(
"profileRootName") )
205 profTree =
new TTree(
"filter",
"filter profile");
206 profTree->Branch(
"run", &
run,
"run/I");
207 profTree->Branch(
"event", &
event,
"event/I");
208 profTree->Branch(
"lumi", &
ls,
"lumi/I");
219 if( makeProfileRoot_ ){
231 std::vector<edm::Provenance const*> provenances;
233 const unsigned int nProvenance = provenances.size();
234 for (
unsigned int ip = 0; ip < nProvenance; ip++) {
236 if( provenance.
moduleLabel().data() == tpDigiCollection_.label() ){ hastpDigiCollection_ = 1; }
237 if( provenance.
moduleLabel().data() == ebReducedRecHitCollection_.label() || provenance.
moduleLabel().data() == eeReducedRecHitCollection_.label() ){
238 hasReducedRecHits_++;
240 if( hastpDigiCollection_ && hasReducedRecHits_>=2 ){
break; }
243 if( debug_ )
std::cout<<
"\nhastpDigiCollection_ : "<<hastpDigiCollection_<<
" hasReducedRecHits_ : "<<hasReducedRecHits_<<std::endl;
246 const unsigned int nHist = history.
size();
248 releaseVersion_ = history[nHist-2].releaseVersion();
249 TString tmpTstr(releaseVersion_);
250 TObjArray *
split = tmpTstr.Tokenize(
"_");
251 int majorV = TString(split->At(1)->GetName()).Atoi();
252 int minorV = TString(split->At(2)->GetName()).Atoi();
254 if( debug_ )
std::cout<<
"processName : "<<history[nHist-2].processName().
data()<<
" releaseVersion : "<<releaseVersion_<<std::endl;
261 if( !hastpDigiCollection_ && !hasReducedRecHits_ ){ useTPmethod_ =
false; useHITmethod_ =
false;
263 std::cout<<
"\nWARNING ... Cannot find either tpDigiCollection_ or reducedRecHitCollecion_ ?!"<<std::endl;
264 std::cout<<
" Will NOT DO ANY FILTERING !"<<std::endl;
267 else if( hastpDigiCollection_ ){ useTPmethod_ =
true; useHITmethod_ =
false; }
269 else if( majorV >=5 || (majorV==4 && minorV >=2) ){ useTPmethod_ =
false; useHITmethod_ =
true; }
270 else{ useTPmethod_ =
false; useHITmethod_ =
false;
272 std::cout<<
"\nWARNING ... TP filter can ONLY be used in AOD after 42X"<<std::endl;
273 std::cout<<
" Will NOT DO ANY FILTERING !"<<std::endl;
277 if( debug_ )
std::cout<<
"useTPmethod_ : "<<useTPmethod_<<
" useHITmethod_ : "<<useHITmethod_<<std::endl;
279 getEventInfoForFilterOnce_ =
true;
293 iEvent.
getByToken(tpDigiCollectionToken_, pTPDigis);
294 if ( !pTPDigis.isValid() ) {
edm::LogWarning(
"EcalDeadCellTriggerPrimitiveFilter") <<
"Can't get the product " << tpDigiCollection_.instance()
295 <<
" with label " << tpDigiCollection_.label();
return; }
300 iEvent.
getByToken(ebReducedRecHitCollectionToken_,barrelReducedRecHitsHandle);
301 iEvent.
getByToken(eeReducedRecHitCollectionToken_,endcapReducedRecHitsHandle);
311 if (debug_ && verbose_ >=2)
std::cout <<
"***envSet***" << std::endl;
313 ecalScale_.setEventSetup( iSetup );
319 if( !ecalStatus.isValid() )
throw "Failed to get ECAL channel status!";
327 std::vector<int> cutFlowFlagTmpVec; std::vector<std::string> cutFlowStrTmpVec;
329 loadEventInfo(iEvent, iSetup);
331 if( !getEventInfoForFilterOnce_ ){ loadEventInfoForFilter(iEvent); }
340 loadEcalDigis(iEvent, iSetup);
341 evtTagged = setEvtTPstatus(etValToBeFlagged_, 13);
345 loadEcalRecHits(iEvent, iSetup);
346 evtTagged = setEvtRecHitstatus(etValToBeFlagged_, 13, 13);
349 if( evtTagged ){ pass =
false; totFilteredCnt++; }
351 if( makeProfileRoot_ ){
353 cutFlowFlagTmpVec.push_back(evtTagged); cutFlowStrTmpVec.push_back(
"TP");
355 cutFlowFlagTmpPtr = &cutFlowFlagTmpVec;
356 cutFlowStrTmpPtr = &cutFlowStrTmpVec;
361 if(debug_ && verbose_ >=2){
362 int evtstatusABS =
abs(evtTagged);
363 printf(
"\nrun : %8u event : %10llu lumi : %4u evtTPstatus ABS : %d 13 : % 2d\n",
run,
event,
ls, evtstatusABS, evtTagged);
366 std::auto_ptr<bool> pOut(
new bool(pass) );
369 if (taggingMode_)
return true;
388 getChannelStatusMaps();
389 if( debug_ && verbose_ >=2)
std::cout<<
"EcalAllDeadChannelsValMap.size() : "<<EcalAllDeadChannelsValMap.size()<<
" EcalAllDeadChannelsBitMap.size() : "<<EcalAllDeadChannelsBitMap.size()<<std::endl;
395 if( debug_ && verbose_ >=2)
std::cout<<
"***begin setEvtTPstatusRecHits***"<<std::endl;
397 accuTTetMap.clear(); accuTTchnMap.clear(); TTzsideMap.clear();
398 accuSCetMap.clear(); accuSCchnMap.clear(); SCzsideMap.clear();
399 avoidDuplicateVec.clear();
411 for (ebrechit = HitecalEB.
begin(); ebrechit != HitecalEB.
end(); ebrechit++) {
415 std::map<DetId, vector<double> >::iterator valItor = EcalAllDeadChannelsValMap.find(det);
416 if( valItor == EcalAllDeadChannelsValMap.end() )
continue;
418 double theta = valItor->second.back();
420 std::map<DetId, vector<int> >::iterator bitItor = EcalAllDeadChannelsBitMap.find(det);
421 if( bitItor == EcalAllDeadChannelsBitMap.end() )
continue;
423 std::map<DetId, EcalTrigTowerDetId>::iterator ttItor = EcalAllDeadChannelsTTMap.find(det);
424 if( ttItor == EcalAllDeadChannelsTTMap.end() )
continue;
426 int status = bitItor->second.back();
429 if( chnStatus >0 && status == chnStatus ) toDo =
true;
430 if( chnStatus <0 && status >=
abs(chnStatus) ) toDo =
true;
433 if( !ebrechit->isRecovered() ) toDo =
false;
439 int ttzside = ttDetId.
zside();
441 std::vector<DetId> vid = ttMap_->constituentsOf(ttDetId);
443 for(std::vector<DetId>::const_iterator dit = vid.begin(); dit != vid.end(); ++dit ) {
444 std::map<DetId, vector<int> >::iterator bit2Itor = EcalAllDeadChannelsBitMap.find( (*dit) );
445 if( bit2Itor == EcalAllDeadChannelsBitMap.end() ){ towerTestCnt ++;
continue; }
446 if( towerTest >0 && bit2Itor->second.back() == towerTest )
continue;
447 if( towerTest <0 && bit2Itor->
second.back() >=
abs(towerTest) )
continue;
450 if( towerTestCnt !=0 && debug_ && verbose_ >=2)
std::cout<<
"towerTestCnt : "<<towerTestCnt<<
" for towerTest : "<<towerTest<<std::endl;
452 std::vector<DetId>::iterator avoidItor; avoidItor =
find( avoidDuplicateVec.begin(), avoidDuplicateVec.end(), det);
453 if( avoidItor == avoidDuplicateVec.end() ){
454 avoidDuplicateVec.push_back(det);
459 std::map<EcalTrigTowerDetId, double>::iterator ttetItor = accuTTetMap.find(ttDetId);
460 if( ttetItor == accuTTetMap.end() ){
461 accuTTetMap[ttDetId] = ebrechit->energy()*
sin(theta);
462 accuTTchnMap[ttDetId] = 1;
463 TTzsideMap[ttDetId] = ttzside;
465 accuTTetMap[ttDetId] += ebrechit->energy()*
sin(theta);
466 accuTTchnMap[ttDetId] ++;
472 for (eerechit = HitecalEE.
begin(); eerechit != HitecalEE.
end(); eerechit++) {
476 std::map<DetId, vector<double> >::iterator valItor = EcalAllDeadChannelsValMap.find(det);
477 if( valItor == EcalAllDeadChannelsValMap.end() )
continue;
479 double theta = valItor->second.back();
481 std::map<DetId, vector<int> >::iterator bitItor = EcalAllDeadChannelsBitMap.find(det);
482 if( bitItor == EcalAllDeadChannelsBitMap.end() )
continue;
484 std::map<DetId, EcalTrigTowerDetId>::iterator ttItor = EcalAllDeadChannelsTTMap.find(det);
485 if( ttItor == EcalAllDeadChannelsTTMap.end() )
continue;
487 int status = bitItor->second.back();
490 if( chnStatus >0 && status == chnStatus ) toDo =
true;
491 if( chnStatus <0 && status >=
abs(chnStatus) ) toDo =
true;
494 if( !eerechit->isRecovered() ) toDo =
false;
503 std::vector<DetId> vid = ttMap_->constituentsOf(ttDetId);
505 for(std::vector<DetId>::const_iterator dit = vid.begin(); dit != vid.end(); ++dit ) {
506 std::map<DetId, vector<int> >::iterator bit2Itor = EcalAllDeadChannelsBitMap.find( (*dit) );
507 if( bit2Itor == EcalAllDeadChannelsBitMap.end() ){ towerTestCnt ++;
continue; }
508 if( towerTest >0 && bit2Itor->second.back() == towerTest )
continue;
509 if( towerTest <0 && bit2Itor->
second.back() >=
abs(towerTest) )
continue;
512 if( towerTestCnt !=0 && debug_ && verbose_ >=2)
std::cout<<
"towerTestCnt : "<<towerTestCnt<<
" for towerTest : "<<towerTest<<std::endl;
517 std::vector<DetId>::iterator avoidItor; avoidItor =
find( avoidDuplicateVec.begin(), avoidDuplicateVec.end(), det);
518 if( avoidItor == avoidDuplicateVec.end() ){
519 avoidDuplicateVec.push_back(det);
524 std::map<EcalScDetId, double>::iterator scetItor = accuSCetMap.find(sc);
525 if( scetItor == accuSCetMap.end() ){
526 accuSCetMap[sc] = eerechit->energy()*
sin(theta);
527 accuSCchnMap[sc] = 1;
528 SCzsideMap[sc] = sc.zside();
530 accuSCetMap[sc] += eerechit->energy()*
sin(theta);
537 std::map<EcalTrigTowerDetId, double>::iterator ttetItor;
538 for( ttetItor = accuTTetMap.begin(); ttetItor != accuTTetMap.end(); ttetItor++){
542 double ttetVal = ttetItor->second;
544 std::map<EcalTrigTowerDetId, int>::iterator ttchnItor = accuTTchnMap.find(ttDetId);
545 if( ttchnItor == accuTTchnMap.end() ){
cout<<
"\nERROR cannot find ttDetId : "<<ttDetId<<
" in accuTTchnMap?!"<<endl<<endl; }
547 std::map<EcalTrigTowerDetId, int>::iterator ttzsideItor = TTzsideMap.find(ttDetId);
548 if( ttzsideItor == TTzsideMap.end() ){
cout<<
"\nERROR cannot find ttDetId : "<<ttDetId<<
" in TTzsideMap?!"<<endl<<endl; }
550 if( ttchnItor->second != 25 && debug_ && verbose_ >=2)
cout<<
"WARNING ... ttchnCnt : "<<ttchnItor->second<<
" NOT equal 25!"<<endl;
552 if( ttetVal >= tpValCut ){ isPassCut = 1; isPassCut *= ttzsideItor->second; }
557 std::map<EcalScDetId, double>::iterator scetItor;
558 for( scetItor = accuSCetMap.begin(); scetItor != accuSCetMap.end(); scetItor++){
562 double scetVal = scetItor->second;
564 std::map<EcalScDetId, int>::iterator scchnItor = accuSCchnMap.find(scDetId);
565 if( scchnItor == accuSCchnMap.end() ){
cout<<
"\nERROR cannot find scDetId : "<<scDetId<<
" in accuSCchnMap?!"<<endl<<endl; }
567 std::map<EcalScDetId, int>::iterator sczsideItor = SCzsideMap.find(scDetId);
568 if( sczsideItor == SCzsideMap.end() ){
cout<<
"\nERROR cannot find scDetId : "<<scDetId<<
" in SCzsideMap?!"<<endl<<endl; }
570 if( scchnItor->second != 25 && debug_ && verbose_ >=2)
cout<<
"WARNING ... scchnCnt : "<<scchnItor->second<<
" NOT equal 25!"<<endl;
572 if( scetVal >= tpValCut ){ isPassCut = 1; isPassCut *= sczsideItor->second; }
576 if( debug_ && verbose_ >=2)
std::cout<<
"***end setEvtTPstatusRecHits***"<<std::endl;
585 if( debug_ && verbose_ >=2)
std::cout<<
"***begin setEvtTPstatus***"<<std::endl;
589 std::map<DetId, std::vector<int> >::iterator bitItor;
590 for(bitItor = EcalAllDeadChannelsBitMap.begin(); bitItor != EcalAllDeadChannelsBitMap.end(); bitItor++){
592 DetId maskedDetId = bitItor->first;
593 int subdet = bitItor->second.front(),
status = bitItor->second.back();
596 if( !doEEfilter_ && subdet != 1 )
continue;
598 std::map<DetId, EcalTrigTowerDetId>::iterator ttItor = EcalAllDeadChannelsTTMap.find(maskedDetId);
599 if( ttItor == EcalAllDeadChannelsTTMap.end() )
continue;
602 if( chnStatus >0 &&
status == chnStatus ) toDo =
true;
603 if( chnStatus <0 && status >=
abs(chnStatus) ) toDo =
true;
608 int ttzside = ttDetId.
zside();
611 tpDigis = pTPDigis.product();
613 if( tp != tpDigis->
end() ){
614 double tpEt = ecalScale_.getTPGInGeV( tp->compressedEt(), tp->id() );
615 if(tpEt >= tpValCut ){ isPassCut = 1; isPassCut *= ttzside; }
620 if( debug_ && verbose_ >=2)
std::cout<<
"***end setEvtTPstatus***"<<std::endl;
628 EcalAllDeadChannelsValMap.clear(); EcalAllDeadChannelsBitMap.clear();
631 for(
int ieta=-85; ieta<=85; ieta++ ){
632 for(
int iphi=0; iphi<=360; iphi++ ){
638 int status = ( chit != ecalStatus->end() ) ? chit->getStatusCode() & 0x1F : -1;
646 if(status >= maskedEcalChannelStatusThreshold_){
647 std::vector<double> valVec; std::vector<int> bitVec;
648 valVec.push_back(eta); valVec.push_back(phi); valVec.push_back(theta);
649 bitVec.push_back(1); bitVec.push_back(ieta); bitVec.push_back(iphi); bitVec.push_back(status);
650 EcalAllDeadChannelsValMap.insert( std::make_pair(detid, valVec) );
651 EcalAllDeadChannelsBitMap.insert( std::make_pair(detid, bitVec) );
657 for(
int ix=0; ix<=100; ix++ ){
658 for(
int iy=0; iy<=100; iy++ ){
659 for(
int iz=-1; iz<=1; iz++ ){
665 int status = ( chit != ecalStatus->end() ) ? chit->getStatusCode() & 0x1F : -1;
673 if(status >= maskedEcalChannelStatusThreshold_){
674 std::vector<double> valVec; std::vector<int> bitVec;
675 valVec.push_back(eta); valVec.push_back(phi); valVec.push_back(theta);
676 bitVec.push_back(2); bitVec.push_back(ix); bitVec.push_back(iy); bitVec.push_back(iz); bitVec.push_back(status);
677 EcalAllDeadChannelsValMap.insert( std::make_pair(detid, valVec) );
678 EcalAllDeadChannelsBitMap.insert( std::make_pair(detid, bitVec) );
684 EcalAllDeadChannelsTTMap.clear();
685 std::map<DetId, std::vector<int> >::iterator bitItor;
686 for(bitItor = EcalAllDeadChannelsBitMap.begin(); bitItor != EcalAllDeadChannelsBitMap.end(); bitItor++){
687 const DetId id = bitItor->first;
689 EcalAllDeadChannelsTTMap.insert(std::make_pair(
id, ttDetId) );
std::string releaseVersion_
EventNumber_t event() const
const bool makeProfileRoot_
edm::EDGetTokenT< EcalRecHitCollection > ebReducedRecHitCollectionToken_
void getAllProvenance(std::vector< Provenance const * > &provenances) const
std::map< DetId, std::vector< int > > EcalAllDeadChannelsBitMap
std::map< EcalScDetId, int > accuSCchnMap
edm::ESHandle< CaloGeometry > geometry
edm::LuminosityBlockNumber_t ls
bool getByToken(EDGetToken token, Handle< PROD > &result) const
const std::string profileRootName_
#define DEFINE_FWK_MODULE(type)
Sin< T >::type sin(const T &t)
Geom::Phi< T > phi() const
std::vector< EcalRecHit >::const_iterator const_iterator
Geom::Theta< T > theta() const
unsigned long long EventNumber_t
edm::LuminosityBlockNumber_t luminosityBlock() const
collection_type const & data() const
std::vector< DetId > avoidDuplicateVec
const int maskedEcalChannelStatusThreshold_
int getChannelStatusMaps()
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
std::map< EcalTrigTowerDetId, int > accuTTchnMap
unsigned int LuminosityBlockNumber_t
std::map< EcalScDetId, double > accuSCetMap
const double etValToBeFlagged_
int zside() const
get the z-side of the tower (1/-1)
virtual const CaloCellGeometry * getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
static bool validDetId(int i, int j)
check if a valid index combination
virtual void endJob() override
Geom::Theta< T > theta() const
U second(std::pair< T, U > const &p)
edm::EDGetTokenT< EcalRecHitCollection > eeReducedRecHitCollectionToken_
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
bool getEventInfoForFilterOnce_
Abs< T >::type abs(const T &t)
void loadEventInfoForFilter(const edm::Event &iEvent)
std::vector< std::string > * cutFlowStrTmpPtr
edm::Handle< EcalRecHitCollection > endcapReducedRecHitsHandle
How EventSelector::AcceptEvent() decides whether to accept an event for output otherwise it is excluding the probing of A single or multiple positive and the trigger will pass if any such matching triggers are PASS or EXCEPTION[A criterion thatmatches no triggers at all is detected and causes a throw.] A single negative with an expectation of appropriate bit checking in the decision and the trigger will pass if any such matching triggers are FAIL or EXCEPTION A wildcarded negative criterion that matches more than one trigger in the trigger but the state exists so we define the behavior If all triggers are the negative crieriion will lead to accepting the event(this again matches the behavior of"!*"before the partial wildcard feature was incorporated).The per-event"cost"of each negative criterion with multiple relevant triggers is about the same as!*was in the past
virtual void beginJob() override
edm::Handle< EcalRecHitCollection > barrelReducedRecHitsHandle
static const int ETAPHIMODE
std::map< EcalScDetId, int > SCzsideMap
edm::ESHandle< EcalTrigTowerConstituentsMap > ttMap_
const_iterator end() const
virtual ProcessHistory const & processHistory() const
std::map< EcalTrigTowerDetId, double > accuTTetMap
std::map< DetId, std::vector< double > > EcalAllDeadChannelsValMap
int setEvtRecHitstatus(const double &tpValCut, const int &chnStatus, const int &towerTest)
static bool validDetId(int crystal_ix, int crystal_iy, int iz)
std::vector< Item >::const_iterator const_iterator
void loadEventInfo(const edm::Event &iEvent, const edm::EventSetup &iSetup)
std::string const & moduleLabel() const
virtual bool filter(edm::Event &, const edm::EventSetup &) override
Geom::Phi< T > phi() const
const edm::InputTag ebReducedRecHitCollection_
virtual void beginRun(const edm::Run &, const edm::EventSetup &) override
tuple EcalDeadCellTriggerPrimitiveFilter
ESHandle< TrackerGeometry > geometry
iterator find(key_type k)
EcalDeadCellTriggerPrimitiveFilter(const edm::ParameterSet &)
edm::EDGetTokenT< EcalTrigPrimDigiCollection > tpDigiCollectionToken_
std::map< EcalTrigTowerDetId, int > TTzsideMap
const edm::InputTag eeReducedRecHitCollection_
void loadEcalRecHits(edm::Event &iEvent, const edm::EventSetup &iSetup)
virtual void envSet(const edm::EventSetup &)
const edm::InputTag tpDigiCollection_
std::map< DetId, EcalTrigTowerDetId > EcalAllDeadChannelsTTMap
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
edm::ESHandle< EcalChannelStatus > ecalStatus
int setEvtTPstatus(const double &tpCntCut, const int &chnStatus)
void loadEcalDigis(edm::Event &iEvent, const edm::EventSetup &iSetup)
std::vector< int > * cutFlowFlagTmpPtr
const_iterator begin() const
~EcalDeadCellTriggerPrimitiveFilter()
edm::Handle< EcalTrigPrimDigiCollection > pTPDigis