]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - STEER/AliSimulation.cxx
29-jun-2007 NvE Unused variable dist2 removed from IceDwalk::Exec to prevent
[u/mrichter/AliRoot.git] / STEER / AliSimulation.cxx
... / ...
CommitLineData
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/* $Id$ */
17
18///////////////////////////////////////////////////////////////////////////////
19// //
20// class for running generation, simulation and digitization //
21// //
22// Hits, sdigits and digits are created for all detectors by typing: //
23// //
24// AliSimulation sim; //
25// sim.Run(); //
26// //
27// The Run method returns kTRUE in case of successful execution. //
28// The number of events can be given as argument to the Run method or it //
29// can be set by //
30// //
31// sim.SetNumberOfEvents(n); //
32// //
33// The name of the configuration file can be passed as argument to the //
34// AliSimulation constructor or can be specified by //
35// //
36// sim.SetConfigFile("..."); //
37// //
38// The generation of particles and the simulation of detector hits can be //
39// switched on or off by //
40// //
41// sim.SetRunGeneration(kTRUE); // generation of primary particles //
42// sim.SetRunSimulation(kFALSE); // but no tracking //
43// //
44// For which detectors sdigits and digits will be created, can be steered //
45// by //
46// //
47// sim.SetMakeSDigits("ALL"); // make sdigits for all detectors //
48// sim.SetMakeDigits("ITS TPC"); // make digits only for ITS and TPC //
49// //
50// The argument is a (case sensitive) string with the names of the //
51// detectors separated by a space. An empty string ("") can be used to //
52// disable the creation of sdigits or digits. The special string "ALL" //
53// selects all available detectors. This is the default. //
54// //
55// The creation of digits from hits instead of from sdigits can be selected //
56// by //
57// //
58// sim.SetMakeDigitsFromHits("TRD"); //
59// //
60// The argument is again a string with the selected detectors. Be aware that //
61// this feature is not available for all detectors and that merging is not //
62// possible, when digits are created directly from hits. //
63// //
64// Background events can be merged by calling //
65// //
66// sim.MergeWith("background/galice.root", 2); //
67// //
68// The first argument is the file name of the background galice file. The //
69// second argument is the number of signal events per background event. //
70// By default this number is calculated from the number of available //
71// background events. MergeWith can be called several times to merge more //
72// than two event streams. It is assumed that the sdigits were already //
73// produced for the background events. //
74// //
75// The output of raw data can be switched on by calling //
76// //
77// sim.SetWriteRawData("MUON"); // write raw data for MUON //
78// //
79// The default output format of the raw data are DDL files. They are //
80// converted to a DATE file, if a file name is given as second argument. //
81// For this conversion the program "dateStream" is required. If the file //
82// name has the extension ".root", the DATE file is converted to a root //
83// file. The program "alimdc" is used for this purpose. For the conversion //
84// to DATE and root format the two conversion programs have to be installed. //
85// Only the raw data in the final format is kept if the third argument is //
86// kTRUE. //
87// //
88// The methods RunSimulation, RunSDigitization, RunDigitization, //
89// RunHitsDigitization and WriteRawData can be used to run only parts of //
90// the full simulation chain. The creation of raw data DDL files and their //
91// conversion to the DATE or root format can be run directly by calling //
92// the methods WriteRawFiles, ConvertRawFilesToDate and ConvertDateToRoot. //
93// //
94// The default number of events per file, which is usually set in the //
95// config file, can be changed for individual detectors and data types //
96// by calling //
97// //
98// sim.SetEventsPerFile("PHOS", "Reconstructed Points", 3); //
99// //
100// The first argument is the detector, the second one the data type and the //
101// last one the number of events per file. Valid data types are "Hits", //
102// "Summable Digits", "Digits", "Reconstructed Points" and "Tracks". //
103// The number of events per file has to be set before the simulation of //
104// hits. Otherwise it has no effect. //
105// //
106///////////////////////////////////////////////////////////////////////////////
107
108#include <TVirtualMCApplication.h>
109#include <TGeoManager.h>
110#include <TObjString.h>
111#include <TStopwatch.h>
112#include <TSystem.h>
113#include <TFile.h>
114
115#include "AliCDBStorage.h"
116#include "AliCDBEntry.h"
117#include "AliCDBManager.h"
118#include "AliGeomManager.h"
119#include "AliAlignObj.h"
120#include "AliCentralTrigger.h"
121#include "AliDAQ.h"
122#include "AliDigitizer.h"
123#include "AliGenerator.h"
124#include "AliLog.h"
125#include "AliModule.h"
126#include "AliRun.h"
127#include "AliRunDigitizer.h"
128#include "AliRunLoader.h"
129#include "AliSimulation.h"
130#include "AliVertexGenFile.h"
131#include "AliCentralTrigger.h"
132#include "AliCTPRawData.h"
133#include "AliRawReaderFile.h"
134#include "AliESD.h"
135#include "AliHeader.h"
136#include "AliGenEventHeader.h"
137#include "AliMC.h"
138
139ClassImp(AliSimulation)
140
141AliSimulation *AliSimulation::fgInstance = 0;
142
143//_____________________________________________________________________________
144AliSimulation::AliSimulation(const char* configFileName, const char* cdbUri,
145 const char* name, const char* title) :
146 TNamed(name, title),
147
148 fRunGeneration(kTRUE),
149 fRunSimulation(kTRUE),
150 fLoadAlignFromCDB(kTRUE),
151 fLoadAlObjsListOfDets("ALL"),
152 fMakeSDigits("ALL"),
153 fMakeDigits("ALL"),
154 fMakeTrigger(""),
155 fMakeDigitsFromHits(""),
156 fWriteRawData(""),
157 fRawDataFileName(""),
158 fDeleteIntermediateFiles(kFALSE),
159 fStopOnError(kFALSE),
160
161 fNEvents(1),
162 fConfigFileName(configFileName),
163 fGAliceFileName("galice.root"),
164 fEventsPerFile(),
165 fBkgrdFileNames(NULL),
166 fAlignObjArray(NULL),
167 fUseBkgrdVertex(kTRUE),
168 fRegionOfInterest(kFALSE),
169 fCDBUri(cdbUri),
170 fSpecCDBUri(),
171 fEmbeddingFlag(kFALSE)
172{
173// create simulation object with default parameters
174 fgInstance = this;
175 SetGAliceFile("galice.root");
176}
177
178//_____________________________________________________________________________
179AliSimulation::AliSimulation(const AliSimulation& sim) :
180 TNamed(sim),
181
182 fRunGeneration(sim.fRunGeneration),
183 fRunSimulation(sim.fRunSimulation),
184 fLoadAlignFromCDB(sim.fLoadAlignFromCDB),
185 fLoadAlObjsListOfDets(sim.fLoadAlObjsListOfDets),
186 fMakeSDigits(sim.fMakeSDigits),
187 fMakeDigits(sim.fMakeDigits),
188 fMakeTrigger(sim.fMakeTrigger),
189 fMakeDigitsFromHits(sim.fMakeDigitsFromHits),
190 fWriteRawData(sim.fWriteRawData),
191 fRawDataFileName(""),
192 fDeleteIntermediateFiles(kFALSE),
193 fStopOnError(sim.fStopOnError),
194
195 fNEvents(sim.fNEvents),
196 fConfigFileName(sim.fConfigFileName),
197 fGAliceFileName(sim.fGAliceFileName),
198 fEventsPerFile(),
199 fBkgrdFileNames(NULL),
200 fAlignObjArray(NULL),
201 fUseBkgrdVertex(sim.fUseBkgrdVertex),
202 fRegionOfInterest(sim.fRegionOfInterest),
203 fCDBUri(sim.fCDBUri),
204 fSpecCDBUri(),
205 fEmbeddingFlag(sim.fEmbeddingFlag)
206{
207// copy constructor
208
209 for (Int_t i = 0; i < sim.fEventsPerFile.GetEntriesFast(); i++) {
210 if (!sim.fEventsPerFile[i]) continue;
211 fEventsPerFile.Add(sim.fEventsPerFile[i]->Clone());
212 }
213
214 fBkgrdFileNames = new TObjArray;
215 for (Int_t i = 0; i < sim.fBkgrdFileNames->GetEntriesFast(); i++) {
216 if (!sim.fBkgrdFileNames->At(i)) continue;
217 fBkgrdFileNames->Add(sim.fBkgrdFileNames->At(i)->Clone());
218 }
219
220 for (Int_t i = 0; i < sim.fSpecCDBUri.GetEntriesFast(); i++) {
221 if (sim.fSpecCDBUri[i]) fSpecCDBUri.Add(sim.fSpecCDBUri[i]->Clone());
222 }
223 fgInstance = this;
224}
225
226//_____________________________________________________________________________
227AliSimulation& AliSimulation::operator = (const AliSimulation& sim)
228{
229// assignment operator
230
231 this->~AliSimulation();
232 new(this) AliSimulation(sim);
233 return *this;
234}
235
236//_____________________________________________________________________________
237AliSimulation::~AliSimulation()
238{
239// clean up
240
241 fEventsPerFile.Delete();
242// if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
243// delete fAlignObjArray; fAlignObjArray=0;
244
245 if (fBkgrdFileNames) {
246 fBkgrdFileNames->Delete();
247 delete fBkgrdFileNames;
248 }
249
250 fSpecCDBUri.Delete();
251 if (fgInstance==this) fgInstance = 0;
252}
253
254
255//_____________________________________________________________________________
256void AliSimulation::SetNumberOfEvents(Int_t nEvents)
257{
258// set the number of events for one run
259
260 fNEvents = nEvents;
261}
262
263//_____________________________________________________________________________
264void AliSimulation::InitCDBStorage()
265{
266// activate a default CDB storage
267// First check if we have any CDB storage set, because it is used
268// to retrieve the calibration and alignment constants
269
270 AliCDBManager* man = AliCDBManager::Instance();
271 if (man->IsDefaultStorageSet())
272 {
273 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
274 AliWarning("Default CDB storage has been already set !");
275 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
276 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
277 fCDBUri = "";
278 }
279 else {
280 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
281 AliDebug(2, Form("Default CDB storage is set to: %s",fCDBUri.Data()));
282 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
283 man->SetDefaultStorage(fCDBUri);
284 }
285
286 // Now activate the detector specific CDB storage locations
287 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
288 TObject* obj = fSpecCDBUri[i];
289 if (!obj) continue;
290 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
291 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
292 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
293 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
294 }
295 man->Print();
296}
297
298//_____________________________________________________________________________
299void AliSimulation::SetDefaultStorage(const char* uri) {
300// Store the desired default CDB storage location
301// Activate it later within the Run() method
302
303 fCDBUri = uri;
304
305}
306
307//_____________________________________________________________________________
308void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
309// Store a detector-specific CDB storage location
310// Activate it later within the Run() method
311
312 AliCDBPath aPath(calibType);
313 if(!aPath.IsValid()){
314 AliError(Form("Not a valid path: %s", calibType));
315 return;
316 }
317
318 TObject* obj = fSpecCDBUri.FindObject(calibType);
319 if (obj) fSpecCDBUri.Remove(obj);
320 fSpecCDBUri.Add(new TNamed(calibType, uri));
321
322}
323
324//_____________________________________________________________________________
325void AliSimulation::SetConfigFile(const char* fileName)
326{
327// set the name of the config file
328
329 fConfigFileName = fileName;
330}
331
332//_____________________________________________________________________________
333void AliSimulation::SetGAliceFile(const char* fileName)
334{
335// set the name of the galice file
336// the path is converted to an absolute one if it is relative
337
338 fGAliceFileName = fileName;
339 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
340 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
341 fGAliceFileName);
342 fGAliceFileName = absFileName;
343 delete[] absFileName;
344 }
345
346 AliDebug(2, Form("galice file name set to %s", fileName));
347}
348
349//_____________________________________________________________________________
350void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
351 Int_t nEvents)
352{
353// set the number of events per file for the given detector and data type
354// ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
355
356 TNamed* obj = new TNamed(detector, type);
357 obj->SetUniqueID(nEvents);
358 fEventsPerFile.Add(obj);
359}
360
361//_____________________________________________________________________________
362Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
363{
364 // Read the alignment objects from CDB.
365 // Each detector is supposed to have the
366 // alignment objects in DET/Align/Data CDB path.
367 // All the detector objects are then collected,
368 // sorted by geometry level (starting from ALIC) and
369 // then applied to the TGeo geometry.
370 // Finally an overlaps check is performed.
371
372 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
373 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
374 return kFALSE;
375 }
376 Bool_t delRunLoader = kFALSE;
377 if (!runLoader) {
378 runLoader = LoadRun("READ");
379 if (!runLoader) return kFALSE;
380 delRunLoader = kTRUE;
381 }
382
383 // Export ideal geometry
384 AliGeomManager::GetGeometry()->Export("geometry.root");
385
386 // Load alignment data from CDB and apply to geometry through AliGeomManager
387 if(fLoadAlignFromCDB){
388
389 TString detStr = fLoadAlObjsListOfDets;
390 TString loadAlObjsListOfDets = "";
391
392 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
393 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
394 AliModule* det = (AliModule*) detArray->At(iDet);
395 if (!det || !det->IsActive()) continue;
396 if (IsSelected(det->GetName(), detStr)) {
397 //add det to list of dets to be aligned from CDB
398 loadAlObjsListOfDets += det->GetName();
399 loadAlObjsListOfDets += " ";
400 }
401 } // end loop over detectors
402 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
403 }else{
404 // Check if the array with alignment objects was
405 // provided by the user. If yes, apply the objects
406 // to the present TGeo geometry
407 if (fAlignObjArray) {
408 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
409 AliError("The misalignment of one or more volumes failed!"
410 "Compare the list of simulated detectors and the list of detector alignment data!");
411 if (delRunLoader) delete runLoader;
412 return kFALSE;
413 }
414 }
415 }
416
417 // Update the internal geometry of modules (ITS needs it)
418 TString detStr = fLoadAlObjsListOfDets;
419 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
420 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
421
422 AliModule* det = (AliModule*) detArray->At(iDet);
423 if (!det || !det->IsActive()) continue;
424 if (IsSelected(det->GetName(), detStr)) {
425 det->UpdateInternalGeometry();
426 }
427 } // end loop over detectors
428
429
430 if (delRunLoader) delete runLoader;
431
432 return kTRUE;
433}
434
435
436//_____________________________________________________________________________
437Bool_t AliSimulation::SetRunNumber()
438{
439 // Set the CDB manager run number
440 // The run number is retrieved from gAlice
441
442 if(AliCDBManager::Instance()->GetRun() < 0) {
443 AliRunLoader* runLoader = LoadRun("READ");
444 if (!runLoader) return kFALSE;
445 else {
446 AliCDBManager::Instance()->SetRun(runLoader->GetAliRun()->GetRunNumber());
447 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
448 delete runLoader;
449 }
450 }
451 return kTRUE;
452}
453
454//_____________________________________________________________________________
455void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
456{
457// add a file with background events for merging
458
459 TObjString* fileNameStr = new TObjString(fileName);
460 fileNameStr->SetUniqueID(nSignalPerBkgrd);
461 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
462 fBkgrdFileNames->Add(fileNameStr);
463}
464
465void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
466{
467// add a file with background events for embeddin
468 MergeWith(fileName, nSignalPerBkgrd);
469 fEmbeddingFlag = kTRUE;
470}
471
472//_____________________________________________________________________________
473Bool_t AliSimulation::Run(Int_t nEvents)
474{
475// run the generation, simulation and digitization
476
477 InitCDBStorage();
478
479 if (nEvents > 0) fNEvents = nEvents;
480
481 // generation and simulation -> hits
482 if (fRunGeneration) {
483 if (!RunSimulation()) if (fStopOnError) return kFALSE;
484 }
485
486 // Set run number in CDBManager (if it is not already set in RunSimulation)
487 if (!SetRunNumber()) if (fStopOnError) return kFALSE;
488
489 // If RunSimulation was not called, load the geometry and misalign it
490 if (!AliGeomManager::GetGeometry()) {
491 // Initialize the geometry manager
492 AliGeomManager::LoadGeometry("geometry.root");
493 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
494 // Misalign geometry
495 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
496 }
497
498 // hits -> summable digits
499 if (!fMakeSDigits.IsNull()) {
500 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
501 }
502
503 // summable digits -> digits
504 if (!fMakeDigits.IsNull()) {
505 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
506 if (fStopOnError) return kFALSE;
507 }
508 }
509
510 // hits -> digits
511 if (!fMakeDigitsFromHits.IsNull()) {
512 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
513 AliWarning(Form("Merging and direct creation of digits from hits "
514 "was selected for some detectors. "
515 "No merging will be done for the following detectors: %s",
516 fMakeDigitsFromHits.Data()));
517 }
518 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
519 if (fStopOnError) return kFALSE;
520 }
521 }
522
523 // digits -> trigger
524 if (!RunTrigger(fMakeTrigger)) {
525 if (fStopOnError) return kFALSE;
526 }
527
528 // digits -> raw data
529 if (!fWriteRawData.IsNull()) {
530 if (!WriteRawData(fWriteRawData, fRawDataFileName,
531 fDeleteIntermediateFiles)) {
532 if (fStopOnError) return kFALSE;
533 }
534 }
535
536 return kTRUE;
537}
538
539//_____________________________________________________________________________
540Bool_t AliSimulation::RunTrigger(const char* descriptors)
541{
542 // run the trigger
543
544 TStopwatch stopwatch;
545 stopwatch.Start();
546
547 AliRunLoader* runLoader = LoadRun("READ");
548 if (!runLoader) return kFALSE;
549 TString des = descriptors;
550
551 if (des.IsNull()) {
552 if (gAlice->GetTriggerDescriptor() != "") {
553 des = gAlice->GetTriggerDescriptor();
554 }
555 else {
556 AliWarning("No trigger descriptor is specified. Skipping the trigger simulation...");
557 return kTRUE;
558 }
559 }
560
561 runLoader->MakeTree( "GG" );
562 AliCentralTrigger* aCTP = runLoader->GetTrigger();
563 // Load Descriptors
564 aCTP->LoadDescriptor( des );
565
566 // digits -> trigger
567 if( !aCTP->RunTrigger( runLoader ) ) {
568 if (fStopOnError) {
569 // delete aCTP;
570 return kFALSE;
571 }
572 }
573
574 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
575 stopwatch.RealTime(),stopwatch.CpuTime()));
576
577 delete runLoader;
578
579 return kTRUE;
580}
581
582//_____________________________________________________________________________
583Bool_t AliSimulation::WriteTriggerRawData()
584{
585 // Writes the CTP (trigger) DDL raw data
586 // Details of the format are given in the
587 // trigger TDR - pages 134 and 135.
588 AliCTPRawData writer;
589 writer.RawData();
590
591 return kTRUE;
592}
593
594//_____________________________________________________________________________
595Bool_t AliSimulation::RunSimulation(Int_t nEvents)
596{
597// run the generation and simulation
598
599 TStopwatch stopwatch;
600 stopwatch.Start();
601
602 if (!gAlice) {
603 AliError("no gAlice object. Restart aliroot and try again.");
604 return kFALSE;
605 }
606 if (gAlice->Modules()->GetEntries() > 0) {
607 AliError("gAlice was already run. Restart aliroot and try again.");
608 return kFALSE;
609 }
610
611 AliInfo(Form("initializing gAlice with config file %s",
612 fConfigFileName.Data()));
613 StdoutToAliInfo(StderrToAliError(
614 gAlice->Init(fConfigFileName.Data());
615 ););
616
617 // Get the trigger descriptor string
618 // Either from AliSimulation or from
619 // gAlice
620 if (fMakeTrigger.IsNull()) {
621 if (gAlice->GetTriggerDescriptor() != "")
622 fMakeTrigger = gAlice->GetTriggerDescriptor();
623 }
624 else
625 gAlice->SetTriggerDescriptor(fMakeTrigger.Data());
626
627 // Set run number in CDBManager
628 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
629
630 AliRunLoader* runLoader = gAlice->GetRunLoader();
631 if (!runLoader) {
632 AliError(Form("gAlice has no run loader object. "
633 "Check your config file: %s", fConfigFileName.Data()));
634 return kFALSE;
635 }
636 SetGAliceFile(runLoader->GetFileName());
637
638 // Misalign geometry
639#if ROOT_VERSION_CODE < 331527
640 AliGeomManager::SetGeometry(gGeoManager);
641 MisalignGeometry(runLoader);
642#endif
643
644// AliRunLoader* runLoader = gAlice->GetRunLoader();
645// if (!runLoader) {
646// AliError(Form("gAlice has no run loader object. "
647// "Check your config file: %s", fConfigFileName.Data()));
648// return kFALSE;
649// }
650// SetGAliceFile(runLoader->GetFileName());
651
652 if (!gAlice->Generator()) {
653 AliError(Form("gAlice has no generator object. "
654 "Check your config file: %s", fConfigFileName.Data()));
655 return kFALSE;
656 }
657 if (nEvents <= 0) nEvents = fNEvents;
658
659 // get vertex from background file in case of merging
660 if (fUseBkgrdVertex &&
661 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
662 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
663 const char* fileName = ((TObjString*)
664 (fBkgrdFileNames->At(0)))->GetName();
665 AliInfo(Form("The vertex will be taken from the background "
666 "file %s with nSignalPerBackground = %d",
667 fileName, signalPerBkgrd));
668 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
669 gAlice->Generator()->SetVertexGenerator(vtxGen);
670 }
671
672 if (!fRunSimulation) {
673 gAlice->Generator()->SetTrackingFlag(0);
674 }
675
676 // set the number of events per file for given detectors and data types
677 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
678 if (!fEventsPerFile[i]) continue;
679 const char* detName = fEventsPerFile[i]->GetName();
680 const char* typeName = fEventsPerFile[i]->GetTitle();
681 TString loaderName(detName);
682 loaderName += "Loader";
683 AliLoader* loader = runLoader->GetLoader(loaderName);
684 if (!loader) {
685 AliError(Form("RunSimulation", "no loader for %s found\n"
686 "Number of events per file not set for %s %s",
687 detName, typeName, detName));
688 continue;
689 }
690 AliDataLoader* dataLoader =
691 loader->GetDataLoader(typeName);
692 if (!dataLoader) {
693 AliError(Form("no data loader for %s found\n"
694 "Number of events per file not set for %s %s",
695 typeName, detName, typeName));
696 continue;
697 }
698 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
699 AliDebug(1, Form("number of events per file set to %d for %s %s",
700 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
701 }
702
703 AliInfo("running gAlice");
704 StdoutToAliInfo(StderrToAliError(
705 gAlice->Run(nEvents);
706 ););
707
708 delete runLoader;
709
710 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
711 stopwatch.RealTime(),stopwatch.CpuTime()));
712
713 return kTRUE;
714}
715
716//_____________________________________________________________________________
717Bool_t AliSimulation::RunSDigitization(const char* detectors)
718{
719// run the digitization and produce summable digits
720
721 TStopwatch stopwatch;
722 stopwatch.Start();
723
724 AliRunLoader* runLoader = LoadRun();
725 if (!runLoader) return kFALSE;
726
727 TString detStr = detectors;
728 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
729 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
730 AliModule* det = (AliModule*) detArray->At(iDet);
731 if (!det || !det->IsActive()) continue;
732 if (IsSelected(det->GetName(), detStr)) {
733 AliInfo(Form("creating summable digits for %s", det->GetName()));
734 TStopwatch stopwatchDet;
735 stopwatchDet.Start();
736 det->Hits2SDigits();
737 AliInfo(Form("Execution time for %s: R:%.2fs C:%.2fs",
738 det->GetName(),stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
739 }
740 }
741
742 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
743 AliError(Form("the following detectors were not found: %s",
744 detStr.Data()));
745 if (fStopOnError) return kFALSE;
746 }
747
748 delete runLoader;
749
750 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
751 stopwatch.RealTime(),stopwatch.CpuTime()));
752
753 return kTRUE;
754}
755
756
757//_____________________________________________________________________________
758Bool_t AliSimulation::RunDigitization(const char* detectors,
759 const char* excludeDetectors)
760{
761// run the digitization and produce digits from sdigits
762
763 TStopwatch stopwatch;
764 stopwatch.Start();
765
766 while (AliRunLoader::GetRunLoader()) delete AliRunLoader::GetRunLoader();
767 if (gAlice) delete gAlice;
768 gAlice = NULL;
769
770 Int_t nStreams = 1;
771 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
772 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
773 AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
774 // manager->SetEmbeddingFlag(fEmbeddingFlag);
775 manager->SetInputStream(0, fGAliceFileName.Data());
776 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
777 const char* fileName = ((TObjString*)
778 (fBkgrdFileNames->At(iStream-1)))->GetName();
779 manager->SetInputStream(iStream, fileName);
780 }
781
782 TString detStr = detectors;
783 TString detExcl = excludeDetectors;
784 manager->GetInputStream(0)->ImportgAlice();
785 AliRunLoader* runLoader =
786 AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
787 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
788 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
789 AliModule* det = (AliModule*) detArray->At(iDet);
790 if (!det || !det->IsActive()) continue;
791 if (IsSelected(det->GetName(), detStr) &&
792 !IsSelected(det->GetName(), detExcl)) {
793 AliDigitizer* digitizer = det->CreateDigitizer(manager);
794
795 if (!digitizer) {
796 AliError(Form("no digitizer for %s", det->GetName()));
797 if (fStopOnError) return kFALSE;
798 } else {
799 digitizer->SetRegionOfInterest(fRegionOfInterest);
800 }
801 }
802 }
803
804 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
805 AliError(Form("the following detectors were not found: %s",
806 detStr.Data()));
807 if (fStopOnError) return kFALSE;
808 }
809
810 if (!manager->GetListOfTasks()->IsEmpty()) {
811 AliInfo("executing digitization");
812 manager->Exec("");
813 }
814
815 delete manager;
816
817 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
818 stopwatch.RealTime(),stopwatch.CpuTime()));
819
820 return kTRUE;
821}
822
823//_____________________________________________________________________________
824Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
825{
826// run the digitization and produce digits from hits
827
828 TStopwatch stopwatch;
829 stopwatch.Start();
830
831 AliRunLoader* runLoader = LoadRun("READ");
832 if (!runLoader) return kFALSE;
833
834 TString detStr = detectors;
835 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
836 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
837 AliModule* det = (AliModule*) detArray->At(iDet);
838 if (!det || !det->IsActive()) continue;
839 if (IsSelected(det->GetName(), detStr)) {
840 AliInfo(Form("creating digits from hits for %s", det->GetName()));
841 det->Hits2Digits();
842 }
843 }
844
845 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
846 AliError(Form("the following detectors were not found: %s",
847 detStr.Data()));
848 if (fStopOnError) return kFALSE;
849 }
850
851 delete runLoader;
852 //PH Temporary fix to avoid interference with the PHOS loder/getter
853 //PH The problem has to be solved in more general way 09/06/05
854
855 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
856 stopwatch.RealTime(),stopwatch.CpuTime()));
857
858 return kTRUE;
859}
860
861//_____________________________________________________________________________
862Bool_t AliSimulation::WriteRawData(const char* detectors,
863 const char* fileName,
864 Bool_t deleteIntermediateFiles)
865{
866// convert the digits to raw data
867// First DDL raw data files for the given detectors are created.
868// If a file name is given, the DDL files are then converted to a DATE file.
869// If deleteIntermediateFiles is true, the DDL raw files are deleted
870// afterwards.
871// If the file name has the extension ".root", the DATE file is converted
872// to a root file.
873// If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
874
875 TStopwatch stopwatch;
876 stopwatch.Start();
877
878 if (!WriteRawFiles(detectors)) {
879 if (fStopOnError) return kFALSE;
880 }
881
882 TString dateFileName(fileName);
883 if (!dateFileName.IsNull()) {
884 Bool_t rootOutput = dateFileName.EndsWith(".root");
885 if (rootOutput) dateFileName += ".date";
886 if (!ConvertRawFilesToDate(dateFileName)) {
887 if (fStopOnError) return kFALSE;
888 }
889 if (deleteIntermediateFiles) {
890 AliRunLoader* runLoader = LoadRun("READ");
891 if (runLoader) for (Int_t iEvent = 0;
892 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
893 char command[256];
894 sprintf(command, "rm -r raw%d", iEvent);
895 gSystem->Exec(command);
896 }
897 }
898
899 if (rootOutput) {
900 if (!ConvertDateToRoot(dateFileName, fileName)) {
901 if (fStopOnError) return kFALSE;
902 }
903 if (deleteIntermediateFiles) {
904 gSystem->Unlink(dateFileName);
905 }
906 }
907 }
908
909 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
910 stopwatch.RealTime(),stopwatch.CpuTime()));
911
912 return kTRUE;
913}
914
915//_____________________________________________________________________________
916Bool_t AliSimulation::WriteRawFiles(const char* detectors)
917{
918// convert the digits to raw data DDL files
919
920 AliRunLoader* runLoader = LoadRun("READ");
921 if (!runLoader) return kFALSE;
922
923 // write raw data to DDL files
924 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
925 AliInfo(Form("processing event %d", iEvent));
926 runLoader->GetEvent(iEvent);
927 TString baseDir = gSystem->WorkingDirectory();
928 char dirName[256];
929 sprintf(dirName, "raw%d", iEvent);
930 gSystem->MakeDirectory(dirName);
931 if (!gSystem->ChangeDirectory(dirName)) {
932 AliError(Form("couldn't change to directory %s", dirName));
933 if (fStopOnError) return kFALSE; else continue;
934 }
935
936 TString detStr = detectors;
937 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
938 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
939 AliModule* det = (AliModule*) detArray->At(iDet);
940 if (!det || !det->IsActive()) continue;
941 if (IsSelected(det->GetName(), detStr)) {
942 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
943 det->Digits2Raw();
944 }
945 }
946
947 if (!WriteTriggerRawData())
948 if (fStopOnError) return kFALSE;
949
950 gSystem->ChangeDirectory(baseDir);
951 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
952 AliError(Form("the following detectors were not found: %s",
953 detStr.Data()));
954 if (fStopOnError) return kFALSE;
955 }
956 }
957
958 delete runLoader;
959 return kTRUE;
960}
961
962//_____________________________________________________________________________
963Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName)
964{
965// convert raw data DDL files to a DATE file with the program "dateStream"
966
967 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
968 if (!path) {
969 AliError("the program dateStream was not found");
970 if (fStopOnError) return kFALSE;
971 } else {
972 delete[] path;
973 }
974
975 AliRunLoader* runLoader = LoadRun("READ");
976 if (!runLoader) return kFALSE;
977
978 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
979 char command[256];
980 // Note the option -s. It is used in order to avoid
981 // the generation of SOR/EOR events.
982 sprintf(command, "dateStream -s -D -o %s -# %d -C",
983 dateFileName, runLoader->GetNumberOfEvents());
984 FILE* pipe = gSystem->OpenPipe(command, "w");
985
986 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
987 fprintf(pipe, "GDC\n");
988 Float_t ldc = 0;
989 Int_t prevLDC = -1;
990
991 // loop over detectors and DDLs
992 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
993 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
994
995 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
996 Int_t ldcID = Int_t(ldc + 0.0001);
997 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
998
999 char rawFileName[256];
1000 sprintf(rawFileName, "raw%d/%s",
1001 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1002
1003 // check existence and size of raw data file
1004 FILE* file = fopen(rawFileName, "rb");
1005 if (!file) continue;
1006 fseek(file, 0, SEEK_END);
1007 unsigned long size = ftell(file);
1008 fclose(file);
1009 if (!size) continue;
1010
1011 if (ldcID != prevLDC) {
1012 fprintf(pipe, " LDC Id %d\n", ldcID);
1013 prevLDC = ldcID;
1014 }
1015 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1016 }
1017 }
1018 }
1019
1020 Int_t result = gSystem->ClosePipe(pipe);
1021
1022 delete runLoader;
1023 return (result == 0);
1024}
1025
1026//_____________________________________________________________________________
1027Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1028 const char* rootFileName)
1029{
1030// convert a DATE file to a root file with the program "alimdc"
1031
1032 // ALIMDC setup
1033 const Int_t kDBSize = 2000000000;
1034 const Int_t kTagDBSize = 1000000000;
1035 const Bool_t kFilter = kFALSE;
1036 const Int_t kCompression = 0;
1037
1038 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1039 if (!path) {
1040 AliError("the program alimdc was not found");
1041 if (fStopOnError) return kFALSE;
1042 } else {
1043 delete[] path;
1044 }
1045
1046 AliInfo(Form("converting DATE file %s to root file %s",
1047 dateFileName, rootFileName));
1048
1049 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1050 const char* tagDBFS = "/tmp/mdc1/tags";
1051
1052 // User defined file system locations
1053 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1054 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1055 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1056 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1057 if (gSystem->Getenv("ALIMDC_TAGDB"))
1058 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1059
1060 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1061 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1062 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1063
1064 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1065 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1066 gSystem->Exec(Form("mkdir %s",tagDBFS));
1067
1068 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1069 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1070 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1071
1072 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1073 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1074 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1075
1076 return (result == 0);
1077}
1078
1079
1080//_____________________________________________________________________________
1081AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1082{
1083// delete existing run loaders, open a new one and load gAlice
1084
1085 while (AliRunLoader::GetRunLoader()) delete AliRunLoader::GetRunLoader();
1086 AliRunLoader* runLoader =
1087 AliRunLoader::Open(fGAliceFileName.Data(),
1088 AliConfig::GetDefaultEventFolderName(), mode);
1089 if (!runLoader) {
1090 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1091 return NULL;
1092 }
1093 runLoader->LoadgAlice();
1094 gAlice = runLoader->GetAliRun();
1095 if (!gAlice) {
1096 AliError(Form("no gAlice object found in file %s",
1097 fGAliceFileName.Data()));
1098 return NULL;
1099 }
1100 return runLoader;
1101}
1102
1103//_____________________________________________________________________________
1104Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1105{
1106// get or calculate the number of signal events per background event
1107
1108 if (!fBkgrdFileNames) return 1;
1109 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1110 if (nBkgrdFiles == 0) return 1;
1111
1112 // get the number of signal events
1113 if (nEvents <= 0) {
1114 AliRunLoader* runLoader =
1115 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1116 if (!runLoader) return 1;
1117
1118 nEvents = runLoader->GetNumberOfEvents();
1119 delete runLoader;
1120 }
1121
1122 Int_t result = 0;
1123 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1124 // get the number of background events
1125 const char* fileName = ((TObjString*)
1126 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1127 AliRunLoader* runLoader =
1128 AliRunLoader::Open(fileName, "BKGRD");
1129 if (!runLoader) continue;
1130 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1131 delete runLoader;
1132
1133 // get or calculate the number of signal per background events
1134 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1135 if (nSignalPerBkgrd <= 0) {
1136 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1137 } else if (result && (result != nSignalPerBkgrd)) {
1138 AliInfo(Form("the number of signal events per background event "
1139 "will be changed from %d to %d for stream %d",
1140 nSignalPerBkgrd, result, iBkgrdFile+1));
1141 nSignalPerBkgrd = result;
1142 }
1143
1144 if (!result) result = nSignalPerBkgrd;
1145 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1146 AliWarning(Form("not enough background events (%d) for %d signal events "
1147 "using %d signal per background events for stream %d",
1148 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1149 }
1150 }
1151
1152 return result;
1153}
1154
1155//_____________________________________________________________________________
1156Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1157{
1158// check whether detName is contained in detectors
1159// if yes, it is removed from detectors
1160
1161 // check if all detectors are selected
1162 if ((detectors.CompareTo("ALL") == 0) ||
1163 detectors.BeginsWith("ALL ") ||
1164 detectors.EndsWith(" ALL") ||
1165 detectors.Contains(" ALL ")) {
1166 detectors = "ALL";
1167 return kTRUE;
1168 }
1169
1170 // search for the given detector
1171 Bool_t result = kFALSE;
1172 if ((detectors.CompareTo(detName) == 0) ||
1173 detectors.BeginsWith(detName+" ") ||
1174 detectors.EndsWith(" "+detName) ||
1175 detectors.Contains(" "+detName+" ")) {
1176 detectors.ReplaceAll(detName, "");
1177 result = kTRUE;
1178 }
1179
1180 // clean up the detectors string
1181 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1182 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1183 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1184
1185 return result;
1186}
1187
1188Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName)
1189{
1190//
1191// Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1192// These can be used for embedding of MC tracks into RAW data using the standard
1193// merging procedure.
1194//
1195// If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1196//
1197 if (!gAlice) {
1198 AliError("no gAlice object. Restart aliroot and try again.");
1199 return kFALSE;
1200 }
1201 if (gAlice->Modules()->GetEntries() > 0) {
1202 AliError("gAlice was already run. Restart aliroot and try again.");
1203 return kFALSE;
1204 }
1205
1206 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1207 StdoutToAliInfo(StderrToAliError(gAlice->Init(fConfigFileName.Data());););
1208//
1209// Initialize CDB
1210 InitCDBStorage();
1211 AliCDBManager* man = AliCDBManager::Instance();
1212 man->SetRun(0); // Should this come from rawdata header ?
1213
1214 Int_t iDet;
1215 //
1216 // Get the runloader
1217 AliRunLoader* runLoader = gAlice->GetRunLoader();
1218 //
1219 // Open esd file if available
1220 TFile* esdFile = TFile::Open(esdFileName);
1221 Bool_t esdOK = (esdFile != 0);
1222 AliESD* esd = new AliESD;
1223 TTree* treeESD = 0;
1224 if (esdOK) {
1225 treeESD = (TTree*) esdFile->Get("esdTree");
1226 if (!treeESD) {
1227 AliWarning("No ESD tree found");
1228 esdOK = kFALSE;
1229 } else {
1230 treeESD->SetBranchAddress("ESD", &esd);
1231 }
1232 }
1233 //
1234 // Create the RawReader
1235 AliRawReaderFile* rawReader = new AliRawReaderFile(rawDirectory);
1236 //
1237 // Get list of detectors
1238 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1239 //
1240 // Get Header
1241 AliHeader* header = runLoader->GetHeader();
1242 //
1243 // Event loop
1244 Int_t nev = 0;
1245 while(kTRUE) {
1246 if (!(rawReader->NextEvent())) break;
1247 //
1248 // Detector loop
1249 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1250 AliModule* det = (AliModule*) detArray->At(iDet);
1251 AliInfo(Form("Calling Raw2SDigits for %s\n", det->GetName()));
1252 det->Raw2SDigits(rawReader);
1253 rawReader->Reset();
1254 } // detectors
1255
1256 //
1257 // If ESD information available obtain reconstructed vertex and store in header.
1258 if (esdOK) {
1259 treeESD->GetEvent(nev);
1260 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
1261 Double_t position[3];
1262 esdVertex->GetXYZ(position);
1263 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
1264 TArrayF mcV;
1265 mcV.Set(3);
1266 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
1267 mcHeader->SetPrimaryVertex(mcV);
1268 header->Reset(0,nev);
1269 header->SetGenEventHeader(mcHeader);
1270 printf("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]);
1271 }
1272 nev++;
1273//
1274// Finish the event
1275 runLoader->TreeE()->Fill();
1276 runLoader->SetNextEvent();
1277 } // events
1278
1279 delete rawReader;
1280//
1281// Finish the run
1282 runLoader->CdGAFile();
1283 runLoader->WriteHeader("OVERWRITE");
1284 runLoader->WriteRunLoader();
1285
1286 return kTRUE;
1287}