]>
Commit | Line | Data |
---|---|---|
68861244 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
16 | ||
17 | //_________________________________________________________________________ | |
18 | // This is a TTask that makes TOF-Digits out of TOF-Hits. | |
19 | // A TOF Digit is essentially a TOF hit with the smearing for strip | |
20 | // time resolution and simulation of the ADC correlation signal. | |
21 | // Digits are written to TreeD in branch "TOF". | |
22 | // AliTOFDigitizer with all current parameters used for digitization | |
23 | // (time resolution and ADC parameter) is written | |
24 | // to TreeD branch "AliTOFDigitizer". | |
25 | // Both branches have the same title. If necessary one can produce | |
26 | // another set of Digits with different parameters. Two versions | |
27 | // can be distunguished using titles of the branches. | |
28 | // User case: | |
29 | // root [0] AliTOFDigitizer * s = new AliTOFDigitizer("galice.root") | |
30 | // Warning in <AliITSgeomSPD425Long::Default Creator>: Detector size may not be write. | |
31 | // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated | |
32 | // root [1] s->ExecuteTask() | |
33 | // // Makes Digits for all events stored in galice.root | |
34 | // root [2] s->SetTimeRes(100.) | |
35 | // // One can change parameters of digitization | |
36 | // root [3] s->SetDigitsBranch("Time Resolution 100. ps") | |
37 | // // and write them into the new branch | |
38 | // root [4] s->ExecuteTask("deb all tim") | |
39 | // // available parameters: | |
40 | // deb - print number of produced Digits | |
41 | // deb all - print number and list of produced Digits | |
42 | // tim - print benchmarking information | |
43 | // | |
44 | // -- Author : F. Pierella (Bologna University) pierella@bo.infn.it | |
45 | ////////////////////////////////////////////////////////////////////////////// | |
46 | ||
47 | #include "TFile.h" | |
48 | #include "TTask.h" | |
49 | #include "TTree.h" | |
50 | #include "TSystem.h" | |
51 | #include "TROOT.h" | |
52 | #include "TFolder.h" | |
53 | #include "TBenchmark.h" | |
0efc163f | 54 | #include "TObjString.h" |
68861244 | 55 | |
56 | #include <iomanip.h> | |
57 | ||
58 | #include "AliRun.h" | |
0efc163f | 59 | #include "AliTOF.h" |
68861244 | 60 | #include "AliTOFdigit.h" |
61 | #include "AliTOFhit.h" | |
62 | #include "AliTOFDigitizer.h" | |
63 | ||
64 | ||
65 | ClassImp(AliTOFDigitizer) | |
66 | ||
67 | ||
68 | //____________________________________________________________________________ | |
69 | AliTOFDigitizer::AliTOFDigitizer():TTask("AliTOFDigitizer","") | |
70 | { | |
71 | // ctor | |
72 | fTimeRes = 100; // ps | |
73 | fChrgRes = 100.; // pC | |
74 | fNevents = 0 ; | |
75 | fDigits = 0 ; | |
76 | fHits = 0 ; | |
77 | fIsInitialized = kFALSE ; | |
78 | ||
79 | } | |
80 | ||
81 | //____________________________________________________________________________ | |
82 | AliTOFDigitizer::AliTOFDigitizer(const char* headerFile, const char *digitsTitle):TTask("AliTOFDigitizer","") | |
83 | { | |
84 | // ctor | |
85 | fTimeRes = 100; // ps | |
86 | fChrgRes = 100.; // pC | |
87 | fNevents = 0 ; | |
88 | fDigitsTitle = digitsTitle ; | |
89 | fHeadersFile = headerFile ; | |
90 | fIsInitialized = kFALSE ; | |
91 | Init(); | |
92 | } | |
93 | ||
94 | //____________________________________________________________________________ | |
95 | AliTOFDigitizer::~AliTOFDigitizer() | |
96 | { | |
97 | // dtor | |
98 | if(fDigits) | |
99 | delete fDigits ; | |
100 | if(fHits) | |
101 | delete fHits ; | |
102 | } | |
103 | //____________________________________________________________________________ | |
104 | void AliTOFDigitizer::Init() | |
105 | { | |
106 | // Initialization: open root-file, allocate arrays for hits and digits, | |
107 | // attach task Digitizer to the list of TOF tasks | |
108 | // | |
109 | // Initialization can not be done in the default constructor | |
110 | ||
111 | if(!fIsInitialized){ | |
112 | ||
113 | if(fHeadersFile.IsNull()) | |
114 | fHeadersFile="galice.root" ; | |
115 | ||
116 | TFile * file = (TFile*) gROOT->GetFile(fHeadersFile.Data() ) ; | |
117 | ||
118 | //if file was not opened yet, read gAlice | |
119 | if(file == 0){ | |
120 | if(fHeadersFile.Contains("rfio")) | |
121 | file = TFile::Open(fHeadersFile,"update") ; | |
122 | else | |
123 | file = new TFile(fHeadersFile.Data(),"update") ; | |
124 | gAlice = (AliRun *) file->Get("gAlice") ; | |
125 | } | |
126 | fHits = new TClonesArray("AliTOFhit" , 405); | |
127 | fDigits = new TClonesArray("AliTOFdigit", 405); | |
128 | ||
129 | //add Task to //FPAlice/tasks/(S)Digitizer/TOF | |
6b77798a | 130 | TList * aList = gROOT->GetListOfBrowsables(); |
131 | if (aList) { | |
132 | TFolder * alice = (TFolder*)aList->FindObject("FPAlice") ; | |
133 | if (alice) { | |
134 | TTask * aliceSD = (TTask*)alice->FindObject("tasks/(S)Digitizer") ; | |
135 | if (aliceSD) { | |
136 | TTask * tofD = (TTask*)aliceSD->GetListOfTasks()->FindObject("TOF") ; | |
137 | tofD->Add(this) ; | |
138 | } | |
139 | } | |
140 | } | |
68861244 | 141 | |
142 | fIsInitialized = kTRUE ; | |
143 | } | |
144 | } | |
145 | //____________________________________________________________________________ | |
146 | void AliTOFDigitizer::Exec(Option_t *option) | |
147 | { | |
148 | Int_t tracks[3]; // track info | |
149 | Int_t vol[5]; // dummy location for digit | |
150 | Float_t digit[2]; // TOF digit variables | |
151 | ||
68861244 | 152 | // Collects all hits in the same active volume into digit |
153 | ||
154 | if(!fIsInitialized) | |
155 | Init() ; | |
156 | ||
157 | if(strstr(option,"tim")) | |
158 | gBenchmark->Start("TOFDigitizer"); | |
159 | ||
160 | fNevents = (Int_t) gAlice->TreeE()->GetEntries() ; | |
161 | ||
162 | Int_t ievent ; | |
163 | // start loop on events | |
164 | for(ievent = 0; ievent < fNevents; ievent++){ | |
165 | gAlice->GetEvent(ievent) ; | |
166 | gAlice->SetEvent(ievent) ; | |
167 | ||
168 | if(gAlice->TreeH()==0){ | |
169 | cout << "AliTOFDigitizer: There is no Hit Tree" << endl; | |
170 | return ; | |
171 | } | |
172 | ||
173 | //set address of the hits | |
174 | TBranch * branch = gAlice->TreeH()->GetBranch("TOF"); | |
175 | if (branch) | |
176 | branch->SetAddress(&fHits); | |
177 | else{ | |
178 | cout << "ERROR in AliTOFDigitizer: "<< endl ; | |
179 | cout << " no branch TOF in TreeH"<< endl ; | |
180 | cout << " do nothing " << endl ; | |
181 | return ; | |
182 | } | |
183 | ||
184 | fDigits->Clear(); | |
185 | Int_t ndigits = 0 ; | |
186 | ||
187 | //Now made Digits from hits, for TOF it is the same a part for the tof smearing | |
188 | // and some missing MC variables | |
189 | Int_t itrack ; | |
190 | for (itrack=0; itrack < gAlice->GetNtrack(); itrack++){ | |
191 | ||
192 | //=========== Get the TOF branch from Hits Tree for the Primary track itrack | |
193 | branch->GetEntry(itrack,0); | |
194 | ||
195 | Int_t i; | |
196 | for ( i = 0 ; i < fHits->GetEntries() ; i++ ) { | |
197 | ||
198 | AliTOFhit * hit = (AliTOFhit*)fHits->At(i) ; | |
199 | ||
200 | vol[0] = hit->GetSector(); | |
201 | vol[1] = hit->GetPlate(); | |
202 | vol[2] = hit->GetPadx(); | |
203 | vol[3] = hit->GetPadz(); | |
204 | vol[4] = hit->GetStrip(); | |
205 | // 95% of efficiency to be inserted here | |
206 | // edge effect to be inserted here | |
207 | // cross talk to be inserted here | |
208 | // simulation of the detector response | |
209 | Float_t idealtime = hit->GetTof(); // unit s | |
210 | idealtime *= 1.E+12; // conversion from s to ps | |
211 | // fTimeRes is given usually in ps | |
0efc163f | 212 | Float_t tdctime = gRandom->Gaus(idealtime, fTimeRes); |
68861244 | 213 | digit[0] = tdctime; |
214 | // typical Landau Distribution to be inserted here | |
215 | Float_t idealcharge = hit->GetEdep(); | |
0efc163f | 216 | Float_t adccharge = gRandom->Gaus(idealcharge, fChrgRes); |
68861244 | 217 | digit[1] = adccharge; |
218 | // to be added a check for overlaps | |
219 | new((*fDigits)[ndigits]) AliTOFdigit(tracks, vol, digit); | |
220 | ndigits++ ; | |
221 | } | |
222 | ||
223 | } // end loop over tracks | |
224 | ||
68861244 | 225 | ndigits = fDigits->GetEntriesFast() ; |
226 | printf("AliTOFDigitizer: Total number of digits %d\n",ndigits); | |
227 | ||
228 | if(gAlice->TreeD() == 0) | |
229 | gAlice->MakeTree("D") ; | |
230 | ||
231 | //check, if this branch already exits | |
232 | TBranch * digitsBranch = 0; | |
233 | TBranch * digitizerBranch = 0; | |
234 | ||
235 | TObjArray * branches = gAlice->TreeD()->GetListOfBranches() ; | |
236 | Int_t ibranch; | |
237 | Bool_t tofNotFound = kTRUE ; | |
238 | Bool_t digitizerNotFound = kTRUE ; | |
239 | ||
240 | for(ibranch = 0;ibranch <branches->GetEntries();ibranch++){ | |
241 | ||
242 | if(tofNotFound){ | |
243 | digitsBranch=(TBranch *) branches->At(ibranch) ; | |
244 | if( (strcmp("TOF",digitsBranch->GetName())==0 ) && | |
245 | (fDigitsTitle.CompareTo(digitsBranch->GetTitle()) == 0) ) | |
246 | tofNotFound = kFALSE ; | |
247 | } | |
248 | if(digitizerNotFound){ | |
249 | digitizerBranch = (TBranch *) branches->At(ibranch) ; | |
250 | if( (strcmp(digitizerBranch->GetName(),"AliTOFDigitizer") == 0)&& | |
251 | (fDigitsTitle.CompareTo(digitizerBranch->GetTitle()) == 0) ) | |
252 | digitizerNotFound = kFALSE ; | |
253 | } | |
254 | } | |
255 | ||
256 | if(!(digitizerNotFound && tofNotFound)){ | |
257 | cout << "AliTOFdigitizer error:" << endl ; | |
258 | cout << "Can not overwrite existing branches: do not write" << endl ; | |
259 | return ; | |
260 | } | |
261 | ||
262 | //Make (if necessary) branches | |
263 | char * file =0; | |
264 | if(gSystem->Getenv("CONFIG_SPLIT_FILE")){ //generating file name | |
265 | file = new char[strlen(gAlice->GetBaseFile())+20] ; | |
266 | sprintf(file,"%s/TOF.Digits.root",gAlice->GetBaseFile()) ; | |
267 | } | |
268 | ||
269 | TDirectory *cwd = gDirectory; | |
270 | ||
271 | //First list of digits | |
272 | Int_t bufferSize = 32000 ; | |
273 | digitsBranch = gAlice->TreeD()->Branch("TOF",&fDigits,bufferSize); | |
274 | digitsBranch->SetTitle(fDigitsTitle.Data()); | |
275 | if (file) { | |
276 | digitsBranch->SetFile(file); | |
277 | TIter next( digitsBranch->GetListOfBranches()); | |
278 | TBranch * subbr; | |
279 | while ((subbr=(TBranch*)next())) { | |
280 | subbr->SetFile(file); | |
281 | } | |
282 | cwd->cd(); | |
283 | } | |
284 | ||
285 | //second - Digitizer | |
286 | Int_t splitlevel = 0 ; | |
287 | AliTOFDigitizer * digtz = this ; | |
288 | digitizerBranch = gAlice->TreeD()->Branch("AliTOFDigitizer","AliTOFDigitizer", | |
289 | &digtz,bufferSize,splitlevel); | |
290 | digitizerBranch->SetTitle(fDigitsTitle.Data()); | |
291 | if (file) { | |
292 | digitizerBranch->SetFile(file); | |
293 | TIter next( digitizerBranch->GetListOfBranches()); | |
294 | TBranch * subbr ; | |
295 | while ((subbr=(TBranch*)next())) { | |
296 | subbr->SetFile(file); | |
297 | } | |
298 | cwd->cd(); | |
299 | delete file; | |
300 | } | |
301 | ||
302 | digitsBranch->Fill(); | |
303 | digitizerBranch->Fill(); | |
304 | gAlice->TreeD()->Write(0,TObject::kOverwrite) ; | |
305 | ||
306 | if(strstr(option,"deb")) | |
307 | PrintDigits(option) ; | |
308 | ||
309 | } | |
310 | ||
311 | if(strstr(option,"tim")){ | |
312 | gBenchmark->Stop("TOFDigitizer"); | |
313 | cout << "AliTOFDigitizer:" << endl ; | |
314 | cout << " took " << gBenchmark->GetCpuTime("TOFDigitizer") << " seconds for Digitizing " | |
315 | << gBenchmark->GetCpuTime("TOFDigitizer")/fNevents << " seconds per event " << endl ; | |
316 | cout << endl ; | |
317 | } | |
318 | ||
319 | ||
320 | } | |
321 | //__________________________________________________________________ | |
322 | void AliTOFDigitizer::SetDigitsBranch(const char* title ) | |
323 | { | |
324 | // Setting title to branch Digits | |
325 | if(!fDigitsTitle.IsNull()) | |
326 | cout << "AliTOFdigitizer: changing Digits file from " <<fDigitsTitle.Data() << " to " << title << endl ; | |
327 | fDigitsTitle=title ; | |
328 | } | |
329 | //__________________________________________________________________ | |
330 | void AliTOFDigitizer::Print(Option_t* option)const | |
331 | { | |
332 | // Prints parameters of Digitizer | |
333 | cout << "------------------- "<< GetName() << " -------------" << endl ; | |
334 | cout << " Writing Digits to branch with title " << fDigitsTitle.Data() << endl ; | |
335 | cout << " with digitization parameters Time resolution = " << fTimeRes << endl ; | |
336 | cout << " Adc smearing parameter = " << fChrgRes << endl ; | |
337 | cout << "---------------------------------------------------"<<endl ; | |
338 | ||
339 | } | |
340 | //__________________________________________________________________ | |
341 | Bool_t AliTOFDigitizer::operator==( AliTOFDigitizer const &digtz )const | |
342 | { | |
343 | // Equal operator. | |
344 | // Dititizers are equal if their time resolution and Adc | |
345 | // smearing parameter are equal | |
346 | ||
347 | if( (fTimeRes==digtz.fTimeRes)&&(fChrgRes==digtz.fChrgRes)) | |
348 | return kTRUE ; | |
349 | else | |
350 | return kFALSE ; | |
351 | } | |
352 | //__________________________________________________________________ | |
353 | void AliTOFDigitizer::PrintDigits(Option_t * option) | |
354 | { | |
355 | // Prints list of digits produced in the current pass of AliTOFDigitizer | |
356 | ||
357 | cout << "AliTOFDigitizer: " << endl ; | |
358 | cout << " Number of entries in Digits list " << fDigits->GetEntriesFast() << endl ; | |
359 | cout << endl ; | |
360 | ||
361 | if(strstr(option,"all")){// print all digits | |
362 | ||
363 | //loop over digits | |
364 | AliTOFdigit * digit; | |
365 | cout << "Digit # " << " Time of flight(ps) " << | |
366 | " ADC(pC)" << " Sector #" << " Plate #" << | |
367 | " Strip #" << " Pad # (x) " << " Pad # (z) " << endl; | |
368 | Int_t index ; | |
369 | for (index = 0 ; index < fDigits->GetEntries() ; index++) { | |
370 | digit = (AliTOFdigit * ) fDigits->At(index) ; | |
371 | ||
372 | // set the width of the output with the setw(int width) manipulator | |
373 | ||
374 | cout << setw(7) << index << " " | |
375 | << setw(13) << digit->GetTdc() << " " | |
376 | << setw(13) << digit->GetAdc() << " " | |
377 | << setw(6) << digit->GetSector() << " " | |
378 | << setw(6) << digit->GetPlate() << " " | |
379 | << setw(6) << digit->GetStrip() << " " | |
380 | << setw(7) << digit->GetPadx() << " " | |
381 | << setw(7) << digit->GetPadz() << " " << endl; | |
382 | ||
383 | } // end loop on digits | |
384 | ||
385 | } // close if "all" option selected | |
386 | } |