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