#include <HcalHaloAlgo.h>

Public Member Functions
reco::HcalHaloData	Calculate (const CaloGeometry &TheCaloGeometry, edm::Handle< HBHERecHitCollection > &TheHBHERecHits, edm::Handle< CaloTowerCollection > &TheCaloTowers)

reco::HcalHaloData	Calculate (const CaloGeometry &TheCaloGeometry, edm::Handle< HBHERecHitCollection > &TheHBHERecHits)

float	GetHBRecHitEnergyThreshold ()

float	GetHERecHitEnergyThreshold ()

float	GetPhiWedgeEnergyThreshold ()

int	GetPhiWedgeNHitsThreshold ()

	HcalHaloAlgo ()

void	SetPhiWedgeEnergyThreshold (float SumE)

void	SetPhiWedgeNHitsThreshold (int nhits)

void	SetPhiWedgeThresholds (float SumE, int nhits)

void	SetRecHitEnergyThresholds (float HB, float HE)

	~HcalHaloAlgo ()

Private Attributes
float	HBRecHitEnergyThreshold

float	HERecHitEnergyThreshold

int	NHitsThreshold

float	SumEnergyThreshold

Detailed Description

Definition at line 30 of file HcalHaloAlgo.h.

Constructor & Destructor Documentation

HcalHaloAlgo::HcalHaloAlgo ( )

Definition at line 21 of file HcalHaloAlgo.cc.

 {
   HBRecHitEnergyThreshold = 0.;
   HERecHitEnergyThreshold = 0.;
   SumEnergyThreshold = 0.;
   NHitsThreshold = 0;
 }

HcalHaloAlgo::~HcalHaloAlgo ( )

inline

Definition at line 35 of file HcalHaloAlgo.h.

35 {}

Member Function Documentation

HcalHaloData HcalHaloAlgo::Calculate	(	const CaloGeometry &	TheCaloGeometry,
		edm::Handle< HBHERecHitCollection > &	TheHBHERecHits,
		edm::Handle< CaloTowerCollection > &	TheCaloTowers
	)

Definition at line 34 of file HcalHaloAlgo.cc.

References funct::abs(), Abs(), CompareTime(), CompareTowers(), CaloTower::emEt(), F(), reco::HcalHaloData::GetPhiWedges(), reco::HcalHaloData::getProblematicStrips(), CaloTower::hadEt(), HcalBarrel, HcalEndcap, i, hit::id, CaloTower::id(), HcalDetId::ieta(), CaloTower::ieta(), CaloTower::iphi(), j, bookConverter::max, CaloTower::numProblematicHcalCells(), and reco::PhiWedge::SetPlusZOriginConfidence().

