- corrected small typo that prevented the AliHLTPluginBase instance from being created
[u/mrichter/AliRoot.git] / PMD / AliPMDtracker.cxx
CommitLineData
7dcaf913 1/***************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//-----------------------------------------------------//
17// //
18// Date : March 25 2004 //
19// This reads the file PMD.RecPoints.root(TreeR), //
20// calls the Clustering algorithm and stores the //
21// clustering output in PMD.RecPoints.root(TreeR) //
22// //
23//-----------------------------------------------------//
24
25#include <Riostream.h>
26#include <TMath.h>
7dcaf913 27#include <TTree.h>
7dcaf913 28#include <TObjArray.h>
29#include <TClonesArray.h>
30#include <TFile.h>
31#include <TBranch.h>
32#include <TNtuple.h>
33#include <TParticle.h>
34
32b0cad8 35#include <TGeoMatrix.h>
36
37#include "AliGeomManager.h"
38
7dcaf913 39#include "AliPMDcluster.h"
40#include "AliPMDclupid.h"
41#include "AliPMDrecpoint1.h"
b0e4d1e1 42#include "AliPMDrecdata.h"
43#include "AliPMDrechit.h"
7dcaf913 44#include "AliPMDUtility.h"
45#include "AliPMDDiscriminator.h"
26f1ae17 46#include "AliPMDEmpDiscriminator.h"
7dcaf913 47#include "AliPMDtracker.h"
48
49#include "AliESDPmdTrack.h"
af885e0f 50#include "AliESDEvent.h"
ecee2a1a 51#include "AliLog.h"
7dcaf913 52
53ClassImp(AliPMDtracker)
54
55AliPMDtracker::AliPMDtracker():
56 fTreeR(0),
722ccc67 57 fRecpoints(new TClonesArray("AliPMDrecpoint1", 10)),
b0e4d1e1 58 fRechits(new TClonesArray("AliPMDrechit", 10)),
7dcaf913 59 fPMDcontin(new TObjArray()),
60 fPMDcontout(new TObjArray()),
7dcaf913 61 fPMDutil(new AliPMDUtility()),
62 fPMDrecpoint(0),
63 fPMDclin(0),
64 fPMDclout(0),
7dcaf913 65 fXvertex(0.),
66 fYvertex(0.),
67 fZvertex(0.),
68 fSigmaX(0.),
69 fSigmaY(0.),
70 fSigmaZ(0.)
71{
72 //
73 // Default Constructor
74 //
75}
76//--------------------------------------------------------------------//
a48edddd 77AliPMDtracker:: AliPMDtracker(const AliPMDtracker & /* tracker */):
2332574a 78 TObject(/* tracker */),
79 fTreeR(0),
80 fRecpoints(NULL),
b0e4d1e1 81 fRechits(NULL),
2332574a 82 fPMDcontin(NULL),
83 fPMDcontout(NULL),
84 fPMDutil(NULL),
85 fPMDrecpoint(0),
86 fPMDclin(0),
87 fPMDclout(0),
88 fXvertex(0.),
89 fYvertex(0.),
90 fZvertex(0.),
91 fSigmaX(0.),
92 fSigmaY(0.),
93 fSigmaZ(0.)
a48edddd 94{
95 // copy constructor
96 AliError("Copy constructor not allowed");
97}
98
99//--------------------------------------------------------------------//
100AliPMDtracker& AliPMDtracker::operator=(const AliPMDtracker & /* tracker */)
101{
102 // assignment operator
103 AliError("Assignment operator not allowed");
104 return *this;
105}
106
107//--------------------------------------------------------------------//
7dcaf913 108AliPMDtracker::~AliPMDtracker()
109{
110 // Destructor
111 if (fRecpoints)
112 {
722ccc67 113 fRecpoints->Clear();
7dcaf913 114 }
b0e4d1e1 115 if (fRechits)
116 {
117 fRechits->Clear();
118 }
119
7dcaf913 120 if (fPMDcontin)
121 {
122 fPMDcontin->Delete();
123 delete fPMDcontin;
124 fPMDcontin=0;
a1e9ac09 125
7dcaf913 126 }
127 if (fPMDcontout)
722ccc67 128 {
7dcaf913 129 fPMDcontout->Delete();
130 delete fPMDcontout;
131 fPMDcontout=0;
a1e9ac09 132
7dcaf913 133 }
316c6cd9 134 delete fPMDutil;
7dcaf913 135}
136//--------------------------------------------------------------------//
137void AliPMDtracker::LoadClusters(TTree *treein)
138{
139 // Load the Reconstructed tree
140 fTreeR = treein;
141}
142//--------------------------------------------------------------------//
af885e0f 143void AliPMDtracker::Clusters2Tracks(AliESDEvent *event)
7dcaf913 144{
145 // Converts digits to recpoints after running clustering
146 // algorithm on CPV plane and PREshower plane
147 //
148
be8b7039 149 Int_t idet = 0;
150 Int_t ismn = 0;
151 Int_t trackno = 1, trackpid = 0;
152 Float_t clusdata[6] = {0.,0.,0.,0.,0.,0.};
b0e4d1e1 153
154 Int_t *irow;
155 Int_t *icol;
156 Int_t *itra;
157 Int_t *ipid;
158 Float_t *cadc;
159
160 AliPMDrechit *rechit = 0x0;
7dcaf913 161
162 TBranch *branch = fTreeR->GetBranch("PMDRecpoint");
ecee2a1a 163 if (!branch)
164 {
165 AliError("PMDRecpoint branch not found");
166 return;
167 }
7dcaf913 168 branch->SetAddress(&fRecpoints);
b0e4d1e1 169
170 TBranch *branch1 = fTreeR->GetBranch("PMDRechit");
171 if (!branch1)
172 {
173 AliError("PMDRechit branch not found");
174 return;
175 }
176 branch1->SetAddress(&fRechits);
d270ca46 177
178 Int_t ncrhit = 0;
26f1ae17 179 Int_t nmodules = (Int_t) branch->GetEntries();
180
ecee2a1a 181 AliDebug(1,Form("Number of modules filled in treeR = %d",nmodules));
7dcaf913 182 for (Int_t imodule = 0; imodule < nmodules; imodule++)
183 {
26f1ae17 184 branch->GetEntry(imodule);
7dcaf913 185 Int_t nentries = fRecpoints->GetLast();
ecee2a1a 186 AliDebug(2,Form("Number of clusters per modules filled in treeR = %d"
187 ,nentries));
b0e4d1e1 188
7dcaf913 189 for(Int_t ient = 0; ient < nentries+1; ient++)
190 {
191 fPMDrecpoint = (AliPMDrecpoint1*)fRecpoints->UncheckedAt(ient);
192 idet = fPMDrecpoint->GetDetector();
193 ismn = fPMDrecpoint->GetSMNumber();
194 clusdata[0] = fPMDrecpoint->GetClusX();
195 clusdata[1] = fPMDrecpoint->GetClusY();
196 clusdata[2] = fPMDrecpoint->GetClusADC();
197 clusdata[3] = fPMDrecpoint->GetClusCells();
26f1ae17 198 clusdata[4] = fPMDrecpoint->GetClusSigmaX();
199 clusdata[5] = fPMDrecpoint->GetClusSigmaY();
4b1e197e 200
e6ba3040 201 if (clusdata[4] >= 0. && clusdata[5] >= 0.)
920e13db 202 {
b0e4d1e1 203 // extract the associated cell information
204 branch1->GetEntry(ncrhit);
205 Int_t nenbr1 = fRechits->GetLast() + 1;
206
207 irow = new Int_t[nenbr1];
208 icol = new Int_t[nenbr1];
209 itra = new Int_t[nenbr1];
210 ipid = new Int_t[nenbr1];
211 cadc = new Float_t[nenbr1];
4b1e197e 212
b0e4d1e1 213 for (Int_t ient1 = 0; ient1 < nenbr1; ient1++)
214 {
215 rechit = (AliPMDrechit*)fRechits->UncheckedAt(ient1);
216 //irow[ient1] = rechit->GetCellX();
217 //icol[ient1] = rechit->GetCellY();
218 itra[ient1] = rechit->GetCellTrack();
219 ipid[ient1] = rechit->GetCellPid();
220 cadc[ient1] = rechit->GetCellAdc();
221 }
4b1e197e 222 if (idet == 0)
223 {
224 AssignTrPidToCluster(nenbr1, itra, ipid, cadc,
225 trackno, trackpid);
226 }
227 else if (idet == 1)
228 {
229 trackno = itra[0];
230 trackpid = ipid[0];
231 }
b0e4d1e1 232
233 delete [] irow;
234 delete [] icol;
235 delete [] itra;
236 delete [] ipid;
237 delete [] cadc;
238
239 fPMDclin = new AliPMDrecdata(idet,ismn,trackno,trackpid,clusdata);
920e13db 240 fPMDcontin->Add(fPMDclin);
b0e4d1e1 241
242 ncrhit++;
920e13db 243 }
7dcaf913 244 }
245 }
246
69c38808 247 AliPMDEmpDiscriminator pmddiscriminator;
248 pmddiscriminator.Discrimination(fPMDcontin,fPMDcontout);
7dcaf913 249
32b0cad8 250 // alignment implemention
251
252 Double_t sectr[4][3] = { {0.,0.,0.},{0.,0.,0.},{0.,0.,0.},{0.,0.,0.}};
253 TString snsector="PMD/Sector";
254 TString symname;
255 TGeoHMatrix gpmdor;
256
257 for(Int_t isector=1; isector<=4; isector++)
258 {
259 symname = snsector;
260 symname += isector;
261 TGeoHMatrix *gpmdal = AliGeomManager::GetMatrix(symname);
262 Double_t *tral = gpmdal->GetTranslation();
263
264 AliGeomManager::GetOrigGlobalMatrix(symname, gpmdor);
265 Double_t *tror = gpmdor.GetTranslation();
266
267 for(Int_t ixyz=0; ixyz<3; ixyz++)
268 {
269 sectr[isector-1][ixyz] = tral[ixyz] - tror[ixyz];
270 }
271 }
272
50555ba1 273 const Float_t kzpos = 361.5; // middle of the PMD
01c4d84a 274
e38d8aee 275 Int_t ix = -1, iy = -1;
be8b7039 276 Int_t det = 0, smn = 0, trno = 1, trpid = 0, mstat = 0;
277 Float_t xpos = 0., ypos = 0.;
278 Float_t adc = 0., ncell = 0., radx = 0., rady = 0.;
01c4d84a 279 Float_t xglobal = 0., yglobal = 0., zglobal = 0;
be8b7039 280 Float_t pid = 0.;
7dcaf913 281
32b0cad8 282 fPMDutil->ApplyAlignment(sectr);
7dcaf913 283
284 Int_t nentries2 = fPMDcontout->GetEntries();
ecee2a1a 285 AliDebug(1,Form("Number of clusters coming after discrimination = %d"
286 ,nentries2));
7dcaf913 287 for (Int_t ient1 = 0; ient1 < nentries2; ient1++)
288 {
289 fPMDclout = (AliPMDclupid*)fPMDcontout->UncheckedAt(ient1);
290
291 det = fPMDclout->GetDetector();
292 smn = fPMDclout->GetSMN();
b0e4d1e1 293 trno = fPMDclout->GetClusTrackNo();
294 trpid = fPMDclout->GetClusTrackPid();
295 mstat = fPMDclout->GetClusMatching();
7dcaf913 296 xpos = fPMDclout->GetClusX();
297 ypos = fPMDclout->GetClusY();
298 adc = fPMDclout->GetClusADC();
299 ncell = fPMDclout->GetClusCells();
b0e4d1e1 300 radx = fPMDclout->GetClusSigmaX();
6db43724 301 // Here in the variable "rady" we are keeping the row and col
302 // of the single isolated cluster having ncell = 1 for offline
303 // calibration
304
305 if ((radx > 999. && radx < 1000.) && ncell == 1)
306 {
307 if (smn < 12)
308 {
309 ix = (Int_t) (ypos +0.5);
310 iy = (Int_t) xpos;
311 }
312 else if (smn > 12 && smn < 24)
313 {
314 ix = (Int_t) xpos;
315 iy = (Int_t) (ypos +0.5);
316 }
317 rady = (Float_t) (ix*100 + iy);
318 }
319 else
320 {
321 rady = fPMDclout->GetClusSigmaY();
322 }
7dcaf913 323 pid = fPMDclout->GetClusPID();
324
325 //
7dcaf913 326 /**********************************************************************
327 * det : Detector, 0: PRE & 1:CPV *
01c4d84a 328 * smn : Serial Module Number 0 to 23 for each plane *
7dcaf913 329 * xpos : x-position of the cluster *
330 * ypos : y-position of the cluster *
331 * THESE xpos & ypos are not the true xpos and ypos *
332 * for some of the unit modules. They are rotated. *
333 * adc : ADC contained in the cluster *
334 * ncell : Number of cells contained in the cluster *
335 * rad : radius of the cluster (1d fit) *
7dcaf913 336 **********************************************************************/
337 //
01c4d84a 338
26f1ae17 339
340 if (det == 0)
341 {
dc461f61 342 zglobal = kzpos + 1.65; // PREshower plane
26f1ae17 343 }
344 else if (det == 1)
345 {
dc461f61 346 zglobal = kzpos - 1.65; // CPV plane
26f1ae17 347 }
348
dc461f61 349 fPMDutil->RectGeomCellPos(smn,xpos,ypos,xglobal,yglobal,zglobal);
350
7dcaf913 351 // Fill ESD
352
353 AliESDPmdTrack *esdpmdtr = new AliESDPmdTrack();
354
355 esdpmdtr->SetDetector(det);
c6f4d28b 356 esdpmdtr->SetSmn(smn);
357 esdpmdtr->SetClusterTrackNo(trno);
358 esdpmdtr->SetClusterTrackPid(trpid);
359 esdpmdtr->SetClusterMatching(mstat);
360
26f1ae17 361 esdpmdtr->SetClusterX(xglobal);
362 esdpmdtr->SetClusterY(yglobal);
363 esdpmdtr->SetClusterZ(zglobal);
7dcaf913 364 esdpmdtr->SetClusterADC(adc);
26f1ae17 365 esdpmdtr->SetClusterCells(ncell);
7dcaf913 366 esdpmdtr->SetClusterPID(pid);
c6f4d28b 367 esdpmdtr->SetClusterSigmaX(radx);
368 esdpmdtr->SetClusterSigmaY(rady);
7dcaf913 369
370 event->AddPmdTrack(esdpmdtr);
69c38808 371 delete esdpmdtr;
7dcaf913 372 }
a1e9ac09 373
374 fPMDcontin->Delete();
375 fPMDcontout->Delete();
376
7dcaf913 377}
378//--------------------------------------------------------------------//
b0e4d1e1 379void AliPMDtracker::AssignTrPidToCluster(Int_t nentry, Int_t *itra,
380 Int_t *ipid, Float_t *cadc,
381 Int_t &trackno, Int_t &trackpid)
382{
383 // assign the track number and the corresponding pid to a cluster
384 // split cluster part will be done at the time of calculating eff/pur
385
386 Int_t *phentry = new Int_t [nentry];
a6621fb4 387 Int_t *hadentry = new Int_t [nentry];
b0e4d1e1 388 Int_t *trenergy = 0x0;
a6621fb4 389 Int_t *trpid = 0x0;
b0e4d1e1 390 Int_t *sortcoord = 0x0;
391
392 Int_t ngtrack = 0;
a6621fb4 393 Int_t nhtrack = 0;
b0e4d1e1 394 for (Int_t i = 0; i < nentry; i++)
395 {
396 phentry[i] = -1;
a6621fb4 397 hadentry[i] = -1;
398
b0e4d1e1 399 if (ipid[i] == 22)
400 {
401 phentry[ngtrack] = i;
402 ngtrack++;
403 }
a6621fb4 404 else if (ipid[i] != 22)
405 {
406 hadentry[nhtrack] = i;
407 nhtrack++;
408 }
b0e4d1e1 409 }
4b1e197e 410
a6621fb4 411 Int_t nghadtrack = ngtrack + nhtrack;
412
b0e4d1e1 413 if (ngtrack == 0)
414 {
415 // hadron track
416 // no need of track number, set to -1
417 trackpid = 8;
418 trackno = -1;
419 }
a6621fb4 420 else if (ngtrack >= 1)
b0e4d1e1 421 {
a6621fb4 422 // one or more than one photon track + charged track
423 // find out which track deposits maximum energy and
424 // assign that track number and track pid
425
426 trenergy = new Int_t [nghadtrack];
427 trpid = new Int_t [nghadtrack];
428 sortcoord = new Int_t [nghadtrack];
b0e4d1e1 429 for (Int_t i = 0; i < ngtrack; i++)
430 {
431 trenergy[i] = 0.;
a6621fb4 432 trpid[i] = -1;
b0e4d1e1 433 for (Int_t j = 0; j < nentry; j++)
434 {
435 if (ipid[j] == 22 && itra[j] == itra[phentry[i]])
436 {
437 trenergy[i] += cadc[j];
a6621fb4 438 trpid[i] = 22;
439 }
440 }
441 }
442 for (Int_t i = ngtrack; i < nghadtrack; i++)
443 {
444 trenergy[i] = 0.;
445 trpid[i] = -1;
446 for (Int_t j = 0; j < nentry; j++)
447 {
448 if (ipid[j] != 22 && itra[j] == itra[hadentry[i-ngtrack]])
449 {
450 trenergy[i] += cadc[j];
451 trpid[i] = ipid[j];
b0e4d1e1 452 }
453 }
454 }
4b1e197e 455
b0e4d1e1 456 Bool_t jsort = true;
a6621fb4 457 TMath::Sort(nghadtrack,trenergy,sortcoord,jsort);
4b1e197e 458
b0e4d1e1 459 Int_t gtr = sortcoord[0];
a6621fb4 460 if (trpid[gtr] == 22)
461 {
462 trackpid = 22;
463 trackno = itra[phentry[gtr]]; // highest adc track
464 }
465 else
466 {
467 trackpid = 8;
468 trackno = -1;
469 }
4b1e197e 470
b0e4d1e1 471 delete [] trenergy;
a6621fb4 472 delete [] trpid;
b0e4d1e1 473 delete [] sortcoord;
4b1e197e 474
a6621fb4 475 } // end of ngtrack >= 1
4b1e197e 476
b0e4d1e1 477}
478//--------------------------------------------------------------------//
7dcaf913 479void AliPMDtracker::SetVertex(Double_t vtx[3], Double_t evtx[3])
480{
481 fXvertex = vtx[0];
482 fYvertex = vtx[1];
483 fZvertex = vtx[2];
484 fSigmaX = evtx[0];
485 fSigmaY = evtx[1];
486 fSigmaZ = evtx[2];
487}
488//--------------------------------------------------------------------//
7dcaf913 489void AliPMDtracker::ResetClusters()
490{
491 if (fRecpoints) fRecpoints->Clear();
492}
493//--------------------------------------------------------------------//