Skip to content

Commit

Permalink
added start/end pos/dir info from kalman output
Browse files Browse the repository at this point in the history
  • Loading branch information
@LiamOS
LiamOS committed Jul 12, 2024
1 parent 64dd7f8 commit 6d947dc
Show file tree
Hide file tree
Showing 5 changed files with 74 additions and 27 deletions.
22 changes: 22 additions & 0 deletions src/TMS_Kalman.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -118,6 +118,17 @@ void TMS_Kalman::RunKalman() {
}

SetMomentum(1./KalmanNodes.back().CurrentState.qp);

//std::cout << "filter start pos : " << KalmanNodes.back().CurrentState.x << ", " << KalmanNodes.back().CurrentState.y << ", " << KalmanNodes.back().CurrentState.z << std::endl;
SetStartPosition(KalmanNodes.back().CurrentState.x, KalmanNodes.back().CurrentState.y, KalmanNodes.back().CurrentState.z);
//std::cout << "filter end pos : " << KalmanNodes.at(1).CurrentState.x << ", " << KalmanNodes.at(1).CurrentState.y << ", " << KalmanNodes.at(1).CurrentState.z << std::endl;
SetEndPosition(KalmanNodes.at(1).CurrentState.x, KalmanNodes.at(1).CurrentState.y, KalmanNodes.at(1).CurrentState.z);

//std::cout << "filter start dir: " << KalmanNodes.back().CurrentState.dxdz << ", " << KalmanNodes.back().CurrentState.dydz << std::endl;
SetStartDirection(KalmanNodes.back().CurrentState.dxdz, KalmanNodes.back().CurrentState.dydz);
//std::cout << "filter end dir: " << KalmanNodes.at(1).CurrentState.dxdz << ", " << KalmanNodes.at(1).CurrentState.dydz << std::endl;
SetEndDirection(KalmanNodes.at(1).CurrentState.dxdz, KalmanNodes.at(1).CurrentState.dydz);

