Belle II Software development
Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | Related Functions | List of all members
SoftGaussMassConstraint Class Reference

Implements constraint 0 = mass1 - mass2 - m. More...

#include <SoftGaussMassConstraint.h>

Inheritance diagram for SoftGaussMassConstraint:
SoftGaussParticleConstraint BaseSoftConstraint BaseConstraint

Public Member Functions

 SoftGaussMassConstraint (double sigma_, double mass_=0.)
 Constructor.
 
virtual ~SoftGaussMassConstraint ()
 Virtual destructor.
 
virtual double getValue () const override
 Returns the value of the constraint function.
 
virtual void getDerivatives (int idim, double der[]) const override
 Get first order derivatives.
 
virtual double getMass (int flag=1)
 Get the actual invariant mass of the fit objects with a given flag.
 
virtual void setMass (double mass_)
 Sets the target mass of the constraint.
 
virtual void setFOList (std::vector< ParticleFitObject * > *fitobjects_)
 Adds several ParticleFitObject objects to the list.
 
virtual void addToFOList (ParticleFitObject &fitobject, int flag=1)
 Adds one ParticleFitObject objects to the list.
 
virtual void resetFOList ()
 Resests ParticleFitObject list.
 
virtual double getChi2 () const override
 Returns the chi2.
 
virtual double getError () const override
 Returns the error on the value of the constraint.
 
virtual double getSigma () const
 Returns the sigma.
 
virtual double setSigma (double sigma_)
 Sets the sigma.
 
virtual void add2ndDerivativesToMatrix (double *M, int idim) const override
 Adds second order derivatives to global covariance matrix M.
 
virtual void addToGlobalChi2DerVector (double *y, int idim) const override
 Add derivatives of chi squared to global derivative matrix.
 
void invalidateCache () const
 Invalidates any cached values for the next event.
 
void test1stDerivatives ()
 
void test2ndDerivatives ()
 
double num1stDerivative (int ifo, int ilocal, double eps)
 Evaluates numerically the 1st derivative w.r.t. a parameter.
 
double num2ndDerivative (int ifo1, int ilocal1, double eps1, int ifo2, int ilocal2, double eps2)
 Evaluates numerically the 2nd derivative w.r.t. 2 parameters.
 
int getVarBasis () const
 
virtual const char * getName () const
 Returns the name of the constraint.
 
void setName (const char *name_)
 Set object's name.
 
virtual std::ostream & print (std::ostream &os) const
 print object to ostream
 

Protected Types

enum  { VAR_BASIS = BaseDefs::VARBASIS_EPXYZ }
 
typedef std::vector< ParticleFitObject * > FitObjectContainer
 Vector of pointers to ParticleFitObjects.
 
typedef FitObjectContainer::iterator FitObjectIterator
 Iterator through vector of pointers to ParticleFitObjects.
 
typedef FitObjectContainer::const_iterator ConstFitObjectIterator
 Constant iterator through vector of pointers to ParticleFitObjects.
 

Protected Member Functions

virtual bool secondDerivatives (int i, int j, double *derivatives) const override
 Second derivatives with respect to the 4-vectors of Fit objects i and j; result false if all derivatives are zero.
 
virtual bool firstDerivatives (int i, double *derivatives) const override
 First derivatives with respect to the 4-vector of Fit objects i; result false if all derivatives are zero.
 

Protected Attributes

double mass
 The mass difference between object sets 1 and 2.
 
FitObjectContainer fitobjects
 The FitObjectContainer.
 
std::vector< double > derivatives
 The derivatives.
 
std::vector< int > flags
 The flags can be used to divide the FitObjectContainer into several subsets used for example to implement an equal mass constraint (see MassConstraint).
 
double sigma
 The sigma of the Gaussian.
 
char * name
 

Related Functions

(Note that these are not member functions.)

std::ostream & operator<< (std::ostream &os, const BaseConstraint &bc)
 Prints out a BaseConstraint, using its print method.
 

Detailed Description

Implements constraint 0 = mass1 - mass2 - m.

This class implements different mass constraints:

  • the invariant mass of several objects should be m
  • the difference of the invariant masses between two sets of objects should be m (normally m=0 in this case).

Author: Jenny List, Benno List Last update:

Date
2008/02/12 16:43:26

by:

Author
blist

Definition at line 45 of file SoftGaussMassConstraint.h.

Member Typedef Documentation

◆ ConstFitObjectIterator

typedef FitObjectContainer::const_iterator ConstFitObjectIterator
protectedinherited

Constant iterator through vector of pointers to ParticleFitObjects.

Definition at line 177 of file SoftGaussParticleConstraint.h.

◆ FitObjectContainer

typedef std::vector<ParticleFitObject*> FitObjectContainer
protectedinherited

Vector of pointers to ParticleFitObjects.

Definition at line 173 of file SoftGaussParticleConstraint.h.

◆ FitObjectIterator

typedef FitObjectContainer::iterator FitObjectIterator
protectedinherited

Iterator through vector of pointers to ParticleFitObjects.

Definition at line 175 of file SoftGaussParticleConstraint.h.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
protectedinherited

