]>
Commit | Line | Data |
---|---|---|
793ff80c | 1 | Int_t AliTRDanalyzeHits() |
2 | { | |
3 | // | |
4 | // Analyzes the hits and fills QA-histograms | |
5 | // | |
6 | ||
7 | Int_t rc = 0; | |
8 | ||
9 | if (!gAlice) { | |
10 | cout << "<AliTRDanalyzeHits> No AliRun object found" << endl; | |
11 | rc = 1; | |
12 | return rc; | |
13 | } | |
14 | gAlice->GetEvent(0); | |
15 | ||
16 | // Get the pointer to the TRD detector | |
f9428ca8 | 17 | AliTRD *TRD = (AliTRD *) gAlice->GetDetector("TRD"); |
793ff80c | 18 | if (!TRD) { |
19 | cout << "<AliTRDanalyzeHits> No TRD detector found" << endl; | |
20 | rc = 2; | |
21 | return rc; | |
22 | } | |
23 | ||
f9428ca8 | 24 | // Get the pointer to the geometry object |
25 | AliTRDgeometryFull *TRDgeometry; | |
26 | if (TRD) { | |
27 | TRDgeometry = (AliTRDgeometryFull *) TRD->GetGeometry(); | |
28 | } | |
29 | else { | |
30 | cout << "<AliTRDanalyzeHits> No TRD geometry found" << endl; | |
31 | rc = 3; | |
32 | return rc; | |
33 | } | |
34 | ||
793ff80c | 35 | // Define the histograms |
f9428ca8 | 36 | TH1F *hQdedx = new TH1F("hQdedx","Charge dedx-hits",100,0.0,1000.0); |
37 | TH1F *hQtr = new TH1F("hQtr" ,"Charge TR-hits" ,100,0.0,1000.0); | |
38 | ||
39 | Float_t rmin = TRDgeometry->Rmin(); | |
40 | Float_t rmax = TRDgeometry->Rmax(); | |
41 | Float_t length = TRDgeometry->GetChamberLengthI(0); | |
42 | Float_t width = TRDgeometry->GetChamberWidth(0); | |
43 | Int_t ncol = TRDgeometry->GetColMax(0); | |
44 | Int_t nrow = TRDgeometry->GetRowMax(0,2,13); | |
45 | Int_t ntime = ((Int_t) (rmax - rmin) / TRDgeometry->GetTimeBinSize()); | |
46 | ||
47 | TH2F *hZY = new TH2F("hZY" ,"Y vs Z (chamber 0)", nrow,-length/2.,length/2. | |
48 | ,ntime, rmin, rmax); | |
49 | TH2F *hXZ = new TH2F("hXZ" ,"Z vs X (plane 0)" , ncol, -width/2., width/2. | |
50 | , nrow,-length/2.,length/2.); | |
51 | ||
52 | TH1F *hNprimP = new TH1F("hNprimP","Number of primary electrons per cm (MIP pion)" | |
53 | ,50,0.0,100.0); | |
54 | TH1F *hNprimE = new TH1F("hNprimE","Number of primary electrons per cm (3GeV electron)" | |
55 | ,50,0.0,100.0); | |
56 | TH1F *hNtotP = new TH1F("hNtotP" ,"Number of electrons per cm (MIP pion)" | |
57 | ,50,0.0,1000.0); | |
58 | TH1F *hNtotE = new TH1F("hNtotE" ,"Number of electrons per cm (3GeV electron)" | |
59 | ,50,0.0,1000.0); | |
793ff80c | 60 | |
61 | // Get the pointer hit tree | |
62 | TTree *HitTree = gAlice->TreeH(); | |
63 | if (!HitTree) { | |
64 | cout << "<AliTRDanalyzeHits> No hit tree found" << endl; | |
f9428ca8 | 65 | rc = 4; |
793ff80c | 66 | return rc; |
67 | } | |
68 | ||
69 | Int_t countHits = 0; | |
70 | Int_t nBytes = 0; | |
71 | ||
72 | // Get the number of entries in the hit tree | |
73 | // (Number of primary particles creating a hit somewhere) | |
74 | Int_t nTrack = (Int_t) HitTree->GetEntries(); | |
75 | cout << "<AliTRDanalyzeHits> Found " << nTrack | |
76 | << " primary particles with hits" << endl; | |
77 | ||
78 | // Loop through all entries in the tree | |
79 | for (Int_t iTrack = 0; iTrack < nTrack; iTrack++) { | |
80 | ||
81 | gAlice->ResetHits(); | |
82 | nBytes += HitTree->GetEvent(iTrack); | |
83 | ||
f9428ca8 | 84 | Int_t nPrimE = 0; |
85 | Int_t nPrimP = 0; | |
86 | Int_t nTotE = 0; | |
87 | Int_t nTotP = 0; | |
88 | ||
793ff80c | 89 | // Get the number of hits in the TRD created by this particle |
90 | Int_t nHit = TRD->Hits()->GetEntriesFast(); | |
793ff80c | 91 | |
92 | // Loop through the TRD hits | |
93 | for (Int_t iHit = 0; iHit < nHit; iHit++) { | |
94 | ||
95 | countHits++; | |
96 | ||
97 | AliTRDhit *hit = (AliTRDhit *) TRD->Hits()->UncheckedAt(iHit); | |
98 | ||
f9428ca8 | 99 | Float_t x = hit->X(); |
100 | Float_t y = hit->Y(); | |
101 | Float_t z = hit->Z(); | |
102 | Float_t q = hit->GetCharge(); | |
103 | Int_t track = hit->Track(); | |
104 | Int_t det = hit->GetDetector(); | |
105 | Int_t plane = TRDgeometry->GetPlane(det); | |
793ff80c | 106 | |
f9428ca8 | 107 | if (hit->FromDrift()) { |
793ff80c | 108 | hQdedx->Fill(q); |
f9428ca8 | 109 | hZY->Fill(z,y); |
110 | if (plane == 0) { | |
111 | hXZ->Fill(x,z); | |
112 | } | |
4a0fe73c | 113 | } |
f9428ca8 | 114 | else if (hit->FromTRphoton()) { |
793ff80c | 115 | hQtr->Fill(TMath::Abs(q)); |
4a0fe73c | 116 | } |
793ff80c | 117 | |
f9428ca8 | 118 | TParticle *part = gAlice->Particle(track); |
119 | ||
120 | if ((plane == 0) && (hit->FromDrift())) { | |
121 | ||
122 | // 3 GeV electrons | |
123 | if ((part->GetPdgCode() == 11) && | |
124 | (part->P() > 2.9)) { | |
125 | nPrimE++; | |
126 | nTotE += ((Int_t) q); | |
127 | } | |
128 | ||
129 | // MIP pions | |
130 | if ((part->GetPdgCode() == -211) && | |
131 | (part->P() > 0.5)) { | |
132 | nPrimP++; | |
133 | nTotP += ((Int_t) q); | |
134 | } | |
135 | ||
136 | } | |
793ff80c | 137 | |
138 | } | |
139 | ||
f9428ca8 | 140 | if (nPrimE > 0) hNprimE->Fill(((Double_t) nPrimE)/3.); |
141 | if (nPrimP > 0) hNprimP->Fill(((Double_t) nPrimP)/3.); | |
142 | if (nTotE > 0) hNtotE->Fill(((Double_t) nTotE)/3.); | |
143 | if (nTotP > 0) hNtotP->Fill(((Double_t) nTotP)/3.); | |
144 | ||
793ff80c | 145 | } |
146 | ||
147 | cout << "<AliTRDanalyzeHits> Found " << countHits << " hits in total" << endl; | |
148 | ||
149 | TCanvas *cHits = new TCanvas("cHits","AliTRDanalyzeHits",50,50,600,600); | |
150 | cHits->Divide(2,2); | |
151 | cHits->cd(1); | |
152 | hXZ->Draw("COL"); | |
153 | cHits->cd(2); | |
154 | hZY->Draw("COL"); | |
155 | cHits->cd(3); | |
156 | gPad->SetLogy(); | |
157 | hQdedx->Draw(); | |
158 | cHits->cd(4); | |
159 | gPad->SetLogy(); | |
160 | hQtr->Draw(); | |
161 | ||
f9428ca8 | 162 | TCanvas *cNel = new TCanvas("cNel","Ionization",50,50,600,600); |
163 | cNel->Divide(2,2); | |
164 | cNel->cd(1); | |
165 | hNprimE->SetStats(); | |
166 | hNprimE->Draw(); | |
167 | cNel->cd(2); | |
168 | hNprimP->SetStats(); | |
169 | hNprimP->Draw(); | |
170 | cNel->cd(3); | |
171 | hNtotE->SetStats(); | |
172 | hNtotE->Draw(); | |
173 | cNel->cd(4); | |
174 | hNtotP->SetStats(); | |
175 | hNtotP->Draw(); | |
176 | ||
177 | TFile *fout = new TFile("TRD_ionization.root","RECREATE"); | |
178 | hNprimE->Write(); | |
179 | hNprimP->Write(); | |
180 | hNtotE->Write(); | |
181 | hNtotP->Write(); | |
182 | fout->Close(); | |
183 | ||
793ff80c | 184 | return rc; |
185 | ||
186 | } |