Skip to content

Dc v1fresh #1524

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Closed
raffaelladevita wants to merge 2 commits into from
Closed

Dc v1fresh #1524

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
3e22298
revert fixes to pathlength to go back to dc v1
raffaelladevita Oct 15, 2026
53fc94d
revert fixes to pathlength and changes to T2D to go back to dc v1
raffaelladevita Oct 16, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
14 changes: 4 additions & 10 deletions ...s/clas-tracking/src/main/java/org/jlab/clas/tracking/kalmanfilter/zReference/KFitter.java
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -718,8 +718,6 @@ private boolean filter(int k, boolean forward, double annealingFactor) {
}
}

// Since no vertex inforamtion, the starting point for path length is the final point at the last layer.
// After vertex information is obtained, transition for the starting point from the final point to vertex will be taken.
private boolean filter(int k, boolean forward) {
StateVec sVec = sv.transported(forward).get(k);
org.jlab.clas.tracking.kalmanfilter.AMeasVecs.MeasVec mVec = mv.measurements.get(k);
Expand Down Expand Up @@ -895,7 +893,7 @@ public void calcFinalChisq(int sector, boolean nofilter) {
sv.transport(sector, k, 0, sVec, mv, this.getSwimmer(), forward);

StateVec svc = sv.transported(forward).get(0);
path += (forward ? 1 : -1) * svc.deltaPath;
path += svc.deltaPath;
svc.setPathLength(path);

double V0 = mv.measurements.get(0).surface.unc[0];
Expand Down Expand Up @@ -941,7 +939,7 @@ public void calcFinalChisq(int sector, boolean nofilter) {

double h = mv.hDoca(line, mv.measurements.get(k1 + 1).surface.wireLine[0]);
svc = sv.transported(forward).get(k1 + 1);
path += (forward ? 1 : -1) * svc.deltaPath;
path += svc.deltaPath;
svc.setPathLength(path);
svc.setProjector(mv.measurements.get(k1 + 1).surface.wireLine[0].origin().x());
svc.setProjectorDoca(h);
Expand Down Expand Up @@ -994,7 +992,7 @@ private void calcFinalChisqDAF(int sector, boolean nofilter) {
sv.transport(sector, k, 0, sVec, mv, this.getSwimmer(), forward);

StateVec svc = sv.transported(forward).get(0);
path += (forward ? 1 : -1) * svc.deltaPath;
path += svc.deltaPath;
svc.setPathLength(path);

Point3D point = new Point3D(svc.x, svc.y, mv.measurements.get(0).surface.measPoint.z());
Expand Down Expand Up @@ -1068,7 +1066,7 @@ private void calcFinalChisqDAF(int sector, boolean nofilter) {
}

svc = sv.transported(forward).get(k1 + 1);
path += (forward ? 1 : -1) * svc.deltaPath;
path += svc.deltaPath;
svc.setPathLength(path);

point = new Point3D(sv.transported(forward).get(k1 + 1).x, sv.transported(forward).get(k1 + 1).y, mv.measurements.get(k1 + 1).surface.measPoint.z());
Expand Down Expand Up @@ -1140,10 +1138,6 @@ public Matrix propagateToVtx(int sector, double Zf) {
return sv.transport(sector, finalStateVec.k, Zf, finalStateVec, mv, this.getSwimmer());
}

public double getDeltaPathToVtx(int sector, double Zf) {
return sv.getDeltaPath(sector, finalStateVec.k, Zf, finalStateVec, mv, this.getSwimmer());
}

//Todo: apply the common funciton to replace current function above
@Override
public void runFitter(AStateVecs sv, AMeasVecs mv) {
Expand Down
10 changes: 2 additions & 8 deletions ...racking/src/main/java/org/jlab/clas/tracking/kalmanfilter/zReference/KFitterStraight.java
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -411,8 +411,6 @@ private void calcFinalChisq(int sector) {
calcFinalChisq(sector, false);
}

// Since no vertex inforamtion, the starting point for path length is the final point at the last layer.
// After vertex information is obtained, transition for the starting point from the final point to vertex will be taken.
private void calcFinalChisq(int sector, boolean nofilter) {
int k = svzLength - 1;
this.chi2 = 0;
Expand All @@ -437,7 +435,7 @@ private void calcFinalChisq(int sector, boolean nofilter) {
sv.transport(sector, k, 0, sVec, mv, this.getSwimmer(), forward);

StateVec svc = sv.transported(forward).get(0);
path += (forward ? 1 : -1) * svc.deltaPath;
path += svc.deltaPath;
svc.setPathLength(path);

double V0 = mv.measurements.get(0).surface.unc[0];
Expand Down Expand Up @@ -484,7 +482,7 @@ private void calcFinalChisq(int sector, boolean nofilter) {

double h = mv.hDoca(line, mv.measurements.get(k1 + 1).surface.wireLine[0]);
svc = sv.transported(forward).get(k1+1);
path += (forward ? 1 : -1) * svc.deltaPath;
path += svc.deltaPath;
svc.setPathLength(path);
svc.setProjector(mv.measurements.get(k1 + 1).surface.wireLine[0].origin().x());
svc.setProjectorDoca(h);
Expand Down Expand Up @@ -515,10 +513,6 @@ public Matrix propagateToVtx(int sector, double Zf) {
return sv.transport(sector, finalStateVec.k, Zf, finalStateVec, mv, this.getSwimmer());
}

public double getDeltaPathToVtx(int sector, double Zf) {
return sv.getDeltaPath(sector, finalStateVec.k, Zf, finalStateVec, mv, this.getSwimmer());
}

@Override
public void runFitter(AStateVecs sv, AMeasVecs mv) {
throw new UnsupportedOperationException("Not supported yet.");
Expand Down
98 changes: 0 additions & 98 deletions ...clas-tracking/src/main/java/org/jlab/clas/tracking/kalmanfilter/zReference/StateVecs.java
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -152,104 +152,6 @@ public Matrix transport(int sector, int i, double Zf, StateVec iVec, AMeasVecs m

}

/**
*
* @param sector
* @param i initial state vector index
* @param Zf
* @param iVec state vector at the initial index
* @param mv measurements
*/
public double getDeltaPath(int sector, int i, double Zf, StateVec iVec, AMeasVecs mv, Swim swimmer) { // s = signed step-size

double stepSize = 1.0;
StateVec fVec = new StateVec(0);
fVec.x = iVec.x;
fVec.y = iVec.y;
fVec.z = iVec.z;
fVec.tx = iVec.tx;
fVec.ty = iVec.ty;
fVec.Q = iVec.Q;
fVec.B = iVec.B;
Matrix5x5.copy(iVec.CM, fVec.CM);

double s = 0;
double zInit = mv.measurements.get(i).surface.measPoint.z();
double BatMeas = iVec.B;

double z = zInit;

while (Math.signum(Zf - zInit) * z < Math.signum(Zf - zInit) * Zf) {
z = fVec.z;
if (z == Zf) {
break;
}

double x = fVec.x;
double y = fVec.y;
double tx = fVec.tx;
double ty = fVec.ty;
double Q = fVec.Q;
double dPath = fVec.deltaPath;
Matrix cMat = new Matrix();
Matrix5x5.copy(fVec.CM, cMat);
s = Math.signum(Zf - zInit) * stepSize;

// LOGGER.log(Level.FINEST, " from "+(float)Z[i]+" to "+(float)Z[f]+" at "+(float)z+" By is "+bf[1]+" B is "+Math.sqrt(bf[0]*bf[0]+bf[1]*bf[1]+bf[2]*bf[2])/Bmax+" stepSize is "+s);
if (Math.signum(Zf - zInit) * (z + s) > Math.signum(Zf - zInit) * Zf) {
s = Math.signum(Zf - zInit) * Math.abs(Zf - z);
}

//rk.RK4transport(sector, Q, x, y, z, tx, ty, s, swimmer, cMat, fVec, dPath);
rk.RK4transport(sector, s, swimmer, cMat, fVec);

// Q process noise matrix estimate
double p = Math.abs(1. / iVec.Q);

double X0 = this.getX0(mv.measurements.get(i).surface, z, Z);
double t_ov_X0 = Math.abs(s) / X0;//path length in radiation length units = t/X0 [true path length/ X0] ; Ar radiation length = 14 cm

double beta = this.beta;
if (beta > 1.0 || beta <= 0) {
beta = 1.0;
}

double sctRMS = 0;

////// Todo: Modify multi-scattering or remove it; After update, some parameters, like iteration termintion chonditions, may need to be updated.
// Speed of light should be 1
// From one measurement site to another, F and Q should be calculated separaetely with multiple steps; and then C' = FTCF + Q
if (Math.abs(s) > 0) {
sctRMS = ((0.0136) / (beta * PhysicsConstants.speedOfLight() * p)) * Math.sqrt(t_ov_X0)
* (1 + 0.038 * Math.log(t_ov_X0));
}

double cov_txtx = (1 + tx * tx) * (1 + tx * tx + ty * ty) * sctRMS * sctRMS;
double cov_tyty = (1 + ty * ty) * (1 + tx * tx + ty * ty) * sctRMS * sctRMS;
double cov_txty = tx * ty * (1 + tx * tx + ty * ty) * sctRMS * sctRMS;

fMS.set(
0, 0, 0, 0, 0,
0, 0, 0, 0, 0,
0, 0, cov_txtx, cov_txty, 0,
0, 0, cov_txty, cov_tyty, 0,
0, 0, 0, 0, 0
);

Matrix5x5.copy(fVec.CM, copyMatrix);
Matrix5x5.add(copyMatrix, fMS, fVec.CM);

if (Math.abs(fVec.B - BatMeas) < 0.0001) {
stepSize *= 2;
}

BatMeas = fVec.B;
}

return fVec.deltaPath;

}

/**
*
* @param sector
Expand Down
26 changes: 13 additions & 13 deletions reconstruction/dc/src/main/java/org/jlab/rec/dc/Constants.java
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -96,7 +96,7 @@ public static Constants getInstance() {
public int SECTORSELECT = 0;
public int NSUPERLAYERTRACKING = 5;
private boolean USETSTART = true;
//private boolean USETIMETBETA = false;
private boolean USETIMETBETA = false;
private boolean CHECKBETA = false;
private int T2D = 1; // 1=polynomial, 0=exponential
private boolean USEDOUBLETS = false;
Expand All @@ -107,11 +107,11 @@ public static Constants getInstance() {
public static final String TT = "/daq/tt/dc";
public static final String DOCARES = "/calibration/dc/signal_generation/doca_resolution";
public static final String TIME2DIST = "/calibration/dc/time_to_distance/time2dist";
public static final String T2DPRESSURE = "/calibration/dc/v2/t2d_pressure";
public static final String T2DPRESSURE = "/calibration/dc/time_to_distance/t2d_pressure";
public static final String PRESSURE = "/hall/weather/pressure";
public static final String T2DPRESSUREREF = "/calibration/dc/v2/ref_pressure";
public static final String T0CORRECTION = "/calibration/dc/v2/t0";
public static final String TDCTCUTS = "/calibration/dc/v2/tdc_cuts";
public static final String T2DPRESSUREREF = "/calibration/dc/time_to_distance/ref_pressure";
public static final String T0CORRECTION = "/calibration/dc/time_corrections/T0Corrections";
public static final String TDCTCUTS = "/calibration/dc/time_corrections/tdctimingcuts";
public static final String WIRESTAT = "/calibration/dc/tracking/wire_status";
public static final String TIMEJITTER = "/calibration/dc/time_jitter";
public static final String BEAMPOS = "/geometry/beam/position";
Expand Down Expand Up @@ -298,13 +298,13 @@ public void setUSETSTART(boolean usetstart) {
USETSTART = usetstart;
}

// public boolean useUSETIMETBETA() {
// return USETIMETBETA;
// }
//
// public void setUSETIMETBETA(boolean usetimebeta) {
// USETIMETBETA = usetimebeta;
// }
public boolean useUSETIMETBETA() {
return USETIMETBETA;
}
public void setUSETIMETBETA(boolean usetimebeta) {
USETIMETBETA = usetimebeta;
}

public double getWIREDIST() {
return WIREDIST;
Expand Down Expand Up @@ -424,7 +424,7 @@ public void printConfig(String engine) {
LOGGER.log(Level.INFO, "["+engine+"] run with wire ministagger = " + MINISTAGGERSTATUS.getName());
LOGGER.log(Level.INFO, "["+engine+"] run with wire feedthroughs = " + FEEDTHROUGHSSTATUS.getName());
LOGGER.log(Level.INFO, "["+engine+"] run with wire distortions = " + ENDPLATESBOWING);
//LOGGER.log(Level.INFO, "["+engine+"] run with with time Beta correction (is false for doca Beta correction) = " + USETIMETBETA);
LOGGER.log(Level.INFO, "["+engine+"] run with with time Beta correction (is false for doca Beta correction) = " + USETIMETBETA);
LOGGER.log(Level.INFO, "["+engine+"] run with with Beta cut = " + CHECKBETA);
LOGGER.log(Level.INFO, "["+engine+"] run with time to distance function set to exponential/polynomial (0/1) = " + T2D);
LOGGER.log(Level.INFO, "["+engine+"] run with with hit doublets recovery = " + USEDOUBLETS);
Expand Down
49 changes: 43 additions & 6 deletions reconstruction/dc/src/main/java/org/jlab/rec/dc/banks/HitReader.java
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,7 @@
import org.jlab.rec.dc.hit.Hit;
import org.jlab.rec.dc.timetodistance.TimeToDistanceEstimator;

import org.jlab.clas.swimtools.Swimmer;
import org.jlab.detector.banks.RawBank.OrderType;
import org.jlab.detector.banks.RawDataBank;
import org.jlab.detector.calib.utils.ConstantsManager;
Expand Down Expand Up @@ -54,6 +55,8 @@ public class HitReader {
private List<FittedHit> _HBHits; //hit-based tracking hit information
private List<FittedHit> _TBHits; //time-based tracking hit information

private final double timeBuf = 25.0;

private static final Logger LOGGER = Logger.getLogger(HitReader.class.getName());

public HitReader(Banks names, DCGeant4Factory detector) {
Expand Down Expand Up @@ -204,15 +207,49 @@ public void fetch_DCHits(DataEvent event) {

double timeCutMin = 0;
double timeCutMax = 0;
double timeCutLC = 0;

int region = ((superlayer + 1) / 2);
long hash = IndexedTable.DEFAULT_GENERATOR.hashCode(0,region,wire);
timeCutMin = tdccuts.getIntValueByHash("min", hash);
timeCutMax = tdccuts.getIntValueByHash("max", hash);

switch (region) {
case 1:
timeCutMin = tdccuts.getIntValue("MinEdge", 0, region, 0);
timeCutMax = tdccuts.getIntValue("MaxEdge", 0, region, 0);
break;
case 2:
if (wire <= 56) {
timeCutLC = tdccuts.getIntValue("LinearCoeff", 0, region, 1);
timeCutMin = tdccuts.getIntValue("MinEdge", 0, region, 1);
timeCutMax = tdccuts.getIntValue("MaxEdge", 0, region, 1);
}
if (wire > 56) {
timeCutLC = tdccuts.getIntValue("LinearCoeff", 0, region, 56);
timeCutMin = tdccuts.getIntValue("MinEdge", 0, region, 56);
timeCutMax = tdccuts.getIntValue("MaxEdge", 0, region, 56);
}
break;
case 3:
timeCutMin = tdccuts.getIntValue("MinEdge", 0, region, 0);
timeCutMax = tdccuts.getIntValue("MaxEdge", 0, region, 0)+timeBuf;
break;
}
boolean passTimingCut = false;

if (tdc > timeCutMin && tdc < timeCutMax)
if (region == 1 && tdc > timeCutMin && tdc < timeCutMax)
passTimingCut = true;
if (region == 2) {
double Bscale = Swimmer.getTorScale() * Swimmer.getTorScale();
if (wire >= 56) {
if (tdc > timeCutMin &&
tdc < timeCutMax + timeCutLC * (double) (112 - wire / 56) * Bscale)
passTimingCut = true;
} else {
if (tdc > timeCutMin &&
tdc < timeCutMax + timeCutLC * (double) (56 - wire / 56) * Bscale)
passTimingCut = true;
}
}
if (region == 3 && tdc > timeCutMin && tdc < timeCutMax)
passTimingCut = true;

if (passTimingCut) { // cut on spurious hits
Expand Down Expand Up @@ -657,8 +694,8 @@ private double[] getT0(int sector, int superlayer,
int cable = this.getCableID1to6(layer, wire);
int slot = this.getSlotID1to7(wire);

double t0 = t0Table.getDoubleValue("t0correction", sector, superlayer, slot, cable);
double t0E = t0Table.getDoubleValue("t0error", sector, superlayer, slot, cable);
double t0 = t0Table.getDoubleValue("T0Correction", sector, superlayer, slot, cable);
double t0E = t0Table.getDoubleValue("T0Error", sector, superlayer, slot, cable);

T0Corr[0] = t0;
T0Corr[1] = t0E;
Expand Down
43 changes: 21 additions & 22 deletions reconstruction/dc/src/main/java/org/jlab/rec/dc/hit/FittedHit.java
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -383,28 +383,27 @@ public void set_TimeToDistance(DataEvent event, double trkAngle, double B, Index
double deltatime_beta = 0;
double deltadoca_beta = 0;

// if(Constants.getInstance().useUSETIMETBETA()==true) {
// deltatime_beta = calcDeltaTimeBetaTFCN(this.get_Time(), tab, beta);
// }
//
// if(event.hasBank("MC::Particle")==false) {
// distance = tde.interpolateOnGrid(B, Math.toDegrees(ralpha),
// this.getCorrectedTime(this.get_Time(), deltatime_beta),
// secIdx, slIdx);
// } else {
// distance = tde.interpolateOnGrid(B, Math.toDegrees(ralpha),
// this.getCorrectedTime(this.get_Time(), 0),
// secIdx, slIdx) ;
// }
// //get deltadoca
// if(Constants.getInstance().useUSETIMETBETA()==false) {
// deltadoca_beta = calcDeltaDocaBeta(distance, tab, beta);
// }
//
// distance -=deltadoca_beta;
//this.set_DeltaTimeBeta(deltatime_beta);
//this.set_DeltaDocaBeta(deltadoca_beta);
distance = tde.interpolateOnGrid(B, Math.toDegrees(ralpha), beta, this.get_Time(), secIdx, slIdx);
if(Constants.getInstance().useUSETIMETBETA()==true) {
deltatime_beta = calcDeltaTimeBetaTFCN(this.get_Time(), tab, beta);
}

if(event.hasBank("MC::Particle")==false) {
distance = tde.interpolateOnGrid(B, Math.toDegrees(ralpha),
this.getCorrectedTime(this.get_Time(), deltatime_beta),
secIdx, slIdx);
} else {
distance = tde.interpolateOnGrid(B, Math.toDegrees(ralpha),
this.getCorrectedTime(this.get_Time(), 0),
secIdx, slIdx) ;
}
//get deltadoca
if(Constants.getInstance().useUSETIMETBETA()==false) {
deltadoca_beta = calcDeltaDocaBeta(distance, tab, beta);
}

distance -=deltadoca_beta;
this.set_DeltaTimeBeta(deltatime_beta);
this.set_DeltaDocaBeta(deltadoca_beta);

}

Expand Down
Loading
Loading