Definition at line 189 of file SoftGaussParticleConstraint.h.

189{ VAR_BASIS = BaseDefs::VARBASIS_EPXYZ }; // this means that the constraint knows about E,px,py,pz

Constructor & Destructor Documentation

◆ SoftGaussMassConstraint()

SoftGaussMassConstraint ( double  sigma_,
double  mass_ = 0. 
)
explicit

Constructor.

Parameters
sigma_The sigma value
mass_The mass difference between object sets 1 and 2

Definition at line 35 of file SoftGaussMassConstraint.cc.

36 : SoftGaussParticleConstraint(sigma_),
37 mass(mass_)
38 {}
Belle2::OrcaKinFit::SoftGaussMassConstraint::mass
double mass
The mass difference between object sets 1 and 2.
Definition: SoftGaussMassConstraint.h:74
Belle2::OrcaKinFit::SoftGaussParticleConstraint::SoftGaussParticleConstraint
SoftGaussParticleConstraint(double sigma_)
Creates an empty SoftGaussParticleConstraint object.
Definition: SoftGaussParticleConstraint.cc:31

Member Function Documentation

◆ add2ndDerivativesToMatrix()

void add2ndDerivativesToMatrix ( double *  M,
int  idim 
) const
overridevirtualinherited

Adds second order derivatives to global covariance matrix M.

Calculates the second derivative of the constraint g w.r.t.

the various parameters and adds it to the global covariance matrix

We denote with P_i the 4-vector of the i-th ParticleFitObject, then

\[ \frac{\partial ^2 g}{\partial a_k \partial a_l} = \sum_i \sum_j \frac{\partial ^2 g}{\partial P_i \partial P_j} \cdot \frac{\partial P_i}{\partial a_k} \cdot \frac{\partial P_j}{\partial a_l} + \sum_i \frac{\partial g}{\partial P_i} \cdot \frac{\partial^2 P_i}{\partial a_k \partial a_l} \]

Here, $\frac{\partial P_i}{\partial a_k}$ is a $4 \times n_i$ Matrix, where $n_i$ is the number of parameters of FitObject i; Correspondingly, $\frac{\partial^2 P_i}{\partial a_k \partial a_l}$ is a $4 \times n_i \times n_i$ matrix. Also, $\frac{\partial ^2 g}{\partial P_i \partial P_j}$ is a $4\times 4$ matrix for a given i and j, and $\frac{\partial g}{\partial P_i}$ is a 4-vector (though not a Lorentz-vector!).

First, treat the part

\[ \frac{\partial ^2 g}{\partial P_i \partial P_j} \cdot \frac{\partial P_i}{\partial a_k} \cdot \frac{\partial P_j}{\partial a_l} \]

Second, treat the parts

\[ \sum_i \frac{\partial g}{\partial P_i} \cdot \frac{\partial^2 P_i}{\partial a_k \partial a_l} \]

and

\[ \frac{\partial^2 h}{\partial g^2} \sum_i \frac{\partial g}{\partial P_i} \cdot \frac{\partial P_i}{\partial a_k} \sum_j \frac{\partial g}{\partial P_j} \cdot \frac{\partial P_j}{\partial a_l} \]

Here, $\frac{\partial g}{\partial P_i}$ is a 4-vector, which we pass on to the FitObject

Parameters
MCovariance matrix, at least idim x idim
idimFirst dimension of the array

Implements BaseSoftConstraint.

Definition at line 92 of file SoftGaussParticleConstraint.cc.

