1 #ifndef ClusterShapes_h
2 #define ClusterShapes_h
11 #include "HelixClass.h"
40 ClusterShapes(
int nhits,
float* a,
float* x,
float* y,
float* z);
51 void setErrors(
float *ex,
float* ey,
float *ez);
79 float* getCentreOfGravityErrors();
84 inline float* getCenterOfGravityErrors() {
return getCentreOfGravityErrors() ; }
92 float* getEigenValInertiaErrors();
102 float* getEigenVecInertiaErrors();
154 int fit3DProfile(
float& chi2,
float& a,
float& b,
float& c,
float& d,
float& xl0,
155 float * xStart,
int& index_xStart,
float* X0,
float* Rm);
179 float* X0,
float* Rm);
224 int FitHelix(
int max_iter,
int status_out,
int parametrisation,
225 double* parameter,
double* dparameter,
double& chi2,
double& distmax,
int direction=1);
228 int FitHelix(
int max_iter,
int status_out,
int parametrisation,
229 float* parameter,
float* dparameter,
float& chi2,
float& distmax,
int direction=1);
233 float getEmax(
float* xStart,
int& index_xStart,
float* X0,
float* Rm);
236 float getsmax(
float* xStart,
int& index_xStart,
float* X0,
float* Rm);
239 float getxt90(
float* xStart,
int& index_xStart,
float* X0,
float* Rm);
242 float getxl20(
float* xStart,
int& index_xStart,
float* X0,
float* Rm);
245 void gethits(
float* xStart,
int& index_xStart,
float* X0,
float* Rm,
float *okxl,
float *okxt,
float *oke);
250 inline float radius() {
return _radius; }
307 float getRhitMean(
float* xStart,
int& index_xStart,
float* X0,
float* Rm);
310 float getRhitRMS(
float* xStart,
int& index_xStart,
float* X0,
float* Rm);
337 float _analogGravity[3];
343 float _ValAnalogInertia[3];
344 float _VecAnalogInertia[9];
346 int _ifNotEigensystem;
355 float _eccentricity ;
363 float findDistance(
int i);
364 float vecProduct(
float * x1,
float * x2);
365 float vecProject(
float * x,
float * axis);
366 double DistanceHelix(
double x,
double y,
double z,
double X0,
double Y0,
double R0,
double bz,
367 double phi0,
double * distRPhiZ);
368 int transformToEigensystem(
float* xStart,
int& index_xStart,
float* X0,
float* Xm);
369 float calculateChi2Fit3DProfileSimple(
float a,
float b,
float c,
float d);
370 float calculateChi2Fit3DProfileAdvanced(
float E0,
float a,
float b,
float d,
float t0);
371 int fit3DProfileSimple(
float& chi2,
float& a,
float& b,
float& c,
float& d);
372 int fit3DProfileAdvanced(
float& chi2,
double* par_init,
double* par,
int npar,
373 float* t,
float* s,
float* E,
float E0);
float getElipsoid_r_ave()
average radius of the ellipsoid (qubic root of volume)
Definition: ClusterShapes.h:281
float * getEigenValInertia()
array of the inertias of mass (i. e. energy) corresponding to the three main axes of inertia...
Definition: ClusterShapes.cc:684
float getElipsoid_r1()
largest spatial axis length of the ellipsoid derived by the inertia tensor (by their eigenvalues and ...
Definition: ClusterShapes.h:257
ClusterShapes(int nhits, float *a, float *x, float *y, float *z)
Constructor.
Definition: ClusterShapes.cc:574
float * getCenterOfGravity()
US spelling of getCentreOfGravity.
Definition: ClusterShapes.h:82
float * getEigenVecInertia()
array of the three main axes of inertia (9 entries) starting with the axis corresponding to the small...
Definition: ClusterShapes.cc:696
int getNumberOfHits()
returns the number of elements of the cluster
Definition: ClusterShapes.cc:659
float getChi2Fit3DProfileAdvanced(float E0, float a, float b, float d, float t0, float *X0, float *Rm)
returns the chi2 of the fit in the method Fit3DProfile (if advanced parametrisation is used) for a gi...
Definition: ClusterShapes.cc:871
int fit3DProfile(float &chi2, float &a, float &b, float &c, float &d, float &xl0, float *xStart, int &index_xStart, float *X0, float *Rm)
performs a least square fit on the shape of an electro- magnetic-shower, which is defined as: A[i] = ...
Definition: ClusterShapes.cc:761
float getElipsoid_vol()
volume of the ellipsoid
Definition: ClusterShapes.h:276
float getTotalAmplitude()
returns the accumulated amplitude for the whole cluster (just the sum of the energy in all the entrie...
Definition: ClusterShapes.cc:665
void setErrors(float *ex, float *ey, float *ez)
Defining errors for Helix Fit.
Definition: ClusterShapes.cc:631
float getWidth()
'mean' width of the cluster perpendicular to the main principal axis, defined as: width := sqrt( 1/To...
Definition: ClusterShapes.cc:708
float getElipsoid_r_back()
distance from centre of gravity to the point nearest to IP projected on the main principal axis ...
Definition: ClusterShapes.h:304
int getEigenSytemCoordinates(float *xlong, float *xtrans)
returns the coordinates of the cluster transformed into the CoG-System.
Definition: ClusterShapes.cc:715
float getElipsoid_r_forw()
distance from centre of gravity to the point most far away from IP projected on the main principal ax...
Definition: ClusterShapes.h:298
float radius()
distance to the centre of gravity measured from IP (absolut value of the vector to the centre of grav...
Definition: ClusterShapes.h:250
float getElipsoid_density()
density of the ellipsoid defined by: totAmpl/vol
Definition: ClusterShapes.h:286
void setHitTypes(int *ityp)
Defining hit types for Helix Fit : type 1 - cyllindrical detector type 2 - Z disk detector...
Definition: ClusterShapes.cc:641
Utility class to derive properties of clusters, such as centre of gravity, axes of inertia...
Definition: ClusterShapes.h:25
float getChi2Fit3DProfileSimple(float a, float b, float c, float d, float *X0, float *Rm)
returns the chi2 of the fit in the method Fit3DProfile (if simple parametrisation is used)for a given...
Definition: ClusterShapes.cc:853
float getElipsoid_eccentricity()
eccentricity of the ellipsoid defined by: Width/r1
Definition: ClusterShapes.h:292
float getElipsoid_r3()
smallest spatial axis length of the ellipsoid derived by the inertia tensor (by their eigenvalues and...
Definition: ClusterShapes.h:271
int FitHelix(int max_iter, int status_out, int parametrisation, double *parameter, double *dparameter, double &chi2, double &distmax, int direction=1)
performs a least square fit on a helix path in space, which which is defined as (Cartesian coordiante...
Definition: ClusterShapes.cc:922
~ClusterShapes()
Destructor.
Definition: ClusterShapes.cc:613
float getElipsoid_r2()
medium spatial axis length of the ellipsoid derived by the inertia tensor (by their eigenvalues and e...
Definition: ClusterShapes.h:264
float * getCentreOfGravity()
returns an array, which represents a vector from the origin of the coordiante system, i. e. IP, to the centre of gravity of the cluster.
Definition: ClusterShapes.cc:672