b52b2b3e477a60fd2b81e1da17e450d524ecfd2b
[u/mrichter/AliRoot.git] / STEER / AliRunDigitizer.cxx
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 $Log$
18 Revision 1.12  2001/12/10 16:40:52  jchudoba
19 Import gAlice from the signal file before InitGlobal() to allow detectors to use it during initialization
20
21 Revision 1.11  2001/12/03 07:10:13  jchudoba
22 Default ctor cannot create new objects, create dummy default ctor which leaves object in not well defined state - to be used only by root for I/O
23
24 Revision 1.10  2001/11/15 11:07:25  jchudoba
25 Set to zero new pointers to TPC and TRD special trees in the default ctor. Add const to all Get functions. Remove unused constant, rename constant according coding rules.
26
27 Revision 1.9  2001/11/15 09:00:11  jchudoba
28 Add special treatment for TPC and TRD, they use different trees than other detectors
29
30 Revision 1.8  2001/10/21 18:38:43  hristov
31 Several pointers were set to zero in the default constructors to avoid memory management problems
32
33 Revision 1.7  2001/10/04 15:56:07  jchudoba
34 TTask inheritance
35
36 Revision 1.4  2001/09/19 06:23:50  jchudoba
37 Move some tasks to AliStream and AliMergeCombi classes
38
39 Revision 1.3  2001/07/30 14:04:18  jchudoba
40 correct bug in the initialization
41
42 Revision 1.2  2001/07/28 10:44:32  hristov
43 Loop variable declared once; typos corrected
44
45 Revision 1.1  2001/07/27 12:59:00  jchudoba
46 Manager class for merging/digitization
47
48 */
49
50 ////////////////////////////////////////////////////////////////////////
51 //
52 // AliRunDigitizer.cxx
53 //
54 // Manager object for merging/digitization
55 //
56 // Instance of this class manages the digitization and/or merging of
57 // Sdigits into Digits. 
58 //
59 // Only one instance of this class is created in the macro:
60 //   AliRunDigitizer * manager = 
61 //      new AliRunDigitizer(nInputStreams,SPERB);
62 // where nInputStreams is number of input streams and SPERB is
63 // signals per background variable, which determines how combinations
64 // of signal and background events are generated.
65 // Then instances of specific detector digitizers are created:
66 //   AliMUONDigitizer *dMUON  = new AliMUONDigitizer(manager)
67 // and the I/O configured (you have to specify input files 
68 // and an output file). The manager connects appropriate trees from 
69 // the input files according a combination returned by AliMergeCombi 
70 // class. It creates TreeD in the output and runs once per 
71 // event Digitize method of all existing AliDetDigitizers 
72 // (without any option). AliDetDigitizers ask manager
73 // for a TTree with input (manager->GetInputTreeS(Int_t i),
74 // merge all inputs, digitize it, and save it in the TreeD 
75 // obtained by manager->GetTreeD(). Output events are stored with 
76 // numbers from 0, this default can be changed by 
77 // manager->SetFirstOutputEventNr(Int_t) method. The particle numbers
78 // in the output are shifted by MASK, which is taken from manager.
79 //
80 // The default output is to the signal file (stream 0). This can be 
81 // changed with the SetOutputFile(TString fn)  method.
82 //
83 // Single input file is permitted. Maximum kMaxStreamsToMerge can be merged.
84 // Input from the memory (on-the-fly merging) is not yet 
85 // supported, as well as access to the input data by invoking methods
86 // on the output data.
87 //
88 // Access to the some data is via gAlice for now (supposing the 
89 // same geometry in all input files), gAlice is taken from the first 
90 // input file on the first stream.
91 //
92 // Example with MUON digitizer, no merging, just digitization
93 //
94 //  AliRunDigitizer * manager = new AliRunDigitizer(1,1);
95 //  manager->SetInputStream(0,"galice.root");
96 //  AliMUONDigitizer *dMUON  = new AliMUONDigitizer(manager);
97 //  manager->Exec("");
98 //
99 // Example with MUON digitizer, merge all events from 
100 //   galice.root (signal) file with events from bgr.root 
101 //   (background) file. Number of merged events is
102 //   min(number of events in galice.root, number of events in bgr.root)
103 //
104 //  AliRunDigitizer * manager = new AliRunDigitizer(2,1);
105 //  manager->SetInputStream(0,"galice.root");
106 //  manager->SetInputStream(1,"bgr.root");
107 //  AliMUONDigitizer *dMUON  = new AliMUONDigitizer(manager);
108 //  manager->Exec("");
109 //
110 // Example with MUON digitizer, save digits in a new file digits.root,
111 //   process only 1 event
112 //
113 //  AliRunDigitizer * manager = new AliRunDigitizer(2,1);
114 //  manager->SetInputStream(0,"galice.root");
115 //  manager->SetInputStream(1,"bgr.root");
116 //  manager->SetOutputFile("digits.root");
117 //  AliMUONDigitizer *dMUON  = new AliMUONDigitizer(manager);
118 //  manager->SetNrOfEventsToWrite(1);
119 //  manager->Exec("");
120 //
121 //////////////////////////////////////////////////////////////////////// 
122
123 // system includes
124
125 #include <iostream.h>
126
127 // ROOT includes
128
129 #include "TFile.h"
130 #include "TTree.h"
131 #include "TList.h"
132
133 // AliROOT includes
134
135 #include "AliRunDigitizer.h"
136 #include "AliDigitizer.h"
137 #include "AliRun.h"
138 #include "AliHeader.h"
139 #include "TParticle.h"
140 #include "AliStream.h"
141 #include "AliMergeCombi.h"
142
143 ClassImp(AliRunDigitizer)
144
145 ////////////////////////////////////////////////////////////////////////
146 AliRunDigitizer::AliRunDigitizer()
147 {
148 // root requires default ctor, where no new objects can be created
149 // do not use this ctor, it is supplied only for root needs
150   
151 // just set all pointers - data members to 0
152   fOutput = 0;
153   fTreeD = 0;
154   fTreeDTPC = 0;
155   fTreeDTRD = 0;
156   fInputStreams = 0;
157   for (Int_t i=0;i<kMaxStreamsToMerge;i++) {
158     fArrayTreeS[i]=fArrayTreeH[i]=fArrayTreeTPCS[i]=fArrayTreeTRDS[i]=NULL;
159     fInputFiles[i]=0;
160   }
161   fCombi = 0;
162
163 }
164
165 ////////////////////////////////////////////////////////////////////////
166 AliRunDigitizer::AliRunDigitizer(Int_t nInputStreams, Int_t sperb) : TTask("AliRunDigitizer","The manager for Merging")
167 {
168 // ctor which should be used to create a manager for merging/digitization
169   if (nInputStreams == 0) {
170     Error("AliRunDigitizer","Specify nr of input streams");
171     return;
172   }
173   Int_t i;
174   fNinputs = nInputStreams;
175   fOutputFileName = "";
176   fOutputDirName = ".";
177   fCombination.Set(kMaxStreamsToMerge);
178   for (i=0;i<kMaxStreamsToMerge;i++) {
179     fArrayTreeS[i]=fArrayTreeH[i]=fArrayTreeTPCS[i]=fArrayTreeTRDS[i]=NULL;
180     fCombination[i]=-1;
181   }
182   fkMASKSTEP = 10000000;
183   fkMASK[0] = 0;
184   for (i=1;i<kMaxStreamsToMerge;i++) {
185     fkMASK[i] = fkMASK[i-1] + fkMASKSTEP;
186   }
187   fInputStreams = new TClonesArray("AliStream",nInputStreams);
188   TClonesArray &lInputStreams = *fInputStreams;
189 // the first Input is open RW to be output as well
190   new(lInputStreams[0]) AliStream("UPDATE");
191   for (i=1;i<nInputStreams;i++) {
192     new(lInputStreams[i]) AliStream("READ");
193   }
194   fOutput = 0;
195   fEvent = 0;
196   fNrOfEventsToWrite = -1;
197   fNrOfEventsWritten = 0;
198   fCopyTreesFromInput = -1;
199   fCombi = new AliMergeCombi(nInputStreams,sperb);
200   fDebug = 0;
201   fTreeD = 0;
202   fTreeDTPC = 0;
203   fTreeDTRD = 0;
204   fTreeDTPCBaseName = "TreeD_75x40_100x60_";
205   fTreeTPCSBaseName = "TreeS_75x40_100x60_";
206
207   for (i=0; i<kMaxStreamsToMerge; i++) fInputFiles[i]=0;
208 }
209
210 ////////////////////////////////////////////////////////////////////////
211
212 AliRunDigitizer::~AliRunDigitizer() {
213 // dtor
214
215   if (fInputStreams) {
216     delete fInputStreams;
217     fInputStreams = 0;
218   }
219   if (fCombi) {
220     delete fCombi;
221     fCombi = 0;
222   }
223
224 }
225 ////////////////////////////////////////////////////////////////////////
226 void AliRunDigitizer::AddDigitizer(AliDigitizer *digitizer)
227 {
228 // add digitizer to the list of active digitizers
229   this->Add(digitizer);
230 }
231 ////////////////////////////////////////////////////////////////////////
232 void AliRunDigitizer::SetInputStream(Int_t i, char *inputFile)
233 {
234   if (i > fInputStreams->GetLast()) {
235     Error("SetInputStream","Input stream number too high");
236     return;
237   }
238   static_cast<AliStream*>(fInputStreams->At(i))->AddFile(inputFile);
239 }
240
241 ////////////////////////////////////////////////////////////////////////
242 void AliRunDigitizer::Digitize(Option_t* option)
243 {
244 // get a new combination of inputs, connect input trees and loop 
245 // over all digitizers
246
247 // take gAlice from the first input file. It is needed to access
248 //  geometry data
249 // If gAlice is already in memory, use it
250   if (!gAlice) {
251     if (!static_cast<AliStream*>(fInputStreams->At(0))->ImportgAlice()) {
252       cerr<<"gAlice object not found in the first file of "
253           <<"the 1st stream"<<endl;
254       return;
255     }
256   }
257   if (!InitGlobal()) {
258     cerr<<"False from InitGlobal"<<endl;
259     return;
260   }
261   Int_t eventsCreated = 0;
262 // loop until there is anything on the input in case fNrOfEventsToWrite < 0
263   while ((eventsCreated++ < fNrOfEventsToWrite) || (fNrOfEventsToWrite < 0)) {
264     if (!ConnectInputTrees()) break;
265     InitEvent();
266 // loop over all registered digitizers and let them do the work
267     ExecuteTasks(option);
268     CleanTasks();
269     FinishEvent();
270   }
271   FinishGlobal();
272 }
273
274 ////////////////////////////////////////////////////////////////////////
275 Bool_t AliRunDigitizer::ConnectInputTrees()
276 {
277 // fill arrays fArrayTreeS, fArrayTreeH and fArrayTreeTPCS with 
278 // pointers to the correct events according fCombination values
279 // null pointers can be in the output, AliDigitizer has to check it
280
281   TTree *tree;
282   char treeName[50];
283   Int_t serialNr;
284   Int_t eventNr[kMaxStreamsToMerge], delta[kMaxStreamsToMerge];
285   fCombi->Combination(eventNr, delta);
286   for (Int_t i=0;i<fNinputs;i++) {
287     if (delta[i] == 1) {
288       AliStream *iStream = static_cast<AliStream*>(fInputStreams->At(i));
289       if (!iStream->NextEventInStream(serialNr)) return kFALSE;
290       fInputFiles[i]=iStream->CurrentFile();
291       sprintf(treeName,"TreeS%d",serialNr);
292       tree = static_cast<TTree*>(iStream->CurrentFile()->Get(treeName));
293       if (fArrayTreeS[i]) {
294         delete fArrayTreeS[i];
295         fArrayTreeS[i] = 0;
296       }
297       fArrayTreeS[i] = tree;
298       sprintf(treeName,"TreeH%d",serialNr);
299       tree = static_cast<TTree*>(iStream->CurrentFile()->Get(treeName));
300       if (fArrayTreeH[i]) {
301         delete fArrayTreeH[i];
302         fArrayTreeH[i] = 0;
303       }
304       fArrayTreeH[i] = tree;
305       sprintf(treeName,"%s%d",fTreeTPCSBaseName,serialNr);
306       tree = static_cast<TTree*>(iStream->CurrentFile()->Get(treeName));
307       if (fArrayTreeTPCS[i]) {
308         delete fArrayTreeTPCS[i];
309         fArrayTreeTPCS[i] = 0;
310       }
311       fArrayTreeTPCS[i] = tree;
312       sprintf(treeName,"TreeS%d_TRD",serialNr);
313       tree = static_cast<TTree*>(iStream->CurrentFile()->Get(treeName));
314       if (fArrayTreeTRDS[i]) {
315         delete fArrayTreeTRDS[i];
316         fArrayTreeTRDS[i] = 0;
317       }
318       fArrayTreeTRDS[i] = tree;
319     } else if (delta[i] != 0) {
320       Error("ConnectInputTrees","Only delta 0 or 1 is implemented");
321       return kFALSE;
322     }
323   }
324   return kTRUE;
325 }
326
327 ////////////////////////////////////////////////////////////////////////
328 Bool_t AliRunDigitizer::InitGlobal()
329 {
330 // called once before Digitize() is called, initialize digitizers and output
331
332   TList* subTasks = this->GetListOfTasks();
333   if (subTasks) {
334     subTasks->ForEach(AliDigitizer,Init)();
335   }  
336   return kTRUE;
337 }
338
339 ////////////////////////////////////////////////////////////////////////
340
341 void AliRunDigitizer::SetOutputFile(TString fn)
342 // the output will be to separate file, not to the signal file
343 {
344   fOutputFileName = fn;
345   (static_cast<AliStream*>(fInputStreams->At(0)))->ChangeMode("READ");
346   InitOutputGlobal();
347 }
348
349 ////////////////////////////////////////////////////////////////////////
350 Bool_t AliRunDigitizer::InitOutputGlobal()
351 {
352 // Creates the output file, called by InitEvent()
353
354   TString fn;
355   fn = fOutputDirName + '/' + fOutputFileName;
356   fOutput = new TFile(fn,"recreate");
357   if (GetDebug()>2) {
358     cerr<<"AliRunDigitizer::InitOutputGlobal(): file "<<fn.Data()<<" was opened"<<endl;
359   }
360   if (fOutput) return kTRUE;
361   Error("InitOutputGlobal","Could not create output file.");
362   return kFALSE;
363 }
364
365
366 ////////////////////////////////////////////////////////////////////////
367 void AliRunDigitizer::InitEvent()
368 {
369 // Creates TreeDxx in the output file, called from Digitize() once for 
370 //  each event. xx = fEvent
371
372   if (GetDebug()>2) 
373     cerr<<"AliRunDigitizer::InitEvent: fEvent = "<<fEvent<<endl;
374
375 // if fOutputFileName was not given, write output to signal file
376   if (fOutputFileName == "") {
377     fOutput = (static_cast<AliStream*>(fInputStreams->At(0)))->CurrentFile();
378   }
379   fOutput->cd();
380   char treeName[30];
381   sprintf(treeName,"TreeD%d",fEvent);
382   fTreeD = static_cast<TTree*>(fOutput->Get(treeName));
383   if (!fTreeD) {
384     fTreeD = new TTree(treeName,"Digits");
385     fTreeD->Write(0,TObject::kOverwrite);
386   }
387
388 // special tree for TPC
389   sprintf(treeName,"%s%d",fTreeDTPCBaseName,fEvent);
390   fTreeDTPC = static_cast<TTree*>(fOutput->Get(treeName));
391   if (!fTreeDTPC) {
392     fTreeDTPC = new TTree(treeName,"TPC_Digits");
393     fTreeDTPC->Write(0,TObject::kOverwrite);
394   }
395
396 // special tree for TRD
397   sprintf(treeName,"TreeD%d_TRD",fEvent);
398   fTreeDTRD = static_cast<TTree*>(fOutput->Get(treeName));
399   if (!fTreeDTRD) {
400     fTreeDTRD = new TTree(treeName,"TRD_Digits");
401     fTreeDTRD->Write(0,TObject::kOverwrite);
402   }
403
404 }
405
406 ////////////////////////////////////////////////////////////////////////
407 void AliRunDigitizer::FinishEvent()
408 {
409 // called at the end of loop over digitizers
410
411   Int_t i;
412   fOutput->cd();
413   if (fCopyTreesFromInput > -1) {
414     char treeName[20];
415     i = fCopyTreesFromInput; 
416     sprintf(treeName,"TreeK%d",fCombination[i]);
417     fInputFiles[i]->Get(treeName)->Clone()->Write();
418     sprintf(treeName,"TreeH%d",fCombination[i]);
419     fInputFiles[i]->Get(treeName)->Clone()->Write();
420   }
421   fEvent++;
422   fNrOfEventsWritten++;
423   if (fTreeD) {
424     delete fTreeD;
425     fTreeD = 0;
426   }
427   if (fTreeDTPC) {
428     delete fTreeDTPC;
429     fTreeDTPC = 0;
430   }
431   if (fTreeDTRD) {
432     delete fTreeDTRD;
433     fTreeDTRD = 0;
434   }
435 }
436 ////////////////////////////////////////////////////////////////////////
437 void AliRunDigitizer::FinishGlobal()
438 {
439 // called at the end of Exec
440 // save unique objects to the output file
441
442   fOutput->cd();
443   this->Write();
444   if (fCopyTreesFromInput > -1) {
445     fInputFiles[fCopyTreesFromInput]->Get("TE")->Clone()->Write();
446     gAlice->Write();
447   }
448   fOutput->Close();
449 }
450
451
452 ////////////////////////////////////////////////////////////////////////
453 Int_t  AliRunDigitizer::GetNParticles(Int_t event) const
454 {
455 // return number of particles in all input files for a given
456 // event (as numbered in the output file)
457 // return -1 if some file cannot be accessed
458
459   Int_t sum = 0;
460   Int_t sumI;
461   for (Int_t i = 0; i < fNinputs; i++) {
462     sumI = GetNParticles(GetInputEventNumber(event,i), i);
463     if (sumI < 0) return -1;
464     sum += sumI;
465   }
466   return sum;
467 }
468
469 ////////////////////////////////////////////////////////////////////////
470 Int_t  AliRunDigitizer::GetNParticles(Int_t event, Int_t input) const
471 {
472 // return number of particles in input file input for a given
473 // event (as numbered in this input file)
474 // return -1 if some error
475
476 // Must be revised in the version with AliStream
477
478   return -1;
479
480 /*
481   TFile *file = ConnectInputFile(input);
482   if (!file) {
483     Error("GetNParticles","Cannot open input file");
484     return -1;
485   }
486
487 // find the header and get Nprimaries and Nsecondaries
488   TTree* tE = (TTree *)file->Get("TE") ;
489   if (!tE) {
490     Error("GetNParticles","input file does not contain TE");
491     return -1;
492   }
493   AliHeader* header;
494   header = 0;
495   tE->SetBranchAddress("Header", &header);
496   if (!tE->GetEntry(event)) {
497     Error("GetNParticles","event %d not found",event);
498     return -1;
499   }
500   if (GetDebug()>2) {
501     cerr<<"Nprimary: "<< header->GetNprimary()<<endl;
502     cerr<<"Nsecondary: "<<header->GetNsecondary()<<endl;
503   }
504   return header->GetNprimary() + header->GetNsecondary();
505 */
506 }
507
508 ////////////////////////////////////////////////////////////////////////
509 Int_t* AliRunDigitizer::GetInputEventNumbers(Int_t event) const
510 {
511 // return pointer to an int array with input event numbers which were
512 // merged in the output event event
513
514 // simplified for now, implement later
515   Int_t * a = new Int_t[kMaxStreamsToMerge];
516   for (Int_t i = 0; i < fNinputs; i++) {
517     a[i] = event;
518   }
519   return a;
520 }
521 ////////////////////////////////////////////////////////////////////////
522 Int_t AliRunDigitizer::GetInputEventNumber(Int_t event, Int_t input) const
523 {
524 // return an event number of an eventInput from input file input
525 // which was merged to create output event event
526
527 // simplified for now, implement later
528   return event;
529 }
530 ////////////////////////////////////////////////////////////////////////
531 TParticle* AliRunDigitizer::GetParticle(Int_t i, Int_t event) const
532 {
533 // return pointer to particle with index i (index with mask)
534
535 // decode the MASK
536   Int_t input = i/fkMASKSTEP;
537   return GetParticle(i,input,GetInputEventNumber(event,input));
538 }
539
540 ////////////////////////////////////////////////////////////////////////
541 TParticle* AliRunDigitizer::GetParticle(Int_t i, Int_t input, Int_t event) const
542 {
543 // return pointer to particle with index i in the input file input
544 // (index without mask)
545 // event is the event number in the file input
546 // return 0 i fit does not exist
547
548 // Must be revised in the version with AliStream
549
550   return 0;
551 /*
552   TFile *file = ConnectInputFile(input);
553   if (!file) {
554     Error("GetParticle","Cannot open input file");
555     return 0;
556   }
557
558 // find the header and get Nprimaries and Nsecondaries
559   TTree* tE = (TTree *)file->Get("TE") ;
560   if (!tE) {
561     Error("GetParticle","input file does not contain TE");
562     return 0;
563   }
564   AliHeader* header;
565   header = 0;
566   tE->SetBranchAddress("Header", &header);
567   if (!tE->GetEntry(event)) {
568     Error("GetParticle","event %d not found",event);
569     return 0;
570   }
571   
572 // connect TreeK  
573   char treeName[30];
574   sprintf(treeName,"TreeK%d",event);  
575   TTree* tK = static_cast<TTree*>(file->Get(treeName));
576   if (!tK) {
577     Error("GetParticle","input file does not contain TreeK%d",event);
578     return 0;
579   }
580   TParticle *particleBuffer;
581   particleBuffer = 0;
582   tK->SetBranchAddress("Particles", &particleBuffer);
583
584
585 // algorithmic way of getting entry index
586 // (primary particles are filled after secondaries)
587   Int_t entry;
588   if (i<header->GetNprimary())
589     entry = i+header->GetNsecondary();
590   else 
591     entry = i-header->GetNprimary();
592   Int_t bytesRead = tK->GetEntry(entry);
593 //  new ((*fParticles)[nentries]) TParticle(*fParticleBuffer);
594   if (bytesRead)
595     return particleBuffer;
596   return  0;
597 */
598 }
599
600 ////////////////////////////////////////////////////////////////////////
601 void AliRunDigitizer::ExecuteTask(Option_t* option)
602 {
603 // overwrite ExecuteTask to do Digitize only
604
605   if (!IsActive()) return;
606   Digitize(option);
607   fHasExecuted = kTRUE;
608   return;
609 }