93 {
94
101 double s = getSigma();
102 double fact = 2 * getValue() / (s * s);
103
104 // Derivatives \f$\frac{\partial ^2 g}{\partial P_i \partial P_j}\f$ at fixed i, j
105 // 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
106 double d2GdPidPj[16];
107 // Derivatives \f$\frac {\partial P_i}{\partial a_k}\f$ for all i;
108 // k is local parameter number
109 // dPidAk[KMAX*4*i + 4*k + ii] is \f$\frac {\partial P_{i,ii}}{\partial a_k}\f$,
110 // with ii=0, 1, 2, 3 for E, px, py, pz
111 const int KMAX = 4;
112 const int n = fitobjects.size();
113 auto* dPidAk = new double[n * KMAX * 4];
114 bool* dPidAkval = new bool[n];
115
116 for (int i = 0; i < n; ++i) dPidAkval[i] = false;
117
118 // Derivatives \f$\frac{\partial ^2 g}{\partial P_i \partial a_l}\f$ at fixed i
119 // d2GdPdAl[4*l + ii] is \f$\frac{\partial ^2 g}{\partial P_{i,ii} \partial a_l}\f$
120 double d2GdPdAl[4 * KMAX];
121 // Derivatives \f$\frac{\partial ^2 g}{\partial a_k \partial a_l}\f$
122 double d2GdAkdAl[KMAX * KMAX] = {0};
123
124 // Global parameter numbers: parglobal[KMAX*i+klocal]
125 // is global parameter number of local parameter klocal of i-th Fit object
126 int* parglobal = new int[KMAX * n];
127
128 for (int i = 0; i < n; ++i) {
129 const ParticleFitObject* foi = fitobjects[i];
130 assert(foi);
131 for (int klocal = 0; klocal < foi->getNPar(); ++klocal) {
132 parglobal [KMAX * i + klocal] = foi->getGlobalParNum(klocal);
133 }
134 }
135
136
137 for (int i = 0; i < n; ++i) {
138 const ParticleFitObject* foi = fitobjects[i];
139 assert(foi);
140 for (int j = 0; j < n; ++j) {
141 const ParticleFitObject* foj = fitobjects[j];
142 assert(foj);
143 if (secondDerivatives(i, j, d2GdPidPj)) {
144 if (!dPidAkval[i]) {
145 foi->getDerivatives(dPidAk + i * (KMAX * 4), KMAX * 4);
146 dPidAkval[i] = true;
147 }
148 if (!dPidAkval[j]) {
149 foj->getDerivatives(dPidAk + j * (KMAX * 4), KMAX * 4);
150 dPidAkval[j] = true;
151 }
152 // Now sum over E/px/Py/Pz for object j:
153 // \f[
154 // \frac{\partial ^2 g}{\partial P_{i,ii} \partial a_l}
155 // = (sum_{j}) sum_{jj} frac{\partial ^2 g}{\partial P_{i,ii} \partial P_{j,jj}}
156 // \cdot \frac{\partial P_{j,jj}}{\partial a_l}
157 // \f]
158 // We're summing over jj here
159 for (int llocal = 0; llocal < foj->getNPar(); ++llocal) {
160 for (int ii = 0; ii < 4; ++ii) {
161 int ind1 = 4 * ii;
162 int ind2 = (KMAX * 4) * j + 4 * llocal;
163 double& r = d2GdPdAl[4 * llocal + ii];
164 r = d2GdPidPj[ ind1] * dPidAk[ ind2]; // E
165 r += d2GdPidPj[++ind1] * dPidAk[++ind2]; // px
166 r += d2GdPidPj[++ind1] * dPidAk[++ind2]; // py
167 r += d2GdPidPj[++ind1] * dPidAk[++ind2]; // pz
168 }
169 }
170 // Now sum over E/px/Py/Pz for object i, i.e. sum over ii:
171 // \f[
172 // \frac{\partial ^2 g}{\partial a_k \partial a_l}
173 // = \sum_{ii} \frac{\partial ^2 g}{\partial P_{i,ii} \partial a_l} \cdot
174 // \frac{\partial P_{i,ii}}{\partial a_k}
175 // \f]
176 for (int klocal = 0; klocal < foi->getNPar(); ++klocal) {
177 for (int llocal = 0; llocal < foj->getNPar(); ++llocal) {
178 int ind1 = 4 * llocal;
179 int ind2 = (KMAX * 4) * i + 4 * klocal;
180 double& r = d2GdAkdAl[KMAX * klocal + llocal];
181 r = d2GdPdAl[ ind1] * dPidAk[ ind2]; //E
182 r += d2GdPdAl[++ind1] * dPidAk[++ind2]; // px
183 r += d2GdPdAl[++ind1] * dPidAk[++ind2]; // py
184 r += d2GdPdAl[++ind1] * dPidAk[++ind2]; // pz
185 }
186 }
187 // Now expand the local parameter numbers to global ones
188 for (int klocal = 0; klocal < foi->getNPar(); ++klocal) {
189 int kglobal = parglobal [KMAX * i + klocal];
190 for (int llocal = 0; llocal < foj->getNPar(); ++llocal) {
191 int lglobal = parglobal [KMAX * j + llocal];
192 M [idim * kglobal + lglobal] += fact * d2GdAkdAl[KMAX * klocal + llocal];
193 }
194 }
195 }
196 }
197 }
216 auto* v = new double[idim];
217 for (int i = 0; i < idim; ++i) v[i] = 0;
218 double sqrtfact2 = sqrt(2.0) / s;
219
220 double dgdpi[4];
221 for (int i = 0; i < n; ++i) {
222 const ParticleFitObject* foi = fitobjects[i];
223 assert(foi);
224 if (firstDerivatives(i, dgdpi)) {
225 foi->addTo2ndDerivatives(M, idim, fact, dgdpi, getVarBasis());
226 foi->addToGlobalChi2DerVector(v, idim, sqrtfact2, dgdpi, getVarBasis());
227 }
228 }
229
230 for (int i = 0; i < idim; ++i) {
231 if (double vi = v[i]) {
232 int ioffs = i * idim;
233 for (double* pvj = v; pvj < v + idim; ++pvj) {
234 M[ioffs++] += vi * (*pvj);
235 }
236 }
237 }
238
239
240 delete[] dPidAk;
241 delete[] dPidAkval;
242 delete[] parglobal;
243 delete[] v;
244 }
Belle2::OrcaKinFit::SoftGaussParticleConstraint::getValue
virtual double getValue() const override=0
Returns the value of the constraint function.
Belle2::OrcaKinFit::SoftGaussParticleConstraint::fitobjects
FitObjectContainer fitobjects
The FitObjectContainer.
Definition: SoftGaussParticleConstraint.h:179
Belle2::OrcaKinFit::SoftGaussParticleConstraint::secondDerivatives
virtual bool secondDerivatives(int i, int j, double *derivatives) const =0
Second derivatives with respect to the 4-vectors of Fit objects i and j; result false if all derivati...
Belle2::OrcaKinFit::SoftGaussParticleConstraint::firstDerivatives
virtual bool firstDerivatives(int i, double *derivatives) const =0
First derivatives with respect to the 4-vector of Fit objects i; result false if all derivatives are ...
Belle2::OrcaKinFit::SoftGaussParticleConstraint::getSigma
virtual double getSigma() const
Returns the sigma.
Definition: SoftGaussParticleConstraint.cc:37
Belle2::sqrt
double sqrt(double a)
sqrt for double
Definition: beamHelpers.h:28

