Multi5x5ClusterAlgo.h

Go to the documentation of this file.

#ifndef RecoECAL_ECALClusters_Multi5x5ClusterAlgo_h
#define RecoECAL_ECALClusters_Multi5x5ClusterAlgo_h
 
#include "Geometry/CaloGeometry/interface/CaloGeometry.h"
 
#include "DataFormats/Math/interface/Point3D.h"
#include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
#include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
#include "DataFormats/EcalDetId/interface/EBDetId.h"
#include "DataFormats/EgammaReco/interface/BasicCluster.h"
#include "DataFormats/Math/interface/RectangularEtaPhiRegion.h"
 
#include "RecoCaloTools/Navigation/interface/CaloNavigator.h"
#include "RecoCaloTools/Navigation/interface/EcalBarrelNavigator.h"
#include "RecoCaloTools/Navigation/interface/EcalEndcapNavigator.h"
#include "Geometry/CaloTopology/interface/EcalBarrelHardcodedTopology.h"
#include "RecoEcal/EgammaCoreTools/interface/PositionCalc.h"
#include "DataFormats/CaloRecHit/interface/CaloID.h"
 
// C/C++ headers
#include <string>
#include <vector>
#include <set>
 
typedef std::map<DetId, EcalRecHit> RecHitsMap;
 
class Multi5x5ClusterAlgo {
public:
  //the 5x5 clustering algo by default makes basic clusters which may not contain their seed crystal if they are close by to other clusters
  //however we would like to post-fix the basic clusters to ensure they always contain their seed crystal
  //so we define a proto basic cluster class which contains all the information which would be in a basic cluster
  //which allows the addition of its seed and the removal of a seed of another cluster easily
  class ProtoBasicCluster {
    float energy_;
    EcalRecHit seed_;
    std::vector<std::pair<DetId, float> > hits_;
    bool containsSeed_;
 
  public:
    ProtoBasicCluster();
    ProtoBasicCluster(float iEnergy, const EcalRecHit &iSeed, std::vector<std::pair<DetId, float> > &iHits)
        : energy_(iEnergy), seed_(iSeed) {
      hits_.swap(iHits);
      containsSeed_ = isSeedCrysInHits_();
    }
 
    float energy() const { return energy_; }
    const EcalRecHit &seed() const { return seed_; }
    const std::vector<std::pair<DetId, float> > &hits() const { return hits_; }
    bool containsSeed() const { return containsSeed_; }
 
    bool removeHit(const EcalRecHit &hitToRM);
    bool addSeed();
 
  private:
    bool isSeedCrysInHits_() const;
  };
 
  Multi5x5ClusterAlgo() {}
 
  Multi5x5ClusterAlgo(double ebst,
                      double ecst,
                      const std::vector<int> &v_chstatus,
                      const PositionCalc &posCalc,
                      bool reassignSeedCrysToClusterItSeeds = false)
      : ecalBarrelSeedThreshold(ebst),
        ecalEndcapSeedThreshold(ecst),
        v_chstatus_(v_chstatus),
        reassignSeedCrysToClusterItSeeds_(reassignSeedCrysToClusterItSeeds) {
    posCalculator_ = posCalc;
    std::sort(v_chstatus_.begin(), v_chstatus_.end());
  }
 
  virtual ~Multi5x5ClusterAlgo() {}
 
  // this is the method that will start the clusterisation
  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<RectangularEtaPhiRegion> &regions = std::vector<RectangularEtaPhiRegion>());
 
  typedef math::XYZPoint Point;
 
private:
  //algo to compute position of clusters
  PositionCalc posCalculator_;
 
  reco::CaloID::Detectors detector_;
 
  // Energy required for a seed:
  double ecalBarrelSeedThreshold;
  double ecalEndcapSeedThreshold;
 
  // collection of all rechits
  const EcalRecHitCollection *recHits_;
 
  // The vector of seeds:
  std::vector<EcalRecHit> seeds;
 
  std::vector<std::pair<DetId, int> > whichClusCrysBelongsTo_;
 
  // The set of used DetID's
  std::set<DetId> used_s;
  std::set<DetId> canSeed_s;  // set of crystals not to be added but which can seed
  // a new 3x3 (e.g. the outer crystals in a 5x5)
 
  // The vector of DetId's in the cluster currently reconstructed
  std::vector<std::pair<DetId, float> > current_v;
 
  // The vector of clusters
  std::vector<reco::BasicCluster> clusters_v;
  std::vector<ProtoBasicCluster> protoClusters_;
  // recHit flag to be excluded from seeding
  std::vector<int> v_chstatus_;
 
  bool reassignSeedCrysToClusterItSeeds_;  //the seed of the 5x5 crystal is sometimes in another basic cluster, however we may want to put it back into the cluster it seeds
 
  void mainSearch(const EcalRecHitCollection *hits,
                  const CaloSubdetectorGeometry *geometry_p,
                  const CaloSubdetectorTopology *topology_p,
                  const CaloSubdetectorGeometry *geometryES_p);
 
  // Is the crystal at the navigator position a
  // local maxiumum in energy?
  bool checkMaxima(CaloNavigator<DetId> &navigator, const EcalRecHitCollection *hits);
 
  // prepare the 5x5 taking care over which crystals
  // are allowed to seed new clusters and which are not
  // after the preparation is complete
  void prepareCluster(CaloNavigator<DetId> &navigator,
                      const EcalRecHitCollection *hits,
                      const CaloSubdetectorGeometry *geometry);
 
  // Add the crystal with DetId det to the current
  // vector of crystals if it meets certain criteria
  void addCrystal(const DetId &det);
 
  // take the crystals in the current_v and build
  // them into a BasicCluster
  void makeCluster(const EcalRecHitCollection *hits,
                   const CaloSubdetectorGeometry *geometry_p,
                   const CaloSubdetectorGeometry *geometryES_p,
                   const EcalRecHitCollection::const_iterator &seedIt,
                   bool seedOutside);
};
 
#endif