#include <ClusterShapeAlgo.h>

Public Member Functions
reco::ClusterShape	Calculate (const reco::BasicCluster &passedCluster, const EcalRecHitCollection hits, const CaloSubdetectorGeometry geometry, const CaloSubdetectorTopology *topology)
	ClusterShapeAlgo (const edm::ParameterSet &par)
	ClusterShapeAlgo ()
Private Types
enum	{ Eta, Phi }
Private Member Functions
double	absZernikeMoment (const reco::BasicCluster &passedCluster, int n, int m, double R0=6.6)
double	calc_AbsZernikeMoment (const reco::BasicCluster &passedCluster, int n, int m, double R0)
void	Calculate_2ndEnergy (const reco::BasicCluster &passedCluster, const EcalRecHitCollection *hits)
void	Calculate_BarrelBasketEnergyFraction (const reco::BasicCluster &passedCluster, const EcalRecHitCollection hits, const int EtaPhi, const CaloSubdetectorGeometry geometry)
void	Calculate_ComplexZernikeMoments (const reco::BasicCluster &passedCluster)
void	Calculate_Covariances (const reco::BasicCluster &passedCluster, const EcalRecHitCollection hits, const CaloSubdetectorGeometry geometry)
void	Calculate_e2x2 ()
void	Calculate_e2x5Bottom ()
void	Calculate_e2x5Left ()
void	Calculate_e2x5Right ()
void	Calculate_e2x5Top ()
void	Calculate_e3x2 ()
void	Calculate_e3x3 ()
void	Calculate_e4x4 ()
void	Calculate_e5x5 ()
void	Calculate_EnergyDepTopology (const reco::BasicCluster &passedCluster, const EcalRecHitCollection hits, const CaloSubdetectorGeometry geometry, bool logW=true)
void	Calculate_lat (const reco::BasicCluster &passedCluster)
void	Calculate_Polynomials (double rho)
void	Calculate_TopEnergy (const reco::BasicCluster &passedCluster, const EcalRecHitCollection *hits)
void	Create_Map (const EcalRecHitCollection hits, const CaloSubdetectorTopology topology)
double	f00 (double r)
double	f11 (double r)
double	f20 (double r)
double	f22 (double r)
double	f31 (double r)
double	f33 (double r)
double	f40 (double r)
double	f42 (double r)
double	f44 (double r)
double	f51 (double r)
double	f53 (double r)
double	f55 (double r)
double	factorial (int n) const
double	fast_AbsZernikeMoment (const reco::BasicCluster &passedCluster, int n, int m, double R0)
Private Attributes
double	A20_
double	A42_
double	covEtaEta_
double	covEtaPhi_
double	covPhiPhi_
double	e2nd_
DetId	e2ndId_
double	e2x2_
int	e2x2_Diagonal_X_
int	e2x2_Diagonal_Y_
double	e2x5Bottom_
double	e2x5Left_
double	e2x5Right_
double	e2x5Top_
double	e3x2_
double	e3x2Ratio_
double	e3x3_
double	e4x4_
double	e5x5_
double	eMax_
DetId	eMaxId_
std::vector< double >	energyBasketFractionEta_
std::vector< double >	energyBasketFractionPhi_
std::vector < EcalClusterEnergyDeposition >	energyDistribution_
std::pair< DetId, double >	energyMap_ [5][5]
double	etaLat_
std::vector< double >	fcn_
double	lat_
edm::ParameterSet	parameterSet_
double	phiLat_

Detailed Description

calculates and creates a ClusterShape object

Author:: Michael A. Balazs, UVa

Version:

Id:: ClusterShapeAlgo.h,v 1.17 2010/11/16 15:09:27 argiro Exp

calculates and creates a ClusterShape object

Author:: Michael A. Balazs, UVa

Version:

Id:: SuperClusterShapeAlgo.h,v 1.1 2008/02/11 13:08:05 kkaadze Exp

Definition at line 37 of file ClusterShapeAlgo.h.

Member Enumeration Documentation

anonymous enum [private]

Enumerator:

Eta
Phi

Definition at line 110 of file ClusterShapeAlgo.h.

{ Eta, Phi };

Constructor & Destructor Documentation

ClusterShapeAlgo::ClusterShapeAlgo ( const edm::ParameterSet & par )

Definition at line 19 of file ClusterShapeAlgo.cc.

                                                             : 
  parameterSet_(par) {}

ClusterShapeAlgo::ClusterShapeAlgo ( ) [inline]

Definition at line 42 of file ClusterShapeAlgo.h.

{ };

Member Function Documentation

double ClusterShapeAlgo::absZernikeMoment	(	const reco::BasicCluster &	passedCluster,
		int	n,
		int	m,
		double	R0 = `6.6`
	)		`[private]`

Definition at line 500 of file ClusterShapeAlgo.cc.

References calc_AbsZernikeMoment(), and fast_AbsZernikeMoment().

