-#include <AliITSVertexerIons.h>
-#include "stdlib.h"
-#include <TMath.h>
-#include <TRandom.h>
-#include <TObjArray.h>
-#include <TROOT.h>
-#include <TFile.h>
-#include <TTree.h>
-#include <Riostream.h>
-
-#include "AliRun.h"
-#include "AliITS.h"
-#include "AliITSgeom.h"
+#include "AliRunLoader.h"
+#include "AliITSDetTypeRec.h"
+#include "AliITSVertexerIons.h"
+#include "AliITSVertexerZ.h"
+#include "AliESDVertex.h"
+#include "AliITSgeomTGeo.h"
#include "AliITSLoader.h"
#include "AliITSRecPoint.h"
-#include "AliGenerator.h"
-#include "AliMagF.h"
-
+#include "AliLog.h"
+#include <TTree.h>
#include <TH1.h>
#include <TF1.h>
-#include <TCanvas.h>
-#include <AliESDVertex.h>
-#include <TObjArray.h>
-#include <TObject.h>
-#include <AliITSVertexerPPZ.h>
+
//////////////////////////////////////////////////////////////////////
// AliITSVertexerIons is a class for full 3D primary vertex //
// finding optimized for Ion-Ion interactions //
// //
//////////////////////////////////////////////////////////////////////
-ClassImp(AliITSVertexerIons)
-
-
+ClassImp(AliITSVertexerIons)
+
//______________________________________________________________________
- AliITSVertexerIons::AliITSVertexerIons():AliITSVertexer() {
+ AliITSVertexerIons::AliITSVertexerIons():AliITSVertexer(),
+fNpThreshold(0),
+fMaxDeltaPhi(0),
+fMaxDeltaZ(0){
// Default Constructor
- fITS = 0;
+ //fITS = 0;
SetNpThreshold();
SetMaxDeltaPhi();
SetMaxDeltaZ();
}
-//______________________________________________________________________
-AliITSVertexerIons::AliITSVertexerIons(TString fn):AliITSVertexer(fn) {
- // Standard constructor
-
- fITS = 0;
- SetNpThreshold();
- SetMaxDeltaPhi();
- SetMaxDeltaZ();}
-
-
//______________________________________________________________________
AliITSVertexerIons::~AliITSVertexerIons() {
// Default Destructor
- fITS = 0;
+ //fITS = 0;
}
//______________________________________________________________________
-AliESDVertex* AliITSVertexerIons::FindVertexForCurrentEvent(Int_t evnumber){
+AliESDVertex* AliITSVertexerIons::FindVertexForCurrentEvent(TTree *itsClusterTree){
+// Defines the AliESDVertex for the current event
fCurrentVertex = 0;
- AliRunLoader *rl = AliRunLoader::GetRunLoader();
- AliITSLoader* itsloader = (AliITSLoader*) rl->GetLoader("ITSLoader");
- if(!fITS) {
- fITS =(AliITS *)gAlice->GetDetector("ITS");
- if(!fITS) {
- Error("FindVertexForCurrentEvent","AliITS object was not found");
- return fCurrentVertex;
- }
- }
- fITS->SetTreeAddress();
- AliITSgeom *g2 = fITS->GetITSgeom();
- TClonesArray *recpoints = fITS->RecPoints();
+ AliITSDetTypeRec detTypeRec;
+
+ detTypeRec.SetTreeAddressR(itsClusterTree);
+
+ TClonesArray *recpoints = detTypeRec.RecPoints();
AliITSRecPoint *pnt;
- TTree *tr = itsloader->TreeR();
+
Int_t npoints=0;
Int_t nopoints1=40000;
Int_t nopoints2=40000;
- Float_t l[3], p[3];
+ Float_t p[3];
Double_t *z1, *z2, *y1, *y2, *x1, *x2, *phi1, *phi2, *r1, *r2;
z1=new Double_t[nopoints1];
Double_t r=0;
Int_t np1=0, np2=0;
- for(Int_t i=g2->GetStartSPD();i<=g2->GetLastSPD();i++) {
- fITS->ResetRecPoints();
- tr->GetEvent(i);
+ for(Int_t i=AliITSgeomTGeo::GetModuleIndex(1,1,1);i<=AliITSgeomTGeo::GetModuleIndex(2,1,1)-1;i++) {
+ detTypeRec.ResetRecPoints();
+ itsClusterTree->GetEvent(i);
npoints = recpoints->GetEntries();
for (Int_t ipoint=0;ipoint<npoints;ipoint++) {
pnt = (AliITSRecPoint*)recpoints->UncheckedAt(ipoint);
- l[0]=pnt->GetX();
+ /*
+ l[0]=pnt->GetDetLocalX();
l[1]=0;
- l[2]=pnt->GetZ();
+ l[2]=pnt->GetDetLocalZ();
g2->LtoG(i, l, p);
+ */
+ pnt->GetGlobalXYZ(p);
r=TMath::Sqrt(TMath::Power(p[0],2)+TMath::Power(p[1],2));
if(i<80 && TMath::Abs(p[2])<14.35) {
}
if(np1<fNpThreshold) {
- Warning("FindVertexForCurrentEvent","AliITSVertexerIons finder is not reliable for low multiplicity events. Switching to AliITSVertexerPPZ with default parameters...\n");
+ Warning("FindVertexForCurrentEvent","AliITSVertexerIons finder is not reliable for low multiplicity events. Switching to AliITSVertexerZ with default parameters...\n");
Warning("FindVertexForCurrentEvent","N rec points = %d - Threshold is %d",np1,fNpThreshold);
- AliITSVertexerPPZ *dovert = new AliITSVertexerPPZ("default");
- fCurrentVertex =dovert->FindVertexForCurrentEvent(rl->GetEventNumber());
+ AliITSVertexerZ *dovert = new AliITSVertexerZ();
+ fCurrentVertex =dovert->FindVertexForCurrentEvent(itsClusterTree);
delete dovert;
return fCurrentVertex;
}
Error("FindVertexForCurrentEvent","No points in the pixer layers");
return fCurrentVertex;
}
- if(fDebug>0) cout << "N. points layer 1 and 2 : " << np1 << " " << np2 << endl;
+ AliDebug(1,Form("N. points layer 1 and 2 : %d %d",np1,np2));
Double_t asparx = (mxpiu-mxmeno)/(mxpiu+mxmeno);
Double_t aspary = (mypiu-mymeno)/(mypiu+mymeno);
Double_t x0 = 2.86*asparx;
Double_t y0 = 2.86*aspary;
- if(fDebug>0) cout << "Rough xy vertex = " << x0 << " " << y0 << endl;
+ AliDebug(1,Form("Rough xy vertex = %f %f",x0,y0));
for(Int_t i=0;i<np1;i++) {
x1[i]-=x0;
// Fitting ...
Double_t max;
- Int_t bin_max;
- Double_t max_center;
+ Int_t binMax;
+ Double_t maxcenter;
max = hzv->GetMaximum();
- bin_max=hzv->GetMaximumBin();
- max_center=hzv->GetBinCenter(bin_max);
+ binMax=hzv->GetMaximumBin();
+ maxcenter=hzv->GetBinCenter(binMax);
Double_t dxy=1.5;
- TF1 *fz = new TF1 ("fz","([0]*exp(-0.5*((x-[1])/[2])*((x-[1])/[2])))+[3]",max_center-dxy,max_center+dxy);
+ TF1 *fz = new TF1 ("fz","([0]*exp(-0.5*((x-[1])/[2])*((x-[1])/[2])))+[3]",maxcenter-dxy,maxcenter+dxy);
fz->SetParameter(0,max);
- fz->SetParameter(1,max_center);
+ fz->SetParameter(1,maxcenter);
fz->SetParameter(2,0.1);
fz->SetLineColor(kRed);
hzv->Fit("fz","RQ0");
Double_t resolution[3]={0,0,0};
Double_t snr[3]={0,0,0};
Char_t name[30];
- if(fDebug>0)Info("FindVertexForCurrentEvent","Vertex found for event %d",evnumber);
sprintf(name,"Vertex");
fCurrentVertex = new AliESDVertex(position,resolution,snr,name);
return fCurrentVertex;
TF1 *fx = new TF1 ("fx","([0]*exp(-0.5*((x-[1])/[2])*((x-[1])/[2])))+[3]+[4]*x+[5]*x*x",x0-dxy,x0+dxy);
fx->SetParameter(0,100);
Double_t dist=0.3;
- Double_t x_approx=FindMaxAround(x0,hxv,dist);
- if(fDebug>0) cout << "x_approx = " << x_approx << endl;
- fx->SetParameter(1,x_approx);
- Double_t dif_centroid=0.07;
- fx->SetParLimits(1,x_approx-dif_centroid,x_approx+dif_centroid);
+ Double_t xapprox=FindMaxAround(x0,hxv,dist);
+ AliDebug(1,Form("xapprox = %f",xapprox));
+ fx->SetParameter(1,xapprox);
+ Double_t difcentroid=0.07;
+ fx->SetParLimits(1,xapprox-difcentroid,xapprox+difcentroid);
fx->SetParameter(2,0.1);
fx->SetLineColor(kRed);
hxv->Fit("fx","RQW0");
TF1 *fy = new TF1 ("fy","([0]*exp(-0.5*((x-[1])/[2])*((x-[1])/[2])))+[3]+[4]*x+[5]*x*x",y0-dxy,y0+dxy);
fy->SetParameter(0,100);
- Double_t y_approx=FindMaxAround(y0,hyv,dist);
- if(fDebug>0) cout << "y_approx = " << y_approx << endl;
- fy->SetParameter(1,y_approx);
- fy->SetParLimits(1,y_approx-dif_centroid,y_approx+dif_centroid);
+ Double_t yapprox=FindMaxAround(y0,hyv,dist);
+ AliDebug(1,Form("yapprox = %f",yapprox));
+ fy->SetParameter(1,yapprox);
+ fy->SetParLimits(1,yapprox-difcentroid,yapprox+difcentroid);
fy->SetParameter(2,0.1);
fy->SetLineColor(kRed);
hyv->Fit("fy","RQW0");
Double_t snr[3]={0,0,0};
Char_t name[30];
- if(fDebug>0)Info("FindVertexForCurrentEvent","Vertex found for event %d",evnumber);
sprintf(name,"Vertex");
fCurrentVertex = new AliESDVertex(position,resolution,snr,name);
+
return fCurrentVertex;
}
//______________________________________________________________________
void AliITSVertexerIons::PhiFunc(Double_t &x,Double_t &y,Double_t &phi) {
+ // Method for the computation of the phi angle given x and y. The result is in degrees.
if(y>0 && x>0) phi=(TMath::ATan((Double_t)(y/x))*57.29578);
if(y>0 && x<0) phi=(TMath::ATan((Double_t)(y/x))*57.29578)+180;
if(y<0 && x<0) phi=(TMath::ATan((Double_t)(y/x))*57.29578)+180;
if(y<0 && x>0) phi=(TMath::ATan((Double_t)(y/x))*57.29578)+360;;
}
-//______________________________________________________________________
-void AliITSVertexerIons::FindVertices(){
- // computes the vertices of the events in the range FirstEvent - LastEvent
- AliRunLoader *rl = AliRunLoader::GetRunLoader();
- AliITSLoader* itsloader = (AliITSLoader*) rl->GetLoader("ITSLoader");
- itsloader->LoadRecPoints("read");
- for(Int_t i=fFirstEvent;i<=fLastEvent;i++){
- rl->GetEvent(i);
- FindVertexForCurrentEvent(i);
- if(fCurrentVertex){
- WriteCurrentVertex();
- }
- else {
- if(fDebug>0){
- cout<<"Vertex not found for event "<<i<<endl;
- }
- }
- }
-}
-
//________________________________________________________
void AliITSVertexerIons::PrintStatus() const {
// Print current status
cout <<"=======================================================\n";
- cout <<" Debug flag: "<<fDebug<<endl;
- cout<<"First event to be processed "<<fFirstEvent;
- cout<<"\n Last event to be processed "<<fLastEvent<<endl;
if(fCurrentVertex)fCurrentVertex->PrintStatus();
}
//________________________________________________________
Double_t AliITSVertexerIons::FindMaxAround(Double_t point, TH1F *h, Double_t distance) {
- Int_t max_content=0;
- Int_t max_bin=0;
+ // It finds the maximum of h within point-distance and distance+point.
+ Int_t maxcontent=0;
+ Int_t maxbin=0;
for(Int_t i=0;i<h->GetNbinsX();i++) {
Int_t content=(Int_t)h->GetBinContent(i);
Double_t center=(Double_t)h->GetBinCenter(i);
- if(fabs(center-point)>distance) continue;
- if(content>max_content) {max_content=content;max_bin=i;}
+ if(TMath::Abs(center-point)>distance) continue;
+ if(content>maxcontent) {maxcontent=content;maxbin=i;}
}
- Double_t max=h->GetBinCenter(max_bin);
+ Double_t max=h->GetBinCenter(maxbin);
return max;
}