Improve HSV/HSL/HSI conversions

This adds the clamping and hue normalization for sources to HSx-->RGB conversions. This also adds the clamping for destinations to HSI-->RGB conversion.
JuliaGraphics · Feb 28, 2020 · 87d1dc4 · 87d1dc4
1 parent d07286d
commit 87d1dc4
Show file tree

Hide file tree

Showing 3 changed files with 132 additions and 101 deletions.
diff --git a/src/conversions.jl b/src/conversions.jl
@@ -71,6 +71,7 @@ end
 # -----------------
 
 correct_gamut(c::CV) where {CV<:AbstractRGB} = CV(clamp01(red(c)), clamp01(green(c)), clamp01(blue(c)))
+correct_gamut(c::CV) where {T<:Union{N0f8,N0f16,N0f32}, CV<:AbstractRGB{T}} = c
 clamp01(v::T) where {T<:Fractional} = ifelse(v < zero(T), zero(T), ifelse(v > one(T), one(T), v))
 
 function srgb_compand(v::Fractional)
@@ -81,58 +82,65 @@ end
 
 cnvt(::Type{CV}, c::AbstractRGB) where {CV<:AbstractRGB} = CV(red(c), green(c), blue(c))
 
-function cnvt(::Type{CV}, c::HSV) where CV<:AbstractRGB
-    h = c.h / 60
-    i = floor(Int, h)
-    f = h - i
-    if i & 1 == 0
-        f = 1 - f
-    end
-    m = c.v * (1 - c.s)
-    n = c.v * (1 - c.s * f)
-    if i == 6 || i == 0; CV(c.v, n, m)
-    elseif i == 1;       CV(n, c.v, m)
-    elseif i == 2;       CV(m, c.v, n)
-    elseif i == 3;       CV(m, n, c.v)
-    elseif i == 4;       CV(n, m, c.v)
-    else;                CV(c.v, m, n)
-    end
-end
-
-function qtrans(u, v, hue)
-    hue = normalize_hue(hue)
-
-    if     hue < 60;  u + (v - u) * hue / 60
-    elseif hue < 180; v
-    elseif hue < 240; u + (v - u) * (240 - hue) / 60
-    else;             u
-    end
-end
-
-function cnvt(::Type{CV}, c::HSL) where CV<:AbstractRGB
-    v = c.l <= 0.5 ? c.l * (1 + c.s) : c.l + c.s - (c.l * c.s)
-    u = 2 * c.l - v
 
