#include <Multi5x5ClusterAlgo.h>

Public Types
typedef math::XYZPoint	Point
	point in the space More...

enum	VerbosityLevel { pDEBUG = 0, pWARNING = 1, pINFO = 2, pERROR = 3 }

Public Member Functions
std::vector< reco::BasicCluster >	makeClusters (const EcalRecHitCollection hits, const CaloSubdetectorGeometry geometry, const CaloSubdetectorTopology topology_p, const CaloSubdetectorGeometry geometryES_p, reco::CaloID::Detectors detector, bool regional=false, const std::vector< EcalEtaPhiRegion > &regions=std::vector< EcalEtaPhiRegion >())

	Multi5x5ClusterAlgo ()

	Multi5x5ClusterAlgo (double ebst, double ecst, std::vector< int > v_chstatus, const PositionCalc &posCalc, VerbosityLevel the_verbosity=pERROR)

void	setVerbosity (VerbosityLevel the_verbosity)

virtual	~Multi5x5ClusterAlgo ()

Private Member Functions
void	addCrystal (const DetId &det)

bool	checkMaxima (CaloNavigator< DetId > &navigator, const EcalRecHitCollection *hits)

void	mainSearch (const EcalRecHitCollection hits, const CaloSubdetectorGeometry geometry_p, const CaloSubdetectorTopology topology_p, const CaloSubdetectorGeometry geometryES_p)

void	makeCluster (const EcalRecHitCollection hits, const CaloSubdetectorGeometry geometry_p, const CaloSubdetectorGeometry *geometryES_p, const EcalRecHitCollection::const_iterator &seedIt)

void	prepareCluster (CaloNavigator< DetId > &navigator, const EcalRecHitCollection hits, const CaloSubdetectorGeometry geometry)

Private Attributes
std::set< DetId >	canSeed_s

std::vector< reco::BasicCluster >	clusters_v

std::vector< std::pair< DetId, float > >	current_v

reco::CaloID::Detectors	detector_
	The ecal region used. More...

double	ecalBarrelSeedThreshold

double	ecalEndcapSeedThreshold

PositionCalc	posCalculator_

const EcalRecHitCollection *	recHits_

std::vector< EcalRecHit >	seeds

std::set< DetId >	used_s

std::vector< int >	v_chstatus_

VerbosityLevel	verbosity

Detailed Description

Definition at line 28 of file Multi5x5ClusterAlgo.h.

Member Typedef Documentation

typedef math::XYZPoint Multi5x5ClusterAlgo::Point

point in the space

Definition at line 62 of file Multi5x5ClusterAlgo.h.

Member Enumeration Documentation

enum Multi5x5ClusterAlgo::VerbosityLevel

Enumerator
pDEBUG
pWARNING
pINFO
pERROR

Definition at line 32 of file Multi5x5ClusterAlgo.h.

32 { pDEBUG = 0, pWARNING = 1, pINFO = 2, pERROR = 3 };

Multi5x5ClusterAlgo::pWARNING

Definition: Multi5x5ClusterAlgo.h:32

Multi5x5ClusterAlgo::pINFO

Definition: Multi5x5ClusterAlgo.h:32

Multi5x5ClusterAlgo::pDEBUG

Definition: Multi5x5ClusterAlgo.h:32

Multi5x5ClusterAlgo::pERROR

Definition: Multi5x5ClusterAlgo.h:32

Constructor & Destructor Documentation

Multi5x5ClusterAlgo::Multi5x5ClusterAlgo ( )

inline

Definition at line 34 of file Multi5x5ClusterAlgo.h.

34 {

35 }

Multi5x5ClusterAlgo::Multi5x5ClusterAlgo	(	double	ebst,
		double	ecst,
		std::vector< int >	v_chstatus,
		const PositionCalc &	posCalc,
		VerbosityLevel	the_verbosity = `pERROR`
	)

inline

Definition at line 37 of file Multi5x5ClusterAlgo.h.

References posCalculator_, python.multivaluedict::sort(), and v_chstatus_.

                                                                                                                                                : 
     ecalBarrelSeedThreshold(ebst), ecalEndcapSeedThreshold(ecst),  v_chstatus_(v_chstatus), verbosity(the_verbosity) {
     posCalculator_ = posCalc;
     std::sort( v_chstatus_.begin(), v_chstatus_.end() );
   }

virtual Multi5x5ClusterAlgo::~Multi5x5ClusterAlgo ( )

inlinevirtual

