Author: [log in to unmask]
Date: Thu Aug 20 11:24:33 2015
New Revision: 3377
Log:
add docs
Modified:
java/trunk/ecal-recon/src/main/java/org/hps/recon/ecal/Ecal3PoleFunction.java
java/trunk/ecal-recon/src/main/java/org/hps/recon/ecal/EcalPulseFitter.java
Modified: java/trunk/ecal-recon/src/main/java/org/hps/recon/ecal/Ecal3PoleFunction.java
=============================================================================
--- java/trunk/ecal-recon/src/main/java/org/hps/recon/ecal/Ecal3PoleFunction.java (original)
+++ java/trunk/ecal-recon/src/main/java/org/hps/recon/ecal/Ecal3PoleFunction.java Thu Aug 20 11:24:33 2015
@@ -2,13 +2,11 @@
import hep.aida.ref.function.AbstractIFunction;
-/*
- *
+/**
* "3-Pole Function:" x**2 * exp(-x/tau)
*
* Here x is time, and we have pedestal and time offsets, and it's normalized:
* PEDESTAL + INTEGRAL / WIDTH^3 / 2 * (TIME-TIME0)**2 * exp(-(TIME-TIME0)/WIDTH)
- *
*/
public class Ecal3PoleFunction extends AbstractIFunction {
Modified: java/trunk/ecal-recon/src/main/java/org/hps/recon/ecal/EcalPulseFitter.java
=============================================================================
--- java/trunk/ecal-recon/src/main/java/org/hps/recon/ecal/EcalPulseFitter.java (original)
+++ java/trunk/ecal-recon/src/main/java/org/hps/recon/ecal/EcalPulseFitter.java Thu Aug 20 11:24:33 2015
@@ -8,8 +8,8 @@
import hep.aida.IFunction;
import hep.aida.IFunctionFactory;
-import java.io.FileWriter;
-import java.io.IOException;
+//import java.io.FileWriter;
+//import java.io.IOException;
import java.util.logging.Level;
import java.util.logging.Logger;
@@ -20,32 +20,53 @@
import org.lcsim.geometry.Detector;
import org.lcsim.util.aida.AIDA;
+/**
+ * Fit ECal FADC Mode-1 waveform with a "3-pole" function for extraction of
+ * pulse integral and time.
+ *
+ * Called from EcalRawConverter. Uses Ecal3PoleFunction.
+ *
+ * Limits are placed on pulse characteristics to avoid wasting time fitting
+ * pulses near the edge of the window. Width parameter can be free or fixed.
+ * Pedestal parameter is initialized from samples before threshold crossing.
+ */
public class EcalPulseFitter {
- public int debug = 0;
- public String fitFileName = null;
+ private EcalConditions ecalConditions = null;
+
+ /**
+ * If this is false, the width will be a free parameter in the fit:
+ */
public boolean fixShapeParameter = true;
- private final static int nsPerSample=4;
- private FileWriter fitFileWriter = null;
- private EcalConditions ecalConditions = null;
-
- // don't bother fitting pulses with threshold crossing outside this sample range:
+ /**
+ * don't bother fitting pulses with threshold crossing outside this sample range:
+ */
public int threshRange[]={7,20}; // (28 ns <--> 80 ns)
- // restrict fit's time0 parameter to this range: (units=samples)
+ /**
+ * restrict fit's time0 parameter to this range: (units=samples)
+ */
public int t0limits[]={1,30};
- // fit sample range relative to threshold crossing:
+ /**
+ * fit sample range relative to threshold crossing
+ */
private final static int fitRange[]={-10,15};
-
- // sample range relative to threshold crossing used to initialize pedestal fit parameter:
+
+ /**
+ * sample range relative to threshold crossing used to initialize pedestal fit parameter:
+ */
private final static int pedRange[]={-10,-5};
- // if this is positive, it will override individual widths:
+ /**
+ * if this is positive, it will override individual widths and use one global width (units = samples)
+ */
public double globalThreePoleWidth=-999;
- // these are 442 channels' widths, measured from 2015 data:
+ /**
+ * 442 channels' widths, measured from 2015 data (units = samples)
+ */
private final static double threePoleWidths[]={
2.44057,2.40057,2.53389,2.32095,2.44267,2.43631,2.40292,2.46911,2.36032,2.42132,2.41473,2.43864,2.37710,2.43278,2.49126,2.41739,2.40560,
2.42042,2.51278,2.46731,2.42647,2.35636,2.55308,2.47517,2.38748,2.48324,2.53726,2.54780,2.50771,2.45188,2.35011,2.35110,2.37575,2.37649,
@@ -75,49 +96,32 @@
2.38524,2.39403,2.42442,2.56829,2.56063,2.58340,2.42897,2.58796,2.56262,2.53480,2.43707,2.38258,2.46584,2.40596,2.42590,2.48281,2.47800
};
+ // Stuff for debugging:
+ //public int debug = 0;
+ //public String fitFileName = null;
+ //private FileWriter fitFileWriter = null;
+
AIDA aida = AIDA.defaultInstance();
IAnalysisFactory analysisFactory = aida.analysisFactory();
IFunctionFactory functionFactory = analysisFactory.createFunctionFactory(null);
IFitFactory fitFactory = analysisFactory.createFitFactory();
IFitter fitter=fitFactory.createFitter();
IFunction fitFcn3Pole=new Ecal3PoleFunction();
- IFunction fitFcnLinear=functionFactory.createFunctionFromScript("fitFcnLinear",1,"p0+x[0]*p1","p0,p1","",null);
IDataPointSet fitData=aida.analysisFactory().createDataPointSetFactory(null).create("ADC DataPointSet", 2);
- public IFitResult fitLeadingEdge(RawTrackerHit hit,int threshCross) {
-
- final short samples[] = hit.getADCValues();
- final double noise = findChannel(hit.getCellID()).getCalibration().getNoise();
-
- int imaxADC=0;
- for (int ii=threshCross; ii<samples.length-1; ii++)
- {
- if (ii<0) continue;
- if (ii>=samples.length) break;
- if (samples[ii+1] < samples[ii]) {
- imaxADC=ii;
- break;
- }
- }
- fitData.clear();
- int nFitPoints=0;
- for (int ii=threshCross-1; ii<imaxADC-1; ii++) {
- fitData.addPoint();
- fitData.point(nFitPoints).coordinate(0).setValue(ii);
- fitData.point(nFitPoints).coordinate(1).setValue(samples[ii]);
- fitData.point(nFitPoints).coordinate(1).setErrorMinus(noise);
- fitData.point(nFitPoints).coordinate(1).setErrorPlus(noise);
- nFitPoints++;
- }
-
- if (nFitPoints<2) return null;
-
- return fitter.fit(fitData,fitFcnLinear);
- }
-
+
+ /**
+ * Perform and return "3-pole" fit of ECAL raw waveform.
+ *
+ * @param hit the RawTrackerHit (Mode-1 FADC ECal readout) to be fit
+ * @param threshCross the sample of threshold crossing, used to initialize fit parameters
+ * @param maxADC the ADC at pulse maximum, used to initialize fit parameters
+ *
+ * @return fit result, which can be null if fit fails or certain conditions are not met.
+ */
public IFitResult fitPulse(RawTrackerHit hit,int threshCross,double maxADC) {
- if (debug>0) System.err.println("FITTING.....................................................");
+ //if (debug>0) System.err.println("FITTING.....................................................");
final short samples[] = hit.getADCValues();
final double noise = findChannel(hit.getCellID()).getCalibration().getNoise();
@@ -151,7 +155,7 @@
{
if (ii<0) continue;
if (ii>=samples.length) break;
- if (debug>0) System.err.print(ii+":"+samples[ii]+" ");
+ //if (debug>0) System.err.print(ii+":"+samples[ii]+" ");
// if (samples[ii] > maxADC) maxADC=samples[ii];
if (iPeakADC<0 && ii<samples.length-1 && samples[ii+1] < samples[ii]) {
iPeakADC=ii;
@@ -164,13 +168,13 @@
fitData.point(nFitPoints).coordinate(1).setErrorPlus(noise);
nFitPoints++;
}
- if (debug>0) System.err.print("\n");
+ //if (debug>0) System.err.print("\n");
// don't bother trying to fit:
if (nFitPoints < 10) return null;
if (maxADC < ped) return null;
- if (debug>0) System.err.println("------- "+ped+" "+threshCross+" "+maxADC);
+ //if (debug>0) System.err.println("------- "+ped+" "+threshCross+" "+maxADC);
final double pulseIntegral = sumADC-ped*nFitPoints;
@@ -188,20 +192,27 @@
fitter.fitParameterSettings("integral").setBounds(0,999999);
if (fixShapeParameter) fitter.fitParameterSettings("width").setFixed(true);
+ /*
if (debug>0) {
System.err.println(String.format("A= %8.2f",fitFcn3Pole.parameter("integral")));
System.err.println(String.format("T= %8.2f",fitFcn3Pole.parameter("time0")*4));
System.err.println(String.format("P= %8.2f",fitFcn3Pole.parameter("pedestal")));
System.err.println(String.format("S= %8.2f",fitFcn3Pole.parameter("width")));
- } else Logger.getLogger("org.freehep.math.minuit").setLevel(Level.OFF);
+ } else
+ */
+
+ // did this because something else was turning it back on every event
+ // once that thing stops doing that, this can be removed.
+ Logger.getLogger("org.freehep.math.minuit").setLevel(Level.OFF);
IFitResult fitResult = fitter.fit(fitData,fitFcn3Pole);
- writeFit(samples,fitResult,cid);
-
+
+/*
if (debug>0) {
+ writeFit(samples,fitResult,cid);
final double P = fitResult.fittedParameter("pedestal");
final double I = fitResult.fittedParameter("integral");
- final double T = fitResult.fittedParameter("time0")*nsPerSample;
+ final double T = fitResult.fittedParameter("time0")*4; // 4 ns per sample
final double S = fitResult.fittedParameter("width");
final double Q = fitResult.quality();
final double E[] = fitResult.errors();
@@ -213,8 +224,9 @@
System.err.println(String.format("M = %8.2f",maxADC-P));
System.err.println(String.format("Q = %8.2f",Q));
}
-
- /*
+*/
+
+/*
// now fit just leading edge:
IFitResult timeResult = null;
if (doTimeFit) {
@@ -225,12 +237,23 @@
fitter.fit(fitData,fitFcn3Pole);
}
IFitResult fitResults[] = {fitResult,timeResult};
- */
+*/
return fitResult;
}
-
+ public void setDetector(Detector detector) {
+ // ECAL combined conditions object.
+ ecalConditions = DatabaseConditionsManager.getInstance().getEcalConditions();
+ }
+
+ public EcalChannelConstants findChannel(long cellID) {
+ return ecalConditions.getChannelConstants(ecalConditions.getChannelCollection().findGeometric(cellID));
+ }
+
+
+
+/*
public void writeFit(short samples[],IFitResult fit,final int cid) {
if (fitFileName == null) return;
if (fitFileWriter == null) {
@@ -262,13 +285,43 @@
}
return ped;
}
- public void setDetector(Detector detector) {
- // ECAL combined conditions object.
- ecalConditions = DatabaseConditionsManager.getInstance().getEcalConditions();
- }
- public EcalChannelConstants findChannel(long cellID) {
- return ecalConditions.getChannelConstants(ecalConditions.getChannelCollection().findGeometric(cellID));
- }
+*/
+
+
+
+/*
+ public IFitResult fitLeadingEdge(RawTrackerHit hit,int threshCross) {
+
+ IFunction fitFcnLinear=functionFactory.createFunctionFromScript("fitFcnLinear",1,"p0+x[0]*p1","p0,p1","",null);
+ final short samples[] = hit.getADCValues();
+ final double noise = findChannel(hit.getCellID()).getCalibration().getNoise();
+
+ int imaxADC=0;
+ for (int ii=threshCross; ii<samples.length-1; ii++)
+ {
+ if (ii<0) continue;
+ if (ii>=samples.length) break;
+ if (samples[ii+1] < samples[ii]) {
+ imaxADC=ii;
+ break;
+ }
+ }
+ fitData.clear();
+ int nFitPoints=0;
+ for (int ii=threshCross-1; ii<imaxADC-1; ii++) {
+ fitData.addPoint();
+ fitData.point(nFitPoints).coordinate(0).setValue(ii);
+ fitData.point(nFitPoints).coordinate(1).setValue(samples[ii]);
+ fitData.point(nFitPoints).coordinate(1).setErrorMinus(noise);
+ fitData.point(nFitPoints).coordinate(1).setErrorPlus(noise);
+ nFitPoints++;
+ }
+
+ if (nFitPoints<2) return null;
+
+ return fitter.fit(fitData,fitFcnLinear);
+ }
+*/
}
|
