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