Suppression of some warnings
[u/mrichter/AliRoot.git] / TOF / AliTOFProb.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 //_________________________________________________________________________
17 // TTask class for TOF Probabilities (ALICE-Offline week June 2002).
18 // Use case: start root and execute the following macro
19 // ntuple.root is assumed to be the ntuple with the reconstructed
20 // data
21 /*
22 {
23   // Dynamically link some shared libs
24   if (gClassTable->GetID("AliRun") < 0) {
25     gROOT->LoadMacro("loadlibs.C");
26     loadlibs();
27   }
28   // create an instance of the AliTOFProb class
29   // You have to pass the ntuple
30  AliTOFProb* tofprob=new AliTOFProb("ntuple.root");
31  tofprob->Exec("");
32 }
33 */
34 // 
35 //-- Author: F. Pierella | pierella@bo.infn.it
36 //////////////////////////////////////////////////////////////////////////////
37
38 #include "TStyle.h"
39 #include "TTask.h"
40 #include "TTree.h"
41 #include "TSystem.h"
42 #include "TFile.h"
43 #include "TH3.h"
44 //#include "TH3.h"
45 #include "TCanvas.h"
46 #include "TFile.h"
47 #include <TF1.h>
48 #include <TF2.h>
49 #include "TROOT.h"
50 #include "TFolder.h"
51 #include "TNtuple.h"
52 #include "TLeaf.h"
53
54 #include "AliConst.h"
55 #include "AliTOFConstants.h"
56 #include "AliTOFProb.h"
57
58 #include <stdlib.h>
59 #include <Riostream.h>
60 #include <Riostream.h>
61
62 ClassImp(AliTOFProb)
63
64 //____________________________________________________________________________ 
65   AliTOFProb::AliTOFProb():TTask("AliTOFProb","") 
66 {
67   // default ctor - set the pointer member vars to zero
68
69   fhfile    = 0;
70   fNtuple   = 0;
71   fgen      = 0;
72   foutfileName  = 0;
73 }
74            
75 //____________________________________________________________________________ 
76   AliTOFProb::AliTOFProb(char* headerFile):TTask("AliTOFProb","") 
77 {
78
79   fhfile = TFile::Open(headerFile); // connect file with ntuple
80   foutfileName=headerFile;
81   fNtuple   = 0;
82   fgen      = 0;
83
84   Init();
85   // add Task to //root/Tasks folder
86   TTask * roottasks = (TTask*)gROOT->GetRootFolder()->FindObject("Tasks") ; 
87   roottasks->Add(this) ; 
88 }
89 //____________________________________________________________________________ 
90 void AliTOFProb::Init()
91 {
92
93   fNtuple= (TNtuple*)fhfile->Get("Ntuple"); // get ntuple from file
94   Int_t nvar = fNtuple->GetNvar(); cout <<"N of var.="<< nvar << endl;
95   fNtuple->GetEvent(0);
96
97 }
98
99 //____________________________________________________________________________ 
100   AliTOFProb::~AliTOFProb()
101 {
102   //
103   // dtor (free used memory)
104   //
105
106   if (fhfile)
107     {
108       delete fhfile;
109       fhfile = 0;
110     }
111
112
113   if (fNtuple)
114     {
115       delete fNtuple;
116       fNtuple = 0;
117     }
118
119   if (fgen)
120     {
121       delete fgen;
122       fgen = 0;
123     }
124
125   if (foutfileName)
126     {
127       delete foutfileName;
128       foutfileName = 0;
129     }
130 }
131
132
133 //____________________________________________________________________________
134 void AliTOFProb::Exec(const Option_t * /*dummyOpt*/) 
135
136   //
137   // Performs Prob for TOF detector
138   // 
139   fTask=1;
140   fhfile->cd();
141   // create the histo
142   TH3* h3D=new TH3F("h3D", "", 120, -0.1, 1.1, 120, -0.1, 1.1, 120, -0.1, 1.1);
143   TH3* h3Ddt=new TH3F("h3Ddt", "", 360, -120., 120., 360, -120., 120., 360, -120., 120.);
144
145   TH1* hpion  =new TH1F("hpion", "", 120, -0.1, 1.1);
146   TH1* hkaon  =new TH1F("hkaon", "", 120, -0.1, 1.1);
147   TH1* hproton=new TH1F("hproton", "", 120, -0.1, 1.1);
148
149   TH1* hdtpion  =new TH1F("hdtpion", "", 600, -120., 120.);
150   TH1* hdtkaon  =new TH1F("hdtkaon", "", 600, -120., 120.);
151   TH1* hdtproton=new TH1F("hdtproton", "", 600, -120., 120.);
152
153   const Float_t timeResolution=0.12;
154
155   Int_t nparticles = (Int_t)fNtuple->GetEntries();
156   cout << " Number of nparticles =" << nparticles << endl;
157   if (nparticles <= 0) return;
158   
159   for (Int_t i=0; i < nparticles; i++) {
160     fNtuple->GetEvent(i);
161     //Int_t event=(Int_t)(fNtuple->GetLeaf("event")->GetValue());
162     Int_t pdgcode=(Int_t)(fNtuple->GetLeaf("ipart")->GetValue());
163     // unused for the time being in the context of probabilities
164     // probably it will be included in AliESD approach
165     //Float_t mass=fNtuple->GetLeaf("mext")->GetValue(0);
166     Int_t matc=(Int_t)(fNtuple->GetLeaf("matc")->GetValue(0));
167     Int_t imam=(Int_t)(fNtuple->GetLeaf("imam")->GetValue(0));
168     Float_t time=fNtuple->GetLeaf("text")->GetValue(0); // [ns]
169     // generate a gaussian around time with sigma timeresolution [ns]
170     Float_t lowTime=time-1000.*timeResolution;
171     Float_t upTime =time+1000.*timeResolution;
172     TF1* gauss=new TF1("tg","gaus",lowTime,upTime);
173
174     Double_t p1=time;
175     Double_t p2=timeResolution;
176     //Float_t twopi=2*TMath::Pi();
177     //Float_t sqrtwopi=TMath::Sqrt(twopi);
178     //Double_t p0=1./(p2*sqrtwopi);
179     // choising a different normalization
180     Double_t p0=1.;
181     gauss->SetParameters(p0,p1,p2);
182
183     Float_t px=fNtuple->GetLeaf("pxvtx")->GetValue(0);
184     Float_t py=fNtuple->GetLeaf("pyvtx")->GetValue(0);
185     Float_t pz=fNtuple->GetLeaf("pzvtx")->GetValue(0);
186     // track length [cm]
187     Float_t toflen=fNtuple->GetLeaf("leng")->GetValue(0);
188
189     Float_t pvtx=TMath::Sqrt(px*px+py*py+pz*pz);
190     Float_t ptvtx=TMath::Sqrt(px*px+py*py);
191
192     Bool_t isSelected=(imam == 0 && pz !=0 && TMath::ATan(TMath::Abs(ptvtx/pz))>TMath::Pi()*45./180.);
193     Bool_t trackWithTime=(matc==2 || matc==3 || matc==4);
194
195     Int_t abspdgcode=TMath::Abs(pdgcode);
196     switch(abspdgcode){
197     case 321:
198       break;
199     case 2212:
200       break;
201     case 11:
202       break;
203     default:
204       break;
205     }
206
207     if (isSelected && trackWithTime && toflen>370.)
208       {//only primary +/-45  
209         Float_t beta[3]={0.,0.,0.};
210         Float_t timeth[3]={0.,0.,0.};
211         Float_t massarray[3]={0.13957,0.493677,0.9382723};
212         for (Int_t j=0; j < 3; j++) {
213           Float_t dummy=(pvtx/massarray[j])*(pvtx/massarray[j]);
214           beta[j]=TMath::Sqrt(dummy/(1.+dummy));
215           timeth[j]=toflen/(29.9792458*beta[j]);
216         }
217         Float_t pPion   =(Float_t)gauss->Eval(timeth[0],0.,0.);
218         Float_t pKaon   =(Float_t)gauss->Eval(timeth[1],0.,0.);
219         Float_t pProton =(Float_t)gauss->Eval(timeth[2],0.,0.);
220
221         Float_t dtpion=time-timeth[0];
222         Float_t dtkaon=time-timeth[1];
223         Float_t dtproton=time-timeth[2];
224         hdtpion->Fill(dtpion);
225         hdtkaon->Fill(dtkaon);
226         hdtproton->Fill(dtproton);
227
228         // normalizing probabilities
229         Float_t sumProb=pPion+pKaon+pProton;
230         pPion=pPion/sumProb;
231         pKaon=pKaon/sumProb;
232         pProton=pProton/sumProb;
233         hpion->Fill(pPion);
234         hkaon->Fill(pKaon);
235         hproton->Fill(pProton);
236         //cout << "pPion " << pPion << " pKaon " << pKaon << " pProton " << pProton << endl;
237         h3D->Fill(pPion,pKaon,pProton);
238         h3Ddt->Fill(dtpion,dtkaon,dtproton);
239         //cout << "pion " << timeth[0] << " kaon " << timeth[1] << " proton " << timeth[2] << endl;
240       }// End of cuts appling
241
242     delete gauss;
243
244   }// End of loop over particles
245
246   char outFileName[100];
247   strcpy(outFileName,"3DProbHisto.root");
248   TFile *houtfile = new TFile(outFileName,"recreate");
249   houtfile->cd();
250   h3D->Write(0,TObject::kOverwrite);
251   h3Ddt->Write(0,TObject::kOverwrite);
252   hpion->Write(0,TObject::kOverwrite);
253   hkaon->Write(0,TObject::kOverwrite);
254   hproton->Write(0,TObject::kOverwrite);
255   hdtpion->Write(0,TObject::kOverwrite);
256   hdtkaon->Write(0,TObject::kOverwrite);
257   hdtproton->Write(0,TObject::kOverwrite);
258
259   houtfile->Close();
260   houtfile->Write(0,TObject::kOverwrite);
261 }
262
263
264 //__________________________________________________________________
265 Bool_t AliTOFProb::operator==( AliTOFProb const & /*tofrec*/)const
266 {
267   // dummy version of Equal operator.
268   // requested by coding conventions
269   return kTRUE;
270
271 }