X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=ITS%2FAliITSVertexer3D.cxx;h=7b9acf4b9245115b98eea09037f5df7753d5652d;hb=960f0a954a0b69dcb4b9ce990a58f97f1ddb0d34;hp=bfa6a1916cea63426cf0eae01cbf4d755c9e1cae;hpb=87104b066948c2f9a241252a88d8f779299d35fb;p=u%2Fmrichter%2FAliRoot.git diff --git a/ITS/AliITSVertexer3D.cxx b/ITS/AliITSVertexer3D.cxx index bfa6a1916ce..7b9acf4b924 100644 --- a/ITS/AliITSVertexer3D.cxx +++ b/ITS/AliITSVertexer3D.cxx @@ -13,12 +13,14 @@ * provided "as is" without express or implied warranty. * **************************************************************************/ #include +#include "AliRunLoader.h" #include "AliESDVertex.h" #include "AliLog.h" #include "AliStrLine.h" #include "AliTracker.h" #include "AliITSDetTypeRec.h" #include "AliITSRecPoint.h" +#include "AliITSRecPointContainer.h" #include "AliITSgeomTGeo.h" #include "AliVertexerTracks.h" #include "AliITSVertexer3D.h" @@ -27,10 +29,12 @@ ///////////////////////////////////////////////////////////////// // this class implements a method to determine // the 3 coordinates of the primary vertex -// for p-p collisions -// It can be used successfully with Pb-Pb collisions +// optimized for +// p-p collisions //////////////////////////////////////////////////////////////// +const Int_t AliITSVertexer3D::fgkMaxNumOfClDefault = 500; + ClassImp(AliITSVertexer3D) /* $Id$ */ @@ -58,7 +62,10 @@ AliITSVertexer3D::AliITSVertexer3D(): fDiffPhiforPileup(0.), fBinSizeR(0.), fBinSizeZ(0.), - fPileupAlgo(0) + fPileupAlgo(0), + fMaxNumOfCl(fgkMaxNumOfClDefault), + fDoDownScale(kFALSE), + fGenerForDownScale(0) { // Default constructor SetCoarseDiffPhiCut(); @@ -81,6 +88,7 @@ AliITSVertexer3D::AliITSVertexer3D(): Double_t binsize=0.02; // default 200 micron Int_t nbins=static_cast(1+2*fZCutDiamond/binsize); fZHisto=new TH1F("hz","",nbins,-fZCutDiamond,-fZCutDiamond+binsize*nbins); + fGenerForDownScale=new TRandom3(987654321); } //______________________________________________________________________ @@ -88,6 +96,7 @@ AliITSVertexer3D::~AliITSVertexer3D() { // Destructor fLines.Clear("C"); if(fZHisto) delete fZHisto; + if(fGenerForDownScale) delete fGenerForDownScale; } //______________________________________________________________________ @@ -115,11 +124,11 @@ AliESDVertex* AliITSVertexer3D::FindVertexForCurrentEvent(TTree *itsClusterTree) rc=Prepare3DVertex(0); if(fVert3D.GetNContributors()>0){ fLines.Clear("C"); - Int_t nolines = FindTracklets(itsClusterTree,1); + nolines = FindTracklets(itsClusterTree,1); if(nolines>=2){ rc=Prepare3DVertex(1); if(fPileupAlgo == 2 && rc == 0) FindVertex3DIterative(); - else if(fPileupAlgo<2 && rc == 0) FindVertex3D(itsClusterTree); + else if(fPileupAlgo!=2 && rc == 0) FindVertex3D(itsClusterTree); if(rc!=0) fVert3D.SetNContributors(0); // exclude this vertex } } @@ -139,13 +148,14 @@ AliESDVertex* AliITSVertexer3D::FindVertexForCurrentEvent(TTree *itsClusterTree) Double_t chi2=99999.; Int_t nContr=vtxz->GetNContributors(); fCurrentVertex = new AliESDVertex(position,covmatrix,chi2,nContr); + fCurrentVertex->SetDispersion(vtxz->GetDispersion()); fCurrentVertex->SetTitle("vertexer: Z"); fCurrentVertex->SetName("SPDVertexZ"); delete vtxz; } - } - FindMultiplicity(itsClusterTree); + } + if(fComputeMultiplicity) FindMultiplicity(itsClusterTree); return fCurrentVertex; } @@ -168,6 +178,7 @@ void AliITSVertexer3D::FindVertex3D(TTree *itsClusterTree){ switch(fPileupAlgo){ case 0: PileupFromZ(); break; case 1: FindOther3DVertices(itsClusterTree); break; + case 3: break; // no pileup algo default: AliError("Wrong pileup algorithm"); break; } if(fNoVertices==1){ @@ -467,11 +478,6 @@ Int_t AliITSVertexer3D::FindTracklets(TTree *itsClusterTree, Int_t optCuts){ // considered. Straight lines (=tracklets)are formed. // The tracklets are processed in Prepare3DVertex - if(!GetDetTypeRec())AliFatal("DetTypeRec pointer has not been set"); - - TTree *tR = itsClusterTree; - fDetTypeRec->ResetRecPoints(); - fDetTypeRec->SetTreeAddressR(tR); TClonesArray *itsRec = 0; if(optCuts==0) fZHisto->Reset(); // gc1 are local and global coordinates for layer 1 @@ -480,12 +486,10 @@ Int_t AliITSVertexer3D::FindTracklets(TTree *itsClusterTree, Int_t optCuts){ // gc2 are local and global coordinates for layer 2 Float_t gc2f[3]={0.,0.,0.}; Double_t gc2[3]={0.,0.,0.}; - - itsRec = fDetTypeRec->RecPoints(); - TBranch *branch = NULL; - branch = tR->GetBranch("ITSRecPoints"); - if(!branch){ - AliWarning("Null pointer for RecPoints branch, vertex not calculated"); + AliITSRecPointContainer* rpcont=AliITSRecPointContainer::Instance(); + itsRec=rpcont->FetchClusters(0,itsClusterTree); + if(!rpcont->IsSPDActive()){ + AliWarning("No SPD rec points found, 3D vertex not calculated"); return -1; } @@ -516,28 +520,31 @@ Int_t AliITSVertexer3D::FindTracklets(TTree *itsClusterTree, Int_t optCuts){ Int_t nrpL1 = 0; // number of rec points on layer 1 Int_t nrpL2 = 0; // number of rec points on layer 2 - - // By default irstL1=0 and lastL1=79 - Int_t firstL1 = AliITSgeomTGeo::GetModuleIndex(1,1,1); - Int_t lastL1 = AliITSgeomTGeo::GetModuleIndex(2,1,1)-1; - for(Int_t module= firstL1; module<=lastL1;module++){ // count number of recopints on layer 1 - branch->GetEvent(module); - nrpL1+= itsRec->GetEntries(); - fDetTypeRec->ResetRecPoints(); - } - //By default firstL2=80 and lastL2=239 - Int_t firstL2 = AliITSgeomTGeo::GetModuleIndex(2,1,1); - Int_t lastL2 = AliITSgeomTGeo::GetModuleIndex(3,1,1)-1; - for(Int_t module= firstL2; module<=lastL2;module++){ // count number of recopints on layer 2 - branch->GetEvent(module); - nrpL2+= itsRec->GetEntries(); - fDetTypeRec->ResetRecPoints(); - } + nrpL1=rpcont->GetNClustersInLayerFast(1); + nrpL2=rpcont->GetNClustersInLayerFast(2); if(nrpL1 == 0 || nrpL2 == 0){ AliDebug(1,Form("No RecPoints in at least one SPD layer (%d %d)",nrpL1,nrpL2)); return -1; } AliDebug(1,Form("RecPoints on Layer 1,2 = %d, %d\n",nrpL1,nrpL2)); + fDoDownScale=kFALSE; + Float_t factDownScal=1.; + Int_t origLaddersOnLayer2=fLadOnLay2; + + if(nrpL1>fMaxNumOfCl || nrpL2>fMaxNumOfCl){ + SetLaddersOnLayer2(2); + fDoDownScale=kTRUE; + factDownScal=(Float_t)fMaxNumOfCl*(Float_t)fMaxNumOfCl/(Float_t)nrpL1/(Float_t)nrpL2; + if(optCuts==1){ + factDownScal*=(fCoarseDiffPhiCut/fDiffPhiMax)*10; + if(factDownScal>1.){ + fDoDownScale=kFALSE; + SetLaddersOnLayer2(origLaddersOnLayer2); + } + } + else if(optCuts==2) return -1; + if(fDoDownScale)AliDebug(1,Form("Too many recpoints on SPD(%d %d ), downscale by %f",nrpL1,nrpL2,factDownScal)); + } Double_t a[3]={xbeam,ybeam,0.}; Double_t b[3]={xbeam,ybeam,10.}; @@ -547,23 +554,20 @@ Int_t AliITSVertexer3D::FindTracklets(TTree *itsClusterTree, Int_t optCuts){ Int_t nolines = 0; // Loop on modules of layer 1 + Int_t firstL1 = AliITSgeomTGeo::GetModuleIndex(1,1,1); + Int_t lastL1 = AliITSgeomTGeo::GetModuleIndex(2,1,1)-1; for(Int_t modul1= firstL1; modul1<=lastL1;modul1++){ // Loop on modules of layer 1 if(!fUseModule[modul1]) continue; - UShort_t ladder=int(modul1/4)+1; // ladders are numbered starting from 1 - branch->GetEvent(modul1); - Int_t nrecp1 = itsRec->GetEntries(); - static TClonesArray prpl1("AliITSRecPoint",nrecp1); - prpl1.SetOwner(); - for(Int_t j=0;jAt(j); - new(prpl1[j])AliITSRecPoint(*recp); - } - fDetTypeRec->ResetRecPoints(); + + UShort_t ladder=modul1/4+1; // ladders are numbered starting from 1 + TClonesArray *prpl1=rpcont->UncheckedGetClusters(modul1); + Int_t nrecp1 = prpl1->GetEntries(); for(Int_t j=0;jkMaxCluPerMod) continue; UShort_t idClu1=modul1*kMaxCluPerMod+j; if(fUsedCluster.TestBitNumber(idClu1)) continue; - AliITSRecPoint *recp1 = (AliITSRecPoint*)prpl1.At(j); + if(fDoDownScale && fGenerForDownScale->Rndm()>factDownScal) continue; + AliITSRecPoint *recp1 = (AliITSRecPoint*)prpl1->At(j); recp1->GetGlobalXYZ(gc1f); for(Int_t ico=0;ico<3;ico++)gc1[ico]=gc1f[ico]; @@ -575,12 +579,13 @@ Int_t AliITSVertexer3D::FindTracklets(TTree *itsClusterTree, Int_t optCuts){ if(ladmod>AliITSgeomTGeo::GetNLadders(2)) ladmod=ladmod-AliITSgeomTGeo::GetNLadders(2); Int_t modul2=AliITSgeomTGeo::GetModuleIndex(2,ladmod,k+1); if(!fUseModule[modul2]) continue; - branch->GetEvent(modul2); + itsRec=rpcont->UncheckedGetClusters(modul2); Int_t nrecp2 = itsRec->GetEntries(); for(Int_t j2=0;j2kMaxCluPerMod) continue; UShort_t idClu2=modul2*kMaxCluPerMod+j2; if(fUsedCluster.TestBitNumber(idClu2)) continue; + AliITSRecPoint *recp2 = (AliITSRecPoint*)itsRec->At(j2); recp2->GetGlobalXYZ(gc2f); for(Int_t ico=0;ico<3;ico++)gc2[ico]=gc2f[ico]; @@ -607,6 +612,7 @@ Int_t AliITSVertexer3D::FindTracklets(TTree *itsClusterTree, Int_t optCuts){ Double_t deltaZ=cp[2]-zvert; if(TMath::Abs(deltaZ)>dZmax)continue; + if(nolines == 0){ if(fLines.GetEntriesFast()>0)fLines.Clear("C"); } @@ -652,12 +658,16 @@ Int_t AliITSVertexer3D::FindTracklets(TTree *itsClusterTree, Int_t optCuts){ new(fLines[nolines++])AliStrLine(gc1,sigmasq,wmat,gc2,kTRUE,idClu1,idClu2); } - fDetTypeRec->ResetRecPoints(); } } } - prpl1.Clear(); } + + if(fDoDownScale){ + SetLaddersOnLayer2(origLaddersOnLayer2); + } + + if(nolines == 0)return -2; return nolines; } @@ -688,6 +698,13 @@ Int_t AliITSVertexer3D::Prepare3DVertex(Int_t optCuts){ deltaR=fMaxRCut; } + Double_t origBinSizeR=fBinSizeR; + Double_t origBinSizeZ=fBinSizeZ; + if(fDoDownScale){ + SetBinSizeR(0.05); + SetBinSizeZ(0.05); + } + Double_t rl=-fCoarseMaxRCut; Double_t rh=fCoarseMaxRCut; Double_t zl=-fZCutDiamond; @@ -770,6 +787,47 @@ Int_t AliITSVertexer3D::Prepare3DVertex(Int_t optCuts){ } delete h3dcs; + // Finer Histo in limited range in case of high mult. + if(fDoDownScale){ + SetBinSizeR(0.01); + SetBinSizeZ(0.01); + Double_t xl=peak[0]-0.3; + Double_t xh=peak[0]+0.3; + Double_t yl=peak[1]-0.3; + Double_t yh=peak[1]+0.3; + zl=peak[2]-0.5; + zh=peak[2]+0.5; + Int_t nbxfs=(Int_t)((xh-xl)/fBinSizeR+0.0001); + Int_t nbyfs=(Int_t)((yh-yl)/fBinSizeR+0.0001); + Int_t nbzfs=(Int_t)((zh-zl)/fBinSizeZ+0.0001); + + TH3F *h3dfs = new TH3F("h3dfs","xyz distribution",nbxfs,xl,xh,nbyfs,yl,yh,nbzfs,zl,zh); + for(Int_t i=0; iCross(l2,point); + if(retc<0)continue; + h3dfs->Fill(point[0],point[1],point[2]); + } + } + ntrkl=0; + ntimes=0; + + Double_t newpeak[3]={0.,0.,0.}; + FindPeaks(h3dfs,newpeak,ntrkl,ntimes); + if(ntimes==1){ + for(Int_t iCoo=0; iCoo<3; iCoo++) peak[iCoo]=newpeak[iCoo]; + binsizer=fBinSizeR; + binsizez=fBinSizeZ; + } + delete h3dfs; + SetBinSizeR(origBinSizeR); + SetBinSizeZ(origBinSizeZ); + } + + // Second selection loop Double_t bs=(binsizer+binsizez)/2.; @@ -787,10 +845,33 @@ Int_t AliITSVertexer3D::Prepare3DVertex(Int_t optCuts){ // make a further selection on tracklets based on this first candidate fVert3D.GetXYZ(peak); AliDebug(1,Form("FIRST V candidate: %f ; %f ; %f",peak[0],peak[1],peak[2])); + Int_t *validate2 = new Int_t [fLines.GetEntriesFast()]; + for(Int_t i=0; iGetDistFromPoint(peak)> fDCAcut)fLines.RemoveAt(i); + if(optCuts==2){ // temporarily only for pileup + for(Int_t j=i+1; jGetDCA(l2)<0.00001){ + Int_t delta1=(Int_t)l1->GetIdPoint(0)-(Int_t)l2->GetIdPoint(0); + Int_t delta2=(Int_t)l1->GetIdPoint(1)-(Int_t)l2->GetIdPoint(1); + Int_t deltamod1=(Int_t)l1->GetIdPoint(0)/kMaxCluPerMod + -(Int_t)l2->GetIdPoint(0)/kMaxCluPerMod; + Int_t deltamod2=(Int_t)l1->GetIdPoint(1)/kMaxCluPerMod + -(Int_t)l2->GetIdPoint(1)/kMaxCluPerMod; + // remove tracklets sharing a point + if( (delta1==0 && deltamod2==0) || + (delta2==0 && deltamod1==0) ) validate2[j]=0; + } + } + } } + for(Int_t i=0; i1){// this new tracklet selection is used @@ -983,5 +1064,7 @@ void AliITSVertexer3D::PrintStatus() const { printf("Pileup algo: %d\n",fPileupAlgo); printf("Min DCA to 1st vertex for pileup (algo 0 and 1): %f\n",fDCAforPileup); printf("Cut on distance between pair-vertices (algo 2): %f\n",fCutOnPairs); + printf("Maximum number of clusters allowed on L1 or L2: %d\n",fMaxNumOfCl); printf("=======================================================\n"); } +