Merge branch 'main' into 493-get_new23-final

src/IceFloeTracker.jl

@@ -72,12 +72,10 @@ include("special_strels.jl")
 include("tilingutils.jl")
 include("histogram_equalization.jl")
 include("new2-new3-final.jl")
+include("brighten.jl")
 include("morph_fill.jl")
 include("imcomplement.jl")
-
-const sk_measure = PyNULL()
-const sk_exposure = PyNULL()
-const getlatlon = PyNULL()
+include("imadjust.jl")
 
 function get_version_from_toml(pth=dirname(dirname(pathof(IceFloeTracker))))::VersionNumber
     toml = TOML.parsefile(joinpath(pth, "Project.toml"))
@@ -86,9 +84,16 @@ end
 
 const IFTVERSION = get_version_from_toml()
 
+const sk_measure = PyNULL()
+const sk_morphology = PyNULL()
+const sk_exposure = PyNULL()
+const getlatlon = PyNULL()
+
 function __init__()
-    copy!(sk_measure, pyimport_conda("skimage.measure", "scikit-image=0.24.0"))
-    copy!(sk_exposure, pyimport_conda("skimage.exposure", "scikit-image=0.24.0"))
+    skimage = "scikit-image=0.24.0"
+    copy!(sk_measure, pyimport_conda("skimage.measure", skimage))
+    copy!(sk_exposure, pyimport_conda("skimage.exposure", skimage))
+    copy!(sk_morphology, pyimport_conda("skimage.morphology", skimage))
     pyimport_conda("pyproj", "pyproj=3.6.0")
     pyimport_conda("rasterio", "rasterio=1.3.7")
     pyimport_conda("jinja2", "jinja2=3.1.2")

src/brighten.jl

@@ -0,0 +1,28 @@
+"""
+    get_brighten_mask(equalized_gray_reconstructed_img, gamma_green)
+# Arguments
+- `equalized_gray_reconstructed_img`: The equalized gray reconstructed image (uint8 in Matlab).
+- `gamma_green`: The gamma value for the green channel (also uint8).
+# Returns
+Difference equalized_gray_reconstructed_img - gamma_green clamped between 0 and 255.
+"""
+function get_brighten_mask(equalized_gray_reconstructed_img, gamma_green)
+    return to_uint8(equalized_gray_reconstructed_img - gamma_green)
+end
+
+"""
+    imbrighten(img, brighten_mask, bright_factor)
+Brighten the image using a mask and a brightening factor.
+# Arguments
+- `img`: The input image.
+- `brighten_mask`: A mask indicating the pixels to brighten.
+- `bright_factor`: The factor by which to brighten the pixels.
+# Returns
+- The brightened image.
+"""
+function imbrighten(img, brighten_mask, bright_factor)
+    img = Float64.(img)
+    brighten_mask = brighten_mask .> 0
+    img[brighten_mask] .= img[brighten_mask] * bright_factor
+    return img = to_uint8(img)
+end

src/histogram_equalization.jl

@@ -4,6 +4,11 @@ function to_uint8(arr::AbstractMatrix{T}) where {T<:AbstractFloat}
     return img
 end
 
+function to_uint8(arr::AbstractMatrix{T}) where {T<:Integer}
+    img = clamp.(arr, 0, 255)
+    return img
+end
+
 function anisotropic_diffusion_3D(I)
     rgbchannels = get_rgb_channels(I)
 

src/imadjust.jl

@@ -0,0 +1,27 @@
+# TODO: make imadjust more general to work with Gray{Float64} and RGB{Float64} types
+"""
+    imadjust(img; low, high)
+
+Adjust the contrast of an image using linear stretching. The image is normalized to [0, 1] and then stretched to the range [low, high].
+
+# Arguments
+- `img`: The input image.
+- `low`: The lower bound of the stretched image. Default is 0.01.
+- `high`: The upper bound of the stretched image. Default is 0.99.
+
+# Returns
+
+The contrast-adjusted image in the range [0, 255].
+"""
+function imadjust(
+    img::AbstractArray{<:Integer}; low::T=0.01, high::T=0.99
+)::Matrix{Int} where {T<:AbstractFloat}
+    img = img ./ 255
+    imgflat = vec(img)
+    plow = StatsBase.percentile(imgflat, low * 100)
+    phigh = StatsBase.percentile(imgflat, high * 100)
+
+    f = LinearStretching((plow, phigh) => (0.0, 1.0))
+
+    return to_uint8(adjust_histogram(img, f) * 255)
+end

src/segmentation_a_direct.jl

@@ -109,7 +109,19 @@ function get_holes(img, min_opening_area=20, se=IceFloeTracker.se_disk4())
     out = IceFloeTracker.MorphSE.fill_holes(out)
 
     return out .!= img
+end
+
+function fillholes!(img)
+    img[get_holes(img)] .= true
+    return nothing
+end
 
+function get_segment_mask(ice_mask, tiled_binmask)
+    Threads.@threads for img in (ice_mask, tiled_binmask)
+        fillholes!(img)
+    end
+    segment_mask = ice_mask .&& tiled_binmask
+    return segment_mask
 end
 
 function branchbridge(img)

src/utils.jl

@@ -120,6 +120,38 @@ function imextendedmin(img::AbstractArray; h::Int=2, conn::Int=2)::BitMatrix
     return Bool.(mask_minima)
 end
 
