Estimate a log-linear formula for the extracted transfer function LUT.

Cargo.lock

@@ -211,6 +211,7 @@ dependencies = [
  "clap",
  "colorbox",
  "exr",
+ "simplers_optimization",
 ]
 
 [[package]]
@@ -228,6 +229,15 @@ dependencies = [
  "getrandom",
 ]
 
+[[package]]
+name = "num-traits"
+version = "0.2.14"
+source = "registry+https://github.com/rust-lang/crates.io-index"
+checksum = "9a64b1ec5cda2586e284722486d802acf1f7dbdc623e2bfc57e65ca1cd099290"
+dependencies = [
+ "autocfg",
+]
+
 [[package]]
 name = "num_cpus"
 version = "1.13.1"
@@ -238,6 +248,15 @@ dependencies = [
  "libc",
 ]
 
+[[package]]
+name = "ordered-float"
+version = "2.10.0"
+source = "registry+https://github.com/rust-lang/crates.io-index"
+checksum = "7940cf2ca942593318d07fcf2596cdca60a85c9e7fab408a5e21a4f9dcd40d87"
+dependencies = [
+ "num-traits",
+]
+
 [[package]]
 name = "os_str_bytes"
 version = "6.0.0"
@@ -267,6 +286,16 @@ dependencies = [
  "syn",
 ]
 
+[[package]]
+name = "priority-queue"
+version = "1.2.1"
+source = "registry+https://github.com/rust-lang/crates.io-index"
+checksum = "00ba480ac08d3cfc40dea10fd466fd2c14dee3ea6fc7873bc4079eda2727caf0"
+dependencies = [
+ "autocfg",
+ "indexmap",
+]
+
 [[package]]
 name = "proc-macro2"
 version = "1.0.37"
@@ -291,6 +320,17 @@ version = "1.1.0"
 source = "registry+https://github.com/rust-lang/crates.io-index"
 checksum = "d29ab0c6d3fc0ee92fe66e2d99f700eab17a8d57d1c1d3b748380fb20baa78cd"
 
+[[package]]
+name = "simplers_optimization"
+version = "0.4.3"
+source = "registry+https://github.com/rust-lang/crates.io-index"
+checksum = "2cd97912bb2a16575a2c632c2a2f2bac8a527827ceaddd73e0ccc12d86adec43"
+dependencies = [
+ "num-traits",
+ "ordered-float",
+ "priority-queue",
+]
+
 [[package]]
 name = "smallvec"
 version = "1.8.0"

Cargo.toml

@@ -7,3 +7,4 @@ edition = "2021"
 exr = "1.4.1"
 clap = { version = "3.1.8", default-features = false, features=["std"] }
 colorbox = { git = "https://github.com/cessen/colorbox", branch = "master" }
+simplers_optimization = "0.4.3"

src/linear_log.rs

@@ -0,0 +1,107 @@
+/// A composite linear-log function.
+///
+/// `slope` is the slope of the linear segment.
+/// `base` is the log base.
+/// The offsets shift the linear and log parts of the curve along
+/// the linear color axis.
+pub fn linear_to_log(x: f64, line_offset: f64, slope: f64, log_offset: f64, base: f64) -> f64 {
+    // Transition point between log and linear.
+    let transition = 1.0 / (slope * base.ln());
+
+    let k = transition + log_offset;
+    let l = (transition - line_offset + log_offset) * slope - transition.log(base);
+
+    if x <= k {
+        (x - line_offset) * slope
+    } else {
+        (x - log_offset).log(base) + l
+    }
+}
+
+/// A composite linear-log function.
+///
+/// `slope` is the slope of the linear segment.
+/// `base` is the log base.
+/// The offsets shift the linear and log parts of the curve along
+/// the linear color axis.
+pub fn log_to_linear(x: f64, line_offset: f64, slope: f64, log_offset: f64, base: f64) -> f64 {
+    // Transition point between log and linear.
+    let transition = 1.0 / (slope * base.ln());
+
+    let k = (transition - line_offset + log_offset) * slope;
+    let l = (transition - line_offset + log_offset) * slope - transition.log(base);
+
+    if x <= k {
+        (x / slope) + line_offset
+    } else {
+        base.powf(x - l) + log_offset
+    }
+}
+
+//-------------------------------------------------------------
+
+/// Generates Rust code for a linear-to-log transfer function with the
+/// given parameters.
+pub fn generate_linear_to_log(line_offset: f64, slope: f64, log_offset: f64, base: f64) -> String {
+    let transition = 1.0 / (slope * base.ln());
+    let k = transition + log_offset;
+    let l = (transition - line_offset + log_offset) * slope - transition.log(base);
+
+    format!(
+        r#"
+pub fn linear_to_log(x: f32) -> f32 {{
+    const A: f32 = {};
+    const B: f32 = {};
+    const C: f32 = {};
+    const D: f32 = {};
+    const E: f32 = {};
+    const F: f32 = {};
+
+    if x <= A {{
+        (x - B) * C
+    }} else {{
+        (x - D).log2() * (1.0 / E) + F
+    }}
+}}
+"#,
+        k,
+        line_offset,
+        slope,
+        log_offset,
+        base.log2(),
+        l,
+    )
+}
+
+/// Generates Rust code for a log-to-linear transfer function with the
+/// given parameters.
+pub fn generate_log_to_linear(line_offset: f64, slope: f64, log_offset: f64, base: f64) -> String {
+    let transition = 1.0 / (slope * base.ln());
+    let k = (transition - line_offset + log_offset) * slope;
+    let l = (transition - line_offset + log_offset) * slope - transition.log(base);
+
+    format!(
+        r#"
+pub fn log_to_linear(x: f32) -> f32 {{
+    const A: f32 = {};
+    const B: f32 = {};
+    const C: f32 = {};
+    const D: f32 = {};
+    const E: f32 = {};
+    const F: f32 = {};
+
+    if x <= A {{
+        (x * (1.0 / C)) + B
+    }} else {{
+        ((x - F) * E).exp2() + D
+    }}
+}}
+"#,
+        k,
+        line_offset,
+        slope,
+        log_offset,
+        base.log2(),
+        l,
+    )
+}

