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