Difference between revisions of "Brief Kinematic Fitting Example"

Revision as of 18:22, 15 March 2012

Brief summary

Kinematic fitting in the HallD sim-recon code is provided by the DKinFit class.

DKinFit kfit;

In the simplest mode, the fitter takes a fully specified reaction (4-momenta and covariance matrices for all initial and final state particles) and fits using the constraint of overall conservation of 4-momentum. Missing final state particles can be specified, and additional mass constraints can be applied. The results of the fit can be used to filter out background that does not satisfy the fit hypothesis and to improve the 4-momenta and covariance matrices.

The kinematic fitter takes as its input two vectors of DKinematicData objects, containing the initial and final state particles, specified by SetInitial() and SetFinal(). Each DKinematicData contains both a 4-momentum and a kinematic covariance matrix.

vector<DKinematicData> initial;
vector<DKinematicData> final;

DKinematicData beam_photon,target,reconstructed_proton,reconstructed_pi1,...;
//fill these objects appropriately using objects (DChargedTrackHypothesis, DNeutralParticleHypothesis, DBeamPhoton) from factories
//...

initial.push_back(beam_photon);
initial.push_back(target);

final.push_back(reconstructed_proton);
final.push_back(reconstructed_pi1);
final.push_back(...);

kfit.SetFinal(final);
kfit.SetInitial(initial);

If there is a missing particle in the final state, use the SetMissingParticle()

const double neutron_mass=.940;
kfit.SetMissingParticle(neutron_mass);

If there is a mass constraint (e.g. two photons from a pi0 decay), use SetMassConstraint()

vector<DKinematicData> final(3);
final[0] = reconstructed_proton;
final[1] = reconstructed_photon1;
final[2] = reconstructed_photon2;

vector<int> constrained_particles_indices(2);
//1 and 2 are the indices of the photons in final
constrained_particles_indices[0]=1;
constrained_particles_indices[1]=2;

const double pi0_mass=0.135;
kfit.SetMassConstraint(pi0_mass,constrained_particles_indices);

kfit.SetFinal(final);

After specifying initial and final state particles, missing particles, and mass constraints, simply use the Fit() function

kfit.Fit();

Useful functions for accessing the results of the fit:

• Prob(): returns the confidence level of fit, based on
• GetPull(int n): returns the n'th pull. There is one pull for each measured quantity, so three pulls for each measured particle (p_x, p_y, and p_z). Pulls 0-2 are the pulls of the first initial particle (0 for p_x, 1 for p_y, 2 for p_z), pulls 3-5 are the pulls of the second initial particle, pulls 6-8 are the pulls of the first final particle, etc.
• vector<DKinematicData> GetInitial_out(), vector<DKinematicData> GetFinal_out(): Returns DKinematicData with post-fit "corrected" 4-momentum and covariance matrices.
• GetMissingError(): returns the error on the missing particle momentum. Unfortunately, it uses the (p,lambda,phi) parameterization at the moment for some reason, so it's somewhat useless.

b1pi example

vector<DKinematicData> initial;
vector<DKinematicData> final;
vector<DKinematicData> initial_post;
vector<DKinematicData> final_post;

//Add the beam particle to the initial state
DKinematicData beam_photon;
DVector3 mom(0.0, 0.0, 9.0);
beam_photon.setMomentum(mom);
beam_photon.setPosition(pos);
beam_photon.setMass(0.0);
beam_photon.setCharge(0.0);
beam_photon.setMassFixed();
beam_photon.clearErrorMatrix();

initial.push_back(beam_photon);

//Add the target to the initial state
DKinematicData kd_targ;
kd_targ.setMass(0.93827);
kd_targ.setCharge(1.0);
kd_targ.setMassFixed();
kd_targ.setMomentum( DVector3(0.0, 0.0, 0.0) );
kd_targ.clearErrorMatrix();
initial.push_back(kd_targ);

//Add all our particles to the final state
final.push_back(proton[k]);
final.push_back(pip[l]);
final.push_back(pip[m]);
final.push_back(pim[n]);
final.push_back(pim[o]);

double cl=kfit.Prob();