-#include "iostream.h"
-#include "TMath.h"
+#if !defined(__CINT__) || defined(__MAKECINT__)
+#include <TClassTable.h>
+#include <TFile.h>
+#include <TTree.h>
+#include <Riostream.h>
+#include <AliRun.h>
+#include <AliHeader.h>
+#include <AliGenEventHeader.h>
+#include <AliITSVertexerIons.h>
+#include <AliRunLoader.h>
+#include <AliITSLoader.h>
-void AliITSFindPrimaryVertex(Int_t evNumber1=0,Int_t evNumber2=0, const char *filename="galice.root") {
+#endif
-
- ///////////////// Dynamically link some shared libs ////////////////////////////////
+void AliITSFindPrimaryVertex(Int_t evNumber1=0,Int_t NumbofEv=1, const char *filename="galice.root") {
+
+ Int_t evNumber2 = evNumber1+NumbofEv;
if (gClassTable->GetID("AliRun") < 0) {
- gROOT->LoadMacro("loadlibs.C");
- loadlibs();
+ gROOT->Macro("loadlibs.C");
} else {
- delete gAlice;
- gAlice=0;
+ if(gAlice){
+ delete gAlice->GetRunLoader();
+ delete gAlice;
+ gAlice=0;
+ }
}
-// Connect the Root Galice file containing Geometry, Kine and Hits
- TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(filename);
- if (!file) file = new TFile(filename,"UPDATE");
-
-
-// Get AliRun object from file or create it if not on file
- if (!gAlice) {
- gAlice = (AliRun*)file->Get("gAlice");
- if (gAlice) printf("AliRun object found on file\n");
- if (!gAlice) gAlice = new AliRun("gAlice","Alice test program");
- }
-
-// Loop over events
-//
- Int_t Nh=0;
- Int_t Nh1=0;
- for (int nev=0; nev<= evNumber2; nev++) {
- gAlice->SetEvent(nev);
- Int_t nparticles = gAlice->GetEvent(nev);
- cout << "nev " << nev <<endl;
- cout << "nparticles " << nparticles <<endl;
- if (nev < evNumber1) continue;
- if (nparticles <= 0) return;
-
- TStopwatch timer;
- timer.Start();
-
- AliITSVertex *V = new AliITSVertex();
- V->Exec();
+ AliRunLoader* rl = AliRunLoader::Open("galice.root");
+ if (rl == 0x0){
+ ::Error("AliITSFindPrimaryVertex.C","Can not open session RL=NULL");
+ return;
+ }
+
+ Int_t retval = rl->LoadHeader();
+ if (retval){
+ cerr<<"AliITSFindPrimaryVertex.C : LoadHeader returned error"<<endl;
+ return;
+ }
+ retval = rl->LoadKinematics();
+ if (retval){
+ cerr<<"AliITSFindPrimaryVertex.C : LoadKinematics returned error"<<endl;
+ return;
+ }
- timer.Stop();
- timer.Print();
+ // Open output file for vertices (default name: ITS.Vertex.root
+ // and Create vertexer
+ AliITSVertexerIons *vertexer = new AliITSVertexerIons("default");
+ AliITSVertex *V;
+ // Loop over events
+ //
+
+ for (int nev=evNumber1; nev< evNumber2; nev++) {
+ cout<<"=============================================================\n";
+ cout<<" Processing event "<<nev<<endl;
+ cout<<"=============================================================\n";
+ rl->GetEvent(nev);
+ cout << "nev " << nev <<endl;
+ // The true Z coord. is fetched for comparison
+ AliHeader *header = rl->GetHeader();
+ AliGenEventHeader* genEventHeader = header->GenEventHeader();
+ TArrayF primaryVertex(3);
+ genEventHeader->PrimaryVertex(primaryVertex);
+
+ TStopwatch timer;
+ timer.Start();
- cout << endl << "Xv = " << V->GetXv() << " cm" << endl;
- cout << "X resolution = " << V->GetXRes()*10000 << " microns" << endl;
- cout << "Signal/Noise for X = " << V->GetXSNR() << endl;
- cout << endl << "Yv = " << V->GetYv() << " cm" << endl;
- cout << "Y resolution = " << V->GetYRes()*10000 << " microns" << endl;
- cout << "Signal/Noise for Y = " << V->GetYSNR() << endl;
- cout << endl << "Zv = " << V->GetZv() << " cm" << endl;
- cout << "Z Resolution = " << V->GetZRes()*10000 << " microns" << endl;
- cout << "Signal/Noise for Z = " << V->GetZSNR() <<endl;
+ V=vertexer->FindVertexForCurrentEvent(nev);
+ if(V){
+ Double_t pos[3];
+ for(Int_t kk=0;kk<3;kk++)pos[kk]=(Double_t)primaryVertex[kk];
+ V->SetTruePos(pos);
+ }
+ timer.Stop();
+ timer.Print();
+ if(!V)continue;
+ cout << endl << "Xv = " << V->GetXv() << " cm" << endl;
+ cout << "X resolution = " << V->GetXRes()*10000 << " microns" << endl;
+ cout << "Signal/Noise for X = " << V->GetXSNR() << endl;
+ cout << endl << "Yv = " << V->GetYv() << " cm" << endl;
+ cout << "Y resolution = " << V->GetYRes()*10000 << " microns" << endl;
+ cout << "Signal/Noise for Y = " << V->GetYSNR() << endl;
+ cout << endl << "Zv = " << V->GetZv() << " cm" << endl;
+ cout << "Z Resolution = " << V->GetZRes()*10000 << " microns" << endl;
+ cout << "Signal/Noise for Z = " << V->GetZSNR() <<endl;
- delete V;
+ vertexer->WriteCurrentVertex();
}
- file->Close();
+
}