51{
53
54 double clusterE = 0;
55 double clusterEta = 0;
56 double cluster_SIGNIFICANCE = 0;
57 double cluster_time = 0;
58 double cluster_SECOND_TIME = 0;
59 double cluster_CENTER_LAMBDA = 0;
60 double cluster_CENTER_MAG = 0;
61 double cluster_ENG_FRAC_EM_INCL = 0;
62 double cluster_FIRST_ENG_DENS = 0;
63 double cluster_LONGITUDINAL = 0;
64 double cluster_LATERAL = 0;
65 double cluster_PTD = 0;
66 double cluster_ISOLATION = 0;
67
68 std::vector<float> transformedFeatures;
69 bool ok{};
70
72 {
77 cluster_SECOND_TIME /= (Gaudi::Units::nanosecond * Gaudi::Units::nanosecond);
79 cluster_CENTER_LAMBDA /= Gaudi::Units::millimeter;
82 cluster_FIRST_ENG_DENS /= (Gaudi::Units::GeV / Gaudi::Units::millimeter3);
87 cluster_time = cluster->time() / Gaudi::Units::nanosecond;
88
89 if (!ok) {
91 return StatusCode::FAILURE;
92 }
93
94 float e_EM = 0.0;
95 for (
size_t s = CaloSampling::PreSamplerB;
s < CaloSampling::Unknown;
s++)
96 {
97 if (s == CaloSampling::EMB1 || s == CaloSampling::EMB2 || s == CaloSampling::EMB3 || s == CaloSampling::EME1 || s == CaloSampling::EME2 || s == CaloSampling::EME3 || s == CaloSampling::FCAL0)
98 {
100 }
101 }
102 cluster_ENG_FRAC_EM_INCL = e_EM / cluster->rawE();
103
104 std::vector<float> rawValues = {
105 static_cast<float>(clusterE),
106 static_cast<float>(clusterEta),
107 static_cast<float>(cluster_SIGNIFICANCE),
108 static_cast<float>(cluster_time),
109 static_cast<float>(cluster_SECOND_TIME),
110 static_cast<float>(cluster_CENTER_LAMBDA),
111 static_cast<float>(cluster_CENTER_MAG),
112 static_cast<float>(cluster_ENG_FRAC_EM_INCL),
113 static_cast<float>(cluster_FIRST_ENG_DENS),
114 static_cast<float>(cluster_LONGITUDINAL),
115 static_cast<float>(cluster_LATERAL),
116 static_cast<float>(cluster_PTD),
117 static_cast<float>(cluster_ISOLATION),
118 static_cast<float>(nPrimVtx),
119 avgMu
120 };
121
123 {
125 const float raw = rawValues.at(i);
127 transformedFeatures.push_back(transformed);
128 }
129 }
130
132 std::vector<int64_t> inputShape = {numClusters, 15};
133
135 inputData["features"] = std::make_pair(
136 inputShape, std::move(transformedFeatures));
137
139
140 outputData["mus"] = std::make_pair(
141 std::vector<int64_t>{numClusters, 3}, std::vector<float>{});
142 outputData["sigmas"] = std::make_pair(
143 std::vector<int64_t>{numClusters, 3}, std::vector<float>{});
144 outputData["alphas"] = std::make_pair(
145 std::vector<int64_t>{numClusters, 3}, std::vector<float>{});
146
148
149 std::vector<float> &onnx_mus = std::get<std::vector<float>>(outputData["mus"].second);
150 std::vector<float> &onnx_sigma2s = std::get<std::vector<float>>(outputData["sigmas"].second);
151 std::vector<float> &onnx_alphas = std::get<std::vector<float>>(outputData["alphas"].second);
152
153 if (msgLvl(MSG::DEBUG)) {
154 int nan_in_mus = 0;
155 int nan_in_sigma2s = 0;
156 int nan_in_alphas = 0;
157
158 for (float val : onnx_mus)
159 if (std::isnan(val))
160 nan_in_mus++;
161
162 for (float val : onnx_sigma2s)
163 if (std::isnan(val))
164 nan_in_sigma2s++;
165
166 for (float val : onnx_alphas)
167 if (std::isnan(val))
168 nan_in_alphas++;
169
170 if (nan_in_mus > 0)
171 ATH_MSG_DEBUG(nan_in_mus <<
" NaN value found in `mus` output layer during ONNX inference");
172
173 if (nan_in_sigma2s)
174 ATH_MSG_DEBUG(nan_in_sigma2s <<
" NaN value found in `sigmas` output layer during ONNX inference");
175
176 if (nan_in_alphas)
177 ATH_MSG_DEBUG(nan_in_alphas <<
" NaN value found in `alphas` output layer during ONNX inference");
178 }
179
180 clusterE_ML_vec.clear();
181 clusterE_ML_Unc_vec.clear();
182 clusterE_ML_vec.reserve(numClusters);
183 clusterE_ML_Unc_vec.reserve(numClusters);
184
185 for (
int i = 0;
i < numClusters; ++
i)
186 {
187 bool calibrateCluster = true;
188 for (size_t j=0; j<3; ++j) {
189 if (std::isnan(onnx_mus[i*3+j]) || std::isnan(onnx_sigma2s[i*3+j]) || std::isnan(onnx_alphas[i*3+j])) {
190 calibrateCluster = false;
191 break;
192 }
193 }
196 if (calibrateCluster) {
197 std::vector<float> current_mus = {onnx_mus[
i * 3], onnx_mus[
i * 3 + 1], onnx_mus[
i * 3 + 2]};
198 std::vector<float> current_sigma2s = {onnx_sigma2s[
i * 3], onnx_sigma2s[
i * 3 + 1], onnx_sigma2s[
i * 3 + 2]};
199 std::vector<float> current_alphas = {onnx_alphas[
i * 3], onnx_alphas[
i * 3 + 1], onnx_alphas[
i * 3 + 2]};
200
202 r = std::pow(10, mode);
204 s = std::abs(std::log(10) *
r) * onnx_s;
205
206 if (!std::isfinite(
r) || std::abs(
r) < 1e-6) {
207 ATH_MSG_WARNING(
"ML-correction factor to cluster energy (used as denominator) is " <<
r <<
"; The ML-correction factor is reset to 1. Uncertainty is set to 0.");
210 }
211 }
212
214
215 clusterE_ML_vec.push_back(cluster_energy);
216 clusterE_ML_Unc_vec.push_back(static_cast<double>(s));
217 }
218
219 return StatusCode::SUCCESS;
220}
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_WARNING(x)
@ PTD
relative spread of pT of constiuent cells = sqrt(n)*RMS/Mean
@ SECOND_TIME
Second moment of cell time distribution in cluster.
@ LATERAL
Normalized lateral moment.
@ LONGITUDINAL
Normalized longitudinal moment.
@ FIRST_ENG_DENS
First Moment in E/V.
@ CENTER_LAMBDA
Shower depth at Cluster Centroid.
@ SIGNIFICANCE
Cluster significance.
@ CENTER_MAG
Cluster Centroid ( )
@ ISOLATION
Energy weighted fraction of non-clustered perimeter cells.
CaloSampling::CaloSample CaloSample
Amg::Vector3D transform(Amg::Vector3D &v, Amg::Transform3D &tr)
Transform a point from a Trasformation3D.
std::map< std::string, InferenceData > OutputDataMap
std::map< std::string, InferenceData > InputDataMap
float modes(const std::vector< float > &mus, const std::vector< float > &log_sigma2s, const std::vector< float > &alphas)
float sigma_stoch(const std::vector< float > &mus, const std::vector< float > &log_sigma2s, const std::vector< float > &alphas)
CaloCluster_v1 CaloCluster
Define the latest version of the calorimeter cluster class.