File indexing completed on 2022-03-09 00:50:55
0001 #include "DataFormats/EgammaCandidates/interface/Photon.h"
0002 #include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
0003
0004 using namespace reco;
0005
0006 Photon::Photon(const LorentzVector& p4, const Point& caloPos, const PhotonCoreRef& core, const Point& vtx)
0007 : RecoCandidate(0, p4, vtx, 22),
0008 caloPosition_(caloPos),
0009 photonCore_(core),
0010 pixelSeed_(false),
0011 haloTaggerMVAVal_(99) {}
0012
0013 Photon::Photon(const Photon& rhs)
0014 : RecoCandidate(rhs),
0015 caloPosition_(rhs.caloPosition_),
0016 photonCore_(rhs.photonCore_),
0017 pixelSeed_(rhs.pixelSeed_),
0018 fiducialFlagBlock_(rhs.fiducialFlagBlock_),
0019 isolationR04_(rhs.isolationR04_),
0020 isolationR03_(rhs.isolationR03_),
0021 showerShapeBlock_(rhs.showerShapeBlock_),
0022 full5x5_showerShapeBlock_(rhs.full5x5_showerShapeBlock_),
0023 saturationInfo_(rhs.saturationInfo_),
0024 eCorrections_(rhs.eCorrections_),
0025 mipVariableBlock_(rhs.mipVariableBlock_),
0026 pfIsolation_(rhs.pfIsolation_),
0027 haloTaggerMVAVal_(rhs.haloTaggerMVAVal_) {}
0028
0029 Photon::~Photon() {}
0030
0031 Photon* Photon::clone() const { return new Photon(*this); }
0032
0033 bool Photon::overlap(const Candidate& c) const {
0034 const RecoCandidate* o = dynamic_cast<const RecoCandidate*>(&c);
0035 return (o != nullptr && (checkOverlap(superCluster(), o->superCluster())));
0036 return false;
0037 }
0038
0039 void Photon::setVertex(const Point& vertex) {
0040 math::XYZVectorF direction = caloPosition() - vertex;
0041 double energy = this->energy();
0042 math::XYZVectorF momentum = direction.unit() * energy;
0043 math::XYZTLorentzVector lv(momentum.x(), momentum.y(), momentum.z(), energy);
0044 setP4(lv);
0045 LeafCandidate::setVertex(vertex);
0046 }
0047
0048 reco::SuperClusterRef Photon::superCluster() const { return this->photonCore()->superCluster(); }
0049
0050 int Photon::conversionTrackProvenance(const edm::RefToBase<reco::Track>& convTrack) const {
0051 const reco::ConversionRefVector& conv2leg = this->photonCore()->conversions();
0052 const reco::ConversionRefVector& conv1leg = this->photonCore()->conversionsOneLeg();
0053
0054 int origin = -1;
0055 bool isEg = false, isPf = false;
0056
0057 for (unsigned iConv = 0; iConv < conv2leg.size(); iConv++) {
0058 std::vector<edm::RefToBase<reco::Track> > convtracks = conv2leg[iConv]->tracks();
0059 for (unsigned itk = 0; itk < convtracks.size(); itk++) {
0060 if (convTrack == convtracks[itk])
0061 isEg = true;
0062 }
0063 }
0064
0065 for (unsigned iConv = 0; iConv < conv1leg.size(); iConv++) {
0066 std::vector<edm::RefToBase<reco::Track> > convtracks = conv1leg[iConv]->tracks();
0067 for (unsigned itk = 0; itk < convtracks.size(); itk++) {
0068 if (convTrack == convtracks[itk])
0069 isPf = true;
0070 }
0071 }
0072
0073 if (isEg)
0074 origin = egamma;
0075 if (isPf)
0076 origin = pflow;
0077 if (isEg && isPf)
0078 origin = both;
0079
0080 return origin;
0081 }
0082
0083 void Photon::setCorrectedEnergy(P4type type, float newEnergy, float delta_e, bool setToRecoCandidate) {
0084 math::XYZTLorentzVectorD newP4 = p4();
0085 newP4 *= newEnergy / newP4.e();
0086 switch (type) {
0087 case ecal_standard:
0088 eCorrections_.scEcalEnergy = newEnergy;
0089 eCorrections_.scEcalEnergyError = delta_e;
0090 break;
0091 case ecal_photons:
0092 eCorrections_.phoEcalEnergy = newEnergy;
0093 eCorrections_.phoEcalEnergyError = delta_e;
0094 break;
0095 case regression1:
0096 eCorrections_.regression1Energy = newEnergy;
0097 eCorrections_.regression1EnergyError = delta_e;
0098 [[fallthrough]];
0099 case regression2:
0100 eCorrections_.regression2Energy = newEnergy;
0101 eCorrections_.regression2EnergyError = delta_e;
0102 break;
0103 default:
0104 throw cms::Exception("reco::Photon") << "unexpected p4 type: " << type;
0105 }
0106 setP4(type, newP4, delta_e, setToRecoCandidate);
0107 }
0108
0109 float Photon::getCorrectedEnergy(P4type type) const {
0110 switch (type) {
0111 case ecal_standard:
0112 return eCorrections_.scEcalEnergy;
0113 break;
0114 case ecal_photons:
0115 return eCorrections_.phoEcalEnergy;
0116 break;
0117 case regression1:
0118 return eCorrections_.regression1Energy;
0119 case regression2:
0120 return eCorrections_.regression2Energy;
0121 break;
0122 default:
0123 throw cms::Exception("reco::Photon") << "unexpected p4 type " << type << " cannot return the energy value: ";
0124 }
0125 }
0126
0127 float Photon::getCorrectedEnergyError(P4type type) const {
0128 switch (type) {
0129 case ecal_standard:
0130 return eCorrections_.scEcalEnergyError;
0131 break;
0132 case ecal_photons:
0133 return eCorrections_.phoEcalEnergyError;
0134 break;
0135 case regression1:
0136 return eCorrections_.regression1EnergyError;
0137 case regression2:
0138 return eCorrections_.regression2EnergyError;
0139 break;
0140 default:
0141 throw cms::Exception("reco::Photon")
0142 << "unexpected p4 type " << type << " cannot return the uncertainty on the energy: ";
0143 }
0144 }
0145
0146 void Photon::setP4(P4type type, const LorentzVector& p4, float error, bool setToRecoCandidate) {
0147 switch (type) {
0148 case ecal_standard:
0149 eCorrections_.scEcalP4 = p4;
0150 eCorrections_.scEcalEnergyError = error;
0151 break;
0152 case ecal_photons:
0153 eCorrections_.phoEcalP4 = p4;
0154 eCorrections_.phoEcalEnergyError = error;
0155 break;
0156 case regression1:
0157 eCorrections_.regression1P4 = p4;
0158 eCorrections_.regression1EnergyError = error;
0159 [[fallthrough]];
0160 case regression2:
0161 eCorrections_.regression2P4 = p4;
0162 eCorrections_.regression2EnergyError = error;
0163 break;
0164 default:
0165 throw cms::Exception("reco::Photon") << "unexpected p4 type: " << type;
0166 }
0167 if (setToRecoCandidate) {
0168 setP4(p4);
0169 eCorrections_.candidateP4type = type;
0170 }
0171 }
0172
0173 const Candidate::LorentzVector& Photon::p4(P4type type) const {
0174 switch (type) {
0175 case ecal_standard:
0176 return eCorrections_.scEcalP4;
0177 case ecal_photons:
0178 return eCorrections_.phoEcalP4;
0179 case regression1:
0180 return eCorrections_.regression1P4;
0181 case regression2:
0182 return eCorrections_.regression2P4;
0183 default:
0184 throw cms::Exception("reco::Photon") << "unexpected p4 type: " << type << " cannot return p4 ";
0185 }
0186 }
0187
0188 void Photon::hcalToRun2EffDepth() {
0189 auto& ss1 = showerShapeBlock_;
0190 auto& ss2 = full5x5_showerShapeBlock_;
0191 auto& iv1 = isolationR03_;
0192 auto& iv2 = isolationR04_;
0193
0194 for (uint id = 2u; id < ss1.hcalOverEcal.size(); ++id) {
0195 ss1.hcalOverEcal[1] += ss1.hcalOverEcal[id];
0196 ss1.hcalOverEcalBc[1] += ss1.hcalOverEcalBc[id];
0197
0198 ss1.hcalOverEcal[id] = 0.f;
0199 ss1.hcalOverEcalBc[id] = 0.f;
0200
0201 ss2.hcalOverEcal[1] += ss2.hcalOverEcal[id];
0202 ss2.hcalOverEcalBc[1] += ss2.hcalOverEcalBc[id];
0203
0204 ss2.hcalOverEcal[id] = 0.f;
0205 ss2.hcalOverEcalBc[id] = 0.f;
0206
0207 iv1.hcalRecHitSumEt[1] += iv1.hcalRecHitSumEt[id];
0208 iv1.hcalRecHitSumEtBc[1] += iv1.hcalRecHitSumEtBc[id];
0209
0210 iv1.hcalRecHitSumEt[id] = 0.f;
0211 iv1.hcalRecHitSumEtBc[id] = 0.f;
0212
0213 iv2.hcalRecHitSumEt[1] += iv2.hcalRecHitSumEt[id];
0214 iv2.hcalRecHitSumEtBc[1] += iv2.hcalRecHitSumEtBc[id];
0215
0216 iv2.hcalRecHitSumEt[id] = 0.f;
0217 iv2.hcalRecHitSumEtBc[id] = 0.f;
0218 }
0219 }