◆ addToFOList()

virtual void addToFOList ( ParticleFitObject &  fitobject,
int  flag = 1 
)
inlinevirtualinherited

Adds one ParticleFitObject objects to the list.

Definition at line 92 of file SoftGaussParticleConstraint.h.

94 {
95 fitobjects.push_back(&fitobject);
96 flags.push_back(flag);
97 };
Belle2::OrcaKinFit::SoftGaussParticleConstraint::flags
std::vector< int > flags
The flags can be used to divide the FitObjectContainer into several subsets used for example to imple...
Definition: SoftGaussParticleConstraint.h:184

◆ addToGlobalChi2DerVector()

void addToGlobalChi2DerVector ( double *  y,
int  idim 
) const
overridevirtualinherited

Add derivatives of chi squared to global derivative matrix.

Parameters
yVector of chi2 derivatives
idimVector size

Implements BaseSoftConstraint.

Definition at line 246 of file SoftGaussParticleConstraint.cc.

247 {
248 double dgdpi[4];
249 double s = getSigma();
250 double r = 2 * getValue() / (s * s);
251 for (unsigned int i = 0; i < fitobjects.size(); ++i) {
252 const ParticleFitObject* foi = fitobjects[i];
253 assert(foi);
254 if (firstDerivatives(i, dgdpi)) {
255 foi->addToGlobalChi2DerVector(y, idim, r, dgdpi, getVarBasis());
256 }
257 }
258 }

◆ firstDerivatives()

bool firstDerivatives ( int  i,
double *  derivatives 
) const
overrideprotectedvirtual

First derivatives with respect to the 4-vector of Fit objects i; result false if all derivatives are zero.

Parameters
inumber of 1st FitObject
derivativesThe result 4-vector

Implements SoftGaussParticleConstraint.

Definition at line 191 of file SoftGaussMassConstraint.cc.

192 {
193 double totE = 0;
194 double totpx = 0;
195 double totpy = 0;
196 double totpz = 0;
197 int index = (flags[i] == 1) ? 0 : 1; // default is 1, but 2 may indicate fitobjects for a second W -> equal mass constraint!
198 for (unsigned int j = 0; j < fitobjects.size(); ++j) {
199 int jndex = (flags[j] == 1) ? 0 : 1; // default is 1, but 2 may indicate fitobjects for a second W -> equal mass constraint!
200 const ParticleFitObject* foj = fitobjects[j];
201 assert(foj);
202 if (index == jndex) {
203 totE += foj->getE();
204 totpx += foj->getPx();
205 totpy += foj->getPy();
206 totpz += foj->getPz();
207 }
208 }
209
210 if (totE <= 0) {
211 B2ERROR("SoftGaussMassConstraint::firstDerivatives: totE = " << totE);
212 }
213
214 double m = std::sqrt(std::abs(totE * totE - totpx * totpx - totpy * totpy - totpz * totpz));
215 if (index) m = -m;
216
217 dderivatives[0] = totE / m;
218 dderivatives[1] = -totpx / m;
219 dderivatives[2] = -totpy / m;
220 dderivatives[3] = -totpz / m;
221 return true;
222 }

◆ getChi2()

double getChi2 ( ) const
overridevirtualinherited

Returns the chi2.

Implements BaseSoftConstraint.

Definition at line 48 of file SoftGaussParticleConstraint.cc.

49 {
50 double r = getValue() / getSigma();
51 return r * r;
52 }

◆ getDerivatives()

void getDerivatives ( int  idim,
double  der[] 
) const
overridevirtual

Get first order derivatives.

Call this with a predefined array "der" with the necessary number of entries!

Parameters
idimFirst dimension of the array
derArray of derivatives, at least idim x idim

Implements SoftGaussParticleConstraint.

Definition at line 68 of file SoftGaussMassConstraint.cc.

