File indexing completed on 2024-07-03 04:18:12
0001 #include "RecoLocalCalo/HGCalRecProducers/plugins/HGCalCLUEAlgo.h"
0002
0003
0004 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0005 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0006 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0007 #include "Geometry/Records/interface/IdealGeometryRecord.h"
0008
0009 #include "RecoEcal/EgammaCoreTools/interface/PositionCalc.h"
0010
0011 #include "DataFormats/CaloRecHit/interface/CaloID.h"
0012 #include "oneapi/tbb/task_arena.h"
0013 #include "oneapi/tbb.h"
0014 #include <limits>
0015 #include "DataFormats/DetId/interface/DetId.h"
0016
0017 using namespace hgcal_clustering;
0018
0019 template <typename T, typename STRATEGY>
0020 void HGCalCLUEAlgoT<T, STRATEGY>::getEventSetupPerAlgorithm(const edm::EventSetup& es) {
0021 cells_.clear();
0022 numberOfClustersPerLayer_.clear();
0023 cells_.resize(2 * (maxlayer_ + 1));
0024 numberOfClustersPerLayer_.resize(2 * (maxlayer_ + 1), 0);
0025 }
0026
0027 template <typename T, typename STRATEGY>
0028 void HGCalCLUEAlgoT<T, STRATEGY>::populate(const HGCRecHitCollection& hits) {
0029
0030 if (dependSensor_) {
0031
0032
0033 computeThreshold();
0034 }
0035
0036 for (unsigned int i = 0; i < hits.size(); ++i) {
0037 const HGCRecHit& hgrh = hits[i];
0038 DetId detid = hgrh.detid();
0039 unsigned int layerOnSide = (rhtools_.getLayerWithOffset(detid) - 1);
0040
0041
0042
0043 float sigmaNoise = 1.f;
0044 if (dependSensor_) {
0045 int thickness_index = rhtools_.getSiThickIndex(detid);
0046 if (thickness_index == -1)
0047 thickness_index = maxNumberOfThickIndices_;
0048
0049 double storedThreshold = thresholds_[layerOnSide][thickness_index];
0050 if (detid.det() == DetId::HGCalHSi || detid.subdetId() == HGCHEF) {
0051 storedThreshold = thresholds_[layerOnSide][thickness_index + deltasi_index_regemfac_];
0052 }
0053 sigmaNoise = v_sigmaNoise_[layerOnSide][thickness_index];
0054
0055 if (hgrh.energy() < storedThreshold)
0056 continue;
0057
0058 }
0059 if (!dependSensor_ && hgrh.energy() < ecut_)
0060 continue;
0061 const GlobalPoint position(rhtools_.getPosition(detid));
0062 int offset = ((rhtools_.zside(detid) + 1) >> 1) * maxlayer_;
0063 int layer = layerOnSide + offset;
0064
0065 if (cells_[layer].layerDim3 == std::numeric_limits<float>::infinity())
0066 cells_[layer].layerDim3 = position.z();
0067
0068 cells_[layer].detid.emplace_back(detid);
0069 if constexpr (std::is_same_v<STRATEGY, HGCalScintillatorStrategy>) {
0070 cells_[layer].dim1.emplace_back(position.eta());
0071 cells_[layer].dim2.emplace_back(position.phi());
0072 }
0073 else {
0074 cells_[layer].dim1.emplace_back(position.x());
0075 cells_[layer].dim2.emplace_back(position.y());
0076 }
0077 cells_[layer].weight.emplace_back(hgrh.energy());
0078 cells_[layer].sigmaNoise.emplace_back(sigmaNoise);
0079 }
0080 }
0081
0082 template <typename T, typename STRATEGY>
0083 void HGCalCLUEAlgoT<T, STRATEGY>::prepareDataStructures(unsigned int l) {
0084 auto cellsSize = cells_[l].detid.size();
0085 cells_[l].rho.resize(cellsSize, 0.f);
0086 cells_[l].delta.resize(cellsSize, 9999999);
0087 cells_[l].nearestHigher.resize(cellsSize, -1);
0088 cells_[l].clusterIndex.resize(cellsSize, -1);
0089 cells_[l].followers.resize(cellsSize);
0090 cells_[l].isSeed.resize(cellsSize, false);
0091 }
0092
0093
0094
0095
0096
0097 template <typename T, typename STRATEGY>
0098 void HGCalCLUEAlgoT<T, STRATEGY>::makeClusters() {
0099
0100 tbb::this_task_arena::isolate([&] {
0101 tbb::parallel_for(size_t(0), size_t(2 * maxlayer_ + 2), [&](size_t i) {
0102 prepareDataStructures(i);
0103 T lt;
0104 lt.clear();
0105 lt.fill(cells_[i].dim1, cells_[i].dim2);
0106
0107 float delta;
0108 if constexpr (std::is_same_v<STRATEGY, HGCalSiliconStrategy>) {
0109
0110 float delta_c;
0111 if (i % maxlayer_ < lastLayerEE_)
0112 delta_c = vecDeltas_[0];
0113 else if (i % maxlayer_ < (firstLayerBH_ - 1))
0114 delta_c = vecDeltas_[1];
0115 else
0116 delta_c = vecDeltas_[2];
0117 delta = delta_c;
0118 } else {
0119 float delta_r = vecDeltas_[3];
0120 delta = delta_r;
0121 }
0122 LogDebug("HGCalCLUEAlgo") << "maxlayer: " << maxlayer_ << " lastLayerEE: " << lastLayerEE_
0123 << " firstLayerBH: " << firstLayerBH_ << "\n";
0124
0125 calculateLocalDensity(lt, i, delta);
0126 calculateDistanceToHigher(lt, i, delta);
0127 numberOfClustersPerLayer_[i] = findAndAssignClusters(i, delta);
0128 });
0129 });
0130 #if DEBUG_CLUSTERS_ALPAKA
0131 hgcalUtils::DumpLegacySoA dumperLegacySoA;
0132 dumperLegacySoA.dumpInfos(cells_, moduleType_);
0133 #endif
0134 }
0135
0136 template <typename T, typename STRATEGY>
0137 std::vector<reco::BasicCluster> HGCalCLUEAlgoT<T, STRATEGY>::getClusters(bool) {
0138 std::vector<int> offsets(numberOfClustersPerLayer_.size(), 0);
0139
0140 int maxClustersOnLayer = numberOfClustersPerLayer_[0];
0141
0142 for (unsigned layerId = 1; layerId < offsets.size(); ++layerId) {
0143 offsets[layerId] = offsets[layerId - 1] + numberOfClustersPerLayer_[layerId - 1];
0144
0145 maxClustersOnLayer = std::max(maxClustersOnLayer, numberOfClustersPerLayer_[layerId]);
0146 }
0147
0148 auto totalNumberOfClusters = offsets.back() + numberOfClustersPerLayer_.back();
0149 clusters_v_.resize(totalNumberOfClusters);
0150 std::vector<std::vector<int>> cellsIdInCluster;
0151 cellsIdInCluster.reserve(maxClustersOnLayer);
0152
0153 for (unsigned int layerId = 0; layerId < 2 * maxlayer_ + 2; ++layerId) {
0154 cellsIdInCluster.resize(numberOfClustersPerLayer_[layerId]);
0155 auto& cellsOnLayer = cells_[layerId];
0156 unsigned int numberOfCells = cellsOnLayer.detid.size();
0157 auto firstClusterIdx = offsets[layerId];
0158
0159 for (unsigned int i = 0; i < numberOfCells; ++i) {
0160 auto clusterIndex = cellsOnLayer.clusterIndex[i];
0161 if (clusterIndex != -1)
0162 cellsIdInCluster[clusterIndex].push_back(i);
0163 }
0164
0165 std::vector<std::pair<DetId, float>> thisCluster;
0166
0167 for (auto& cl : cellsIdInCluster) {
0168 float maxEnergyValue = std::numeric_limits<float>::min();
0169 int maxEnergyCellIndex = -1;
0170 DetId maxEnergyDetId;
0171 float energy = 0.f;
0172 int seedDetId = -1;
0173 float x = 0.f;
0174 float y = 0.f;
0175 float z = cellsOnLayer.layerDim3;
0176
0177 for (auto cellIdx : cl) {
0178 energy += cellsOnLayer.weight[cellIdx];
0179 if (cellsOnLayer.weight[cellIdx] > maxEnergyValue) {
0180 maxEnergyValue = cellsOnLayer.weight[cellIdx];
0181 maxEnergyCellIndex = cellIdx;
0182 maxEnergyDetId = cellsOnLayer.detid[cellIdx];
0183 }
0184 thisCluster.emplace_back(cellsOnLayer.detid[cellIdx], 1.f);
0185 if (cellsOnLayer.isSeed[cellIdx]) {
0186 seedDetId = cellsOnLayer.detid[cellIdx];
0187 }
0188 }
0189
0190 float total_weight_log = 0.f;
0191 float total_weight = energy;
0192
0193 if constexpr (std::is_same_v<STRATEGY, HGCalSiliconStrategy>) {
0194 auto thick = rhtools_.getSiThickIndex(maxEnergyDetId);
0195 for (auto cellIdx : cl) {
0196 const float d1 = cellsOnLayer.dim1[cellIdx] - cellsOnLayer.dim1[maxEnergyCellIndex];
0197 const float d2 = cellsOnLayer.dim2[cellIdx] - cellsOnLayer.dim2[maxEnergyCellIndex];
0198 if ((d1 * d1 + d2 * d2) < positionDeltaRho2_) {
0199 float Wi = std::max(thresholdW0_[thick] + std::log(cellsOnLayer.weight[cellIdx] / energy), 0.);
0200 x += cellsOnLayer.dim1[cellIdx] * Wi;
0201 y += cellsOnLayer.dim2[cellIdx] * Wi;
0202 total_weight_log += Wi;
0203 }
0204 }
0205 } else {
0206 for (auto cellIdx : cl) {
0207 auto position = rhtools_.getPosition(cellsOnLayer.detid[cellIdx]);
0208 x += position.x() * cellsOnLayer.weight[cellIdx];
0209 y += position.y() * cellsOnLayer.weight[cellIdx];
0210 }
0211 }
0212
0213 if constexpr (std::is_same_v<STRATEGY, HGCalSiliconStrategy>) {
0214 total_weight = total_weight_log;
0215 }
0216
0217 if (total_weight != 0.) {
0218 float inv_tot_weight = 1.f / total_weight;
0219 x *= inv_tot_weight;
0220 y *= inv_tot_weight;
0221 } else {
0222 x = cellsOnLayer.dim1[maxEnergyCellIndex];
0223 y = cellsOnLayer.dim2[maxEnergyCellIndex];
0224 }
0225 math::XYZPoint position = math::XYZPoint(x, y, z);
0226
0227 auto globalClusterIndex = cellsOnLayer.clusterIndex[cl[0]] + firstClusterIdx;
0228
0229 clusters_v_[globalClusterIndex] =
0230 reco::BasicCluster(energy, position, reco::CaloID::DET_HGCAL_ENDCAP, thisCluster, algoId_);
0231 clusters_v_[globalClusterIndex].setSeed(seedDetId);
0232 thisCluster.clear();
0233 }
0234
0235 cellsIdInCluster.clear();
0236 }
0237 return clusters_v_;
0238 }
0239 template <typename T, typename STRATEGY>
0240 void HGCalCLUEAlgoT<T, STRATEGY>::calculateLocalDensity(const T& lt,
0241 const unsigned int layerId,
0242 float delta,
0243 HGCalSiliconStrategy strategy) {
0244 auto& cellsOnLayer = cells_[layerId];
0245 unsigned int numberOfCells = cellsOnLayer.detid.size();
0246 for (unsigned int i = 0; i < numberOfCells; i++) {
0247 std::array<int, 4> search_box = lt.searchBox(cellsOnLayer.dim1[i] - delta,
0248 cellsOnLayer.dim1[i] + delta,
0249 cellsOnLayer.dim2[i] - delta,
0250 cellsOnLayer.dim2[i] + delta);
0251
0252 for (int xBin = search_box[0]; xBin < search_box[1] + 1; ++xBin) {
0253 for (int yBin = search_box[2]; yBin < search_box[3] + 1; ++yBin) {
0254 int binId = lt.getGlobalBinByBin(xBin, yBin);
0255 size_t binSize = lt[binId].size();
0256
0257 for (unsigned int j = 0; j < binSize; j++) {
0258 unsigned int otherId = lt[binId][j];
0259 if (distance(lt, i, otherId, layerId) < delta) {
0260 cellsOnLayer.rho[i] += (i == otherId ? 1.f : 0.5f) * cellsOnLayer.weight[otherId];
0261 }
0262 }
0263 }
0264 }
0265 LogDebug("HGCalCLUEAlgo") << "Debugging calculateLocalDensity: \n"
0266 << " cell: " << i << " eta: " << cellsOnLayer.dim1[i] << " phi: " << cellsOnLayer.dim2[i]
0267 << " energy: " << cellsOnLayer.weight[i] << " density: " << cellsOnLayer.rho[i] << "\n";
0268 }
0269 }
0270 template <typename T, typename STRATEGY>
0271 void HGCalCLUEAlgoT<T, STRATEGY>::calculateLocalDensity(const T& lt,
0272 const unsigned int layerId,
0273 float delta,
0274 HGCalScintillatorStrategy strategy) {
0275 auto& cellsOnLayer = cells_[layerId];
0276 unsigned int numberOfCells = cellsOnLayer.detid.size();
0277 for (unsigned int i = 0; i < numberOfCells; i++) {
0278 std::array<int, 4> search_box = lt.searchBox(cellsOnLayer.dim1[i] - delta,
0279 cellsOnLayer.dim1[i] + delta,
0280 cellsOnLayer.dim2[i] - delta,
0281 cellsOnLayer.dim2[i] + delta);
0282 cellsOnLayer.rho[i] += cellsOnLayer.weight[i];
0283 float northeast(0), northwest(0), southeast(0), southwest(0), all(0);
0284 for (int etaBin = search_box[0]; etaBin < search_box[1] + 1; ++etaBin) {
0285 for (int phiBin = search_box[2]; phiBin < search_box[3] + 1; ++phiBin) {
0286 int phi = (phiBin % T::type::nRows);
0287 int binId = lt.getGlobalBinByBin(etaBin, phi);
0288 size_t binSize = lt[binId].size();
0289 for (unsigned int j = 0; j < binSize; j++) {
0290 unsigned int otherId = lt[binId][j];
0291 if (distance(lt, i, otherId, layerId) < delta) {
0292 int iPhi = HGCScintillatorDetId(cellsOnLayer.detid[i]).iphi();
0293 int otherIPhi = HGCScintillatorDetId(cellsOnLayer.detid[otherId]).iphi();
0294 int iEta = HGCScintillatorDetId(cellsOnLayer.detid[i]).ieta();
0295 int otherIEta = HGCScintillatorDetId(cellsOnLayer.detid[otherId]).ieta();
0296 int dIPhi = otherIPhi - iPhi;
0297 dIPhi += abs(dIPhi) < 2 ? 0
0298 : dIPhi < 0 ? scintMaxIphi_
0299 : -scintMaxIphi_;
0300 int dIEta = otherIEta - iEta;
0301 LogDebug("HGCalCLUEAlgo") << " Debugging calculateLocalDensity for Scintillator: \n"
0302 << " cell: " << otherId << " energy: " << cellsOnLayer.weight[otherId]
0303 << " otherIPhi: " << otherIPhi << " iPhi: " << iPhi << " otherIEta: " << otherIEta
0304 << " iEta: " << iEta << "\n";
0305
0306 if (otherId != i) {
0307 auto neighborCellContribution = 0.5f * cellsOnLayer.weight[otherId];
0308 all += neighborCellContribution;
0309 if (dIPhi >= 0 && dIEta >= 0)
0310 northeast += neighborCellContribution;
0311 if (dIPhi <= 0 && dIEta >= 0)
0312 southeast += neighborCellContribution;
0313 if (dIPhi >= 0 && dIEta <= 0)
0314 northwest += neighborCellContribution;
0315 if (dIPhi <= 0 && dIEta <= 0)
0316 southwest += neighborCellContribution;
0317 }
0318 LogDebug("HGCalCLUEAlgo") << " Debugging calculateLocalDensity for Scintillator: \n"
0319 << " northeast: " << northeast << " southeast: " << southeast
0320 << " northwest: " << northwest << " southwest: " << southwest << "\n";
0321 }
0322 }
0323 }
0324 }
0325 float neighborsval = (std::max(northeast, northwest) > std::max(southeast, southwest))
0326 ? std::max(northeast, northwest)
0327 : std::max(southeast, southwest);
0328 if (use2x2_)
0329 cellsOnLayer.rho[i] += neighborsval;
0330 else
0331 cellsOnLayer.rho[i] += all;
0332 LogDebug("HGCalCLUEAlgo") << "Debugging calculateLocalDensity: \n"
0333 << " cell: " << i << " eta: " << cellsOnLayer.dim1[i] << " phi: " << cellsOnLayer.dim2[i]
0334 << " energy: " << cellsOnLayer.weight[i] << " density: " << cellsOnLayer.rho[i] << "\n";
0335 }
0336 }
0337 template <typename T, typename STRATEGY>
0338 void HGCalCLUEAlgoT<T, STRATEGY>::calculateLocalDensity(const T& lt, const unsigned int layerId, float delta) {
0339 if constexpr (std::is_same_v<STRATEGY, HGCalSiliconStrategy>) {
0340 calculateLocalDensity(lt, layerId, delta, HGCalSiliconStrategy());
0341 } else {
0342 calculateLocalDensity(lt, layerId, delta, HGCalScintillatorStrategy());
0343 }
0344 }
0345
0346 template <typename T, typename STRATEGY>
0347 void HGCalCLUEAlgoT<T, STRATEGY>::calculateDistanceToHigher(const T& lt, const unsigned int layerId, float delta) {
0348 auto& cellsOnLayer = cells_[layerId];
0349 unsigned int numberOfCells = cellsOnLayer.detid.size();
0350
0351 for (unsigned int i = 0; i < numberOfCells; i++) {
0352
0353 float maxDelta = std::numeric_limits<float>::max();
0354 float i_delta = maxDelta;
0355 int i_nearestHigher = -1;
0356 float rho_max = 0.f;
0357 auto range = outlierDeltaFactor_ * delta;
0358 std::array<int, 4> search_box = lt.searchBox(cellsOnLayer.dim1[i] - range,
0359 cellsOnLayer.dim1[i] + range,
0360 cellsOnLayer.dim2[i] - range,
0361 cellsOnLayer.dim2[i] + range);
0362
0363 for (int dim1Bin = search_box[0]; dim1Bin < search_box[1] + 1; ++dim1Bin) {
0364 for (int dim2Bin = search_box[2]; dim2Bin < search_box[3] + 1; ++dim2Bin) {
0365
0366 size_t binId = lt.getGlobalBinByBin(dim1Bin, dim2Bin);
0367 if constexpr (std::is_same_v<STRATEGY, HGCalScintillatorStrategy>)
0368 binId = lt.getGlobalBinByBin(dim1Bin, (dim2Bin % T::type::nRows));
0369
0370 size_t binSize = lt[binId].size();
0371
0372
0373 for (unsigned int j = 0; j < binSize; j++) {
0374 unsigned int otherId = lt[binId][j];
0375 float dist = distance2(lt, i, otherId, layerId);
0376 bool foundHigher =
0377 (cellsOnLayer.rho[otherId] > cellsOnLayer.rho[i]) ||
0378 (cellsOnLayer.rho[otherId] == cellsOnLayer.rho[i] && cellsOnLayer.detid[otherId] > cellsOnLayer.detid[i]);
0379 if (!foundHigher) {
0380 continue;
0381 }
0382 if ((dist < i_delta) || ((dist == i_delta) && (cellsOnLayer.rho[otherId] > rho_max)) ||
0383 ((dist == i_delta) && (cellsOnLayer.rho[otherId] == rho_max) &&
0384 (cellsOnLayer.detid[otherId] > cellsOnLayer.detid[i]))) {
0385 rho_max = cellsOnLayer.rho[otherId];
0386 i_delta = dist;
0387 i_nearestHigher = otherId;
0388 }
0389 }
0390 }
0391 }
0392 bool foundNearestHigherInSearchBox = (i_delta != maxDelta);
0393 if (foundNearestHigherInSearchBox) {
0394 cellsOnLayer.delta[i] = std::sqrt(i_delta);
0395 cellsOnLayer.nearestHigher[i] = i_nearestHigher;
0396 } else {
0397
0398
0399 cellsOnLayer.delta[i] = maxDelta;
0400 cellsOnLayer.nearestHigher[i] = -1;
0401 }
0402
0403 LogDebug("HGCalCLUEAlgo") << "Debugging calculateDistanceToHigher: \n"
0404 << " cell: " << i << " eta: " << cellsOnLayer.dim1[i] << " phi: " << cellsOnLayer.dim2[i]
0405 << " energy: " << cellsOnLayer.weight[i] << " density: " << cellsOnLayer.rho[i]
0406 << " nearest higher: " << cellsOnLayer.nearestHigher[i]
0407 << " distance: " << cellsOnLayer.delta[i] << "\n";
0408 }
0409 }
0410
0411 template <typename T, typename STRATEGY>
0412 int HGCalCLUEAlgoT<T, STRATEGY>::findAndAssignClusters(const unsigned int layerId, float delta) {
0413
0414
0415
0416
0417 unsigned int nClustersOnLayer = 0;
0418 auto& cellsOnLayer = cells_[layerId];
0419 unsigned int numberOfCells = cellsOnLayer.detid.size();
0420 std::vector<int> localStack;
0421
0422 for (unsigned int i = 0; i < numberOfCells; i++) {
0423 float rho_c = kappa_ * cellsOnLayer.sigmaNoise[i];
0424
0425 cellsOnLayer.clusterIndex[i] = -1;
0426 bool isSeed = (cellsOnLayer.delta[i] > delta) && (cellsOnLayer.rho[i] >= rho_c);
0427 bool isOutlier = (cellsOnLayer.delta[i] > outlierDeltaFactor_ * delta) && (cellsOnLayer.rho[i] < rho_c);
0428 if (isSeed) {
0429 cellsOnLayer.clusterIndex[i] = nClustersOnLayer;
0430 cellsOnLayer.isSeed[i] = true;
0431 nClustersOnLayer++;
0432 localStack.push_back(i);
0433
0434 } else if (!isOutlier) {
0435 cellsOnLayer.followers[cellsOnLayer.nearestHigher[i]].push_back(i);
0436 }
0437 }
0438
0439
0440 while (!localStack.empty()) {
0441 int endStack = localStack.back();
0442 auto& thisSeed = cellsOnLayer.followers[endStack];
0443 localStack.pop_back();
0444
0445
0446 for (int j : thisSeed) {
0447
0448 cellsOnLayer.clusterIndex[j] = cellsOnLayer.clusterIndex[endStack];
0449
0450 localStack.push_back(j);
0451 }
0452 }
0453 return nClustersOnLayer;
0454 }
0455
0456 template <typename T, typename STRATEGY>
0457 void HGCalCLUEAlgoT<T, STRATEGY>::computeThreshold() {
0458
0459
0460
0461
0462
0463
0464
0465
0466 if (initialized_)
0467 return;
0468
0469 initialized_ = true;
0470
0471 std::vector<double> dummy;
0472
0473 dummy.resize(maxNumberOfThickIndices_ + !isNose_, 0);
0474 thresholds_.resize(maxlayer_, dummy);
0475 v_sigmaNoise_.resize(maxlayer_, dummy);
0476
0477 for (unsigned ilayer = 1; ilayer <= maxlayer_; ++ilayer) {
0478 for (unsigned ithick = 0; ithick < maxNumberOfThickIndices_; ++ithick) {
0479 float sigmaNoise = 0.001f * fcPerEle_ * nonAgedNoises_[ithick] * dEdXweights_[ilayer] /
0480 (fcPerMip_[ithick] * thicknessCorrection_[ithick]);
0481 thresholds_[ilayer - 1][ithick] = sigmaNoise * ecut_;
0482 v_sigmaNoise_[ilayer - 1][ithick] = sigmaNoise;
0483 LogDebug("HGCalCLUEAlgo") << "ilayer: " << ilayer << " nonAgedNoises: " << nonAgedNoises_[ithick]
0484 << " fcPerEle: " << fcPerEle_ << " fcPerMip: " << fcPerMip_[ithick]
0485 << " noiseMip: " << fcPerEle_ * nonAgedNoises_[ithick] / fcPerMip_[ithick]
0486 << " sigmaNoise: " << sigmaNoise << "\n";
0487 }
0488
0489 if (!isNose_) {
0490 float scintillators_sigmaNoise = 0.001f * noiseMip_ * dEdXweights_[ilayer] / sciThicknessCorrection_;
0491 thresholds_[ilayer - 1][maxNumberOfThickIndices_] = ecut_ * scintillators_sigmaNoise;
0492 v_sigmaNoise_[ilayer - 1][maxNumberOfThickIndices_] = scintillators_sigmaNoise;
0493 LogDebug("HGCalCLUEAlgo") << "ilayer: " << ilayer << " noiseMip: " << noiseMip_
0494 << " scintillators_sigmaNoise: " << scintillators_sigmaNoise << "\n";
0495 }
0496 }
0497 }
0498
0499
0500 template class HGCalCLUEAlgoT<HGCalSiliconLayerTiles, HGCalSiliconStrategy>;
0501 template class HGCalCLUEAlgoT<HGCalScintillatorLayerTiles, HGCalScintillatorStrategy>;
0502 template class HGCalCLUEAlgoT<HFNoseLayerTiles, HGCalSiliconStrategy>;