Line.h

Go to the documentation of this file.

#pragma once
 
#include <iostream>
#include <utility>
#include <cmath>
 
#include <TF3.h>
 
#include "LinearAlgebra.h"
 
 
namespace Line
{
    //____________________________________________________________________________
    struct FitLine2D
    {
        double m;              // slope
        double c;              // y-intercept
        double err_m;          // uncertainty in slope
        double err_c;          // uncertainty in y-intercept
        double fit_chi2ndf;    // goodness of linear fit
        int trkID;             // AMS trkID
        int g3PID;             // Geant3 PID of track
 
        FitLine2D(double m, double c, double err_m, double err_c, double fit_chi2ndf, int trkID, int g3PID)
            : m(m), c(c), err_m(err_m), err_c(err_c), fit_chi2ndf(fit_chi2ndf), trkID(trkID), g3PID(g3PID) {}
    };
    
    
    //____________________________________________________________________________
    struct FitLine3DMinuit
    {
        int trkID;
        int g3PID;
 
        bool fitSuccess;  // Flag to indicate if the fit was successful
        std::vector<double> fitParameters;  // Fit parameters
        std::vector<double> fitErrors;      // Errors on fit parameters
 
        double chi2;      // Chi-square of the fit
        int ndf;          // Degrees of freedom
 
        FitLine3DMinuit() : trkID(-99999999), g3PID(-99999999), fitSuccess(false), chi2(0.0), ndf(0) {}  // default constructor
    };
    
    
    //__________________________________________________________________________
    struct FitLine3D
    {
        TFitResultPtr fit;  // Store TFitResultPtr instead of TF3
        int trkID;
        int g3PID;
 
        FitLine3D() : trkID(-99999999), g3PID(-99999999) {}  // default constructor
    };
    
    
    //___________________________________________________________________________
    struct Line3DParametricForm
    {
        int trkID;
        int g3PID;
        std::vector<double> point;      
        std::vector<double> direction;  
 
        Line3DParametricForm() : trkID(-99999999), g3PID(-99999999), point(), direction() {} 
 
        Line3DParametricForm(int trkID, int g3PID, const std::vector<double>& pt, const std::vector<double>& dir)
            : trkID(trkID), g3PID(g3PID), point(pt), direction(dir) {} 
    };
    
    
    //___________________________________________________________________________
    Line3DParametricForm Getline3DParametricForm(const double& m1, const double& b1,  
                                                 const double& m2, const double& b2,  
                                                 const double& m3, const double& b3,  
                                                 const int& trkID,
                                                 const int& g3PID,
                                                 bool debugMode) 
    {
        
        // const std::vector<double> vX0Y0Z0 = LA::SolveSystemOfThreeEquations(m1, b1, m2, b2, m3, b3, debugMode);
        
        const std::vector<double> vX0Y0Z0 = LA::ROOTEquationSolver(m1, b1, m2, b2, m3, b3, debugMode);         
        
        if (vX0Y0Z0.size() != 3)
            return Line3DParametricForm(); 
 
        double x0 = vX0Y0Z0.at(0);
        double y0 = vX0Y0Z0.at(1);
        double z0 = vX0Y0Z0.at(2);
 
        double a = 1;  // x does not change if parallel to x-axis
        double b = m3; // y does not change if parallel to y-axis
        double c = m1; // Use m2 if line is parallel to x or y-axis
 
        if (debugMode)
        {
            std::cout << "\nThe parametric form of the 3D line is:" << std::endl;
            std::cout << "x = " << x0 << " + " << a << ".t" << std::endl;
            std::cout << "y = " << y0 << " + " << b << ".t" << std::endl;
            std::cout << "z = " << z0 << " + " << c << ".t" << std::endl;
        }
        
        const std::vector<double> vDir = {a, b, c};
        return Line3DParametricForm(trkID, g3PID, vX0Y0Z0, vDir);
    }
    
    
    
    
    //__________________________________________________________________________
    double distanceFromLine3D(double *coords, double *params)
    {
        //______________________________________
        
        double x = coords[0];
        double y = coords[1];
        double z = coords[2];
 
        double px = params[0];  
        double py = params[1];  
        double pz = params[2];  
        double vx = params[3];  
        double vy = params[4];  
        double vz = params[5];  
        
        //______________________________________
        
        double u = ((x - px) * vx + (y - py) * vy + (z - pz) * vz) / (vx * vx + vy * vy + vz * vz);
        
        
        double dx = px + u * vx - x;
        double dy = py + u * vy - y;
        double dz = pz + u * vz - z;
        
        
        return std::sqrt(dx*dx + dy*dy + dz*dz);
    }
}