Trying to figure out why the solver diverges...

4 files changed, 174 insertions, 0 deletions

diff --git a/test_compare_with_tables.py b/test_compare_with_tables.py
new file mode 100644
index 0000000..a13be00
--- /dev/null
+++ b/test_compare_with_tables.py
@@ -0,0 +1,68 @@
+import numpy as np
+from colour import *
+from colour.difference import delta_E_CIE1976
+from colour.models import RGB_COLOURSPACES
+from colour.recovery import RGB_to_sd_Jakob2019, Jakob2019Interpolator
+from gsoc_common import *
+
+colourspace = RGB_COLOURSPACES["ACES2065-1"]
+
+
+def stats(table, prefix):
+    table = np.array(table)
+    print("%s min/max/avg/std: %g/%g/%g/%g" % (prefix, np.min(table),
+           np.max(table), np.mean(table), np.std(table)))
+
+
+if __name__ == "__main__":
+    interp = Jakob2019Interpolator()
+    interp.from_file("data/aces2065_1.coeff")
+
+    errors = []
+    errors_interp = []
+    rel_deltas = []
+    for i in range(100):
+        RGB = np.random.random(3)
+
+        for use_feedback in [False, True]:
+            recovered_sd, error = RGB_to_sd_Jakob2019(
+                RGB,
+                colourspace,
+                return_error=True,
+                use_feedback=use_feedback
+            )
+
+            if error < 1e-3:
+                break
+
+        XYZ = RGB_to_XYZ(
+            RGB,
+            colourspace.whitepoint,
+            D65_xy,
+            colourspace.RGB_to_XYZ_matrix,
+        )
+
+        Lab = XYZ_to_Lab(XYZ, D65_xy)
+
+        sd_interp = interp.RGB_to_sd(RGB)
+        XYZ_interp = sd_to_XYZ(sd_interp, illuminant=D65) / 100
+        Lab_interp = XYZ_to_Lab(XYZ_interp, D65_xy)
+        error_interp = delta_E_CIE1976(Lab, Lab_interp)
+
+        XYZ_recovered = sd_to_XYZ(recovered_sd, illuminant=D65) / 100
+        Lab_recovered = XYZ_to_Lab(XYZ_recovered, D65_xy)
+        rel_delta = error / error_interp
+
+        print("%.5f %.5f %.5f  E%10.6f  I%7.3f  RΔ%.6f"
+              % (*RGB, error, error_interp, rel_delta))
+
+        errors.append(error)
+        errors_interp.append(error_interp)
+        rel_deltas.append(rel_delta)
+
+    print(colourspace.name)
+    stats(errors, "Colour")
+    stats(errors_interp, "Interpolator")
+    stats(rel_deltas, "Rel. delta")
+
+
diff --git a/test_diff.py b/test_diff.py
new file mode 100644
index 0000000..775752f
--- /dev/null
+++ b/test_diff.py
@@ -0,0 +1,51 @@
+import numpy as np
+from scipy.optimize import minimize
+from colour import *
+from colour.recovery.jakob2019 import error_function
+from matplotlib import pyplot as plt
+from gsoc_common import plot_comparison
+
+
+# This test checks if derivatives are calculated correctly by comparing them
+# to finite differences.
+if __name__ == "__main__":
+    shape = SpectralShape(360, 830, 1)
+    cmfs = STANDARD_OBSERVER_CMFS["CIE 1931 2 Degree Standard Observer"].align(shape)
+
+    illuminant = SpectralDistribution(ILLUMINANT_SDS["D65"]).align(shape)
+    illuminant_XYZ = sd_to_XYZ(illuminant) / 100
+
+    target = np.array([22.955, 46.374, 3.5981]) # Some arbitrary Lab colour
+    xs = np.linspace(-100, 100, 500)
+    h = xs[1] - xs[0]
+
+    # Vary one coefficient at a time
+    for c_index in range(3):
+        errors = np.empty(len(xs))
+        derrors = np.empty(len(xs))
+
+        for i, x in enumerate(xs):
+            c = np.array([0.0, 0, 0])
+            c[c_index] = x
+
+            error, derror_dc = error_function(
+                c, target, shape, cmfs, illuminant, illuminant_XYZ
+            )
+
+            errors[i] = error
+            derrors[i] = derror_dc[c_index]
+
+
+        plt.subplot(2, 3, 1 + c_index)
+        plt.xlabel("c%d" % c_index)
+        plt.ylabel("ΔE")
+        plt.plot(xs, errors)
+
+        plt.subplot(2, 3, 4 + c_index)
+        plt.xlabel("c%d" % c_index)
+        plt.ylabel("dΔE/dc%d" % c_index)
+
+        plt.plot(xs, derrors, "k-")
+        plt.plot(xs[:-1] + h / 2, np.diff(errors) / h, "r:")
+
+    plt.show()
diff --git a/test_generate.py b/test_generate.py
new file mode 100644
index 0000000..20d5b14
--- /dev/null
+++ b/test_generate.py
@@ -0,0 +1,19 @@
+import numpy as np
+from colour import *
+from colour.difference import delta_E_CIE1976
+from colour.models import RGB_COLOURSPACES
+from colour.recovery import RGB_to_sd_Jakob2019, Jakob2019Interpolator
+from gsoc_common import *
+
+colourspace = RGB_COLOURSPACES["sRGB"]
+
+
+if __name__ == "__main__":
+    interp = Jakob2019Interpolator()
+    interp.generate(
+        RGB_COLOURSPACES["sRGB"],
+        STANDARD_OBSERVER_CMFS['CIE 1931 2 Degree Standard Observer'],
+        D65,
+        16,
+        64
+    )
diff --git a/test_what.py b/test_what.py
new file mode 100644
index 0000000..92b98e8
--- /dev/null
+++ b/test_what.py
@@ -0,0 +1,36 @@
+import numpy as np
+from colour import *
+from colour.difference import delta_E_CIE1976
+from colour.models import RGB_COLOURSPACES
+from colour.recovery import RGB_to_sd_Jakob2019, Jakob2019Interpolator
+from gsoc_common import *
+
+colourspace = RGB_COLOURSPACES["sRGB"]
+
+
+if __name__ == "__main__":
+    RGB = np.array([0.1668, 0.000, 0.0334])
+
+
+    matched_sd, error = RGB_to_sd_Jakob2019(
+        RGB,
+        colourspace,
+        #coefficients_0=(57.94046366, -39.65092867, 2.48839434),
+        return_error=True,
+        use_feedback=True,
+    )
+
+    print("Error is %g" % error)
+
+    XYZ = RGB_to_XYZ(
+        RGB,
+        colourspace.whitepoint,
+        D65_xy,
+        colourspace.RGB_to_XYZ_matrix,
+    )
+
+    print("Target was %s %s" % (colourspace.name, RGB))
+    print("...XYZ = %s" % XYZ)
+    print("...Lab = %s" % XYZ_to_Lab(XYZ, D65_xy))
+
+    plot_comparison(XYZ, matched_sd, "target", error, D65)