Line data Source code
1 : #include "TObjArray.h"
2 : #include "TMath.h"
3 :
4 : #include "AliLog.h"
5 : #include "AliMathBase.h"
6 :
7 : #include "AliTRDseedV1.h"
8 : #include "AliTRDcluster.h"
9 : #include "AliTRDpadPlane.h"
10 : #include "AliTRDtrackletOflHelper.h"
11 :
12 48 : ClassImp(AliTRDtrackletOflHelper)
13 : //___________________________________________________________________
14 : AliTRDtrackletOflHelper::AliTRDtrackletOflHelper()
15 458 : :TObject()
16 458 : ,fRow(-1)
17 458 : ,fClusters(NULL)
18 458 : ,fPadPlane(NULL)
19 2290 : {
20 : // Default constructor
21 458 : fCol[0]=144; fCol[1]=0;
22 458 : fTBrange[0]=1; fTBrange[1]=21;
23 916 : }
24 :
25 : //___________________________________________________________________
26 : AliTRDtrackletOflHelper::AliTRDtrackletOflHelper(const AliTRDtrackletOflHelper &ref)
27 0 : :TObject(ref)
28 0 : ,fRow(-1)
29 0 : ,fClusters(NULL)
30 0 : ,fPadPlane(ref.fPadPlane)
31 0 : {
32 : // Copy constructor
33 0 : fRow = ref.fRow; fCol[0]=ref.fCol[0]; fCol[1]=ref.fCol[1];
34 0 : fTBrange[0] = ref.fTBrange[0];fTBrange[1] = ref.fTBrange[1];
35 : Int_t n(0);
36 0 : if(ref.fClusters){
37 0 : if((n = ref.fClusters->GetEntriesFast())) {
38 0 : fClusters = new TObjArray(n);
39 0 : for(Int_t ic(n); ic--;) fClusters->AddAt(ref.fClusters->At(ic), ic);
40 0 : }
41 : }
42 0 : }
43 :
44 : //___________________________________________________________________
45 : AliTRDtrackletOflHelper& AliTRDtrackletOflHelper::operator=(const AliTRDtrackletOflHelper &rhs)
46 : {
47 80 : if(this != &rhs){
48 40 : TObject::operator=(rhs);
49 40 : fPadPlane = rhs.fPadPlane;
50 40 : fRow = rhs.fRow;
51 40 : fCol[0] = rhs.fCol[0];
52 40 : fCol[1] = rhs.fCol[1];
53 40 : fTBrange[0]= rhs.fTBrange[0];
54 40 : fTBrange[1]= rhs.fTBrange[1];
55 40 : if(rhs.fClusters){
56 40 : Int_t n(rhs.fClusters->GetEntriesFast());
57 96 : if(!fClusters) fClusters = new TObjArray(n);
58 40 : if(fClusters->GetEntriesFast() != n){
59 30 : fClusters->Clear();
60 30 : fClusters->Expand(n);
61 30 : }
62 1668 : for(Int_t ic(n); ic--;) fClusters->AddAt(rhs.fClusters->At(ic), ic);
63 40 : } else {
64 0 : if(fClusters) delete fClusters;
65 : }
66 : }
67 40 : return *this;
68 0 : }
69 :
70 : //___________________________________________________________________
71 : AliTRDtrackletOflHelper::~AliTRDtrackletOflHelper()
72 1832 : {
73 : // Clean helper
74 914 : if(fClusters) fClusters->Clear();
75 716 : delete fClusters;
76 916 : }
77 :
78 : //___________________________________________________________________
79 : Int_t AliTRDtrackletOflHelper::Expand(TObjArray *cls, Int_t *mark, Int_t groupId)
80 : {
81 : // Allocate new clusters for this helper. If a subset of clusters is to be allocated
82 : // this can be specified via the identifier "groupId" which have to be compared
83 : // with individual elements in the array "mark".
84 : //
85 : // Return total no. of clusters in helper
86 :
87 216 : if(!fPadPlane || !fClusters ){
88 0 : AliError("Helper not initialized.");
89 0 : return 0;
90 : }
91 72 : Int_t ncl(fClusters->GetEntriesFast());
92 72 : Int_t mcl(cls->GetEntriesFast());
93 72 : fClusters->Expand(mcl + ncl);
94 72 : fCol[0]=144; fCol[1]=0;
95 : AliTRDcluster *c(NULL);
96 4612 : for(Int_t icl(0), jcl(ncl); icl<mcl; icl++){
97 2234 : if(!(c=(AliTRDcluster*)cls->At(icl))) continue;
98 4468 : if(mark && mark[icl]!=groupId) continue;
99 796 : fClusters->AddAt(c, jcl++);
100 796 : Int_t col(c->GetPadCol());
101 796 : Short_t *sig(c->GetSignals());
102 12736 : for(Int_t icol(0), jcol(col-3); icol<7; icol++, jcol++){
103 5572 : if(sig[icol]==0) continue;
104 4812 : if(jcol<fCol[0]) fCol[0]=jcol;
105 5112 : if(jcol>fCol[1]) fCol[1]=jcol;
106 : }
107 796 : }
108 :
109 216 : AliDebug(1, Form("Segment[%d] Clusters[%2d] pad[%3d %3d].", groupId, fClusters->GetEntriesFast(), fCol[0], fCol[1]));
110 72 : return fClusters->GetEntriesFast();
111 72 : }
112 :
113 :
114 : //___________________________________________________________________
115 : Int_t AliTRDtrackletOflHelper::Init(AliTRDpadPlane *p, TObjArray *cls, Int_t *mark, Int_t groupId)
116 : {
117 : // Allocate clusters for this helper. If a subset of clusters is to be allocated
118 : // this can be specified via the identifier "groupId" which have to be compared
119 : // with individual elements in the array "mark".
120 : //
121 : // Return no. of clusters allocated
122 :
123 994 : if(!p){
124 0 : AliError("PadPlane not initialized."); return 0;
125 : }
126 497 : if(!cls){
127 0 : AliError("Cluster array not initialized."); return 0;
128 : }
129 497 : Int_t ncl(cls->GetEntriesFast());
130 497 : if(ncl<2){
131 0 : AliDebug(1, Form("Segment[%d] failed n[%d].", groupId, ncl)); return 0;
132 : }
133 : Int_t mcl(ncl);
134 497 : if(mark){
135 : mcl = 0;
136 6061 : for(Int_t icl(ncl); icl--;)
137 16189 : if(cls->At(icl) && mark[icl]==groupId) mcl++;
138 267 : }
139 497 : if(mcl<2){
140 0 : AliDebug(1, Form("Segment[%d] failed n[%d] in group.", groupId, mcl)); return 0;
141 : }
142 957 : if(!fClusters) fClusters = new TObjArray(mcl);
143 : else{
144 267 : fClusters->Clear();
145 267 : fClusters->Expand(mcl);
146 : }
147 497 : fCol[0]=144; fCol[1]=0; fRow=-1;
148 : AliTRDcluster *c(NULL);
149 22202 : for(Int_t icl(0), jcl(0); icl<ncl; icl++){
150 10604 : if(!(c=(AliTRDcluster*)cls->At(icl))) continue;
151 16131 : if(mark && mark[icl]!=groupId) continue;
152 10212 : fClusters->AddAt(c, jcl++);
153 10709 : if(fRow<0) fRow = c->GetPadRow();
154 10212 : Int_t col(c->GetPadCol());
155 10212 : Short_t *sig(c->GetSignals());
156 163392 : for(Int_t icol(0), jcol(col-3); icol<7; icol++, jcol++){
157 71484 : if(sig[icol]==0) continue;
158 58336 : if(jcol<fCol[0]) fCol[0]=jcol;
159 60478 : if(jcol>fCol[1]) fCol[1]=jcol;
160 : }
161 10212 : }
162 497 : fPadPlane = p;
163 :
164 1491 : AliDebug(1, Form("Segment[%d] Clusters[%2d] pad[%3d %3d].", groupId, fClusters->GetEntriesFast(), fCol[0], fCol[1]));
165 497 : return fClusters->GetEntriesFast();
166 497 : }
167 :
168 : //___________________________________________________________________
169 : Int_t AliTRDtrackletOflHelper::ClassifyTopology()
170 : {
171 : // Classify topology and return classification code
172 : // 0 - normal tracklet
173 : // 1 - delta ray candidate
174 : // 2 - secondary candidate
175 : // 3 - "elephant" candidate
176 : // 4 - unknown topology
177 :
178 460 : if(!fClusters){
179 0 : AliError("Helper not initialized. Missing clusters.");
180 0 : return kUnknown;
181 : }
182 230 : Int_t ncl(fClusters->GetEntries());
183 230 : if(!ncl){
184 0 : AliError("Helper not initialized. No cluster allocated.");
185 0 : return kUnknown;
186 : }
187 : const Int_t kRange(22); // DEFINE based on vd
188 :
189 : // compute local occupancy
190 230 : Int_t localOccupancy[AliTRDseedV1::kNtb]; memset(localOccupancy, 0, AliTRDseedV1::kNtb*sizeof(Int_t));
191 : AliTRDcluster *c(NULL);
192 230 : Double_t sy[kNcls]; Int_t mcl(0);
193 10024 : for(Int_t icl(ncl), time; icl--;){
194 5077 : c = (AliTRDcluster*)fClusters->At(icl);
195 5077 : time = c->GetPadTime();
196 5077 : if(time==0 || time>=kRange){
197 590 : sy[icl] = -1; // mark clusters outsde drift volume
198 590 : continue;
199 : }
200 : // protect against wrong error param.
201 4487 : if(c->GetSigmaY2() < 1.e-5) sy[icl] = 0.02;
202 4487 : else sy[icl] = TMath::Min(TMath::Sqrt(c->GetSigmaY2()), 0.04);
203 4487 : localOccupancy[time]++;
204 4487 : mcl++;
205 : }
206 230 : if(!mcl){
207 0 : AliWarning("No clusters in the active area.");
208 0 : return kUnknown;
209 : }
210 : //compute number of 0 bins and high occupancy
211 : Int_t goodOccupancy(0), highOccupancy(0), lowOccupancy(0);
212 10120 : for(Int_t itb(1); itb<kRange; itb++){
213 4830 : switch(localOccupancy[itb]){
214 559 : case 0: lowOccupancy++; break;
215 4083 : case 1: goodOccupancy++; break;
216 188 : default: highOccupancy++; break;
217 : }
218 : }
219 690 : AliDebug(2, Form("H[%2d | %5.2f%%] L[%2d | %5.2f%%] N[%2d | %5.2f%% | %5.2f%%]",
220 : highOccupancy, 100.*highOccupancy/kRange,
221 : lowOccupancy, 100.*lowOccupancy/kRange,
222 : goodOccupancy, 100.*goodOccupancy/kRange, 100.*goodOccupancy/mcl));
223 :
224 : // filter
225 230 : Double_t dy[kNcls];
226 446 : if(goodOccupancy==mcl) return kNormal;
227 14 : else if(Double_t(goodOccupancy)/kRange > 0.9){
228 0 : if(highOccupancy == 0 ) {
229 0 : if(!FitPSR(dy)) return kUnknown;
230 : Int_t nin(0);
231 0 : for(Int_t idy(ncl); idy--;){
232 0 : if(sy[idy]<0.) continue;
233 0 : if(dy[idy] > 3*sy[idy]) continue;
234 0 : nin ++;
235 : }
236 0 : if(Double_t(nin)/mcl > 0.8) return kNormal;
237 0 : else return kUnknown;
238 0 : } else return kDeltaRay;
239 14 : } else return kUnknown;
240 690 : }
241 :
242 :
243 : //___________________________________________________________________
244 : void AliTRDtrackletOflHelper::FindSolidCls(Bool_t *mark, Int_t *q)
245 : {
246 : // Find clusters produced by large fluctuations of energy deposits
247 : // Largest charge and well separation from neighbors
248 :
249 : Int_t ntb(AliTRDseedV1::kNtb);
250 0 : Int_t idx[ntb+1];
251 0 : TMath::Sort(ntb, q, idx, kTRUE);
252 0 : Int_t qmax = Int_t(0.3*q[idx[0]]);
253 0 : mark[0] = kFALSE;
254 0 : for(Int_t icl(ntb-5); icl<ntb; icl++) mark[icl] = kFALSE;
255 0 : for(Int_t icl(0); icl<ntb; icl++){
256 0 : Int_t jcl(idx[icl]);
257 0 : if(!mark[jcl]) continue;
258 0 : if(q[jcl-1]>q[jcl] || q[jcl+1]>q[jcl]){
259 0 : mark[jcl] = kFALSE;
260 0 : continue;
261 : }
262 0 : if(q[jcl] < qmax){
263 0 : mark[jcl] = kFALSE;
264 0 : continue;
265 : }
266 0 : for(Int_t kcl=TMath::Max(0, jcl-2); kcl<jcl+3; kcl++){
267 0 : if(kcl==jcl) continue;
268 0 : mark[kcl] = kFALSE;
269 0 : }
270 0 : }
271 0 : }
272 :
273 : //___________________________________________________________________
274 : Bool_t AliTRDtrackletOflHelper::FitPSR(Double_t dy[200], Bool_t useSolid)
275 : {
276 : // Fit tracklet in Pad System of Reference [PSR] to avoid uncertainty related to
277 : // Lorentz angle correction
278 :
279 0 : Bool_t mark[200]; memset(mark, 0, 200*sizeof(Bool_t));
280 0 : Int_t q[200]; memset(q, 0, 200*sizeof(Int_t));
281 0 : Int_t ncl(fClusters->GetEntries());
282 : AliTRDcluster *c(NULL);
283 0 : for(Int_t icl(ncl); icl--;){
284 0 : c = (AliTRDcluster*)fClusters->At(icl);
285 0 : if(c->GetPadRow() != fRow) continue;
286 0 : Int_t time(c->GetPadTime());
287 0 : mark[time] = kTRUE; q[time] = Int_t(c->GetQ());
288 : }
289 0 : if(useSolid) FindSolidCls(mark, q);
290 :
291 0 : Double_t
292 : x[AliTRDseedV1::kNtb], y[AliTRDseedV1::kNtb], sy[AliTRDseedV1::kNtb],
293 : xf[AliTRDseedV1::kNtb], yf[AliTRDseedV1::kNtb], syf[AliTRDseedV1::kNtb],
294 : par[3];
295 : Int_t jcl(0), kcl(0);
296 0 : for(Int_t icl(0); icl<ncl; icl++){
297 0 : c = (AliTRDcluster*)fClusters->At(icl);
298 0 : if(c->GetPadRow() != fRow) continue;
299 0 : Int_t col(c->GetPadCol());
300 0 : Int_t time(c->GetPadTime());
301 0 : Double_t center(c->GetCenter());
302 0 : Double_t cw(fPadPlane->GetColSize(col));
303 : //Double_t corr = AliTRDcluster::GetYcorr(AliTRDgeometry::GetLayer(det), center);
304 0 : y[jcl] = fPadPlane->GetColPos(col) + (.5 + center)*cw /*+ corr*/;
305 0 : x[jcl] = c->GetX();
306 0 : sy[jcl]= TMath::Sqrt(c->GetSigmaY2());
307 0 : if(mark[time]){
308 0 : yf[kcl] = y[jcl];
309 0 : xf[kcl] = x[jcl];
310 0 : syf[kcl]= useSolid ? 0.5/time:sy[jcl];
311 0 : kcl++;
312 0 : }
313 0 : jcl++;
314 0 : }
315 0 : Fit(kcl, xf, yf, syf, par);
316 :
317 0 : for(Int_t icl(0); icl<jcl; icl++){
318 0 : Double_t dx(x[icl] - par[2]);
319 0 : dy[icl] = y[icl] - (par[0] + par[1]*dx);
320 : }
321 0 : return kTRUE;
322 0 : }
323 :
324 : //___________________________________________________________________
325 : Bool_t AliTRDtrackletOflHelper::Fit(Int_t n, Double_t *x, Double_t *y, Double_t *sy, Double_t *par, Double_t sCut, Double_t *cov)
326 : {
327 : // Iterative robust tracklet fit
328 1120 : if(n<3) return kFALSE;
329 : //select reference radial position
330 557 : if(par[2]<0.){ //compute reference radial position as <x>
331 0 : par[2] = 0.; for(Int_t ic(n); ic--;) par[2] += x[ic]; par[2] /= n;
332 0 : }
333 557 : AliTRDtrackerV1::AliTRDLeastSquare &f=Fitter();
334 4456 : for(Int_t iter(0); iter<3; iter++){
335 1671 : f.Reset();
336 : Int_t jp(0);
337 74730 : for(Int_t ip(0); ip<n; ip++){
338 35694 : Double_t dx(x[ip]-par[2]);
339 35694 : Double_t dy(y[ip] - (par[0] + par[1]*dx));
340 68277 : if(iter && TMath::Abs(dy)>sCut*sy[ip]) continue;
341 26907 : f.AddPoint(&dx, y[ip], sy[ip]);
342 26907 : jp++;
343 62601 : }
344 1787 : if(jp<3) continue;
345 1555 : if(!f.Eval()) continue;
346 1555 : par[0]=f.GetFunctionParameter(0);
347 1555 : par[1]=f.GetFunctionParameter(1);
348 2231 : if(cov) f.GetCovarianceMatrix(cov);
349 4665 : AliDebugGeneral("AliTRDtrackletOflHelper::Fit()", 2, Form("Iter[%d] Ncl[%2d/%2d] par[%f %f %f]", iter, jp, n, par[0], par[1], par[2]));
350 1555 : }
351 : return kTRUE;
352 559 : }
353 :
354 : //___________________________________________________________________
355 : Bool_t AliTRDtrackletOflHelper::Fit(Double_t *par, Double_t sCut) const
356 : {
357 : // Wrapper for clusters attach to this for static Fit function
358 : //
359 662 : Int_t n(fClusters->GetEntriesFast()), jc(0), dr(0);
360 662 : Double_t corr = TMath::Tan(TMath::DegToRad()*fPadPlane->GetTiltingAngle())*
361 331 : fPadPlane->GetLengthIPad();
362 331 : Double_t x[kNcls], y[kNcls], sy[kNcls];
363 : AliTRDcluster *c(NULL);
364 16224 : for(Int_t ic(0); ic<n; ic++){
365 7781 : c = (AliTRDcluster*)fClusters->At(ic);
366 7781 : if(!c->IsInChamber()) continue;
367 7258 : dr = c->GetPadRow() - fRow;
368 7258 : x[jc] = c->GetX();
369 7258 : y[jc] = c->GetY()+corr*dr;
370 21774 : sy[jc]= c->GetSigmaY2()>0?(TMath::Min(TMath::Sqrt(c->GetSigmaY2()), 0.08)):0.08;
371 7258 : jc++;
372 7258 : }
373 662 : return Fit(jc, x, y, sy, par, sCut);
374 331 : }
375 :
376 : //___________________________________________________________________
377 : void AliTRDtrackletOflHelper::GetColSignals(Int_t col, Int_t adc[32], Bool_t mainRow) const
378 : {
379 0 : memset(adc, 0, 32*sizeof(Int_t));
380 0 : if(col<fCol[0] || col>fCol[1]) return;
381 :
382 : AliTRDcluster *cc(NULL);
383 0 : for(Int_t ic(fClusters->GetEntriesFast()); ic--;){
384 0 : cc = (AliTRDcluster*)fClusters->At(ic);
385 0 : if((mainRow && cc->GetPadRow()!=fRow) ||
386 0 : (!mainRow && cc->GetPadRow()==fRow)) continue;
387 0 : Short_t *sig = cc->GetSignals();
388 0 : Int_t padcol(cc->GetPadCol());
389 0 : Int_t time(cc->GetPadTime());
390 0 : for(Int_t icol(0), jcol(padcol-3); icol<7; icol++, jcol++){
391 0 : if(jcol!=col) continue;
392 0 : adc[time]+=sig[icol];
393 0 : }
394 : }
395 0 : }
396 :
397 : //___________________________________________________________________
398 : Int_t AliTRDtrackletOflHelper::GetRMS(Double_t &r, Double_t &m, Double_t &s, Double_t xm) const
399 : {
400 : // Calculate Rotation[r], Mean y[m] (at radial position [xm]) and Sigma[s] (of a gaussian distribution in the tracklet [SR])
401 : // for clusters attach to this helper. The Rotation and Mean are calculated without tilt correction option.
402 : // It returns the number of clusters in 1 sigma cut.
403 :
404 662 : Int_t n(fClusters->GetEntriesFast());
405 331 : if(n<=2) return n;
406 331 : Double_t par[3] = {0., 0., -1.};
407 662 : if(xm>0.) par[2] = xm; // select reference radial position from outside
408 331 : Fit(par);
409 331 : m = par[0];
410 331 : r = par[1];
411 331 : xm= par[2];
412 :
413 662 : Double_t corr = TMath::Tan(TMath::DegToRad()*fPadPlane->GetTiltingAngle())*
414 331 : fPadPlane->GetLengthIPad();
415 331 : Double_t y[kNcls];
416 : AliTRDcluster *c(NULL);
417 16224 : for(Int_t ic(n); ic--;){
418 7781 : c = (AliTRDcluster*)fClusters->At(ic);
419 7781 : Double_t x(c->GetX() - xm);
420 7781 : Int_t dr(c->GetPadRow() - fRow);
421 7781 : y[ic] = c->GetY()+corr*dr - (par[0] + par[1]*x);
422 : }
423 331 : Double_t m1(0.);
424 331 : AliMathBase::EvaluateUni(n, y, m1, s, 0);
425 : Int_t n0(0);
426 16224 : for(Int_t ic(n); ic--;){
427 7781 : c = (AliTRDcluster*)fClusters->At(ic);
428 23343 : Double_t sy = c->GetSigmaY2()>0?(TMath::Min(TMath::Sqrt(c->GetSigmaY2()), 0.08)):0.08;
429 11076 : if(TMath::Abs(y[ic]-m1) <= sy) n0++;
430 : }
431 : return n0;
432 662 : }
433 :
434 : //___________________________________________________________________
435 : Double_t AliTRDtrackletOflHelper::GetQ() const
436 : {
437 : // Calculate total charge NOT normalized to inclination
438 :
439 : Double_t q(0.);
440 662 : Int_t n(fClusters->GetEntriesFast());
441 : AliTRDcluster *c(NULL);
442 16224 : for(Int_t ic(n); ic--;){
443 7781 : c = (AliTRDcluster*)fClusters->At(ic);
444 7781 : q += TMath::Abs(c->GetQ());
445 : }
446 331 : return q;
447 : }
448 :
449 : //___________________________________________________________________
450 : Double_t AliTRDtrackletOflHelper::GetSyMean() const
451 : {
452 : // Calculate mean uncertainty of clusters
453 :
454 : Double_t sym(0.);
455 662 : Int_t n(fClusters->GetEntriesFast());
456 : AliTRDcluster *c(NULL);
457 16224 : for(Int_t ic(n); ic--;){
458 7781 : c = (AliTRDcluster*)fClusters->At(ic);
459 23343 : sym += c->GetSigmaY2()>0?(TMath::Min(TMath::Sqrt(c->GetSigmaY2()), 0.08)):0.08;
460 : }
461 331 : return sym/=n;
462 : }
463 :
464 : //___________________________________________________________________
465 : Int_t AliTRDtrackletOflHelper::Segmentation(Int_t n, Double_t *x, Double_t *y, Int_t *Index)
466 : {
467 : // Segmentation of clusters in the tracklet roads
468 : //
469 : // The user supply the coordinates of the clusters (x,y) and their number "n". Also
470 : // The array "Index" has to be allocated by the user with a size equal or larger than "n".
471 : // On return the function returns the number of segments found and the array "Index" i-th element
472 : // is filled with the index of the tracklet segment to which the i-th cluster was assigned too.
473 : //
474 : // Observation
475 : // The parameter which controls the segmentation is set inside the function "kGapSize" and it is used
476 : // for both "x" and "y" segmentations. An improvement can be a parametrization as function of angle of
477 : // incidence.
478 : //
479 : // author:
480 : // Alex Bercuci <a.bercuci@gsi.de>
481 :
482 94 : if(!n || !x || !y || !Index){
483 0 : AliErrorGeneral("AliTRDtrackletOflHelper::Segmentation()", "One of the input arrays non initialized.");
484 0 : return 0;
485 : }
486 : const Int_t kBuffer = 200;
487 47 : if(n>kBuffer){
488 0 : AliWarningGeneral("AliTRDtrackletOflHelper::Segmentation()", Form("Input array size %d exceed buffer %d. Truncate.", n, kBuffer));
489 : n = kBuffer;
490 0 : }
491 : const Double_t kGapSize(0.2); // cm
492 : Int_t ng(0),
493 : nc(0);
494 47 : Double_t xx[kBuffer], dy;
495 47 : Int_t idx[kBuffer+1], jdx[kBuffer], kdx[kBuffer];
496 47 : TMath::Sort(n, y, idx);
497 1978 : for(Int_t iy(0); iy<n; iy++){
498 2779 : dy = iy>0?(TMath::Abs(y[idx[iy-1]]-y[idx[iy]])):0.;
499 942 : if(dy>kGapSize){
500 30 : TMath::Sort(nc, xx, jdx);
501 864 : for(Int_t ic(0), jc0, jc1; ic<nc; ic++){
502 1176 : jc0 = ic>0?kdx[jdx[ic-1]]:0;
503 402 : jc1 = kdx[jdx[ic]];
504 402 : dy = TMath::Abs(y[jc0] - y[jc1]);
505 412 : if(ic && dy>kGapSize) ng++;
506 402 : Index[jc1] = ng;
507 1206 : AliDebugGeneral("AliTRDtrackletOflHelper::Segmentation()", 4, Form(" y[%2d]=%+f x[%+f] %2d -> %2d -> %2d ng[%2d]", jc1, y[jc1], x[jc1], ic, jdx[ic], jc1, ng));
508 : }
509 30 : ng++;
510 : nc=0;
511 30 : }
512 942 : xx[nc] = x[idx[iy]];
513 942 : kdx[nc]= idx[iy];
514 2826 : AliDebugGeneral("AliTRDtrackletOflHelper::Segmentation()", 4, Form("y[%2d]=%+f -> %2d", idx[iy], y[idx[iy]], nc));
515 942 : nc++;
516 : }
517 47 : if(nc){
518 47 : TMath::Sort(nc, xx, jdx);
519 1174 : for(Int_t ic(0), jc0, jc1; ic<nc; ic++){
520 1573 : jc0 = ic>0?kdx[jdx[ic-1]]:0;
521 540 : jc1 = kdx[jdx[ic]];
522 540 : dy = TMath::Abs(y[jc0] - y[jc1]);
523 548 : if(ic && dy>kGapSize) ng++;
524 540 : Index[jc1] = ng;
525 1620 : AliDebugGeneral("AliTRDtrackletOflHelper::Segmentation()", 4, Form(" y[%2d]=%+f x[%+f|%+f] %2d -> %2d -> %2d ng[%2d]\n", jc1, y[jc1], xx[jdx[ic]], x[jc1], ic, jdx[ic], jc1, ng));
526 : }
527 47 : ng++;
528 47 : }
529 : return ng;
530 94 : }
531 :
532 : //___________________________________________________________________
533 : void AliTRDtrackletOflHelper::SetTbRange(Float_t t0, Float_t vd)
534 : {
535 : // Set first time bin and total number of time bins corresponding to clusters
536 : // in chamber based on the calibrated info "t0" and drift velocity "vd"
537 :
538 : // TO CHECK
539 0 : fTBrange[0] = Int_t(t0*0.1);
540 0 : fTBrange[1] = Int_t(vd*0.1);
541 0 : }
542 :
543 : #include "TH2.h"
544 : #include "TGraph.h"
545 : #include "TGraphErrors.h"
546 : #include "TCanvas.h"
547 : #include "TLegend.h"
548 : #include "TGaxis.h"
549 : //___________________________________________________________________
550 : void AliTRDtrackletOflHelper::View(TVirtualPad *vpad)
551 : {
552 : // Visualization support. Draw this tracklet segment
553 :
554 :
555 : Int_t row(-1), det(-1);
556 0 : Int_t n(fClusters->GetEntriesFast());
557 0 : AliInfo(Form("Processing Clusters[%2d] Class[%d].", n ,ClassifyTopology()));
558 :
559 : // prepare drawing objects
560 0 : Int_t ncols(fPadPlane->GetNcols());
561 0 : Double_t cw(fPadPlane->GetColSize(1));
562 0 : TH2 *h2 = new TH2I("h2pm", ";y_{Local} [cm]; pad time;Charge", ncols, fPadPlane->GetColPos(1)-cw, fPadPlane->GetColPos(ncols-1)+cw, 31, -30.5, 0.5);
563 0 : h2->SetMarkerColor(kWhite);h2->SetMarkerSize(2.);
564 0 : h2->SetFillColor(9);
565 :
566 0 : TGraph *gcls = new TGraph(n);
567 0 : gcls->SetMarkerColor(kBlack);gcls->SetMarkerStyle(20);
568 0 : TGraph *gclsLC = new TGraph(n);
569 0 : gclsLC->SetMarkerColor(kBlack);gclsLC->SetMarkerStyle(28);
570 0 : TGraph *gclsSC = new TGraph(n);
571 0 : gclsSC->SetMarkerColor(kBlack);gclsSC->SetMarkerStyle(4);gclsSC->SetMarkerSize(1.5);
572 0 : TGraph *gFP = new TGraph(n);
573 0 : gFP->SetMarkerColor(kBlack);gFP->SetLineWidth(2);
574 :
575 : // fill signal data
576 0 : Bool_t map[200]; memset(map, 0, 200*sizeof(Bool_t));
577 0 : Int_t qa[200]; memset(qa, 0, 200*sizeof(Int_t));
578 : AliTRDcluster *cc(NULL); Int_t tm = 30; Double_t ym=0.;
579 0 : for(Int_t ic(n); ic--;){
580 0 : cc = (AliTRDcluster*)fClusters->At(ic);
581 0 : if(cc->GetPadRow() != fRow) continue; //TODO extend for 2 pad rows
582 0 : Short_t *sig = cc->GetSignals();
583 0 : Int_t col(cc->GetPadCol());
584 0 : det = cc->GetDetector();
585 0 : row = cc->GetPadRow();
586 0 : Int_t time(cc->GetPadTime());
587 0 : map[time] = kTRUE; qa[time] = (Int_t)cc->GetQ();
588 0 : for(Int_t ipad(0), jpad(col-2); ipad<7; ipad++, jpad++){
589 0 : Int_t q = (Int_t)h2->GetBinContent(jpad, 31-time);
590 0 : h2->SetBinContent(jpad, 31-time, q+sig[ipad]);
591 : }
592 0 : Double_t y0 = fPadPlane->GetColPos(col) + (.5 + cc->GetCenter())*cw;
593 : //h2->GetXaxis()->GetBinCenter(pad)+cc->GetCenter();
594 0 : gcls->SetPoint(ic, y0, -time);
595 0 : if(time <= tm) {ym = cc->GetY() - y0; tm = time;}
596 : }
597 :
598 : // draw special clusters (solid and Lorentz corrected and fit)
599 0 : FindSolidCls(map, qa);
600 0 : Double_t ddy[200]; FitPSR(ddy, kTRUE);
601 0 : Double_t dt, dy;
602 0 : for(Int_t ic(0), jc(0); ic<n; ic++){
603 0 : cc = (AliTRDcluster*)fClusters->At(ic);
604 0 : if(cc->GetPadRow() != fRow) continue; //TODO extend for 2 pad rows
605 0 : gcls->GetPoint(ic, dy, dt);
606 0 : gclsLC->SetPoint(ic, cc->GetY()-ym, dt);
607 0 : gFP->SetPoint(ic, dy-ddy[ic], dt);
608 0 : if(map[cc->GetPadTime()]) gclsSC->SetPoint(jc++, dy, dt);
609 : }
610 :
611 :
612 : // prepare frame histo
613 0 : h2->SetName(Form("h2s%03d%02d", det, row));
614 0 : Int_t binSrt(fCol[0]+1), binSop(fCol[1]+1);
615 : TH1 *h1(NULL);
616 :
617 : // show everything
618 0 : if(!vpad){
619 0 : vpad = new TCanvas(Form("c%03d%02d", det, row), Form("D %03d [%02d_%d_%d] R[%02d]", det, AliTRDgeometry::GetSector(det), AliTRDgeometry::GetStack(det), AliTRDgeometry::GetLayer(det), row), 700, 500);
620 0 : }
621 : TVirtualPad *pp(NULL);
622 0 : vpad->Divide(2,1,2.e-5,2.e-5);
623 0 : pp = vpad->cd(1); pp->SetRightMargin(0.0001);pp->SetTopMargin(0.1);pp->SetBorderMode(0); pp->SetFillColor(kWhite);
624 0 : h1 = h2->ProjectionY(); h1->GetYaxis()->SetTitle("Total Charge");
625 0 : h1->SetTitle(Form("D%03d[%02d_%d_%d] R[%02d]",
626 0 : det,AliTRDgeometry::GetSector(det), AliTRDgeometry::GetStack(det), AliTRDgeometry::GetLayer(det), row));
627 0 : TAxis *ax(h1->GetXaxis()); for(Int_t ib(0), jb(1); ib<ax->GetNbins(); ib++, jb++) ax->SetBinLabel(jb, Form("%d", -Int_t(ax->GetBinCenter(jb))));
628 0 : h1->Draw("hbar3");
629 :
630 0 : pp = vpad->cd(2);
631 0 : pp->SetRightMargin(0.15);pp->SetLeftMargin(0.0001);pp->SetTopMargin(0.1);
632 0 : pp->SetBorderMode(0); pp->SetFillColor(kWhite);
633 0 : ax=h2->GetXaxis();
634 0 : ax->SetRange(TMath::Max(1, binSrt-1), TMath::Min(ncols, binSop+1));
635 0 : h2->Draw("coltextz");
636 0 : gPad->Update();
637 0 : TGaxis *axis = new TGaxis(gPad->GetUxmin(),
638 0 : gPad->GetUymax(),
639 0 : gPad->GetUxmax(),
640 0 : gPad->GetUymax(),
641 0 : TMath::Max(0, binSrt-2)-0.5, TMath::Min(ncols-1, binSop)+0.5, 510,"-L");
642 :
643 0 : axis->SetNdivisions(103+binSop-binSrt);
644 0 : axis->SetTitle("pad col");
645 0 : axis->Draw();
646 :
647 0 : TLegend *leg = new TLegend(0.01, 0.1, 0.84, 0.21);
648 0 : leg->SetBorderSize(1); leg->SetFillColor(kYellow-9);leg->SetTextSize(0.04);
649 0 : gcls->Draw("p"); leg->AddEntry(gcls, "Cls. in Pad [SR]", "p");
650 0 : gclsLC->Draw("p"); leg->AddEntry(gclsLC, "Lorentz Corr. Cls.", "p");
651 0 : gclsSC->Draw("p"); leg->AddEntry(gclsSC, "Solid Cls.", "p");
652 0 : gFP->Draw("l"); leg->AddEntry(gFP, "Fit in Pad [SR]", "l");
653 0 : leg->Draw();
654 0 : }
655 :
656 : //___________________________________________________________________
657 : AliTRDtrackerV1::AliTRDLeastSquare& AliTRDtrackletOflHelper::Fitter()
658 : {
659 1120 : static AliTRDtrackerV1::AliTRDLeastSquare f;
660 557 : return f;
661 0 : }
|
