Merge pull request #458 from mronian/add-hist-eq

Histogram Equalisation

docs/src/function_reference.md

@@ -282,6 +282,7 @@ findlocalminima
 
 ```@docs
 imhist
+histeq
 ```
 
 # Filtering kernels

src/Images.jl

@@ -217,6 +217,7 @@ export # types
     imgaussiannoise,
     imgradients,
     imhist,
+    histeq,
     imlaplacian,
     imlineardiffusion,
     imlog,
@@ -280,7 +281,7 @@ Algorithms:
     - Resizing: `restrict`, `imresize` (not yet exported)
     - Filtering: `imfilter`, `imfilter_fft`, `imfilter_gaussian`, `imfilter_LoG`, `imROF`, `ncc`, `padarray`
     - Filtering kernels: `ando[345]`, `guassian2d`, `imaverage`, `imdog`, `imlaplacian`, `prewitt`, `sobel`
-    - Exposure : `imhist`
+    - Exposure : `imhist`, `histeq`
     - Gradients: `backdiffx`, `backdiffy`, `forwarddiffx`, `forwarddiffy`, `imgradients`
     - Edge detection: `imedge`, `imgradients`, `thin_edges`, `magnitude`, `phase`, `magnitudephase`, `orientation`
     - Corner detection: `imcorner`

src/algorithms.jl

@@ -1239,6 +1239,105 @@ end
 
 imhist{T<:Union{Gray,Number}}(img::AbstractArray{T}, nbins, minval, maxval) = imhist(img, nbins, convert(T, minval), convert(T, maxval))
 
+"""
+```
+hist_equalised_img = histeq(img, nbins, dtype = "8bit")
+hist_equalised_img = histeq(img, nbins, minval, maxval, dtype = "8bit")
+```
+
+Returns a histogram equalised image with a granularity of approximately `nbins` 
+number of bins. An optional `dtype` argument (defaulting to 8bit) can be specified
+to choose the number of bits of the returned image. 
+
+The `histeq` function can handle a variety of input types. The returned image depends 
+on the input type. If the input is an `Image` then the resulting image is of the same type
+and has the same properties. 
+
+For coloured images, the input is converted to YCbCr type and the Y channel is equalised. This
+is the combined with the Cb and Cr channels and the resulting image converted to the same type 
+as the input.
+
+If minval and maxval are specified then only the intensities in the range 
+(minval, maxval) are equalised.
+
+"""
+
+RET_TYPE = Dict("8bit" => Gray{U8}, "16bit" => Gray{FixedPointNumbers.UFixed{UInt16,16}})
+Y_MIN = 16
+Y_MAX = 235
+
+function _prep_image_for_histeq(img::AbstractArray, dtype)
+    img_shape = size(img)
+    if dtype == "8bit"
+        img = [convert(base_colorant_type(c){FixedPointNumbers.UFixed{UInt8,8}}, c) for c in img]
+    elseif dtype == "16bit"
+        img = [convert(base_colorant_type(c){FixedPointNumbers.UFixed{UInt16,16}}, c) for c in img]
+    end
+    img = reshape(img, img_shape)
+    img
+end
+
+histeq{T<:Colorant}(img::AbstractArray{T}, nbins, dtype = "8bit") = histeq(_prep_image_for_histeq(img, dtype), nbins, zero(YCbCr), zero(YCbCr))
+histeq{T<:Colorant}(img::AbstractArray{T}, nbins, minval, maxval, dtype = "8bit") = histeq(_prep_image_for_histeq(img, dtype), nbins, minval, maxval)
+
+function recompose_y(c::Color, eq_Y::Float32)
+    c_ycbcr = convert(YCbCr, color(c))
+    ret_c = YCbCr(eq_Y, c_ycbcr.cb, c_ycbcr.cr)
+    ret_c
+end
+
+function recompose_y(c::TransparentColor, eq_Y::Float32)
+    c_ycbcr = convert(YCbCrA, color(c))
+    ret_c = YCbCrA(eq_Y, c_ycbcr.cb, c_ycbcr.cr, c_ycbcr.alpha)
+    ret_c
+end
+
+function histeq{T<:Colorant}(img::AbstractArray{T}, nbins, minval, maxval)
+    Y = [convert(YCbCr, color(c)).y for c in img]
+    if maxval == zero(YCbCr{Float32})
+        eq_Y = histeq(Y, nbins, Y_MIN, Y_MAX)
+    else
+        eq_Y = histeq(Y, nbins, minval, maxval)
+    end
+    hist_equalised_img = reshape([convert(T, recompose_y(c, eq_Y[i])) for (i, c) in enumerate(img)], size(img))
+    hist_equalised_img
+end
+
+function histeq(img::AbstractImage, nbins, dtype = "8bit")
+    hist_equalised_img = histeq(data(img), nbins, dtype)
+    hist_equalised_img = shareproperties(img, hist_equalised_img)
+    hist_equalised_img
+end
+
+function histeq(img::AbstractImage, nbins, minval, maxval, dtype = "8bit")
+    hist_equalised_img = histeq(data(img), nbins, minval, maxval, dtype)
+    hist_equalised_img = shareproperties(img, hist_equalised_img)
+    hist_equalised_img
+end
+
+function histeq{T<:TransparentGray}(img::AbstractArray{T}, args...)
+    opaque_img = [color(c) for c in img]
+    hist_equalised_img = histeq(opaque_img, args...)
+    hist_equalised_img = reshape([convert(T, AGray(hist_equalised_img[i], alpha(c))) for (i, c) in enumerate(img)], size(img))
+    hist_equalised_img    
+end
+
+histeq{T<:Gray}(img::AbstractArray{T}, nbins, dtype = "8bit") = histeq(_prep_image_for_histeq(img, dtype), nbins, ColorVectorSpace.typemin(RET_TYPE[dtype]), ColorVectorSpace.typemax(RET_TYPE[dtype]))
+
+function histeq{T<:Union{Gray,Number}}(img::AbstractArray{T}, nbins::Int, minval::T, maxval::T)
+    bins, histogram = imhist(img, nbins, minval, maxval)
+    cdf = cumsum(histogram[2:end-1])
+    img_shape = size(img)
+    cdf_length = size(cdf)[1]
+    hist_equalised_img = [max(1,Int(ceil((x-minval)*cdf_length/(maxval-minval)))) for x in img]
+    hist_equalised_img = [minval + ((cdf[x]-cdf[1])*(maxval-minval)/(cdf[end]-cdf[1])) for x in hist_equalised_img]
+    hist_equalised_img = reshape(hist_equalised_img, img_shape)
+    hist_equalised_img
+end
+
+histeq{T<:Number}(img::AbstractArray{T}, nbins) = histeq(img, nbins, minfinite(img), maxfinite(img))
+
+histeq{T<:Union{Gray,Number}}(img::AbstractArray{T}, nbins, minval, maxval) = histeq(img, Int(nbins), convert(T, minval), convert(T, maxval))
 
 """
 `imgr = restrict(img[, region])` performs two-fold reduction in size