CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Classes | Public Member Functions | Static Public Member Functions | Static Private Member Functions
EcalClusterToolsT< noZS > Class Template Reference

#include <EcalClusterTools.h>

Classes

struct  EcalClusterEnergyDeposition
 

Public Member Functions

 EcalClusterToolsT ()
 
 ~EcalClusterToolsT ()
 

Static Public Member Functions

static Cluster2ndMoments cluster2ndMoments (const reco::BasicCluster &basicCluster, const EcalRecHitCollection &recHits, double phiCorrectionFactor=0.8, double w0=4.7, bool useLogWeights=true)
 
static Cluster2ndMoments cluster2ndMoments (const reco::SuperCluster &superCluster, const EcalRecHitCollection &recHits, double phiCorrectionFactor=0.8, double w0=4.7, bool useLogWeights=true)
 
static Cluster2ndMoments cluster2ndMoments (const std::vector< std::pair< const EcalRecHit *, float > > &RH_ptrs_fracs, double phiCorrectionFactor=0.8, double w0=4.7, bool useLogWeights=true)
 
static std::vector< float > covariances (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, const CaloGeometry *geometry, float w0=4.7)
 
static float e1x3 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e1x5 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e2nd (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
 
static float e2x2 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e2x5Bottom (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e2x5Left (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e2x5Max (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e2x5Right (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e2x5Top (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e3x1 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e3x2 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e3x3 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e4x4 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e5x1 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float e5x5 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float eBottom (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float eLeft (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float eMax (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
 
static std::vector< float > energyBasketFractionEta (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
 
static std::vector< float > energyBasketFractionPhi (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
 
static float eRight (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float eTop (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static float getFraction (const std::vector< std::pair< DetId, float > > &v_id, DetId id)
 
static std::pair< DetId, float > getMaximum (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
 
static std::pair< DetId, float > getMaximum (const std::vector< std::pair< DetId, float > > &v_id, const EcalRecHitCollection *recHits)
 
static std::vector< float > lat (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, bool logW=true, float w0=4.7)
 
static std::vector< float > localCovariances (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, float w0=4.7)
 
static std::vector< DetIdmatrixDetId (const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
 
static float matrixEnergy (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
 
static float recHitEnergy (DetId id, const EcalRecHitCollection *recHits)
 
static std::vector< float > roundnessBarrelSuperClusters (const reco::SuperCluster &superCluster, const EcalRecHitCollection &recHits, int weightedPositionMethod=0, float energyThreshold=0.0)
 
static std::vector< float > roundnessBarrelSuperClustersUserExtended (const reco::SuperCluster &superCluster, const EcalRecHitCollection &recHits, int ieta_delta=0, int iphi_delta=0, float energyRHThresh=0.00000, int weightedPositionMethod=0)
 
static std::vector< float > roundnessSelectedBarrelRecHits (const std::vector< std::pair< const EcalRecHit *, float > > &rhVector, int weightedPositionMethod=0)
 
static std::vector< float > scLocalCovariances (const reco::SuperCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, float w0=4.7)
 
static double zernike20 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, double R0=6.6, bool logW=true, float w0=4.7)
 
static double zernike42 (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, double R0=6.6, bool logW=true, float w0=4.7)
 

Static Private Member Functions

static double absZernikeMoment (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, int n, int m, double R0, bool logW, float w0)
 
static double calc_AbsZernikeMoment (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, int n, int m, double R0, bool logW, float w0)
 
static float computeWeight (float eRH, float energyTotal, int weightedPositionMethod)
 
static int deltaIEta (int seed_ieta, int rh_ieta)
 
static int deltaIPhi (int seed_iphi, int rh_iphi)
 
static double f00 (double r)
 
static double f11 (double r)
 
static double f20 (double r)
 
static double f22 (double r)
 
static double f31 (double r)
 
static double f33 (double r)
 
static double f40 (double r)
 
static double f42 (double r)
 
static double f44 (double r)
 
static double f51 (double r)
 
static double f53 (double r)
 
static double f55 (double r)
 
static double factorial (int n)
 
static double fast_AbsZernikeMoment (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, int n, int m, double R0, bool logW, float w0)
 
static float getDPhiEndcap (const DetId &crysId, float meanX, float meanY)
 
static std::vector
< EcalClusterEnergyDeposition
getEnergyDepTopology (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, bool logW, float w0)
 
static float getIEta (const DetId &id)
 
static float getIPhi (const DetId &id)
 
static float getNormedIX (const DetId &id)
 
static float getNormedIY (const DetId &id)
 
static float getNrCrysDiffInEta (const DetId &crysId, const DetId &orginId)
 
static float getNrCrysDiffInPhi (const DetId &crysId, const DetId &orginId)
 
static std::vector< int > getSeedPosition (const std::vector< std::pair< const EcalRecHit *, float > > &RH_ptrs)
 
static float getSumEnergy (const std::vector< std::pair< const EcalRecHit *, float > > &RH_ptrs_fracs)
 
static std::pair< float, float > mean5x5PositionInLocalCrysCoord (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static std::pair< float, float > mean5x5PositionInXY (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
 
static math::XYZVector meanClusterPosition (const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, const CaloGeometry *geometry)
 

Detailed Description

template<bool noZS>
class EcalClusterToolsT< noZS >

Definition at line 76 of file EcalClusterTools.h.

Constructor & Destructor Documentation

template<bool noZS>
EcalClusterToolsT< noZS >::EcalClusterToolsT ( )
inline

Definition at line 78 of file EcalClusterTools.h.

78 {};
template<bool noZS>
EcalClusterToolsT< noZS >::~EcalClusterToolsT ( )
inline

Definition at line 79 of file EcalClusterTools.h.

79 {};

Member Function Documentation

template<bool noZS>
double EcalClusterToolsT< noZS >::absZernikeMoment ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloGeometry geometry,
int  n,
int  m,
double  R0,
bool  logW,
float  w0 
)
staticprivate

Definition at line 951 of file EcalClusterTools.h.

952 {
953  // 1. Check if n,m are correctly
954  if ((m>n) || ((n-m)%2 != 0) || (n<0) || (m<0)) return -1;
955 
956  // 2. Check if n,R0 are within validity Range :
957  // n>20 or R0<2.19cm just makes no sense !
958  if ((n>20) || (R0<=2.19)) return -1;
959  if (n<=5) return fast_AbsZernikeMoment(cluster, recHits, geometry, n, m, R0, logW, w0 );
960  else return calc_AbsZernikeMoment(cluster, recHits, geometry, n, m, R0, logW, w0 );
961 }
static double calc_AbsZernikeMoment(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, int n, int m, double R0, bool logW, float w0)
static double fast_AbsZernikeMoment(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, int n, int m, double R0, bool logW, float w0)
template<bool noZS>
double EcalClusterToolsT< noZS >::calc_AbsZernikeMoment ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloGeometry geometry,
int  n,
int  m,
double  R0,
bool  logW,
float  w0 
)
staticprivate

Definition at line 1000 of file EcalClusterTools.h.

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

1001 {
1002  double r, ph, e, Re=0, Im=0, f_nm;
1003  double TotalEnergy = cluster.energy();
1004  std::vector<EcalClusterEnergyDeposition> energyDistribution = getEnergyDepTopology( cluster, recHits, geometry, logW, w0 );
1005  int clusterSize=energyDistribution.size();
1006  if(clusterSize<3) return 0.0;
1007 
1008  for (int i = 0; i < clusterSize; ++i)
1009  {
1010  r = energyDistribution[i].r / R0;
1011  if (r < 1) {
1012  ph = energyDistribution[i].phi;
1013  e = energyDistribution[i].deposited_energy;
1014  f_nm = 0;
1015  for (int s=0; s<=(n-m)/2; s++) {
1016  if (s%2==0) {
1017  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));
1018  } else {
1019  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));
1020  }
1021  }
1022  Re = Re + e/TotalEnergy * f_nm * cos( (double) m*ph);
1023  Im = Im - e/TotalEnergy * f_nm * sin( (double) m*ph);
1024  }
1025  }
1026  return sqrt(Re*Re+Im*Im);
1027 }
int i
Definition: DBlmapReader.cc:9
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
static std::vector< EcalClusterEnergyDeposition > getEnergyDepTopology(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, bool logW, float w0)
T sqrt(T t)
Definition: SSEVec.h:48
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
static double factorial(int n)
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
template<bool noZS>
Cluster2ndMoments EcalClusterToolsT< noZS >::cluster2ndMoments ( const reco::BasicCluster basicCluster,
const EcalRecHitCollection recHits,
double  phiCorrectionFactor = 0.8,
double  w0 = 4.7,
bool  useLogWeights = true 
)
static

Definition at line 1236 of file EcalClusterTools.h.

References edm::SortedCollection< T, SORT >::find(), first, i, and edm::second().

Referenced by EcalClusterToolsT< noZS >::cluster2ndMoments().

1236  {
1237 
1238  // for now implemented only for EB:
1239  // if( fabs( basicCluster.eta() ) < 1.479 ) {
1240 
1241  std::vector<std::pair<const EcalRecHit*, float> > RH_ptrs_fracs;
1242 
1243  const std::vector< std::pair<DetId, float> >& myHitsPair = basicCluster.hitsAndFractions();
1244 
1245  for(unsigned int i=0; i<myHitsPair.size(); i++){
1246  //get pointer to recHit object
1247  EcalRecHitCollection::const_iterator myRH = recHits.find(myHitsPair[i].first);
1248  RH_ptrs_fracs.push_back( std::make_pair(&(*myRH) , myHitsPair[i].second) );
1249  }
1250 
1251  return EcalClusterToolsT<noZS>::cluster2ndMoments(RH_ptrs_fracs, phiCorrectionFactor, w0, useLogWeights);
1252 }
int i
Definition: DBlmapReader.cc:9
std::vector< EcalRecHit >::const_iterator const_iterator
U second(std::pair< T, U > const &p)
bool first
Definition: L1TdeRCT.cc:75
static Cluster2ndMoments cluster2ndMoments(const reco::BasicCluster &basicCluster, const EcalRecHitCollection &recHits, double phiCorrectionFactor=0.8, double w0=4.7, bool useLogWeights=true)
iterator find(key_type k)
template<bool noZS>
Cluster2ndMoments EcalClusterToolsT< noZS >::cluster2ndMoments ( const reco::SuperCluster superCluster,
const EcalRecHitCollection recHits,
double  phiCorrectionFactor = 0.8,
double  w0 = 4.7,
bool  useLogWeights = true 
)
static

Definition at line 1255 of file EcalClusterTools.h.

References Cluster2ndMoments::alpha, EcalClusterToolsT< noZS >::cluster2ndMoments(), reco::SuperCluster::seed(), Cluster2ndMoments::sMaj, and Cluster2ndMoments::sMin.

1255  {
1256 
1257  // for now returns second moments of supercluster seed cluster:
1258  Cluster2ndMoments returnMoments;
1259  returnMoments.sMaj = -1.;
1260  returnMoments.sMin = -1.;
1261  returnMoments.alpha = 0.;
1262 
1263  // for now implemented only for EB:
1264  // if( fabs( superCluster.eta() ) < 1.479 ) {
1265  returnMoments = EcalClusterToolsT<noZS>::cluster2ndMoments( *(superCluster.seed()), recHits, phiCorrectionFactor, w0, useLogWeights);
1266  // }
1267 
1268  return returnMoments;
1269 
1270 }
static Cluster2ndMoments cluster2ndMoments(const reco::BasicCluster &basicCluster, const EcalRecHitCollection &recHits, double phiCorrectionFactor=0.8, double w0=4.7, bool useLogWeights=true)
const CaloClusterPtr & seed() const
seed BasicCluster
Definition: SuperCluster.h:66
template<bool noZS>
Cluster2ndMoments EcalClusterToolsT< noZS >::cluster2ndMoments ( const std::vector< std::pair< const EcalRecHit *, float > > &  RH_ptrs_fracs,
double  phiCorrectionFactor = 0.8,
double  w0 = 4.7,
bool  useLogWeights = true 
)
static

Definition at line 1273 of file EcalClusterTools.h.

References funct::abs(), Cluster2ndMoments::alpha, cuy::denominator, EcalRecHit::detid(), EcalBarrel, EcalRecHit::energy(), EcalClusterToolsT< noZS >::getSumEnergy(), i, fff_deleter::log, max(), phi, Cluster2ndMoments::sMaj, Cluster2ndMoments::sMin, mathSSE::sqrt(), and DetId::subdetId().

1273  {
1274 
1275  double mid_eta(0),mid_phi(0),mid_x(0),mid_y(0);
1276 
1277  double Etot = EcalClusterToolsT<noZS>::getSumEnergy( RH_ptrs_fracs );
1278 
1279  double max_phi=-10.;
1280  double min_phi=100.;
1281 
1282 
1283  std::vector<double> etaDetId;
1284  std::vector<double> phiDetId;
1285  std::vector<double> xDetId;
1286  std::vector<double> yDetId;
1287  std::vector<double> wiDetId;
1288 
1289  unsigned int nCry=0;
1290  double denominator=0.;
1291  bool isBarrel(1);
1292 
1293  // loop over rechits and compute weights:
1294  for(std::vector<std::pair<const EcalRecHit*, float> >::const_iterator rhf_ptr = RH_ptrs_fracs.begin(); rhf_ptr != RH_ptrs_fracs.end(); rhf_ptr++){
1295  const EcalRecHit* rh_ptr = rhf_ptr->first;
1296 
1297 
1298  //get iEta, iPhi
1299  double temp_eta(0),temp_phi(0),temp_x(0),temp_y(0);
1300  isBarrel = rh_ptr->detid().subdetId()==EcalBarrel;
1301 
1302  if(isBarrel) {
1303  temp_eta = (getIEta(rh_ptr->detid()) > 0. ? getIEta(rh_ptr->detid()) + 84.5 : getIEta(rh_ptr->detid()) + 85.5);
1304  temp_phi= getIPhi(rh_ptr->detid()) - 0.5;
1305  }
1306  else {
1307  temp_eta = getIEta(rh_ptr->detid());
1308  temp_x = getNormedIX(rh_ptr->detid());
1309  temp_y = getNormedIY(rh_ptr->detid());
1310  }
1311 
1312  double temp_ene=rh_ptr->energy() * (noZS ? 1.0 : rhf_ptr->second);
1313 
1314  double temp_wi=((useLogWeights) ?
1315  std::max(0.0, w0 + std::log( std::abs(temp_ene)/Etot ))
1316  : temp_ene);
1317 
1318 
1319  if(temp_phi>max_phi) max_phi=temp_phi;
1320  if(temp_phi<min_phi) min_phi=temp_phi;
1321  etaDetId.push_back(temp_eta);
1322  phiDetId.push_back(temp_phi);
1323  xDetId.push_back(temp_x);
1324  yDetId.push_back(temp_y);
1325  wiDetId.push_back(temp_wi);
1326  denominator+=temp_wi;
1327  nCry++;
1328  }
1329 
1330  if(isBarrel){
1331  // correct phi wrap-around:
1332  if(max_phi==359.5 && min_phi==0.5){
1333  for(unsigned int i=0; i<nCry; i++){
1334  if(phiDetId[i] - 179. > 0.) phiDetId[i]-=360.;
1335  mid_phi+=phiDetId[i]*wiDetId[i];
1336  mid_eta+=etaDetId[i]*wiDetId[i];
1337  }
1338  } else{
1339  for(unsigned int i=0; i<nCry; i++){
1340  mid_phi+=phiDetId[i]*wiDetId[i];
1341  mid_eta+=etaDetId[i]*wiDetId[i];
1342  }
1343  }
1344  }else{
1345  for(unsigned int i=0; i<nCry; i++){
1346  mid_eta+=etaDetId[i]*wiDetId[i];
1347  mid_x+=xDetId[i]*wiDetId[i];
1348  mid_y+=yDetId[i]*wiDetId[i];
1349  }
1350  }
1351 
1352  mid_eta/=denominator;
1353  mid_phi/=denominator;
1354  mid_x/=denominator;
1355  mid_y/=denominator;
1356 
1357 
1358  // See = sigma eta eta
1359  // Spp = (B field corrected) sigma phi phi
1360  // See = (B field corrected) sigma eta phi
1361  double See=0.;
1362  double Spp=0.;
1363  double Sep=0.;
1364  double deta(0),dphi(0);
1365  // compute (phi-corrected) covariance matrix:
1366  for(unsigned int i=0; i<nCry; i++) {
1367  if(isBarrel) {
1368  deta = etaDetId[i]-mid_eta;
1369  dphi = phiDetId[i]-mid_phi;
1370  } else {
1371  deta = etaDetId[i]-mid_eta;
1372  float hitLocalR2 = (xDetId[i]-mid_x)*(xDetId[i]-mid_x)+(yDetId[i]-mid_y)*(yDetId[i]-mid_y);
1373  float hitR2 = xDetId[i]*xDetId[i]+yDetId[i]*yDetId[i];
1374  float meanR2 = mid_x*mid_x+mid_y*mid_y;
1375  float hitR = sqrt(hitR2);
1376  float meanR = sqrt(meanR2);
1377  float phi = acos((hitR2+meanR2-hitLocalR2)/(2*hitR*meanR));
1378  dphi = hitR*phi;
1379 
1380  }
1381  See += (wiDetId[i]* deta * deta) / denominator;
1382  Spp += phiCorrectionFactor*(wiDetId[i]* dphi * dphi) / denominator;
1383  Sep += sqrt(phiCorrectionFactor)*(wiDetId[i]*deta*dphi) / denominator;
1384  }
1385 
1386  Cluster2ndMoments returnMoments;
1387 
1388  // compute matrix eigenvalues:
1389  returnMoments.sMaj = ((See + Spp) + sqrt((See - Spp)*(See - Spp) + 4.*Sep*Sep)) / 2.;
1390  returnMoments.sMin = ((See + Spp) - sqrt((See - Spp)*(See - Spp) + 4.*Sep*Sep)) / 2.;
1391 
1392  returnMoments.alpha = atan( (See - Spp + sqrt( (Spp - See)*(Spp - See) + 4.*Sep*Sep )) / (2.*Sep));
1393 
1394  return returnMoments;
1395 
1396 }
int i
Definition: DBlmapReader.cc:9
const DetId & detid() const
Definition: EcalRecHit.h:71
list denominator
Definition: cuy.py:484
const T & max(const T &a, const T &b)
static float getIEta(const DetId &id)
static float getIPhi(const DetId &id)
T sqrt(T t)
Definition: SSEVec.h:48
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
float energy() const
Definition: EcalRecHit.h:68
static float getNormedIY(const DetId &id)
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
static float getSumEnergy(const std::vector< std::pair< const EcalRecHit *, float > > &RH_ptrs_fracs)
static float getNormedIX(const DetId &id)
Definition: DDAxes.h:10
template<bool noZS>
static float EcalClusterToolsT< noZS >::computeWeight ( float  eRH,
float  energyTotal,
int  weightedPositionMethod 
)
staticprivate
template<bool noZS>
std::vector< float > EcalClusterToolsT< noZS >::covariances ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology,
const CaloGeometry geometry,
float  w0 = 4.7 
)
static

Definition at line 785 of file EcalClusterTools.h.

References cuy::denominator, dPhi(), relval_parameters_module::energy, PV3DBase< T, PVType, FrameType >::eta(), f, cropTnPTrees::frac, CaloSubdetectorGeometry::getGeometry(), CaloCellGeometry::getPosition(), CaloGeometry::getSubdetectorGeometry(), CaloTopology::getSubdetectorTopology(), i, j, fff_deleter::log, max(), Geom::Phi< T >::phi(), PV3DBase< T, PVType, FrameType >::phi(), Geom::pi(), position, Geom::twoPi(), findQualityFiles::v, and w().

Referenced by GsfElectronFull5x5Filler::calculateShowerShape_full5x5(), GsfElectronAlgo::calculateShowerShape_full5x5(), and GEDPhotonProducer::fillPhotonCollection().

786 {
787  float e_5x5 = e5x5( cluster, recHits, topology );
788  float covEtaEta, covEtaPhi, covPhiPhi;
789  if (e_5x5 >= 0.) {
790  //double w0_ = parameterMap_.find("W0")->second;
791  const std::vector< std::pair<DetId, float>>& v_id =cluster.hitsAndFractions();
792  math::XYZVector meanPosition = meanClusterPosition( cluster, recHits, topology, geometry );
793 
794  // now we can calculate the covariances
795  double numeratorEtaEta = 0;
796  double numeratorEtaPhi = 0;
797  double numeratorPhiPhi = 0;
798  double denominator = 0;
799 
800  DetId id = getMaximum( v_id, recHits ).first;
801  CaloNavigator<DetId> cursor = CaloNavigator<DetId>( id, topology->getSubdetectorTopology( id ) );
802  for ( int i = -2; i <= 2; ++i ) {
803  for ( int j = -2; j <= 2; ++j ) {
804  cursor.home();
805  cursor.offsetBy( i, j );
806  float frac=getFraction(v_id,*cursor);
807  float energy = recHitEnergy( *cursor, recHits )*frac;
808 
809  if ( energy <= 0 ) continue;
810 
811  GlobalPoint position = geometry->getSubdetectorGeometry(*cursor)->getGeometry(*cursor)->getPosition();
812 
813  double dPhi = position.phi() - meanPosition.phi();
814  if (dPhi > + Geom::pi()) { dPhi = Geom::twoPi() - dPhi; }
815  if (dPhi < - Geom::pi()) { dPhi = Geom::twoPi() + dPhi; }
816 
817  double dEta = position.eta() - meanPosition.eta();
818  double w = 0.;
819  w = std::max(0.0f, w0 + std::log( energy / e_5x5 ));
820 
821  denominator += w;
822  numeratorEtaEta += w * dEta * dEta;
823  numeratorEtaPhi += w * dEta * dPhi;
824  numeratorPhiPhi += w * dPhi * dPhi;
825  }
826  }
827 
828  if (denominator != 0.0) {
829  covEtaEta = numeratorEtaEta / denominator;
830  covEtaPhi = numeratorEtaPhi / denominator;
831  covPhiPhi = numeratorPhiPhi / denominator;
832  } else {
833  covEtaEta = 999.9;
834  covEtaPhi = 999.9;
835  covPhiPhi = 999.9;
836  }
837 
838  } else {
839  // Warn the user if there was no energy in the cells and return zeroes.
840  // std::cout << "\ClusterShapeAlgo::Calculate_Covariances: no energy in supplied cells.\n";
841  covEtaEta = 0;
842  covEtaPhi = 0;
843  covPhiPhi = 0;
844  }
845  std::vector<float> v;
846  v.push_back( covEtaEta );
847  v.push_back( covEtaPhi );
848  v.push_back( covPhiPhi );
849  return v;
850 }
const CaloSubdetectorGeometry * getSubdetectorGeometry(const DetId &id) const
access the subdetector geometry for the given subdetector directly
Definition: CaloGeometry.cc:43
int i
Definition: DBlmapReader.cc:9
static math::XYZVector meanClusterPosition(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, const CaloGeometry *geometry)
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
static int position[TOTALCHAMBERS][3]
Definition: ReadPGInfo.cc:509
virtual const CaloCellGeometry * getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
list denominator
Definition: cuy.py:484
double dPhi(double phi1, double phi2)
Definition: JetUtil.h:30
const T & max(const T &a, const T &b)
T phi() const
Definition: Phi.h:41
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
int j
Definition: DBlmapReader.cc:9
double f[11][100]
static float getFraction(const std::vector< std::pair< DetId, float > > &v_id, DetId id)
Definition: DetId.h:18
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
const CaloSubdetectorTopology * getSubdetectorTopology(const DetId &id) const
access the subdetector Topology for the given subdetector directly
Definition: CaloTopology.cc:26
T eta() const
Definition: PV3DBase.h:76
double pi()
Definition: Pi.h:31
double twoPi()
Definition: Pi.h:32
T w() const
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
static float e5x5(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
template<bool noZS>
int EcalClusterToolsT< noZS >::deltaIEta ( int  seed_ieta,
int  rh_ieta 
)
staticprivate

Definition at line 1618 of file EcalClusterTools.h.

1618  {
1619  // get rid of the fact that there is no ieta=0
1620  if(seed_ieta < 0) seed_ieta++;
1621  if(rh_ieta < 0) rh_ieta++;
1622  int rel_ieta = rh_ieta - seed_ieta;
1623  return rel_ieta;
1624 }
template<bool noZS>
int EcalClusterToolsT< noZS >::deltaIPhi ( int  seed_iphi,
int  rh_iphi 
)
staticprivate

Definition at line 1606 of file EcalClusterTools.h.

1606  {
1607  int rel_iphi = rh_iphi - seed_iphi;
1608  // take care of cyclic variable iphi [1,360]
1609  if(rel_iphi > 180) rel_iphi = rel_iphi - 360;
1610  if(rel_iphi < -180) rel_iphi = rel_iphi + 360;
1611  return rel_iphi;
1612 }
template<bool noZS>
float EcalClusterToolsT< noZS >::e1x3 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 501 of file EcalClusterTools.h.

References first.

502 {
503  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
504  return matrixEnergy( cluster, recHits, topology, id, 0, 0, -1, 1 );
505 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e1x5 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 487 of file EcalClusterTools.h.

References first.

Referenced by GsfElectronFull5x5Filler::calculateShowerShape_full5x5(), GsfElectronAlgo::calculateShowerShape_full5x5(), and GEDPhotonProducer::fillPhotonCollection().

488 {
489  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
490  return matrixEnergy( cluster, recHits, topology, id, 0, 0, -2, 2 );
491 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e2nd ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits 
)
static

Definition at line 425 of file EcalClusterTools.h.

References first, i, and edm::second().

426 {
427  std::vector<float> energies;
428  const std::vector< std::pair<DetId, float> >& v_id = cluster.hitsAndFractions();
429  energies.reserve( v_id.size() );
430  if ( v_id.size() < 2 ) return 0;
431  for ( size_t i = 0; i < v_id.size(); ++i ) {
432  energies.push_back( recHitEnergy( v_id[i].first, recHits ) * (noZS ? 1.0 : v_id[i].second) );
433  }
434  std::partial_sort( energies.begin(), energies.begin()+2, energies.end(), std::greater<float>() );
435  return energies[1];
436 
437 
438 }
int i
Definition: DBlmapReader.cc:9
U second(std::pair< T, U > const &p)
bool first
Definition: L1TdeRCT.cc:75
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
template<bool noZS>
float EcalClusterToolsT< noZS >::e2x2 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 371 of file EcalClusterTools.h.

References first, and max().

372 {
373  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
374  std::list<float> energies;
375  float max_E = matrixEnergy( cluster, recHits, topology, id, -1, 0, -1, 0 );
376  max_E = std::max( max_E, matrixEnergy( cluster, recHits, topology, id, -1, 0, 0, 1 ) );
377  max_E = std::max( max_E, matrixEnergy( cluster, recHits, topology, id, 0, 1, 0, 1 ) );
378  max_E = std::max( max_E, matrixEnergy( cluster, recHits, topology, id, 0, 1, -1, 0 ) );
379  return max_E;
380 }
const T & max(const T &a, const T &b)
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e2x5Bottom ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 462 of file EcalClusterTools.h.

References first.

463 {
464  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
465  return matrixEnergy( cluster, recHits, topology, id, -2, 2, -2, -1 );
466 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e2x5Left ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 448 of file EcalClusterTools.h.

References first.

449 {
450  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
451  return matrixEnergy( cluster, recHits, topology, id, -2, -1, -2, 2 );
452 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e2x5Max ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 471 of file EcalClusterTools.h.

References first.

Referenced by GsfElectronFull5x5Filler::calculateShowerShape_full5x5(), GsfElectronAlgo::calculateShowerShape_full5x5(), and GEDPhotonProducer::fillPhotonCollection().

472 {
473  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
474 
475  // 1x5 strip left of seed
476  float left = matrixEnergy( cluster, recHits, topology, id, -1, -1, -2, 2 );
477  // 1x5 strip right of seed
478  float right = matrixEnergy( cluster, recHits, topology, id, 1, 1, -2, 2 );
479  // 1x5 strip containing seed
480  float centre = matrixEnergy( cluster, recHits, topology, id, 0, 0, -2, 2 );
481 
482  // Return the maximum of (left+center) or (right+center) strip
483  return left > right ? left+centre : right+centre;
484 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e2x5Right ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 441 of file EcalClusterTools.h.

References first.

442 {
443  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
444  return matrixEnergy( cluster, recHits, topology, id, 1, 2, -2, 2 );
445 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e2x5Top ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 455 of file EcalClusterTools.h.

References first.

456 {
457  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
458  return matrixEnergy( cluster, recHits, topology, id, -2, 2, 1, 2 );
459 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e3x1 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 508 of file EcalClusterTools.h.

References first.

509 {
510  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
511  return matrixEnergy( cluster, recHits, topology, id, -1, 1, 0, 0 );
512 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e3x2 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 383 of file EcalClusterTools.h.

References first, and max().

384 {
385  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
386  float max_E = matrixEnergy( cluster, recHits, topology, id, -1, 1, -1, 0 );
387  max_E = std::max( max_E, matrixEnergy( cluster, recHits, topology, id, 0, 1, -1, 1 ) );
388  max_E = std::max( max_E, matrixEnergy( cluster, recHits, topology, id, -1, 1, 0, 1 ) );
389  max_E = std::max( max_E, matrixEnergy( cluster, recHits, topology, id, -1, 0, -1, 1 ) );
390  return max_E;
391 }
const T & max(const T &a, const T &b)
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e3x3 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 394 of file EcalClusterTools.h.

References first.

Referenced by GsfElectronFull5x5Filler::calculateShowerShape_full5x5(), GsfElectronAlgo::calculateShowerShape_full5x5(), and GEDPhotonProducer::fillPhotonCollection().

395 {
396  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
397  return matrixEnergy( cluster, recHits, topology, id, -1, 1, -1, 1 );
398 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e4x4 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 401 of file EcalClusterTools.h.

References first, and max().

402 {
403  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
404  float max_E = matrixEnergy( cluster, recHits, topology, id, -1, 2, -2, 1 );
405  max_E = std::max( max_E, matrixEnergy( cluster, recHits, topology, id, -2, 1, -2, 1 ) );
406  max_E = std::max( max_E, matrixEnergy( cluster, recHits, topology, id, -2, 1, -1, 2 ) );
407  max_E = std::max( max_E, matrixEnergy( cluster, recHits, topology, id, -1, 2, -1, 2 ) );
408  return max_E;
409 }
const T & max(const T &a, const T &b)
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e5x1 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 494 of file EcalClusterTools.h.

References first.

495 {
496  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
497  return matrixEnergy( cluster, recHits, topology, id, -2, 2, 0, 0 );
498 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::e5x5 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 412 of file EcalClusterTools.h.

References first.

Referenced by GsfElectronFull5x5Filler::calculateShowerShape_full5x5(), GsfElectronAlgo::calculateShowerShape_full5x5(), and GEDPhotonProducer::fillPhotonCollection().

413 {
414  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
415  return matrixEnergy( cluster, recHits, topology, id, -2, 2, -2, 2 );
416 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::eBottom ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 536 of file EcalClusterTools.h.

References first.

537 {
538  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
539  return matrixEnergy( cluster, recHits, topology, id, 0, 0, -1, -1 );
540 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::eLeft ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 515 of file EcalClusterTools.h.

References first.

516 {
517  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
518  return matrixEnergy( cluster, recHits, topology, id, -1, -1, 0, 0 );
519 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::eMax ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits 
)
static

Definition at line 419 of file EcalClusterTools.h.

References edm::second().

Referenced by GEDPhotonProducer::fillPhotonCollection().

420 {
421  return getMaximum( cluster.hitsAndFractions(), recHits ).second;
422 }
U second(std::pair< T, U > const &p)
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
template<bool noZS>
std::vector< float > EcalClusterToolsT< noZS >::energyBasketFractionEta ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits 
)
static

Definition at line 543 of file EcalClusterTools.h.

References EcalBarrel, first, i, EBDetId::im(), EBDetId::kModulesPerSM, EBDetId::positiveZ(), and python.multivaluedict::sort().

544 {
545  std::vector<float> basketFraction( 2 * EBDetId::kModulesPerSM );
546  float clusterEnergy = cluster.energy();
547  const std::vector< std::pair<DetId, float> >& v_id = cluster.hitsAndFractions();
548  if ( v_id[0].first.subdetId() != EcalBarrel ) {
549  edm::LogWarning("EcalClusterToolsT<noZS>::energyBasketFractionEta") << "Trying to get basket fraction for endcap basic-clusters. Basket fractions can be obtained ONLY for barrel basic-clusters. Returning empty vector.";
550  return basketFraction;
551  }
552  for ( size_t i = 0; i < v_id.size(); ++i ) {
553  basketFraction[ EBDetId(v_id[i].first).im()-1 + EBDetId(v_id[i].first).positiveZ()*EBDetId::kModulesPerSM ] += recHitEnergy( v_id[i].first, recHits ) * v_id[i].second / clusterEnergy;
554  }
555  std::sort( basketFraction.rbegin(), basketFraction.rend() );
556  return basketFraction;
557 }
int i
Definition: DBlmapReader.cc:9
int im() const
get the number of module inside the SM (1-4)
Definition: EBDetId.h:66
static const int kModulesPerSM
Definition: EBDetId.h:147
bool first
Definition: L1TdeRCT.cc:75
bool positiveZ() const
Definition: EBDetId.h:78
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
template<bool noZS>
std::vector< float > EcalClusterToolsT< noZS >::energyBasketFractionPhi ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits 
)
static

Definition at line 560 of file EcalClusterTools.h.

References EcalBarrel, first, i, EBDetId::iphi(), EBDetId::kCrystalsInPhi, EBDetId::kTowersInPhi, EBDetId::MAX_IPHI, EBDetId::positiveZ(), and python.multivaluedict::sort().

561 {
562  std::vector<float> basketFraction( 2 * (EBDetId::MAX_IPHI / EBDetId::kCrystalsInPhi) );
563  float clusterEnergy = cluster.energy();
564  const std::vector< std::pair<DetId, float> >& v_id = cluster.hitsAndFractions();
565  if ( v_id[0].first.subdetId() != EcalBarrel ) {
566  edm::LogWarning("EcalClusterToolsT<noZS>::energyBasketFractionPhi") << "Trying to get basket fraction for endcap basic-clusters. Basket fractions can be obtained ONLY for barrel basic-clusters. Returning empty vector.";
567  return basketFraction;
568  }
569  for ( size_t i = 0; i < v_id.size(); ++i ) {
570  basketFraction[ (EBDetId(v_id[i].first).iphi()-1)/EBDetId::kCrystalsInPhi + EBDetId(v_id[i].first).positiveZ()*EBDetId::kTowersInPhi] += recHitEnergy( v_id[i].first, recHits ) * (noZS ? 1.0 : v_id[i].second) / clusterEnergy;
571  }
572  std::sort( basketFraction.rbegin(), basketFraction.rend() );
573  return basketFraction;
574 }
int i
Definition: DBlmapReader.cc:9
static const int kTowersInPhi
Definition: EBDetId.h:146
int iphi() const
get the crystal iphi
Definition: EBDetId.h:53
static const int kCrystalsInPhi
Definition: EBDetId.h:149
bool first
Definition: L1TdeRCT.cc:75
bool positiveZ() const
Definition: EBDetId.h:78
static const int MAX_IPHI
Definition: EBDetId.h:144
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
template<bool noZS>
float EcalClusterToolsT< noZS >::eRight ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 522 of file EcalClusterTools.h.

References first.

523 {
524  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
525  return matrixEnergy( cluster, recHits, topology, id, 1, 1, 0, 0 );
526 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
float EcalClusterToolsT< noZS >::eTop ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
static

Definition at line 529 of file EcalClusterTools.h.

References first.

530 {
531  DetId id = getMaximum( cluster.hitsAndFractions(), recHits ).first;
532  return matrixEnergy( cluster, recHits, topology, id, 0, 0, 1, 1 );
533 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
static float matrixEnergy(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
Definition: DetId.h:18
template<bool noZS>
static double EcalClusterToolsT< noZS >::f00 ( double  r)
inlinestaticprivate

Definition at line 210 of file EcalClusterTools.h.

210 { return 1; }
template<bool noZS>
static double EcalClusterToolsT< noZS >::f11 ( double  r)
inlinestaticprivate

Definition at line 211 of file EcalClusterTools.h.

References alignCSCRings::r.

211 { return r; }
template<bool noZS>
static double EcalClusterToolsT< noZS >::f20 ( double  r)
inlinestaticprivate

Definition at line 212 of file EcalClusterTools.h.

212 { return 2.0*r*r-1.0; }
template<bool noZS>
static double EcalClusterToolsT< noZS >::f22 ( double  r)
inlinestaticprivate

Definition at line 213 of file EcalClusterTools.h.

References alignCSCRings::r.

213 { return r*r; }
template<bool noZS>
static double EcalClusterToolsT< noZS >::f31 ( double  r)
inlinestaticprivate

Definition at line 214 of file EcalClusterTools.h.

References alignCSCRings::r.

214 { return 3.0*r*r*r - 2.0*r; }
template<bool noZS>
static double EcalClusterToolsT< noZS >::f33 ( double  r)
inlinestaticprivate

Definition at line 215 of file EcalClusterTools.h.

References alignCSCRings::r.

215 { return r*r*r; }
template<bool noZS>
static double EcalClusterToolsT< noZS >::f40 ( double  r)
inlinestaticprivate

Definition at line 216 of file EcalClusterTools.h.

216 { return 6.0*r*r*r*r-6.0*r*r+1.0; }
template<bool noZS>
static double EcalClusterToolsT< noZS >::f42 ( double  r)
inlinestaticprivate

Definition at line 217 of file EcalClusterTools.h.

References alignCSCRings::r.

217 { return 4.0*r*r*r*r-3.0*r*r; }
template<bool noZS>
static double EcalClusterToolsT< noZS >::f44 ( double  r)
inlinestaticprivate

Definition at line 218 of file EcalClusterTools.h.

References alignCSCRings::r.

218 { return r*r*r*r; }
template<bool noZS>
static double EcalClusterToolsT< noZS >::f51 ( double  r)
inlinestaticprivate

Definition at line 219 of file EcalClusterTools.h.

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

219 { return 10.0*pow(r,5)-12.0*pow(r,3)+3.0*r; }
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
template<bool noZS>
static double EcalClusterToolsT< noZS >::f53 ( double  r)
inlinestaticprivate

Definition at line 220 of file EcalClusterTools.h.

References funct::pow().

220 { return 5.0*pow(r,5) - 4.0*pow(r,3); }
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
template<bool noZS>
static double EcalClusterToolsT< noZS >::f55 ( double  r)
inlinestaticprivate

Definition at line 221 of file EcalClusterTools.h.

References funct::pow().

221 { return pow(r,5); }
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
template<bool noZS>
static double EcalClusterToolsT< noZS >::factorial ( int  n)
inlinestaticprivate

Definition at line 227 of file EcalClusterTools.h.

References i, and n.

227  {
228  double res = 1.;
229  for (int i = 2; i <= n; ++i) res *= i;
230  return res;
231  }
int i
Definition: DBlmapReader.cc:9
template<bool noZS>
double EcalClusterToolsT< noZS >::fast_AbsZernikeMoment ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloGeometry geometry,
int  n,
int  m,
double  R0,
bool  logW,
float  w0 
)
staticprivate

Definition at line 964 of file EcalClusterTools.h.

References funct::cos(), alignCSCRings::e, i, cmsHarvester::index, phi, alignCSCRings::r, funct::sin(), and mathSSE::sqrt().

965 {
966  double r,ph,e,Re=0,Im=0;
967  double TotalEnergy = cluster.energy();
968  int index = (n/2)*(n/2)+(n/2)+m;
969  std::vector<EcalClusterEnergyDeposition> energyDistribution = getEnergyDepTopology( cluster, recHits, geometry, logW, w0 );
970  int clusterSize = energyDistribution.size();
971  if(clusterSize < 3) return 0.0;
972 
973  for (int i=0; i<clusterSize; i++)
974  {
975  r = energyDistribution[i].r / R0;
976  if (r<1) {
977  std::vector<double> pol;
978  pol.push_back( f00(r) );
979  pol.push_back( f11(r) );
980  pol.push_back( f20(r) );
981  pol.push_back( f22(r) );
982  pol.push_back( f31(r) );
983  pol.push_back( f33(r) );
984  pol.push_back( f40(r) );
985  pol.push_back( f42(r) );
986  pol.push_back( f44(r) );
987  pol.push_back( f51(r) );
988  pol.push_back( f53(r) );
989  pol.push_back( f55(r) );
990  ph = (energyDistribution[i]).phi;
991  e = energyDistribution[i].deposited_energy;
992  Re = Re + e/TotalEnergy * pol[index] * cos( (double) m * ph);
993  Im = Im - e/TotalEnergy * pol[index] * sin( (double) m * ph);
994  }
995  }
996  return sqrt(Re*Re+Im*Im);
997 }
int i
Definition: DBlmapReader.cc:9
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
static double f33(double r)
static double f11(double r)
static std::vector< EcalClusterEnergyDeposition > getEnergyDepTopology(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, bool logW, float w0)
static double f55(double r)
T sqrt(T t)
Definition: SSEVec.h:48
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
static double f40(double r)
static double f31(double r)
static double f44(double r)
static double f00(double r)
static double f53(double r)
static double f42(double r)
static double f51(double r)
static double f20(double r)
static double f22(double r)
Definition: DDAxes.h:10
template<bool noZS>
float EcalClusterToolsT< noZS >::getDPhiEndcap ( const DetId crysId,
float  meanX,
float  meanY 
)
staticprivate

Definition at line 1133 of file EcalClusterTools.h.

References dPhi(), phi, mathSSE::sqrt(), and tmp.

1134 {
1135  float iXNorm = getNormedIX(crysId);
1136  float iYNorm = getNormedIY(crysId);
1137 
1138  float hitLocalR2 = (iXNorm-meanX)*(iXNorm-meanX)+(iYNorm-meanY)*(iYNorm-meanY);
1139  float hitR2 = iXNorm*iXNorm+iYNorm*iYNorm;
1140  float meanR2 = meanX*meanX+meanY*meanY;
1141  float hitR = sqrt(hitR2);
1142  float meanR = sqrt(meanR2);
1143 
1144  float tmp = (hitR2+meanR2-hitLocalR2)/(2*hitR*meanR);
1145  if (tmp<-1) tmp =-1;
1146  if (tmp>1) tmp=1;
1147  float phi = acos(tmp);
1148  float dPhi = hitR*phi;
1149 
1150  return dPhi;
1151 }
double dPhi(double phi1, double phi2)
Definition: JetUtil.h:30
T sqrt(T t)
Definition: SSEVec.h:48
static float getNormedIY(const DetId &id)
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
static float getNormedIX(const DetId &id)
Definition: DDAxes.h:10
template<bool noZS>
std::vector< typename EcalClusterToolsT< noZS >::EcalClusterEnergyDeposition > EcalClusterToolsT< noZS >::getEnergyDepTopology ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloGeometry geometry,
bool  logW,
float  w0 
)
staticprivate

Definition at line 577 of file EcalClusterTools.h.

References funct::abs(), EcalClusterToolsT< noZS >::EcalClusterEnergyDeposition::deposited_energy, diffTreeTool::diff, edm::SortedCollection< T, SORT >::end(), EcalRecHit::energy(), edm::SortedCollection< T, SORT >::find(), CaloSubdetectorGeometry::getGeometry(), CaloCellGeometry::getPosition(), CaloGeometry::getSubdetectorGeometry(), fff_deleter::log, LogDebug, M_PI, max(), EcalClusterToolsT< noZS >::EcalClusterEnergyDeposition::phi, EcalClusterToolsT< noZS >::EcalClusterEnergyDeposition::r, histoStyle::weight, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

578 {
579  std::vector<typename EcalClusterToolsT<noZS>::EcalClusterEnergyDeposition> energyDistribution;
580  // init a map of the energy deposition centered on the
581  // cluster centroid. This is for momenta calculation only.
582  CLHEP::Hep3Vector clVect(cluster.position().x(), cluster.position().y(), cluster.position().z());
583  CLHEP::Hep3Vector clDir(clVect);
584  clDir*=1.0/clDir.mag();
585  // in the transverse plane, axis perpendicular to clusterDir
586  CLHEP::Hep3Vector theta_axis(clDir.y(),-clDir.x(),0.0);
587  theta_axis *= 1.0/theta_axis.mag();
588  CLHEP::Hep3Vector phi_axis = theta_axis.cross(clDir);
589 
590  const std::vector< std::pair<DetId, float> >& clusterDetIds = cluster.hitsAndFractions();
591 
593  EcalRecHit testEcalRecHit;
594  std::vector< std::pair<DetId, float> >::const_iterator posCurrent;
595  // loop over crystals
596  for(posCurrent=clusterDetIds.begin(); posCurrent!=clusterDetIds.end(); ++posCurrent) {
597  EcalRecHitCollection::const_iterator itt = recHits->find( (*posCurrent).first );
598  testEcalRecHit=*itt;
599 
600  if(( (*posCurrent).first != DetId(0)) && (recHits->find( (*posCurrent).first ) != recHits->end())) {
601  clEdep.deposited_energy = testEcalRecHit.energy() * (noZS ? 1.0 : (*posCurrent).second);
602  // if logarithmic weight is requested, apply cut on minimum energy of the recHit
603  if(logW) {
604  //double w0 = parameterMap_.find("W0")->second;
605 
606  double weight = std::max(0.0, w0 + log(std::abs(clEdep.deposited_energy)/cluster.energy()) );
607  if(weight==0) {
608  LogDebug("ClusterShapeAlgo") << "Crystal has insufficient energy: E = "
609  << clEdep.deposited_energy << " GeV; skipping... ";
610  continue;
611  }
612  else LogDebug("ClusterShapeAlgo") << "===> got crystal. Energy = " << clEdep.deposited_energy << " GeV. ";
613  }
614  DetId id_ = (*posCurrent).first;
615  const CaloCellGeometry *this_cell = geometry->getSubdetectorGeometry(id_)->getGeometry(id_);
616  GlobalPoint cellPos = this_cell->getPosition();
617  CLHEP::Hep3Vector gblPos (cellPos.x(),cellPos.y(),cellPos.z()); //surface position?
618  // Evaluate the distance from the cluster centroid
619  CLHEP::Hep3Vector diff = gblPos - clVect;
620  // Important: for the moment calculation, only the "lateral distance" is important
621  // "lateral distance" r_i = distance of the digi position from the axis Origin-Cluster Center
622  // ---> subtract the projection on clDir
623  CLHEP::Hep3Vector DigiVect = diff - diff.dot(clDir)*clDir;
624  clEdep.r = DigiVect.mag();
625  LogDebug("ClusterShapeAlgo") << "E = " << clEdep.deposited_energy
626  << "\tdiff = " << diff.mag()
627  << "\tr = " << clEdep.r;
628  clEdep.phi = DigiVect.angle(theta_axis);
629  if(DigiVect.dot(phi_axis)<0) clEdep.phi = 2 * M_PI - clEdep.phi;
630  energyDistribution.push_back(clEdep);
631  }
632  }
633  return energyDistribution;
634 }
#define LogDebug(id)
const CaloSubdetectorGeometry * getSubdetectorGeometry(const DetId &id) const
access the subdetector geometry for the given subdetector directly
Definition: CaloGeometry.cc:43
std::vector< EcalRecHit >::const_iterator const_iterator
T y() const
Definition: PV3DBase.h:63
virtual const CaloCellGeometry * getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
const T & max(const T &a, const T &b)
T z() const
Definition: PV3DBase.h:64
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
float energy() const
Definition: EcalRecHit.h:68
#define M_PI
const_iterator end() const
Definition: DetId.h:18
iterator find(key_type k)
int weight
Definition: histoStyle.py:50
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
T x() const
Definition: PV3DBase.h:62
template<bool noZS>
float EcalClusterToolsT< noZS >::getFraction ( const std::vector< std::pair< DetId, float > > &  v_id,
DetId  id 
)
static

Definition at line 253 of file EcalClusterTools.h.

References first, cropTnPTrees::frac, and i.

254  {
255  if(noZS) return 1.0;
256  float frac = 0.0;
257  for ( size_t i = 0; i < v_id.size(); ++i ) {
258  if(v_id[i].first.rawId()==id.rawId()){
259  frac= v_id[i].second;
260  break;
261  }
262  }
263  return frac;
264 }
int i
Definition: DBlmapReader.cc:9
bool first
Definition: L1TdeRCT.cc:75
template<bool noZS>
float EcalClusterToolsT< noZS >::getIEta ( const DetId id)
staticprivate

Definition at line 1034 of file EcalClusterTools.h.

References DetId::Ecal, EcalBarrel, EcalEndcap, EBDetId::ieta(), and mathSSE::sqrt().

1035 {
1036  if(id.det()==DetId::Ecal){
1037  if(id.subdetId()==EcalBarrel){
1038  EBDetId ebId(id);
1039  return ebId.ieta();
1040  }else if(id.subdetId()==EcalEndcap){
1041  float iXNorm = getNormedIX(id);
1042  float iYNorm = getNormedIY(id);
1043 
1044  return std::sqrt(iXNorm*iXNorm+iYNorm*iYNorm);
1045  }
1046  }
1047  return 0.;
1048 }
T sqrt(T t)
Definition: SSEVec.h:48
static float getNormedIY(const DetId &id)
static float getNormedIX(const DetId &id)
template<bool noZS>
float EcalClusterToolsT< noZS >::getIPhi ( const DetId id)
staticprivate

Definition at line 1056 of file EcalClusterTools.h.

References DetId::Ecal, EcalBarrel, and EBDetId::iphi().

1057 {
1058  if(id.det()==DetId::Ecal){
1059  if(id.subdetId()==EcalBarrel){
1060  EBDetId ebId(id);
1061  return ebId.iphi();
1062  }
1063  }
1064  return 0.;
1065 }
template<bool noZS>
std::pair< DetId, float > EcalClusterToolsT< noZS >::getMaximum ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits 
)
static

Definition at line 282 of file EcalClusterTools.h.

283 {
284  return getMaximum( cluster.hitsAndFractions(), recHits );
285 }
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
template<bool noZS>
std::pair< DetId, float > EcalClusterToolsT< noZS >::getMaximum ( const std::vector< std::pair< DetId, float > > &  v_id,
const EcalRecHitCollection recHits 
)
static

Definition at line 267 of file EcalClusterTools.h.

References relval_parameters_module::energy, first, i, and max().

268 {
269  float max = 0;
270  DetId id(0);
271  for ( size_t i = 0; i < v_id.size(); ++i ) {
272  float energy = recHitEnergy( v_id[i].first, recHits ) * (noZS ? 1.0 : v_id[i].second);
273  if ( energy > max ) {
274  max = energy;
275  id = v_id[i].first;
276  }
277  }
278  return std::pair<DetId, float>(id, max);
279 }
int i
Definition: DBlmapReader.cc:9
const T & max(const T &a, const T &b)
bool first
Definition: L1TdeRCT.cc:75
Definition: DetId.h:18
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
template<bool noZS>
float EcalClusterToolsT< noZS >::getNormedIX ( const DetId id)
staticprivate

Definition at line 1069 of file EcalClusterTools.h.

References DetId::Ecal, EcalEndcap, and EEDetId::ix().

1070 {
1071  if(id.det()==DetId::Ecal && id.subdetId()==EcalEndcap){
1072  EEDetId eeId(id);
1073  int iXNorm = eeId.ix()-50;
1074  if(iXNorm<=0) iXNorm--;
1075  return iXNorm;
1076  }
1077  return 0;
1078 }
template<bool noZS>
float EcalClusterToolsT< noZS >::getNormedIY ( const DetId id)
staticprivate

Definition at line 1082 of file EcalClusterTools.h.

References DetId::Ecal, EcalEndcap, and EEDetId::iy().

1083 {
1084  if(id.det()==DetId::Ecal && id.subdetId()==EcalEndcap){
1085  EEDetId eeId(id);
1086  int iYNorm = eeId.iy()-50;
1087  if(iYNorm<=0) iYNorm--;
1088  return iYNorm;
1089  }
1090  return 0;
1091 }
template<bool noZS>
float EcalClusterToolsT< noZS >::getNrCrysDiffInEta ( const DetId crysId,
const DetId orginId 
)
staticprivate

Definition at line 1095 of file EcalClusterTools.h.

References EcalBarrel, and DetId::subdetId().

1096 {
1097  float crysIEta = getIEta(crysId);
1098  float orginIEta = getIEta(orginId);
1099  bool isBarrel = orginId.subdetId()==EcalBarrel;
1100 
1101  float nrCrysDiff = crysIEta-orginIEta;
1102 
1103  //no iEta=0 in barrel, so if go from positive to negative
1104  //need to reduce abs(detEta) by 1
1105  if(isBarrel){
1106  if(crysIEta*orginIEta<0){ // -1 to 1 transition
1107  if(crysIEta>0) nrCrysDiff--;
1108  else nrCrysDiff++;
1109  }
1110  }
1111  return nrCrysDiff;
1112 }
static float getIEta(const DetId &id)
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
template<bool noZS>
float EcalClusterToolsT< noZS >::getNrCrysDiffInPhi ( const DetId crysId,
const DetId orginId 
)
staticprivate

Definition at line 1116 of file EcalClusterTools.h.

References EcalBarrel, and DetId::subdetId().

1117 {
1118  float crysIPhi = getIPhi(crysId);
1119  float orginIPhi = getIPhi(orginId);
1120  bool isBarrel = orginId.subdetId()==EcalBarrel;
1121 
1122  float nrCrysDiff = crysIPhi-orginIPhi;
1123 
1124  if(isBarrel){ //if barrel, need to map into 0-180
1125  if (nrCrysDiff > + 180) { nrCrysDiff = nrCrysDiff - 360; }
1126  if (nrCrysDiff < - 180) { nrCrysDiff = nrCrysDiff + 360; }
1127  }
1128  return nrCrysDiff;
1129 }
static float getIPhi(const DetId &id)
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
template<bool noZS>
std::vector< int > EcalClusterToolsT< noZS >::getSeedPosition ( const std::vector< std::pair< const EcalRecHit *, float > > &  RH_ptrs)
staticprivate

Definition at line 1628 of file EcalClusterTools.h.

References EcalRecHit::detid(), and EcalRecHit::energy().

Referenced by EcalClusterToolsT< noZS >::roundnessBarrelSuperClustersUserExtended(), and EcalClusterToolsT< noZS >::roundnessSelectedBarrelRecHits().

1628  {
1629  std::vector<int> seedPosition;
1630  float eSeedRH = 0;
1631  int iEtaSeedRH = 0;
1632  int iPhiSeedRH = 0;
1633 
1634  for(std::vector<std::pair<const EcalRecHit*,float> >::const_iterator rhf_ptr = RH_ptrs_fracs.begin(); rhf_ptr != RH_ptrs_fracs.end(); rhf_ptr++){
1635  const EcalRecHit* rh_ptr = rhf_ptr->first;
1636  //get iEta, iPhi
1637  EBDetId EBdetIdi( rh_ptr->detid() );
1638  float rh_energy = rh_ptr->energy() * (noZS ? 1.0 : rhf_ptr->second);
1639 
1640  if(eSeedRH < rh_energy){
1641  eSeedRH = rh_energy;
1642  iEtaSeedRH = EBdetIdi.ieta();
1643  iPhiSeedRH = EBdetIdi.iphi();
1644  }
1645 
1646  }// for loop
1647 
1648  seedPosition.push_back(iEtaSeedRH);
1649  seedPosition.push_back(iPhiSeedRH);
1650  return seedPosition;
1651 }
const DetId & detid() const
Definition: EcalRecHit.h:71
float energy() const
Definition: EcalRecHit.h:68
template<bool noZS>
float EcalClusterToolsT< noZS >::getSumEnergy ( const std::vector< std::pair< const EcalRecHit *, float > > &  RH_ptrs_fracs)
staticprivate

Definition at line 1655 of file EcalClusterTools.h.

Referenced by EcalClusterToolsT< noZS >::cluster2ndMoments(), and EcalClusterToolsT< noZS >::roundnessSelectedBarrelRecHits().

1655  {
1656  float sumE = 0.;
1657  for( const auto& hAndF : RH_ptrs_fracs ) {
1658  sumE += hAndF.first->energy() * (noZS ? 1.0 : hAndF.second);
1659  }
1660  return sumE;
1661 }
template<bool noZS>
std::vector< float > EcalClusterToolsT< noZS >::lat ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloGeometry geometry,
bool  logW = true,
float  w0 = 4.7 
)
static

Definition at line 637 of file EcalClusterTools.h.

References funct::cos(), reco::e1, reco::e2, i, n, phi, alignCSCRings::r, funct::sin(), and tmp.

638 {
639  std::vector<EcalClusterToolsT::EcalClusterEnergyDeposition> energyDistribution = getEnergyDepTopology( cluster, recHits, geometry, logW, w0 );
640 
641  std::vector<float> lat;
642  double r, redmoment=0;
643  double phiRedmoment = 0 ;
644  double etaRedmoment = 0 ;
645  int n,n1,n2,tmp;
646  int clusterSize=energyDistribution.size();
647  float etaLat_, phiLat_, lat_;
648  if (clusterSize<3) {
649  etaLat_ = 0.0 ;
650  lat_ = 0.0;
651  lat.push_back(0.);
652  lat.push_back(0.);
653  lat.push_back(0.);
654  return lat;
655  }
656 
657  n1=0; n2=1;
658  if (energyDistribution[1].deposited_energy >
659  energyDistribution[0].deposited_energy)
660  {
661  tmp=n2; n2=n1; n1=tmp;
662  }
663  for (int i=2; i<clusterSize; i++) {
664  n=i;
665  if (energyDistribution[i].deposited_energy >
666  energyDistribution[n1].deposited_energy)
667  {
668  tmp = n2;
669  n2 = n1; n1 = i; n=tmp;
670  } else {
671  if (energyDistribution[i].deposited_energy >
672  energyDistribution[n2].deposited_energy)
673  {
674  tmp=n2; n2=i; n=tmp;
675  }
676  }
677 
678  r = energyDistribution[n].r;
679  redmoment += r*r* energyDistribution[n].deposited_energy;
680  double rphi = r * cos (energyDistribution[n].phi) ;
681  phiRedmoment += rphi * rphi * energyDistribution[n].deposited_energy;
682  double reta = r * sin (energyDistribution[n].phi) ;
683  etaRedmoment += reta * reta * energyDistribution[n].deposited_energy;
684  }
685  double e1 = energyDistribution[n1].deposited_energy;
686  double e2 = energyDistribution[n2].deposited_energy;
687 
688  lat_ = redmoment/(redmoment+2.19*2.19*(e1+e2));
689  phiLat_ = phiRedmoment/(phiRedmoment+2.19*2.19*(e1+e2));
690  etaLat_ = etaRedmoment/(etaRedmoment+2.19*2.19*(e1+e2));
691 
692  lat.push_back(etaLat_);
693  lat.push_back(phiLat_);
694  lat.push_back(lat_);
695  return lat;
696 }
int i
Definition: DBlmapReader.cc:9
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
static std::vector< EcalClusterEnergyDeposition > getEnergyDepTopology(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, bool logW, float w0)
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
Float e1
Definition: deltaR.h:22
static std::vector< float > lat(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, bool logW=true, float w0=4.7)
Float e2
Definition: deltaR.h:23
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
Definition: DDAxes.h:10
template<bool noZS>
std::vector< float > EcalClusterToolsT< noZS >::localCovariances ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology,
float  w0 = 4.7 
)
static

Definition at line 857 of file EcalClusterTools.h.

References cuy::denominator, dPhi(), EcalBarrel, relval_parameters_module::energy, f, cropTnPTrees::frac, CaloTopology::getSubdetectorTopology(), fff_deleter::log, max(), DetId::subdetId(), findQualityFiles::v, and w().

Referenced by GsfElectronFull5x5Filler::calculateShowerShape_full5x5(), GsfElectronAlgo::calculateShowerShape_full5x5(), and GEDPhotonProducer::fillPhotonCollection().

858 {
859 
860  float e_5x5 = e5x5( cluster, recHits, topology );
861  float covEtaEta, covEtaPhi, covPhiPhi;
862 
863  if (e_5x5 >= 0.) {
864  //double w0_ = parameterMap_.find("W0")->second;
865  const std::vector< std::pair<DetId, float> >& v_id = cluster.hitsAndFractions();
866  std::pair<float,float> mean5x5PosInNrCrysFromSeed = mean5x5PositionInLocalCrysCoord( cluster, recHits, topology );
867  std::pair<float,float> mean5x5XYPos = mean5x5PositionInXY(cluster,recHits,topology);
868 
869  // now we can calculate the covariances
870  double numeratorEtaEta = 0;
871  double numeratorEtaPhi = 0;
872  double numeratorPhiPhi = 0;
873  double denominator = 0;
874 
875  //these allow us to scale the localCov by the crystal size
876  //so that the localCovs have the same average value as the normal covs
877  const double barrelCrysSize = 0.01745; //approximate size of crystal in eta,phi in barrel
878  const double endcapCrysSize = 0.0447; //the approximate crystal size sigmaEtaEta was corrected to in the endcap
879 
880  DetId seedId = getMaximum( v_id, recHits ).first;
881 
882  bool isBarrel=seedId.subdetId()==EcalBarrel;
883  const double crysSize = isBarrel ? barrelCrysSize : endcapCrysSize;
884 
885  CaloNavigator<DetId> cursor = CaloNavigator<DetId>( seedId, topology->getSubdetectorTopology( seedId ) );
886 
887  for ( int eastNr = -2; eastNr <= 2; ++eastNr ) { //east is eta in barrel
888  for ( int northNr = -2; northNr <= 2; ++northNr ) { //north is phi in barrel
889  cursor.home();
890  cursor.offsetBy( eastNr, northNr);
891  float frac = getFraction(v_id,*cursor);
892  float energy = recHitEnergy( *cursor, recHits )*frac;
893  if ( energy <= 0 ) continue;
894 
895  float dEta = getNrCrysDiffInEta(*cursor,seedId) - mean5x5PosInNrCrysFromSeed.first;
896  float dPhi = 0;
897 
898  if(isBarrel) dPhi = getNrCrysDiffInPhi(*cursor,seedId) - mean5x5PosInNrCrysFromSeed.second;
899  else dPhi = getDPhiEndcap(*cursor,mean5x5XYPos.first,mean5x5XYPos.second);
900 
901 
902  double w = std::max(0.0f,w0 + std::log( energy / e_5x5 ));
903 
904  denominator += w;
905  numeratorEtaEta += w * dEta * dEta;
906  numeratorEtaPhi += w * dEta * dPhi;
907  numeratorPhiPhi += w * dPhi * dPhi;
908  } //end east loop
909  }//end north loop
910 
911 
912  //multiplying by crysSize to make the values compariable to normal covariances
913  if (denominator != 0.0) {
914  covEtaEta = crysSize*crysSize* numeratorEtaEta / denominator;
915  covEtaPhi = crysSize*crysSize* numeratorEtaPhi / denominator;
916  covPhiPhi = crysSize*crysSize* numeratorPhiPhi / denominator;
917  } else {
918  covEtaEta = 999.9;
919  covEtaPhi = 999.9;
920  covPhiPhi = 999.9;
921  }
922 
923 
924  } else {
925  // Warn the user if there was no energy in the cells and return zeroes.
926  // std::cout << "\ClusterShapeAlgo::Calculate_Covariances: no energy in supplied cells.\n";
927  covEtaEta = 0;
928  covEtaPhi = 0;
929  covPhiPhi = 0;
930  }
931  std::vector<float> v;
932  v.push_back( covEtaEta );
933  v.push_back( covEtaPhi );
934  v.push_back( covPhiPhi );
935  return v;
936 }
list denominator
Definition: cuy.py:484
double dPhi(double phi1, double phi2)
Definition: JetUtil.h:30
const T & max(const T &a, const T &b)
static float getDPhiEndcap(const DetId &crysId, float meanX, float meanY)
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
double f[11][100]
static float getNrCrysDiffInEta(const DetId &crysId, const DetId &orginId)
static float getFraction(const std::vector< std::pair< DetId, float > > &v_id, DetId id)
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
static float getNrCrysDiffInPhi(const DetId &crysId, const DetId &orginId)
static std::pair< float, float > mean5x5PositionInXY(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
Definition: DetId.h:18
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
const CaloSubdetectorTopology * getSubdetectorTopology(const DetId &id) const
access the subdetector Topology for the given subdetector directly
Definition: CaloTopology.cc:26
static std::pair< float, float > mean5x5PositionInLocalCrysCoord(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
T w() const
static float e5x5(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
template<bool noZS>
std::vector< DetId > EcalClusterToolsT< noZS >::matrixDetId ( const CaloTopology topology,
DetId  id,
int  ixMin,
int  ixMax,
int  iyMin,
int  iyMax 
)
static

Definition at line 355 of file EcalClusterTools.h.

References CaloTopology::getSubdetectorTopology(), i, j, and findQualityFiles::v.

356 {
357  CaloNavigator<DetId> cursor = CaloNavigator<DetId>( id, topology->getSubdetectorTopology( id ) );
358  std::vector<DetId> v;
359  for ( int i = ixMin; i <= ixMax; ++i ) {
360  for ( int j = iyMin; j <= iyMax; ++j ) {
361  cursor.home();
362  cursor.offsetBy( i, j );
363  if ( *cursor != DetId(0) ) v.push_back( *cursor );
364  }
365  }
366  return v;
367 }
int i
Definition: DBlmapReader.cc:9
int j
Definition: DBlmapReader.cc:9
Definition: DetId.h:18
const CaloSubdetectorTopology * getSubdetectorTopology(const DetId &id) const
access the subdetector Topology for the given subdetector directly
Definition: CaloTopology.cc:26
template<bool noZS>
float EcalClusterToolsT< noZS >::matrixEnergy ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology,
DetId  id,
int  ixMin,
int  ixMax,
int  iyMin,
int  iyMax 
)
static

Definition at line 329 of file EcalClusterTools.h.

References relval_parameters_module::energy, cropTnPTrees::frac, CaloTopology::getSubdetectorTopology(), i, and j.

330 {
331  //take into account fractions
332  // fast version
333  CaloNavigator<DetId> cursor = CaloNavigator<DetId>( id, topology->getSubdetectorTopology( id ) );
334  float energy = 0;
335  const std::vector< std::pair<DetId, float> >& v_id = cluster.hitsAndFractions();
336  for ( int i = ixMin; i <= ixMax; ++i ) {
337  for ( int j = iyMin; j <= iyMax; ++j ) {
338  cursor.home();
339  cursor.offsetBy( i, j );
340  float frac=getFraction(v_id,*cursor);
341  energy += recHitEnergy( *cursor, recHits )*frac;
342  }
343  }
344  // slow elegant version
345  //float energy = 0;
346  //std::vector<DetId> v_id = matrixDetId( topology, id, ixMin, ixMax, iyMin, iyMax );
347  //for ( std::vector<DetId>::const_iterator it = v_id.begin(); it != v_id.end(); ++it ) {
348  // energy += recHitEnergy( *it, recHits );
349  //}
350  return energy;
351 }
int i
Definition: DBlmapReader.cc:9
int j
Definition: DBlmapReader.cc:9
static float getFraction(const std::vector< std::pair< DetId, float > > &v_id, DetId id)
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
const CaloSubdetectorTopology * getSubdetectorTopology(const DetId &id) const
access the subdetector Topology for the given subdetector directly
Definition: CaloTopology.cc:26
template<bool noZS>
std::pair< float, float > EcalClusterToolsT< noZS >::mean5x5PositionInLocalCrysCoord ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
staticprivate

Definition at line 725 of file EcalClusterTools.h.

References relval_parameters_module::energy, and CastorDataFrameFilter_impl::energySum().

726 {
727  DetId seedId = getMaximum( cluster, recHits ).first;
728  float meanDEta=0.;
729  float meanDPhi=0.;
730  float energySum=0.;
731 
732  const std::vector<std::pair<DetId,float> >& hsAndFs = cluster.hitsAndFractions();
733  std::vector<DetId> v_id = matrixDetId( topology,seedId, -2, 2, -2, 2 );
734  for( const std::pair<DetId,float>& hAndF : hsAndFs ) {
735  for ( std::vector<DetId>::const_iterator it = v_id.begin(); it != v_id.end(); ++it ) {
736  if( hAndF.first != *it && !noZS ) continue;
737  float energy = recHitEnergy(*it,recHits) * hAndF.second;
738  if(energy<0.) continue;//skipping negative energy crystals
739  meanDEta += energy * getNrCrysDiffInEta(*it,seedId);
740  meanDPhi += energy * getNrCrysDiffInPhi(*it,seedId);
741  energySum +=energy;
742  }
743  if(noZS) break;
744  }
745  meanDEta /=energySum;
746  meanDPhi /=energySum;
747  return std::pair<float,float>(meanDEta,meanDPhi);
748 }
static std::vector< DetId > matrixDetId(const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
static float getNrCrysDiffInEta(const DetId &crysId, const DetId &orginId)
static float getNrCrysDiffInPhi(const DetId &crysId, const DetId &orginId)
Definition: DetId.h:18
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
double energySum(const DataFrame &df, int fs, int ls)
template<bool noZS>
std::pair< float, float > EcalClusterToolsT< noZS >::mean5x5PositionInXY ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology 
)
staticprivate

Definition at line 757 of file EcalClusterTools.h.

References EcalBarrel, relval_parameters_module::energy, CastorDataFrameFilter_impl::energySum(), and DetId::subdetId().

758 {
759  DetId seedId = getMaximum( cluster, recHits ).first;
760 
761  std::pair<float,float> meanXY(0.,0.);
762  if(seedId.subdetId()==EcalBarrel) return meanXY;
763 
764  float energySum=0.;
765 
766  const std::vector<std::pair<DetId,float> >& hsAndFs = cluster.hitsAndFractions();
767  std::vector<DetId> v_id = matrixDetId( topology,seedId, -2, 2, -2, 2 );
768  for( const std::pair<DetId,float>& hAndF : hsAndFs ) {
769  for ( std::vector<DetId>::const_iterator it = v_id.begin(); it != v_id.end(); ++it ) {
770  if( hAndF.first != *it && !noZS) continue;
771  float energy = recHitEnergy(*it,recHits) * hAndF.second;
772  if(energy<0.) continue;//skipping negative energy crystals
773  meanXY.first += energy * getNormedIX(*it);
774  meanXY.second += energy * getNormedIY(*it);
775  energySum +=energy;
776  }
777  if(noZS) break;
778  }
779  meanXY.first/=energySum;
780  meanXY.second/=energySum;
781  return meanXY;
782 }
static std::vector< DetId > matrixDetId(const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
static float getNormedIY(const DetId &id)
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
Definition: DetId.h:18
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
double energySum(const DataFrame &df, int fs, int ls)
static float getNormedIX(const DetId &id)
template<bool noZS>
math::XYZVector EcalClusterToolsT< noZS >::meanClusterPosition ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology,
const CaloGeometry geometry 
)
staticprivate

Definition at line 699 of file EcalClusterTools.h.

References first, CaloSubdetectorGeometry::getGeometry(), CaloCellGeometry::getPosition(), CaloGeometry::getSubdetectorGeometry(), position, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

700 {
701  // find mean energy position of a 5x5 cluster around the maximum
702  math::XYZVector meanPosition(0.0, 0.0, 0.0);
703  const std::vector<std::pair<DetId,float> >& hsAndFs = cluster.hitsAndFractions();
704  std::vector<DetId> v_id = matrixDetId( topology, getMaximum( cluster, recHits ).first, -2, 2, -2, 2 );
705  for( const std::pair<DetId,float>& hitAndFrac : hsAndFs ) {
706  for( std::vector<DetId>::const_iterator it = v_id.begin(); it != v_id.end(); ++it ) {
707  if( hitAndFrac.first != *it && !noZS) continue;
708  GlobalPoint positionGP = geometry->getSubdetectorGeometry( *it )->getGeometry( *it )->getPosition();
709  math::XYZVector position(positionGP.x(),positionGP.y(),positionGP.z());
710  meanPosition = meanPosition + recHitEnergy( *it, recHits ) * position * hitAndFrac.second;
711  }
712  if(noZS) break;
713  }
714  return meanPosition / e5x5( cluster, recHits, topology );
715 }
const CaloSubdetectorGeometry * getSubdetectorGeometry(const DetId &id) const
access the subdetector geometry for the given subdetector directly
Definition: CaloGeometry.cc:43
static std::vector< DetId > matrixDetId(const CaloTopology *topology, DetId id, int ixMin, int ixMax, int iyMin, int iyMax)
T y() const
Definition: PV3DBase.h:63
static int position[TOTALCHAMBERS][3]
Definition: ReadPGInfo.cc:509
virtual const CaloCellGeometry * getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
T z() const
Definition: PV3DBase.h:64
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
bool first
Definition: L1TdeRCT.cc:75
XYZVectorD XYZVector
spatial vector with cartesian internal representation
Definition: Vector3D.h:30
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
T x() const
Definition: PV3DBase.h:62
static float e5x5(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
template<bool noZS>
float EcalClusterToolsT< noZS >::recHitEnergy ( DetId  id,
const EcalRecHitCollection recHits 
)
static

Definition at line 289 of file EcalClusterTools.h.

References EcalBarrel, edm::SortedCollection< T, SORT >::end(), edm::SortedCollection< T, SORT >::find(), EcalRecHit::kDiWeird, EcalRecHit::kTowerRecovered, and EcalRecHit::kWeird.

290 {
291  if ( id == DetId(0) ) {
292  return 0;
293  } else {
294  EcalRecHitCollection::const_iterator it = recHits->find( id );
295  if ( it != recHits->end() ) {
296  if( noZS && ( it->checkFlag(EcalRecHit::kTowerRecovered) ||
297  it->checkFlag(EcalRecHit::kWeird) ||
298  (it->detid().subdetId() == EcalBarrel &&
299  it->checkFlag(EcalRecHit::kDiWeird) )
300  )
301  ) {
302  return 0.0;
303  } else {
304  return (*it).energy();
305  }
306  } else {
307  //throw cms::Exception("EcalRecHitNotFound") << "The recHit corresponding to the DetId" << id.rawId() << " not found in the EcalRecHitCollection";
308  // the recHit is not in the collection (hopefully zero suppressed)
309  return 0;
310  }
311  }
312  return 0;
313 }
std::vector< EcalRecHit >::const_iterator const_iterator
const_iterator end() const
Definition: DetId.h:18
iterator find(key_type k)
template<bool noZS>
std::vector< float > EcalClusterToolsT< noZS >::roundnessBarrelSuperClusters ( const reco::SuperCluster superCluster,
const EcalRecHitCollection recHits,
int  weightedPositionMethod = 0,
float  energyThreshold = 0.0 
)
static

Definition at line 1406 of file EcalClusterTools.h.

References edm::SortedCollection< T, SORT >::end(), edm::SortedCollection< T, SORT >::find(), first, reco::CaloCluster::hitsAndFractions(), i, EcalClusterToolsT< noZS >::roundnessSelectedBarrelRecHits(), edm::second(), and groupFilesInBlocks::temp.

1406  {//int positionWeightingMethod=0){
1407  std::vector<std::pair<const EcalRecHit*, float> > RH_ptrs_fracs;
1408  const std::vector< std::pair<DetId, float> >& myHitsPair = superCluster.hitsAndFractions();
1409  for(unsigned int i=0; i< myHitsPair.size(); ++i){
1410  //get pointer to recHit object
1411  EcalRecHitCollection::const_iterator myRH = recHits.find(myHitsPair[i].first);
1412  if( myRH != recHits.end() && myRH->energy()*(noZS ? 1.0 : myHitsPair[i].second) > energyThreshold){
1413  //require rec hit to have positive energy
1414  RH_ptrs_fracs.push_back( std::make_pair(&(*myRH) , myHitsPair[i].second) );
1415  }
1416  }
1417  std::vector<float> temp = EcalClusterToolsT<noZS>::roundnessSelectedBarrelRecHits(RH_ptrs_fracs,weightedPositionMethod);
1418  return temp;
1419 }
int i
Definition: DBlmapReader.cc:9
std::vector< EcalRecHit >::const_iterator const_iterator
const std::vector< std::pair< DetId, float > > & hitsAndFractions() const
Definition: CaloCluster.h:190
U second(std::pair< T, U > const &p)
bool first
Definition: L1TdeRCT.cc:75
const_iterator end() const
static std::vector< float > roundnessSelectedBarrelRecHits(const std::vector< std::pair< const EcalRecHit *, float > > &rhVector, int weightedPositionMethod=0)
iterator find(key_type k)
template<bool noZS>
std::vector< float > EcalClusterToolsT< noZS >::roundnessBarrelSuperClustersUserExtended ( const reco::SuperCluster superCluster,
const EcalRecHitCollection recHits,
int  ieta_delta = 0,
int  iphi_delta = 0,
float  energyRHThresh = 0.00000,
int  weightedPositionMethod = 0 
)
static

Definition at line 1426 of file EcalClusterTools.h.

References funct::abs(), edm::SortedCollection< T, SORT >::begin(), edm::SortedCollection< T, SORT >::end(), edm::SortedCollection< T, SORT >::find(), first, EcalClusterToolsT< noZS >::getSeedPosition(), reco::CaloCluster::hitsAndFractions(), i, EcalClusterToolsT< noZS >::roundnessSelectedBarrelRecHits(), and edm::second().

1426  {
1427 
1428  std::vector<std::pair<const EcalRecHit*, float> > RH_ptrs_fracs;
1429  const std::vector< std::pair<DetId, float> >& myHitsPair = superCluster.hitsAndFractions();
1430  for(unsigned int i=0; i<myHitsPair.size(); ++i){
1431  //get pointer to recHit object
1432  EcalRecHitCollection::const_iterator myRH = recHits.find(myHitsPair[i].first);
1433  if(myRH != recHits.end() && myRH->energy()*(noZS ? 1.0 : myHitsPair[i].second) > energyRHThresh)
1434  RH_ptrs_fracs.push_back( std::make_pair(&(*myRH) , myHitsPair[i].second) );
1435  }
1436 
1437 
1438  std::vector<int> seedPosition = EcalClusterToolsT<noZS>::getSeedPosition( RH_ptrs_fracs );
1439 
1440  for(EcalRecHitCollection::const_iterator rh = recHits.begin(); rh != recHits.end(); rh++){
1441  EBDetId EBdetIdi( rh->detid() );
1442  float the_fraction = 0;
1443  //if(rh != recHits.end())
1444  bool inEtaWindow = ( abs( deltaIEta(seedPosition[0],EBdetIdi.ieta()) ) <= ieta_delta );
1445  bool inPhiWindow = ( abs( deltaIPhi(seedPosition[1],EBdetIdi.iphi()) ) <= iphi_delta );
1446  bool passEThresh = ( rh->energy() > energyRHThresh );
1447  bool alreadyCounted = false;
1448 
1449  // figure out if the rechit considered now is already inside the SC
1450  bool is_SCrh_inside_recHits = false;
1451  for(unsigned int i=0; i<myHitsPair.size(); i++){
1452  EcalRecHitCollection::const_iterator SCrh = recHits.find(myHitsPair[i].first);
1453  if(SCrh != recHits.end()){
1454  the_fraction = myHitsPair[i].second;
1455  is_SCrh_inside_recHits = true;
1456  if( rh->detid() == SCrh->detid() ) alreadyCounted = true;
1457  }
1458  }//for loop over SC's recHits
1459 
1460  if( is_SCrh_inside_recHits && !alreadyCounted && passEThresh && inEtaWindow && inPhiWindow){
1461  RH_ptrs_fracs.push_back( std::make_pair(&(*rh),the_fraction) );
1462  }
1463 
1464  }//for loop over rh
1465  return EcalClusterToolsT<noZS>::roundnessSelectedBarrelRecHits(RH_ptrs_fracs,weightedPositionMethod);
1466 }
int i
Definition: DBlmapReader.cc:9
std::vector< EcalRecHit >::const_iterator const_iterator
const std::vector< std::pair< DetId, float > > & hitsAndFractions() const
Definition: CaloCluster.h:190
static int deltaIEta(int seed_ieta, int rh_ieta)
static std::vector< int > getSeedPosition(const std::vector< std::pair< const EcalRecHit *, float > > &RH_ptrs)
static int deltaIPhi(int seed_iphi, int rh_iphi)
U second(std::pair< T, U > const &p)
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
bool first
Definition: L1TdeRCT.cc:75
const_iterator end() const
static std::vector< float > roundnessSelectedBarrelRecHits(const std::vector< std::pair< const EcalRecHit *, float > > &rhVector, int weightedPositionMethod=0)
iterator find(key_type k)
const_iterator begin() const
template<bool noZS>
std::vector< float > EcalClusterToolsT< noZS >::roundnessSelectedBarrelRecHits ( const std::vector< std::pair< const EcalRecHit *, float > > &  rhVector,
int  weightedPositionMethod = 0 
)
static

Definition at line 1472 of file EcalClusterTools.h.

References cuy::denominator, EcalRecHit::detid(), EcalRecHit::energy(), EcalClusterToolsT< noZS >::getSeedPosition(), EcalClusterToolsT< noZS >::getSumEnergy(), i, fff_deleter::log, max(), groupFilesInBlocks::temp, and histoStyle::weight.

Referenced by EcalClusterToolsT< noZS >::roundnessBarrelSuperClusters(), and EcalClusterToolsT< noZS >::roundnessBarrelSuperClustersUserExtended().

1472  {//int weightedPositionMethod = 0){
1473  //positionWeightingMethod = 0 linear weighting, 1 log energy weighting
1474 
1475  std::vector<float> shapes; // this is the returning vector
1476 
1477  //make sure photon has more than one crystal; else roundness and angle suck
1478  if(RH_ptrs_fracs.size()<2){
1479  shapes.push_back( -3 );
1480  shapes.push_back( -3 );
1481  return shapes;
1482  }
1483 
1484  //Find highest E RH (Seed) and save info, compute sum total energy used
1485  std::vector<int> seedPosition = EcalClusterToolsT<noZS>::getSeedPosition( RH_ptrs_fracs );// *recHits);
1486  int tempInt = seedPosition[0];
1487  if(tempInt <0) tempInt++;
1488  float energyTotal = EcalClusterToolsT<noZS>::getSumEnergy( RH_ptrs_fracs );
1489 
1490  //1st loop over rechits: compute new weighted center position in coordinates centered on seed
1491  float centerIEta = 0.;
1492  float centerIPhi = 0.;
1493  float denominator = 0.;
1494 
1495  for(std::vector<std::pair<const EcalRecHit*,float> >::const_iterator rhf_ptr = RH_ptrs_fracs.begin(); rhf_ptr != RH_ptrs_fracs.end(); rhf_ptr++){
1496  const EcalRecHit* rh_ptr = rhf_ptr->first;
1497  //get iEta, iPhi
1498  EBDetId EBdetIdi( rh_ptr->detid() );
1499  if(fabs(energyTotal) < 0.0001){
1500  // don't /0, bad!
1501  shapes.push_back( -2 );
1502  shapes.push_back( -2 );
1503  return shapes;
1504  }
1505  float rh_energy = rh_ptr->energy() * (noZS ? 1.0 : rhf_ptr->second);
1506  float weight = 0;
1507  if(fabs(weightedPositionMethod)<0.0001){ //linear
1508  weight = rh_energy/energyTotal;
1509  }else{ //logrithmic
1510  weight = std::max(0.0, 4.2 + log(rh_energy/energyTotal));
1511  }
1512  denominator += weight;
1513  centerIEta += weight*deltaIEta(seedPosition[0],EBdetIdi.ieta());
1514  centerIPhi += weight*deltaIPhi(seedPosition[1],EBdetIdi.iphi());
1515  }
1516  if(fabs(denominator) < 0.0001){
1517  // don't /0, bad!
1518  shapes.push_back( -2 );
1519  shapes.push_back( -2 );
1520  return shapes;
1521  }
1522  centerIEta = centerIEta / denominator;
1523  centerIPhi = centerIPhi / denominator;
1524 
1525 
1526  //2nd loop over rechits: compute inertia tensor
1527  TMatrixDSym inertia(2); //initialize 2d inertia tensor
1528  double inertia00 = 0.;
1529  double inertia01 = 0.;// = inertia10 b/c matrix should be symmetric
1530  double inertia11 = 0.;
1531  int i = 0;
1532  for(std::vector<std::pair<const EcalRecHit*,float> >::const_iterator rhf_ptr = RH_ptrs_fracs.begin(); rhf_ptr != RH_ptrs_fracs.end(); rhf_ptr++){
1533  const EcalRecHit* rh_ptr = rhf_ptr->first;
1534  //get iEta, iPhi
1535  EBDetId EBdetIdi( rh_ptr->detid() );
1536 
1537  if(fabs(energyTotal) < 0.0001){
1538  // don't /0, bad!
1539  shapes.push_back( -2 );
1540  shapes.push_back( -2 );
1541  return shapes;
1542  }
1543  float rh_energy = rh_ptr->energy() * (noZS ? 1.0 : rhf_ptr->second);
1544  float weight = 0;
1545  if(fabs(weightedPositionMethod) < 0.0001){ //linear
1546  weight = rh_energy/energyTotal;
1547  }else{ //logrithmic
1548  weight = std::max(0.0, 4.2 + log(rh_energy/energyTotal));
1549  }
1550 
1551  float ieta_rh_to_center = deltaIEta(seedPosition[0],EBdetIdi.ieta()) - centerIEta;
1552  float iphi_rh_to_center = deltaIPhi(seedPosition[1],EBdetIdi.iphi()) - centerIPhi;
1553 
1554  inertia00 += weight*iphi_rh_to_center*iphi_rh_to_center;
1555  inertia01 -= weight*iphi_rh_to_center*ieta_rh_to_center;
1556  inertia11 += weight*ieta_rh_to_center*ieta_rh_to_center;
1557  i++;
1558  }
1559 
1560  inertia[0][0] = inertia00;
1561  inertia[0][1] = inertia01; // use same number here
1562  inertia[1][0] = inertia01; // and here to insure symmetry
1563  inertia[1][1] = inertia11;
1564 
1565 
1566  //step 1 find principal axes of inertia
1567  TMatrixD eVectors(2,2);
1568  TVectorD eValues(2);
1569  //std::cout<<"EcalClusterToolsT<noZS>::showerRoundness- about to compute eVectors"<<std::endl;
1570  eVectors=inertia.EigenVectors(eValues); //ordered highest eV to lowest eV (I checked!)
1571  //and eVectors are in columns of matrix! I checked!
1572  //and they are normalized to 1
1573 
1574 
1575 
1576  //step 2 select eta component of smaller eVal's eVector
1577  TVectorD smallerAxis(2);//easiest to spin SC on this axis (smallest eVal)
1578  smallerAxis[0]=eVectors[0][1];//row,col //eta component
1579  smallerAxis[1]=eVectors[1][1]; //phi component
1580 
1581  //step 3 compute interesting quatities
1582  Double_t temp = fabs(smallerAxis[0]);// closer to 1 ->beamhalo, closer to 0 something else
1583  if(fabs(eValues[0]) < 0.0001){
1584  // don't /0, bad!
1585  shapes.push_back( -2 );
1586  shapes.push_back( -2 );
1587  return shapes;
1588  }
1589 
1590  float Roundness = eValues[1]/eValues[0];
1591  float Angle=acos(temp);
1592 
1593  if( -0.00001 < Roundness && Roundness < 0) Roundness = 0.;
1594  if( -0.00001 < Angle && Angle < 0 ) Angle = 0.;
1595 
1596  shapes.push_back( Roundness );
1597  shapes.push_back( Angle );
1598  return shapes;
1599 
1600 }
int i
Definition: DBlmapReader.cc:9
const DetId & detid() const
Definition: EcalRecHit.h:71
static int deltaIEta(int seed_ieta, int rh_ieta)
static std::vector< int > getSeedPosition(const std::vector< std::pair< const EcalRecHit *, float > > &RH_ptrs)
static int deltaIPhi(int seed_iphi, int rh_iphi)
list denominator
Definition: cuy.py:484
const T & max(const T &a, const T &b)
float energy() const
Definition: EcalRecHit.h:68
static float getSumEnergy(const std::vector< std::pair< const EcalRecHit *, float > > &RH_ptrs_fracs)
Definition: Angle.h:17
int weight
Definition: histoStyle.py:50
template<bool noZS>
std::vector< float > EcalClusterToolsT< noZS >::scLocalCovariances ( const reco::SuperCluster cluster,
const EcalRecHitCollection recHits,
const CaloTopology topology,
float  w0 = 4.7 
)
static

Definition at line 1154 of file EcalClusterTools.h.

References cuy::denominator, dPhi(), EcalBarrel, relval_parameters_module::energy, f, first, cropTnPTrees::frac, CaloTopology::getSubdetectorTopology(), reco::CaloCluster::hitsAndFractions(), i, fff_deleter::log, max(), reco::SuperCluster::seed(), DetId::subdetId(), findQualityFiles::v, and w().

1155 {
1156  const reco::BasicCluster bcluster = *(cluster.seed());
1157 
1158  float e_5x5 = e5x5(bcluster, recHits, topology);
1159  float covEtaEta, covEtaPhi, covPhiPhi;
1160 
1161  if (e_5x5 >= 0.) {
1162  const std::vector<std::pair<DetId, float> >& v_id = cluster.hitsAndFractions();
1163  std::pair<float,float> mean5x5PosInNrCrysFromSeed = mean5x5PositionInLocalCrysCoord(bcluster, recHits, topology);
1164  std::pair<float,float> mean5x5XYPos = mean5x5PositionInXY(cluster,recHits,topology);
1165  // now we can calculate the covariances
1166  double numeratorEtaEta = 0;
1167  double numeratorEtaPhi = 0;
1168  double numeratorPhiPhi = 0;
1169  double denominator = 0;
1170 
1171  const double barrelCrysSize = 0.01745; //approximate size of crystal in eta,phi in barrel
1172  const double endcapCrysSize = 0.0447; //the approximate crystal size sigmaEtaEta was corrected to in the endcap
1173 
1174  DetId seedId = getMaximum(v_id, recHits).first;
1175  bool isBarrel=seedId.subdetId()==EcalBarrel;
1176 
1177  const double crysSize = isBarrel ? barrelCrysSize : endcapCrysSize;
1178 
1179  for (size_t i = 0; i < v_id.size(); ++i) {
1180  CaloNavigator<DetId> cursor = CaloNavigator<DetId>(v_id[i].first, topology->getSubdetectorTopology(v_id[i].first));
1181  float frac = getFraction(v_id,*cursor);
1182  float energy = recHitEnergy(*cursor, recHits)*frac;
1183 
1184  if (energy <= 0) continue;
1185 
1186  float dEta = getNrCrysDiffInEta(*cursor,seedId) - mean5x5PosInNrCrysFromSeed.first;
1187  float dPhi = 0;
1188  if(isBarrel) dPhi = getNrCrysDiffInPhi(*cursor,seedId) - mean5x5PosInNrCrysFromSeed.second;
1189  else dPhi = getDPhiEndcap(*cursor,mean5x5XYPos.first,mean5x5XYPos.second);
1190 
1191 
1192 
1193  double w = 0.;
1194  w = std::max(0.0f, w0 + std::log( energy / e_5x5 ));
1195 
1196  denominator += w;
1197  numeratorEtaEta += w * dEta * dEta;
1198  numeratorEtaPhi += w * dEta * dPhi;
1199  numeratorPhiPhi += w * dPhi * dPhi;
1200  }
1201 
1202  //multiplying by crysSize to make the values compariable to normal covariances
1203  if (denominator != 0.0) {
1204  covEtaEta = crysSize*crysSize* numeratorEtaEta / denominator;
1205  covEtaPhi = crysSize*crysSize* numeratorEtaPhi / denominator;
1206  covPhiPhi = crysSize*crysSize* numeratorPhiPhi / denominator;
1207  } else {
1208  covEtaEta = 999.9;
1209  covEtaPhi = 999.9;
1210  covPhiPhi = 999.9;
1211  }
1212 
1213  } else {
1214  // Warn the user if there was no energy in the cells and return zeroes.
1215  // std::cout << "\ClusterShapeAlgo::Calculate_Covariances: no energy in supplied cells.\n";
1216  covEtaEta = 0;
1217  covEtaPhi = 0;
1218  covPhiPhi = 0;
1219  }
1220 
1221  std::vector<float> v;
1222  v.push_back( covEtaEta );
1223  v.push_back( covEtaPhi );
1224  v.push_back( covPhiPhi );
1225 
1226  return v;
1227 }
int i
Definition: DBlmapReader.cc:9
const std::vector< std::pair< DetId, float > > & hitsAndFractions() const
Definition: CaloCluster.h:190
list denominator
Definition: cuy.py:484
double dPhi(double phi1, double phi2)
Definition: JetUtil.h:30
const T & max(const T &a, const T &b)
static float getDPhiEndcap(const DetId &crysId, float meanX, float meanY)
static std::pair< DetId, float > getMaximum(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits)
double f[11][100]
static float getNrCrysDiffInEta(const DetId &crysId, const DetId &orginId)
static float getFraction(const std::vector< std::pair< DetId, float > > &v_id, DetId id)
bool first
Definition: L1TdeRCT.cc:75
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
static float getNrCrysDiffInPhi(const DetId &crysId, const DetId &orginId)
static std::pair< float, float > mean5x5PositionInXY(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
Definition: DetId.h:18
static float recHitEnergy(DetId id, const EcalRecHitCollection *recHits)
const CaloSubdetectorTopology * getSubdetectorTopology(const DetId &id) const
access the subdetector Topology for the given subdetector directly
Definition: CaloTopology.cc:26
static std::pair< float, float > mean5x5PositionInLocalCrysCoord(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
const CaloClusterPtr & seed() const
seed BasicCluster
Definition: SuperCluster.h:66
T w() const
static float e5x5(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloTopology *topology)
template<bool noZS>
double EcalClusterToolsT< noZS >::zernike20 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloGeometry geometry,
double  R0 = 6.6,
bool  logW = true,
float  w0 = 4.7 
)
static

Definition at line 939 of file EcalClusterTools.h.

940 {
941  return absZernikeMoment( cluster, recHits, geometry, 2, 0, R0, logW, w0 );
942 }
static double absZernikeMoment(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, int n, int m, double R0, bool logW, float w0)
template<bool noZS>
double EcalClusterToolsT< noZS >::zernike42 ( const reco::BasicCluster cluster,
const EcalRecHitCollection recHits,
const CaloGeometry geometry,
double  R0 = 6.6,
bool  logW = true,
float  w0 = 4.7 
)
static

Definition at line 945 of file EcalClusterTools.h.

946 {
947  return absZernikeMoment( cluster, recHits, geometry, 4, 2, R0, logW, w0 );
948 }
static double absZernikeMoment(const reco::BasicCluster &cluster, const EcalRecHitCollection *recHits, const CaloGeometry *geometry, int n, int m, double R0, bool logW, float w0)