VertexerZ updated: no histograms for Z finding, new definition of errors
[u/mrichter/AliRoot.git] / ITS / AliITSVertexerZTest.C
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include <TFile.h>
3 #include <TGeoManager.h>
4 #include <TH1F.h>
5 #include <TH2F.h>
6 #include <TTree.h>
7 #include <Riostream.h>
8 #include <AliRun.h>
9 #include <AliHeader.h>
10 #include <AliGenEventHeader.h>
11 #include <AliITSVertexerZ.h>
12 #include <AliRunLoader.h>
13 #include <AliITSLoader.h>
14
15 #endif
16
17 void AliITSVertexerZTest(Float_t delphi=0.05,Float_t window=3.,Float_t initx=0., Float_t inity=0.,TString FileIn="galice.root"){
18   // delphi ---> azimuthal range to accept tracklets
19   // window ---> window in Z around the peak of tracklets proj. in mm
20   Int_t kDebug = 50;
21   TH2F *diff2 = new TH2F("diff2","Zfound vs Ztrue",100,-6,6,100,-6,6); 
22   TH1F *diff1 = new TH1F("diff1","Zfound - Ztrue(#mu m)",100,-500,500);
23   delete gAlice;
24   gAlice = 0;
25   AliRunLoader *rl = AliRunLoader::Open(FileIn.Data());
26   Int_t retval = rl->LoadgAlice();
27   if(retval){
28     cerr<<"AliITSVertexerZTest.C: AliRun object not found"<<endl;
29     return;
30   }
31   retval = rl->LoadHeader();
32   if (retval){
33     cerr<<"AliITSVertexerZTest.C : LoadHeader returned error"<<endl;
34     return;
35   }
36   retval = rl->LoadKinematics();
37   if (retval){
38     cerr<<"AliITSVertexerZTest.C : LoadKinematics returned error"<<endl;
39     return;
40   }
41   TGeoManager::Import("geometry.root");
42   AliITSLoader* ITSloader =  (AliITSLoader*) rl->GetLoader("ITSLoader");
43
44   if(!ITSloader){
45     cerr<<"AliITSVertexerZTest.C :  ITS loader not found"<<endl;
46     return;
47   }
48   //  ITSloader->LoadRecPoints("read");
49   //  TFile *fo = new TFile("vertici.root","recreate");
50   //  AliITSVertexerPPZ *dovert = new AliITSVertexerPPZ("default",initx,inity);
51   AliITSVertexerZ *dovert = new AliITSVertexerZ("default",initx,inity);
52   dovert->SetDebug(0);
53   //  dovert->SetDiffPhiMax(delphi);
54   //  dovert->SetWindow(window);
55   dovert->PrintStatus();
56   Int_t sigmazero=0;
57   AliESDVertex *vert = 0;
58   for(Int_t i=0; i<rl->TreeE()->GetEntries(); i++){
59     rl->GetEvent(i);
60     // The true Z coord. is fetched for comparison
61     AliHeader *header = rl->GetHeader();
62     AliGenEventHeader* genEventHeader = header->GenEventHeader();
63     TArrayF primaryVertex(3);
64     genEventHeader->PrimaryVertex(primaryVertex);
65     vert = dovert->FindVertexForCurrentEvent(i);
66     if(kDebug>0){
67       // Prints the results
68       cout <<"========================================================\n";
69       cout << "Event number: "<<i<<")  Z Vertex:"<<endl;
70       if(vert){
71         cout<<"FOUND: "<<vert->GetZv()<<"; ";
72         cout<<vert->GetZRes()<<"; "<<vert->GetNContributors()<<endl;
73         cout <<" True Z position "<<primaryVertex[2]<<", diff= ";
74         cout<<(primaryVertex[2]-vert->GetZv())*10000.<<endl;
75       } else {
76         cout<<"NOT FOUND"<<endl;
77       }
78     }
79     if(vert){
80       Double_t pos[3];
81       for(Int_t kk=0;kk<3;kk++)pos[kk]=(Double_t)primaryVertex[kk];
82       vert->SetTruePos(pos);
83       Float_t found = vert->GetZv();
84       diff2->Fill(primaryVertex[2],found);
85       found = 10000.*(found-primaryVertex[2]);
86       if(vert->GetZRes()!=0){
87         diff1->Fill(found);
88       } else {
89         sigmazero++;
90       }
91       dovert->WriteCurrentVertex();
92     }
93   }
94   if(kDebug>0){
95     cout<<"Only one tracklet (sigma = 0) "<<sigmazero<<endl;
96   }
97   
98 }