Definition at line 43 of file Multi5x5ClusterAlgo.h.

44 {

45 }

Member Function Documentation

void Multi5x5ClusterAlgo::addCrystal ( const DetId & det )

private

Definition at line 346 of file Multi5x5ClusterAlgo.cc.

References current_v, edm::SortedCollection< T, SORT >::end(), edm::SortedCollection< T, SORT >::find(), recHits_, and used_s.

Referenced by prepareCluster().

 {   
 
    EcalRecHitCollection::const_iterator thisIt =  recHits_->find(det);
    if ((thisIt != recHits_->end()) && (thisIt->id() != DetId(0)))
    { 
       if ((used_s.find(thisIt->id()) == used_s.end())) 
       {
      //std::cout << "   ... this is a good crystal and will be added" << std::endl;
          current_v.push_back( std::pair<DetId, float>(det, 1.) ); // by default hit energy fractions are set at 1.
          used_s.insert(det);
       }
    } 
   
 }

bool Multi5x5ClusterAlgo::checkMaxima	(	CaloNavigator< DetId > &	navigator,
		const EcalRecHitCollection *	hits
	)

private

Definition at line 237 of file Multi5x5ClusterAlgo.cc.

References CaloNavigator< T >::east(), edm::SortedCollection< T, SORT >::end(), spr::find(), edm::SortedCollection< T, SORT >::find(), CaloNavigator< T >::home(), i, CaloNavigator< T >::north(), CaloNavigator< T >::pos(), recHits_, CaloNavigator< T >::south(), v_chstatus_, and CaloNavigator< T >::west().

Referenced by mainSearch().

 {
 
    bool maxima = true;
    EcalRecHitCollection::const_iterator thisHit;
    EcalRecHitCollection::const_iterator seedHit = hits->find(navigator.pos());
    double seedEnergy = seedHit->energy();
 
    std::vector<DetId> swissCrossVec;
    swissCrossVec.clear();
 
    swissCrossVec.push_back(navigator.west());
    navigator.home();
    swissCrossVec.push_back(navigator.east());
    navigator.home();
    swissCrossVec.push_back(navigator.north());
    navigator.home();
    swissCrossVec.push_back(navigator.south());
    navigator.home();
 
    std::vector<DetId>::const_iterator detItr;
    for (unsigned int i = 0; i < swissCrossVec.size(); ++i)
    {
 
       // look for this hit
       thisHit = recHits_->find(swissCrossVec[i]);
 
       // continue if this hit was not found
       if  ((swissCrossVec[i] == DetId(0)) || thisHit == recHits_->end()) continue; 
 
       // the recHit has to be skipped in the local maximum search if it was found
       // in the map of channels to be excluded 
       uint32_t rhFlag = thisHit->recoFlag();
       std::vector<int>::const_iterator vit = std::find(v_chstatus_.begin(), v_chstatus_.end(), rhFlag);
       if (vit != v_chstatus_.end()) continue;
 
       // if this crystal has more energy than the seed then we do 
       // not have a local maxima
       if (thisHit->energy() > seedEnergy)
       {
          maxima = false;
          break;
       }
    }
 
    return maxima;
 
 }

void Multi5x5ClusterAlgo::mainSearch	(	const EcalRecHitCollection *	hits,
		const CaloSubdetectorGeometry *	geometry_p,
		const CaloSubdetectorTopology *	topology_p,
		const CaloSubdetectorGeometry *	geometryES_p
	)

private

Definition at line 112 of file Multi5x5ClusterAlgo.cc.

References canSeed_s, checkMaxima(), gather_cfg::cout, current_v, spr::find(), edm::SortedCollection< T, SORT >::find(), makeCluster(), pINFO, prepareCluster(), seeds, used_s, v_chstatus_, and verbosity.

