File indexing completed on 2021-02-14 14:23:04
0001 #include "L1Trigger/L1THGCal/interface/backend/HGCalShowerShape.h"
0002 #include "DataFormats/L1THGCal/interface/HGCalMulticluster.h"
0003 #include "DataFormats/L1THGCal/interface/HGCalCluster.h"
0004 #include "DataFormats/Math/interface/deltaPhi.h"
0005
0006 #include <unordered_map>
0007 #include <numeric>
0008
0009 HGCalShowerShape::HGCalShowerShape(const edm::ParameterSet& conf)
0010 : threshold_(conf.getParameter<double>("shape_threshold")),
0011 distance_(conf.getParameter<double>("shape_distance")) {}
0012
0013
0014
0015 float HGCalShowerShape::meanX(const std::vector<pair<float, float>>& energy_X_tc) const {
0016 float Etot = 0;
0017 float X_sum = 0;
0018
0019 for (const auto& energy_X : energy_X_tc) {
0020 X_sum += energy_X.first * energy_X.second;
0021 Etot += energy_X.first;
0022 }
0023
0024 float X_mean = 0;
0025 if (Etot > 0)
0026 X_mean = X_sum / Etot;
0027 return X_mean;
0028 }
0029
0030 int HGCalShowerShape::firstLayer(const l1t::HGCalMulticluster& c3d) const {
0031 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0032
0033 int firstLayer = 999;
0034
0035 for (const auto& id_clu : clustersPtrs) {
0036 if (!pass(*id_clu.second, c3d))
0037 continue;
0038 int layer = triggerTools_.layerWithOffset(id_clu.second->detId());
0039 if (layer < firstLayer)
0040 firstLayer = layer;
0041 }
0042
0043 return firstLayer;
0044 }
0045
0046 int HGCalShowerShape::maxLayer(const l1t::HGCalMulticluster& c3d) const {
0047 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0048 std::unordered_map<int, float> layers_pt;
0049 float max_pt = 0.;
0050 int max_layer = 0;
0051 for (const auto& id_cluster : clustersPtrs) {
0052 if (!pass(*id_cluster.second, c3d))
0053 continue;
0054 unsigned layer = triggerTools_.layerWithOffset(id_cluster.second->detId());
0055 auto itr_insert = layers_pt.emplace(layer, 0.);
0056 itr_insert.first->second += id_cluster.second->pt();
0057 if (itr_insert.first->second > max_pt) {
0058 max_pt = itr_insert.first->second;
0059 max_layer = layer;
0060 }
0061 }
0062 return max_layer;
0063 }
0064
0065 int HGCalShowerShape::lastLayer(const l1t::HGCalMulticluster& c3d) const {
0066 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0067
0068 int lastLayer = -999;
0069
0070 for (const auto& id_clu : clustersPtrs) {
0071 if (!pass(*id_clu.second, c3d))
0072 continue;
0073 int layer = triggerTools_.layerWithOffset(id_clu.second->detId());
0074 if (layer > lastLayer)
0075 lastLayer = layer;
0076 }
0077
0078 return lastLayer;
0079 }
0080
0081 int HGCalShowerShape::coreShowerLength(const l1t::HGCalMulticluster& c3d,
0082 const HGCalTriggerGeometryBase& triggerGeometry) const {
0083 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0084 unsigned nlayers = triggerTools_.layers(ForwardSubdetector::ForwardEmpty);
0085 std::vector<bool> layers(nlayers);
0086 for (const auto& id_cluster : clustersPtrs) {
0087 if (!pass(*id_cluster.second, c3d))
0088 continue;
0089 unsigned layer = triggerGeometry.triggerLayer(id_cluster.second->detId());
0090 if (triggerTools_.isNose(id_cluster.second->detId()))
0091 nlayers = triggerTools_.layers(ForwardSubdetector::HFNose);
0092 else {
0093 nlayers = triggerTools_.layers(ForwardSubdetector::ForwardEmpty);
0094 }
0095 if (layer == 0 || layer > nlayers)
0096 continue;
0097 layers[layer - 1] = true;
0098 }
0099 int length = 0;
0100 int maxlength = 0;
0101 for (bool layer : layers) {
0102 if (layer)
0103 length++;
0104 else
0105 length = 0;
0106 if (length > maxlength)
0107 maxlength = length;
0108 }
0109 return maxlength;
0110 }
0111
0112 float HGCalShowerShape::percentileLayer(const l1t::HGCalMulticluster& c3d,
0113 const HGCalTriggerGeometryBase& triggerGeometry,
0114 float quantile) const {
0115 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0116 unsigned nlayers = triggerTools_.layers(ForwardSubdetector::ForwardEmpty);
0117 std::vector<double> layers(nlayers, 0);
0118 for (const auto& id_clu : clustersPtrs) {
0119 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0120
0121 for (const auto& id_tc : triggerCells) {
0122 if (!pass(*id_tc.second, c3d))
0123 continue;
0124 unsigned layer = triggerGeometry.triggerLayer(id_tc.second->detId());
0125 if (triggerTools_.isNose(id_tc.second->detId()))
0126 nlayers = triggerTools_.layers(ForwardSubdetector::HFNose);
0127 else {
0128 nlayers = triggerTools_.layers(ForwardSubdetector::ForwardEmpty);
0129 }
0130 if (layer == 0 || layer > nlayers)
0131 continue;
0132 layers[layer - 1] += id_tc.second->pt();
0133 }
0134 }
0135 std::partial_sum(layers.begin(), layers.end(), layers.begin());
0136 double pt_threshold = layers.back() * quantile;
0137 unsigned percentile = 0;
0138 for (double pt : layers) {
0139 if (pt > pt_threshold) {
0140 break;
0141 }
0142 percentile++;
0143 }
0144
0145 double pt0 = (percentile > 0 ? layers[percentile - 1] : 0.);
0146 double pt1 = (percentile < layers.size() ? layers[percentile] : layers.back());
0147 return percentile + (pt1 - pt0 > 0. ? (pt_threshold - pt0) / (pt1 - pt0) : 0.);
0148 }
0149
0150 float HGCalShowerShape::percentileTriggerCells(const l1t::HGCalMulticluster& c3d, float quantile) const {
0151 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0152 std::set<double> ordered_tcs;
0153 double pt_sum = 0.;
0154 for (const auto& id_clu : clustersPtrs) {
0155 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0156
0157 for (const auto& id_tc : triggerCells) {
0158 if (!pass(*id_tc.second, c3d))
0159 continue;
0160 ordered_tcs.emplace(id_tc.second->pt());
0161 pt_sum += id_tc.second->pt();
0162 }
0163 }
0164 double pt_threshold = pt_sum * quantile;
0165 double partial_sum = 0.;
0166 double partial_sum_prev = 0.;
0167 int ntc = 0;
0168 for (auto itr = ordered_tcs.rbegin(); itr != ordered_tcs.rend(); ++itr) {
0169 partial_sum_prev = partial_sum;
0170 partial_sum += *itr;
0171 ntc++;
0172 if (partial_sum > pt_threshold) {
0173 break;
0174 }
0175 }
0176
0177 return ntc - 1 +
0178 (partial_sum - partial_sum_prev > 0. ? (pt_threshold - partial_sum_prev) / (partial_sum - partial_sum_prev)
0179 : 0.);
0180 }
0181
0182 float HGCalShowerShape::sigmaEtaEtaTot(const l1t::HGCalMulticluster& c3d) const {
0183 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0184
0185 std::vector<std::pair<float, float>> tc_energy_eta;
0186
0187 for (const auto& id_clu : clustersPtrs) {
0188 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0189
0190 for (const auto& id_tc : triggerCells) {
0191 if (!pass(*id_tc.second, c3d))
0192 continue;
0193 tc_energy_eta.emplace_back(std::make_pair(id_tc.second->energy(), id_tc.second->eta()));
0194 }
0195 }
0196
0197 float SeeTot = sigmaXX(tc_energy_eta, c3d.eta());
0198
0199 return SeeTot;
0200 }
0201
0202 float HGCalShowerShape::sigmaPhiPhiTot(const l1t::HGCalMulticluster& c3d) const {
0203 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0204
0205 std::vector<std::pair<float, float>> tc_energy_phi;
0206
0207 for (const auto& id_clu : clustersPtrs) {
0208 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0209
0210 for (const auto& id_tc : triggerCells) {
0211 if (!pass(*id_tc.second, c3d))
0212 continue;
0213 tc_energy_phi.emplace_back(std::make_pair(id_tc.second->energy(), id_tc.second->phi()));
0214 }
0215 }
0216
0217 float SppTot = sigmaPhiPhi(tc_energy_phi, c3d.phi());
0218
0219 return SppTot;
0220 }
0221
0222 float HGCalShowerShape::sigmaRRTot(const l1t::HGCalMulticluster& c3d) const {
0223 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0224
0225 std::vector<std::pair<float, float>> tc_energy_r;
0226
0227 for (const auto& id_clu : clustersPtrs) {
0228 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0229
0230 for (const auto& id_tc : triggerCells) {
0231 if (!pass(*id_tc.second, c3d))
0232 continue;
0233 float r = (id_tc.second->position().z() != 0.
0234 ? std::sqrt(pow(id_tc.second->position().x(), 2) + pow(id_tc.second->position().y(), 2)) /
0235 std::abs(id_tc.second->position().z())
0236 : 0.);
0237 tc_energy_r.emplace_back(std::make_pair(id_tc.second->energy(), r));
0238 }
0239 }
0240
0241 float r_mean = meanX(tc_energy_r);
0242 float Szz = sigmaXX(tc_energy_r, r_mean);
0243
0244 return Szz;
0245 }
0246
0247 float HGCalShowerShape::sigmaEtaEtaMax(const l1t::HGCalMulticluster& c3d) const {
0248 std::unordered_map<int, std::vector<std::pair<float, float>>> tc_layer_energy_eta;
0249 std::unordered_map<int, LorentzVector> layer_LV;
0250
0251 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0252
0253 for (const auto& id_clu : clustersPtrs) {
0254 unsigned layer = triggerTools_.layerWithOffset(id_clu.second->detId());
0255
0256 layer_LV[layer] += id_clu.second->p4();
0257
0258 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0259
0260 for (const auto& id_tc : triggerCells) {
0261 if (!pass(*id_tc.second, c3d))
0262 continue;
0263 tc_layer_energy_eta[layer].emplace_back(std::make_pair(id_tc.second->energy(), id_tc.second->eta()));
0264 }
0265 }
0266
0267 float SigmaEtaEtaMax = 0;
0268
0269 for (auto& tc_iter : tc_layer_energy_eta) {
0270 const std::vector<std::pair<float, float>>& energy_eta_layer = tc_iter.second;
0271 const LorentzVector& LV_layer = layer_LV[tc_iter.first];
0272 float SigmaEtaEtaLayer = sigmaXX(energy_eta_layer, LV_layer.eta());
0273 if (SigmaEtaEtaLayer > SigmaEtaEtaMax)
0274 SigmaEtaEtaMax = SigmaEtaEtaLayer;
0275 }
0276
0277 return SigmaEtaEtaMax;
0278 }
0279
0280 float HGCalShowerShape::sigmaPhiPhiMax(const l1t::HGCalMulticluster& c3d) const {
0281 std::unordered_map<int, std::vector<std::pair<float, float>>> tc_layer_energy_phi;
0282 std::unordered_map<int, LorentzVector> layer_LV;
0283
0284 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0285
0286 for (const auto& id_clu : clustersPtrs) {
0287 unsigned layer = triggerTools_.layerWithOffset(id_clu.second->detId());
0288
0289 layer_LV[layer] += id_clu.second->p4();
0290
0291 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0292
0293 for (const auto& id_tc : triggerCells) {
0294 if (!pass(*id_tc.second, c3d))
0295 continue;
0296 tc_layer_energy_phi[layer].emplace_back(std::make_pair(id_tc.second->energy(), id_tc.second->phi()));
0297 }
0298 }
0299
0300 float SigmaPhiPhiMax = 0;
0301
0302 for (auto& tc_iter : tc_layer_energy_phi) {
0303 const std::vector<std::pair<float, float>>& energy_phi_layer = tc_iter.second;
0304 const LorentzVector& LV_layer = layer_LV[tc_iter.first];
0305 float SigmaPhiPhiLayer = sigmaPhiPhi(energy_phi_layer, LV_layer.phi());
0306 if (SigmaPhiPhiLayer > SigmaPhiPhiMax)
0307 SigmaPhiPhiMax = SigmaPhiPhiLayer;
0308 }
0309
0310 return SigmaPhiPhiMax;
0311 }
0312
0313 float HGCalShowerShape::sigmaRRMax(const l1t::HGCalMulticluster& c3d) const {
0314 std::unordered_map<int, std::vector<std::pair<float, float>>> tc_layer_energy_r;
0315
0316 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0317
0318 for (const auto& id_clu : clustersPtrs) {
0319 unsigned layer = triggerTools_.layerWithOffset(id_clu.second->detId());
0320
0321 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0322
0323 for (const auto& id_tc : triggerCells) {
0324 if (!pass(*id_tc.second, c3d))
0325 continue;
0326 float r = (id_tc.second->position().z() != 0.
0327 ? std::sqrt(pow(id_tc.second->position().x(), 2) + pow(id_tc.second->position().y(), 2)) /
0328 std::abs(id_tc.second->position().z())
0329 : 0.);
0330 tc_layer_energy_r[layer].emplace_back(std::make_pair(id_tc.second->energy(), r));
0331 }
0332 }
0333
0334 float SigmaRRMax = 0;
0335
0336 for (auto& tc_iter : tc_layer_energy_r) {
0337 const std::vector<std::pair<float, float>>& energy_r_layer = tc_iter.second;
0338 float r_mean_layer = meanX(energy_r_layer);
0339 float SigmaRRLayer = sigmaXX(energy_r_layer, r_mean_layer);
0340 if (SigmaRRLayer > SigmaRRMax)
0341 SigmaRRMax = SigmaRRLayer;
0342 }
0343
0344 return SigmaRRMax;
0345 }
0346
0347 float HGCalShowerShape::sigmaRRMean(const l1t::HGCalMulticluster& c3d, float radius) const {
0348 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0349
0350 std::unordered_map<int, std::vector<edm::Ptr<l1t::HGCalTriggerCell>>> layers_tcs;
0351 for (const auto& id_clu : clustersPtrs) {
0352 unsigned layer = triggerTools_.layerWithOffset(id_clu.second->detId());
0353 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0354 for (const auto& id_tc : triggerCells) {
0355 if (!pass(*id_tc.second, c3d))
0356 continue;
0357 layers_tcs[layer].emplace_back(id_tc.second);
0358 }
0359 }
0360
0361
0362 std::unordered_map<int, std::vector<std::pair<float, float>>> tc_layers_energy_r;
0363 for (const auto& layer_tcs : layers_tcs) {
0364 int layer = layer_tcs.first;
0365 edm::Ptr<l1t::HGCalTriggerCell> max_tc = layer_tcs.second.front();
0366 for (const auto& tc : layer_tcs.second) {
0367 if (tc->energy() > max_tc->energy())
0368 max_tc = tc;
0369 }
0370 for (const auto& tc : layer_tcs.second) {
0371 double dx = tc->position().x() - max_tc->position().x();
0372 double dy = tc->position().y() - max_tc->position().y();
0373 double distance_to_max = std::sqrt(dx * dx + dy * dy);
0374 if (distance_to_max < radius) {
0375 float r = (tc->position().z() != 0.
0376 ? std::sqrt(tc->position().x() * tc->position().x() + tc->position().y() * tc->position().y()) /
0377 std::abs(tc->position().z())
0378 : 0.);
0379 tc_layers_energy_r[layer].emplace_back(std::make_pair(tc->energy(), r));
0380 }
0381 }
0382 }
0383
0384
0385 std::vector<std::pair<float, float>> layers_energy_srr2;
0386 for (const auto& layer_energy_r : tc_layers_energy_r) {
0387 const auto& energies_r = layer_energy_r.second;
0388 float r_mean_layer = meanX(energies_r);
0389 float srr = sigmaXX(energies_r, r_mean_layer);
0390 double energy_sum = 0.;
0391 for (const auto& energy_r : energies_r) {
0392 energy_sum += energy_r.first;
0393 }
0394 layers_energy_srr2.emplace_back(std::make_pair(energy_sum, srr * srr));
0395 }
0396
0397 float srr2_mean = meanX(layers_energy_srr2);
0398 return std::sqrt(srr2_mean);
0399 }
0400
0401 float HGCalShowerShape::eMax(const l1t::HGCalMulticluster& c3d) const {
0402 std::unordered_map<int, float> layer_energy;
0403
0404 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0405
0406 for (const auto& id_clu : clustersPtrs) {
0407 if (!pass(*id_clu.second, c3d))
0408 continue;
0409 unsigned layer = triggerTools_.layerWithOffset(id_clu.second->detId());
0410 layer_energy[layer] += id_clu.second->energy();
0411 }
0412
0413 float EMax = 0;
0414
0415 for (const auto& layer : layer_energy) {
0416 if (layer.second > EMax)
0417 EMax = layer.second;
0418 }
0419
0420 return EMax;
0421 }
0422
0423 float HGCalShowerShape::meanZ(const l1t::HGCalMulticluster& c3d) const {
0424 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0425
0426 std::vector<std::pair<float, float>> tc_energy_z;
0427
0428 for (const auto& id_clu : clustersPtrs) {
0429 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0430
0431 for (const auto& id_tc : triggerCells) {
0432 if (!pass(*id_tc.second, c3d))
0433 continue;
0434 tc_energy_z.emplace_back(id_tc.second->energy(), id_tc.second->position().z());
0435 }
0436 }
0437
0438 return meanX(tc_energy_z);
0439 }
0440
0441 float HGCalShowerShape::sigmaZZ(const l1t::HGCalMulticluster& c3d) const {
0442 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalCluster>>& clustersPtrs = c3d.constituents();
0443
0444 std::vector<std::pair<float, float>> tc_energy_z;
0445
0446 for (const auto& id_clu : clustersPtrs) {
0447 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& triggerCells = id_clu.second->constituents();
0448
0449 for (const auto& id_tc : triggerCells) {
0450 if (!pass(*id_tc.second, c3d))
0451 continue;
0452 tc_energy_z.emplace_back(std::make_pair(id_tc.second->energy(), id_tc.second->position().z()));
0453 }
0454 }
0455
0456 float z_mean = meanX(tc_energy_z);
0457 float Szz = sigmaXX(tc_energy_z, z_mean);
0458
0459 return Szz;
0460 }
0461
0462 float HGCalShowerShape::sigmaEtaEtaTot(const l1t::HGCalCluster& c2d) const {
0463 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& cellsPtrs = c2d.constituents();
0464
0465 std::vector<std::pair<float, float>> tc_energy_eta;
0466
0467 for (const auto& id_cell : cellsPtrs) {
0468 if (!pass(*id_cell.second, c2d))
0469 continue;
0470 tc_energy_eta.emplace_back(std::make_pair(id_cell.second->energy(), id_cell.second->eta()));
0471 }
0472
0473 float See = sigmaXX(tc_energy_eta, c2d.eta());
0474
0475 return See;
0476 }
0477
0478 float HGCalShowerShape::sigmaPhiPhiTot(const l1t::HGCalCluster& c2d) const {
0479 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& cellsPtrs = c2d.constituents();
0480
0481 std::vector<std::pair<float, float>> tc_energy_phi;
0482
0483 for (const auto& id_cell : cellsPtrs) {
0484 if (!pass(*id_cell.second, c2d))
0485 continue;
0486 tc_energy_phi.emplace_back(std::make_pair(id_cell.second->energy(), id_cell.second->phi()));
0487 }
0488
0489 float Spp = sigmaPhiPhi(tc_energy_phi, c2d.phi());
0490
0491 return Spp;
0492 }
0493
0494 float HGCalShowerShape::sigmaRRTot(const l1t::HGCalCluster& c2d) const {
0495 const std::unordered_map<uint32_t, edm::Ptr<l1t::HGCalTriggerCell>>& cellsPtrs = c2d.constituents();
0496
0497 std::vector<std::pair<float, float>> tc_energy_r;
0498
0499 for (const auto& id_cell : cellsPtrs) {
0500 if (!pass(*id_cell.second, c2d))
0501 continue;
0502 float r = (id_cell.second->position().z() != 0.
0503 ? std::sqrt(pow(id_cell.second->position().x(), 2) + pow(id_cell.second->position().y(), 2)) /
0504 std::abs(id_cell.second->position().z())
0505 : 0.);
0506 tc_energy_r.emplace_back(std::make_pair(id_cell.second->energy(), r));
0507 }
0508
0509 float r_mean = meanX(tc_energy_r);
0510 float Srr = sigmaXX(tc_energy_r, r_mean);
0511
0512 return Srr;
0513 }
0514
0515 void HGCalShowerShape::fillShapes(l1t::HGCalMulticluster& c3d, const HGCalTriggerGeometryBase& triggerGeometry) const {
0516 c3d.showerLength(showerLength(c3d));
0517 c3d.coreShowerLength(coreShowerLength(c3d, triggerGeometry));
0518 c3d.firstLayer(firstLayer(c3d));
0519 c3d.maxLayer(maxLayer(c3d));
0520 c3d.sigmaEtaEtaTot(sigmaEtaEtaTot(c3d));
0521 c3d.sigmaEtaEtaMax(sigmaEtaEtaMax(c3d));
0522 c3d.sigmaPhiPhiTot(sigmaPhiPhiTot(c3d));
0523 c3d.sigmaPhiPhiMax(sigmaPhiPhiMax(c3d));
0524 c3d.sigmaZZ(sigmaZZ(c3d));
0525 c3d.sigmaRRTot(sigmaRRTot(c3d));
0526 c3d.sigmaRRMax(sigmaRRMax(c3d));
0527 c3d.sigmaRRMean(sigmaRRMean(c3d));
0528 c3d.eMax(eMax(c3d));
0529 c3d.zBarycenter(meanZ(c3d));
0530 c3d.layer10percent(percentileLayer(c3d, triggerGeometry, 0.10));
0531 c3d.layer50percent(percentileLayer(c3d, triggerGeometry, 0.50));
0532 c3d.layer90percent(percentileLayer(c3d, triggerGeometry, 0.90));
0533 c3d.triggerCells67percent(percentileTriggerCells(c3d, 0.67));
0534 c3d.triggerCells90percent(percentileTriggerCells(c3d, 0.90));
0535 }