Referenced by Calculate_ComplexZernikeMoments().

                                                                   {
  // 1. Check if n,m are correctly
  if ((m>n) || ((n-m)%2 != 0) || (n<0) || (m<0)) return -1;

  // 2. Check if n,R0 are within validity Range :
  // n>20 or R0<2.19cm  just makes no sense !
  if ((n>20) || (R0<=2.19)) return -1;
  if (n<=5) return fast_AbsZernikeMoment(passedCluster,n,m,R0);
  else return calc_AbsZernikeMoment(passedCluster,n,m,R0);
}

double ClusterShapeAlgo::calc_AbsZernikeMoment	(	const reco::BasicCluster &	passedCluster,
		int	n,
		int	m,
		double	R0
	)		`[private]`

Definition at line 561 of file ClusterShapeAlgo.cc.

References funct::cos(), alignCSCRings::e, energyDistribution_, factorial(), i, m, phi, funct::pow(), alignCSCRings::r, query::result, alignCSCRings::s, funct::sin(), and mathSSE::sqrt().

Referenced by absZernikeMoment().

                                                                        {
  double r,ph,e,Re=0,Im=0,f_nm,result;
  double TotalEnergy = passedCluster.energy();
  std::vector< std::pair<DetId, float> > clusterDetIds = passedCluster.hitsAndFractions();
  int clusterSize=energyDistribution_.size();
  if(clusterSize<3) return 0.0;

  for (int i=0; i<clusterSize; i++)
    { 
      r = energyDistribution_[i].r / R0;
      if (r<1) {
        ph = (energyDistribution_[i]).phi;
        e = energyDistribution_[i].deposited_energy;
        f_nm=0;
        for (int s=0; s<=(n-m)/2; s++) {
          if (s%2==0)
            { 
              f_nm = f_nm + factorial(n-s)*pow(r,(double) (n-2*s))/(factorial(s)*factorial((n+m)/2-s)*factorial((n-m)/2-s));
            }else {
              f_nm = f_nm - factorial(n-s)*pow(r,(double) (n-2*s))/(factorial(s)*factorial((n+m)/2-s)*factorial((n-m)/2-s));
            }
        }
        Re = Re + e/TotalEnergy * f_nm * cos( (double) m*ph);
        Im = Im - e/TotalEnergy * f_nm * sin( (double) m*ph);
      }
    }
  result = sqrt(Re*Re+Im*Im);

  return result;
}

reco::ClusterShape ClusterShapeAlgo::Calculate	(	const reco::BasicCluster &	passedCluster,
		const EcalRecHitCollection *	hits,
		const CaloSubdetectorGeometry *	geometry,
		const CaloSubdetectorTopology *	topology
	)

Definition at line 22 of file ClusterShapeAlgo.cc.

References A20_, A42_, Calculate_2ndEnergy(), Calculate_BarrelBasketEnergyFraction(), Calculate_ComplexZernikeMoments(), Calculate_Covariances(), Calculate_e2x2(), Calculate_e2x5Bottom(), Calculate_e2x5Left(), Calculate_e2x5Right(), Calculate_e2x5Top(), Calculate_e3x2(), Calculate_e3x3(), Calculate_e4x4(), Calculate_e5x5(), Calculate_EnergyDepTopology(), Calculate_lat(), Calculate_TopEnergy(), covEtaEta_, covEtaPhi_, covPhiPhi_, Create_Map(), e2nd_, e2ndId_, e2x2_, e2x5Bottom_, e2x5Left_, e2x5Right_, e2x5Top_, e3x2_, e3x2Ratio_, e3x3_, e4x4_, e5x5_, eMax_, eMaxId_, energyBasketFractionEta_, energyBasketFractionPhi_, Eta, etaLat_, lat_, Phi, and phiLat_.

Referenced by CosmicClusterProducer::clusterizeECALPart(), IslandClusterProducer::clusterizeECALPart(), and Pi0FixedMassWindowCalibration::duringLoop().

{
  Calculate_TopEnergy(passedCluster,hits);
  Calculate_2ndEnergy(passedCluster,hits);
  Create_Map(hits,topology);
  Calculate_e2x2();
  Calculate_e3x2();
  Calculate_e3x3();
  Calculate_e4x4();
  Calculate_e5x5();
  Calculate_e2x5Right();
  Calculate_e2x5Left();
  Calculate_e2x5Top();
  Calculate_e2x5Bottom();
  Calculate_Covariances(passedCluster,hits,geometry);
  Calculate_BarrelBasketEnergyFraction(passedCluster,hits, Eta, geometry);
  Calculate_BarrelBasketEnergyFraction(passedCluster,hits, Phi, geometry);
  Calculate_EnergyDepTopology (passedCluster,hits,geometry,true) ;
  Calculate_lat(passedCluster);
  Calculate_ComplexZernikeMoments(passedCluster);

  return reco::ClusterShape(covEtaEta_, covEtaPhi_, covPhiPhi_, eMax_, eMaxId_,
                            e2nd_, e2ndId_, e2x2_, e3x2_, e3x3_,e4x4_, e5x5_,
                            e2x5Right_, e2x5Left_, e2x5Top_, e2x5Bottom_,
                            e3x2Ratio_, lat_, etaLat_, phiLat_, A20_, A42_,
                            energyBasketFractionEta_, energyBasketFractionPhi_);
}

void ClusterShapeAlgo::Calculate_2ndEnergy	(	const reco::BasicCluster &	passedCluster,
		const EcalRecHitCollection *	hits
	)		`[private]`