+"""
+    impose_minima(I::AbstractArray{T}, BW::AbstractArray{Bool}) where {T<:Integer}
+
+Use morphological reconstruction to enforce minima on the input image `I` at the positions where the binary mask `BW` is non-zero.
+
+It supports both integer and grayscale images using different implementations for each.
+"""
+function impose_minima(I::AbstractArray{T}, BW::AbstractArray{Bool}) where {T<:Integer}
+    marker = 255 .* BW
+    mask = imcomplement(min.(I .+ 1, 255 .- marker))
+    reconstructed = sk_morphology.reconstruction(marker, mask)
+    return IceFloeTracker.imcomplement(Int.(reconstructed))
+end
+
+function impose_minima(
+    I::AbstractArray{T}, BW::AbstractMatrix{Bool}
+) where {T<:AbstractFloat}
+    # compute shift
+    a, b = extrema(I)
+    rng = b - a
+    h = rng == 0 ? 0.1 : rng / 1000
+
+    marker = -Inf * BW .+ (Inf * .!BW)
+    mask = min.(I .+ h, marker)
+
+    return 1 .- sk_morphology.reconstruction(1 .- marker, 1 .- mask)
+end
+
+function imregionalmin(A, conn=2)
+    return ImageMorphology.local_minima(A; connectivity=conn) .> 0
+end
+
 """
     bwdist(bwimg)
 
@@ -142,11 +174,11 @@ function padnhood(img, I, nhood)
     tofill = SizedMatrix{3,3}(zeros(Int, 3, 3))
     @views if I == CartesianIndex(1, 1) # top left corner`
         tofill[2:3, 2:3] = img[nhood]
-    elseif I == CartesianIndex(maxr, 1) # bottom left corner 
+    elseif I == CartesianIndex(maxr, 1) # bottom left corner
         tofill[1:2, 2:3] = img[nhood]
-    elseif I == CartesianIndex(1, maxc) # top right corner 
+    elseif I == CartesianIndex(1, maxc) # top right corner
         tofill[2:3, 1:2] = img[nhood]
-    elseif I == CartesianIndex(maxr, maxc) # bottom right corner 
+    elseif I == CartesianIndex(maxr, maxc) # bottom right corner
         tofill[1:2, 1:2] = img[nhood]
     elseif I[1] == 1 # top edge (first row)
         tofill[2:3, 1:3] = img[nhood]
@@ -169,10 +201,10 @@ Get decimal representation of a bit vector `v` with the leading bit at its leftm
 
 Example
 ```
-julia> _bin9todec([0 0 0 0 0 0 0 0 0])    
+julia> _bin9todec([0 0 0 0 0 0 0 0 0])
 0
 
-julia> _bin9todec([1 1 1 1 1 1 1 1 1])    
+julia> _bin9todec([1 1 1 1 1 1 1 1 1])
 511
 ```
 """

test/runtests.jl

@@ -9,6 +9,7 @@ using Random
 using Test
 using TestImages
 using TiledIteration
+using ZipFile
 include("test_error_rate.jl")
 include("config.jl")
 

test/test-brighten.jl

@@ -0,0 +1,22 @@
+using IceFloeTracker: get_brighten_mask, imbrighten
+
+@testset "brighten tests" begin
+    @testset "get_brighten_mask" begin
+        img = rand(0:255, 5, 5)
+        bumped_img = img .+ 1
+        mask = get_brighten_mask(img, bumped_img)
+        @test all(mask .== 0)
+    end
+
+    @testset "imbrighten tests" begin
+        img = [1 2; 3 4]
+        brighten_mask = [1 0; 1 0]
+
+        test_cases = [(1.25, [1 2; 4 4]), (0.1, [0 2; 0 4]), (0.9, img)]
+
+        for (bright_factor, expected_result) in test_cases
+            result = imbrighten(img, brighten_mask, bright_factor)
+            @test result == expected_result
+        end
+    end
+end

test/test-imadjust.jl

@@ -0,0 +1,5 @@
+@testset "imadjust" begin
+    Random.seed!(123)
+    img = rand(0:255, 100, 100)
+    @test sum(IceFloeTracker.imadjust(img)) == 1291155
+end

test/test-img-processing.jl

@@ -1,6 +1,10 @@
-using ZipFile
-using DelimitedFiles: readdlm
-using IceFloeTracker: imgradientmag, to_uint8, imbinarize, adjustgamma, get_holes, se_disk4
+using IceFloeTracker:
+    imgradientmag,
+    to_uint8,
+    imbinarize,
+    adjustgamma,
+    get_holes,
+    impose_minima
 
 @testset "misc. image processing" begin
     r = ZipFile.Reader("test_inputs/coins.zip")
@@ -24,4 +28,11 @@ using IceFloeTracker: imgradientmag, to_uint8, imbinarize, adjustgamma, get_hole
         bw = coins .> 100
         @test sum(get_holes(bw)) == 2536
     end
+
+    @testset "impose_minima" begin
+        img = readdlm("test_inputs/imposemin.csv", ',', Int)
+        marker = falses(size(img))
+        marker[65:70, 65:70] .= true
+        @test sum(impose_minima(img, marker)) == 7675653
+    end
 end