Revert "Revert "#103626: Commit DCal geometry to master" since the files are broken."
[u/mrichter/AliRoot.git] / TPC / AliTPCCmpNG.C
CommitLineData
e8558587 1////////////////////////////////////////////////////////////////////////
2//
3// name: AliTPCCmpNG.C
4//
5// date: 24.09.2002
6// author: Jiri Chudoba
7// version: 1.0
8// description:
9// define a class TPCGenTrack
10// save TPC related properties of tracks into a single tree
11//
12// input:
13// Int_t nEvents ... nr of events to process
14// Int_t firstEventNr ... first event number (starts from 0)
15// char* fnRecTracks .. name of file with reconstructed tracks
16// char* fnHits ... name of file with hits and Kine Tree
17// char* fnDigits ... name of file with digits
18// char* fnTracks .. output file name, default genTracks.root
19//
20// How to use:
21// Typical usage:
22// .L AliTPCCmpNG.C+
23// TPCFindGenTracks *t = new TPCFindGenTracks("galice.root","tpc.digits.root")
24// t->Exec();
25// .q
26// aliroot
27// .L AliTPCCmpNG.C+
28// TPCCmpTr *t2 = new TPCCmpTr("tpc.tracks.root","genTracks.root","cmpTracks.root");
29// t2->Exec();
30//
31// Details:
32// ========
33//
34// Step 1 - summurize information from simulation
35// ===============================================
36// Compile macro with ACLIC:
37// .L AliTPCCmpNG.C+
38// create an object TPCFindGenTracks, which processes information
39// from simulations. As input it needs:
40// object gAlice: to get magnetic field
41// TreeK: to get parameters of generated particles
42// TreeTR: to get track parameters at the TPC entry
43// TreeD: to get number of digits and digits pattern
44// for a given track in TPC
45// These input objects can be in different files, gAlice, TreeK and
46// TreeTR are in the file fnHits, TreeD in the file fnDigits (can be
47// the same as fnHits. Output is written to the file fnRes
48// ("genTracks.root" by default). Use can specify number of
49// events to process and the first event number:
50// TPCFindGenTracks *t = new TPCFindGenTracks("galice.root","tpc.digits.root","genTracks.root",1,0)
51// The filenames in the example on previous line are defaults, user can
52// specify just the file name with gAlice object (and TreeTR and TreeK),
53// so equivalent shorter initialization is:
54// TPCFindGenTracks *t = new TPCFindGenTracks("galice.root")
55// The task is done by calling Exec() method:
56// t->Exec();
57// User can set different debug levels by invoking:
58// t->SetDebug(debugLevel)
59// Number of events to process and the first event number can be
60// specified as parameters to Exec:
61// t->Exec(nEvents, firstEvent)
62// Then you have to quit root to get rid of problems with deleting gAlice
63// object (it is not deleted, but read again in the following step):
64//
65// Step 2 - compare reconstructed tracks with simulated
66// ====================================================
67//
68// Load (and compile) the macro:
69// .L AliTPCCmpNG.C+
70// Create object TPCCmpTr, which does the comparison. As input it requires
71// name of the file with reconstructed TPC tracks. You can specify
72// name of the file with genTrack tree (summarized info about simulation),
73// file with gAlice object, output file name, number of events to process
74// and first event number:
75// TPCCmpTr *t2 = new TPCCmpTr("tpc.tracks.root","genTracks.root","cmpTracks.root","galice.root",1,0);
76// The interface is quite similar to the TPCFindGenTracks class.
77// Then just invoke Exec() method:
78// t2->Exec();
79//
80// Step 3 - study the results
81// ==========================
82// Load the outoput TTree and you can do Draw(), Scan() or other
83// usual things to do with TTree:
84// TFile *f = new TFile("cmpTracks.root")
85// TTree *t = (TTree*)f->Get("TPCcmpTracks")
86// t->Draw("fVDist[3]","fReconstructed")
87//
88// History:
89//
90// 24.09.02 - first version
91// 24.01.03 - v7, before change from TPC Special Hits to TrackReferences
92// 26.01.03 - change from TPC Special Hits to TrackReferences
93// (loop over TreeTR instead of TreeH)
94// 28.01.03 - v8 last version before removing TPC special point
95// 28.01.03 - remove TPC special point, loop over TreeH
96// store TParticle and AliTrack
97// 29.01.03 - v9 last version before moving the loop over rec. tracks
98// into separate step
99// 03.02.03 - rename to AliTPCCmpNG.C, remove the part with rec. tracks
100// (will be addded in a macro AliTPCCmpTr.C
101//
102//
103////////////////////////////////////////////////////////////////////////
104
105#if !defined(__CINT__) || defined(__MAKECINT__)
106#include "iostream.h"
107#include <stdio.h>
108#include <string.h>
109#include "Rtypes.h"
110#include "TFile.h"
111#include "TTree.h"
112#include "TString.h"
113#include "TBenchmark.h"
114#include "TStopwatch.h"
115#include "TParticle.h"
116#include "AliRun.h"
117#include "AliStack.h"
118#include "AliTPCtrack.h"
119#include "AliSimDigits.h"
120#include "AliTPCParam.h"
121#include "TParticle.h"
122#include "AliTPC.h"
123#include "AliDetector.h"
124#include "AliTrackReference.h"
125#include "TSystem.h"
126#include "TTimer.h"
127// #include "AliConst.h"
128#endif
129
130#include "AliTPCTracking.C"
131
132// include ML.C:
133////////////////////////////////////////////////////////////////////////
134//
135// name: ML.C (Macro Library - collection of functions used in diff macros
136// date: 08.05.2002
137// last update: 08.05.2002
138// author: Jiri Chudoba
139// version: 1.0
140// description:
141//
142// History:
143//
144// 08.05.02 - first version
145//
146////////////////////////////////////////////////////////////////////////
147
148#if !defined(__CINT__) || defined(__MAKECINT__)
149#include "iostream.h"
150#include "Rtypes.h"
151#include "TSystem.h"
152#include "TTimer.h"
153#include "Getline.h"
154#include "TChain.h"
155#include "TString.h"
156#include "TFile.h"
157#include "AliRun.h"
158
159void WaitForReturn();
160Bool_t ImportgAlice(TFile *file);
161TFile* OpenAliceFile(char *fn, Bool_t importgAlice = kFALSE, char *mode = "read");
162Int_t Chain(TString baseDir, TString subDirNameMask, TString fn, TChain& chain);
163
164#endif
165
166////////////////////////////////////////////////////////////////////////
167void WaitForReturn()
168{
169//
170// wait until user press return;
171//
172 char *input;
173 Bool_t done = kFALSE;
174 TTimer *timer = new TTimer("gSystem->ProcessEvents();", 50, kFALSE);
175
176 do {
177 timer->TurnOn();
178 timer->Reset();
179 // Now let's read the input, we can use here any
180 // stdio or iostream reading methods. like std::cin >> myinputl;
181 input = Getline("Type <return> to continue: ");
182 timer->TurnOff();
183 if (input) done = kTRUE;
184 } while (!done);
185}
186
187////////////////////////////////////////////////////////////////////////
188
189Int_t Chain(TString baseDir, TString subDirNameMask, TString fn, TChain& chain) {
190// chain all files fn in the subdirectories which match subDirName
191// return number of chained files
192
193// open baseDir, loop over subdirs
194
195 void *dirp = gSystem->OpenDirectory(baseDir);
196 if (!dirp) {
197 cerr<<"Could not open base directory "<<baseDir.Data()<<endl;
198 return 0;
199 }
200 const char *afile;
201 Long_t id, size, flag, modTime;
202 Int_t rc = 0;
203 char relName[100];
204 while((afile = gSystem->GetDirEntry(dirp))) {
205// printf("file: %s\n",afile);
206 if (strstr(afile,subDirNameMask.Data())) {
207 sprintf(relName,"%s/%s/%s",baseDir.Data(),afile,fn.Data());
208// cerr<<"relName = "<<relName<<endl;
209 if(!gSystem->GetPathInfo(relName, &id, &size, &flag, &modTime)) {
210 rc += chain.Add(relName);
211 }
212 }
213 }
214 gSystem->FreeDirectory(dirp);
215// cerr<<"rc = "<<rc<<endl;
216 return rc;
217}
218////////////////////////////////////////////////////////////////////////
219Bool_t ImportgAlice(TFile *file) {
220// read in gAlice object from the file
221 gAlice = (AliRun*)file->Get("gAlice");
222 if (!gAlice) return kFALSE;
223 return kTRUE;
224}
225////////////////////////////////////////////////////////////////////////
226TFile* OpenAliceFile(char *fn, Bool_t importgAlice, char *mode) {
227 TFile *file = TFile::Open(fn,mode);
228 if (!file->IsOpen()) {
229 cerr<<"OpenAliceFile: cannot open file "<<fn<<" in mode "
230 <<mode<<endl;
231 return 0;
232 }
233 if (!importgAlice) return file;
234 if (ImportgAlice(file)) return file;
235 return 0;
236}
237////////////////////////////////////////////////////////////////////////
238
239
240////////////////////////////////////////////////////////////////////////
241//
242// Start of implementation of the class TPCGenInfo
243//
244////////////////////////////////////////////////////////////////////////
245class TPCGenInfo: public TObject {
246
247public:
248 TPCGenInfo();
249 AliTrackReference *fTrackRef; // track reference saved in the output tree
250 TParticle *fParticle; // generated particle
251 Int_t fLabel; // track label
252
253 Int_t fRowsWithDigitsInn; // number of rows with digits in the inner sectors
254 Int_t fRowsWithDigits; // number of rows with digits in the outer sectors
255 Int_t fRowsTrackLength; // last - first row with digit
256 Int_t fDigitsInSeed; // digits in the default seed rows
257
258 ClassDef(TPCGenInfo,1) // container for
259};
260ClassImp(TPCGenInfo)
261////////////////////////////////////////////////////////////////////////
262TPCGenInfo::TPCGenInfo()
263{
264 fTrackRef = 0;
265 fParticle = 0;
266}
267////////////////////////////////////////////////////////////////////////
268
269////////////////////////////////////////////////////////////////////////
270//
271// End of implementation of the class TPCGenInfo
272//
273////////////////////////////////////////////////////////////////////////
274
275
276
277////////////////////////////////////////////////////////////////////////
278//
279// Start of implementation of the class digitRow
280//
281////////////////////////////////////////////////////////////////////////
282const Int_t kgRowBytes = 32;
283
284class digitRow: public TObject {
285
286public:
287 digitRow();
288// digitRow(){;}
289 virtual ~digitRow(){;}
290 void SetRow(Int_t row);
291 Bool_t TestRow(Int_t row);
292 digitRow & operator=(const digitRow &digOld);
293 Int_t RowsOn(Int_t upto=8*kgRowBytes);
294 Int_t Last();
295 Int_t First();
296 void Reset();
297
298//private:
299 UChar_t fDig[kgRowBytes];
300
301 ClassDef(digitRow,1) // container for digit pattern
302};
303ClassImp(digitRow)
304////////////////////////////////////////////////////////////////////////
305digitRow::digitRow()
306{
307 Reset();
308}
309////////////////////////////////////////////////////////////////////////
310digitRow & digitRow::operator=(const digitRow &digOld)
311{
312 for (Int_t i = 0; i<kgRowBytes; i++) fDig[i] = digOld.fDig[i];
313 return (*this);
314}
315////////////////////////////////////////////////////////////////////////
316void digitRow::SetRow(Int_t row)
317{
318 if (row >= 8*kgRowBytes) {
319 cerr<<"digitRow::SetRow: index "<<row<<" out of bounds."<<endl;
320 return;
321 }
322 Int_t iC = row/8;
323 Int_t iB = row%8;
324 SETBIT(fDig[iC],iB);
325}
326
327////////////////////////////////////////////////////////////////////////
328Bool_t digitRow::TestRow(Int_t row)
329{
330//
331// return kTRUE if row is on
332//
333 Int_t iC = row/8;
334 Int_t iB = row%8;
335 return TESTBIT(fDig[iC],iB);
336}
337////////////////////////////////////////////////////////////////////////
338Int_t digitRow::RowsOn(Int_t upto)
339{
340//
341// returns number of rows with a digit
342// count only rows less equal row number upto
343//
344 Int_t total = 0;
345 for (Int_t i = 0; i<kgRowBytes; i++) {
346 for (Int_t j = 0; j < 8; j++) {
347 if (i*8+j > upto) return total;
348 if (TESTBIT(fDig[i],j)) total++;
349 }
350 }
351 return total;
352}
353////////////////////////////////////////////////////////////////////////
354void digitRow::Reset()
355{
356//
357// resets all rows to zero
358//
359 for (Int_t i = 0; i<kgRowBytes; i++) {
360 fDig[i] <<= 8;
361 }
362}
363////////////////////////////////////////////////////////////////////////
364Int_t digitRow::Last()
365{
366//
367// returns the last row number with a digit
368// returns -1 if now digits
369//
370 for (Int_t i = kgRowBytes-1; i>=0; i--) {
371 for (Int_t j = 7; j >= 0; j--) {
372 if TESTBIT(fDig[i],j) return i*8+j;
373 }
374 }
375 return -1;
376}
377////////////////////////////////////////////////////////////////////////
378Int_t digitRow::First()
379{
380//
381// returns the first row number with a digit
382// returns -1 if now digits
383//
384 for (Int_t i = 0; i<kgRowBytes; i++) {
385 for (Int_t j = 0; j < 8; j++) {
386 if (TESTBIT(fDig[i],j)) return i*8+j;
387 }
388 }
389 return -1;
390}
391
392////////////////////////////////////////////////////////////////////////
393//
394// end of implementation of a class digitRow
395//
396////////////////////////////////////////////////////////////////////////
397////////////////////////////////////////////////////////////////////////
398//
399// Start of implementation of the class TPCFindGenTracks
400//
401////////////////////////////////////////////////////////////////////////
402
403class TPCFindGenTracks {
404
405public:
406 TPCFindGenTracks();
407 TPCFindGenTracks(char* fnHits,
408 char* fnDigits ="tpc.digits.root",
409 char* fnRes ="genTracks.root",
410 Int_t nEvents=1, Int_t firstEvent=0);
411 virtual ~TPCFindGenTracks();
412 void Reset();
413 Int_t Exec();
414 Int_t Exec(Int_t nEvents, Int_t firstEventNr);
415 void CreateTreeGenTracks();
416 void CloseOutputFile();
417 Int_t TreeKLoop();
418 Int_t TreeTRLoop();
419 Int_t TreeDLoop();
420 void SetFirstEventNr(Int_t i) {fFirstEventNr = i;}
421 void SetNEvents(Int_t i) {fNEvents = i;}
422 void SetDebug(Int_t level) {fDebug = level;}
423
424 Float_t TR2LocalX(AliTrackReference *trackRef,
425 AliTPCParam *paramTPC);
426
427public:
428 Int_t fDebug; //! debug flag
429 Int_t fEventNr; //! current event number
430 Int_t fLabel; //! track label
431 Int_t fNEvents; //! number of events to process
432 Int_t fFirstEventNr; //! first event to process
433 Int_t fNParticles; //! number of particles in TreeK
434 TTree *fTreeGenTracks; //! output tree with generated tracks
435 char *fFnRes; //! output file name with stored tracks
436 char *fFnHits; //! input file name with hits
437 char *fFnDigits; //! input file name with digits
438 TFile *fFileGenTracks; //! output file with stored fTreeGenTracks
439 TFile *fFileHits; //! input file with hits
440 TFile *fFileTreeD; //! input file with digits
441 digitRow *fDigitRow; //! pointer to the object saved in Branch
442 digitRow *fContainerDigitRow; //! big container for partial information
443 AliTrackReference *fTrackRef; //! track reference saved in the output tree
444 AliTrackReference *fContainerTR;//! big container for partial information
445 Int_t *fIndexTR; //! index of particle label in the fContainerTR
446 Int_t fLastIndexTR; //! last used index in fIndexTR
447
448 AliTPCParam* fParamTPC; //! AliTPCParam
449
450 Double_t fVPrim[3]; //! primary vertex position
451 Double_t fVDist[4]; //! distance of the particle vertex from primary vertex
452 // the fVDist[3] contains size of the 3-vector
453 TParticle *fParticle; //! generated particle
454
455 Int_t fRowsWithDigitsInn; //! number of rows with digits in the inner sectors
456 Int_t fRowsWithDigits; //! number of rows with digits in the outer sectors
457 Int_t fRowsTrackLength; //! last - first row with digit
458 Int_t fDigitsInSeed; //! digits in the default seed rows
459
460private:
461
462// some constants for the original non-pareller tracking (by Y.Belikov)
463 static const Int_t seedRow11 = 158; // nRowUp - 1
464 static const Int_t seedRow12 = 139; // nRowUp - 1 - (Int_t) 0.125*nRowUp
465 static const Int_t seedRow21 = 149; // seedRow11 - shift
466 static const Int_t seedRow22 = 130; // seedRow12 - shift
d8d3b5b8 467 static const Double_t kRaddeg = 180./TMath::Pi();
e8558587 468
469 static const Int_t fgMaxIndexTR = 50000; // maximum number of tracks with a track ref
470 static const Int_t fgMaxParticles = 2000000; // maximum number of generated particles
471 static const Double_t fgPtCut = .001; // do not store particles with generated pT less than this
472 static const Float_t fgTrackRefLocalXMax = 82.95;
473 static const Float_t fgTrackRefLocalXMaxDelta = 500.;
474
475 ClassDef(TPCFindGenTracks,1) // class which creates and fills tree with TPCGenTrack objects
476};
477ClassImp(TPCFindGenTracks)
478
479////////////////////////////////////////////////////////////////////////
480TPCFindGenTracks::TPCFindGenTracks()
481{
482 Reset();
483}
484
485////////////////////////////////////////////////////////////////////////
486TPCFindGenTracks::TPCFindGenTracks(char* fnHits, char* fnDigits,
487 char* fnRes,
488 Int_t nEvents, Int_t firstEvent)
489{
490 Reset();
491 fFirstEventNr = firstEvent;
492 fEventNr = firstEvent;
493 fNEvents = nEvents;
494 fFnRes = fnRes;
495 fFnHits = fnHits;
496 fFnDigits = fnDigits;
497}
498////////////////////////////////////////////////////////////////////////
499void TPCFindGenTracks::Reset()
500{
501 fDigitRow = 0;
502 fEventNr = 0;
503 fNEvents = 0;
504 fTreeGenTracks = 0;
505 fFnRes = "genTracks.root";
506 fFnHits = "rfio:galice.root";
507 fFnDigits = "rfio:its.tpc.trd.digits.root";
508 fFileGenTracks = 0;
509 fFileHits =0;
510 fFileTreeD =0;
511 fContainerDigitRow = 0;
512 fContainerTR = 0;
513 fIndexTR = 0;
514 fLastIndexTR = -1;
515 fParticle = 0;
516 fTrackRef = 0;
517 fDebug = 0;
518 fVPrim[0] = -1000.;
519 fVPrim[1] = -1000.;
520 fVPrim[2] = -1000.;
521 fParamTPC = 0;
522
523}
524////////////////////////////////////////////////////////////////////////
525TPCFindGenTracks::~TPCFindGenTracks()
526{
527 ;
528}
529////////////////////////////////////////////////////////////////////////
530Int_t TPCFindGenTracks::Exec(Int_t nEvents, Int_t firstEventNr)
531{
532 fNEvents = nEvents;
533 fFirstEventNr = firstEventNr;
534 return Exec();
535}
536
537////////////////////////////////////////////////////////////////////////
538Int_t TPCFindGenTracks::Exec()
539{
540 TStopwatch timer;
541 timer.Start();
542
543 fDigitRow = new digitRow();
544 CreateTreeGenTracks();
545 if (!fTreeGenTracks) return 1;
546 fFileHits = OpenAliceFile(fFnHits,kTRUE,"read"); //gAlice is read here
547 if (!fFileHits) return 1;
548 fFileHits->cd();
549 SetFieldFactor();
550
551 fFileTreeD = TFile::Open(fFnDigits,"read");
552 if (!fFileTreeD->IsOpen()) {
553 cerr<<"Cannot open file "<<fFnDigits<<endl;
554 return 1;
555 }
556
557 fParamTPC = LoadTPCParam(fFileTreeD);
558 if (!fParamTPC) {
559 cerr<<"TPC parameters not found and could not be created"<<endl;
560 return 1;
561 }
562
563 for (fEventNr = fFirstEventNr; fEventNr < fFirstEventNr+fNEvents;
564 fEventNr++) {
565 fNParticles = gAlice->GetEvent(fEventNr);
566 fContainerDigitRow = new digitRow[fNParticles];
567 fContainerTR = new AliTrackReference[fgMaxIndexTR];
568 fIndexTR = new Int_t[fNParticles];
569 for (Int_t i = 0; i<fNParticles; i++) {
570 fIndexTR[i] = -1;
571 }
572
573 cout<<"Start to process event "<<fEventNr<<endl;
574 cout<<"\tfNParticles = "<<fNParticles<<endl;
575 if (fDebug>2) cout<<"\tStart loop over TreeD"<<endl;
576 if (TreeDLoop()>0) return 1;
577 if (fDebug>2) cout<<"\tStart loop over TreeTR"<<endl;
578 if (TreeTRLoop()>0) return 1;
579 if (fDebug>2) cout<<"\tStart loop over TreeK"<<endl;
580 if (TreeKLoop()>0) return 1;
581 if (fDebug>2) cout<<"\tEnd loop over TreeK"<<endl;
582
583 delete [] fContainerDigitRow;
584 delete [] fContainerTR;
585 delete [] fIndexTR;
586 }
587
588 CloseOutputFile();
589
590 cerr<<"Exec finished"<<endl;
591 fFileHits->cd();
592 delete gAlice;
593 fFileHits->Close();
594 delete fFileHits;
595
596 timer.Stop();
597 timer.Print();
598 return 0;
599}
600////////////////////////////////////////////////////////////////////////
601void TPCFindGenTracks::CreateTreeGenTracks()
602{
603 fFileGenTracks = TFile::Open(fFnRes,"RECREATE");
604 if (!fFileGenTracks) {
605 cerr<<"Error in CreateTreeGenTracks: cannot open file "<<fFnRes<<endl;
606 return;
607 }
608 fTreeGenTracks = new TTree("genTracksTree","genTracksTree");
609 TBranch *branchBits = fTreeGenTracks->Branch("bitsBranch", "digitRow", &fDigitRow, 4000, 0);
610 if (!branchBits) {
611 cerr<<"Error in CreateTreeGenTracks: cannot create branch."<<endl;
612 return;
613 }
614 fTreeGenTracks->Branch("fEventNr",&fEventNr,"fEventNr/I");
615 fTreeGenTracks->Branch("fLabel",&fLabel,"fLabel/I");
616 fTreeGenTracks->Branch("fRowsWithDigitsInn",&fRowsWithDigitsInn,"fRowsWithDigitsInn/I");
617 fTreeGenTracks->Branch("fRowsWithDigits",&fRowsWithDigits,"fRowsWithDigits/I");
618 fTreeGenTracks->Branch("fRowsTrackLength",&fRowsTrackLength,"fRowsTrackLength/I");
619 fTreeGenTracks->Branch("fDigitsInSeed",&fDigitsInSeed,"fDigitsInSeed/I");
620
621 fTreeGenTracks->Branch("Particle","TParticle",&fParticle);
622 fTreeGenTracks->Branch("fVDist",&fVDist,"fVDist[4]/D");
623 fTreeGenTracks->Branch("TR","AliTrackReference",&fTrackRef);
624
625 fTreeGenTracks->AutoSave();
626}
627////////////////////////////////////////////////////////////////////////
628void TPCFindGenTracks::CloseOutputFile()
629{
630 if (!fFileGenTracks) {
631 cerr<<"File "<<fFnRes<<" not found as an open file."<<endl;
632 return;
633 }
634 fFileGenTracks->cd();
635 fTreeGenTracks->Write();
636 delete fTreeGenTracks;
637 fFileGenTracks->Close();
638 delete fFileGenTracks;
639 return;
640}
641
642////////////////////////////////////////////////////////////////////////
643Int_t TPCFindGenTracks::TreeKLoop()
644{
645//
646// open the file with treeK
647// loop over all entries there and save information about some tracks
648//
649
650 fFileHits->cd();
651 if (fDebug > 0) {
652 cout<<"There are "<<fNParticles<<" primary and secondary particles in event "
653 <<fEventNr<<endl;
654 }
655 AliStack * stack = gAlice->Stack();
656 if (!stack) {cerr<<"Stack was not found!\n"; return 1;}
657
658// not all generators give primary vertex position. Take the vertex
659// of the particle 0 as primary vertex.
660 fParticle = stack->ParticleFromTreeK(0);
661 fVPrim[0] = fParticle->Vx();
662 fVPrim[1] = fParticle->Vy();
663 fVPrim[2] = fParticle->Vz();
664
665 for (Int_t iParticle = 0; iParticle < fNParticles; iParticle++) {
666// for (Int_t iParticle = 0; iParticle < fDebug; iParticle++) {
667
668// load only particles with TR
669 if (fIndexTR[iParticle] < 0) continue;
670 fParticle = stack->ParticleFromTreeK(iParticle);
671 if (fDebug > 3 && iParticle < 10) {
672 cout<<"processing particle "<<iParticle<<" ";
673 fParticle->Print();
674 }
675
676// fill the tree
677
678 fLabel = iParticle;
679 fVDist[0] = fParticle->Vx()-fVPrim[0];
680 fVDist[1] = fParticle->Vy()-fVPrim[1];
681 fVDist[2] = fParticle->Vz()-fVPrim[2];
682 fVDist[3] = TMath::Sqrt(fVDist[0]*fVDist[0]+fVDist[1]*fVDist[1]+fVDist[2]*fVDist[2]);
683 fDigitRow = &(fContainerDigitRow[iParticle]);
684 fRowsWithDigitsInn = fDigitRow->RowsOn(63); // 63 = number of inner rows
685 fRowsWithDigits = fDigitRow->RowsOn();
686 fRowsTrackLength = fDigitRow->Last() - fDigitRow->First();
687 if (fDebug > 2 && iParticle < 10) {
688 cerr<<"Fill track with a label "<<iParticle<<endl;
689 }
690 fDigitsInSeed = 0;
691 if (fDigitRow->TestRow(seedRow11) && fDigitRow->TestRow(seedRow12))
692 fDigitsInSeed = 1;
693 if (fDigitRow->TestRow(seedRow21) && fDigitRow->TestRow(seedRow22))
694 fDigitsInSeed += 10;
695
696 if (fIndexTR[iParticle] >= 0) {
697 fTrackRef = &(fContainerTR[fIndexTR[iParticle]]);
698// cerr<<"Debug: fTrackRef->X() = "<<fTrackRef->X()<<endl;
699 } else {
700 fTrackRef->SetTrack(-1);
701 }
702
703 fTreeGenTracks->Fill();
704
705 }
706 fTreeGenTracks->AutoSave();
707// delete gAlice; gAlice = 0;
708// fFileHits->Close();
709
710 if (fDebug > 2) cerr<<"end of TreeKLoop"<<endl;
711
712 return 0;
713}
714////////////////////////////////////////////////////////////////////////
715Int_t TPCFindGenTracks::TreeDLoop()
716{
717//
718// open the file with treeD
719// loop over all entries there and save information about some tracks
720//
721
722
723// Int_t nrow_up=fParamTPC->GetNRowUp();
724// Int_t nrows=fParamTPC->GetNRowLow()+nrow_up;
725 Int_t nInnerSector = fParamTPC->GetNInnerSector();
726 Int_t rowShift = 0;
727 Int_t zero=fParamTPC->GetZeroSup();
728// Int_t gap=Int_t(0.125*nrows), shift=Int_t(0.5*gap);
729
730 char treeDName[100];
731 sprintf(treeDName,"TreeD_75x40_100x60_150x60_%d",fEventNr);
732 TTree *treeD=(TTree*)fFileTreeD->Get(treeDName);
733 AliSimDigits digitsAddress, *digits=&digitsAddress;
734 treeD->GetBranch("Segment")->SetAddress(&digits);
735
736 Int_t sectorsByRows=(Int_t)treeD->GetEntries();
737 if (fDebug > 1) cout<<"\tsectorsByRows = "<<sectorsByRows<<endl;
738 for (Int_t i=0; i<sectorsByRows; i++) {
739// for (Int_t i=5720; i<sectorsByRows; i++) {
740 if (!treeD->GetEvent(i)) continue;
741 Int_t sec,row;
742 fParamTPC->AdjustSectorRow(digits->GetID(),sec,row);
743 if (fDebug > 1) cout<<sec<<' '<<row<<" \r";
744// cerr<<sec<<' '<<row<<endl;
745
746// here I expect that upper sectors follow lower sectors
747 if (sec > nInnerSector) rowShift = fParamTPC->GetNRowLow();
748 digits->First();
749 do {
750 Int_t iRow=digits->CurrentRow();
751 Int_t iColumn=digits->CurrentColumn();
752 Short_t digitValue = digits->CurrentDigit();
753// cout<<"Inner loop: sector, iRow, iColumn "
754// <<sec<<" "<<iRow<<" "<<iColumn<<endl;
755 if (digitValue >= zero) {
756 Int_t label;
757 for (Int_t j = 0; j<3; j++) {
758 label = digits->GetTrackID(iRow,iColumn,j);
759 if (label >= fNParticles) {
760 cerr<<"particle label too big: fNParticles, label "
761 <<fNParticles<<" "<<label<<endl;
762 }
763 if (label >= 0 && label <= fNParticles) {
764// if (label >= 0 && label <= fDebug) {
765 if (fDebug > 6 ) {
766 cout<<"Inner loop: sector, iRow, iColumn, label, value, row "
767 <<sec<<" "
768 <<iRow<<" "<<iColumn<<" "<<label<<" "<<digitValue
769 <<" "<<row<<endl;
770 }
771 fContainerDigitRow[label].SetRow(row+rowShift);
772 }
773 }
774 }
775 } while (digits->Next());
776 }
777
778 if (fDebug > 2) cerr<<"end of TreeDLoop"<<endl;
779
780 return 0;
781}
782
783////////////////////////////////////////////////////////////////////////
784Int_t TPCFindGenTracks::TreeTRLoop()
785{
786//
787// loop over TrackReferences and store the first one for each track
788//
789
790 TTree *treeTR=gAlice->TreeTR();
791 if (!treeTR) {
792 cerr<<"TreeTR not found"<<endl;
793 return 1;
794 }
795 Int_t nPrimaries = (Int_t) treeTR->GetEntries();
796 if (fDebug > 1) cout<<"There are "<<nPrimaries<<" entries in TreeTR"<<endl;
797 TBranch *TPCBranchTR = treeTR->GetBranch("TPC");
798 if (!TPCBranchTR) {
799 cerr<<"TPC branch in TR not found"<<endl;
800 return 1;
801 }
802 TClonesArray* TPCArrayTR = new TClonesArray("AliTrackReference");
803 TPCBranchTR->SetAddress(&TPCArrayTR);
804 for (Int_t iPrimPart = 0; iPrimPart<nPrimaries; iPrimPart++) {
805 TPCBranchTR->GetEntry(iPrimPart);
806 for (Int_t iTrackRef = 0; iTrackRef < TPCArrayTR->GetEntriesFast(); iTrackRef++) {
807 AliTrackReference *trackRef = (AliTrackReference*)TPCArrayTR->At(iTrackRef);
808 Int_t label = trackRef->GetTrack();
809
810// save info in the fContainerTR
811 if (label<0 || label > fNParticles) {
812 cerr<<"Wrong label: "<<label<<endl;
813 continue;
814// return 1;
815 }
816
817// store only if we do not have any track reference yet for this label
818 if (fIndexTR[label] >= 0) continue;
819
820// store only references with localX < fgTrackRefLocalXMax +- fgTrackRefLocalXMaxDelta
821// and the pT > fgPtCut
822 Float_t localX = TR2LocalX(trackRef,fParamTPC);
823 if (TMath::Abs(localX-fgTrackRefLocalXMax)>fgTrackRefLocalXMaxDelta) continue;
824 if (trackRef->Pt() < fgPtCut) continue;
825
826
827// cout<<"label, xg "<<label<<" "<<xg<<endl;
828 if (fLastIndexTR >= fgMaxIndexTR-1) {
829 cerr<<"Too many tracks with track reference. Increase the constant"
830 <<" fgMaxIndexTR"<<endl;
831 return 1;
832 }
833 fIndexTR[label] = ++fLastIndexTR;
834
835// someone was lazy to implement copy ctor in AliTrackReference, so we have to do
836// it here "manually"
837 fContainerTR[fIndexTR[label]].SetPosition(trackRef->X(),trackRef->Y(),trackRef->Z());
838
839 fContainerTR[fIndexTR[label]].SetMomentum(trackRef->Px(),trackRef->Py(),trackRef->Pz());
840 fContainerTR[fIndexTR[label]].SetTrack(trackRef->GetTrack());
841 fContainerTR[fIndexTR[label]].SetLength(trackRef->GetLength());
842// cerr<<"Debug: trackRef->X(), stored: "<<trackRef->X()<<" "
843// << fContainerTR[fIndexTR[label]].X()<<endl;
844
845 }
846 }
847 return 0;
848}
849////////////////////////////////////////////////////////////////////////
850Float_t TPCFindGenTracks::TR2LocalX(AliTrackReference *trackRef,
851 AliTPCParam *paramTPC) {
852
853 Float_t x[3] = { trackRef->X(),trackRef->Y(),trackRef->Z()};
854 Int_t index[4];
855 paramTPC->Transform0to1(x,index);
856 paramTPC->Transform1to2(x,index);
857 return x[0];
858}
859////////////////////////////////////////////////////////////////////////
860
861
862
863////////////////////////////////////////////////////////////////////////
864//
865// Start of implementation of the class TPCCmpTr
866//
867////////////////////////////////////////////////////////////////////////
868
869class TPCCmpTr {
870
871public:
872 TPCCmpTr();
873 TPCCmpTr(char* fnRecTracks,
874 char* fnGenTracks ="genTracks.root",
875 char* fnCmpRes ="cmpTracks.root",
876 char* fnGalice ="galice.root",
877 Int_t nEvents=1, Int_t firstEvent=0);
878 virtual ~TPCCmpTr();
879 void Reset();
880 Int_t Exec();
881 Int_t Exec(Int_t nEvents, Int_t firstEventNr);
882 void CreateTreeCmp();
883 void CloseOutputFile();
884 Bool_t ConnectGenTree();
885 Int_t TreeGenLoop(Int_t eventNr);
886 Int_t TreeTLoop(Int_t eventNr);
887 void SetFirstEventNr(Int_t i) {fFirstEventNr = i;}
888 void SetNEvents(Int_t i) {fNEvents = i;}
889 void SetDebug(Int_t level) {fDebug = level;}
890
891// tmp method, should go to TrackReferenceTPC
892 Float_t TR2LocalX(AliTrackReference *trackRef,
893 AliTPCParam *paramTPC);
894
895private:
896 digitRow *fDigitRow; //! pointer to the object saved in Branch
897 Int_t fEventNr; //! current event number
898 Int_t fLabel; //! track label
899 Int_t fNEvents; //! number of events to process
900 Int_t fFirstEventNr; //! first event to process
901
902 char *fFnCmp; //! output file name with cmp tracks
903 TFile *fFileCmp; //! output file with cmp tracks
904 TTree *fTreeCmp; //! output tree with cmp tracks
905
906 char *fFnGenTracks; //! input file name with gen tracks
907 TFile *fFileGenTracks;
908 TTree *fTreeGenTracks;
909
910 char *fFnHits; //! input file name with gAlice object (needed for B)
911 TFile *fFileHits; //! input file with gAlice
912
913 char *fFnRecTracks; //! input file name with tpc rec. tracks
914 TFile *fFileRecTracks; //! input file with reconstructed tracks
915 TTree *fTreeRecTracks; //! tree with reconstructed tracks
916
917 AliTPCtrack *fTPCTrack; //! pointer to TPC track to connect branch
918 Int_t *fIndexRecTracks; //! index of particle label in the TreeT_TPC
919
920 Int_t fRowsWithDigitsInn; //! number of rows with digits in the inner sectors
921 Int_t fRowsWithDigits; //! number of rows with digits in the outer sectors
922 Int_t fRowsTrackLength; //! last - first row with digit
923 Int_t fDigitsInSeed; //! digits in the default seed rows
924 TParticle *fParticle; //! generated particle
925 Double_t fVDist[4]; //! distance of the particle vertex from primary vertex
926 // the fVDist[3] contains size of the 3-vector
927 AliTrackReference *fTrackRef; //! track reference saved in the output tree
928 Int_t fReconstructed; //! flag if track was reconstructed
929 AliTPCParam* fParamTPC; //! AliTPCParam
930
931 Int_t fNParticles; //! number of particles in the input tree genTracks
932 Int_t fDebug; //! debug flag
933
934 Int_t fNextTreeGenEntryToRead; //! last entry already read from genTracks tree
935 TPCGenInfo *fGenInfo; //! container for all the details
936
937 Double_t fRecPhi; // reconstructed phi angle (0;2*kPI)
938 Double_t fLambda; // reconstructed
939 Double_t fRecPt_1; // reconstructed
940 Float_t fdEdx; // reconstructed dEdx
941
942
943 ClassDef(TPCCmpTr,1) // class which creates and fills tree with TPCGenTrack objects
944};
945ClassImp(TPCCmpTr)
946
947////////////////////////////////////////////////////////////////////////
948TPCCmpTr::TPCCmpTr()
949{
950 Reset();
951}
952
953////////////////////////////////////////////////////////////////////////
954TPCCmpTr::TPCCmpTr(char* fnRecTracks,
955 char* fnGenTracks,
956 char* fnCmp,
957 char* fnGalice,
958 Int_t nEvents, Int_t firstEvent)
959{
960 Reset();
961 fFnRecTracks = fnRecTracks;
962 fFnGenTracks = fnGenTracks;
963 fFnCmp = fnCmp;
964 fFnHits = fnGalice;
965 fFirstEventNr = firstEvent;
966 fEventNr = firstEvent;
967 fNEvents = nEvents;
968}
969////////////////////////////////////////////////////////////////////////
970void TPCCmpTr::Reset()
971{
972 fDigitRow = 0;
973 fEventNr = 0;
974 fNEvents = 0;
975 fTreeCmp = 0;
976 fFnCmp = "cmpTracks.root";
977 fFnHits = "galice.root";
978 fFileGenTracks = 0;
979 fFileHits =0;
980 fParticle = 0;
981 fTrackRef = 0;
982
983 fRowsWithDigitsInn = 0;
984 fRowsWithDigits = 0;
985 fRowsTrackLength = 0;
986 fDigitsInSeed = 0;
987
988 fDebug = 0;
989
990 fParamTPC = 0;
991 fFnRecTracks = "tpc.tracks.root";
992 fTreeRecTracks = 0;
993 fFileRecTracks = 0;
994 fTPCTrack = 0;
995 fGenInfo = 0;
996}
997////////////////////////////////////////////////////////////////////////
998TPCCmpTr::~TPCCmpTr()
999{
1000 ;
1001}
1002////////////////////////////////////////////////////////////////////////
1003Int_t TPCCmpTr::Exec(Int_t nEvents, Int_t firstEventNr)
1004{
1005 fNEvents = nEvents;
1006 fFirstEventNr = firstEventNr;
1007 return Exec();
1008}
1009
1010////////////////////////////////////////////////////////////////////////
1011Int_t TPCCmpTr::Exec()
1012{
1013 TStopwatch timer;
1014 timer.Start();
1015
1016 fDigitRow = new digitRow();
1017 CreateTreeCmp();
1018 if (!fTreeCmp) {
1019 cerr<<"output tree not created"<<endl;
1020 return 1;
1021 }
1022 fFileHits = OpenAliceFile(fFnHits,kTRUE,"read"); //gAlice is read here
1023 if (!fFileHits) {
1024 cerr<<"Cannot open file with gAlice object "<<fFnHits<<endl;
1025 return 1;
1026 }
1027 fFileHits->cd();
1028 SetFieldFactor();
1029
1030 fParamTPC = LoadTPCParam(fFileHits);
1031 if (!fParamTPC) {
1032 cerr<<"TPC parameters not found and could not be created"<<endl;
1033 return 1;
1034 }
1035
1036 if (!ConnectGenTree()) {
1037 cerr<<"Cannot connect tree with generated tracks"<<endl;
1038 return 1;
1039 }
1040 fFileHits->cd();
1041 fNextTreeGenEntryToRead = 0;
1042 cerr<<"fFirstEventNr, fNEvents: "<<fFirstEventNr<<" "<<fNEvents<<endl;
1043 for (Int_t eventNr = fFirstEventNr; eventNr < fFirstEventNr+fNEvents;
1044 eventNr++) {
1045 fNParticles = gAlice->GetEvent(fEventNr);
1046 fIndexRecTracks = new Int_t[fNParticles];
1047 for (Int_t i = 0; i<fNParticles; i++) {
1048 fIndexRecTracks[i] = -1;
1049 }
1050
1051 cout<<"Start to process event "<<fEventNr<<endl;
1052 cout<<"\tfNParticles = "<<fNParticles<<endl;
1053 if (fDebug>2) cout<<"\tStart loop over TreeT"<<endl;
1054 if (TreeTLoop(eventNr)>0) return 1;
1055
1056 if (fDebug>2) cout<<"\tStart loop over tree genTracks"<<endl;
1057 if (TreeGenLoop(eventNr)>0) return 1;
1058 if (fDebug>2) cout<<"\tEnd loop over tree genTracks"<<endl;
1059
1060 delete fTreeRecTracks;
1061
1062 delete [] fIndexRecTracks;
1063 }
1064
1065 CloseOutputFile();
1066
1067 cerr<<"Exec finished"<<endl;
1068 fFileHits->cd();
1069 delete gAlice;
1070 fFileHits->Close();
1071 delete fFileHits;
1072
1073 timer.Stop();
1074 timer.Print();
1075 return 0;
1076}
1077////////////////////////////////////////////////////////////////////////
1078Bool_t TPCCmpTr::ConnectGenTree()
1079{
1080//
1081// connect all branches from the genTracksTree
1082// use the same variables as for the new cmp tree, it may work
1083//
1084 fFileGenTracks = TFile::Open(fFnGenTracks,"READ");
1085 if (!fFileGenTracks) {
1086 cerr<<"Error in ConnectGenTree: cannot open file "<<fFnGenTracks<<endl;
1087 return kFALSE;
1088 }
1089 fTreeGenTracks = (TTree*)fFileGenTracks->Get("genTracksTree");
1090 if (!fTreeGenTracks) {
1091 cerr<<"Error in ConnectGenTree: cannot find genTracksTree in the file "
1092 <<fFnGenTracks<<endl;
1093 return kFALSE;
1094 }
1095 fTreeGenTracks->SetBranchAddress("fEventNr",&fEventNr);
1096 fTreeGenTracks->SetBranchAddress("fLabel",&fLabel);
1097 fTreeGenTracks->SetBranchAddress("fRowsWithDigitsInn",&fRowsWithDigitsInn);
1098 fTreeGenTracks->SetBranchAddress("fRowsWithDigits",&fRowsWithDigits);
1099 fTreeGenTracks->SetBranchAddress("fRowsTrackLength",&fRowsTrackLength);
1100 fTreeGenTracks->SetBranchAddress("fDigitsInSeed",&fDigitsInSeed);
1101 fTreeGenTracks->SetBranchAddress("Particle",&fParticle);
1102 fTreeGenTracks->SetBranchAddress("fVDist",fVDist);
1103 fTreeGenTracks->SetBranchAddress("TR",&fTrackRef);
1104
1105 if (fDebug > 1) {
1106 cout<<"Number of gen. tracks with TR: "<<fTreeGenTracks->GetEntries()<<endl;
1107 }
1108 return kTRUE;
1109}
1110
1111
1112////////////////////////////////////////////////////////////////////////
1113void TPCCmpTr::CreateTreeCmp()
1114{
1115 fFileCmp = TFile::Open(fFnCmp,"RECREATE");
1116 if (!fFileCmp) {
1117 cerr<<"Error in CreateTreeCmp: cannot open file "<<fFnCmp<<endl;
1118 return;
1119 }
1120 fTreeCmp = new TTree("TPCcmpTracks","TPCcmpTracks");
1121 TBranch *branchBits = fTreeCmp->Branch("bitsBranch", "digitRow", &fDigitRow, 4000, 0);
1122 if (!branchBits) {
1123 cerr<<"Error in CreateTreeCmp: cannot create branch."<<endl;
1124 return;
1125 }
1126 fTreeCmp->Branch("fEventNr",&fEventNr,"fEventNr/I");
1127 fTreeCmp->Branch("fLabel",&fLabel,"fLabel/I");
1128 fTreeCmp->Branch("fRowsWithDigitsInn",&fRowsWithDigitsInn,"fRowsWithDigitsInn/I");
1129 fTreeCmp->Branch("fRowsWithDigits",&fRowsWithDigits,"fRowsWithDigits/I");
1130 fTreeCmp->Branch("fRowsTrackLength",&fRowsTrackLength,"fRowsTrackLength/I");
1131 fTreeCmp->Branch("fDigitsInSeed",&fDigitsInSeed,"fDigitsInSeed/I");
1132
1133 fTreeCmp->Branch("fReconstructed",&fReconstructed,"fReconstructed/I");
1134 fTreeCmp->Branch("fTPCTrack","AliTPCtrack",&fTPCTrack);
1135
1136 fTreeCmp->Branch("Particle","TParticle",&fParticle);
1137 fTreeCmp->Branch("fVDist",&fVDist,"fVDist[4]/D");
1138 fTreeCmp->Branch("TR","AliTrackReference",&fTrackRef);
1139
1140 fTreeCmp->AutoSave();
1141}
1142////////////////////////////////////////////////////////////////////////
1143void TPCCmpTr::CloseOutputFile()
1144{
1145 if (!fFileCmp) {
1146 cerr<<"File "<<fFnCmp<<" not found as an open file."<<endl;
1147 return;
1148 }
1149 fFileCmp->cd();
1150 fTreeCmp->Write();
1151 delete fTreeCmp;
1152 fFileCmp->Close();
1153 delete fFileCmp;
1154 return;
1155}
1156////////////////////////////////////////////////////////////////////////
1157
1158Float_t TPCCmpTr::TR2LocalX(AliTrackReference *trackRef,
1159 AliTPCParam *paramTPC) {
1160
1161 Float_t x[3] = { trackRef->X(),trackRef->Y(),trackRef->Z()};
1162 Int_t index[4];
1163 paramTPC->Transform0to1(x,index);
1164 paramTPC->Transform1to2(x,index);
1165 return x[0];
1166}
1167////////////////////////////////////////////////////////////////////////
1168
1169Int_t TPCCmpTr::TreeTLoop(Int_t eventNr)
1170{
1171//
1172// loop over all TPC reconstructed tracks and store info in memory
1173//
1174
1175
1176 if (!fFileRecTracks) fFileRecTracks = TFile::Open(fFnRecTracks,"read");
1177 if (!fFileRecTracks->IsOpen()) {
1178 cerr<<"Cannot open file "<<fFnRecTracks<<endl;
1179 return 1;
1180 }
1181
1182 char treeNameBase[11] = "TreeT_TPC_";
1183 char treeName[20];
1184 sprintf(treeName,"%s%d",treeNameBase,eventNr);
1185
1186 fTreeRecTracks=(TTree*)fFileRecTracks->Get(treeName);
1187 if (!fTreeRecTracks) {
1188 cerr<<"Can't get a tree with TPC rec. tracks named "<<treeName<<endl;
1189 return 1;
1190 }
1191
1192 Int_t nEntries = (Int_t) fTreeRecTracks->GetEntries();
1193 if (fDebug > 2) cout<<"Event, rec. tracks: "<<eventNr<<" "
1194 <<nEntries<<endl;
1195 TBranch * br= fTreeRecTracks->GetBranch("tracks");
1196 br->SetAddress(&fTPCTrack);
1197
1198 for (Int_t iEntry=0; iEntry<nEntries;iEntry++){
1199 br->GetEntry(iEntry);
1200 Int_t label = fTPCTrack->GetLabel();
1201 fIndexRecTracks[label] = iEntry;
1202 }
1203
1204 if (fDebug > 2) cerr<<"end of TreeTLoop"<<endl;
1205
1206 return 0;
1207}
1208////////////////////////////////////////////////////////////////////////
1209Int_t TPCCmpTr::TreeGenLoop(Int_t eventNr)
1210{
1211//
1212// loop over all entries for a given event, find corresponding
1213// rec. track and store in the fTreeCmp
1214//
1215
1216 fFileGenTracks->cd();
1217 Int_t entry = fNextTreeGenEntryToRead;
1218 Double_t nParticlesTR = fTreeGenTracks->GetEntriesFast();
1219 cerr<<"fNParticles, nParticlesTR, fNextTreeGenEntryToRead: "<<fNParticles<<" "
1220 <<nParticlesTR<<" "<<fNextTreeGenEntryToRead<<endl;
1221 while (entry < nParticlesTR) {
1222 fTreeGenTracks->GetEntry(entry);
1223 entry++;
1224 if (fEventNr < eventNr) continue;
1225 if (fEventNr > eventNr) break;
1226 fNextTreeGenEntryToRead = entry-1;
1227 if (fDebug > 2 && fLabel < 10) {
1228 cerr<<"Fill track with a label "<<fLabel<<endl;
1229 }
1230
1231 fReconstructed = 0;
1232 fdEdx = 0.;
1233 fRecPhi = fLambda = fRecPt_1 = 0.;
1234
1235 if (fDebug > 2) {
1236 cerr<<"fLabel, fIndexRecTracks[fLabel] "<<fLabel<<" "<<fIndexRecTracks[fLabel]<<endl;
1237 }
1238 if (fIndexRecTracks[fLabel] >= 0) {
1239 Int_t nBytes = fTreeRecTracks->GetEvent(fIndexRecTracks[fLabel]);
1240 if (nBytes > 0) {
1241 fReconstructed = 1;
1242 fdEdx = fTPCTrack->GetdEdx();
1243 Double_t Localx = TR2LocalX(fTrackRef,fParamTPC);
1244 if (fDebug > 3) {
1245 cerr<<"Track local X before prolongation: "<<fTPCTrack->GetX()<<endl;
1246 }
1247 fTPCTrack->PropagateTo(Localx);
1248 Double_t par[5];
1249 if (fDebug > 3) {
1250 cerr<<"Track local X after prolongation: "<<fTPCTrack->GetX()<<endl;
1251 }
1252 fTPCTrack->GetExternalParameters(Localx,par);
1253 fRecPhi=TMath::ASin(par[2]) + fTPCTrack->GetAlpha();
d8d3b5b8 1254 if (fRecPhi<0) fRecPhi+=2*TMath::Pi();
1255 if (fRecPhi>=2*TMath::Pi()) fRecPhi-=2*TMath::Pi();
e8558587 1256// fRecPhi = (fRecPhi)*kRaddeg;
1257 fLambda = TMath::ATan(par[3]);
1258 fRecPt_1 = TMath::Abs(par[4]);
1259 }
1260 }
1261
1262
1263 fTreeCmp->Fill();
1264
1265 }
1266 fTreeCmp->AutoSave();
1267// delete gAlice; gAlice = 0;
1268// fFileHits->Close();
1269
1270 if (fDebug > 2) cerr<<"end of TreeKLoop"<<endl;
1271
1272 return 0;
1273}
1274////////////////////////////////////////////////////////////////////////