Definition at line 82 of file ClusterShapeAlgo.cc.

References e2nd_, e2ndId_, eMaxId_, edm::SortedCollection< T, SORT >::end(), CaloRecHit::energy(), edm::SortedCollection< T, SORT >::find(), and EcalRecHit::id().

Referenced by Calculate().

{
  double e2nd=0;
  DetId e2ndId(0);

  std::vector< std::pair<DetId, float> > clusterDetIds = passedCluster.hitsAndFractions();
  std::vector< std::pair<DetId, float> >::iterator posCurrent;

  EcalRecHit testEcalRecHit;

  for(posCurrent = clusterDetIds.begin(); posCurrent != clusterDetIds.end(); posCurrent++)
  { 
    if (( (*posCurrent).first != DetId(0)) && (hits->find( (*posCurrent).first ) != hits->end()))
    {
      EcalRecHitCollection::const_iterator itt = hits->find( (*posCurrent).first );
      testEcalRecHit = *itt;

      if(testEcalRecHit.energy() * (*posCurrent).second > e2nd && testEcalRecHit.id() != eMaxId_)
      {
        e2nd = testEcalRecHit.energy() * (*posCurrent).second;
        e2ndId = testEcalRecHit.id();
      }
    }
  }

  e2nd_ = e2nd;
  e2ndId_ = e2ndId;
}

void ClusterShapeAlgo::Calculate_BarrelBasketEnergyFraction	(	const reco::BasicCluster &	passedCluster,
		const EcalRecHitCollection *	hits,
		const int	EtaPhi,
		const CaloSubdetectorGeometry *	geometry
	)		`[private]`

Definition at line 384 of file ClusterShapeAlgo.cc.

References abs, EcalBarrel, relval_parameters_module::energy, energyBasketFractionEta_, energyBasketFractionPhi_, Eta, edm::SortedCollection< T, SORT >::find(), first, i, EBDetId::ieta(), EBDetId::iphi(), EBDetId::MAX_IPHI, EBDetId::MIN_IPHI, and Phi.

Referenced by Calculate().

{
  if(  (hits!=0) && ( ((*hits)[0]).id().subdetId() != EcalBarrel )  ) {
     //std::cout << "No basket correction for endacap!" << std::endl;
     return;
  }

  std::map<int,double> indexedBasketEnergy;
  std::vector< std::pair<DetId, float> > clusterDetIds = passedCluster.hitsAndFractions();
  const EcalBarrelGeometry* subDetGeometry = (const EcalBarrelGeometry*) geometry;

  for(std::vector< std::pair<DetId, float> >::iterator posCurrent = clusterDetIds.begin(); posCurrent != clusterDetIds.end(); posCurrent++)
  {
    int basketIndex = 999;

    if(EtaPhi == Eta) {
      int unsignedIEta = abs(EBDetId( (*posCurrent).first ).ieta());
      std::vector<int> etaBasketSize = subDetGeometry->getEtaBaskets();

      for(unsigned int i = 0; i < etaBasketSize.size(); i++) {
        unsignedIEta -= etaBasketSize[i];
        if(unsignedIEta - 1 < 0)
        {
          basketIndex = i;
          break;
        }
      }
      basketIndex = (basketIndex+1)*(EBDetId( (*posCurrent).first ).ieta() > 0 ? 1 : -1);

    } else if(EtaPhi == Phi) {
      int halfNumBasketInPhi = (EBDetId::MAX_IPHI - EBDetId::MIN_IPHI + 1)/subDetGeometry->getBasketSizeInPhi()/2;

      basketIndex = (EBDetId( (*posCurrent).first ).iphi() - 1)/subDetGeometry->getBasketSizeInPhi()
                  - (EBDetId( (clusterDetIds[0]).first ).iphi() - 1)/subDetGeometry->getBasketSizeInPhi();

      if(basketIndex >= halfNumBasketInPhi)             basketIndex -= 2*halfNumBasketInPhi;
      else if(basketIndex <  -1*halfNumBasketInPhi)     basketIndex += 2*halfNumBasketInPhi;

    } else throw(std::runtime_error("\n\nOh No! Calculate_BarrelBasketEnergyFraction called on invalid index.\n\n"));

    indexedBasketEnergy[basketIndex] += (hits->find( (*posCurrent).first ))->energy();
  }

  std::vector<double> energyFraction;
  for(std::map<int,double>::iterator posCurrent = indexedBasketEnergy.begin(); posCurrent != indexedBasketEnergy.end(); posCurrent++)
  {
    energyFraction.push_back(indexedBasketEnergy[posCurrent->first]/passedCluster.energy());
  }

  switch(EtaPhi)
  {
    case Eta: energyBasketFractionEta_ = energyFraction; break;
    case Phi: energyBasketFractionPhi_ = energyFraction; break;
  }

}

void ClusterShapeAlgo::Calculate_ComplexZernikeMoments ( const reco::BasicCluster & passedCluster ) [private]

Definition at line 493 of file ClusterShapeAlgo.cc.

References A20_, A42_, and absZernikeMoment().

