#include "AliITStrackerSA.h"
#include "AliITSVertexerIons.h"
#include "AliITSVertexerFast.h"
-#include "AliITSVertexerPPZ.h"
+#include "AliITSVertexer3D.h"
#include "AliITSVertexerZ.h"
+#include "AliITSVertexerCosmics.h"
#include "AliESD.h"
#include "AliITSpidESD.h"
#include "AliITSpidESD1.h"
Info("CreateVertexer","a AliITSVertexerFast object has been selected\n");
return new AliITSVertexerFast(smear);
}
- if(selectedVertexer.Contains("ppz") || selectedVertexer.Contains("PPZ")){
- Info("CreateVertexer","a AliITSVertexerPPZ object has been selected\n");
- return new AliITSVertexerPPZ("null");
+ if(selectedVertexer.Contains("3d") || selectedVertexer.Contains("3D")){
+ Info("CreateVertexer","a AliITSVertexer3D object has been selected\n");
+ return new AliITSVertexer3D("null");
+ }
+ if(selectedVertexer.Contains("cosmics") || selectedVertexer.Contains("COSMICS")){
+ Info("CreateVertexer","a AliITSVertexerCosmics object has been selected\n");
+ return new AliITSVertexerCosmics();
}
// by default an AliITSVertexerZ object is instatiated
Info("CreateVertexer","a AliITSVertexerZ object has been selected\n");
//////////////////////////////////////////////////////////////////////
// Base class for primary vertex reconstruction //
// AliESDVertexer is a class for full 3D primary vertex finding //
-// derived classes: AliITSVertexerIons AliITSvertexerPPZ //
-// AliITSVertexer3D //
+// derived classes: AliITSVertexerIons AliITSvertexer3D //
+// AliITSVertexerCosmics //
//////////////////////////////////////////////////////////////////////
//______________________________________________________________________
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+
+#include <TClonesArray.h>
+#include "AliLog.h"
+#include "AliESDVertex.h"
+#include "AliRunLoader.h"
+#include "AliITSLoader.h"
+#include "AliITSgeom.h"
+#include "AliITSRecPoint.h"
+#include "AliITSVertexerCosmics.h"
+#include "AliStrLine.h"
+
+//------------------------------------------------------------------------
+// This class implements a method to construct a "fake" primary
+// vertex for cosmic events in which the muon crosses the SPD inner
+// layer (SPD1). A fake primary vertex is needed for the reconstruction,
+// with e.g. AliITStrackerSA, of the two tracks produced by the muon
+// in the ITS.
+// We build pairs of clusters on SPD1 and define the fake vertex as
+// the mid-point of the straight line joining the two clusters.
+// We reject the backgroung by requiring at least one clusters on SPD2
+// closer than fMaxDistOnSPD2 to the tracklet prolongation.
+// We can reject (potentially pathological) events with the muon track
+// tangential to the SPD1 layer by the requiring the radial position of
+// the vertex to be smaller than fMaxVtxRadius.
+// Due to background clusters, more than one vertex per event can
+// be found. We consider the first found.
+// The number of contributors set in the AliESDVertex object is the
+// number of vertices found in the event; if this number is <1,
+// the procedure could not find a vertex position and by default
+// the vertex coordinates are set to (0,0,0) with large errors (100,100,100)
+// Number of contributors = 0 --> No SPD1 tracklets matching criteria
+// Number of contributors = -1 --> No SPD1 recpoints
+//
+// Origin: A.Dainese, andrea.dainese@lnl.infn.it
+//-------------------------------------------------------------------------
+
+ClassImp(AliITSVertexerCosmics)
+
+//-------------------------------------------------------------------------
+AliITSVertexerCosmics::AliITSVertexerCosmics():AliITSVertexer(),
+fFirstSPD1(0),
+fLastSPD1(0),
+fFirstSPD2(0),
+fLastSPD2(0),
+fMaxDistOnSPD2(0),
+fMaxVtxRadius(0),
+fMinDist2Vtxs(0)
+{
+ // Default constructor
+ SetSPD1Modules();
+ SetSPD2Modules();
+ SetMaxDistOnSPD2();
+ SetMaxVtxRadius();
+ SetMinDist2Vtxs();
+}
+//--------------------------------------------------------------------------
+AliESDVertex* AliITSVertexerCosmics::FindVertexForCurrentEvent(Int_t evnumber)
+{
+ // Defines the AliESDVertex for the current event
+
+ fCurrentVertex = 0;
+ AliRunLoader *rl =AliRunLoader::GetRunLoader();
+ AliITSLoader* itsLoader = (AliITSLoader*)rl->GetLoader("ITSLoader");
+ AliITSgeom* geom = itsLoader->GetITSgeom();
+ itsLoader->LoadRecPoints();
+ rl->GetEvent(evnumber);
+
+ TTree *rpTree = itsLoader->TreeR();
+
+ TClonesArray *recpoints=new TClonesArray("AliITSRecPoint",10000);
+ rpTree->SetBranchAddress("ITSRecPoints",&recpoints);
+
+ Double_t xclspd1[100],yclspd1[100],zclspd1[100],modclspd1[100];
+ Int_t nclspd1stored=0;
+ Double_t xclspd2[100],yclspd2[100],zclspd2[100],modclspd2[100];
+ Int_t nclspd2stored=0;
+ Int_t nrecpoints,nrecpointsSPD1=0;
+
+ Double_t gc[3]={0.,0.,0.};
+ Double_t lc[3]={0.,0.,0.};
+ Int_t lay,lad,det;
+
+ Double_t x[100],y[100],z[100],p1[3],p2[3],p3[3];
+ Int_t nvtxs;
+ Bool_t good,matchtospd2;
+ Double_t xvtx,yvtx,zvtx,rvtx;
+
+ for(Int_t imodule=fFirstSPD1; imodule<fLastSPD2; imodule++) { // SPD
+ rpTree->GetEvent(imodule);
+ geom->GetModuleId(imodule,lay,lad,det);
+ nrecpoints=recpoints->GetEntriesFast();
+ if(imodule<fLastSPD1) nrecpointsSPD1 += nrecpoints;
+ //printf("cosmics: module %d clusters %d\n",imodule,nrecpoints);
+ for(Int_t irp=0; irp<nrecpoints; irp++) {
+ AliITSRecPoint *rp=(AliITSRecPoint*)recpoints->UncheckedAt(irp);
+ // Local coordinates of this recpoint
+ lc[0]=rp->GetDetLocalX();
+ lc[2]=rp->GetDetLocalZ();
+ geom->LtoG(imodule,lc,gc); // global coordinates
+ if(lay==1) { // store SPD1 clusters
+ xclspd1[nclspd1stored]=gc[0];
+ yclspd1[nclspd1stored]=gc[1];
+ zclspd1[nclspd1stored]=gc[2];
+ modclspd1[nclspd1stored]=imodule;
+ nclspd1stored++;
+ }
+ if(lay==2) { // store SPD2 clusters
+ xclspd2[nclspd2stored]=gc[0];
+ yclspd2[nclspd2stored]=gc[1];
+ zclspd2[nclspd2stored]=gc[2];
+ modclspd2[nclspd2stored]=imodule;
+ nclspd2stored++;
+ }
+ if(nclspd1stored>100 || nclspd2stored>100)
+ AliFatal("More than 100 clusters per layer in SPD");
+ }// end clusters in a module
+ }// end SPD modules for a given event
+
+ // build fake vertices
+ nvtxs=0;
+ // SPD1 - first cluster
+ for(Int_t i1spd1=0; i1spd1<nclspd1stored; i1spd1++) {
+ p1[0]=xclspd1[i1spd1]; p1[1]=yclspd1[i1spd1]; p1[2]=zclspd1[i1spd1];
+ // SPD1 - second cluster
+ for(Int_t i2spd1=i1spd1+1; i2spd1<nclspd1stored; i2spd1++) {
+ if(modclspd1[i1spd1]==modclspd1[i2spd1]) continue;
+ p2[0]=xclspd1[i2spd1]; p2[1]=yclspd1[i2spd1]; p2[2]=zclspd1[i2spd1];
+ // look for point on SPD2
+ AliStrLine spd1line(p1,p2,kTRUE);
+ matchtospd2 = kFALSE;
+ for(Int_t ispd2=0; ispd2<nclspd2stored; ispd2++) {
+ p3[0]=xclspd2[ispd2]; p3[1]=yclspd2[ispd2]; p3[2]=zclspd2[ispd2];
+ //printf(" %f\n",spd1line.GetDistFromPoint(p3));
+ if(spd1line.GetDistFromPoint(p3)<fMaxDistOnSPD2)
+ { matchtospd2 = kTRUE; break; }
+ }
+ if(!matchtospd2) continue;
+ xvtx = 0.5*(xclspd1[i1spd1]+xclspd1[i2spd1]);
+ yvtx = 0.5*(yclspd1[i1spd1]+yclspd1[i2spd1]);
+ zvtx = 0.5*(zclspd1[i1spd1]+zclspd1[i2spd1]);
+ rvtx = TMath::Sqrt(xvtx*xvtx+yvtx*yvtx);
+ if(rvtx>fMaxVtxRadius) continue;
+ good = kTRUE;
+ for(Int_t iv=0; iv<nvtxs; iv++) {
+ if(TMath::Sqrt((xvtx- x[iv])*(xvtx- x[iv])+
+ (yvtx- y[iv])*(yvtx- y[iv])+
+ (zvtx- z[iv])*(zvtx- z[iv])) < fMinDist2Vtxs)
+ good = kFALSE;
+ }
+ if(good) {
+ x[nvtxs]=xvtx;
+ y[nvtxs]=yvtx;
+ z[nvtxs]=zvtx;
+ nvtxs++;
+ }
+ } // SPD1 - second cluster
+ } // SPD1 - first cluster
+
+
+ Double_t pos[3]={0.,0.,0.};
+ Double_t err[3]={100.,100.,100.};
+ if(nvtxs) {
+ pos[0]=x[0];
+ pos[1]=y[0];
+ pos[2]=z[0];
+ err[0]=0.1;
+ err[1]=0.1;
+ err[2]=0.1;
+ }
+ if(!nrecpointsSPD1) nvtxs=-1;
+ fCurrentVertex = new AliESDVertex(pos,err,"cosmics");
+ fCurrentVertex->SetNContributors(nvtxs);
+ fCurrentVertex->SetTitle("cosmics fake vertex");
+
+ //if(nvtxs>0) fCurrentVertex->Print();
+
+ delete recpoints;
+ itsLoader->UnloadRecPoints();
+
+ return fCurrentVertex;
+}
+//-------------------------------------------------------------------------
+void AliITSVertexerCosmics::FindVertices()
+{
+ // computes the vertices of the events in the range FirstEvent - LastEvent
+ AliRunLoader *rl = AliRunLoader::GetRunLoader();
+ AliITSLoader* itsLoader = (AliITSLoader*) rl->GetLoader("ITSLoader");
+ itsLoader->ReloadRecPoints();
+ for(Int_t i=fFirstEvent;i<=fLastEvent;i++){
+ // printf("Processing event %d\n",i);
+ rl->GetEvent(i);
+ FindVertexForCurrentEvent(i);
+ if(fCurrentVertex){
+ WriteCurrentVertex();
+ }
+ }
+}
+//-------------------------------------------------------------------------
+void AliITSVertexerCosmics::PrintStatus() const
+{
+ // Print current status
+ printf("=======================================================\n");
+ printf(" fMaxDistOnSPD2: %f\n",fMaxDistOnSPD2);
+ printf(" fMaxVtxRadius: %f\n",fMaxVtxRadius);
+ printf(" fMinDist2Vtxs: %f\n",fMinDist2Vtxs);
+ printf(" First layer first and last modules: %d, %d\n",fFirstSPD1,fLastSPD1);
+ printf(" Second layer first and last modules: %d, %d\n",fFirstSPD2,fLastSPD2);
+ printf("=======================================================\n");
+}
+//-------------------------------------------------------------------------
--- /dev/null
+#ifndef ALIITSVERTEXERCOSMICS_H
+#define ALIITSVERTEXERCOSMICS_H
+/* Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//#include "AliESDVertex.h"
+#include "AliITSVertexer.h"
+
+//-----------------------------------------------------------------------
+//
+// Class for constructing a fake primary vertex for cosmics events
+//
+// Origin: A.Dainese andrea.dainese@lnl.infn.it
+//-----------------------------------------------------------------------
+
+class AliESDVertex;
+
+class AliITSVertexerCosmics : public AliITSVertexer {
+
+ public:
+
+ AliITSVertexerCosmics();
+ virtual ~AliITSVertexerCosmics() {}
+ AliESDVertex* FindVertexForCurrentEvent(Int_t evnumb);
+ void FindVertices();
+ void PrintStatus() const;
+ void SetSPD1Modules(Int_t m1=0,Int_t m2=79) {fFirstSPD1 = m1; fLastSPD1 = m2;}
+ void SetSPD2Modules(Int_t m1=80, Int_t m2=239) {fFirstSPD2 = m1; fLastSPD2 = m2;}
+ void SetMaxDistOnSPD2(Double_t max=0.1) {fMaxDistOnSPD2=max;}
+ Double_t GetMaxDistOnSPD2() const {return fMaxDistOnSPD2;}
+ void SetMaxVtxRadius(Double_t maxr=3.5) {fMaxVtxRadius=maxr;}
+ Double_t GetMaVtxRadius() const {return fMaxVtxRadius;}
+ void SetMinDist2Vtxs(Double_t mind=0.1) {fMinDist2Vtxs=mind;}
+ Double_t GetMinDist2Vtxs() const {return fMinDist2Vtxs;}
+
+ private:
+
+ Int_t fFirstSPD1; // first module of the first pixel layer used
+ Int_t fLastSPD1; // last module of the first pixel layer used
+ Int_t fFirstSPD2; // first module of the second pixel layer used
+ Int_t fLastSPD2; // last module of the second pixel layer used
+ Double_t fMaxDistOnSPD2; // max dca between SPD1 tracklet and SPD2 cls
+ Double_t fMaxVtxRadius; // maximum radial pos of vertex
+ Double_t fMinDist2Vtxs; // minimum distance between two vertices
+
+ ClassDef(AliITSVertexerCosmics,1);
+};
+
+#endif
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
- * *
- * Author: The ALICE Off-line Project. *
- * Contributors are mentioned in the code where appropriate. *
- * *
- * Permission to use, copy, modify and distribute this software and its *
- * documentation strictly for non-commercial purposes is hereby granted *
- * without fee, provided that the above copyright notice appears in all *
- * copies and that both the copyright notice and this permission notice *
- * appear in the supporting documentation. The authors make no claims *
- * about the suitability of this software for any purpose. It is *
- * provided "as is" without express or implied warranty. *
- **************************************************************************/
-#include <AliITSVertexerPPZ.h>
-#include <TArrayF.h>
-#include <TH1F.h>
-#include <TObjArray.h>
-#include <TTree.h>
-#include <TClonesArray.h>
-#include "AliITSDetTypeRec.h"
-#include "AliITSgeom.h"
-#include "AliITSLoader.h"
-#include "AliITSRecPoint.h"
-/////////////////////////////////////////////////////////////////////////
-// //
-// This class is intended to compute the Z coordinate of the //
-// primary vertex for p-p interactions, or in general when the //
-// number of secondaries is too slow to use the class //
-// AliITSVertexerIons, which in turn is optimized for A-A collsions //
-// Origin: masera@to.infn.it 9/12/2002 //
-// //
-/////////////////////////////////////////////////////////////////////////
-ClassImp(AliITSVertexerPPZ)
-
-
-
-//______________________________________________________________________
-AliITSVertexerPPZ::AliITSVertexerPPZ():AliITSVertexer(),
-fFirstL1(0),
-fLastL1(0),
-fFirstL2(0),
-fLastL2(0),
-fDiffPhiMax(0),
-fX0(0),
-fY0(0),
-fZFound(0),
-fZsig(0),
-fWindow(0){
- // Default Constructor
-
- SetDiffPhiMax(0);
- SetFirstLayerModules(0);
- SetSecondLayerModules(0);
- //fITS = 0;
- SetWindow(0);
-}
-
-AliITSVertexerPPZ::AliITSVertexerPPZ(TString fn, Float_t x0, Float_t y0):AliITSVertexer(fn),
-fFirstL1(0),
-fLastL1(0),
-fFirstL2(0),
-fLastL2(0),
-fDiffPhiMax(0),
-fX0(x0),
-fY0(y0),
-fZFound(0),
-fZsig(0),
-fWindow(0) {
- // Standard constructor
- SetDiffPhiMax();
- fX0 = x0;
- fY0 = y0;
- SetFirstLayerModules();
- SetSecondLayerModules();
- //fITS = 0;
- SetWindow();
-}
-
-//______________________________________________________________________
-AliITSVertexerPPZ::~AliITSVertexerPPZ() {
- // Default Destructor
- //fITS = 0;
-}
-
-//________________________________________________________
-void AliITSVertexerPPZ::EvalZ(TH1F *hist,Int_t sepa, Int_t ncoinc, TArrayF *zval) {
-
- //Evaluation of Z
- Float_t deltaVal = hist->GetBinWidth(1)*fWindow; // max window in Z for searching
- fZFound=0;
- fZsig=0;
- Int_t nN=0;
- Int_t nbinNotZero=0;
- Float_t totst = 0.;
- Float_t totst2 = 0.;
- Float_t curz = 0.;
- for(Int_t i=1;i<=sepa;i++){
- Float_t cont=hist->GetBinContent(i);
- if(cont!=0)curz = hist->GetBinLowEdge(i)+0.5*hist->GetBinWidth(i);
- totst+=cont;
- totst2+=cont*cont;
- nN++;
- if(cont!=0)nbinNotZero++;
- }
- if(nbinNotZero==0){fZFound=-100; fZsig=-100; return;}
- if(nbinNotZero==1){
- fZFound = curz;
- fZsig=0;
- fCurrentVertex = new AliESDVertex(fZFound,fZsig,nbinNotZero);
- return;
- }
- Float_t errsq = totst2/(nN-1)-totst*totst/nN/(nN-1);
- if(errsq>=0){
- totst2=TMath::Sqrt(errsq);
- }
- else {
- Error("EvalZ","Negative variance: %d - %12.7f - %12.7f",nN,totst2,totst);
- fZFound=-100;
- fZsig=-100;
- return;
- }
- totst /= nN;
- Float_t cut = totst+totst2*2.;
- // if(fDebug>1)cout<<"totst, totst2, cut: "<<totst<<", "<<totst2<<", "<<cut<<endl;
- Float_t val1=hist->GetBinLowEdge(sepa);
- Float_t val2=hist->GetBinLowEdge(1);
- Float_t valm = 0.;
- Float_t zmax = 0.;
- for(Int_t i=1;i<=sepa;i++){
- Float_t cont=hist->GetBinContent(i);
- if(cont>valm){
- valm = cont;
- zmax = hist->GetBinLowEdge(i)+0.5*hist->GetBinWidth(1);
- }
- if(cont>cut){
- curz=hist->GetBinLowEdge(i);
- if(curz<val1)val1=curz;
- if(curz>val2)val2=curz;
- }
- }
- val2+=hist->GetBinWidth(1);
- if((val2-val1)>deltaVal){
- val1 = zmax-deltaVal/2.;
- val2 = zmax+deltaVal/2.;
- // if(fDebug>0)cout<<"val1 and val2 recomputed\n";
- }
- // if(fDebug>0)cout<<"Values for Z finding: "<<val1<<" "<<val2<<" "<<val2-val1<<endl;
- fZFound=0;
- fZsig=0;
- nN=0;
- for(Int_t i=0; i<ncoinc; i++){
- Float_t z=zval->At(i);
- if(z<val1)continue;
- if(z>val2)continue;
-
- fZFound+=z;
- fZsig+=z*z;
- /* weights defined by the curvature
- Float_t wei = 1./curv->At(i);
- fZFound+=z*wei;
- fZsig+=wei;
- */
- nN++;
- }
- if(nN<1){fZFound=-110; fZsig=-110; return;}
- if(nN==1){
- fZsig = 0;
- fCurrentVertex = new AliESDVertex(fZFound,fZsig,nN);
- return;
- }
- errsq = (fZsig/(nN-1)-fZFound*fZFound/nN/(nN-1))/nN;
- if(errsq>=0.){
- fZsig=TMath::Sqrt(errsq);
- }
- else {
- Error("evalZ","Negative variance: %d - %12.7f %12.7f",nN,fZsig,fZFound);
- fZsig=0;
- }
- fZFound=fZFound/nN;
- /* weights defined by the curvature
- fZsig=1./fZsig;
- fZFound*=fZsig;
- fZsig = TMath::Sqrt(fZsig);
- */
- fCurrentVertex = new AliESDVertex(fZFound,fZsig,nN);
-}
-
-//______________________________________________________________________
-AliESDVertex* AliITSVertexerPPZ::FindVertexForCurrentEvent(Int_t evnumber){
- // Defines the AliESDVertex for the current event
- fCurrentVertex = 0;
- fZFound = -999;
- fZsig = -999;
- AliRunLoader *rl =AliRunLoader::GetRunLoader();
- AliITSLoader* iTSloader = (AliITSLoader*)rl->GetLoader("ITSLoader");
- TDirectory * olddir = gDirectory;
- rl->CdGAFile();
- AliITSgeom* geom = (AliITSgeom*)gDirectory->Get("AliITSgeom");
- olddir->cd();
-
- AliITSDetTypeRec detTypeRec;
-
- if(!geom) {
- Error("FindVertexForCurrentEvent","ITS geometry is not defined");
- return fCurrentVertex;
- }
- TTree *tR=0;
- TClonesArray *itsRec = 0;
- Float_t lc[3]; for(Int_t ii=0; ii<3; ii++) lc[ii]=0.;
- Float_t gc[3]; for(Int_t ii=0; ii<3; ii++) gc[ii]=0.;
- Float_t lc2[3]; for(Int_t ii=0; ii<3; ii++) lc2[ii]=0.;
- Float_t gc2[3]; for(Int_t ii=0; ii<3; ii++) gc2[ii]=0.;
-
- tR = iTSloader->TreeR();
- if(!tR){
- Error("FindVertexForCurrentEvent","TreeR not found");
- return fCurrentVertex;
- }
- detTypeRec.SetTreeAddressR(tR);
- itsRec = detTypeRec.RecPoints();
- // missing
- // TClonesArray dummy("AliITSRecPoint",10000), *clusters=&dummy;
- TBranch *branch;
- branch = tR->GetBranch("ITSRecPoints");
- if(!branch){
- branch = tR->GetBranch("ITSRecPointsF");
- }
- //}
- if(!branch){
- Error("FindVertexForCurrentEvent","branch for ITS rec points not found");
- return fCurrentVertex;
- }
- Float_t zave=0;
- Float_t rmszav=0;
- Float_t zave2=0;
- Int_t firipixe=0;
- for(Int_t module= fFirstL1; module<=fLastL1;module++){
- branch->GetEvent(module);
- Int_t nrecp1 = itsRec->GetEntries();
- for(Int_t i=0; i<nrecp1;i++){
- AliITSRecPoint *current = (AliITSRecPoint*)itsRec->At(i);
- lc[0]=current->GetDetLocalX();
- lc[2]=current->GetDetLocalZ();
- geom->LtoG(module,lc,gc);
- zave+=gc[2];
- zave2+=gc[2]*gc[2];
- firipixe++;
- }
- detTypeRec.ResetRecPoints();
- }
- if(firipixe>1){
- rmszav=TMath::Sqrt(zave2/(firipixe-1)-zave*zave/firipixe/(firipixe-1));
- zave=zave/firipixe;
- // if(fDebug>1)cout<<"++++++++++++++++++++++++++++++++++++++++++++++++++++++\n Number of firing pixels: "<<firipixe<<endl;
- }
- else {
- fZFound = -200;
- fZsig = -200;
- Warning("FindVertexForCurrentEvent","No rec points on first layer for this event");
- return fCurrentVertex;
- }
- Float_t zlim1=zave-rmszav;
- Float_t zlim2=zave+rmszav;
- Int_t sepa=(Int_t)((zlim2-zlim1)*10.+1.);
- zlim2=zlim1 + sepa/10.;
- TH1F *zvdis = new TH1F("z_ev","zv distr",sepa,zlim1,zlim2);
- /*
- if(fDebug>0){
- cout<<"Z limits: "<<zlim1<<" "<<zlim2<<"; Bins= "<<sepa<<endl;
- cout<<"Bin width: "<<zvdis->GetBinWidth(1)*10000.<<" micron\n";
- }
- */
- Int_t sizarr=100;
- TArrayF *zval = new TArrayF(sizarr);
- // TArrayF *curv = new TArrayF(sizarr);
- Int_t ncoinc=0;
- for(Int_t module= fFirstL1; module<=fLastL1;module++){
- // if(fDebug>0)cout<<"processing module "<<module<<" \r";
- branch->GetEvent(module);
- Int_t nrecp1 = itsRec->GetEntries();
- TObjArray *poiL1 = new TObjArray(nrecp1);
- for(Int_t i=0; i<nrecp1;i++)poiL1->AddAt(itsRec->At(i),i);
- detTypeRec.ResetRecPoints();
- for(Int_t i=0; i<nrecp1;i++){
- AliITSRecPoint *current = (AliITSRecPoint*)poiL1->At(i);
- lc[0]=current->GetDetLocalX();
- lc[2]=current->GetDetLocalZ();
- geom->LtoG(module,lc,gc);
- gc[0]-=fX0;
- gc[1]-=fY0;
- Float_t r1=TMath::Sqrt(gc[0]*gc[0]+gc[1]*gc[1]);
- Float_t phi1 = TMath::ATan2(gc[1],gc[0]);
- if(phi1<0)phi1=2*TMath::Pi()+phi1;
- // if(fDebug>1)cout<<"module "<<module<<" "<<gc[0]<<" "<<gc[1]<<" "<<gc[2]<<" "<<phi1<<" \n";
- for(Int_t modul2=fFirstL2; modul2<=fLastL2; modul2++){
- branch->GetEvent(modul2);
- Int_t nrecp2 = itsRec->GetEntries();
- for(Int_t j=0; j<nrecp2;j++){
- AliITSRecPoint *recp = (AliITSRecPoint*)itsRec->At(j);
- lc2[0]=recp->GetDetLocalX();
- lc2[2]=recp->GetDetLocalZ();
- geom->LtoG(modul2,lc2,gc2);
- gc2[0]-=fX0;
- gc2[1]-=fY0;
- Float_t r2=TMath::Sqrt(gc2[0]*gc2[0]+gc2[1]*gc2[1]);
- Float_t zr0=(r2*gc[2]-r1*gc2[2])/(r2-r1);
- Float_t phi2 = TMath::ATan2(gc2[1],gc2[0]);
- if(phi2<0)phi2=2.*TMath::Pi()+phi2;
- Float_t diff = TMath::Abs(phi2-phi1);
- if(diff>TMath::Pi())diff=2.*TMath::Pi()-diff;
- if(zr0>zlim1 && zr0<zlim2){
- if(diff<fDiffPhiMax ){
- zvdis->Fill(zr0);
- zval->AddAt(zr0,ncoinc);
- /* uncomment these lines to use curvature as a weight
- Float_t cu = Curv(0.,0.,gc[0],gc[1],gc2[0],gc[1]);
- curv->AddAt(cu,ncoinc);
- */
- ncoinc++;
- if(ncoinc==(sizarr-1)){
- sizarr+=100;
- zval->Set(sizarr);
- //uncomment next line to use curvature as weight
- // curv->Set(sizarr);
- }
- }
- }
- }
- detTypeRec.ResetRecPoints();
- }
- }
- delete poiL1;
- } // loop on modules
- /*
- if(fDebug>0){
- cout<<endl<<"Number of coincidences = "<<ncoinc<<endl;
- }
- */
- // EvalZ(zvdis,sepa,ncoinc,zval,curv);
- EvalZ(zvdis,sepa,ncoinc,zval);
- delete zvdis;
- delete zval;
- // delete curv;
- if(fCurrentVertex){
- char name[30];
- sprintf(name,"Vertex_%d",evnumber);
- // fCurrentVertex->SetName(name);
- fCurrentVertex->SetTitle("vertexer: PPZ");
- }
-
- return fCurrentVertex;
-}
-
-//______________________________________________________________________
-void AliITSVertexerPPZ::FindVertices(){
- // computes the vertices of the events in the range FirstEvent - LastEvent
- AliRunLoader *rl = AliRunLoader::GetRunLoader();
- AliITSLoader* iTSloader = (AliITSLoader*) rl->GetLoader("ITSLoader");
- iTSloader->ReloadRecPoints();
- for(Int_t i=fFirstEvent;i<=fLastEvent;i++){
- cout<<"Processing event "<<i<<endl;
- rl->GetEvent(i);
- FindVertexForCurrentEvent(i);
- if(fCurrentVertex){
- WriteCurrentVertex();
- }
- else {
- /*
- if(fDebug>0){
- cout<<"Vertex not found for event "<<i<<endl;
- cout<<"fZFound = "<<fZFound<<", fZsig= "<<fZsig<<endl;
- }
- */
- }
- }
-}
-
-//________________________________________________________
-void AliITSVertexerPPZ::PrintStatus() const {
- // Print current status
- cout <<"=======================================================\n";
- cout <<" First layer first and last modules: "<<fFirstL1<<", ";
- cout <<fLastL1<<endl;
- cout <<" Second layer first and last modules: "<<fFirstL2<<", ";
- cout <<fLastL2<<endl;
- cout <<" Max Phi difference: "<<fDiffPhiMax<<endl;
- cout <<" Window for Z search: "<<fWindow<<endl;
- cout <<" Current Z "<<fZFound<<"; Z sig "<<fZsig<<endl;
- cout<<"First event to be processed "<<fFirstEvent;
- cout<<"\n Last event to be processed "<<fLastEvent<<endl;
-}
-
-//________________________________________________________
-Float_t AliITSVertexerPPZ::Curv(Double_t x1,Double_t y1,
- Double_t x2,Double_t y2,
- Double_t x3,Double_t y3)
-{
- //-----------------------------------------------------------------
- // Initial approximation of the track curvature (Y. Belikov) squared
- //-----------------------------------------------------------------
- Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1);
- Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)-
- (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2));
- Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)-
- (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1));
-
- Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b);
-
- Float_t val = (Float_t)(( -xr*yr/TMath::Sqrt(xr*xr+yr*yr))
- *( -xr*yr/TMath::Sqrt(xr*xr+yr*yr)));
- return val;
-}
+++ /dev/null
-#ifndef ALIITSVERTEXERPPZ_H
-#define ALIITSVERTEXERPPZ_H
-/* Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice */
-
-
-#include <AliITSVertexer.h>
-
-/////////////////////////////////////////////////////////////////////
-// //
-// Class for primary vertex Z coordinate reconstruction //
-// Optimized for p-p events (in general: low multiplicity events) //
-// //
-/////////////////////////////////////////////////////////////////////
-
-class TH1F;
-class TArrayF;
-
-class AliITSVertexerPPZ : public AliITSVertexer {
-
- public:
- AliITSVertexerPPZ();
- AliITSVertexerPPZ(TString fn, Float_t x0=0., Float_t y0=0.); // standard constructor
- virtual ~AliITSVertexerPPZ(); // destructor
- virtual AliESDVertex* FindVertexForCurrentEvent(Int_t event);
- virtual void FindVertices();
- virtual Float_t GetZFound() const {return fZFound;}
- virtual Float_t GetZsig() const {return fZsig;}
- virtual void PrintStatus() const;
- virtual void SetDiffPhiMax(Float_t pm = 0.05){fDiffPhiMax = pm;}
- virtual void SetFirstLayerModules(Int_t m1 = 0, Int_t m2 = 79){fFirstL1 = m1; fLastL1 = m2;}
- virtual void SetSecondLayerModules(Int_t m1 = 80, Int_t m2 = 239){fFirstL2 = m1; fLastL2 = m2;}
- virtual void SetWindow(Float_t w=3.){fWindow = w;}
- static Float_t Curv(Double_t x1,Double_t y1, Double_t x2,Double_t y2,
- Double_t x3,Double_t y3);
-
-
- protected:
- Int_t fFirstL1; // first module of the first pixel layer
- Int_t fLastL1; // last module of the first pixel layer
- Int_t fFirstL2; // first module of the second pixel layer
- Int_t fLastL2; // last module of the second pixel layer
- Float_t fDiffPhiMax; // Maximum delta phi allowed among corr. pixels
- Float_t fX0; // Nominal x coordinate of the vertex
- Float_t fY0; // Nominal y coordinate of the vertex
- //AliITS *fITS; //! pointer to the AliITS object
- Float_t fZFound; //! found value for the current event
- Float_t fZsig; //! RMS of Z
- Float_t fWindow; // window width for Z search in mm (3 mm by def.)
-
- private:
- void EvalZ(TH1F *hist,Int_t sepa, Int_t ncoinc, TArrayF *zval);
-
- ClassDef(AliITSVertexerPPZ,3);
-};
-
-#endif
#pragma link C++ class AliITSVertexer+;
#pragma link C++ class AliITSVertexerIons+;
-#pragma link C++ class AliITSVertexerPPZ+;
+#pragma link C++ class AliITSVertexerCosmics+;
#pragma link C++ class AliITSVertexerZ+;
#pragma link C++ class AliITSVertexer3D+;
#pragma link C++ class AliITSZPoint+;
AliITStrackerMI.cxx \
AliITSVertexer.cxx \
AliITSVertexerIons.cxx \
- AliITSVertexerPPZ.cxx \
+ AliITSVertexerCosmics.cxx \
AliITSVertexer3D.cxx \
AliITSVertexerZ.cxx \
AliITSVertexerFast.cxx \