Referenced by makeClusters().

 {
 
    if (verbosity < pINFO)
    {
       std::cout << "Building clusters............" << std::endl;
    }
 
    // Loop over seeds:
    std::vector<EcalRecHit>::iterator it;
    for (it = seeds.begin(); it != seeds.end(); it++)
    {
 
       // check if this crystal is able to seed
       // (event though it is already used)
       bool usedButCanSeed = false;
       if (canSeed_s.find(it->id()) != canSeed_s.end()) usedButCanSeed = true;
 
       // avoid seeding for anomalous channels (recoFlag based)
       uint32_t rhFlag = (*it).recoFlag();
       std::vector<int>::const_iterator vit = std::find( v_chstatus_.begin(), v_chstatus_.end(), rhFlag );
       if ( vit != v_chstatus_.end() ) continue; // the recHit has to be excluded from seeding
 
       // make sure the current seed does not belong to a cluster already.
       if ((used_s.find(it->id()) != used_s.end()) && (usedButCanSeed == false))
       {
          if (it == seeds.begin())
      {
         if (verbosity < pINFO)
             {
                std::cout << "##############################################################" << std::endl;
            std::cout << "DEBUG ALERT: Highest energy seed already belongs to a cluster!" << std::endl;
            std::cout << "##############################################################" << std::endl;
         }
      }
 
       // seed crystal is used or is used and cannot seed a cluster
           // so continue to the next seed crystal...
       continue;
       }
 
       // clear the vector of hits in current cluster
       current_v.clear();
 
       // Create a navigator at the seed and get seed
       // energy
       CaloNavigator<DetId> navigator(it->id(), topology_p);
       DetId seedId = navigator.pos();
       EcalRecHitCollection::const_iterator seedIt = hits->find(seedId);
       navigator.setHome(seedId);
 
       // Is the seed a local maximum?
       bool localMaxima = checkMaxima(navigator, hits);
 
       if (localMaxima)
       {
          // build the 5x5 taking care over which crystals
          // can seed new clusters and which can't
          prepareCluster(navigator, hits, geometry_p);
       }
 
       // If some crystals in the current vector then 
       // make them into a cluster 
       if (current_v.size() > 0) 
       {
          makeCluster(hits, geometry_p, geometryES_p, seedIt);
       }
 
    }  // End loop on seed crystals
 
 }

void Multi5x5ClusterAlgo::makeCluster	(	const EcalRecHitCollection *	hits,
		const CaloSubdetectorGeometry *	geometry_p,
		const CaloSubdetectorGeometry *	geometryES_p,
		const EcalRecHitCollection::const_iterator &	seedIt
	)

private

Definition at line 188 of file Multi5x5ClusterAlgo.cc.

References PositionCalc::Calculate_Location(), clusters_v, gather_cfg::cout, current_v, reco::CaloID::DET_ECAL_BARREL, reco::CaloID::DET_ECAL_ENDCAP, detector_, EcalBarrel, CaloRecHit::energy(), relval_parameters_module::energy, edm::SortedCollection< T, SORT >::find(), reco::CaloCluster::multi5x5, pINFO, posCalculator_, position, and verbosity.

Referenced by mainSearch().

 {
 
    double energy = 0;
    //double chi2   = 0;
    reco::CaloID caloID;
    Point position;
    position = posCalculator_.Calculate_Location(current_v, hits,geometry, geometryES);
   
    std::vector<std::pair<DetId, float> >::iterator it;
    for (it = current_v.begin(); it != current_v.end(); it++)
    {
       EcalRecHitCollection::const_iterator itt = hits->find( (*it).first );
       EcalRecHit hit_p = *itt;
       energy += hit_p.energy();
       //chi2 += 0;
       if ( (*it).first.subdetId() == EcalBarrel ) {
               caloID = reco::CaloID::DET_ECAL_BARREL;
       } else {
               caloID = reco::CaloID::DET_ECAL_ENDCAP;
       }
 
    }
    //chi2 /= energy;
 
    if (verbosity < pINFO)
    { 
       std::cout << "******** NEW CLUSTER ********" << std::endl;
       std::cout << "No. of crystals = " << current_v.size() << std::endl;
       std::cout << "     Energy     = " << energy << std::endl;
       std::cout << "     Phi        = " << position.phi() << std::endl;
       std::cout << "     Eta        = " << position.eta() << std::endl;
       std::cout << "*****************************" << std::endl;
     }
 
    // to be a valid cluster the cluster energy
    // must be at least the seed energy
    double seedEnergy = seedIt->energy();
    if (energy >= seedEnergy)
    {
       //clusters_v.push_back(reco::BasicCluster(energy, position, chi2, current_v, reco::CaloCluster::island));
       clusters_v.push_back(reco::BasicCluster(energy, position, reco::CaloID(detector_), current_v, reco::CaloCluster::multi5x5, seedIt->id()));
    }
 
 }

std::vector< reco::BasicCluster > Multi5x5ClusterAlgo::makeClusters	(	const EcalRecHitCollection *	hits,
		const CaloSubdetectorGeometry *	geometry,
		const CaloSubdetectorTopology *	topology_p,
		const CaloSubdetectorGeometry *	geometryES_p,
		reco::CaloID::Detectors	detector,
		bool	regional = `false`,
		const std::vector< EcalEtaPhiRegion > &	regions = `std::vector<EcalEtaPhiRegion>()`
	)

Definition at line 16 of file Multi5x5ClusterAlgo.cc.

References edm::SortedCollection< T, SORT >::begin(), canSeed_s, clusters_v, gather_cfg::cout, reco::CaloID::DET_ECAL_BARREL, reco::CaloID::DET_ECAL_ENDCAP, reco::CaloID::DET_NONE, detector_, ecalBarrelSeedThreshold, ecalEndcapSeedThreshold, edm::SortedCollection< T, SORT >::end(), relval_parameters_module::energy, ET, CaloSubdetectorGeometry::getGeometry(), CaloCellGeometry::getPosition(), mainSearch(), pINFO, position, recHits_, seeds, funct::sin(), python.multivaluedict::sort(), PV3DBase< T, PVType, FrameType >::theta(), dtT0WireCalibration_cfg::threshold, used_s, and verbosity.

Referenced by Multi5x5ClusterProducer::clusterizeECALPart(), and EgammaHLTMulti5x5ClusterProducer::clusterizeECALPart().

 {
   seeds.clear();
   used_s.clear();
   canSeed_s.clear();
   clusters_v.clear();
 
   recHits_ = hits;
 
   double threshold = 0;
   std::string ecalPart_string;
   detector_ = reco::CaloID::DET_NONE;
   if (detector == reco::CaloID::DET_ECAL_ENDCAP) 
     {
       detector_ = reco::CaloID::DET_ECAL_ENDCAP;
       threshold = ecalEndcapSeedThreshold;
       ecalPart_string = "EndCap";
     }
   if (detector == reco::CaloID::DET_ECAL_BARREL) 
     {
       detector_ = reco::CaloID::DET_ECAL_BARREL;
       threshold = ecalBarrelSeedThreshold;
       ecalPart_string = "Barrel";
     }
 
   if (verbosity < pINFO)
     {
       std::cout << "-------------------------------------------------------------" << std::endl;
       std::cout << "Island algorithm invoked for ECAL" << ecalPart_string << std::endl;
       std::cout << "Looking for seeds, energy threshold used = " << threshold << " GeV" <<std::endl;
     }
 
   int nregions=0;
   if(regional) nregions=regions.size();
 
   if(!regional || nregions) {
 
     EcalRecHitCollection::const_iterator it;
     for(it = hits->begin(); it != hits->end(); it++)
       {
     double energy = it->energy();
     if (energy < threshold) continue; // need to check to see if this line is useful!
 
     const CaloCellGeometry *thisCell = geometry_p->getGeometry(it->id());
     GlobalPoint position = thisCell->getPosition();
 
     // Require that RecHit is within clustering region in case
     // of regional reconstruction
     bool withinRegion = false;
     if (regional) {
       std::vector<EcalEtaPhiRegion>::const_iterator region;
       for (region=regions.begin(); region!=regions.end(); region++) {
         if (region->inRegion(position)) {
           withinRegion =  true;
           break;
         }
       }
     }
 
     if (!regional || withinRegion) {
       float ET = it->energy() * sin(position.theta());
       if (ET > threshold) seeds.push_back(*it);
     }
       }
     
   }
   
    sort(seeds.begin(), seeds.end(), EcalRecHitLess());
 
    if (verbosity < pINFO)
    {
       std::cout << "Total number of seeds found in event = " << seeds.size() << std::endl;
    }
 
    mainSearch(hits, geometry_p, topology_p, geometryES_p);
    sort(clusters_v.rbegin(), clusters_v.rend(), ClusterEtLess());
 
    if (verbosity < pINFO)
    {
       std::cout << "---------- end of main search. clusters have been sorted ----" << std::endl;
    }
   
    return clusters_v;
  
 }

void Multi5x5ClusterAlgo::prepareCluster	(	CaloNavigator< DetId > &	navigator,
		const EcalRecHitCollection *	hits,
		const CaloSubdetectorGeometry *	geometry
	)

private

Definition at line 287 of file Multi5x5ClusterAlgo.cc.

References abs, addCrystal(), canSeed_s, CaloNavigator< T >::home(), and CaloNavigator< T >::offsetBy().

