diff options
Diffstat (limited to 'matrices.py')
| -rwxr-xr-x | matrices.py | 184 |
1 files changed, 184 insertions, 0 deletions
diff --git a/matrices.py b/matrices.py new file mode 100755 index 0000000..485374c --- /dev/null +++ b/matrices.py @@ -0,0 +1,184 @@ +#!/usr/bin/env python +# Both Python 2 and Python 3 will work + +import sys, math + +FORMAT = '%.50f' +DELIMITER = ', ' + +def _divrow(R, d): + return [c / d for c in R] + +def _subrow(A, B, m): + return [a - b * m for a, b in zip(A, B)] + +def transpose(M): + r, c = len(M), len(M[0]) + return [[M[y][x] for y in range(r)] for x in range(c)] + +def to_vector(M): + if len(M[0]) == 1: + return tuple(R[0] for R in M) + elif len(M) == 1: + return M[0] + raise Exception('Invalid arguments') + +def _multiply(A, B): + if isinstance(A[0], int) or isinstance(A[0], float): + A = [A] + if isinstance(B[0], int) or isinstance(B[0], float): + B = [B] + ar, ac, br, bc = len(A), len(A[0]), len(B), len(B[0]) + if ac != br: + if ac == bc: + B, br, bc = transpose(B), bc, br + elif ar == br: + A, ar, ac = transpose(A), ac, ar + else: + raise Exception('Invalid arguments') + return [[sum(A[r][i] * B[i][c] for i in range(ac)) for c in range(bc)] for r in range(ar)] + +def multiply(*Ms): + R, Ms = Ms[-1], list(Ms)[-2::-1] + for M in Ms: + R = _multiply(M, R) + return R + +def invert(M): + r, c = len(M), len(M[0]) + if r != c: + raise Exception('Invalid arguments') + I = [[1 if x == y else 0 for x in range(c)] for y in range(r)] + M = [MR + IR for MR, IR in zip(M, I)] + for r1 in range(r): + if M[r1][r1] == 0: + for r2 in range(r1 + 1, r): + if M[r2][r1] != 0: + break + if r2 == r: + raise Exception('Not invertable') + M[r1], M[r2] = M[r2], M[r1] + M[r1] = _divrow(M[r1], M[r1][r1]) + for r2 in range(r1 + 1, r): + M[r2] = _subrow(M[r2], M[r1], M[r2][r1]) + for r1 in reversed(range(r)): + for r2 in reversed(range(r1)): + M[r2] = _subrow(M[r2], M[r1], M[r2][r1]) + return [R[c:] for R in M] + +def str_matrix(M): + rfmt = '%s%s%s%s%s' % (FORMAT, DELIMITER, FORMAT, DELIMITER, FORMAT) + fmt = '%s\n%s\n%s' % (rfmt, rfmt, rfmt) + return fmt % tuple(M[0] + M[1] +M[2]) + +def print_matrix(M): + print(str_matrix(M)) + +def matrix(m11, m12, m13, m21, m22, m23, m31, m32, m33): + return [[m11, m12, m13], [m21, m22, m23], [m31, m32, m33]] + +def _cos(deg): + return math.cos(deg * math.pi / 180.) + +def _sin(deg): + return math.sin(deg * math.pi / 180.) + +def _yiq_i(k, a, b): + return 0.877 * _cos(33.) * (a - k) - (0.492 * _sin(33.)) * (b - k) + +def _yiq_q(k, a, b): + return 0.877 * _sin(33.) * (a - k) + (0.492 * _cos(33.)) * (b - k) + + +matrices = { + ('srgb', 'yes') : matrix(0.253, 0.684, 0.063, + 0.50, -0.50, 0., + 0.25, 0.25, -0.50), + + ('ycgco', 'srgb') : matrix(1., -1., 1., + 1., 1., 0., + 1., -1., -1.), + + ('srgb', 'yiq') : matrix(0.299, 0.587, 0.114, + _yiq_i(0.299, 1, 0), _yiq_i(0.587, 0, 0), _yiq_i(0.114, 0, 1), + _yiq_q(0.299, 1, 0), _yiq_q(0.587, 0, 0), _yiq_q(0.114, 0, 1)), + + ('srgb', 'ydbdr') : matrix( 0.299, 0.587, 0.114, + -0.450, -0.883, 1.333, + -1.333, 1.116, 0.217), + + ('srgb', 'ypbpr') : matrix( 0.2126, 0.7152, 0.0722, + -0.2126, -0.7152, 1 - 0.0722, + 1 - 0.2126, -0.7152, -0.0722), + + ('srgb', 'ciexyz') : matrix(0.412457445582367576708548995157, + 0.357575865245515878143578447634, + 0.180437247826399665973085006954, + 0.212673370378408277403536885686, + 0.715151730491031756287156895269, + 0.072174899130559869164791564344, + 0.019333942761673460208893260415, + 0.119191955081838593666354597644, + 0.950302838552371742508739771438), + + ('yuv', 'ydbdr') : matrix(1., 0., 0., + 0., 3.069, 0., + 0., 0., -2.169), +} + +if __name__ == '__main__': + models = set() + for (f, t) in list(matrices.keys()): + models.add(f) + models.add(t) + if (t, f) not in matrices: + matrices[(t, f)] = invert(matrices[(f, t)]) + + while True: + ms = list(matrices.keys()) + added = 0 + for m1 in models: + for m2 in models: + if m1 == m2: + continue + for m3 in models: + if m1 == m3 or m2 == m3: + continue + if (m1, m3) in ms: + continue + if (m1, m2) not in ms: + continue + if (m2, m3) not in ms: + continue + added += 1 + m = multiply(matrices[(m2, m3)], matrices[(m1, m2)]) + matrices[(m1, m3)] = m + if not added: + break + + for (f, t) in sorted(matrices.keys()): + M = multiply(matrices[(t, f)], matrices[(f, t)]) + for r in range(len(M)): + for c in range(len(M[r])): + x = M[r][c] + x -= 1 if r == c else 0 + x *= x + if x > 0.0001: + sys.stderr.write('(%s, %s)(%s, %s) != I\n' % (t, f, f, t)) + for r in M: + for c in r: + sys.stderr.write('\t%lg' % c) + sys.stderr.write('\n') + sys.exit(1) + + for (f, t) in sorted(matrices.keys()): + M = matrices[(f, t)] + macro = '%s_TO_%s' % (f.upper(), t.upper()) + print('') + print('#ifndef %s' % macro) + print('# define %s MATRIX(\\' % macro) + for row in M: + a = ', '.join([FORMAT] * len(row)) + b = ')' if row is M[-1] else ',\\' + print(('\t' + a + b) % tuple(row)) + print('#endif') |