Update of TRD code
[u/mrichter/AliRoot.git] / TRD / AliTRDanalyzeDigits.C
CommitLineData
793ff80c 1Int_t AliTRDanalyzeDigits()
2{
3 //
4 // Analyzes the digits
5 //
6
7 Int_t rc = 0;
8
e2c86a4a 9 if (!gAlice) {
10 cout << "<AliTRDanalyzeDigits> No AliRun object found" << endl;
11 rc = 1;
12 return rc;
13 }
4a0fe73c 14 Int_t nPart = gAlice->GetEvent(0);
e2c86a4a 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
4a0fe73c 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
e2c86a4a 33 // Define the histograms
4a0fe73c 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);
e2c86a4a 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;
4a0fe73c 51 rc = 4;
e2c86a4a 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;
4a0fe73c 61 rc = 5;
e2c86a4a 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
4a0fe73c 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
e2c86a4a 99 if (amp > 0) {
4a0fe73c 100
e2c86a4a 101 countDigits++;
e2c86a4a 102 TRDmatrix->SetSignal(row,col,time,amp);
4a0fe73c 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);
e2c86a4a 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);
4a0fe73c 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);
e2c86a4a 151 gPad->SetLogy();
4a0fe73c 152 hAmpPi->Draw();
e2c86a4a 153
793ff80c 154 return rc;
155
156}