New version of the test procedure
[u/mrichter/AliRoot.git] / TRD / AliTRDanalyzeDigits.C
1 Int_t AliTRDanalyzeDigits()
2 {
3   //
4   // Analyzes the digits
5   //
6
7   Int_t rc = 0;
8
9   if (!gAlice) {
10     cout << "<AliTRDanalyzeDigits> No AliRun object found" << endl;
11     rc = 1;
12     return rc;
13   }
14   Int_t nPart = gAlice->GetEvent(0);
15
16   // Get the pointer to the TRD detector 
17   AliTRD *TRD = (AliTRD *) gAlice->GetDetector("TRD");
18   if (!TRD) {
19     cout << "<AliTRDanalyzeDigits> No TRD detector found" << endl;
20     rc = 2;
21     return rc;
22   }
23
24   // Get the digitizer object
25   TFile *file = (TFile *) gROOT->GetListOfFiles()->FindObject("TRD_test.root");
26   AliTRDdigitizer *Digitizer = (AliTRDdigitizer *) file->Get("digitizer");
27   if (!Digitizer) {
28     cout << "<AliTRDanalyzeDigits> No digitizer object found" << endl;
29     rc = 3;
30     return rc;
31   }
32
33   // Define the histograms
34   Int_t adcRange = ((Int_t) Digitizer->GetADCoutRange());
35   TH1F *hAmpAll   = new TH1F("hAmpAll"  ,"Amplitude of the digits (all)"
36                             ,adcRange+1,-0.5,((Float_t) adcRange)+0.5);
37   TH1F *hAmpEl    = new TH1F("hAmpEl"   ,"Amplitude of the digits (electrons)"
38                             ,adcRange+1,-0.5,((Float_t) adcRange)+0.5);
39   TH1F *hAmpPi    = new TH1F("hAmpPi"   ,"Amplitude of the digits (pions)"
40                             ,adcRange+1,-0.5,((Float_t) adcRange)+0.5);
41   TH1F *hAmpNoise = new TH1F("hAmpNoise","Amplitude of the digits (noise)"
42                             ,adcRange+1,-0.5,((Float_t) adcRange)+0.5);
43
44   // Get the pointer to the geometry object
45   AliTRDgeometry *TRDgeometry;
46   if (TRD) {
47     TRDgeometry = TRD->GetGeometry();
48   }
49   else {
50     cout << "<AliTRDanalyzeDigits> No TRD geometry found" << endl;
51     rc = 4;
52     return rc;
53   }
54
55   // Create the digits manager
56   AliTRDdigitsManager *DigitsManager = new AliTRDdigitsManager();
57
58   // Read the digits from the file
59   if (!(DigitsManager->ReadDigits())) {
60     cout << "<AliTRDanalyzeDigits> Cannot read the digits" << endl;
61     rc = 5;
62     return rc;
63   }
64   
65   // Define the detector matrix for one chamber
66   Int_t iSec = TRD->GetSensSector();
67   Int_t iCha = TRD->GetSensChamber();
68   Int_t iPla = 0;
69   Int_t  rowMax = TRDgeometry->GetRowMax(iPla,iCha,iSec);
70   Int_t  colMax = TRDgeometry->GetColMax(iPla);
71   Int_t timeMax = TRDgeometry->GetTimeMax();
72   cout << "<AliTRDanalyzeDigits> Geometry: rowMax = "  <<  rowMax
73                                       << " colMax = "  <<  colMax
74                                       << " timeMax = " << timeMax << endl;
75   AliTRDmatrix *TRDmatrix = new AliTRDmatrix(rowMax,colMax,timeMax,iSec,iCha,iPla);
76   // Get the detector number
77   Int_t iDet = TRDgeometry->GetDetector(iPla,iCha,iSec);
78   cout << "<AliTRDanalyzeDigits> iSec = " << iSec
79                             << " iCha = " << iCha 
80                             << " iPla = " << iPla 
81                             << " iDet = " << iDet << endl;
82
83   // Loop through the detector pixel
84   Int_t countDigits = 0;
85   for (Int_t time = 0; time < timeMax; time++) {
86     for (Int_t  col = 0;  col <  colMax;  col++) {
87       for (Int_t  row = 0;  row <  rowMax;  row++) {
88
89         AliTRDdigit *Digit = DigitsManager->GetDigit(row,col,time,iDet);
90         Int_t amp = Digit->GetAmp();
91
92         Int_t track0 = DigitsManager->GetTrack(0,row,col,time,iDet);
93         Int_t track1 = DigitsManager->GetTrack(1,row,col,time,iDet);
94         TParticle *Part = 0;
95         if (track0 > -1) {
96           Part = gAlice->Particle(track0);
97         }
98
99         if (amp > 0) {
100
101           countDigits++;
102           TRDmatrix->SetSignal(row,col,time,amp);
103
104         }
105
106         // Total spectrum
107         hAmpAll->Fill(amp);
108
109         // Noise spectrum
110         if (track0 < 0) {
111           hAmpNoise->Fill(amp);
112         }          
113
114         // Electron digit
115         if ((Part) && (Part->GetPdgCode() ==   11) && (track1 < 0)) {
116           hAmpEl->Fill(amp);
117         }
118
119         // Pion digit
120         if ((Part) && (Part->GetPdgCode() == -211) && (track1 < 0)) {
121           hAmpPi->Fill(amp);
122         }
123
124         delete Digit;
125
126       }
127     }
128   }
129
130   cout << "<AliTRDanalyzeDigits> Found " << countDigits << " digits in total" << endl;
131
132   // Display the detector matrix
133   TRDmatrix->Draw();
134   TRDmatrix->ProjRow();
135   TRDmatrix->ProjCol();
136   TRDmatrix->ProjTime();
137
138   TCanvas *cDigits = new TCanvas("cDigits","AliTRDanalyzeDigits",50,50,600,600);
139   cDigits->Divide(2,2);
140   cDigits->cd(1);
141   gPad->SetLogy();
142   hAmpAll->Draw();
143   cDigits->cd(2);
144   gPad->SetLogy();
145   gPad->SetLogx();
146   hAmpNoise->Draw();
147   cDigits->cd(3);
148   gPad->SetLogy();
149   hAmpEl->Draw();
150   cDigits->cd(4);
151   gPad->SetLogy();
152   hAmpPi->Draw();
153
154   return rc;
155
156 }