69 {
70 double totE[2] = {0, 0};
71 double totpx[2] = {0, 0};
72 double totpy[2] = {0, 0};
73 double totpz[2] = {0, 0};
74 bool valid[2] = {false, false};
75 for (unsigned int i = 0; i < fitobjects.size(); i++) {
76 int index = (flags[i] == 1) ? 0 : 1; // default is 1, but 2 may indicate fitobjects for a second W -> equal mass constraint!
77 valid[index] = true;
78 totE[index] += fitobjects[i]->getE();
79 totpx[index] += fitobjects[i]->getPx();
80 totpy[index] += fitobjects[i]->getPy();
81 totpz[index] += fitobjects[i]->getPz();
82 }
83 double m2[2];
84 double m_inv[2] = {0, 0};
85 for (int index = 0; index < 2; ++index) {
86 m2[index] = totE[index] * totE[index] - totpx[index] * totpx[index]
87 - totpy[index] * totpy[index] - totpz[index] * totpz[index];
88 if (m2[index] < 0 && m2[index] > -1E-9) m2[index] = 0;
89 if (m2[index] < 0 && valid[index]) {
90 B2ERROR("SoftGaussMassConstraint::getDerivatives: m2<0!");
91 for (unsigned int j = 0; j < fitobjects.size(); j++) {
92 int jndex = (flags[j] == 1) ? 0 : 1;
93 if (jndex == index) {
94 B2ERROR(fitobjects[j]->getName() << ": E=" << fitobjects[j]->getE() << ", px=" << fitobjects[j]->getPx()
95 << ", py=" << fitobjects[j]->getPy() << ", pz=" << fitobjects[j]->getPz());
96 }
97 }
98 B2ERROR("sum: E=" << totE[index] << ", px=" << totpx[index]
99 << ", py=" << totpy[index] << ", pz=" << totpz[index] << ", m2=" << m2[index]);
100 }
101 if (m2[index] != 0) m_inv[index] = 1 / std::sqrt(std::abs(m2[index]));
102 }
103
104 for (unsigned int i = 0; i < fitobjects.size(); i++) {
105 int index = (flags[i] == 1) ? 0 : 1; // default is 1, but 2 may indicate fitobjects for a second W -> equal mass constraint!
106 for (int ilocal = 0; ilocal < fitobjects[i]->getNPar(); ilocal++) {
107 if (!fitobjects[i]->isParamFixed(ilocal)) {
108 int iglobal = fitobjects[i]->getGlobalParNum(ilocal);
109 assert(iglobal >= 0 && iglobal < idim);
110 if (m2[index] != 0) {
111 der[iglobal] = totE[index] * fitobjects[i]->getDE(ilocal)
112 - totpx[index] * fitobjects[i]->getDPx(ilocal)
113 - totpy[index] * fitobjects[i]->getDPy(ilocal)
114 - totpz[index] * fitobjects[i]->getDPz(ilocal);
115 der[iglobal] *= m_inv[index];
116 } else der[iglobal] = 1;
117 if (index == 1) der[iglobal] *= -1.;
118 }
119 }
120 }
121 }
E
R E
internal precision of FFTW codelets
Definition: DiscreteCosineTransform_31points.cc:36
Belle2::OrcaKinFit::BaseConstraint::getName
virtual const char * getName() const
Returns the name of the constraint.
Definition: BaseConstraint.cc:56

◆ getError()

double getError ( ) const
overridevirtualinherited

Returns the error on the value of the constraint.

Reimplemented from BaseConstraint.

Definition at line 54 of file SoftGaussParticleConstraint.cc.

55 {
56 double dgdpi[4];
57 double error2 = 0;
58 for (unsigned int i = 0; i < fitobjects.size(); ++i) {
59 const ParticleFitObject* foi = fitobjects[i];
60 assert(foi);
61 if (firstDerivatives(i, dgdpi)) {
62 error2 += foi->getError2(dgdpi, getVarBasis());
63 }
64 }
65 return std::sqrt(std::abs(error2));
66 }

◆ getMass()

double getMass ( int  flag = 1)
virtual

Get the actual invariant mass of the fit objects with a given flag.

Parameters
flagThe flag

Definition at line 123 of file SoftGaussMassConstraint.cc.

124 {
125 double totE = 0;
126 double totpx = 0;
127 double totpy = 0;
128 double totpz = 0;
129 for (unsigned int i = 0; i < fitobjects.size(); i++) {
130 if (flags[i] == flag) {
131 totE += fitobjects[i]->getE();
132 totpx += fitobjects[i]->getPx();
133 totpy += fitobjects[i]->getPy();
134 totpz += fitobjects[i]->getPz();
135 }
136 }
137 return std::sqrt(std::abs(totE * totE - totpx * totpx - totpy * totpy - totpz * totpz));
138 }

◆ getName()

const char * getName ( ) const
virtualinherited

Returns the name of the constraint.

Definition at line 56 of file BaseConstraint.cc.

57 {
58 return name;
59 }

◆ getSigma()

double getSigma ( ) const
virtualinherited

Returns the sigma.

Definition at line 37 of file SoftGaussParticleConstraint.cc.

38 {
39 return sigma;
40 }

◆ getValue()

double getValue ( ) const
overridevirtual

Returns the value of the constraint function.

Implements SoftGaussParticleConstraint.

Definition at line 44 of file SoftGaussMassConstraint.cc.

