1 /**************************************************************************
2 * Author: Panos Christakoglou. *
3 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 //-----------------------------------------------------------------
17 // AliAODTagCreator 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 //-----------------------------------------------------------------
23 #include <Riostream.h>
29 #include <TLorentzVector.h>
33 #include <TGridResult.h>
36 #include "AliRunTag.h"
37 #include "AliEventTag.h"
40 #include "AliAODEvent.h"
41 #include "AliAODVertex.h"
44 #include "AliAODTagCreator.h"
47 ClassImp(AliAODTagCreator)
50 //______________________________________________________________________________
51 AliAODTagCreator::AliAODTagCreator() :
53 //==============Default constructor for a AliAODTagCreator================
56 //______________________________________________________________________________
57 AliAODTagCreator::~AliAODTagCreator() {
58 //================Default destructor for a AliAODTagCreator===================
61 //______________________________________________________________________________
62 Bool_t AliAODTagCreator::ReadGridCollection(TGridResult *fresult) {
63 // Reads the entry of the TGridResult and creates the tags
64 Int_t nEntries = fresult->GetEntries();
73 for(Int_t i = 0; i < nEntries; i++) {
74 alienUrl = fresult->GetKey(i,"turl");
75 guid = fresult->GetKey(i,"guid");
76 if(fresult->GetKey(i,"size")) size = atol (fresult->GetKey(i,"size"));
77 md5 = fresult->GetKey(i,"md5");
78 turl = fresult->GetKey(i,"turl");
79 if(md5 && !strlen(guid)) md5 = 0;
80 if(guid && !strlen(guid)) guid = 0;
82 TFile *f = TFile::Open(alienUrl,"READ");
83 //CreateTag(f,guid,md5,turl,size,counter);
92 //______________________________________________________________________________
93 Bool_t AliAODTagCreator::ReadLocalCollection(const char *localpath) {
94 // Checks the different subdirs of the given local path and in the
95 // case where it finds an AliAODs.root file it creates the tags
97 void *dira = gSystem->OpenDirectory(localpath);
99 const char * dirname = 0x0;
100 const char * filename = 0x0;
101 const char * pattern = "AliAODs.root";
104 while((dirname = gSystem->GetDirEntry(dira))) {
105 sprintf(fPath,"%s/%s",localpath,dirname);
106 void *dirb = gSystem->OpenDirectory(fPath);
107 while((filename = gSystem->GetDirEntry(dirb))) {
108 if(strstr(filename,pattern)) {
109 TString fAODFileName;
110 fAODFileName = fPath;
112 fAODFileName += pattern;
113 TFile *f = TFile::Open(fAODFileName,"READ");
114 //CreateTag(f,fAODFileName,counter);
120 }//child directory's entry loop
121 }//parent directory's entry loop
126 //______________________________________________________________________________
127 Bool_t AliAODTagCreator::ReadCAFCollection(const char *filename) {
128 // Temporary solution for CAF: Takes as an input the ascii file that
129 // lists the AODs stored in the SE of the CAF and creates the tags.
131 // Open the input stream
137 // Read the input list of files and add them to the chain
140 if (!esdfile.Contains("root")) continue; // protection
141 TFile *f = TFile::Open(esdfile,"READ");
142 //CreateTag(f,esdfile,counter);
152 //__________________________________________________________________________
153 void AliAODTagCreator::CreateAODTags(Int_t fFirstEvent, Int_t fLastEvent) {
154 //creates tag files for AODs
156 Float_t fLowPtCut = 1.0;
157 Float_t fHighPtCut = 3.0;
158 Float_t fVeryHighPtCut = 10.0;
161 Double_t partFrac[10] = {0.01, 0.01, 0.85, 0.10, 0.05, 0., 0., 0., 0., 0.};
163 // Creates the tags for all the events in a given AOD file
165 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
166 Int_t nPos, nNeg, nNeutr;
167 Int_t nKinks, nV0s, nCascades;
168 Int_t nK0s, nNeutrons, nPi0s, nGamas;
169 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
170 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
171 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
172 Float_t maxPt = .0, meanPt = .0, totalP = .0;
174 AliRunTag *tag = new AliRunTag();
175 TTree ttag("T","A Tree with event tags");
176 TBranch * btag = ttag.Branch("AliTAG", &tag);
177 btag->SetCompressionLevel(9);
179 //reading the esd tag file
180 TChain *oldTagTree = new TChain("T");
181 const char * tagPattern = "ESD.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)) oldTagTree->Add(name);
189 AliInfo(Form("Chained tag files: %d",oldTagTree->GetEntries()));
191 //reading the esd tag file
192 AliRunTag *oldtag = new AliRunTag();
194 oldTagTree->SetBranchAddress("AliTAG",&oldtag);
195 oldTagTree->GetEntry(0);
196 tag->CopyStandardContent(oldtag);
197 const TClonesArray *evTagList = oldtag->GetEventTags();
199 AliInfo(Form("Creating the AOD tags......."));
201 TFile *file = TFile::Open("AliAOD.root");
202 if (!file || !file->IsOpen()) {
203 AliError(Form("opening failed"));
207 TTree *aodTree = (TTree*)file->Get("aodTree");
208 AliAODEvent *aod = new AliAODEvent();
209 aod->ReadFromTree(aodTree);
212 if(fLastEvent == -1) lastEvent = (Int_t)aodTree->GetEntries();
213 else lastEvent = fLastEvent;
216 Int_t nEvents = aodTree->GetEntries();
217 for (Int_t iEventNumber = 0; iEventNumber < nEvents; iEventNumber++) {
218 AliEventTag *evTag = (AliEventTag *)evTagList->At(iEventNumber);
220 nPos = 0; nNeg = 0; nNeutr =0;
221 nKinks = 0; nV0s = 0; nCascades = 0;
222 nK0s = 0; nNeutrons = 0; nPi0s = 0; nGamas = 0;
223 nProtons = 0; nKaons = 0; nPions = 0; nMuons = 0; nElectrons = 0;
224 nCh1GeV = 0; nCh3GeV = 0; nCh10GeV = 0;
225 nMu1GeV = 0; nMu3GeV = 0; nMu10GeV = 0;
226 nEl1GeV = 0; nEl3GeV = 0; nEl10GeV = 0;
227 maxPt = .0; meanPt = .0; totalP = .0;
230 aodTree->GetEvent(iEventNumber);
233 aod->GetStdContent();
235 Int_t nTracks = aod->GetNTracks();
236 // loop over vertices
237 Int_t nVtxs = aod->GetNVertices();
238 for (Int_t nVtx = 0; nVtx < nVtxs; nVtx++) {
240 AliAODVertex *vertex = aod->GetVertex(nVtx);
241 if(vertex->GetType() == 1) nKinks += 1;
242 if(vertex->GetType() == 2) nV0s += 1;
243 if(vertex->GetType() == 3) nCascades += 1;
245 for (Int_t nTr = 0; nTr < nTracks; nTr++) {
246 AliAODTrack *track = aod->GetTrack(nTr);
248 Double_t fPt = track->Pt();
249 if(fPt > maxPt) maxPt = fPt;
250 if(track->Charge() > 0) {
252 if(fPt > fLowPtCut) nCh1GeV++;
253 if(fPt > fHighPtCut) nCh3GeV++;
254 if(fPt > fVeryHighPtCut) nCh10GeV++;
256 if(track->Charge() < 0) {
258 if(fPt > fLowPtCut) nCh1GeV++;
259 if(fPt > fHighPtCut) nCh3GeV++;
260 if(fPt > fVeryHighPtCut) nCh10GeV++;
262 if(track->Charge() == 0) nNeutr++;
265 const Double32_t *prob = track->PID();
267 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
268 if(rcc == 0.0) continue;
271 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
274 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
276 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
278 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
280 if ((w[1]>w[4])&&(w[1]>w[3])&&(w[1]>w[2])&&(w[1]>w[0])) {
282 if(fPt > fLowPtCut) nMu1GeV++;
283 if(fPt > fHighPtCut) nMu3GeV++;
284 if(fPt > fVeryHighPtCut) nMu10GeV++;
287 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
289 if(fPt > fLowPtCut) nEl1GeV++;
290 if(fPt > fHighPtCut) nEl3GeV++;
291 if(fPt > fVeryHighPtCut) nEl10GeV++;
294 totalP += track->P();
298 // Fill the event tags
300 meanPt = meanPt/ntrack;
302 evTag->SetEventId(iEventNumber+1);
304 evTag->SetNumOfTracks(nTracks);
305 evTag->SetNumOfPosTracks(nPos);
306 evTag->SetNumOfNegTracks(nNeg);
307 evTag->SetNumOfNeutrTracks(nNeutr);
309 evTag->SetNumOfV0s(nV0s);
310 evTag->SetNumOfCascades(nCascades);
311 evTag->SetNumOfKinks(nKinks);
313 evTag->SetNumOfProtons(nProtons);
314 evTag->SetNumOfKaons(nKaons);
315 evTag->SetNumOfPions(nPions);
316 evTag->SetNumOfMuons(nMuons);
317 evTag->SetNumOfElectrons(nElectrons);
318 evTag->SetNumOfPhotons(nGamas);
319 evTag->SetNumOfPi0s(nPi0s);
320 evTag->SetNumOfNeutrons(nNeutrons);
321 evTag->SetNumOfKaon0s(nK0s);
323 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
324 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
325 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
326 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
327 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
328 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
329 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
330 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
331 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
333 evTag->SetTotalMomentum(totalP);
334 evTag->SetMeanPt(meanPt);
335 evTag->SetMaxPt(maxPt);
336 tag->AddEventTag(*evTag);
338 if(fLastEvent == -1) lastEvent = (Int_t)aodTree->GetEntries();
339 else lastEvent = fLastEvent;
345 sprintf(fileName, "Run%d.Event%d_%d.AOD.tag.root",
346 tag->GetRunId(),fFirstEvent,lastEvent );
347 AliInfo(Form("writing tags to file %s", fileName));
348 AliDebug(1, Form("writing tags to file %s", fileName));
350 TFile* ftag = TFile::Open(fileName, "recreate");
358 //_____________________________________________________________________________
359 /*void AliAODTagCreator::CreateTag(TFile* file, const char *guid, const char *md5, const char *turl, Long64_t size, Int_t Counter) {
360 //private method that creates tag files
361 TString fguid = guid;
363 TString fturl = turl;
367 Double_t fMUONMASS = 0.105658369;
370 Double_t fThetaX, fThetaY, fPyz, fChisquare;
371 Double_t fPxRec, fPyRec, fPzRec, fEnergy;
373 TLorentzVector fEPvector;
375 Float_t fZVertexCut = 10.0;
376 Float_t fRhoVertexCut = 2.0;
378 Float_t fLowPtCut = 1.0;
379 Float_t fHighPtCut = 3.0;
380 Float_t fVeryHighPtCut = 10.0;
383 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
385 // Creates the tags for all the events in a given AOD file
386 Bool_t fIsSim = kTRUE;
388 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
389 Int_t nPos, nNeg, nNeutr;
390 Int_t nK0s, nNeutrons, nPi0s, nGamas;
391 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
392 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
393 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
394 Float_t maxPt = .0, meanPt = .0, totalP = .0;
396 Int_t iRunNumber = 0;
399 AliRunTag *tag = new AliRunTag();
400 AliEventTag *evTag = new AliEventTag();
401 TTree ttag("T","A Tree with event tags");
402 TBranch * btag = ttag.Branch("AliTAG", &tag);
403 btag->SetCompressionLevel(9);
405 AliInfo(Form("Creating the tags......."));
407 Int_t firstEvent = 0,lastEvent = 0;
408 TTree *t = (TTree*) file->Get("esdTree");
409 AliESDEvent *esd = new AliESDEvent();
410 esd->ReadFromTree(t);
413 Int_t iInitRunNumber = esd->GetRunNumber();
415 Int_t iNumberOfEvents = (Int_t)t->GetEntries();
416 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
444 t->GetEntry(iEventNumber);
445 iRunNumber = esd->GetRunNumber();
446 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD!!!");
447 const AliESDVertex * vertexIn = esd->GetVertex();
448 fVertexName = vertexIn->GetName();
449 if(fVertexName == "default") fVertexflag = 0;
451 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
452 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
453 if(esdTrack->GetLabel() != 0) fIsSim = kTRUE;
454 else if(esdTrack->GetLabel() == 0) fIsSim = kFALSE;
455 UInt_t status = esdTrack->GetStatus();
457 //select only tracks with ITS refit
458 if ((status&AliESDtrack::kITSrefit)==0) continue;
459 //select only tracks with TPC refit
460 if ((status&AliESDtrack::kTPCrefit)==0) continue;
462 //select only tracks with the "combined PID"
463 if ((status&AliESDtrack::kESDpid)==0) continue;
465 esdTrack->GetPxPyPz(p);
466 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
467 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
470 if(fPt > maxPt) maxPt = fPt;
472 if(esdTrack->GetSign() > 0) {
474 if(fPt > fLowPtCut) nCh1GeV++;
475 if(fPt > fHighPtCut) nCh3GeV++;
476 if(fPt > fVeryHighPtCut) nCh10GeV++;
478 if(esdTrack->GetSign() < 0) {
480 if(fPt > fLowPtCut) nCh1GeV++;
481 if(fPt > fHighPtCut) nCh3GeV++;
482 if(fPt > fVeryHighPtCut) nCh10GeV++;
484 if(esdTrack->GetSign() == 0) nNeutr++;
488 esdTrack->GetESDpid(prob);
491 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
492 if(rcc == 0.0) continue;
495 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
498 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
500 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
502 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
504 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
506 if(fPt > fLowPtCut) nEl1GeV++;
507 if(fPt > fHighPtCut) nEl3GeV++;
508 if(fPt > fVeryHighPtCut) nEl10GeV++;
516 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
517 // loop over all reconstructed tracks (also first track of combination)
518 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
519 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
520 if (muonTrack == 0x0) continue;
522 // Coordinates at vertex
523 fZ = muonTrack->GetZ();
524 fY = muonTrack->GetBendingCoor();
525 fX = muonTrack->GetNonBendingCoor();
527 fThetaX = muonTrack->GetThetaX();
528 fThetaY = muonTrack->GetThetaY();
530 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
531 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
532 fPxRec = fPzRec * TMath::Tan(fThetaX);
533 fPyRec = fPzRec * TMath::Tan(fThetaY);
534 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
536 //ChiSquare of the track if needed
537 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
538 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
539 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
541 // total number of muons inside a vertex cut
542 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
544 if(fEPvector.Pt() > fLowPtCut) {
546 if(fEPvector.Pt() > fHighPtCut) {
548 if (fEPvector.Pt() > fVeryHighPtCut) {
556 // Fill the event tags
557 if(ntrack != 0) meanPt = meanPt/ntrack;
559 evTag->SetEventId(iEventNumber+1);
560 evTag->SetGUID(fguid);
562 evTag->SetTURL(fturl);
563 evTag->SetSize(size);
564 evTag->SetVertexX(vertexIn->GetXv());
565 evTag->SetVertexY(vertexIn->GetYv());
566 evTag->SetVertexZ(vertexIn->GetZv());
567 evTag->SetVertexZError(vertexIn->GetZRes());
568 evTag->SetVertexFlag(fVertexflag);
570 evTag->SetT0VertexZ(esd->GetT0zVertex());
572 evTag->SetTriggerMask(esd->GetTriggerMask());
573 evTag->SetTriggerCluster(esd->GetTriggerCluster());
575 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
576 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
577 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
578 evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
579 evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
580 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
583 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
584 evTag->SetNumOfPosTracks(nPos);
585 evTag->SetNumOfNegTracks(nNeg);
586 evTag->SetNumOfNeutrTracks(nNeutr);
588 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
589 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
590 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
591 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
593 evTag->SetNumOfProtons(nProtons);
594 evTag->SetNumOfKaons(nKaons);
595 evTag->SetNumOfPions(nPions);
596 evTag->SetNumOfMuons(nMuons);
597 evTag->SetNumOfElectrons(nElectrons);
598 evTag->SetNumOfPhotons(nGamas);
599 evTag->SetNumOfPi0s(nPi0s);
600 evTag->SetNumOfNeutrons(nNeutrons);
601 evTag->SetNumOfKaon0s(nK0s);
603 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
604 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
605 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
606 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
607 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
608 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
609 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
610 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
611 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
613 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
614 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
616 evTag->SetTotalMomentum(totalP);
617 evTag->SetMeanPt(meanPt);
618 evTag->SetMaxPt(maxPt);
620 tag->SetRunId(iInitRunNumber);
621 if(fIsSim) tag->SetDataType(0);
622 else tag->SetDataType(1);
623 tag->AddEventTag(*evTag);
625 lastEvent = iNumberOfEvents;
632 TString localFileName = "Run"; localFileName += tag->GetRunId();
633 localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
634 localFileName += ".ESD.tag.root";
639 fileName = localFileName.Data();
640 AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
642 else if(fStorage == 1) {
643 TString alienLocation = "/alien";
644 alienLocation += gGrid->Pwd();
645 alienLocation += fgridpath.Data();
646 alienLocation += "/";
647 alienLocation += localFileName;
648 alienLocation += "?se=";
649 alienLocation += fSE.Data();
650 fileName = alienLocation.Data();
651 AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
654 TFile* ftag = TFile::Open(fileName, "recreate");
665 //_____________________________________________________________________________
666 /*void AliESDTagCreator::CreateTag(TFile* file, const char *filepath, Int_t Counter) {
667 //private method that creates tag files
671 Double_t fMUONMASS = 0.105658369;
674 Double_t fThetaX, fThetaY, fPyz, fChisquare;
675 Double_t fPxRec, fPyRec, fPzRec, fEnergy;
677 TLorentzVector fEPvector;
679 Float_t fZVertexCut = 10.0;
680 Float_t fRhoVertexCut = 2.0;
682 Float_t fLowPtCut = 1.0;
683 Float_t fHighPtCut = 3.0;
684 Float_t fVeryHighPtCut = 10.0;
687 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
689 // Creates the tags for all the events in a given ESD file
690 Bool_t fIsSim = kTRUE;
692 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
693 Int_t nPos, nNeg, nNeutr;
694 Int_t nK0s, nNeutrons, nPi0s, nGamas;
695 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
696 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
697 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
698 Float_t maxPt = .0, meanPt = .0, totalP = .0;
700 Int_t iRunNumber = 0;
703 AliRunTag *tag = new AliRunTag();
704 AliEventTag *evTag = new AliEventTag();
705 TTree ttag("T","A Tree with event tags");
706 TBranch * btag = ttag.Branch("AliTAG", &tag);
707 btag->SetCompressionLevel(9);
709 AliInfo(Form("Creating the tags......."));
711 Int_t firstEvent = 0,lastEvent = 0;
713 TTree *t = (TTree*) file->Get("esdTree");
714 AliESDEvent *esd = new AliESDEvent();
715 esd->ReadFromTree(t);
718 Int_t iInitRunNumber = esd->GetRunNumber();
720 Int_t iNumberOfEvents = (Int_t)t->GetEntries();
721 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
749 t->GetEntry(iEventNumber);
750 iRunNumber = esd->GetRunNumber();
751 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD!!!");
752 const AliESDVertex * vertexIn = esd->GetVertex();
753 fVertexName = vertexIn->GetName();
754 if(fVertexName == "default") fVertexflag = 0;
756 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
757 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
758 if(esdTrack->GetLabel() != 0) fIsSim = kTRUE;
759 else if(esdTrack->GetLabel() == 0) fIsSim = kFALSE;
760 UInt_t status = esdTrack->GetStatus();
762 //select only tracks with ITS refit
763 if ((status&AliESDtrack::kITSrefit)==0) continue;
764 //select only tracks with TPC refit
765 if ((status&AliESDtrack::kTPCrefit)==0) continue;
767 //select only tracks with the "combined PID"
768 if ((status&AliESDtrack::kESDpid)==0) continue;
770 esdTrack->GetPxPyPz(p);
771 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
772 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
775 if(fPt > maxPt) maxPt = fPt;
777 if(esdTrack->GetSign() > 0) {
779 if(fPt > fLowPtCut) nCh1GeV++;
780 if(fPt > fHighPtCut) nCh3GeV++;
781 if(fPt > fVeryHighPtCut) nCh10GeV++;
783 if(esdTrack->GetSign() < 0) {
785 if(fPt > fLowPtCut) nCh1GeV++;
786 if(fPt > fHighPtCut) nCh3GeV++;
787 if(fPt > fVeryHighPtCut) nCh10GeV++;
789 if(esdTrack->GetSign() == 0) nNeutr++;
793 esdTrack->GetESDpid(prob);
796 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
797 if(rcc == 0.0) continue;
800 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
803 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
805 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
807 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
809 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
811 if(fPt > fLowPtCut) nEl1GeV++;
812 if(fPt > fHighPtCut) nEl3GeV++;
813 if(fPt > fVeryHighPtCut) nEl10GeV++;
821 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
822 // loop over all reconstructed tracks (also first track of combination)
823 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
824 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
825 if (muonTrack == 0x0) continue;
827 // Coordinates at vertex
828 fZ = muonTrack->GetZ();
829 fY = muonTrack->GetBendingCoor();
830 fX = muonTrack->GetNonBendingCoor();
832 fThetaX = muonTrack->GetThetaX();
833 fThetaY = muonTrack->GetThetaY();
835 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
836 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
837 fPxRec = fPzRec * TMath::Tan(fThetaX);
838 fPyRec = fPzRec * TMath::Tan(fThetaY);
839 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
841 //ChiSquare of the track if needed
842 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
843 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
844 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
846 // total number of muons inside a vertex cut
847 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
849 if(fEPvector.Pt() > fLowPtCut) {
851 if(fEPvector.Pt() > fHighPtCut) {
853 if (fEPvector.Pt() > fVeryHighPtCut) {
861 // Fill the event tags
862 if(ntrack != 0) meanPt = meanPt/ntrack;
864 evTag->SetEventId(iEventNumber+1);
865 evTag->SetPath(filepath);
867 evTag->SetVertexX(vertexIn->GetXv());
868 evTag->SetVertexY(vertexIn->GetYv());
869 evTag->SetVertexZ(vertexIn->GetZv());
870 evTag->SetVertexZError(vertexIn->GetZRes());
871 evTag->SetVertexFlag(fVertexflag);
873 evTag->SetT0VertexZ(esd->GetT0zVertex());
875 evTag->SetTriggerMask(esd->GetTriggerMask());
876 evTag->SetTriggerCluster(esd->GetTriggerCluster());
878 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
879 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
880 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
881 evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
882 evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
883 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
886 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
887 evTag->SetNumOfPosTracks(nPos);
888 evTag->SetNumOfNegTracks(nNeg);
889 evTag->SetNumOfNeutrTracks(nNeutr);
891 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
892 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
893 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
894 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
896 evTag->SetNumOfProtons(nProtons);
897 evTag->SetNumOfKaons(nKaons);
898 evTag->SetNumOfPions(nPions);
899 evTag->SetNumOfMuons(nMuons);
900 evTag->SetNumOfElectrons(nElectrons);
901 evTag->SetNumOfPhotons(nGamas);
902 evTag->SetNumOfPi0s(nPi0s);
903 evTag->SetNumOfNeutrons(nNeutrons);
904 evTag->SetNumOfKaon0s(nK0s);
906 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
907 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
908 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
909 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
910 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
911 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
912 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
913 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
914 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
916 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
917 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
919 evTag->SetTotalMomentum(totalP);
920 evTag->SetMeanPt(meanPt);
921 evTag->SetMaxPt(maxPt);
923 tag->SetRunId(iInitRunNumber);
924 if(fIsSim) tag->SetDataType(0);
925 else tag->SetDataType(1);
926 tag->AddEventTag(*evTag);
928 lastEvent = iNumberOfEvents;
935 TString localFileName = "Run"; localFileName += tag->GetRunId();
936 localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
937 localFileName += ".ESD.tag.root";
942 fileName = localFileName.Data();
943 AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
945 else if(fStorage == 1) {
946 TString alienLocation = "/alien";
947 alienLocation += gGrid->Pwd();
948 alienLocation += fgridpath.Data();
949 alienLocation += "/";
950 alienLocation += localFileName;
951 alienLocation += "?se=";
952 alienLocation += fSE.Data();
953 fileName = alienLocation.Data();
954 AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
957 TFile* ftag = TFile::Open(fileName, "recreate");