EcalBasicClusterLocalContCorrection Class Reference

#include <EcalBasicClusterLocalContCorrection.h>

Public Member Functions
void	checkInit () const
	EcalBasicClusterLocalContCorrection (edm::ConsumesCollector &&cc)
void	init (const edm::EventSetup &es)
float	operator() (const reco::BasicCluster &, const EcalRecHitCollection &) const

Private Member Functions
int	getEcalModule (DetId id) const

Private Attributes
CaloGeometry const *	caloGeometry_
const edm::ESGetToken < CaloGeometry, CaloGeometryRecord >	caloGeometryToken_
EcalClusterLocalContCorrParameters const *	params_
const edm::ESGetToken < EcalClusterLocalContCorrParameters, EcalClusterLocalContCorrParametersRcd >	paramsToken_

Detailed Description

Function to correct em object energy for energy not contained in a 5x5 crystal area in the calorimeter

$Id: EcalBasicClusterLocalContCorrection.h $Date: $Revision:

Author

Federico Ferri, CEA Saclay, November 2008

Definition at line 22 of file EcalBasicClusterLocalContCorrection.h.

Constructor & Destructor Documentation

EcalBasicClusterLocalContCorrection::EcalBasicClusterLocalContCorrection

(

edm::ConsumesCollector &&

cc

)

Definition at line 8 of file EcalBasicClusterLocalContCorrection.cc.

    : paramsToken_{cc.esConsumes()}, caloGeometryToken_{cc.esConsumes()} {}

Member Function Documentation

void EcalBasicClusterLocalContCorrection::checkInit

(

)

const

Definition at line 16 of file EcalBasicClusterLocalContCorrection.cc.

References Exception, and params_.

                                                          {
  if (!params_) {
    // non-initialized function parameters: throw exception
    throw cms::Exception("EcalBasicClusterLocalContCorrection::checkInit()")
        << "Trying to access an uninitialized crack correction function.\n"
           "Please call `init( edm::EventSetup &)' before any use of the function.\n";
  }
}

int EcalBasicClusterLocalContCorrection::getEcalModule ( DetId  id ) const

private

Definition at line 157 of file EcalBasicClusterLocalContCorrection.cc.

References mod().

                                                                     {
  int mod = 0;
  int ieta = (EBDetId(id)).ieta();

  if (fabs(ieta) <= 25)
    mod = 0;
  if (fabs(ieta) > 25 && fabs(ieta) <= 45)
    mod = 1;
  if (fabs(ieta) > 45 && fabs(ieta) <= 65)
    mod = 2;
  if (fabs(ieta) > 65 && fabs(ieta) <= 85)
    mod = 3;

  return (mod);
}

void EcalBasicClusterLocalContCorrection::init

(

const edm::EventSetup &

es

)

Definition at line 11 of file EcalBasicClusterLocalContCorrection.cc.

References caloGeometry_, caloGeometryToken_, edm::EventSetup::getData(), params_, and paramsToken_.

                                                                      {
  params_ = &es.getData(paramsToken_);
  caloGeometry_ = &es.getData(caloGeometryToken_);
}

float EcalBasicClusterLocalContCorrection::operator()	(	const reco::BasicCluster &	basicCluster,
		const EcalRecHitCollection &	recHit
	)		const

Definition at line 28 of file EcalBasicClusterLocalContCorrection.cc.