src/main.rs

@@ -1,9 +1,11 @@
+mod linear_log;
+mod optimize_log;
 mod test_image;
 
 use std::{fs::File, io::BufWriter, path::Path};
 
 use clap::{Arg, Command};
-use colorbox::transfer_functions::srgb;
+// use colorbox::transfer_functions::srgb;
 
 use test_image::{GRADIENT_LEN, RES_X, RES_Y};
 
@@ -46,7 +48,7 @@ fn main() {
     } else {
         let input_path = args.value_of("input").unwrap();
 
-        let base_image = test_image::build();
+        // let base_image = test_image::build();
         let mut input_image = {
             let (image, res_x, res_y) = read_rgb_exr(input_path).unwrap();
             assert_eq!(
@@ -87,6 +89,8 @@ fn main() {
             gray_b.push(rgb[2]);
         }
 
+        optimize_log::find_parameters(&gray_r);
+
         // Write the LUT.
         colorbox::formats::cube::write_1d(
             BufWriter::new(File::create("test.cube").unwrap()),

src/optimize_log.rs

@@ -0,0 +1,57 @@
+use crate::linear_log::log_to_linear as log_to_lin;
+
+pub fn find_parameters(lut: &[f32]) {
+    let lin_norm = 1.0 / (lut.len() - 1) as f64;
+
+    // Compute the stuff that we can without estimation.
+    let offset = lut[0] as f64;
+    let slope = lin_norm / (lut[1] as f64 - lut[0] as f64);
+
+    // Select a range of points from the lookup table to fit to.
+    let idxs: Vec<_> = (0..lut.len()).step_by(lut.len() / 256).collect();
+    let coords: Vec<(f64, f64)> = idxs
+        .iter()
+        .map(|i| (*i as f64 * lin_norm, lut[*i] as f64))
+        .collect();
+
+    // Do the fitting.
+    let f = |v: &[f64]| {
+        let mut avg_sqr_err = 0.0f64;
+        for (x, y) in coords.iter().copied() {
+            let e = (log_to_lin(x, offset, slope, v[0], v[1]) - y).abs() / y.abs();
+            avg_sqr_err += e * e;
+        }
+        let last_y = lut[lut.len() - 1] as f64;
+        let e = (log_to_lin(1.0, offset, slope, v[0], v[1]) - last_y).abs() / last_y.abs();
+        avg_sqr_err += e * e;
+        avg_sqr_err
+    };
+    let input_interval = vec![(-0.2, 0.2), (1.1, 1000.0)];
+    let (_, params) = simplers_optimization::Optimizer::minimize(&f, &input_interval, 1000000);
+
+    // Calculate the error of our model.
+    let mut max_err = 0.0f64;
+    let mut avg_err = 0.0f64;
+    let mut avg_samples = 0usize;
+    for (i, y) in lut.iter().map(|y| *y as f64).enumerate() {
+        // We only record error for values that aren't crazy tiny, since
+        // their relative error isn't representative.
+        if y.abs() > 0.0001 {
+            let x = i as f64 * lin_norm;
+            let e = (log_to_lin(x, offset, slope, params[0], params[1]) - y).abs()
+            / y.abs();
+            max_err = max_err.max(e);
+            avg_err += e;
+            avg_samples += 1;
+        }
+    }
+    avg_err /= avg_samples as f64;
+
+    println!("Max Err: {:.4}%\nAvg Err: {:.4}%", max_err * 100.0, avg_err * 100.0);
+
+    println!(
+        "{}{}",
+        crate::linear_log::generate_linear_to_log(offset, slope, params[0], params[1],),
+        crate::linear_log::generate_log_to_linear(offset, slope, params[0], params[1],),
+    );
+}