Referenced by reco::HcalHaloDataProducer::produce().

 {
   HcalHaloData TheHcalHaloData;
   
   // Store Energy sum of rechits as a function of iPhi (iPhi goes from 1 to 72)
   float SumE[73];
   // Store Number of rechits as a function of iPhi 
   int NumHits[73];
   // Store minimum time of rechit as a function of iPhi
   float MinTimeHits[73];
   // Store maximum time of rechit as a function of iPhi
   float MaxTimeHits[73];
   for(unsigned int i = 0 ; i < 73 ; i++ ) 
     {
       SumE[i] = 0;
       NumHits[i]= 0;
       MinTimeHits[i] = 0.;
       MaxTimeHits[i] = 0.;
     }
   
   for( HBHERecHitCollection::const_iterator hit = TheHBHERecHits->begin() ; hit != TheHBHERecHits->end() ; hit++ )
     {
       HcalDetId id = HcalDetId(hit->id());                                                                                                    
       switch ( id.subdet() )                                                                                         
     {      
     case HcalBarrel:                                                                           
       if(hit->energy() < HBRecHitEnergyThreshold )continue;
       break;                                                                                                                  
     case HcalEndcap:                                                                                          
       if(hit->energy() < HERecHitEnergyThreshold ) continue;
       break;
     default:
       continue;
     }
       
       int iEta = id.ieta();
       int iPhi = id.iphi();
       if(iPhi < 73 && TMath::Abs(iEta) < 23 )
     { 
       SumE[iPhi]+= hit->energy();
       NumHits[iPhi] ++;
       
       float time = hit->time();
       MinTimeHits[iPhi] = time < MinTimeHits[iPhi] ? time : MinTimeHits[iPhi];
       MaxTimeHits[iPhi] = time > MaxTimeHits[iPhi] ? time : MaxTimeHits[iPhi];
     }
     }
   
   for( int iPhi = 1 ; iPhi < 73 ; iPhi++ )
     {
       if( SumE[iPhi] >= SumEnergyThreshold && NumHits[iPhi] > NHitsThreshold )
     {
       // Build PhiWedge and store to HcalHaloData if energy or #hits pass thresholds
       PhiWedge wedge(SumE[iPhi], iPhi, NumHits[iPhi], MinTimeHits[iPhi], MaxTimeHits[iPhi]);
       
       // Loop over rechits again to calculate direction based on timing info
       std::vector<const HBHERecHit*> Hits;
       for( HBHERecHitCollection::const_iterator hit = TheHBHERecHits->begin() ; hit != TheHBHERecHits->end() ; hit++ )
         {
 
           HcalDetId id = HcalDetId(hit->id());
           if( id.iphi() != iPhi ) continue;
           if( TMath::Abs(id.ieta() ) > 22 ) continue;  // has to overlap geometrically w/ HB
           switch ( id.subdet() )
         {
         case HcalBarrel:
           if(hit->energy() < HBRecHitEnergyThreshold )continue;
           break;
         case HcalEndcap:
           if(hit->energy() < HERecHitEnergyThreshold ) continue;
           break;
         default:
           continue;
         }
           Hits.push_back(&(*hit));
         }
           
       std::sort( Hits.begin() , Hits.end() , CompareTime);
       float MinusToPlus = 0.;
       float PlusToMinus = 0.;
       for( unsigned int i = 0 ; i < Hits.size() ; i++ )
         {
           HcalDetId id_i = HcalDetId(Hits[i]->id() );
           int ieta_i = id_i.ieta();
           for( unsigned int j = (i+1) ; j < Hits.size() ; j++ )
         {
           HcalDetId id_j = HcalDetId(Hits[j]->id() );
           int ieta_j = id_j.ieta();
           if( ieta_i > ieta_j ) PlusToMinus += TMath::Abs(ieta_i - ieta_j ) ;
           else MinusToPlus += TMath::Abs(ieta_i - ieta_j);
         }
         }
       float PlusZOriginConfidence = (PlusToMinus + MinusToPlus )? PlusToMinus / ( PlusToMinus + MinusToPlus ) : -1. ;
       wedge.SetPlusZOriginConfidence( PlusZOriginConfidence );
       TheHcalHaloData.GetPhiWedges().push_back( wedge );
     }
     }
 
 
   // Don't use HF.
   int maxAbsIEta = 29;
 
 
   std::map<int, float> iPhiHadEtMap;
   std::vector<const CaloTower*> sortedCaloTowers;
   for(CaloTowerCollection::const_iterator tower = TheCaloTowers->begin(); tower != TheCaloTowers->end(); tower++) {
     if(abs(tower->ieta()) > maxAbsIEta) continue;
 
     int iPhi = tower->iphi();
     if(!iPhiHadEtMap.count(iPhi)) iPhiHadEtMap[iPhi] = 0.0;
     iPhiHadEtMap[iPhi] += tower->hadEt();
 
     if(tower->numProblematicHcalCells() > 0) sortedCaloTowers.push_back(&(*tower));
 
   }
 
 
   // Sort towers such that lowest iphi and ieta are first, highest last, and towers
   // with same iphi value are consecutive. Then we can do everything else in one loop.
   std::sort(sortedCaloTowers.begin(), sortedCaloTowers.end(), CompareTowers);
 
   HaloTowerStrip strip;
 
 
   int prevIEta = -99, prevIPhi = -99;
   float prevHadEt = 0.;
   float prevEmEt = 0.;
   std::pair<uint8_t, CaloTowerDetId> prevPair, towerPair;
   bool wasContiguous = true;
 
   // Loop through and store a vector of pairs (problematicCells, DetId) for each contiguous strip we find
   for(unsigned int i = 0; i < sortedCaloTowers.size(); i++) {
     const CaloTower* tower = sortedCaloTowers[i];
 
     towerPair = std::make_pair((uint8_t)tower->numProblematicHcalCells(), tower->id());
 
     bool newIPhi = tower->iphi() != prevIPhi;
     bool isContiguous = tower->ieta() == 1 ? tower->ieta() - 2 == prevIEta : tower->ieta() - 1 == prevIEta;
 
     isContiguous = isContiguous || (tower->ieta() == -maxAbsIEta);
     if(newIPhi) isContiguous = false;
 
     if(!wasContiguous && isContiguous) {
       strip.cellTowerIds.push_back(prevPair);
       strip.cellTowerIds.push_back(towerPair);
       strip.hadEt += prevHadEt + tower->hadEt();
       strip.emEt += prevEmEt + tower->emEt();
     }
 
     if(wasContiguous && isContiguous) {
       strip.cellTowerIds.push_back(towerPair);
       strip.hadEt += tower->hadEt();
       strip.emEt += tower->emEt();
     }
 
     if((wasContiguous && !isContiguous) || i == sortedCaloTowers.size()-1) { //ended the strip, so flush it
 
       if(strip.cellTowerIds.size() > 3) {
 
         int iPhi = strip.cellTowerIds.at(0).second.iphi();
         int iPhiLower = (iPhi == 1) ? 72 : iPhi - 1;
         int iPhiUpper = (iPhi == 72) ? 1 : iPhi + 1;
 
         float energyRatio = 0.0;
         if(iPhiHadEtMap.count(iPhiLower)) energyRatio += iPhiHadEtMap[iPhiLower];
         if(iPhiHadEtMap.count(iPhiUpper)) energyRatio += iPhiHadEtMap[iPhiUpper];
         iPhiHadEtMap[iPhi] = max(iPhiHadEtMap[iPhi], 0.001F);
 
         energyRatio /= iPhiHadEtMap[iPhi];
         strip.energyRatio = energyRatio;
 
         TheHcalHaloData.getProblematicStrips().push_back( strip );
 
       }
       strip = HaloTowerStrip();
     }
 
     wasContiguous = isContiguous;
     prevPair = towerPair;
     prevEmEt = tower->emEt();
     prevIPhi = tower->iphi();
     prevIEta = tower->ieta();
     prevHadEt = tower->hadEt();
   }
 
   return TheHcalHaloData;
   
 }

HcalHaloData HcalHaloAlgo::Calculate	(	const CaloGeometry &	TheCaloGeometry,
		edm::Handle< HBHERecHitCollection > &	TheHBHERecHits
	)

Definition at line 29 of file HcalHaloAlgo.cc.

                                                                                                                          {
     edm::Handle<CaloTowerCollection> TheCaloTowers;
     return Calculate(TheCaloGeometry, TheHBHERecHits, TheCaloTowers);
 }