References funct::abs(), funct::cos(), HLT_FULL_cff::depth, DetId::Ecal, EcalBarrel, reco::CaloCluster::energy(), PV3DBase< T, PVType, FrameType >::eta(), first, relativeConstraints::geom, CaloSubdetectorGeometry::getGeometry(), reco::CaloCluster::hitsAndFractions(), EBDetId::ieta(), EBDetId::iphi(), log, M_PI, colinearityKinematic::Phi, PV3DBase< T, PVType, FrameType >::phi(), reco::CaloCluster::position(), createTree::pp, reco::CaloCluster::seed(), PV3DBase< T, PVType, FrameType >::theta(), and X0.

                                                                                                {
  checkInit();

  // number of parameters needed by this parametrization
  constexpr size_t nparams = 24;

  //correction factor to be returned, and to be calculated in this present function:
  double correction_factor = 1.;
  double fetacor = 1.;  //eta dependent part of the correction factor
  double fphicor = 1.;  //phi dependent part of the correction factor

  //--------------if barrel calculate local position wrt xtal center -------------------
  const CaloSubdetectorGeometry *geom =
      caloGeometry_->getSubdetectorGeometry(DetId::Ecal, EcalBarrel);  //EcalBarrel = 1

  const math::XYZPoint &position_ = basicCluster.position();
  double Theta = -position_.theta() + 0.5 * M_PI;
  double Eta = position_.eta();
  double Phi = TVector2::Phi_mpi_pi(position_.phi());

  //Calculate expected depth of the maximum shower from energy (like in PositionCalc::Calculate_Location()):
  // The parameters X0 and T0 are hardcoded here because these values were used to calculate the corrections:
  constexpr float X0 = 0.89;
  constexpr float T0 = 7.4;
  double depth = X0 * (T0 + log(basicCluster.energy()));

  //search which crystal is closest to the cluster position and call it crystalseed:
  //std::vector<DetId> crystals_vector = *scRef.getHitsByDetId();   //deprecated
  std::vector<std::pair<DetId, float> > crystals_vector = basicCluster.hitsAndFractions();
  float dphimin = 999.;
  float detamin = 999.;
  int ietaclosest = 0;
  int iphiclosest = 0;

  for (unsigned int icry = 0; icry != crystals_vector.size(); ++icry) {
    EBDetId crystal(crystals_vector[icry].first);
    auto cell = geom->getGeometry(crystal);  // problema qui
    GlobalPoint center_pos = cell->getPosition(depth);
    double EtaCentr = center_pos.eta();
    double PhiCentr = TVector2::Phi_mpi_pi(center_pos.phi());
    if (std::abs(EtaCentr - Eta) < detamin) {
      detamin = std::abs(EtaCentr - Eta);
      ietaclosest = crystal.ieta();
    }
    if (std::abs(TVector2::Phi_mpi_pi(PhiCentr - Phi)) < dphimin) {
      dphimin = std::abs(TVector2::Phi_mpi_pi(PhiCentr - Phi));
      iphiclosest = crystal.iphi();
    }
  }

  EBDetId crystalseed(ietaclosest, iphiclosest);

  // Get center cell position from shower depth
  auto cell = geom->getGeometry(crystalseed);
  GlobalPoint center_pos = cell->getPosition(depth);

  //PHI
  double PhiCentr = TVector2::Phi_mpi_pi(center_pos.phi());
  double PhiWidth = (M_PI / 180.);
  double PhiCry = (TVector2::Phi_mpi_pi(Phi - PhiCentr)) / PhiWidth;
  if (PhiCry > 0.5)
    PhiCry = 0.5;
  if (PhiCry < -0.5)
    PhiCry = -0.5;
  //flip to take into account ECAL barrel symmetries:
  if (ietaclosest < 0)
    PhiCry *= -1.;

  //ETA
  double ThetaCentr = -center_pos.theta() + 0.5 * M_PI;
  double ThetaWidth = (M_PI / 180.) * std::cos(ThetaCentr);
  double EtaCry = (Theta - ThetaCentr) / ThetaWidth;
  if (EtaCry > 0.5)
    EtaCry = 0.5;
  if (EtaCry < -0.5)
    EtaCry = -0.5;
  //flip to take into account ECAL barrel symmetries:
  if (ietaclosest < 0)
    EtaCry *= -1.;

  //-------------- end calculate local position -------------

  size_t payloadsize = params_->params().size();

  if (payloadsize < nparams)
    edm::LogError("Invalid Payload") << "Parametrization requires " << nparams << " parameters but only " << payloadsize
                                     << " are found in DB. Perhaps incompatible Global Tag" << std::endl;

  if (payloadsize > nparams)
    edm::LogWarning("Size mismatch ") << "Parametrization requires " << nparams << " parameters but " << payloadsize
                                      << " are found in DB. Perhaps incompatible Global Tag" << std::endl;

  std::pair<double, double> localPosition(EtaCry, PhiCry);

  //--- local cluster coordinates
  float localEta = localPosition.first;
  float localPhi = localPosition.second;

  //--- ecal module
  int imod = getEcalModule(basicCluster.seed());

  //-- corrections parameters
  float pe[3], pp[3];
  pe[0] = (params_->params())[0 + imod * 3];
  pe[1] = (params_->params())[1 + imod * 3];
  pe[2] = (params_->params())[2 + imod * 3];
  pp[0] = (params_->params())[12 + imod * 3];
  pp[1] = (params_->params())[13 + imod * 3];
  pp[2] = (params_->params())[14 + imod * 3];

  //--- correction vs local eta
  fetacor = pe[0] + pe[1] * localEta + pe[2] * localEta * localEta;

  //--- correction vs local phi
  fphicor = pp[0] + pp[1] * localPhi + pp[2] * localPhi * localPhi;

  //if the seed crystal is neighbourgh of a supermodule border, don't apply the phi dependent  containment corrections, but use the larger crack corrections instead.
  int iphimod20 = std::abs(iphiclosest % 20);
  if (iphimod20 <= 1)
    fphicor = 1.;

  correction_factor = (1. / fetacor) * (1. / fphicor);

  //return the correction factor. Use it to multiply the cluster energy.
  return correction_factor;
}

Member Data Documentation

CaloGeometry const* EcalBasicClusterLocalContCorrection::caloGeometry_

private

Definition at line 42 of file EcalBasicClusterLocalContCorrection.h.

Referenced by init().

const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> EcalBasicClusterLocalContCorrection::caloGeometryToken_

private

Definition at line 39 of file EcalBasicClusterLocalContCorrection.h.

Referenced by init().

EcalClusterLocalContCorrParameters const* EcalBasicClusterLocalContCorrection::params_

private

Definition at line 41 of file EcalBasicClusterLocalContCorrection.h.

Referenced by checkInit(), and init().

const edm::ESGetToken<EcalClusterLocalContCorrParameters, EcalClusterLocalContCorrParametersRcd> EcalBasicClusterLocalContCorrection::paramsToken_

private

Definition at line 38 of file EcalBasicClusterLocalContCorrection.h.

Referenced by init().