AutoCalibrator always returns a zero offset

[jvavrek]

As far as I can tell, there's no way to generate a LinearEnergyCal object via AutoCalibrator without getting a calibration that has a zero offset. I don't know if there's a nice way to incorporate a non-zero offset into the minimization of fit_gain() but a useful workaround is to use AutoCalibrator to find the best channel/energy pairs and then pass them to a full linear fit:

if np.isclose(cal.ch2kev(0), 0.):
    # These are returned sorted by becquerel
    xdata = cal.fit_channels
    ydata = cal.fit_energies

    def lin(x, p0, p1):
        return p0 + p1 * x

    popt, _ = curve_fit(lin, xdata, ydata)
    print('Found new gain: %.6f keV/channel.' % popt[1])
    print('Found new incp: %.2f keV.' % popt[0])
    cal = bq.LinearEnergyCal().from_coeffs({'m': popt[1], 'b': popt[0]})

[jccurtis]

Hey @jvavrek , thanks for this! I agree that the functionality of the autocalibrator should be expanded, especially for non-linear calibrations (i.e. NaI + PMT). @markbandstra how difficult is it to adjust multiple params in the brute force calibration search? Could we use arbitrary functions? We need to have a more general conversation about the ecal workflow since a full calibration should include (peak) fits to refine the adc_value for each feature.

[jccurtis]

Also general note:

You don't need to initialize bq.LinearEnergyCal since from_coeffs is a class method which acts as an alternate constructor:

cal = bq.LinearEnergyCal.from_coeffs({'m': popt[1], 'b': popt[0]})

In python3 the preferred string formatting is the newer .format(...) style. This method was introduced in python2.7 and is much more powerful and efficient than the older C-style % format. While it is not deprecated, bq follows the newer standard. link

print('Found new gain: {:.6f} keV/channel.'.format(popt[1]))
print('Found new incp: {:.2f} keV.'.format(popt[0]))

Or even:

print('Found new gain: {1:.6f} keV/channel.\nFound new incp: {0:.2f} keV.'.format(*popt))

😄

[markbandstra]

The autocalibrator was designed to be pretty simple, with only a single gain parameter. I agree that it should be expanded to include an offset at least. This effort will require reworking the math used (which is basically a chi-squared minimization with a brute force search over calibration line attributions) but it might not be that bad. The caveat is that the more parameters are included, the worse the attributions might get.

@jccurtis points out that restructuring this to use RANSAC might be a good approach. https://en.wikipedia.org/wiki/Random_sample_consensus