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