Referenced by Calculate().

                                                                                            {
  // Calculate only the moments which go into the default cluster shape
  // (moments with m>=2 are the only sensitive to azimuthal shape)
  A20_ = absZernikeMoment(passedCluster,2,0);
  A42_ = absZernikeMoment(passedCluster,4,2);
}

void ClusterShapeAlgo::Calculate_Covariances	(	const reco::BasicCluster &	passedCluster,
		const EcalRecHitCollection *	hits,
		const CaloSubdetectorGeometry *	geometry
	)		`[private]`

Definition at line 312 of file ClusterShapeAlgo.cc.

References covEtaEta_, covEtaPhi_, covPhiPhi_, RecoTauValidation_cfi::denominator, dPhi(), e5x5_, energyMap_, PV3DBase< T, PVType, FrameType >::eta(), CaloSubdetectorGeometry::getGeometry(), edm::ParameterSet::getParameter(), CaloCellGeometry::getPosition(), i, j, funct::log(), max(), parameterSet_, PV3DBase< T, PVType, FrameType >::phi(), Geom::Phi< T >::phi(), Geom::pi(), position, edm::second(), Geom::twoPi(), w(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by Calculate().

{
  if (e5x5_ > 0.)
    {
      double w0_ = parameterSet_.getParameter<double>("W0");
      
      
      // first find energy-weighted mean position - doing it when filling the energy map might save time
      math::XYZVector meanPosition(0.0, 0.0, 0.0);
      for (int i = 0; i < 5; ++i)
        {
          for (int j = 0; j < 5; ++j)
            {
              DetId id = energyMap_[i][j].first;
              if (id != DetId(0))
                {
                  GlobalPoint positionGP = geometry->getGeometry(id)->getPosition();
                  math::XYZVector position(positionGP.x(),positionGP.y(),positionGP.z());
                  meanPosition = meanPosition + energyMap_[i][j].second * position;
                }
            }
        }
      
      meanPosition /= e5x5_;
      
      // now we can calculate the covariances
      double numeratorEtaEta = 0;
      double numeratorEtaPhi = 0;
      double numeratorPhiPhi = 0;
      double denominator     = 0;
      
      for (int i = 0; i < 5; ++i)
        {
          for (int j = 0; j < 5; ++j)
            {
              DetId id = energyMap_[i][j].first;
              if (id != DetId(0))
                {
                  GlobalPoint position = geometry->getGeometry(id)->getPosition();
                  
                  double dPhi = position.phi() - meanPosition.phi();
                  if (dPhi > + Geom::pi()) { dPhi = Geom::twoPi() - dPhi; }
                  if (dPhi < - Geom::pi()) { dPhi = Geom::twoPi() + dPhi; }
                  
                  double dEta = position.eta() - meanPosition.eta();
                  double w = 0.;
                  if ( energyMap_[i][j].second > 0.)
                    w = std::max(0.0, w0_ + log( energyMap_[i][j].second / e5x5_));
                  
                  denominator += w;
                  numeratorEtaEta += w * dEta * dEta;
                  numeratorEtaPhi += w * dEta * dPhi;
                  numeratorPhiPhi += w * dPhi * dPhi;
                }
            }
        }
      
      covEtaEta_ = numeratorEtaEta / denominator;
      covEtaPhi_ = numeratorEtaPhi / denominator;
      covPhiPhi_ = numeratorPhiPhi / denominator;
    }
  else 
    {
      // Warn the user if there was no energy in the cells and return zeroes.
      //       std::cout << "\ClusterShapeAlgo::Calculate_Covariances:  no energy in supplied cells.\n";
      covEtaEta_ = 0;
      covEtaPhi_ = 0;
      covPhiPhi_ = 0;
    }
}

void ClusterShapeAlgo::Calculate_e2x2 ( ) [private]

Definition at line 133 of file ClusterShapeAlgo.cc.

References e2x2_, e2x2_Diagonal_X_, e2x2_Diagonal_Y_, and energyMap_.

Referenced by Calculate().

{
  double e2x2Max = 0;
  double e2x2Test = 0;

  int deltaX=0, deltaY=0;

  for(int corner = 0; corner < 4; corner++)
  {
    switch(corner)
    {
      case 0: deltaX = -1; deltaY = -1; break;
      case 1: deltaX = -1; deltaY =  1; break;
      case 2: deltaX =  1; deltaY = -1; break;
      case 3: deltaX =  1; deltaY =  1; break;
    }

    e2x2Test  = energyMap_[2][2].second;
    e2x2Test += energyMap_[2+deltaY][2].second;
    e2x2Test += energyMap_[2][2+deltaX].second;
    e2x2Test += energyMap_[2+deltaY][2+deltaX].second;

    if(e2x2Test > e2x2Max)
    {
                        e2x2Max = e2x2Test;
                        e2x2_Diagonal_X_ = 2+deltaX;
                        e2x2_Diagonal_Y_ = 2+deltaY;
    }
  }

  e2x2_ = e2x2Max;

}

void ClusterShapeAlgo::Calculate_e2x5Bottom ( ) [private]

Definition at line 288 of file ClusterShapeAlgo.cc.

References e2x5Bottom_, energyMap_, i, and j.

Referenced by Calculate().

{
double e2x5B=0.0;
  for(int i = 0; i <= 4; i++){
    for(int j = 0; j <= 4; j++){

      if(i>2){e2x5B +=energyMap_[i][j].second;}
    }
  }
  e2x5Bottom_=e2x5B;
}

void ClusterShapeAlgo::Calculate_e2x5Left ( ) [private]

Definition at line 277 of file ClusterShapeAlgo.cc.

References e2x5Left_, energyMap_, i, and j.

Referenced by Calculate().

{
double e2x5L=0.0;
  for(int i = 0; i <= 4; i++){
    for(int j = 0; j <= 4; j++){
      if(j<2){e2x5L +=energyMap_[i][j].second;}
    }
  }
  e2x5Left_=e2x5L;
}

void ClusterShapeAlgo::Calculate_e2x5Right ( ) [private]

Definition at line 266 of file ClusterShapeAlgo.cc.

References e2x5Right_, energyMap_, i, and j.

Referenced by Calculate().

{
double e2x5R=0.0;
  for(int i = 0; i <= 4; i++){
    for(int j = 0; j <= 4; j++){
      if(j>2){e2x5R +=energyMap_[i][j].second;}
    }
  }
  e2x5Right_=e2x5R;
}

void ClusterShapeAlgo::Calculate_e2x5Top ( ) [private]

Definition at line 300 of file ClusterShapeAlgo.cc.

References e2x5Top_, energyMap_, i, and j.

Referenced by Calculate().

{
double e2x5T=0.0;
  for(int i = 0; i <= 4; i++){
    for(int j = 0; j <= 4; j++){

      if(i<2){e2x5T +=energyMap_[i][j].second;}
    }
  }
  e2x5Top_=e2x5T;
}

void ClusterShapeAlgo::Calculate_e3x2 ( ) [private]

Definition at line 167 of file ClusterShapeAlgo.cc.

References e3x2_, e3x2Ratio_, energyMap_, and edm::second().

Referenced by Calculate().

{
  double e3x2 = 0.0;
  double e3x2Ratio = 0.0, e3x2RatioNumerator = 0.0, e3x2RatioDenominator = 0.0;

  int e2ndX = 2, e2ndY = 2;
  int deltaY = 0, deltaX = 0;

  double nextEnergy = -999;
  int nextEneryDirection = -1;

  for(int cardinalDirection = 0; cardinalDirection < 4; cardinalDirection++)
  {
    switch(cardinalDirection)
    {
      case 0: deltaX = -1; deltaY =  0; break;
      case 1: deltaX =  1; deltaY =  0; break;
      case 2: deltaX =  0; deltaY = -1; break;
      case 3: deltaX =  0; deltaY =  1; break;
    }
   
    if(energyMap_[2+deltaY][2+deltaX].second >= nextEnergy)
    {
        nextEnergy = energyMap_[2+deltaY][2+deltaX].second;
        nextEneryDirection = cardinalDirection;
       
        e2ndX = 2+deltaX;
        e2ndY = 2+deltaY;
    }
  }
 
  switch(nextEneryDirection)
  {
    case 0: ;
    case 1: deltaX = 0; deltaY = 1; break;
    case 2: ;
    case 3: deltaX = 1; deltaY = 0; break;
  }

  for(int sign = -1; sign <= 1; sign++)
      e3x2 += (energyMap_[2+deltaY*sign][2+deltaX*sign].second + energyMap_[e2ndY+deltaY*sign][e2ndX+deltaX*sign].second);
 
  e3x2RatioNumerator   = energyMap_[e2ndY+deltaY][e2ndX+deltaX].second + energyMap_[e2ndY-deltaY][e2ndX-deltaX].second;
  e3x2RatioDenominator = 0.5 + energyMap_[2+deltaY][2+deltaX].second + energyMap_[2-deltaY][2-deltaX].second;
  e3x2Ratio = e3x2RatioNumerator / e3x2RatioDenominator;

  e3x2_ = e3x2;
  e3x2Ratio_ = e3x2Ratio;
}

void ClusterShapeAlgo::Calculate_e3x3 ( ) [private]

Definition at line 217 of file ClusterShapeAlgo.cc.

References e3x3_, energyMap_, i, j, and edm::second().

Referenced by Calculate().

{
  double e3x3=0;

  for(int i = 1; i <= 3; i++)
    for(int j = 1; j <= 3; j++)
      e3x3 += energyMap_[j][i].second;

  e3x3_ = e3x3;

}

void ClusterShapeAlgo::Calculate_e4x4 ( ) [private]

Definition at line 229 of file ClusterShapeAlgo.cc.

References e2x2_Diagonal_X_, e2x2_Diagonal_Y_, e4x4_, energyMap_, i, j, and edm::second().

Referenced by Calculate().

{
  double e4x4=0;
        
  int startX=-1, startY=-1;

        switch(e2x2_Diagonal_X_)
        {
                case 1: startX = 0; break;
                case 3: startX = 1; break;
        }

        switch(e2x2_Diagonal_Y_)
        {
                case 1: startY = 0; break;
                case 3: startY = 1; break;
        }

  for(int i = startX; i <= startX+3; i++)
          for(int j = startY; j <= startY+3; j++)
          e4x4 += energyMap_[j][i].second;

  e4x4_ = e4x4;
}

void ClusterShapeAlgo::Calculate_e5x5 ( ) [private]

Definition at line 254 of file ClusterShapeAlgo.cc.

References e5x5_, energyMap_, i, j, and edm::second().

Referenced by Calculate().

{
  double e5x5=0;

  for(int i = 0; i <= 4; i++)
    for(int j = 0; j <= 4; j++)
      e5x5 += energyMap_[j][i].second;

  e5x5_ = e5x5;

}

void ClusterShapeAlgo::Calculate_EnergyDepTopology	(	const reco::BasicCluster &	passedCluster,
		const EcalRecHitCollection *	hits,
		const CaloSubdetectorGeometry *	geometry,
		bool	logW = `true`
	)		`[private]`

Definition at line 593 of file ClusterShapeAlgo.cc.

References EcalClusterEnergyDeposition::deposited_energy, diffTreeTool::diff, edm::SortedCollection< T, SORT >::end(), CaloRecHit::energy(), energyDistribution_, edm::SortedCollection< T, SORT >::find(), CaloSubdetectorGeometry::getGeometry(), edm::ParameterSet::getParameter(), CaloCellGeometry::getPosition(), funct::log(), LogDebug, M_PI, max(), parameterSet_, EcalClusterEnergyDeposition::phi, EcalClusterEnergyDeposition::r, CommonMethods::weight(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by Calculate().

                                                               {
  // resets the energy distribution
  energyDistribution_.clear();

  // init a map of the energy deposition centered on the
  // cluster centroid. This is for momenta calculation only.
  CLHEP::Hep3Vector clVect(passedCluster.position().x(),
                           passedCluster.position().y(),
                           passedCluster.position().z());
  CLHEP::Hep3Vector clDir(clVect);
  clDir*=1.0/clDir.mag();
  // in the transverse plane, axis perpendicular to clusterDir
  CLHEP::Hep3Vector theta_axis(clDir.y(),-clDir.x(),0.0);
  theta_axis *= 1.0/theta_axis.mag();
  CLHEP::Hep3Vector phi_axis = theta_axis.cross(clDir);

  std::vector< std::pair<DetId, float> > clusterDetIds = passedCluster.hitsAndFractions();

  EcalClusterEnergyDeposition clEdep;
  EcalRecHit testEcalRecHit;
  std::vector< std::pair<DetId, float> >::iterator posCurrent;
  // loop over crystals
  for(posCurrent=clusterDetIds.begin(); posCurrent!=clusterDetIds.end(); ++posCurrent) {
    EcalRecHitCollection::const_iterator itt = hits->find( (*posCurrent).first );
    testEcalRecHit=*itt;

    if(( (*posCurrent).first != DetId(0)) && (hits->find( (*posCurrent).first ) != hits->end())) {
      clEdep.deposited_energy = testEcalRecHit.energy();

      // if logarithmic weight is requested, apply cut on minimum energy of the recHit
      if(logW) {
        double w0_ = parameterSet_.getParameter<double>("W0");

        if ( clEdep.deposited_energy == 0 ) {
          LogDebug("ClusterShapeAlgo") << "Crystal has zero energy; skipping... ";
          continue;
        }
        double weight = std::max(0.0, w0_ + log(fabs(clEdep.deposited_energy)/passedCluster.energy()) );
        if(weight==0) {
          LogDebug("ClusterShapeAlgo") << "Crystal has insufficient energy: E = " 
                                       << clEdep.deposited_energy << " GeV; skipping... ";
          continue;
        }
        else LogDebug("ClusterShapeAlgo") << "===> got crystal. Energy = " << clEdep.deposited_energy << " GeV. ";
      }
      DetId id_ = (*posCurrent).first;
      const CaloCellGeometry *this_cell = geometry->getGeometry(id_);
      GlobalPoint cellPos = this_cell->getPosition();
      CLHEP::Hep3Vector gblPos (cellPos.x(),cellPos.y(),cellPos.z()); //surface position?
      // Evaluate the distance from the cluster centroid
      CLHEP::Hep3Vector diff = gblPos - clVect;
      // Important: for the moment calculation, only the "lateral distance" is important
      // "lateral distance" r_i = distance of the digi position from the axis Origin-Cluster Center
      // ---> subtract the projection on clDir
      CLHEP::Hep3Vector DigiVect = diff - diff.dot(clDir)*clDir;
      clEdep.r = DigiVect.mag();
      LogDebug("ClusterShapeAlgo") << "E = " << clEdep.deposited_energy
                                   << "\tdiff = " << diff.mag()
                                   << "\tr = " << clEdep.r;
      clEdep.phi = DigiVect.angle(theta_axis);
      if(DigiVect.dot(phi_axis)<0) clEdep.phi = 2*M_PI - clEdep.phi;
      energyDistribution_.push_back(clEdep);
    }
  } 
}

void ClusterShapeAlgo::Calculate_lat ( const reco::BasicCluster & passedCluster ) [private]

Definition at line 444 of file ClusterShapeAlgo.cc.

References funct::cos(), energyDistribution_, etaLat_, i, lat_, n, phi, phiLat_, alignCSCRings::r, funct::sin(), and tmp.

Referenced by Calculate().

                                                                          {

  double r,redmoment=0;
  double phiRedmoment = 0 ;
  double etaRedmoment = 0 ;
  int n,n1,n2,tmp;
  int clusterSize=energyDistribution_.size();
  if (clusterSize<3) {
    etaLat_ = 0.0 ; 
    lat_ = 0.0;
    return; 
  }
  
  n1=0; n2=1;
  if (energyDistribution_[1].deposited_energy > 
      energyDistribution_[0].deposited_energy) 
    {
      tmp=n2; n2=n1; n1=tmp;
    }
  for (int i=2; i<clusterSize; i++) {
    n=i;
    if (energyDistribution_[i].deposited_energy > 
        energyDistribution_[n1].deposited_energy) 
      {
        tmp = n2;
        n2 = n1; n1 = i; n=tmp;
      } else {
        if (energyDistribution_[i].deposited_energy > 
            energyDistribution_[n2].deposited_energy) 
          {
            tmp=n2; n2=i; n=tmp;
          }
      }

    r = energyDistribution_[n].r;
    redmoment += r*r* energyDistribution_[n].deposited_energy;
    double rphi = r * cos (energyDistribution_[n].phi) ;
    phiRedmoment += rphi * rphi * energyDistribution_[n].deposited_energy;
    double reta = r * sin (energyDistribution_[n].phi) ;
    etaRedmoment += reta * reta * energyDistribution_[n].deposited_energy;
  } 
  double e1 = energyDistribution_[n1].deposited_energy;
  double e2 = energyDistribution_[n2].deposited_energy;
  
  lat_ = redmoment/(redmoment+2.19*2.19*(e1+e2));
  phiLat_ = phiRedmoment/(phiRedmoment+2.19*2.19*(e1+e2));
  etaLat_ = etaRedmoment/(etaRedmoment+2.19*2.19*(e1+e2));
}

void ClusterShapeAlgo::Calculate_Polynomials ( double rho ) [private]

Definition at line 662 of file ClusterShapeAlgo.cc.

References f00(), f11(), f20(), f22(), f31(), f33(), f40(), f42(), f44(), f51(), f53(), f55(), and fcn_.

Referenced by fast_AbsZernikeMoment().

                                                       {
  fcn_.push_back(f00(rho));
  fcn_.push_back(f11(rho));
  fcn_.push_back(f20(rho));
  fcn_.push_back(f31(rho));
  fcn_.push_back(f22(rho));
  fcn_.push_back(f33(rho));
  fcn_.push_back(f40(rho));
  fcn_.push_back(f51(rho));
  fcn_.push_back(f42(rho));
  fcn_.push_back(f53(rho));
  fcn_.push_back(f44(rho));
  fcn_.push_back(f55(rho));
}

void ClusterShapeAlgo::Calculate_TopEnergy	(	const reco::BasicCluster &	passedCluster,
		const EcalRecHitCollection *	hits
	)		`[private]`

Definition at line 53 of file ClusterShapeAlgo.cc.

References jptDQMConfig_cff::eMax, eMax_, eMaxId_, edm::SortedCollection< T, SORT >::end(), CaloRecHit::energy(), edm::SortedCollection< T, SORT >::find(), and EcalRecHit::id().

Referenced by Calculate().

{
  double eMax=0;
  DetId eMaxId(0);

  std::vector< std::pair<DetId, float> > clusterDetIds = passedCluster.hitsAndFractions();
  std::vector< std::pair<DetId, float> >::iterator posCurrent;

  EcalRecHit testEcalRecHit;

  for(posCurrent = clusterDetIds.begin(); posCurrent != clusterDetIds.end(); posCurrent++)
  {
    if (((*posCurrent).first != DetId(0)) && (hits->find((*posCurrent).first) != hits->end()))
    {
      EcalRecHitCollection::const_iterator itt = hits->find((*posCurrent).first);
      testEcalRecHit = *itt;

      if(testEcalRecHit.energy() * (*posCurrent).second > eMax)
      {
        eMax = testEcalRecHit.energy() * (*posCurrent).second;
        eMaxId = testEcalRecHit.id();
      }
    }
  }

  eMax_ = eMax;
  eMaxId_ = eMaxId;
}

void ClusterShapeAlgo::Create_Map	(	const EcalRecHitCollection *	hits,
		const CaloSubdetectorTopology *	topology
	)		`[private]`

Definition at line 111 of file ClusterShapeAlgo.cc.

References eMaxId_, edm::SortedCollection< T, SORT >::end(), CaloRecHit::energy(), energyMap_, edm::SortedCollection< T, SORT >::find(), CaloNavigator< T >::home(), EcalRecHit::id(), CaloNavigator< T >::offsetBy(), x, and detailsBasic3DVector::y.

Referenced by Calculate().

{
  EcalRecHit tempEcalRecHit;
  CaloNavigator<DetId> posCurrent = CaloNavigator<DetId>(eMaxId_,topology );

  for(int x = 0; x < 5; x++)
    for(int y = 0; y < 5; y++)
    {
      posCurrent.home();
      posCurrent.offsetBy(-2+x,-2+y);

      if((*posCurrent != DetId(0)) && (hits->find(*posCurrent) != hits->end()))
      {
                                EcalRecHitCollection::const_iterator itt = hits->find(*posCurrent);
                                tempEcalRecHit = *itt;
                                energyMap_[y][x] = std::make_pair(tempEcalRecHit.id(),tempEcalRecHit.energy()); 
      }
      else
                                energyMap_[y][x] = std::make_pair(DetId(0), 0);  
    }
}

double ClusterShapeAlgo::f00 ( double r ) [private]

Definition at line 512 of file ClusterShapeAlgo.cc.

Referenced by Calculate_Polynomials().

{ return 1; }

double ClusterShapeAlgo::f11 ( double r ) [private]

Definition at line 514 of file ClusterShapeAlgo.cc.

References alignCSCRings::r.

Referenced by Calculate_Polynomials().

{ return r; }

double ClusterShapeAlgo::f20 ( double r ) [private]

Definition at line 516 of file ClusterShapeAlgo.cc.

Referenced by Calculate_Polynomials().

{ return 2.0*r*r-1.0; }

double ClusterShapeAlgo::f22 ( double r ) [private]

Definition at line 518 of file ClusterShapeAlgo.cc.

References alignCSCRings::r.

Referenced by Calculate_Polynomials().

{ return r*r; }

double ClusterShapeAlgo::f31 ( double r ) [private]

Definition at line 520 of file ClusterShapeAlgo.cc.

References alignCSCRings::r.

Referenced by Calculate_Polynomials().

{ return 3.0*r*r*r - 2.0*r; }

double ClusterShapeAlgo::f33 ( double r ) [private]

Definition at line 522 of file ClusterShapeAlgo.cc.

References alignCSCRings::r.

Referenced by Calculate_Polynomials().

{ return r*r*r; }

double ClusterShapeAlgo::f40 ( double r ) [private]

Definition at line 524 of file ClusterShapeAlgo.cc.

Referenced by Calculate_Polynomials().

{ return 6.0*r*r*r*r-6.0*r*r+1.0; }

double ClusterShapeAlgo::f42 ( double r ) [private]

Definition at line 526 of file ClusterShapeAlgo.cc.

References alignCSCRings::r.

Referenced by Calculate_Polynomials().

{ return 4.0*r*r*r*r-3.0*r*r; }

double ClusterShapeAlgo::f44 ( double r ) [private]

Definition at line 528 of file ClusterShapeAlgo.cc.

References alignCSCRings::r.

Referenced by Calculate_Polynomials().

{ return r*r*r*r; }

double ClusterShapeAlgo::f51 ( double r ) [private]

Definition at line 530 of file ClusterShapeAlgo.cc.

References funct::pow(), and alignCSCRings::r.

Referenced by Calculate_Polynomials().

{ return 10.0*pow(r,5)-12.0*pow(r,3)+3.0*r; }

double ClusterShapeAlgo::f53 ( double r ) [private]

Definition at line 532 of file ClusterShapeAlgo.cc.

References funct::pow().

Referenced by Calculate_Polynomials().

{ return 5.0*pow(r,5) - 4.0*pow(r,3); }

double ClusterShapeAlgo::f55 ( double r ) [private]

Definition at line 534 of file ClusterShapeAlgo.cc.

References funct::pow().

Referenced by Calculate_Polynomials().

{ return pow(r,5); }

double ClusterShapeAlgo::factorial ( int n ) const [private]

Definition at line 677 of file ClusterShapeAlgo.cc.

References i, and n.

Referenced by calc_AbsZernikeMoment().

                                              {
  double res=1.0;
  for(int i=2; i<=n; i++) res*=(double) i;
  return res;
}

double ClusterShapeAlgo::fast_AbsZernikeMoment	(	const reco::BasicCluster &	passedCluster,
		int	n,
		int	m,
		double	R0
	)		`[private]`

Definition at line 536 of file ClusterShapeAlgo.cc.

References Calculate_Polynomials(), funct::cos(), alignCSCRings::e, energyDistribution_, fcn_, i, getHLTprescales::index, phi, alignCSCRings::r, query::result, funct::sin(), and mathSSE::sqrt().

Referenced by absZernikeMoment().

                                                                        {
  double r,ph,e,Re=0,Im=0,result;
  double TotalEnergy = passedCluster.energy();
  int index = (n/2)*(n/2)+(n/2)+m;
  int clusterSize=energyDistribution_.size();
  if(clusterSize<3) return 0.0;

  for (int i=0; i<clusterSize; i++)
    { 
      r = energyDistribution_[i].r / R0;
      if (r<1) {
        fcn_.clear();
        Calculate_Polynomials(r);
        ph = (energyDistribution_[i]).phi;
        e = energyDistribution_[i].deposited_energy;
        Re = Re + e/TotalEnergy * fcn_[index] * cos( (double) m * ph);
        Im = Im - e/TotalEnergy * fcn_[index] * sin( (double) m * ph);
      }
    }
  result = sqrt(Re*Re+Im*Im);

  return result;
}