// The objects pointed by the following pointer are not owned
// by this class and are not deleted
fCurrentVertex = 0;
- if(fTrksToSkip) { delete [] fTrksToSkip; fTrksToSkip=NULL; }
- if(fIdSel) { delete [] fIdSel; fIdSel=NULL; }
+ if(fTrksToSkip) { delete fTrksToSkip; fTrksToSkip=NULL; }
+ if(fIdSel) { delete fIdSel; fIdSel=NULL; }
}
//----------------------------------------------------------------------------
AliESDVertex* AliVertexerTracks::FindPrimaryVertex(AliVEvent *vEvent)
FindPrimaryVertex(&trkArrayOrig,idOrig);
if(fMode==0) trkArrayOrig.Delete();
- delete [] idOrig; idOrig=NULL;
+ delete idOrig; idOrig=NULL;
if(f) {
f->Close(); delete f; f = NULL;
// fill fTrkArraySel, for VertexFinder()
fIdSel = new UShort_t[nTrksOrig];
PrepareTracks(*trkArrayOrig,idOrig,0);
- if(fIdSel) { delete [] fIdSel; fIdSel=NULL; }
+ if(fIdSel) { delete fIdSel; fIdSel=NULL; }
Double_t cutsave = fDCAcut;
fDCAcut = fDCAcutIter0;
// vertex finder
// between initVertex and fCurrentVertex)
for(Int_t iter=1; iter<=2; iter++) {
if(fOnlyFitter && iter==1) continue;
- if(fIdSel) { delete [] fIdSel; fIdSel=NULL; }
+ if(fIdSel) { delete fIdSel; fIdSel=NULL; }
fIdSel = new UShort_t[nTrksOrig];
Int_t nTrksSel = PrepareTracks(*trkArrayOrig,idOrig,iter);
AliDebug(1,Form("N tracks selected in iteration %d: %d",iter,nTrksSel));
indices[jj] = fIdSel[jj];
fCurrentVertex->SetIndices(nIndices,indices);
}
- delete [] indices; indices=NULL;
+ if (indices) {delete indices; indices=NULL;}
//
// set vertex title
AliDebug(1,Form("xyz: %f %f %f; nc %d",fCurrentVertex->GetXv(),fCurrentVertex->GetYv(),fCurrentVertex->GetZv(),fCurrentVertex->GetNContributors()));
// clean up
- delete [] fIdSel; fIdSel=NULL;
+ delete fIdSel; fIdSel=NULL;
fTrkArraySel.Delete();
- if(fTrksToSkip) { delete [] fTrksToSkip; fTrksToSkip=NULL; }
+ if(fTrksToSkip) { delete fTrksToSkip; fTrksToSkip=NULL; }
//
return fCurrentVertex;
delete track; continue;
}
- Bool_t propagateOK = kFALSE;
+ Bool_t propagateOK = kFALSE, cutond0z0 = kTRUE;
// propagate track to vertex
if(optImpParCut<2 || fOnlyFitter) { // optImpParCut==1 or 0
propagateOK = track->PropagateToDCA(initVertex,GetFieldkG(),100.,d0z0,covd0z0);
fCurrentVertex->GetSigmaXYZ(sigmaCurr);
normdistx = TMath::Abs(fCurrentVertex->GetXv()-fNominalPos[0])/TMath::Sqrt(sigmaCurr[0]*sigmaCurr[0]+fNominalCov[0]);
normdisty = TMath::Abs(fCurrentVertex->GetYv()-fNominalPos[1])/TMath::Sqrt(sigmaCurr[1]*sigmaCurr[1]+fNominalCov[2]);
+ AliDebug(1,Form("normdistx %f %f %f",fCurrentVertex->GetXv(),fNominalPos[0],TMath::Sqrt(sigmaCurr[0]*sigmaCurr[0]+fNominalCov[0])));
+ AliDebug(1,Form("normdisty %f %f %f",fCurrentVertex->GetYv(),fNominalPos[1],TMath::Sqrt(sigmaCurr[1]*sigmaCurr[1]+fNominalCov[2])));
+ AliDebug(1,Form("sigmaCurr %f %f %f",sigmaCurr[0],sigmaCurr[1],TMath::Sqrt(fNominalCov[0])+TMath::Sqrt(fNominalCov[2])));
if(normdistx < 3. && normdisty < 3. &&
(sigmaCurr[0]+sigmaCurr[1])<(TMath::Sqrt(fNominalCov[0])+TMath::Sqrt(fNominalCov[2]))) {
propagateOK = track->PropagateToDCA(fCurrentVertex,GetFieldkG(),100.,d0z0,covd0z0);
} else {
propagateOK = track->PropagateToDCA(initVertex,GetFieldkG(),100.,d0z0,covd0z0);
+ if(fConstraint) cutond0z0=kFALSE;
}
}
maxd0rphi = fNSigma*sigmad0;
if(optImpParCut==1) maxd0rphi *= 5.;
maxd0rphi = TMath::Min(maxd0rphi,fFiducialR);
- //sigmad0z0 = TMath::Sqrt(covd0z0[0]+covd0z0[2]); // for future improvement
- //maxd0z0 = 10.*fNSigma*sigmad0z0;
+ //sigmad0z0 = TMath::Sqrt(covd0z0[0]+covd0z0[2]);
AliDebug(1,Form("trk %d; id %d; |d0| = %f; d0 cut = %f; |z0| = %f; |d0|oplus|z0| = %f; d0z0 cut = %f",i,(Int_t)idOrig[i],TMath::Abs(d0z0[0]),maxd0rphi,TMath::Abs(d0z0[1]),TMath::Sqrt(d0z0[0]*d0z0[0]+d0z0[1]*d0z0[1]),fMaxd0z0));
// if fConstraint=kTRUE, during iteration 2,
// select tracks with d0oplusz0 < fMaxd0z0
if((!fConstraint && optImpParCut>0 && fVert.GetNContributors()>0) ||
- ( fConstraint && optImpParCut==2)) {
+ ( fConstraint && optImpParCut==2 && cutond0z0)) {
if(nTrksOrig>=3 &&
TMath::Sqrt(d0z0[0]*d0z0[0]+d0z0[1]*d0z0[1])>fMaxd0z0) {
AliDebug(1," rejected");
AliESDVertex* AliVertexerTracks::RemoveTracksFromVertex(AliESDVertex *inVtx,
TObjArray *trkArray,
UShort_t *id,
- Float_t *diamondxy)
+ Float_t *diamondxy) const
{
//
// Removes tracks in trksTree from fit of inVtx
sumWi -= wWi;
sumWiri -= wWiri;
- // track chi2
+ // track contribution to chi2
TMatrixD deltar = rv; deltar -= ri;
TMatrixD wWideltar(wWi,TMatrixD::kMult,deltar);
Double_t chi2i = deltar(0,0)*wWideltar(0,0)+
cov[5] = vVnew(2,2);
// store data in the vertex object
- AliESDVertex *outVtx = new AliESDVertex(position,cov,chi2,nUsedTrks);
+ AliESDVertex *outVtx = new AliESDVertex(position,cov,chi2,nUsedTrks+1); // the +1 is for the constraint
outVtx->SetTitle(inVtx->GetTitle());
UShort_t *inindices = inVtx->GetIndices();
Int_t nIndices = nUsedTrks;
}
}
outVtx->SetIndices(nIndices,outindices);
- delete [] outindices;
+ if (outindices) delete outindices;
/*
printf("Vertex before removing tracks:");
return outVtx;
}
//---------------------------------------------------------------------------
+AliESDVertex* AliVertexerTracks::RemoveConstraintFromVertex(AliESDVertex *inVtx,
+ Float_t *diamondxyz,
+ Float_t *diamondcov) const
+{
+//
+// Removes diamond constraint from fit of inVtx
+//
+
+ if(!strstr(inVtx->GetTitle(),"VertexerTracksWithConstraint")) {
+ printf("ERROR: primary vertex has no constraint: cannot remove it\n");
+ return 0x0;
+ }
+ if(inVtx->GetNContributors()<3) {
+ printf("ERROR: primary vertex has less than 2 tracks: cannot remove contraint\n");
+ return 0x0;
+ }
+
+ // diamond constraint
+ TMatrixD vVb(3,3);
+ vVb(0,0) = diamondcov[0];
+ vVb(0,1) = diamondcov[1];
+ vVb(0,2) = 0.;
+ vVb(1,0) = diamondcov[1];
+ vVb(1,1) = diamondcov[2];
+ vVb(1,2) = 0.;
+ vVb(2,0) = 0.;
+ vVb(2,1) = 0.;
+ vVb(2,2) = diamondcov[5];
+ TMatrixD vVbInv(TMatrixD::kInverted,vVb);
+ TMatrixD rb(3,1);
+ rb(0,0) = diamondxyz[0];
+ rb(1,0) = diamondxyz[1];
+ rb(2,0) = diamondxyz[2];
+ TMatrixD vVbInvrb(vVbInv,TMatrixD::kMult,rb);
+
+ // input vertex
+ TMatrixD rv(3,1);
+ rv(0,0) = inVtx->GetXv();
+ rv(1,0) = inVtx->GetYv();
+ rv(2,0) = inVtx->GetZv();
+ TMatrixD vV(3,3);
+ Double_t cov[6];
+ inVtx->GetCovMatrix(cov);
+ vV(0,0) = cov[0];
+ vV(0,1) = cov[1]; vV(1,0) = cov[1];
+ vV(1,1) = cov[2];
+ vV(0,2) = cov[3]; vV(2,0) = cov[3];
+ vV(1,2) = cov[4]; vV(2,1) = cov[4];
+ vV(2,2) = cov[5];
+ TMatrixD vVInv(TMatrixD::kInverted,vV);
+ TMatrixD vVInvrv(vVInv,TMatrixD::kMult,rv);
+
+
+ TMatrixD sumWi = vVInv - vVbInv;
+
+
+ TMatrixD sumWiri = vVInvrv - vVbInvrb;
+
+ TMatrixD rvnew(3,1);
+ TMatrixD vVnew(3,3);
+
+ // new inverted of weights matrix
+ TMatrixD invsumWi(TMatrixD::kInverted,sumWi);
+ vVnew = invsumWi;
+ // new position of primary vertex
+ rvnew.Mult(vVnew,sumWiri);
+
+ Double_t position[3];
+ position[0] = rvnew(0,0);
+ position[1] = rvnew(1,0);
+ position[2] = rvnew(2,0);
+ cov[0] = vVnew(0,0);
+ cov[1] = vVnew(0,1);
+ cov[2] = vVnew(1,1);
+ cov[3] = vVnew(0,2);
+ cov[4] = vVnew(1,2);
+ cov[5] = vVnew(2,2);
+
+
+ Double_t chi2 = inVtx->GetChi2();
+
+ // diamond constribution to chi2
+ TMatrixD deltar = rv; deltar -= rb;
+ TMatrixD vVbInvdeltar(vVbInv,TMatrixD::kMult,deltar);
+ Double_t chi2b = deltar(0,0)*vVbInvdeltar(0,0)+
+ deltar(1,0)*vVbInvdeltar(1,0)+
+ deltar(2,0)*vVbInvdeltar(2,0);
+ // remove from total chi2
+ chi2 -= chi2b;
+
+ // store data in the vertex object
+ AliESDVertex *outVtx = new AliESDVertex(position,cov,chi2,inVtx->GetNContributors()-1);
+ outVtx->SetTitle("VertexerTracksNoConstraint");
+ UShort_t *inindices = inVtx->GetIndices();
+ Int_t nIndices = inVtx->GetNIndices();
+ outVtx->SetIndices(nIndices,inindices);
+
+ return outVtx;
+}
+//---------------------------------------------------------------------------
void AliVertexerTracks::SetCuts(Double_t *cuts)
{
//
}
if(!fTrkArraySel.IsEmpty()) fTrkArraySel.Delete();
- if(fIdSel) {delete [] fIdSel; fIdSel=NULL;}
- if(fTrksToSkip) {delete [] fTrksToSkip; fTrksToSkip=NULL;}
+ if(fIdSel) {delete fIdSel; fIdSel=NULL;}
+ if(fTrksToSkip) {delete fTrksToSkip; fTrksToSkip=NULL;}
return;
}
AliESDVertex* AliVertexerTracks::VertexForSelectedTracks(TObjArray *trkArray,
UShort_t *id,
Bool_t optUseFitter,
- Bool_t optPropagate)
+ Bool_t optPropagate,
+ Bool_t optUseDiamondConstraint)
{
//
// Return vertex from tracks (AliExternalTrackParam) in array
//
fCurrentVertex = 0;
- SetConstraintOff();
+ // set optUseDiamondConstraint=TRUE only if you are reconstructing the
+ // primary vertex!
+ if(optUseDiamondConstraint) {
+ SetConstraintOn();
+ } else {
+ SetConstraintOff();
+ }
// get tracks and propagate them to initial vertex position
fIdSel = new UShort_t[(Int_t)trkArray->GetEntriesFast()];
Int_t nTrksSel = PrepareTracks(*trkArray,id,0);
- if(nTrksSel < TMath::Max(2,fMinTracks)) {
+ if((!optUseDiamondConstraint && nTrksSel<TMath::Max(2,fMinTracks)) ||
+ (optUseDiamondConstraint && nTrksSel<1)) {
TooFewTracks();
return fCurrentVertex;
}
// vertex finder
- switch (fAlgo) {
- case 1: StrLinVertexFinderMinDist(1); break;
- case 2: StrLinVertexFinderMinDist(0); break;
- case 3: HelixVertexFinder(); break;
- case 4: VertexFinder(1); break;
- case 5: VertexFinder(0); break;
- default: printf("Wrong algorithm\n"); break;
+ if(nTrksSel==1) {
+ AliDebug(1,"Just one track");
+ OneTrackVertFinder();
+ } else {
+ switch (fAlgo) {
+ case 1: StrLinVertexFinderMinDist(1); break;
+ case 2: StrLinVertexFinderMinDist(0); break;
+ case 3: HelixVertexFinder(); break;
+ case 4: VertexFinder(1); break;
+ case 5: VertexFinder(0); break;
+ default: printf("Wrong algorithm\n"); break;
+ }
}
AliDebug(1," Vertex finding completed\n");
Double_t d0z0[2],covd0z0[3];
AliExternalTrackParam *t = 0;
if(fCurrentVertex->GetNContributors()>0) {
- indices = new UShort_t[fCurrentVertex->GetNContributors()];
+ indices = new UShort_t[fTrkArraySel.GetEntriesFast()];
for(Int_t jj=0; jj<(Int_t)fTrkArraySel.GetEntriesFast(); jj++) {
indices[jj] = fIdSel[jj];
t = (AliExternalTrackParam*)fTrkArraySel.At(jj);
}
// clean up
- if (indices) {delete [] indices; indices=NULL;}
- delete [] fIdSel; fIdSel=NULL;
+ if (indices) {delete indices; indices=NULL;}
+ delete fIdSel; fIdSel=NULL;
fTrkArraySel.Delete();
return fCurrentVertex;
}
//----------------------------------------------------------------------------
AliESDVertex* AliVertexerTracks::VertexForSelectedESDTracks(TObjArray *trkArray,
- Bool_t optUseFitter,
- Bool_t optPropagate)
+ Bool_t optUseFitter,
+ Bool_t optPropagate,
+ Bool_t optUseDiamondConstraint)
+
{
//
// Return vertex from array of ESD tracks
id[i] = (UShort_t)esdt->GetID();
}
- VertexForSelectedTracks(trkArray,id,optUseFitter,optPropagate);
+ VertexForSelectedTracks(trkArray,id,optUseFitter,optPropagate,optUseDiamondConstraint);
- delete [] id; id=NULL;
+ delete id; id=NULL;
return fCurrentVertex;
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff