MarlinUtil  1.12.1
ClusterShapes.h
1 #ifndef ClusterShapes_h
2 #define ClusterShapes_h
3 
4 
5 #include <stdio.h>
6 #include <iostream>
7 #include <iomanip>
8 #include <string>
9 #include <sstream>
10 #include <cstdlib>
11 #include "HelixClass.h"
12 #include <math.h>
13 
14 
26 
27 public:
28 
40  ClusterShapes(int nhits, float* a, float* x, float* y, float* z);
41 
42 
47 
51  void setErrors(float *ex, float* ey, float *ez);
52 
58  void setHitTypes(int *ityp);
59 
60 
64  int getNumberOfHits();
65 
70  float getTotalAmplitude();
71 
77  float* getCentreOfGravity();
78  // this is (for now) a pure dummy to allow MarlinPandora development!
79  float* getCentreOfGravityErrors();
80 
82  inline float* getCenterOfGravity() { return getCentreOfGravity() ; }
83  // this is (for now) a pure dummy to allow MarlinPandora development!
84  inline float* getCenterOfGravityErrors() { return getCentreOfGravityErrors() ; }
85 
90  float* getEigenValInertia();
91  // this is (for now) a pure dummy to allow MarlinPandora development!
92  float* getEigenValInertiaErrors();
93 
100  float* getEigenVecInertia();
101  // this is (for now) a pure dummy to allow MarlinPandora development!
102  float* getEigenVecInertiaErrors();
103 
111  float getWidth();
112 
121  int getEigenSytemCoordinates(float* xlong, float* xtrans);
122 
133  int getEigenSytemCoordinates(float* xlong, float* xtrans, float* a);
134 
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);
156 
166  float getChi2Fit3DProfileSimple(float a, float b, float c, float d, float* X0,
167  float* Rm);
168 
178  float getChi2Fit3DProfileAdvanced(float E0, float a, float b, float d, float t0,
179  float* X0, float* Rm);
180 
224  int FitHelix(int max_iter, int status_out, int parametrisation,
225  double* parameter, double* dparameter, double& chi2, double& distmax, int direction=1);
226 
227 
228  int FitHelix(int max_iter, int status_out, int parametrisation,
229  float* parameter, float* dparameter, float& chi2, float& distmax, int direction=1);
230 
231  //here add my functions(variables estimated with detector base)
232  //maximum deposit energy of hits
233  float getEmax(float* xStart, int& index_xStart, float* X0, float* Rm);
234 
235  //shower max of the hits from the shower start hit
236  float getsmax(float* xStart, int& index_xStart, float* X0, float* Rm);
237 
238  //radius where 90% of the cluster energy exists
239  float getxt90(float* xStart, int& index_xStart, float* X0, float* Rm);
240 
241  //length where less than 20% of the cluster energy exists
242  float getxl20(float* xStart, int& index_xStart, float* X0, float* Rm);
243 
244  //for cluster study
245  void gethits(float* xStart, int& index_xStart, float* X0, float* Rm, float *okxl, float *okxt, float *oke);
250  inline float radius() { return _radius; }
251 
257  inline float getElipsoid_r1() { return _r1; }
258 
264  inline float getElipsoid_r2() { return _r2; }
265 
271  inline float getElipsoid_r3() { return _r3; }
272 
276  inline float getElipsoid_vol() { return _vol; }
277 
281  inline float getElipsoid_r_ave() { return _r_ave; }
282 
286  inline float getElipsoid_density() { return _density; }
287 
292  inline float getElipsoid_eccentricity() { return _eccentricity; }
293 
298  inline float getElipsoid_r_forw() { return _r1_forw; }
299 
304  inline float getElipsoid_r_back() { return _r1_back; }
305 
306  //Mean of the radius of the hits
307  float getRhitMean(float* xStart, int& index_xStart, float* X0, float* Rm);
308 
309  //RMS of the radius of the hits
310  float getRhitRMS(float* xStart, int& index_xStart, float* X0, float* Rm);
311 
312 
313 
314 private:
315 
316  int _nHits;
317 
318  float* _aHit;
319  float* _xHit;
320  float* _yHit;
321  float* _zHit;
322  float* _exHit;
323  float* _eyHit;
324  float* _ezHit;
325  int* _types;
326  float* _xl;
327  float* _xt;
328  float* _t;
329  float* _s;
330 
331  int _ifNotGravity;
332  float _totAmpl;
333  float _radius;
334  float _xgr;
335  float _ygr;
336  float _zgr;
337  float _analogGravity[3];
338 
339  int _ifNotWidth;
340  float _analogWidth;
341 
342  int _ifNotInertia;
343  float _ValAnalogInertia[3];
344  float _VecAnalogInertia[9];
345 
346  int _ifNotEigensystem;
347 
348  int _ifNotElipsoid;
349  float _r1 ; // Cluster spatial axis length -- the largest
350  float _r2 ; // Cluster spatial axis length -- less
351  float _r3 ; // Cluster spatial axis length -- less
352  float _vol ; // Cluster ellipsoid volume
353  float _r_ave ; // Cluster average radius (qubic root)
354  float _density ; // Cluster density
355  float _eccentricity ; // Cluster Eccentricity
356  float _r1_forw ;
357  float _r1_back ;
358 
359  void findElipsoid();
360  void findGravity();
361  void findInertia();
362  void findWidth();
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);
374 
375  // private methods for non-linear, multidim. fitting (helix)
376  // static int functParametrisation1(const gsl_vector* par, void* data, gsl_vector* f);
377  // static int dfunctParametrisation1(const gsl_vector* par, void* d, gsl_matrix* J);
378  // static int fdfParametrisation1(const gsl_vector* par, void* d, gsl_vector* f, gsl_matrix* J);
379 
380 
381 };
382 
383 
384 #endif
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