]>
Commit | Line | Data |
---|---|---|
fa148e6c | 1 | |
2 | #if !defined(__CINT__) || defined(__MAKECINT__) | |
3 | ||
4 | #include "TH1.h" | |
5 | #include "TFile.h" | |
6 | #include "TCanvas.h" | |
7 | #include "TProfile.h" | |
8 | #include "TClonesArray.h" | |
9 | ||
10 | #include "AliModule.h" | |
11 | #include "AliFRAMEv1.h" | |
12 | ||
13 | #include "AliTRDv1.h" | |
14 | #include "AliTRDhit.h" | |
15 | #include "AliTRDsim.h" | |
16 | #include "AliTRDsimple.h" | |
17 | #include "AliTRDsimpleGen.h" | |
18 | #include "AliTRDsimpleMC.h" | |
19 | #include "AliTRDgeometry.h" | |
20 | #include "AliTRDdigitizer.h" | |
21 | #include "AliTRDdigitsManager.h" | |
22 | #include "AliTRDdataArrayI.h" | |
23 | ||
24 | #endif | |
25 | ||
26 | void AliTRDrunSimple() | |
27 | { | |
28 | ||
29 | /////////////////////////////////////////////////////////////////////////////// | |
30 | // | |
31 | // Macro to run the simplified simulator | |
32 | // | |
33 | /////////////////////////////////////////////////////////////////////////////// | |
34 | ||
35 | //_____________________________________________________________________________ | |
36 | // | |
37 | // Initialization part | |
38 | //_____________________________________________________________________________ | |
39 | // | |
40 | ||
41 | // The simple simulator | |
42 | AliTRDsimple *simple = new AliTRDsimple(); | |
43 | simple->Init(); | |
44 | ||
45 | // Initialize a dummy frame so that the TRD works | |
46 | new AliFRAMEv1("FRAME","Space Frame"); | |
47 | ||
48 | // Initialize the TRD detector | |
49 | AliTRDv1 *trd = new AliTRDv1("TRD","TRD slow simulator"); | |
ee7e3ede | 50 | trd->SetHitTypeStandard(); |
fa148e6c | 51 | |
52 | // Needed for some material properties | |
53 | trd->CreateMaterials(); | |
54 | ||
55 | // Select the gas mixture (0: 97% Xe + 3% isobutane, 1: 85% Xe + 15% CO2) | |
56 | trd->SetGasMix(1); | |
57 | ||
ee7e3ede | 58 | // Define the parameter object |
59 | // If no external parameter object is defined, | |
60 | // default parameter will be used | |
61 | AliTRDparameter *parameter = new AliTRDparameter("TRDparameter" | |
62 | ,"TRD parameter class"); | |
63 | ||
64 | parameter->SetADCthreshold(0); | |
65 | parameter->SetNTimeBin(30); // The number of timebins | |
66 | parameter->SetExpandTimeBin(5,15); // The additional timebins | |
fa148e6c | 67 | |
68 | // Switch on TR production | |
69 | AliTRDsim *trdsim = trd->CreateTR(); | |
70 | trdsim->SetNFoils(100); | |
71 | trdsim->SetFoilThick(0.0013); | |
72 | trdsim->SetGapThick(0.0500); | |
73 | ||
74 | // Initialize the TRD object | |
75 | trd->Init(); | |
76 | ||
77 | // The digitizer | |
78 | AliTRDdigitizer *digitizer = new AliTRDdigitizer("TRDdigitizer","Digitizer class"); | |
ee7e3ede | 79 | |
80 | digitizer->SetParameter(parameter); | |
fa148e6c | 81 | digitizer->SetSimple(); |
fa148e6c | 82 | digitizer->InitDetector(); |
fa148e6c | 83 | |
84 | // The event generator | |
85 | AliTRDsimpleGen *generator = simple->GetGenerator(); | |
86 | generator->SetMomentum(3.0,3.0); | |
87 | generator->SetPdg(11); // Electrons | |
fa148e6c | 88 | |
89 | //_____________________________________________________________________________ | |
90 | // | |
91 | // Histograms | |
92 | //_____________________________________________________________________________ | |
93 | // | |
94 | ||
ee7e3ede | 95 | Int_t timeMax = parameter->GetTimeTotal(); |
96 | Int_t adcRange = ((Int_t) parameter->GetADCoutRange()); | |
fa148e6c | 97 | |
98 | TH1F *hQ = new TH1F("hQ" ,"Charge per hit (all)" ,100,0.0,1000.0); | |
99 | TH1F *hQdedx = new TH1F("hQdedx","Charge per hit (dedx)",100,0.0,1000.0); | |
100 | TH1F *hQtr = new TH1F("hQtr ","Charge per hit (tr)" ,100,0.0,1000.0); | |
101 | ||
102 | TProfile *hQX = new TProfile("hQX" ,"Charge per hit vs X (all)" ,35,0.0,3.5); | |
103 | TProfile *hQXdedx = new TProfile("hQXdedx","Charge per hit vs X (dedx)",35,0.0,3.5); | |
104 | TProfile *hQXtr = new TProfile("hQXtr ","Charge per hit vs X (tr)" ,35,0.0,3.5); | |
105 | ||
106 | TH1F *hNstep = new TH1F("hNstep","Number of steps / cm" , 151,-0.5, 150.5); | |
107 | TH1F *hNel = new TH1F("hNel" ,"Number of electrons / cm",1001,-0.5,1000.5); | |
108 | ||
109 | TH1F *hAmp = new TH1F("hAmp","Amplitude of the digits" | |
110 | ,adcRange+1,-0.5,((Float_t) adcRange)+0.5); | |
111 | TProfile *hAmpTime = new TProfile("hAmpTime","Amplitude vs timebin" | |
112 | ,timeMax,-0.5,((Float_t) timeMax)-0.5); | |
113 | ||
114 | //_____________________________________________________________________________ | |
115 | // | |
116 | // Event loop | |
117 | //_____________________________________________________________________________ | |
118 | // | |
119 | ||
120 | // Number of events (tracks) | |
121 | Int_t nEvent = 10000; | |
122 | ||
ee7e3ede | 123 | Float_t x0 = parameter->GetTime0(0) - AliTRDgeometry::DrThick(); |
fa148e6c | 124 | |
125 | TClonesArray *hitsArray = trd->Hits(); | |
126 | ||
127 | for (Int_t iEvent = 0; iEvent < nEvent; iEvent++) { | |
128 | ||
129 | if (!(iEvent % 100) && (iEvent)) { | |
130 | printf("Event no. %d\n",iEvent); | |
131 | } | |
132 | ||
133 | // Generate the hits for one track | |
134 | simple->ProcessEvent(iEvent); | |
135 | ||
136 | // Analyze the hits | |
137 | Float_t nElectrons = 0.0; | |
138 | for (Int_t iHit = 0; iHit < hitsArray->GetEntries(); iHit++) { | |
139 | ||
140 | AliTRDhit *hit = (AliTRDhit *) hitsArray->UncheckedAt(iHit); | |
141 | Int_t charge = TMath::Abs(hit->GetCharge()); | |
142 | Float_t x = hit->X() - x0; | |
143 | hQ->Fill(charge); | |
144 | hQX->Fill(x,charge); | |
145 | if (hit->FromDrift() || | |
146 | hit->FromAmplification()) { | |
147 | hQdedx->Fill(charge); | |
148 | hQXdedx->Fill(x,charge); | |
149 | nElectrons += charge; | |
150 | } | |
151 | if (hit->FromTRphoton()) { | |
152 | hQtr->Fill(charge); | |
153 | hQXtr->Fill(x,charge); | |
154 | } | |
155 | ||
156 | } | |
157 | ||
158 | hNstep->Fill(((AliTRDsimpleMC *) gMC)->GetNStep() / 3.5); | |
159 | hNel->Fill(nElectrons / 3.5); | |
160 | ||
161 | // Digitize the hits | |
162 | digitizer->MakeDigits(); | |
163 | ||
164 | // Analyze the digits | |
165 | Int_t row = 8; | |
166 | Int_t col = 48; | |
167 | AliTRDdigitsManager *digitsManager = digitizer->Digits(); | |
168 | AliTRDdataArrayI *digitsArray = digitsManager->GetDigits(12); | |
169 | for (Int_t time = 0; time < timeMax; time++) { | |
170 | ||
171 | Int_t amp = digitsArray->GetDataUnchecked(row,col,time); | |
172 | if (amp > 0) { | |
173 | hAmp->Fill((Double_t) amp); | |
174 | hAmpTime->Fill((Double_t) time, (Double_t) amp); | |
175 | } | |
176 | ||
177 | } | |
178 | ||
179 | } | |
180 | ||
181 | //_____________________________________________________________________________ | |
182 | // | |
183 | // Plot the spectra | |
184 | //_____________________________________________________________________________ | |
185 | // | |
186 | ||
187 | TCanvas *cHit = new TCanvas("cHit","Hits",50,50,800,600); | |
188 | cHit->Divide(3,2); | |
189 | ||
190 | cHit->cd(1); | |
191 | gPad->SetLogy(); | |
192 | hQ->SetLineColor(4); | |
193 | hQ->SetXTitle("Q (a.u.)"); | |
194 | hQ->SetYTitle("Entries"); | |
195 | hQ->Draw(); | |
196 | hQdedx->SetLineColor(3); | |
197 | hQdedx->Draw("SAME"); | |
198 | hQtr->SetLineColor(2); | |
199 | hQtr->Draw("SAME"); | |
200 | ||
201 | cHit->cd(2); | |
202 | gPad->SetLogy(); | |
203 | hQX->SetLineColor(4); | |
204 | hQX->SetXTitle("x (cm)"); | |
205 | hQX->SetYTitle("<Q> (a.u.)"); | |
206 | hQX->Draw("HIST"); | |
207 | hQXdedx->SetLineColor(3); | |
208 | hQXdedx->Draw("SAMEHIST"); | |
209 | ||
210 | cHit->cd(3); | |
211 | gPad->SetLogy(); | |
212 | hQXtr->SetLineColor(2); | |
213 | hQXtr->SetXTitle("x (cm)"); | |
214 | hQXtr->SetYTitle("<Q> (a.u.)"); | |
215 | hQXtr->Draw("HIST"); | |
216 | ||
217 | cHit->cd(4); | |
218 | hNstep->SetLineColor(4); | |
219 | hNstep->SetXTitle("N_{step} / cm"); | |
220 | hNstep->SetYTitle("Entries"); | |
221 | hNstep->Draw(); | |
222 | ||
223 | cHit->cd(5); | |
224 | hNel->SetLineColor(4); | |
225 | hNel->SetXTitle("N_{electron} / cm"); | |
226 | hNel->SetYTitle("Entries"); | |
227 | hNel->Draw(); | |
228 | ||
229 | TCanvas *cDigit = new TCanvas("cDigit","Digits",50,50,600,400); | |
230 | cDigit->Divide(2,1); | |
231 | ||
232 | cDigit->cd(1); | |
233 | gPad->SetLogy(); | |
234 | hAmp->SetLineColor(2); | |
235 | hAmp->SetXTitle("ADC channels"); | |
236 | hAmp->SetYTitle("Entries"); | |
237 | hAmp->Draw(); | |
238 | ||
239 | cDigit->cd(2); | |
240 | hAmpTime->SetLineColor(2); | |
241 | hAmpTime->SetXTitle("Timebin"); | |
242 | hAmpTime->SetYTitle("<ADC channels>"); | |
243 | hAmpTime->Draw("HIST"); | |
244 | ||
245 | //_____________________________________________________________________________ | |
246 | // | |
247 | // Save the histograms | |
248 | //_____________________________________________________________________________ | |
249 | // | |
250 | ||
251 | TFile *fileOut = new TFile("simple.root","RECREATE"); | |
252 | ||
253 | hQ->Write(); | |
254 | hQdedx->Write(); | |
255 | hQtr->Write(); | |
256 | ||
257 | hQX->Write(); | |
258 | hQXdedx->Write(); | |
259 | hQXtr->Write(); | |
260 | ||
261 | hNstep->Write(); | |
262 | hNel->Write(); | |
263 | ||
264 | hAmp->Write(); | |
265 | hAmpTime->Write(); | |
266 | ||
267 | fileOut->Close(); | |
268 | ||
269 | } |