test/minimalFittingExample/main2.cc

0001 #include <ConstField.h>
0002 #include <Exception.h>
0003 #include <FieldManager.h>
0004 #include <KalmanFitterRefTrack.h>
0005 #include <StateOnPlane.h>
0006 #include <Track.h>
0007 #include <TrackPoint.h>
0008
0009 #include <MaterialEffects.h>
0010 #include <RKTrackRep.h>
0011 #include <TGeoMaterialInterface.h>
0012
0013 #include <EventDisplay.h>
0014
0015 #include <PlanarMeasurement.h>
0016
0017 #include <TEveManager.h>
0018 #include <TGeoManager.h>
0019 #include <TVector3.h>
0020 #include <vector>
0021
0022 #include "TDatabasePDG.h"
0023 #include <TMath.h>
0024
0025
0026
0027
0028 int main() {
0029
0030   // init geometry and mag. field
0031   new TGeoManager("Geometry", "Geane geometry");
0032   TGeoManager::Import("genfitGeom.root");
0033   genfit::MaterialEffects::getInstance()->init(new genfit::TGeoMaterialInterface());
0034   genfit::FieldManager::getInstance()->init(new genfit::ConstField(0. ,10., 0.)); // 1 T
0035
0036
0037   // init event display
0038   genfit::EventDisplay* display = genfit::EventDisplay::getInstance();
0039
0040
0041   // init fitter
0042   genfit::AbsKalmanFitter* fitter = new genfit::KalmanFitterRefTrack();
0043
0044
0045   // particle pdg code; pion hypothesis
0046   const int pdg = 211;
0047
0048   // start values for the fit, e.g. from pattern recognition
0049   TVector3 pos(0, 0, 0);
0050   TVector3 mom(0, 0, 3);
0051
0052
0053   // trackrep
0054   genfit::AbsTrackRep* rep = new genfit::RKTrackRep(pdg);
0055
0056   // create track
0057   genfit::Track fitTrack(rep, pos, mom);
0058
0059
0060   const int detId(0); // detector ID
0061   int planeId(0); // detector plane ID
0062   int hitId(0); // hit ID
0063
0064   double detectorResolution(0.001); // resolution of planar detectors
0065   TMatrixDSym hitCov(2);
0066   hitCov.UnitMatrix();
0067   hitCov *= detectorResolution*detectorResolution;
0068
0069
0070   // add some planar hits to track with coordinates I just made up
0071   TVectorD hitCoords(2);
0072   hitCoords[0] = 0;
0073   hitCoords[1] = 0;
0074   genfit::PlanarMeasurement* measurement = new genfit::PlanarMeasurement(hitCoords, hitCov, detId, ++hitId, nullptr);
0075   measurement->setPlane(genfit::SharedPlanePtr(new genfit::DetPlane(TVector3(0,0,0), TVector3(1,0,0), TVector3(0,1,0))), ++planeId);
0076   fitTrack.insertPoint(new genfit::TrackPoint(measurement, &fitTrack));
0077
0078   hitCoords[0] = -0.15;
0079   hitCoords[1] = 0;
0080   measurement = new genfit::PlanarMeasurement(hitCoords, hitCov, detId, ++hitId, nullptr);
0081   measurement->setPlane(genfit::SharedPlanePtr(new genfit::DetPlane(TVector3(0,0,10), TVector3(1,0,0), TVector3(0,1,0))), ++planeId);
0082   fitTrack.insertPoint(new genfit::TrackPoint(measurement, &fitTrack));
0083
0084   hitCoords[0] = -0.4;
0085   hitCoords[1] = 0;
0086   measurement = new genfit::PlanarMeasurement(hitCoords, hitCov, detId, ++hitId, nullptr);
0087   measurement->setPlane(genfit::SharedPlanePtr(new genfit::DetPlane(TVector3(0,0,20), TVector3(1,0,0), TVector3(0,1,0))), ++planeId);
0088   fitTrack.insertPoint(new genfit::TrackPoint(measurement, &fitTrack));
0089
0090
0091
0092   //check
0093   fitTrack.checkConsistency();
0094
0095   // do the fit
0096   fitter->processTrack(&fitTrack);
0097
0098   // print fit result
0099   fitTrack.getFittedState().Print();
0100
0101   //check
0102   fitTrack.checkConsistency();
0103
0104
0105   display->addEvent(&fitTrack);
0106
0107
0108   delete fitter;
0109
0110   // open event display
0111   display->open();
0112
0113 }
0114
0115