if (std::isnan(momentum) || std::isinf(momentum)){
std::cerr << "[TMS_Kalmann.cpp] Weirdness -- Momentum from fitter isn't a sane number: " << momentum << std::endl;
}
Expand Down Expand Up @@ -328,6 +339,17 @@ void TMS_Kalman::Predict(TMS_KalmanNode &Node) {
CurrentState.dxdz = FilteredVec[2];
CurrentState.dydz = FilteredVec[3];

//'CovarianceMatrix.Print();
//GainMatrix.Print();
//CurrentState.Print();
if ( (CurrentState.x < TMS_Const::TMS_Start[0]) || (CurrentState.x > TMS_Const::TMS_End[0]) ) // point outside x region
{
std::cerr << "lol x value outside TMS: " << CurrentState.y << "\tTMS: [" << TMS_Const::TMS_Start[0] << ", "<< TMS_Const::TMS_End[0] << "]" << std::endl;
}
if ( (CurrentState.y < TMS_Const::TMS_Start[1]) || (CurrentState.y > TMS_Const::TMS_End[1]) ) // point outside y region
{
std::cerr << "lol y value outside TMS: " << CurrentState.y << "\tTMS: [" << TMS_Const::TMS_Start[1] << ", "<< TMS_Const::TMS_End[1] << "]" << std::endl;
}

Node.SetRecoXY(CurrentState);
if (Talk) Node.PrintTrueReco();
Expand Down
12 changes: 12 additions & 0 deletions src/TMS_Kalman.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -205,9 +205,21 @@ class TMS_Kalman {
TMS_Kalman();
TMS_Kalman(std::vector<TMS_Hit> &Candidates);

double Start[3];
double End[3];
double StartDirection[3];
double EndDirection[3];

double GetKEEstimateFromLength(double startx, double endx, double startz, double endz);

void SetMomentum(double mom) {momentum = mom;}
// Set direction unit vectors from only x and y slope
void SetStartDirection(double ax, double ay) {StartDirection[0]=ax; StartDirection[1]=ay; StartDirection[2]=sqrt(1 - ax*ax - ay*ay);};
void SetEndDirection (double ax, double ay) {EndDirection[0]=ax; EndDirection[1]=ay; EndDirection[2]=sqrt(1 - ax*ax - ay*ay);};

// Set position unit vectors
void SetStartPosition(double ax, double ay, double az) {Start[0]=ax; Start[1]=ay; Start[2]=az;};
void SetEndPosition (double ax, double ay, double az) {End[0]=ax; End[1]=ay; End[2]=az;};

double GetMomentum() {return momentum;}

Expand Down
17 changes: 13 additions & 4 deletions src/TMS_Track.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,8 @@ class TMS_Track {
int Charge;
double Start[3]; // Start point in x,y,z
double End[3]; // End point in x,y,z
double Direction[3]; // Unit vector in track direction
double StartDirection[3]; // Unit vector in track direction at start
double EndDirection[3]; // Unit vector in track direction at end
double Length;
double Occupancy;
double EnergyDeposit;
Expand All @@ -34,9 +35,17 @@ class TMS_Track {
TMS_TrueParticle GetTrueParticle() {return fTrueParticle;};

// Manually set variables
void SetEnergyDeposit(double val) {EnergyDeposit = val;};
void SetEnergyRange (double val) {EnergyRange = val;};
void SetMomentum (double val) {Momentum = val;};
void SetEnergyDeposit (double val) {EnergyDeposit = val;};
void SetEnergyRange (double val) {EnergyRange = val;};
void SetMomentum (double val) {Momentum = val;};

// Set direction unit vectors from only x and y slope
void SetStartDirection(double ax, double ay, double az) {StartDirection[0]=ax; StartDirection[1]=ay; StartDirection[2]=az;};
void SetEndDirection (double ax, double ay, double az) {EndDirection[0]=ax; EndDirection[1]=ay; EndDirection[2]=az;};

// Set position unit vectors
void SetStartPosition(double ax, double ay, double az) {Start[0]=ax; Start[1]=ay; Start[2]=az;};
void SetEndPosition (double ax, double ay, double az) {End[0]=ax; End[1]=ay; End[2]=az;};

int nHits;
std::vector<TMS_Hit> Hits;
Expand Down
47 changes: 25 additions & 22 deletions src/TMS_TreeWriter.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -227,19 +227,20 @@ void TMS_TreeWriter::MakeBranches() {
Reco_Tree->Branch("SliceNo", &SliceNo, "SliceNo/I");
Reco_Tree->Branch("SpillNo", &SpillNo, "SpillNo/I");

Reco_Tree->Branch("nTracks", &nTracks, "nTracks/I");
Reco_Tree->Branch("nHits", nHitsIn3DTrack, "nHits[nTracks]/I");
Reco_Tree->Branch("TrackHitPos", RecoTrackHitPos, "TrackHitPos[nTracks][200][3]/F");
Reco_Tree->Branch("nKalmanNodes", nKalmanNodes, "nKalmanNodes[nTracks]/I");
Reco_Tree->Branch("KalmanPos", RecoTrackKalmanPos, "TrackHitPos[nTracks][200][3]/F");
Reco_Tree->Branch("KalmanTruePos",RecoTrackKalmanTruePos, "TrackHitTruePos[nTracks][200][3]/F");
Reco_Tree->Branch("StartPos", RecoTrackStartPos, "StartPos[nTracks][3]/F");
Reco_Tree->Branch("Direction", RecoTrackDirection, "Direction[nTracks][3]/F");
Reco_Tree->Branch("EndPos", RecoTrackEndPos, "EndPos[nTracks][3]/F");
Reco_Tree->Branch("EnergyRange", RecoTrackEnergyRange, "EnergyRange[nTracks]/F");
Reco_Tree->Branch("EnergyDeposit",RecoTrackEnergyDeposit, "EnergyDeposit[nTracks]/F");
Reco_Tree->Branch("Momentum", RecoTrackMomentum, "Momentum[nTracks]/F");
Reco_Tree->Branch("Length", RecoTrackLength, "Length[nTracks]/F");
Reco_Tree->Branch("nTracks", &nTracks, "nTracks/I");
Reco_Tree->Branch("nHits", nHitsIn3DTrack, "nHits[nTracks]/I");
Reco_Tree->Branch("TrackHitPos", RecoTrackHitPos, "TrackHitPos[nTracks][200][3]/F");
Reco_Tree->Branch("nKalmanNodes", nKalmanNodes, "nKalmanNodes[nTracks]/I");
Reco_Tree->Branch("KalmanPos", RecoTrackKalmanPos, "TrackHitPos[nTracks][200][3]/F");
Reco_Tree->Branch("KalmanTruePos", RecoTrackKalmanTruePos, "TrackHitTruePos[nTracks][200][3]/F");
Reco_Tree->Branch("StartPos", RecoTrackStartPos, "StartPos[nTracks][3]/F");
Reco_Tree->Branch("StartDirection",RecoTrackStartDirection,"StartDirection[nTracks][3]/F");
Reco_Tree->Branch("EndDirection", RecoTrackEndDirection, "StartDirection[nTracks][3]/F");
Reco_Tree->Branch("EndPos", RecoTrackEndPos, "EndPos[nTracks][3]/F");
Reco_Tree->Branch("EnergyRange", RecoTrackEnergyRange, "EnergyRange[nTracks]/F");
Reco_Tree->Branch("EnergyDeposit", RecoTrackEnergyDeposit, "EnergyDeposit[nTracks]/F");
Reco_Tree->Branch("Momentum", RecoTrackMomentum, "Momentum[nTracks]/F");
Reco_Tree->Branch("Length", RecoTrackLength, "Length[nTracks]/F");


MakeTruthBranches(Truth_Info);
Expand Down Expand Up @@ -1214,7 +1215,8 @@ void TMS_TreeWriter::Fill(TMS_Event &event) {
for (int j = 0; j < 3; j++) {
RecoTrackStartPos[itTrack][j] = RecoTrack->Start[j];
RecoTrackEndPos[itTrack][j] = RecoTrack->End[j];
RecoTrackDirection[itTrack][j] = RecoTrack->Direction[j];
RecoTrackStartDirection[itTrack][j] = RecoTrack->StartDirection[j];
RecoTrackEndDirection[itTrack][j] = RecoTrack->EndDirection[j];
}

for (unsigned int j = 0; j < RecoTrack->KalmanNodes.size(); ++j) {
Expand Down Expand Up @@ -1248,13 +1250,13 @@ void TMS_TreeWriter::Fill(TMS_Event &event) {
// std::cout << "TreeWriter hit position: " << RecoTrackHitPos[itTrack][j][0] << " " << RecoTrackHitPos[itTrack][j][1] << " " << RecoTrackHitPos[itTrack][j][2] << std::endl;
}
// Can manually compute direction if it hasn't been set
if ( (RecoTrackDirection[itTrack][0] == 0) && (RecoTrackDirection[itTrack][1] == 0) && (RecoTrackDirection[itTrack][2] == 0) )
{ // If true it seems the direction hasn't been set
for (int j = 0; j < 3; j++)
{ // Right now no need to make sure this is a unit vector
RecoTrackDirection[itTrack][j] = RecoTrack->End[j] - RecoTrack->Start[j];
}
}
// if ( (RecoTrackDirection[itTrack][0] == 0) && (RecoTrackDirection[itTrack][1] == 0) && (RecoTrackDirection[itTrack][2] == 0) )
// { // If true it seems the direction hasn't been set
// for (int j = 0; j < 3; j++)
// { // Right now no need to make sure this is a unit vector
// RecoTrackDirection[itTrack][j] = RecoTrack->End[j] - RecoTrack->Start[j];
// }
// }

// Now fill truth info
if (itTrack >= __TMS_MAX_LINES__) {
Expand Down Expand Up @@ -1666,8 +1668,9 @@ void TMS_TreeWriter::Clear() {
for (int i = 0; i < __TMS_MAX_TRACKS__; ++i) {
for (int j = 0; j < 3; ++j) {
RecoTrackStartPos[i][j] = DEFAULT_CLEARING_FLOAT;
RecoTrackDirection[i][j] = DEFAULT_CLEARING_FLOAT;
RecoTrackStartDirection[i][j] = DEFAULT_CLEARING_FLOAT;
RecoTrackEndPos[i][j] = DEFAULT_CLEARING_FLOAT;
RecoTrackEndDirection[i][j] = DEFAULT_CLEARING_FLOAT;
}
for (int k = 0; k < __TMS_MAX_LINE_HITS__; ++k) {
RecoTrackHitPos[i][k][0] = DEFAULT_CLEARING_FLOAT;
Expand Down
3 changes: 2 additions & 1 deletion src/TMS_TreeWriter.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -50,8 +50,9 @@ class TMS_TreeWriter {
float RecoTrackKalmanTruePos[__TMS_MAX_TRACKS__][__TMS_MAX_LINE_HITS__][3];
float RecoTrackHitPos[__TMS_MAX_TRACKS__][__TMS_MAX_LINE_HITS__][3]; // Due to a lack of variables, but as this is taken from line hits, it would make sense (maybe times 2?)
float RecoTrackStartPos[__TMS_MAX_TRACKS__][3];
float RecoTrackDirection[__TMS_MAX_TRACKS__][3];
float RecoTrackStartDirection[__TMS_MAX_TRACKS__][3];
float RecoTrackEndPos[__TMS_MAX_TRACKS__][3];
float RecoTrackEndDirection[__TMS_MAX_TRACKS__][3];
float RecoTrackEnergyRange[__TMS_MAX_TRACKS__];
float RecoTrackEnergyDeposit[__TMS_MAX_TRACKS__];
float RecoTrackMomentum[__TMS_MAX_TRACKS__];
Expand Down

0 comments on commit 6d947dc

Please sign in to comment.