Fixes for building of DA (Anshul)
[u/mrichter/AliRoot.git] / TPC / AliTPCSavePID.C
1 ////////////////////////////////////////////////
2 // Dubna, 22 Jan 2003                         //
3 ////////////////////////////////////////////////
4 Int_t AliTPCSavePID(Int_t emax=3) {
5 ///////////////////////////////////////
6 #include<fstream.h>
7 TFile *inkin = TFile::Open("galice.root");
8 if (!inkin->IsOpen()) {
9   cerr<<"Can't open galice.root !\n";
10 }                                                                                        
11 gAlice = (AliRun*)inkin->Get("gAlice");
12 cout<<"AliRun object found on file "<<gAlice<<endl;
13 cout<<"!!!! field ="<<gAlice->Field()->SolenoidField()<<endl;
14 AliKalmanTrack::SetFieldMap(gAlice->Field());
15 inkin->Close();                                  
16  cout<<" Convconst="<<AliKalmanTrack::GetConvConst()<<endl;
17 /////////////////////////////////////// 
18 AliTPCtrack *iotrack=0; 
19 TFile *tf=TFile::Open("AliTPCtracks.root");
20    if (!tf->IsOpen()) {cerr<<"Can't open AliTPCtracks.root !\n"; return 3;}
21                                         
22 AliTPCPid pid;
23 AliTPCtrackPid pidtmp;                                                                
24 AliTPCtrackPid *outpid=&pidtmp;
25 TFile fpid("AliTPCtracksPid.root","recreate");
26
27  Int_t nentr=0,i=0;
28  for(int je=0;je<emax;je++){
29    char tname[100]; sprintf(tname,"TreeT_TPC_%d",je);
30    TTree *tracktree=(TTree*)tf->Get(tname);
31    if (!tracktree) {cerr<<"Can't get a tree with TPC tracks !\n"; return 4;}
32
33    TBranch *tbranch=tracktree->GetBranch("tracks");
34    nentr=(Int_t)tracktree->GetEntries();
35    cout<<nentr<<" "<<" tracks in track tree "<<tname<<"."<<endl;
36
37     char tpname[100]; sprintf(tpname,"TreeT%d",je);
38     TTree *ptree = new  TTree(tpname,"Tree with PID"); 
39     ptree->Branch("pids","AliTPCtrackPid",&outpid,32000,1); 
40
41    for (i=0; i<nentr; i++) {
42        iotrack=new AliTPCtrack;
43        tbranch->SetAddress(&iotrack);
44        tracktree->GetEvent(i);
45       Double_t par[5],xx;xx=85.1951;
46       iotrack->PropagateTo(xx); 
47       iotrack->GetExternalParameters(xx,par);
48       //cout<<" par="<<par[0]<<" "<<par[1]<<" "<<par[2]<<" "<<par[3]<<" "<<par[4]<<endl;
49       Float_t phi=TMath::ASin(par[2]) + iotrack->GetAlpha();
50       if (phi<-TMath::Pi()) phi+=2*TMath::Pi();
51       if (phi>=TMath::Pi()) phi-=2*TMath::Pi();
52       Float_t lam=TMath::ATan(par[3]); 
53       Float_t pt_1=TMath::Abs(par[4]);
54       //cout<<"pt_1,lam="<<pt_1<<" "<<lam<<endl;
55     if(pt_1!=0.){
56       Float_t mom=1./(pt_1*TMath::Cos(lam));
57       Float_t dedx=iotrack->GetdEdx();
58       Int_t pcode=pid.GetPcode(dedx/50.,mom);
59        pidtmp.fPcode=pcode;
60        pidtmp.fSignal=dedx;
61        pidtmp.fMom=mom;
62        pidtmp.fPhi=phi;
63        pidtmp.fLam=lam;
64        pidtmp.fWpi=pid.fWpi;
65        pidtmp.fWk=pid.fWk;
66        pidtmp.fWp=pid.fWp;
67        pidtmp.fLabel=TMath::Abs(iotrack->GetLabel());
68        cout<<"tlab,dedx,mom,pcode="<<pidtmp.fLabel<<" "<<dedx<<" "<<mom<<" "<<pcode<<endl;      
69        ptree->Fill();
70      }//if(pt...
71    }//for(i...
72   cout<<"Event "<<je+1<<" ok."<<endl;      
73    delete tracktree; 
74  }//End for(je...)   
75 tf->Close();
76 fpid.Write();
77 cout<<"File AliTPCtracksPid.root written"<<endl;
78 return 0;
79 }