New version of FMD code (A.Maevskaia)
[u/mrichter/AliRoot.git] / FMD / TestReconstruction.C
CommitLineData
37c55dc0 1void TestReconstruction (Int_t vol=1, const Int_t nRings=128, const Int_t nSectors=20)
dc8af42e 2{
37c55dc0 3 // Dynamically link some shared libs
dc8af42e 4 if (gClassTable->GetID("AliRun") < 0)
5 {
37c55dc0 6 gROOT->LoadMacro("loadlibs.C");
dc8af42e 7 loadlibs();
8 }
37c55dc0 9 // Connect the Root Galice file containing Geometry, Kine and Hits
10 char filename[]="galice.root";
11 TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
12 if (!file) file = new TFile(filename);
13 // Get AliRun object from file or create it if not on file
14 if (!gAlice)
15 {
16 gAlice = (AliRun*)file->Get("gAlice");
17 if (gAlice) printf("\nAliRun object found on file\n");
18 if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
19 }
20 Int_t nbytes = 0;
21 Int_t ipart;
22 Int_t volume;
23 Int_t np[nRings][nSectors];
24
dc8af42e 25
37c55dc0 26 for (int i=0;i<nRings;i++)
27 for(int j=0;j<nSectors;j++)
28 np[i][j]=0;
29
30 TH1F *hNReal = new TH1F("hNReal","Real number of particles",50,0,50);
31 TH1F *hNRec = new TH1F("hNRec ","Reconst. number of particles",50,0,50);
32
33 Int_t nparticles = gAlice->GetEvent(0);
34 if (nparticles <= 0) return;
35 printf("\nnparticles=%d\n",nparticles);
36
37 gAlice->TreeR()->GetEvent(0);
38 AliFMD *FMD = (AliFMD*)gAlice->GetDetector("FMD");
39 TClonesArray *Particles = gAlice->Particles();
40 TParticle *particle;
41 AliFMDhit *fmdHit;
42 AliFMDReconstParticles *fmdRP;
43 if (FMD)
44 {
45 TClonesArray *FMDhits = FMD->Hits();
46 TClonesArray *FMDrec = FMD->ReconParticles();
47 TTree *TH = gAlice->TreeH();
48 Int_t ntracks = TH->GetEntries();
49 if (ntracks<=0) return;
50
51 Int_t nPads=FMDrec->GetEntries();
52
53
54#ifdef DEBUG
55 cout<<"\n(AliFMDReconstParticles*)FMDrec->UncheckedAt(0)="<<(AliFMDReconstParticles*)FMDrec->UncheckedAt(0);
56 cout<<"\nFMDrec->UncheckedAt(0)="<<FMDrec->UncheckedAt(0);
57#endif
58
59 for (Int_t track=0; track<ntracks;track++)
60 {
61 gAlice->ResetHits();
62 nbytes += TH->GetEvent(track);//?
63 particle=(TParticle*)Particles->UncheckedAt(track);
64 // Int_t numpart=particle->GetKF();
65 //Float_t eta=particle->GetEta();
66
67 Int_t nhits=FMDhits->GetEntriesFast();
68 for (Int_t hit=0;hit<nhits;hit++)
69 {
70 fmdHit = (AliFMDhit*)FMDhits->UncheckedAt(hit);
71 volume=fmdHit->Volume();
72 if(volume==vol)
73 {
74 np[fmdHit->NumberOfRing()-1][fmdHit->NumberOfSector()-1]++;
75 }
76 }
77 }
78 //Int_t nRecPart=FMDrec->GetEntriesFast();
79 Int_t nDeterm=0; Int_t nReal=0;
80 for (Int_t pad=0;pad<nPads;pad++)
81 {
82 fmdRP = (AliFMDReconstParticles*)FMDrec->UncheckedAt(pad);
83 volume=fmdRP->GetVolume();
84 if(volume==vol)
85 {
86;
87#ifdef DEBUG
88 fmdDigit = (AliFMDdigit*)FMDdig->UncheckedAt(pad);
89 cout<<"\nfmdDigit->ADCsignal()="<<fmdDigit->ADCsignal();
90 cout<<"\nfmdDigit->NumberOfRing()="<<fmdDigit->NumberOfRing();
91 cout<<"\nfmdDigit->NumberOfSector()="<<fmdDigit->NumberOfSector();
92#endif
93 nDeterm+=fmdRP->GetNumberOfReconstParticles();
94 nReal+=np[fmdRP->GetNumberOfRing()-1][fmdRP->GetNumberOfSector()-1]; //-1=?
95 Int_t RecRing=fmdRP->GetNumberOfRing()-1;
96 Int_t RecSector=fmdRP->GetNumberOfSector()-1;
97 hNReal->Fill(np[RecRing][RecSector]);
98 hNRec->Fill(fmdRP->GetNumberOfReconstParticles());
99 }
100 }
101 }
102 cout<<"\nReal="<<nReal<<
103 " nDeterm="<<nDeterm<<
104 "\nerror="<<float(nDeterm-nReal)/float(nReal)<<endl;
105
106
107 TCanvas *c1 = new TCanvas("c1","Alice FMD ",400,10,800,800);
108 hNReal->SetFillColor(2);
109 hNReal->Draw();
110 hNRec->SetFillColor(4);
111 hNRec->Draw("same");
112
113}
114
115
116
117
118
119