Bug Summary

File:analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc
Warning:line 177, column 39
The left operand of '*' is a garbage value

Annotated Source Code

Press '?' to see keyboard shortcuts

clang -cc1 -cc1 -triple x86_64-unknown-linux-gnu -O3 -analyze -disable-free -clear-ast-before-backend -disable-llvm-verifier -discard-value-names -main-file-name SoftGaussParticleConstraint.cc -analyzer-checker=core -analyzer-checker=apiModeling -analyzer-checker=unix -analyzer-checker=deadcode -analyzer-checker=cplusplus -analyzer-checker=security.insecureAPI.UncheckedReturn -analyzer-checker=security.insecureAPI.getpw -analyzer-checker=security.insecureAPI.gets -analyzer-checker=security.insecureAPI.mktemp -analyzer-checker=security.insecureAPI.mkstemp -analyzer-checker=security.insecureAPI.vfork -analyzer-checker=nullability.NullPassedToNonnull -analyzer-checker=nullability.NullReturnedFromNonnull -analyzer-output plist -w -setup-static-analyzer -mrelocation-model pic -pic-level 2 -fhalf-no-semantic-interposition -mframe-pointer=none -fmath-errno -ffp-contract=on -fno-rounding-math -mconstructor-aliases -funwind-tables=2 -target-cpu x86-64 -tune-cpu generic -debugger-tuning=gdb -fdebug-compilation-dir=/data/b2soft/buildbot/development/build -fcoverage-compilation-dir=/data/b2soft/buildbot/development/build -resource-dir /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/lib/clang/21 -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/c++ -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/c++/x86_64-redhat-linux -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/c++/backward -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/include -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/python3.12 -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/include/CLHEP -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/Geant4 -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/include/root -isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/include/belle_legacy -I include/ -D _PACKAGE_="analysis" -D G4UI_USE_TCSH -D RaveDllExport= -D HAS_SQLITE -D HAS_CALLGRIND -I include -I /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/include/libxml2 -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/bin/../lib64/gcc/x86_64-redhat-linux/15.2.0/../../../../include/c++ -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/bin/../lib64/gcc/x86_64-redhat-linux/15.2.0/../../../../include/c++/x86_64-redhat-linux -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/bin/../lib64/gcc/x86_64-redhat-linux/15.2.0/../../../../include/c++/backward -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/lib/clang/21/include -internal-isystem /usr/local/include -internal-isystem /cvmfs/belle.cern.ch/el9/externals/v02-04-00/Linux_x86_64/common/bin/../lib64/gcc/x86_64-redhat-linux/15.2.0/../../../../x86_64-redhat-linux/include -internal-externc-isystem /include -internal-externc-isystem /usr/include -Wno-missing-braces -Wno-unused-command-line-argument -std=c++20 -fdeprecated-macro -ferror-limit 19 -fgnuc-version=4.2.1 -fno-implicit-modules -fskip-odr-check-in-gmf -fcxx-exceptions -fexceptions -vectorize-loops -vectorize-slp -analyzer-output=html -faddrsig -D__GCC_HAVE_DWARF2_CFI_ASM=1 -o /scan_build/2026-05-31-004316-385593-1 -x c++ analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * Forked from https://github.com/iLCSoft/MarlinKinfit *
6 * *
7 * Further information about the fit engine and the user interface *
8 * provided in MarlinKinfit can be found at *
9 * https://www.desy.de/~blist/kinfit/doc/html/ *
10 * and in the LCNotes LC-TOOL-2009-001 and LC-TOOL-2009-004 available *
11 * from http://www-flc.desy.de/lcnotes/ *
12 * *
13 * See git log for contributors and copyright holders. *
14 * This file is licensed under LGPL-3.0, see LICENSE.md. *
15 **************************************************************************/
16
17#include "analysis/OrcaKinFit/SoftGaussParticleConstraint.h"
18#include "analysis/OrcaKinFit/ParticleFitObject.h"
19#include <framework/logging/Logger.h>
20#include <iostream>
21#include <cmath>
22using namespace std;
23
24namespace Belle2 {
25 namespace OrcaKinFit {
26
27 SoftGaussParticleConstraint::SoftGaussParticleConstraint(double sigma_)
28 : sigma(sigma_)
29 {
30 invalidateCache();
31 }
32
33 double SoftGaussParticleConstraint::getSigma() const
34 {
35 return sigma;
36 }
37
38 double SoftGaussParticleConstraint::setSigma(double sigma_)
39 {
40 if (sigma_ != 0) sigma = std::abs(sigma_);
41 return sigma;
42 }
43
44 double SoftGaussParticleConstraint::getChi2() const
45 {
46 double r = getValue() / getSigma();
47 return r * r;
48 }
49
50 double SoftGaussParticleConstraint::getError() const
51 {
52 double dgdpi[4];
53 double error2 = 0;
54 for (unsigned int i = 0; i < fitobjects.size(); ++i) {
55 const ParticleFitObject* foi = fitobjects[i];
56 assert(foi)(static_cast <bool> (foi) ? void (0) : __assert_fail ("foi"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 56
, __extension__ __PRETTY_FUNCTION__))
;
57 if (firstDerivatives(i, dgdpi)) {
58 error2 += foi->getError2(dgdpi, getVarBasis());
59 }
60 }
61 return std::sqrt(std::abs(error2));
62 }
63
64 /**
65 * Calculates the second derivative of the constraint g w.r.t. the various parameters
66 * and adds it to the global covariance matrix
67 *
68 * We denote with P_i the 4-vector of the i-th ParticleFitObject,
69 * then
70 * \f[
71 * \frac{\partial ^2 g}{\partial a_k \partial a_l}
72 * = \sum_i \sum_j \frac{\partial ^2 g}{\partial P_i \partial P_j} \cdot
73 * \frac{\partial P_i}{\partial a_k} \cdot \frac{\partial P_j}{\partial a_l}
74 * + \sum_i \frac{\partial g}{\partial P_i} \cdot
75 * \frac{\partial^2 P_i}{\partial a_k \partial a_l}
76 * \f]
77 * Here, \f$\frac{\partial P_i}{\partial a_k}\f$ is a \f$4 \times n_i\f$ Matrix, where
78 * \f$n_i\f$ is the number of parameters of FitObject i;
79 * Correspondingly, \f$\frac{\partial^2 P_i}{\partial a_k \partial a_l}\f$ is a
80 * \f$4 \times n_i \times n_i\f$ matrix.
81 * Also, \f$\frac{\partial ^2 g}{\partial P_i \partial P_j}\f$ is a \f$4\times 4\f$ matrix
82 * for a given i and j, and \f$\frac{\partial g}{\partial P_i}\f$ is a 4-vector
83 * (though not a Lorentz-vector!).
84 *
85 */
86
87
88 void SoftGaussParticleConstraint::add2ndDerivativesToMatrix(double* M, int idim) const
89 {
90
91 /** First, treat the part
92 * \f[
93 * \frac{\partial ^2 g}{\partial P_i \partial P_j} \cdot
94 * \frac{\partial P_i}{\partial a_k} \cdot \frac{\partial P_j}{\partial a_l}
95 * \f]
96 */
97 double s = getSigma();
98 double fact = 2 * getValue() / (s * s);
99
100 // Derivatives \f$\frac{\partial ^2 g}{\partial P_i \partial P_j}\f$ at fixed i, j
101 // d2GdPidPj[4*ii+jj] is derivative w.r.t. P_i,ii and P_j,jj, where ii=0,1,2,3 for E,px,py,pz
102 double d2GdPidPj[16];
103 // Derivatives \f$\frac {\partial P_i}{\partial a_k}\f$ for all i;
104 // k is local parameter number
105 // dPidAk[KMAX*4*i + 4*k + ii] is \f$\frac {\partial P_{i,ii}}{\partial a_k}\f$,
106 // with ii=0, 1, 2, 3 for E, px, py, pz
107 const int KMAX = 4;
108 const int n = fitobjects.size();
109 auto* dPidAk = new double[n * KMAX * 4];
110 bool* dPidAkval = new bool[n];
111
112 for (int i = 0; i
0.1
'i' is < 'n'
< n
; ++i) dPidAkval[i] = false;
1
Loop condition is true. Entering loop body
2
Assuming 'i' is >= 'n'
3
Loop condition is false. Execution continues on line 116
113
114 // Derivatives \f$\frac{\partial ^2 g}{\partial P_i \partial a_l}\f$ at fixed i
115 // d2GdPdAl[4*l + ii] is \f$\frac{\partial ^2 g}{\partial P_{i,ii} \partial a_l}\f$
116 double d2GdPdAl[4 * KMAX];
117 // Derivatives \f$\frac{\partial ^2 g}{\partial a_k \partial a_l}\f$
118 double d2GdAkdAl[KMAX * KMAX] = {0};
119
120 // Global parameter numbers: parglobal[KMAX*i+klocal]
121 // is global parameter number of local parameter klocal of i-th Fit object
122 int* parglobal = new int[KMAX * n];
123
124 for (int i = 0; i < n; ++i) {
4
Loop condition is true. Entering loop body
9
Loop condition is false. Execution continues on line 133
125 const ParticleFitObject* foi = fitobjects[i];
126 assert(foi)(static_cast <bool> (foi) ? void (0) : __assert_fail ("foi"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 126
, __extension__ __PRETTY_FUNCTION__))
;
5
Assuming 'foi' is non-null
6
'?' condition is true
127 for (int klocal = 0; klocal < foi->getNPar(); ++klocal) {
7
Assuming the condition is false
8
Loop condition is false. Execution continues on line 124
128 parglobal [KMAX * i + klocal] = foi->getGlobalParNum(klocal);
129 }
130 }
131
132
133 for (int i = 0; i < n; ++i) {
10
Loop condition is true. Entering loop body
134 const ParticleFitObject* foi = fitobjects[i];
135 assert(foi)(static_cast <bool> (foi) ? void (0) : __assert_fail ("foi"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 135
, __extension__ __PRETTY_FUNCTION__))
;
11
Assuming 'foi' is non-null
12
'?' condition is true
136 for (int j = 0; j < n; ++j) {
13
Loop condition is true. Entering loop body
137 const ParticleFitObject* foj = fitobjects[j];
138 assert(foj)(static_cast <bool> (foj) ? void (0) : __assert_fail ("foj"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 138
, __extension__ __PRETTY_FUNCTION__))
;
14
Assuming 'foj' is non-null
15
'?' condition is true
139 if (secondDerivatives(i, j, d2GdPidPj)) {
16
Assuming the condition is true
17
Taking true branch
140 if (!dPidAkval[i]) {
18
Taking true branch
141 foi->getDerivatives(dPidAk + i * (KMAX * 4), KMAX * 4);
142 dPidAkval[i] = true;
143 }
144 if (!dPidAkval[j]) {
19
Taking false branch
145 foj->getDerivatives(dPidAk + j * (KMAX * 4), KMAX * 4);
146 dPidAkval[j] = true;
147 }
148 // Now sum over E/px/Py/Pz for object j:
149 // \f[
150 // \frac{\partial ^2 g}{\partial P_{i,ii} \partial a_l}
151 // = (sum_{j}) sum_{jj} frac{\partial ^2 g}{\partial P_{i,ii} \partial P_{j,jj}}
152 // \cdot \frac{\partial P_{j,jj}}{\partial a_l}
153 // \f]
154 // We're summing over jj here
155 for (int llocal = 0; llocal < foj->getNPar(); ++llocal) {
20
Assuming the condition is false
21
Loop condition is false. Execution continues on line 172
156 for (int ii = 0; ii < 4; ++ii) {
157 int ind1 = 4 * ii;
158 int ind2 = (KMAX * 4) * j + 4 * llocal;
159 double& r = d2GdPdAl[4 * llocal + ii];
160 r = d2GdPidPj[ ind1] * dPidAk[ ind2]; // E
161 r += d2GdPidPj[++ind1] * dPidAk[++ind2]; // px
162 r += d2GdPidPj[++ind1] * dPidAk[++ind2]; // py
163 r += d2GdPidPj[++ind1] * dPidAk[++ind2]; // pz
164 }
165 }
166 // Now sum over E/px/Py/Pz for object i, i.e. sum over ii:
167 // \f[
168 // \frac{\partial ^2 g}{\partial a_k \partial a_l}
169 // = \sum_{ii} \frac{\partial ^2 g}{\partial P_{i,ii} \partial a_l} \cdot
170 // \frac{\partial P_{i,ii}}{\partial a_k}
171 // \f]
172 for (int klocal = 0; klocal < foi->getNPar(); ++klocal) {
22
Assuming the condition is true
23
Loop condition is true. Entering loop body
173 for (int llocal = 0; llocal < foj->getNPar(); ++llocal) {
24
Assuming the condition is true
25
Loop condition is true. Entering loop body
174 int ind1 = 4 * llocal;
26
'ind1' initialized to 0
175 int ind2 = (KMAX * 4) * i + 4 * klocal;
176 double& r = d2GdAkdAl[KMAX * klocal + llocal];
177 r = d2GdPdAl[ ind1] * dPidAk[ ind2]; //E
27
The left operand of '*' is a garbage value
178 r += d2GdPdAl[++ind1] * dPidAk[++ind2]; // px
179 r += d2GdPdAl[++ind1] * dPidAk[++ind2]; // py
180 r += d2GdPdAl[++ind1] * dPidAk[++ind2]; // pz
181 }
182 }
183 // Now expand the local parameter numbers to global ones
184 for (int klocal = 0; klocal < foi->getNPar(); ++klocal) {
185 int kglobal = parglobal [KMAX * i + klocal];
186 for (int llocal = 0; llocal < foj->getNPar(); ++llocal) {
187 int lglobal = parglobal [KMAX * j + llocal];
188 M [idim * kglobal + lglobal] += fact * d2GdAkdAl[KMAX * klocal + llocal];
189 }
190 }
191 }
192 }
193 }
194 /** Second, treat the parts
195 * \f[
196 * \sum_i \frac{\partial g}{\partial P_i} \cdot
197 * \frac{\partial^2 P_i}{\partial a_k \partial a_l}
198 * \f]
199 * and
200 * \f[
201 * \frac{\partial^2 h}{\partial g^2}
202 * \sum_i \frac{\partial g}{\partial P_i} \cdot
203 * \frac{\partial P_i}{\partial a_k}
204 * \sum_j \frac{\partial g}{\partial P_j} \cdot
205 * \frac{\partial P_j}{\partial a_l}
206 * \f]
207 *
208 * Here, \f$\frac{\partial g}{\partial P_i}\f$ is a 4-vector, which we pass on to
209 * the FitObject
210 */
211
212 auto* v = new double[idim];
213 for (int i = 0; i < idim; ++i) v[i] = 0;
214 double sqrtfact2 = sqrt(2.0) / s;
215
216 double dgdpi[4];
217 for (int i = 0; i < n; ++i) {
218 const ParticleFitObject* foi = fitobjects[i];
219 assert(foi)(static_cast <bool> (foi) ? void (0) : __assert_fail ("foi"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 219
, __extension__ __PRETTY_FUNCTION__))
;
220 if (firstDerivatives(i, dgdpi)) {
221 foi->addTo2ndDerivatives(M, idim, fact, dgdpi, getVarBasis());
222 foi->addToGlobalChi2DerVector(v, idim, sqrtfact2, dgdpi, getVarBasis());
223 }
224 }
225
226 for (int i = 0; i < idim; ++i) {
227 if (double vi = v[i]) {
228 int ioffs = i * idim;
229 for (double* pvj = v; pvj < v + idim; ++pvj) {
230 M[ioffs++] += vi * (*pvj);
231 }
232 }
233 }
234
235
236 delete[] dPidAk;
237 delete[] dPidAkval;
238 delete[] parglobal;
239 delete[] v;
240 }
241
242 void SoftGaussParticleConstraint::addToGlobalChi2DerVector(double* y, int idim) const
243 {
244 double dgdpi[4];
245 double s = getSigma();
246 double r = 2 * getValue() / (s * s);
247 for (unsigned int i = 0; i < fitobjects.size(); ++i) {
248 const ParticleFitObject* foi = fitobjects[i];
249 assert(foi)(static_cast <bool> (foi) ? void (0) : __assert_fail ("foi"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 249
, __extension__ __PRETTY_FUNCTION__))
;
250 if (firstDerivatives(i, dgdpi)) {
251 foi->addToGlobalChi2DerVector(y, idim, r, dgdpi, getVarBasis());
252 }
253 }
254 }
255
256 void SoftGaussParticleConstraint::test1stDerivatives()
257 {
258 B2INFO("SoftGaussParticleConstraint::test1stDerivatives for " << getName())do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "SoftGaussParticleConstraint::test1stDerivatives for "
<< getName(); Belle2::LogSystem::Instance().sendMessage
(Belle2::LogMessage(Belle2::LogConfig::c_Info, std::move(varStream
), "analysis", __PRETTY_FUNCTION__, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc"
, 258, 0)); }; } } while(false)
;
259 double y[100];
260 for (double& i : y) i = 0;
261 addToGlobalChi2DerVector(y, 100);
262 double eps = 0.00001;
263 for (unsigned int ifo = 0; ifo < fitobjects.size(); ++ifo) {
264 ParticleFitObject* fo = fitobjects[ifo];
265 assert(fo)(static_cast <bool> (fo) ? void (0) : __assert_fail ("fo"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 265
, __extension__ __PRETTY_FUNCTION__))
;
266 for (int ilocal = 0; ilocal < fo->getNPar(); ++ilocal) {
267 int iglobal = fo->getGlobalParNum(ilocal);
268 double calc = y[iglobal];
269 double num = num1stDerivative(ifo, ilocal, eps);
270 B2INFO("fo: " << fo->getName() << " par " << ilocal << "/"do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "fo: " << fo->getName() <<
" par " << ilocal << "/" << iglobal <<
" (" << fo->getParamName(ilocal) << ") calc: "
<< calc << " - num: " << num << " = "
<< calc - num; Belle2::LogSystem::Instance().sendMessage
(Belle2::LogMessage(Belle2::LogConfig::c_Info, std::move(varStream
), "analysis", __PRETTY_FUNCTION__, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc"
, 272, 0)); }; } } while(false)
271 << iglobal << " (" << fo->getParamName(ilocal)do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "fo: " << fo->getName() <<
" par " << ilocal << "/" << iglobal <<
" (" << fo->getParamName(ilocal) << ") calc: "
<< calc << " - num: " << num << " = "
<< calc - num; Belle2::LogSystem::Instance().sendMessage
(Belle2::LogMessage(Belle2::LogConfig::c_Info, std::move(varStream
), "analysis", __PRETTY_FUNCTION__, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc"
, 272, 0)); }; } } while(false)
272 << ") calc: " << calc << " - num: " << num << " = " << calc - num)do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "fo: " << fo->getName() <<
" par " << ilocal << "/" << iglobal <<
" (" << fo->getParamName(ilocal) << ") calc: "
<< calc << " - num: " << num << " = "
<< calc - num; Belle2::LogSystem::Instance().sendMessage
(Belle2::LogMessage(Belle2::LogConfig::c_Info, std::move(varStream
), "analysis", __PRETTY_FUNCTION__, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc"
, 272, 0)); }; } } while(false)
;
273 }
274 }
275 }
276 void SoftGaussParticleConstraint::test2ndDerivatives()
277 {
278 B2INFO("SoftGaussParticleConstraint::test2ndDerivatives for " << getName())do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "SoftGaussParticleConstraint::test2ndDerivatives for "
<< getName(); Belle2::LogSystem::Instance().sendMessage
(Belle2::LogMessage(Belle2::LogConfig::c_Info, std::move(varStream
), "analysis", __PRETTY_FUNCTION__, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc"
, 278, 0)); }; } } while(false)
;
279 const int idim = 100;
280 auto* M = new double[idim * idim];
281 for (int i = 0; i < idim * idim; ++i) M[i] = 0;
282 add2ndDerivativesToMatrix(M, idim);
283 double eps = 0.0001;
284 B2INFO("eps=" << eps)do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "eps=" << eps; Belle2::LogSystem::
Instance().sendMessage(Belle2::LogMessage(Belle2::LogConfig::
c_Info, std::move(varStream), "analysis", __PRETTY_FUNCTION__
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 284
, 0)); }; } } while(false)
;
285
286 for (unsigned int ifo1 = 0; ifo1 < fitobjects.size(); ++ifo1) {
287 ParticleFitObject* fo1 = fitobjects[ifo1];
288 assert(fo1)(static_cast <bool> (fo1) ? void (0) : __assert_fail ("fo1"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 288
, __extension__ __PRETTY_FUNCTION__))
;
289 for (unsigned int ifo2 = ifo1; ifo2 < fitobjects.size(); ++ifo2) {
290 ParticleFitObject* fo2 = fitobjects[ifo2];
291 assert(fo2)(static_cast <bool> (fo2) ? void (0) : __assert_fail ("fo2"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 291
, __extension__ __PRETTY_FUNCTION__))
;
292 for (int ilocal1 = 0; ilocal1 < fo1->getNPar(); ++ilocal1) {
293 int iglobal1 = fo1->getGlobalParNum(ilocal1);
294 for (int ilocal2 = (ifo1 == ifo2 ? ilocal1 : 0); ilocal2 < fo2->getNPar(); ++ilocal2) {
295 int iglobal2 = fo2->getGlobalParNum(ilocal2);
296 double calc = M[idim * iglobal1 + iglobal2];
297 double num = num2ndDerivative(ifo1, ilocal1, eps, ifo2, ilocal2, eps);
298 B2INFO("fo1: " << fo1->getName() << " par " << ilocal1 << "/"do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "fo1: " << fo1->getName() <<
" par " << ilocal1 << "/" << iglobal1 <<
" (" << fo1->getParamName(ilocal1) << "), fo2: "
<< fo2->getName() << " par " << ilocal2
<< "/" << iglobal2 << " (" << fo2->
getParamName(ilocal2) << ") calc: " << calc <<
" - num: " << num << " = " << calc - num; Belle2
::LogSystem::Instance().sendMessage(Belle2::LogMessage(Belle2
::LogConfig::c_Info, std::move(varStream), "analysis", __PRETTY_FUNCTION__
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 302
, 0)); }; } } while(false)
299 << iglobal1 << " (" << fo1->getParamName(ilocal1)do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "fo1: " << fo1->getName() <<
" par " << ilocal1 << "/" << iglobal1 <<
" (" << fo1->getParamName(ilocal1) << "), fo2: "
<< fo2->getName() << " par " << ilocal2
<< "/" << iglobal2 << " (" << fo2->
getParamName(ilocal2) << ") calc: " << calc <<
" - num: " << num << " = " << calc - num; Belle2
::LogSystem::Instance().sendMessage(Belle2::LogMessage(Belle2
::LogConfig::c_Info, std::move(varStream), "analysis", __PRETTY_FUNCTION__
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 302
, 0)); }; } } while(false)
300 << "), fo2: " << fo2->getName() << " par " << ilocal2 << "/"do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "fo1: " << fo1->getName() <<
" par " << ilocal1 << "/" << iglobal1 <<
" (" << fo1->getParamName(ilocal1) << "), fo2: "
<< fo2->getName() << " par " << ilocal2
<< "/" << iglobal2 << " (" << fo2->
getParamName(ilocal2) << ") calc: " << calc <<
" - num: " << num << " = " << calc - num; Belle2
::LogSystem::Instance().sendMessage(Belle2::LogMessage(Belle2
::LogConfig::c_Info, std::move(varStream), "analysis", __PRETTY_FUNCTION__
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 302
, 0)); }; } } while(false)
301 << iglobal2 << " (" << fo2->getParamName(ilocal2)do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "fo1: " << fo1->getName() <<
" par " << ilocal1 << "/" << iglobal1 <<
" (" << fo1->getParamName(ilocal1) << "), fo2: "
<< fo2->getName() << " par " << ilocal2
<< "/" << iglobal2 << " (" << fo2->
getParamName(ilocal2) << ") calc: " << calc <<
" - num: " << num << " = " << calc - num; Belle2
::LogSystem::Instance().sendMessage(Belle2::LogMessage(Belle2
::LogConfig::c_Info, std::move(varStream), "analysis", __PRETTY_FUNCTION__
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 302
, 0)); }; } } while(false)
302 << ") calc: " << calc << " - num: " << num << " = " << calc - num)do { if (Belle2::LogSystem::Instance().isLevelEnabled(Belle2::
LogConfig::c_Info, 0, "analysis")) { { LogVariableStream varStream
; varStream << "fo1: " << fo1->getName() <<
" par " << ilocal1 << "/" << iglobal1 <<
" (" << fo1->getParamName(ilocal1) << "), fo2: "
<< fo2->getName() << " par " << ilocal2
<< "/" << iglobal2 << " (" << fo2->
getParamName(ilocal2) << ") calc: " << calc <<
" - num: " << num << " = " << calc - num; Belle2
::LogSystem::Instance().sendMessage(Belle2::LogMessage(Belle2
::LogConfig::c_Info, std::move(varStream), "analysis", __PRETTY_FUNCTION__
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 302
, 0)); }; } } while(false)
;
303 }
304 }
305 }
306 }
307 delete[] M;
308 }
309
310
311 double SoftGaussParticleConstraint::num1stDerivative(int ifo, int ilocal, double eps)
312 {
313 ParticleFitObject* fo = fitobjects[ifo];
314 assert(fo)(static_cast <bool> (fo) ? void (0) : __assert_fail ("fo"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 314
, __extension__ __PRETTY_FUNCTION__))
;
315 double save = fo->getParam(ilocal);
316 fo->setParam(ilocal, save + eps);
317 double v1 = getChi2();
318 fo->setParam(ilocal, save - eps);
319 double v2 = getChi2();
320 double result = (v1 - v2) / (2 * eps);
321 fo->setParam(ilocal, save);
322 return result;
323 }
324
325 double SoftGaussParticleConstraint::num2ndDerivative(int ifo1, int ilocal1, double eps1,
326 int ifo2, int ilocal2, double eps2)
327 {
328 double result;
329
330 if (ifo1 == ifo2 && ilocal1 == ilocal2) {
331 ParticleFitObject* fo = fitobjects[ifo1];
332 assert(fo)(static_cast <bool> (fo) ? void (0) : __assert_fail ("fo"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 332
, __extension__ __PRETTY_FUNCTION__))
;
333 double save = fo->getParam(ilocal1);
334 double v0 = getChi2();
335 fo->setParam(ilocal1, save + eps1);
336 double v1 = getChi2();
337 fo->setParam(ilocal1, save - eps1);
338 double v2 = getChi2();
339 result = (v1 + v2 - 2 * v0) / (eps1 * eps1);
340 fo->setParam(ilocal1, save);
341 } else {
342 ParticleFitObject* fo1 = fitobjects[ifo1];
343 assert(fo1)(static_cast <bool> (fo1) ? void (0) : __assert_fail ("fo1"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 343
, __extension__ __PRETTY_FUNCTION__))
;
344 ParticleFitObject* fo2 = fitobjects[ifo2];
345 assert(fo2)(static_cast <bool> (fo2) ? void (0) : __assert_fail ("fo2"
, "analysis/OrcaKinFit/src/SoftGaussParticleConstraint.cc", 345
, __extension__ __PRETTY_FUNCTION__))
;
346 double save1 = fo1->getParam(ilocal1);
347 double save2 = fo2->getParam(ilocal2);
348 fo1->setParam(ilocal1, save1 + eps1);
349 fo2->setParam(ilocal2, save2 + eps2);
350 double v11 = getChi2();
351 fo2->setParam(ilocal2, save2 - eps2);
352 double v12 = getChi2();
353 fo1->setParam(ilocal1, save1 - eps1);
354 double v22 = getChi2();
355 fo2->setParam(ilocal2, save2 + eps2);
356 double v21 = getChi2();
357 result = (v11 + v22 - v12 - v21) / (4 * eps1 * eps2);
358 fo1->setParam(ilocal1, save1);
359 fo2->setParam(ilocal2, save2);
360 }
361 return result;
362 }
363
364 int SoftGaussParticleConstraint::getVarBasis()
365 {
366 return VAR_BASIS;
367 }
368 }// end OrcaKinFit namespace
369} // end Belle2 namespace
370