File indexing completed on 2024-04-06 12:10:19
0001 #include "EgammaAnalysis/ElectronTools/interface/ElectronEnergyRegressionEvaluate.h"
0002 #include <cmath>
0003 #include <cassert>
0004
0005 #ifndef STANDALONE
0006 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0007 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0008 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0009 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0010 #include "DataFormats/EcalDetId/interface/EEDetId.h"
0011 #include "FWCore/ParameterSet/interface/FileInPath.h"
0012 #endif
0013
0014 ElectronEnergyRegressionEvaluate::ElectronEnergyRegressionEvaluate()
0015 : fIsInitialized(kFALSE),
0016 fVersionType(kNoTrkVar),
0017 forestCorrection_eb(nullptr),
0018 forestCorrection_ee(nullptr),
0019 forestUncertainty_eb(nullptr),
0020 forestUncertainty_ee(nullptr) {}
0021
0022 ElectronEnergyRegressionEvaluate::~ElectronEnergyRegressionEvaluate() {}
0023
0024
0025 void ElectronEnergyRegressionEvaluate::initialize(std::string weightsFile,
0026 ElectronEnergyRegressionEvaluate::ElectronEnergyRegressionType type) {
0027
0028 TFile file(edm::FileInPath(weightsFile.c_str()).fullPath().c_str());
0029
0030 forestCorrection_eb = (GBRForest *)file.Get("EBCorrection");
0031 forestCorrection_ee = (GBRForest *)file.Get("EECorrection");
0032 forestUncertainty_eb = (GBRForest *)file.Get("EBUncertainty");
0033 forestUncertainty_ee = (GBRForest *)file.Get("EEUncertainty");
0034
0035
0036 assert(forestCorrection_eb);
0037 assert(forestCorrection_ee);
0038 assert(forestUncertainty_eb);
0039 assert(forestUncertainty_ee);
0040
0041
0042 fVersionType = type;
0043 fIsInitialized = kTRUE;
0044 }
0045
0046 #ifndef STANDALONE
0047 double ElectronEnergyRegressionEvaluate::calculateRegressionEnergy(
0048 const reco::GsfElectron *ele, SuperClusterHelper &mySCHelper, double rho, double nvertices, bool printDebug) {
0049 if (!fIsInitialized) {
0050 std::cout << "Error: Electron Energy Regression has not been initialized yet. return 0. \n";
0051 return 0;
0052 }
0053
0054 if (printDebug) {
0055 std::cout << "Regression Type: " << fVersionType << std::endl;
0056 std::cout << "Electron : " << ele->pt() << " " << ele->eta() << " " << ele->phi() << "\n";
0057 }
0058
0059 if (fVersionType == kNoTrkVar) {
0060 return regressionValueNoTrkVar(mySCHelper.rawEnergy(),
0061 mySCHelper.eta(),
0062 mySCHelper.phi(),
0063 mySCHelper.r9(),
0064 mySCHelper.etaWidth(),
0065 mySCHelper.phiWidth(),
0066 mySCHelper.clustersSize(),
0067 mySCHelper.hadronicOverEm(),
0068 rho,
0069 nvertices,
0070 mySCHelper.seedEta(),
0071 mySCHelper.seedPhi(),
0072 mySCHelper.seedEnergy(),
0073 mySCHelper.e3x3(),
0074 mySCHelper.e5x5(),
0075 mySCHelper.sigmaIetaIeta(),
0076 mySCHelper.spp(),
0077 mySCHelper.sep(),
0078 mySCHelper.eMax(),
0079 mySCHelper.e2nd(),
0080 mySCHelper.eTop(),
0081 mySCHelper.eBottom(),
0082 mySCHelper.eLeft(),
0083 mySCHelper.eRight(),
0084 mySCHelper.e2x5Max(),
0085 mySCHelper.e2x5Top(),
0086 mySCHelper.e2x5Bottom(),
0087 mySCHelper.e2x5Left(),
0088 mySCHelper.e2x5Right(),
0089 mySCHelper.ietaSeed(),
0090 mySCHelper.iphiSeed(),
0091 mySCHelper.etaCrySeed(),
0092 mySCHelper.phiCrySeed(),
0093 mySCHelper.preshowerEnergyOverRaw(),
0094 printDebug);
0095 }
0096 if (fVersionType == kWithSubCluVar) {
0097 return regressionValueWithSubClusters(mySCHelper.rawEnergy(),
0098 mySCHelper.eta(),
0099 mySCHelper.phi(),
0100 mySCHelper.r9(),
0101 mySCHelper.etaWidth(),
0102 mySCHelper.phiWidth(),
0103 mySCHelper.clustersSize(),
0104 mySCHelper.hadronicOverEm(),
0105 rho,
0106 nvertices,
0107 mySCHelper.seedEta(),
0108 mySCHelper.seedPhi(),
0109 mySCHelper.seedEnergy(),
0110 mySCHelper.e3x3(),
0111 mySCHelper.e5x5(),
0112 mySCHelper.sigmaIetaIeta(),
0113 mySCHelper.spp(),
0114 mySCHelper.sep(),
0115 mySCHelper.eMax(),
0116 mySCHelper.e2nd(),
0117 mySCHelper.eTop(),
0118 mySCHelper.eBottom(),
0119 mySCHelper.eLeft(),
0120 mySCHelper.eRight(),
0121 mySCHelper.e2x5Max(),
0122 mySCHelper.e2x5Top(),
0123 mySCHelper.e2x5Bottom(),
0124 mySCHelper.e2x5Left(),
0125 mySCHelper.e2x5Right(),
0126 mySCHelper.ietaSeed(),
0127 mySCHelper.iphiSeed(),
0128 mySCHelper.etaCrySeed(),
0129 mySCHelper.phiCrySeed(),
0130 mySCHelper.preshowerEnergyOverRaw(),
0131 ele->ecalDrivenSeed(),
0132 ele->isEBEtaGap(),
0133 ele->isEBPhiGap(),
0134 ele->isEEDeeGap(),
0135 mySCHelper.eSubClusters(),
0136 mySCHelper.subClusterEnergy(1),
0137 mySCHelper.subClusterEta(1),
0138 mySCHelper.subClusterPhi(1),
0139 mySCHelper.subClusterEmax(1),
0140 mySCHelper.subClusterE3x3(1),
0141 mySCHelper.subClusterEnergy(2),
0142 mySCHelper.subClusterEta(2),
0143 mySCHelper.subClusterPhi(2),
0144 mySCHelper.subClusterEmax(2),
0145 mySCHelper.subClusterE3x3(2),
0146 mySCHelper.subClusterEnergy(3),
0147 mySCHelper.subClusterEta(3),
0148 mySCHelper.subClusterPhi(3),
0149 mySCHelper.subClusterEmax(3),
0150 mySCHelper.subClusterE3x3(3),
0151 mySCHelper.nPreshowerClusters(),
0152 mySCHelper.eESClusters(),
0153 mySCHelper.esClusterEnergy(0),
0154 mySCHelper.esClusterEta(0),
0155 mySCHelper.esClusterPhi(0),
0156 mySCHelper.esClusterEnergy(1),
0157 mySCHelper.esClusterEta(1),
0158 mySCHelper.esClusterPhi(1),
0159 mySCHelper.esClusterEnergy(2),
0160 mySCHelper.esClusterEta(2),
0161 mySCHelper.esClusterPhi(2),
0162 ele->isEB(),
0163 printDebug);
0164 } else if (fVersionType == kNoTrkVarV1) {
0165 return regressionValueNoTrkVarV1(mySCHelper.rawEnergy(),
0166 mySCHelper.eta(),
0167 mySCHelper.phi(),
0168 mySCHelper.r9(),
0169 mySCHelper.etaWidth(),
0170 mySCHelper.phiWidth(),
0171 mySCHelper.clustersSize(),
0172 mySCHelper.hadronicOverEm(),
0173 rho,
0174 nvertices,
0175 mySCHelper.seedEta(),
0176 mySCHelper.seedPhi(),
0177 mySCHelper.seedEnergy(),
0178 mySCHelper.e3x3(),
0179 mySCHelper.e5x5(),
0180 mySCHelper.sigmaIetaIeta(),
0181 mySCHelper.spp(),
0182 mySCHelper.sep(),
0183 mySCHelper.eMax(),
0184 mySCHelper.e2nd(),
0185 mySCHelper.eTop(),
0186 mySCHelper.eBottom(),
0187 mySCHelper.eLeft(),
0188 mySCHelper.eRight(),
0189 mySCHelper.e2x5Max(),
0190 mySCHelper.e2x5Top(),
0191 mySCHelper.e2x5Bottom(),
0192 mySCHelper.e2x5Left(),
0193 mySCHelper.e2x5Right(),
0194 mySCHelper.ietaSeed(),
0195 mySCHelper.iphiSeed(),
0196 mySCHelper.etaCrySeed(),
0197 mySCHelper.phiCrySeed(),
0198 mySCHelper.preshowerEnergyOverRaw(),
0199 ele->ecalDrivenSeed(),
0200 printDebug);
0201 } else if (fVersionType == kWithTrkVarV1) {
0202 return regressionValueWithTrkVarV1(mySCHelper.rawEnergy(),
0203 mySCHelper.eta(),
0204 mySCHelper.phi(),
0205 mySCHelper.r9(),
0206 mySCHelper.etaWidth(),
0207 mySCHelper.phiWidth(),
0208 mySCHelper.clustersSize(),
0209 mySCHelper.hadronicOverEm(),
0210 rho,
0211 nvertices,
0212 mySCHelper.seedEta(),
0213 mySCHelper.seedPhi(),
0214 mySCHelper.seedEnergy(),
0215 mySCHelper.e3x3(),
0216 mySCHelper.e5x5(),
0217 mySCHelper.sigmaIetaIeta(),
0218 mySCHelper.spp(),
0219 mySCHelper.sep(),
0220 mySCHelper.eMax(),
0221 mySCHelper.e2nd(),
0222 mySCHelper.eTop(),
0223 mySCHelper.eBottom(),
0224 mySCHelper.eLeft(),
0225 mySCHelper.eRight(),
0226 mySCHelper.e2x5Max(),
0227 mySCHelper.e2x5Top(),
0228 mySCHelper.e2x5Bottom(),
0229 mySCHelper.e2x5Left(),
0230 mySCHelper.e2x5Right(),
0231 mySCHelper.ietaSeed(),
0232 mySCHelper.iphiSeed(),
0233 mySCHelper.etaCrySeed(),
0234 mySCHelper.phiCrySeed(),
0235 mySCHelper.preshowerEnergyOverRaw(),
0236 ele->ecalDrivenSeed(),
0237 ele->trackMomentumAtVtx().R(),
0238 fmax(ele->fbrem(), -1.0),
0239 ele->charge(),
0240 fmin(ele->eSuperClusterOverP(), 20.0),
0241 ele->trackMomentumError(),
0242 ele->correctedEcalEnergyError(),
0243 ele->classification(),
0244 printDebug);
0245 } else if (fVersionType == kWithTrkVarV2) {
0246 return regressionValueWithTrkVarV2(
0247 mySCHelper.rawEnergy(),
0248 mySCHelper.eta(),
0249 mySCHelper.phi(),
0250 mySCHelper.r9(),
0251 mySCHelper.etaWidth(),
0252 mySCHelper.phiWidth(),
0253 mySCHelper.clustersSize(),
0254 mySCHelper.hadronicOverEm(),
0255 rho,
0256 nvertices,
0257 mySCHelper.seedEta(),
0258 mySCHelper.seedPhi(),
0259 mySCHelper.seedEnergy(),
0260 mySCHelper.e3x3(),
0261 mySCHelper.e5x5(),
0262 mySCHelper.sigmaIetaIeta(),
0263 mySCHelper.spp(),
0264 mySCHelper.sep(),
0265 mySCHelper.eMax(),
0266 mySCHelper.e2nd(),
0267 mySCHelper.eTop(),
0268 mySCHelper.eBottom(),
0269 mySCHelper.eLeft(),
0270 mySCHelper.eRight(),
0271 mySCHelper.e2x5Max(),
0272 mySCHelper.e2x5Top(),
0273 mySCHelper.e2x5Bottom(),
0274 mySCHelper.e2x5Left(),
0275 mySCHelper.e2x5Right(),
0276 mySCHelper.ietaSeed(),
0277 mySCHelper.iphiSeed(),
0278 mySCHelper.etaCrySeed(),
0279 mySCHelper.phiCrySeed(),
0280 mySCHelper.preshowerEnergyOverRaw(),
0281 ele->ecalDrivenSeed(),
0282 ele->trackMomentumAtVtx().R(),
0283 fmax(ele->fbrem(), -1.0),
0284 ele->charge(),
0285 fmin(ele->eSuperClusterOverP(), 20.0),
0286 ele->trackMomentumError(),
0287 ele->correctedEcalEnergyError(),
0288 ele->classification(),
0289 fmin(std::abs(ele->deltaEtaSuperClusterTrackAtVtx()), 0.6),
0290 ele->deltaPhiSuperClusterTrackAtVtx(),
0291 ele->deltaEtaSeedClusterTrackAtCalo(),
0292 ele->deltaPhiSeedClusterTrackAtCalo(),
0293 ele->gsfTrack()->chi2() / ele->gsfTrack()->ndof(),
0294 (ele->closestCtfTrackRef().isNonnull() ? ele->closestCtfTrackRef()->hitPattern().trackerLayersWithMeasurement()
0295 : -1),
0296 fmin(ele->eEleClusterOverPout(), 20.0),
0297 printDebug);
0298 } else {
0299 std::cout << "Warning: Electron Regression Type " << fVersionType
0300 << " is not supported. Reverting to default electron momentum.\n";
0301 return ele->p();
0302 }
0303 }
0304
0305 double ElectronEnergyRegressionEvaluate::calculateRegressionEnergyUncertainty(
0306 const reco::GsfElectron *ele, SuperClusterHelper &mySCHelper, double rho, double nvertices, bool printDebug) {
0307 if (!fIsInitialized) {
0308 std::cout << "Error: Electron Energy Regression has not been initialized yet. return 0. \n";
0309 return 0;
0310 }
0311
0312 if (printDebug) {
0313 std::cout << "Regression Type: " << fVersionType << std::endl;
0314 std::cout << "Electron : " << ele->pt() << " " << ele->eta() << " " << ele->phi() << "\n";
0315 }
0316
0317 if (fVersionType == kNoTrkVar) {
0318 return regressionUncertaintyNoTrkVar(mySCHelper.rawEnergy(),
0319 mySCHelper.eta(),
0320 mySCHelper.phi(),
0321 mySCHelper.r9(),
0322 mySCHelper.etaWidth(),
0323 mySCHelper.phiWidth(),
0324 mySCHelper.clustersSize(),
0325 mySCHelper.hadronicOverEm(),
0326 rho,
0327 nvertices,
0328 mySCHelper.seedEta(),
0329 mySCHelper.seedPhi(),
0330 mySCHelper.seedEnergy(),
0331 mySCHelper.e3x3(),
0332 mySCHelper.e5x5(),
0333 mySCHelper.sigmaIetaIeta(),
0334 mySCHelper.spp(),
0335 mySCHelper.sep(),
0336 mySCHelper.eMax(),
0337 mySCHelper.e2nd(),
0338 mySCHelper.eTop(),
0339 mySCHelper.eBottom(),
0340 mySCHelper.eLeft(),
0341 mySCHelper.eRight(),
0342 mySCHelper.e2x5Max(),
0343 mySCHelper.e2x5Top(),
0344 mySCHelper.e2x5Bottom(),
0345 mySCHelper.e2x5Left(),
0346 mySCHelper.e2x5Right(),
0347 mySCHelper.ietaSeed(),
0348 mySCHelper.iphiSeed(),
0349 mySCHelper.etaCrySeed(),
0350 mySCHelper.phiCrySeed(),
0351 mySCHelper.preshowerEnergyOverRaw(),
0352 printDebug);
0353 }
0354 if (fVersionType == kWithSubCluVar) {
0355 return regressionUncertaintyWithSubClusters(mySCHelper.rawEnergy(),
0356 mySCHelper.eta(),
0357 mySCHelper.phi(),
0358 mySCHelper.r9(),
0359 mySCHelper.etaWidth(),
0360 mySCHelper.phiWidth(),
0361 mySCHelper.clustersSize(),
0362 mySCHelper.hadronicOverEm(),
0363 rho,
0364 nvertices,
0365 mySCHelper.seedEta(),
0366 mySCHelper.seedPhi(),
0367 mySCHelper.seedEnergy(),
0368 mySCHelper.e3x3(),
0369 mySCHelper.e5x5(),
0370 mySCHelper.sigmaIetaIeta(),
0371 mySCHelper.spp(),
0372 mySCHelper.sep(),
0373 mySCHelper.eMax(),
0374 mySCHelper.e2nd(),
0375 mySCHelper.eTop(),
0376 mySCHelper.eBottom(),
0377 mySCHelper.eLeft(),
0378 mySCHelper.eRight(),
0379 mySCHelper.e2x5Max(),
0380 mySCHelper.e2x5Top(),
0381 mySCHelper.e2x5Bottom(),
0382 mySCHelper.e2x5Left(),
0383 mySCHelper.e2x5Right(),
0384 mySCHelper.ietaSeed(),
0385 mySCHelper.iphiSeed(),
0386 mySCHelper.etaCrySeed(),
0387 mySCHelper.phiCrySeed(),
0388 mySCHelper.preshowerEnergyOverRaw(),
0389 ele->ecalDrivenSeed(),
0390 ele->isEBEtaGap(),
0391 ele->isEBPhiGap(),
0392 ele->isEEDeeGap(),
0393 mySCHelper.eSubClusters(),
0394 mySCHelper.subClusterEnergy(1),
0395 mySCHelper.subClusterEta(1),
0396 mySCHelper.subClusterPhi(1),
0397 mySCHelper.subClusterEmax(1),
0398 mySCHelper.subClusterE3x3(1),
0399 mySCHelper.subClusterEnergy(2),
0400 mySCHelper.subClusterEta(2),
0401 mySCHelper.subClusterPhi(2),
0402 mySCHelper.subClusterEmax(2),
0403 mySCHelper.subClusterE3x3(2),
0404 mySCHelper.subClusterEnergy(3),
0405 mySCHelper.subClusterEta(3),
0406 mySCHelper.subClusterPhi(3),
0407 mySCHelper.subClusterEmax(3),
0408 mySCHelper.subClusterE3x3(3),
0409 mySCHelper.nPreshowerClusters(),
0410 mySCHelper.eESClusters(),
0411 mySCHelper.esClusterEnergy(0),
0412 mySCHelper.esClusterEta(0),
0413 mySCHelper.esClusterPhi(0),
0414 mySCHelper.esClusterEnergy(1),
0415 mySCHelper.esClusterEta(1),
0416 mySCHelper.esClusterPhi(1),
0417 mySCHelper.esClusterEnergy(2),
0418 mySCHelper.esClusterEta(2),
0419 mySCHelper.esClusterPhi(2),
0420 ele->isEB(),
0421 printDebug);
0422 } else if (fVersionType == kNoTrkVarV1) {
0423 return regressionUncertaintyNoTrkVarV1(mySCHelper.rawEnergy(),
0424 mySCHelper.eta(),
0425 mySCHelper.phi(),
0426 mySCHelper.r9(),
0427 mySCHelper.etaWidth(),
0428 mySCHelper.phiWidth(),
0429 mySCHelper.clustersSize(),
0430 mySCHelper.hadronicOverEm(),
0431 rho,
0432 nvertices,
0433 mySCHelper.seedEta(),
0434 mySCHelper.seedPhi(),
0435 mySCHelper.seedEnergy(),
0436 mySCHelper.e3x3(),
0437 mySCHelper.e5x5(),
0438 mySCHelper.sigmaIetaIeta(),
0439 mySCHelper.spp(),
0440 mySCHelper.sep(),
0441 mySCHelper.eMax(),
0442 mySCHelper.e2nd(),
0443 mySCHelper.eTop(),
0444 mySCHelper.eBottom(),
0445 mySCHelper.eLeft(),
0446 mySCHelper.eRight(),
0447 mySCHelper.e2x5Max(),
0448 mySCHelper.e2x5Top(),
0449 mySCHelper.e2x5Bottom(),
0450 mySCHelper.e2x5Left(),
0451 mySCHelper.e2x5Right(),
0452 mySCHelper.ietaSeed(),
0453 mySCHelper.iphiSeed(),
0454 mySCHelper.etaCrySeed(),
0455 mySCHelper.phiCrySeed(),
0456 mySCHelper.preshowerEnergyOverRaw(),
0457 ele->ecalDrivenSeed(),
0458 printDebug);
0459 } else if (fVersionType == kWithTrkVarV1) {
0460 return regressionUncertaintyWithTrkVarV1(mySCHelper.rawEnergy(),
0461 mySCHelper.eta(),
0462 mySCHelper.phi(),
0463 mySCHelper.r9(),
0464 mySCHelper.etaWidth(),
0465 mySCHelper.phiWidth(),
0466 mySCHelper.clustersSize(),
0467 mySCHelper.hadronicOverEm(),
0468 rho,
0469 nvertices,
0470 mySCHelper.seedEta(),
0471 mySCHelper.seedPhi(),
0472 mySCHelper.seedEnergy(),
0473 mySCHelper.e3x3(),
0474 mySCHelper.e5x5(),
0475 mySCHelper.sigmaIetaIeta(),
0476 mySCHelper.spp(),
0477 mySCHelper.sep(),
0478 mySCHelper.eMax(),
0479 mySCHelper.e2nd(),
0480 mySCHelper.eTop(),
0481 mySCHelper.eBottom(),
0482 mySCHelper.eLeft(),
0483 mySCHelper.eRight(),
0484 mySCHelper.e2x5Max(),
0485 mySCHelper.e2x5Top(),
0486 mySCHelper.e2x5Bottom(),
0487 mySCHelper.e2x5Left(),
0488 mySCHelper.e2x5Right(),
0489 mySCHelper.ietaSeed(),
0490 mySCHelper.iphiSeed(),
0491 mySCHelper.etaCrySeed(),
0492 mySCHelper.phiCrySeed(),
0493 mySCHelper.preshowerEnergyOverRaw(),
0494 ele->ecalDrivenSeed(),
0495 ele->trackMomentumAtVtx().R(),
0496 fmax(ele->fbrem(), -1.0),
0497 ele->charge(),
0498 fmin(ele->eSuperClusterOverP(), 20.0),
0499 ele->trackMomentumError(),
0500 ele->correctedEcalEnergyError(),
0501 ele->classification(),
0502 printDebug);
0503 } else if (fVersionType == kWithTrkVarV2) {
0504 return regressionUncertaintyWithTrkVarV2(
0505 mySCHelper.rawEnergy(),
0506 mySCHelper.eta(),
0507 mySCHelper.phi(),
0508 mySCHelper.r9(),
0509 mySCHelper.etaWidth(),
0510 mySCHelper.phiWidth(),
0511 mySCHelper.clustersSize(),
0512 mySCHelper.hadronicOverEm(),
0513 rho,
0514 nvertices,
0515 mySCHelper.seedEta(),
0516 mySCHelper.seedPhi(),
0517 mySCHelper.seedEnergy(),
0518 mySCHelper.e3x3(),
0519 mySCHelper.e5x5(),
0520 mySCHelper.sigmaIetaIeta(),
0521 mySCHelper.spp(),
0522 mySCHelper.sep(),
0523 mySCHelper.eMax(),
0524 mySCHelper.e2nd(),
0525 mySCHelper.eTop(),
0526 mySCHelper.eBottom(),
0527 mySCHelper.eLeft(),
0528 mySCHelper.eRight(),
0529 mySCHelper.e2x5Max(),
0530 mySCHelper.e2x5Top(),
0531 mySCHelper.e2x5Bottom(),
0532 mySCHelper.e2x5Left(),
0533 mySCHelper.e2x5Right(),
0534 mySCHelper.ietaSeed(),
0535 mySCHelper.iphiSeed(),
0536 mySCHelper.etaCrySeed(),
0537 mySCHelper.phiCrySeed(),
0538 mySCHelper.preshowerEnergyOverRaw(),
0539 ele->ecalDrivenSeed(),
0540 ele->trackMomentumAtVtx().R(),
0541 fmax(ele->fbrem(), -1.0),
0542 ele->charge(),
0543 fmin(ele->eSuperClusterOverP(), 20.0),
0544 ele->trackMomentumError(),
0545 ele->correctedEcalEnergyError(),
0546 ele->classification(),
0547 fmin(std::abs(ele->deltaEtaSuperClusterTrackAtVtx()), 0.6),
0548 ele->deltaPhiSuperClusterTrackAtVtx(),
0549 ele->deltaEtaSeedClusterTrackAtCalo(),
0550 ele->deltaPhiSeedClusterTrackAtCalo(),
0551 ele->gsfTrack()->chi2() / ele->gsfTrack()->ndof(),
0552 (ele->closestCtfTrackRef().isNonnull() ? ele->closestCtfTrackRef()->hitPattern().trackerLayersWithMeasurement()
0553 : -1),
0554 fmin(ele->eEleClusterOverPout(), 20.0),
0555 printDebug);
0556 } else {
0557 std::cout << "Warning: Electron Regression Type " << fVersionType
0558 << " is not supported. Reverting to default electron momentum.\n";
0559 return ele->p();
0560 }
0561 }
0562 #endif
0563
0564 double ElectronEnergyRegressionEvaluate::regressionValueNoTrkVar(double SCRawEnergy,
0565 double scEta,
0566 double scPhi,
0567 double R9,
0568 double etawidth,
0569 double phiwidth,
0570 double NClusters,
0571 double HoE,
0572 double rho,
0573 double vertices,
0574 double EtaSeed,
0575 double PhiSeed,
0576 double ESeed,
0577 double E3x3Seed,
0578 double E5x5Seed,
0579 double see,
0580 double spp,
0581 double sep,
0582 double EMaxSeed,
0583 double E2ndSeed,
0584 double ETopSeed,
0585 double EBottomSeed,
0586 double ELeftSeed,
0587 double ERightSeed,
0588 double E2x5MaxSeed,
0589 double E2x5TopSeed,
0590 double E2x5BottomSeed,
0591 double E2x5LeftSeed,
0592 double E2x5RightSeed,
0593 double IEtaSeed,
0594 double IPhiSeed,
0595 double EtaCrySeed,
0596 double PhiCrySeed,
0597 double PreShowerOverRaw,
0598 bool printDebug) {
0599
0600 if (fIsInitialized == kFALSE) {
0601 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
0602 return 0;
0603 }
0604
0605
0606 if (!(fVersionType == kNoTrkVar)) {
0607 std::cout << "Error: Regression VersionType " << fVersionType
0608 << " is not supported to use function regressionValueNoTrkVar.\n";
0609 return 0;
0610 }
0611 assert(forestCorrection_ee);
0612
0613
0614 float *vals = (std::abs(scEta) <= 1.479) ? new float[38] : new float[31];
0615 if (std::abs(scEta) <= 1.479) {
0616 vals[0] = SCRawEnergy;
0617 vals[1] = scEta;
0618 vals[2] = scPhi;
0619 vals[3] = R9;
0620 vals[4] = E5x5Seed / SCRawEnergy;
0621 vals[5] = etawidth;
0622 vals[6] = phiwidth;
0623 vals[7] = NClusters;
0624 vals[8] = HoE;
0625 vals[9] = rho;
0626 vals[10] = vertices;
0627 vals[11] = EtaSeed - scEta;
0628 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
0629 vals[13] = ESeed / SCRawEnergy;
0630 vals[14] = E3x3Seed / ESeed;
0631 vals[15] = E5x5Seed / ESeed;
0632 vals[16] = see;
0633 vals[17] = spp;
0634 vals[18] = sep;
0635 vals[19] = EMaxSeed / ESeed;
0636 vals[20] = E2ndSeed / ESeed;
0637 vals[21] = ETopSeed / ESeed;
0638 vals[22] = EBottomSeed / ESeed;
0639 vals[23] = ELeftSeed / ESeed;
0640 vals[24] = ERightSeed / ESeed;
0641 vals[25] = E2x5MaxSeed / ESeed;
0642 vals[26] = E2x5TopSeed / ESeed;
0643 vals[27] = E2x5BottomSeed / ESeed;
0644 vals[28] = E2x5LeftSeed / ESeed;
0645 vals[29] = E2x5RightSeed / ESeed;
0646 vals[30] = IEtaSeed;
0647 vals[31] = IPhiSeed;
0648 vals[32] = ((int)IEtaSeed) % 5;
0649 vals[33] = ((int)IPhiSeed) % 2;
0650 vals[34] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
0651 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
0652 vals[35] = ((int)IPhiSeed) % 20;
0653 vals[36] = EtaCrySeed;
0654 vals[37] = PhiCrySeed;
0655 } else {
0656 vals[0] = SCRawEnergy;
0657 vals[1] = scEta;
0658 vals[2] = scPhi;
0659 vals[3] = R9;
0660 vals[4] = E5x5Seed / SCRawEnergy;
0661 vals[5] = etawidth;
0662 vals[6] = phiwidth;
0663 vals[7] = NClusters;
0664 vals[8] = HoE;
0665 vals[9] = rho;
0666 vals[10] = vertices;
0667 vals[11] = EtaSeed - scEta;
0668 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
0669 vals[13] = ESeed / SCRawEnergy;
0670 vals[14] = E3x3Seed / ESeed;
0671 vals[15] = E5x5Seed / ESeed;
0672 vals[16] = see;
0673 vals[17] = spp;
0674 vals[18] = sep;
0675 vals[19] = EMaxSeed / ESeed;
0676 vals[20] = E2ndSeed / ESeed;
0677 vals[21] = ETopSeed / ESeed;
0678 vals[22] = EBottomSeed / ESeed;
0679 vals[23] = ELeftSeed / ESeed;
0680 vals[24] = ERightSeed / ESeed;
0681 vals[25] = E2x5MaxSeed / ESeed;
0682 vals[26] = E2x5TopSeed / ESeed;
0683 vals[27] = E2x5BottomSeed / ESeed;
0684 vals[28] = E2x5LeftSeed / ESeed;
0685 vals[29] = E2x5RightSeed / ESeed;
0686 vals[30] = PreShowerOverRaw;
0687 }
0688
0689
0690 double regressionResult = 0;
0691 Int_t BinIndex = -1;
0692
0693 if (fVersionType == kNoTrkVar) {
0694 if (std::abs(scEta) <= 1.479) {
0695 regressionResult = SCRawEnergy * forestCorrection_eb->GetResponse(vals);
0696 BinIndex = 0;
0697 } else {
0698 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestCorrection_ee->GetResponse(vals);
0699 BinIndex = 1;
0700 }
0701 }
0702
0703
0704 if (printDebug) {
0705 if (std::abs(scEta) <= 1.479) {
0706 std::cout << "Barrel :";
0707 for (unsigned int v = 0; v < 38; ++v)
0708 std::cout << vals[v] << ", ";
0709 std::cout << "\n";
0710 } else {
0711 std::cout << "Endcap :";
0712 for (unsigned int v = 0; v < 31; ++v)
0713 std::cout << vals[v] << ", ";
0714 std::cout << "\n";
0715 }
0716 std::cout << "BinIndex : " << BinIndex << "\n";
0717 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
0718 std::cout << "regression energy = " << regressionResult << std::endl;
0719 }
0720
0721
0722 delete[] vals;
0723 return regressionResult;
0724 }
0725
0726 double ElectronEnergyRegressionEvaluate::regressionUncertaintyNoTrkVar(double SCRawEnergy,
0727 double scEta,
0728 double scPhi,
0729 double R9,
0730 double etawidth,
0731 double phiwidth,
0732 double NClusters,
0733 double HoE,
0734 double rho,
0735 double vertices,
0736 double EtaSeed,
0737 double PhiSeed,
0738 double ESeed,
0739 double E3x3Seed,
0740 double E5x5Seed,
0741 double see,
0742 double spp,
0743 double sep,
0744 double EMaxSeed,
0745 double E2ndSeed,
0746 double ETopSeed,
0747 double EBottomSeed,
0748 double ELeftSeed,
0749 double ERightSeed,
0750 double E2x5MaxSeed,
0751 double E2x5TopSeed,
0752 double E2x5BottomSeed,
0753 double E2x5LeftSeed,
0754 double E2x5RightSeed,
0755 double IEtaSeed,
0756 double IPhiSeed,
0757 double EtaCrySeed,
0758 double PhiCrySeed,
0759 double PreShowerOverRaw,
0760 bool printDebug) {
0761
0762 if (fIsInitialized == kFALSE) {
0763 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
0764 return 0;
0765 }
0766
0767
0768 if (!(fVersionType == kNoTrkVar)) {
0769 std::cout << "Error: Regression VersionType " << fVersionType
0770 << " is not supported to use function regressionValueNoTrkVar.\n";
0771 return 0;
0772 }
0773
0774
0775 float *vals = (std::abs(scEta) <= 1.479) ? new float[38] : new float[31];
0776 if (std::abs(scEta) <= 1.479) {
0777 vals[0] = SCRawEnergy;
0778 vals[1] = scEta;
0779 vals[2] = scPhi;
0780 vals[3] = R9;
0781 vals[4] = E5x5Seed / SCRawEnergy;
0782 vals[5] = etawidth;
0783 vals[6] = phiwidth;
0784 vals[7] = NClusters;
0785 vals[8] = HoE;
0786 vals[9] = rho;
0787 vals[10] = vertices;
0788 vals[11] = EtaSeed - scEta;
0789 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
0790 vals[13] = ESeed / SCRawEnergy;
0791 vals[14] = E3x3Seed / ESeed;
0792 vals[15] = E5x5Seed / ESeed;
0793 vals[16] = see;
0794 vals[17] = spp;
0795 vals[18] = sep;
0796 vals[19] = EMaxSeed / ESeed;
0797 vals[20] = E2ndSeed / ESeed;
0798 vals[21] = ETopSeed / ESeed;
0799 vals[22] = EBottomSeed / ESeed;
0800 vals[23] = ELeftSeed / ESeed;
0801 vals[24] = ERightSeed / ESeed;
0802 vals[25] = E2x5MaxSeed / ESeed;
0803 vals[26] = E2x5TopSeed / ESeed;
0804 vals[27] = E2x5BottomSeed / ESeed;
0805 vals[28] = E2x5LeftSeed / ESeed;
0806 vals[29] = E2x5RightSeed / ESeed;
0807 vals[30] = IEtaSeed;
0808 vals[31] = IPhiSeed;
0809 vals[32] = ((int)IEtaSeed) % 5;
0810 vals[33] = ((int)IPhiSeed) % 2;
0811 vals[34] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
0812 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
0813 vals[35] = ((int)IPhiSeed) % 20;
0814 vals[36] = EtaCrySeed;
0815 vals[37] = PhiCrySeed;
0816 } else {
0817 vals[0] = SCRawEnergy;
0818 vals[1] = scEta;
0819 vals[2] = scPhi;
0820 vals[3] = R9;
0821 vals[4] = E5x5Seed / SCRawEnergy;
0822 vals[5] = etawidth;
0823 vals[6] = phiwidth;
0824 vals[7] = NClusters;
0825 vals[8] = HoE;
0826 vals[9] = rho;
0827 vals[10] = vertices;
0828 vals[11] = EtaSeed - scEta;
0829 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
0830 vals[13] = ESeed / SCRawEnergy;
0831 vals[14] = E3x3Seed / ESeed;
0832 vals[15] = E5x5Seed / ESeed;
0833 vals[16] = see;
0834 vals[17] = spp;
0835 vals[18] = sep;
0836 vals[19] = EMaxSeed / ESeed;
0837 vals[20] = E2ndSeed / ESeed;
0838 vals[21] = ETopSeed / ESeed;
0839 vals[22] = EBottomSeed / ESeed;
0840 vals[23] = ELeftSeed / ESeed;
0841 vals[24] = ERightSeed / ESeed;
0842 vals[25] = E2x5MaxSeed / ESeed;
0843 vals[26] = E2x5TopSeed / ESeed;
0844 vals[27] = E2x5BottomSeed / ESeed;
0845 vals[28] = E2x5LeftSeed / ESeed;
0846 vals[29] = E2x5RightSeed / ESeed;
0847 vals[30] = PreShowerOverRaw;
0848 }
0849
0850
0851 double regressionResult = 0;
0852 Int_t BinIndex = -1;
0853
0854 if (fVersionType == kNoTrkVar) {
0855 if (std::abs(scEta) <= 1.479) {
0856 regressionResult = SCRawEnergy * forestUncertainty_eb->GetResponse(vals);
0857 BinIndex = 0;
0858 } else {
0859 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestUncertainty_ee->GetResponse(vals);
0860 BinIndex = 1;
0861 }
0862 }
0863
0864
0865 if (printDebug) {
0866 if (std::abs(scEta) <= 1.479) {
0867 std::cout << "Barrel :";
0868 for (unsigned int v = 0; v < 38; ++v)
0869 std::cout << vals[v] << ", ";
0870 std::cout << "\n";
0871 } else {
0872 std::cout << "Endcap :";
0873 for (unsigned int v = 0; v < 31; ++v)
0874 std::cout << vals[v] << ", ";
0875 std::cout << "\n";
0876 }
0877 std::cout << "BinIndex : " << BinIndex << "\n";
0878 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
0879 std::cout << "regression energy uncertainty = " << regressionResult << std::endl;
0880 }
0881
0882
0883 delete[] vals;
0884 return regressionResult;
0885 }
0886
0887 double ElectronEnergyRegressionEvaluate::regressionValueNoTrkVarV1(double SCRawEnergy,
0888 double scEta,
0889 double scPhi,
0890 double R9,
0891 double etawidth,
0892 double phiwidth,
0893 double NClusters,
0894 double HoE,
0895 double rho,
0896 double vertices,
0897 double EtaSeed,
0898 double PhiSeed,
0899 double ESeed,
0900 double E3x3Seed,
0901 double E5x5Seed,
0902 double see,
0903 double spp,
0904 double sep,
0905 double EMaxSeed,
0906 double E2ndSeed,
0907 double ETopSeed,
0908 double EBottomSeed,
0909 double ELeftSeed,
0910 double ERightSeed,
0911 double E2x5MaxSeed,
0912 double E2x5TopSeed,
0913 double E2x5BottomSeed,
0914 double E2x5LeftSeed,
0915 double E2x5RightSeed,
0916 double IEtaSeed,
0917 double IPhiSeed,
0918 double EtaCrySeed,
0919 double PhiCrySeed,
0920 double PreShowerOverRaw,
0921 int IsEcalDriven,
0922 bool printDebug) {
0923
0924 if (fIsInitialized == kFALSE) {
0925 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
0926 return 0;
0927 }
0928
0929
0930 if (!(fVersionType == kNoTrkVarV1)) {
0931 std::cout << "Error: Regression VersionType " << fVersionType
0932 << " is not supported to use function regressionValueNoTrkVar.\n";
0933 return 0;
0934 }
0935
0936
0937 float *vals = (std::abs(scEta) <= 1.479) ? new float[39] : new float[32];
0938 if (std::abs(scEta) <= 1.479) {
0939 vals[0] = SCRawEnergy;
0940 vals[1] = scEta;
0941 vals[2] = scPhi;
0942 vals[3] = R9;
0943 vals[4] = E5x5Seed / SCRawEnergy;
0944 vals[5] = etawidth;
0945 vals[6] = phiwidth;
0946 vals[7] = NClusters;
0947 vals[8] = HoE;
0948 vals[9] = rho;
0949 vals[10] = vertices;
0950 vals[11] = EtaSeed - scEta;
0951 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
0952 vals[13] = ESeed / SCRawEnergy;
0953 vals[14] = E3x3Seed / ESeed;
0954 vals[15] = E5x5Seed / ESeed;
0955 vals[16] = see;
0956 vals[17] = spp;
0957 vals[18] = sep;
0958 vals[19] = EMaxSeed / ESeed;
0959 vals[20] = E2ndSeed / ESeed;
0960 vals[21] = ETopSeed / ESeed;
0961 vals[22] = EBottomSeed / ESeed;
0962 vals[23] = ELeftSeed / ESeed;
0963 vals[24] = ERightSeed / ESeed;
0964 vals[25] = E2x5MaxSeed / ESeed;
0965 vals[26] = E2x5TopSeed / ESeed;
0966 vals[27] = E2x5BottomSeed / ESeed;
0967 vals[28] = E2x5LeftSeed / ESeed;
0968 vals[29] = E2x5RightSeed / ESeed;
0969 vals[30] = IsEcalDriven;
0970 vals[31] = IEtaSeed;
0971 vals[32] = IPhiSeed;
0972 vals[33] = ((int)IEtaSeed) % 5;
0973 vals[34] = ((int)IPhiSeed) % 2;
0974 vals[35] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
0975 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
0976 vals[36] = ((int)IPhiSeed) % 20;
0977 vals[37] = EtaCrySeed;
0978 vals[38] = PhiCrySeed;
0979 } else {
0980 vals[0] = SCRawEnergy;
0981 vals[1] = scEta;
0982 vals[2] = scPhi;
0983 vals[3] = R9;
0984 vals[4] = E5x5Seed / SCRawEnergy;
0985 vals[5] = etawidth;
0986 vals[6] = phiwidth;
0987 vals[7] = NClusters;
0988 vals[8] = HoE;
0989 vals[9] = rho;
0990 vals[10] = vertices;
0991 vals[11] = EtaSeed - scEta;
0992 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
0993 vals[13] = ESeed / SCRawEnergy;
0994 vals[14] = E3x3Seed / ESeed;
0995 vals[15] = E5x5Seed / ESeed;
0996 vals[16] = see;
0997 vals[17] = spp;
0998 vals[18] = sep;
0999 vals[19] = EMaxSeed / ESeed;
1000 vals[20] = E2ndSeed / ESeed;
1001 vals[21] = ETopSeed / ESeed;
1002 vals[22] = EBottomSeed / ESeed;
1003 vals[23] = ELeftSeed / ESeed;
1004 vals[24] = ERightSeed / ESeed;
1005 vals[25] = E2x5MaxSeed / ESeed;
1006 vals[26] = E2x5TopSeed / ESeed;
1007 vals[27] = E2x5BottomSeed / ESeed;
1008 vals[28] = E2x5LeftSeed / ESeed;
1009 vals[29] = E2x5RightSeed / ESeed;
1010 vals[30] = IsEcalDriven;
1011 vals[31] = PreShowerOverRaw;
1012 }
1013
1014
1015 double regressionResult = 0;
1016 Int_t BinIndex = -1;
1017
1018 if (fVersionType == kNoTrkVarV1) {
1019 if (std::abs(scEta) <= 1.479) {
1020 regressionResult = SCRawEnergy * forestCorrection_eb->GetResponse(vals);
1021 BinIndex = 0;
1022 } else {
1023 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestCorrection_ee->GetResponse(vals);
1024 BinIndex = 1;
1025 }
1026 }
1027
1028
1029 if (printDebug) {
1030 if (std::abs(scEta) <= 1.479) {
1031 std::cout << "Barrel :";
1032 for (unsigned int v = 0; v < 39; ++v)
1033 std::cout << vals[v] << ", ";
1034 std::cout << "\n";
1035 } else {
1036 std::cout << "Endcap :";
1037 for (unsigned int v = 0; v < 32; ++v)
1038 std::cout << vals[v] << ", ";
1039 std::cout << "\n";
1040 }
1041 std::cout << "BinIndex : " << BinIndex << "\n";
1042 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
1043 std::cout << "regression energy = " << regressionResult << std::endl;
1044 }
1045
1046
1047 delete[] vals;
1048 return regressionResult;
1049 }
1050
1051 double ElectronEnergyRegressionEvaluate::regressionUncertaintyNoTrkVarV1(double SCRawEnergy,
1052 double scEta,
1053 double scPhi,
1054 double R9,
1055 double etawidth,
1056 double phiwidth,
1057 double NClusters,
1058 double HoE,
1059 double rho,
1060 double vertices,
1061 double EtaSeed,
1062 double PhiSeed,
1063 double ESeed,
1064 double E3x3Seed,
1065 double E5x5Seed,
1066 double see,
1067 double spp,
1068 double sep,
1069 double EMaxSeed,
1070 double E2ndSeed,
1071 double ETopSeed,
1072 double EBottomSeed,
1073 double ELeftSeed,
1074 double ERightSeed,
1075 double E2x5MaxSeed,
1076 double E2x5TopSeed,
1077 double E2x5BottomSeed,
1078 double E2x5LeftSeed,
1079 double E2x5RightSeed,
1080 double IEtaSeed,
1081 double IPhiSeed,
1082 double EtaCrySeed,
1083 double PhiCrySeed,
1084 double PreShowerOverRaw,
1085 int IsEcalDriven,
1086 bool printDebug) {
1087
1088 if (fIsInitialized == kFALSE) {
1089 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
1090 return 0;
1091 }
1092
1093
1094 if (!(fVersionType == kNoTrkVarV1)) {
1095 std::cout << "Error: Regression VersionType " << fVersionType
1096 << " is not supported to use function regressionValueNoTrkVar.\n";
1097 return 0;
1098 }
1099
1100
1101 float *vals = (std::abs(scEta) <= 1.479) ? new float[39] : new float[32];
1102 if (std::abs(scEta) <= 1.479) {
1103 vals[0] = SCRawEnergy;
1104 vals[1] = scEta;
1105 vals[2] = scPhi;
1106 vals[3] = R9;
1107 vals[4] = E5x5Seed / SCRawEnergy;
1108 vals[5] = etawidth;
1109 vals[6] = phiwidth;
1110 vals[7] = NClusters;
1111 vals[8] = HoE;
1112 vals[9] = rho;
1113 vals[10] = vertices;
1114 vals[11] = EtaSeed - scEta;
1115 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1116 vals[13] = ESeed / SCRawEnergy;
1117 vals[14] = E3x3Seed / ESeed;
1118 vals[15] = E5x5Seed / ESeed;
1119 vals[16] = see;
1120 vals[17] = spp;
1121 vals[18] = sep;
1122 vals[19] = EMaxSeed / ESeed;
1123 vals[20] = E2ndSeed / ESeed;
1124 vals[21] = ETopSeed / ESeed;
1125 vals[22] = EBottomSeed / ESeed;
1126 vals[23] = ELeftSeed / ESeed;
1127 vals[24] = ERightSeed / ESeed;
1128 vals[25] = E2x5MaxSeed / ESeed;
1129 vals[26] = E2x5TopSeed / ESeed;
1130 vals[27] = E2x5BottomSeed / ESeed;
1131 vals[28] = E2x5LeftSeed / ESeed;
1132 vals[29] = E2x5RightSeed / ESeed;
1133 vals[30] = IsEcalDriven;
1134 vals[31] = IEtaSeed;
1135 vals[32] = IPhiSeed;
1136 vals[33] = ((int)IEtaSeed) % 5;
1137 vals[34] = ((int)IPhiSeed) % 2;
1138 vals[35] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
1139 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
1140 vals[36] = ((int)IPhiSeed) % 20;
1141 vals[37] = EtaCrySeed;
1142 vals[38] = PhiCrySeed;
1143 } else {
1144 vals[0] = SCRawEnergy;
1145 vals[1] = scEta;
1146 vals[2] = scPhi;
1147 vals[3] = R9;
1148 vals[4] = E5x5Seed / SCRawEnergy;
1149 vals[5] = etawidth;
1150 vals[6] = phiwidth;
1151 vals[7] = NClusters;
1152 vals[8] = HoE;
1153 vals[9] = rho;
1154 vals[10] = vertices;
1155 vals[11] = EtaSeed - scEta;
1156 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1157 vals[13] = ESeed / SCRawEnergy;
1158 vals[14] = E3x3Seed / ESeed;
1159 vals[15] = E5x5Seed / ESeed;
1160 vals[16] = see;
1161 vals[17] = spp;
1162 vals[18] = sep;
1163 vals[19] = EMaxSeed / ESeed;
1164 vals[20] = E2ndSeed / ESeed;
1165 vals[21] = ETopSeed / ESeed;
1166 vals[22] = EBottomSeed / ESeed;
1167 vals[23] = ELeftSeed / ESeed;
1168 vals[24] = ERightSeed / ESeed;
1169 vals[25] = E2x5MaxSeed / ESeed;
1170 vals[26] = E2x5TopSeed / ESeed;
1171 vals[27] = E2x5BottomSeed / ESeed;
1172 vals[28] = E2x5LeftSeed / ESeed;
1173 vals[29] = E2x5RightSeed / ESeed;
1174 vals[30] = IsEcalDriven;
1175 vals[31] = PreShowerOverRaw;
1176 }
1177
1178
1179 double regressionResult = 0;
1180 Int_t BinIndex = -1;
1181
1182 if (fVersionType == kNoTrkVarV1) {
1183 if (std::abs(scEta) <= 1.479) {
1184 regressionResult = SCRawEnergy * forestUncertainty_eb->GetResponse(vals);
1185 BinIndex = 0;
1186 } else {
1187 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestUncertainty_ee->GetResponse(vals);
1188 BinIndex = 1;
1189 }
1190 }
1191
1192
1193 if (printDebug) {
1194 if (std::abs(scEta) <= 1.479) {
1195 std::cout << "Barrel :";
1196 for (unsigned int v = 0; v < 39; ++v)
1197 std::cout << vals[v] << ", ";
1198 std::cout << "\n";
1199 } else {
1200 std::cout << "Endcap :";
1201 for (unsigned int v = 0; v < 32; ++v)
1202 std::cout << vals[v] << ", ";
1203 std::cout << "\n";
1204 }
1205 std::cout << "BinIndex : " << BinIndex << "\n";
1206 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
1207 std::cout << "regression energy uncertainty = " << regressionResult << std::endl;
1208 }
1209
1210
1211 delete[] vals;
1212 return regressionResult;
1213 }
1214
1215
1216
1217 double ElectronEnergyRegressionEvaluate::regressionValueWithTrkVar(double electronP,
1218 double SCRawEnergy,
1219 double scEta,
1220 double scPhi,
1221 double R9,
1222 double etawidth,
1223 double phiwidth,
1224 double NClusters,
1225 double HoE,
1226 double rho,
1227 double vertices,
1228 double EtaSeed,
1229 double PhiSeed,
1230 double ESeed,
1231 double E3x3Seed,
1232 double E5x5Seed,
1233 double see,
1234 double spp,
1235 double sep,
1236 double EMaxSeed,
1237 double E2ndSeed,
1238 double ETopSeed,
1239 double EBottomSeed,
1240 double ELeftSeed,
1241 double ERightSeed,
1242 double E2x5MaxSeed,
1243 double E2x5TopSeed,
1244 double E2x5BottomSeed,
1245 double E2x5LeftSeed,
1246 double E2x5RightSeed,
1247 double pt,
1248 double GsfTrackPIn,
1249 double fbrem,
1250 double Charge,
1251 double EoP,
1252 double IEtaSeed,
1253 double IPhiSeed,
1254 double EtaCrySeed,
1255 double PhiCrySeed,
1256 double PreShowerOverRaw,
1257 bool printDebug) {
1258
1259 if (fIsInitialized == kFALSE) {
1260 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
1261 return 0;
1262 }
1263
1264
1265 assert(fVersionType == kWithTrkVar);
1266
1267 float *vals = (std::abs(scEta) <= 1.479) ? new float[43] : new float[36];
1268 if (std::abs(scEta) <= 1.479) {
1269 vals[0] = SCRawEnergy;
1270 vals[1] = scEta;
1271 vals[2] = scPhi;
1272 vals[3] = R9;
1273 vals[4] = E5x5Seed / SCRawEnergy;
1274 vals[5] = etawidth;
1275 vals[6] = phiwidth;
1276 vals[7] = NClusters;
1277 vals[8] = HoE;
1278 vals[9] = rho;
1279 vals[10] = vertices;
1280 vals[11] = EtaSeed - scEta;
1281 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1282 vals[13] = ESeed / SCRawEnergy;
1283 vals[14] = E3x3Seed / ESeed;
1284 vals[15] = E5x5Seed / ESeed;
1285 vals[16] = see;
1286 vals[17] = spp;
1287 vals[18] = sep;
1288 vals[19] = EMaxSeed / ESeed;
1289 vals[20] = E2ndSeed / ESeed;
1290 vals[21] = ETopSeed / ESeed;
1291 vals[22] = EBottomSeed / ESeed;
1292 vals[23] = ELeftSeed / ESeed;
1293 vals[24] = ERightSeed / ESeed;
1294 vals[25] = E2x5MaxSeed / ESeed;
1295 vals[26] = E2x5TopSeed / ESeed;
1296 vals[27] = E2x5BottomSeed / ESeed;
1297 vals[28] = E2x5LeftSeed / ESeed;
1298 vals[29] = E2x5RightSeed / ESeed;
1299 vals[30] = pt;
1300 vals[31] = GsfTrackPIn;
1301 vals[32] = fbrem;
1302 vals[33] = Charge;
1303 vals[34] = EoP;
1304 vals[35] = IEtaSeed;
1305 vals[36] = IPhiSeed;
1306 vals[37] = ((int)IEtaSeed) % 5;
1307 vals[38] = ((int)IPhiSeed) % 2;
1308 vals[39] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
1309 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
1310 vals[40] = ((int)IPhiSeed) % 20;
1311 vals[41] = EtaCrySeed;
1312 vals[42] = PhiCrySeed;
1313 }
1314
1315 else {
1316 vals[0] = SCRawEnergy;
1317 vals[1] = scEta;
1318 vals[2] = scPhi;
1319 vals[3] = R9;
1320 vals[4] = E5x5Seed / SCRawEnergy;
1321 vals[5] = etawidth;
1322 vals[6] = phiwidth;
1323 vals[7] = NClusters;
1324 vals[8] = HoE;
1325 vals[9] = rho;
1326 vals[10] = vertices;
1327 vals[11] = EtaSeed - scEta;
1328 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1329 vals[13] = ESeed / SCRawEnergy;
1330 vals[14] = E3x3Seed / ESeed;
1331 vals[15] = E5x5Seed / ESeed;
1332 vals[16] = see;
1333 vals[17] = spp;
1334 vals[18] = sep;
1335 vals[19] = EMaxSeed / ESeed;
1336 vals[20] = E2ndSeed / ESeed;
1337 vals[21] = ETopSeed / ESeed;
1338 vals[22] = EBottomSeed / ESeed;
1339 vals[23] = ELeftSeed / ESeed;
1340 vals[24] = ERightSeed / ESeed;
1341 vals[25] = E2x5MaxSeed / ESeed;
1342 vals[26] = E2x5TopSeed / ESeed;
1343 vals[27] = E2x5BottomSeed / ESeed;
1344 vals[28] = E2x5LeftSeed / ESeed;
1345 vals[29] = E2x5RightSeed / ESeed;
1346 vals[30] = pt;
1347 vals[31] = GsfTrackPIn;
1348 vals[32] = fbrem;
1349 vals[33] = Charge;
1350 vals[34] = EoP;
1351 vals[35] = PreShowerOverRaw;
1352 }
1353
1354
1355 double regressionResult = 0;
1356
1357 if (fVersionType == kWithTrkVar) {
1358 if (std::abs(scEta) <= 1.479)
1359 regressionResult = SCRawEnergy * forestCorrection_eb->GetResponse(vals);
1360 else
1361 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestCorrection_ee->GetResponse(vals);
1362 }
1363
1364
1365 if (printDebug) {
1366 if (scEta <= 1.479) {
1367 std::cout << "Barrel :";
1368 for (unsigned int v = 0; v < 43; ++v)
1369 std::cout << vals[v] << ", ";
1370 std::cout << "\n";
1371 } else {
1372 std::cout << "Endcap :";
1373 for (unsigned int v = 0; v < 36; ++v)
1374 std::cout << vals[v] << ", ";
1375 std::cout << "\n";
1376 }
1377 std::cout << "pt = " << pt << " : SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw
1378 << std::endl;
1379 std::cout << "regression energy = " << regressionResult << std::endl;
1380 }
1381
1382
1383 delete[] vals;
1384 return regressionResult;
1385 }
1386
1387
1388
1389 double ElectronEnergyRegressionEvaluate::regressionUncertaintyWithTrkVar(double electronP,
1390 double SCRawEnergy,
1391 double scEta,
1392 double scPhi,
1393 double R9,
1394 double etawidth,
1395 double phiwidth,
1396 double NClusters,
1397 double HoE,
1398 double rho,
1399 double vertices,
1400 double EtaSeed,
1401 double PhiSeed,
1402 double ESeed,
1403 double E3x3Seed,
1404 double E5x5Seed,
1405 double see,
1406 double spp,
1407 double sep,
1408 double EMaxSeed,
1409 double E2ndSeed,
1410 double ETopSeed,
1411 double EBottomSeed,
1412 double ELeftSeed,
1413 double ERightSeed,
1414 double E2x5MaxSeed,
1415 double E2x5TopSeed,
1416 double E2x5BottomSeed,
1417 double E2x5LeftSeed,
1418 double E2x5RightSeed,
1419 double pt,
1420 double GsfTrackPIn,
1421 double fbrem,
1422 double Charge,
1423 double EoP,
1424 double IEtaSeed,
1425 double IPhiSeed,
1426 double EtaCrySeed,
1427 double PhiCrySeed,
1428 double PreShowerOverRaw,
1429 bool printDebug) {
1430
1431 if (fIsInitialized == kFALSE) {
1432 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
1433 return 0;
1434 }
1435
1436
1437 assert(fVersionType == kWithTrkVar);
1438
1439 float *vals = (std::abs(scEta) <= 1.479) ? new float[43] : new float[36];
1440 if (std::abs(scEta) <= 1.479) {
1441 vals[0] = SCRawEnergy;
1442 vals[1] = scEta;
1443 vals[2] = scPhi;
1444 vals[3] = R9;
1445 vals[4] = E5x5Seed / SCRawEnergy;
1446 vals[5] = etawidth;
1447 vals[6] = phiwidth;
1448 vals[7] = NClusters;
1449 vals[8] = HoE;
1450 vals[9] = rho;
1451 vals[10] = vertices;
1452 vals[11] = EtaSeed - scEta;
1453 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1454 vals[13] = ESeed / SCRawEnergy;
1455 vals[14] = E3x3Seed / ESeed;
1456 vals[15] = E5x5Seed / ESeed;
1457 vals[16] = see;
1458 vals[17] = spp;
1459 vals[18] = sep;
1460 vals[19] = EMaxSeed / ESeed;
1461 vals[20] = E2ndSeed / ESeed;
1462 vals[21] = ETopSeed / ESeed;
1463 vals[22] = EBottomSeed / ESeed;
1464 vals[23] = ELeftSeed / ESeed;
1465 vals[24] = ERightSeed / ESeed;
1466 vals[25] = E2x5MaxSeed / ESeed;
1467 vals[26] = E2x5TopSeed / ESeed;
1468 vals[27] = E2x5BottomSeed / ESeed;
1469 vals[28] = E2x5LeftSeed / ESeed;
1470 vals[29] = E2x5RightSeed / ESeed;
1471 vals[30] = pt;
1472 vals[31] = GsfTrackPIn;
1473 vals[32] = fbrem;
1474 vals[33] = Charge;
1475 vals[34] = EoP;
1476 vals[35] = IEtaSeed;
1477 vals[36] = IPhiSeed;
1478 vals[37] = ((int)IEtaSeed) % 5;
1479 vals[38] = ((int)IPhiSeed) % 2;
1480 vals[39] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
1481 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
1482 vals[40] = ((int)IPhiSeed) % 20;
1483 vals[41] = EtaCrySeed;
1484 vals[42] = PhiCrySeed;
1485 }
1486
1487 else {
1488 vals[0] = SCRawEnergy;
1489 vals[1] = scEta;
1490 vals[2] = scPhi;
1491 vals[3] = R9;
1492 vals[4] = E5x5Seed / SCRawEnergy;
1493 vals[5] = etawidth;
1494 vals[6] = phiwidth;
1495 vals[7] = NClusters;
1496 vals[8] = HoE;
1497 vals[9] = rho;
1498 vals[10] = vertices;
1499 vals[11] = EtaSeed - scEta;
1500 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1501 vals[13] = ESeed / SCRawEnergy;
1502 vals[14] = E3x3Seed / ESeed;
1503 vals[15] = E5x5Seed / ESeed;
1504 vals[16] = see;
1505 vals[17] = spp;
1506 vals[18] = sep;
1507 vals[19] = EMaxSeed / ESeed;
1508 vals[20] = E2ndSeed / ESeed;
1509 vals[21] = ETopSeed / ESeed;
1510 vals[22] = EBottomSeed / ESeed;
1511 vals[23] = ELeftSeed / ESeed;
1512 vals[24] = ERightSeed / ESeed;
1513 vals[25] = E2x5MaxSeed / ESeed;
1514 vals[26] = E2x5TopSeed / ESeed;
1515 vals[27] = E2x5BottomSeed / ESeed;
1516 vals[28] = E2x5LeftSeed / ESeed;
1517 vals[29] = E2x5RightSeed / ESeed;
1518 vals[30] = pt;
1519 vals[31] = GsfTrackPIn;
1520 vals[32] = fbrem;
1521 vals[33] = Charge;
1522 vals[34] = EoP;
1523 vals[35] = PreShowerOverRaw;
1524 }
1525
1526
1527 double regressionResult = 0;
1528
1529 if (fVersionType == kWithTrkVar) {
1530 if (std::abs(scEta) <= 1.479)
1531 regressionResult = SCRawEnergy * forestUncertainty_eb->GetResponse(vals);
1532 else
1533 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestUncertainty_ee->GetResponse(vals);
1534 }
1535
1536
1537 if (printDebug) {
1538 if (scEta <= 1.479) {
1539 std::cout << "Barrel :";
1540 for (unsigned int v = 0; v < 43; ++v)
1541 std::cout << vals[v] << ", ";
1542 std::cout << "\n";
1543 } else {
1544 std::cout << "Endcap :";
1545 for (unsigned int v = 0; v < 36; ++v)
1546 std::cout << vals[v] << ", ";
1547 std::cout << "\n";
1548 }
1549 std::cout << "pt = " << pt << " : SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw
1550 << std::endl;
1551 std::cout << "regression energy uncertainty = " << regressionResult << std::endl;
1552 }
1553
1554
1555 delete[] vals;
1556 return regressionResult;
1557 }
1558
1559 double ElectronEnergyRegressionEvaluate::regressionValueWithTrkVarV1(double SCRawEnergy,
1560 double scEta,
1561 double scPhi,
1562 double R9,
1563 double etawidth,
1564 double phiwidth,
1565 double NClusters,
1566 double HoE,
1567 double rho,
1568 double vertices,
1569 double EtaSeed,
1570 double PhiSeed,
1571 double ESeed,
1572 double E3x3Seed,
1573 double E5x5Seed,
1574 double see,
1575 double spp,
1576 double sep,
1577 double EMaxSeed,
1578 double E2ndSeed,
1579 double ETopSeed,
1580 double EBottomSeed,
1581 double ELeftSeed,
1582 double ERightSeed,
1583 double E2x5MaxSeed,
1584 double E2x5TopSeed,
1585 double E2x5BottomSeed,
1586 double E2x5LeftSeed,
1587 double E2x5RightSeed,
1588 double IEtaSeed,
1589 double IPhiSeed,
1590 double EtaCrySeed,
1591 double PhiCrySeed,
1592 double PreShowerOverRaw,
1593 int IsEcalDriven,
1594 double GsfTrackPIn,
1595 double fbrem,
1596 double Charge,
1597 double EoP,
1598 double TrackMomentumError,
1599 double EcalEnergyError,
1600 int Classification,
1601 bool printDebug) {
1602
1603 if (fIsInitialized == kFALSE) {
1604 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
1605 return 0;
1606 }
1607
1608
1609 assert(fVersionType == kWithTrkVarV1);
1610
1611 float *vals = (std::abs(scEta) <= 1.479) ? new float[46] : new float[39];
1612 if (std::abs(scEta) <= 1.479) {
1613 vals[0] = SCRawEnergy;
1614 vals[1] = scEta;
1615 vals[2] = scPhi;
1616 vals[3] = R9;
1617 vals[4] = E5x5Seed / SCRawEnergy;
1618 vals[5] = etawidth;
1619 vals[6] = phiwidth;
1620 vals[7] = NClusters;
1621 vals[8] = HoE;
1622 vals[9] = rho;
1623 vals[10] = vertices;
1624 vals[11] = EtaSeed - scEta;
1625 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1626 vals[13] = ESeed / SCRawEnergy;
1627 vals[14] = E3x3Seed / ESeed;
1628 vals[15] = E5x5Seed / ESeed;
1629 vals[16] = see;
1630 vals[17] = spp;
1631 vals[18] = sep;
1632 vals[19] = EMaxSeed / ESeed;
1633 vals[20] = E2ndSeed / ESeed;
1634 vals[21] = ETopSeed / ESeed;
1635 vals[22] = EBottomSeed / ESeed;
1636 vals[23] = ELeftSeed / ESeed;
1637 vals[24] = ERightSeed / ESeed;
1638 vals[25] = E2x5MaxSeed / ESeed;
1639 vals[26] = E2x5TopSeed / ESeed;
1640 vals[27] = E2x5BottomSeed / ESeed;
1641 vals[28] = E2x5LeftSeed / ESeed;
1642 vals[29] = E2x5RightSeed / ESeed;
1643 vals[30] = IsEcalDriven;
1644 vals[31] = GsfTrackPIn;
1645 vals[32] = fbrem;
1646 vals[33] = Charge;
1647 vals[34] = EoP;
1648 vals[35] = TrackMomentumError / GsfTrackPIn;
1649 vals[36] = EcalEnergyError / SCRawEnergy;
1650 vals[37] = Classification;
1651 vals[38] = IEtaSeed;
1652 vals[39] = IPhiSeed;
1653 vals[40] = ((int)IEtaSeed) % 5;
1654 vals[41] = ((int)IPhiSeed) % 2;
1655 vals[42] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
1656 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
1657 vals[43] = ((int)IPhiSeed) % 20;
1658 vals[44] = EtaCrySeed;
1659 vals[45] = PhiCrySeed;
1660 }
1661
1662 else {
1663 vals[0] = SCRawEnergy;
1664 vals[1] = scEta;
1665 vals[2] = scPhi;
1666 vals[3] = R9;
1667 vals[4] = E5x5Seed / SCRawEnergy;
1668 vals[5] = etawidth;
1669 vals[6] = phiwidth;
1670 vals[7] = NClusters;
1671 vals[8] = HoE;
1672 vals[9] = rho;
1673 vals[10] = vertices;
1674 vals[11] = EtaSeed - scEta;
1675 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1676 vals[13] = ESeed / SCRawEnergy;
1677 vals[14] = E3x3Seed / ESeed;
1678 vals[15] = E5x5Seed / ESeed;
1679 vals[16] = see;
1680 vals[17] = spp;
1681 vals[18] = sep;
1682 vals[19] = EMaxSeed / ESeed;
1683 vals[20] = E2ndSeed / ESeed;
1684 vals[21] = ETopSeed / ESeed;
1685 vals[22] = EBottomSeed / ESeed;
1686 vals[23] = ELeftSeed / ESeed;
1687 vals[24] = ERightSeed / ESeed;
1688 vals[25] = E2x5MaxSeed / ESeed;
1689 vals[26] = E2x5TopSeed / ESeed;
1690 vals[27] = E2x5BottomSeed / ESeed;
1691 vals[28] = E2x5LeftSeed / ESeed;
1692 vals[29] = E2x5RightSeed / ESeed;
1693 vals[30] = IsEcalDriven;
1694 vals[31] = GsfTrackPIn;
1695 vals[32] = fbrem;
1696 vals[33] = Charge;
1697 vals[34] = EoP;
1698 vals[35] = TrackMomentumError / GsfTrackPIn;
1699 vals[36] = EcalEnergyError / SCRawEnergy;
1700 vals[37] = Classification;
1701 vals[38] = PreShowerOverRaw;
1702 }
1703
1704
1705 double regressionResult = 0;
1706
1707 if (fVersionType == kWithTrkVarV1) {
1708 if (std::abs(scEta) <= 1.479)
1709 regressionResult = SCRawEnergy * forestCorrection_eb->GetResponse(vals);
1710 else
1711 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestCorrection_ee->GetResponse(vals);
1712 }
1713
1714
1715 if (printDebug) {
1716 if (std::abs(scEta) <= 1.479) {
1717 std::cout << "Barrel :";
1718 for (unsigned int v = 0; v < 46; ++v)
1719 std::cout << vals[v] << ", ";
1720 std::cout << "\n";
1721 } else {
1722 std::cout << "Endcap :";
1723 for (unsigned int v = 0; v < 39; ++v)
1724 std::cout << vals[v] << ", ";
1725 std::cout << "\n";
1726 }
1727 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
1728 std::cout << "regression energy = " << regressionResult << std::endl;
1729 }
1730
1731
1732 delete[] vals;
1733 return regressionResult;
1734 }
1735
1736 double ElectronEnergyRegressionEvaluate::regressionUncertaintyWithTrkVarV1(double SCRawEnergy,
1737 double scEta,
1738 double scPhi,
1739 double R9,
1740 double etawidth,
1741 double phiwidth,
1742 double NClusters,
1743 double HoE,
1744 double rho,
1745 double vertices,
1746 double EtaSeed,
1747 double PhiSeed,
1748 double ESeed,
1749 double E3x3Seed,
1750 double E5x5Seed,
1751 double see,
1752 double spp,
1753 double sep,
1754 double EMaxSeed,
1755 double E2ndSeed,
1756 double ETopSeed,
1757 double EBottomSeed,
1758 double ELeftSeed,
1759 double ERightSeed,
1760 double E2x5MaxSeed,
1761 double E2x5TopSeed,
1762 double E2x5BottomSeed,
1763 double E2x5LeftSeed,
1764 double E2x5RightSeed,
1765 double IEtaSeed,
1766 double IPhiSeed,
1767 double EtaCrySeed,
1768 double PhiCrySeed,
1769 double PreShowerOverRaw,
1770 int IsEcalDriven,
1771 double GsfTrackPIn,
1772 double fbrem,
1773 double Charge,
1774 double EoP,
1775 double TrackMomentumError,
1776 double EcalEnergyError,
1777 int Classification,
1778 bool printDebug) {
1779
1780 if (fIsInitialized == kFALSE) {
1781 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
1782 return 0;
1783 }
1784
1785
1786 assert(fVersionType == kWithTrkVarV1);
1787
1788 float *vals = (std::abs(scEta) <= 1.479) ? new float[46] : new float[39];
1789 if (std::abs(scEta) <= 1.479) {
1790 vals[0] = SCRawEnergy;
1791 vals[1] = scEta;
1792 vals[2] = scPhi;
1793 vals[3] = R9;
1794 vals[4] = E5x5Seed / SCRawEnergy;
1795 vals[5] = etawidth;
1796 vals[6] = phiwidth;
1797 vals[7] = NClusters;
1798 vals[8] = HoE;
1799 vals[9] = rho;
1800 vals[10] = vertices;
1801 vals[11] = EtaSeed - scEta;
1802 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1803 vals[13] = ESeed / SCRawEnergy;
1804 vals[14] = E3x3Seed / ESeed;
1805 vals[15] = E5x5Seed / ESeed;
1806 vals[16] = see;
1807 vals[17] = spp;
1808 vals[18] = sep;
1809 vals[19] = EMaxSeed / ESeed;
1810 vals[20] = E2ndSeed / ESeed;
1811 vals[21] = ETopSeed / ESeed;
1812 vals[22] = EBottomSeed / ESeed;
1813 vals[23] = ELeftSeed / ESeed;
1814 vals[24] = ERightSeed / ESeed;
1815 vals[25] = E2x5MaxSeed / ESeed;
1816 vals[26] = E2x5TopSeed / ESeed;
1817 vals[27] = E2x5BottomSeed / ESeed;
1818 vals[28] = E2x5LeftSeed / ESeed;
1819 vals[29] = E2x5RightSeed / ESeed;
1820 vals[30] = IsEcalDriven;
1821 vals[31] = GsfTrackPIn;
1822 vals[32] = fbrem;
1823 vals[33] = Charge;
1824 vals[34] = EoP;
1825 vals[35] = TrackMomentumError / GsfTrackPIn;
1826 vals[36] = EcalEnergyError / SCRawEnergy;
1827 vals[37] = Classification;
1828 vals[38] = IEtaSeed;
1829 vals[39] = IPhiSeed;
1830 vals[40] = ((int)IEtaSeed) % 5;
1831 vals[41] = ((int)IPhiSeed) % 2;
1832 vals[42] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
1833 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
1834 vals[43] = ((int)IPhiSeed) % 20;
1835 vals[44] = EtaCrySeed;
1836 vals[45] = PhiCrySeed;
1837 }
1838
1839 else {
1840 vals[0] = SCRawEnergy;
1841 vals[1] = scEta;
1842 vals[2] = scPhi;
1843 vals[3] = R9;
1844 vals[4] = E5x5Seed / SCRawEnergy;
1845 vals[5] = etawidth;
1846 vals[6] = phiwidth;
1847 vals[7] = NClusters;
1848 vals[8] = HoE;
1849 vals[9] = rho;
1850 vals[10] = vertices;
1851 vals[11] = EtaSeed - scEta;
1852 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1853 vals[13] = ESeed / SCRawEnergy;
1854 vals[14] = E3x3Seed / ESeed;
1855 vals[15] = E5x5Seed / ESeed;
1856 vals[16] = see;
1857 vals[17] = spp;
1858 vals[18] = sep;
1859 vals[19] = EMaxSeed / ESeed;
1860 vals[20] = E2ndSeed / ESeed;
1861 vals[21] = ETopSeed / ESeed;
1862 vals[22] = EBottomSeed / ESeed;
1863 vals[23] = ELeftSeed / ESeed;
1864 vals[24] = ERightSeed / ESeed;
1865 vals[25] = E2x5MaxSeed / ESeed;
1866 vals[26] = E2x5TopSeed / ESeed;
1867 vals[27] = E2x5BottomSeed / ESeed;
1868 vals[28] = E2x5LeftSeed / ESeed;
1869 vals[29] = E2x5RightSeed / ESeed;
1870 vals[30] = IsEcalDriven;
1871 vals[31] = GsfTrackPIn;
1872 vals[32] = fbrem;
1873 vals[33] = Charge;
1874 vals[34] = EoP;
1875 vals[35] = TrackMomentumError / GsfTrackPIn;
1876 vals[36] = EcalEnergyError / SCRawEnergy;
1877 vals[37] = Classification;
1878 vals[38] = PreShowerOverRaw;
1879 }
1880
1881
1882 double regressionResult = 0;
1883
1884 if (fVersionType == kWithTrkVarV1) {
1885 if (std::abs(scEta) <= 1.479)
1886 regressionResult = SCRawEnergy * forestUncertainty_eb->GetResponse(vals);
1887 else
1888 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestUncertainty_ee->GetResponse(vals);
1889 }
1890
1891
1892 if (printDebug) {
1893 if (std::abs(scEta) <= 1.479) {
1894 std::cout << "Barrel :";
1895 for (unsigned int v = 0; v < 46; ++v)
1896 std::cout << vals[v] << ", ";
1897 std::cout << "\n";
1898 } else {
1899 std::cout << "Endcap :";
1900 for (unsigned int v = 0; v < 39; ++v)
1901 std::cout << vals[v] << ", ";
1902 std::cout << "\n";
1903 }
1904 std::cout << " SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
1905 std::cout << "regression energy uncertainty = " << regressionResult << std::endl;
1906 }
1907
1908
1909 delete[] vals;
1910 return regressionResult;
1911 }
1912
1913 double ElectronEnergyRegressionEvaluate::regressionValueWithTrkVarV1(std::vector<double> &inputvars, bool printDebug) {
1914
1915 if (fIsInitialized == kFALSE) {
1916 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
1917 return 0;
1918 }
1919
1920
1921 assert(fVersionType == kWithTrkVarV1);
1922
1923
1924 assert(inputvars.size() == 42);
1925
1926 double SCRawEnergy = inputvars[0];
1927 double scEta = inputvars[1];
1928 double scPhi = inputvars[2];
1929 double R9 = inputvars[3];
1930 double etawidth = inputvars[4];
1931 double phiwidth = inputvars[5];
1932 double NClusters = inputvars[6];
1933 double HoE = inputvars[7];
1934 double rho = inputvars[8];
1935 double vertices = inputvars[9];
1936 double EtaSeed = inputvars[10];
1937 double PhiSeed = inputvars[11];
1938 double ESeed = inputvars[12];
1939 double E3x3Seed = inputvars[13];
1940 double E5x5Seed = inputvars[14];
1941 double see = inputvars[15];
1942 double spp = inputvars[16];
1943 double sep = inputvars[17];
1944 double EMaxSeed = inputvars[18];
1945 double E2ndSeed = inputvars[19];
1946 double ETopSeed = inputvars[20];
1947 double EBottomSeed = inputvars[21];
1948 double ELeftSeed = inputvars[22];
1949 double ERightSeed = inputvars[23];
1950 double E2x5MaxSeed = inputvars[24];
1951 double E2x5TopSeed = inputvars[25];
1952 double E2x5BottomSeed = inputvars[26];
1953 double E2x5LeftSeed = inputvars[27];
1954 double E2x5RightSeed = inputvars[28];
1955 double IEtaSeed = inputvars[29];
1956 double IPhiSeed = inputvars[30];
1957 double EtaCrySeed = inputvars[31];
1958 double PhiCrySeed = inputvars[32];
1959 double PreShowerOverRaw = inputvars[33];
1960 int IsEcalDriven = inputvars[34];
1961 double GsfTrackPIn = inputvars[35];
1962 double fbrem = inputvars[36];
1963 double Charge = inputvars[37];
1964 double EoP = inputvars[38];
1965 double TrackMomentumError = inputvars[39];
1966 double EcalEnergyError = inputvars[40];
1967 int Classification = inputvars[41];
1968
1969 float *vals = (std::abs(scEta) <= 1.479) ? new float[46] : new float[39];
1970 if (std::abs(scEta) <= 1.479) {
1971 vals[0] = SCRawEnergy;
1972 vals[1] = scEta;
1973 vals[2] = scPhi;
1974 vals[3] = R9;
1975 vals[4] = E5x5Seed / SCRawEnergy;
1976 vals[5] = etawidth;
1977 vals[6] = phiwidth;
1978 vals[7] = NClusters;
1979 vals[8] = HoE;
1980 vals[9] = rho;
1981 vals[10] = vertices;
1982 vals[11] = EtaSeed - scEta;
1983 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
1984 vals[13] = ESeed / SCRawEnergy;
1985 vals[14] = E3x3Seed / ESeed;
1986 vals[15] = E5x5Seed / ESeed;
1987 vals[16] = see;
1988 vals[17] = spp;
1989 vals[18] = sep;
1990 vals[19] = EMaxSeed / ESeed;
1991 vals[20] = E2ndSeed / ESeed;
1992 vals[21] = ETopSeed / ESeed;
1993 vals[22] = EBottomSeed / ESeed;
1994 vals[23] = ELeftSeed / ESeed;
1995 vals[24] = ERightSeed / ESeed;
1996 vals[25] = E2x5MaxSeed / ESeed;
1997 vals[26] = E2x5TopSeed / ESeed;
1998 vals[27] = E2x5BottomSeed / ESeed;
1999 vals[28] = E2x5LeftSeed / ESeed;
2000 vals[29] = E2x5RightSeed / ESeed;
2001 vals[30] = IsEcalDriven;
2002 vals[31] = GsfTrackPIn;
2003 vals[32] = fbrem;
2004 vals[33] = Charge;
2005 vals[34] = EoP;
2006 vals[35] = TrackMomentumError / GsfTrackPIn;
2007 vals[36] = EcalEnergyError / SCRawEnergy;
2008 vals[37] = Classification;
2009 vals[38] = IEtaSeed;
2010 vals[39] = IPhiSeed;
2011 vals[40] = ((int)IEtaSeed) % 5;
2012 vals[41] = ((int)IPhiSeed) % 2;
2013 vals[42] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
2014 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
2015 vals[43] = ((int)IPhiSeed) % 20;
2016 vals[44] = EtaCrySeed;
2017 vals[45] = PhiCrySeed;
2018 }
2019
2020 else {
2021 vals[0] = SCRawEnergy;
2022 vals[1] = scEta;
2023 vals[2] = scPhi;
2024 vals[3] = R9;
2025 vals[4] = E5x5Seed / SCRawEnergy;
2026 vals[5] = etawidth;
2027 vals[6] = phiwidth;
2028 vals[7] = NClusters;
2029 vals[8] = HoE;
2030 vals[9] = rho;
2031 vals[10] = vertices;
2032 vals[11] = EtaSeed - scEta;
2033 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2034 vals[13] = ESeed / SCRawEnergy;
2035 vals[14] = E3x3Seed / ESeed;
2036 vals[15] = E5x5Seed / ESeed;
2037 vals[16] = see;
2038 vals[17] = spp;
2039 vals[18] = sep;
2040 vals[19] = EMaxSeed / ESeed;
2041 vals[20] = E2ndSeed / ESeed;
2042 vals[21] = ETopSeed / ESeed;
2043 vals[22] = EBottomSeed / ESeed;
2044 vals[23] = ELeftSeed / ESeed;
2045 vals[24] = ERightSeed / ESeed;
2046 vals[25] = E2x5MaxSeed / ESeed;
2047 vals[26] = E2x5TopSeed / ESeed;
2048 vals[27] = E2x5BottomSeed / ESeed;
2049 vals[28] = E2x5LeftSeed / ESeed;
2050 vals[29] = E2x5RightSeed / ESeed;
2051 vals[30] = IsEcalDriven;
2052 vals[31] = GsfTrackPIn;
2053 vals[32] = fbrem;
2054 vals[33] = Charge;
2055 vals[34] = EoP;
2056 vals[35] = TrackMomentumError / GsfTrackPIn;
2057 vals[36] = EcalEnergyError / SCRawEnergy;
2058 vals[37] = Classification;
2059 vals[38] = PreShowerOverRaw;
2060 }
2061
2062
2063 double regressionResult = 0;
2064
2065 if (fVersionType == kWithTrkVarV1) {
2066 if (std::abs(scEta) <= 1.479)
2067 regressionResult = SCRawEnergy * forestCorrection_eb->GetResponse(vals);
2068 else
2069 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestCorrection_ee->GetResponse(vals);
2070 }
2071
2072
2073 if (printDebug) {
2074 if (std::abs(scEta) <= 1.479) {
2075 std::cout << "Barrel :";
2076 for (unsigned int v = 0; v < 46; ++v)
2077 std::cout << vals[v] << ", ";
2078 std::cout << "\n";
2079 } else {
2080 std::cout << "Endcap :";
2081 for (unsigned int v = 0; v < 39; ++v)
2082 std::cout << vals[v] << ", ";
2083 std::cout << "\n";
2084 }
2085 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
2086 std::cout << "regression energy = " << regressionResult << std::endl;
2087 }
2088
2089
2090 delete[] vals;
2091 return regressionResult;
2092 }
2093
2094 double ElectronEnergyRegressionEvaluate::regressionUncertaintyWithTrkVarV1(std::vector<double> &inputvars,
2095 bool printDebug) {
2096
2097 if (fIsInitialized == kFALSE) {
2098 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
2099 return 0;
2100 }
2101
2102
2103 assert(fVersionType == kWithTrkVarV1);
2104
2105
2106 assert(inputvars.size() == 42);
2107
2108 double SCRawEnergy = inputvars[0];
2109 double scEta = inputvars[1];
2110 double scPhi = inputvars[2];
2111 double R9 = inputvars[3];
2112 double etawidth = inputvars[4];
2113 double phiwidth = inputvars[5];
2114 double NClusters = inputvars[6];
2115 double HoE = inputvars[7];
2116 double rho = inputvars[8];
2117 double vertices = inputvars[9];
2118 double EtaSeed = inputvars[10];
2119 double PhiSeed = inputvars[11];
2120 double ESeed = inputvars[12];
2121 double E3x3Seed = inputvars[13];
2122 double E5x5Seed = inputvars[14];
2123 double see = inputvars[15];
2124 double spp = inputvars[16];
2125 double sep = inputvars[17];
2126 double EMaxSeed = inputvars[18];
2127 double E2ndSeed = inputvars[19];
2128 double ETopSeed = inputvars[20];
2129 double EBottomSeed = inputvars[21];
2130 double ELeftSeed = inputvars[22];
2131 double ERightSeed = inputvars[23];
2132 double E2x5MaxSeed = inputvars[24];
2133 double E2x5TopSeed = inputvars[25];
2134 double E2x5BottomSeed = inputvars[26];
2135 double E2x5LeftSeed = inputvars[27];
2136 double E2x5RightSeed = inputvars[28];
2137 double IEtaSeed = inputvars[29];
2138 double IPhiSeed = inputvars[30];
2139 double EtaCrySeed = inputvars[31];
2140 double PhiCrySeed = inputvars[32];
2141 double PreShowerOverRaw = inputvars[33];
2142 int IsEcalDriven = inputvars[34];
2143 double GsfTrackPIn = inputvars[35];
2144 double fbrem = inputvars[36];
2145 double Charge = inputvars[37];
2146 double EoP = inputvars[38];
2147 double TrackMomentumError = inputvars[39];
2148 double EcalEnergyError = inputvars[40];
2149 int Classification = inputvars[41];
2150
2151 float *vals = (std::abs(scEta) <= 1.479) ? new float[46] : new float[39];
2152 if (std::abs(scEta) <= 1.479) {
2153 vals[0] = SCRawEnergy;
2154 vals[1] = scEta;
2155 vals[2] = scPhi;
2156 vals[3] = R9;
2157 vals[4] = E5x5Seed / SCRawEnergy;
2158 vals[5] = etawidth;
2159 vals[6] = phiwidth;
2160 vals[7] = NClusters;
2161 vals[8] = HoE;
2162 vals[9] = rho;
2163 vals[10] = vertices;
2164 vals[11] = EtaSeed - scEta;
2165 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2166 vals[13] = ESeed / SCRawEnergy;
2167 vals[14] = E3x3Seed / ESeed;
2168 vals[15] = E5x5Seed / ESeed;
2169 vals[16] = see;
2170 vals[17] = spp;
2171 vals[18] = sep;
2172 vals[19] = EMaxSeed / ESeed;
2173 vals[20] = E2ndSeed / ESeed;
2174 vals[21] = ETopSeed / ESeed;
2175 vals[22] = EBottomSeed / ESeed;
2176 vals[23] = ELeftSeed / ESeed;
2177 vals[24] = ERightSeed / ESeed;
2178 vals[25] = E2x5MaxSeed / ESeed;
2179 vals[26] = E2x5TopSeed / ESeed;
2180 vals[27] = E2x5BottomSeed / ESeed;
2181 vals[28] = E2x5LeftSeed / ESeed;
2182 vals[29] = E2x5RightSeed / ESeed;
2183 vals[30] = IsEcalDriven;
2184 vals[31] = GsfTrackPIn;
2185 vals[32] = fbrem;
2186 vals[33] = Charge;
2187 vals[34] = EoP;
2188 vals[35] = TrackMomentumError / GsfTrackPIn;
2189 vals[36] = EcalEnergyError / SCRawEnergy;
2190 vals[37] = Classification;
2191 vals[38] = IEtaSeed;
2192 vals[39] = IPhiSeed;
2193 vals[40] = ((int)IEtaSeed) % 5;
2194 vals[41] = ((int)IPhiSeed) % 2;
2195 vals[42] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
2196 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
2197 vals[43] = ((int)IPhiSeed) % 20;
2198 vals[44] = EtaCrySeed;
2199 vals[45] = PhiCrySeed;
2200 }
2201
2202 else {
2203 vals[0] = SCRawEnergy;
2204 vals[1] = scEta;
2205 vals[2] = scPhi;
2206 vals[3] = R9;
2207 vals[4] = E5x5Seed / SCRawEnergy;
2208 vals[5] = etawidth;
2209 vals[6] = phiwidth;
2210 vals[7] = NClusters;
2211 vals[8] = HoE;
2212 vals[9] = rho;
2213 vals[10] = vertices;
2214 vals[11] = EtaSeed - scEta;
2215 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2216 vals[13] = ESeed / SCRawEnergy;
2217 vals[14] = E3x3Seed / ESeed;
2218 vals[15] = E5x5Seed / ESeed;
2219 vals[16] = see;
2220 vals[17] = spp;
2221 vals[18] = sep;
2222 vals[19] = EMaxSeed / ESeed;
2223 vals[20] = E2ndSeed / ESeed;
2224 vals[21] = ETopSeed / ESeed;
2225 vals[22] = EBottomSeed / ESeed;
2226 vals[23] = ELeftSeed / ESeed;
2227 vals[24] = ERightSeed / ESeed;
2228 vals[25] = E2x5MaxSeed / ESeed;
2229 vals[26] = E2x5TopSeed / ESeed;
2230 vals[27] = E2x5BottomSeed / ESeed;
2231 vals[28] = E2x5LeftSeed / ESeed;
2232 vals[29] = E2x5RightSeed / ESeed;
2233 vals[30] = IsEcalDriven;
2234 vals[31] = GsfTrackPIn;
2235 vals[32] = fbrem;
2236 vals[33] = Charge;
2237 vals[34] = EoP;
2238 vals[35] = TrackMomentumError / GsfTrackPIn;
2239 vals[36] = EcalEnergyError / SCRawEnergy;
2240 vals[37] = Classification;
2241 vals[38] = PreShowerOverRaw;
2242 }
2243
2244
2245 double regressionResult = 0;
2246
2247 if (fVersionType == kWithTrkVarV1) {
2248 if (std::abs(scEta) <= 1.479)
2249 regressionResult = SCRawEnergy * forestUncertainty_eb->GetResponse(vals);
2250 else
2251 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestUncertainty_ee->GetResponse(vals);
2252 }
2253
2254
2255 if (printDebug) {
2256 if (std::abs(scEta) <= 1.479) {
2257 std::cout << "Barrel :";
2258 for (unsigned int v = 0; v < 46; ++v)
2259 std::cout << vals[v] << ", ";
2260 std::cout << "\n";
2261 } else {
2262 std::cout << "Endcap :";
2263 for (unsigned int v = 0; v < 39; ++v)
2264 std::cout << vals[v] << ", ";
2265 std::cout << "\n";
2266 }
2267 std::cout << " SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
2268 std::cout << "regression energy uncertainty = " << regressionResult << std::endl;
2269 }
2270
2271
2272 delete[] vals;
2273 return regressionResult;
2274 }
2275
2276 double ElectronEnergyRegressionEvaluate::regressionValueWithTrkVarV2(double SCRawEnergy,
2277 double scEta,
2278 double scPhi,
2279 double R9,
2280 double etawidth,
2281 double phiwidth,
2282 double NClusters,
2283 double HoE,
2284 double rho,
2285 double vertices,
2286 double EtaSeed,
2287 double PhiSeed,
2288 double ESeed,
2289 double E3x3Seed,
2290 double E5x5Seed,
2291 double see,
2292 double spp,
2293 double sep,
2294 double EMaxSeed,
2295 double E2ndSeed,
2296 double ETopSeed,
2297 double EBottomSeed,
2298 double ELeftSeed,
2299 double ERightSeed,
2300 double E2x5MaxSeed,
2301 double E2x5TopSeed,
2302 double E2x5BottomSeed,
2303 double E2x5LeftSeed,
2304 double E2x5RightSeed,
2305 double IEtaSeed,
2306 double IPhiSeed,
2307 double EtaCrySeed,
2308 double PhiCrySeed,
2309 double PreShowerOverRaw,
2310 int IsEcalDriven,
2311 double GsfTrackPIn,
2312 double fbrem,
2313 double Charge,
2314 double EoP,
2315 double TrackMomentumError,
2316 double EcalEnergyError,
2317 int Classification,
2318 double detaIn,
2319 double dphiIn,
2320 double detaCalo,
2321 double dphiCalo,
2322 double GsfTrackChiSqr,
2323 double KFTrackNLayers,
2324 double ElectronEnergyOverPout,
2325 bool printDebug) {
2326
2327 if (fIsInitialized == kFALSE) {
2328 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
2329 return 0;
2330 }
2331
2332
2333 assert(fVersionType == kWithTrkVarV2);
2334
2335 float *vals = (std::abs(scEta) <= 1.479) ? new float[53] : new float[46];
2336 if (std::abs(scEta) <= 1.479) {
2337 vals[0] = SCRawEnergy;
2338 vals[1] = scEta;
2339 vals[2] = scPhi;
2340 vals[3] = R9;
2341 vals[4] = E5x5Seed / SCRawEnergy;
2342 vals[5] = etawidth;
2343 vals[6] = phiwidth;
2344 vals[7] = NClusters;
2345 vals[8] = HoE;
2346 vals[9] = rho;
2347 vals[10] = vertices;
2348 vals[11] = EtaSeed - scEta;
2349 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2350 vals[13] = ESeed / SCRawEnergy;
2351 vals[14] = E3x3Seed / ESeed;
2352 vals[15] = E5x5Seed / ESeed;
2353 vals[16] = see;
2354 vals[17] = spp;
2355 vals[18] = sep;
2356 vals[19] = EMaxSeed / ESeed;
2357 vals[20] = E2ndSeed / ESeed;
2358 vals[21] = ETopSeed / ESeed;
2359 vals[22] = EBottomSeed / ESeed;
2360 vals[23] = ELeftSeed / ESeed;
2361 vals[24] = ERightSeed / ESeed;
2362 vals[25] = E2x5MaxSeed / ESeed;
2363 vals[26] = E2x5TopSeed / ESeed;
2364 vals[27] = E2x5BottomSeed / ESeed;
2365 vals[28] = E2x5LeftSeed / ESeed;
2366 vals[29] = E2x5RightSeed / ESeed;
2367 vals[30] = IsEcalDriven;
2368 vals[31] = GsfTrackPIn;
2369 vals[32] = fbrem;
2370 vals[33] = Charge;
2371 vals[34] = EoP;
2372 vals[35] = TrackMomentumError / GsfTrackPIn;
2373 vals[36] = EcalEnergyError / SCRawEnergy;
2374 vals[37] = Classification;
2375 vals[38] = detaIn;
2376 vals[39] = dphiIn;
2377 vals[40] = detaCalo;
2378 vals[41] = dphiCalo;
2379 vals[42] = GsfTrackChiSqr;
2380 vals[43] = KFTrackNLayers;
2381 vals[44] = ElectronEnergyOverPout;
2382 vals[45] = IEtaSeed;
2383 vals[46] = IPhiSeed;
2384 vals[47] = ((int)IEtaSeed) % 5;
2385 vals[48] = ((int)IPhiSeed) % 2;
2386 vals[49] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
2387 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
2388 vals[50] = ((int)IPhiSeed) % 20;
2389 vals[51] = EtaCrySeed;
2390 vals[52] = PhiCrySeed;
2391 }
2392
2393 else {
2394 vals[0] = SCRawEnergy;
2395 vals[1] = scEta;
2396 vals[2] = scPhi;
2397 vals[3] = R9;
2398 vals[4] = E5x5Seed / SCRawEnergy;
2399 vals[5] = etawidth;
2400 vals[6] = phiwidth;
2401 vals[7] = NClusters;
2402 vals[8] = HoE;
2403 vals[9] = rho;
2404 vals[10] = vertices;
2405 vals[11] = EtaSeed - scEta;
2406 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2407 vals[13] = ESeed / SCRawEnergy;
2408 vals[14] = E3x3Seed / ESeed;
2409 vals[15] = E5x5Seed / ESeed;
2410 vals[16] = see;
2411 vals[17] = spp;
2412 vals[18] = sep;
2413 vals[19] = EMaxSeed / ESeed;
2414 vals[20] = E2ndSeed / ESeed;
2415 vals[21] = ETopSeed / ESeed;
2416 vals[22] = EBottomSeed / ESeed;
2417 vals[23] = ELeftSeed / ESeed;
2418 vals[24] = ERightSeed / ESeed;
2419 vals[25] = E2x5MaxSeed / ESeed;
2420 vals[26] = E2x5TopSeed / ESeed;
2421 vals[27] = E2x5BottomSeed / ESeed;
2422 vals[28] = E2x5LeftSeed / ESeed;
2423 vals[29] = E2x5RightSeed / ESeed;
2424 vals[30] = IsEcalDriven;
2425 vals[31] = GsfTrackPIn;
2426 vals[32] = fbrem;
2427 vals[33] = Charge;
2428 vals[34] = EoP;
2429 vals[35] = TrackMomentumError / GsfTrackPIn;
2430 vals[36] = EcalEnergyError / SCRawEnergy;
2431 vals[37] = Classification;
2432 vals[38] = detaIn;
2433 vals[39] = dphiIn;
2434 vals[40] = detaCalo;
2435 vals[41] = dphiCalo;
2436 vals[42] = GsfTrackChiSqr;
2437 vals[43] = KFTrackNLayers;
2438 vals[44] = ElectronEnergyOverPout;
2439 vals[45] = PreShowerOverRaw;
2440 }
2441
2442
2443 double regressionResult = 0;
2444
2445 if (fVersionType == kWithTrkVarV2) {
2446 if (std::abs(scEta) <= 1.479)
2447 regressionResult = SCRawEnergy * forestCorrection_eb->GetResponse(vals);
2448 else
2449 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestCorrection_ee->GetResponse(vals);
2450 }
2451
2452
2453 if (printDebug) {
2454 if (std::abs(scEta) <= 1.479) {
2455 std::cout << "Barrel :";
2456 for (unsigned int v = 0; v < 53; ++v)
2457 std::cout << vals[v] << ", ";
2458 std::cout << "\n";
2459 } else {
2460 std::cout << "Endcap :";
2461 for (unsigned int v = 0; v < 46; ++v)
2462 std::cout << vals[v] << ", ";
2463 std::cout << "\n";
2464 }
2465 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
2466 std::cout << "regression energy = " << regressionResult << std::endl;
2467 }
2468
2469
2470 delete[] vals;
2471 return regressionResult;
2472 }
2473
2474 double ElectronEnergyRegressionEvaluate::regressionUncertaintyWithTrkVarV2(double SCRawEnergy,
2475 double scEta,
2476 double scPhi,
2477 double R9,
2478 double etawidth,
2479 double phiwidth,
2480 double NClusters,
2481 double HoE,
2482 double rho,
2483 double vertices,
2484 double EtaSeed,
2485 double PhiSeed,
2486 double ESeed,
2487 double E3x3Seed,
2488 double E5x5Seed,
2489 double see,
2490 double spp,
2491 double sep,
2492 double EMaxSeed,
2493 double E2ndSeed,
2494 double ETopSeed,
2495 double EBottomSeed,
2496 double ELeftSeed,
2497 double ERightSeed,
2498 double E2x5MaxSeed,
2499 double E2x5TopSeed,
2500 double E2x5BottomSeed,
2501 double E2x5LeftSeed,
2502 double E2x5RightSeed,
2503 double IEtaSeed,
2504 double IPhiSeed,
2505 double EtaCrySeed,
2506 double PhiCrySeed,
2507 double PreShowerOverRaw,
2508 int IsEcalDriven,
2509 double GsfTrackPIn,
2510 double fbrem,
2511 double Charge,
2512 double EoP,
2513 double TrackMomentumError,
2514 double EcalEnergyError,
2515 int Classification,
2516 double detaIn,
2517 double dphiIn,
2518 double detaCalo,
2519 double dphiCalo,
2520 double GsfTrackChiSqr,
2521 double KFTrackNLayers,
2522 double ElectronEnergyOverPout,
2523 bool printDebug) {
2524
2525 if (fIsInitialized == kFALSE) {
2526 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
2527 return 0;
2528 }
2529
2530
2531 assert(fVersionType == kWithTrkVarV2);
2532
2533 float *vals = (std::abs(scEta) <= 1.479) ? new float[53] : new float[46];
2534 if (std::abs(scEta) <= 1.479) {
2535 vals[0] = SCRawEnergy;
2536 vals[1] = scEta;
2537 vals[2] = scPhi;
2538 vals[3] = R9;
2539 vals[4] = E5x5Seed / SCRawEnergy;
2540 vals[5] = etawidth;
2541 vals[6] = phiwidth;
2542 vals[7] = NClusters;
2543 vals[8] = HoE;
2544 vals[9] = rho;
2545 vals[10] = vertices;
2546 vals[11] = EtaSeed - scEta;
2547 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2548 vals[13] = ESeed / SCRawEnergy;
2549 vals[14] = E3x3Seed / ESeed;
2550 vals[15] = E5x5Seed / ESeed;
2551 vals[16] = see;
2552 vals[17] = spp;
2553 vals[18] = sep;
2554 vals[19] = EMaxSeed / ESeed;
2555 vals[20] = E2ndSeed / ESeed;
2556 vals[21] = ETopSeed / ESeed;
2557 vals[22] = EBottomSeed / ESeed;
2558 vals[23] = ELeftSeed / ESeed;
2559 vals[24] = ERightSeed / ESeed;
2560 vals[25] = E2x5MaxSeed / ESeed;
2561 vals[26] = E2x5TopSeed / ESeed;
2562 vals[27] = E2x5BottomSeed / ESeed;
2563 vals[28] = E2x5LeftSeed / ESeed;
2564 vals[29] = E2x5RightSeed / ESeed;
2565 vals[30] = IsEcalDriven;
2566 vals[31] = GsfTrackPIn;
2567 vals[32] = fbrem;
2568 vals[33] = Charge;
2569 vals[34] = EoP;
2570 vals[35] = TrackMomentumError / GsfTrackPIn;
2571 vals[36] = EcalEnergyError / SCRawEnergy;
2572 vals[37] = Classification;
2573 vals[38] = detaIn;
2574 vals[39] = dphiIn;
2575 vals[40] = detaCalo;
2576 vals[41] = dphiCalo;
2577 vals[42] = GsfTrackChiSqr;
2578 vals[43] = KFTrackNLayers;
2579 vals[44] = ElectronEnergyOverPout;
2580 vals[45] = IEtaSeed;
2581 vals[46] = IPhiSeed;
2582 vals[47] = ((int)IEtaSeed) % 5;
2583 vals[48] = ((int)IPhiSeed) % 2;
2584 vals[49] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
2585 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
2586 vals[50] = ((int)IPhiSeed) % 20;
2587 vals[51] = EtaCrySeed;
2588 vals[52] = PhiCrySeed;
2589 }
2590
2591 else {
2592 vals[0] = SCRawEnergy;
2593 vals[1] = scEta;
2594 vals[2] = scPhi;
2595 vals[3] = R9;
2596 vals[4] = E5x5Seed / SCRawEnergy;
2597 vals[5] = etawidth;
2598 vals[6] = phiwidth;
2599 vals[7] = NClusters;
2600 vals[8] = HoE;
2601 vals[9] = rho;
2602 vals[10] = vertices;
2603 vals[11] = EtaSeed - scEta;
2604 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2605 vals[13] = ESeed / SCRawEnergy;
2606 vals[14] = E3x3Seed / ESeed;
2607 vals[15] = E5x5Seed / ESeed;
2608 vals[16] = see;
2609 vals[17] = spp;
2610 vals[18] = sep;
2611 vals[19] = EMaxSeed / ESeed;
2612 vals[20] = E2ndSeed / ESeed;
2613 vals[21] = ETopSeed / ESeed;
2614 vals[22] = EBottomSeed / ESeed;
2615 vals[23] = ELeftSeed / ESeed;
2616 vals[24] = ERightSeed / ESeed;
2617 vals[25] = E2x5MaxSeed / ESeed;
2618 vals[26] = E2x5TopSeed / ESeed;
2619 vals[27] = E2x5BottomSeed / ESeed;
2620 vals[28] = E2x5LeftSeed / ESeed;
2621 vals[29] = E2x5RightSeed / ESeed;
2622 vals[30] = IsEcalDriven;
2623 vals[31] = GsfTrackPIn;
2624 vals[32] = fbrem;
2625 vals[33] = Charge;
2626 vals[34] = EoP;
2627 vals[35] = TrackMomentumError / GsfTrackPIn;
2628 vals[36] = EcalEnergyError / SCRawEnergy;
2629 vals[37] = Classification;
2630 vals[38] = detaIn;
2631 vals[39] = dphiIn;
2632 vals[40] = detaCalo;
2633 vals[41] = dphiCalo;
2634 vals[42] = GsfTrackChiSqr;
2635 vals[43] = KFTrackNLayers;
2636 vals[44] = ElectronEnergyOverPout;
2637 vals[45] = PreShowerOverRaw;
2638 }
2639
2640
2641 double regressionResult = 0;
2642
2643 if (fVersionType == kWithTrkVarV2) {
2644 if (std::abs(scEta) <= 1.479)
2645 regressionResult = SCRawEnergy * forestUncertainty_eb->GetResponse(vals);
2646 else
2647 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestUncertainty_ee->GetResponse(vals);
2648 }
2649
2650
2651 if (printDebug) {
2652 if (std::abs(scEta) <= 1.479) {
2653 std::cout << "Barrel :";
2654 for (unsigned int v = 0; v < 53; ++v)
2655 std::cout << vals[v] << ", ";
2656 std::cout << "\n";
2657 } else {
2658 std::cout << "Endcap :";
2659 for (unsigned int v = 0; v < 46; ++v)
2660 std::cout << vals[v] << ", ";
2661 std::cout << "\n";
2662 }
2663 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
2664 std::cout << "regression energy uncertainty = " << regressionResult << std::endl;
2665 }
2666
2667
2668 delete[] vals;
2669 return regressionResult;
2670 }
2671
2672 double ElectronEnergyRegressionEvaluate::regressionValueWithTrkVarV2(std::vector<double> &inputvars, bool printDebug) {
2673
2674 if (fIsInitialized == kFALSE) {
2675 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
2676 return 0;
2677 }
2678
2679
2680 assert(fVersionType == kWithTrkVarV2);
2681
2682
2683 assert(inputvars.size() == 49);
2684
2685 double SCRawEnergy = inputvars[0];
2686 double scEta = inputvars[1];
2687 double scPhi = inputvars[2];
2688 double R9 = inputvars[3];
2689 double etawidth = inputvars[4];
2690 double phiwidth = inputvars[5];
2691 double NClusters = inputvars[6];
2692 double HoE = inputvars[7];
2693 double rho = inputvars[8];
2694 double vertices = inputvars[9];
2695 double EtaSeed = inputvars[10];
2696 double PhiSeed = inputvars[11];
2697 double ESeed = inputvars[12];
2698 double E3x3Seed = inputvars[13];
2699 double E5x5Seed = inputvars[14];
2700 double see = inputvars[15];
2701 double spp = inputvars[16];
2702 double sep = inputvars[17];
2703 double EMaxSeed = inputvars[18];
2704 double E2ndSeed = inputvars[19];
2705 double ETopSeed = inputvars[20];
2706 double EBottomSeed = inputvars[21];
2707 double ELeftSeed = inputvars[22];
2708 double ERightSeed = inputvars[23];
2709 double E2x5MaxSeed = inputvars[24];
2710 double E2x5TopSeed = inputvars[25];
2711 double E2x5BottomSeed = inputvars[26];
2712 double E2x5LeftSeed = inputvars[27];
2713 double E2x5RightSeed = inputvars[28];
2714 double IEtaSeed = inputvars[29];
2715 double IPhiSeed = inputvars[30];
2716 double EtaCrySeed = inputvars[31];
2717 double PhiCrySeed = inputvars[32];
2718 double PreShowerOverRaw = inputvars[33];
2719 int IsEcalDriven = inputvars[34];
2720 double GsfTrackPIn = inputvars[35];
2721 double fbrem = inputvars[36];
2722 double Charge = inputvars[37];
2723 double EoP = inputvars[38];
2724 double TrackMomentumError = inputvars[39];
2725 double EcalEnergyError = inputvars[40];
2726 int Classification = inputvars[41];
2727 double detaIn = inputvars[42];
2728 double dphiIn = inputvars[43];
2729 double detaCalo = inputvars[44];
2730 double dphiCalo = inputvars[45];
2731 double GsfTrackChiSqr = inputvars[46];
2732 double KFTrackNLayers = inputvars[47];
2733 double ElectronEnergyOverPout = inputvars[48];
2734
2735 float *vals = (std::abs(scEta) <= 1.479) ? new float[53] : new float[46];
2736 if (std::abs(scEta) <= 1.479) {
2737 vals[0] = SCRawEnergy;
2738 vals[1] = scEta;
2739 vals[2] = scPhi;
2740 vals[3] = R9;
2741 vals[4] = E5x5Seed / SCRawEnergy;
2742 vals[5] = etawidth;
2743 vals[6] = phiwidth;
2744 vals[7] = NClusters;
2745 vals[8] = HoE;
2746 vals[9] = rho;
2747 vals[10] = vertices;
2748 vals[11] = EtaSeed - scEta;
2749 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2750 vals[13] = ESeed / SCRawEnergy;
2751 vals[14] = E3x3Seed / ESeed;
2752 vals[15] = E5x5Seed / ESeed;
2753 vals[16] = see;
2754 vals[17] = spp;
2755 vals[18] = sep;
2756 vals[19] = EMaxSeed / ESeed;
2757 vals[20] = E2ndSeed / ESeed;
2758 vals[21] = ETopSeed / ESeed;
2759 vals[22] = EBottomSeed / ESeed;
2760 vals[23] = ELeftSeed / ESeed;
2761 vals[24] = ERightSeed / ESeed;
2762 vals[25] = E2x5MaxSeed / ESeed;
2763 vals[26] = E2x5TopSeed / ESeed;
2764 vals[27] = E2x5BottomSeed / ESeed;
2765 vals[28] = E2x5LeftSeed / ESeed;
2766 vals[29] = E2x5RightSeed / ESeed;
2767 vals[30] = IsEcalDriven;
2768 vals[31] = GsfTrackPIn;
2769 vals[32] = fbrem;
2770 vals[33] = Charge;
2771 vals[34] = EoP;
2772 vals[35] = TrackMomentumError / GsfTrackPIn;
2773 vals[36] = EcalEnergyError / SCRawEnergy;
2774 vals[37] = Classification;
2775 vals[38] = detaIn;
2776 vals[39] = dphiIn;
2777 vals[40] = detaCalo;
2778 vals[41] = dphiCalo;
2779 vals[42] = GsfTrackChiSqr;
2780 vals[43] = KFTrackNLayers;
2781 vals[44] = ElectronEnergyOverPout;
2782 vals[45] = IEtaSeed;
2783 vals[46] = IPhiSeed;
2784 vals[47] = ((int)IEtaSeed) % 5;
2785 vals[48] = ((int)IPhiSeed) % 2;
2786 vals[49] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
2787 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
2788 vals[50] = ((int)IPhiSeed) % 20;
2789 vals[51] = EtaCrySeed;
2790 vals[52] = PhiCrySeed;
2791 }
2792
2793 else {
2794 vals[0] = SCRawEnergy;
2795 vals[1] = scEta;
2796 vals[2] = scPhi;
2797 vals[3] = R9;
2798 vals[4] = E5x5Seed / SCRawEnergy;
2799 vals[5] = etawidth;
2800 vals[6] = phiwidth;
2801 vals[7] = NClusters;
2802 vals[8] = HoE;
2803 vals[9] = rho;
2804 vals[10] = vertices;
2805 vals[11] = EtaSeed - scEta;
2806 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2807 vals[13] = ESeed / SCRawEnergy;
2808 vals[14] = E3x3Seed / ESeed;
2809 vals[15] = E5x5Seed / ESeed;
2810 vals[16] = see;
2811 vals[17] = spp;
2812 vals[18] = sep;
2813 vals[19] = EMaxSeed / ESeed;
2814 vals[20] = E2ndSeed / ESeed;
2815 vals[21] = ETopSeed / ESeed;
2816 vals[22] = EBottomSeed / ESeed;
2817 vals[23] = ELeftSeed / ESeed;
2818 vals[24] = ERightSeed / ESeed;
2819 vals[25] = E2x5MaxSeed / ESeed;
2820 vals[26] = E2x5TopSeed / ESeed;
2821 vals[27] = E2x5BottomSeed / ESeed;
2822 vals[28] = E2x5LeftSeed / ESeed;
2823 vals[29] = E2x5RightSeed / ESeed;
2824 vals[30] = IsEcalDriven;
2825 vals[31] = GsfTrackPIn;
2826 vals[32] = fbrem;
2827 vals[33] = Charge;
2828 vals[34] = EoP;
2829 vals[35] = TrackMomentumError / GsfTrackPIn;
2830 vals[36] = EcalEnergyError / SCRawEnergy;
2831 vals[37] = Classification;
2832 vals[38] = detaIn;
2833 vals[39] = dphiIn;
2834 vals[40] = detaCalo;
2835 vals[41] = dphiCalo;
2836 vals[42] = GsfTrackChiSqr;
2837 vals[43] = KFTrackNLayers;
2838 vals[44] = ElectronEnergyOverPout;
2839 vals[45] = PreShowerOverRaw;
2840 }
2841
2842
2843 double regressionResult = 0;
2844
2845 if (fVersionType == kWithTrkVarV2) {
2846 if (std::abs(scEta) <= 1.479)
2847 regressionResult = SCRawEnergy * forestCorrection_eb->GetResponse(vals);
2848 else
2849 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestCorrection_ee->GetResponse(vals);
2850 }
2851
2852
2853 if (printDebug) {
2854 if (std::abs(scEta) <= 1.479) {
2855 std::cout << "Barrel :";
2856 for (unsigned int v = 0; v < 53; ++v)
2857 std::cout << vals[v] << ", ";
2858 std::cout << "\n";
2859 } else {
2860 std::cout << "Endcap :";
2861 for (unsigned int v = 0; v < 46; ++v)
2862 std::cout << vals[v] << ", ";
2863 std::cout << "\n";
2864 }
2865 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
2866 std::cout << "regression energy = " << regressionResult << std::endl;
2867 }
2868
2869
2870 delete[] vals;
2871 return regressionResult;
2872 }
2873
2874 double ElectronEnergyRegressionEvaluate::regressionUncertaintyWithTrkVarV2(std::vector<double> &inputvars,
2875 bool printDebug) {
2876
2877 if (fIsInitialized == kFALSE) {
2878 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
2879 return 0;
2880 }
2881
2882
2883 assert(fVersionType == kWithTrkVarV2);
2884
2885
2886 assert(inputvars.size() == 49);
2887
2888 double SCRawEnergy = inputvars[0];
2889 double scEta = inputvars[1];
2890 double scPhi = inputvars[2];
2891 double R9 = inputvars[3];
2892 double etawidth = inputvars[4];
2893 double phiwidth = inputvars[5];
2894 double NClusters = inputvars[6];
2895 double HoE = inputvars[7];
2896 double rho = inputvars[8];
2897 double vertices = inputvars[9];
2898 double EtaSeed = inputvars[10];
2899 double PhiSeed = inputvars[11];
2900 double ESeed = inputvars[12];
2901 double E3x3Seed = inputvars[13];
2902 double E5x5Seed = inputvars[14];
2903 double see = inputvars[15];
2904 double spp = inputvars[16];
2905 double sep = inputvars[17];
2906 double EMaxSeed = inputvars[18];
2907 double E2ndSeed = inputvars[19];
2908 double ETopSeed = inputvars[20];
2909 double EBottomSeed = inputvars[21];
2910 double ELeftSeed = inputvars[22];
2911 double ERightSeed = inputvars[23];
2912 double E2x5MaxSeed = inputvars[24];
2913 double E2x5TopSeed = inputvars[25];
2914 double E2x5BottomSeed = inputvars[26];
2915 double E2x5LeftSeed = inputvars[27];
2916 double E2x5RightSeed = inputvars[28];
2917 double IEtaSeed = inputvars[29];
2918 double IPhiSeed = inputvars[30];
2919 double EtaCrySeed = inputvars[31];
2920 double PhiCrySeed = inputvars[32];
2921 double PreShowerOverRaw = inputvars[33];
2922 int IsEcalDriven = inputvars[34];
2923 double GsfTrackPIn = inputvars[35];
2924 double fbrem = inputvars[36];
2925 double Charge = inputvars[37];
2926 double EoP = inputvars[38];
2927 double TrackMomentumError = inputvars[39];
2928 double EcalEnergyError = inputvars[40];
2929 int Classification = inputvars[41];
2930 double detaIn = inputvars[42];
2931 double dphiIn = inputvars[43];
2932 double detaCalo = inputvars[44];
2933 double dphiCalo = inputvars[45];
2934 double GsfTrackChiSqr = inputvars[46];
2935 double KFTrackNLayers = inputvars[47];
2936 double ElectronEnergyOverPout = inputvars[48];
2937
2938 float *vals = (std::abs(scEta) <= 1.479) ? new float[53] : new float[46];
2939 if (std::abs(scEta) <= 1.479) {
2940 vals[0] = SCRawEnergy;
2941 vals[1] = scEta;
2942 vals[2] = scPhi;
2943 vals[3] = R9;
2944 vals[4] = E5x5Seed / SCRawEnergy;
2945 vals[5] = etawidth;
2946 vals[6] = phiwidth;
2947 vals[7] = NClusters;
2948 vals[8] = HoE;
2949 vals[9] = rho;
2950 vals[10] = vertices;
2951 vals[11] = EtaSeed - scEta;
2952 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
2953 vals[13] = ESeed / SCRawEnergy;
2954 vals[14] = E3x3Seed / ESeed;
2955 vals[15] = E5x5Seed / ESeed;
2956 vals[16] = see;
2957 vals[17] = spp;
2958 vals[18] = sep;
2959 vals[19] = EMaxSeed / ESeed;
2960 vals[20] = E2ndSeed / ESeed;
2961 vals[21] = ETopSeed / ESeed;
2962 vals[22] = EBottomSeed / ESeed;
2963 vals[23] = ELeftSeed / ESeed;
2964 vals[24] = ERightSeed / ESeed;
2965 vals[25] = E2x5MaxSeed / ESeed;
2966 vals[26] = E2x5TopSeed / ESeed;
2967 vals[27] = E2x5BottomSeed / ESeed;
2968 vals[28] = E2x5LeftSeed / ESeed;
2969 vals[29] = E2x5RightSeed / ESeed;
2970 vals[30] = IsEcalDriven;
2971 vals[31] = GsfTrackPIn;
2972 vals[32] = fbrem;
2973 vals[33] = Charge;
2974 vals[34] = EoP;
2975 vals[35] = TrackMomentumError / GsfTrackPIn;
2976 vals[36] = EcalEnergyError / SCRawEnergy;
2977 vals[37] = Classification;
2978 vals[38] = detaIn;
2979 vals[39] = dphiIn;
2980 vals[40] = detaCalo;
2981 vals[41] = dphiCalo;
2982 vals[42] = GsfTrackChiSqr;
2983 vals[43] = KFTrackNLayers;
2984 vals[44] = ElectronEnergyOverPout;
2985 vals[45] = IEtaSeed;
2986 vals[46] = IPhiSeed;
2987 vals[47] = ((int)IEtaSeed) % 5;
2988 vals[48] = ((int)IPhiSeed) % 2;
2989 vals[49] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
2990 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
2991 vals[50] = ((int)IPhiSeed) % 20;
2992 vals[51] = EtaCrySeed;
2993 vals[52] = PhiCrySeed;
2994 }
2995
2996 else {
2997 vals[0] = SCRawEnergy;
2998 vals[1] = scEta;
2999 vals[2] = scPhi;
3000 vals[3] = R9;
3001 vals[4] = E5x5Seed / SCRawEnergy;
3002 vals[5] = etawidth;
3003 vals[6] = phiwidth;
3004 vals[7] = NClusters;
3005 vals[8] = HoE;
3006 vals[9] = rho;
3007 vals[10] = vertices;
3008 vals[11] = EtaSeed - scEta;
3009 vals[12] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
3010 vals[13] = ESeed / SCRawEnergy;
3011 vals[14] = E3x3Seed / ESeed;
3012 vals[15] = E5x5Seed / ESeed;
3013 vals[16] = see;
3014 vals[17] = spp;
3015 vals[18] = sep;
3016 vals[19] = EMaxSeed / ESeed;
3017 vals[20] = E2ndSeed / ESeed;
3018 vals[21] = ETopSeed / ESeed;
3019 vals[22] = EBottomSeed / ESeed;
3020 vals[23] = ELeftSeed / ESeed;
3021 vals[24] = ERightSeed / ESeed;
3022 vals[25] = E2x5MaxSeed / ESeed;
3023 vals[26] = E2x5TopSeed / ESeed;
3024 vals[27] = E2x5BottomSeed / ESeed;
3025 vals[28] = E2x5LeftSeed / ESeed;
3026 vals[29] = E2x5RightSeed / ESeed;
3027 vals[30] = IsEcalDriven;
3028 vals[31] = GsfTrackPIn;
3029 vals[32] = fbrem;
3030 vals[33] = Charge;
3031 vals[34] = EoP;
3032 vals[35] = TrackMomentumError / GsfTrackPIn;
3033 vals[36] = EcalEnergyError / SCRawEnergy;
3034 vals[37] = Classification;
3035 vals[38] = detaIn;
3036 vals[39] = dphiIn;
3037 vals[40] = detaCalo;
3038 vals[41] = dphiCalo;
3039 vals[42] = GsfTrackChiSqr;
3040 vals[43] = KFTrackNLayers;
3041 vals[44] = ElectronEnergyOverPout;
3042 vals[45] = PreShowerOverRaw;
3043 }
3044
3045
3046 double regressionResult = 0;
3047
3048 if (fVersionType == kWithTrkVarV2) {
3049 if (std::abs(scEta) <= 1.479)
3050 regressionResult = SCRawEnergy * forestUncertainty_eb->GetResponse(vals);
3051 else
3052 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestUncertainty_ee->GetResponse(vals);
3053 }
3054
3055
3056 if (printDebug) {
3057 if (std::abs(scEta) <= 1.479) {
3058 std::cout << "Barrel :";
3059 for (unsigned int v = 0; v < 53; ++v)
3060 std::cout << vals[v] << ", ";
3061 std::cout << "\n";
3062 } else {
3063 std::cout << "Endcap :";
3064 for (unsigned int v = 0; v < 46; ++v)
3065 std::cout << vals[v] << ", ";
3066 std::cout << "\n";
3067 }
3068 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
3069 std::cout << "regression energy uncertainty = " << regressionResult << std::endl;
3070 }
3071
3072
3073 delete[] vals;
3074 return regressionResult;
3075 }
3076
3077 double ElectronEnergyRegressionEvaluate::regressionValueWithSubClusters(double SCRawEnergy,
3078 double scEta,
3079 double scPhi,
3080 double R9,
3081 double etawidth,
3082 double phiwidth,
3083 double NClusters,
3084 double HoE,
3085 double rho,
3086 double vertices,
3087 double EtaSeed,
3088 double PhiSeed,
3089 double ESeed,
3090 double E3x3Seed,
3091 double E5x5Seed,
3092 double see,
3093 double spp,
3094 double sep,
3095 double EMaxSeed,
3096 double E2ndSeed,
3097 double ETopSeed,
3098 double EBottomSeed,
3099 double ELeftSeed,
3100 double ERightSeed,
3101 double E2x5MaxSeed,
3102 double E2x5TopSeed,
3103 double E2x5BottomSeed,
3104 double E2x5LeftSeed,
3105 double E2x5RightSeed,
3106 double IEtaSeed,
3107 double IPhiSeed,
3108 double EtaCrySeed,
3109 double PhiCrySeed,
3110 double PreShowerOverRaw,
3111 double isEcalDriven,
3112 double isEtaGap,
3113 double isPhiGap,
3114 double isDeeGap,
3115 double ESubs,
3116 double ESub1,
3117 double EtaSub1,
3118 double PhiSub1,
3119 double EMaxSub1,
3120 double E3x3Sub1,
3121 double ESub2,
3122 double EtaSub2,
3123 double PhiSub2,
3124 double EMaxSub2,
3125 double E3x3Sub2,
3126 double ESub3,
3127 double EtaSub3,
3128 double PhiSub3,
3129 double EMaxSub3,
3130 double E3x3Sub3,
3131 double NPshwClusters,
3132 double EPshwSubs,
3133 double EPshwSub1,
3134 double EtaPshwSub1,
3135 double PhiPshwSub1,
3136 double EPshwSub2,
3137 double EtaPshwSub2,
3138 double PhiPshwSub2,
3139 double EPshwSub3,
3140 double EtaPshwSub3,
3141 double PhiPshwSub3,
3142 bool isEB,
3143 bool printDebug) {
3144
3145 if (fIsInitialized == kFALSE) {
3146 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
3147 return 0;
3148 }
3149
3150
3151 if (!(fVersionType == kWithSubCluVar)) {
3152 std::cout << "Error: Regression VersionType " << fVersionType
3153 << " is not supported to use function regressionValueWithSubClusters.\n";
3154 return 0;
3155 }
3156
3157
3158 float *vals = (isEB) ? new float[61] : new float[65];
3159 if (isEB) {
3160 vals[0] = rho;
3161 vals[1] = vertices;
3162 vals[2] = isEcalDriven;
3163 vals[3] = isEtaGap;
3164 vals[4] = isPhiGap;
3165 vals[5] = isDeeGap;
3166 vals[6] = SCRawEnergy;
3167 vals[7] = scEta;
3168 vals[8] = scPhi;
3169 vals[9] = R9;
3170 vals[10] = etawidth;
3171 vals[11] = phiwidth;
3172 vals[12] = NClusters;
3173 vals[13] = HoE;
3174 vals[14] = EtaSeed - scEta;
3175 vals[15] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
3176 vals[16] = ESeed / SCRawEnergy;
3177 vals[17] = E3x3Seed / ESeed;
3178 vals[18] = E5x5Seed / SCRawEnergy;
3179 vals[19] = E5x5Seed / ESeed;
3180 vals[20] = EMaxSeed / ESeed;
3181 vals[21] = E2ndSeed / ESeed;
3182 vals[22] = ETopSeed / ESeed;
3183 vals[23] = EBottomSeed / ESeed;
3184 vals[24] = ELeftSeed / ESeed;
3185 vals[25] = ERightSeed / ESeed;
3186 vals[26] = E2x5MaxSeed / ESeed;
3187 vals[27] = E2x5TopSeed / ESeed;
3188 vals[28] = E2x5BottomSeed / ESeed;
3189 vals[29] = E2x5LeftSeed / ESeed;
3190 vals[30] = E2x5RightSeed / ESeed;
3191 vals[31] = see;
3192 vals[32] = spp;
3193 vals[33] = sep;
3194 vals[34] = phiwidth / etawidth;
3195 vals[35] = (ELeftSeed + ERightSeed == 0. ? 0. : (ELeftSeed - ERightSeed) / (ELeftSeed + ERightSeed));
3196 vals[36] = (ETopSeed + EBottomSeed == 0. ? 0. : (ETopSeed - EBottomSeed) / (ETopSeed + EBottomSeed));
3197 vals[37] = ESubs / SCRawEnergy;
3198 vals[38] = ESub1 / SCRawEnergy;
3199 vals[39] = (NClusters <= 1 ? 999. : EtaSub1 - EtaSeed);
3200 vals[40] = (NClusters <= 1 ? 999. : atan2(sin(PhiSub1 - PhiSeed), cos(PhiSub1 - PhiSeed)));
3201 vals[41] = (NClusters <= 1 ? 0. : EMaxSub1 / ESub1);
3202 vals[42] = (NClusters <= 1 ? 0. : E3x3Sub1 / ESub1);
3203 vals[43] = ESub2 / SCRawEnergy;
3204 vals[44] = (NClusters <= 2 ? 999. : EtaSub2 - EtaSeed);
3205 vals[45] = (NClusters <= 2 ? 999. : atan2(sin(PhiSub2 - PhiSeed), cos(PhiSub2 - PhiSeed)));
3206 vals[46] = (NClusters <= 2 ? 0. : EMaxSub2 / ESub2);
3207 vals[47] = (NClusters <= 2 ? 0. : E3x3Sub2 / ESub2);
3208 vals[48] = ESub3 / SCRawEnergy;
3209 vals[49] = (NClusters <= 3 ? 999. : EtaSub3 - EtaSeed);
3210 vals[50] = (NClusters <= 3 ? 999. : atan2(sin(PhiSub3 - PhiSeed), cos(PhiSub3 - PhiSeed)));
3211 vals[51] = (NClusters <= 3 ? 0. : EMaxSub3 / ESub3);
3212 vals[52] = (NClusters <= 3 ? 0. : E3x3Sub3 / ESub3);
3213 vals[53] = IEtaSeed;
3214 vals[54] = ((int)IEtaSeed) % 5;
3215 vals[55] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
3216 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
3217 vals[56] = IPhiSeed;
3218 vals[57] = ((int)IPhiSeed) % 2;
3219 vals[58] = ((int)IPhiSeed) % 20;
3220 vals[59] = EtaCrySeed;
3221 vals[60] = PhiCrySeed;
3222 } else {
3223 vals[0] = rho;
3224 vals[1] = vertices;
3225 vals[2] = isEcalDriven;
3226 vals[3] = isEtaGap;
3227 vals[4] = isPhiGap;
3228 vals[5] = isDeeGap;
3229 vals[6] = SCRawEnergy;
3230 vals[7] = scEta;
3231 vals[8] = scPhi;
3232 vals[9] = R9;
3233 vals[10] = etawidth;
3234 vals[11] = phiwidth;
3235 vals[12] = NClusters;
3236 vals[13] = HoE;
3237 vals[14] = EtaSeed - scEta;
3238 vals[15] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
3239 vals[16] = ESeed / SCRawEnergy;
3240 vals[17] = E3x3Seed / ESeed;
3241 vals[18] = E5x5Seed / SCRawEnergy;
3242 vals[19] = E5x5Seed / ESeed;
3243 vals[20] = EMaxSeed / ESeed;
3244 vals[21] = E2ndSeed / ESeed;
3245 vals[22] = ETopSeed / ESeed;
3246 vals[23] = EBottomSeed / ESeed;
3247 vals[24] = ELeftSeed / ESeed;
3248 vals[25] = ERightSeed / ESeed;
3249 vals[26] = E2x5MaxSeed / ESeed;
3250 vals[27] = E2x5TopSeed / ESeed;
3251 vals[28] = E2x5BottomSeed / ESeed;
3252 vals[29] = E2x5LeftSeed / ESeed;
3253 vals[30] = E2x5RightSeed / ESeed;
3254 vals[31] = see;
3255 vals[32] = spp;
3256 vals[33] = sep;
3257 vals[34] = phiwidth / etawidth;
3258 vals[35] = (ELeftSeed + ERightSeed == 0. ? 0. : (ELeftSeed - ERightSeed) / (ELeftSeed + ERightSeed));
3259 vals[36] = (ETopSeed + EBottomSeed == 0. ? 0. : (ETopSeed - EBottomSeed) / (ETopSeed + EBottomSeed));
3260 vals[37] = ESubs / SCRawEnergy;
3261 vals[38] = ESub1 / SCRawEnergy;
3262 vals[39] = (NClusters <= 1 ? 999. : EtaSub1 - EtaSeed);
3263 vals[40] = (NClusters <= 1 ? 999. : atan2(sin(PhiSub1 - PhiSeed), cos(PhiSub1 - PhiSeed)));
3264 vals[41] = (NClusters <= 1 ? 0. : EMaxSub1 / ESub1);
3265 vals[42] = (NClusters <= 1 ? 0. : E3x3Sub1 / ESub1);
3266 vals[43] = ESub2 / SCRawEnergy;
3267 vals[44] = (NClusters <= 2 ? 999. : EtaSub2 - EtaSeed);
3268 vals[45] = (NClusters <= 2 ? 999. : atan2(sin(PhiSub2 - PhiSeed), cos(PhiSub2 - PhiSeed)));
3269 vals[46] = (NClusters <= 2 ? 0. : EMaxSub2 / ESub2);
3270 vals[47] = (NClusters <= 2 ? 0. : E3x3Sub2 / ESub2);
3271 vals[48] = ESub3 / SCRawEnergy;
3272 vals[49] = (NClusters <= 3 ? 999. : EtaSub3 - EtaSeed);
3273 vals[50] = (NClusters <= 3 ? 999. : atan2(sin(PhiSub3 - PhiSeed), cos(PhiSub3 - PhiSeed)));
3274 vals[51] = (NClusters <= 3 ? 0. : EMaxSub3 / ESub3);
3275 vals[52] = (NClusters <= 3 ? 0. : E3x3Sub3 / ESub3);
3276 vals[53] = PreShowerOverRaw;
3277 vals[54] = NPshwClusters;
3278 vals[55] = EPshwSubs / SCRawEnergy;
3279 vals[56] = EPshwSub1 / SCRawEnergy;
3280 vals[57] = (NPshwClusters == 0 ? 999. : EtaPshwSub1 - EtaSeed);
3281 vals[58] = (NPshwClusters == 0 ? 999. : atan2(sin(PhiPshwSub1 - PhiSeed), cos(PhiPshwSub1 - PhiSeed)));
3282 vals[59] = EPshwSub2 / SCRawEnergy;
3283 vals[60] = (NPshwClusters <= 1 ? 999. : EtaPshwSub2 - EtaSeed);
3284 vals[61] = (NPshwClusters <= 1 ? 999. : atan2(sin(PhiPshwSub2 - PhiSeed), cos(PhiPshwSub2 - PhiSeed)));
3285 vals[62] = EPshwSub3 / SCRawEnergy;
3286 vals[63] = (NPshwClusters <= 2 ? 999. : EtaPshwSub3 - EtaSeed);
3287 vals[64] = (NPshwClusters <= 2 ? 999. : atan2(sin(PhiPshwSub3 - PhiSeed), cos(PhiPshwSub3 - PhiSeed)));
3288 }
3289
3290
3291 double regressionResult = 0;
3292 Int_t BinIndex = -1;
3293
3294 if (fVersionType == kWithSubCluVar) {
3295 if (isEB) {
3296 regressionResult = SCRawEnergy * forestCorrection_eb->GetResponse(vals);
3297 BinIndex = 0;
3298 } else {
3299 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestCorrection_ee->GetResponse(vals);
3300 BinIndex = 1;
3301 }
3302 }
3303
3304
3305 if (printDebug) {
3306 if (isEB) {
3307 std::cout << "Barrel :";
3308 for (unsigned int v = 0; v < 61; ++v)
3309 std::cout << v << "=" << vals[v] << ", ";
3310 std::cout << "\n";
3311 } else {
3312 std::cout << "Endcap :";
3313 for (unsigned int v = 0; v < 65; ++v)
3314 std::cout << v << "=" << vals[v] << ", ";
3315 std::cout << "\n";
3316 }
3317 std::cout << "BinIndex : " << BinIndex << "\n";
3318 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
3319 std::cout << "regression energy = " << regressionResult << std::endl;
3320 }
3321
3322
3323 delete[] vals;
3324 return regressionResult;
3325 }
3326
3327 double ElectronEnergyRegressionEvaluate::regressionUncertaintyWithSubClusters(double SCRawEnergy,
3328 double scEta,
3329 double scPhi,
3330 double R9,
3331 double etawidth,
3332 double phiwidth,
3333 double NClusters,
3334 double HoE,
3335 double rho,
3336 double vertices,
3337 double EtaSeed,
3338 double PhiSeed,
3339 double ESeed,
3340 double E3x3Seed,
3341 double E5x5Seed,
3342 double see,
3343 double spp,
3344 double sep,
3345 double EMaxSeed,
3346 double E2ndSeed,
3347 double ETopSeed,
3348 double EBottomSeed,
3349 double ELeftSeed,
3350 double ERightSeed,
3351 double E2x5MaxSeed,
3352 double E2x5TopSeed,
3353 double E2x5BottomSeed,
3354 double E2x5LeftSeed,
3355 double E2x5RightSeed,
3356 double IEtaSeed,
3357 double IPhiSeed,
3358 double EtaCrySeed,
3359 double PhiCrySeed,
3360 double PreShowerOverRaw,
3361 double isEcalDriven,
3362 double isEtaGap,
3363 double isPhiGap,
3364 double isDeeGap,
3365 double ESubs,
3366 double ESub1,
3367 double EtaSub1,
3368 double PhiSub1,
3369 double EMaxSub1,
3370 double E3x3Sub1,
3371 double ESub2,
3372 double EtaSub2,
3373 double PhiSub2,
3374 double EMaxSub2,
3375 double E3x3Sub2,
3376 double ESub3,
3377 double EtaSub3,
3378 double PhiSub3,
3379 double EMaxSub3,
3380 double E3x3Sub3,
3381 double NPshwClusters,
3382 double EPshwSubs,
3383 double EPshwSub1,
3384 double EtaPshwSub1,
3385 double PhiPshwSub1,
3386 double EPshwSub2,
3387 double EtaPshwSub2,
3388 double PhiPshwSub2,
3389 double EPshwSub3,
3390 double EtaPshwSub3,
3391 double PhiPshwSub3,
3392 bool isEB,
3393 bool printDebug) {
3394
3395 if (fIsInitialized == kFALSE) {
3396 printf("ElectronEnergyRegressionEvaluate instance not initialized !!!");
3397 return 0;
3398 }
3399
3400
3401 if (!(fVersionType == kWithSubCluVar)) {
3402 std::cout << "Error: Regression VersionType " << fVersionType
3403 << " is not supported to use function regressionValueWithSubClusters.\n";
3404 return 0;
3405 }
3406
3407
3408 float *vals = (isEB) ? new float[61] : new float[65];
3409 if (isEB) {
3410 vals[0] = rho;
3411 vals[1] = vertices;
3412 vals[2] = isEcalDriven;
3413 vals[3] = isEtaGap;
3414 vals[4] = isPhiGap;
3415 vals[5] = isDeeGap;
3416 vals[6] = SCRawEnergy;
3417 vals[7] = scEta;
3418 vals[8] = scPhi;
3419 vals[9] = R9;
3420 vals[10] = etawidth;
3421 vals[11] = phiwidth;
3422 vals[12] = NClusters;
3423 vals[13] = HoE;
3424 vals[14] = EtaSeed - scEta;
3425 vals[15] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
3426 vals[16] = ESeed / SCRawEnergy;
3427 vals[17] = E3x3Seed / ESeed;
3428 vals[18] = E5x5Seed / SCRawEnergy;
3429 vals[19] = E5x5Seed / ESeed;
3430 vals[20] = EMaxSeed / ESeed;
3431 vals[21] = E2ndSeed / ESeed;
3432 vals[22] = ETopSeed / ESeed;
3433 vals[23] = EBottomSeed / ESeed;
3434 vals[24] = ELeftSeed / ESeed;
3435 vals[25] = ERightSeed / ESeed;
3436 vals[26] = E2x5MaxSeed / ESeed;
3437 vals[27] = E2x5TopSeed / ESeed;
3438 vals[28] = E2x5BottomSeed / ESeed;
3439 vals[29] = E2x5LeftSeed / ESeed;
3440 vals[30] = E2x5RightSeed / ESeed;
3441 vals[31] = see;
3442 vals[32] = spp;
3443 vals[33] = sep;
3444 vals[34] = phiwidth / etawidth;
3445 vals[35] = (ELeftSeed + ERightSeed == 0. ? 0. : (ELeftSeed - ERightSeed) / (ELeftSeed + ERightSeed));
3446 vals[36] = (ETopSeed + EBottomSeed == 0. ? 0. : (ETopSeed - EBottomSeed) / (ETopSeed + EBottomSeed));
3447 vals[37] = ESubs / SCRawEnergy;
3448 vals[38] = ESub1 / SCRawEnergy;
3449 vals[39] = (NClusters <= 1 ? 999. : EtaSub1 - EtaSeed);
3450 vals[40] = (NClusters <= 1 ? 999. : atan2(sin(PhiSub1 - PhiSeed), cos(PhiSub1 - PhiSeed)));
3451 vals[41] = (NClusters <= 1 ? 0. : EMaxSub1 / ESub1);
3452 vals[42] = (NClusters <= 1 ? 0. : E3x3Sub1 / ESub1);
3453 vals[43] = ESub2 / SCRawEnergy;
3454 vals[44] = (NClusters <= 2 ? 999. : EtaSub2 - EtaSeed);
3455 vals[45] = (NClusters <= 2 ? 999. : atan2(sin(PhiSub2 - PhiSeed), cos(PhiSub2 - PhiSeed)));
3456 vals[46] = (NClusters <= 2 ? 0. : EMaxSub2 / ESub2);
3457 vals[47] = (NClusters <= 2 ? 0. : E3x3Sub2 / ESub2);
3458 vals[48] = ESub3 / SCRawEnergy;
3459 vals[49] = (NClusters <= 3 ? 999. : EtaSub3 - EtaSeed);
3460 vals[50] = (NClusters <= 3 ? 999. : atan2(sin(PhiSub3 - PhiSeed), cos(PhiSub3 - PhiSeed)));
3461 vals[51] = (NClusters <= 3 ? 0. : EMaxSub3 / ESub3);
3462 vals[52] = (NClusters <= 3 ? 0. : E3x3Sub3 / ESub3);
3463 vals[53] = IEtaSeed;
3464 vals[54] = ((int)IEtaSeed) % 5;
3465 vals[55] = (std::abs(IEtaSeed) <= 25) * (((int)IEtaSeed) % 25) +
3466 (std::abs(IEtaSeed) > 25) * (((int)(IEtaSeed - 25 * std::abs(IEtaSeed) / IEtaSeed)) % 20);
3467 vals[56] = IPhiSeed;
3468 vals[57] = ((int)IPhiSeed) % 2;
3469 vals[58] = ((int)IPhiSeed) % 20;
3470 vals[59] = EtaCrySeed;
3471 vals[60] = PhiCrySeed;
3472 } else {
3473 vals[0] = rho;
3474 vals[1] = vertices;
3475 vals[2] = isEcalDriven;
3476 vals[3] = isEtaGap;
3477 vals[4] = isPhiGap;
3478 vals[5] = isDeeGap;
3479 vals[6] = SCRawEnergy;
3480 vals[7] = scEta;
3481 vals[8] = scPhi;
3482 vals[9] = R9;
3483 vals[10] = etawidth;
3484 vals[11] = phiwidth;
3485 vals[12] = NClusters;
3486 vals[13] = HoE;
3487 vals[14] = EtaSeed - scEta;
3488 vals[15] = atan2(sin(PhiSeed - scPhi), cos(PhiSeed - scPhi));
3489 vals[16] = ESeed / SCRawEnergy;
3490 vals[17] = E3x3Seed / ESeed;
3491 vals[18] = E5x5Seed / SCRawEnergy;
3492 vals[19] = E5x5Seed / ESeed;
3493 vals[20] = EMaxSeed / ESeed;
3494 vals[21] = E2ndSeed / ESeed;
3495 vals[22] = ETopSeed / ESeed;
3496 vals[23] = EBottomSeed / ESeed;
3497 vals[24] = ELeftSeed / ESeed;
3498 vals[25] = ERightSeed / ESeed;
3499 vals[26] = E2x5MaxSeed / ESeed;
3500 vals[27] = E2x5TopSeed / ESeed;
3501 vals[28] = E2x5BottomSeed / ESeed;
3502 vals[29] = E2x5LeftSeed / ESeed;
3503 vals[30] = E2x5RightSeed / ESeed;
3504 vals[31] = see;
3505 vals[32] = spp;
3506 vals[33] = sep;
3507 vals[34] = phiwidth / etawidth;
3508 vals[35] = (ELeftSeed + ERightSeed == 0. ? 0. : (ELeftSeed - ERightSeed) / (ELeftSeed + ERightSeed));
3509 vals[36] = (ETopSeed + EBottomSeed == 0. ? 0. : (ETopSeed - EBottomSeed) / (ETopSeed + EBottomSeed));
3510 vals[37] = ESubs / SCRawEnergy;
3511 vals[38] = ESub1 / SCRawEnergy;
3512 vals[39] = (NClusters <= 1 ? 999. : EtaSub1 - EtaSeed);
3513 vals[40] = (NClusters <= 1 ? 999. : atan2(sin(PhiSub1 - PhiSeed), cos(PhiSub1 - PhiSeed)));
3514 vals[41] = (NClusters <= 1 ? 0. : EMaxSub1 / ESub1);
3515 vals[42] = (NClusters <= 1 ? 0. : E3x3Sub1 / ESub1);
3516 vals[43] = ESub2 / SCRawEnergy;
3517 vals[44] = (NClusters <= 2 ? 999. : EtaSub2 - EtaSeed);
3518 vals[45] = (NClusters <= 2 ? 999. : atan2(sin(PhiSub2 - PhiSeed), cos(PhiSub2 - PhiSeed)));
3519 vals[46] = (NClusters <= 2 ? 0. : EMaxSub2 / ESub2);
3520 vals[47] = (NClusters <= 2 ? 0. : E3x3Sub2 / ESub2);
3521 vals[48] = ESub3 / SCRawEnergy;
3522 vals[49] = (NClusters <= 3 ? 999. : EtaSub3 - EtaSeed);
3523 vals[50] = (NClusters <= 3 ? 999. : atan2(sin(PhiSub3 - PhiSeed), cos(PhiSub3 - PhiSeed)));
3524 vals[51] = (NClusters <= 3 ? 0. : EMaxSub3 / ESub3);
3525 vals[52] = (NClusters <= 3 ? 0. : E3x3Sub3 / ESub3);
3526 vals[53] = PreShowerOverRaw;
3527 vals[54] = NPshwClusters;
3528 vals[55] = EPshwSubs / SCRawEnergy;
3529 vals[56] = EPshwSub1 / SCRawEnergy;
3530 vals[57] = (NPshwClusters <= 0 ? 999. : EtaPshwSub1 - EtaSeed);
3531 vals[58] = (NPshwClusters <= 0 ? 999. : atan2(sin(PhiPshwSub1 - PhiSeed), cos(PhiPshwSub1 - PhiSeed)));
3532 vals[59] = EPshwSub2 / SCRawEnergy;
3533 vals[60] = (NPshwClusters <= 1 ? 999. : EtaPshwSub2 - EtaSeed);
3534 vals[61] = (NPshwClusters <= 1 ? 999. : atan2(sin(PhiPshwSub2 - PhiSeed), cos(PhiPshwSub2 - PhiSeed)));
3535 vals[62] = EPshwSub3 / SCRawEnergy;
3536 vals[63] = (NPshwClusters <= 2 ? 999. : EtaPshwSub3 - EtaSeed);
3537 vals[64] = (NPshwClusters <= 2 ? 999. : atan2(sin(PhiPshwSub3 - PhiSeed), cos(PhiPshwSub3 - PhiSeed)));
3538 }
3539
3540
3541 double regressionResult = 0;
3542 Int_t BinIndex = -1;
3543
3544 if (fVersionType == kWithSubCluVar) {
3545 if (isEB) {
3546 regressionResult = SCRawEnergy * forestUncertainty_eb->GetResponse(vals);
3547 BinIndex = 0;
3548 } else {
3549 regressionResult = (SCRawEnergy * (1 + PreShowerOverRaw)) * forestUncertainty_ee->GetResponse(vals);
3550 BinIndex = 1;
3551 }
3552 }
3553
3554
3555 if (printDebug) {
3556 if (isEB) {
3557 std::cout << "Barrel :";
3558 for (unsigned int v = 0; v < 38; ++v)
3559 std::cout << vals[v] << ", ";
3560 std::cout << "\n";
3561 } else {
3562 std::cout << "Endcap :";
3563 for (unsigned int v = 0; v < 31; ++v)
3564 std::cout << vals[v] << ", ";
3565 std::cout << "\n";
3566 }
3567 std::cout << "BinIndex : " << BinIndex << "\n";
3568 std::cout << "SCRawEnergy = " << SCRawEnergy << " : PreShowerOverRaw = " << PreShowerOverRaw << std::endl;
3569 std::cout << "regression energy uncertainty = " << regressionResult << std::endl;
3570 }
3571
3572
3573 delete[] vals;
3574 return regressionResult;
3575 }
