diff --git a/src/calng/correction_addons/peakfinder9.py b/src/calng/correction_addons/peakfinder9.py
index 4a9d199f4c7dd77cadb579d827160d40bee685be..d70353757c4dbf667529b657cfd46f4de8d3dbe2 100644
--- a/src/calng/correction_addons/peakfinder9.py
+++ b/src/calng/correction_addons/peakfinder9.py
@@ -17,6 +17,17 @@ class Peakfinder9(BaseCorrectionAddon):
(
UINT32_ELEMENT(schema)
.key(f"{prefix}.windowRadius")
+ .description(
+ "The algorithm works by first deciding (in parallel, looking at each "
+ "pixel) if 'current pixel' is a candidate as start of a peak. Then, it "
+ "moves out in rings, expanding peak up until at most windowRadius. "
+ "The border of the window is used for background estimation. The other "
+ "parameters adjust rejection criteria for the 'current pixel' "
+ "(minPeakValueOverNeighbors, minSnrMaxPixel), criteria for adding "
+ "pixels to the peak (minSnrPeakPixels), and rejection criteria for the "
+ "full peak after expansion (minSnrWholePeak, minPeakPixels, "
+ "maxPeakPixels)."
+ )
.tags("managed")
.assignmentOptional()
.defaultValue(2)
@@ -25,6 +36,10 @@ class Peakfinder9(BaseCorrectionAddon):
UINT32_ELEMENT(schema)
.key(f"{prefix}.maxPeaks")
+ .description(
+ "For GPU buffer reasons, we need an upper limit on how many peaks "
+ "there will be per frame. Just leave this reasonably high."
+ )
.tags("managed")
.assignmentOptional()
.defaultValue(500)
@@ -33,6 +48,11 @@ class Peakfinder9(BaseCorrectionAddon):
FLOAT_ELEMENT(schema)
.key(f"{prefix}.minPeakValueOverNeighbors")
+ .description(
+ "Current pixel is rejected as start of peak if it does not exceed the "
+ "maximum pixel value on the window border by at least this amount; "
+ "note, this is in absolute terms unlike other parameters."
+ )
.tags("managed")
.assignmentOptional()
.defaultValue(10)
@@ -41,6 +61,11 @@ class Peakfinder9(BaseCorrectionAddon):
FLOAT_ELEMENT(schema)
.key(f"{prefix}.minSnrMaxPixel")
+ .description(
+ "Current pixel is rejected as start of peak if it does not exceed mean "
+ "by at least mean by this amount times sigma. Both mean and sigma are "
+ "based on window border pixels."
+ )
.tags("managed")
.assignmentOptional()
.defaultValue(5)
@@ -49,6 +74,12 @@ class Peakfinder9(BaseCorrectionAddon):
FLOAT_ELEMENT(schema)
.key(f"{prefix}.minSnrPeakPixels")
+ .description(
+ "If current pixel was not rejected until now, we start expanding peak "
+ "layer by layer. For a pixel to be added to the peak, it must (be in "
+ "the next layer, and) exceed the mean by at least this number times "
+ "sigma."
+ )
.tags("managed")
.assignmentOptional()
.defaultValue(4)
@@ -57,14 +88,47 @@ class Peakfinder9(BaseCorrectionAddon):
FLOAT_ELEMENT(schema)
.key(f"{prefix}.minSnrWholePeak")
+ .description(
+ "After expanding peak with additional peak pixels, the whole peak is "
+ "rejected if its total mass does not exceed the mean by at least this "
+ "number times sigma."
+ )
.tags("managed")
.assignmentOptional()
.defaultValue(6)
.reconfigurable()
.commit(),
+ UINT32_ELEMENT(schema)
+ .key(f"{prefix}.minPeakPixels")
+ .description(
+ "After expanding peak, it is rejected if the number of pixels is not "
+ "at least this."
+ )
+ .tags("managed")
+ .assignmentOptional()
+ .defaultValue(1)
+ .reconfigurable()
+ .commit(),
+
+ UINT32_ELEMENT(schema)
+ .key(f"{prefix}.maxPeakPixels")
+ .description(
+ "After expanding peak, it is rejected if the number of pixels exceeds "
+ "this."
+ )
+ .tags("managed")
+ .assignmentOptional()
+ .defaultValue(100)
+ .reconfigurable()
+ .commit(),
+
FLOAT_ELEMENT(schema)
.key(f"{prefix}.minSigma")
+ .description(
+ "When background estimate for stdev is very low, this is used to set "
+ "minimum sigma value used for the thresholds above."
+ )
.tags("managed")
.assignmentOptional()
.defaultValue(5)
@@ -73,6 +137,7 @@ class Peakfinder9(BaseCorrectionAddon):
UINT32_ELEMENT(schema)
.key(f"{prefix}.blockX")
+ .description("GPU execution tuning parameter you probably want to ignore.")
.tags("managed")
.assignmentOptional()
.defaultValue(1)
@@ -81,6 +146,7 @@ class Peakfinder9(BaseCorrectionAddon):
UINT32_ELEMENT(schema)
.key(f"{prefix}.blockY")
+ .description("GPU execution tuning parameter you probably want to ignore.")
.tags("managed")
.assignmentOptional()
.defaultValue(1)
@@ -89,6 +155,7 @@ class Peakfinder9(BaseCorrectionAddon):
UINT32_ELEMENT(schema)
.key(f"{prefix}.blockZ")
+ .description("GPU execution tuning parameter you probably want to ignore.")
.tags("managed")
.assignmentOptional()
.defaultValue(64)
@@ -214,6 +281,8 @@ class Peakfinder9(BaseCorrectionAddon):
cupy.float32(self._config["minSnrPeakPixels"]),
cupy.float32(self._config["minSnrWholePeak"]),
cupy.float32(self._config["minSigma"]),
+ cupy.uint32(self._config["minPeakPixels"]),
+ cupy.uint32(self._config["maxPeakPixels"]),
cupy.uint32(self._config["maxPeaks"]),
)
diff --git a/src/calng/kernels/peakfinder9_gpu.cu b/src/calng/kernels/peakfinder9_gpu.cu
index a07bdfb3c27140dd9c807cc6836abdd1480f8932..159e9584c7c02018b1869673a0aa8e4d664c434d 100644
--- a/src/calng/kernels/peakfinder9_gpu.cu
+++ b/src/calng/kernels/peakfinder9_gpu.cu
@@ -58,14 +58,19 @@ public:
extern "C" __global__ void pf9(const unsigned short num_frames,
const unsigned short ss_dim,
const unsigned short fs_dim,
- const unsigned short window_radius,
- const float min_peak_over_border,
- const float min_snr_biggest_pixel,
- const float min_snr_peak_pixels,
- const float min_snr_whole_peak,
- const float min_sigma,
- const unsigned int max_peaks,
- const float* image,
+ const unsigned short window_radius, // radius to border: border is for background estimate and bounds peak size
+ const float min_peak_over_border, // for "current pixel" to be candidate, it must be this much (absolute value) over max value on border
+ // for the next parameters, snr is based on stdev of border pixels
+ // note: hardcoded limit says there must be at least 10 non-masked border pixels for this
+ const float min_snr_biggest_pixel, // "current pixel" is compared against this
+ // we "expand" peak in rings
+ const float min_snr_peak_pixels, // additional peak pixels get compared against this
+ const float min_snr_whole_peak, // total mass of peak after expansion is compared against this
+ const float min_sigma, // when background has suspiciously low stdev, this sigma kicks in
+ const unsigned min_peak_pixels, // if less than this number of pixels are considered part of peak, reject
+ const unsigned max_peak_pixels, // ditto but upper bound (count must be <=)
+ const unsigned max_peaks, // for buffer reasons
+ const float* image, // input data
unsigned int* output_counts,
float* output_x,
float* output_y,
@@ -150,29 +155,28 @@ extern "C" __global__ void pf9(const unsigned short num_frames,
float peak_weighted_ss;
float peak_weighted_fs;
float peak_total_mass = 0;
+ unsigned num_peak_pixels = 1; // will expand from center pixel, count it already
{
- /* TODO: more compact form */
float peak_weighted_ss_nom = 0;
float peak_weighted_fs_nom = 0;
// expand peak in rings
- {
- const float peak_pixel_threshold = mean + min_snr_peak_pixels * sigma;
- for (unsigned short layer=1; layer<=window_radius; ++layer) {
- float total_mass_before = peak_total_mass;
- masked_frame.fun_ring(ss_index, fs_index, layer, [&] (unsigned short i, unsigned short j, float val) {
- if (val > peak_pixel_threshold) {
- float val_over_mean = val - mean;
- peak_total_mass += val_over_mean;
- peak_weighted_ss_nom += val_over_mean * static_cast<float>(i);
- peak_weighted_fs_nom += val_over_mean * static_cast<float>(j);
- }
- });
- // in case nothing was added, stop expanding
- if (peak_total_mass == total_mass_before) {
- break;
- }
- }
- }
+ const float peak_pixel_threshold = mean + min_snr_peak_pixels * sigma;
+ for (unsigned short layer=1; layer<=window_radius; ++layer) {
+ unsigned num_peak_pixels_before = num_peak_pixels;
+ masked_frame.fun_ring(ss_index, fs_index, layer, [&] (unsigned short i, unsigned short j, float val) {
+ if (val > peak_pixel_threshold) {
+ float val_over_mean = val - mean;
+ peak_total_mass += val_over_mean;
+ peak_weighted_ss_nom += val_over_mean * static_cast<float>(i);
+ peak_weighted_fs_nom += val_over_mean * static_cast<float>(j);
+ ++num_peak_pixels;
+ }
+ });
+ // in case nothing was added, stop expanding
+ if (num_peak_pixels == num_peak_pixels_before) {
+ break;
+ }
+ }
if (peak_total_mass <= mean + min_snr_whole_peak * sigma) {
return;
}