diff --git a/common-tools/clas-detector/src/main/java/org/jlab/detector/pulse/ModeAHDC.java b/common-tools/clas-detector/src/main/java/org/jlab/detector/pulse/ModeAHDC.java index a1beb05060..862915e17c 100644 --- a/common-tools/clas-detector/src/main/java/org/jlab/detector/pulse/ModeAHDC.java +++ b/common-tools/clas-detector/src/main/java/org/jlab/detector/pulse/ModeAHDC.java @@ -26,7 +26,7 @@ public class ModeAHDC extends HipoExtractor { //Saturation threshold should be 4095 (2^12-1) //But in practice waveforms saturate below it - private final short ADC_LIMIT = 3500; + private final short ADC_LIMIT = 3600; //number of consecutive samples exceeding threshold to consider a saturation plateau private final int consecutiveSaturatingSamples = 3; //Sampling time in ns @@ -113,8 +113,11 @@ public void baselineComputation(AHDCPulse p) * @param type an int that is the new wf type to be applied */ public void assignValidType(int type, AHDCPulse p){ - //We only assign a new type if a more restrictive type was not defined before - if(type>p.wftype) p.wftype = (short) type; + // if a signal has already been classified as type 1 (saturated), do nothing + if (p.wftype != 1) { + //We only assign a new type if a more restrictive type was not defined before + if(type>p.wftype) p.wftype = (short) type; + } } /** @@ -186,18 +189,26 @@ public int waveformADCProcessing(AHDCPulse p) //Define the pulse time as the peak time p.time = (p.binMax + p.time_ZS)*this.samplingTime; + + // Compute adcMax + if (p.adcMax == 0) { + assignValidType(5, p); + return 0; + } // classification before the fit + + //If there are 2*npts+1 points around the peak //(if the peak is not at an edge of the window) //Then the peak ADC value is revisited to be the average of these //TO DO: just use the adcmax??? int npts = 1; - if (p.adcMax == 0) { assignValidType(5, p);} // classification before the fit if ((p.binMax - npts >= 0) && (p.binMax + npts <= p.numberOfBins - 1)){ p.adcMax = 0; for (int bin = p.binMax - npts; bin <= p.binMax + npts; bin++) p.adcMax += p.samples[bin]; p.adcMax = p.adcMax/(2*npts+1); } + return 0; } @@ -215,64 +226,101 @@ public int waveformCFAprocessing(AHDCPulse p){ p.leadingEdgeTime = -9999; p.trailingEdgeTime = -9999; p.timeOverThreshold = -9999; - - //Set the CFA threshold - float threshold = this.amplitudeFractionCFA*p.adcMax; - - //Crossing the threshold before the peak - int binRise = -99; - for (int bin = 0; bin < p.binMax - 1; bin++){ - if (p.samples[bin] < threshold && p.samples[bin+1] >= threshold) - binRise = bin; //Here we keep only the last time the signal crosses the threshold - } - //If the waveform does not cross the ths before the peak - //it was early or remant from a previous event and we cannot define a leading time - //we set the time at zero - if(binRise==-99) { - this.assignValidType(5, p); - p.leadingEdgeTime = (0 + p.time_ZS)*this.samplingTime; + + if (p.binMax < 0) return 0; + + // update adcMax for saturated signals + float slope_max = 0; + int bin_slope_max = -1; + if (p.wftype == 1) { + for (int bin = 0; bin < p.binMax; bin++) { + float slope = (p.samples[bin+1] - p.samples[bin])/samplingTime; + if (slope > slope_max) { + slope_max = slope; + bin_slope_max = bin; + } + } + // update adcMax: empirical formula + p.adcMax = -16.371574f + 191.277306f*slope_max; } - else - { + + // threshold + float threshold = this.amplitudeFractionCFA*p.adcMax; + + if (threshold < ADC_LIMIT) { + + //Crossing the threshold before the peak + int binRise = -99; float slopeRise = 0; - //linear interpolation - //threshold = leadingtime*(ADC1-ADC0)+ADC0 - if (binRise + 1 <= p.numberOfBins-1) - slopeRise = p.samples[binRise+1] - p.samples[binRise]; - float fittedBinRise = (slopeRise == 0) ? binRise : binRise + (threshold - p.samples[binRise])/slopeRise; - p.leadingEdgeTime = (fittedBinRise + p.time_ZS)*this.samplingTime; - } - - //Crossing the threshold back down after the peak - int binFall = -99; - for (int bin = p.binMax; bin < p.numberOfBins-1; bin++){ - if (p.samples[bin] > threshold - && p.samples[bin+1] <= threshold - && p.samples[bin]>0 - && p.samples[bin+1]>0){ - binFall = bin+1; - break; //We keep only the first time the signal crosses back the threshold + for (int bin = 0; bin < p.binMax - 1; bin++){ + if (p.samples[bin] < threshold && p.samples[bin+1] >= threshold && p.samples[bin+1] < ADC_LIMIT) { + binRise = bin; //Here we keep only the last time the signal crosses the threshold + slopeRise = p.samples[binRise+1] - p.samples[binRise]; + } + else if (bin > 0 && p.samples[bin] < threshold && p.samples[bin+1] >= threshold) { // we enter here if p.samples[bin+1] >= ADC_LIMIT; in that scenario, the slope is underestimated + binRise = bin; + slopeRise = p.samples[binRise] - p.samples[binRise-1]; // we extrapolate using the previous slope + } + } + //If the waveform does not cross the ths before the peak + //it was early or remant from a previous event and we cannot define a leading time + //we set the time at zero + if(binRise==-99) { + this.assignValidType(5, p); + p.leadingEdgeTime = (0 + p.time_ZS)*this.samplingTime; + } + else + { + //linear interpolation + //threshold = leadingtime*(ADC1-ADC0)+ADC0 + float fittedBinRise = (slopeRise == 0) ? binRise : binRise + (threshold - p.samples[binRise])/slopeRise; + p.leadingEdgeTime = (fittedBinRise + p.time_ZS)*this.samplingTime; + } + + //Crossing the threshold back down after the peak + int binFall = -99; + for (int bin = p.binMax; bin < p.numberOfBins-1; bin++){ + if (p.samples[bin] > threshold + && p.samples[bin+1] <= threshold + && p.samples[bin]>0 + && p.samples[bin+1]>0){ + binFall = bin+1; + break; //We keep only the first time the signal crosses back the threshold + } + } + //If the waveform does not cross the ths again + //it was late and falling edge cannot be defined + //the time correspond to the end of the signal + if(binFall==-99) { + this.assignValidType(2, p); + p.trailingEdgeTime = (p.effectiveNumberOfBins -1 + p.time_ZS)*this.samplingTime; + } + else + { + float slopeFall = 0; + if (binFall - 1 >= 0) + slopeFall = p.samples[binFall] - p.samples[binFall-1]; + float fittedBinFall = (slopeFall == 0) ? binFall : binFall-1 + (threshold - p.samples[binFall-1])/slopeFall; + p.trailingEdgeTime = (fittedBinFall + p.time_ZS)*this.samplingTime; + } + + p.timeOverThreshold = p.trailingEdgeTime - p.leadingEdgeTime; + //if ((this.pulse.timeOverThreshold < 300) || (this.pulse.timeOverThreshold > 750)) this.assignValidType(4); // bad tot + if (p.timeOverThreshold < 300) this.assignValidType(4, p); // tot too small + } + else { // there are necessary type 1, so we use the slope_max + if (p.wftype == 1) { + // extraplote the leading edge of the signal using the slope: caution, slope_max is in ADC/ns + float extrapolatedTime = (threshold-p.samples[bin_slope_max])/slope_max + samplingTime*bin_slope_max; + p.leadingEdgeTime = extrapolatedTime + p.time_ZS*samplingTime; + p.timeOverThreshold = p.effectiveNumberOfBins*samplingTime - extrapolatedTime; + } else { // unexpected behavior + p.leadingEdgeTime = 0; + p.timeOverThreshold = 0; } } - //If the waveform does not cross the ths again - //it was late and falling edge cannot be defined - //the time correspond to the end of the signal - if(binFall==-99) { - this.assignValidType(2, p); - p.trailingEdgeTime = (p.effectiveNumberOfBins -1 + p.time_ZS)*this.samplingTime; - } - else - { - float slopeFall = 0; - if (binFall - 1 >= 0) - slopeFall = p.samples[binFall] - p.samples[binFall-1]; - float fittedBinFall = (slopeFall == 0) ? binFall : binFall-1 + (threshold - p.samples[binFall-1])/slopeFall; - p.trailingEdgeTime = (fittedBinFall + p.time_ZS)*this.samplingTime; - } - p.timeOverThreshold = p.trailingEdgeTime - p.leadingEdgeTime; - //if ((this.pulse.timeOverThreshold < 300) || (this.pulse.timeOverThreshold > 750)) this.assignValidType(4); // bad tot - if (p.timeOverThreshold < 300) this.assignValidType(4, p); // tot too small + return 0; }