• comparison TCFIT with KFparticle : here
• summary : here
• plot : red = 'S' parameter from kfparticle (s =l/p) , black = slength_tcfit/pt (D0)
• pdf file to summarize : here
• update :

April, 23th :

a bug was found when plotting the decay length for real data (the same behavior appears for simulation with fixed vertex) : this pdf

The solution was to set the pVertex (used for KFParticle to (PrimVertexX[l], PrimVertexY[l], PrimVertexZ[l]) , instead of (0,0, PrimVertexZ[l]) as originally)

And there is no move of StDcaGeometry to get the DCA of tracks to PV.

The correction is explained here

April, 27th :

I used the KFParticle code over K0short (since the decay length : ~2cm) is much larger than for D0/D+

The summary is in this webpage

Another issue appears : on slide 7, one can see that the invariant mass seems broader ; it is confirmed by this plot

that shows :

1. in red   : invariant mass of K0 calculated with the fitted momenta from KFParticle
2. in blue : invariant mass of K0 calculated with the global momenta from GlobalTracks branch

We see that the inv. mass with the global momenta has a lower resolution

sigma_mass (global) ~ 5.9 MeV
sigma_mass (fitted) ~ 9 MeV

As explained in this presentation by S. Gorbunov, that gave an example of how-to-apply KFParticle, I may have missed to set the daughter particl
to the production vertex
--> currently investigating

Another possibility is to look at the error of the invariant mass and the error of the decay length : maybe some cuts has to be applied.

Looking at KFParticleBase.cxx::GetMass() :

Int_t KFParticleBase::GetMass( Double_t &m, Double_t &error ) const
{
//* Calculate particle mass

// s = sigma^2 of m2/2

Double_t s = ( fP[3]*fP[3]*fC[9] + fP[4]*fP[4]*fC[14] + fP[5]*fP[5]*fC[20]
+ fP[6]*fP[6]*fC[27]
+2*( + fP[3]*fP[4]*fC[13] + fP[5]*(fP[3]*fC[18] + fP[4]*fC[19])
- fP[6]*( fP[3]*fC[24] + fP[4]*fC[25] + fP[5]*fC[26] ) )
);
Double_t m2 = TMath::Abs(fP[6]*fP[6] - fP[3]*fP[3] - fP[4]*fP[4] - fP[5]*fP[5]);
m = TMath::Sqrt(m2);
if( m>1.e-10 ){
if( s>=0 ){
error = TMath::Sqrt(s)/m;
return 0;
}
}
error = 1.e20;
return 1;
}

It looks like the error on the mass is error ~ sqrt(decaylength) / mass;

Reminder : KFParticleBase.h :

double GetX () const { return fP[0]; }
double GetY () const { return fP[1]; }
double GetZ () const { return fP[2]; }
double GetPx () const { return fP[3]; }
double GetPy () const { return fP[4]; }
double GetPz () const { return fP[5]; }
double GetE () const { return fP[6]; }
double GetS () const { return fP[7]; }
int GetQ () const { return fQ; }
double GetChi2 () const { return fChi2; }
int GetNDF () const { return fNDF; }

double GetParameter ( int i ) const { return fP[i]; }
double GetCovariance( int i ) const { return fC[i]; }
double GetCovariance( int i, int j ) const { return fC[IJ(i,j)]; }

 2 interesting plots are :

• error of decaylength vs error of mass
• mass vs its error

From the first plot, it seems to have a correlation btw both, ie error on the mass of the order of ~0.003 is obtained for low error on dcay length.

This is directly related to the second plot : dm/m, where we see that the expected K0 mass (0.49) (the plateau on the plot) has an error of ~0.003

April, 28 th

Test  on real data : day 110

Test 1 : set the productionVertex to both daughters via :

          particle[1].SetProductionVertex(DP);
          particle[2].SetProductionVertex(DP);

plot is here

Test 2 : did the same + substract daughters from PV :

        //this is a test
        particle[1].SubtractFromVertex(particle[0]);
        particle[2].SubtractFromVertex(particle[0]);
        particle[0].AddDaughter(DP);
       
        DP.SetProductionVertex(pVertex); // it is already set in KFParticleBase.cxx::Construct()
       
        particle[1].SetProductionVertex(DP);
        particle[2].SetProductionVertex(DP);