Referenced by mainSearch().

 {
 
    DetId thisDet;
    std::set<DetId>::iterator setItr;
 
    // now add the 5x5 taking care to mark the edges
    // as able to seed and where overlapping in the central
    // region with crystals that were previously able to seed
    // change their status so they are not able to seed
    //std::cout << std::endl;
    for (int dx = -2; dx < 3; ++dx)
    {
       for (int dy = -2; dy < 3; ++ dy)
       {
 
       // navigate in free steps forming
           // a full 5x5
           thisDet = navigator.offsetBy(dx, dy);
           navigator.home();
 
           // add the current crystal
       //std::cout << "adding " << dx << ", " << dy << std::endl;
       addCrystal(thisDet);
 
       // now consider if we are in an edge (outer 16)
           // or central (inner 9) region
           if ((abs(dx) > 1) || (abs(dy) > 1))
           {    
              // this is an "edge" so should be allowed to seed
              // provided it is not already used
              //std::cout << "   setting can seed" << std::endl;
              canSeed_s.insert(thisDet);
           }  // end if "edge"
           else 
           {
              // or else we are in the central 3x3
              // and must remove any of these crystals from the canSeed set
              setItr = canSeed_s.find(thisDet);
              if (setItr != canSeed_s.end())
              {
                 //std::cout << "   unsetting can seed" << std::endl;
                 canSeed_s.erase(setItr);
              }
           }  // end if "centre"
 
 
       } // end loop on dy
 
    } // end loop on dx
 
    //std::cout << "*** " << std::endl;
    //std::cout << " current_v contains " << current_v.size() << std::endl;
    //std::cout << "*** " << std::endl;
 }

void Multi5x5ClusterAlgo::setVerbosity ( VerbosityLevel the_verbosity )

inline

Definition at line 47 of file Multi5x5ClusterAlgo.h.

References verbosity.

     {
       verbosity = the_verbosity;
     }

Member Data Documentation

std::set<DetId> Multi5x5ClusterAlgo::canSeed_s

private

Definition at line 84 of file Multi5x5ClusterAlgo.h.

Referenced by mainSearch(), makeClusters(), and prepareCluster().

std::vector<reco::BasicCluster> Multi5x5ClusterAlgo::clusters_v

private

Definition at line 92 of file Multi5x5ClusterAlgo.h.

Referenced by makeCluster(), and makeClusters().

std::vector<std::pair<DetId, float> > Multi5x5ClusterAlgo::current_v

private

Definition at line 89 of file Multi5x5ClusterAlgo.h.

Referenced by addCrystal(), mainSearch(), and makeCluster().

reco::CaloID::Detectors Multi5x5ClusterAlgo::detector_

private

The ecal region used.

Definition at line 70 of file Multi5x5ClusterAlgo.h.

Referenced by makeCluster(), and makeClusters().

double Multi5x5ClusterAlgo::ecalBarrelSeedThreshold

private

Definition at line 73 of file Multi5x5ClusterAlgo.h.

Referenced by makeClusters().

double Multi5x5ClusterAlgo::ecalEndcapSeedThreshold

private

Definition at line 74 of file Multi5x5ClusterAlgo.h.

Referenced by makeClusters().

PositionCalc Multi5x5ClusterAlgo::posCalculator_

private

Definition at line 67 of file Multi5x5ClusterAlgo.h.

Referenced by makeCluster(), and Multi5x5ClusterAlgo().

const EcalRecHitCollection* Multi5x5ClusterAlgo::recHits_

private

Definition at line 77 of file Multi5x5ClusterAlgo.h.

Referenced by addCrystal(), checkMaxima(), and makeClusters().

std::vector<EcalRecHit> Multi5x5ClusterAlgo::seeds

private

Definition at line 80 of file Multi5x5ClusterAlgo.h.

Referenced by mainSearch(), and makeClusters().

std::set<DetId> Multi5x5ClusterAlgo::used_s

private

Definition at line 83 of file Multi5x5ClusterAlgo.h.

Referenced by addCrystal(), mainSearch(), and makeClusters().

std::vector<int> Multi5x5ClusterAlgo::v_chstatus_

private

Definition at line 95 of file Multi5x5ClusterAlgo.h.

Referenced by checkMaxima(), mainSearch(), and Multi5x5ClusterAlgo().

VerbosityLevel Multi5x5ClusterAlgo::verbosity

private

Definition at line 98 of file Multi5x5ClusterAlgo.h.

Referenced by mainSearch(), makeCluster(), makeClusters(), and setVerbosity().

Public Types

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Member Typedef Documentation

Member Enumeration Documentation

Constructor & Destructor Documentation

Member Function Documentation

Member Data Documentation