-    if c.s == 0; CV(c.l, c.l, c.l)
-    else;        CV(qtrans(u, v, c.h + 120),
-                    qtrans(u, v, c.h),
-                    qtrans(u, v, c.h - 120))
-    end
-end
-
-function cnvt(::Type{CV}, c::HSI) where CV<:AbstractRGB
-    h, s, i = normalize_hue(c.h), c.s, c.i
-    is = i*s
-    if h < 120
-        cosr = cosd(h) / cosd(60-h)
-        CV(i+is*cosr, i+is*(1-cosr), i-is)
-    elseif h < 240
-        cosr = cosd(h-120) / cosd(180-h)
-        CV(i-is, i+is*cosr, i+is*(1-cosr))
+function _hsx_to_rgb(im::UInt8, v, n, m)
+    r = ifelse((im & 0b100001) != 0x0, v, ifelse((im & 0b010010) != 0x0, n, m))
+    g = ifelse((im & 0b000110) != 0x0, v, ifelse((im & 0b001001) != 0x0, n, m))
+    b = ifelse((im & 0b011000) != 0x0, v, ifelse((im & 0b100100) != 0x0, n, m))
+    return (r, g, b)
+end
+function _hsx_to_rgb(im::UInt8, v::T, n::T, m::T) where T <:Union{Float16, Float32, Float64}
+    vu, nu, mu = reinterpret.(Unsigned, (v, n, m)) # prompt the compiler to use conditional moves
+    r = ifelse((im & 0b100001) != 0x0, vu, ifelse((im & 0b010010) != 0x0, nu, mu))
+    g = ifelse((im & 0b000110) != 0x0, vu, ifelse((im & 0b001001) != 0x0, nu, mu))
+    b = ifelse((im & 0b011000) != 0x0, vu, ifelse((im & 0b100100) != 0x0, nu, mu))
+    return reinterpret.(T, (r, g, b))
+end
+
+function cnvt(::Type{CV}, c::HSV) where {T, CV<:AbstractRGB{T}}
+    F = promote_type(T, eltype(c))
+    h, s, v = div60(F(c.h)), clamp01(F(c.s)), clamp01(F(c.v))
+    hi = unsafe_trunc(Int32, h) # instead of floor
+    i = h >= 0 ? hi : hi - 0x1
+    f = (i & 0x1 != 0) ? h - i : 1 - (h - i)
+    im = 0x1 << (mod6(UInt8, i) & 0x07)
+    m = v * (1 - s)
+    n = v * (1 - s * f)
+    r, g, b = _hsx_to_rgb(im, v, n, m)
+    T <: FixedPoint && typemax(T) >= 1 ? CV(r % T, g % T, b % T) : CV(r, g, b)
+end
+
+function cnvt(::Type{CV}, c::HSL) where {T, CV<:AbstractRGB{T}}
+    F = promote_type(T, eltype(c))
+    h, s, l = div60(F(c.h)), clamp01(F(c.s)), clamp01(F(c.l))
+    a = @fastmath min(l, 1 - l) * s
+    v = a + l
+    hi = unsafe_trunc(Int32, h) # instead of floor
+    i = h >= 0 ? hi : hi - 0x1
+    f = (i & 0x1 != 0) ? h - i : 1 - (h - i)
+    im = 0x1 << (mod6(UInt8, i) & 0x07)
+    m = v - 2 * a
+    n = v - 2 * a * f
+    r, g, b = _hsx_to_rgb(im, v, n, m)
+    T <: FixedPoint && typemax(T) >= 1 ? CV(r % T, g % T, b % T) : CV(r, g, b)
+end
+
+function cnvt(::Type{CV}, c::HSI) where {T, CV<:AbstractRGB{T}}
+    F = promote_type(T, eltype(c))
+    h, s, i = deg2rad(normalize_hue(F(c.h))), clamp01(F(c.s)), clamp01(F(c.i))
+    is = i * s
+    if h < F(2π/3)
+        @fastmath cosr = cos(h) / cos(F(π/3)-h)
+        r0, g0, b0 = muladd(is, cosr, i), muladd(is, 1-cosr, i), i - is
+    elseif h < F(4π/3)
+        @fastmath cosr = cos(h-F(2π/3)) / cos(F(π)-h)
+        r0, g0, b0 = i - is, muladd(is, cosr, i), muladd(is, 1-cosr, i)
     else
-        cosr = cosd(h-240) / cosd(300-h)
-        CV(i+is*(1-cosr), i-is, i+is*cosr)
+        @fastmath cosr = cos(h-F(4π/3)) / cos(F(5π/3)-h)
+        r0, g0, b0 = muladd(is, 1-cosr, i), i - is, muladd(is, cosr, i)
     end
+    r, g, b = min(r0, oneunit(F)), min(g0, oneunit(F)), min(b0, oneunit(F))
+    T <: FixedPoint && typemax(T) >= 1 ? CV(r % T, g % T, b % T) : CV(r, g, b)
 end
 
 function cnvt(::Type{CV}, c::XYZ) where CV<:AbstractRGB
@@ -175,26 +183,23 @@ end
 # -----------------
 
 function cnvt(::Type{HSV{T}}, c::AbstractRGB) where T
-    c_min = Float64(min(red(c), green(c), blue(c)))
-    c_max = Float64(max(red(c), green(c), blue(c)))
-    if c_min == c_max
-        return HSV{T}(zero(T), zero(T), c_max)
-    end
-
-    if c_min == red(c)
-        f = Float64(green(c)) - Float64(blue(c))
-        i = 3
-    elseif c_min == green(c)
-        f = Float64(blue(c)) - Float64(red(c))
-        i = 5
+    F = promote_type(T, eltype(c))
+    r, g, b = F.((red(c), green(c), blue(c)))
+    @fastmath c_min = min(min(r, g), b)
+    @fastmath c_max = max(max(r, g), b)
+    s0 = c_max - c_min
+    s0 == 0 && return HSV{T}(zero(T), zero(T), T(c_max))
+
+    if c_max == r
+        h0 = (g - b) / s0
+        h = (h0 < 0) * 6 + h0
+    elseif c_max == g
+        h = 2 + (b - r) / s0
     else
-        f = Float64(red(c)) - Float64(green(c))
-        i = 1
+        h = 4 + (r - g) / s0
     end
 
-    HSV{T}(60 * (i - f / (c_max - c_min)),
-        (c_max - c_min) / c_max,
-        c_max)
+    HSV{T}(h * 60, s0 / c_max, c_max)
 end
 
 
@@ -205,28 +210,26 @@ cnvt(::Type{HSV{T}}, c::Color3) where {T} = cnvt(HSV{T}, convert(RGB{T}, c))
 # -----------------
 
 function cnvt(::Type{HSL{T}}, c::AbstractRGB) where T
-    r, g, b = T(red(c)), T(green(c)), T(blue(c))
-    c_min = min(r, g, b)
-    c_max = max(r, g, b)
-    l = (c_max + c_min) / 2
-
-    if c_max == c_min
-        return HSL(zero(T), zero(T), l)
-    end
-
-    if l < 0.5; s = (c_max - c_min) / (c_max + c_min)
-    else;       s = (c_max - c_min) / (convert(T, 2) - c_max - c_min)
-    end
-
-    if c_max == red(c)
-        h = (g - b) / (c_max - c_min)
-    elseif c_max == green(c)
-        h = convert(T, 2) + (b - r) / (c_max - c_min)
+    F = promote_type(T, eltype(c))
+    r, g, b = F(red(c)), F(green(c)), F(blue(c))
+    c_min = @fastmath min(min(r, g), b)
+    c_max = @fastmath max(max(r, g), b)
+    l0 = c_max + c_min
+    s0 = c_max - c_min
+    s0 == 0 && return HSL{T}(zero(T), zero(T), l0 * T(0.5))
+
+    s = @fastmath s0 / min(l0, 2 - l0)
+
+    if c_max == r
+        h0 = (g - b) / s0
+        h = (h0 < 0) * 6 + h0
+    elseif c_max == g
+        h = 2 + (b - r) / s0
     else
-        h = convert(T, 4) + (r - g) / (c_max - c_min)
+        h = 4 + (r - g) / s0
     end
 
-    HSL{T}(normalize_hue(h * 60), s, l)
+    HSL{T}(h * 60, s, l0 * F(0.5))
 end
 
 
@@ -236,22 +239,20 @@ cnvt(::Type{HSL{T}}, c::Color3) where {T} = cnvt(HSL{T}, convert(RGB{T}, c))
 # Everything to HSI
 # -----------------
 
-function cnvt(::Type{HSI{T}}, c::AbstractRGB) where T
+# Since acosd() is slow, the following is "inline-worthy".
+@inline function cnvt(::Type{HSI{T}}, c::AbstractRGB) where T
     rgb = correct_gamut(c)
-    r, g, b = float(red(rgb)), float(green(rgb)), float(blue(rgb))
-    isum = r+g+b
-    dnorm = sqrt(((r-g)^2 + (r-b)^2 + (g-b)^2)/2)
-    dnorm = dnorm == 0 ? oftype(dnorm, 1) : dnorm
-    i = isum/3
-    m = min(r, g, b)
-    s = i > 0 ? 1-m/i : zero(1 - m/i)
-    val = (r-(g+b)/2)/dnorm
-    val = clamp(val, -oneunit(val), oneunit(val))
-    h = acosd(val)
-    if b > g
-        h = 360-h
-    end
-    HSI{T}(h, s, i)
+    F = promote_type(T, eltype(c))
+    r, g, b = F(red(rgb)), F(green(rgb)), F(blue(rgb))
+    dnorm = @fastmath sqrt(((r-g)^2 + (r-b)^2 + (g-b)^2) * F(0.5))
+    isum = r + g + b
+    i = isum / 3
+    dnorm == 0 && return HSI{T}(T(90), zero(T), T(i))
+    val = muladd(g + b, F(-0.5), r) / dnorm
+    h = @fastmath acosd(clamp(val, -oneunit(F), oneunit(F)))
+    m = @fastmath min(min(r, g), b)
+    s = 1 - m/i
+    HSI{T}(b > g ? 360 - h : h, s, i)
 end
 
 cnvt(::Type{HSI{T}}, c::Color3) where {T} = cnvt(HSI{T}, convert(RGB{T}, c))

