CMS 3D CMS Logo

Public Types | Public Member Functions | Private Member Functions | Private Attributes

Multi5x5ClusterAlgo Class Reference

#include <Multi5x5ClusterAlgo.h>

List of all members.

Public Types

typedef math::XYZPoint Point
 point in the space
enum  VerbosityLevel { pDEBUG = 0, pWARNING = 1, pINFO = 2, pERROR = 3 }

Public Member Functions

std::vector< reco::BasicClustermakeClusters (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 (double ebst, double ecst, std::vector< int > v_chstatus, const PositionCalc &posCalc, VerbosityLevel the_verbosity=pERROR)
 Multi5x5ClusterAlgo ()
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< DetIdcanSeed_s
std::vector< reco::BasicClusterclusters_v
std::vector< std::pair< DetId,
float > > 
current_v
reco::CaloID::Detectors detector_
 The ecal region used.
double ecalBarrelSeedThreshold
double ecalEndcapSeedThreshold
PositionCalc posCalculator_
const EcalRecHitCollectionrecHits_
std::vector< EcalRecHitseeds
std::set< DetIdused_s
std::vector< int > v_chstatus_
VerbosityLevel verbosity

Detailed Description

Definition at line 28 of file Multi5x5ClusterAlgo.h.


Member Typedef Documentation

point in the space

Definition at line 62 of file Multi5x5ClusterAlgo.h.


Member Enumeration Documentation

Enumerator:
pDEBUG 
pWARNING 
pINFO 
pERROR 

Definition at line 32 of file Multi5x5ClusterAlgo.h.

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

Constructor & Destructor Documentation

Multi5x5ClusterAlgo::Multi5x5ClusterAlgo ( ) [inline]

Definition at line 34 of file Multi5x5ClusterAlgo.h.

                        {
  }
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 ( ) [inline, virtual]

Definition at line 43 of file Multi5x5ClusterAlgo.h.

    {
    }

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(), edm::SortedCollection< T, SORT >::find(), spr::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, edm::SortedCollection< T, SORT >::find(), spr::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(), crabWrap::threshold, used_s, and verbosity.

Referenced by EgammaHLTMulti5x5ClusterProducer::clusterizeECALPart(), and Multi5x5ClusterProducer::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

Definition at line 84 of file Multi5x5ClusterAlgo.h.

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

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().

The ecal region used.

Definition at line 70 of file Multi5x5ClusterAlgo.h.

Referenced by makeCluster(), and makeClusters().

Definition at line 73 of file Multi5x5ClusterAlgo.h.

Referenced by makeClusters().

Definition at line 74 of file Multi5x5ClusterAlgo.h.

Referenced by makeClusters().

Definition at line 67 of file Multi5x5ClusterAlgo.h.

Referenced by makeCluster(), and Multi5x5ClusterAlgo().

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().

Definition at line 98 of file Multi5x5ClusterAlgo.h.

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