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