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