45 {
46 double totE[2] = {0, 0};
47 double totpx[2] = {0, 0};
48 double totpy[2] = {0, 0};
49 double totpz[2] = {0, 0};
50 for (unsigned int i = 0; i < fitobjects.size(); i++) {
51 int index = (flags[i] == 1) ? 0 : 1; // default is 1, but 2 may indicate fitobjects for a second W -> equal mass constraint!
52 totE[index] += fitobjects[i]->getE();
53 totpx[index] += fitobjects[i]->getPx();
54 totpy[index] += fitobjects[i]->getPy();
55 totpz[index] += fitobjects[i]->getPz();
56 }
57 double result = -mass;
58 result += std::sqrt(std::abs(totE[0] * totE[0] - totpx[0] * totpx[0] - totpy[0] * totpy[0] - totpz[0] * totpz[0]));
59 result -= std::sqrt(std::abs(totE[1] * totE[1] - totpx[1] * totpx[1] - totpy[1] * totpy[1] - totpz[1] * totpz[1]));
60 return result;
61 }

◆ getVarBasis()

int getVarBasis ( ) const
inherited

Definition at line 368 of file SoftGaussParticleConstraint.cc.

369 {
370 return VAR_BASIS;
371 }

◆ invalidateCache()

void invalidateCache ( ) const
inlineinherited

Invalidates any cached values for the next event.

Definition at line 138 of file SoftGaussParticleConstraint.h.

138{}

◆ num1stDerivative()

double num1stDerivative ( int  ifo,
int  ilocal,
double  eps 
)
inherited

Evaluates numerically the 1st derivative w.r.t. a parameter.

Parameters
ifoNumber of FitObject
ilocalLocal parameter number
epsvariation of local parameter

Definition at line 315 of file SoftGaussParticleConstraint.cc.

316 {
317 ParticleFitObject* fo = fitobjects[ifo];
318 assert(fo);
319 double save = fo->getParam(ilocal);
320 fo->setParam(ilocal, save + eps);
321 double v1 = getChi2();
322 fo->setParam(ilocal, save - eps);
323 double v2 = getChi2();
324 double result = (v1 - v2) / (2 * eps);
325 fo->setParam(ilocal, save);
326 return result;
327 }
Belle2::OrcaKinFit::SoftGaussParticleConstraint::getChi2
virtual double getChi2() const override
Returns the chi2.
Definition: SoftGaussParticleConstraint.cc:48
VXDTFFilterTest::v2
const std::vector< double > v2
MATLAB generated random vector.
Definition: decorrelationMatrix.cc:42
VXDTFFilterTest::v1
const std::vector< double > v1
MATLAB generated random vector.
Definition: decorrelationMatrix.cc:28

◆ num2ndDerivative()

double num2ndDerivative ( int  ifo1,
int  ilocal1,
double  eps1,
int  ifo2,
int  ilocal2,
double  eps2 
)
inherited

Evaluates numerically the 2nd derivative w.r.t. 2 parameters.

Parameters
ifo1Number of 1st FitObject
ilocal11st local parameter number
eps1variation of 1st local parameter
ifo2Number of 1st FitObject
ilocal21st local parameter number
eps2variation of 2nd local parameter

Definition at line 329 of file SoftGaussParticleConstraint.cc.

331 {
332 double result;
333
334 if (ifo1 == ifo2 && ilocal1 == ilocal2) {
335 ParticleFitObject* fo = fitobjects[ifo1];
336 assert(fo);
337 double save = fo->getParam(ilocal1);
338 double v0 = getChi2();
339 fo->setParam(ilocal1, save + eps1);
340 double v1 = getChi2();
341 fo->setParam(ilocal1, save - eps1);
342 double v2 = getChi2();
343 result = (v1 + v2 - 2 * v0) / (eps1 * eps1);
344 fo->setParam(ilocal1, save);
345 } else {
346 ParticleFitObject* fo1 = fitobjects[ifo1];
347 assert(fo1);
348 ParticleFitObject* fo2 = fitobjects[ifo2];
349 assert(fo2);
350 double save1 = fo1->getParam(ilocal1);
351 double save2 = fo2->getParam(ilocal2);
352 fo1->setParam(ilocal1, save1 + eps1);
353 fo2->setParam(ilocal2, save2 + eps2);
354 double v11 = getChi2();
355 fo2->setParam(ilocal2, save2 - eps2);
356 double v12 = getChi2();
357 fo1->setParam(ilocal1, save1 - eps1);
358 double v22 = getChi2();
359 fo2->setParam(ilocal2, save2 + eps2);
360 double v21 = getChi2();
361 result = (v11 + v22 - v12 - v21) / (4 * eps1 * eps2);
362 fo1->setParam(ilocal1, save1);
363 fo2->setParam(ilocal2, save2);
364 }
365 return result;
366 }

◆ print()

std::ostream & print ( std::ostream &  os) const
virtualinherited

print object to ostream

Parameters
osThe output stream

Definition at line 76 of file BaseConstraint.cc.

77 {
78 os << getName() << "=" << getValue();
79 return os;
80 }
Belle2::OrcaKinFit::BaseConstraint::getValue
virtual double getValue() const =0
Returns the value of the constraint function.

◆ resetFOList()

virtual void resetFOList ( )
inlinevirtualinherited

Resests ParticleFitObject list.

Definition at line 99 of file SoftGaussParticleConstraint.h.

100 {
101 fitobjects.resize(0);
102 flags.resize(0);
103 };