diff --git a/src/utilities.jl b/src/utilities.jl
@@ -1,6 +1,12 @@
 # Helper data for CIE observer functions
 include("cie_data.jl")
 
+# for optimization
+div60(x) = x / 60
+div60(x::T) where T <: Union{Float32, Float64} = fma(x, T(1/960), x * T(0x1p-6))
+
+# mod6 supports the input `x` in [-2^28, 2^29]
+mod6(::Type{T}, x::Int32) where T = unsafe_trunc(T, x - 6 * ((x * Int64(0x2aaaaaaa) + 0x20000000) >> 0x20))
 
 # Linear interpolation in [a, b] where x is in [0,1],
 # or coerced to be if not.

diff --git a/test/conversion.jl b/test/conversion.jl
@@ -117,7 +117,7 @@ using ColorTypes: eltype_default, parametric3
 
     @test convert(RGB{N0f8}, red24) == RGB{N0f8}(1,0,0)
     @test convert(RGBA{N0f8}, red32) == RGBA{N0f8}(1,0,0,1)
-    @test convert(HSVA{Float64}, red32) == HSVA{Float64}(360, 1, 1, 1)
+    @test convert(HSVA{Float64}, red32) == HSVA{Float64}(0, 1, 1, 1)
 
     @test_throws MethodError AlphaColor(RGB(1,0,0), r8(0xff))
 
@@ -202,6 +202,30 @@ using ColorTypes: eltype_default, parametric3
     hsi = convert(HSI, c)
     @test hsi.i > 0.96 && hsi.h ≈ 210
 
+    # {HSV, HSL, HSI} --> RGB (issue #379)
+    @testset "HSx --> RGB" begin
+        @test convert(RGB, HSV{Float32}( 780, 1, 1)) === RGB{Float32}(1,1,0)
+        @test convert(RGB, HSV{Float32}(   0, 1, 1)) === RGB{Float32}(1,0,0)
+        @test convert(RGB, HSV{Float32}(-780, 1, 1)) === RGB{Float32}(1,0,1)
+        @test convert(RGB, HSV{Float32}(30,  2, .5)) === RGB{Float32}(.5,.25,0)
+        @test convert(RGB, HSV{Float32}(30, .5, -1)) === RGB{Float32}(0,0,0)
+        @test convert(RGB{Float64}, HSV{BigFloat}(-360120, 2, 1)) === RGB{Float64}(0,0,1)
+
+        @test convert(RGB, HSL{Float32}( 780, 1, .5)) === RGB{Float32}(1,1,0)
+        @test convert(RGB, HSL{Float32}(   0, 1, .5)) === RGB{Float32}(1,0,0)
+        @test convert(RGB, HSL{Float32}(-780, 1, .5)) === RGB{Float32}(1,0,1)
+        @test convert(RGB, HSL{Float32}(30,  2, .25)) === RGB{Float32}(.5,.25,0)
+        @test convert(RGB, HSL{Float32}(30, .5,  -1)) === RGB{Float32}(0,0,0)
+        @test convert(RGB{Float64}, HSL{BigFloat}(-360120, 2, .5)) === RGB{Float64}(0,0,1)
+
+        @test convert(RGB, HSI{Float32}( 780, .5, .5)) ≈ RGB{Float32}(.625,.625,.25)
+        @test convert(RGB, HSI{Float32}(   0, .5, .5)) ≈ RGB{Float32}(1,.25,.25)
+        @test convert(RGB, HSI{Float32}(-780, .5, .5)) ≈ RGB{Float32}(.625,.25,.625)
+        @test convert(RGB, HSI{Float32}(30,  2, .25)) ≈ RGB{Float32}(.5,.25,0)
+        @test convert(RGB, HSI{Float32}(30, .5,  -1)) ≈ RGB{Float32}(0,0,0)
+        @test convert(RGB{Float64}, HSI{BigFloat}(-360120, .5, .5)) ≈ RGB{Float64}(.25,.25,1)
+    end
+
     # Test accuracy of conversion
     include("test_conversions.jl")