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