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