diff options
Diffstat (limited to 'crl/cri.py')
-rw-r--r-- | crl/cri.py | 117 |
1 files changed, 117 insertions, 0 deletions
diff --git a/crl/cri.py b/crl/cri.py new file mode 100644 index 0000000..af62b1f --- /dev/null +++ b/crl/cri.py @@ -0,0 +1,117 @@ +import numpy as np +import scipy.spatial +import crl.color as color +import crl.tables as tables + +class CRIError(Exception): + pass + +class CRIResult: + def __init__(self, method_name): + self.method_name = method_name + self.SCRIs = [] + self.tcs = [] + self.tcs_ref = [] + + def __str__(self): + ret = f"CRI result for {self.method_name}\n" + ret += f"\tXYZ = {self.XYZ}\n" + ret += f"\tCCT = {self.CCT} K (delta uv = {self.duv})\n" + SCRIs = ", ".join("%.1f" % SCRI for SCRI in self.SCRIs) + ret += f"\tSCRIs = [{SCRIs}]\n" + ret += f"\tCRI = {self.CRI}" + return ret + +def CAT(sample, ref, test): + def c(u, v, _): + return (4 - u - 10 * v) / v + def d(u, v, _): + return (1.708 * v - 1.481 * u + 0.404) / v + + cs, ds = c(*sample), d(*sample) + cr, dr = c(*ref), d(*ref) + ct, dt = c(*test), d(*test) + + u = (10.872 + 0.404 * cr / ct * cs - 4 * dr / dt * ds) \ + / (16.518 + 1.481 * cr / ct * cs - dr / dt * ds) + v = 5.52 / (16.518 + 1.481 * cr / ct * cs - dr / dt * ds) + + return color.uvY(u, v, sample[2]) + +def cri_Ra(S, observer=color.Observer_2deg, ignore_errors=False): + res = CRIResult("CRI Ra") + res.S = S.normed() + res.XYZ = color.XYZ(res.S, observer=observer) + res.CCT, res.duv = color.cct(res.XYZ) + + if res.duv > 5e-2 and not ignore_errors: + raise CRIError(f"spectrum not white enough with delta uv = {res.duv}") + + if res.CCT < 5000: + white = color.BlackBodySpectrum(res.CCT).normed() + else: + if res.CCT > 25000 and ignore_errors: + res.CCT = 25000 + white = color.CIEDaylightSpectrum(res.CCT).normed() + + res.XYZ_ref = color.XYZ(white, observer=observer) + + class UVWstarRef(color.UVWstar): + white_uvY = color.uvY(white, observer=observer) + + class UVWstar(color.UVWstar): + white_uvY = color.uvY(S, observer=observer) + + for i in range(14): + tcs = color.Spectrum(tables.tcs_ra[:, 0], tables.tcs_ra[:, 1 + i]) + + lit = color.uvY(color.CombinedSpectrum(tcs, res.S)) + lit = CAT(lit, UVWstarRef.white_uvY, UVWstar.white_uvY) + lit = UVWstarRef(lit) + lit_ref = UVWstarRef(color.CombinedSpectrum(tcs, white), observer=observer) + + res.tcs.append(lit) + res.tcs_ref.append(lit_ref) + + dE = np.linalg.norm(lit.array - lit_ref.array) + res.SCRIs.append(100 - 4.6 * dE) + + res.CRI = np.mean(res.SCRIs[:8]) + res.ECRI = np.mean(res.SCRIs) + return res + +def cri_GAI(spec, observer=color.Observer_2deg): + points = [] + for i in range(8): + tcs = color.Spectrum(color.tables.tcs_ra[:, 0], color.tables.tcs_ra[:, 1 + i]) + u, v, _ = color.uvY76(color.CombinedSpectrum(tcs, spec), observer=observer) + points.append([u, v]) + points = np.array(points) + hull = scipy.spatial.ConvexHull(points) + + # The magic number is the hull.volume of the E illuminant + GAI = hull.volume * 100 / 0.007351716795637625 + + return GAI, hull.points[hull.vertices] + +def cri_FSI(spec): + wvl = np.arange(380, 730, 1) + Y = spec.Yi(wvl) + dwvl = wvl[1] - wvl[0] + Y = spec.Yi(wvl) + Y /= np.sum(Y * dwvl) + + CEE = (wvl - 380) / (730 - 380) + + Is = [] + for i, _ in enumerate(wvl): + C = np.cumsum(Y * dwvl) + D = (C - CEE) ** 2 + Is.append(sum(D * dwvl)) + Y = np.roll(Y, 1) + + FSI = np.mean(Is) + return FSI + +def cri_FSCI(spec): + return 100 - 5.1 * cri_FSI(spec) |