Corrected version (Jochen Klein).
[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// //
973388c2 50// In case of raw-data reconstruction the user can modify the default //
51// number of events per digits/clusters/tracks file. In case the option //
52// is not used the number is set 1. In case the user provides 0, than //
53// the number of events is equal to the number of events inside the //
54// raw-data file (i.e. one digits/clusters/tracks file): //
55// //
56// rec.SetNumberOfEventsPerFile(...); //
57// //
58// //
596a855f 59// The name of the galice file can be changed from the default //
e583c30d 60// "galice.root" by passing it as argument to the AliReconstruction //
61// constructor or by //
596a855f 62// //
63// rec.SetGAliceFile("..."); //
64// //
59697224 65// The local reconstruction can be switched on or off for individual //
66// detectors by //
596a855f 67// //
59697224 68// rec.SetRunLocalReconstruction("..."); //
596a855f 69// //
70// The argument is a (case sensitive) string with the names of the //
71// detectors separated by a space. The special string "ALL" selects all //
72// available detectors. This is the default. //
73// //
c71de921 74// The reconstruction of the primary vertex position can be switched off by //
75// //
76// rec.SetRunVertexFinder(kFALSE); //
77// //
b8cd5251 78// The tracking and the creation of ESD tracks can be switched on for //
79// selected detectors by //
596a855f 80// //
b8cd5251 81// rec.SetRunTracking("..."); //
596a855f 82// //
c84a5e9e 83// Uniform/nonuniform field tracking switches (default: uniform field) //
84// //
1d99986f 85// rec.SetUniformFieldTracking(); ( rec.SetUniformFieldTracking(kFALSE); ) //
c84a5e9e 86// //
596a855f 87// The filling of additional ESD information can be steered by //
88// //
89// rec.SetFillESD("..."); //
90// //
b8cd5251 91// Again, for both methods the string specifies the list of detectors. //
92// The default is "ALL". //
93// //
94// The call of the shortcut method //
95// //
96// rec.SetRunReconstruction("..."); //
97// //
98// is equivalent to calling SetRunLocalReconstruction, SetRunTracking and //
99// SetFillESD with the same detector selecting string as argument. //
596a855f 100// //
c71de921 101// The reconstruction requires digits or raw data as input. For the creation //
102// of digits and raw data have a look at the class AliSimulation. //
596a855f 103// //
24f7a148 104// For debug purposes the method SetCheckPointLevel can be used. If the //
105// argument is greater than 0, files with ESD events will be written after //
106// selected steps of the reconstruction for each event: //
107// level 1: after tracking and after filling of ESD (final) //
108// level 2: in addition after each tracking step //
109// level 3: in addition after the filling of ESD for each detector //
110// If a final check point file exists for an event, this event will be //
111// skipped in the reconstruction. The tracking and the filling of ESD for //
112// a detector will be skipped as well, if the corresponding check point //
113// file exists. The ESD event will then be loaded from the file instead. //
114// //
596a855f 115///////////////////////////////////////////////////////////////////////////////
116
024a7e64 117#include <TArrayF.h>
118#include <TFile.h>
119#include <TSystem.h>
120#include <TROOT.h>
121#include <TPluginManager.h>
3103d196 122#include <TGeoManager.h>
2bdb9d38 123#include <TLorentzVector.h>
596a855f 124
125#include "AliReconstruction.h"
87932dab 126#include "AliCodeTimer.h"
b8cd5251 127#include "AliReconstructor.h"
815c2b38 128#include "AliLog.h"
596a855f 129#include "AliRunLoader.h"
130#include "AliRun.h"
b649205a 131#include "AliRawReaderFile.h"
132#include "AliRawReaderDate.h"
133#include "AliRawReaderRoot.h"
001397cd 134#include "AliRawEventHeaderBase.h"
af885e0f 135#include "AliESDEvent.h"
faffd83e 136#include "AliESDMuonTrack.h"
1d99986f 137#include "AliESDfriend.h"
2257f27e 138#include "AliESDVertex.h"
faffd83e 139#include "AliESDcascade.h"
140#include "AliESDkink.h"
141#include "AliESDtrack.h"
142#include "AliESDCaloCluster.h"
32e449be 143#include "AliMultiplicity.h"
c84a5e9e 144#include "AliTracker.h"
2257f27e 145#include "AliVertexer.h"
c5e3e5d1 146#include "AliVertexerTracks.h"
5e4ff34d 147#include "AliV0vertexer.h"
148#include "AliCascadeVertexer.h"
596a855f 149#include "AliHeader.h"
150#include "AliGenEventHeader.h"
b26c3770 151#include "AliPID.h"
596a855f 152#include "AliESDpid.h"
ff8bb5ae 153#include "AliESDtrack.h"
f3a97c86 154
08e1a23e 155#include "AliESDTagCreator.h"
a1069ee1 156#include "AliAODTagCreator.h"
f3a97c86 157
25be1e5c 158#include "AliGeomManager.h"
98937d93 159#include "AliTrackPointArray.h"
b0314964 160#include "AliCDBManager.h"
6bae477a 161#include "AliCDBEntry.h"
162#include "AliAlignObj.h"
f3a97c86 163
b647652d 164#include "AliCentralTrigger.h"
165#include "AliCTPRawStream.h"
166
f29f1726 167#include "AliAODEvent.h"
168#include "AliAODHeader.h"
169#include "AliAODTrack.h"
170#include "AliAODVertex.h"
171#include "AliAODCluster.h"
172
173
596a855f 174ClassImp(AliReconstruction)
175
176
177//_____________________________________________________________________________
b384f8a4 178const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
c757bafd 179
180//_____________________________________________________________________________
024cf675 181AliReconstruction::AliReconstruction(const char* gAliceFilename, const char* cdbUri,
e583c30d 182 const char* name, const char* title) :
183 TNamed(name, title),
184
c84a5e9e 185 fUniformField(kTRUE),
2257f27e 186 fRunVertexFinder(kTRUE),
1f46a9ae 187 fRunHLTTracking(kFALSE),
e66fbafb 188 fRunMuonTracking(kFALSE),
1d99986f 189 fStopOnError(kFALSE),
190 fWriteAlignmentData(kFALSE),
d64bd07d 191 fCleanESD(kTRUE),
1d99986f 192 fWriteESDfriend(kFALSE),
a7807689 193 fWriteAOD(kFALSE),
b647652d 194 fFillTriggerESD(kTRUE),
1d99986f 195
196 fRunLocalReconstruction("ALL"),
b8cd5251 197 fRunTracking("ALL"),
e583c30d 198 fFillESD("ALL"),
199 fGAliceFileName(gAliceFilename),
b649205a 200 fInput(""),
35042093 201 fEquipIdMap(""),
b26c3770 202 fFirstEvent(0),
203 fLastEvent(-1),
973388c2 204 fNumberOfEventsPerFile(1),
24f7a148 205 fCheckPointLevel(0),
b8cd5251 206 fOptions(),
6bae477a 207 fLoadAlignFromCDB(kTRUE),
208 fLoadAlignData("ALL"),
46698ae4 209 fESDPar(""),
e583c30d 210
211 fRunLoader(NULL),
b649205a 212 fRawReader(NULL),
b8cd5251 213
98937d93 214 fVertexer(NULL),
9178838a 215 fDiamondProfile(NULL),
98937d93 216
6bae477a 217 fAlignObjArray(NULL),
ec92bee0 218 fCDBUri(cdbUri),
219 fSpecCDBUri()
596a855f 220{
221// create reconstruction object with default parameters
b8cd5251 222
223 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
224 fReconstructor[iDet] = NULL;
225 fLoader[iDet] = NULL;
226 fTracker[iDet] = NULL;
227 }
e47c4c2e 228 AliPID pid;
596a855f 229}
230
231//_____________________________________________________________________________
232AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
e583c30d 233 TNamed(rec),
234
c84a5e9e 235 fUniformField(rec.fUniformField),
2257f27e 236 fRunVertexFinder(rec.fRunVertexFinder),
1f46a9ae 237 fRunHLTTracking(rec.fRunHLTTracking),
e66fbafb 238 fRunMuonTracking(rec.fRunMuonTracking),
1d99986f 239 fStopOnError(rec.fStopOnError),
240 fWriteAlignmentData(rec.fWriteAlignmentData),
d64bd07d 241 fCleanESD(rec.fCleanESD),
1d99986f 242 fWriteESDfriend(rec.fWriteESDfriend),
a7807689 243 fWriteAOD(rec.fWriteAOD),
b647652d 244 fFillTriggerESD(rec.fFillTriggerESD),
1d99986f 245
246 fRunLocalReconstruction(rec.fRunLocalReconstruction),
e583c30d 247 fRunTracking(rec.fRunTracking),
248 fFillESD(rec.fFillESD),
249 fGAliceFileName(rec.fGAliceFileName),
b649205a 250 fInput(rec.fInput),
35042093 251 fEquipIdMap(rec.fEquipIdMap),
b26c3770 252 fFirstEvent(rec.fFirstEvent),
253 fLastEvent(rec.fLastEvent),
973388c2 254 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
24f7a148 255 fCheckPointLevel(0),
b8cd5251 256 fOptions(),
6bae477a 257 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
258 fLoadAlignData(rec.fLoadAlignData),
46698ae4 259 fESDPar(rec.fESDPar),
e583c30d 260
261 fRunLoader(NULL),
b649205a 262 fRawReader(NULL),
b8cd5251 263
98937d93 264 fVertexer(NULL),
9178838a 265 fDiamondProfile(NULL),
98937d93 266
6bae477a 267 fAlignObjArray(rec.fAlignObjArray),
ec92bee0 268 fCDBUri(rec.fCDBUri),
269 fSpecCDBUri()
596a855f 270{
271// copy constructor
272
ec92bee0 273 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
efd2085e 274 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
275 }
b8cd5251 276 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
277 fReconstructor[iDet] = NULL;
278 fLoader[iDet] = NULL;
279 fTracker[iDet] = NULL;
280 }
ec92bee0 281 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
282 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
283 }
596a855f 284}
285
286//_____________________________________________________________________________
287AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
288{
289// assignment operator
290
291 this->~AliReconstruction();
292 new(this) AliReconstruction(rec);
293 return *this;
294}
295
296//_____________________________________________________________________________
297AliReconstruction::~AliReconstruction()
298{
299// clean up
300
e583c30d 301 CleanUp();
efd2085e 302 fOptions.Delete();
ec92bee0 303 fSpecCDBUri.Delete();
87932dab 304
305 AliCodeTimer::Instance()->Print();
596a855f 306}
307
024cf675 308//_____________________________________________________________________________
309void AliReconstruction::InitCDBStorage()
310{
311// activate a default CDB storage
312// First check if we have any CDB storage set, because it is used
313// to retrieve the calibration and alignment constants
314
315 AliCDBManager* man = AliCDBManager::Instance();
ec92bee0 316 if (man->IsDefaultStorageSet())
024cf675 317 {
318 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
ec92bee0 319 AliWarning("Default CDB storage has been already set !");
320 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
024cf675 321 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
ec92bee0 322 fCDBUri = "";
323 }
324 else {
b8ec52f6 325 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
326 AliDebug(2, Form("Default CDB storage is set to: %s",fCDBUri.Data()));
327 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
ec92bee0 328 man->SetDefaultStorage(fCDBUri);
329 }
330
331 // Now activate the detector specific CDB storage locations
c3a7b59a 332 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
333 TObject* obj = fSpecCDBUri[i];
334 if (!obj) continue;
b8ec52f6 335 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
336 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
337 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
c3a7b59a 338 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
ec92bee0 339 }
727d0766 340 man->Print();
024cf675 341}
342
343//_____________________________________________________________________________
344void AliReconstruction::SetDefaultStorage(const char* uri) {
ec92bee0 345// Store the desired default CDB storage location
346// Activate it later within the Run() method
024cf675 347
ec92bee0 348 fCDBUri = uri;
024cf675 349
350}
351
352//_____________________________________________________________________________
c3a7b59a 353void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
ec92bee0 354// Store a detector-specific CDB storage location
355// Activate it later within the Run() method
024cf675 356
c3a7b59a 357 AliCDBPath aPath(calibType);
358 if(!aPath.IsValid()){
359 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
360 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
361 if(!strcmp(calibType, fgkDetectorName[iDet])) {
362 aPath.SetPath(Form("%s/*", calibType));
363 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
364 break;
365 }
366 }
367 if(!aPath.IsValid()){
368 AliError(Form("Not a valid path or detector: %s", calibType));
369 return;
370 }
371 }
372
373 // check that calibType refers to a "valid" detector name
374 Bool_t isDetector = kFALSE;
375 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
376 TString detName = fgkDetectorName[iDet];
377 if(aPath.GetLevel0() == detName) {
378 isDetector = kTRUE;
379 break;
380 }
381 }
382
383 if(!isDetector) {
384 AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
385 return;
386 }
387
388 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
ec92bee0 389 if (obj) fSpecCDBUri.Remove(obj);
c3a7b59a 390 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
024cf675 391
392}
393
46698ae4 394
395
396
6bae477a 397//_____________________________________________________________________________
398Bool_t AliReconstruction::SetRunNumber()
399{
400 // The method is called in Run() in order
401 // to set a correct run number.
402 // In case of raw data reconstruction the
403 // run number is taken from the raw data header
404
405 if(AliCDBManager::Instance()->GetRun() < 0) {
406 if (!fRunLoader) {
407 AliError("No run loader is found !");
408 return kFALSE;
409 }
410 // read run number from gAlice
ec92bee0 411 if(fRunLoader->GetAliRun())
412 AliCDBManager::Instance()->SetRun(fRunLoader->GetAliRun()->GetRunNumber());
413 else {
414 if(fRawReader) {
415 if(fRawReader->NextEvent()) {
416 AliCDBManager::Instance()->SetRun(fRawReader->GetRunNumber());
417 fRawReader->RewindEvents();
418 }
419 else {
420 AliError("No raw-data events found !");
421 return kFALSE;
422 }
423 }
424 else {
425 AliError("Neither gAlice nor RawReader objects are found !");
426 return kFALSE;
427 }
428 }
429 AliInfo(Form("CDB Run number: %d",AliCDBManager::Instance()->GetRun()));
6bae477a 430 }
431 return kTRUE;
432}
433
434//_____________________________________________________________________________
6bae477a 435Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
436{
437 // Read the alignment objects from CDB.
438 // Each detector is supposed to have the
439 // alignment objects in DET/Align/Data CDB path.
440 // All the detector objects are then collected,
441 // sorted by geometry level (starting from ALIC) and
442 // then applied to the TGeo geometry.
443 // Finally an overlaps check is performed.
444
445 // Load alignment data from CDB and fill fAlignObjArray
446 if(fLoadAlignFromCDB){
6bae477a 447
25be1e5c 448 TString detStr = detectors;
449 TString loadAlObjsListOfDets = "";
450
451 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
452 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
453 loadAlObjsListOfDets += fgkDetectorName[iDet];
454 loadAlObjsListOfDets += " ";
455 } // end loop over detectors
98e303d9 456 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
25be1e5c 457 }else{
458 // Check if the array with alignment objects was
459 // provided by the user. If yes, apply the objects
460 // to the present TGeo geometry
461 if (fAlignObjArray) {
462 if (gGeoManager && gGeoManager->IsClosed()) {
98e303d9 463 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
25be1e5c 464 AliError("The misalignment of one or more volumes failed!"
465 "Compare the list of simulated detectors and the list of detector alignment data!");
466 return kFALSE;
467 }
468 }
469 else {
470 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
6bae477a 471 return kFALSE;
472 }
473 }
6bae477a 474 }
25be1e5c 475
8e245d15 476 delete fAlignObjArray; fAlignObjArray=0;
a03b0371 477
6bae477a 478 return kTRUE;
479}
596a855f 480
481//_____________________________________________________________________________
482void AliReconstruction::SetGAliceFile(const char* fileName)
483{
484// set the name of the galice file
485
486 fGAliceFileName = fileName;
487}
488
efd2085e 489//_____________________________________________________________________________
490void AliReconstruction::SetOption(const char* detector, const char* option)
491{
492// set options for the reconstruction of a detector
493
494 TObject* obj = fOptions.FindObject(detector);
495 if (obj) fOptions.Remove(obj);
496 fOptions.Add(new TNamed(detector, option));
497}
498
596a855f 499
500//_____________________________________________________________________________
4a33489c 501Bool_t AliReconstruction::Run(const char* input)
596a855f 502{
503// run the reconstruction
504
87932dab 505 AliCodeTimerAuto("")
506
b649205a 507 // set the input
508 if (!input) input = fInput.Data();
509 TString fileName(input);
510 if (fileName.EndsWith("/")) {
511 fRawReader = new AliRawReaderFile(fileName);
512 } else if (fileName.EndsWith(".root")) {
513 fRawReader = new AliRawReaderRoot(fileName);
514 } else if (!fileName.IsNull()) {
515 fRawReader = new AliRawReaderDate(fileName);
516 fRawReader->SelectEvents(7);
517 }
35042093 518 if (!fEquipIdMap.IsNull() && fRawReader)
519 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
520
b649205a 521
f08fc9f5 522 // get the run loader
523 if (!InitRunLoader()) return kFALSE;
596a855f 524
ec92bee0 525 // Initialize the CDB storage
526 InitCDBStorage();
527
6bae477a 528 // Set run number in CDBManager (if it is not already set by the user)
529 if (!SetRunNumber()) if (fStopOnError) return kFALSE;
530
531 // Import ideal TGeo geometry and apply misalignment
532 if (!gGeoManager) {
533 TString geom(gSystem->DirName(fGAliceFileName));
534 geom += "/geometry.root";
98e303d9 535 AliGeomManager::LoadGeometry(geom.Data());
6bae477a 536 if (!gGeoManager) if (fStopOnError) return kFALSE;
537 }
8e245d15 538
539 AliCDBManager* man = AliCDBManager::Instance();
6bae477a 540 if (!MisalignGeometry(fLoadAlignData)) if (fStopOnError) return kFALSE;
541
596a855f 542 // local reconstruction
59697224 543 if (!fRunLocalReconstruction.IsNull()) {
544 if (!RunLocalReconstruction(fRunLocalReconstruction)) {
e583c30d 545 if (fStopOnError) {CleanUp(); return kFALSE;}
596a855f 546 }
547 }
b26c3770 548// if (!fRunVertexFinder && fRunTracking.IsNull() &&
549// fFillESD.IsNull()) return kTRUE;
2257f27e 550
551 // get vertexer
552 if (fRunVertexFinder && !CreateVertexer()) {
553 if (fStopOnError) {
554 CleanUp();
555 return kFALSE;
556 }
557 }
596a855f 558
f08fc9f5 559 // get trackers
b8cd5251 560 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
24f7a148 561 if (fStopOnError) {
562 CleanUp();
563 return kFALSE;
564 }
596a855f 565 }
24f7a148 566
b26c3770 567 // get the possibly already existing ESD file and tree
af885e0f 568 AliESDEvent* esd = new AliESDEvent(); AliESDEvent* hltesd = new AliESDEvent();
b26c3770 569 TFile* fileOld = NULL;
1f46a9ae 570 TTree* treeOld = NULL; TTree *hlttreeOld = NULL;
b26c3770 571 if (!gSystem->AccessPathName("AliESDs.root")){
572 gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
573 fileOld = TFile::Open("AliESDs.old.root");
574 if (fileOld && fileOld->IsOpen()) {
575 treeOld = (TTree*) fileOld->Get("esdTree");
46698ae4 576 if (treeOld)esd->ReadFromTree(treeOld);
1f46a9ae 577 hlttreeOld = (TTree*) fileOld->Get("HLTesdTree");
46698ae4 578 if (hlttreeOld) hltesd->ReadFromTree(hlttreeOld);
b26c3770 579 }
580 }
581
36711aa4 582 // create the ESD output file and tree
596a855f 583 TFile* file = TFile::Open("AliESDs.root", "RECREATE");
46698ae4 584 file->SetCompressionLevel(2);
596a855f 585 if (!file->IsOpen()) {
815c2b38 586 AliError("opening AliESDs.root failed");
b26c3770 587 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
596a855f 588 }
46698ae4 589
36711aa4 590 TTree* tree = new TTree("esdTree", "Tree with ESD objects");
af885e0f 591 esd = new AliESDEvent();
46698ae4 592 esd->CreateStdContent();
593 esd->WriteToTree(tree);
594
1f46a9ae 595 TTree* hlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
af885e0f 596 hltesd = new AliESDEvent();
46698ae4 597 hltesd->CreateStdContent();
598 hltesd->WriteToTree(hlttree);
599
600 /* CKB Why?
1f46a9ae 601 delete esd; delete hltesd;
602 esd = NULL; hltesd = NULL;
46698ae4 603 */
500d54ab 604 // create the branch with ESD additions
5728d3d5 605
606
607
46698ae4 608 AliESDfriend *esdf = 0;
1d99986f 609 if (fWriteESDfriend) {
46698ae4 610 esdf = new AliESDfriend();
611 TBranch *br=tree->Branch("ESDfriend.","AliESDfriend", &esdf);
612 br->SetFile("AliESDfriends.root");
613 esd->AddObject(esdf);
1d99986f 614 }
5728d3d5 615
46698ae4 616
17c86e90 617 // Get the diamond profile from OCDB
618 AliCDBEntry* entry = AliCDBManager::Instance()
619 ->Get("GRP/Calib/MeanVertex");
620
621 if(entry) {
622 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
623 } else {
624 AliError("No diamond profile found in OCDB!");
625 }
626
20e5681c 627 AliVertexerTracks tVertexer(AliTracker::GetBz());
9178838a 628 if(fDiamondProfile) tVertexer.SetVtxStart(fDiamondProfile);
c5e3e5d1 629
596a855f 630 // loop over events
b649205a 631 if (fRawReader) fRawReader->RewindEvents();
f08fc9f5 632
596a855f 633 for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
b26c3770 634 if (fRawReader) fRawReader->NextEvent();
4a33489c 635 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
b26c3770 636 // copy old ESD to the new one
637 if (treeOld) {
46698ae4 638 esd->ReadFromTree(treeOld);
b26c3770 639 treeOld->GetEntry(iEvent);
640 }
641 tree->Fill();
1f46a9ae 642 if (hlttreeOld) {
46698ae4 643 esd->ReadFromTree(hlttreeOld);
1f46a9ae 644 hlttreeOld->GetEntry(iEvent);
645 }
646 hlttree->Fill();
b26c3770 647 continue;
648 }
46698ae4 649
815c2b38 650 AliInfo(Form("processing event %d", iEvent));
596a855f 651 fRunLoader->GetEvent(iEvent);
24f7a148 652
bb0901a4 653 char aFileName[256];
654 sprintf(aFileName, "ESD_%d.%d_final.root",
f08fc9f5 655 fRunLoader->GetHeader()->GetRun(),
656 fRunLoader->GetHeader()->GetEventNrInRun());
bb0901a4 657 if (!gSystem->AccessPathName(aFileName)) continue;
24f7a148 658
b26c3770 659 // local reconstruction
660 if (!fRunLocalReconstruction.IsNull()) {
661 if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
662 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
663 }
664 }
665
46698ae4 666
f08fc9f5 667 esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1f46a9ae 668 hltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
31fd97b2 669 esd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
670 hltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
46698ae4 671
d6ee376f 672 // Set magnetic field from the tracker
673 esd->SetMagneticField(AliTracker::GetBz());
674 hltesd->SetMagneticField(AliTracker::GetBz());
596a855f 675
46698ae4 676
677
2e3550da 678 // Fill raw-data error log into the ESD
679 if (fRawReader) FillRawDataErrorLog(iEvent,esd);
680
2257f27e 681 // vertex finder
682 if (fRunVertexFinder) {
683 if (!ReadESD(esd, "vertex")) {
684 if (!RunVertexFinder(esd)) {
b26c3770 685 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
2257f27e 686 }
687 if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
688 }
689 }
690
1f46a9ae 691 // HLT tracking
692 if (!fRunTracking.IsNull()) {
693 if (fRunHLTTracking) {
694 hltesd->SetVertex(esd->GetVertex());
695 if (!RunHLTTracking(hltesd)) {
696 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
697 }
698 }
699 }
700
e66fbafb 701 // Muon tracking
b8cd5251 702 if (!fRunTracking.IsNull()) {
e66fbafb 703 if (fRunMuonTracking) {
761350a6 704 if (!RunMuonTracking(esd)) {
b26c3770 705 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
24f7a148 706 }
596a855f 707 }
708 }
709
e66fbafb 710 // barrel tracking
711 if (!fRunTracking.IsNull()) {
21c573b7 712 if (!ReadESD(esd, "tracking")) {
713 if (!RunTracking(esd)) {
714 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
e66fbafb 715 }
21c573b7 716 if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
e66fbafb 717 }
718 }
21c573b7 719
596a855f 720 // fill ESD
721 if (!fFillESD.IsNull()) {
722 if (!FillESD(esd, fFillESD)) {
b26c3770 723 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
596a855f 724 }
725 }
001397cd 726 // fill Event header information from the RawEventHeader
727 if (fRawReader){FillRawEventHeaderESD(esd);}
596a855f 728
729 // combined PID
730 AliESDpid::MakePID(esd);
24f7a148 731 if (fCheckPointLevel > 1) WriteESD(esd, "PID");
596a855f 732
b647652d 733 if (fFillTriggerESD) {
734 if (!ReadESD(esd, "trigger")) {
735 if (!FillTriggerESD(esd)) {
736 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
737 }
738 if (fCheckPointLevel > 1) WriteESD(esd, "trigger");
739 }
740 }
741
a03cd2e0 742 //Try to improve the reconstructed primary vertex position using the tracks
c060d7fe 743 AliESDVertex *pvtx=0;
744 Bool_t dovertex=kTRUE;
745 TObject* obj = fOptions.FindObject("ITS");
746 if (obj) {
747 TString optITS = obj->GetTitle();
748 if (optITS.Contains("cosmics") || optITS.Contains("COSMICS"))
749 dovertex=kFALSE;
750 }
751 if(dovertex) pvtx=tVertexer.FindPrimaryVertex(esd);
17c86e90 752 if(fDiamondProfile) esd->SetDiamond(fDiamondProfile);
753
a03cd2e0 754 if (pvtx)
755 if (pvtx->GetStatus()) {
756 // Store the improved primary vertex
757 esd->SetPrimaryVertex(pvtx);
758 // Propagate the tracks to the DCA to the improved primary vertex
759 Double_t somethingbig = 777.;
760 Double_t bz = esd->GetMagneticField();
761 Int_t nt=esd->GetNumberOfTracks();
762 while (nt--) {
763 AliESDtrack *t = esd->GetTrack(nt);
764 t->RelateToVertex(pvtx, bz, somethingbig);
765 }
766 }
c5e3e5d1 767
5e4ff34d 768 {
769 // V0 finding
770 AliV0vertexer vtxer;
771 vtxer.Tracks2V0vertices(esd);
772
773 // Cascade finding
774 AliCascadeVertexer cvtxer;
775 cvtxer.V0sTracks2CascadeVertices(esd);
776 }
777
596a855f 778 // write ESD
d64bd07d 779 if (fCleanESD) CleanESD(esd);
1d99986f 780 if (fWriteESDfriend) {
46698ae4 781 new (esdf) AliESDfriend(); // Reset...
782 esd->GetESDfriend(esdf);
1d99986f 783 }
500d54ab 784 tree->Fill();
785
786 // write HLT ESD
787 hlttree->Fill();
1d99986f 788
f3a97c86 789 if (fCheckPointLevel > 0) WriteESD(esd, "final");
46698ae4 790 esd->Reset();
791 hltesd->Reset();
5728d3d5 792 if (fWriteESDfriend) {
793 new (esdf) AliESDfriend(); // Reset...
794 }
46698ae4 795 // esdf->Reset();
796 // delete esdf; esdf = 0;
797 }
798
799
5728d3d5 800
801
46698ae4 802 tree->GetUserInfo()->Add(esd);
803 hlttree->GetUserInfo()->Add(hltesd);
804
805
806
807 if(fESDPar.Contains("ESD.par")){
808 AliInfo("Attaching ESD.par to Tree");
809 TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
810 tree->GetUserInfo()->Add(fn);
596a855f 811 }
812
46698ae4 813
36711aa4 814 file->cd();
a9c0e6db 815 if (fWriteESDfriend)
816 tree->SetBranchStatus("ESDfriend*",0);
36711aa4 817 tree->Write();
1f46a9ae 818 hlttree->Write();
f3a97c86 819
a7807689 820 if (fWriteAOD) {
f29f1726 821 TFile *aodFile = TFile::Open("AliAOD.root", "RECREATE");
822 ESDFile2AODFile(file, aodFile);
823 aodFile->Close();
a7807689 824 }
825
46698ae4 826 gROOT->cd();
151e0c96 827 CleanUp(file, fileOld);
828
f3a97c86 829 // Create tags for the events in the ESD tree (the ESD tree is always present)
830 // In case of empty events the tags will contain dummy values
08e1a23e 831 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
832 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent);
a1069ee1 833 if (fWriteAOD) {
834 AliAODTagCreator *aodtagCreator = new AliAODTagCreator();
835 aodtagCreator->CreateAODTags(fFirstEvent,fLastEvent);
836 }
596a855f 837
838 return kTRUE;
839}
840
841
842//_____________________________________________________________________________
59697224 843Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
596a855f 844{
59697224 845// run the local reconstruction
596a855f 846
87932dab 847 AliCodeTimerAuto("")
030b532d 848
8e245d15 849 AliCDBManager* man = AliCDBManager::Instance();
850 Bool_t origCache = man->GetCacheFlag();
851
596a855f 852 TString detStr = detectors;
b8cd5251 853 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
854 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
855 AliReconstructor* reconstructor = GetReconstructor(iDet);
856 if (!reconstructor) continue;
b26c3770 857 if (reconstructor->HasLocalReconstruction()) continue;
b8cd5251 858
87932dab 859 AliCodeTimerStart(Form("running reconstruction for %s", fgkDetectorName[iDet]));
b8cd5251 860 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
87932dab 861
862 AliCodeTimerStart(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
8e245d15 863 AliInfo(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
864
865 man->SetCacheFlag(kTRUE);
866 TString calibPath = Form("%s/Calib/*", fgkDetectorName[iDet]);
867 man->GetAll(calibPath); // entries are cached!
868
87932dab 869 AliCodeTimerStop(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
870
b8cd5251 871 if (fRawReader) {
872 fRawReader->RewindEvents();
873 reconstructor->Reconstruct(fRunLoader, fRawReader);
874 } else {
875 reconstructor->Reconstruct(fRunLoader);
596a855f 876 }
87932dab 877
878 AliCodeTimerStop(Form("running reconstruction for %s", fgkDetectorName[iDet]));
8e245d15 879
880 // unload calibration data
98e303d9 881 man->UnloadFromCache(calibPath);
882 //man->ClearCache();
596a855f 883 }
884
8e245d15 885 man->SetCacheFlag(origCache);
886
596a855f 887 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
815c2b38 888 AliError(Form("the following detectors were not found: %s",
889 detStr.Data()));
596a855f 890 if (fStopOnError) return kFALSE;
891 }
892
893 return kTRUE;
894}
895
896//_____________________________________________________________________________
b26c3770 897Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
898{
899// run the local reconstruction
900
87932dab 901 AliCodeTimerAuto("")
b26c3770 902
903 TString detStr = detectors;
904 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
905 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
906 AliReconstructor* reconstructor = GetReconstructor(iDet);
907 if (!reconstructor) continue;
908 AliLoader* loader = fLoader[iDet];
909
910 // conversion of digits
911 if (fRawReader && reconstructor->HasDigitConversion()) {
912 AliInfo(Form("converting raw data digits into root objects for %s",
913 fgkDetectorName[iDet]));
87932dab 914 AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
915 fgkDetectorName[iDet]));
b26c3770 916 loader->LoadDigits("update");
917 loader->CleanDigits();
918 loader->MakeDigitsContainer();
919 TTree* digitsTree = loader->TreeD();
920 reconstructor->ConvertDigits(fRawReader, digitsTree);
921 loader->WriteDigits("OVERWRITE");
922 loader->UnloadDigits();
b26c3770 923 }
924
925 // local reconstruction
926 if (!reconstructor->HasLocalReconstruction()) continue;
927 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
87932dab 928 AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
b26c3770 929 loader->LoadRecPoints("update");
930 loader->CleanRecPoints();
931 loader->MakeRecPointsContainer();
932 TTree* clustersTree = loader->TreeR();
933 if (fRawReader && !reconstructor->HasDigitConversion()) {
934 reconstructor->Reconstruct(fRawReader, clustersTree);
935 } else {
936 loader->LoadDigits("read");
937 TTree* digitsTree = loader->TreeD();
938 if (!digitsTree) {
939 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
940 if (fStopOnError) return kFALSE;
941 } else {
942 reconstructor->Reconstruct(digitsTree, clustersTree);
943 }
944 loader->UnloadDigits();
945 }
946 loader->WriteRecPoints("OVERWRITE");
947 loader->UnloadRecPoints();
b26c3770 948 }
949
950 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
951 AliError(Form("the following detectors were not found: %s",
952 detStr.Data()));
953 if (fStopOnError) return kFALSE;
954 }
5f8272e1 955
b26c3770 956 return kTRUE;
957}
958
959//_____________________________________________________________________________
af885e0f 960Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
596a855f 961{
962// run the barrel tracking
963
87932dab 964 AliCodeTimerAuto("")
030b532d 965
2257f27e 966 AliESDVertex* vertex = NULL;
967 Double_t vtxPos[3] = {0, 0, 0};
968 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
969 TArrayF mcVertex(3);
a6b0b91b 970 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
971 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
972 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
973 }
2257f27e 974
b8cd5251 975 if (fVertexer) {
17c86e90 976 if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
815c2b38 977 AliInfo("running the ITS vertex finder");
b26c3770 978 if (fLoader[0]) fLoader[0]->LoadRecPoints();
b8cd5251 979 vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
b26c3770 980 if (fLoader[0]) fLoader[0]->UnloadRecPoints();
2257f27e 981 if(!vertex){
815c2b38 982 AliWarning("Vertex not found");
c710f220 983 vertex = new AliESDVertex();
d1a50cb5 984 vertex->SetName("default");
2257f27e 985 }
986 else {
d1a50cb5 987 vertex->SetName("reconstructed");
2257f27e 988 }
989
990 } else {
815c2b38 991 AliInfo("getting the primary vertex from MC");
2257f27e 992 vertex = new AliESDVertex(vtxPos, vtxErr);
993 }
994
995 if (vertex) {
996 vertex->GetXYZ(vtxPos);
997 vertex->GetSigmaXYZ(vtxErr);
998 } else {
815c2b38 999 AliWarning("no vertex reconstructed");
2257f27e 1000 vertex = new AliESDVertex(vtxPos, vtxErr);
1001 }
1002 esd->SetVertex(vertex);
32e449be 1003 // if SPD multiplicity has been determined, it is stored in the ESD
25be1e5c 1004 AliMultiplicity *mult = fVertexer->GetMultiplicity();
32e449be 1005 if(mult)esd->SetMultiplicity(mult);
1006
b8cd5251 1007 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1008 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1009 }
2257f27e 1010 delete vertex;
1011
2257f27e 1012 return kTRUE;
1013}
1014
1015//_____________________________________________________________________________
af885e0f 1016Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1f46a9ae 1017{
1018// run the HLT barrel tracking
1019
87932dab 1020 AliCodeTimerAuto("")
1f46a9ae 1021
1022 if (!fRunLoader) {
1023 AliError("Missing runLoader!");
1024 return kFALSE;
1025 }
1026
1027 AliInfo("running HLT tracking");
1028
1029 // Get a pointer to the HLT reconstructor
1030 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1031 if (!reconstructor) return kFALSE;
1032
1033 // TPC + ITS
1034 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1035 TString detName = fgkDetectorName[iDet];
1036 AliDebug(1, Form("%s HLT tracking", detName.Data()));
1037 reconstructor->SetOption(detName.Data());
1038 AliTracker *tracker = reconstructor->CreateTracker(fRunLoader);
1039 if (!tracker) {
1040 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1041 if (fStopOnError) return kFALSE;
9dcc06e1 1042 continue;
1f46a9ae 1043 }
1044 Double_t vtxPos[3];
1045 Double_t vtxErr[3]={0.005,0.005,0.010};
1046 const AliESDVertex *vertex = esd->GetVertex();
1047 vertex->GetXYZ(vtxPos);
1048 tracker->SetVertex(vtxPos,vtxErr);
1049 if(iDet != 1) {
1050 fLoader[iDet]->LoadRecPoints("read");
1051 TTree* tree = fLoader[iDet]->TreeR();
1052 if (!tree) {
1053 AliError(Form("Can't get the %s cluster tree", detName.Data()));
1054 return kFALSE;
1055 }
1056 tracker->LoadClusters(tree);
1057 }
1058 if (tracker->Clusters2Tracks(esd) != 0) {
1059 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1060 return kFALSE;
1061 }
1062 if(iDet != 1) {
1063 tracker->UnloadClusters();
1064 }
1065 delete tracker;
1066 }
1067
1f46a9ae 1068 return kTRUE;
1069}
1070
1071//_____________________________________________________________________________
af885e0f 1072Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
e66fbafb 1073{
1074// run the muon spectrometer tracking
1075
87932dab 1076 AliCodeTimerAuto("")
e66fbafb 1077
1078 if (!fRunLoader) {
1079 AliError("Missing runLoader!");
1080 return kFALSE;
1081 }
1082 Int_t iDet = 7; // for MUON
1083
1084 AliInfo("is running...");
1085
1086 // Get a pointer to the MUON reconstructor
1087 AliReconstructor *reconstructor = GetReconstructor(iDet);
1088 if (!reconstructor) return kFALSE;
1089
1090
1091 TString detName = fgkDetectorName[iDet];
1092 AliDebug(1, Form("%s tracking", detName.Data()));
1093 AliTracker *tracker = reconstructor->CreateTracker(fRunLoader);
1094 if (!tracker) {
1095 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1096 return kFALSE;
1097 }
1098
1099 // create Tracks
1100 fLoader[iDet]->LoadTracks("update");
1101 fLoader[iDet]->CleanTracks();
1102 fLoader[iDet]->MakeTracksContainer();
1103
1104 // read RecPoints
761350a6 1105 fLoader[iDet]->LoadRecPoints("read");
1106 tracker->LoadClusters(fLoader[iDet]->TreeR());
1107
1108 Int_t rv = tracker->Clusters2Tracks(esd);
1109
1110 fLoader[iDet]->UnloadRecPoints();
1111
1112 if ( rv )
1113 {
e66fbafb 1114 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1115 return kFALSE;
1116 }
761350a6 1117
1118 tracker->UnloadClusters();
1119
e66fbafb 1120 fLoader[iDet]->UnloadRecPoints();
1121
1122 fLoader[iDet]->WriteTracks("OVERWRITE");
1123 fLoader[iDet]->UnloadTracks();
1124
1125 delete tracker;
1126
e66fbafb 1127 return kTRUE;
1128}
1129
1130
1131//_____________________________________________________________________________
af885e0f 1132Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
2257f27e 1133{
1134// run the barrel tracking
1135
87932dab 1136 AliCodeTimerAuto("")
24f7a148 1137
815c2b38 1138 AliInfo("running tracking");
596a855f 1139
91b876d1 1140 //Fill the ESD with the T0 info (will be used by the TOF)
1f4331b3 1141 if (fReconstructor[11])
1142 GetReconstructor(11)->FillESD(fRunLoader, esd);
91b876d1 1143
b8cd5251 1144 // pass 1: TPC + ITS inwards
1145 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1146 if (!fTracker[iDet]) continue;
1147 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
24f7a148 1148
b8cd5251 1149 // load clusters
1150 fLoader[iDet]->LoadRecPoints("read");
1151 TTree* tree = fLoader[iDet]->TreeR();
1152 if (!tree) {
1153 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
24f7a148 1154 return kFALSE;
1155 }
b8cd5251 1156 fTracker[iDet]->LoadClusters(tree);
1157
1158 // run tracking
1159 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
1160 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
24f7a148 1161 return kFALSE;
1162 }
b8cd5251 1163 if (fCheckPointLevel > 1) {
1164 WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
1165 }
878e1fe1 1166 // preliminary PID in TPC needed by the ITS tracker
1167 if (iDet == 1) {
1168 GetReconstructor(1)->FillESD(fRunLoader, esd);
b26c3770 1169 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
878e1fe1 1170 AliESDpid::MakePID(esd);
1171 }
b8cd5251 1172 }
596a855f 1173
b8cd5251 1174 // pass 2: ALL backwards
1175 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1176 if (!fTracker[iDet]) continue;
1177 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
1178
1179 // load clusters
1180 if (iDet > 1) { // all except ITS, TPC
1181 TTree* tree = NULL;
7b61cd9c 1182 fLoader[iDet]->LoadRecPoints("read");
1183 tree = fLoader[iDet]->TreeR();
b8cd5251 1184 if (!tree) {
1185 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
24f7a148 1186 return kFALSE;
1187 }
b8cd5251 1188 fTracker[iDet]->LoadClusters(tree);
1189 }
24f7a148 1190
b8cd5251 1191 // run tracking
1192 if (fTracker[iDet]->PropagateBack(esd) != 0) {
1193 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
49dfd67a 1194 // return kFALSE;
b8cd5251 1195 }
1196 if (fCheckPointLevel > 1) {
1197 WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
1198 }
24f7a148 1199
b8cd5251 1200 // unload clusters
1201 if (iDet > 2) { // all except ITS, TPC, TRD
1202 fTracker[iDet]->UnloadClusters();
7b61cd9c 1203 fLoader[iDet]->UnloadRecPoints();
b8cd5251 1204 }
8f37df88 1205 // updated PID in TPC needed by the ITS tracker -MI
1206 if (iDet == 1) {
1207 GetReconstructor(1)->FillESD(fRunLoader, esd);
1208 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1209 AliESDpid::MakePID(esd);
1210 }
b8cd5251 1211 }
596a855f 1212
98937d93 1213 // write space-points to the ESD in case alignment data output
1214 // is switched on
1215 if (fWriteAlignmentData)
1216 WriteAlignmentData(esd);
1217
b8cd5251 1218 // pass 3: TRD + TPC + ITS refit inwards
1219 for (Int_t iDet = 2; iDet >= 0; iDet--) {
1220 if (!fTracker[iDet]) continue;
1221 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
596a855f 1222
b8cd5251 1223 // run tracking
1224 if (fTracker[iDet]->RefitInward(esd) != 0) {
1225 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
49dfd67a 1226 // return kFALSE;
b8cd5251 1227 }
1228 if (fCheckPointLevel > 1) {
1229 WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
1230 }
596a855f 1231
b8cd5251 1232 // unload clusters
1233 fTracker[iDet]->UnloadClusters();
1234 fLoader[iDet]->UnloadRecPoints();
1235 }
ff8bb5ae 1236 //
1237 // Propagate track to the vertex - if not done by ITS
1238 //
1239 Int_t ntracks = esd->GetNumberOfTracks();
1240 for (Int_t itrack=0; itrack<ntracks; itrack++){
1241 const Double_t kRadius = 3; // beam pipe radius
1242 const Double_t kMaxStep = 5; // max step
1243 const Double_t kMaxD = 123456; // max distance to prim vertex
1244 Double_t fieldZ = AliTracker::GetBz(); //
1245 AliESDtrack * track = esd->GetTrack(itrack);
1246 if (!track) continue;
1247 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
a7265806 1248 AliTracker::PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
ff8bb5ae 1249 track->RelateToVertex(esd->GetVertex(),fieldZ, kMaxD);
1250 }
1251
596a855f 1252 return kTRUE;
1253}
1254
1255//_____________________________________________________________________________
d64bd07d 1256Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
1257 //
1258 // Remove the data which are not needed for the physics analysis.
1259 //
1260
1261 AliInfo("Cleaning the ESD...");
1262
1263 const AliESDVertex *vertex=esd->GetVertex();
1264 Double_t vz=vertex->GetZv();
1265
1266 Int_t nTracks=esd->GetNumberOfTracks();
1267 for (Int_t i=0; i<nTracks; i++) {
1268 AliESDtrack *track=esd->GetTrack(i);
1269
1270 Float_t xy,z; track->GetImpactParameters(xy,z);
1271 if (TMath::Abs(xy) < 50.) continue;
1272 if (vertex->GetStatus())
1273 if (TMath::Abs(vz-z) < 5.) continue;
1274
1275 esd->RemoveTrack(i);
1276 }
1277
1278 return kTRUE;
1279}
1280
1281//_____________________________________________________________________________
af885e0f 1282Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
596a855f 1283{
1284// fill the event summary data
1285
87932dab 1286 AliCodeTimerAuto("")
030b532d 1287
596a855f 1288 TString detStr = detectors;
b8cd5251 1289 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1290 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1291 AliReconstructor* reconstructor = GetReconstructor(iDet);
1292 if (!reconstructor) continue;
1293
1294 if (!ReadESD(esd, fgkDetectorName[iDet])) {
1295 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
b26c3770 1296 TTree* clustersTree = NULL;
1297 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
1298 fLoader[iDet]->LoadRecPoints("read");
1299 clustersTree = fLoader[iDet]->TreeR();
1300 if (!clustersTree) {
1301 AliError(Form("Can't get the %s clusters tree",
1302 fgkDetectorName[iDet]));
1303 if (fStopOnError) return kFALSE;
1304 }
1305 }
1306 if (fRawReader && !reconstructor->HasDigitConversion()) {
1307 reconstructor->FillESD(fRawReader, clustersTree, esd);
1308 } else {
1309 TTree* digitsTree = NULL;
1310 if (fLoader[iDet]) {
1311 fLoader[iDet]->LoadDigits("read");
1312 digitsTree = fLoader[iDet]->TreeD();
1313 if (!digitsTree) {
1314 AliError(Form("Can't get the %s digits tree",
1315 fgkDetectorName[iDet]));
1316 if (fStopOnError) return kFALSE;
1317 }
1318 }
1319 reconstructor->FillESD(digitsTree, clustersTree, esd);
1320 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1321 }
1322 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
1323 fLoader[iDet]->UnloadRecPoints();
1324 }
1325
b8cd5251 1326 if (fRawReader) {
1327 reconstructor->FillESD(fRunLoader, fRawReader, esd);
1328 } else {
1329 reconstructor->FillESD(fRunLoader, esd);
24f7a148 1330 }
b8cd5251 1331 if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
596a855f 1332 }
1333 }
1334
1335 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
815c2b38 1336 AliError(Form("the following detectors were not found: %s",
1337 detStr.Data()));
596a855f 1338 if (fStopOnError) return kFALSE;
1339 }
1340
1341 return kTRUE;
1342}
1343
b647652d 1344//_____________________________________________________________________________
af885e0f 1345Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
b647652d 1346{
1347 // Reads the trigger decision which is
1348 // stored in Trigger.root file and fills
1349 // the corresponding esd entries
1350
87932dab 1351 AliCodeTimerAuto("")
1352
b647652d 1353 AliInfo("Filling trigger information into the ESD");
1354
1355 if (fRawReader) {
1356 AliCTPRawStream input(fRawReader);
1357 if (!input.Next()) {
1358 AliError("No valid CTP (trigger) DDL raw data is found ! The trigger information is not stored in the ESD !");
1359 return kFALSE;
1360 }
1361 esd->SetTriggerMask(input.GetClassMask());
1362 esd->SetTriggerCluster(input.GetClusterMask());
1363 }
1364 else {
1365 AliRunLoader *runloader = AliRunLoader::GetRunLoader();
1366 if (runloader) {
1367 if (!runloader->LoadTrigger()) {
1368 AliCentralTrigger *aCTP = runloader->GetTrigger();
1369 esd->SetTriggerMask(aCTP->GetClassMask());
1370 esd->SetTriggerCluster(aCTP->GetClusterMask());
1371 }
1372 else {
1373 AliWarning("No trigger can be loaded! The trigger information is not stored in the ESD !");
1374 return kFALSE;
1375 }
1376 }
1377 else {
1378 AliError("No run loader is available! The trigger information is not stored in the ESD !");
1379 return kFALSE;
1380 }
1381 }
1382
1383 return kTRUE;
1384}
596a855f 1385
001397cd 1386
1387
1388
1389
1390//_____________________________________________________________________________
af885e0f 1391Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
001397cd 1392{
1393 //
1394 // Filling information from RawReader Header
1395 //
1396
1397 AliInfo("Filling information from RawReader Header");
31fd97b2 1398 esd->SetBunchCrossNumber(0);
1399 esd->SetOrbitNumber(0);
9bcc1e45 1400 esd->SetPeriodNumber(0);
001397cd 1401 esd->SetTimeStamp(0);
1402 esd->SetEventType(0);
1403 const AliRawEventHeaderBase * eventHeader = fRawReader->GetEventHeader();
1404 if (eventHeader){
9bcc1e45 1405
1406 const UInt_t *id = eventHeader->GetP("Id");
1407 esd->SetBunchCrossNumber((id)[1]&0x00000fff);
1408 esd->SetOrbitNumber((((id)[0]<<20)&0xf00000)|(((id)[1]>>12)&0xfffff));
1409 esd->SetPeriodNumber(((id)[0]>>4)&0x0fffffff);
1410
001397cd 1411 esd->SetTimeStamp((eventHeader->Get("Timestamp")));
31fd97b2 1412 esd->SetEventType((eventHeader->Get("Type")));
001397cd 1413 }
1414
1415 return kTRUE;
1416}
1417
1418
596a855f 1419//_____________________________________________________________________________
1420Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
1421{
1422// check whether detName is contained in detectors
1423// if yes, it is removed from detectors
1424
1425 // check if all detectors are selected
1426 if ((detectors.CompareTo("ALL") == 0) ||
1427 detectors.BeginsWith("ALL ") ||
1428 detectors.EndsWith(" ALL") ||
1429 detectors.Contains(" ALL ")) {
1430 detectors = "ALL";
1431 return kTRUE;
1432 }
1433
1434 // search for the given detector
1435 Bool_t result = kFALSE;
1436 if ((detectors.CompareTo(detName) == 0) ||
1437 detectors.BeginsWith(detName+" ") ||
1438 detectors.EndsWith(" "+detName) ||
1439 detectors.Contains(" "+detName+" ")) {
1440 detectors.ReplaceAll(detName, "");
1441 result = kTRUE;
1442 }
1443
1444 // clean up the detectors string
1445 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1446 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1447 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1448
1449 return result;
1450}
e583c30d 1451
1452//_____________________________________________________________________________
f08fc9f5 1453Bool_t AliReconstruction::InitRunLoader()
1454{
1455// get or create the run loader
1456
1457 if (gAlice) delete gAlice;
1458 gAlice = NULL;
1459
b26c3770 1460 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
1461 // load all base libraries to get the loader classes
1462 TString libs = gSystem->GetLibraries();
1463 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1464 TString detName = fgkDetectorName[iDet];
1465 if (detName == "HLT") continue;
1466 if (libs.Contains("lib" + detName + "base.so")) continue;
1467 gSystem->Load("lib" + detName + "base.so");
1468 }
f08fc9f5 1469 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
1470 if (!fRunLoader) {
1471 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1472 CleanUp();
1473 return kFALSE;
1474 }
b26c3770 1475 fRunLoader->CdGAFile();
1476 if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
1477 if (fRunLoader->LoadgAlice() == 0) {
1478 gAlice = fRunLoader->GetAliRun();
c84a5e9e 1479 AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
b26c3770 1480 }
f08fc9f5 1481 }
1482 if (!gAlice && !fRawReader) {
1483 AliError(Form("no gAlice object found in file %s",
1484 fGAliceFileName.Data()));
1485 CleanUp();
1486 return kFALSE;
1487 }
1488
6cae184e 1489 //PH This is a temporary fix to give access to the kinematics
1490 //PH that is needed for the labels of ITS clusters
1491 fRunLoader->LoadKinematics();
1492
f08fc9f5 1493 } else { // galice.root does not exist
1494 if (!fRawReader) {
1495 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
1496 CleanUp();
1497 return kFALSE;
1498 }
1499 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
1500 AliConfig::GetDefaultEventFolderName(),
1501 "recreate");
1502 if (!fRunLoader) {
1503 AliError(Form("could not create run loader in file %s",
1504 fGAliceFileName.Data()));
1505 CleanUp();
1506 return kFALSE;
1507 }
1508 fRunLoader->MakeTree("E");
1509 Int_t iEvent = 0;
1510 while (fRawReader->NextEvent()) {
1511 fRunLoader->SetEventNumber(iEvent);
1512 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1513 iEvent, iEvent);
1514 fRunLoader->MakeTree("H");
1515 fRunLoader->TreeE()->Fill();
1516 iEvent++;
1517 }
1518 fRawReader->RewindEvents();
973388c2 1519 if (fNumberOfEventsPerFile > 0)
1520 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
1521 else
1522 fRunLoader->SetNumberOfEventsPerFile(iEvent);
f08fc9f5 1523 fRunLoader->WriteHeader("OVERWRITE");
1524 fRunLoader->CdGAFile();
1525 fRunLoader->Write(0, TObject::kOverwrite);
1526// AliTracker::SetFieldMap(???);
1527 }
1528
1529 return kTRUE;
1530}
1531
1532//_____________________________________________________________________________
b8cd5251 1533AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
c757bafd 1534{
f08fc9f5 1535// get the reconstructor object and the loader for a detector
c757bafd 1536
b8cd5251 1537 if (fReconstructor[iDet]) return fReconstructor[iDet];
1538
1539 // load the reconstructor object
1540 TPluginManager* pluginManager = gROOT->GetPluginManager();
1541 TString detName = fgkDetectorName[iDet];
1542 TString recName = "Ali" + detName + "Reconstructor";
f08fc9f5 1543 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
b8cd5251 1544
1545 if (detName == "HLT") {
1546 if (!gROOT->GetClass("AliLevel3")) {
4aa41877 1547 gSystem->Load("libAliHLTSrc.so");
1548 gSystem->Load("libAliHLTMisc.so");
1549 gSystem->Load("libAliHLTHough.so");
1550 gSystem->Load("libAliHLTComp.so");
b8cd5251 1551 }
1552 }
1553
1554 AliReconstructor* reconstructor = NULL;
1555 // first check if a plugin is defined for the reconstructor
1556 TPluginHandler* pluginHandler =
1557 pluginManager->FindHandler("AliReconstructor", detName);
f08fc9f5 1558 // if not, add a plugin for it
1559 if (!pluginHandler) {
b8cd5251 1560 AliDebug(1, Form("defining plugin for %s", recName.Data()));
b26c3770 1561 TString libs = gSystem->GetLibraries();
1562 if (libs.Contains("lib" + detName + "base.so") ||
1563 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
b8cd5251 1564 pluginManager->AddHandler("AliReconstructor", detName,
1565 recName, detName + "rec", recName + "()");
1566 } else {
1567 pluginManager->AddHandler("AliReconstructor", detName,
1568 recName, detName, recName + "()");
c757bafd 1569 }
b8cd5251 1570 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
1571 }
1572 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1573 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
c757bafd 1574 }
b8cd5251 1575 if (reconstructor) {
1576 TObject* obj = fOptions.FindObject(detName.Data());
1577 if (obj) reconstructor->SetOption(obj->GetTitle());
b26c3770 1578 reconstructor->Init(fRunLoader);
b8cd5251 1579 fReconstructor[iDet] = reconstructor;
1580 }
1581
f08fc9f5 1582 // get or create the loader
1583 if (detName != "HLT") {
1584 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
1585 if (!fLoader[iDet]) {
1586 AliConfig::Instance()
1587 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
1588 detName, detName);
1589 // first check if a plugin is defined for the loader
bb0901a4 1590 pluginHandler =
f08fc9f5 1591 pluginManager->FindHandler("AliLoader", detName);
1592 // if not, add a plugin for it
1593 if (!pluginHandler) {
1594 TString loaderName = "Ali" + detName + "Loader";
1595 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
1596 pluginManager->AddHandler("AliLoader", detName,
1597 loaderName, detName + "base",
1598 loaderName + "(const char*, TFolder*)");
1599 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
1600 }
1601 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1602 fLoader[iDet] =
1603 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
1604 fRunLoader->GetEventFolder());
1605 }
1606 if (!fLoader[iDet]) { // use default loader
1607 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
1608 }
1609 if (!fLoader[iDet]) {
1610 AliWarning(Form("couldn't get loader for %s", detName.Data()));
6667b602 1611 if (fStopOnError) return NULL;
f08fc9f5 1612 } else {
1613 fRunLoader->AddLoader(fLoader[iDet]);
1614 fRunLoader->CdGAFile();
1615 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
1616 fRunLoader->Write(0, TObject::kOverwrite);
1617 }
1618 }
1619 }
1620
b8cd5251 1621 return reconstructor;
c757bafd 1622}
1623
1624//_____________________________________________________________________________
2257f27e 1625Bool_t AliReconstruction::CreateVertexer()
1626{
1627// create the vertexer
1628
b8cd5251 1629 fVertexer = NULL;
1630 AliReconstructor* itsReconstructor = GetReconstructor(0);
59697224 1631 if (itsReconstructor) {
b8cd5251 1632 fVertexer = itsReconstructor->CreateVertexer(fRunLoader);
2257f27e 1633 }
b8cd5251 1634 if (!fVertexer) {
815c2b38 1635 AliWarning("couldn't create a vertexer for ITS");
2257f27e 1636 if (fStopOnError) return kFALSE;
1637 }
1638
1639 return kTRUE;
1640}
1641
1642//_____________________________________________________________________________
b8cd5251 1643Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
24f7a148 1644{
f08fc9f5 1645// create the trackers
24f7a148 1646
b8cd5251 1647 TString detStr = detectors;
1648 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1649 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1650 AliReconstructor* reconstructor = GetReconstructor(iDet);
1651 if (!reconstructor) continue;
1652 TString detName = fgkDetectorName[iDet];
1f46a9ae 1653 if (detName == "HLT") {
1654 fRunHLTTracking = kTRUE;
1655 continue;
1656 }
e66fbafb 1657 if (detName == "MUON") {
1658 fRunMuonTracking = kTRUE;
1659 continue;
1660 }
1661
f08fc9f5 1662
1663 fTracker[iDet] = reconstructor->CreateTracker(fRunLoader);
1664 if (!fTracker[iDet] && (iDet < 7)) {
1665 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
8250d5f5 1666 if (fStopOnError) return kFALSE;
1667 }
1668 }
1669
24f7a148 1670 return kTRUE;
1671}
1672
1673//_____________________________________________________________________________
b26c3770 1674void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
e583c30d 1675{
1676// delete trackers and the run loader and close and delete the file
1677
b8cd5251 1678 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1679 delete fReconstructor[iDet];
1680 fReconstructor[iDet] = NULL;
1681 fLoader[iDet] = NULL;
1682 delete fTracker[iDet];
1683 fTracker[iDet] = NULL;
1684 }
1685 delete fVertexer;
1686 fVertexer = NULL;
9178838a 1687 delete fDiamondProfile;
1688 fDiamondProfile = NULL;
e583c30d 1689
1690 delete fRunLoader;
1691 fRunLoader = NULL;
b649205a 1692 delete fRawReader;
1693 fRawReader = NULL;
e583c30d 1694
1695 if (file) {
1696 file->Close();
1697 delete file;
1698 }
b26c3770 1699
1700 if (fileOld) {
1701 fileOld->Close();
1702 delete fileOld;
1703 gSystem->Unlink("AliESDs.old.root");
1704 }
e583c30d 1705}
24f7a148 1706
1707
1708//_____________________________________________________________________________
af885e0f 1709
1710Bool_t AliReconstruction::ReadESD(AliESDEvent*& esd, const char* recStep) const
24f7a148 1711{
1712// read the ESD event from a file
1713
1714 if (!esd) return kFALSE;
1715 char fileName[256];
1716 sprintf(fileName, "ESD_%d.%d_%s.root",
31fd97b2 1717 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
24f7a148 1718 if (gSystem->AccessPathName(fileName)) return kFALSE;
1719
f3a97c86 1720 AliInfo(Form("reading ESD from file %s", fileName));
815c2b38 1721 AliDebug(1, Form("reading ESD from file %s", fileName));
24f7a148 1722 TFile* file = TFile::Open(fileName);
1723 if (!file || !file->IsOpen()) {
815c2b38 1724 AliError(Form("opening %s failed", fileName));
24f7a148 1725 delete file;
1726 return kFALSE;
1727 }
1728
1729 gROOT->cd();
1730 delete esd;
af885e0f 1731 esd = (AliESDEvent*) file->Get("ESD");
24f7a148 1732 file->Close();
1733 delete file;
1734 return kTRUE;
af885e0f 1735
24f7a148 1736}
1737
af885e0f 1738
1739
24f7a148 1740//_____________________________________________________________________________
af885e0f 1741void AliReconstruction::WriteESD(AliESDEvent* esd, const char* recStep) const
24f7a148 1742{
1743// write the ESD event to a file
1744
1745 if (!esd) return;
1746 char fileName[256];
1747 sprintf(fileName, "ESD_%d.%d_%s.root",
31fd97b2 1748 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
24f7a148 1749
815c2b38 1750 AliDebug(1, Form("writing ESD to file %s", fileName));
24f7a148 1751 TFile* file = TFile::Open(fileName, "recreate");
1752 if (!file || !file->IsOpen()) {
815c2b38 1753 AliError(Form("opening %s failed", fileName));
24f7a148 1754 } else {
1755 esd->Write("ESD");
1756 file->Close();
1757 }
1758 delete file;
1759}
f3a97c86 1760
1761
1762
1763
f3a97c86 1764
a7807689 1765//_____________________________________________________________________________
f29f1726 1766void AliReconstruction::ESDFile2AODFile(TFile* esdFile, TFile* aodFile)
a7807689 1767{
f29f1726 1768 // write all files from the given esd file to an aod file
85ba66b8 1769
f29f1726 1770 // create an AliAOD object
1771 AliAODEvent *aod = new AliAODEvent();
1772 aod->CreateStdContent();
1773
1774 // go to the file
1775 aodFile->cd();
1776
1777 // create the tree
b97637d4 1778 TTree *aodTree = new TTree("aodTree", "AliAOD tree");
f29f1726 1779 aodTree->Branch(aod->GetList());
1780
1781 // connect to ESD
1782 TTree *t = (TTree*) esdFile->Get("esdTree");
af885e0f 1783 AliESDEvent *esd = new AliESDEvent();
53ec9628 1784 esd->ReadFromTree(t);
f29f1726 1785
53ec9628 1786 Int_t nEvents = t->GetEntries();
f29f1726 1787
1788 // set arrays and pointers
1789 Float_t posF[3];
1790 Double_t pos[3];
1791 Double_t p[3];
1792 Double_t covVtx[6];
1793 Double_t covTr[21];
1794 Double_t pid[10];
1795
1796 // loop over events and fill them
1797 for (Int_t iEvent = 0; iEvent < nEvents; ++iEvent) {
53ec9628 1798 t->GetEntry(iEvent);
f29f1726 1799
1800 // Multiplicity information needed by the header (to be revised!)
1801 Int_t nTracks = esd->GetNumberOfTracks();
1802 Int_t nPosTracks = 0;
1803 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack)
b97637d4 1804 if (esd->GetTrack(iTrack)->Charge()> 0) nPosTracks++;
f29f1726 1805
85ba66b8 1806 // Access the header
1807 AliAODHeader *header = aod->GetHeader();
1808
1809 // fill the header
ade23daf 1810 header->SetRunNumber (esd->GetRunNumber() );
1811 header->SetBunchCrossNumber(esd->GetBunchCrossNumber());
1812 header->SetOrbitNumber (esd->GetOrbitNumber() );
1813 header->SetPeriodNumber (esd->GetPeriodNumber() );
1814 header->SetTriggerMask (esd->GetTriggerMask() );
1815 header->SetTriggerCluster (esd->GetTriggerCluster() );
1816 header->SetEventType (esd->GetEventType() );
1817 header->SetMagneticField (esd->GetMagneticField() );
1818 header->SetZDCN1Energy (esd->GetZDCN1Energy() );
1819 header->SetZDCP1Energy (esd->GetZDCP1Energy() );
1820 header->SetZDCN2Energy (esd->GetZDCN2Energy() );
1821 header->SetZDCP2Energy (esd->GetZDCP2Energy() );
1822 header->SetZDCEMEnergy (esd->GetZDCEMEnergy() );
a1d4139d 1823 header->SetRefMultiplicity (nTracks);
1824 header->SetRefMultiplicityPos(nPosTracks);
1825 header->SetRefMultiplicityNeg(nTracks - nPosTracks);
1826 header->SetMuonMagFieldScale(-999.); // FIXME
1827 header->SetCentrality(-999.); // FIXME
f29f1726 1828
1829 Int_t nV0s = esd->GetNumberOfV0s();
1830 Int_t nCascades = esd->GetNumberOfCascades();
1831 Int_t nKinks = esd->GetNumberOfKinks();
1832 Int_t nVertices = nV0s + nCascades + nKinks;
1833
1834 aod->ResetStd(nTracks, nVertices);
1835 AliAODTrack *aodTrack;
1836
f29f1726 1837 // Array to take into account the tracks already added to the AOD
1838 Bool_t * usedTrack = NULL;
1839 if (nTracks>0) {
1840 usedTrack = new Bool_t[nTracks];
1841 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
1842 }
1843 // Array to take into account the V0s already added to the AOD
1844 Bool_t * usedV0 = NULL;
1845 if (nV0s>0) {
1846 usedV0 = new Bool_t[nV0s];
1847 for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
1848 }
1849 // Array to take into account the kinks already added to the AOD
1850 Bool_t * usedKink = NULL;
1851 if (nKinks>0) {
1852 usedKink = new Bool_t[nKinks];
1853 for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
1854 }
1855
1856 // Access to the AOD container of vertices
1857 TClonesArray &vertices = *(aod->GetVertices());
1858 Int_t jVertices=0;
1859
1860 // Access to the AOD container of tracks
1861 TClonesArray &tracks = *(aod->GetTracks());
1862 Int_t jTracks=0;
1863
1864 // Add primary vertex. The primary tracks will be defined
1865 // after the loops on the composite objects (V0, cascades, kinks)
1866 const AliESDVertex *vtx = esd->GetPrimaryVertex();
1867
1868 vtx->GetXYZ(pos); // position
1869 vtx->GetCovMatrix(covVtx); //covariance matrix
1870
1871 AliAODVertex * primary = new(vertices[jVertices++])
02153d58 1872 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
f29f1726 1873
1874 // Create vertices starting from the most complex objects
1875
1876 // Cascades
1877 for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
1878 AliESDcascade *cascade = esd->GetCascade(nCascade);
1879
1880 cascade->GetXYZ(pos[0], pos[1], pos[2]);
1881 cascade->GetPosCovXi(covVtx);
1882
1883 // Add the cascade vertex
1884 AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(pos,
1885 covVtx,
1886 cascade->GetChi2Xi(), // = chi2/NDF since NDF = 2*2-3
1887 primary,
02153d58 1888 nCascade,
f29f1726 1889 AliAODVertex::kCascade);
1890
1891 primary->AddDaughter(vcascade);
1892
1893 // Add the V0 from the cascade. The ESD class have to be optimized...
85ba66b8 1894 // Now we have to search for the corresponding V0 in the list of V0s
f29f1726 1895 // using the indeces of the positive and negative tracks
1896
1897 Int_t posFromV0 = cascade->GetPindex();
1898 Int_t negFromV0 = cascade->GetNindex();
1899
1900
1901 AliESDv0 * v0 = 0x0;
1902 Int_t indV0 = -1;
1903
1904 for (Int_t iV0=0; iV0<nV0s; ++iV0) {
1905
1906 v0 = esd->GetV0(iV0);
1907 Int_t posV0 = v0->GetPindex();
1908 Int_t negV0 = v0->GetNindex();
1909
1910 if (posV0==posFromV0 && negV0==negFromV0) {
1911 indV0 = iV0;
1912 break;
1913 }
1914 }
1915
1916 AliAODVertex * vV0FromCascade = 0x0;
1917
1918 if (indV0>-1 && !usedV0[indV0] ) {
1919
1920 // the V0 exists in the array of V0s and is not used
1921
1922 usedV0[indV0] = kTRUE;
1923
1924 v0->GetXYZ(pos[0], pos[1], pos[2]);
1925 v0->GetPosCov(covVtx);
1926
1927 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
1928 covVtx,
1929 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
1930 vcascade,
02153d58 1931 indV0,
f29f1726 1932 AliAODVertex::kV0);
1933 } else {
1934
1935 // the V0 doesn't exist in the array of V0s or was used
1936 cerr << "Error: event " << iEvent << " cascade " << nCascade
1937 << " The V0 " << indV0
1938 << " doesn't exist in the array of V0s or was used!" << endl;
1939
1940 cascade->GetXYZ(pos[0], pos[1], pos[2]);
1941 cascade->GetPosCov(covVtx);
1942
1943 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
1944 covVtx,
1945 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
1946 vcascade,
02153d58 1947 indV0,
f29f1726 1948 AliAODVertex::kV0);
1949 vcascade->AddDaughter(vV0FromCascade);
1950 }
1951
1952 // Add the positive tracks from the V0
1953
1954 if (! usedTrack[posFromV0]) {
1955
1956 usedTrack[posFromV0] = kTRUE;
1957
1958 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
1959 esdTrack->GetPxPyPz(p);
1960 esdTrack->GetXYZ(pos);
1961 esdTrack->GetCovarianceXYZPxPyPz(covTr);
1962 esdTrack->GetESDpid(pid);
1963
1964 vV0FromCascade->AddDaughter(aodTrack =
1965 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
1966 esdTrack->GetLabel(),
1967 p,
1968 kTRUE,
1969 pos,
1970 kFALSE,
1971 covTr,
b97637d4 1972 (Short_t)esdTrack->Charge(),
f29f1726 1973 esdTrack->GetITSClusterMap(),
1974 pid,
1975 vV0FromCascade,
1976 kTRUE, // check if this is right
1977 kFALSE, // check if this is right
1978 AliAODTrack::kSecondary)
1979 );
1980 aodTrack->ConvertAliPIDtoAODPID();
1981 }
1982 else {
1983 cerr << "Error: event " << iEvent << " cascade " << nCascade
1984 << " track " << posFromV0 << " has already been used!" << endl;
1985 }
1986
1987 // Add the negative tracks from the V0
1988
1989 if (!usedTrack[negFromV0]) {
1990
1991 usedTrack[negFromV0] = kTRUE;
1992
1993 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
1994 esdTrack->GetPxPyPz(p);
1995 esdTrack->GetXYZ(pos);
1996 esdTrack->GetCovarianceXYZPxPyPz(covTr);
1997 esdTrack->GetESDpid(pid);
1998
1999 vV0FromCascade->AddDaughter(aodTrack =
2000 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2001 esdTrack->GetLabel(),
2002 p,
2003 kTRUE,
2004 pos,
2005 kFALSE,
2006 covTr,
b97637d4 2007 (Short_t)esdTrack->Charge(),
f29f1726 2008 esdTrack->GetITSClusterMap(),
2009 pid,
2010 vV0FromCascade,
2011 kTRUE, // check if this is right
2012 kFALSE, // check if this is right
2013 AliAODTrack::kSecondary)
2014 );
2015 aodTrack->ConvertAliPIDtoAODPID();
2016 }
2017 else {
2018 cerr << "Error: event " << iEvent << " cascade " << nCascade
2019 << " track " << negFromV0 << " has already been used!" << endl;
2020 }
2021
2022 // Add the bachelor track from the cascade
2023
2024 Int_t bachelor = cascade->GetBindex();
2025
2026 if(!usedTrack[bachelor]) {
2027
2028 usedTrack[bachelor] = kTRUE;
2029
2030 AliESDtrack *esdTrack = esd->GetTrack(bachelor);
2031 esdTrack->GetPxPyPz(p);
2032 esdTrack->GetXYZ(pos);
2033 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2034 esdTrack->GetESDpid(pid);
2035
2036 vcascade->AddDaughter(aodTrack =
2037 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2038 esdTrack->GetLabel(),
2039 p,
2040 kTRUE,
2041 pos,
2042 kFALSE,
2043 covTr,
b97637d4 2044 (Short_t)esdTrack->Charge(),
f29f1726 2045 esdTrack->GetITSClusterMap(),
2046 pid,
2047 vcascade,
2048 kTRUE, // check if this is right
2049 kFALSE, // check if this is right
2050 AliAODTrack::kSecondary)
2051 );
2052 aodTrack->ConvertAliPIDtoAODPID();
2053 }
2054 else {
2055 cerr << "Error: event " << iEvent << " cascade " << nCascade
2056 << " track " << bachelor << " has already been used!" << endl;
2057 }
2058
2059 // Add the primary track of the cascade (if any)
2060
2061 } // end of the loop on cascades
2062
2063 // V0s
2064
2065 for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
2066
2067 if (usedV0[nV0]) continue; // skip if aready added to the AOD
2068
2069 AliESDv0 *v0 = esd->GetV0(nV0);
2070
2071 v0->GetXYZ(pos[0], pos[1], pos[2]);
2072 v0->GetPosCov(covVtx);
2073
2074 AliAODVertex * vV0 =
2075 new(vertices[jVertices++]) AliAODVertex(pos,
2076 covVtx,
2077 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2078 primary,
02153d58 2079 nV0,
f29f1726 2080 AliAODVertex::kV0);
2081 primary->AddDaughter(vV0);
2082
2083 Int_t posFromV0 = v0->GetPindex();
2084 Int_t negFromV0 = v0->GetNindex();
2085
2086 // Add the positive tracks from the V0
2087
2088 if (!usedTrack[posFromV0]) {
2089
2090 usedTrack[posFromV0] = kTRUE;
2091
2092 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2093 esdTrack->GetPxPyPz(p);
2094 esdTrack->GetXYZ(pos);
2095 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2096 esdTrack->GetESDpid(pid);
2097
2098 vV0->AddDaughter(aodTrack =
2099 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2100 esdTrack->GetLabel(),
2101 p,
2102 kTRUE,
2103 pos,
2104 kFALSE,
2105 covTr,
b97637d4 2106 (Short_t)esdTrack->Charge(),
f29f1726 2107 esdTrack->GetITSClusterMap(),
2108 pid,
2109 vV0,
2110 kTRUE, // check if this is right
2111 kFALSE, // check if this is right
2112 AliAODTrack::kSecondary)
2113 );
2114 aodTrack->ConvertAliPIDtoAODPID();
2115 }
2116 else {
2117 cerr << "Error: event " << iEvent << " V0 " << nV0
2118 << " track " << posFromV0 << " has already been used!" << endl;
2119 }
a7807689 2120
f29f1726 2121 // Add the negative tracks from the V0
2122
2123 if (!usedTrack[negFromV0]) {
2124
2125 usedTrack[negFromV0] = kTRUE;
2126
2127 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2128 esdTrack->GetPxPyPz(p);
2129 esdTrack->GetXYZ(pos);
2130 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2131 esdTrack->GetESDpid(pid);
2132
2133 vV0->AddDaughter(aodTrack =
2134 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2135 esdTrack->GetLabel(),
2136 p,
2137 kTRUE,
2138 pos,
2139 kFALSE,
2140 covTr,
b97637d4 2141 (Short_t)esdTrack->Charge(),
f29f1726 2142 esdTrack->GetITSClusterMap(),
2143 pid,
2144 vV0,
2145 kTRUE, // check if this is right
2146 kFALSE, // check if this is right
2147 AliAODTrack::kSecondary)
2148 );
2149 aodTrack->ConvertAliPIDtoAODPID();
2150 }
2151 else {
2152 cerr << "Error: event " << iEvent << " V0 " << nV0
2153 << " track " << negFromV0 << " has already been used!" << endl;
2154 }
2155
2156 } // end of the loop on V0s
2157
2158 // Kinks: it is a big mess the access to the information in the kinks
2159 // The loop is on the tracks in order to find the mother and daugther of each kink
2160
2161
2162 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
2163
2164
2165 AliESDtrack * esdTrack = esd->GetTrack(iTrack);
2166
2167 Int_t ikink = esdTrack->GetKinkIndex(0);
2168
2169 if (ikink) {
2170 // Negative kink index: mother, positive: daughter
2171
2172 // Search for the second track of the kink
2173
2174 for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
2175
2176 AliESDtrack * esdTrack1 = esd->GetTrack(jTrack);
2177
2178 Int_t jkink = esdTrack1->GetKinkIndex(0);
2179
2180 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
2181
2182 // The two tracks are from the same kink
2183
2184 if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
2185
2186 Int_t imother = -1;
2187 Int_t idaughter = -1;
2188
2189 if (ikink<0 && jkink>0) {
2190
2191 imother = iTrack;
2192 idaughter = jTrack;
2193 }
2194 else if (ikink>0 && jkink<0) {
2195
2196 imother = jTrack;
2197 idaughter = iTrack;
2198 }
2199 else {
2200 cerr << "Error: Wrong combination of kink indexes: "
2201 << ikink << " " << jkink << endl;
2202 continue;
2203 }
2204
2205 // Add the mother track
2206
2207 AliAODTrack * mother = NULL;
2208
2209 if (!usedTrack[imother]) {
2210
2211 usedTrack[imother] = kTRUE;
2212
2213 AliESDtrack *esdTrack = esd->GetTrack(imother);
2214 esdTrack->GetPxPyPz(p);
2215 esdTrack->GetXYZ(pos);
2216 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2217 esdTrack->GetESDpid(pid);
2218
2219 mother =
2220 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2221 esdTrack->GetLabel(),
2222 p,
2223 kTRUE,
2224 pos,
2225 kFALSE,
2226 covTr,
b97637d4 2227 (Short_t)esdTrack->Charge(),
f29f1726 2228 esdTrack->GetITSClusterMap(),
2229 pid,
2230 primary,
2231 kTRUE, // check if this is right
2232 kTRUE, // check if this is right
2233 AliAODTrack::kPrimary);
2234 primary->AddDaughter(mother);
2235 mother->ConvertAliPIDtoAODPID();
2236 }
2237 else {
2238 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2239 << " track " << imother << " has already been used!" << endl;
2240 }
2241
2242 // Add the kink vertex
2243 AliESDkink * kink = esd->GetKink(TMath::Abs(ikink)-1);
2244
2245 AliAODVertex * vkink =
2246 new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
2247 NULL,
2248 0.,
2249 mother,
02153d58 2250 esdTrack->GetID(), // This is the track ID of the mother's track!
f29f1726 2251 AliAODVertex::kKink);
2252 // Add the daughter track
2253
2254 AliAODTrack * daughter = NULL;
2255
2256 if (!usedTrack[idaughter]) {
2257
2258 usedTrack[idaughter] = kTRUE;
2259
2260 AliESDtrack *esdTrack = esd->GetTrack(idaughter);
2261 esdTrack->GetPxPyPz(p);
2262 esdTrack->GetXYZ(pos);
2263 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2264 esdTrack->GetESDpid(pid);
2265
2266 daughter =
2267 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2268 esdTrack->GetLabel(),
2269 p,
2270 kTRUE,
2271 pos,
2272 kFALSE,
2273 covTr,
b97637d4 2274 (Short_t)esdTrack->Charge(),
f29f1726 2275 esdTrack->GetITSClusterMap(),
2276 pid,
2277 vkink,
2278 kTRUE, // check if this is right
2279 kTRUE, // check if this is right
2280 AliAODTrack::kPrimary);
2281 vkink->AddDaughter(daughter);
2282 daughter->ConvertAliPIDtoAODPID();
2283 }
2284 else {
2285 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2286 << " track " << idaughter << " has already been used!" << endl;
2287 }
2288
2289
2290 }
2291 }
2292
2293 }
2294
2295 }
2296
2297
2298 // Tracks (primary and orphan)
2299
2300 for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
2301
2302
2303 if (usedTrack[nTrack]) continue;
2304
2305 AliESDtrack *esdTrack = esd->GetTrack(nTrack);
2306 esdTrack->GetPxPyPz(p);
2307 esdTrack->GetXYZ(pos);
2308 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2309 esdTrack->GetESDpid(pid);
2310
2311 Float_t impactXY, impactZ;
2312
2313 esdTrack->GetImpactParameters(impactXY,impactZ);
2314
2315 if (impactXY<3) {
2316 // track inside the beam pipe
2317
2318 primary->AddDaughter(aodTrack =
2319 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2320 esdTrack->GetLabel(),
2321 p,
2322 kTRUE,
2323 pos,
2324 kFALSE,
2325 covTr,
b97637d4 2326 (Short_t)esdTrack->Charge(),
f29f1726 2327 esdTrack->GetITSClusterMap(),
2328 pid,
2329 primary,
2330 kTRUE, // check if this is right
2331 kTRUE, // check if this is right
2332 AliAODTrack::kPrimary)
2333 );
2334 aodTrack->ConvertAliPIDtoAODPID();
2335 }
2336 else {
2337 // outside the beam pipe: orphan track
2338 aodTrack =
2339 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2340 esdTrack->GetLabel(),
2341 p,
2342 kTRUE,
2343 pos,
2344 kFALSE,
2345 covTr,
b97637d4 2346 (Short_t)esdTrack->Charge(),
f29f1726 2347 esdTrack->GetITSClusterMap(),
2348 pid,
2349 NULL,
2350 kFALSE, // check if this is right
2351 kFALSE, // check if this is right
2352 AliAODTrack::kOrphan);
2353 aodTrack->ConvertAliPIDtoAODPID();
2354 }
2355 } // end of loop on tracks
2356
2357 // muon tracks
2358 Int_t nMuTracks = esd->GetNumberOfMuonTracks();
2359 for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) {
2360
2361 AliESDMuonTrack *esdMuTrack = esd->GetMuonTrack(nMuTrack);
2362 p[0] = esdMuTrack->Px();
2363 p[1] = esdMuTrack->Py();
2364 p[2] = esdMuTrack->Pz();
2365 pos[0] = primary->GetX();
2366 pos[1] = primary->GetY();
2367 pos[2] = primary->GetZ();
2368
2369 // has to be changed once the muon pid is provided by the ESD
2370 for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.;
2371
85ba66b8 2372 primary->AddDaughter(aodTrack =
f29f1726 2373 new(tracks[jTracks++]) AliAODTrack(0, // no ID provided
2374 0, // no label provided
2375 p,
2376 kTRUE,
2377 pos,
2378 kFALSE,
2379 NULL, // no covariance matrix provided
b97637d4 2380 esdMuTrack->Charge(),
e704c7d4 2381 0, // ITSClusterMap is set below
f29f1726 2382 pid,
2383 primary,
85ba66b8 2384 kFALSE, // muon tracks are not used to fit the primary vtx
2385 kFALSE, // not used for vertex fit
f29f1726 2386 AliAODTrack::kPrimary)
2387 );
85ba66b8 2388
2389 aodTrack->SetHitsPatternInTrigCh(esdMuTrack->GetHitsPatternInTrigCh());
2390 Int_t track2Trigger = esdMuTrack->GetMatchTrigger();
2391 aodTrack->SetMatchTrigger(track2Trigger);
2392 if (track2Trigger)
2393 aodTrack->SetChi2MatchTrigger(esdMuTrack->GetChi2MatchTrigger());
2394 else
2395 aodTrack->SetChi2MatchTrigger(0.);
f29f1726 2396 }
2397
2398 // Access to the AOD container of clusters
2399 TClonesArray &clusters = *(aod->GetClusters());
2400 Int_t jClusters=0;
2401
2402 // Calo Clusters
2403 Int_t nClusters = esd->GetNumberOfCaloClusters();
2404
2405 for (Int_t iClust=0; iClust<nClusters; ++iClust) {
2406
2407 AliESDCaloCluster * cluster = esd->GetCaloCluster(iClust);
2408
2409 Int_t id = cluster->GetID();
2410 Int_t label = -1;
53ec9628 2411 Float_t energy = cluster->E();
2412 cluster->GetPosition(posF);
f29f1726 2413 AliAODVertex *prodVertex = primary;
2414 AliAODTrack *primTrack = NULL;
2415 Char_t ttype=AliAODCluster::kUndef;
85ba66b8 2416
f29f1726 2417 if (cluster->IsPHOS()) ttype=AliAODCluster::kPHOSNeutral;
2418 else if (cluster->IsEMCAL()) {
85ba66b8 2419
f29f1726 2420 if (cluster->GetClusterType() == AliESDCaloCluster::kPseudoCluster)
2421 ttype = AliAODCluster::kEMCALPseudoCluster;
2422 else
2423 ttype = AliAODCluster::kEMCALClusterv1;
85ba66b8 2424
f29f1726 2425 }
85ba66b8 2426
f29f1726 2427 new(clusters[jClusters++]) AliAODCluster(id,
2428 label,
2429 energy,
2430 pos,
2431 NULL, // no covariance matrix provided
2432 NULL, // no pid for clusters provided
2433 prodVertex,
2434 primTrack,
2435 ttype);
85ba66b8 2436
f29f1726 2437 } // end of loop on calo clusters
85ba66b8 2438
2439 // tracklets
2440 const AliMultiplicity *mult = esd->GetMultiplicity();
2441 if (mult) {
2442 if (mult->GetNumberOfTracklets()>0) {
2443 aod->GetTracklets()->CreateContainer(mult->GetNumberOfTracklets());
2444
2445 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
2446 aod->GetTracklets()->SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n));
2447 }
2448 }
2449 } else {
2450 Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
2451 }
2452
f29f1726 2453 delete [] usedTrack;
2454 delete [] usedV0;
2455 delete [] usedKink;
85ba66b8 2456
f29f1726 2457 // fill the tree for this event
2458 aodTree->Fill();
2459 } // end of event loop
85ba66b8 2460
f29f1726 2461 aodTree->GetUserInfo()->Add(aod);
85ba66b8 2462
f29f1726 2463 // close ESD file
2464 esdFile->Close();
85ba66b8 2465
f29f1726 2466 // write the tree to the specified file
2467 aodFile = aodTree->GetCurrentFile();
2468 aodFile->cd();
2469 aodTree->Write();
85ba66b8 2470
a7807689 2471 return;
2472}
2473
af885e0f 2474void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
98937d93 2475{
2476 // Write space-points which are then used in the alignment procedures
2477 // For the moment only ITS, TRD and TPC
2478
2479 // Load TOF clusters
d528ee75 2480 if (fTracker[3]){
2481 fLoader[3]->LoadRecPoints("read");
2482 TTree* tree = fLoader[3]->TreeR();
2483 if (!tree) {
2484 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
2485 return;
2486 }
2487 fTracker[3]->LoadClusters(tree);
98937d93 2488 }
98937d93 2489 Int_t ntracks = esd->GetNumberOfTracks();
2490 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2491 {
2492 AliESDtrack *track = esd->GetTrack(itrack);
2493 Int_t nsp = 0;
ef7253ac 2494 Int_t idx[200];
98937d93 2495 for (Int_t iDet = 3; iDet >= 0; iDet--)
2496 nsp += track->GetNcls(iDet);
2497 if (nsp) {
2498 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2499 track->SetTrackPointArray(sp);
2500 Int_t isptrack = 0;
2501 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2502 AliTracker *tracker = fTracker[iDet];
2503 if (!tracker) continue;
2504 Int_t nspdet = track->GetNcls(iDet);
98937d93 2505 if (nspdet <= 0) continue;
2506 track->GetClusters(iDet,idx);
2507 AliTrackPoint p;
2508 Int_t isp = 0;
2509 Int_t isp2 = 0;
2510 while (isp < nspdet) {
2511 Bool_t isvalid = tracker->GetTrackPoint(idx[isp2],p); isp2++;
160db090 2512 const Int_t kNTPCmax = 159;
2513 if (iDet==1 && isp2>kNTPCmax) break; // to be fixed
98937d93 2514 if (!isvalid) continue;
2515 sp->AddPoint(isptrack,&p); isptrack++; isp++;
2516 }
98937d93 2517 }
2518 }
2519 }
d528ee75 2520 if (fTracker[3]){
2521 fTracker[3]->UnloadClusters();
2522 fLoader[3]->UnloadRecPoints();
2523 }
98937d93 2524}
2e3550da 2525
2526//_____________________________________________________________________________
af885e0f 2527void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2e3550da 2528{
2529 // The method reads the raw-data error log
2530 // accumulated within the rawReader.
2531 // It extracts the raw-data errors related to
2532 // the current event and stores them into
2533 // a TClonesArray inside the esd object.
2534
2535 if (!fRawReader) return;
2536
2537 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2538
2539 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2540 if (!log) continue;
2541 if (iEvent != log->GetEventNumber()) continue;
2542
2543 esd->AddRawDataErrorLog(log);
2544 }
2545
2546}
46698ae4 2547
2548TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString fName){
b545009a 2549 // Dump a file content into a char in TNamed
46698ae4 2550 ifstream in;
2551 in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2552 Int_t kBytes = (Int_t)in.tellg();
2553 printf("Size: %d \n",kBytes);
2554 TNamed *fn = 0;
2555 if(in.good()){
2556 char* memblock = new char [kBytes];
2557 in.seekg (0, ios::beg);
2558 in.read (memblock, kBytes);
2559 in.close();
2560 TString fData(memblock,kBytes);
2561 fn = new TNamed(fName,fData);
2562 printf("fData Size: %d \n",fData.Sizeof());
2563 printf("fName Size: %d \n",fName.Sizeof());
2564 printf("fn Size: %d \n",fn->Sizeof());
2565 delete[] memblock;
2566 }
2567 else{
2568 AliInfo(Form("Could not Open %s\n",fPath.Data()));
2569 }
2570
2571 return fn;
2572}
2573
2574void AliReconstruction::TNamedToFile(TTree* fTree, TString fName){
46698ae4 2575 // This is not really needed in AliReconstruction at the moment
2576 // but can serve as a template
2577
2578 TList *fList = fTree->GetUserInfo();
2579 TNamed *fn = (TNamed*)fList->FindObject(fName.Data());
2580 printf("fn Size: %d \n",fn->Sizeof());
2581
2582 TString fTmp(fn->GetName()); // to be 100% sure in principle fName also works
2583 const char* cdata = fn->GetTitle();
2584 printf("fTmp Size %d\n",fTmp.Sizeof());
2585
2586 int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2587 printf("calculated size %d\n",size);
2588 ofstream out(fName.Data(),ios::out | ios::binary);
2589 out.write(cdata,size);
2590 out.close();
2591
2592}
2593