File indexing completed on 2025-03-29 02:43:21
0001 #include "RecoParticleFlow/PFClusterProducer/interface/PFClusterEMEnergyCorrector.h"
0002
0003 #include "vdt/vdtMath.h"
0004 #include <array>
0005
0006 namespace {
0007 typedef reco::PFCluster::EEtoPSAssociation::value_type EEPSPair;
0008 bool sortByKey(const EEPSPair &a, const EEPSPair &b) { return a.first < b.first; }
0009
0010 double getOffset(const double lo, const double hi) { return lo + 0.5 * (hi - lo); }
0011 double getScale(const double lo, const double hi) { return 0.5 * (hi - lo); }
0012 }
0013
0014 PFClusterEMEnergyCorrector::PFClusterEMEnergyCorrector(const edm::ParameterSet &conf, edm::ConsumesCollector &&cc)
0015 : ecalClusterToolsESGetTokens_{cc}, ecalReadoutToolsESGetTokens_{conf, cc}, calibrator_(new PFEnergyCalibration) {
0016 applyCrackCorrections_ = conf.getParameter<bool>("applyCrackCorrections");
0017 applyMVACorrections_ = conf.getParameter<bool>("applyMVACorrections");
0018 srfAwareCorrection_ = conf.getParameter<bool>("srfAwareCorrection");
0019 setEnergyUncertainty_ = conf.getParameter<bool>("setEnergyUncertainty");
0020 maxPtForMVAEvaluation_ = conf.getParameter<double>("maxPtForMVAEvaluation");
0021
0022 if (applyMVACorrections_) {
0023 meanlimlowEB_ = -0.336;
0024 meanlimhighEB_ = 0.916;
0025 meanoffsetEB_ = getOffset(meanlimlowEB_, meanlimhighEB_);
0026 meanscaleEB_ = getScale(meanlimlowEB_, meanlimhighEB_);
0027
0028 meanlimlowEE_ = -0.336;
0029 meanlimhighEE_ = 0.916;
0030 meanoffsetEE_ = getOffset(meanlimlowEE_, meanlimhighEE_);
0031 meanscaleEE_ = getScale(meanlimlowEE_, meanlimhighEE_);
0032
0033 sigmalimlowEB_ = 0.001;
0034 sigmalimhighEB_ = 0.4;
0035 sigmaoffsetEB_ = getOffset(sigmalimlowEB_, sigmalimhighEB_);
0036 sigmascaleEB_ = getScale(sigmalimlowEB_, sigmalimhighEB_);
0037
0038 sigmalimlowEE_ = 0.001;
0039 sigmalimhighEE_ = 0.4;
0040 sigmaoffsetEE_ = getOffset(sigmalimlowEE_, sigmalimhighEE_);
0041 sigmascaleEE_ = getScale(sigmalimlowEE_, sigmalimhighEE_);
0042
0043 recHitsEB_ = cc.consumes<EcalRecHitCollection>(conf.getParameter<edm::InputTag>("recHitsEBLabel"));
0044 recHitsEE_ = cc.consumes<EcalRecHitCollection>(conf.getParameter<edm::InputTag>("recHitsEELabel"));
0045 autoDetectBunchSpacing_ = conf.getParameter<bool>("autoDetectBunchSpacing");
0046
0047 if (autoDetectBunchSpacing_) {
0048 bunchSpacing_ = cc.consumes<unsigned int>(edm::InputTag("bunchSpacingProducer"));
0049 bunchSpacingManual_ = 0;
0050 } else {
0051 bunchSpacingManual_ = conf.getParameter<int>("bunchSpacing");
0052 }
0053
0054 condnames_mean_25ns_.insert(condnames_mean_25ns_.end(),
0055 {"ecalPFClusterCorV2_EB_pfSize1_mean_25ns",
0056 "ecalPFClusterCorV2_EB_pfSize2_mean_25ns",
0057 "ecalPFClusterCorV2_EB_pfSize3_ptbin1_mean_25ns",
0058 "ecalPFClusterCorV2_EB_pfSize3_ptbin2_mean_25ns",
0059 "ecalPFClusterCorV2_EB_pfSize3_ptbin3_mean_25ns",
0060 "ecalPFClusterCorV2_EE_pfSize1_mean_25ns",
0061 "ecalPFClusterCorV2_EE_pfSize2_mean_25ns",
0062 "ecalPFClusterCorV2_EE_pfSize3_ptbin1_mean_25ns",
0063 "ecalPFClusterCorV2_EE_pfSize3_ptbin2_mean_25ns",
0064 "ecalPFClusterCorV2_EE_pfSize3_ptbin3_mean_25ns"});
0065 condnames_sigma_25ns_.insert(condnames_sigma_25ns_.end(),
0066 {"ecalPFClusterCorV2_EB_pfSize1_sigma_25ns",
0067 "ecalPFClusterCorV2_EB_pfSize2_sigma_25ns",
0068 "ecalPFClusterCorV2_EB_pfSize3_ptbin1_sigma_25ns",
0069 "ecalPFClusterCorV2_EB_pfSize3_ptbin2_sigma_25ns",
0070 "ecalPFClusterCorV2_EB_pfSize3_ptbin3_sigma_25ns",
0071 "ecalPFClusterCorV2_EE_pfSize1_sigma_25ns",
0072 "ecalPFClusterCorV2_EE_pfSize2_sigma_25ns",
0073 "ecalPFClusterCorV2_EE_pfSize3_ptbin1_sigma_25ns",
0074 "ecalPFClusterCorV2_EE_pfSize3_ptbin2_sigma_25ns",
0075 "ecalPFClusterCorV2_EE_pfSize3_ptbin3_sigma_25ns"});
0076 condnames_mean_50ns_.insert(condnames_mean_50ns_.end(),
0077 {"ecalPFClusterCorV2_EB_pfSize1_mean_50ns",
0078 "ecalPFClusterCorV2_EB_pfSize2_mean_50ns",
0079 "ecalPFClusterCorV2_EB_pfSize3_ptbin1_mean_50ns",
0080 "ecalPFClusterCorV2_EB_pfSize3_ptbin2_mean_50ns",
0081 "ecalPFClusterCorV2_EB_pfSize3_ptbin3_mean_50ns",
0082 "ecalPFClusterCorV2_EE_pfSize1_mean_50ns",
0083 "ecalPFClusterCorV2_EE_pfSize2_mean_50ns",
0084 "ecalPFClusterCorV2_EE_pfSize3_ptbin1_mean_50ns",
0085 "ecalPFClusterCorV2_EE_pfSize3_ptbin2_mean_50ns",
0086 "ecalPFClusterCorV2_EE_pfSize3_ptbin3_mean_50ns"});
0087 condnames_sigma_50ns_.insert(condnames_sigma_50ns_.end(),
0088 {"ecalPFClusterCorV2_EB_pfSize1_sigma_50ns",
0089 "ecalPFClusterCorV2_EB_pfSize2_sigma_50ns",
0090 "ecalPFClusterCorV2_EB_pfSize3_ptbin1_sigma_50ns",
0091 "ecalPFClusterCorV2_EB_pfSize3_ptbin2_sigma_50ns",
0092 "ecalPFClusterCorV2_EB_pfSize3_ptbin3_sigma_50ns",
0093 "ecalPFClusterCorV2_EE_pfSize1_sigma_50ns",
0094 "ecalPFClusterCorV2_EE_pfSize2_sigma_50ns",
0095 "ecalPFClusterCorV2_EE_pfSize3_ptbin1_sigma_50ns",
0096 "ecalPFClusterCorV2_EE_pfSize3_ptbin2_sigma_50ns",
0097 "ecalPFClusterCorV2_EE_pfSize3_ptbin3_sigma_50ns"});
0098
0099 if (srfAwareCorrection_) {
0100 sigmalimlowEE_ = 0.001;
0101 sigmalimhighEE_ = 0.1;
0102 sigmaoffsetEE_ = getOffset(sigmalimlowEE_, sigmalimhighEE_);
0103 sigmascaleEE_ = getScale(sigmalimlowEE_, sigmalimhighEE_);
0104
0105 ebSrFlagToken_ = cc.consumes<EBSrFlagCollection>(conf.getParameter<edm::InputTag>("ebSrFlagLabel"));
0106 eeSrFlagToken_ = cc.consumes<EESrFlagCollection>(conf.getParameter<edm::InputTag>("eeSrFlagLabel"));
0107
0108 condnames_mean_.insert(condnames_mean_.end(),
0109 {"ecalPFClusterCor2017V2_EB_ZS_mean_25ns",
0110 "ecalPFClusterCor2017V2_EB_Full_ptbin1_mean_25ns",
0111 "ecalPFClusterCor2017V2_EB_Full_ptbin2_mean_25ns",
0112 "ecalPFClusterCor2017V2_EB_Full_ptbin3_mean_25ns",
0113 "ecalPFClusterCor2017V2_EE_ZS_mean_25ns",
0114 "ecalPFClusterCor2017V2_EE_Full_ptbin1_mean_25ns",
0115 "ecalPFClusterCor2017V2_EE_Full_ptbin2_mean_25ns",
0116 "ecalPFClusterCor2017V2_EE_Full_ptbin3_mean_25ns"});
0117
0118 condnames_sigma_.insert(condnames_sigma_.end(),
0119 {"ecalPFClusterCor2017V2_EB_ZS_sigma_25ns",
0120 "ecalPFClusterCor2017V2_EB_Full_ptbin1_sigma_25ns",
0121 "ecalPFClusterCor2017V2_EB_Full_ptbin2_sigma_25ns",
0122 "ecalPFClusterCor2017V2_EB_Full_ptbin3_sigma_25ns",
0123 "ecalPFClusterCor2017V2_EE_ZS_sigma_25ns",
0124 "ecalPFClusterCor2017V2_EE_Full_ptbin1_sigma_25ns",
0125 "ecalPFClusterCor2017V2_EE_Full_ptbin2_sigma_25ns",
0126 "ecalPFClusterCor2017V2_EE_Full_ptbin3_sigma_25ns"});
0127
0128 for (short i = 0; i < (short)condnames_mean_.size(); i++) {
0129 forestMeanTokens_25ns_.emplace_back(cc.esConsumes(edm::ESInputTag("", condnames_mean_[i])));
0130 forestSigmaTokens_25ns_.emplace_back(cc.esConsumes(edm::ESInputTag("", condnames_sigma_[i])));
0131 }
0132 } else {
0133 for (short i = 0; i < (short)condnames_mean_25ns_.size(); i++) {
0134 forestMeanTokens_25ns_.emplace_back(cc.esConsumes(edm::ESInputTag("", condnames_mean_25ns_[i])));
0135 forestSigmaTokens_25ns_.emplace_back(cc.esConsumes(edm::ESInputTag("", condnames_sigma_25ns_[i])));
0136 }
0137 for (short i = 0; i < (short)condnames_mean_50ns_.size(); i++) {
0138 forestMeanTokens_50ns_.emplace_back(cc.esConsumes(edm::ESInputTag("", condnames_mean_50ns_[i])));
0139 forestSigmaTokens_50ns_.emplace_back(cc.esConsumes(edm::ESInputTag("", condnames_sigma_50ns_[i])));
0140 }
0141 }
0142 }
0143 }
0144
0145 void PFClusterEMEnergyCorrector::correctEnergies(const edm::Event &evt,
0146 const edm::EventSetup &es,
0147 const reco::PFCluster::EEtoPSAssociation &assoc,
0148 reco::PFClusterCollection &cs) {
0149
0150
0151 if (!applyMVACorrections_) {
0152 for (unsigned int idx = 0; idx < cs.size(); ++idx) {
0153 reco::PFCluster &cluster = cs[idx];
0154 bool iseb = cluster.layer() == PFLayer::ECAL_BARREL;
0155 float ePS1 = 0., ePS2 = 0.;
0156 if (!iseb)
0157 getAssociatedPSEnergy(idx, assoc, ePS1, ePS2);
0158 double correctedEnergy = calibrator_->energyEm(cluster, ePS1, ePS2, applyCrackCorrections_);
0159 cluster.setCorrectedEnergy(correctedEnergy);
0160 }
0161 return;
0162 }
0163
0164
0165 EcalClusterLazyTools lazyTool(evt, ecalClusterToolsESGetTokens_.get(es), recHitsEB_, recHitsEE_);
0166 EcalReadoutTools readoutTool(evt, es, ecalReadoutToolsESGetTokens_);
0167
0168 if (!srfAwareCorrection_) {
0169 int bunchspacing = 450;
0170 if (autoDetectBunchSpacing_) {
0171 edm::Handle<unsigned int> bunchSpacingH;
0172 evt.getByToken(bunchSpacing_, bunchSpacingH);
0173 bunchspacing = *bunchSpacingH;
0174 } else {
0175 bunchspacing = bunchSpacingManual_;
0176 }
0177
0178 const unsigned int ncor = (bunchspacing == 25) ? condnames_mean_25ns_.size() : condnames_mean_50ns_.size();
0179
0180 std::vector<edm::ESHandle<GBRForestD> > forestH_mean(ncor);
0181 std::vector<edm::ESHandle<GBRForestD> > forestH_sigma(ncor);
0182
0183 if (bunchspacing == 25) {
0184 for (unsigned int icor = 0; icor < ncor; ++icor) {
0185 forestH_mean[icor] = es.getHandle(forestMeanTokens_25ns_[icor]);
0186 forestH_sigma[icor] = es.getHandle(forestSigmaTokens_25ns_[icor]);
0187 }
0188 } else {
0189 for (unsigned int icor = 0; icor < ncor; ++icor) {
0190 forestH_mean[icor] = es.getHandle(forestMeanTokens_50ns_[icor]);
0191 forestH_sigma[icor] = es.getHandle(forestSigmaTokens_50ns_[icor]);
0192 }
0193 }
0194
0195 std::array<float, 5> eval;
0196 for (unsigned int idx = 0; idx < cs.size(); ++idx) {
0197 reco::PFCluster &cluster = cs[idx];
0198 bool iseb = cluster.layer() == PFLayer::ECAL_BARREL;
0199 float ePS1 = 0., ePS2 = 0.;
0200 if (!iseb)
0201 getAssociatedPSEnergy(idx, assoc, ePS1, ePS2);
0202
0203 double e = cluster.energy();
0204 double pt = cluster.pt();
0205 double evale = e;
0206 if (maxPtForMVAEvaluation_ > 0. && pt > maxPtForMVAEvaluation_) {
0207 evale *= maxPtForMVAEvaluation_ / pt;
0208 }
0209 double invE = (e == 0.) ? 0. : 1. / e;
0210 int size = lazyTool.n5x5(cluster);
0211
0212 int ietaix = 0;
0213 int iphiiy = 0;
0214 if (iseb) {
0215 EBDetId ebseed(cluster.seed());
0216 ietaix = ebseed.ieta();
0217 iphiiy = ebseed.iphi();
0218 } else {
0219 EEDetId eeseed(cluster.seed());
0220 ietaix = eeseed.ix();
0221 iphiiy = eeseed.iy();
0222 }
0223
0224
0225 int coridx = std::min(size, 3) - 1;
0226 if (coridx == 2) {
0227 if (pt > 4.5) {
0228 coridx += 1;
0229 }
0230 if (pt > 18.) {
0231 coridx += 1;
0232 }
0233 }
0234 if (!iseb) {
0235 coridx += 5;
0236 }
0237
0238 const GBRForestD &meanforest = *forestH_mean[coridx].product();
0239 const GBRForestD &sigmaforest = *forestH_sigma[coridx].product();
0240
0241
0242 eval[0] = evale;
0243 eval[1] = ietaix;
0244 eval[2] = iphiiy;
0245 if (!iseb) {
0246 eval[3] = ePS1 * invE;
0247 eval[4] = ePS2 * invE;
0248 }
0249
0250
0251 double rawmean = meanforest.GetResponse(eval.data());
0252 double rawsigma = sigmaforest.GetResponse(eval.data());
0253
0254
0255 double mean = iseb ? meanoffsetEB_ + meanscaleEB_ * vdt::fast_sin(rawmean)
0256 : meanoffsetEE_ + meanscaleEE_ * vdt::fast_sin(rawmean);
0257 double sigma = iseb ? sigmaoffsetEB_ + sigmascaleEB_ * vdt::fast_sin(rawsigma)
0258 : sigmaoffsetEE_ + sigmascaleEE_ * vdt::fast_sin(rawsigma);
0259
0260
0261
0262 double ecor = vdt::fast_exp(mean) * e;
0263 double sigmacor = sigma * ecor;
0264
0265 cluster.setCorrectedEnergy(ecor);
0266 if (setEnergyUncertainty_)
0267 cluster.setCorrectedEnergyUncertainty(sigmacor);
0268 else
0269 cluster.setCorrectedEnergyUncertainty(0.);
0270 }
0271 return;
0272 }
0273
0274
0275 edm::Handle<EBSrFlagCollection> ebSrFlags;
0276 edm::Handle<EESrFlagCollection> eeSrFlags;
0277 evt.getByToken(ebSrFlagToken_, ebSrFlags);
0278 evt.getByToken(eeSrFlagToken_, eeSrFlags);
0279 if (not ebSrFlags.isValid() or not eeSrFlags.isValid())
0280 edm::LogInfo("PFClusterEMEnergyCorrector") << "SrFlagCollection information is not available. The ECAL PFCluster "
0281 "corrections will assume \"full readout\" for all hits.";
0282
0283 const unsigned int ncor = forestMeanTokens_25ns_.size();
0284 std::vector<edm::ESHandle<GBRForestD> > forestH_mean(ncor);
0285 std::vector<edm::ESHandle<GBRForestD> > forestH_sigma(ncor);
0286
0287 for (unsigned int icor = 0; icor < ncor; ++icor) {
0288 forestH_mean[icor] = es.getHandle(forestMeanTokens_25ns_[icor]);
0289 forestH_sigma[icor] = es.getHandle(forestSigmaTokens_25ns_[icor]);
0290 }
0291
0292 std::array<float, 6> evalEB;
0293 std::array<float, 5> evalEE;
0294
0295 for (unsigned int idx = 0; idx < cs.size(); ++idx) {
0296 reco::PFCluster &cluster = cs[idx];
0297 bool iseb = cluster.layer() == PFLayer::ECAL_BARREL;
0298 float ePS1 = 0., ePS2 = 0.;
0299 if (!iseb)
0300 getAssociatedPSEnergy(idx, assoc, ePS1, ePS2);
0301
0302 double e = cluster.energy();
0303 double pt = cluster.pt();
0304 double evale = e;
0305 if (maxPtForMVAEvaluation_ > 0. && pt > maxPtForMVAEvaluation_) {
0306 evale *= maxPtForMVAEvaluation_ / pt;
0307 }
0308 double invE = (e == 0.) ? 0. : 1. / e;
0309 int size = lazyTool.n5x5(cluster);
0310 int reducedHits = size;
0311 if (size >= 3)
0312 reducedHits = 3;
0313
0314 int ietaix = 0;
0315 int iphiiy = 0;
0316 if (iseb) {
0317 EBDetId ebseed(cluster.seed());
0318 ietaix = ebseed.ieta();
0319 iphiiy = ebseed.iphi();
0320 } else {
0321 EEDetId eeseed(cluster.seed());
0322 ietaix = eeseed.ix();
0323 iphiiy = eeseed.iy();
0324 }
0325
0326
0327 int signeta = (ietaix > 0) ? 1 : -1;
0328 int ietamod20 = (std::abs(ietaix) < 26) ? ietaix - signeta : (ietaix - 26 * signeta) % 20;
0329 int iphimod20 = (iphiiy - 1) % 20;
0330
0331
0332 int clusFlag = 3;
0333 if (iseb) {
0334 if (ebSrFlags.isValid()) {
0335 auto ecalUnit = readoutTool.readOutUnitOf(static_cast<EBDetId>(cluster.seed()));
0336 EBSrFlagCollection::const_iterator srf = ebSrFlags->find(ecalUnit);
0337 if (srf != ebSrFlags->end())
0338 clusFlag = srf->value();
0339 }
0340 } else {
0341 if (eeSrFlags.isValid()) {
0342 auto ecalUnit = readoutTool.readOutUnitOf(static_cast<EEDetId>(cluster.seed()));
0343 EESrFlagCollection::const_iterator srf = eeSrFlags->find(ecalUnit);
0344 if (srf != eeSrFlags->end())
0345 clusFlag = srf->value();
0346 }
0347 }
0348
0349
0350 int coridx = 0;
0351 int regind = 0;
0352 if (!iseb)
0353 regind = 4;
0354
0355
0356
0357
0358
0359
0360
0361
0362
0363
0364 int ZS_bit = clusFlag >> 0 & 1;
0365 int FR_bit = clusFlag >> 1 & 1;
0366
0367 if (ZS_bit != 0 && FR_bit == 0)
0368 coridx = 0 + regind;
0369 else {
0370 if (pt < 2.5)
0371 coridx = 1 + regind;
0372 else if (pt >= 2.5 && pt < 6.)
0373 coridx = 2 + regind;
0374 else if (pt >= 6.)
0375 coridx = 3 + regind;
0376 }
0377 if (ZS_bit == 0 || clusFlag > 7) {
0378 edm::LogWarning("PFClusterEMEnergyCorrector")
0379 << "We can only correct regions readout in ZS (flag 1,5) or FULL readout (flag 3,7). Flag " << clusFlag
0380 << " is not recognized."
0381 << "\n"
0382 << "Assuming FULL readout and continuing";
0383 }
0384
0385 const GBRForestD &meanforest = *forestH_mean[coridx].product();
0386 const GBRForestD &sigmaforest = *forestH_sigma[coridx].product();
0387
0388
0389 if (iseb) {
0390 evalEB[0] = evale;
0391 evalEB[1] = ietaix;
0392 evalEB[2] = iphiiy;
0393 evalEB[3] = ietamod20;
0394 evalEB[4] = iphimod20;
0395 evalEB[5] = reducedHits;
0396 } else {
0397 evalEE[0] = evale;
0398 evalEE[1] = ietaix;
0399 evalEE[2] = iphiiy;
0400 evalEE[3] = (ePS1 + ePS2) * invE;
0401 evalEE[4] = reducedHits;
0402 }
0403
0404
0405 double rawmean = 1;
0406 double rawsigma = 0;
0407
0408 if (iseb) {
0409 rawmean = meanforest.GetResponse(evalEB.data());
0410 rawsigma = sigmaforest.GetResponse(evalEB.data());
0411 } else {
0412 rawmean = meanforest.GetResponse(evalEE.data());
0413 rawsigma = sigmaforest.GetResponse(evalEE.data());
0414 }
0415
0416
0417
0418
0419 double mean = iseb ? meanoffsetEB_ + meanscaleEB_ * vdt::fast_sin(rawmean)
0420 : meanoffsetEE_ + meanscaleEE_ * vdt::fast_sin(rawmean);
0421 double sigma = iseb ? sigmaoffsetEB_ + sigmascaleEB_ * vdt::fast_sin(rawsigma)
0422 : sigmaoffsetEE_ + sigmascaleEE_ * vdt::fast_sin(rawsigma);
0423
0424
0425
0426 double ecor = iseb ? vdt::fast_exp(mean) * e : vdt::fast_exp(mean) * (e + ePS1 + ePS2);
0427 double sigmacor = sigma * ecor;
0428
0429 LogDebug("PFClusterEMEnergyCorrector")
0430 << "ieta : iphi : ietamod20 : iphimod20 : size : reducedHits = " << ietaix << " " << iphiiy << " " << ietamod20
0431 << " " << iphimod20 << " " << size << " " << reducedHits << "\n"
0432 << "isEB : eraw : ePS1 : ePS2 : (eps1+eps2)/raw : Flag = " << iseb << " " << evale << " " << ePS1 << " " << ePS2
0433 << " " << (ePS1 + ePS2) / evale << " " << clusFlag << "\n"
0434 << "response : correction = " << exp(mean) << " " << ecor;
0435
0436 cluster.setCorrectedEnergy(ecor);
0437 if (setEnergyUncertainty_)
0438 cluster.setCorrectedEnergyUncertainty(sigmacor);
0439 else
0440 cluster.setCorrectedEnergyUncertainty(0.);
0441 }
0442 }
0443
0444 void PFClusterEMEnergyCorrector::getAssociatedPSEnergy(const size_t clusIdx,
0445 const reco::PFCluster::EEtoPSAssociation &assoc,
0446 float &e1,
0447 float &e2) {
0448 e1 = 0;
0449 e2 = 0;
0450 auto ee_key_val = std::make_pair(clusIdx, edm::Ptr<reco::PFCluster>());
0451 const auto clustops = std::equal_range(assoc.begin(), assoc.end(), ee_key_val, sortByKey);
0452 for (auto i_ps = clustops.first; i_ps != clustops.second; ++i_ps) {
0453 edm::Ptr<reco::PFCluster> psclus(i_ps->second);
0454 switch (psclus->layer()) {
0455 case PFLayer::PS1:
0456 e1 += psclus->energy();
0457 break;
0458 case PFLayer::PS2:
0459 e2 += psclus->energy();
0460 break;
0461 default:
0462 break;
0463 }
0464 }
0465 }