From fab8d7b292eaea2afdce7c15c722a5f5ead15e6d Mon Sep 17 00:00:00 2001 From: Michael Paolone Date: Wed, 1 Apr 2026 15:37:39 -0600 Subject: [PATCH] Added new T2D functional form and new wire-by-wire T2D database table call: time_to_distance_wire --- .../java/org/jlab/rec/ahdc/Hit/HitReader.java | 377 +++++++++--------- .../constants/CalibrationConstantsLoader.java | 59 ++- .../org/jlab/service/ahdc/AHDCEngine.java | 1 + 3 files changed, 236 insertions(+), 201 deletions(-) diff --git a/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Hit/HitReader.java b/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Hit/HitReader.java index f7d2e3d444..2ab61b5367 100644 --- a/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Hit/HitReader.java +++ b/reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Hit/HitReader.java @@ -19,194 +19,195 @@ public HitReader(DataEvent event, AlertDCDetector detector, boolean simulation) fetch_AHDCHits(event, detector); if (simulation) fetch_TrueAHDCHits(event); } - - public final void fetch_AHDCHits(DataEvent event, AlertDCDetector detector) { - - ArrayList hits = new ArrayList<>(); - - if (!event.hasBank("AHDC::adc")) { - this.set_AHDCHits(hits); - return; - } - - double startTime = 0.0; - if (event.hasBank("REC::Event") && !sim) { - DataBank bankRecEvent = event.getBank("REC::Event"); - startTime = bankRecEvent.getFloat("startTime", 0); - } - - RawDataBank bankDGTZ = new RawDataBank("AHDC::adc"); - bankDGTZ.read(event); - - for (int i = 0; i < bankDGTZ.rows(); i++) { - - int id = bankDGTZ.trueIndex(i) + 1; - int number = bankDGTZ.getByte("layer", i); // e.g. 11,12,21,... (this matches CCDB "layer") - int layer = number % 10; - int superlayer = (number % 100) / 10; - int sector = bankDGTZ.getInt("sector", i); - int wire = bankDGTZ.getShort("component", i); - - // RAW quantities from bank - double adcRaw = bankDGTZ.getInt("ADC", i); - double leadingEdgeTime = bankDGTZ.getFloat("leadingEdgeTime", i); - double timeOverThreshold = bankDGTZ.getFloat("timeOverThreshold", i); - double adcOffset = bankDGTZ.getFloat("ped", i); - int wfType = bankDGTZ.getShort("wfType", i); - - // CCDB key - int key_value = sector * 10000 + number * 100 + wire; - - // ----------------------------- - // Raw hit cuts - // ----------------------------- - // double[] rawHitCuts = CalibrationConstantsLoader.AHDC_RAW_HIT_CUTS.get(key_value); - //if (rawHitCuts == null) continue; - - - double[] rawHitCuts = CalibrationConstantsLoader.AHDC_RAW_HIT_CUTS.get(key_value); - if (rawHitCuts == null) {throw new IllegalStateException("Missing CCDB table /calibration/alert/ahdc/raw_hit_cuts for key=" + key_value+ " (check run/variation + key mapping)"); - } + public double T2Dfunction(int key, double time){ + double[] time2distance = CalibrationConstantsLoader.AHDC_TIME_TO_DISTANCE_WIRE.get(key); + if (time2distance == null){ + throw new IllegalStateException("Missing CCDB /calibration/alert/ahdc/time_to_distance for key=" + key + " (check run/variation + key mapping)"); + } + + // T2D function consists of three 1st order polynomials (p1, p2, p3) and two transition functions (t1, t2). + + double p1 = (time2distance[0] + time2distance[1]*time); + double p2 = (time2distance[2] + time2distance[3]*time); + double p3 = (time2distance[4] + time2distance[5]*time); + + double t1 = 1.0/(1.0 + Math.exp(-(time - time2distance[6])/time2distance[7])); + double t2 = 1.0/(1.0 + Math.exp(-(time - time2distance[8])/time2distance[9])); + + double retval = (p1)*(1.0 - t1) + (t1)*(p2)*(1.0 - t2) + (t2)*(p3); + return retval; + } + + public final void fetch_AHDCHits(DataEvent event, AlertDCDetector detector) { + + ArrayList hits = new ArrayList<>(); + + if (!event.hasBank("AHDC::adc")) { + this.set_AHDCHits(hits); + return; + } + + double startTime = 0.0; + if (event.hasBank("REC::Event") && !sim) { + DataBank bankRecEvent = event.getBank("REC::Event"); + startTime = bankRecEvent.getFloat("startTime", 0); + } + + RawDataBank bankDGTZ = new RawDataBank("AHDC::adc"); + bankDGTZ.read(event); + + for (int i = 0; i < bankDGTZ.rows(); i++) { + + int id = bankDGTZ.trueIndex(i) + 1; + int number = bankDGTZ.getByte("layer", i); // e.g. 11,12,21,... (this matches CCDB "layer") + int layer = number % 10; + int superlayer = (number % 100) / 10; + int sector = bankDGTZ.getInt("sector", i); + int wire = bankDGTZ.getShort("component", i); + + // RAW quantities from bank + double adcRaw = bankDGTZ.getInt("ADC", i); + double leadingEdgeTime = bankDGTZ.getFloat("leadingEdgeTime", i); + double timeOverThreshold = bankDGTZ.getFloat("timeOverThreshold", i); + double adcOffset = bankDGTZ.getFloat("ped", i); + int wfType = bankDGTZ.getShort("wfType", i); + + // CCDB key + int key_value = sector * 10000 + number * 100 + wire; + + // ----------------------------- + // Raw hit cuts + // ----------------------------- + // double[] rawHitCuts = CalibrationConstantsLoader.AHDC_RAW_HIT_CUTS.get(key_value); + //if (rawHitCuts == null) continue; + + + double[] rawHitCuts = CalibrationConstantsLoader.AHDC_RAW_HIT_CUTS.get(key_value); + if (rawHitCuts == null) {throw new IllegalStateException("Missing CCDB table /calibration/alert/ahdc/raw_hit_cuts for key=" + key_value+ " (check run/variation + key mapping)"); + } - double t_min = rawHitCuts[0]; - double t_max = rawHitCuts[1]; - double tot_min = rawHitCuts[2]; - double tot_max = rawHitCuts[3]; - double adc_min = rawHitCuts[4]; - double adc_max = rawHitCuts[5]; - double ped_min = rawHitCuts[6]; - double ped_max = rawHitCuts[7]; - - // ----------------------------- - // Time calibration + t->d - // ----------------------------- - //double[] timeOffsets = CalibrationConstantsLoader.AHDC_TIME_OFFSETS.get(key_value); - //if (timeOffsets == null) continue; - - // double[] timeOffsets = CalibrationConstantsLoader.AHDC_TIME_OFFSETS.get(key_value); - //if (timeOffsets == null) { - //throw new IllegalStateException("Missing AHDC time_offsets for key=" + key_value); - //} - - double[] timeOffsets = CalibrationConstantsLoader.AHDC_TIME_OFFSETS.get(key_value); + double t_min = rawHitCuts[0]; + double t_max = rawHitCuts[1]; + double tot_min = rawHitCuts[2]; + double tot_max = rawHitCuts[3]; + double adc_min = rawHitCuts[4]; + double adc_max = rawHitCuts[5]; + double ped_min = rawHitCuts[6]; + double ped_max = rawHitCuts[7]; + + // ----------------------------- + // Time calibration + t->d + // ----------------------------- + //double[] timeOffsets = CalibrationConstantsLoader.AHDC_TIME_OFFSETS.get(key_value); + //if (timeOffsets == null) continue; + + // double[] timeOffsets = CalibrationConstantsLoader.AHDC_TIME_OFFSETS.get(key_value); + //if (timeOffsets == null) { + //throw new IllegalStateException("Missing AHDC time_offsets for key=" + key_value); + //} + + double[] timeOffsets = CalibrationConstantsLoader.AHDC_TIME_OFFSETS.get(key_value); - if (timeOffsets == null) { - throw new IllegalStateException("Missing CCDB /calibration/alert/ahdc/time_offsets for key=" + key_value + " (check run/variation + key mapping)"); - } - - - - - double[] time2distance = CalibrationConstantsLoader.AHDC_TIME_TO_DISTANCE.get(10101); - if (time2distance == null) continue; - - double t0 = timeOffsets[0]; - - double p0 = time2distance[0]; - double p1 = time2distance[1]; - double p2 = time2distance[2]; - double p3 = time2distance[3]; - double p4 = time2distance[4]; - double p5 = time2distance[5]; - - double time = leadingEdgeTime - t0 - startTime; - - // ----------------------------- - // ToT correction (new CCDB) - // convention: ToT_corr = ToT_raw * totCorr - // ----------------------------- - double totUsed = timeOverThreshold; - if (!sim) { - double[] totArr = CalibrationConstantsLoader.AHDC_TIME_OVER_THRESHOLD.get(key_value); - if (totArr != null && totArr.length > 0) { - double totCorr = totArr[0]; - totUsed = timeOverThreshold * totCorr; - } - } - - // ----------------------------- - // Hit selection (cuts) - // NOTE: we cut on totUsed (corrected ToT). If you want RAW-ToT cuts, - // replace totUsed with timeOverThreshold in the condition. - // ----------------------------- - boolean passCuts = - (wfType <= 2) && - (adcRaw >= adc_min) && (adcRaw <= adc_max) && - (time >= t_min) && (time <= t_max) && - (timeOverThreshold >= tot_min) && (timeOverThreshold <= tot_max)&& - //(totUsed >= tot_min) && (totUsed <= tot_max) && - (adcOffset >= ped_min) && (adcOffset <= ped_max); - - if (!passCuts && !sim) continue; - - // ----------------------------- - // DOCA from calibrated time - // ----------------------------- - double doca = p0 - + p1*Math.pow(time, 1.0) - + p2*Math.pow(time, 2.0) - + p3*Math.pow(time, 3.0) - + p4*Math.pow(time, 4.0) - + p5*Math.pow(time, 5.0); - - if (time < 0) doca = 0.0; - - // ----------------------------- - // ADC gain calibration (new gains schema: gainCorr is index 0) - // convention: ADC_cal = ADC_raw * gainCorr - // ----------------------------- - double adcCal = adcRaw; - if (!sim) { - double[] gainArr = CalibrationConstantsLoader.AHDC_ADC_GAINS.get(key_value); - if (gainArr != null && gainArr.length > 0) { - double gainCorr = gainArr[0]; - adcCal = adcRaw * gainCorr; - } - } - - Hit h = new Hit(id, superlayer, layer, wire, doca, adcRaw, time); - h.setWirePosition(detector); - h.setADC(adcCal); // place to store calibrated ADC - h.setToT(totUsed); // place to store caibrated ToT - hits.add(h); - } - - this.set_AHDCHits(hits); - } - - public final void fetch_TrueAHDCHits(DataEvent event) { - - ArrayList truehits = new ArrayList<>(); - - if (event.hasBank("MC::True")) { - DataBank bankSIMU = event.getBank("MC::True"); - for (int i = 0; i < bankSIMU.rows(); i++) { - int pid = bankSIMU.getInt("pid", i); - double x_true = bankSIMU.getFloat("avgX", i); - double y_true = bankSIMU.getFloat("avgY", i); - double z_true = bankSIMU.getFloat("avgZ", i); - double trackE = bankSIMU.getFloat("trackE", i); - - truehits.add(new TrueHit(pid, x_true, y_true, z_true, trackE)); - } - } - - this.set_TrueAHDCHits(truehits); - } - - public ArrayList get_AHDCHits() { - return _AHDCHits; - } - - public void set_AHDCHits(ArrayList hits) { - this._AHDCHits = hits; - } - - public ArrayList get_TrueAHDCHits() { - return _TrueAHDCHits; - } - - public void set_TrueAHDCHits(ArrayList trueHits) { - this._TrueAHDCHits = trueHits; - } + if (timeOffsets == null) { + throw new IllegalStateException("Missing CCDB /calibration/alert/ahdc/time_offsets for key=" + key_value + " (check run/variation + key mapping)"); + } + + + + double t0 = timeOffsets[0]; + double time = leadingEdgeTime - t0 - startTime; + + // ----------------------------- + // ToT correction (new CCDB) + // convention: ToT_corr = ToT_raw * totCorr + // ----------------------------- + double totUsed = timeOverThreshold; + if (!sim) { + double[] totArr = CalibrationConstantsLoader.AHDC_TIME_OVER_THRESHOLD.get(key_value); + if (totArr != null && totArr.length > 0) { + double totCorr = totArr[0]; + totUsed = timeOverThreshold * totCorr; + } + } + + // ----------------------------- + // Hit selection (cuts) + // NOTE: we cut on totUsed (corrected ToT). If you want RAW-ToT cuts, + // replace totUsed with timeOverThreshold in the condition. + // ----------------------------- + boolean passCuts = + (wfType <= 2) && + (adcRaw >= adc_min) && (adcRaw <= adc_max) && + (time >= t_min) && (time <= t_max) && + (timeOverThreshold >= tot_min) && (timeOverThreshold <= tot_max)&& + //(totUsed >= tot_min) && (totUsed <= tot_max) && + (adcOffset >= ped_min) && (adcOffset <= ped_max); + + if (!passCuts && !sim) continue; + + // ----------------------------- + // DOCA from calibrated time + // ----------------------------- + double doca = T2Dfunction(key_value,time); + + if (time < 0) doca = 0.0; + + // ----------------------------- + // ADC gain calibration (new gains schema: gainCorr is index 0) + // convention: ADC_cal = ADC_raw * gainCorr + // ----------------------------- + double adcCal = adcRaw; + if (!sim) { + double[] gainArr = CalibrationConstantsLoader.AHDC_ADC_GAINS.get(key_value); + if (gainArr != null && gainArr.length > 0) { + double gainCorr = gainArr[0]; + adcCal = adcRaw * gainCorr; + } + } + + Hit h = new Hit(id, superlayer, layer, wire, doca, adcRaw, time); + h.setWirePosition(detector); + h.setADC(adcCal); // place to store calibrated ADC + h.setToT(totUsed); // place to store caibrated ToT + hits.add(h); + } + + this.set_AHDCHits(hits); + } + + public final void fetch_TrueAHDCHits(DataEvent event) { + + ArrayList truehits = new ArrayList<>(); + + if (event.hasBank("MC::True")) { + DataBank bankSIMU = event.getBank("MC::True"); + for (int i = 0; i < bankSIMU.rows(); i++) { + int pid = bankSIMU.getInt("pid", i); + double x_true = bankSIMU.getFloat("avgX", i); + double y_true = bankSIMU.getFloat("avgY", i); + double z_true = bankSIMU.getFloat("avgZ", i); + double trackE = bankSIMU.getFloat("trackE", i); + + truehits.add(new TrueHit(pid, x_true, y_true, z_true, trackE)); + } + } + + this.set_TrueAHDCHits(truehits); + } + + public ArrayList get_AHDCHits() { + return _AHDCHits; + } + + public void set_AHDCHits(ArrayList hits) { + this._AHDCHits = hits; + } + + public ArrayList get_TrueAHDCHits() { + return _TrueAHDCHits; + } + + public void set_TrueAHDCHits(ArrayList trueHits) { + this._TrueAHDCHits = trueHits; + } } diff --git a/reconstruction/alert/src/main/java/org/jlab/rec/alert/constants/CalibrationConstantsLoader.java b/reconstruction/alert/src/main/java/org/jlab/rec/alert/constants/CalibrationConstantsLoader.java index 0da9fef09c..be45c2de2e 100644 --- a/reconstruction/alert/src/main/java/org/jlab/rec/alert/constants/CalibrationConstantsLoader.java +++ b/reconstruction/alert/src/main/java/org/jlab/rec/alert/constants/CalibrationConstantsLoader.java @@ -24,21 +24,22 @@ public CalibrationConstantsLoader() { // Maps for constants from database // AHDC - public static Map AHDC_TIME_OFFSETS = new HashMap<>(); ///< {t0,dt0,extra1,extra2,chi2ndf} - public static Map AHDC_TIME_TO_DISTANCE = new HashMap<>(); ///< {p0..p5, dp0..dp5, chi2ndf} - public static Map AHDC_RAW_HIT_CUTS = new HashMap<>(); ///< {t_min,t_max,tot_min,tot_max,adc_min,adc_max,ped_min,ped_max} - + public static Map AHDC_TIME_OFFSETS = new HashMap<>(); ///< {t0,dt0,extra1,extra2,chi2ndf} + public static Map AHDC_TIME_TO_DISTANCE = new HashMap<>(); ///< {p0..p5, dp0..dp5, chi2ndf} + public static Map AHDC_TIME_TO_DISTANCE_WIRE = new HashMap<>(); ///< T2D function for every wire + public static Map AHDC_RAW_HIT_CUTS = new HashMap<>(); ///< {t_min,t_max,tot_min,tot_max,adc_min,adc_max,ped_min,ped_max} + // UPDATED SCHEMA: keys (sector,layer,component), columns: gainCorr,dgainCorr,extra1,extra2,extra3 - public static Map AHDC_ADC_GAINS = new HashMap<>(); ///< {gainCorr, dgainCorr, extra1, extra2, extra3} - + public static Map AHDC_ADC_GAINS = new HashMap<>(); ///< {gainCorr, dgainCorr, extra1, extra2, extra3} + // NEW ToT TABLE: keys (sector,layer,component), columns: totCorr,dtotCorr,extra1,extra2,extra3 - public static Map AHDC_TIME_OVER_THRESHOLD = new HashMap<>(); ///< {totCorr, dtotCorr, extra1, extra2, extra3} - - // ATOF - public static Map ATOF_EFFECTIVE_VELOCITY = new HashMap<>(); ///< {veff,dveff,extra1,extra2} - public static Map ATOF_TIME_WALK = new HashMap<>(); ///< {tw0..tw3, dtw0..dtw3, chi2ndf} - public static Map ATOF_ATTENUATION_LENGTH = new HashMap<>(); ///< {attlen,dattlen,extra1,extra2} - public static Map ATOF_TIME_OFFSETS = new HashMap<>(); ///< {t0,upstream_downstream,wedge_bar,extra1,extra2} + public static Map AHDC_TIME_OVER_THRESHOLD = new HashMap<>(); ///< {totCorr, dtotCorr, extra1, extra2, extra3} + + // ATOF + public static Map ATOF_EFFECTIVE_VELOCITY = new HashMap<>(); ///< {veff,dveff,extra1,extra2} + public static Map ATOF_TIME_WALK = new HashMap<>(); ///< {tw0..tw3, dtw0..dtw3, chi2ndf} + public static Map ATOF_ATTENUATION_LENGTH = new HashMap<>(); ///< {attlen,dattlen,extra1,extra2} + public static Map ATOF_TIME_OFFSETS = new HashMap<>(); ///< {t0,upstream_downstream,wedge_bar,extra1,extra2} public static synchronized void Load(int runno, ConstantsManager manager) { @@ -46,6 +47,7 @@ public static synchronized void Load(int runno, ConstantsManager manager) { IndexedTable ahdc_timeOffsets = manager.getConstants(runno, "/calibration/alert/ahdc/time_offsets"); IndexedTable ahdc_time2distance = manager.getConstants(runno, "/calibration/alert/ahdc/time_to_distance"); + IndexedTable ahdc_time2distanceWire = manager.getConstants(runno, "/calibration/alert/ahdc/time_to_distance_wire"); IndexedTable ahdc_rawHitCuts = manager.getConstants(runno, "/calibration/alert/ahdc/raw_hit_cuts"); // Gains table (UPDATED SCHEMA) @@ -105,6 +107,37 @@ public static synchronized void Load(int runno, ConstantsManager manager) { double params[] = { p0, p1, p2, p3, p4, p5, dp0, dp1, dp2, dp3, dp4, dp5, chi2ndf }; AHDC_TIME_TO_DISTANCE.put(Integer.valueOf(key), params); } + + //////////////////////////////////////////////////////////////////////////////////////////////////////////////// + // AHDC time to distance per wire + // See implimentation and functional form in reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Hit/HitReader.java + for (int i = 0; i < ahdc_time2distanceWire.getRowCount(); i++) { + int sector = Integer.parseInt((String) ahdc_time2distanceWire.getValueAt(i, 0)); + int layer = Integer.parseInt((String) ahdc_time2distanceWire.getValueAt(i, 1)); + int component = Integer.parseInt((String) ahdc_time2distanceWire.getValueAt(i, 2)); + + double p1_int = ahdc_time2distanceWire.getDoubleValue("p1_int", sector, layer, component); + double p1_slope = ahdc_time2distanceWire.getDoubleValue("p1_slope", sector, layer, component); + double p2_int = ahdc_time2distanceWire.getDoubleValue("p2_int", sector, layer, component); + double p2_slope = ahdc_time2distanceWire.getDoubleValue("p2_slope", sector, layer, component); + double p3_int = ahdc_time2distanceWire.getDoubleValue("p3_int", sector, layer, component); + double p3_slope = ahdc_time2distanceWire.getDoubleValue("p3_slope", sector, layer, component); + double t1_x0 = ahdc_time2distanceWire.getDoubleValue("t1_x0", sector, layer, component); + double t1_width = ahdc_time2distanceWire.getDoubleValue("t1_width", sector, layer, component); + double t2_x0 = ahdc_time2distanceWire.getDoubleValue("t2_x0", sector, layer, component); + double t2_width = ahdc_time2distanceWire.getDoubleValue("t2_width", sector, layer, component); + double z0 = ahdc_time2distanceWire.getDoubleValue("z0", sector, layer, component); + double z1 = ahdc_time2distanceWire.getDoubleValue("z1", sector, layer, component); + double z2 = ahdc_time2distanceWire.getDoubleValue("z2", sector, layer, component); + double extra1 = ahdc_time2distanceWire.getDoubleValue("extra1", sector, layer, component); + double extra2 = ahdc_time2distanceWire.getDoubleValue("extra2", sector, layer, component); + double chi2ndf = ahdc_time2distanceWire.getDoubleValue("chi2ndf", sector, layer, component); + + int key = sector * 10000 + layer * 100 + component; + double params[] = { p1_int, p1_slope, p2_int, p2_slope, p3_int, p3_slope, t1_x0, t1_width, t2_x0, t2_width, z0, z1, z2, extra1, extra2, chi2ndf }; + AHDC_TIME_TO_DISTANCE_WIRE.put(Integer.valueOf(key), params); + } + //////////////////////////////////////////////////////////////////////////////////////////////////////////////// // AHDC raw hit cuts 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 eb3dbdca24..4bd3ab1115 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 @@ -88,6 +88,7 @@ else if (Objects.equals(this.getEngineConfigString("Mode"), ModeTrackFinding.CV_ String[] alertTables = new String[] { "/calibration/alert/ahdc/time_offsets", "/calibration/alert/ahdc/time_to_distance", + "/calibration/alert/ahdc/time_to_distance_wire", "/calibration/alert/ahdc/raw_hit_cuts", "/calibration/alert/atof/effective_velocity", "/calibration/alert/atof/time_walk",