]> git.uio.no Git - u/mrichter/AliRoot.git/blame - STEER/AliReconstruction.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
The CTP (trigger) infor is written into the ESD. The corresponding entries are fTrigg...
[u/mrichter/AliRoot.git] / STEER / AliReconstruction.cxx
CommitLineData
596a855f 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 the reconstruction //
21// //
22// Clusters and tracks are created for all detectors and all events by //
23// typing: //
24// //
25// AliReconstruction rec; //
26// rec.Run(); //
27// //
28// The Run method returns kTRUE in case of successful execution. //
c71de921 29// //
30// If the input to the reconstruction are not simulated digits but raw data, //
31// this can be specified by an argument of the Run method or by the method //
32// //
33// rec.SetInput("..."); //
34// //
35// The input formats and the corresponding argument are: //
36// - DDL raw data files: directory name, ends with "/" //
37// - raw data root file: root file name, extension ".root" //
38// - raw data DATE file: DATE file name, any other non-empty string //
39// - MC root files : empty string, default //
40// //
b26c3770 41// By default all events are reconstructed. The reconstruction can be //
42// limited to a range of events by giving the index of the first and the //
43// last event as an argument to the Run method or by calling //
44// //
45// rec.SetEventRange(..., ...); //
46// //
47// The index -1 (default) can be used for the last event to indicate no //
48// upper limit of the event range. //
49// //
596a855f 50// The name of the galice file can be changed from the default //
e583c30d 51// "galice.root" by passing it as argument to the AliReconstruction //
52// constructor or by //
596a855f 53// //
54// rec.SetGAliceFile("..."); //
55// //
59697224 56// The local reconstruction can be switched on or off for individual //
57// detectors by //
596a855f 58// //
59697224 59// rec.SetRunLocalReconstruction("..."); //
596a855f 60// //
61// The argument is a (case sensitive) string with the names of the //
62// detectors separated by a space. The special string "ALL" selects all //
63// available detectors. This is the default. //
64// //
c71de921 65// The reconstruction of the primary vertex position can be switched off by //
66// //
67// rec.SetRunVertexFinder(kFALSE); //
68// //
b8cd5251 69// The tracking and the creation of ESD tracks can be switched on for //
70// selected detectors by //
596a855f 71// //
b8cd5251 72// rec.SetRunTracking("..."); //
596a855f 73// //
c84a5e9e 74// Uniform/nonuniform field tracking switches (default: uniform field) //
75// //
1d99986f 76// rec.SetUniformFieldTracking(); ( rec.SetUniformFieldTracking(kFALSE); ) //
c84a5e9e 77// //
596a855f 78// The filling of additional ESD information can be steered by //
79// //
80// rec.SetFillESD("..."); //
81// //
b8cd5251 82// Again, for both methods the string specifies the list of detectors. //
83// The default is "ALL". //
84// //
85// The call of the shortcut method //
86// //
87// rec.SetRunReconstruction("..."); //
88// //
89// is equivalent to calling SetRunLocalReconstruction, SetRunTracking and //
90// SetFillESD with the same detector selecting string as argument. //
596a855f 91// //
c71de921 92// The reconstruction requires digits or raw data as input. For the creation //
93// of digits and raw data have a look at the class AliSimulation. //
596a855f 94// //
24f7a148 95// For debug purposes the method SetCheckPointLevel can be used. If the //
96// argument is greater than 0, files with ESD events will be written after //
97// selected steps of the reconstruction for each event: //
98// level 1: after tracking and after filling of ESD (final) //
99// level 2: in addition after each tracking step //
100// level 3: in addition after the filling of ESD for each detector //
101// If a final check point file exists for an event, this event will be //
102// skipped in the reconstruction. The tracking and the filling of ESD for //
103// a detector will be skipped as well, if the corresponding check point //
104// file exists. The ESD event will then be loaded from the file instead. //
105// //
596a855f 106///////////////////////////////////////////////////////////////////////////////
107
024a7e64 108#include <TArrayF.h>
109#include <TFile.h>
110#include <TSystem.h>
111#include <TROOT.h>
112#include <TPluginManager.h>
fd46e2d2 113#include <TStopwatch.h>
3103d196 114#include <TGeoManager.h>
2bdb9d38 115#include <TLorentzVector.h>
596a855f 116
117#include "AliReconstruction.h"
b8cd5251 118#include "AliReconstructor.h"
815c2b38 119#include "AliLog.h"
596a855f 120#include "AliRunLoader.h"
121#include "AliRun.h"
b649205a 122#include "AliRawReaderFile.h"
123#include "AliRawReaderDate.h"
124#include "AliRawReaderRoot.h"
596a855f 125#include "AliESD.h"
1d99986f 126#include "AliESDfriend.h"
2257f27e 127#include "AliESDVertex.h"
c84a5e9e 128#include "AliTracker.h"
2257f27e 129#include "AliVertexer.h"
596a855f 130#include "AliHeader.h"
131#include "AliGenEventHeader.h"
b26c3770 132#include "AliPID.h"
596a855f 133#include "AliESDpid.h"
ff8bb5ae 134#include "AliESDtrack.h"
f3a97c86 135
136#include "AliRunTag.h"
f3a97c86 137#include "AliDetectorTag.h"
138#include "AliEventTag.h"
139
98937d93 140#include "AliTrackPointArray.h"
b0314964 141#include "AliCDBManager.h"
6bae477a 142#include "AliCDBEntry.h"
143#include "AliAlignObj.h"
f3a97c86 144
b647652d 145#include "AliCentralTrigger.h"
146#include "AliCTPRawStream.h"
147
596a855f 148ClassImp(AliReconstruction)
149
150
c757bafd 151//_____________________________________________________________________________
482070f2 152const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "RICH", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "START", "VZERO", "CRT", "HLT"};
c757bafd 153
596a855f 154//_____________________________________________________________________________
024cf675 155AliReconstruction::AliReconstruction(const char* gAliceFilename, const char* cdbUri,
e583c30d 156 const char* name, const char* title) :
157 TNamed(name, title),
158
c84a5e9e 159 fUniformField(kTRUE),
2257f27e 160 fRunVertexFinder(kTRUE),
1f46a9ae 161 fRunHLTTracking(kFALSE),
1d99986f 162 fStopOnError(kFALSE),
163 fWriteAlignmentData(kFALSE),
164 fWriteESDfriend(kFALSE),
b647652d 165 fFillTriggerESD(kTRUE),
1d99986f 166
167 fRunLocalReconstruction("ALL"),
b8cd5251 168 fRunTracking("ALL"),
e583c30d 169 fFillESD("ALL"),
170 fGAliceFileName(gAliceFilename),
b649205a 171 fInput(""),
b26c3770 172 fFirstEvent(0),
173 fLastEvent(-1),
24f7a148 174 fCheckPointLevel(0),
b8cd5251 175 fOptions(),
6bae477a 176 fLoadAlignFromCDB(kTRUE),
177 fLoadAlignData("ALL"),
e583c30d 178
179 fRunLoader(NULL),
b649205a 180 fRawReader(NULL),
b8cd5251 181
98937d93 182 fVertexer(NULL),
183
6bae477a 184 fAlignObjArray(NULL),
024cf675 185 fCDBUri(cdbUri)
596a855f 186{
187// create reconstruction object with default parameters
b8cd5251 188
189 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
190 fReconstructor[iDet] = NULL;
191 fLoader[iDet] = NULL;
192 fTracker[iDet] = NULL;
193 }
e47c4c2e 194 AliPID pid;
596a855f 195}
196
197//_____________________________________________________________________________
198AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
e583c30d 199 TNamed(rec),
200
c84a5e9e 201 fUniformField(rec.fUniformField),
2257f27e 202 fRunVertexFinder(rec.fRunVertexFinder),
1f46a9ae 203 fRunHLTTracking(rec.fRunHLTTracking),
1d99986f 204 fStopOnError(rec.fStopOnError),
205 fWriteAlignmentData(rec.fWriteAlignmentData),
206 fWriteESDfriend(rec.fWriteESDfriend),
b647652d 207 fFillTriggerESD(rec.fFillTriggerESD),
1d99986f 208
209 fRunLocalReconstruction(rec.fRunLocalReconstruction),
e583c30d 210 fRunTracking(rec.fRunTracking),
211 fFillESD(rec.fFillESD),
212 fGAliceFileName(rec.fGAliceFileName),
b649205a 213 fInput(rec.fInput),
b26c3770 214 fFirstEvent(rec.fFirstEvent),
215 fLastEvent(rec.fLastEvent),
24f7a148 216 fCheckPointLevel(0),
b8cd5251 217 fOptions(),
6bae477a 218 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
219 fLoadAlignData(rec.fLoadAlignData),
e583c30d 220
221 fRunLoader(NULL),
b649205a 222 fRawReader(NULL),
b8cd5251 223
98937d93 224 fVertexer(NULL),
225
6bae477a 226 fAlignObjArray(rec.fAlignObjArray),
024cf675 227 fCDBUri(rec.fCDBUri)
596a855f 228{
229// copy constructor
230
efd2085e 231 for (Int_t i = 0; i < fOptions.GetEntriesFast(); i++) {
232 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
233 }
b8cd5251 234 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
235 fReconstructor[iDet] = NULL;
236 fLoader[iDet] = NULL;
237 fTracker[iDet] = NULL;
238 }
596a855f 239}
240
241//_____________________________________________________________________________
242AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
243{
244// assignment operator
245
246 this->~AliReconstruction();
247 new(this) AliReconstruction(rec);
248 return *this;
249}
250
251//_____________________________________________________________________________
252AliReconstruction::~AliReconstruction()
253{
254// clean up
255
e583c30d 256 CleanUp();
efd2085e 257 fOptions.Delete();
596a855f 258}
259
024cf675 260//_____________________________________________________________________________
261void AliReconstruction::InitCDBStorage()
262{
263// activate a default CDB storage
264// First check if we have any CDB storage set, because it is used
265// to retrieve the calibration and alignment constants
266
267 AliCDBManager* man = AliCDBManager::Instance();
268 if (!man->IsDefaultStorageSet())
269 {
270 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
271 AliWarning("Default CDB storage not yet set");
272 AliWarning(Form("Using default storage declared in AliSimulation: %s",fCDBUri.Data()));
273 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
274 SetDefaultStorage(fCDBUri);
e08887cd 275
276 Int_t cdbRun = AliCDBManager::Instance()->GetRun();
277 if(cdbRun == -1){
278 AliWarning("AliCDBManager's run number temporarily set to 0!!");
279 AliCDBManager::Instance()->SetRun(0);
280 }
024cf675 281 }
282
283}
284
285//_____________________________________________________________________________
286void AliReconstruction::SetDefaultStorage(const char* uri) {
287// activate a default CDB storage
288
289 AliCDBManager::Instance()->SetDefaultStorage(uri);
290
291}
292
293//_____________________________________________________________________________
294void AliReconstruction::SetSpecificStorage(const char* detName, const char* uri) {
295// activate a detector-specific CDB storage
296
297 AliCDBManager::Instance()->SetSpecificStorage(detName, uri);
298
299}
300
6bae477a 301//_____________________________________________________________________________
302Bool_t AliReconstruction::SetRunNumber()
303{
304 // The method is called in Run() in order
305 // to set a correct run number.
306 // In case of raw data reconstruction the
307 // run number is taken from the raw data header
308
309 if(AliCDBManager::Instance()->GetRun() < 0) {
310 if (!fRunLoader) {
311 AliError("No run loader is found !");
312 return kFALSE;
313 }
314 // read run number from gAlice
315 AliCDBManager::Instance()->SetRun(fRunLoader->GetAliRun()->GetRunNumber());
316 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
317 }
318 return kTRUE;
319}
320
321//_____________________________________________________________________________
322Bool_t AliReconstruction::ApplyAlignObjsToGeom(TObjArray* alObjArray)
323{
324 // Read collection of alignment objects (AliAlignObj derived) saved
325 // in the TClonesArray ClArrayName and apply them to the geometry
326 // manager singleton.
327 //
328 alObjArray->Sort();
329 Int_t nvols = alObjArray->GetEntriesFast();
330
dc0984f8 331 Bool_t flag = kTRUE;
332
6bae477a 333 for(Int_t j=0; j<nvols; j++)
334 {
335 AliAlignObj* alobj = (AliAlignObj*) alObjArray->UncheckedAt(j);
dc0984f8 336 if (alobj->ApplyToGeometry() == kFALSE) flag = kFALSE;
6bae477a 337 }
338
339 if (AliDebugLevelClass() >= 1) {
340 gGeoManager->CheckOverlaps(20);
341 TObjArray* ovexlist = gGeoManager->GetListOfOverlaps();
342 if(ovexlist->GetEntriesFast()){
343 AliError("The application of alignment objects to the geometry caused huge overlaps/extrusions!");
344 }
345 }
346
dc0984f8 347 return flag;
6bae477a 348
349}
350
351//_____________________________________________________________________________
352Bool_t AliReconstruction::SetAlignObjArraySingleDet(const char* detName)
353{
354 // Fills array of single detector's alignable objects from CDB
355
356 AliDebug(2, Form("Loading alignment data for detector: %s",detName));
357
358 AliCDBEntry *entry;
359
360 AliCDBPath path(detName,"Align","Data");
361
362 entry=AliCDBManager::Instance()->Get(path.GetPath());
363 if(!entry){
364 AliDebug(2,Form("Couldn't load alignment data for detector %s",detName));
365 return kFALSE;
366 }
367 entry->SetOwner(1);
368 TClonesArray *alignArray = (TClonesArray*) entry->GetObject();
369 alignArray->SetOwner(0);
370 AliDebug(2,Form("Found %d alignment objects for %s",
371 alignArray->GetEntries(),detName));
372
373 AliAlignObj *alignObj=0;
374 TIter iter(alignArray);
375
376 // loop over align objects in detector
377 while( ( alignObj=(AliAlignObj *) iter.Next() ) ){
378 fAlignObjArray->Add(alignObj);
379 }
380 // delete entry --- Don't delete, it is cached!
381
382 AliDebug(2, Form("fAlignObjArray entries: %d",fAlignObjArray->GetEntries() ));
383 return kTRUE;
384
385}
386
387//_____________________________________________________________________________
388Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
389{
390 // Read the alignment objects from CDB.
391 // Each detector is supposed to have the
392 // alignment objects in DET/Align/Data CDB path.
393 // All the detector objects are then collected,
394 // sorted by geometry level (starting from ALIC) and
395 // then applied to the TGeo geometry.
396 // Finally an overlaps check is performed.
397
398 // Load alignment data from CDB and fill fAlignObjArray
399 if(fLoadAlignFromCDB){
400 if(!fAlignObjArray) fAlignObjArray = new TObjArray();
401
402 //fAlignObjArray->RemoveAll();
403 fAlignObjArray->Clear();
404 fAlignObjArray->SetOwner(0);
405
406 TString detStr = detectors;
407 TString dataNotLoaded="";
408 TString dataLoaded="";
409
410 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
411 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
412 if(!SetAlignObjArraySingleDet(fgkDetectorName[iDet])){
413 dataNotLoaded += fgkDetectorName[iDet];
414 dataNotLoaded += " ";
415 } else {
416 dataLoaded += fgkDetectorName[iDet];
417 dataLoaded += " ";
418 }
419 } // end loop over detectors
420
421 if ((detStr.CompareTo("ALL") == 0)) detStr = "";
422 dataNotLoaded += detStr;
423 AliInfo(Form("Alignment data loaded for: %s",
424 dataLoaded.Data()));
425 AliInfo(Form("Didn't/couldn't load alignment data for: %s",
426 dataNotLoaded.Data()));
427 } // fLoadAlignFromCDB flag
428
429 // Check if the array with alignment objects was
430 // provided by the user. If yes, apply the objects
431 // to the present TGeo geometry
432 if (fAlignObjArray) {
433 if (gGeoManager && gGeoManager->IsClosed()) {
434 if (ApplyAlignObjsToGeom(fAlignObjArray) == kFALSE) {
dc0984f8 435 AliError("The misalignment of one or more volumes failed!"
436 "Compare the list of simulated detectors and the list of detector alignment data!");
6bae477a 437 return kFALSE;
438 }
439 }
440 else {
441 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
442 return kFALSE;
443 }
444 }
445
446 return kTRUE;
447}
596a855f 448
449//_____________________________________________________________________________
450void AliReconstruction::SetGAliceFile(const char* fileName)
451{
452// set the name of the galice file
453
454 fGAliceFileName = fileName;
455}
456
efd2085e 457//_____________________________________________________________________________
458void AliReconstruction::SetOption(const char* detector, const char* option)
459{
460// set options for the reconstruction of a detector
461
462 TObject* obj = fOptions.FindObject(detector);
463 if (obj) fOptions.Remove(obj);
464 fOptions.Add(new TNamed(detector, option));
465}
466
596a855f 467
468//_____________________________________________________________________________
b26c3770 469Bool_t AliReconstruction::Run(const char* input,
470 Int_t firstEvent, Int_t lastEvent)
596a855f 471{
472// run the reconstruction
473
024cf675 474 InitCDBStorage();
b0314964 475
b649205a 476 // set the input
477 if (!input) input = fInput.Data();
478 TString fileName(input);
479 if (fileName.EndsWith("/")) {
480 fRawReader = new AliRawReaderFile(fileName);
481 } else if (fileName.EndsWith(".root")) {
482 fRawReader = new AliRawReaderRoot(fileName);
483 } else if (!fileName.IsNull()) {
484 fRawReader = new AliRawReaderDate(fileName);
485 fRawReader->SelectEvents(7);
486 }
487
f08fc9f5 488 // get the run loader
489 if (!InitRunLoader()) return kFALSE;
596a855f 490
6bae477a 491 // Set run number in CDBManager (if it is not already set by the user)
492 if (!SetRunNumber()) if (fStopOnError) return kFALSE;
493
494 // Import ideal TGeo geometry and apply misalignment
495 if (!gGeoManager) {
496 TString geom(gSystem->DirName(fGAliceFileName));
497 geom += "/geometry.root";
498 TGeoManager::Import(geom.Data());
499 if (!gGeoManager) if (fStopOnError) return kFALSE;
500 }
501 if (!MisalignGeometry(fLoadAlignData)) if (fStopOnError) return kFALSE;
502
503 // Temporary fix by A.Gheata
504 // Could be removed with the next Root version (>5.11)
505 if (gGeoManager) {
506 TIter next(gGeoManager->GetListOfVolumes());
507 TGeoVolume *vol;
508 while ((vol = (TGeoVolume *)next())) {
509 if (vol->GetVoxels()) {
510 if (vol->GetVoxels()->NeedRebuild()) {
511 vol->GetVoxels()->Voxelize();
512 vol->FindOverlaps();
513 }
514 }
515 }
516 }
517
596a855f 518 // local reconstruction
59697224 519 if (!fRunLocalReconstruction.IsNull()) {
520 if (!RunLocalReconstruction(fRunLocalReconstruction)) {
e583c30d 521 if (fStopOnError) {CleanUp(); return kFALSE;}
596a855f 522 }
523 }
b26c3770 524// if (!fRunVertexFinder && fRunTracking.IsNull() &&
525// fFillESD.IsNull()) return kTRUE;
2257f27e 526
527 // get vertexer
528 if (fRunVertexFinder && !CreateVertexer()) {
529 if (fStopOnError) {
530 CleanUp();
531 return kFALSE;
532 }
533 }
596a855f 534
f08fc9f5 535 // get trackers
b8cd5251 536 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
24f7a148 537 if (fStopOnError) {
538 CleanUp();
539 return kFALSE;
540 }
596a855f 541 }
24f7a148 542
9db3a215 543
544 TStopwatch stopwatch;
545 stopwatch.Start();
546
b26c3770 547 // get the possibly already existing ESD file and tree
1f46a9ae 548 AliESD* esd = new AliESD; AliESD* hltesd = new AliESD;
b26c3770 549 TFile* fileOld = NULL;
1f46a9ae 550 TTree* treeOld = NULL; TTree *hlttreeOld = NULL;
b26c3770 551 if (!gSystem->AccessPathName("AliESDs.root")){
552 gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
553 fileOld = TFile::Open("AliESDs.old.root");
554 if (fileOld && fileOld->IsOpen()) {
555 treeOld = (TTree*) fileOld->Get("esdTree");
556 if (treeOld) treeOld->SetBranchAddress("ESD", &esd);
1f46a9ae 557 hlttreeOld = (TTree*) fileOld->Get("HLTesdTree");
558 if (hlttreeOld) hlttreeOld->SetBranchAddress("ESD", &hltesd);
b26c3770 559 }
560 }
561
36711aa4 562 // create the ESD output file and tree
596a855f 563 TFile* file = TFile::Open("AliESDs.root", "RECREATE");
564 if (!file->IsOpen()) {
815c2b38 565 AliError("opening AliESDs.root failed");
b26c3770 566 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
596a855f 567 }
36711aa4 568 TTree* tree = new TTree("esdTree", "Tree with ESD objects");
569 tree->Branch("ESD", "AliESD", &esd);
1f46a9ae 570 TTree* hlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
571 hlttree->Branch("ESD", "AliESD", &hltesd);
572 delete esd; delete hltesd;
573 esd = NULL; hltesd = NULL;
1d99986f 574
575 // create the file and tree with ESD additions
576 TFile *filef=0; TTree *treef=0; AliESDfriend *esdf=0;
577 if (fWriteESDfriend) {
578 filef = TFile::Open("AliESDfriends.root", "RECREATE");
579 if (!filef->IsOpen()) {
580 AliError("opening AliESDfriends.root failed");
581 }
582 treef = new TTree("esdFriendTree", "Tree with ESD friends");
583 treef->Branch("ESDfriend", "AliESDfriend", &esdf);
584 }
585
36711aa4 586 gROOT->cd();
596a855f 587
588 // loop over events
b649205a 589 if (fRawReader) fRawReader->RewindEvents();
f08fc9f5 590
596a855f 591 for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
b26c3770 592 if (fRawReader) fRawReader->NextEvent();
593 if ((iEvent < firstEvent) || ((lastEvent >= 0) && (iEvent > lastEvent))) {
594 // copy old ESD to the new one
595 if (treeOld) {
596 treeOld->SetBranchAddress("ESD", &esd);
597 treeOld->GetEntry(iEvent);
598 }
599 tree->Fill();
1f46a9ae 600 if (hlttreeOld) {
601 hlttreeOld->SetBranchAddress("ESD", &hltesd);
602 hlttreeOld->GetEntry(iEvent);
603 }
604 hlttree->Fill();
b26c3770 605 continue;
606 }
607
815c2b38 608 AliInfo(Form("processing event %d", iEvent));
596a855f 609 fRunLoader->GetEvent(iEvent);
24f7a148 610
611 char fileName[256];
612 sprintf(fileName, "ESD_%d.%d_final.root",
f08fc9f5 613 fRunLoader->GetHeader()->GetRun(),
614 fRunLoader->GetHeader()->GetEventNrInRun());
24f7a148 615 if (!gSystem->AccessPathName(fileName)) continue;
616
b26c3770 617 // local reconstruction
618 if (!fRunLocalReconstruction.IsNull()) {
619 if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
620 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
621 }
622 }
623
1f46a9ae 624 esd = new AliESD; hltesd = new AliESD;
f08fc9f5 625 esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1f46a9ae 626 hltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
f08fc9f5 627 esd->SetEventNumber(fRunLoader->GetHeader()->GetEventNrInRun());
1f46a9ae 628 hltesd->SetEventNumber(fRunLoader->GetHeader()->GetEventNrInRun());
f08fc9f5 629 if (gAlice) {
7ebb06c3 630 esd->SetMagneticField(AliTracker::GetBz());
631 hltesd->SetMagneticField(AliTracker::GetBz());
f08fc9f5 632 } else {
633 // ???
634 }
596a855f 635
2257f27e 636 // vertex finder
637 if (fRunVertexFinder) {
638 if (!ReadESD(esd, "vertex")) {
639 if (!RunVertexFinder(esd)) {
b26c3770 640 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
2257f27e 641 }
642 if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
643 }
644 }
645
1f46a9ae 646 // HLT tracking
647 if (!fRunTracking.IsNull()) {
648 if (fRunHLTTracking) {
649 hltesd->SetVertex(esd->GetVertex());
650 if (!RunHLTTracking(hltesd)) {
651 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
652 }
653 }
654 }
655
596a855f 656 // barrel tracking
b8cd5251 657 if (!fRunTracking.IsNull()) {
24f7a148 658 if (!ReadESD(esd, "tracking")) {
659 if (!RunTracking(esd)) {
b26c3770 660 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
24f7a148 661 }
662 if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
596a855f 663 }
664 }
665
666 // fill ESD
667 if (!fFillESD.IsNull()) {
668 if (!FillESD(esd, fFillESD)) {
b26c3770 669 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
596a855f 670 }
671 }
672
673 // combined PID
674 AliESDpid::MakePID(esd);
24f7a148 675 if (fCheckPointLevel > 1) WriteESD(esd, "PID");
596a855f 676
b647652d 677 if (fFillTriggerESD) {
678 if (!ReadESD(esd, "trigger")) {
679 if (!FillTriggerESD(esd)) {
680 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
681 }
682 if (fCheckPointLevel > 1) WriteESD(esd, "trigger");
683 }
684 }
685
596a855f 686 // write ESD
36711aa4 687 tree->Fill();
1f46a9ae 688 // write HLT ESD
689 hlttree->Fill();
24f7a148 690
1d99986f 691 // write ESD friend
692 if (fWriteESDfriend) {
d75007f6 693 esdf=new AliESDfriend();
694 esd->GetESDfriend(esdf);
1d99986f 695 treef->Fill();
696 }
697
f3a97c86 698 if (fCheckPointLevel > 0) WriteESD(esd, "final");
699
1f46a9ae 700 delete esd; delete hltesd;
701 esd = NULL; hltesd = NULL;
596a855f 702 }
703
9db3a215 704 AliInfo(Form("Execution time for filling ESD : R:%.2fs C:%.2fs",
705 stopwatch.RealTime(),stopwatch.CpuTime()));
706
36711aa4 707 file->cd();
708 tree->Write();
1f46a9ae 709 hlttree->Write();
f3a97c86 710
1d99986f 711 if (fWriteESDfriend) {
712 filef->cd();
713 treef->Write(); delete treef; filef->Close(); delete filef;
714 }
715
f3a97c86 716 // Create tags for the events in the ESD tree (the ESD tree is always present)
717 // In case of empty events the tags will contain dummy values
718 CreateTag(file);
b26c3770 719 CleanUp(file, fileOld);
596a855f 720
721 return kTRUE;
722}
723
724
725//_____________________________________________________________________________
59697224 726Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
596a855f 727{
59697224 728// run the local reconstruction
596a855f 729
030b532d 730 TStopwatch stopwatch;
731 stopwatch.Start();
732
596a855f 733 TString detStr = detectors;
b8cd5251 734 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
735 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
736 AliReconstructor* reconstructor = GetReconstructor(iDet);
737 if (!reconstructor) continue;
b26c3770 738 if (reconstructor->HasLocalReconstruction()) continue;
b8cd5251 739
740 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
741 TStopwatch stopwatchDet;
742 stopwatchDet.Start();
743 if (fRawReader) {
744 fRawReader->RewindEvents();
745 reconstructor->Reconstruct(fRunLoader, fRawReader);
746 } else {
747 reconstructor->Reconstruct(fRunLoader);
596a855f 748 }
5f8272e1 749 AliInfo(Form("Execution time for %s: R:%.2fs C:%.2fs",
750 fgkDetectorName[iDet],
751 stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
596a855f 752 }
753
754 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
815c2b38 755 AliError(Form("the following detectors were not found: %s",
756 detStr.Data()));
596a855f 757 if (fStopOnError) return kFALSE;
758 }
759
5f8272e1 760 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
761 stopwatch.RealTime(),stopwatch.CpuTime()));
030b532d 762
596a855f 763 return kTRUE;
764}
765
b26c3770 766//_____________________________________________________________________________
767Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
768{
769// run the local reconstruction
770
771 TStopwatch stopwatch;
772 stopwatch.Start();
773
774 TString detStr = detectors;
775 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
776 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
777 AliReconstructor* reconstructor = GetReconstructor(iDet);
778 if (!reconstructor) continue;
779 AliLoader* loader = fLoader[iDet];
780
781 // conversion of digits
782 if (fRawReader && reconstructor->HasDigitConversion()) {
783 AliInfo(Form("converting raw data digits into root objects for %s",
784 fgkDetectorName[iDet]));
785 TStopwatch stopwatchDet;
786 stopwatchDet.Start();
787 loader->LoadDigits("update");
788 loader->CleanDigits();
789 loader->MakeDigitsContainer();
790 TTree* digitsTree = loader->TreeD();
791 reconstructor->ConvertDigits(fRawReader, digitsTree);
792 loader->WriteDigits("OVERWRITE");
793 loader->UnloadDigits();
5f8272e1 794 AliInfo(Form("Execution time for %s: R:%.2fs C:%.2fs",
795 fgkDetectorName[iDet],
796 stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
b26c3770 797 }
798
799 // local reconstruction
800 if (!reconstructor->HasLocalReconstruction()) continue;
801 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
802 TStopwatch stopwatchDet;
803 stopwatchDet.Start();
804 loader->LoadRecPoints("update");
805 loader->CleanRecPoints();
806 loader->MakeRecPointsContainer();
807 TTree* clustersTree = loader->TreeR();
808 if (fRawReader && !reconstructor->HasDigitConversion()) {
809 reconstructor->Reconstruct(fRawReader, clustersTree);
810 } else {
811 loader->LoadDigits("read");
812 TTree* digitsTree = loader->TreeD();
813 if (!digitsTree) {
814 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
815 if (fStopOnError) return kFALSE;
816 } else {
817 reconstructor->Reconstruct(digitsTree, clustersTree);
818 }
819 loader->UnloadDigits();
820 }
821 loader->WriteRecPoints("OVERWRITE");
822 loader->UnloadRecPoints();
5f8272e1 823 AliDebug(1,Form("Execution time for %s: R:%.2fs C:%.2fs",
824 fgkDetectorName[iDet],
825 stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
b26c3770 826 }
827
828 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
829 AliError(Form("the following detectors were not found: %s",
830 detStr.Data()));
831 if (fStopOnError) return kFALSE;
832 }
5f8272e1 833
834 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
835 stopwatch.RealTime(),stopwatch.CpuTime()));
b26c3770 836
837 return kTRUE;
838}
839
596a855f 840//_____________________________________________________________________________
2257f27e 841Bool_t AliReconstruction::RunVertexFinder(AliESD*& esd)
596a855f 842{
843// run the barrel tracking
844
030b532d 845 TStopwatch stopwatch;
846 stopwatch.Start();
847
2257f27e 848 AliESDVertex* vertex = NULL;
849 Double_t vtxPos[3] = {0, 0, 0};
850 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
851 TArrayF mcVertex(3);
a6b0b91b 852 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
853 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
854 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
855 }
2257f27e 856
b8cd5251 857 if (fVertexer) {
815c2b38 858 AliInfo("running the ITS vertex finder");
b26c3770 859 if (fLoader[0]) fLoader[0]->LoadRecPoints();
b8cd5251 860 vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
b26c3770 861 if (fLoader[0]) fLoader[0]->UnloadRecPoints();
2257f27e 862 if(!vertex){
815c2b38 863 AliWarning("Vertex not found");
c710f220 864 vertex = new AliESDVertex();
d1a50cb5 865 vertex->SetName("default");
2257f27e 866 }
867 else {
868 vertex->SetTruePos(vtxPos); // store also the vertex from MC
d1a50cb5 869 vertex->SetName("reconstructed");
2257f27e 870 }
871
872 } else {
815c2b38 873 AliInfo("getting the primary vertex from MC");
2257f27e 874 vertex = new AliESDVertex(vtxPos, vtxErr);
875 }
876
877 if (vertex) {
878 vertex->GetXYZ(vtxPos);
879 vertex->GetSigmaXYZ(vtxErr);
880 } else {
815c2b38 881 AliWarning("no vertex reconstructed");
2257f27e 882 vertex = new AliESDVertex(vtxPos, vtxErr);
883 }
884 esd->SetVertex(vertex);
b8cd5251 885 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
886 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
887 }
2257f27e 888 delete vertex;
889
5f8272e1 890 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
891 stopwatch.RealTime(),stopwatch.CpuTime()));
2257f27e 892
893 return kTRUE;
894}
895
1f46a9ae 896//_____________________________________________________________________________
897Bool_t AliReconstruction::RunHLTTracking(AliESD*& esd)
898{
899// run the HLT barrel tracking
900
901 TStopwatch stopwatch;
902 stopwatch.Start();
903
904 if (!fRunLoader) {
905 AliError("Missing runLoader!");
906 return kFALSE;
907 }
908
909 AliInfo("running HLT tracking");
910
911 // Get a pointer to the HLT reconstructor
912 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
913 if (!reconstructor) return kFALSE;
914
915 // TPC + ITS
916 for (Int_t iDet = 1; iDet >= 0; iDet--) {
917 TString detName = fgkDetectorName[iDet];
918 AliDebug(1, Form("%s HLT tracking", detName.Data()));
919 reconstructor->SetOption(detName.Data());
920 AliTracker *tracker = reconstructor->CreateTracker(fRunLoader);
921 if (!tracker) {
922 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
923 if (fStopOnError) return kFALSE;
9dcc06e1 924 continue;
1f46a9ae 925 }
926 Double_t vtxPos[3];
927 Double_t vtxErr[3]={0.005,0.005,0.010};
928 const AliESDVertex *vertex = esd->GetVertex();
929 vertex->GetXYZ(vtxPos);
930 tracker->SetVertex(vtxPos,vtxErr);
931 if(iDet != 1) {
932 fLoader[iDet]->LoadRecPoints("read");
933 TTree* tree = fLoader[iDet]->TreeR();
934 if (!tree) {
935 AliError(Form("Can't get the %s cluster tree", detName.Data()));
936 return kFALSE;
937 }
938 tracker->LoadClusters(tree);
939 }
940 if (tracker->Clusters2Tracks(esd) != 0) {
941 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
942 return kFALSE;
943 }
944 if(iDet != 1) {
945 tracker->UnloadClusters();
946 }
947 delete tracker;
948 }
949
5f8272e1 950 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
951 stopwatch.RealTime(),stopwatch.CpuTime()));
1f46a9ae 952
953 return kTRUE;
954}
955
2257f27e 956//_____________________________________________________________________________
957Bool_t AliReconstruction::RunTracking(AliESD*& esd)
958{
959// run the barrel tracking
960
961 TStopwatch stopwatch;
962 stopwatch.Start();
24f7a148 963
815c2b38 964 AliInfo("running tracking");
596a855f 965
b8cd5251 966 // pass 1: TPC + ITS inwards
967 for (Int_t iDet = 1; iDet >= 0; iDet--) {
968 if (!fTracker[iDet]) continue;
969 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
24f7a148 970
b8cd5251 971 // load clusters
972 fLoader[iDet]->LoadRecPoints("read");
973 TTree* tree = fLoader[iDet]->TreeR();
974 if (!tree) {
975 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
24f7a148 976 return kFALSE;
977 }
b8cd5251 978 fTracker[iDet]->LoadClusters(tree);
979
980 // run tracking
981 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
982 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
24f7a148 983 return kFALSE;
984 }
b8cd5251 985 if (fCheckPointLevel > 1) {
986 WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
987 }
878e1fe1 988 // preliminary PID in TPC needed by the ITS tracker
989 if (iDet == 1) {
990 GetReconstructor(1)->FillESD(fRunLoader, esd);
b26c3770 991 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
878e1fe1 992 AliESDpid::MakePID(esd);
993 }
b8cd5251 994 }
596a855f 995
b8cd5251 996 // pass 2: ALL backwards
997 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
998 if (!fTracker[iDet]) continue;
999 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
1000
1001 // load clusters
1002 if (iDet > 1) { // all except ITS, TPC
1003 TTree* tree = NULL;
7b61cd9c 1004 fLoader[iDet]->LoadRecPoints("read");
1005 tree = fLoader[iDet]->TreeR();
b8cd5251 1006 if (!tree) {
1007 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
24f7a148 1008 return kFALSE;
1009 }
b8cd5251 1010 fTracker[iDet]->LoadClusters(tree);
1011 }
24f7a148 1012
b8cd5251 1013 // run tracking
1014 if (fTracker[iDet]->PropagateBack(esd) != 0) {
1015 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
1016 return kFALSE;
1017 }
1018 if (fCheckPointLevel > 1) {
1019 WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
1020 }
24f7a148 1021
b8cd5251 1022 // unload clusters
1023 if (iDet > 2) { // all except ITS, TPC, TRD
1024 fTracker[iDet]->UnloadClusters();
7b61cd9c 1025 fLoader[iDet]->UnloadRecPoints();
b8cd5251 1026 }
8f37df88 1027 // updated PID in TPC needed by the ITS tracker -MI
1028 if (iDet == 1) {
1029 GetReconstructor(1)->FillESD(fRunLoader, esd);
1030 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1031 AliESDpid::MakePID(esd);
1032 }
b8cd5251 1033 }
596a855f 1034
98937d93 1035 // write space-points to the ESD in case alignment data output
1036 // is switched on
1037 if (fWriteAlignmentData)
1038 WriteAlignmentData(esd);
1039
b8cd5251 1040 // pass 3: TRD + TPC + ITS refit inwards
1041 for (Int_t iDet = 2; iDet >= 0; iDet--) {
1042 if (!fTracker[iDet]) continue;
1043 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
596a855f 1044
b8cd5251 1045 // run tracking
1046 if (fTracker[iDet]->RefitInward(esd) != 0) {
1047 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
1048 return kFALSE;
1049 }
1050 if (fCheckPointLevel > 1) {
1051 WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
1052 }
596a855f 1053
b8cd5251 1054 // unload clusters
1055 fTracker[iDet]->UnloadClusters();
1056 fLoader[iDet]->UnloadRecPoints();
1057 }
ff8bb5ae 1058 //
1059 // Propagate track to the vertex - if not done by ITS
1060 //
1061 Int_t ntracks = esd->GetNumberOfTracks();
1062 for (Int_t itrack=0; itrack<ntracks; itrack++){
1063 const Double_t kRadius = 3; // beam pipe radius
1064 const Double_t kMaxStep = 5; // max step
1065 const Double_t kMaxD = 123456; // max distance to prim vertex
1066 Double_t fieldZ = AliTracker::GetBz(); //
1067 AliESDtrack * track = esd->GetTrack(itrack);
1068 if (!track) continue;
1069 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1070 track->PropagateTo(kRadius, track->GetMass(),kMaxStep,kTRUE);
1071 track->RelateToVertex(esd->GetVertex(),fieldZ, kMaxD);
1072 }
1073
5f8272e1 1074 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
1075 stopwatch.RealTime(),stopwatch.CpuTime()));
030b532d 1076
596a855f 1077 return kTRUE;
1078}
1079
1080//_____________________________________________________________________________
24f7a148 1081Bool_t AliReconstruction::FillESD(AliESD*& esd, const TString& detectors)
596a855f 1082{
1083// fill the event summary data
1084
030b532d 1085 TStopwatch stopwatch;
1086 stopwatch.Start();
815c2b38 1087 AliInfo("filling ESD");
030b532d 1088
596a855f 1089 TString detStr = detectors;
b8cd5251 1090 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1091 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1092 AliReconstructor* reconstructor = GetReconstructor(iDet);
1093 if (!reconstructor) continue;
1094
1095 if (!ReadESD(esd, fgkDetectorName[iDet])) {
1096 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
b26c3770 1097 TTree* clustersTree = NULL;
1098 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
1099 fLoader[iDet]->LoadRecPoints("read");
1100 clustersTree = fLoader[iDet]->TreeR();
1101 if (!clustersTree) {
1102 AliError(Form("Can't get the %s clusters tree",
1103 fgkDetectorName[iDet]));
1104 if (fStopOnError) return kFALSE;
1105 }
1106 }
1107 if (fRawReader && !reconstructor->HasDigitConversion()) {
1108 reconstructor->FillESD(fRawReader, clustersTree, esd);
1109 } else {
1110 TTree* digitsTree = NULL;
1111 if (fLoader[iDet]) {
1112 fLoader[iDet]->LoadDigits("read");
1113 digitsTree = fLoader[iDet]->TreeD();
1114 if (!digitsTree) {
1115 AliError(Form("Can't get the %s digits tree",
1116 fgkDetectorName[iDet]));
1117 if (fStopOnError) return kFALSE;
1118 }
1119 }
1120 reconstructor->FillESD(digitsTree, clustersTree, esd);
1121 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1122 }
1123 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
1124 fLoader[iDet]->UnloadRecPoints();
1125 }
1126
b8cd5251 1127 if (fRawReader) {
1128 reconstructor->FillESD(fRunLoader, fRawReader, esd);
1129 } else {
1130 reconstructor->FillESD(fRunLoader, esd);
24f7a148 1131 }
b8cd5251 1132 if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
596a855f 1133 }
1134 }
1135
1136 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
815c2b38 1137 AliError(Form("the following detectors were not found: %s",
1138 detStr.Data()));
596a855f 1139 if (fStopOnError) return kFALSE;
1140 }
1141
5f8272e1 1142 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
1143 stopwatch.RealTime(),stopwatch.CpuTime()));
030b532d 1144
596a855f 1145 return kTRUE;
1146}
1147
b647652d 1148//_____________________________________________________________________________
1149Bool_t AliReconstruction::FillTriggerESD(AliESD*& esd)
1150{
1151 // Reads the trigger decision which is
1152 // stored in Trigger.root file and fills
1153 // the corresponding esd entries
1154
1155 AliInfo("Filling trigger information into the ESD");
1156
1157 if (fRawReader) {
1158 AliCTPRawStream input(fRawReader);
1159 if (!input.Next()) {
1160 AliError("No valid CTP (trigger) DDL raw data is found ! The trigger information is not stored in the ESD !");
1161 return kFALSE;
1162 }
1163 esd->SetTriggerMask(input.GetClassMask());
1164 esd->SetTriggerCluster(input.GetClusterMask());
1165 }
1166 else {
1167 AliRunLoader *runloader = AliRunLoader::GetRunLoader();
1168 if (runloader) {
1169 if (!runloader->LoadTrigger()) {
1170 AliCentralTrigger *aCTP = runloader->GetTrigger();
1171 esd->SetTriggerMask(aCTP->GetClassMask());
1172 esd->SetTriggerCluster(aCTP->GetClusterMask());
1173 }
1174 else {
1175 AliWarning("No trigger can be loaded! The trigger information is not stored in the ESD !");
1176 return kFALSE;
1177 }
1178 }
1179 else {
1180 AliError("No run loader is available! The trigger information is not stored in the ESD !");
1181 return kFALSE;
1182 }
1183 }
1184
1185 return kTRUE;
1186}
596a855f 1187
1188//_____________________________________________________________________________
1189Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
1190{
1191// check whether detName is contained in detectors
1192// if yes, it is removed from detectors
1193
1194 // check if all detectors are selected
1195 if ((detectors.CompareTo("ALL") == 0) ||
1196 detectors.BeginsWith("ALL ") ||
1197 detectors.EndsWith(" ALL") ||
1198 detectors.Contains(" ALL ")) {
1199 detectors = "ALL";
1200 return kTRUE;
1201 }
1202
1203 // search for the given detector
1204 Bool_t result = kFALSE;
1205 if ((detectors.CompareTo(detName) == 0) ||
1206 detectors.BeginsWith(detName+" ") ||
1207 detectors.EndsWith(" "+detName) ||
1208 detectors.Contains(" "+detName+" ")) {
1209 detectors.ReplaceAll(detName, "");
1210 result = kTRUE;
1211 }
1212
1213 // clean up the detectors string
1214 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1215 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1216 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1217
1218 return result;
1219}
e583c30d 1220
f08fc9f5 1221//_____________________________________________________________________________
1222Bool_t AliReconstruction::InitRunLoader()
1223{
1224// get or create the run loader
1225
1226 if (gAlice) delete gAlice;
1227 gAlice = NULL;
1228
b26c3770 1229 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
1230 // load all base libraries to get the loader classes
1231 TString libs = gSystem->GetLibraries();
1232 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1233 TString detName = fgkDetectorName[iDet];
1234 if (detName == "HLT") continue;
1235 if (libs.Contains("lib" + detName + "base.so")) continue;
1236 gSystem->Load("lib" + detName + "base.so");
1237 }
f08fc9f5 1238 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
1239 if (!fRunLoader) {
1240 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1241 CleanUp();
1242 return kFALSE;
1243 }
b26c3770 1244 fRunLoader->CdGAFile();
1245 if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
1246 if (fRunLoader->LoadgAlice() == 0) {
1247 gAlice = fRunLoader->GetAliRun();
c84a5e9e 1248 AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
b26c3770 1249 }
f08fc9f5 1250 }
1251 if (!gAlice && !fRawReader) {
1252 AliError(Form("no gAlice object found in file %s",
1253 fGAliceFileName.Data()));
1254 CleanUp();
1255 return kFALSE;
1256 }
1257
1258 } else { // galice.root does not exist
1259 if (!fRawReader) {
1260 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
1261 CleanUp();
1262 return kFALSE;
1263 }
1264 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
1265 AliConfig::GetDefaultEventFolderName(),
1266 "recreate");
1267 if (!fRunLoader) {
1268 AliError(Form("could not create run loader in file %s",
1269 fGAliceFileName.Data()));
1270 CleanUp();
1271 return kFALSE;
1272 }
1273 fRunLoader->MakeTree("E");
1274 Int_t iEvent = 0;
1275 while (fRawReader->NextEvent()) {
1276 fRunLoader->SetEventNumber(iEvent);
1277 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1278 iEvent, iEvent);
1279 fRunLoader->MakeTree("H");
1280 fRunLoader->TreeE()->Fill();
1281 iEvent++;
1282 }
1283 fRawReader->RewindEvents();
1284 fRunLoader->WriteHeader("OVERWRITE");
1285 fRunLoader->CdGAFile();
1286 fRunLoader->Write(0, TObject::kOverwrite);
1287// AliTracker::SetFieldMap(???);
1288 }
1289
1290 return kTRUE;
1291}
1292
c757bafd 1293//_____________________________________________________________________________
b8cd5251 1294AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
c757bafd 1295{
f08fc9f5 1296// get the reconstructor object and the loader for a detector
c757bafd 1297
b8cd5251 1298 if (fReconstructor[iDet]) return fReconstructor[iDet];
1299
1300 // load the reconstructor object
1301 TPluginManager* pluginManager = gROOT->GetPluginManager();
1302 TString detName = fgkDetectorName[iDet];
1303 TString recName = "Ali" + detName + "Reconstructor";
f08fc9f5 1304 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
b8cd5251 1305
1306 if (detName == "HLT") {
1307 if (!gROOT->GetClass("AliLevel3")) {
1308 gSystem->Load("libAliL3Src.so");
1309 gSystem->Load("libAliL3Misc.so");
1310 gSystem->Load("libAliL3Hough.so");
1311 gSystem->Load("libAliL3Comp.so");
1312 }
1313 }
1314
1315 AliReconstructor* reconstructor = NULL;
1316 // first check if a plugin is defined for the reconstructor
1317 TPluginHandler* pluginHandler =
1318 pluginManager->FindHandler("AliReconstructor", detName);
f08fc9f5 1319 // if not, add a plugin for it
1320 if (!pluginHandler) {
b8cd5251 1321 AliDebug(1, Form("defining plugin for %s", recName.Data()));
b26c3770 1322 TString libs = gSystem->GetLibraries();
1323 if (libs.Contains("lib" + detName + "base.so") ||
1324 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
b8cd5251 1325 pluginManager->AddHandler("AliReconstructor", detName,
1326 recName, detName + "rec", recName + "()");
1327 } else {
1328 pluginManager->AddHandler("AliReconstructor", detName,
1329 recName, detName, recName + "()");
c757bafd 1330 }
b8cd5251 1331 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
1332 }
1333 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1334 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
c757bafd 1335 }
b8cd5251 1336 if (reconstructor) {
1337 TObject* obj = fOptions.FindObject(detName.Data());
1338 if (obj) reconstructor->SetOption(obj->GetTitle());
b26c3770 1339 reconstructor->Init(fRunLoader);
b8cd5251 1340 fReconstructor[iDet] = reconstructor;
1341 }
1342
f08fc9f5 1343 // get or create the loader
1344 if (detName != "HLT") {
1345 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
1346 if (!fLoader[iDet]) {
1347 AliConfig::Instance()
1348 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
1349 detName, detName);
1350 // first check if a plugin is defined for the loader
1351 TPluginHandler* pluginHandler =
1352 pluginManager->FindHandler("AliLoader", detName);
1353 // if not, add a plugin for it
1354 if (!pluginHandler) {
1355 TString loaderName = "Ali" + detName + "Loader";
1356 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
1357 pluginManager->AddHandler("AliLoader", detName,
1358 loaderName, detName + "base",
1359 loaderName + "(const char*, TFolder*)");
1360 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
1361 }
1362 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1363 fLoader[iDet] =
1364 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
1365 fRunLoader->GetEventFolder());
1366 }
1367 if (!fLoader[iDet]) { // use default loader
1368 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
1369 }
1370 if (!fLoader[iDet]) {
1371 AliWarning(Form("couldn't get loader for %s", detName.Data()));
6667b602 1372 if (fStopOnError) return NULL;
f08fc9f5 1373 } else {
1374 fRunLoader->AddLoader(fLoader[iDet]);
1375 fRunLoader->CdGAFile();
1376 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
1377 fRunLoader->Write(0, TObject::kOverwrite);
1378 }
1379 }
1380 }
1381
b8cd5251 1382 return reconstructor;
c757bafd 1383}
1384
2257f27e 1385//_____________________________________________________________________________
1386Bool_t AliReconstruction::CreateVertexer()
1387{
1388// create the vertexer
1389
b8cd5251 1390 fVertexer = NULL;
1391 AliReconstructor* itsReconstructor = GetReconstructor(0);
59697224 1392 if (itsReconstructor) {
b8cd5251 1393 fVertexer = itsReconstructor->CreateVertexer(fRunLoader);
2257f27e 1394 }
b8cd5251 1395 if (!fVertexer) {
815c2b38 1396 AliWarning("couldn't create a vertexer for ITS");
2257f27e 1397 if (fStopOnError) return kFALSE;
1398 }
1399
1400 return kTRUE;
1401}
1402
24f7a148 1403//_____________________________________________________________________________
b8cd5251 1404Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
24f7a148 1405{
f08fc9f5 1406// create the trackers
24f7a148 1407
b8cd5251 1408 TString detStr = detectors;
1409 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1410 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1411 AliReconstructor* reconstructor = GetReconstructor(iDet);
1412 if (!reconstructor) continue;
1413 TString detName = fgkDetectorName[iDet];
1f46a9ae 1414 if (detName == "HLT") {
1415 fRunHLTTracking = kTRUE;
1416 continue;
1417 }
f08fc9f5 1418
1419 fTracker[iDet] = reconstructor->CreateTracker(fRunLoader);
1420 if (!fTracker[iDet] && (iDet < 7)) {
1421 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
8250d5f5 1422 if (fStopOnError) return kFALSE;
1423 }
1424 }
1425
24f7a148 1426 return kTRUE;
1427}
1428
e583c30d 1429//_____________________________________________________________________________
b26c3770 1430void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
e583c30d 1431{
1432// delete trackers and the run loader and close and delete the file
1433
b8cd5251 1434 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1435 delete fReconstructor[iDet];
1436 fReconstructor[iDet] = NULL;
1437 fLoader[iDet] = NULL;
1438 delete fTracker[iDet];
1439 fTracker[iDet] = NULL;
1440 }
1441 delete fVertexer;
1442 fVertexer = NULL;
e583c30d 1443
1444 delete fRunLoader;
1445 fRunLoader = NULL;
b649205a 1446 delete fRawReader;
1447 fRawReader = NULL;
e583c30d 1448
1449 if (file) {
1450 file->Close();
1451 delete file;
1452 }
b26c3770 1453
1454 if (fileOld) {
1455 fileOld->Close();
1456 delete fileOld;
1457 gSystem->Unlink("AliESDs.old.root");
1458 }
e583c30d 1459}
24f7a148 1460
1461
1462//_____________________________________________________________________________
1463Bool_t AliReconstruction::ReadESD(AliESD*& esd, const char* recStep) const
1464{
1465// read the ESD event from a file
1466
1467 if (!esd) return kFALSE;
1468 char fileName[256];
1469 sprintf(fileName, "ESD_%d.%d_%s.root",
1470 esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1471 if (gSystem->AccessPathName(fileName)) return kFALSE;
1472
f3a97c86 1473 AliInfo(Form("reading ESD from file %s", fileName));
815c2b38 1474 AliDebug(1, Form("reading ESD from file %s", fileName));
24f7a148 1475 TFile* file = TFile::Open(fileName);
1476 if (!file || !file->IsOpen()) {
815c2b38 1477 AliError(Form("opening %s failed", fileName));
24f7a148 1478 delete file;
1479 return kFALSE;
1480 }
1481
1482 gROOT->cd();
1483 delete esd;
1484 esd = (AliESD*) file->Get("ESD");
1485 file->Close();
1486 delete file;
1487 return kTRUE;
1488}
1489
1490//_____________________________________________________________________________
1491void AliReconstruction::WriteESD(AliESD* esd, const char* recStep) const
1492{
1493// write the ESD event to a file
1494
1495 if (!esd) return;
1496 char fileName[256];
1497 sprintf(fileName, "ESD_%d.%d_%s.root",
1498 esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1499
815c2b38 1500 AliDebug(1, Form("writing ESD to file %s", fileName));
24f7a148 1501 TFile* file = TFile::Open(fileName, "recreate");
1502 if (!file || !file->IsOpen()) {
815c2b38 1503 AliError(Form("opening %s failed", fileName));
24f7a148 1504 } else {
1505 esd->Write("ESD");
1506 file->Close();
1507 }
1508 delete file;
1509}
f3a97c86 1510
1511
1512
1513
1514//_____________________________________________________________________________
1515void AliReconstruction::CreateTag(TFile* file)
1516{
2bdb9d38 1517 /////////////
1518 //muon code//
1519 ////////////
56982dd1 1520 Double_t fMUONMASS = 0.105658369;
2bdb9d38 1521 //Variables
56982dd1 1522 Double_t fX,fY,fZ ;
1523 Double_t fThetaX, fThetaY, fPyz, fChisquare;
1524 Double_t fPxRec,fPyRec, fPzRec, fEnergy;
1525 Int_t fCharge;
1526 TLorentzVector fEPvector;
1527
1528 Float_t fZVertexCut = 10.0;
1529 Float_t fRhoVertexCut = 2.0;
1530
1531 Float_t fLowPtCut = 1.0;
1532 Float_t fHighPtCut = 3.0;
1533 Float_t fVeryHighPtCut = 10.0;
2bdb9d38 1534 ////////////
1535
1536 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
1537
089bf903 1538 // Creates the tags for all the events in a given ESD file
4302e20f 1539 Int_t ntrack;
089bf903 1540 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
1541 Int_t nPos, nNeg, nNeutr;
1542 Int_t nK0s, nNeutrons, nPi0s, nGamas;
1543 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
1544 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
1545 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
1546 Float_t maxPt = .0, meanPt = .0, totalP = .0;
d1a50cb5 1547 Int_t fVertexflag;
1548 TString fVertexName;
4302e20f 1549
f3a97c86 1550 AliRunTag *tag = new AliRunTag();
f3a97c86 1551 AliEventTag *evTag = new AliEventTag();
1552 TTree ttag("T","A Tree with event tags");
2bdb9d38 1553 TBranch * btag = ttag.Branch("AliTAG", &tag);
f3a97c86 1554 btag->SetCompressionLevel(9);
2bdb9d38 1555
f3a97c86 1556 AliInfo(Form("Creating the tags......."));
1557
1558 if (!file || !file->IsOpen()) {
1559 AliError(Form("opening failed"));
1560 delete file;
1561 return ;
2bdb9d38 1562 }
1563 Int_t firstEvent = 0,lastEvent = 0;
f3a97c86 1564 TTree *t = (TTree*) file->Get("esdTree");
1565 TBranch * b = t->GetBranch("ESD");
1566 AliESD *esd = 0;
1567 b->SetAddress(&esd);
2bdb9d38 1568
f3a97c86 1569 tag->SetRunId(esd->GetRunNumber());
2bdb9d38 1570
089bf903 1571 Int_t iNumberOfEvents = b->GetEntries();
2bdb9d38 1572 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
1573 ntrack = 0;
1574 nPos = 0;
1575 nNeg = 0;
1576 nNeutr =0;
1577 nK0s = 0;
1578 nNeutrons = 0;
1579 nPi0s = 0;
1580 nGamas = 0;
1581 nProtons = 0;
1582 nKaons = 0;
1583 nPions = 0;
1584 nMuons = 0;
1585 nElectrons = 0;
1586 nCh1GeV = 0;
1587 nCh3GeV = 0;
1588 nCh10GeV = 0;
1589 nMu1GeV = 0;
1590 nMu3GeV = 0;
1591 nMu10GeV = 0;
1592 nEl1GeV = 0;
1593 nEl3GeV = 0;
1594 nEl10GeV = 0;
1595 maxPt = .0;
1596 meanPt = .0;
1597 totalP = .0;
d1a50cb5 1598 fVertexflag = 1;
1599
2bdb9d38 1600 b->GetEntry(iEventNumber);
1601 const AliESDVertex * vertexIn = esd->GetVertex();
d1a50cb5 1602 fVertexName = vertexIn->GetName();
1603 if(fVertexName == "default") fVertexflag = 0;
2bdb9d38 1604 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
1605 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
1606 UInt_t status = esdTrack->GetStatus();
f3a97c86 1607
2bdb9d38 1608 //select only tracks with ITS refit
1609 if ((status&AliESDtrack::kITSrefit)==0) continue;
1610 //select only tracks with TPC refit
1611 if ((status&AliESDtrack::kTPCrefit)==0) continue;
4302e20f 1612
2bdb9d38 1613 //select only tracks with the "combined PID"
1614 if ((status&AliESDtrack::kESDpid)==0) continue;
1615 Double_t p[3];
1616 esdTrack->GetPxPyPz(p);
1617 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
1618 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
1619 totalP += momentum;
1620 meanPt += fPt;
1621 if(fPt > maxPt) maxPt = fPt;
4302e20f 1622
2bdb9d38 1623 if(esdTrack->GetSign() > 0) {
1624 nPos++;
56982dd1 1625 if(fPt > fLowPtCut) nCh1GeV++;
1626 if(fPt > fHighPtCut) nCh3GeV++;
1627 if(fPt > fVeryHighPtCut) nCh10GeV++;
2bdb9d38 1628 }
1629 if(esdTrack->GetSign() < 0) {
1630 nNeg++;
56982dd1 1631 if(fPt > fLowPtCut) nCh1GeV++;
1632 if(fPt > fHighPtCut) nCh3GeV++;
1633 if(fPt > fVeryHighPtCut) nCh10GeV++;
2bdb9d38 1634 }
1635 if(esdTrack->GetSign() == 0) nNeutr++;
4302e20f 1636
2bdb9d38 1637 //PID
1638 Double_t prob[5];
1639 esdTrack->GetESDpid(prob);
4302e20f 1640
2bdb9d38 1641 Double_t rcc = 0.0;
1642 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
1643 if(rcc == 0.0) continue;
1644 //Bayes' formula
1645 Double_t w[5];
1646 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
4302e20f 1647
2bdb9d38 1648 //protons
1649 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
1650 //kaons
1651 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
1652 //pions
1653 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
1654 //electrons
1655 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
1656 nElectrons++;
56982dd1 1657 if(fPt > fLowPtCut) nEl1GeV++;
1658 if(fPt > fHighPtCut) nEl3GeV++;
1659 if(fPt > fVeryHighPtCut) nEl10GeV++;
2bdb9d38 1660 }
1661 ntrack++;
1662 }//track loop
1663
1664 /////////////
1665 //muon code//
1666 ////////////
1667 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
1668 // loop over all reconstructed tracks (also first track of combination)
1669 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
1670 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
1671 if (muonTrack == 0x0) continue;
4302e20f 1672
2bdb9d38 1673 // Coordinates at vertex
56982dd1 1674 fZ = muonTrack->GetZ();
1675 fY = muonTrack->GetBendingCoor();
1676 fX = muonTrack->GetNonBendingCoor();
4302e20f 1677
56982dd1 1678 fThetaX = muonTrack->GetThetaX();
1679 fThetaY = muonTrack->GetThetaY();
4302e20f 1680
56982dd1 1681 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
1682 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
1683 fPxRec = fPzRec * TMath::Tan(fThetaX);
1684 fPyRec = fPzRec * TMath::Tan(fThetaY);
1685 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
f3a97c86 1686
2bdb9d38 1687 //ChiSquare of the track if needed
56982dd1 1688 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
1689 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
1690 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
4302e20f 1691
2bdb9d38 1692 // total number of muons inside a vertex cut
56982dd1 1693 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
2bdb9d38 1694 nMuons++;
56982dd1 1695 if(fEPvector.Pt() > fLowPtCut) {
2bdb9d38 1696 nMu1GeV++;
56982dd1 1697 if(fEPvector.Pt() > fHighPtCut) {
2bdb9d38 1698 nMu3GeV++;
56982dd1 1699 if (fEPvector.Pt() > fVeryHighPtCut) {
2bdb9d38 1700 nMu10GeV++;
1701 }
1702 }
1703 }
1704 }
1705 }//muon track loop
1706
1707 // Fill the event tags
1708 if(ntrack != 0)
1709 meanPt = meanPt/ntrack;
1710
1711 evTag->SetEventId(iEventNumber+1);
1712 evTag->SetVertexX(vertexIn->GetXv());
1713 evTag->SetVertexY(vertexIn->GetYv());
1714 evTag->SetVertexZ(vertexIn->GetZv());
d1a50cb5 1715 evTag->SetVertexZError(vertexIn->GetZRes());
1716 evTag->SetVertexFlag(fVertexflag);
1717
2bdb9d38 1718 evTag->SetT0VertexZ(esd->GetT0zVertex());
1719
b647652d 1720 evTag->SetTrigger(esd->GetTriggerMask());
2bdb9d38 1721
32a5cab4 1722 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
1723 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
1724 evTag->SetZDCNeutron2Energy(esd->GetZDCN2Energy());
1725 evTag->SetZDCProton2Energy(esd->GetZDCP2Energy());
2bdb9d38 1726 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
1727 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
1728
1729
1730 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
1731 evTag->SetNumOfPosTracks(nPos);
1732 evTag->SetNumOfNegTracks(nNeg);
1733 evTag->SetNumOfNeutrTracks(nNeutr);
1734
1735 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
1736 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
1737 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
1738 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
1739
1740 evTag->SetNumOfProtons(nProtons);
1741 evTag->SetNumOfKaons(nKaons);
1742 evTag->SetNumOfPions(nPions);
1743 evTag->SetNumOfMuons(nMuons);
1744 evTag->SetNumOfElectrons(nElectrons);
1745 evTag->SetNumOfPhotons(nGamas);
1746 evTag->SetNumOfPi0s(nPi0s);
1747 evTag->SetNumOfNeutrons(nNeutrons);
1748 evTag->SetNumOfKaon0s(nK0s);
1749
1750 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
1751 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
1752 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
1753 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
1754 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
1755 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
1756 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
1757 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
1758 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
1759
85c60a8e 1760 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
1761 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
2bdb9d38 1762
1763 evTag->SetTotalMomentum(totalP);
1764 evTag->SetMeanPt(meanPt);
1765 evTag->SetMaxPt(maxPt);
1766
1767 tag->AddEventTag(*evTag);
1768 }
089bf903 1769 lastEvent = iNumberOfEvents;
f3a97c86 1770
1771 ttag.Fill();
1772 tag->Clear();
1773
1774 char fileName[256];
1775 sprintf(fileName, "Run%d.Event%d_%d.ESD.tag.root",
1776 tag->GetRunId(),firstEvent,lastEvent );
1777 AliInfo(Form("writing tags to file %s", fileName));
1778 AliDebug(1, Form("writing tags to file %s", fileName));
1779
1780 TFile* ftag = TFile::Open(fileName, "recreate");
1781 ftag->cd();
1782 ttag.Write();
1783 ftag->Close();
1784 file->cd();
1785 delete tag;
f3a97c86 1786 delete evTag;
1787}
1788
98937d93 1789void AliReconstruction::WriteAlignmentData(AliESD* esd)
1790{
1791 // Write space-points which are then used in the alignment procedures
1792 // For the moment only ITS, TRD and TPC
1793
1794 // Load TOF clusters
d528ee75 1795 if (fTracker[3]){
1796 fLoader[3]->LoadRecPoints("read");
1797 TTree* tree = fLoader[3]->TreeR();
1798 if (!tree) {
1799 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
1800 return;
1801 }
1802 fTracker[3]->LoadClusters(tree);
98937d93 1803 }
98937d93 1804 Int_t ntracks = esd->GetNumberOfTracks();
1805 for (Int_t itrack = 0; itrack < ntracks; itrack++)
1806 {
1807 AliESDtrack *track = esd->GetTrack(itrack);
1808 Int_t nsp = 0;
ef7253ac 1809 Int_t idx[200];
98937d93 1810 for (Int_t iDet = 3; iDet >= 0; iDet--)
1811 nsp += track->GetNcls(iDet);
1812 if (nsp) {
1813 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
1814 track->SetTrackPointArray(sp);
1815 Int_t isptrack = 0;
1816 for (Int_t iDet = 3; iDet >= 0; iDet--) {
1817 AliTracker *tracker = fTracker[iDet];
1818 if (!tracker) continue;
1819 Int_t nspdet = track->GetNcls(iDet);
98937d93 1820 if (nspdet <= 0) continue;
1821 track->GetClusters(iDet,idx);
1822 AliTrackPoint p;
1823 Int_t isp = 0;
1824 Int_t isp2 = 0;
1825 while (isp < nspdet) {
1826 Bool_t isvalid = tracker->GetTrackPoint(idx[isp2],p); isp2++;
160db090 1827 const Int_t kNTPCmax = 159;
1828 if (iDet==1 && isp2>kNTPCmax) break; // to be fixed
98937d93 1829 if (!isvalid) continue;
1830 sp->AddPoint(isptrack,&p); isptrack++; isp++;
1831 }
98937d93 1832 }
1833 }
1834 }
d528ee75 1835 if (fTracker[3]){
1836 fTracker[3]->UnloadClusters();
1837 fLoader[3]->UnloadRecPoints();
1838 }
98937d93 1839}