Using symbolic names instead of volume paths (Raffaele)
[u/mrichter/AliRoot.git] / STEER / AliTagCreator.cxx
CommitLineData
b45e5084 1/**************************************************************************
2 * Author: Panos Christakoglou. *
3 * Contributors are mentioned in the code where appropriate. *
4 * *
5 * Permission to use, copy, modify and distribute this software and its *
6 * documentation strictly for non-commercial purposes is hereby granted *
7 * without fee, provided that the above copyright notice appears in all *
8 * copies and that both the copyright notice and this permission notice *
9 * appear in the supporting documentation. The authors make no claims *
10 * about the suitability of this software for any purpose. It is *
11 * provided "as is" without express or implied warranty. *
12 **************************************************************************/
13
14/* $Id$ */
15
16//-----------------------------------------------------------------
17// AliTagCreator class
18// This is the class to deal with the tag creation (post process)
19// Origin: Panos Christakoglou, UOA-CERN, Panos.Christakoglou@cern.ch
20//-----------------------------------------------------------------
21
22//ROOT
c336835b 23#include <Riostream.h>
b45e5084 24#include <TFile.h>
25#include <TString.h>
cb1645b7 26#include <TTree.h>
4fd84456 27#include <TSystem.h>
c336835b 28#include <TSystemDirectory.h>
4fd84456 29#include <TChain.h>
2bdb9d38 30#include <TLorentzVector.h>
b45e5084 31
32//ROOT-AliEn
33#include <TGrid.h>
b45e5084 34#include <TGridResult.h>
b45e5084 35
36//AliRoot
37#include "AliRunTag.h"
38#include "AliEventTag.h"
39#include "AliESD.h"
40#include "AliESDVertex.h"
b45e5084 41#include "AliLog.h"
42
43
44#include "AliTagCreator.h"
45
46
47ClassImp(AliTagCreator)
48
49
50//______________________________________________________________________________
442646f2 51 AliTagCreator::AliTagCreator() :
52 TObject(),
53 fSE("ALICE::CERN::se01"),
54 fgridpath(""),
55 fStorage(0)
56{
b45e5084 57 //==============Default constructor for a AliTagCreator==================
b45e5084 58}
59
60//______________________________________________________________________________
c336835b 61AliTagCreator::~AliTagCreator() {
b45e5084 62//================Default destructor for a AliTagCreator=======================
63}
64
65//______________________________________________________________________________
c336835b 66void AliTagCreator::SetStorage(Int_t storage) {
cb1645b7 67 // Sets correctly the storage: 0 for local, 1 for GRID
28afeb2e 68 fStorage = storage;
69 if(fStorage == 0)
70 AliInfo(Form("Tags will be stored locally...."));
71 if(fStorage == 1)
72 AliInfo(Form("Tags will be stored in the grid...."));
73 if((fStorage != 0)&&(fStorage != 1))
74 {
75 AliInfo(Form("Storage was not properly set!!!"));
76 abort();
77 }
78}
79
b45e5084 80
81//______________________________________________________________________________
c336835b 82Bool_t AliTagCreator::ReadGridCollection(TGridResult *fresult) {
cb1645b7 83 // Reads the entry of the TGridResult and creates the tags
84 Int_t nEntries = fresult->GetEntries();
b45e5084 85
cb1645b7 86 TString alienUrl;
b45e5084 87 const char *guid;
e16601cf 88 const char *md5;
89 const char *turl;
90 Long64_t size = -1;
91
b45e5084 92 Int_t counter = 0;
c336835b 93 for(Int_t i = 0; i < nEntries; i++) {
94 alienUrl = fresult->GetKey(i,"turl");
95 guid = fresult->GetKey(i,"guid");
96 if(fresult->GetKey(i,"size")) size = atol (fresult->GetKey(i,"size"));
97 md5 = fresult->GetKey(i,"md5");
98 turl = fresult->GetKey(i,"turl");
99 if(md5 && !strlen(guid)) md5 = 0;
100 if(guid && !strlen(guid)) guid = 0;
101
102 TFile *f = TFile::Open(alienUrl,"READ");
103 CreateTag(f,guid,md5,turl,size,counter);
104 f->Close();
105 delete f;
106 counter += 1;
107 }//grid result loop
108
109 return kTRUE;
110}
111
112//______________________________________________________________________________
113Bool_t AliTagCreator::ReadLocalCollection(const char *localpath) {
114 // Checks the different subdirs of the given local path and in the
115 // case where it finds an AliESDs.root file it creates the tags
116
117 void *dira = gSystem->OpenDirectory(localpath);
118 Char_t fPath[256];
119 const char * dirname = 0x0;
120 const char * filename = 0x0;
121 const char * pattern = "AliESDs.root";
122
123 Int_t counter = 0;
2005260b 124 while((dirname = gSystem->GetDirEntry(dira))) {
c336835b 125 sprintf(fPath,"%s/%s",localpath,dirname);
126 void *dirb = gSystem->OpenDirectory(fPath);
2005260b 127 while((filename = gSystem->GetDirEntry(dirb))) {
c336835b 128 if(strstr(filename,pattern)) {
129 TString fESDFileName;
130 fESDFileName = fPath;
131 fESDFileName += "/";
132 fESDFileName += pattern;
133 TFile *f = TFile::Open(fESDFileName,"READ");
134 CreateTag(f,fESDFileName,counter);
135 f->Close();
136 delete f;
137
138 counter += 1;
139 }//pattern check
140 }//child directory's entry loop
141 }//parent directory's entry loop
142
143 return kTRUE;
144}
145
146//______________________________________________________________________________
147Bool_t AliTagCreator::ReadCAFCollection(const char *filename) {
148 // Temporary solution for CAF: Takes as an input the ascii file that
149 // lists the ESDs stored in the SE of the CAF and creates the tags.
150
151 // Open the input stream
152 ifstream in;
153 in.open(filename);
154
155 Int_t counter = 0;
156 TString esdfile;
157 // Read the input list of files and add them to the chain
158 while(in.good()) {
159 in >> esdfile;
160 if (!esdfile.Contains("root")) continue; // protection
161 TFile *f = TFile::Open(esdfile,"READ");
162 CreateTag(f,esdfile,counter);
163 f->Close();
164 delete f;
165
166 counter += 1;
167 }
b45e5084 168
169 return kTRUE;
170}
171
c336835b 172
4fd84456 173//__________________________________________________________________________
c336835b 174Bool_t AliTagCreator::MergeTags() {
4fd84456 175 //Merges the tags and stores the merged tag file
176 //locally if fStorage=0 or in the grid if fStorage=1
177 AliInfo(Form("Merging tags....."));
178 TChain *fgChain = new TChain("T");
179
180 if(fStorage == 0) {
181 const char * tagPattern = "tag";
182 // Open the working directory
183 void * dirp = gSystem->OpenDirectory(gSystem->pwd());
184 const char * name = 0x0;
185 // Add all files matching *pattern* to the chain
186 while((name = gSystem->GetDirEntry(dirp))) {
187 if (strstr(name,tagPattern))
188 fgChain->Add(name);
189 }//directory loop
190 AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
191 }//local mode
192
193 else if(fStorage == 1) {
194 TString alienLocation = gGrid->Pwd();
195 alienLocation += fgridpath.Data();
196 alienLocation += "/";
197
198 TGridResult *tagresult = gGrid->Query(alienLocation,"*tag.root","","");
199 Int_t nEntries = tagresult->GetEntries();
200 for(Int_t i = 0; i < nEntries; i++) {
201 TString alienUrl = tagresult->GetKey(i,"turl");
202 fgChain->Add(alienUrl);
203 }//grid result loop
204 AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
205 }//grid mode
206
207 AliRunTag *tag = new AliRunTag;
208 AliEventTag *evTag = new AliEventTag;
209 fgChain->SetBranchAddress("AliTAG",&tag);
210
211 //Defining new tag objects
212 AliRunTag *newTag = new AliRunTag();
213 TTree ttag("T","A Tree with event tags");
214 TBranch * btag = ttag.Branch("AliTAG", &newTag);
215 btag->SetCompressionLevel(9);
216 for(Int_t iTagFiles = 0; iTagFiles < fgChain->GetEntries(); iTagFiles++) {
217 fgChain->GetEntry(iTagFiles);
218 newTag->SetRunId(tag->GetRunId());
219 const TClonesArray *tagList = tag->GetEventTags();
442646f2 220 for(Int_t j = 0; j < tagList->GetEntries(); j++) {
221 evTag = (AliEventTag *) tagList->At(j);
222 newTag->AddEventTag(*evTag);
223 }
224 ttag.Fill();
225 newTag->Clear();
226 }//tag file loop
227
228 TString localFileName = "Run"; localFileName += tag->GetRunId();
229 localFileName += ".Merged"; localFileName += ".ESD.tag.root";
230
231 TString filename = 0x0;
232
233 if(fStorage == 0) {
234 filename = localFileName.Data();
235 AliInfo(Form("Writing merged tags to local file: %s",filename.Data()));
236 }
237 else if(fStorage == 1) {
238 TString alienFileName = "/alien";
239 alienFileName += gGrid->Pwd();
240 alienFileName += fgridpath.Data();
241 alienFileName += "/";
242 alienFileName += localFileName;
243 alienFileName += "?se=";
244 alienFileName += fSE.Data();
245 filename = alienFileName.Data();
246 AliInfo(Form("Writing merged tags to grid file: %s",filename.Data()));
247 }
248
249 TFile* ftag = TFile::Open(filename, "recreate");
250 ftag->cd();
251 ttag.Write();
252 ftag->Close();
253
254 delete tag;
255 delete newTag;
256
257 return kTRUE;
258}
259
260//__________________________________________________________________________
261Bool_t AliTagCreator::MergeTags(TGridResult *result) {
262 //Merges the tags that are listed in the TGridResult
263 AliInfo(Form("Merging tags....."));
264 TChain *fgChain = new TChain("T");
265
266 Int_t nEntries = result->GetEntries();
267
268 TString alienUrl;
269 for(Int_t i = 0; i < nEntries; i++) {
270 alienUrl = result->GetKey(i,"turl");
271 fgChain->Add(alienUrl);
272 }
273 AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
274 AliRunTag *tag = new AliRunTag;
275 AliEventTag *evTag = new AliEventTag;
276 fgChain->SetBranchAddress("AliTAG",&tag);
277
278 //Defining new tag objects
279 AliRunTag *newTag = new AliRunTag();
280 TTree ttag("T","A Tree with event tags");
281 TBranch * btag = ttag.Branch("AliTAG", &newTag);
282 btag->SetCompressionLevel(9);
283 for(Int_t iTagFiles = 0; iTagFiles < fgChain->GetEntries(); iTagFiles++) {
284 fgChain->GetEntry(iTagFiles);
285 newTag->SetRunId(tag->GetRunId());
286 const TClonesArray *tagList = tag->GetEventTags();
4fd84456 287 for(Int_t j = 0; j < tagList->GetEntries(); j++) {
288 evTag = (AliEventTag *) tagList->At(j);
289 newTag->AddEventTag(*evTag);
290 }
291 ttag.Fill();
292 newTag->Clear();
293 }//tag file loop
294
295 TString localFileName = "Run"; localFileName += tag->GetRunId();
296 localFileName += ".Merged"; localFileName += ".ESD.tag.root";
297
4fd84456 298 TString filename = 0x0;
299
300 if(fStorage == 0) {
301 filename = localFileName.Data();
302 AliInfo(Form("Writing merged tags to local file: %s",filename.Data()));
303 }
c336835b 304 else if(fStorage == 1) {
305 TString alienFileName = "/alien";
306 alienFileName += gGrid->Pwd();
307 alienFileName += fgridpath.Data();
308 alienFileName += "/";
309 alienFileName += localFileName;
310 alienFileName += "?se=";
311 alienFileName += fSE.Data();
4fd84456 312 filename = alienFileName.Data();
313 AliInfo(Form("Writing merged tags to grid file: %s",filename.Data()));
314 }
315
316 TFile* ftag = TFile::Open(filename, "recreate");
317 ftag->cd();
318 ttag.Write();
319 ftag->Close();
320
321 delete tag;
4fd84456 322 delete newTag;
323
324 return kTRUE;
325}
b45e5084 326
327//_____________________________________________________________________________
c336835b 328void AliTagCreator::CreateTag(TFile* file, const char *guid, const char *md5, const char *turl, Long64_t size, Int_t Counter) {
329 //private method that creates tag files
330 TString fguid = guid;
331 TString fmd5 = md5;
332 TString fturl = turl;
2bdb9d38 333 /////////////
334 //muon code//
335 ////////////
56982dd1 336 Double_t fMUONMASS = 0.105658369;
2bdb9d38 337 //Variables
56982dd1 338 Double_t fX,fY,fZ ;
339 Double_t fThetaX, fThetaY, fPyz, fChisquare;
340 Double_t fPxRec, fPyRec, fPzRec, fEnergy;
341 Int_t fCharge;
342 TLorentzVector fEPvector;
343
344 Float_t fZVertexCut = 10.0;
345 Float_t fRhoVertexCut = 2.0;
346
347 Float_t fLowPtCut = 1.0;
348 Float_t fHighPtCut = 3.0;
349 Float_t fVeryHighPtCut = 10.0;
2bdb9d38 350 ////////////
351
352 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
353
cb1645b7 354 // Creates the tags for all the events in a given ESD file
c336835b 355 Bool_t fIsSim = kTRUE;
b45e5084 356 Int_t ntrack;
cb1645b7 357 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
358 Int_t nPos, nNeg, nNeutr;
359 Int_t nK0s, nNeutrons, nPi0s, nGamas;
360 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
361 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
362 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
363 Float_t maxPt = .0, meanPt = .0, totalP = .0;
d1a50cb5 364 Int_t fVertexflag;
ea865c85 365 Int_t iRunNumber = 0;
d1a50cb5 366 TString fVertexName;
b45e5084 367
368 AliRunTag *tag = new AliRunTag();
b45e5084 369 AliEventTag *evTag = new AliEventTag();
370 TTree ttag("T","A Tree with event tags");
4fd84456 371 TBranch * btag = ttag.Branch("AliTAG", &tag);
b45e5084 372 btag->SetCompressionLevel(9);
373
374 AliInfo(Form("Creating the tags......."));
375
376 Int_t firstEvent = 0,lastEvent = 0;
377 TTree *t = (TTree*) file->Get("esdTree");
378 TBranch * b = t->GetBranch("ESD");
379 AliESD *esd = 0;
380 b->SetAddress(&esd);
381
ea865c85 382 b->GetEntry(0);
383 Int_t iInitRunNumber = esd->GetRunNumber();
384
cb1645b7 385 Int_t iNumberOfEvents = b->GetEntries();
2bdb9d38 386 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
387 ntrack = 0;
388 nPos = 0;
389 nNeg = 0;
390 nNeutr =0;
391 nK0s = 0;
392 nNeutrons = 0;
393 nPi0s = 0;
394 nGamas = 0;
395 nProtons = 0;
396 nKaons = 0;
397 nPions = 0;
398 nMuons = 0;
399 nElectrons = 0;
400 nCh1GeV = 0;
401 nCh3GeV = 0;
402 nCh10GeV = 0;
403 nMu1GeV = 0;
404 nMu3GeV = 0;
405 nMu10GeV = 0;
406 nEl1GeV = 0;
407 nEl3GeV = 0;
408 nEl10GeV = 0;
409 maxPt = .0;
410 meanPt = .0;
411 totalP = .0;
d1a50cb5 412 fVertexflag = 1;
2bdb9d38 413
414 b->GetEntry(iEventNumber);
ea865c85 415 iRunNumber = esd->GetRunNumber();
416 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD!!!");
2bdb9d38 417 const AliESDVertex * vertexIn = esd->GetVertex();
d1a50cb5 418 fVertexName = vertexIn->GetName();
419 if(fVertexName == "default") fVertexflag = 0;
420
2bdb9d38 421 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
422 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
c336835b 423 if(esdTrack->GetLabel() != 0) fIsSim = kTRUE;
424 else if(esdTrack->GetLabel() == 0) fIsSim = kFALSE;
2bdb9d38 425 UInt_t status = esdTrack->GetStatus();
b45e5084 426
2bdb9d38 427 //select only tracks with ITS refit
428 if ((status&AliESDtrack::kITSrefit)==0) continue;
429 //select only tracks with TPC refit
430 if ((status&AliESDtrack::kTPCrefit)==0) continue;
b45e5084 431
2bdb9d38 432 //select only tracks with the "combined PID"
433 if ((status&AliESDtrack::kESDpid)==0) continue;
434 Double_t p[3];
435 esdTrack->GetPxPyPz(p);
436 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
437 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
438 totalP += momentum;
439 meanPt += fPt;
440 if(fPt > maxPt) maxPt = fPt;
b45e5084 441
2bdb9d38 442 if(esdTrack->GetSign() > 0) {
443 nPos++;
56982dd1 444 if(fPt > fLowPtCut) nCh1GeV++;
445 if(fPt > fHighPtCut) nCh3GeV++;
446 if(fPt > fVeryHighPtCut) nCh10GeV++;
2bdb9d38 447 }
448 if(esdTrack->GetSign() < 0) {
449 nNeg++;
56982dd1 450 if(fPt > fLowPtCut) nCh1GeV++;
451 if(fPt > fHighPtCut) nCh3GeV++;
452 if(fPt > fVeryHighPtCut) nCh10GeV++;
2bdb9d38 453 }
454 if(esdTrack->GetSign() == 0) nNeutr++;
b45e5084 455
2bdb9d38 456 //PID
457 Double_t prob[5];
458 esdTrack->GetESDpid(prob);
b45e5084 459
2bdb9d38 460 Double_t rcc = 0.0;
461 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
462 if(rcc == 0.0) continue;
463 //Bayes' formula
464 Double_t w[5];
465 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
b45e5084 466
2bdb9d38 467 //protons
468 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
469 //kaons
470 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
471 //pions
472 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
473 //electrons
474 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
475 nElectrons++;
56982dd1 476 if(fPt > fLowPtCut) nEl1GeV++;
477 if(fPt > fHighPtCut) nEl3GeV++;
478 if(fPt > fVeryHighPtCut) nEl10GeV++;
2bdb9d38 479 }
480 ntrack++;
481 }//esd track loop
482
483 /////////////
484 //muon code//
485 ////////////
486 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
487 // loop over all reconstructed tracks (also first track of combination)
488 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
489 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
490 if (muonTrack == 0x0) continue;
b45e5084 491
2bdb9d38 492 // Coordinates at vertex
56982dd1 493 fZ = muonTrack->GetZ();
494 fY = muonTrack->GetBendingCoor();
495 fX = muonTrack->GetNonBendingCoor();
b45e5084 496
56982dd1 497 fThetaX = muonTrack->GetThetaX();
498 fThetaY = muonTrack->GetThetaY();
b45e5084 499
56982dd1 500 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
501 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
502 fPxRec = fPzRec * TMath::Tan(fThetaX);
503 fPyRec = fPzRec * TMath::Tan(fThetaY);
504 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
b45e5084 505
2bdb9d38 506 //ChiSquare of the track if needed
56982dd1 507 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
508 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
509 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
b45e5084 510
2bdb9d38 511 // total number of muons inside a vertex cut
56982dd1 512 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
2bdb9d38 513 nMuons++;
56982dd1 514 if(fEPvector.Pt() > fLowPtCut) {
2bdb9d38 515 nMu1GeV++;
56982dd1 516 if(fEPvector.Pt() > fHighPtCut) {
2bdb9d38 517 nMu3GeV++;
56982dd1 518 if (fEPvector.Pt() > fVeryHighPtCut) {
2bdb9d38 519 nMu10GeV++;
520 }
521 }
522 }
523 }
524 }//muon track loop
525
526 // Fill the event tags
c336835b 527 if(ntrack != 0) meanPt = meanPt/ntrack;
2bdb9d38 528
529 evTag->SetEventId(iEventNumber+1);
c336835b 530 evTag->SetGUID(fguid);
531 evTag->SetMD5(fmd5);
532 evTag->SetTURL(fturl);
2bdb9d38 533 evTag->SetSize(size);
534 evTag->SetVertexX(vertexIn->GetXv());
535 evTag->SetVertexY(vertexIn->GetYv());
536 evTag->SetVertexZ(vertexIn->GetZv());
d1a50cb5 537 evTag->SetVertexZError(vertexIn->GetZRes());
538 evTag->SetVertexFlag(fVertexflag);
2bdb9d38 539
540 evTag->SetT0VertexZ(esd->GetT0zVertex());
541
8bd8ac26 542 evTag->SetTriggerMask(esd->GetTriggerMask());
543 evTag->SetTriggerCluster(esd->GetTriggerCluster());
2bdb9d38 544
32a5cab4 545 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
546 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
2bdb9d38 547 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
32a5cab4 548 evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
549 evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
2bdb9d38 550 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
551
552
553 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
554 evTag->SetNumOfPosTracks(nPos);
555 evTag->SetNumOfNegTracks(nNeg);
556 evTag->SetNumOfNeutrTracks(nNeutr);
557
558 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
559 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
560 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
561 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
562
563 evTag->SetNumOfProtons(nProtons);
564 evTag->SetNumOfKaons(nKaons);
565 evTag->SetNumOfPions(nPions);
566 evTag->SetNumOfMuons(nMuons);
567 evTag->SetNumOfElectrons(nElectrons);
568 evTag->SetNumOfPhotons(nGamas);
569 evTag->SetNumOfPi0s(nPi0s);
570 evTag->SetNumOfNeutrons(nNeutrons);
571 evTag->SetNumOfKaon0s(nK0s);
572
573 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
574 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
575 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
576 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
577 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
578 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
579 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
580 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
581 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
582
85c60a8e 583 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
584 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
2bdb9d38 585
586 evTag->SetTotalMomentum(totalP);
587 evTag->SetMeanPt(meanPt);
588 evTag->SetMaxPt(maxPt);
589
ea865c85 590 tag->SetRunId(iInitRunNumber);
c336835b 591 if(fIsSim) tag->SetDataType(0);
592 else tag->SetDataType(1);
2bdb9d38 593 tag->AddEventTag(*evTag);
594 }//event loop
cb1645b7 595 lastEvent = iNumberOfEvents;
b45e5084 596
597 t->Delete("");
2bdb9d38 598
b45e5084 599 ttag.Fill();
600 tag->Clear();
601
cb1645b7 602 TString localFileName = "Run"; localFileName += tag->GetRunId();
603 localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
604 localFileName += ".ESD.tag.root";
b45e5084 605
cb1645b7 606 TString fileName;
28afeb2e 607
c336835b 608 if(fStorage == 0) {
609 fileName = localFileName.Data();
610 AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
611 }
612 else if(fStorage == 1) {
613 TString alienLocation = "/alien";
614 alienLocation += gGrid->Pwd();
615 alienLocation += fgridpath.Data();
616 alienLocation += "/";
617 alienLocation += localFileName;
618 alienLocation += "?se=";
619 alienLocation += fSE.Data();
620 fileName = alienLocation.Data();
621 AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
622 }
b45e5084 623
cb1645b7 624 TFile* ftag = TFile::Open(fileName, "recreate");
b45e5084 625 ftag->cd();
626 ttag.Write();
627 ftag->Close();
628
629 delete ftag;
630 delete esd;
631
632 delete tag;
b45e5084 633}
634
c336835b 635//_____________________________________________________________________________
636void AliTagCreator::CreateTag(TFile* file, const char *filepath, Int_t Counter) {
637 //private method that creates tag files
638 /////////////
639 //muon code//
640 ////////////
641 Double_t fMUONMASS = 0.105658369;
642 //Variables
643 Double_t fX,fY,fZ ;
644 Double_t fThetaX, fThetaY, fPyz, fChisquare;
645 Double_t fPxRec, fPyRec, fPzRec, fEnergy;
646 Int_t fCharge;
647 TLorentzVector fEPvector;
648
649 Float_t fZVertexCut = 10.0;
650 Float_t fRhoVertexCut = 2.0;
651
652 Float_t fLowPtCut = 1.0;
653 Float_t fHighPtCut = 3.0;
654 Float_t fVeryHighPtCut = 10.0;
655 ////////////
656
657 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
658
659 // Creates the tags for all the events in a given ESD file
660 Bool_t fIsSim = kTRUE;
661 Int_t ntrack;
662 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
663 Int_t nPos, nNeg, nNeutr;
664 Int_t nK0s, nNeutrons, nPi0s, nGamas;
665 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
666 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
667 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
668 Float_t maxPt = .0, meanPt = .0, totalP = .0;
669 Int_t fVertexflag;
670 Int_t iRunNumber = 0;
671 TString fVertexName;
672
673 AliRunTag *tag = new AliRunTag();
674 AliEventTag *evTag = new AliEventTag();
675 TTree ttag("T","A Tree with event tags");
676 TBranch * btag = ttag.Branch("AliTAG", &tag);
677 btag->SetCompressionLevel(9);
678
679 AliInfo(Form("Creating the tags......."));
680
681 Int_t firstEvent = 0,lastEvent = 0;
682
683 TTree *t = (TTree*) file->Get("esdTree");
684 TBranch * b = t->GetBranch("ESD");
685 AliESD *esd = 0;
686 b->SetAddress(&esd);
687
688 b->GetEntry(0);
689 Int_t iInitRunNumber = esd->GetRunNumber();
690
691 Int_t iNumberOfEvents = b->GetEntries();
692 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
693 ntrack = 0;
694 nPos = 0;
695 nNeg = 0;
696 nNeutr =0;
697 nK0s = 0;
698 nNeutrons = 0;
699 nPi0s = 0;
700 nGamas = 0;
701 nProtons = 0;
702 nKaons = 0;
703 nPions = 0;
704 nMuons = 0;
705 nElectrons = 0;
706 nCh1GeV = 0;
707 nCh3GeV = 0;
708 nCh10GeV = 0;
709 nMu1GeV = 0;
710 nMu3GeV = 0;
711 nMu10GeV = 0;
712 nEl1GeV = 0;
713 nEl3GeV = 0;
714 nEl10GeV = 0;
715 maxPt = .0;
716 meanPt = .0;
717 totalP = .0;
718 fVertexflag = 1;
719
720 b->GetEntry(iEventNumber);
721 iRunNumber = esd->GetRunNumber();
722 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD!!!");
723 const AliESDVertex * vertexIn = esd->GetVertex();
724 fVertexName = vertexIn->GetName();
725 if(fVertexName == "default") fVertexflag = 0;
726
727 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
728 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
729 if(esdTrack->GetLabel() != 0) fIsSim = kTRUE;
730 else if(esdTrack->GetLabel() == 0) fIsSim = kFALSE;
731 UInt_t status = esdTrack->GetStatus();
732
733 //select only tracks with ITS refit
734 if ((status&AliESDtrack::kITSrefit)==0) continue;
735 //select only tracks with TPC refit
736 if ((status&AliESDtrack::kTPCrefit)==0) continue;
737
738 //select only tracks with the "combined PID"
739 if ((status&AliESDtrack::kESDpid)==0) continue;
740 Double_t p[3];
741 esdTrack->GetPxPyPz(p);
742 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
743 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
744 totalP += momentum;
745 meanPt += fPt;
746 if(fPt > maxPt) maxPt = fPt;
747
748 if(esdTrack->GetSign() > 0) {
749 nPos++;
750 if(fPt > fLowPtCut) nCh1GeV++;
751 if(fPt > fHighPtCut) nCh3GeV++;
752 if(fPt > fVeryHighPtCut) nCh10GeV++;
753 }
754 if(esdTrack->GetSign() < 0) {
755 nNeg++;
756 if(fPt > fLowPtCut) nCh1GeV++;
757 if(fPt > fHighPtCut) nCh3GeV++;
758 if(fPt > fVeryHighPtCut) nCh10GeV++;
759 }
760 if(esdTrack->GetSign() == 0) nNeutr++;
761
762 //PID
763 Double_t prob[5];
764 esdTrack->GetESDpid(prob);
765
766 Double_t rcc = 0.0;
767 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
768 if(rcc == 0.0) continue;
769 //Bayes' formula
770 Double_t w[5];
771 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
772
773 //protons
774 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
775 //kaons
776 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
777 //pions
778 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
779 //electrons
780 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
781 nElectrons++;
782 if(fPt > fLowPtCut) nEl1GeV++;
783 if(fPt > fHighPtCut) nEl3GeV++;
784 if(fPt > fVeryHighPtCut) nEl10GeV++;
785 }
786 ntrack++;
787 }//esd track loop
788
789 /////////////
790 //muon code//
791 ////////////
792 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
793 // loop over all reconstructed tracks (also first track of combination)
794 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
795 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
796 if (muonTrack == 0x0) continue;
797
798 // Coordinates at vertex
799 fZ = muonTrack->GetZ();
800 fY = muonTrack->GetBendingCoor();
801 fX = muonTrack->GetNonBendingCoor();
802
803 fThetaX = muonTrack->GetThetaX();
804 fThetaY = muonTrack->GetThetaY();
805
806 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
807 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
808 fPxRec = fPzRec * TMath::Tan(fThetaX);
809 fPyRec = fPzRec * TMath::Tan(fThetaY);
810 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
811
812 //ChiSquare of the track if needed
813 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
814 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
815 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
816
817 // total number of muons inside a vertex cut
818 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
819 nMuons++;
820 if(fEPvector.Pt() > fLowPtCut) {
821 nMu1GeV++;
822 if(fEPvector.Pt() > fHighPtCut) {
823 nMu3GeV++;
824 if (fEPvector.Pt() > fVeryHighPtCut) {
825 nMu10GeV++;
826 }
827 }
828 }
829 }
830 }//muon track loop
831
832 // Fill the event tags
833 if(ntrack != 0) meanPt = meanPt/ntrack;
834
835 evTag->SetEventId(iEventNumber+1);
836 evTag->SetPath(filepath);
837
838 evTag->SetVertexX(vertexIn->GetXv());
839 evTag->SetVertexY(vertexIn->GetYv());
840 evTag->SetVertexZ(vertexIn->GetZv());
841 evTag->SetVertexZError(vertexIn->GetZRes());
842 evTag->SetVertexFlag(fVertexflag);
843
844 evTag->SetT0VertexZ(esd->GetT0zVertex());
845
846 evTag->SetTriggerMask(esd->GetTriggerMask());
847 evTag->SetTriggerCluster(esd->GetTriggerCluster());
848
849 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
850 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
851 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
852 evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
853 evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
854 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
855
856
857 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
858 evTag->SetNumOfPosTracks(nPos);
859 evTag->SetNumOfNegTracks(nNeg);
860 evTag->SetNumOfNeutrTracks(nNeutr);
861
862 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
863 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
864 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
865 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
866
867 evTag->SetNumOfProtons(nProtons);
868 evTag->SetNumOfKaons(nKaons);
869 evTag->SetNumOfPions(nPions);
870 evTag->SetNumOfMuons(nMuons);
871 evTag->SetNumOfElectrons(nElectrons);
872 evTag->SetNumOfPhotons(nGamas);
873 evTag->SetNumOfPi0s(nPi0s);
874 evTag->SetNumOfNeutrons(nNeutrons);
875 evTag->SetNumOfKaon0s(nK0s);
876
877 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
878 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
879 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
880 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
881 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
882 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
883 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
884 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
885 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
886
887 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
888 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
889
890 evTag->SetTotalMomentum(totalP);
891 evTag->SetMeanPt(meanPt);
892 evTag->SetMaxPt(maxPt);
893
894 tag->SetRunId(iInitRunNumber);
895 if(fIsSim) tag->SetDataType(0);
896 else tag->SetDataType(1);
897 tag->AddEventTag(*evTag);
898 }//event loop
899 lastEvent = iNumberOfEvents;
900
901 t->Delete("");
902
903 ttag.Fill();
904 tag->Clear();
905
906 TString localFileName = "Run"; localFileName += tag->GetRunId();
907 localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
908 localFileName += ".ESD.tag.root";
909
910 TString fileName;
911
912 if(fStorage == 0) {
913 fileName = localFileName.Data();
914 AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
915 }
916 else if(fStorage == 1) {
917 TString alienLocation = "/alien";
918 alienLocation += gGrid->Pwd();
919 alienLocation += fgridpath.Data();
920 alienLocation += "/";
921 alienLocation += localFileName;
922 alienLocation += "?se=";
923 alienLocation += fSE.Data();
924 fileName = alienLocation.Data();
925 AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
926 }
927
928 TFile* ftag = TFile::Open(fileName, "recreate");
929 ftag->cd();
930 ttag.Write();
931 ftag->Close();
932
933 delete ftag;
934 delete esd;
935
936 delete tag;
937}