RGB to YCbCr using SIMD vectors lose some data

Question

I'm writing JPEG decoder/encoder in Rust and I have some problem with RGB ↔ YCbCr conversion.

My code:

use std::simd::f32x4;

fn clamp<T>(val: T, min: T, max: T) -> T
where T: PartialOrd {
    if val < min { min }
    else if max < val { max }
    else { val }
}

// in oryginal code there are 2 methods, one for processors with SSE3 and for rest
// both do the same and give the same results
pub fn sum_f32x4(f32x4(a, b, c, d): f32x4) -> f32 {
    a + b + c + d
}

pub fn rgb_to_ycbcr(r: u8, g: u8, b: u8) -> (u8, u8, u8) {
    let rgb = f32x4(r as f32, g as f32, b as f32, 1.0);
    let y  = sum_f32x4(rgb * f32x4( 0.2990,  0.5870,  0.1140,   0.0));
    let cb = sum_f32x4(rgb * f32x4(-0.1687, -0.3313,  0.5000, 128.0));
    let cr = sum_f32x4(rgb * f32x4( 0.5000, -0.4187, -0.0813, 128.0));

    (y as u8, cb as u8, cr as u8)
}

pub fn ycbcr_to_rgb(y: u8, cb: u8, cr: u8) -> (u8, u8, u8) {
    let ycbcr = f32x4(y as f32, cb as f32 - 128.0f32, cr as f32 - 128.0f32, 0.0);
    let r = sum_f32x4(ycbcr * f32x4(1.0,  0.00000,  1.40200, 0.0));
    let g = sum_f32x4(ycbcr * f32x4(1.0, -0.34414, -0.71414, 0.0));
    let b = sum_f32x4(ycbcr * f32x4(1.0,  1.77200,  0.00000, 0.0));

    (clamp(r, 0., 255.) as u8, clamp(g, 0., 255.) as u8, clamp(b, 0., 255.) as u8)
}

fn main() {
    assert_eq!(rgb_to_ycbcr(  0,  71, 171), ( 61, 189,  84));
    // assert_eq!(rgb_to_ycbcr(  0,  71, 169), ( 61, 189,  84)); // will fail
    // for some reason we always lose data on blue channel
    assert_eq!(ycbcr_to_rgb( 61, 189,  84), (  0,  71, 169));
}

For some reason booth tests (in comments) passes. I would rather expect that at least one of them will fail. Am I wrong? At least it should stop at some point, but when I change jpeg::color::utils::rgb_to_ycbcr(0, 71, 171) to jpeg::color::utils::rgb_to_ycbcr(0, 71, 169) then test fails as YCbCr value has changed, so I will lose my blue channel forever.

What do you mean by "lose your blue channel"? And why are you expecting changing the input to give the same output for `rgb_to_ycbcr`? — huon, Jan 22 '15 at 00:05
When I convert from RGB `( 0, 71, 171)` to YCbCr `( 61, 189, 84)` and again to RGB `( 0, 71, 169)` the blue channel value is different, and no matter what values I set they are always different (I haven't found "stop" value at which it doesn't happen). — Hauleth, Jan 22 '15 at 00:09
I think this is either an algorithmic or a floating point problem. If I just use normal floats and no SIMD, I [get the same answers](http://is.gd/B8fuj3). — Shepmaster, Jan 22 '15 at 00:09
Also, why do you believe that that RGB ->YCbCr -> RGB will have no losses or changes? Different color spaces can represent different colors, and we have a fixed amount of precision. — Shepmaster, Jan 22 '15 at 00:12
I do not believe that values will always be the same. I believe that at some point they should stop changing. — Hauleth, Jan 22 '15 at 00:14
Have you seen [this SO question about invertible YUV / RGB transforms](http://stackoverflow.com/questions/17892346/how-to-convert-rgb-yuv-rgb-both-ways) ? — Shepmaster, Jan 22 '15 at 00:26
I suspect this may be because `as u8` (like the equivalent cast in C) truncates, it doesn't round to the nearest. You may be interested in calling [`.round`](http://doc.rust-lang.org/nightly/std/num/trait.Float.html#tymethod.round) first. — huon, Jan 22 '15 at 00:36

score 4 · Accepted Answer · answered Jan 22 '15 at 00:42

@dbaupp put the nail in the coffin with the suggestion to use round:

#![allow(unstable)]

use std::simd::{f32x4};
use std::num::Float;

fn clamp(val: f32) -> u8 {
    if val < 0.0 { 0 }
    else if val > 255.0 { 255 }
    else { val.round() as u8 }
}

fn sum_f32x4(v: f32x4) -> f32 {
    v.0 + v.1 + v.2 + v.3
}

pub fn rgb_to_ycbcr((r, g, b): (u8, u8, u8)) -> (u8, u8, u8) {
    let rgb = f32x4(r as f32, g as f32, b as f32, 1.0);
    let y  = sum_f32x4(rgb * f32x4( 0.299000,  0.587000,  0.114000,   0.0));
    let cb = sum_f32x4(rgb * f32x4(-0.168736, -0.331264,  0.500000, 128.0));
    let cr = sum_f32x4(rgb * f32x4( 0.500000, -0.418688, -0.081312, 128.0));

    (clamp(y), clamp(cb), clamp(cr))
}

pub fn ycbcr_to_rgb((y, cb, cr): (u8, u8, u8)) -> (u8, u8, u8) {
    let ycbcr = f32x4(y as f32, cb as f32 - 128.0f32, cr as f32 - 128.0f32, 0.0);
    let r = sum_f32x4(ycbcr * f32x4(1.0,  0.00000,  1.40200, 0.0));
    let g = sum_f32x4(ycbcr * f32x4(1.0, -0.34414, -0.71414, 0.0));
    let b = sum_f32x4(ycbcr * f32x4(1.0,  1.77200,  0.00000, 0.0));

    (clamp(r), clamp(g), clamp(b))
}

fn main() {
    let mut rgb = (0, 71, 16);
    println!("{:?}", rgb);

    for _ in 0..100 {
        let yuv = rgb_to_ycbcr(rgb);
        rgb = ycbcr_to_rgb(yuv);

        println!("{:?}", rgb);
        }
}

Note that I also increased the precision of your values in rgb_to_ycbcr from the Wikipedia page. I also clamp in both functions, as well as calling round. Now the output is:

(0u8, 71u8, 16u8)
(1u8, 72u8, 16u8)
(1u8, 72u8, 16u8)

With the last value repeating for the entire loop.

RGB to YCbCr using SIMD vectors lose some data

1 Answers1