◆ secondDerivatives()

bool secondDerivatives ( int  i,
int  j,
double *  derivatives 
) const
overrideprotectedvirtual

Second derivatives with respect to the 4-vectors of Fit objects i and j; result false if all derivatives are zero.

Parameters
inumber of 1st FitObject
jnumber of 2nd FitObject
derivativesThe result 4x4 matrix

Implements SoftGaussParticleConstraint.

Definition at line 145 of file SoftGaussMassConstraint.cc.

146 {
147 B2DEBUG(14, "SoftGaussMassConstraint::secondDerivatives: i=" << i << ", j=" << j);
148 int index = (flags[i] == 1) ? 0 : 1; // default is 1, but 2 may indicate fitobjects for a second W -> equal mass constraint!
149 int jndex = (flags[j] == 1) ? 0 : 1; // default is 1, but 2 may indicate fitobjects for a second W -> equal mass constraint!
150 if (index != jndex) return false;
151 double totE = 0;
152 double totpx = 0;
153 double totpy = 0;
154 double totpz = 0;
155 for (unsigned int k = 0; k < fitobjects.size(); ++k) {
156 int kndex = (flags[k] == 1) ? 0 : 1; // default is 1, but 2 may indicate fitobjects for a second W -> equal mass constraint!
157 const ParticleFitObject* fok = fitobjects[k];
158 assert(fok);
159 if (index == kndex) {
160 totE += fok->getE();
161 totpx += fok->getPx();
162 totpy += fok->getPy();
163 totpz += fok->getPz();
164 }
165 }
166
167 if (totE <= 0) {
168 B2ERROR("SoftGaussMassConstraint::secondDerivatives: totE = " << totE);
169 }
170
171 double m2 = std::abs(totE * totE - totpx * totpx - totpy * totpy - totpz * totpz);
172 double m = std::sqrt(m2);
173 if (index) m = -m;
174 double minv3 = 1 / (m * m * m);
175
176 assert(dderivatives);
177 for (int k = 0; k < 16; ++k) dderivatives[k] = 0;
178 dderivatives[4 * 0 + 0] = (m2 - totE * totE) * minv3;
179 dderivatives[4 * 0 + 1] = dderivatives[4 * 1 + 0] = totE * totpx * minv3;
180 dderivatives[4 * 0 + 2] = dderivatives[4 * 2 + 0] = totE * totpy * minv3;
181 dderivatives[4 * 0 + 3] = dderivatives[4 * 3 + 0] = totE * totpz * minv3;
182 dderivatives[4 * 1 + 1] = -(m2 + totpx * totpx) * minv3;
183 dderivatives[4 * 1 + 2] = dderivatives[4 * 2 + 1] = -totpx * totpy * minv3;
184 dderivatives[4 * 1 + 3] = dderivatives[4 * 3 + 1] = -totpx * totpz * minv3;
185 dderivatives[4 * 2 + 2] = -(m2 + totpy * totpy) * minv3;
186 dderivatives[4 * 2 + 3] = dderivatives[4 * 3 + 2] = -totpy * totpz * minv3;
187 dderivatives[4 * 3 + 3] = -(m2 + totpz * totpz) * minv3;
188 return true;
189 }

◆ setFOList()

virtual void setFOList ( std::vector< ParticleFitObject * > *  fitobjects_)
inlinevirtualinherited

Adds several ParticleFitObject objects to the list.

Parameters
fitobjects_A list of BaseFitObject objects

Definition at line 83 of file SoftGaussParticleConstraint.h.

85 {
86 for (int i = 0; i < (int) fitobjects_->size(); i++) {
87 fitobjects.push_back((*fitobjects_)[i]);
88 flags.push_back(1);
89 }
90 };

◆ setMass()

void setMass ( double  mass_)
virtual

Sets the target mass of the constraint.

Parameters
mass_The new mass

Definition at line 140 of file SoftGaussMassConstraint.cc.

141 {
142 mass = mass_;
143 }

◆ setName()

void setName ( const char *  name_)
inherited

Set object's name.

Definition at line 61 of file BaseConstraint.cc.

62 {
63 if (name_ == nullptr) return;
64 size_t l = strlen(name_);
65 if (name) delete[] name;
66 name = new char[l + 1];
67 strcpy(name, name_);
68 }

◆ setSigma()

double setSigma ( double  sigma_)
virtualinherited

Sets the sigma.

Parameters
sigma_The new sigma value

Definition at line 42 of file SoftGaussParticleConstraint.cc.

43 {
44 if (sigma_ != 0) sigma = std::abs(sigma_);
45 return sigma;
46 }

◆ test1stDerivatives()

void test1stDerivatives ( )
inherited

Definition at line 260 of file SoftGaussParticleConstraint.cc.

