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 // 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 //-----------------------------------------------------------------
23 #include <Riostream.h>
29 #include <TLorentzVector.h>
33 #include <TGridResult.h>
36 #include "AliRunTag.h"
37 #include "AliEventTag.h"
39 #include "AliESDEvent.h"
40 #include "AliESDVertex.h"
43 #include "AliAODEvent.h"
44 #include "AliAODVertex.h"
45 #include "AliAODTrack.h"
47 #include "AliTagCreator.h"
50 ClassImp(AliTagCreator)
53 //______________________________________________________________________________
54 AliTagCreator::AliTagCreator() :
56 fSE("ALICE::CERN::se"),
60 //==============Default constructor for a AliTagCreator==================
63 //______________________________________________________________________________
64 AliTagCreator::~AliTagCreator() {
65 //================Default destructor for a AliTagCreator=======================
68 //______________________________________________________________________________
69 void AliTagCreator::SetStorage(Int_t storage) {
70 // Sets correctly the storage: 0 for local, 1 for GRID
73 AliInfo(Form("Tags will be stored locally...."));
75 AliInfo(Form("Tags will be stored in the grid...."));
76 if((fStorage != 0)&&(fStorage != 1))
78 AliInfo(Form("Storage was not properly set!!!"));
84 //______________________________________________________________________________
85 Bool_t AliTagCreator::ReadGridCollection(TGridResult *fresult) {
86 // Reads the entry of the TGridResult and creates the tags
87 Int_t nEntries = fresult->GetEntries();
96 for(Int_t i = 0; i < nEntries; i++) {
97 alienUrl = fresult->GetKey(i,"turl");
98 guid = fresult->GetKey(i,"guid");
99 if(fresult->GetKey(i,"size")) size = atol (fresult->GetKey(i,"size"));
100 md5 = fresult->GetKey(i,"md5");
101 turl = fresult->GetKey(i,"turl");
102 if(md5 && !strlen(guid)) md5 = 0;
103 if(guid && !strlen(guid)) guid = 0;
105 TFile *f = TFile::Open(alienUrl,"READ");
106 CreateTag(f,guid,md5,turl,size,counter);
115 //______________________________________________________________________________
116 Bool_t AliTagCreator::ReadLocalCollection(const char *localpath) {
117 // Checks the different subdirs of the given local path and in the
118 // case where it finds an AliESDs.root file it creates the tags
120 void *dira = gSystem->OpenDirectory(localpath);
122 const char * dirname = 0x0;
123 const char * filename = 0x0;
124 const char * pattern = "AliESDs.root";
127 while((dirname = gSystem->GetDirEntry(dira))) {
128 sprintf(fPath,"%s/%s",localpath,dirname);
129 void *dirb = gSystem->OpenDirectory(fPath);
130 while((filename = gSystem->GetDirEntry(dirb))) {
131 if(strstr(filename,pattern)) {
132 TString fESDFileName;
133 fESDFileName = fPath;
135 fESDFileName += pattern;
136 TFile *f = TFile::Open(fESDFileName,"READ");
137 CreateTag(f,fESDFileName,counter);
143 }//child directory's entry loop
144 }//parent directory's entry loop
149 //______________________________________________________________________________
150 Bool_t AliTagCreator::ReadCAFCollection(const char *filename) {
151 // Temporary solution for CAF: Takes as an input the ascii file that
152 // lists the ESDs stored in the SE of the CAF and creates the tags.
154 // Open the input stream
160 // Read the input list of files and add them to the chain
163 if (!esdfile.Contains("root")) continue; // protection
164 TFile *f = TFile::Open(esdfile,"READ");
165 CreateTag(f,esdfile,counter);
176 //__________________________________________________________________________
177 Bool_t AliTagCreator::MergeTags() {
178 //Merges the tags and stores the merged tag file
179 //locally if fStorage=0 or in the grid if fStorage=1
180 AliInfo(Form("Merging tags....."));
181 TChain *fgChain = new TChain("T");
184 const char * tagPattern = "tag";
185 // Open the working directory
186 void * dirp = gSystem->OpenDirectory(gSystem->pwd());
187 const char * name = 0x0;
188 // Add all files matching *pattern* to the chain
189 while((name = gSystem->GetDirEntry(dirp))) {
190 if (strstr(name,tagPattern)) fgChain->Add(name);
192 AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
195 else if(fStorage == 1) {
196 TString alienLocation = gGrid->Pwd();
197 alienLocation += fgridpath.Data();
198 alienLocation += "/";
200 TGridResult *tagresult = gGrid->Query(alienLocation,"*tag.root","","");
201 Int_t nEntries = tagresult->GetEntries();
202 for(Int_t i = 0; i < nEntries; i++) {
203 TString alienUrl = tagresult->GetKey(i,"turl");
204 fgChain->Add(alienUrl);
206 AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
209 AliRunTag *tag = new AliRunTag;
210 fgChain->SetBranchAddress("AliTAG",&tag);
211 fgChain->GetEntry(0);
212 TString localFileName = "Run"; localFileName += tag->GetRunId();
213 localFileName += ".Merged"; localFileName += ".ESD.tag.root";
215 TString filename = 0x0;
218 filename = localFileName.Data();
219 AliInfo(Form("Writing merged tags to local file: %s",filename.Data()));
221 else if(fStorage == 1) {
222 TString alienFileName = "/alien";
223 alienFileName += gGrid->Pwd();
224 alienFileName += fgridpath.Data();
225 alienFileName += "/";
226 alienFileName += localFileName;
227 alienFileName += "?se=";
228 alienFileName += fSE.Data();
229 filename = alienFileName.Data();
230 AliInfo(Form("Writing merged tags to grid file: %s",filename.Data()));
233 fgChain->Merge(filename);
238 //__________________________________________________________________________
239 Bool_t AliTagCreator::MergeTags(TGridResult *result) {
240 //Merges the tags that are listed in the TGridResult
241 AliInfo(Form("Merging tags....."));
242 TChain *fgChain = new TChain("T");
244 Int_t nEntries = result->GetEntries();
247 for(Int_t i = 0; i < nEntries; i++) {
248 alienUrl = result->GetKey(i,"turl");
249 fgChain->Add(alienUrl);
251 AliInfo(Form("Chained tag files: %d",fgChain->GetEntries()));
252 AliRunTag *tag = new AliRunTag;
253 fgChain->SetBranchAddress("AliTAG",&tag);
254 fgChain->GetEntry(0);
256 TString localFileName = "Run"; localFileName += tag->GetRunId();
257 localFileName += ".Merged"; localFileName += ".ESD.tag.root";
259 TString filename = 0x0;
262 filename = localFileName.Data();
263 AliInfo(Form("Writing merged tags to local file: %s",filename.Data()));
265 else if(fStorage == 1) {
266 TString alienFileName = "/alien";
267 alienFileName += gGrid->Pwd();
268 alienFileName += fgridpath.Data();
269 alienFileName += "/";
270 alienFileName += localFileName;
271 alienFileName += "?se=";
272 alienFileName += fSE.Data();
273 filename = alienFileName.Data();
274 AliInfo(Form("Writing merged tags to grid file: %s",filename.Data()));
277 fgChain->Merge(filename);
282 //_____________________________________________________________________________
283 void AliTagCreator::CreateTag(TFile* file, const char *guid, const char *md5, const char *turl, Long64_t size, Int_t Counter) {
284 //private method that creates tag files
285 TString fguid = guid;
287 TString fturl = turl;
291 Double_t fMUONMASS = 0.105658369;
294 Double_t fThetaX, fThetaY, fPyz, fChisquare;
295 Double_t fPxRec, fPyRec, fPzRec, fEnergy;
297 TLorentzVector fEPvector;
299 Float_t fZVertexCut = 10.0;
300 Float_t fRhoVertexCut = 2.0;
302 Float_t fLowPtCut = 1.0;
303 Float_t fHighPtCut = 3.0;
304 Float_t fVeryHighPtCut = 10.0;
307 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
309 // Creates the tags for all the events in a given ESD file
310 Bool_t fIsSim = kTRUE;
312 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
313 Int_t nPos, nNeg, nNeutr;
314 Int_t nK0s, nNeutrons, nPi0s, nGamas;
315 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
316 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
317 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
318 Float_t maxPt = .0, meanPt = .0, totalP = .0;
320 Int_t iRunNumber = 0;
323 AliRunTag *tag = new AliRunTag();
324 AliEventTag *evTag = new AliEventTag();
325 TTree ttag("T","A Tree with event tags");
326 TBranch * btag = ttag.Branch("AliTAG", &tag);
327 btag->SetCompressionLevel(9);
329 AliInfo(Form("Creating the tags......."));
331 Int_t firstEvent = 0,lastEvent = 0;
332 TTree *t = (TTree*) file->Get("esdTree");
333 AliESDEvent *esd = new AliESDEvent();
334 esd->ReadFromTree(t);
337 Int_t iInitRunNumber = esd->GetRunNumber();
339 Int_t iNumberOfEvents = (Int_t)t->GetEntries();
340 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
368 t->GetEntry(iEventNumber);
369 iRunNumber = esd->GetRunNumber();
370 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD!!!");
371 const AliESDVertex * vertexIn = esd->GetVertex();
372 fVertexName = vertexIn->GetName();
373 if(fVertexName == "default") fVertexflag = 0;
375 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
376 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
377 if(esdTrack->GetLabel() != 0) fIsSim = kTRUE;
378 else if(esdTrack->GetLabel() == 0) fIsSim = kFALSE;
379 UInt_t status = esdTrack->GetStatus();
381 //select only tracks with ITS refit
382 if ((status&AliESDtrack::kITSrefit)==0) continue;
383 //select only tracks with TPC refit
384 if ((status&AliESDtrack::kTPCrefit)==0) continue;
386 //select only tracks with the "combined PID"
387 if ((status&AliESDtrack::kESDpid)==0) continue;
389 esdTrack->GetPxPyPz(p);
390 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
391 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
394 if(fPt > maxPt) maxPt = fPt;
396 if(esdTrack->GetSign() > 0) {
398 if(fPt > fLowPtCut) nCh1GeV++;
399 if(fPt > fHighPtCut) nCh3GeV++;
400 if(fPt > fVeryHighPtCut) nCh10GeV++;
402 if(esdTrack->GetSign() < 0) {
404 if(fPt > fLowPtCut) nCh1GeV++;
405 if(fPt > fHighPtCut) nCh3GeV++;
406 if(fPt > fVeryHighPtCut) nCh10GeV++;
408 if(esdTrack->GetSign() == 0) nNeutr++;
412 esdTrack->GetESDpid(prob);
415 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
416 if(rcc == 0.0) continue;
419 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
422 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
424 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
426 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
428 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
430 if(fPt > fLowPtCut) nEl1GeV++;
431 if(fPt > fHighPtCut) nEl3GeV++;
432 if(fPt > fVeryHighPtCut) nEl10GeV++;
440 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
441 // loop over all reconstructed tracks (also first track of combination)
442 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
443 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
444 if (muonTrack == 0x0) continue;
446 // Coordinates at vertex
447 fZ = muonTrack->GetZ();
448 fY = muonTrack->GetBendingCoor();
449 fX = muonTrack->GetNonBendingCoor();
451 fThetaX = muonTrack->GetThetaX();
452 fThetaY = muonTrack->GetThetaY();
454 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
455 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
456 fPxRec = fPzRec * TMath::Tan(fThetaX);
457 fPyRec = fPzRec * TMath::Tan(fThetaY);
458 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
460 //ChiSquare of the track if needed
461 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
462 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
463 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
465 // total number of muons inside a vertex cut
466 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
468 if(fEPvector.Pt() > fLowPtCut) {
470 if(fEPvector.Pt() > fHighPtCut) {
472 if (fEPvector.Pt() > fVeryHighPtCut) {
480 // Fill the event tags
481 if(ntrack != 0) meanPt = meanPt/ntrack;
483 evTag->SetEventId(iEventNumber+1);
484 evTag->SetGUID(fguid);
486 evTag->SetTURL(fturl);
487 evTag->SetSize(size);
488 evTag->SetVertexX(vertexIn->GetXv());
489 evTag->SetVertexY(vertexIn->GetYv());
490 evTag->SetVertexZ(vertexIn->GetZv());
491 evTag->SetVertexZError(vertexIn->GetZRes());
492 evTag->SetVertexFlag(fVertexflag);
494 evTag->SetT0VertexZ(esd->GetT0zVertex());
496 evTag->SetTriggerMask(esd->GetTriggerMask());
497 evTag->SetTriggerCluster(esd->GetTriggerCluster());
499 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
500 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
501 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
502 evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
503 evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
504 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
507 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
508 evTag->SetNumOfPosTracks(nPos);
509 evTag->SetNumOfNegTracks(nNeg);
510 evTag->SetNumOfNeutrTracks(nNeutr);
512 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
513 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
514 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
515 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
517 evTag->SetNumOfProtons(nProtons);
518 evTag->SetNumOfKaons(nKaons);
519 evTag->SetNumOfPions(nPions);
520 evTag->SetNumOfMuons(nMuons);
521 evTag->SetNumOfElectrons(nElectrons);
522 evTag->SetNumOfPhotons(nGamas);
523 evTag->SetNumOfPi0s(nPi0s);
524 evTag->SetNumOfNeutrons(nNeutrons);
525 evTag->SetNumOfKaon0s(nK0s);
527 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
528 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
529 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
530 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
531 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
532 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
533 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
534 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
535 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
537 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
538 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
540 evTag->SetTotalMomentum(totalP);
541 evTag->SetMeanPt(meanPt);
542 evTag->SetMaxPt(maxPt);
544 tag->SetRunId(iInitRunNumber);
545 if(fIsSim) tag->SetDataType(0);
546 else tag->SetDataType(1);
547 tag->AddEventTag(*evTag);
549 lastEvent = iNumberOfEvents;
556 TString localFileName = "Run"; localFileName += tag->GetRunId();
557 localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
558 localFileName += ".ESD.tag.root";
563 fileName = localFileName.Data();
564 AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
566 else if(fStorage == 1) {
567 TString alienLocation = "/alien";
568 alienLocation += gGrid->Pwd();
569 alienLocation += fgridpath.Data();
570 alienLocation += "/";
571 alienLocation += localFileName;
572 alienLocation += "?se=";
573 alienLocation += fSE.Data();
574 fileName = alienLocation.Data();
575 AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
578 TFile* ftag = TFile::Open(fileName, "recreate");
589 //_____________________________________________________________________________
590 void AliTagCreator::CreateTag(TFile* file, const char *filepath, Int_t Counter) {
591 //private method that creates tag files
595 Double_t fMUONMASS = 0.105658369;
598 Double_t fThetaX, fThetaY, fPyz, fChisquare;
599 Double_t fPxRec, fPyRec, fPzRec, fEnergy;
601 TLorentzVector fEPvector;
603 Float_t fZVertexCut = 10.0;
604 Float_t fRhoVertexCut = 2.0;
606 Float_t fLowPtCut = 1.0;
607 Float_t fHighPtCut = 3.0;
608 Float_t fVeryHighPtCut = 10.0;
611 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
613 // Creates the tags for all the events in a given ESD file
614 Bool_t fIsSim = kTRUE;
616 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
617 Int_t nPos, nNeg, nNeutr;
618 Int_t nK0s, nNeutrons, nPi0s, nGamas;
619 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
620 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
621 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
622 Float_t maxPt = .0, meanPt = .0, totalP = .0;
624 Int_t iRunNumber = 0;
627 AliRunTag *tag = new AliRunTag();
628 AliEventTag *evTag = new AliEventTag();
629 TTree ttag("T","A Tree with event tags");
630 TBranch * btag = ttag.Branch("AliTAG", &tag);
631 btag->SetCompressionLevel(9);
633 AliInfo(Form("Creating the tags......."));
635 Int_t firstEvent = 0,lastEvent = 0;
637 TTree *t = (TTree*) file->Get("esdTree");
638 AliESDEvent *esd = new AliESDEvent();
639 esd->ReadFromTree(t);
642 Int_t iInitRunNumber = esd->GetRunNumber();
644 Int_t iNumberOfEvents = (Int_t)t->GetEntries();
645 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
673 t->GetEntry(iEventNumber);
674 iRunNumber = esd->GetRunNumber();
675 if(iRunNumber != iInitRunNumber) AliFatal("Inconsistency of run numbers in the AliESD!!!");
676 const AliESDVertex * vertexIn = esd->GetVertex();
677 fVertexName = vertexIn->GetName();
678 if(fVertexName == "default") fVertexflag = 0;
680 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
681 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
682 if(esdTrack->GetLabel() != 0) fIsSim = kTRUE;
683 else if(esdTrack->GetLabel() == 0) fIsSim = kFALSE;
684 UInt_t status = esdTrack->GetStatus();
686 //select only tracks with ITS refit
687 if ((status&AliESDtrack::kITSrefit)==0) continue;
688 //select only tracks with TPC refit
689 if ((status&AliESDtrack::kTPCrefit)==0) continue;
691 //select only tracks with the "combined PID"
692 if ((status&AliESDtrack::kESDpid)==0) continue;
694 esdTrack->GetPxPyPz(p);
695 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
696 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
699 if(fPt > maxPt) maxPt = fPt;
701 if(esdTrack->GetSign() > 0) {
703 if(fPt > fLowPtCut) nCh1GeV++;
704 if(fPt > fHighPtCut) nCh3GeV++;
705 if(fPt > fVeryHighPtCut) nCh10GeV++;
707 if(esdTrack->GetSign() < 0) {
709 if(fPt > fLowPtCut) nCh1GeV++;
710 if(fPt > fHighPtCut) nCh3GeV++;
711 if(fPt > fVeryHighPtCut) nCh10GeV++;
713 if(esdTrack->GetSign() == 0) nNeutr++;
717 esdTrack->GetESDpid(prob);
720 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
721 if(rcc == 0.0) continue;
724 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
727 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
729 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
731 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
733 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
735 if(fPt > fLowPtCut) nEl1GeV++;
736 if(fPt > fHighPtCut) nEl3GeV++;
737 if(fPt > fVeryHighPtCut) nEl10GeV++;
745 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
746 // loop over all reconstructed tracks (also first track of combination)
747 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
748 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
749 if (muonTrack == 0x0) continue;
751 // Coordinates at vertex
752 fZ = muonTrack->GetZ();
753 fY = muonTrack->GetBendingCoor();
754 fX = muonTrack->GetNonBendingCoor();
756 fThetaX = muonTrack->GetThetaX();
757 fThetaY = muonTrack->GetThetaY();
759 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
760 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
761 fPxRec = fPzRec * TMath::Tan(fThetaX);
762 fPyRec = fPzRec * TMath::Tan(fThetaY);
763 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
765 //ChiSquare of the track if needed
766 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
767 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
768 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
770 // total number of muons inside a vertex cut
771 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
773 if(fEPvector.Pt() > fLowPtCut) {
775 if(fEPvector.Pt() > fHighPtCut) {
777 if (fEPvector.Pt() > fVeryHighPtCut) {
785 // Fill the event tags
786 if(ntrack != 0) meanPt = meanPt/ntrack;
788 evTag->SetEventId(iEventNumber+1);
789 evTag->SetPath(filepath);
791 evTag->SetVertexX(vertexIn->GetXv());
792 evTag->SetVertexY(vertexIn->GetYv());
793 evTag->SetVertexZ(vertexIn->GetZv());
794 evTag->SetVertexZError(vertexIn->GetZRes());
795 evTag->SetVertexFlag(fVertexflag);
797 evTag->SetT0VertexZ(esd->GetT0zVertex());
799 evTag->SetTriggerMask(esd->GetTriggerMask());
800 evTag->SetTriggerCluster(esd->GetTriggerCluster());
802 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
803 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
804 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
805 evTag->SetZDCNeutron1Energy(esd->GetZDCN2Energy());
806 evTag->SetZDCProton1Energy(esd->GetZDCP2Energy());
807 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
810 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
811 evTag->SetNumOfPosTracks(nPos);
812 evTag->SetNumOfNegTracks(nNeg);
813 evTag->SetNumOfNeutrTracks(nNeutr);
815 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
816 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
817 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
818 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
820 evTag->SetNumOfProtons(nProtons);
821 evTag->SetNumOfKaons(nKaons);
822 evTag->SetNumOfPions(nPions);
823 evTag->SetNumOfMuons(nMuons);
824 evTag->SetNumOfElectrons(nElectrons);
825 evTag->SetNumOfPhotons(nGamas);
826 evTag->SetNumOfPi0s(nPi0s);
827 evTag->SetNumOfNeutrons(nNeutrons);
828 evTag->SetNumOfKaon0s(nK0s);
830 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
831 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
832 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
833 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
834 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
835 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
836 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
837 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
838 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
840 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
841 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
843 evTag->SetTotalMomentum(totalP);
844 evTag->SetMeanPt(meanPt);
845 evTag->SetMaxPt(maxPt);
847 tag->SetRunId(iInitRunNumber);
848 if(fIsSim) tag->SetDataType(0);
849 else tag->SetDataType(1);
850 tag->AddEventTag(*evTag);
852 lastEvent = iNumberOfEvents;
859 TString localFileName = "Run"; localFileName += tag->GetRunId();
860 localFileName += ".Event"; localFileName += firstEvent; localFileName += "_"; localFileName += lastEvent; localFileName += "."; localFileName += Counter;
861 localFileName += ".ESD.tag.root";
866 fileName = localFileName.Data();
867 AliInfo(Form("Writing tags to local file: %s",fileName.Data()));
869 else if(fStorage == 1) {
870 TString alienLocation = "/alien";
871 alienLocation += gGrid->Pwd();
872 alienLocation += fgridpath.Data();
873 alienLocation += "/";
874 alienLocation += localFileName;
875 alienLocation += "?se=";
876 alienLocation += fSE.Data();
877 fileName = alienLocation.Data();
878 AliInfo(Form("Writing tags to grid file: %s",fileName.Data()));
881 TFile* ftag = TFile::Open(fileName, "recreate");
892 //__________________________________________________________________________
893 void AliTagCreator::CreateAODTags(Int_t fFirstEvent, Int_t fLastEvent) {
894 //creates tag files for AODs
896 Float_t fLowPtCut = 1.0;
897 Float_t fHighPtCut = 3.0;
898 Float_t fVeryHighPtCut = 10.0;
901 Double_t partFrac[10] = {0.01, 0.01, 0.85, 0.10, 0.05, 0., 0., 0., 0., 0.};
903 // Creates the tags for all the events in a given ESD file
905 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
906 Int_t nPos, nNeg, nNeutr;
907 Int_t nKinks, nV0s, nCascades;
908 Int_t nK0s, nNeutrons, nPi0s, nGamas;
909 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
910 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
911 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
912 Float_t maxPt = .0, meanPt = .0, totalP = .0;
914 AliRunTag *tag = new AliRunTag();
915 TTree ttag("T","A Tree with event tags");
916 TBranch * btag = ttag.Branch("AliTAG", &tag);
917 btag->SetCompressionLevel(9);
919 //reading the esd tag file
920 TChain *oldTagTree = new TChain("T");
921 const char * tagPattern = "ESD.tag";
922 // Open the working directory
923 void * dirp = gSystem->OpenDirectory(gSystem->pwd());
924 const char * name = 0x0;
925 // Add all files matching *pattern* to the chain
926 while((name = gSystem->GetDirEntry(dirp))) {
927 if (strstr(name,tagPattern)) oldTagTree->Add(name);
929 AliInfo(Form("Chained tag files: %d",oldTagTree->GetEntries()));
931 //reading the esd tag file
932 AliRunTag *oldtag = new AliRunTag();
934 oldTagTree->SetBranchAddress("AliTAG",&oldtag);
935 oldTagTree->GetEntry(0);
936 tag->CopyStandardContent(oldtag);
937 const TClonesArray *evTagList = oldtag->GetEventTags();
939 cout<<"Creating the tags......."<<endl;
941 TFile *file = TFile::Open("AliAOD.root");
942 if (!file || !file->IsOpen()) {
943 AliError(Form("opening failed"));
947 TTree *aodTree = (TTree*)file->Get("AOD");
948 AliAODEvent *aod = (AliAODEvent*)aodTree->GetUserInfo()->FindObject("AliAODEvent");
949 TIter next(aod->GetList());
953 if(fLastEvent == -1) lastEvent = (Int_t)aodTree->GetEntries();
954 else lastEvent = fLastEvent;
956 while((el=(TNamed*)next()))
957 aodTree->SetBranchAddress(el->GetName(),aod->GetList()->GetObjectRef(el));
960 Int_t nEvents = aodTree->GetEntries();
961 for (Int_t iEventNumber = 0; iEventNumber < nEvents; iEventNumber++) {
962 AliEventTag *evTag = (AliEventTag *)evTagList->At(iEventNumber);
964 nPos = 0; nNeg = 0; nNeutr =0;
965 nKinks = 0; nV0s = 0; nCascades = 0;
966 nK0s = 0; nNeutrons = 0; nPi0s = 0; nGamas = 0;
967 nProtons = 0; nKaons = 0; nPions = 0; nMuons = 0; nElectrons = 0;
968 nCh1GeV = 0; nCh3GeV = 0; nCh10GeV = 0;
969 nMu1GeV = 0; nMu3GeV = 0; nMu10GeV = 0;
970 nEl1GeV = 0; nEl3GeV = 0; nEl10GeV = 0;
971 maxPt = .0; meanPt = .0; totalP = .0;
974 aodTree->GetEvent(iEventNumber);
977 aod->GetStdContent();
979 Int_t nTracks = aod->GetNTracks();
980 cout << "Event: " << iEventNumber+1 << "/" << nEvents <<" - Tracks: "<<nTracks<< endl;
981 // loop over vertices
982 Int_t nVtxs = aod->GetNVertices();
983 for (Int_t nVtx = 0; nVtx < nVtxs; nVtx++) {
985 AliAODVertex *vertex = aod->GetVertex(nVtx);
986 if(vertex->GetType() == 1) nKinks += 1;
987 if(vertex->GetType() == 2) nV0s += 1;
988 if(vertex->GetType() == 3) nCascades += 1;
990 for (Int_t nTr = 0; nTr < nTracks; nTr++) {
991 AliAODTrack *track = aod->GetTrack(nTr);
993 Double_t fPt = track->Pt();
994 if(fPt > maxPt) maxPt = fPt;
995 if(track->Charge() > 0) {
997 if(fPt > fLowPtCut) nCh1GeV++;
998 if(fPt > fHighPtCut) nCh3GeV++;
999 if(fPt > fVeryHighPtCut) nCh10GeV++;
1001 if(track->Charge() < 0) {
1003 if(fPt > fLowPtCut) nCh1GeV++;
1004 if(fPt > fHighPtCut) nCh3GeV++;
1005 if(fPt > fVeryHighPtCut) nCh10GeV++;
1007 if(track->Charge() == 0) nNeutr++;
1010 const Double32_t *prob = track->PID();
1012 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
1013 if(rcc == 0.0) continue;
1016 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
1019 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
1021 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
1023 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
1025 if ((w[1]>w[4])&&(w[1]>w[3])&&(w[1]>w[2])&&(w[1]>w[0])) {
1027 if(fPt > fLowPtCut) nMu1GeV++;
1028 if(fPt > fHighPtCut) nMu3GeV++;
1029 if(fPt > fVeryHighPtCut) nMu10GeV++;
1032 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
1034 if(fPt > fLowPtCut) nEl1GeV++;
1035 if(fPt > fHighPtCut) nEl3GeV++;
1036 if(fPt > fVeryHighPtCut) nEl10GeV++;
1039 totalP += track->P();
1043 // Fill the event tags
1045 meanPt = meanPt/ntrack;
1047 evTag->SetEventId(iEventNumber+1);
1049 evTag->SetNumOfTracks(nTracks);
1050 evTag->SetNumOfPosTracks(nPos);
1051 evTag->SetNumOfNegTracks(nNeg);
1052 evTag->SetNumOfNeutrTracks(nNeutr);
1054 evTag->SetNumOfV0s(nV0s);
1055 evTag->SetNumOfCascades(nCascades);
1056 evTag->SetNumOfKinks(nKinks);
1058 evTag->SetNumOfProtons(nProtons);
1059 evTag->SetNumOfKaons(nKaons);
1060 evTag->SetNumOfPions(nPions);
1061 evTag->SetNumOfMuons(nMuons);
1062 evTag->SetNumOfElectrons(nElectrons);
1063 evTag->SetNumOfPhotons(nGamas);
1064 evTag->SetNumOfPi0s(nPi0s);
1065 evTag->SetNumOfNeutrons(nNeutrons);
1066 evTag->SetNumOfKaon0s(nK0s);
1068 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
1069 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
1070 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
1071 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
1072 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
1073 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
1074 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
1075 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
1076 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
1078 evTag->SetTotalMomentum(totalP);
1079 evTag->SetMeanPt(meanPt);
1080 evTag->SetMaxPt(maxPt);
1081 tag->AddEventTag(*evTag);
1083 if(fLastEvent == -1) lastEvent = (Int_t)aodTree->GetEntries();
1084 else lastEvent = fLastEvent;
1090 sprintf(fileName, "Run%d.Event%d_%d.AOD.tag.root",
1091 tag->GetRunId(),fFirstEvent,lastEvent );
1092 cout<<"writing tags to file "<<fileName<<endl;
1093 AliDebug(1, Form("writing tags to file %s", fileName));
1095 TFile* ftag = TFile::Open(fileName, "recreate");