L3 becomes HLT
[u/mrichter/AliRoot.git] / HLT / TPCLib / AliHLTTPCVertexArray.cxx
CommitLineData
a6c02c85 1// @(#) $Id$
4aa41877 2// Original: AliHLTVertexArray.cxx,v 1.7 2005/06/14 10:55:21 cvetan
a6c02c85 3
4// Author: Uli Frankenfeld <mailto:franken@fi.uib.no>
5//*-- Copyright &copy ALICE HLT Group
6
a6c02c85 7#include "AliHLTTPCLogging.h"
8#include "AliHLTTPCVertexArray.h"
9
10
11/** \class AliHLTTPCVertexArray
12<pre>
13//_____________________________________________________________
14// AliHLTTPCVertexArray
15//
16// The HLTTPC Fast Vertex Finder Base Class
17//
18// usage:
19//
20//for(Int_t sec=0;sec<NSEC;sec++){
21// ResetSector();
22// FillSectorSeed3D(x,y,z);
23// FillSector3D(x,y,z);
24// FindSectorVertex();
25// Double_t z = GetZSector();
26// Double_t zerr = GetZSectorErr();
27// // do somethink with z, zerr
28//}
29</pre>
30*/
31
32#if __GNUC__ >= 3
33using namespace std;
34#endif
35
36ClassImp(AliHLTTPCVertexArray)
37
38void AliHLTTPCVertexArray::AnalyzeSector(Float_t *vertex, Int_t *array, Int_t len)
39{
40 //loop over all seeds and all vertex position
41 LOG(AliHLTTPCLog::kInformational,"AliHLTTPCVertexArray::AnalyzeSector","Analyze")
42 <<AliHLTTPCLog::kDec<<"Number of Seeds: "<<fNSeed<<ENDLOG;
43 for(Int_t i =0;i<fNSeed;i++)
44 for(Int_t bin = 0;bin<len;bin++)
45 array[bin] += Trace(fZSeed[i],fRSeed[i],fSecSeed[i],vertex[bin]);
46}
47
48void AliHLTTPCVertexArray::FindSectorVertex(Double_t pos, Double_t range, Int_t nbin)
49{
50 //define position and range for search and
51 //loop over all seeds
52 //and find position of vertex and error
53 const Int_t klen = nbin;
54 const Double_t kwidth = range;
55 const Double_t kxmin = pos - kwidth/2;
56 const Double_t kstep = kwidth/klen;
57 const Double_t kstart = kxmin + kstep/2.;
58 Int_t * array = new Int_t[klen];
59 Float_t * ver = new Float_t[klen];
60 for(Int_t i=0;i<klen;i++){
61 ver[i] = kstart + kstep * i;
62 array[i] = 0;
63 }
64 AnalyzeSector(ver,array,klen);
65 FindMean(ver,array,klen);
66 delete[] array;
67 delete[] ver;
68}
69
70void AliHLTTPCVertexArray::FindMean(Float_t *vertex,Int_t *array, Int_t len)
71{
72 //find mean and error of array and store it in
73 //fZSector and fZSectorErr
74 const Int_t knbin = len;
75 Int_t xbin =0;
76 Int_t max=0;
77 for(Int_t i = 0;i<knbin;i++){
78 if(array[i]>max){
79 max = array[i];
80 xbin =i;
81 }
82 }
83 Int_t hmax = max/2;
84 Int_t xmin,xmax;
85 Int_t ops = 0;
86 xmin = xbin;
87 while(xmin--){
88 if(xmin<0) {ops++;break;}
89 if(array[xmin]<hmax) {
90 break;
91 }
92 }
93 xmax = xbin;
94 while(xmax++){
95 if(xmax>=knbin) {ops++;break;}
96 if(array[xmax]<hmax){
97 break;
98 }
99 }
100 if(ops){
101 if(xbin >= knbin/2){xmin = 2 * xbin - knbin +1;xmax = knbin-1;}
102 else{xmin = 0;xmax = 2 * xbin;}
103 }
104 Double_t sumw=0;
105 Double_t sumw2=0;
106 Double_t sumwx=0;
107 Double_t sumwx2=0;
108 for(Int_t bin = xmin;bin<=xmax;bin++){
109 sumw += array[bin];
110 sumw2 += array[bin] * array[bin];
111 sumwx += array[bin] * vertex[bin];
112 sumwx2 += array[bin] * vertex[bin] * vertex[bin];
113 }
114 if(sumw){
115 Double_t mean = sumwx/sumw;
116 Double_t rms2 = fabs(sumwx2/sumw - mean*mean);
117 fZSectorErr = sqrt(rms2/sumw);
118 fZSector = mean;
119 }
120 else{fZSectorErr = fZSector = 0;}
121 sumw=sumw2=sumwx=sumwx2=0;
122 for(Int_t bin = 0;bin<knbin;bin++){
123 sumw += array[bin];
124 sumw2 += array[bin] * array[bin];
125 sumwx += array[bin] * vertex[bin];
126 sumwx2 += array[bin] * vertex[bin] * vertex[bin];
127 }
128 if(sumw){
129 Double_t mean = fZSector;
130 Double_t rms2 = fabs(sumwx2/sumw - mean*mean);
131 fZSectorErr = sqrt(rms2/sumw);
132 }
133}
134