colorutils_rs/
xyb.rs

1/*
2 * // Copyright 2024 (c) the Radzivon Bartoshyk. All rights reserved.
3 * //
4 * // Use of this source code is governed by a BSD-style
5 * // license that can be found in the LICENSE file.
6 */
7use crate::utils::mlaf;
8use crate::{EuclideanDistance, Rgb, TaxicabDistance, TransferFunction};
9use num_traits::Pow;
10use std::ops::{Add, AddAssign, Div, DivAssign, Mul, MulAssign, Neg, Sub, SubAssign};
11
12/// XYB is a color space that was designed for use with the JPEG XL Image Coding System.
13///
14/// It is an LMS-based color model inspired by the human visual system, facilitating perceptually uniform quantization.
15/// It uses a gamma of 3 for computationally efficient decoding.
16#[repr(C)]
17#[derive(Copy, Clone, Debug, PartialOrd, PartialEq)]
18pub struct Xyb {
19    pub x: f32,
20    pub y: f32,
21    pub b: f32,
22}
23
24impl Xyb {
25    #[inline]
26    pub fn new(x: f32, y: f32, b: f32) -> Xyb {
27        Xyb { x, y, b }
28    }
29
30    #[inline]
31    /// Converts [Rgb] to [Xyb] using provided [TransferFunction]
32    pub fn from_rgb(rgb: Rgb<u8>, transfer_function: TransferFunction) -> Xyb {
33        let linear_rgb = rgb.to_linear(transfer_function);
34        Self::from_linear_rgb(linear_rgb)
35    }
36
37    #[inline]
38    /// Converts linear [Rgb] to [Xyb]
39    pub fn from_linear_rgb(rgb: Rgb<f32>) -> Xyb {
40        const BIAS_CBRT: f32 = 0.155954200549248620f32;
41        const BIAS: f32 = 0.00379307325527544933;
42        let lgamma = mlaf(
43            0.3f32,
44            rgb.r,
45            mlaf(0.622f32, rgb.g, mlaf(0.078f32, rgb.b, BIAS)),
46        )
47        .cbrt()
48            - BIAS_CBRT;
49        let mgamma = mlaf(
50            0.23f32,
51            rgb.r,
52            mlaf(0.692f32, rgb.g, mlaf(0.078f32, rgb.b, BIAS)),
53        )
54        .cbrt()
55            - BIAS_CBRT;
56        let sgamma = mlaf(
57            0.24342268924547819f32,
58            rgb.r,
59            mlaf(
60                0.20476744424496821f32,
61                rgb.g,
62                mlaf(0.55180986650955360f32, rgb.b, BIAS),
63            ),
64        )
65        .cbrt()
66            - BIAS_CBRT;
67        let x = (lgamma - mgamma) * 0.5f32;
68        let y = (lgamma + mgamma) * 0.5f32;
69        let b = sgamma - mgamma;
70        Xyb::new(x, y, b)
71    }
72
73    #[inline]
74    /// Converts [Xyb] to linear [Rgb]
75    pub fn to_linear_rgb(&self) -> Rgb<f32> {
76        const BIAS_CBRT: f32 = 0.155954200549248620f32;
77        const BIAS: f32 = 0.00379307325527544933;
78        let x_lms = (self.x + self.y) + BIAS_CBRT;
79        let y_lms = (-self.x + self.y) + BIAS_CBRT;
80        let b_lms = (-self.x + self.y + self.b) + BIAS_CBRT;
81        let x_c_lms = (x_lms * x_lms * x_lms) - BIAS;
82        let y_c_lms = (y_lms * y_lms * y_lms) - BIAS;
83        let b_c_lms = (b_lms * b_lms * b_lms) - BIAS;
84        let r = mlaf(
85            11.031566901960783,
86            x_c_lms,
87            mlaf(-9.866943921568629, y_c_lms, -0.16462299647058826 * b_c_lms),
88        );
89        let g = mlaf(
90            -3.254147380392157,
91            x_c_lms,
92            mlaf(4.418770392156863, y_c_lms, -0.16462299647058826 * b_c_lms),
93        );
94        let b = mlaf(
95            -3.6588512862745097,
96            x_c_lms,
97            mlaf(2.7129230470588235, y_c_lms, 1.9459282392156863 * b_c_lms),
98        );
99        Rgb::new(r, g, b)
100    }
101
102    #[inline]
103    /// Converts [Xyb] to [Rgb] using provided [TransferFunction]
104    pub fn to_rgb(&self, transfer_function: TransferFunction) -> Rgb<u8> {
105        let linear_rgb = self.to_linear_rgb();
106        linear_rgb.gamma(transfer_function).to_u8()
107    }
108}
109
110impl Add<f32> for Xyb {
111    type Output = Xyb;
112
113    #[inline]
114    fn add(self, rhs: f32) -> Self::Output {
115        Xyb::new(self.x + rhs, self.y + rhs, self.b + rhs)
116    }
117}
118
119impl Add<Xyb> for Xyb {
120    type Output = Xyb;
121
122    #[inline]
123    fn add(self, rhs: Xyb) -> Self::Output {
124        Xyb::new(self.x + rhs.x, self.y + rhs.y, self.b + rhs.b)
125    }
126}
127
128impl Sub<f32> for Xyb {
129    type Output = Xyb;
130
131    #[inline]
132    fn sub(self, rhs: f32) -> Self::Output {
133        Xyb::new(self.x - rhs, self.y - rhs, self.b - rhs)
134    }
135}
136
137impl Sub<Xyb> for Xyb {
138    type Output = Xyb;
139
140    #[inline]
141    fn sub(self, rhs: Xyb) -> Self::Output {
142        Xyb::new(self.x - rhs.x, self.y - rhs.y, self.b - rhs.b)
143    }
144}
145
146impl Mul<f32> for Xyb {
147    type Output = Xyb;
148
149    #[inline]
150    fn mul(self, rhs: f32) -> Self::Output {
151        Xyb::new(self.x * rhs, self.y * rhs, self.b * rhs)
152    }
153}
154
155impl Mul<Xyb> for Xyb {
156    type Output = Xyb;
157
158    #[inline]
159    fn mul(self, rhs: Xyb) -> Self::Output {
160        Xyb::new(self.x * rhs.x, self.y * rhs.y, self.b * rhs.b)
161    }
162}
163
164impl Div<f32> for Xyb {
165    type Output = Xyb;
166
167    #[inline]
168    fn div(self, rhs: f32) -> Self::Output {
169        Xyb::new(self.x / rhs, self.y / rhs, self.b / rhs)
170    }
171}
172
173impl Div<Xyb> for Xyb {
174    type Output = Xyb;
175
176    #[inline]
177    fn div(self, rhs: Xyb) -> Self::Output {
178        Xyb::new(self.x / rhs.x, self.y / rhs.y, self.b / rhs.b)
179    }
180}
181
182impl Neg for Xyb {
183    type Output = Xyb;
184
185    #[inline]
186    fn neg(self) -> Self::Output {
187        Xyb::new(-self.x, -self.y, -self.b)
188    }
189}
190
191impl Pow<f32> for Xyb {
192    type Output = Xyb;
193
194    #[inline]
195    fn pow(self, rhs: f32) -> Self::Output {
196        Xyb::new(self.x.powf(rhs), self.y.powf(rhs), self.b.powf(rhs))
197    }
198}
199
200impl Pow<Xyb> for Xyb {
201    type Output = Xyb;
202
203    #[inline]
204    fn pow(self, rhs: Xyb) -> Self::Output {
205        Xyb::new(self.x.powf(rhs.x), self.y.powf(rhs.y), self.b.powf(rhs.b))
206    }
207}
208
209impl MulAssign<f32> for Xyb {
210    #[inline]
211    fn mul_assign(&mut self, rhs: f32) {
212        self.x *= rhs;
213        self.y *= rhs;
214        self.b *= rhs;
215    }
216}
217
218impl MulAssign<Xyb> for Xyb {
219    #[inline]
220    fn mul_assign(&mut self, rhs: Xyb) {
221        self.x *= rhs.x;
222        self.y *= rhs.y;
223        self.b *= rhs.b;
224    }
225}
226
227impl AddAssign<f32> for Xyb {
228    #[inline]
229    fn add_assign(&mut self, rhs: f32) {
230        self.x += rhs;
231        self.y += rhs;
232        self.b += rhs;
233    }
234}
235
236impl AddAssign<Xyb> for Xyb {
237    #[inline]
238    fn add_assign(&mut self, rhs: Xyb) {
239        self.x += rhs.x;
240        self.y += rhs.y;
241        self.b += rhs.b;
242    }
243}
244
245impl SubAssign<f32> for Xyb {
246    #[inline]
247    fn sub_assign(&mut self, rhs: f32) {
248        self.x -= rhs;
249        self.y -= rhs;
250        self.b -= rhs;
251    }
252}
253
254impl SubAssign<Xyb> for Xyb {
255    #[inline]
256    fn sub_assign(&mut self, rhs: Xyb) {
257        self.x -= rhs.x;
258        self.y -= rhs.y;
259        self.b -= rhs.b;
260    }
261}
262
263impl DivAssign<f32> for Xyb {
264    #[inline]
265    fn div_assign(&mut self, rhs: f32) {
266        self.x /= rhs;
267        self.y /= rhs;
268        self.b /= rhs;
269    }
270}
271
272impl DivAssign<Xyb> for Xyb {
273    #[inline]
274    fn div_assign(&mut self, rhs: Xyb) {
275        self.x /= rhs.x;
276        self.y /= rhs.y;
277        self.b /= rhs.b;
278    }
279}
280
281impl Xyb {
282    #[inline]
283    pub fn sqrt(&self) -> Xyb {
284        Xyb::new(
285            if self.x < 0. { 0. } else { self.x.sqrt() },
286            if self.y < 0. { 0. } else { self.y.sqrt() },
287            if self.b < 0. { 0. } else { self.b.sqrt() },
288        )
289    }
290
291    #[inline]
292    pub fn cbrt(&self) -> Xyb {
293        Xyb::new(self.x.cbrt(), self.y.cbrt(), self.b.cbrt())
294    }
295}
296
297impl EuclideanDistance for Xyb {
298    fn euclidean_distance(&self, other: Self) -> f32 {
299        let dx = self.x - other.x;
300        let dy = self.y - other.y;
301        let db = self.b - other.b;
302        (dx * dx + dy * dy + db * db).sqrt()
303    }
304}
305
306impl TaxicabDistance for Xyb {
307    fn taxicab_distance(&self, other: Self) -> f32 {
308        let dx = self.x - other.x;
309        let dy = self.y - other.y;
310        let db = self.b - other.b;
311        dx.abs() + dy.abs() + db.abs()
312    }
313}