From a98c78e57bdccb2bbe4b5cf77bfe3999e7b303c6 Mon Sep 17 00:00:00 2001 From: morsch Date: Fri, 9 Jul 2010 10:57:29 +0000 Subject: [PATCH] Pile-up vertices added to AODEvent. (F. Prino) --- ANALYSIS/AliAnalysisTaskESDfilter.cxx | 40 +++++++++ STEER/AliAODEvent.cxx | 114 +++++++++++++++++++++++++- STEER/AliAODEvent.h | 16 +++- STEER/AliESDEvent.cxx | 80 ++++++++---------- STEER/AliESDEvent.h | 8 +- STEER/AliVEvent.h | 9 +- 6 files changed, 217 insertions(+), 50 deletions(-) diff --git a/ANALYSIS/AliAnalysisTaskESDfilter.cxx b/ANALYSIS/AliAnalysisTaskESDfilter.cxx index b061b04066b..323316f461c 100644 --- a/ANALYSIS/AliAnalysisTaskESDfilter.cxx +++ b/ANALYSIS/AliAnalysisTaskESDfilter.cxx @@ -187,12 +187,17 @@ void AliAnalysisTaskESDfilter::ConvertESDtoAOD() { Float_t diamxy[2]={esd->GetDiamondX(),esd->GetDiamondY()}; Float_t diamcov[3]; esd->GetDiamondCovXY(diamcov); header->SetDiamond(diamxy,diamcov); + header->SetDiamondZ(esd->GetDiamondZ(),esd->GetSigma2DiamondZ()); // // Int_t nV0s = esd->GetNumberOfV0s(); Int_t nCascades = esd->GetNumberOfCascades(); Int_t nKinks = esd->GetNumberOfKinks(); Int_t nVertices = nV0s + nCascades /*V0 wihtin cascade already counted*/+ nKinks + 1 /* = prim. vtx*/; + Int_t nPileSPDVertices=1+esd->GetNumberOfPileupVerticesSPD(); // also SPD main vertex + Int_t nPileTrkVertices=esd->GetNumberOfPileupVerticesTracks(); + nVertices+=nPileSPDVertices; + nVertices+=nPileTrkVertices; Int_t nJets = 0; Int_t nCaloClus = esd->GetNumberOfCaloClusters(); Int_t nFmdClus = 0; @@ -280,6 +285,41 @@ void AliAnalysisTaskESDfilter::ConvertESDtoAOD() { if (fDebug > 0) primary->Print(); + + // Add SPD "main" vertex + const AliESDVertex *vtxS = esd->GetPrimaryVertexSPD(); + vtxS->GetXYZ(pos); // position + vtxS->GetCovMatrix(covVtx); //covariance matrix + AliAODVertex * mVSPD = new(vertices[jVertices++]) + AliAODVertex(pos, covVtx, vtxS->GetChi2toNDF(), NULL, -1, AliAODVertex::kMainSPD); + mVSPD->SetName(vtxS->GetName()); + mVSPD->SetTitle(vtxS->GetTitle()); + mVSPD->SetNContributors(vtxS->GetNContributors()); + + // Add SPD pileup vertices + for(Int_t iV=0; iVGetPileupVertexSPD(iV); + vtxP->GetXYZ(pos); // position + vtxP->GetCovMatrix(covVtx); //covariance matrix + AliAODVertex * pVSPD = new(vertices[jVertices++]) + AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupSPD); + pVSPD->SetName(vtxP->GetName()); + pVSPD->SetTitle(vtxP->GetTitle()); + pVSPD->SetNContributors(vtxP->GetNContributors()); + } + + // Add TRK pileup vertices + for(Int_t iV=0; iVGetPileupVertexTracks(iV); + vtxP->GetXYZ(pos); // position + vtxP->GetCovMatrix(covVtx); //covariance matrix + AliAODVertex * pVTRK = new(vertices[jVertices++]) + AliAODVertex(pos, covVtx, vtxP->GetChi2toNDF(), NULL, -1, AliAODVertex::kPileupTracks); + pVTRK->SetName(vtxP->GetName()); + pVTRK->SetTitle(vtxP->GetTitle()); + pVTRK->SetNContributors(vtxP->GetNContributors()); + } + // Create vertices starting from the most complex objects Double_t chi2 = 0.; diff --git a/STEER/AliAODEvent.cxx b/STEER/AliAODEvent.cxx index 4c4110c5d47..4dbb5dd6c8b 100644 --- a/STEER/AliAODEvent.cxx +++ b/STEER/AliAODEvent.cxx @@ -373,7 +373,7 @@ void AliAODEvent::ResetStd(Int_t trkArrSize, fVertices->Delete(); if (vtxArrSize > fVertices->GetSize()) fVertices->Expand(vtxArrSize); - + fV0s->Delete(); if (v0ArrSize > fV0s->GetSize()) fV0s->Expand(v0ArrSize); @@ -618,6 +618,118 @@ void AliAODEvent::ReadFromTree(TTree *tree, Option_t* opt /*= ""*/) fAODObjects->SetOwner(kTRUE); } } +//______________________________________________________________________________ +Int_t AliAODEvent::GetNumberOfPileupVerticesSPD() const{ + // count number of SPD pileup vertices + Int_t nVertices=GetNumberOfVertices(); + Int_t nPileupVertices=0; + for(Int_t iVert=0; iVertGetType()==AliAODVertex::kPileupSPD) nPileupVertices++; + } + return nPileupVertices; +} +//______________________________________________________________________________ +Int_t AliAODEvent::GetNumberOfPileupVerticesTracks() const{ + // count number of track pileup vertices + Int_t nVertices=GetNumberOfVertices(); + Int_t nPileupVertices=0; + for(Int_t iVert=0; iVertGetType()==AliAODVertex::kPileupTracks) nPileupVertices++; + } + return nPileupVertices; +} +//______________________________________________________________________________ +AliAODVertex* AliAODEvent::GetPrimaryVertexSPD() const{ + // + Int_t nVertices=GetNumberOfVertices(); + for(Int_t iVert=0; iVertGetType()==AliAODVertex::kMainSPD) return v; + } + return 0; +} +//______________________________________________________________________________ +AliAODVertex* AliAODEvent::GetPileupVertexSPD(Int_t iV) const{ + // + Int_t nVertices=GetNumberOfVertices(); + Int_t counter=0; + for(Int_t iVert=0; iVertGetType()==AliAODVertex::kPileupSPD){ + if(counter==iV) return v; + ++counter; + } + } + return 0; +} +//______________________________________________________________________________ +AliAODVertex* AliAODEvent::GetPileupVertexTracks(Int_t iV) const{ + // + Int_t nVertices=GetNumberOfVertices(); + Int_t counter=0; + for(Int_t iVert=0; iVertGetType()==AliAODVertex::kPileupTracks){ + if(counter==iV) return v; + ++counter; + } + } + return 0; +} +//______________________________________________________________________________ +Bool_t AliAODEvent::IsPileupFromSPD(Int_t minContributors, + Double_t minZdist, + Double_t nSigmaZdist, + Double_t nSigmaDiamXY, + Double_t nSigmaDiamZ) const{ + // + // This function checks if there was a pile up + // reconstructed with SPD + // + AliAODVertex *mainV=GetPrimaryVertexSPD(); + if(!mainV) return kFALSE; + Int_t nc1=mainV->GetNContributors(); + if(nc1<1) return kFALSE; + Int_t nPileVert=GetNumberOfPileupVerticesSPD(); + if(nPileVert==0) return kFALSE; + Int_t nVertices=GetNumberOfVertices(); + + for(Int_t iVert=0; iVertGetType()!=AliAODVertex::kPileupSPD) continue; + Int_t nc2=pv->GetNContributors(); + if(nc2>=minContributors){ + Double_t z1=mainV->GetZ(); + Double_t z2=pv->GetZ(); + Double_t distZ=TMath::Abs(z2-z1); + Double_t distZdiam=TMath::Abs(z2-GetDiamondZ()); + Double_t cutZdiam=nSigmaDiamZ*GetSigma2DiamondZ(); + if(GetSigma2DiamondZ()<0.0001)cutZdiam=99999.; //protection for missing z diamond information + if(distZ>minZdist && distZdiamGetX(); + Double_t y2=pv->GetY(); + Double_t distXdiam=TMath::Abs(x2-GetDiamondX()); + Double_t distYdiam=TMath::Abs(y2-GetDiamondY()); + Double_t cov1[6],cov2[6]; + mainV->GetCovarianceMatrix(cov1); + pv->GetCovarianceMatrix(cov2); + Double_t errxDist=TMath::Sqrt(cov2[0]+GetSigma2DiamondX()); + Double_t erryDist=TMath::Sqrt(cov2[2]+GetSigma2DiamondY()); + Double_t errzDist=TMath::Sqrt(cov1[5]+cov2[5]); + Double_t cutXdiam=nSigmaDiamXY*errxDist; + if(GetSigma2DiamondX()<0.0001)cutXdiam=99999.; //protection for missing diamond information + Double_t cutYdiam=nSigmaDiamXY*erryDist; + if(GetSigma2DiamondY()<0.0001)cutYdiam=99999.; //protection for missing diamond information + if( (distXdiamnSigmaZdist*errzDist) ){ + return kTRUE; + } + } + } + } + return kFALSE; +} //______________________________________________________________________________ void AliAODEvent::Print(Option_t *) const diff --git a/STEER/AliAODEvent.h b/STEER/AliAODEvent.h index 3b99fb15715..f54cdd38e76 100644 --- a/STEER/AliAODEvent.h +++ b/STEER/AliAODEvent.h @@ -80,7 +80,7 @@ class AliAODEvent : public AliVEvent { void SetMagneticField(Double_t mf){if (fHeader) fHeader->SetMagneticField(mf);} void SetMuonMagFieldScale(Double_t mf){if (fHeader) fHeader->SetMuonMagFieldScale(mf);} void SetDiamond(Float_t xy[2],Float_t cov[3]){if (fHeader) fHeader->SetDiamond(xy,cov);} - + void SetDiamondZ(Float_t z, Float_t sig2z){if (fHeader) fHeader->SetDiamondZ(z,sig2z);} Int_t GetRunNumber() const {return fHeader ? fHeader->GetRunNumber() : -999;} UInt_t GetPeriodNumber() const {return fHeader ? fHeader->GetPeriodNumber() : 0;} UInt_t GetOrbitNumber() const {return fHeader ? fHeader->GetOrbitNumber() : 0;} @@ -89,7 +89,11 @@ class AliAODEvent : public AliVEvent { Double_t GetMuonMagFieldScale() const {return fHeader ? fHeader->GetMuonMagFieldScale() : -999.;} Double_t GetDiamondX() const {return fHeader ? fHeader->GetDiamondX() : -999.;} Double_t GetDiamondY() const {return fHeader ? fHeader->GetDiamondY() : -999.;} + Double_t GetDiamondZ() const {return fHeader ? fHeader->GetDiamondZ() : -999.;} void GetDiamondCovXY(Float_t cov[3]) const {cov[0]=-999.; if(fHeader) fHeader->GetDiamondCovXY(cov);} + Double_t GetSigma2DiamondX() const {return fHeader ? fHeader->GetSigma2DiamondX() : -999.;} + Double_t GetSigma2DiamondY() const {return fHeader ? fHeader->GetSigma2DiamondY() : -999.;} + Double_t GetSigma2DiamondZ() const {return fHeader ? fHeader->GetSigma2DiamondZ() : -999.;} void SetEventType(UInt_t eventType){fHeader->SetEventType(eventType);} void SetTriggerMask(ULong64_t n) {fHeader->SetTriggerMask(n);} @@ -124,6 +128,16 @@ class AliAODEvent : public AliVEvent { // primary vertex virtual AliAODVertex *GetPrimaryVertex() const { return GetVertex(0); } + virtual AliAODVertex *GetPrimaryVertexSPD() const; + + // -- Pileup vertices + Int_t GetNumberOfPileupVerticesTracks() const; + Int_t GetNumberOfPileupVerticesSPD() const; + virtual AliAODVertex *GetPileupVertexSPD(Int_t iV=0) const; + virtual AliAODVertex *GetPileupVertexTracks(Int_t iV=0) const; + virtual Bool_t IsPileupFromSPD(Int_t minContributors=3, Double_t minZdist=0.8, Double_t nSigmaZdist=3., Double_t nSigmaDiamXY=2., Double_t nSigmaDiamZ=5.) const; + + // V0 TClonesArray *GetV0s() const { return fV0s; } diff --git a/STEER/AliESDEvent.cxx b/STEER/AliESDEvent.cxx index 4f8f55e05bf..5ef88937475 100644 --- a/STEER/AliESDEvent.cxx +++ b/STEER/AliESDEvent.cxx @@ -1603,64 +1603,54 @@ Bool_t AliESDEvent::IsHLTTriggerFired(const char* name) const return kTRUE; } -Bool_t AliESDEvent::IsPileupFromSPD(Int_t ncont, Double_t distz, Double_t nSigmaDeltaZ, Double_t nSigmaXY, Int_t option) const { +//______________________________________________________________________________ +Bool_t AliESDEvent::IsPileupFromSPD(Int_t minContributors, + Double_t minZdist, + Double_t nSigmaZdist, + Double_t nSigmaDiamXY, + Double_t nSigmaDiamZ) const{ // // This function checks if there was a pile up // reconstructed with SPD // - Double_t diamx= GetDiamondX(); - Double_t diamsigma2x= GetSigma2DiamondX(); - Double_t diamy= GetDiamondY(); - Double_t diamsigma2y= GetSigma2DiamondY(); - - Double_t sigmax= TMath::Sqrt(diamsigma2x); - Double_t sigmay= TMath::Sqrt(diamsigma2y); - - Double_t z1=fSPDVertex->GetZ(); Int_t nc1=fSPDVertex->GetNContributors(); if(nc1<1) return kFALSE; Int_t nPileVert=GetNumberOfPileupVerticesSPD(); if(nPileVert==0) return kFALSE; + for(Int_t i=0; iGetZ(); - Double_t x2=pv->GetX(); - Double_t y2=pv->GetY(); Int_t nc2=pv->GetNContributors(); - Double_t distanceZ=TMath::Abs(z2-z1); - Double_t distanceX=TMath::Abs(x2-diamx); - Double_t distanceY=TMath::Abs(y2-diamy); - Double_t errzDist=0.; - Double_t errxDist=0.; - Double_t erryDist=0.; - if(option==0){ - Double_t ez1=fSPDVertex->GetZRes(); - Double_t ez2=pv->GetZRes(); - errzDist=TMath::Sqrt(ez1*ez1+ez2*ez2); - }else{ - Double_t resol1=-75.6+834.6/TMath::Sqrt((Double_t)nc1); - resol1/=10000.; - Double_t resol2=-75.6+834.6/TMath::Sqrt((Double_t)nc2); - resol2/=10000.; - errzDist=TMath::Sqrt(resol1*resol1+resol2*resol2); + if(nc2>=minContributors){ + Double_t z1=fSPDVertex->GetZ(); + Double_t z2=pv->GetZ(); + Double_t distZ=TMath::Abs(z2-z1); + Double_t distZdiam=TMath::Abs(z2-GetDiamondZ()); + Double_t cutZdiam=nSigmaDiamZ*GetSigma2DiamondZ(); + if(GetSigma2DiamondZ()<0.0001)cutZdiam=99999.; //protection for missing z diamond information + if(distZ>minZdist && distZdiamGetX(); + Double_t y2=pv->GetY(); + Double_t distXdiam=TMath::Abs(x2-GetDiamondX()); + Double_t distYdiam=TMath::Abs(y2-GetDiamondY()); + Double_t cov1[6],cov2[6]; + fSPDVertex->GetCovarianceMatrix(cov1); + pv->GetCovarianceMatrix(cov2); + Double_t errxDist=TMath::Sqrt(cov2[0]+GetSigma2DiamondX()); + Double_t erryDist=TMath::Sqrt(cov2[2]+GetSigma2DiamondY()); + Double_t errzDist=TMath::Sqrt(cov1[5]+cov2[5]); + Double_t cutXdiam=nSigmaDiamXY*errxDist; + if(GetSigma2DiamondX()<0.0001)cutXdiam=99999.; //protection for missing diamond information + Double_t cutYdiam=nSigmaDiamXY*erryDist; + if(GetSigma2DiamondY()<0.0001)cutYdiam=99999.; //protection for missing diamond information + if( (distXdiamnSigmaZdist*errzDist) ){ + return kTRUE; + } + } } - - - Double_t ex2 = pv->GetXRes(); - Double_t ey2= pv->GetYRes(); - errxDist=TMath::Sqrt(ex2*ex2+sigmax*sigmax); - erryDist=TMath::Sqrt(ey2*ey2+sigmay*sigmay); - - if(nc2>=ncont && distanceZ>nSigmaDeltaZ*errzDist && distanceXdistz) - - return kTRUE; - } - - - - return kFALSE; } - + return kFALSE; +} //______________________________________________________________________________ void AliESDEvent::EstimateMultiplicity(Int_t &tracklets, Int_t &trITSTPC, Int_t &trITSSApure, Double_t eta, Bool_t useDCAFlag,Bool_t useV0Flag) const diff --git a/STEER/AliESDEvent.h b/STEER/AliESDEvent.h index a2122708cda..abda314e1bb 100644 --- a/STEER/AliESDEvent.h +++ b/STEER/AliESDEvent.h @@ -254,13 +254,17 @@ public: return (const AliESDVertex *)(fSPDPileupVertices?fSPDPileupVertices->UncheckedAt(i):0x0); } Char_t AddPileupVertexSPD(const AliESDVertex *vtx); - Bool_t IsPileupFromSPD(Int_t ncont=2, Double_t distz=1., Double_t nSigmaDeltaZ=3., Double_t nSigmaXY=2., Int_t option=0) const; - const AliESDVertex *GetPileupVertexTracks(Int_t i) const { return (const AliESDVertex *)(fTrkPileupVertices?fTrkPileupVertices->UncheckedAt(i):0x0); } Char_t AddPileupVertexTracks(const AliESDVertex *vtx); + virtual Bool_t IsPileupFromSPD(Int_t minContributors=3, + Double_t minZdist=0.8, + Double_t nSigmaZdist=3., + Double_t nSigmaDiamXY=2., + Double_t nSigmaDiamZ=5.) const; + AliESDtrack *GetTrack(Int_t i) const { return (AliESDtrack *)(fTracks?fTracks->UncheckedAt(i):0x0); } diff --git a/STEER/AliVEvent.h b/STEER/AliVEvent.h index 0e79e70bf2f..a1707ba9d4b 100644 --- a/STEER/AliVEvent.h +++ b/STEER/AliVEvent.h @@ -91,7 +91,14 @@ public: // Primary vertex virtual const AliVVertex *GetPrimaryVertex() const {return 0x0;} - + virtual Bool_t IsPileupFromSPD(Int_t /*minContributors*/, + Double_t /*minZdist*/, + Double_t /*nSigmaZdist*/, + Double_t /*nSigmaDiamXY*/, + Double_t /*nSigmaDiamZ*/) + const{ + return kFALSE; + } //---------- end of new stuff -- 2.43.0