]> git.uio.no Git - u/mrichter/AliRoot.git/blame - STEER/AliESDTagCreator.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Removing osolete variable. Adding title to trigger configuration obj also
[u/mrichter/AliRoot.git] / STEER / AliESDTagCreator.cxx
CommitLineData
08e1a23e 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// AliESDTagCreator 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
23#include <Riostream.h>
24#include <TFile.h>
25#include <TString.h>
26#include <TTree.h>
27#include <TSystem.h>
28#include <TChain.h>
103d3fba 29#include <TList.h>
30#include <TObjString.h>
08e1a23e 31#include <TLorentzVector.h>
6b6e4472 32#include <TMap.h>
08e1a23e 33
34//ROOT-AliEn
35#include <TGrid.h>
36#include <TGridResult.h>
37
38//AliRoot
39#include "AliRunTag.h"
40#include "AliEventTag.h"
41#include "AliESD.h"
42#include "AliESDEvent.h"
43#include "AliESDVertex.h"
44#include "AliLog.h"
45
46#include "AliESDTagCreator.h"
47
48
49ClassImp(AliESDTagCreator)
50
51
52//______________________________________________________________________________
53 AliESDTagCreator::AliESDTagCreator() :
103d3fba 54 AliTagCreator(),
55 fChain(new TChain("esdTree")), fGUIDList(new TList()),
517aef4e 56 fMD5List(new TList()), fTURLList(new TList()), fBranches(""),
103d3fba 57 meminfo(new MemInfo_t) {
08e1a23e 58 //==============Default constructor for a AliESDTagCreator================
59}
60
61//______________________________________________________________________________
62AliESDTagCreator::~AliESDTagCreator() {
63//================Default destructor for a AliESDTagCreator===================
103d3fba 64 delete fChain;
65 delete fGUIDList;
66 delete fMD5List;
67 delete fTURLList;
68 delete meminfo;
08e1a23e 69}
70
71//______________________________________________________________________________
72Bool_t AliESDTagCreator::ReadGridCollection(TGridResult *fresult) {
73 // Reads the entry of the TGridResult and creates the tags
74 Int_t nEntries = fresult->GetEntries();
75
76 TString alienUrl;
103d3fba 77 const char* guid;
78 const char* md5;
79 const char* turl;
08e1a23e 80 Long64_t size = -1;
81
82 Int_t counter = 0;
83 for(Int_t i = 0; i < nEntries; i++) {
84 alienUrl = fresult->GetKey(i,"turl");
85 guid = fresult->GetKey(i,"guid");
86 if(fresult->GetKey(i,"size")) size = atol (fresult->GetKey(i,"size"));
87 md5 = fresult->GetKey(i,"md5");
88 turl = fresult->GetKey(i,"turl");
89 if(md5 && !strlen(guid)) md5 = 0;
90 if(guid && !strlen(guid)) guid = 0;
103d3fba 91
92 fChain->Add(alienUrl);
93 //fGUIDList->Add(new TObjString(guid));
94 //fMD5List->Add(new TObjString(md5));
95 //fTURLList->Add(new TObjString(turl));
96
97 //TFile *f = TFile::Open(alienUrl,"READ");
98 //CreateTag(f,guid,md5,turl,size,counter);
99 //f->Close();
100 //delete f;
08e1a23e 101 counter += 1;
102 }//grid result loop
103d3fba 103
104 AliInfo(Form("ESD chain created......."));
105 AliInfo(Form("Chain entries: %d",fChain->GetEntries()));
2856e38b 106 // Switch of branches on user request
107 SwitchOffBranches();
103d3fba 108 CreateTag(fChain,"grid");
109
08e1a23e 110 return kTRUE;
111}
112
113//______________________________________________________________________________
114Bool_t AliESDTagCreator::ReadLocalCollection(const char *localpath) {
115 // Checks the different subdirs of the given local path and in the
116 // case where it finds an AliESDs.root file it creates the tags
117
118 void *dira = gSystem->OpenDirectory(localpath);
119 Char_t fPath[256];
120 const char * dirname = 0x0;
121 const char * filename = 0x0;
122 const char * pattern = "AliESDs.root";
123
124 Int_t counter = 0;
125 while((dirname = gSystem->GetDirEntry(dira))) {
126 sprintf(fPath,"%s/%s",localpath,dirname);
127 void *dirb = gSystem->OpenDirectory(fPath);
128 while((filename = gSystem->GetDirEntry(dirb))) {
129 if(strstr(filename,pattern)) {
130 TString fESDFileName;
131 fESDFileName = fPath;
132 fESDFileName += "/";
133 fESDFileName += pattern;
103d3fba 134
135 fChain->Add(fESDFileName);
136
137 //TFile *f = TFile::Open(fESDFileName,"READ");
138 //CreateTag(f,fESDFileName,counter);
139 //f->Close();
140 //delete f;
08e1a23e 141
142 counter += 1;
143 }//pattern check
144 }//child directory's entry loop
145 }//parent directory's entry loop
146
103d3fba 147 AliInfo(Form("ESD chain created......."));
148 AliInfo(Form("Chain entries: %d",fChain->GetEntries()));
2856e38b 149 // Switch of branches on user request
150 SwitchOffBranches();
103d3fba 151 CreateTag(fChain,"local");
152
08e1a23e 153 return kTRUE;
154}
155
156//______________________________________________________________________________
157Bool_t AliESDTagCreator::ReadCAFCollection(const char *filename) {
158 // Temporary solution for CAF: Takes as an input the ascii file that
159 // lists the ESDs stored in the SE of the CAF and creates the tags.
160
161 // Open the input stream
162 ifstream in;
163 in.open(filename);
164
165 Int_t counter = 0;
166 TString esdfile;
167 // Read the input list of files and add them to the chain
168 while(in.good()) {
169 in >> esdfile;
170 if (!esdfile.Contains("root")) continue; // protection
103d3fba 171
172 fChain->Add(esdfile);
173
174 //TFile *f = TFile::Open(esdfile,"READ");
175 //CreateTag(f,esdfile,counter);
176 //f->Close();
177 //delete f;
08e1a23e 178
179 counter += 1;
180 }
181
103d3fba 182 AliInfo(Form("ESD chain created......."));
183 AliInfo(Form("Chain entries: %d",fChain->GetEntries()));
2856e38b 184 // Switch of branches on user request
185 SwitchOffBranches();
103d3fba 186 CreateTag(fChain,"proof");
187
08e1a23e 188 return kTRUE;
189}
190
103d3fba 191//_____________________________________________________________________________
192void AliESDTagCreator::CreateTag(TChain* fChain, const char *type) {
193 //private method that creates tag files
194 TString fSession = type;
103d3fba 195 TString fguid, fmd5, fturl;
196 TString fTempGuid = 0;
197
198 /////////////
199 //muon code//
200 ////////////
201 Double_t fMUONMASS = 0.105658369;
202 //Variables
203 Double_t fX,fY,fZ ;
204 Double_t fThetaX, fThetaY, fPyz, fChisquare;
205 Double_t fPxRec, fPyRec, fPzRec, fEnergy;
206 Int_t fCharge;
207 TLorentzVector fEPvector;
208
209 Float_t fZVertexCut = 10.0;
210 Float_t fRhoVertexCut = 2.0;
211
212 Float_t fLowPtCut = 1.0;
213 Float_t fHighPtCut = 3.0;
214 Float_t fVeryHighPtCut = 10.0;
215 ////////////
216
217 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
218
219 // Creates the tags for all the events in a given ESD file
220 Bool_t fIsSim = kTRUE;
221 Int_t ntrack;
222 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
223 Int_t nPos, nNeg, nNeutr;
224 Int_t nK0s, nNeutrons, nPi0s, nGamas;
225 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
226 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
227 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
228 Float_t maxPt = .0, meanPt = .0, totalP = .0;
229 Int_t fVertexflag;
230 Int_t iRunNumber = 0;
231 TString fVertexName;
232
103d3fba 233 //gSystem->GetMemInfo(meminfo);
234 //AliInfo(Form("After the tag initialization - Memory used: %d MB",meminfo->fMemUsed));
235 //Int_t tempmem = meminfo->fMemUsed;
236
237 AliInfo(Form("Creating the ESD tags......."));
238
239 Int_t firstEvent = 0,lastEvent = 0;
240 AliESDEvent *esd = new AliESDEvent();
241 esd->ReadFromTree(fChain);
fcc6b05f 242 AliESD *esdold = 0x0;
103d3fba 243
244 //gSystem->GetMemInfo(meminfo);
245 //AliInfo(Form("After the esd initialization - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
246 //tempmem = meminfo->fMemUsed;
247
248 Int_t iInitRunNumber = -1;
fcc6b05f 249 fChain->GetEntry(0);
250 TFile *f = fChain->GetFile();
251 fTempGuid = f->GetUUID().AsString();
252
253 TString localFileName = "Run"; localFileName += esd->GetRunNumber();
254 localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += fChain->GetEntries(); //localFileName += "."; localFileName += Counter;
255 localFileName += ".ESD.tag.root";
256
257 TString fileName;
258
259 if(fStorage == 0) {
260 fileName = localFileName.Data();
261 AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
262 }
263 else if(fStorage == 1) {
264 TString alienLocation = "/alien";
265 alienLocation += gGrid->Pwd();
266 alienLocation += fgridpath.Data();
267 alienLocation += "/";
268 alienLocation += localFileName;
269 alienLocation += "?se=";
270 alienLocation += fSE.Data();
271 fileName = alienLocation.Data();
272 AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
273 }
274
275 TFile* ftag = TFile::Open(fileName, "recreate");
276
277 AliRunTag *tag = new AliRunTag();
278 AliEventTag *evTag = new AliEventTag();
279 TTree ttag("T","A Tree with event tags");
280 TBranch * btag = ttag.Branch("AliTAG", &tag);
281 btag->SetCompressionLevel(9);
282
103d3fba 283 for(Int_t iEventNumber = 0; iEventNumber < fChain->GetEntries(); iEventNumber++) {
284 ntrack = 0; nPos = 0; nNeg = 0; nNeutr =0;
285 nK0s = 0; nNeutrons = 0; nPi0s = 0;
286 nGamas = 0; nProtons = 0; nKaons = 0;
287 nPions = 0; nMuons = 0; nElectrons = 0;
288 nCh1GeV = 0; nCh3GeV = 0; nCh10GeV = 0;
289 nMu1GeV = 0; nMu3GeV = 0; nMu10GeV = 0;
290 nEl1GeV = 0; nEl3GeV = 0; nEl10GeV = 0;
291 maxPt = .0; meanPt = .0; totalP = .0;
292 fVertexflag = 1;
293
fcc6b05f 294 fChain->GetEntry(iEventNumber);
295 esdold = esd->GetAliESDOld();
296 if(esdold) esd->CopyFromOldESD();
297
103d3fba 298 TFile *f = fChain->GetFile();
299 const TUrl *url = f->GetEndpointUrl();
300 fguid = f->GetUUID().AsString();
301 if(fSession == "grid") {
302 TString fturltemp = "alien://"; fturltemp += url->GetFile();
303 fturl = fturltemp(0,fturltemp.Index(".root",5,0,TString::kExact)+5);
304 }
305 else fturl = url->GetFile();
306
307 if(iEventNumber == 0) iInitRunNumber = esd->GetRunNumber();
308 iRunNumber = esd->GetRunNumber();
309 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD - You are trying to merge different runs!!!");
310
311 const AliESDVertex * vertexIn = esd->GetVertex();
312 fVertexName = vertexIn->GetName();
313 if(fVertexName == "default") fVertexflag = 0;
314
315 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
316 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
317 if(esdTrack->GetLabel() != 0) fIsSim = kTRUE;
318 else if(esdTrack->GetLabel() == 0) fIsSim = kFALSE;
319 UInt_t status = esdTrack->GetStatus();
320
321 //select only tracks with ITS refit
322 if ((status&AliESDtrack::kITSrefit)==0) continue;
323 //select only tracks with TPC refit
324 if ((status&AliESDtrack::kTPCrefit)==0) continue;
325
326 //select only tracks with the "combined PID"
327 if ((status&AliESDtrack::kESDpid)==0) continue;
328 Double_t p[3];
329 esdTrack->GetPxPyPz(p);
330 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
331 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
332 totalP += momentum;
333 meanPt += fPt;
334 if(fPt > maxPt) maxPt = fPt;
335
336 if(esdTrack->GetSign() > 0) {
337 nPos++;
338 if(fPt > fLowPtCut) nCh1GeV++;
339 if(fPt > fHighPtCut) nCh3GeV++;
340 if(fPt > fVeryHighPtCut) nCh10GeV++;
341 }
342 if(esdTrack->GetSign() < 0) {
343 nNeg++;
344 if(fPt > fLowPtCut) nCh1GeV++;
345 if(fPt > fHighPtCut) nCh3GeV++;
346 if(fPt > fVeryHighPtCut) nCh10GeV++;
347 }
348 if(esdTrack->GetSign() == 0) nNeutr++;
349
350 //PID
351 Double_t prob[5];
352 esdTrack->GetESDpid(prob);
353
354 Double_t rcc = 0.0;
355 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
356 if(rcc == 0.0) continue;
357 //Bayes' formula
358 Double_t w[5];
359 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
360
361 //protons
362 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
363 //kaons
364 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
365 //pions
366 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
367 //electrons
368 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
369 nElectrons++;
370 if(fPt > fLowPtCut) nEl1GeV++;
371 if(fPt > fHighPtCut) nEl3GeV++;
372 if(fPt > fVeryHighPtCut) nEl10GeV++;
373 }
374 ntrack++;
375 }//esd track loop
376
377 /////////////
378 //muon code//
379 ////////////
380 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
381 // loop over all reconstructed tracks (also first track of combination)
382 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
383 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
384 if (muonTrack == 0x0) continue;
385
386 // Coordinates at vertex
387 fZ = muonTrack->GetZ();
388 fY = muonTrack->GetBendingCoor();
389 fX = muonTrack->GetNonBendingCoor();
390
391 fThetaX = muonTrack->GetThetaX();
392 fThetaY = muonTrack->GetThetaY();
393
394 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
395 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
396 fPxRec = fPzRec * TMath::Tan(fThetaX);
397 fPyRec = fPzRec * TMath::Tan(fThetaY);
398 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
399
400 //ChiSquare of the track if needed
401 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
402 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
403 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
404
405 // total number of muons inside a vertex cut
406 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
407 nMuons++;
408 if(fEPvector.Pt() > fLowPtCut) {
409 nMu1GeV++;
410 if(fEPvector.Pt() > fHighPtCut) {
411 nMu3GeV++;
412 if (fEPvector.Pt() > fVeryHighPtCut) {
413 nMu10GeV++;
414 }
415 }
416 }
417 }
418 }//muon track loop
419
420 // Fill the event tags
421 if(ntrack != 0) meanPt = meanPt/ntrack;
422
fcc6b05f 423 //AliInfo(Form("====================================="));
424 //AliInfo(Form("URL: %s - GUID: %s",fturl.Data(),fguid.Data()));
425 //AliInfo(Form("====================================="));
103d3fba 426
427 evTag->SetEventId(iEventNumber+1);
428 evTag->SetGUID(fguid);
429 if(fSession == "grid") {
430 evTag->SetMD5(0);
431 evTag->SetTURL(fturl);
432 evTag->SetSize(0);
433 }
434 else evTag->SetPath(fturl);
435
436 evTag->SetVertexX(vertexIn->GetXv());
437 evTag->SetVertexY(vertexIn->GetYv());
438 evTag->SetVertexZ(vertexIn->GetZv());
439 evTag->SetVertexZError(vertexIn->GetZRes());
440 evTag->SetVertexFlag(fVertexflag);
441
442 evTag->SetT0VertexZ(esd->GetT0zVertex());
443
444 evTag->SetTriggerMask(esd->GetTriggerMask());
445 evTag->SetTriggerCluster(esd->GetTriggerCluster());
446
447 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
448 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
a85132e7 449 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy(0),esd->GetZDCEMEnergy(1));
103d3fba 450 evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
451 evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
452 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
453
454
455 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
456 evTag->SetNumOfPosTracks(nPos);
457 evTag->SetNumOfNegTracks(nNeg);
458 evTag->SetNumOfNeutrTracks(nNeutr);
459
460 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
461 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
462 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
463 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
464
465 evTag->SetNumOfProtons(nProtons);
466 evTag->SetNumOfKaons(nKaons);
467 evTag->SetNumOfPions(nPions);
468 evTag->SetNumOfMuons(nMuons);
469 evTag->SetNumOfElectrons(nElectrons);
470 evTag->SetNumOfPhotons(nGamas);
471 evTag->SetNumOfPi0s(nPi0s);
472 evTag->SetNumOfNeutrons(nNeutrons);
473 evTag->SetNumOfKaon0s(nK0s);
474
475 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
476 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
477 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
478 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
479 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
480 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
481 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
482 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
483 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
484
485 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
486 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
487
488 evTag->SetTotalMomentum(totalP);
489 evTag->SetMeanPt(meanPt);
490 evTag->SetMaxPt(maxPt);
491
492 tag->SetRunId(iInitRunNumber);
493 if(fIsSim) tag->SetDataType(0);
494 else tag->SetDataType(1);
103d3fba 495
496 if(fguid != fTempGuid) {
497 fTempGuid = fguid;
498 ttag.Fill();
499 tag->Clear("");
500 }
fcc6b05f 501 tag->AddEventTag(*evTag);
502 if(iEventNumber+1 == fChain->GetEntries()) {
503 //AliInfo(Form("File: %s",fturl.Data()));
504 ttag.Fill();
505 tag->Clear("");
506 }
103d3fba 507 }//event loop
508 lastEvent = fChain->GetEntries();
509
510 //gSystem->GetMemInfo(meminfo);
511 //AliInfo(Form("After the event and track loop - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
512 //tempmem = meminfo->fMemUsed;
513
514 //fChain->Delete("");
515
516 //gSystem->GetMemInfo(meminfo);
517 //AliInfo(Form("After the t->Delete - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
518 //tempmem = meminfo->fMemUsed;
519
103d3fba 520 ftag->cd();
103d3fba 521 tag->Clear();
522 ttag.Write();
523 ftag->Close();
524
525 //gSystem->GetMemInfo(meminfo);
526 //AliInfo(Form("After the file closing - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
527 //tempmem = meminfo->fMemUsed;
528
103d3fba 529 delete esd;
103d3fba 530 delete tag;
fcc6b05f 531
103d3fba 532 //gSystem->GetMemInfo(meminfo);
533 //AliInfo(Form("After the delete objects - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
534}
535
08e1a23e 536//_____________________________________________________________________________
537void AliESDTagCreator::CreateTag(TFile* file, const char *guid, const char *md5, const char *turl, Long64_t size, Int_t Counter) {
538 //private method that creates tag files
539 TString fguid = guid;
540 TString fmd5 = md5;
541 TString fturl = turl;
542 /////////////
543 //muon code//
544 ////////////
545 Double_t fMUONMASS = 0.105658369;
546 //Variables
547 Double_t fX,fY,fZ ;
548 Double_t fThetaX, fThetaY, fPyz, fChisquare;
549 Double_t fPxRec, fPyRec, fPzRec, fEnergy;
550 Int_t fCharge;
551 TLorentzVector fEPvector;
552
553 Float_t fZVertexCut = 10.0;
554 Float_t fRhoVertexCut = 2.0;
555
556 Float_t fLowPtCut = 1.0;
557 Float_t fHighPtCut = 3.0;
558 Float_t fVeryHighPtCut = 10.0;
559 ////////////
560
561 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
562
563 // Creates the tags for all the events in a given ESD file
564 Bool_t fIsSim = kTRUE;
565 Int_t ntrack;
566 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
567 Int_t nPos, nNeg, nNeutr;
568 Int_t nK0s, nNeutrons, nPi0s, nGamas;
569 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
570 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
571 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
572 Float_t maxPt = .0, meanPt = .0, totalP = .0;
573 Int_t fVertexflag;
574 Int_t iRunNumber = 0;
575 TString fVertexName;
576
577 AliRunTag *tag = new AliRunTag();
578 AliEventTag *evTag = new AliEventTag();
579 TTree ttag("T","A Tree with event tags");
580 TBranch * btag = ttag.Branch("AliTAG", &tag);
581 btag->SetCompressionLevel(9);
103d3fba 582 gSystem->GetMemInfo(meminfo);
583 AliInfo(Form("After the tag initialization - Memory used: %d MB",meminfo->fMemUsed));
584 Int_t tempmem = meminfo->fMemUsed;
08e1a23e 585
a1069ee1 586 AliInfo(Form("Creating the ESD tags......."));
08e1a23e 587
588 Int_t firstEvent = 0,lastEvent = 0;
589 TTree *t = (TTree*) file->Get("esdTree");
590 AliESDEvent *esd = new AliESDEvent();
591 esd->ReadFromTree(t);
592
103d3fba 593 gSystem->GetMemInfo(meminfo);
594 AliInfo(Form("After the esd initialization - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
595 tempmem = meminfo->fMemUsed;
596
08e1a23e 597 t->GetEntry(0);
598 Int_t iInitRunNumber = esd->GetRunNumber();
599
600 Int_t iNumberOfEvents = (Int_t)t->GetEntries();
601 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
602 ntrack = 0;
603 nPos = 0;
604 nNeg = 0;
605 nNeutr =0;
606 nK0s = 0;
607 nNeutrons = 0;
608 nPi0s = 0;
609 nGamas = 0;
610 nProtons = 0;
611 nKaons = 0;
612 nPions = 0;
613 nMuons = 0;
614 nElectrons = 0;
615 nCh1GeV = 0;
616 nCh3GeV = 0;
617 nCh10GeV = 0;
618 nMu1GeV = 0;
619 nMu3GeV = 0;
620 nMu10GeV = 0;
621 nEl1GeV = 0;
622 nEl3GeV = 0;
623 nEl10GeV = 0;
624 maxPt = .0;
625 meanPt = .0;
626 totalP = .0;
627 fVertexflag = 1;
628
629 t->GetEntry(iEventNumber);
630 iRunNumber = esd->GetRunNumber();
631 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD!!!");
632 const AliESDVertex * vertexIn = esd->GetVertex();
633 fVertexName = vertexIn->GetName();
634 if(fVertexName == "default") fVertexflag = 0;
635
636 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
637 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
638 if(esdTrack->GetLabel() != 0) fIsSim = kTRUE;
639 else if(esdTrack->GetLabel() == 0) fIsSim = kFALSE;
640 UInt_t status = esdTrack->GetStatus();
641
642 //select only tracks with ITS refit
643 if ((status&AliESDtrack::kITSrefit)==0) continue;
644 //select only tracks with TPC refit
645 if ((status&AliESDtrack::kTPCrefit)==0) continue;
646
647 //select only tracks with the "combined PID"
648 if ((status&AliESDtrack::kESDpid)==0) continue;
649 Double_t p[3];
650 esdTrack->GetPxPyPz(p);
651 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
652 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
653 totalP += momentum;
654 meanPt += fPt;
655 if(fPt > maxPt) maxPt = fPt;
656
657 if(esdTrack->GetSign() > 0) {
658 nPos++;
659 if(fPt > fLowPtCut) nCh1GeV++;
660 if(fPt > fHighPtCut) nCh3GeV++;
661 if(fPt > fVeryHighPtCut) nCh10GeV++;
662 }
663 if(esdTrack->GetSign() < 0) {
664 nNeg++;
665 if(fPt > fLowPtCut) nCh1GeV++;
666 if(fPt > fHighPtCut) nCh3GeV++;
667 if(fPt > fVeryHighPtCut) nCh10GeV++;
668 }
669 if(esdTrack->GetSign() == 0) nNeutr++;
670
671 //PID
672 Double_t prob[5];
673 esdTrack->GetESDpid(prob);
674
675 Double_t rcc = 0.0;
676 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
677 if(rcc == 0.0) continue;
678 //Bayes' formula
679 Double_t w[5];
680 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
681
682 //protons
683 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
684 //kaons
685 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
686 //pions
687 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
688 //electrons
689 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
690 nElectrons++;
691 if(fPt > fLowPtCut) nEl1GeV++;
692 if(fPt > fHighPtCut) nEl3GeV++;
693 if(fPt > fVeryHighPtCut) nEl10GeV++;
694 }
695 ntrack++;
696 }//esd track loop
697
698 /////////////
699 //muon code//
700 ////////////
701 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
702 // loop over all reconstructed tracks (also first track of combination)
703 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
704 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
705 if (muonTrack == 0x0) continue;
706
707 // Coordinates at vertex
708 fZ = muonTrack->GetZ();
709 fY = muonTrack->GetBendingCoor();
710 fX = muonTrack->GetNonBendingCoor();
711
712 fThetaX = muonTrack->GetThetaX();
713 fThetaY = muonTrack->GetThetaY();
714
715 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
716 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
717 fPxRec = fPzRec * TMath::Tan(fThetaX);
718 fPyRec = fPzRec * TMath::Tan(fThetaY);
719 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
720
721 //ChiSquare of the track if needed
722 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
723 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
724 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
725
726 // total number of muons inside a vertex cut
727 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
728 nMuons++;
729 if(fEPvector.Pt() > fLowPtCut) {
730 nMu1GeV++;
731 if(fEPvector.Pt() > fHighPtCut) {
732 nMu3GeV++;
733 if (fEPvector.Pt() > fVeryHighPtCut) {
734 nMu10GeV++;
735 }
736 }
737 }
738 }
739 }//muon track loop
740
741 // Fill the event tags
742 if(ntrack != 0) meanPt = meanPt/ntrack;
743
744 evTag->SetEventId(iEventNumber+1);
745 evTag->SetGUID(fguid);
746 evTag->SetMD5(fmd5);
747 evTag->SetTURL(fturl);
748 evTag->SetSize(size);
749 evTag->SetVertexX(vertexIn->GetXv());
750 evTag->SetVertexY(vertexIn->GetYv());
751 evTag->SetVertexZ(vertexIn->GetZv());
752 evTag->SetVertexZError(vertexIn->GetZRes());
753 evTag->SetVertexFlag(fVertexflag);
754
755 evTag->SetT0VertexZ(esd->GetT0zVertex());
756
757 evTag->SetTriggerMask(esd->GetTriggerMask());
758 evTag->SetTriggerCluster(esd->GetTriggerCluster());
759
760 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
761 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
a85132e7 762 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy(0),esd->GetZDCEMEnergy(1));
08e1a23e 763 evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
764 evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
765 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
766
767
768 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
769 evTag->SetNumOfPosTracks(nPos);
770 evTag->SetNumOfNegTracks(nNeg);
771 evTag->SetNumOfNeutrTracks(nNeutr);
772
773 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
774 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
775 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
776 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
777
778 evTag->SetNumOfProtons(nProtons);
779 evTag->SetNumOfKaons(nKaons);
780 evTag->SetNumOfPions(nPions);
781 evTag->SetNumOfMuons(nMuons);
782 evTag->SetNumOfElectrons(nElectrons);
783 evTag->SetNumOfPhotons(nGamas);
784 evTag->SetNumOfPi0s(nPi0s);
785 evTag->SetNumOfNeutrons(nNeutrons);
786 evTag->SetNumOfKaon0s(nK0s);
787
788 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
789 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
790 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
791 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
792 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
793 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
794 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
795 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
796 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
797
798 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
799 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
800
801 evTag->SetTotalMomentum(totalP);
802 evTag->SetMeanPt(meanPt);
803 evTag->SetMaxPt(maxPt);
804
805 tag->SetRunId(iInitRunNumber);
806 if(fIsSim) tag->SetDataType(0);
807 else tag->SetDataType(1);
808 tag->AddEventTag(*evTag);
809 }//event loop
810 lastEvent = iNumberOfEvents;
811
103d3fba 812 gSystem->GetMemInfo(meminfo);
813 AliInfo(Form("After the event and track loop - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
814 tempmem = meminfo->fMemUsed;
08e1a23e 815 t->Delete("");
816
103d3fba 817 gSystem->GetMemInfo(meminfo);
818 AliInfo(Form("After the t->Delete - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
819 tempmem = meminfo->fMemUsed;
820
08e1a23e 821 TString localFileName = "Run"; localFileName += tag->GetRunId();
822 localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
823 localFileName += ".ESD.tag.root";
824
825 TString fileName;
826
827 if(fStorage == 0) {
828 fileName = localFileName.Data();
829 AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
830 }
831 else if(fStorage == 1) {
832 TString alienLocation = "/alien";
833 alienLocation += gGrid->Pwd();
834 alienLocation += fgridpath.Data();
835 alienLocation += "/";
836 alienLocation += localFileName;
837 alienLocation += "?se=";
838 alienLocation += fSE.Data();
839 fileName = alienLocation.Data();
840 AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
841 }
842
843 TFile* ftag = TFile::Open(fileName, "recreate");
844 ftag->cd();
7ee67d2a 845 ttag.Fill();
846 tag->Clear();
08e1a23e 847 ttag.Write();
848 ftag->Close();
849
103d3fba 850 gSystem->GetMemInfo(meminfo);
851 AliInfo(Form("After the file closing - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
852 tempmem = meminfo->fMemUsed;
853
08e1a23e 854 delete ftag;
855 delete esd;
856
857 delete tag;
103d3fba 858 gSystem->GetMemInfo(meminfo);
859 AliInfo(Form("After the delete objects - Memory used: %d MB - Increase: %d MB",meminfo->fMemUsed,meminfo->fMemUsed - tempmem));
08e1a23e 860}
861
862//_____________________________________________________________________________
863void AliESDTagCreator::CreateTag(TFile* file, const char *filepath, Int_t Counter) {
864 //private method that creates tag files
865 /////////////
866 //muon code//
867 ////////////
868 Double_t fMUONMASS = 0.105658369;
869 //Variables
870 Double_t fX,fY,fZ ;
871 Double_t fThetaX, fThetaY, fPyz, fChisquare;
872 Double_t fPxRec, fPyRec, fPzRec, fEnergy;
873 Int_t fCharge;
874 TLorentzVector fEPvector;
875
876 Float_t fZVertexCut = 10.0;
877 Float_t fRhoVertexCut = 2.0;
878
879 Float_t fLowPtCut = 1.0;
880 Float_t fHighPtCut = 3.0;
881 Float_t fVeryHighPtCut = 10.0;
882 ////////////
883
884 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
885
886 // Creates the tags for all the events in a given ESD file
887 Bool_t fIsSim = kTRUE;
888 Int_t ntrack;
889 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
890 Int_t nPos, nNeg, nNeutr;
891 Int_t nK0s, nNeutrons, nPi0s, nGamas;
892 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
893 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
894 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
895 Float_t maxPt = .0, meanPt = .0, totalP = .0;
896 Int_t fVertexflag;
897 Int_t iRunNumber = 0;
898 TString fVertexName;
899
900 AliRunTag *tag = new AliRunTag();
901 AliEventTag *evTag = new AliEventTag();
902 TTree ttag("T","A Tree with event tags");
903 TBranch * btag = ttag.Branch("AliTAG", &tag);
904 btag->SetCompressionLevel(9);
905
a1069ee1 906 AliInfo(Form("Creating the ESD tags......."));
08e1a23e 907
908 Int_t firstEvent = 0,lastEvent = 0;
909
910 TTree *t = (TTree*) file->Get("esdTree");
911 AliESDEvent *esd = new AliESDEvent();
912 esd->ReadFromTree(t);
913
914 t->GetEntry(0);
915 Int_t iInitRunNumber = esd->GetRunNumber();
916
917 Int_t iNumberOfEvents = (Int_t)t->GetEntries();
918 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
919 ntrack = 0;
920 nPos = 0;
921 nNeg = 0;
922 nNeutr =0;
923 nK0s = 0;
924 nNeutrons = 0;
925 nPi0s = 0;
926 nGamas = 0;
927 nProtons = 0;
928 nKaons = 0;
929 nPions = 0;
930 nMuons = 0;
931 nElectrons = 0;
932 nCh1GeV = 0;
933 nCh3GeV = 0;
934 nCh10GeV = 0;
935 nMu1GeV = 0;
936 nMu3GeV = 0;
937 nMu10GeV = 0;
938 nEl1GeV = 0;
939 nEl3GeV = 0;
940 nEl10GeV = 0;
941 maxPt = .0;
942 meanPt = .0;
943 totalP = .0;
944 fVertexflag = 1;
945
946 t->GetEntry(iEventNumber);
947 iRunNumber = esd->GetRunNumber();
948 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD!!!");
949 const AliESDVertex * vertexIn = esd->GetVertex();
950 fVertexName = vertexIn->GetName();
951 if(fVertexName == "default") fVertexflag = 0;
952
953 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
954 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
955 if(esdTrack->GetLabel() != 0) fIsSim = kTRUE;
956 else if(esdTrack->GetLabel() == 0) fIsSim = kFALSE;
957 UInt_t status = esdTrack->GetStatus();
958
959 //select only tracks with ITS refit
960 if ((status&AliESDtrack::kITSrefit)==0) continue;
961 //select only tracks with TPC refit
962 if ((status&AliESDtrack::kTPCrefit)==0) continue;
963
964 //select only tracks with the "combined PID"
965 if ((status&AliESDtrack::kESDpid)==0) continue;
966 Double_t p[3];
967 esdTrack->GetPxPyPz(p);
968 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
969 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
970 totalP += momentum;
971 meanPt += fPt;
972 if(fPt > maxPt) maxPt = fPt;
973
974 if(esdTrack->GetSign() > 0) {
975 nPos++;
976 if(fPt > fLowPtCut) nCh1GeV++;
977 if(fPt > fHighPtCut) nCh3GeV++;
978 if(fPt > fVeryHighPtCut) nCh10GeV++;
979 }
980 if(esdTrack->GetSign() < 0) {
981 nNeg++;
982 if(fPt > fLowPtCut) nCh1GeV++;
983 if(fPt > fHighPtCut) nCh3GeV++;
984 if(fPt > fVeryHighPtCut) nCh10GeV++;
985 }
986 if(esdTrack->GetSign() == 0) nNeutr++;
987
988 //PID
989 Double_t prob[5];
990 esdTrack->GetESDpid(prob);
991
992 Double_t rcc = 0.0;
993 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
994 if(rcc == 0.0) continue;
995 //Bayes' formula
996 Double_t w[5];
997 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
998
999 //protons
1000 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
1001 //kaons
1002 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
1003 //pions
1004 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
1005 //electrons
1006 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
1007 nElectrons++;
1008 if(fPt > fLowPtCut) nEl1GeV++;
1009 if(fPt > fHighPtCut) nEl3GeV++;
1010 if(fPt > fVeryHighPtCut) nEl10GeV++;
1011 }
1012 ntrack++;
1013 }//esd track loop
1014
1015 /////////////
1016 //muon code//
1017 ////////////
1018 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
1019 // loop over all reconstructed tracks (also first track of combination)
1020 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
1021 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
1022 if (muonTrack == 0x0) continue;
1023
1024 // Coordinates at vertex
1025 fZ = muonTrack->GetZ();
1026 fY = muonTrack->GetBendingCoor();
1027 fX = muonTrack->GetNonBendingCoor();
1028
1029 fThetaX = muonTrack->GetThetaX();
1030 fThetaY = muonTrack->GetThetaY();
1031
1032 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
1033 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
1034 fPxRec = fPzRec * TMath::Tan(fThetaX);
1035 fPyRec = fPzRec * TMath::Tan(fThetaY);
1036 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
1037
1038 //ChiSquare of the track if needed
1039 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
1040 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
1041 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
1042
1043 // total number of muons inside a vertex cut
1044 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
1045 nMuons++;
1046 if(fEPvector.Pt() > fLowPtCut) {
1047 nMu1GeV++;
1048 if(fEPvector.Pt() > fHighPtCut) {
1049 nMu3GeV++;
1050 if (fEPvector.Pt() > fVeryHighPtCut) {
1051 nMu10GeV++;
1052 }
1053 }
1054 }
1055 }
1056 }//muon track loop
1057
1058 // Fill the event tags
1059 if(ntrack != 0) meanPt = meanPt/ntrack;
1060
1061 evTag->SetEventId(iEventNumber+1);
1062 evTag->SetPath(filepath);
1063
1064 evTag->SetVertexX(vertexIn->GetXv());
1065 evTag->SetVertexY(vertexIn->GetYv());
1066 evTag->SetVertexZ(vertexIn->GetZv());
1067 evTag->SetVertexZError(vertexIn->GetZRes());
1068 evTag->SetVertexFlag(fVertexflag);
1069
1070 evTag->SetT0VertexZ(esd->GetT0zVertex());
1071
1072 evTag->SetTriggerMask(esd->GetTriggerMask());
1073 evTag->SetTriggerCluster(esd->GetTriggerCluster());
1074
1075 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
1076 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
a85132e7 1077 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy(0),esd->GetZDCEMEnergy(1));
08e1a23e 1078 evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
1079 evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
1080 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
1081
1082
1083 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
1084 evTag->SetNumOfPosTracks(nPos);
1085 evTag->SetNumOfNegTracks(nNeg);
1086 evTag->SetNumOfNeutrTracks(nNeutr);
1087
1088 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
1089 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
1090 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
1091 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
1092
1093 evTag->SetNumOfProtons(nProtons);
1094 evTag->SetNumOfKaons(nKaons);
1095 evTag->SetNumOfPions(nPions);
1096 evTag->SetNumOfMuons(nMuons);
1097 evTag->SetNumOfElectrons(nElectrons);
1098 evTag->SetNumOfPhotons(nGamas);
1099 evTag->SetNumOfPi0s(nPi0s);
1100 evTag->SetNumOfNeutrons(nNeutrons);
1101 evTag->SetNumOfKaon0s(nK0s);
1102
1103 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
1104 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
1105 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
1106 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
1107 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
1108 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
1109 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
1110 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
1111 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
1112
1113 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
1114 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
1115
1116 evTag->SetTotalMomentum(totalP);
1117 evTag->SetMeanPt(meanPt);
1118 evTag->SetMaxPt(maxPt);
1119
1120 tag->SetRunId(iInitRunNumber);
1121 if(fIsSim) tag->SetDataType(0);
1122 else tag->SetDataType(1);
1123 tag->AddEventTag(*evTag);
1124 }//event loop
1125 lastEvent = iNumberOfEvents;
1126
1127 t->Delete("");
1128
08e1a23e 1129 TString localFileName = "Run"; localFileName += tag->GetRunId();
1130 localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
1131 localFileName += ".ESD.tag.root";
1132
1133 TString fileName;
1134
1135 if(fStorage == 0) {
1136 fileName = localFileName.Data();
1137 AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
1138 }
1139 else if(fStorage == 1) {
1140 TString alienLocation = "/alien";
1141 alienLocation += gGrid->Pwd();
1142 alienLocation += fgridpath.Data();
1143 alienLocation += "/";
1144 alienLocation += localFileName;
1145 alienLocation += "?se=";
1146 alienLocation += fSE.Data();
1147 fileName = alienLocation.Data();
1148 AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
1149 }
1150
1151 TFile* ftag = TFile::Open(fileName, "recreate");
1152 ftag->cd();
7ee67d2a 1153 ttag.Fill();
1154 tag->Clear();
08e1a23e 1155 ttag.Write();
1156 ftag->Close();
1157
1158 delete ftag;
1159 delete esd;
1160
1161 delete tag;
1162}
1163
1164//_____________________________________________________________________________
6b6e4472 1165void AliESDTagCreator::CreateESDTags(Int_t fFirstEvent, Int_t fLastEvent, TMap */* grpData */) {
08e1a23e 1166 //GRP
1167 Float_t lhcLuminosity = 0.0;
1168 TString lhcState = "test";
1169 UInt_t detectorMask = 0;
1170
1171 /////////////
1172 //muon code//
1173 ////////////
1174 Double_t fMUONMASS = 0.105658369;
1175 //Variables
1176 Double_t fX,fY,fZ ;
1177 Double_t fThetaX, fThetaY, fPyz, fChisquare;
1178 Double_t fPxRec,fPyRec, fPzRec, fEnergy;
1179 Int_t fCharge;
1180 TLorentzVector fEPvector;
1181
1182 Float_t fZVertexCut = 10.0;
1183 Float_t fRhoVertexCut = 2.0;
1184
1185 Float_t fLowPtCut = 1.0;
1186 Float_t fHighPtCut = 3.0;
1187 Float_t fVeryHighPtCut = 10.0;
1188 ////////////
1189
1190 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
1191
1192 // Creates the tags for all the events in a given ESD file
1193 Int_t ntrack;
1194 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
1195 Int_t nPos, nNeg, nNeutr;
1196 Int_t nK0s, nNeutrons, nPi0s, nGamas;
1197 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
1198 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
1199 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
1200 Float_t maxPt = .0, meanPt = .0, totalP = .0;
1201 Int_t fVertexflag;
1202 Int_t iRunNumber = 0;
1203 TString fVertexName("default");
08e1a23e 1204
a1069ee1 1205 AliInfo(Form("Creating the ESD tags......."));
08e1a23e 1206
1207 TFile *file = TFile::Open("AliESDs.root");
1208 if (!file || !file->IsOpen()) {
1209 AliError(Form("opening failed"));
1210 delete file;
1211 return ;
1212 }
1213 Int_t lastEvent = 0;
1214 TTree *b = (TTree*) file->Get("esdTree");
1215 AliESDEvent *esd = new AliESDEvent();
1216 esd->ReadFromTree(b);
1217
1218 b->GetEntry(fFirstEvent);
1219 Int_t iInitRunNumber = esd->GetRunNumber();
1220
1221 Int_t iNumberOfEvents = (Int_t)b->GetEntries();
fcc6b05f 1222 if(fLastEvent == -1) lastEvent = (Int_t)b->GetEntries();
1223 else lastEvent = fLastEvent;
1224
1225 char fileName[256];
1226 sprintf(fileName, "Run%d.Event%d_%d.ESD.tag.root",
517aef4e 1227 iInitRunNumber,fFirstEvent,lastEvent);
fcc6b05f 1228 AliInfo(Form("writing tags to file %s", fileName));
1229 AliDebug(1, Form("writing tags to file %s", fileName));
1230
1231 TFile* ftag = TFile::Open(fileName, "recreate");
1232
517aef4e 1233 AliRunTag *tag = new AliRunTag();
1234 AliEventTag *evTag = new AliEventTag();
1235 TTree ttag("T","A Tree with event tags");
1236 TBranch * btag = ttag.Branch("AliTAG", &tag);
1237 btag->SetCompressionLevel(9);
1238
08e1a23e 1239 if(fLastEvent != -1) iNumberOfEvents = fLastEvent + 1;
1240 for (Int_t iEventNumber = fFirstEvent; iEventNumber < iNumberOfEvents; iEventNumber++) {
1241 ntrack = 0;
1242 nPos = 0;
1243 nNeg = 0;
1244 nNeutr =0;
1245 nK0s = 0;
1246 nNeutrons = 0;
1247 nPi0s = 0;
1248 nGamas = 0;
1249 nProtons = 0;
1250 nKaons = 0;
1251 nPions = 0;
1252 nMuons = 0;
1253 nElectrons = 0;
1254 nCh1GeV = 0;
1255 nCh3GeV = 0;
1256 nCh10GeV = 0;
1257 nMu1GeV = 0;
1258 nMu3GeV = 0;
1259 nMu10GeV = 0;
1260 nEl1GeV = 0;
1261 nEl3GeV = 0;
1262 nEl10GeV = 0;
1263 maxPt = .0;
1264 meanPt = .0;
1265 totalP = .0;
1266 fVertexflag = 0;
1267
1268 b->GetEntry(iEventNumber);
1269 iRunNumber = esd->GetRunNumber();
1270 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD!!!");
1271 const AliESDVertex * vertexIn = esd->GetVertex();
1272 if (!vertexIn) AliError("ESD has not defined vertex.");
1273 if (vertexIn) fVertexName = vertexIn->GetName();
1274 if(fVertexName != "default") fVertexflag = 1;
1275 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
1276 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
1277 UInt_t status = esdTrack->GetStatus();
1278
1279 //select only tracks with ITS refit
1280 if ((status&AliESDtrack::kITSrefit)==0) continue;
1281 //select only tracks with TPC refit
1282 if ((status&AliESDtrack::kTPCrefit)==0) continue;
1283
1284 //select only tracks with the "combined PID"
1285 if ((status&AliESDtrack::kESDpid)==0) continue;
1286 Double_t p[3];
1287 esdTrack->GetPxPyPz(p);
1288 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
1289 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
1290 totalP += momentum;
1291 meanPt += fPt;
1292 if(fPt > maxPt) maxPt = fPt;
1293
1294 if(esdTrack->GetSign() > 0) {
1295 nPos++;
1296 if(fPt > fLowPtCut) nCh1GeV++;
1297 if(fPt > fHighPtCut) nCh3GeV++;
1298 if(fPt > fVeryHighPtCut) nCh10GeV++;
1299 }
1300 if(esdTrack->GetSign() < 0) {
1301 nNeg++;
1302 if(fPt > fLowPtCut) nCh1GeV++;
1303 if(fPt > fHighPtCut) nCh3GeV++;
1304 if(fPt > fVeryHighPtCut) nCh10GeV++;
1305 }
1306 if(esdTrack->GetSign() == 0) nNeutr++;
1307
1308 //PID
1309 Double_t prob[5];
1310 esdTrack->GetESDpid(prob);
1311
1312 Double_t rcc = 0.0;
1313 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
1314 if(rcc == 0.0) continue;
1315 //Bayes' formula
1316 Double_t w[5];
1317 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
1318
1319 //protons
1320 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
1321 //kaons
1322 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
1323 //pions
1324 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
1325 //electrons
1326 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
1327 nElectrons++;
1328 if(fPt > fLowPtCut) nEl1GeV++;
1329 if(fPt > fHighPtCut) nEl3GeV++;
1330 if(fPt > fVeryHighPtCut) nEl10GeV++;
1331 }
1332 ntrack++;
1333 }//track loop
1334
1335 /////////////
1336 //muon code//
1337 ////////////
1338 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
1339 // loop over all reconstructed tracks (also first track of combination)
1340 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
1341 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
1342 if (muonTrack == 0x0) continue;
1343
1344 // Coordinates at vertex
1345 fZ = muonTrack->GetZ();
1346 fY = muonTrack->GetBendingCoor();
1347 fX = muonTrack->GetNonBendingCoor();
1348
1349 fThetaX = muonTrack->GetThetaX();
1350 fThetaY = muonTrack->GetThetaY();
1351
1352 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
1353 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
1354 fPxRec = fPzRec * TMath::Tan(fThetaX);
1355 fPyRec = fPzRec * TMath::Tan(fThetaY);
1356 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
1357
1358 //ChiSquare of the track if needed
1359 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
1360 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
1361 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
1362
1363 // total number of muons inside a vertex cut
1364 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
1365 nMuons++;
1366 if(fEPvector.Pt() > fLowPtCut) {
1367 nMu1GeV++;
1368 if(fEPvector.Pt() > fHighPtCut) {
1369 nMu3GeV++;
1370 if (fEPvector.Pt() > fVeryHighPtCut) {
1371 nMu10GeV++;
1372 }
1373 }
1374 }
1375 }
1376 }//muon track loop
1377
1378 // Fill the event tags
1379 if(ntrack != 0)
1380 meanPt = meanPt/ntrack;
1381
1382 evTag->SetEventId(iEventNumber+1);
1383 if (vertexIn) {
1384 evTag->SetVertexX(vertexIn->GetXv());
1385 evTag->SetVertexY(vertexIn->GetYv());
1386 evTag->SetVertexZ(vertexIn->GetZv());
1387 evTag->SetVertexZError(vertexIn->GetZRes());
1388 }
1389 evTag->SetVertexFlag(fVertexflag);
1390
1391 evTag->SetT0VertexZ(esd->GetT0zVertex());
1392
1393 evTag->SetTriggerMask(esd->GetTriggerMask());
1394 evTag->SetTriggerCluster(esd->GetTriggerCluster());
1395
1396 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
1397 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
1398 evTag->SetZDCNeutron2Energy(esd->GetZDCN2Energy());
1399 evTag->SetZDCProton2Energy(esd->GetZDCP2Energy());
a85132e7 1400 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy(0),esd->GetZDCEMEnergy(1));
08e1a23e 1401 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
1402
1403
1404 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
1405 evTag->SetNumOfPosTracks(nPos);
1406 evTag->SetNumOfNegTracks(nNeg);
1407 evTag->SetNumOfNeutrTracks(nNeutr);
1408
1409 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
1410 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
1411 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
1412 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
1413
1414 evTag->SetNumOfProtons(nProtons);
1415 evTag->SetNumOfKaons(nKaons);
1416 evTag->SetNumOfPions(nPions);
1417 evTag->SetNumOfMuons(nMuons);
1418 evTag->SetNumOfElectrons(nElectrons);
1419 evTag->SetNumOfPhotons(nGamas);
1420 evTag->SetNumOfPi0s(nPi0s);
1421 evTag->SetNumOfNeutrons(nNeutrons);
1422 evTag->SetNumOfKaon0s(nK0s);
1423
1424 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
1425 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
1426 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
1427 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
1428 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
1429 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
1430 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
1431 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
1432 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
1433
1434 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
1435 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
1436
1437 evTag->SetTotalMomentum(totalP);
1438 evTag->SetMeanPt(meanPt);
1439 evTag->SetMaxPt(maxPt);
1440
1441 tag->SetLHCTag(lhcLuminosity,lhcState);
1442 tag->SetDetectorTag(detectorMask);
1443
1444 tag->SetRunId(iInitRunNumber);
1445 tag->AddEventTag(*evTag);
1446 }
08e1a23e 1447
08e1a23e 1448 ftag->cd();
7ee67d2a 1449 ttag.Fill();
1450 tag->Clear();
08e1a23e 1451 ttag.Write();
1452 ftag->Close();
1453 file->cd();
1454 delete tag;
1455 delete evTag;
1456}
1457
517aef4e 1458//_____________________________________________________________________________
1459void AliESDTagCreator::SwitchOffBranches() const {
2856e38b 1460 //
1461 // Switch of branches on user request
517aef4e 1462 TObjArray * tokens = fBranches.Tokenize(" ");
1463 Int_t ntok = tokens->GetEntries();
1464 for (Int_t i = 0; i < ntok; i++) {
1465 TString str = ((TObjString*) tokens->At(i))->GetString();
1466 fChain->SetBranchStatus(Form("%s%s%s","*", str.Data(), "*"), 0);
1467 AliInfo(Form("Branch %s switched off \n", str.Data()));
1468 }
2856e38b 1469}