0d627872a293422a647d7f85563835c0a70ecbd3
[u/mrichter/AliRoot.git] / STEER / AliTagCreator.cxx
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
23 #include <TObject.h>
24 #include <Riostream.h>
25 #include <TSystem.h>
26 #include <TChain.h>
27 #include <TFile.h>
28 #include <TString.h>
29
30 //ROOT-AliEn
31 #include <TGrid.h>
32 #include <TGridResult.h>
33
34 //AliRoot
35 #include "AliRunTag.h"
36 #include "AliEventTag.h"
37 #include "AliESD.h"
38 #include "AliESDVertex.h"
39 #include "AliPID.h"
40 #include "AliESDpid.h"
41 #include "AliLog.h"
42
43
44 #include "AliTagCreator.h"
45
46
47 ClassImp(AliTagCreator)
48
49
50 //______________________________________________________________________________
51 AliTagCreator::AliTagCreator() //local mode
52 {
53   //==============Default constructor for a AliTagCreator==================
54   fgridpath = "";
55   fUser = "";
56   fPasswd = "";  
57   fSE = "";   
58   fHost = "";
59   fPort = 0; 
60   fStorage = 0; 
61   fresult = 0;
62 }
63
64 //______________________________________________________________________________
65 AliTagCreator::AliTagCreator(const char *host, Int_t port, const char *username)
66 {
67   //==============Default constructor for a AliTagCreator==================
68   fStorage = 0; 
69   fresult = 0;
70   fgridpath = "";
71   fHost = host;
72   fUser = username;
73   fPort = port;
74
75   TString grid = "alien://";
76   grid += fHost;
77   grid += ":";
78   grid += port;
79   
80   //connect to the grid
81   TGrid::Connect(grid.Data(),fUser.Data());
82   if(!gGrid)
83     {
84       AliError("Connection failed!!!");
85       return;
86     }
87 }
88
89
90 //______________________________________________________________________________
91 AliTagCreator::AliTagCreator(const char *host, Int_t port, const char *username, const char *passwd)
92 {
93   //==============Default constructor for a AliTagCreator==================
94   fStorage = 0; 
95   fresult = 0;
96   fgridpath = "";
97   fHost = host;
98   fUser = username;
99   fPasswd = passwd;
100   fPort = port;
101
102   TString grid = "alien://";
103   grid += fHost;
104   grid += ":";
105   grid += port;
106   
107   //connect to the grid
108   TGrid::Connect(grid.Data(),fUser.Data(),fPasswd.Data());
109    if(!gGrid)
110     {
111       AliError("Connection failed!!!");
112       return;
113     }
114 }
115
116
117 //______________________________________________________________________________
118 AliTagCreator::~AliTagCreator()
119 {
120 //================Default destructor for a AliTagCreator=======================
121 }
122
123 //______________________________________________________________________________
124 void AliTagCreator::SetSE(const char *se)
125 {
126   fSE = se;
127 }
128
129 //______________________________________________________________________________
130 void AliTagCreator::SetGridPath(const char *gridpath)
131 {
132   fgridpath = gridpath;
133 }
134
135 //______________________________________________________________________________
136 void AliTagCreator::SetStorage(Int_t storage)
137 {
138   fStorage = storage;
139   if(fStorage == 0)
140     AliInfo(Form("Tags will be stored locally...."));
141   if(fStorage == 1)
142     AliInfo(Form("Tags will be stored in the grid...."));
143   if((fStorage != 0)&&(fStorage != 1))
144     {
145       AliInfo(Form("Storage was not properly set!!!"));
146       abort();
147     }  
148 }
149
150 //______________________________________________________________________________
151 Bool_t AliTagCreator::ConnectToGrid(const char *host, Int_t port, const char *username)
152 {
153   fHost = host;
154   fUser = username;
155   fPort = port;
156
157   TString grid = "alien://";
158   grid += fHost;
159   grid += ":";
160   grid += port;
161   
162   //connect to the grid
163   TGrid::Connect(grid.Data(),fUser.Data());
164   if(!gGrid)
165     {
166       AliError("Connection failed!!!");
167       return kFALSE;
168     } //connect to the grid
169
170   return kTRUE;
171 }
172
173 //______________________________________________________________________________
174 Bool_t AliTagCreator::ReadESDCollection(TGridResult *res)
175 {
176   fresult = res;
177   Int_t NEntries = fresult->GetEntries();
178
179   TString gridname = "alien://";
180   TString AlienUrl;
181   const char *guid;
182  
183   Int_t counter = 0;
184   for(Int_t i = 0; i < NEntries; i++)
185     {
186       AlienUrl = gridname;
187       AlienUrl += fresult->GetKey(i,"lfn");
188       guid = fresult->GetKey(i,"guid");
189       TFile *f = TFile::Open(AlienUrl,"READ");
190       CreateTag(f,guid,counter);
191       f->Close();
192       delete f;  
193       counter += 1;
194     }//grid result loop
195
196   return kTRUE;
197 }
198
199
200 //_____________________________________________________________________________
201 void AliTagCreator::CreateTag(TFile* file, const char *guid, Int_t Counter)
202 {
203   Int_t ntrack;
204   Int_t NProtons, NKaons, NPions, NMuons, NElectrons;
205   Int_t Npos, Nneg, Nneutr;
206   Int_t NK0s, Nneutrons, Npi0s, Ngamas;
207   Int_t Nch1GeV, Nch3GeV, Nch10GeV;
208   Int_t Nmu1GeV, Nmu3GeV, Nmu10GeV;
209   Int_t Nel1GeV, Nel3GeV, Nel10GeV;
210   Float_t MaxPt = .0, MeanPt = .0, TotalP = .0;
211
212   AliRunTag *tag = new AliRunTag();
213   AliDetectorTag *detTag = new AliDetectorTag();
214   AliEventTag *evTag = new AliEventTag();
215   TTree ttag("T","A Tree with event tags");
216   TBranch * btag = ttag.Branch("AliTAG", "AliRunTag", &tag);
217   btag->SetCompressionLevel(9);
218   
219   AliInfo(Form("Creating the tags......."));    
220   
221   Int_t firstEvent = 0,lastEvent = 0;
222   TTree *t = (TTree*) file->Get("esdTree");
223   TBranch * b = t->GetBranch("ESD");
224   AliESD *esd = 0;
225   b->SetAddress(&esd);
226   
227   tag->SetRunId(esd->GetRunNumber());
228   
229   Int_t i_NumberOfEvents = b->GetEntries();
230   for (Int_t i_EventNumber = 0; i_EventNumber < i_NumberOfEvents; i_EventNumber++)
231     {
232       ntrack = 0;
233       Npos = 0;
234       Nneg = 0;
235       Nneutr =0;
236       NK0s = 0;
237       Nneutrons = 0;
238       Npi0s = 0;
239       Ngamas = 0;
240       NProtons = 0;
241       NKaons = 0;
242       NPions = 0;
243       NMuons = 0;
244       NElectrons = 0;     
245       Nch1GeV = 0;
246       Nch3GeV = 0;
247       Nch10GeV = 0;
248       Nmu1GeV = 0;
249       Nmu3GeV = 0;
250       Nmu10GeV = 0;
251       Nel1GeV = 0;
252       Nel3GeV = 0;
253       Nel10GeV = 0;
254       MaxPt = .0;
255       MeanPt = .0;
256       TotalP = .0;
257       
258       b->GetEntry(i_EventNumber);
259       const AliESDVertex * VertexIn = esd->GetVertex();
260       
261       for (Int_t i_TrackNumber = 0; i_TrackNumber < esd->GetNumberOfTracks(); i_TrackNumber++)
262         {
263           AliESDtrack * ESDTrack = esd->GetTrack(i_TrackNumber);
264           UInt_t status = ESDTrack->GetStatus();
265           
266           //select only tracks with ITS refit
267           if ((status&AliESDtrack::kITSrefit)==0) continue;
268           
269           //select only tracks with TPC refit-->remove extremely high Pt tracks
270           if ((status&AliESDtrack::kTPCrefit)==0) continue;
271           
272           //select only tracks with the "combined PID"
273           if ((status&AliESDtrack::kESDpid)==0) continue;
274           Double_t p[3];
275           ESDTrack->GetPxPyPz(p);
276           Double_t P = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
277           Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
278           TotalP += P;
279           MeanPt += fPt;
280           if(fPt > MaxPt)
281             MaxPt = fPt;
282           
283           if(ESDTrack->GetSign() > 0)
284             {
285               Npos++;
286               if(fPt > 1.0)
287                 Nch1GeV++;
288               if(fPt > 3.0)
289                 Nch3GeV++;
290               if(fPt > 10.0)
291                 Nch10GeV++;
292             }
293           if(ESDTrack->GetSign() < 0)
294             {
295               Nneg++;
296               if(fPt > 1.0)
297                 Nch1GeV++;
298               if(fPt > 3.0)
299                 Nch3GeV++;
300               if(fPt > 10.0)
301                 Nch10GeV++;
302             }
303           if(ESDTrack->GetSign() == 0)
304             Nneutr++;
305           
306           //PID
307           Double_t prob[10];
308           ESDTrack->GetESDpid(prob);
309           
310           //K0s
311           if ((prob[8]>prob[7])&&(prob[8]>prob[6])&&(prob[8]>prob[5])&&(prob[8]>prob[4])&&(prob[8]>prob[3])&&(prob[8]>prob[2])&&(prob[8]>prob[1])&&(prob[8]>prob[0]))
312             NK0s++;
313           //neutrons
314           if ((prob[7]>prob[8])&&(prob[7]>prob[6])&&(prob[7]>prob[5])&&(prob[7]>prob[4])&&(prob[7]>prob[3])&&(prob[7]>prob[2])&&(prob[7]>prob[1])&&(prob[7]>prob[0]))
315             Nneutrons++; 
316           //pi0s
317           if ((prob[6]>prob[8])&&(prob[6]>prob[7])&&(prob[6]>prob[5])&&(prob[6]>prob[4])&&(prob[6]>prob[3])&&(prob[6]>prob[2])&&(prob[6]>prob[1])&&(prob[6]>prob[0]))
318             Npi0s++;
319           //gamas
320           if ((prob[5]>prob[8])&&(prob[5]>prob[7])&&(prob[5]>prob[6])&&(prob[5]>prob[4])&&(prob[5]>prob[3])&&(prob[5]>prob[2])&&(prob[5]>prob[1])&&(prob[5]>prob[0]))
321             Ngamas++;
322           //protons
323           if ((prob[4]>prob[8])&&(prob[4]>prob[7])&&(prob[4]>prob[6])&&(prob[4]>prob[5])&&(prob[4]>prob[3])&&(prob[4]>prob[2])&&(prob[4]>prob[1])&&(prob[4]>prob[0]))
324             NProtons++;
325           //kaons
326           if ((prob[3]>prob[8])&&(prob[3]>prob[7])&&(prob[3]>prob[6])&&(prob[3]>prob[5])&&(prob[3]>prob[4])&&(prob[3]>prob[2])&&(prob[3]>prob[1])&&(prob[3]>prob[0]))
327             NKaons++;
328           //pions
329           if ((prob[2]>prob[8])&&(prob[2]>prob[7])&&(prob[2]>prob[6])&&(prob[2]>prob[5])&&(prob[2]>prob[4])&&(prob[2]>prob[3])&&(prob[2]>prob[1])&&(prob[2]>prob[0]))
330             NPions++; 
331           //muons
332           if ((prob[1]>prob[8])&&(prob[1]>prob[7])&&(prob[1]>prob[6])&&(prob[1]>prob[5])&&(prob[1]>prob[4])&&(prob[1]>prob[3])&&(prob[1]>prob[2])&&(prob[1]>prob[0]))
333             {
334               NMuons++;
335               if(fPt > 1.0)
336                 Nmu1GeV++;
337               if(fPt > 3.0)
338                 Nmu3GeV++;
339               if(fPt > 10.0)
340                 Nmu10GeV++;
341             }
342           //electrons
343           if ((prob[0]>prob[8])&&(prob[0]>prob[7])&&(prob[0]>prob[6])&&(prob[0]>prob[5])&&(prob[0]>prob[4])&&(prob[0]>prob[3])&&(prob[0]>prob[2])&&(prob[0]>prob[1]))
344             {
345               NElectrons++;
346               if(fPt > 1.0)
347                 Nel1GeV++;
348               if(fPt > 3.0)
349                 Nel3GeV++;
350               if(fPt > 10.0)
351                 Nel10GeV++;
352             }
353           
354           
355           
356           ntrack++;
357         }//track loop
358       // Fill the event tags 
359       if(ntrack != 0)
360         MeanPt = MeanPt/ntrack;
361       
362       evTag->SetEventId(i_EventNumber+1);
363       evTag->SetGUID(guid);
364       evTag->SetVertexX(VertexIn->GetXv());
365       evTag->SetVertexY(VertexIn->GetYv());
366       evTag->SetVertexZ(VertexIn->GetZv());
367       
368       evTag->SetT0VertexZ(esd->GetT0zVertex());
369       
370       evTag->SetTrigger(esd->GetTrigger());
371       
372       evTag->SetZDCNeutronEnergy(esd->GetZDCNEnergy());
373       evTag->SetZDCProtonEnergy(esd->GetZDCPEnergy());
374       evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
375       evTag->SetNumOfParticipants(esd->GetZDCParticipants());
376       
377       
378       evTag->SetNumOfTracks(esd->GetNumberOfTracks());
379       evTag->SetNumOfPosTracks(Npos);
380       evTag->SetNumOfNegTracks(Nneg);
381       evTag->SetNumOfNeutrTracks(Nneutr);
382       
383       evTag->SetNumOfV0s(esd->GetNumberOfV0s());
384       evTag->SetNumOfCascades(esd->GetNumberOfCascades());
385       evTag->SetNumOfKinks(esd->GetNumberOfKinks());
386       evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
387       
388       evTag->SetNumOfProtons(NProtons);
389       evTag->SetNumOfKaons(NKaons);
390       evTag->SetNumOfPions(NPions);
391       evTag->SetNumOfMuons(NMuons);
392       evTag->SetNumOfElectrons(NElectrons);
393       evTag->SetNumOfPhotons(Ngamas);
394       evTag->SetNumOfPi0s(Npi0s);
395       evTag->SetNumOfNeutrons(Nneutrons);
396       evTag->SetNumOfKaon0s(NK0s);
397       
398       evTag->SetNumOfChargedAbove1GeV(Nch1GeV);
399       evTag->SetNumOfChargedAbove3GeV(Nch3GeV);
400       evTag->SetNumOfChargedAbove10GeV(Nch10GeV);
401       evTag->SetNumOfMuonsAbove1GeV(Nmu1GeV);
402       evTag->SetNumOfMuonsAbove3GeV(Nmu3GeV);
403       evTag->SetNumOfMuonsAbove10GeV(Nmu10GeV);
404       evTag->SetNumOfElectronsAbove1GeV(Nel1GeV);
405       evTag->SetNumOfElectronsAbove3GeV(Nel3GeV);
406       evTag->SetNumOfElectronsAbove10GeV(Nel10GeV);
407       
408       evTag->SetNumOfPHOSTracks(esd->GetNumberOfPHOSParticles());
409       evTag->SetNumOfEMCALTracks(esd->GetNumberOfEMCALParticles());
410       
411       evTag->SetTotalMomentum(TotalP);
412       evTag->SetMeanPt(MeanPt);
413       evTag->SetMaxPt(MaxPt);
414   
415       tag->AddEventTag(evTag);
416     }
417   lastEvent = i_NumberOfEvents;
418   
419   t->Delete("");
420         
421   ttag.Fill();
422   tag->Clear();
423
424   TString LocalfileName = "Run"; LocalfileName += tag->GetRunId(); 
425   LocalfileName += ".Event"; LocalfileName += firstEvent; LocalfileName += "_"; LocalfileName += lastEvent; LocalfileName += "."; LocalfileName += Counter;
426   LocalfileName += ".ESD.tag.root";
427
428   TString AlienLocation = "/alien";
429   AlienLocation += gGrid->Pwd();
430   AlienLocation += fgridpath.Data();
431   AlienLocation += "/";
432   AlienLocation +=  LocalfileName;
433   AlienLocation += "?se=";
434   AlienLocation += fSE.Data();
435
436   TString FileName;
437   
438   if(fStorage == 0)
439     {
440       FileName = LocalfileName.Data();      
441       cout<<"Writing tags to local file: "<<FileName.Data()<<endl;
442     }
443   if(fStorage == 1)
444     {
445       FileName = AlienLocation.Data();
446       cout<<"Writing tags to grid file: "<<FileName.Data()<<endl;
447     }
448
449   TFile* ftag = TFile::Open(FileName, "recreate");
450   ftag->cd();
451   ttag.Write();
452   ftag->Close();
453
454   delete ftag;
455   delete esd;
456
457   delete tag;
458   delete detTag;
459   delete evTag;
460 }
461