261 {
262 B2INFO("SoftGaussParticleConstraint::test1stDerivatives for " << getName());
263 double y[100];
264 for (double& i : y) i = 0;
265 addToGlobalChi2DerVector(y, 100);
266 double eps = 0.00001;
267 for (unsigned int ifo = 0; ifo < fitobjects.size(); ++ifo) {
268 ParticleFitObject* fo = fitobjects[ifo];
269 assert(fo);
270 for (int ilocal = 0; ilocal < fo->getNPar(); ++ilocal) {
271 int iglobal = fo->getGlobalParNum(ilocal);
272 double calc = y[iglobal];
273 double num = num1stDerivative(ifo, ilocal, eps);
274 B2INFO("fo: " << fo->getName() << " par " << ilocal << "/"
275 << iglobal << " (" << fo->getParamName(ilocal)
276 << ") calc: " << calc << " - num: " << num << " = " << calc - num);
277 }
278 }
279 }
Belle2::OrcaKinFit::SoftGaussParticleConstraint::addToGlobalChi2DerVector
virtual void addToGlobalChi2DerVector(double *y, int idim) const override
Add derivatives of chi squared to global derivative matrix.
Definition: SoftGaussParticleConstraint.cc:246
Belle2::OrcaKinFit::SoftGaussParticleConstraint::num1stDerivative
double num1stDerivative(int ifo, int ilocal, double eps)
Evaluates numerically the 1st derivative w.r.t. a parameter.
Definition: SoftGaussParticleConstraint.cc:315

◆ test2ndDerivatives()

void test2ndDerivatives ( )
inherited

Definition at line 280 of file SoftGaussParticleConstraint.cc.

281 {
282 B2INFO("SoftGaussParticleConstraint::test2ndDerivatives for " << getName());
283 const int idim = 100;
284 auto* M = new double[idim * idim];
285 for (int i = 0; i < idim * idim; ++i) M[i] = 0;
286 add2ndDerivativesToMatrix(M, idim);
287 double eps = 0.0001;
288 B2INFO("eps=" << eps);
289
290 for (unsigned int ifo1 = 0; ifo1 < fitobjects.size(); ++ifo1) {
291 ParticleFitObject* fo1 = fitobjects[ifo1];
292 assert(fo1);
293 for (unsigned int ifo2 = ifo1; ifo2 < fitobjects.size(); ++ifo2) {
294 ParticleFitObject* fo2 = fitobjects[ifo2];
295 assert(fo2);
296 for (int ilocal1 = 0; ilocal1 < fo1->getNPar(); ++ilocal1) {
297 int iglobal1 = fo1->getGlobalParNum(ilocal1);
298 for (int ilocal2 = (ifo1 == ifo2 ? ilocal1 : 0); ilocal2 < fo2->getNPar(); ++ilocal2) {
299 int iglobal2 = fo2->getGlobalParNum(ilocal2);
300 double calc = M[idim * iglobal1 + iglobal2];
301 double num = num2ndDerivative(ifo1, ilocal1, eps, ifo2, ilocal2, eps);
302 B2INFO("fo1: " << fo1->getName() << " par " << ilocal1 << "/"
303 << iglobal1 << " (" << fo1->getParamName(ilocal1)
304 << "), fo2: " << fo2->getName() << " par " << ilocal2 << "/"
305 << iglobal2 << " (" << fo2->getParamName(ilocal2)
306 << ") calc: " << calc << " - num: " << num << " = " << calc - num);
307 }
308 }
309 }
310 }
311 delete[] M;
312 }
Belle2::OrcaKinFit::SoftGaussParticleConstraint::add2ndDerivativesToMatrix
virtual void add2ndDerivativesToMatrix(double *M, int idim) const override
Adds second order derivatives to global covariance matrix M.
Definition: SoftGaussParticleConstraint.cc:92
Belle2::OrcaKinFit::SoftGaussParticleConstraint::num2ndDerivative
double num2ndDerivative(int ifo1, int ilocal1, double eps1, int ifo2, int ilocal2, double eps2)
Evaluates numerically the 2nd derivative w.r.t. 2 parameters.
Definition: SoftGaussParticleConstraint.cc:329

Friends And Related Function Documentation

◆ operator<<()

std::ostream & operator<< ( std::ostream &  os,
const BaseConstraint &  bc 
)
related

Prints out a BaseConstraint, using its print method.

Parameters
osThe output stream
bcThe object to print

Definition at line 114 of file BaseConstraint.h.

117 {
118 return bc.print(os);
119 }

Member Data Documentation

◆ derivatives

std::vector<double> derivatives
protectedinherited

The derivatives.

Definition at line 181 of file SoftGaussParticleConstraint.h.

◆ fitobjects

FitObjectContainer fitobjects
protectedinherited

The FitObjectContainer.

Definition at line 179 of file SoftGaussParticleConstraint.h.

◆ flags

std::vector<int> flags
protectedinherited

The flags can be used to divide the FitObjectContainer into several subsets used for example to implement an equal mass constraint (see MassConstraint).

Definition at line 184 of file SoftGaussParticleConstraint.h.

◆ mass

double mass
protected

The mass difference between object sets 1 and 2.

Definition at line 74 of file SoftGaussMassConstraint.h.

◆ name

char* name
protectedinherited

Definition at line 108 of file BaseConstraint.h.

◆ sigma

double sigma
protectedinherited

The sigma of the Gaussian.

Definition at line 187 of file SoftGaussParticleConstraint.h.

The documentation for this class was generated from the following files: