code cleanup

This commit is contained in:
Matt Terry
2013-08-05 12:18:32 -07:00
parent a8403d26a1
commit dc692dddcc
2 changed files with 89 additions and 96 deletions
+11 -4
View File
@@ -1075,13 +1075,20 @@ def lab2lch(lab):
lch = _prepare_lab_array(lab)
a, b = lch[..., 1], lch[..., 2]
lch[..., 1], lch[..., 2] = np.sqrt(a ** 2 + b ** 2), np.arctan2(b, a)
H = lch[..., 2]
H += np.where(H < 0, 2*np.pi, 0) # (-pi, pi) -> (0, 2*pi)
lch[..., 1], lch[..., 2] = _cart2polar_2pi(a, b)
return lch
def _cart2polar_2pi(x, y):
"""convert cartesian coordiantes to polar (uses non-standard theta range!)
NON-STANDARD RANGE! Maps to (0, 2*pi) rather than usual (-pi, +pi)
"""
r, t = np.hypot(x, y), np.arctan2(y, x)
t += np.where(t < 0., 2 * np.pi, 0)
return r, t
def lch2lab(lch):
"""CIE-LCH to CIE-LAB color space conversion.
+78 -92
View File
@@ -21,12 +21,7 @@ from __future__ import division
import numpy as np
def _arctan2pi(b, a):
"""np.arctan2 mapped to (0, 2 * pi)"""
ans = np.arctan2(b, a)
ans += np.where(ans < 0, 2 * np.pi, 0.)
return ans
from skimage.color.colorconv import lab2lch, _cart2polar_2pi
def deltaE_cie76(lab1, lab2):
@@ -50,9 +45,9 @@ def deltaE_cie76(lab1, lab2):
.. [2] A. R. Robertson, "The CIE 1976 color-difference formulae,"
Color Res. Appl. 2, 7-11 (1977).
"""
l1, a1, b1 = np.rollaxis(lab1, -1)[:3]
l2, a2, b2 = np.rollaxis(lab2, -1)[:3]
return np.sqrt((l2 - l1) ** 2 + (a2 - a1) ** 2 + (b2 - b1) ** 2)
L1, a1, b1 = np.rollaxis(lab1, -1)[:3]
L2, a2, b2 = np.rollaxis(lab2, -1)[:3]
return np.sqrt((L2 - L1) ** 2 + (a2 - a1) ** 2 + (b2 - b1) ** 2)
def deltaE_ciede94(lab1, lab2, kH=1, kC=1, kL=1, k1=0.045, k2=0.015):
@@ -104,22 +99,20 @@ def deltaE_ciede94(lab1, lab2, kH=1, kC=1, kL=1, k1=0.045, k2=0.015):
.. [1] http://en.wikipedia.org/wiki/Color_difference
.. [2] http://www.brucelindbloom.com/index.html?Eqn_DeltaE_CIE94.html
"""
l1, a1, b1 = np.rollaxis(lab1, -1)[:3]
l2, a2, b2 = np.rollaxis(lab2, -1)[:3]
L1, C1, h1 = np.rollaxis(lab2lch(lab1), -1)[:3]
L2, C2, h2 = np.rollaxis(lab2lch(lab2), -1)[:3]
dl = l1 - l2
c1 = np.hypot(a1, b1)
c2 = np.hypot(a2, b2)
dc = c1 - c2
dh_ab = np.sqrt(deltaE_cie76(lab1, lab2) ** 2 - dl ** 2 - dc ** 2)
dL = L1 - L2
dC = C1 - C2
dH_ab = np.sqrt(deltaE_cie76(lab1, lab2) ** 2 - dL ** 2 - dC ** 2)
SL = 1
SC = 1 + k1 * c1
SH = 1 + k2 * c1
SC = 1 + k1 * C1
SH = 1 + k2 * C1
ans = ((dl / (kL * SL)) ** 2 +
(dc / (kC * SC)) ** 2 +
(dh_ab / (kH * SH)) ** 2
ans = ((dL / (kL * SL)) ** 2 +
(dC / (kC * SC)) ** 2 +
(dH_ab / (kH * SH)) ** 2
)
return np.sqrt(ans)
@@ -166,76 +159,74 @@ def deltaE_ciede2000(lab1, lab2, kL=1, kC=1, kH=1):
L1, a1, b1 = np.rollaxis(lab1, -1)[:3]
L2, a2, b2 = np.rollaxis(lab2, -1)[:3]
c1 = np.hypot(a1, b1)
c2 = np.hypot(a2, b2)
cbar = 0.5 * (c1 + c2)
c7 = cbar ** 7
# distort `a` based on average chroma
# then convert to lch coordines from distorted `a`
# all subsequence calculations are in the new coordiantes
# (often denoted "prime" in the literature)
Cbar = 0.5 * (np.hypot(a1, b1) + np.hypot(a2, b2))
c7 = Cbar ** 7
G = 0.5 * (1 - np.sqrt(c7 / (c7 + 25 ** 7)))
scale = 1 + G
C1, h1 = _cart2polar_2pi(a1 * scale, b1)
C2, h2 = _cart2polar_2pi(a2 * scale, b2)
# recall that c, h are polar coordiantes. c==r, h==theta
dL_prime = L2 - L1
# cide2000 has four terms to delta_e:
# 1) Luminance term
# 2) Hue term
# 3) Chroma term
# 4) hue Rotation term
# luminance term
Lbar = 0.5 * (L1 + L2)
tmp = (Lbar - 50) ** 2
SL = 1 + 0.015 * tmp / np.sqrt(20 + tmp)
L_term = (L2 - L1) / (kL * SL)
a1_prime = a1 * (1 + G)
a2_prime = a2 * (1 + G)
# chroma term
Cbar = 0.5 * (C1 + C2) # new coordiantes
SC = 1 + 0.045 * Cbar
C_term = (C2 - C1) / (kC * SC)
c1_prime = np.hypot(a1_prime, b1)
c2_prime = np.hypot(a2_prime, b2)
cbar_prime = 0.5 * (c1_prime + c2_prime)
dC_prime = c2_prime - c1_prime
# hue term
h_diff = h2 - h1
h_sum = h1 + h2
CC = C1 * C2
h1_prime = _arctan2pi(b1, a1_prime)
h2_prime = _arctan2pi(b2, a2_prime)
dH = h_diff.copy()
dH[h_diff > np.pi] -= 2 * np.pi
dH[h_diff < -np.pi] += 2 * np.pi
dH[CC == 0.] = 0. # if r == 0, dtheta == 0
dH_term = 2 * np.sqrt(CC) * np.sin(dH / 2)
dh_prime = h2_prime - h1_prime
cc = c1_prime * c2_prime
mask1 = cc == 0.
mask2 = np.logical_and(-mask1, dh_prime > np.pi)
mask3 = np.logical_and(-mask1, dh_prime < -np.pi)
dh_prime = np.where(mask1, 0., dh_prime)
dh_prime += np.where(mask2, 2 * np.pi, 0)
dh_prime -= np.where(mask3, 2 * np.pi, 0)
dH_prime = 2 * np.sqrt(cc) * np.sin(dh_prime / 2)
Hbar_prime = h1_prime + h2_prime
mask0 = np.logical_and(np.abs(h1_prime - h2_prime) > np.pi, cc != 0.)
mask1 = np.logical_and(mask0, Hbar_prime < 2 * np.pi)
mask2 = np.logical_and(mask0, Hbar_prime >= 2 * np.pi)
Hbar_prime += np.where(mask1, 2 * np.pi, 0)
Hbar_prime -= np.where(mask2, 2 * np.pi, 0)
Hbar_prime *= np.where(cc == 0., 2, 1)
Hbar_prime *= 0.5
Hbar = h_sum.copy()
mask = np.logical_and(CC != 0., np.abs(h_diff) > np.pi)
Hbar[mask * (h_sum < 2 * np.pi)] += 2 * np.pi
Hbar[mask * (h_sum >= 2 * np.pi)] -= 2 * np.pi
Hbar[CC == 0.] *= 2
Hbar *= 0.5
T = (1 -
0.17 * np.cos(Hbar_prime - np.deg2rad(30)) +
0.24 * np.cos(2 * Hbar_prime) +
0.32 * np.cos(3 * Hbar_prime + np.deg2rad(6)) -
0.20 * np.cos(4 * Hbar_prime - np.deg2rad(63))
0.17 * np.cos(Hbar - np.deg2rad(30)) +
0.24 * np.cos(2 * Hbar) +
0.32 * np.cos(3 * Hbar + np.deg2rad(6)) -
0.20 * np.cos(4 * Hbar - np.deg2rad(63))
)
dTheta = (np.deg2rad(30) *
np.exp(-((np.rad2deg(Hbar_prime) - 275) / 25) ** 2)
)
c7 = cbar_prime ** 7
SH = 1 + 0.015 * Cbar * T
H_term = dH_term / (kH * SH)
# hue rotation
c7 = Cbar ** 7
Rc = 2 * np.sqrt(c7 / (c7 + 25 ** 7))
dtheta = np.deg2rad(30) * np.exp(-((np.rad2deg(Hbar) - 275) / 25) ** 2)
R_term = -np.sin(2 * dtheta) * Rc * C_term * H_term
term = (Lbar - 50) ** 2
SL = 1 + 0.015 * term / np.sqrt(20 + term)
SC = 1 + 0.045 * cbar_prime
SH = 1 + 0.015 * cbar_prime * T
RT = -np.sin(2 * dTheta) * Rc
l_term = dL_prime / (kL * SL)
c_term = dC_prime / (kC * SC)
h_term = dH_prime / (kH * SH)
r_term = RT * c_term * h_term
dE2 = l_term ** 2
dE2 += c_term ** 2
dE2 += h_term ** 2
dE2 += r_term
# put it all together
dE2 = L_term ** 2
dE2 += C_term ** 2
dE2 += H_term ** 2
dE2 += R_term
return np.sqrt(dE2)
@@ -277,27 +268,22 @@ def deltaE_cmc(lab1, lab2, kL=1, kC=1):
JPC79 colour-difference formula," J. Soc. Dyers Colour. 100, 128-132
(1984).
"""
l1, a1, b1 = np.rollaxis(lab1, -1)[:3]
l2, a2, b2 = np.rollaxis(lab2, -1)[:3]
L1, C1, h1 = np.rollaxis(lab2lch(lab1), -1)[:3]
L2, C2, h2 = np.rollaxis(lab2lch(lab2), -1)[:3]
c1 = np.hypot(a1, b1)
c2 = np.hypot(a2, b2)
dC = c1 - c2
dl = l1 - l2
dH = np.sqrt(deltaE_cie76(lab1, lab2) ** 2 - dl ** 2 - dC ** 2)
dC = C1 - C2
dL = L1 - L2
dH = np.sqrt(deltaE_cie76(lab1, lab2) ** 2 - dL ** 2 - dC ** 2)
dL = l1 - l2
h1 = _arctan2pi(b1, a1)
T = np.where(np.logical_and(np.rad2deg(h1) >= 164, np.rad2deg(h1) <= 345),
0.56 + 0.2 * np.abs(np.cos(h1 + np.deg2rad(168))),
0.36 + 0.4 * np.abs(np.cos(h1 + np.deg2rad(35)))
)
c1_4 = c1 ** 4
c1_4 = C1 ** 4
F = np.sqrt(c1_4 / (c1_4 + 1900))
SL = np.where(l1 < 16, 0.511, 0.040975 * l1 / (1. + 0.01765 * l1))
SC = 0.638 + 0.0638 * c1 / (1. + 0.0131 * c1)
SL = np.where(L1 < 16, 0.511, 0.040975 * L1 / (1. + 0.01765 * L1))
SC = 0.638 + 0.0638 * C1 / (1. + 0.0131 * C1)
SH = SC * (F * T + 1 - F)
dE2 = (dL / (kL * SL)) ** 2