plot is here

for reference, w/o any of these changes, the plots is here

reference : 21k events

test1 : 12k events

test2 : 15k events

 

Note how we use KFParticle :

In MuKpi, we use KFParticle by calling Construct() method :

        KFParticle DP;
        DP.Construct(vDaughters,NDaughters,&pVertex,-1,0);

which refers to :

void KFParticleBase::Construct( const KFParticleBase* vDaughters[], Int_t NDaughters,
                   const KFParticleBase *Parent,  Double_t Mass, Bool_t IsConstrained         )
{
  //* Full reconstruction in one go

  Int_t maxIter = 1;
  bool wasLinearized = fIsLinearized;
  if( !fIsLinearized || IsConstrained ){
    //fVtxGuess[0] = fVtxGuess[1] = fVtxGuess[2] = 0;  //!!!!
    fVtxGuess[0] = GetX();
    fVtxGuess[1] = GetY();
    fVtxGuess[2] = GetZ();
    fIsLinearized = 1;
    maxIter = 3;
  }

  Double_t constraintC[6];

  if( IsConstrained ){
    for(Int_t i=0;i<6;++i) constraintC[i]=fC[i];
  } else {
    for(Int_t i=0;i<6;++i) constraintC[i]=0.;
    constraintC[0] = constraintC[2] = constraintC[5] = 100.;   
  }

  for( Int_t iter=0; iter<maxIter; iter++ ){
    fAtProductionVertex = 0;
    fSFromDecay = 0;
    fP[0] = fVtxGuess[0];
    fP[1] = fVtxGuess[1];
    fP[2] = fVtxGuess[2];
    fP[3] = 0;
    fP[4] = 0;
    fP[5] = 0;
    fP[6] = 0;
    fP[7] = 0;
  
    for(Int_t i=0;i<6; ++i) fC[i]=constraintC[i];
    for(Int_t i=6;i<36;++i) fC[i]=0.;
    fC[35] = 1.;
   
    fNDF  = IsConstrained ?0 :-3;
    fChi2 =  0.;
    fQ = 0;

    for( Int_t itr =0; itr<NDaughters; itr++ ){
      AddDaughter( *vDaughters[itr] );   
    }
    if( iter<maxIter-1){
      for( Int_t i=0; i<3; i++ ) fVtxGuess[i] = fP[i]; 
    }
  }
  fIsLinearized = wasLinearized;   

  if( Mass>=0 ) SetMassConstraint( Mass );
  if( Parent ) SetProductionVertex( *Parent );
}
 

--> So we do not constrained the mass of the mother particle

 May, 27th

The issue is still the resolution of K0's worse with KFParticle than getting the mass using TLorentzVector with global momenta.

The mass (and its error) are calculated in KFParticle with the following method :

Int_t KFParticleBase::GetMass( Double_t &m, Double_t &error ) const
{
  //* Calculate particle mass
 
  // s = sigma^2 of m2/2

  Double_t s = (  fP[3]*fP[3]*fC[9] + fP[4]*fP[4]*fC[14] + fP[5]*fP[5]*fC[20]
          + fP[6]*fP[6]*fC[27]
        +2*( + fP[3]*fP[4]*fC[13] + fP[5]*(fP[3]*fC[18] + fP[4]*fC[19])
             - fP[6]*( fP[3]*fC[24] + fP[4]*fC[25] + fP[5]*fC[26] )   )
         );
  Double_t m2 = TMath::Abs(fP[6]*fP[6] - fP[3]*fP[3] - fP[4]*fP[4] - fP[5]*fP[5]);
  m  = TMath::Sqrt(m2);
  if( m>1.e-10 ){
    if( s>=0 ){
      error = TMath::Sqrt(s)/m;
      return 0;
    }
  }
  error = 1.e20;
  return 1;
}

, where fP[i] is  : double fP[8];  //* Main particle parameters {X,Y,Z,Px,Py,Pz,E,S[=DecayLength/P]}

--> so the formulation of the mass is correct

issue : fP[3,4,5] are the updated momentum of the mother particle ; in principle, at the end of the kalman Fitter, these parameters should be the best estimation

 

debug:

mass TCFIT = 1.88495         // mass from TCFIT
 before refit of mother
 D0.pT = 2.76089 D0.P = 2.92729  D0.Mass = 1.88454
 after refit of mother
 D0.pT = 2.76089 D0.P = 2.92729  D0.Mass = 1.88454  // we see that SetProductionVertex for the mother has no effect because it's done in KFParticleBase.cxx
 before refit
 particle[1].px = -2.38415 particle[1].py = 1.27415 particle[1].pz = 0.851399
 particle[2].px = -0.270076 particle[2].py = -0.509804 particle[2].pz = 0.122535
 mass = 1.88496
 after refit
 particle[1].px = -2.38236 particle[1].py = 1.27383 particle[1].pz = 0.850883
 particle[2].px = -0.270567 particle[2].py = -0.509868 particle[2].pz = 0.121962
 mass = 1.88452

it seems that the mass before SetProductionVertex() for both daughters is in better agreement with mass from KFParticle

code tested :

                KFParticle DP;
        gBenchmark->Start("KFParticle");
        DP.Construct(vDaughters,NDaughters,&pVertex,-1,0);
        //DP.Construct(vDaughters,NDaughters,particle[0],-1,0); this is not working      
        gBenchmark->Stop("KFParticle");
        particle[1].SubtractFromVertex(particle[0]);
        particle[2].SubtractFromVertex(particle[0]);
        particle[0].AddDaughter(DP);
        cout << " before refit of mother " << endl;
        cout << " D0.pT = " <<DP.GetPt() <<" D0.P = " << DP.GetMomentum() <<"  D0.Mass = " << DP.GetMass() << endl;
        DP.SetProductionVertex(pVertex); // it is already set in KFParticleBase.cxx::Construct()
        cout << " after refit of mother " << endl;
        cout << " D0.pT = " <<DP.GetPt() <<" D0.P = " << DP.GetMomentum() <<"  D0.Mass = " << DP.GetMass() << endl;
        cout <<" before refit " << endl;
        cout <<" particle[1].px = " << particle[1].GetPx() << " particle[1].py = " << particle[1].GetPy()<<" particle[1].pz = " << particle[1].GetPz() <<endl;
        cout <<" particle[2].px = " << particle[2].GetPx() << " particle[2].py = " << particle[2].GetPy()<<" particle[2].pz = " << particle[2].GetPz() <<endl;
        TLorentzVector tempo1;
        TLorentzVector tempo2;
        TLorentzVector tempo11;
        TLorentzVector tempo22;
        TVector3 kfp1,kfp2,akfp1,akfp2;
        kfp1 = TVector3(particle[1].GetPx(),particle[1].GetPy(),particle[1].GetPz());
        kfp2 = TVector3(particle[2].GetPx(),particle[2].GetPy(),particle[2].GetPz());
        tempo1.SetVectMag(kfp1,amK);
        tempo2.SetVectMag(kfp2,amPi);
        tempo1+=tempo2;
        cout <<" mass = " << tempo1.M() << endl;
        particle[1].SetProductionVertex(DP);
        particle[2].SetProductionVertex(DP);
        cout <<" after refit" << endl;
        cout <<" particle[1].px = " << particle[1].GetPx() << " particle[1].py = " << particle[1].GetPy()<<" particle[1].pz = " << particle[1].GetPz() <<endl;
        cout <<" particle[2].px = " << particle[2].GetPx() << " particle[2].py = " << particle[2].GetPy()<<" particle[2].pz = " << particle[2].GetPz() <<endl;
        akfp1 = TVector3(particle[1].GetPx(),particle[1].GetPy(),particle[1].GetPz());
        akfp2 = TVector3(particle[2].GetPx(),particle[2].GetPy(),particle[2].GetPz());
        tempo11.SetVectMag(akfp1,amK);
        tempo22.SetVectMag(akfp2,amPi);
        tempo11+=tempo22;
        cout <<" mass = " << tempo11.M() << endl;

****possible solution****

                DP.Construct(vDaughters,NDaughters,&pVertex,-1,0);
        //DP.Construct(vDaughters,NDaughters,particle[0],-1,0); this is not working      
        gBenchmark->Stop("KFParticle");
        particle[1].SubtractFromVertex(particle[0]);
        particle[2].SubtractFromVertex(particle[0]);
        particle[0].AddDaughter(DP);

1. Substract daughters from PV
2. add mother to PV
3. do not SetProductionVertex for daughters