Fixes #278: LoG scale normalization

[gschivre]

Simply replace -σ by -σ^2 in multiLoG and update the docstring of BlobLoG.

[codecov]

Codecov Report

All modified and coverable lines are covered by tests ✅

Project coverage is 91.81%. Comparing base (3b3d88a) to head (c6e854c).
Report is 4 commits behind head on master.

Additional details and impacted files

@@           Coverage Diff           @@
##           master     #279   +/-   ##
=======================================
  Coverage   91.81%   91.81%           
=======================================
  Files          12       12           
  Lines        1662     1662           
=======================================
  Hits         1526     1526           
  Misses        136      136           

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[ViralBShah]

Would be nice for someone else in the project to chime in to get this merged.

[timholy]

I didn't work through this in detail, but are you sure this is correct? Gaussian densities have dimensions of the inverse of volume, and this seems like you might have changed this to inverse volume-squared?

[gschivre]

I didn't work through this in detail, but are you sure this is correct? Gaussian densities have dimensions of the inverse of volume, and this seems like you might have changed this to inverse volume-squared?

From my understanding, the amplitude response of the n-th derivative of a Gaussian scales with σ^n, so for the Laplacian of Gaussian (sum of second order derivative) it should be σ^2.

This can be checked:

using ImageFiltering, LinearAlgebra

# fake image with a single ball of radius r
function NDball(r, N = 2; size = ntuple(d -> 64, N), center = CartesianIndex(ntuple(d -> 32, N)))
    return [1.0 - (norm(c.I .- center.I) >= r) for c in CartesianIndices(size)]
end

# get the blob of highest amplitude
function getTopBlob(img, scales)
    blobs = blob_LoG(img, scales)
    return first(sort(blobs, by = b -> b.amplitude, rev = true))
end

getTopBlob(NDball(10, 1), 1:20)
#= With σ scalling
BlobLoG(location=CartesianIndex(32,), σ=(7,), amplitude=0.06165173598954085)
=#
#= With σ^2 scalling
BlobLoG(location=CartesianIndex(32,), σ=(10,), amplitude=0.4831344763518372)
=#

getTopBlob(NDball(10, 2), 1:20)
#= With σ scalling
BlobLoG(location=CartesianIndex(36, 24), σ=(1, 1), amplitude=0.3046499888904334)
=#
#= With σ^2 scalling
BlobLoG(location=CartesianIndex(32, 32), σ=(7, 7), amplitude=0.7355903572521065)
=#

getTopBlob(NDball(10, 3), 1:20)
#= With σ scalling
BlobLoG(location=CartesianIndex(40, 36, 32), σ=(1, 1, 1), amplitude=0.33406542371161835)
=#
#= With σ^2 scalling
BlobLoG(location=CartesianIndex(32, 32, 32), σ=(6, 6, 6), amplitude=0.9201386216134283)
=#

[getTopBlob(NDball(r, 1), 1:20).amplitude for r in 2:2:14]
#= With σ scalling
7-element Vector{Float64}:
 0.3989422804014327
 0.15830978051758404
 0.10687862635162675
 0.07780732143311218
 0.06165173598954085
 0.051051348521593064
 0.043324744637534995
=#
#= With σ^2 scalling
7-element Vector{Float64}:
 0.46352013527394104
 0.47888456067216945
 0.48169786073672155
 0.48268019904506626
 0.4831344763518372
 0.48331328709165605
 0.4816113967778672
=#

In all cases, the selected scale is wrong. Additionally, for a 1-D image, the amplitude decays with the radius with the σ scaling, but not with the σ^2 scaling.

This is also the same scaling used by scikit-image.