Skip to content

Commit c39bcf9

Browse files
authored
use a new waveform classification (#849)
1 parent 21407b6 commit c39bcf9

File tree

2 files changed

+67
-60
lines changed
  • common-tools/clas-detector/src/main/java/org/jlab/detector/pulse
  • reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Hit

2 files changed

+67
-60
lines changed

common-tools/clas-detector/src/main/java/org/jlab/detector/pulse/ModeAHDC.java

Lines changed: 62 additions & 55 deletions
Original file line numberDiff line numberDiff line change
@@ -37,7 +37,7 @@ public class ModeAHDC extends HipoExtractor {
3737
private final int binDelayCFD = 5;
3838
private final float fractionCFD = 0.5f;
3939
//Number of samples defining the baseline
40-
private final int nBaselineSamples = 5;
40+
private final int nBaselineSamples = 4;
4141
//threshold for the slope to consider part of a wf flat
4242
private final int flateness = 200;
4343
//ADC offset to be considered as the default baseline
@@ -53,7 +53,9 @@ public class ModeAHDC extends HipoExtractor {
5353
private Pulse pulse;
5454
private long time_ZS;
5555
private int binMax;
56-
56+
private int numberOfBins;
57+
private int effectiveNumberOfBins;
58+
5759
//Waveform types:
5860
//0 is good,
5961
//-1 is invalid,
@@ -64,13 +66,24 @@ public class ModeAHDC extends HipoExtractor {
6466
//5 has low ADC ("flat")
6567

6668
/**
67-
* This method checks if the wf starts with a flat baseline
68-
*/
69-
public void baselineSlope()
70-
{
71-
int slope = this.samples[this.nBaselineSamples-1] - this.samples[0];
72-
if(slope<-this.flateness){this.pulse.wftype = 4;}//remnant from a previous sample, will be checked further later
73-
else if(slope>this.flateness){this.pulse.wftype = 2;}//early waveform
69+
*
70+
*/
71+
public void preProcess() {
72+
assignValidType(0);
73+
this.numberOfBins = samples.length;
74+
// check if the signal is too short
75+
int countNotZeros = 0;
76+
boolean FirstNonZeroFromTheTailReached = false;
77+
for (int i = samples.length-1; i >= 0; i--) {
78+
if (samples[i] > 0) {
79+
countNotZeros++;
80+
if (!FirstNonZeroFromTheTailReached) {
81+
FirstNonZeroFromTheTailReached = true;
82+
this.effectiveNumberOfBins = i+1; // this value can be identified as the effective number of bins
83+
}
84+
}
85+
}
86+
if (countNotZeros <= 10) assignInvalidType();
7487
}
7588

7689
/**
@@ -83,20 +96,13 @@ public void baselineComputation()
8396
float adcOffset = -99;
8497

8598
//Check if the wf has sufficient length, if not it is invalid
86-
if (this.samples.length >= this.nBaselineSamples+1){
87-
88-
//Compute the slope of the baseline to know how to adress that wf shape
89-
this.baselineSlope();
99+
if (numberOfBins >= this.nBaselineSamples+1){
90100

91-
//If it is an early pulse, pedestal should be read from DB default
92-
if(this.pulse.wftype == 2) {adcOffset = this.defaultBaseline;} //TO DO: Here we should read ccdb!
93-
else{
94-
//The baseline is the average of the first few samples
95-
for(int i=0; i<this.nBaselineSamples;i++){
96-
adcOffset += this.samples[i];
97-
}
98-
adcOffset = adcOffset/this.nBaselineSamples;
101+
//The baseline is the average of the first few samples
102+
for(int i=0; i<this.nBaselineSamples;i++){
103+
adcOffset += this.samples[i];
99104
}
105+
adcOffset = adcOffset/this.nBaselineSamples;
100106
}
101107
else {this.assignInvalidType();}//invalid, too short wf
102108
this.pulse.pedestal = adcOffset;
@@ -131,14 +137,14 @@ public void assignInvalidType(){
131137
*/
132138
public int waveformADCProcessing()
133139
{
134-
int binNumber = this.samples.length;
140+
//int numberOfBins = this.samples.length;
135141

136142
//Initialize to dummy values
137143
this.pulse.adcMax = -99;
138144
this.binMax = -99;
139145

140146
//For invalid wf, dummy values kept
141-
if(this.pulse.wftype == -1) return 0;
147+
//if(this.pulse.wftype == 6) return 0;
142148

143149
//Checking for saturation
144150
int nsaturating = 0;
@@ -148,11 +154,11 @@ public int waveformADCProcessing()
148154
int endPlateau = -99;
149155

150156
//Looping through samples
151-
for (int bin = 0; bin < binNumber; bin++){
157+
for (int bin = 0; bin < numberOfBins; bin++){
152158
//Baseline subtraction
153159
this.samples[bin] = (short) Math.max(this.samples[bin] - this.pulse.pedestal, 0);
154160
//Look for maximum
155-
if (this.pulse.adcMax < this.samples[bin]){
161+
if ((this.pulse.adcMax < this.samples[bin]) && (bin >= nBaselineSamples)) { // the max should not be in the baseline
156162
this.pulse.adcMax = this.samples[bin];
157163
this.binMax = bin;
158164
}
@@ -173,14 +179,6 @@ public int waveformADCProcessing()
173179
else nsaturating = 0;//this sample is not in a row of saturating samples
174180
}
175181

176-
//If the signal is flat i.e if the ADC max is the same height as the baseline
177-
//within the "flatness" threshold
178-
if(this.pulse.adcMax < this.flateness) this.assignValidType(5);
179-
180-
//Signal if the waveform is at the beginning of the window
181-
//This helps again to get rid of remnants from previous waveforms
182-
if(this.binMax<this.nBaselineSamples) this.assignValidType(4);
183-
184182
//Saturating samples: if the number of consecutive saturating samples
185183
//is greater than the limit, we consider the wf saturated
186184
if(maxNsaturating>=this.consecutiveSaturatingSamples){
@@ -191,15 +189,18 @@ public int waveformADCProcessing()
191189

192190
//Define the pulse time as the peak time
193191
this.pulse.time = (this.binMax + this.time_ZS)*this.samplingTime;
194-
//If there are five points around the peak
192+
//If there are 2*npts+1 points around the peak
195193
//(if the peak is not at an edge of the window)
196194
//Then the peak ADC value is revisited to be the average of these
197195
//TO DO: just use the adcmax???
198-
if ((this.binMax - 2 >= 0) && (this.binMax + 2 <= binNumber - 1)){
196+
int npts = 1;
197+
if (this.pulse.adcMax == 0) { assignValidType(5);} // classification before the fit
198+
if ((this.binMax - npts >= 0) && (this.binMax + npts <= numberOfBins - 1)){
199199
this.pulse.adcMax = 0;
200-
for (int bin = this.binMax - 2; bin <= this.binMax + 2; bin++) this.pulse.adcMax += this.samples[bin];
201-
this.pulse.adcMax = this.pulse.adcMax/5;
200+
for (int bin = this.binMax - npts; bin <= this.binMax + npts; bin++) this.pulse.adcMax += this.samples[bin];
201+
this.pulse.adcMax = this.pulse.adcMax/(2*npts+1);
202202
}
203+
203204
return 0;
204205
}
205206

@@ -213,17 +214,11 @@ public int waveformADCProcessing()
213214
*/
214215
public int waveformCFAprocessing(){
215216

216-
int binNumber = this.samples.length;
217-
218217
//Initialization
219218
this.pulse.leadingEdgeTime = -9999;
220219
this.pulse.trailingEdgeTime = -9999;
221220
this.pulse.timeOverThreshold = -9999;
222221

223-
//Waveforms for which timing cannot be defined
224-
//compute the values even if the type may seem bad
225-
//if(this.pulse.wftype == -1 || this.pulse.wftype>=4) return 0;
226-
227222
//Set the CFA threshold
228223
float threshold = this.amplitudeFractionCFA*this.pulse.adcMax;
229224

@@ -234,22 +229,26 @@ public int waveformCFAprocessing(){
234229
binRise = bin; //Here we keep only the last time the signal crosses the threshold
235230
}
236231
//If the waveform does not cross the ths before the peak
237-
//it was early and we cannot define a leading time
238-
if(binRise==-99) this.assignValidType(4);
232+
//it was early or remant from a previous event and we cannot define a leading time
233+
//we set the time at zero
234+
if(binRise==-99) {
235+
this.assignValidType(5);
236+
this.pulse.leadingEdgeTime = (0 + this.time_ZS)*this.samplingTime;
237+
}
239238
else
240239
{
241240
float slopeRise = 0;
242241
//linear interpolation
243242
//threshold = leadingtime*(ADC1-ADC0)+ADC0
244-
if (binRise + 1 <= binNumber-1)
243+
if (binRise + 1 <= numberOfBins-1)
245244
slopeRise = this.samples[binRise+1] - this.samples[binRise];
246245
float fittedBinRise = (slopeRise == 0) ? binRise : binRise + (threshold - this.samples[binRise])/slopeRise;
247246
this.pulse.leadingEdgeTime = (fittedBinRise + this.time_ZS)*this.samplingTime;
248247
}
249248

250249
//Crossing the threshold back down after the peak
251250
int binFall = -99;
252-
for (int bin = binMax; bin < binNumber-1; bin++){
251+
for (int bin = binMax; bin < numberOfBins-1; bin++){
253252
if (this.samples[bin] > threshold
254253
&& this.samples[bin+1] <= threshold
255254
&& this.samples[bin]>0
@@ -260,7 +259,11 @@ public int waveformCFAprocessing(){
260259
}
261260
//If the waveform does not cross the ths again
262261
//it was late and falling edge cannot be defined
263-
if(binFall==-99) {this.assignValidType(3); return 0;}
262+
//the time correspond to the end of the signal
263+
if(binFall==-99) {
264+
this.assignValidType(2);
265+
this.pulse.trailingEdgeTime = (this.effectiveNumberOfBins -1 + this.time_ZS)*this.samplingTime;
266+
}
264267
else
265268
{
266269
float slopeFall = 0;
@@ -271,6 +274,8 @@ public int waveformCFAprocessing(){
271274
}
272275

273276
this.pulse.timeOverThreshold = this.pulse.trailingEdgeTime - this.pulse.leadingEdgeTime;
277+
//if ((this.pulse.timeOverThreshold < 300) || (this.pulse.timeOverThreshold > 750)) this.assignValidType(4); // bad tot
278+
if (this.pulse.timeOverThreshold < 300) this.assignValidType(4); // tot too small
274279

275280
return 0;
276281
}
@@ -287,19 +292,19 @@ public int computeTimeUsingConstantFractionDiscriminator(){
287292
this.pulse.constantFractionTime = -99;
288293
if(this.pulse.wftype>=4) return 0;
289294

290-
int binNumber = this.samples.length;
291-
float[] signal = new float[binNumber];
295+
//int numberOfBins = this.samples.length;
296+
float[] signal = new float[numberOfBins];
292297

293298
// signal generation
294-
for (int bin = 0; bin < binNumber; bin++){
299+
for (int bin = 0; bin < numberOfBins; bin++){
295300
signal[bin] = (1 - this.fractionCFD)*this.samples[bin]; // we fill it with a fraction of the original signal
296-
if (bin < binNumber - this.binDelayCFD)
301+
if (bin < numberOfBins - this.binDelayCFD)
297302
signal[bin] += -1*this.fractionCFD*this.samples[bin + this.binDelayCFD]; // we advance and invert a complementary fraction of the original signal and superimpose it to the previous signal
298303
}
299304
// determine the two humps
300305
int binHumpSup = 0;
301306
int binHumpInf = 0;
302-
for (int bin = 0; bin < binNumber; bin++){
307+
for (int bin = 0; bin < numberOfBins; bin++){
303308
if (signal[bin] > signal[binHumpSup])
304309
binHumpSup = bin;
305310
}
@@ -314,7 +319,7 @@ public int computeTimeUsingConstantFractionDiscriminator(){
314319
binZero = bin; // last pass below zero
315320
} // at this stage : binZero < constantFractionTime/samplingTime <= binZero + 1 // constantFractionTime is determined by assuming a linear fit between binZero and binZero + 1
316321
float slopeCFD = 0;
317-
if (binZero + 1 <= binNumber)
322+
if (binZero + 1 <= numberOfBins)
318323
slopeCFD = signal[binZero+1] - signal[binZero];
319324
float fittedBinZero = (slopeCFD == 0) ? binZero : binZero + (0 - signal[binZero])/slopeCFD;
320325
this.pulse.constantFractionTime = (fittedBinZero + this.time_ZS)*this.samplingTime;
@@ -362,7 +367,9 @@ public List<Pulse> extract(NamedEntry pars, int id, long timestamp, long time_ZS
362367
this.binMax = -99;
363368

364369
List<Pulse> output = new ArrayList<>();
365-
370+
371+
// Pre Process
372+
this.preProcess();
366373
//Baseline computation
367374
this.baselineComputation();
368375

reconstruction/alert/src/main/java/org/jlab/rec/ahdc/Hit/HitReader.java

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -46,8 +46,9 @@ public final void fetch_AHDCHits(DataEvent event, AlertDCDetector detector) {
4646
int wire = bankDGTZ.getShort("component", i);
4747
double adc = bankDGTZ.getInt("ADC", i);
4848
double leadingEdgeTime = bankDGTZ.getFloat("leadingEdgeTime", i);
49-
//double timeOverThreshold = bankDGTZ.getFloat("timeOverThreshold", i);
50-
//double adcOffset = bankDGTZ.getFloat("ped", i);
49+
double timeOverThreshold = bankDGTZ.getFloat("timeOverThreshold", i);
50+
double adcOffset = bankDGTZ.getFloat("ped", i);
51+
int wfType = bankDGTZ.getShort("wfType", i);
5152
// Retrieve raw hit cuts from CCDB
5253
int key_value = sector*10000 + number*100 + wire;
5354
double[] rawHitCuts = CalibrationConstantsLoader.AHDC_RAW_HIT_CUTS.get( key_value );
@@ -76,9 +77,8 @@ public final void fetch_AHDCHits(DataEvent event, AlertDCDetector detector) {
7677
// Remark: leadingEdgeTime already has the fine timestamp correction
7778
double time = leadingEdgeTime - t0 - startTime;
7879

79-
if ((bankDGTZ.getShort("wfType", i) <= 1) || sim) {
80-
// Apply raw hit cuts
81-
//if (((adc >= adc_min) && (adc <= adc_max) && (time >= t_min) && (time <= t_max) && (timeOverThreshold >= tot_min) && (timeOverThreshold <= tot_max) && (adcOffset >= ped_min) && (adcOffset <= ped_max)) || sim) {
80+
// Hit selection : wfType and additional cuts
81+
if (((wfType <= 2) && (adc >= adc_min) && (adc <= adc_max) && (time >= t_min) && (time <= t_max) && (timeOverThreshold >= tot_min) && (timeOverThreshold <= tot_max) && (adcOffset >= ped_min) && (adcOffset <= ped_max)) || sim) {
8282
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);
8383
if (sim) {
8484
time += 5; // correction from mctime - systematic error of the decoding

0 commit comments

Comments
 (0)