File indexing completed on 2024-04-06 11:58:40
0001 #include <memory>
0002 #include <cmath>
0003 #include "FWCore/Framework/interface/Frameworkfwd.h"
0004 #include "FWCore/Framework/interface/MakerMacros.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "Calibration/EcalCalibAlgos/interface/IMACalibBlock.h"
0007 #include "Calibration/EcalCalibAlgos/interface/L3CalibBlock.h"
0008 #include "Calibration/EcalCalibAlgos/interface/EcalEleCalibLooper.h"
0009 #include "DataFormats/Common/interface/Handle.h"
0010 #include "DataFormats/EgammaReco/interface/ClusterShape.h"
0011 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0012 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0013 #include "Calibration/Tools/interface/calibXMLwriter.h"
0014 #include "CalibCalorimetry/CaloMiscalibTools/interface/CaloMiscalibTools.h"
0015 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterTools.h"
0016 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
0017 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0018 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
0019 #include "DataFormats/EgammaReco/interface/ClusterShapeFwd.h"
0020 #include "Calibration/EcalCalibAlgos/interface/MatrixFillMap.h"
0021 #include "Calibration/EcalCalibAlgos/interface/ClusterFillMap.h"
0022 #include "TH1.h"
0023 #include "TH2.h"
0024 #include "TFile.h"
0025
0026
0027
0028
0029 EcalEleCalibLooper::EcalEleCalibLooper(const edm::ParameterSet& iConfig)
0030 : m_barrelAlCa(iConfig.getParameter<edm::InputTag>("alcaBarrelHitCollection")),
0031 m_endcapAlCa(iConfig.getParameter<edm::InputTag>("alcaEndcapHitCollection")),
0032 m_ElectronLabel(iConfig.getParameter<edm::InputTag>("electronLabel")),
0033 m_recoWindowSidex(iConfig.getParameter<int>("recoWindowSidex")),
0034 m_recoWindowSidey(iConfig.getParameter<int>("recoWindowSidey")),
0035 m_etaWidth(iConfig.getParameter<int>("etaWidth")),
0036 m_phiWidthEB(iConfig.getParameter<int>("phiWidthEB")),
0037 m_etaStart(etaShifter(iConfig.getParameter<int>("etaStart"))),
0038 m_etaEnd(etaShifter(iConfig.getParameter<int>("etaEnd"))),
0039 m_phiStartEB(iConfig.getParameter<int>("phiStartEB")),
0040 m_phiEndEB(iConfig.getParameter<int>("phiEndEB")),
0041 m_radStart(iConfig.getParameter<int>("radStart")),
0042 m_radEnd(iConfig.getParameter<int>("radEnd")),
0043 m_radWidth(iConfig.getParameter<int>("radWidth")),
0044 m_phiStartEE(iConfig.getParameter<int>("phiStartEE")),
0045 m_phiEndEE(iConfig.getParameter<int>("phiEndEE")),
0046 m_phiWidthEE(iConfig.getParameter<int>("phiWidthEE")),
0047 m_maxSelectedNumPerXtal(iConfig.getParameter<int>("maxSelectedNumPerCrystal")),
0048 m_minEnergyPerCrystal(iConfig.getParameter<double>("minEnergyPerCrystal")),
0049 m_maxEnergyPerCrystal(iConfig.getParameter<double>("maxEnergyPerCrystal")),
0050 m_minCoeff(iConfig.getParameter<double>("minCoeff")),
0051 m_maxCoeff(iConfig.getParameter<double>("maxCoeff")),
0052 m_usingBlockSolver(iConfig.getParameter<int>("usingBlockSolver")),
0053 m_loops(iConfig.getParameter<int>("loops")),
0054 m_ebRecHitToken(consumes<EBRecHitCollection>(m_barrelAlCa)),
0055 m_eeRecHitToken(consumes<EERecHitCollection>(m_endcapAlCa)),
0056 m_gsfElectronToken(consumes<reco::GsfElectronCollection>(m_ElectronLabel)),
0057 m_geometryToken(esConsumes()) {
0058 edm::LogInfo("IML") << "[EcalEleCalibLooper][ctor] asserts";
0059 assert(!((m_etaEnd - m_etaStart) % m_etaWidth));
0060
0061 assert(m_etaStart >= 0 && m_etaStart <= 171);
0062 assert(m_etaEnd >= m_etaStart && m_etaEnd <= 171);
0063 assert((m_radEnd - m_radStart) % m_radWidth == 0);
0064 assert(m_radStart >= 0 && m_radStart <= 50);
0065 assert(m_radEnd >= m_radStart && m_radEnd <= 50);
0066 edm::LogInfo("IML") << "[EcalEleCalibLooper][ctor] entering ";
0067 edm::LogInfo("IML") << "[EcalEleCalibLooper][ctor] region definition";
0068 EBRegionDefinition();
0069 EERegionDefinition();
0070
0071 TH2F* EBRegion = new TH2F("EBRegion", "EBRegion", 170, 0, 170, 360, 0, 360);
0072 for (int eta = 0; eta < 170; ++eta)
0073 for (int phi = 0; phi < 360; ++phi) {
0074 EBRegion->Fill(eta, phi, m_xtalRegionId[EBDetId::unhashIndex(eta * 360 + phi).rawId()]);
0075 }
0076 TH2F* EERegion = new TH2F("EERegion", "EERegion", 100, 0, 100, 100, 0, 100);
0077 for (int x = 0; x < 100; ++x)
0078 for (int y = 0; y < 100; ++y) {
0079 if (EEDetId::validDetId(x + 1, y + 1, 1))
0080 EERegion->Fill(x, y, m_xtalRegionId[EEDetId(x + 1, y + 1, -1).rawId()]);
0081 }
0082
0083 TFile out("EBZone.root", "recreate");
0084 EBRegion->Write();
0085 EERegion->Write();
0086 out.Close();
0087 delete EERegion;
0088 delete EBRegion;
0089
0090
0091
0092
0093
0094 edm::LogInfo("IML") << "[EcalEleCalibLooper][ctor] Calib Block";
0095 std::string algorithm = iConfig.getParameter<std::string>("algorithm");
0096 int eventWeight = iConfig.getUntrackedParameter<int>("L3EventWeight", 1);
0097
0098
0099 for (int region = 0; region < EBregionsNum() + 2 * EEregionsNum(); ++region) {
0100 if (algorithm == "IMA")
0101 m_EcalCalibBlocks.push_back(new IMACalibBlock(m_regions.at(region)));
0102 else if (algorithm == "L3")
0103 m_EcalCalibBlocks.push_back(new L3CalibBlock(m_regions.at(region), eventWeight));
0104 else {
0105 edm::LogError("building") << algorithm << " is not a valid calibration algorithm";
0106 exit(1);
0107 }
0108 }
0109 std::string mapFiller = iConfig.getParameter<std::string>("FillType");
0110 if (mapFiller == "Cluster")
0111 m_MapFiller = new ClusterFillMap(m_recoWindowSidex,
0112 m_recoWindowSidey,
0113 m_xtalRegionId,
0114 m_minEnergyPerCrystal,
0115 m_maxEnergyPerCrystal,
0116 m_xtalPositionInRegion,
0117 &m_barrelMap,
0118 &m_endcapMap);
0119 if (mapFiller == "Matrix")
0120 m_MapFiller = new MatrixFillMap(m_recoWindowSidex,
0121 m_recoWindowSidey,
0122 m_xtalRegionId,
0123 m_minEnergyPerCrystal,
0124 m_maxEnergyPerCrystal,
0125 m_xtalPositionInRegion,
0126 &m_barrelMap,
0127 &m_endcapMap);
0128 }
0129
0130
0131
0132
0133 EcalEleCalibLooper::~EcalEleCalibLooper() {
0134 edm::LogInfo("IML") << "[EcalEleCalibLooper][dtor]";
0135 for (std::vector<VEcalCalibBlock*>::iterator calibBlock = m_EcalCalibBlocks.begin();
0136 calibBlock != m_EcalCalibBlocks.end();
0137 ++calibBlock)
0138 delete (*calibBlock);
0139 }
0140
0141
0142
0143
0144 void EcalEleCalibLooper::beginOfJob() { isfirstcall_ = true; }
0145
0146
0147
0148
0149
0150 void EcalEleCalibLooper::startingNewLoop(unsigned int ciclo) {
0151 edm::LogInfo("IML") << "[InvMatrixCalibLooper][Start] entering loop " << ciclo;
0152
0153 for (std::vector<VEcalCalibBlock*>::iterator calibBlock = m_EcalCalibBlocks.begin();
0154 calibBlock != m_EcalCalibBlocks.end();
0155 ++calibBlock)
0156 (*calibBlock)->reset();
0157 for (std::map<int, int>::iterator it = m_xtalNumOfHits.begin(); it != m_xtalNumOfHits.end(); ++it)
0158 it->second = 0;
0159 return;
0160 }
0161
0162
0163
0164
0165
0166 edm::EDLooper::Status EcalEleCalibLooper::duringLoop(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0167
0168
0169
0170 if (isfirstcall_) {
0171 const auto& geometry = iSetup.getData(m_geometryToken);
0172 m_barrelCells = geometry.getValidDetIds(DetId::Ecal, EcalBarrel);
0173 m_endcapCells = geometry.getValidDetIds(DetId::Ecal, EcalEndcap);
0174 for (std::vector<DetId>::const_iterator barrelIt = m_barrelCells.begin(); barrelIt != m_barrelCells.end();
0175 ++barrelIt) {
0176 m_barrelMap[*barrelIt] = 1;
0177 m_xtalNumOfHits[barrelIt->rawId()] = 0;
0178 }
0179 for (std::vector<DetId>::const_iterator endcapIt = m_endcapCells.begin(); endcapIt != m_endcapCells.end();
0180 ++endcapIt) {
0181 m_endcapMap[*endcapIt] = 1;
0182 m_xtalNumOfHits[endcapIt->rawId()] = 0;
0183 }
0184
0185 isfirstcall_ = false;
0186 }
0187
0188
0189 const EBRecHitCollection* barrelHitsCollection = nullptr;
0190 edm::Handle<EBRecHitCollection> barrelRecHitsHandle;
0191 iEvent.getByToken(m_ebRecHitToken, barrelRecHitsHandle);
0192 barrelHitsCollection = barrelRecHitsHandle.product();
0193 if (!barrelRecHitsHandle.isValid()) {
0194 edm::LogError("reading") << "[EcalEleCalibLooper] barrel rec hits not found";
0195 return kContinue;
0196 }
0197
0198
0199 const EERecHitCollection* endcapHitsCollection = nullptr;
0200 edm::Handle<EERecHitCollection> endcapRecHitsHandle;
0201 iEvent.getByToken(m_eeRecHitToken, endcapRecHitsHandle);
0202 endcapHitsCollection = endcapRecHitsHandle.product();
0203 if (!endcapRecHitsHandle.isValid()) {
0204 edm::LogError("reading") << "[EcalEleCalibLooper] endcap rec hits not found";
0205 return kContinue;
0206 }
0207
0208
0209 edm::Handle<reco::GsfElectronCollection> pElectrons;
0210 iEvent.getByToken(m_gsfElectronToken, pElectrons);
0211 if (!pElectrons.isValid()) {
0212 edm::LogError("reading") << "[EcalEleCalibLooper] electrons not found";
0213 return kContinue;
0214 }
0215
0216
0217 for (reco::GsfElectronCollection::const_iterator eleIt = pElectrons->begin(); eleIt != pElectrons->end(); ++eleIt) {
0218 double pSubtract = 0;
0219 double pTk = 0;
0220 std::map<int, double> xtlMap;
0221 DetId Max = 0;
0222 if (std::abs(eleIt->eta()) < 1.49)
0223 Max = EcalClusterTools::getMaximum(eleIt->superCluster()->hitsAndFractions(), barrelHitsCollection).first;
0224 else
0225 Max = EcalClusterTools::getMaximum(eleIt->superCluster()->hitsAndFractions(), endcapHitsCollection).first;
0226 if (Max.det() == 0)
0227 continue;
0228 m_MapFiller->fillMap(
0229 eleIt->superCluster()->hitsAndFractions(), Max, barrelHitsCollection, endcapHitsCollection, xtlMap, pSubtract);
0230 if (m_maxSelectedNumPerXtal > 0 && m_xtalNumOfHits[Max.rawId()] > m_maxSelectedNumPerXtal)
0231 continue;
0232 ++m_xtalNumOfHits[Max.rawId()];
0233 if (m_xtalRegionId[Max.rawId()] == -1)
0234 continue;
0235 pTk = eleIt->trackMomentumAtVtx().R();
0236 m_EcalCalibBlocks.at(m_xtalRegionId[Max.rawId()])->Fill(xtlMap.begin(), xtlMap.end(), pTk, pSubtract);
0237 }
0238
0239 return kContinue;
0240 }
0241
0242
0243
0244
0245
0246
0247 edm::EDLooper::Status EcalEleCalibLooper::endOfLoop(const edm::EventSetup& dumb, unsigned int iCounter) {
0248 edm::LogInfo("IML") << "[InvMatrixCalibLooper][endOfLoop] entering...";
0249 for (std::vector<VEcalCalibBlock*>::iterator calibBlock = m_EcalCalibBlocks.begin();
0250 calibBlock != m_EcalCalibBlocks.end();
0251 ++calibBlock)
0252 (*calibBlock)->solve(m_usingBlockSolver, m_minCoeff, m_maxCoeff);
0253
0254 TH1F* EBcoeffEnd = new TH1F("EBRegion", "EBRegion", 100, 0.5, 2.1);
0255 TH2F* EBcoeffMap = new TH2F("EBcoeff", "EBcoeff", 171, -85, 85, 360, 1, 361);
0256 TH1F* EEPcoeffEnd = new TH1F("EEPRegion", "EEPRegion", 100, 0.5, 2.1);
0257 TH1F* EEMcoeffEnd = new TH1F("EEMRegion", "EEMRegion", 100, 0.5, 2.1);
0258 TH2F* EEPcoeffMap = new TH2F("EEPcoeffMap", "EEPcoeffMap", 101, 1, 101, 101, 0, 101);
0259 TH2F* EEMcoeffMap = new TH2F("EEMcoeffMap", "EEMcoeffMap", 101, 1, 101, 101, 0, 101);
0260
0261 for (std::vector<DetId>::const_iterator barrelIt = m_barrelCells.begin(); barrelIt != m_barrelCells.end();
0262 ++barrelIt) {
0263 EBDetId ee(*barrelIt);
0264 int index = barrelIt->rawId();
0265 if (m_xtalRegionId[index] == -1)
0266 continue;
0267 m_barrelMap[*barrelIt] *= m_EcalCalibBlocks.at(m_xtalRegionId[index])->at(m_xtalPositionInRegion[index]);
0268 EBcoeffEnd->Fill(m_barrelMap[*barrelIt]);
0269 EBcoeffMap->Fill(ee.ieta(), ee.iphi(), m_barrelMap[*barrelIt]);
0270 }
0271
0272
0273 for (std::vector<DetId>::const_iterator endcapIt = m_endcapCells.begin(); endcapIt != m_endcapCells.end();
0274 ++endcapIt) {
0275 EEDetId ee(*endcapIt);
0276 int index = endcapIt->rawId();
0277 if (ee.zside() > 0) {
0278 if (m_xtalRegionId[index] == -1)
0279 continue;
0280 m_endcapMap[*endcapIt] *= m_EcalCalibBlocks.at(m_xtalRegionId[index])->at(m_xtalPositionInRegion[index]);
0281 EEPcoeffEnd->Fill(m_endcapMap[*endcapIt]);
0282 EEPcoeffMap->Fill(ee.ix(), ee.iy(), m_endcapMap[*endcapIt]);
0283 } else {
0284 m_endcapMap[*endcapIt] *= m_EcalCalibBlocks.at(m_xtalRegionId[index])->at(m_xtalPositionInRegion[index]);
0285 EEMcoeffEnd->Fill(m_endcapMap[*endcapIt]);
0286 EEMcoeffMap->Fill(ee.ix(), ee.iy(), m_endcapMap[*endcapIt]);
0287 }
0288 }
0289
0290 edm::LogInfo("IML") << "[InvMatrixCalibLooper][endOfLoop] End of endOfLoop";
0291
0292 char filename[80];
0293 sprintf(filename, "coeffs%d.root", iCounter);
0294 TFile zout(filename, "recreate");
0295 EBcoeffEnd->Write();
0296 EBcoeffMap->Write();
0297 EEPcoeffEnd->Write();
0298 EEPcoeffMap->Write();
0299 EEMcoeffEnd->Write();
0300 EEMcoeffMap->Write();
0301 zout.Close();
0302 delete EBcoeffEnd;
0303 delete EBcoeffMap;
0304 delete EEPcoeffEnd;
0305 delete EEMcoeffEnd;
0306 delete EEPcoeffMap;
0307 delete EEMcoeffMap;
0308 if (iCounter < m_loops - 1)
0309 return kContinue;
0310 else
0311 return kStop;
0312 }
0313
0314
0315
0316
0317
0318 void EcalEleCalibLooper::endOfJob() {
0319 edm::LogInfo("IML") << "[InvMatrixCalibLooper][endOfJob] saving calib coeffs";
0320
0321
0322 calibXMLwriter barrelWriter(EcalBarrel);
0323 calibXMLwriter endcapWriter(EcalEndcap);
0324 for (std::vector<DetId>::const_iterator barrelIt = m_barrelCells.begin(); barrelIt != m_barrelCells.end();
0325 ++barrelIt) {
0326 EBDetId eb(*barrelIt);
0327 barrelWriter.writeLine(eb, m_barrelMap[*barrelIt]);
0328 }
0329 for (std::vector<DetId>::const_iterator endcapIt = m_endcapCells.begin(); endcapIt != m_endcapCells.end();
0330 ++endcapIt) {
0331 EEDetId ee(*endcapIt);
0332 endcapWriter.writeLine(ee, m_endcapMap[*endcapIt]);
0333 }
0334 edm::LogInfo("IML") << "[InvMatrixCalibLooper][endOfJob] Exiting";
0335 }
0336
0337
0338
0339
0340
0341
0342 int EcalEleCalibLooper::EBregionCheck(const int eta, const int phi) const {
0343 if (eta < m_etaStart)
0344 return 1;
0345 if (eta >= m_etaEnd)
0346 return 2;
0347 if (phi < m_phiStartEB)
0348 return 3;
0349 if (phi >= m_phiEndEB)
0350 return 4;
0351 return 0;
0352 }
0353
0354
0355
0356
0357 inline double degrees(double radiants) { return radiants * 180 * (1 / M_PI); }
0358
0359
0360
0361
0362 inline double radiants(int degrees) { return degrees * M_PI * (1 / 180); }
0363
0364
0365
0366
0367 double EcalEleCalibLooper::giveLimit(int degrees) {
0368
0369 if (degrees == 90)
0370 return 90;
0371 return tan(radiants(degrees));
0372 }
0373
0374
0375
0376
0377 inline double Mod(double phi) {
0378 if (phi >= 360 && phi < 720)
0379 return phi - 360;
0380 if (phi >= 720)
0381 return phi - 720;
0382 return phi;
0383 }
0384
0385
0386
0387
0388 int EcalEleCalibLooper::EBRegionId(const int etaXtl, const int phiXtl) const {
0389 if (EBregionCheck(etaXtl, phiXtl))
0390 return -1;
0391 int phifake = m_phiStartEB;
0392 if (m_phiStartEB > m_phiEndEB)
0393 phifake = m_phiStartEB - 360;
0394 int Nphi = (m_phiEndEB - phifake) / m_phiWidthEB;
0395 int etaI = (etaXtl - m_etaStart) / m_etaWidth;
0396 int phiI = (phiXtl - m_phiStartEB) / m_phiWidthEB;
0397 int regionNumEB = phiI + Nphi * etaI;
0398 return (int)regionNumEB;
0399 }
0400
0401
0402
0403
0404 int EcalEleCalibLooper::EERegionId(const int ics, const int ips) const {
0405 if (EEregionCheck(ics, ips))
0406 return -1;
0407 int phifake = m_phiStartEE;
0408 if (m_phiStartEE > m_phiEndEE)
0409 phifake = m_phiStartEE - 360;
0410 double radius = (ics - 50) * (ics - 50) + (ips - 50) * (ips - 50);
0411 radius = sqrt(radius);
0412 int Nphi = (m_phiEndEE - phifake) / m_phiWidthEE;
0413 double phi = atan2(static_cast<double>(ips - 50), static_cast<double>(ics - 50));
0414 phi = degrees(phi);
0415 if (phi < 0)
0416 phi += 360;
0417 int radI = static_cast<int>((radius - m_radStart) / m_radWidth);
0418 int phiI = static_cast<int>((m_phiEndEE - phi) / m_phiWidthEE);
0419 int regionNumEE = phiI + Nphi * radI;
0420 return regionNumEE;
0421 }
0422
0423
0424
0425
0426 inline int EcalEleCalibLooper::EEregionsNum() const {
0427 int phifake = m_phiStartEE;
0428 if (m_phiStartEE > m_phiEndEE)
0429 phifake = m_phiStartEE - 360;
0430 return ((m_radEnd - m_radStart) / m_radWidth) * ((m_phiEndEE - phifake) / m_phiWidthEE);
0431 }
0432
0433
0434
0435
0436 inline int EcalEleCalibLooper::EBregionsNum() const {
0437 int phi = m_phiStartEB;
0438 if (m_phiStartEB > m_phiEndEB)
0439 phi = m_phiStartEB - 360;
0440 return ((m_etaEnd - m_etaStart) / m_etaWidth) * ((m_phiEndEB - phi) / m_phiWidthEB);
0441 }
0442
0443
0444
0445
0446 void EcalEleCalibLooper::EBRegionDefinition() {
0447 int reg = -1;
0448 for (int it = 0; it < EBregionsNum(); ++it)
0449 m_regions.push_back(0);
0450 for (int eta = 0; eta < 170; ++eta)
0451 for (int phi = 0; phi < 360; ++phi) {
0452 reg = EBRegionId(eta, phi);
0453 m_xtalRegionId[EBDetId::unhashIndex(eta * 360 + phi).rawId()] = reg;
0454 if (reg == -1)
0455 continue;
0456 m_xtalPositionInRegion[EBDetId::unhashIndex(eta * 360 + phi).rawId()] = m_regions.at(reg);
0457 ++m_regions.at(reg);
0458 }
0459 }
0460
0461
0462
0463
0464 void EcalEleCalibLooper::EERegionDefinition() {
0465
0466 int EBnum = EBregionsNum();
0467 int EEnum = EEregionsNum();
0468 for (int it = 0; it < 2 * EEnum; ++it)
0469 m_regions.push_back(0);
0470
0471 int reg = -1;
0472 for (int ics = 0; ics < 100; ++ics)
0473 for (int ips = 0; ips < 100; ++ips) {
0474 int ireg = EERegionId(ics, ips);
0475 if (ireg == -1)
0476 reg = -1;
0477 else
0478 reg = EBnum + ireg;
0479 if (EEDetId::validDetId(ics + 1, ips + 1, 1)) {
0480 m_xtalRegionId[EEDetId(ics + 1, ips + 1, 1).rawId()] = reg;
0481 if (reg == -1)
0482 continue;
0483 m_xtalPositionInRegion[EEDetId(ics + 1, ips + 1, 1).rawId()] = m_regions.at(reg);
0484 ++m_regions.at(reg);
0485 }
0486 if (reg != -1)
0487 reg += EEnum;
0488 if (EEDetId::validDetId(ics + 1, ips + 1, -1)) {
0489 m_xtalRegionId[EEDetId(ics + 1, ips + 1, -1).rawId()] = reg;
0490 if (reg == -1)
0491 continue;
0492 m_xtalPositionInRegion[EEDetId(ics + 1, ips + 1, -1).rawId()] = m_regions.at(reg);
0493 ++m_regions.at(reg);
0494 }
0495 }
0496 }
0497
0498
0499
0500
0501 int EcalEleCalibLooper::EEregionCheck(const int ics, const int ips) const {
0502 int x = ics - 50;
0503 int y = ips - 50;
0504 double radius2 = x * x + y * y;
0505 if (radius2 < 10 * 10)
0506 return 1;
0507 if (radius2 > 50 * 50)
0508 return 1;
0509 if (radius2 < m_radStart * m_radStart)
0510 return 2;
0511 if (radius2 >= m_radEnd * m_radEnd)
0512 return 2;
0513 double phi = atan2(static_cast<double>(y), static_cast<double>(x));
0514 phi = degrees(phi);
0515 if (phi < 0)
0516 phi += 360;
0517 if (m_phiStartEE < m_phiEndEE && phi > m_phiStartEE && phi < m_phiEndEE)
0518 return 0;
0519 if (m_phiStartEE > m_phiEndEE && (phi > m_phiStartEE || phi < m_phiEndEE))
0520 return 0;
0521 return 3;
0522 }
0523
0524
0525
0526
0527 int EcalEleCalibLooper::etaShifter(const int etaOld) const {
0528 if (etaOld < 0)
0529 return etaOld + 85;
0530 else if (etaOld > 0)
0531 return etaOld + 84;
0532 assert(0 != etaOld);
0533 return 999;
0534 }
0535
0536