1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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
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
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.
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
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
33 Revision 1.7 2001/10/04 15:56:07 jchudoba
36 Revision 1.4 2001/09/19 06:23:50 jchudoba
37 Move some tasks to AliStream and AliMergeCombi classes
39 Revision 1.3 2001/07/30 14:04:18 jchudoba
40 correct bug in the initialization
42 Revision 1.2 2001/07/28 10:44:32 hristov
43 Loop variable declared once; typos corrected
45 Revision 1.1 2001/07/27 12:59:00 jchudoba
46 Manager class for merging/digitization
50 ////////////////////////////////////////////////////////////////////////
52 // AliRunDigitizer.cxx
54 // Manager object for merging/digitization
56 // Instance of this class manages the digitization and/or merging of
57 // Sdigits into Digits.
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.
80 // The default output is to the signal file (stream 0). This can be
81 // changed with the SetOutputFile(TString fn) method.
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.
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.
92 // Example with MUON digitizer, no merging, just digitization
94 // AliRunDigitizer * manager = new AliRunDigitizer(1,1);
95 // manager->SetInputStream(0,"galice.root");
96 // AliMUONDigitizer *dMUON = new AliMUONDigitizer(manager);
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)
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("");
110 // Example with MUON digitizer, save digits in a new file digits.root,
111 // process only 1 event
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("");
121 ////////////////////////////////////////////////////////////////////////
125 #include <iostream.h>
135 #include "AliRunDigitizer.h"
136 #include "AliDigitizer.h"
138 #include "AliHeader.h"
139 #include "TParticle.h"
140 #include "AliStream.h"
141 #include "AliMergeCombi.h"
143 ClassImp(AliRunDigitizer)
145 ////////////////////////////////////////////////////////////////////////
146 AliRunDigitizer::AliRunDigitizer()
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
151 // just set all pointers - data members to 0
157 for (Int_t i=0;i<kMaxStreamsToMerge;i++) {
158 fArrayTreeS[i]=fArrayTreeH[i]=fArrayTreeTPCS[i]=fArrayTreeTRDS[i]=NULL;
165 ////////////////////////////////////////////////////////////////////////
166 AliRunDigitizer::AliRunDigitizer(Int_t nInputStreams, Int_t sperb) : TTask("AliRunDigitizer","The manager for Merging")
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");
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;
182 fkMASKSTEP = 10000000;
184 for (i=1;i<kMaxStreamsToMerge;i++) {
185 fkMASK[i] = fkMASK[i-1] + fkMASKSTEP;
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");
196 fNrOfEventsToWrite = -1;
197 fNrOfEventsWritten = 0;
198 fCopyTreesFromInput = -1;
199 fCombi = new AliMergeCombi(nInputStreams,sperb);
204 fTreeDTPCBaseName = "TreeD_75x40_100x60_";
205 fTreeTPCSBaseName = "TreeS_75x40_100x60_";
207 for (i=0; i<kMaxStreamsToMerge; i++) fInputFiles[i]=0;
210 ////////////////////////////////////////////////////////////////////////
212 AliRunDigitizer::~AliRunDigitizer() {
216 delete fInputStreams;
225 ////////////////////////////////////////////////////////////////////////
226 void AliRunDigitizer::AddDigitizer(AliDigitizer *digitizer)
228 // add digitizer to the list of active digitizers
229 this->Add(digitizer);
231 ////////////////////////////////////////////////////////////////////////
232 void AliRunDigitizer::SetInputStream(Int_t i, char *inputFile)
234 if (i > fInputStreams->GetLast()) {
235 Error("SetInputStream","Input stream number too high");
238 static_cast<AliStream*>(fInputStreams->At(i))->AddFile(inputFile);
241 ////////////////////////////////////////////////////////////////////////
242 void AliRunDigitizer::Digitize(Option_t* option)
244 // get a new combination of inputs, connect input trees and loop
245 // over all digitizers
247 // take gAlice from the first input file. It is needed to access
249 // If gAlice is already in memory, use it
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;
258 cerr<<"False from InitGlobal"<<endl;
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;
266 // loop over all registered digitizers and let them do the work
267 ExecuteTasks(option);
274 ////////////////////////////////////////////////////////////////////////
275 Bool_t AliRunDigitizer::ConnectInputTrees()
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
284 Int_t eventNr[kMaxStreamsToMerge], delta[kMaxStreamsToMerge];
285 fCombi->Combination(eventNr, delta);
286 for (Int_t i=0;i<fNinputs;i++) {
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];
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];
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;
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;
318 fArrayTreeTRDS[i] = tree;
319 } else if (delta[i] != 0) {
320 Error("ConnectInputTrees","Only delta 0 or 1 is implemented");
327 ////////////////////////////////////////////////////////////////////////
328 Bool_t AliRunDigitizer::InitGlobal()
330 // called once before Digitize() is called, initialize digitizers and output
332 TList* subTasks = this->GetListOfTasks();
334 subTasks->ForEach(AliDigitizer,Init)();
339 ////////////////////////////////////////////////////////////////////////
341 void AliRunDigitizer::SetOutputFile(TString fn)
342 // the output will be to separate file, not to the signal file
344 fOutputFileName = fn;
345 (static_cast<AliStream*>(fInputStreams->At(0)))->ChangeMode("READ");
349 ////////////////////////////////////////////////////////////////////////
350 Bool_t AliRunDigitizer::InitOutputGlobal()
352 // Creates the output file, called by InitEvent()
355 fn = fOutputDirName + '/' + fOutputFileName;
356 fOutput = new TFile(fn,"recreate");
358 cerr<<"AliRunDigitizer::InitOutputGlobal(): file "<<fn.Data()<<" was opened"<<endl;
360 if (fOutput) return kTRUE;
361 Error("InitOutputGlobal","Could not create output file.");
366 ////////////////////////////////////////////////////////////////////////
367 void AliRunDigitizer::InitEvent()
369 // Creates TreeDxx in the output file, called from Digitize() once for
370 // each event. xx = fEvent
373 cerr<<"AliRunDigitizer::InitEvent: fEvent = "<<fEvent<<endl;
375 // if fOutputFileName was not given, write output to signal file
376 if (fOutputFileName == "") {
377 fOutput = (static_cast<AliStream*>(fInputStreams->At(0)))->CurrentFile();
381 sprintf(treeName,"TreeD%d",fEvent);
382 fTreeD = static_cast<TTree*>(fOutput->Get(treeName));
384 fTreeD = new TTree(treeName,"Digits");
385 fTreeD->Write(0,TObject::kOverwrite);
388 // special tree for TPC
389 sprintf(treeName,"%s%d",fTreeDTPCBaseName,fEvent);
390 fTreeDTPC = static_cast<TTree*>(fOutput->Get(treeName));
392 fTreeDTPC = new TTree(treeName,"TPC_Digits");
393 fTreeDTPC->Write(0,TObject::kOverwrite);
396 // special tree for TRD
397 sprintf(treeName,"TreeD%d_TRD",fEvent);
398 fTreeDTRD = static_cast<TTree*>(fOutput->Get(treeName));
400 fTreeDTRD = new TTree(treeName,"TRD_Digits");
401 fTreeDTRD->Write(0,TObject::kOverwrite);
406 ////////////////////////////////////////////////////////////////////////
407 void AliRunDigitizer::FinishEvent()
409 // called at the end of loop over digitizers
413 if (fCopyTreesFromInput > -1) {
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();
422 fNrOfEventsWritten++;
436 ////////////////////////////////////////////////////////////////////////
437 void AliRunDigitizer::FinishGlobal()
439 // called at the end of Exec
440 // save unique objects to the output file
444 if (fCopyTreesFromInput > -1) {
445 fInputFiles[fCopyTreesFromInput]->Get("TE")->Clone()->Write();
452 ////////////////////////////////////////////////////////////////////////
453 Int_t AliRunDigitizer::GetNParticles(Int_t event) const
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
461 for (Int_t i = 0; i < fNinputs; i++) {
462 sumI = GetNParticles(GetInputEventNumber(event,i), i);
463 if (sumI < 0) return -1;
469 ////////////////////////////////////////////////////////////////////////
470 Int_t AliRunDigitizer::GetNParticles(Int_t event, Int_t input) const
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
476 // Must be revised in the version with AliStream
481 TFile *file = ConnectInputFile(input);
483 Error("GetNParticles","Cannot open input file");
487 // find the header and get Nprimaries and Nsecondaries
488 TTree* tE = (TTree *)file->Get("TE") ;
490 Error("GetNParticles","input file does not contain TE");
495 tE->SetBranchAddress("Header", &header);
496 if (!tE->GetEntry(event)) {
497 Error("GetNParticles","event %d not found",event);
501 cerr<<"Nprimary: "<< header->GetNprimary()<<endl;
502 cerr<<"Nsecondary: "<<header->GetNsecondary()<<endl;
504 return header->GetNprimary() + header->GetNsecondary();
508 ////////////////////////////////////////////////////////////////////////
509 Int_t* AliRunDigitizer::GetInputEventNumbers(Int_t event) const
511 // return pointer to an int array with input event numbers which were
512 // merged in the output event event
514 // simplified for now, implement later
515 Int_t * a = new Int_t[kMaxStreamsToMerge];
516 for (Int_t i = 0; i < fNinputs; i++) {
521 ////////////////////////////////////////////////////////////////////////
522 Int_t AliRunDigitizer::GetInputEventNumber(Int_t event, Int_t input) const
524 // return an event number of an eventInput from input file input
525 // which was merged to create output event event
527 // simplified for now, implement later
530 ////////////////////////////////////////////////////////////////////////
531 TParticle* AliRunDigitizer::GetParticle(Int_t i, Int_t event) const
533 // return pointer to particle with index i (index with mask)
536 Int_t input = i/fkMASKSTEP;
537 return GetParticle(i,input,GetInputEventNumber(event,input));
540 ////////////////////////////////////////////////////////////////////////
541 TParticle* AliRunDigitizer::GetParticle(Int_t i, Int_t input, Int_t event) const
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
548 // Must be revised in the version with AliStream
552 TFile *file = ConnectInputFile(input);
554 Error("GetParticle","Cannot open input file");
558 // find the header and get Nprimaries and Nsecondaries
559 TTree* tE = (TTree *)file->Get("TE") ;
561 Error("GetParticle","input file does not contain TE");
566 tE->SetBranchAddress("Header", &header);
567 if (!tE->GetEntry(event)) {
568 Error("GetParticle","event %d not found",event);
574 sprintf(treeName,"TreeK%d",event);
575 TTree* tK = static_cast<TTree*>(file->Get(treeName));
577 Error("GetParticle","input file does not contain TreeK%d",event);
580 TParticle *particleBuffer;
582 tK->SetBranchAddress("Particles", &particleBuffer);
585 // algorithmic way of getting entry index
586 // (primary particles are filled after secondaries)
588 if (i<header->GetNprimary())
589 entry = i+header->GetNsecondary();
591 entry = i-header->GetNprimary();
592 Int_t bytesRead = tK->GetEntry(entry);
593 // new ((*fParticles)[nentries]) TParticle(*fParticleBuffer);
595 return particleBuffer;
600 ////////////////////////////////////////////////////////////////////////
601 void AliRunDigitizer::ExecuteTask(Option_t* option)
603 // overwrite ExecuteTask to do Digitize only
605 if (!IsActive()) return;
607 fHasExecuted = kTRUE;