Problems with conversion testing

[kimikage]

While working on the issue #354, I reviewed the conversion tests and found two problems which are not essential in the discussion #354.

One is that there are duplicate test cases. For example, @test typeof(c) == Cfrom{Tfrom} is evaluated for each Cto{Tto}.

Colors.jl/test/conversion.jl

Lines 84 to 97 in d64fb1a

	for Cto in ColorTypes.parametric3
	for Cfrom in ColorTypes.parametric3
	for Tto in (Float32, Float64)
	for Tfrom in (Float32, Float64)
	c = convert(Cfrom{Tfrom}, redF64)
	@test typeof(c) == Cfrom{Tfrom}
	c1 = convert(Cto, c)
	@test eltype(c1) == Tfrom
	c2 = convert(Cto{Tto}, c)
	@test typeof(c2) == Cto{Tto}
	end
	end
	end
	end

Unfortunately, the modification have almost no speedup effect. However, I think we should fix it because the number of duplicates is too large.

The other is that the accuracy tests are not fair.

Colors.jl/test/conversion.jl

Lines 232 to 238 in d64fb1a

	diff = abs(comp1(t)-comp1(f)) + abs(comp2(t)-comp2(f)) + abs(comp3(t)-comp3(f))
	mag = abs(comp1(t)+comp1(f)) + abs(comp2(t)+comp2(f)) + abs(comp3(t)+comp3(f))
	if showfailure && diff>eps*mag
	original = from[i]
	@show original f t
	end
	errmax = max(errmax, diff/mag)

For example, the HSV has three components: Hue in [0,360), Saturation in [0,1] and Value in [0,1].

julia> function diff(t, f)
           diff = abs(comp1(t)-comp1(f)) + abs(comp2(t)-comp2(f)) + abs(comp3(t)-comp3(f))
           mag = abs(comp1(t)+comp1(f)) + abs(comp2(t)+comp2(f)) + abs(comp3(t)+comp3(f))
           @show diff
           @show mag
           @show diff/mag
       end;

julia> diff(HSV(1e-3,1e-3,0.0), HSV(0,1e-3,0));
diff = 0.001
mag = 0.003
diff / mag = 0.3333333333333333

julia> diff(HSV(1e-3,1e-3,0.0), HSV(1e-3,0,0));
diff = 0.001
mag = 0.003
diff / mag = 0.3333333333333333

julia> diff(HSV(359.998,1,0), HSV(359.999,0,0));
diff = 1.0010000000000332
mag = 720.9970000000001
diff / mag = 0.0013883552913535467

So, although the difference of weights between the components may be ignorable, the normalization with mag is unreasonable. Of course, the magnitude of the mantissa LSB depends on the magnitude of the value (i.e. exponent), but the tolerance 1e-3 is enough larger for Float64.

Then, I wrote new methods for the evaluation.:

# Since `colordiff`(e.g. `DE_2000`) involves a color space conversion, it is
# not suitable for evaluating the conversion itself. On the other hand,
# since the tolerance varies from component to component, a homogeneous
# error evaluation function (e.g. a simple sum of differences) is also not
# appropriate. Therefore, a series of `diffnorm`, which returns the
# normalized Euclidian distance, is defined as follows. They are just for
# testing purposes as the cyclicity of hue is ignored.
sqd(a, b, s=1.0) = ((float(a) - float(b))/s)^2
function diffnorm(a::T, b::T) where {T<:Color3} # RGB,XYZ,xyY,LMS
    sqrt(sqd(comp1(a), comp1(b)) + sqd(comp2(a), comp2(b)) + sqd(comp3(a),comp3(b)))/sqrt(3)
end
function diffnorm(a::T, b::T) where {T<:Union{HSV,HSL,HSI}}
    sqrt(sqd(a.h, b.h, 360) + sqd(a.s, b.s) + sqd(comp3(a), comp3(b)))/sqrt(3)
end
function diffnorm(a::T, b::T) where {T<:Union{Lab,Luv}}
    sqrt(sqd(a.l, b.l, 100) + sqd(comp2(a), comp2(b), 200) + sqd(comp3(a), comp3(b), 200))/sqrt(3)
end
function diffnorm(a::T, b::T) where {T<:Union{LCHab,LCHuv}}
    sqrt(sqd(a.l, b.l, 100) + sqd(a.c, b.c, 100) + sqd(a.h, b.h, 360))/sqrt(3)
end
function diffnorm(a::T, b::T) where {T<:Union{DIN99,DIN99d,DIN99o}}
    sqrt(sqd(a.l, b.l, 100) + sqd(a.a, b.a, 100) + sqd(a.b, b.b, 100))/sqrt(3)
end
function diffnorm(a::T, b::T) where {T<:YIQ}
    sqrt(sqd(a.y, b.y) + sqd(a.i, b.i, 1.2) + sqd(a.q, b.q, 1.2))/sqrt(3)
end
function diffnorm(a::T, b::T) where {T<:YCbCr}
    sqrt(sqd(a.y, b.y, 219) + sqd(a.cb, b.cb, 224) + sqd(a.cr, b.cr, 224))/sqrt(3)
end

(To my regret, they make the tests slightly slower.)

However, in some cases of the conversion to HSI or HSL, the tests failed where the tolerance is 1e-3. (Needless to say, there is not any deep meaning in the number 1e-3)

original = Luv{Float64}(52.41763760855498,152.259323372186,15.212166995935663)
f = HSI{Float64}(343.88180452813356,0.8692773480658311,0.44676933798519713)
t = HSI{Float64}(343.93235717256283,0.8671370103299603,0.44710191307453834)
original = Luv{Float64}(50.02974879505459,157.6413223995605,26.989718351905513)
f = HSI{Float64}(349.5903318751271,0.9823710753987525,0.3749258760832548)
t = HSI{Float64}(349.6518042121588,0.9791646598774862,0.37533793657077363)
┌ Warning: Error on conversions from Luv{Float64} to HSI{Float64}, relative error = 0.00186905075014632

I suppose this is a problem with the transit (i.e. XYZ-->RGB), not HSL or HSI.

[kimikage]

I think the values in "test_conversions.jld2" have much greater errors than eps(Float64).
Where did the values ​​come from? If they were from the independent tools, this suggests that the tools used different coefficients (e.g. different fractional truncation, different white point definition).

By the way, is there any reason why "test_conversions.jld2" holds the from-to pairs even though it deals with the specific 20 colors only?

[kimikage]

Incidentally, although it does not matter in local environment, in the CI (i.e. Travis or AppVeyor), JLD2 (or CodecZlib) may lead to an increase of the build time. Is there any previous study on this?

Edit:
Using JLD (JLD2) is nice idea and it is a potentially powerful tool. However, the parser of Julia is not so slow now.
kimikage@7af3cb0
https://travis-ci.org/kimikage/Colors.jl/builds/590690586

Edit2:
If we stop using "test/test_conversions.jld2", I would recommend re-selecting the 20 colors. I prefer
intentionally fixed colors over random colors. For example:

Color3{Float64}[
    RGB(0,0,0) RGB(1,0,0) RGB(0,1,0) RGB(0,0,1)
    RGB(1,1,1) RGB(1,1,0) RGB(0,1,1) RGB(1,0,1)
    XYZ(0.5,0.5,0.5) LMS(0.4,0.5,0.6) RGB(0.6,0.6,0.6) YIQ(0.5,0,0)
    Lab(60,0,-20) Lab(50,-30,0) Lab(40,0,40) Lab(30,50,0)   
    Luv(25,-10,-10) Luv(75,60,60) YIQ(0.5,0.5,0) YCbCr(100,128,240)
]

[kimikage]

If they were from the independent tools, this suggests that the tools used different coefficients (e.g. different fractional truncation, different white point definition).

XYZ<-->(linear)sRGB transformation matrices can be calculated with the following method:

function mat(wp::XYZ)
    # Primaries in xyz (not XYZ)
    prim_r = BigFloat[6400, 3300,  300] / 10000
    prim_g = BigFloat[3000, 6000, 1000] / 10000
    prim_b = BigFloat[1500,  600, 7900] / 10000
    
    mat_prim = hcat(prim_r, prim_g, prim_b)
    
    s = mat_prim^-1 * [wp.x, wp.y, wp.z]
    
    mat_RGB2XYZ = mat_prim * Diagonal(s)
    mat_XYZ2RGB = mat_RGB2XYZ^-1
    
    (mat_RGB2XYZ, mat_XYZ2RGB)
end

The sRGB's primaries are based on ITU-R Rec. BT.709. The original publication of sRGB (i.e. working draft of IEC 61966-2-1) uses the BT.709's D65 white point, which is defined (rounded) in xy space. However, IEC 61966-2-1 uses the D65 white point which is defined in XYZ space. The two white point definitions are slightly different.

julia> wp_rec709 = convert(XYZ, xyY(BigFloat("0.3127"), BigFloat("0.3290"), 1));

julia> XYZ{Float64}(wp_rec709)
XYZ{Float64}(0.9504559270516717,1.0,1.0890577507598784)

julia> wp_d65 = XYZ{BigFloat}(BigFloat("0.95047"), 1, BigFloat("1.08883"));

julia> Colors.WP_D65
XYZ{Float64}(0.95047,1.0,1.08883)

julia> convert(xyY, Colors.WP_D65)
xyY{Float64}(0.3127266146810121,0.32902313032606195,1.0)

Therefore, there are roughly two types of the transformation matrices.

julia> rgb2xyz, xyz2rgb = mat(wp_rec709);

julia> show(IOContext(stdout, :compact=>false), "text/plain", Float64.(xyz2rgb))
3×3 Array{Float64,2}:
  3.2409699419045213   -1.5373831775700935   -0.4986107602930033
 -0.9692436362808798    1.8759675015077206    0.04155505740717561
  0.05563007969699361  -0.20397695888897657   1.0569715142428786
julia> show(IOContext(stdout, :compact=>false), "text/plain", Float64.(rgb2xyz))
3×3 Array{Float64,2}:
 0.4123907992659595   0.35758433938387796  0.1804807884018343
 0.21263900587151036  0.7151686787677559   0.07219231536073371
 0.01933081871559185  0.11919477979462599  0.9505321522496606

julia> rgb2xyz, xyz2rgb = mat(wp_d65);

julia> show(IOContext(stdout, :compact=>false), "text/plain", Float64.(xyz2rgb))
3×3 Array{Float64,2}:
  3.2404541621141054   -1.5371385127977166   -0.4985314095560162
 -0.9692660305051868    1.8760108454466942    0.04155601753034984
  0.05564343095911469  -0.20402591351675387   1.0572251882231791
julia> show(IOContext(stdout, :compact=>false), "text/plain", Float64.(rgb2xyz))
3×3 Array{Float64,2}:
 0.4124564390896921    0.357576077643909  0.18043748326639894
 0.21267285140562248   0.715152155287818  0.07217499330655958
 0.019333895582329317  0.119192025881303  0.9503040785363677

The conversion matrices in IEC 61966-2-1 are the truncated version of the latter. (To reduce errors, the truncation seems not so simple.)

And now, let's see the Colors.jl's implementation.

Colors.jl/src/conversions.jl

Lines 142 to 145 in 43d93ff

	function cnvt(::Type{CV}, c::XYZ) where CV<:AbstractRGB
	r = 3.2404542*c.x - 1.5371385*c.y - 0.4985314*c.z
	g = -0.9692660*c.x + 1.8760108*c.y + 0.0415560*c.z
	b = 0.0556434*c.x - 0.2040259*c.y + 1.0572252*c.z

Colors.jl/src/conversions.jl

Lines 287 to 291 in 43d93ff

	function cnvt(::Type{XYZ{T}}, c::AbstractRGB) where T
	r, g, b = invert_srgb_compand(red(c)), invert_srgb_compand(green(c)), invert_srgb_compand(blue(c))
	XYZ{T}(0.4124564*r + 0.3575761*g + 0.1804375*b,
	0.2126729*r + 0.7151522*g + 0.0721750*b,
	0.0193339*r + 0.1191920*g + 0.9503041*b)

So, Colors.jl does not follow the IEC 61966-2-1's truncation. I think the digits come from Float32 precision (≈7). If so, I think there is no problem with changing the element type of the matrices to Float32 (e.g. 3.2404542f0). Alternatively, we can also round the elements to 16 significant digits as Float64, since currently (d64fb1a) there is not much difference in speed between Float32 and Float64 on x86_64 systems.

[timholy]

Where did the values ​​come from? If they were from the independent tools, this suggests that the tools used different coefficients (e.g. different fractional truncation, different white point definition).

I think that's correct. It's been so long I don't remember for sure, but I was able to find a comment here: #259 (comment). I would support transitioning to other strategies, I was just worried (partly because of the whole DIN* issue) that we could be wildly off base and wanted to make sure we were at least in the right ballpark.

[kimikage]

The two problems mentioned above were fixed. (Although Appveyor failed the build, it is not a coding mistake but a network problem.)
I would like to submit a new issue to discuss the new test suite for conversions. However, I don't have the definitive solution, and I have some unclosed issues. So, I'll do it later, and until then, I leave this issue open.
I will continue to discuss the sRGB matrices in issue #351 (or its PR).

[kimikage]

The issue #368 took over this.