diff --git a/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/HelixFit/HelixFitJava.java b/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/HelixFit/HelixFitJava.java index 0e0ac9226e..a82914b892 100644 --- a/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/HelixFit/HelixFitJava.java +++ b/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/HelixFit/HelixFitJava.java @@ -997,7 +997,81 @@ public HelixFitObject HelixFit(int PointNum, double szPos[][], int fit_track_to_ return h; } -// Public weighted wrapper (similar to the original HelixFit) + // Compute the path length along the helix described by h, + // for the portion where the radial distance to the beamline + // (sqrt(x^2 + y^2)) is between rMin and rMax [mm]. + private double computePathFromHelix(HelixFitObject h, double rMin, double rMax) { + + // Helix parameters from the fit + double A = h.get_A(); // circle center X in XY-plane + double B = h.get_B(); // circle center Y in XY-plane + double X0 = h.get_X0(); // point of closest approach in XY + double Y0 = h.get_Y0(); + double thDeg = h.get_Theta(); // already stored in degrees in the wrapper + double theta = Math.toRadians(thDeg); // convert to radians for sin() + + // Circle radius in XY-plane + double R = Math.sqrt((X0 - A)*(X0 - A) + (Y0 - B)*(Y0 - B)); + if (R < 1e-6) { + // Degenerate circle + return 0.0; + } + + // Distance from origin to circle center + double C = Math.sqrt(A*A + B*B); + + // Minimal / maximal distance from the circle to the origin + double Dmin = Math.abs(C - R); + double Dmax = C + R; + + // If there is no overlap between [rMin, rMax] and [Dmin, Dmax], path = 0 + if (rMax <= Dmin || rMin >= Dmax) { + return 0.0; + } + + // Clip the requested radial window to what the circle actually covers + double rIn = Math.max(rMin, Dmin); + double rOut = Math.min(rMax, Dmax); + if (rOut <= rIn) { + return 0.0; + } + + // Helper constants for r^2(psi) = C0 + 2 R C cos(psi) + double C0 = A*A + B*B + R*R; + + // Helper to compute psi(r) in [0,π] from r: + // cos(psi) = (r^2 - C0)/(2 R C) + java.util.function.DoubleFunction psiOfR = (double r) -> { + double num = r*r - C0; + double den = 2.0 * R * C; + if (Math.abs(den) < 1e-12) { + return 0.0; + } + double u = num / den; + // Clamp to [-1,1] to avoid NaN from acos + if (u > 1.0) u = 1.0; + if (u < -1.0) u = -1.0; + return Math.acos(u); // in [0, π] + }; + + double psiIn = psiOfR.apply(rIn); + double psiOut = psiOfR.apply(rOut); + + double dPsi = Math.abs(psiOut - psiIn); + + // ds/dpsi = R / sin(theta) + double sinTheta = Math.sin(theta); + if (Math.abs(sinTheta) < 1e-6) { + // Track almost parallel to the beam; this formula is ill-defined. + // Return 0 for safety or some large sentinel. + return 0.0; + } + + double path = (R / sinTheta) * dPsi; + + return path; + } + // Public weighted wrapper (similar to the original HelixFit) public HelixFitObject HelixFitWeighted(int PointNum, double[][] szPos, double[] weights, int fit_track_to_beamline) { @@ -1188,14 +1262,17 @@ public HelixFitObject helix_fit_with_doca_selection(List docaCluste } // weights = null => all weights = 1.0 inside helix_fit_weighted + // Use the DOCa-based bestChi2 as helix chi2 double PI=Math.acos(0.0)*2; //double Rho=0,Phi=0,Theta=0,X0=0,Y0=0,DCA=0,Chi2=0; double Phi_deg; double Theta_deg; - HelixFitObject h = helix_fit_weighted(nPoints, szPosSel, null, fit_track_to_beamline);; - + HelixFitObject h = helix_fit_weighted(nPoints, szPosSel, null, fit_track_to_beamline); + h.set_Chi2(bestChi2); + double path = computePathFromHelix(h, 30.0, 70.0); + h.set_path(path); Phi_deg=Math.toDegrees(h.get_Phi()); if(Phi_deg >= 180){ Phi_deg -= 360; @@ -1211,6 +1288,4 @@ public HelixFitObject helix_fit_with_doca_selection(List docaCluste } -} - - +} \ No newline at end of file diff --git a/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/HelixFit/HelixFitObject.java b/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/HelixFit/HelixFitObject.java index 3e8591a13b..cee8ad91d2 100644 --- a/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/HelixFit/HelixFitObject.java +++ b/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/HelixFit/HelixFitObject.java @@ -16,7 +16,8 @@ public class HelixFitObject { private double _DCA; private double _Chi2; private double _magfield; - + private double _path; + public HelixFitObject(){ //default constructor } @@ -33,6 +34,7 @@ public HelixFitObject(double Rho, double A, double B, double Phi, double Theta, _DCA = DCA; _Chi2 = Chi2; _magfield = 50; + _path = 0; } public double get_Rho(){ return _Rho; @@ -82,6 +84,9 @@ public double get_DCA(){ public double get_Chi2(){ return _Chi2; } + public void set_Chi2(double Chi2){ + _Chi2 = Chi2; + } public double get_Mom(){ return 0.3*_magfield*Math.abs(_Rho)/(10*Math.sin(Math.toRadians(_Theta))); } @@ -103,5 +108,10 @@ public double get_dEdx(){ public void set_magfield(double magfield){ _magfield = magfield; } - + public double get_path(){ + return _path; + } + public void set_path(double path){ + _path = path; + } } \ No newline at end of file diff --git a/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Track/Track.java b/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Track/Track.java index e9944b7846..7efc0cc37e 100644 --- a/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Track/Track.java +++ b/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Track/Track.java @@ -73,6 +73,8 @@ public void setPositionAndMomentum(HelixFitObject helixFitObject) { this.px0 = helixFitObject.get_px(); this.py0 = helixFitObject.get_py(); this.pz0 = helixFitObject.get_pz(); + this.chi2 = helixFitObject.get_Chi2(); + this.path = helixFitObject.get_path(); } public void setPositionAndMomentumVec(double[] x) { @@ -176,15 +178,40 @@ public void set_trackId(int _trackId) { public void set_chi2(double _chi2) { chi2 = _chi2;} public void set_sum_residuals(double _sum_residuals) { sum_residuals = _sum_residuals;} public int get_trackId() {return trackId;} - public int get_n_hits() {return n_hits;} - public int get_sum_adc() {return sum_adc;} + public int get_n_hits() { + if (hits == null) { + return 0; + } + return hits.size(); + } + public int get_sum_adc() { + if (hits == null || hits.isEmpty()) { + return 0; + } + int sum = 0; + for (Hit h : hits) { + sum += (int) Math.round(h.getADC()); + } + return sum; + } public double get_chi2() {return chi2;} public double get_sum_residuals() {return sum_residuals;} // AHDC::track public void set_dEdx(double _dEdx) { dEdx = _dEdx;} public void set_p_drift(double _p_drift) { p_drift = _p_drift;} public void set_path(double _path) { path = _path;} - public double get_dEdx() {return dEdx;} + public double get_dEdx() { + if (path <= 0) { + dEdx = 0; + }else { + int sum = 0; + for (Hit h : hits) { + sum += (int) Math.round(h.getADC()); + } + dEdx = sum/path; + } + return dEdx; + } public double get_p_drift() {return p_drift;} public double get_path() {return path;} @@ -196,4 +223,4 @@ public void set_trackId(int _trackId) { public int get_predicted_ATOF_layer() {return predicted_ATOF_layer;} public int get_predicted_ATOF_wedge() {return predicted_ATOF_wedge;} -} +} \ No newline at end of file diff --git a/reconstruction/alert/src/main/java/org/jlab/service/ahdc/AHDCEngine.java b/reconstruction/alert/src/main/java/org/jlab/service/ahdc/AHDCEngine.java index b6c703bfc1..eb3dbdca24 100644 --- a/reconstruction/alert/src/main/java/org/jlab/service/ahdc/AHDCEngine.java +++ b/reconstruction/alert/src/main/java/org/jlab/service/ahdc/AHDCEngine.java @@ -218,13 +218,14 @@ else if (modeTrackFinding == ModeTrackFinding.CV_Hough) { // V) Global fit int trackid = 0; ArrayList all_docaClusters = new ArrayList<>(); + AHDC_Tracks.removeIf(track -> track.get_Clusters().size() < 3); for (Track track : AHDC_Tracks) { trackid++; track.set_trackId(trackid); List originalClusters = track.get_Clusters(); ArrayList docaClusters = DocaClusterRefiner.buildRefinedClusters(originalClusters); all_docaClusters.addAll(docaClusters); - if (docaClusters == null || docaClusters.size() < 3) { + if (docaClusters == null || docaClusters.size() < 3 || originalClusters == null || originalClusters.size() < 3) { // not enough points, skip helix fit continue; } @@ -314,4 +315,4 @@ public static void main(String[] args) { System.out.println("finished " + (System.nanoTime() - starttime) * Math.pow(10, -9)); } -} +} \ No newline at end of file