cmsmtrx.c revision 12677:a4299d47bd00
1/* 2 * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER. 3 * 4 * This code is free software; you can redistribute it and/or modify it 5 * under the terms of the GNU General Public License version 2 only, as 6 * published by the Free Software Foundation. Oracle designates this 7 * particular file as subject to the "Classpath" exception as provided 8 * by Oracle in the LICENSE file that accompanied this code. 9 * 10 * This code is distributed in the hope that it will be useful, but WITHOUT 11 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 12 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 13 * version 2 for more details (a copy is included in the LICENSE file that 14 * accompanied this code). 15 * 16 * You should have received a copy of the GNU General Public License version 17 * 2 along with this work; if not, write to the Free Software Foundation, 18 * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. 19 * 20 * Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA 21 * or visit www.oracle.com if you need additional information or have any 22 * questions. 23 */ 24 25// This file is available under and governed by the GNU General Public 26// License version 2 only, as published by the Free Software Foundation. 27// However, the following notice accompanied the original version of this 28// file: 29// 30//--------------------------------------------------------------------------------- 31// 32// Little Color Management System 33// Copyright (c) 1998-2012 Marti Maria Saguer 34// 35// Permission is hereby granted, free of charge, to any person obtaining 36// a copy of this software and associated documentation files (the "Software"), 37// to deal in the Software without restriction, including without limitation 38// the rights to use, copy, modify, merge, publish, distribute, sublicense, 39// and/or sell copies of the Software, and to permit persons to whom the Software 40// is furnished to do so, subject to the following conditions: 41// 42// The above copyright notice and this permission notice shall be included in 43// all copies or substantial portions of the Software. 44// 45// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, 46// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO 47// THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND 48// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE 49// LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION 50// OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION 51// WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. 52// 53//--------------------------------------------------------------------------------- 54// 55 56#include "lcms2_internal.h" 57 58 59#define DSWAP(x, y) {cmsFloat64Number tmp = (x); (x)=(y); (y)=tmp;} 60 61 62// Initiate a vector 63void CMSEXPORT _cmsVEC3init(cmsVEC3* r, cmsFloat64Number x, cmsFloat64Number y, cmsFloat64Number z) 64{ 65 r -> n[VX] = x; 66 r -> n[VY] = y; 67 r -> n[VZ] = z; 68} 69 70// Vector substraction 71void CMSEXPORT _cmsVEC3minus(cmsVEC3* r, const cmsVEC3* a, const cmsVEC3* b) 72{ 73 r -> n[VX] = a -> n[VX] - b -> n[VX]; 74 r -> n[VY] = a -> n[VY] - b -> n[VY]; 75 r -> n[VZ] = a -> n[VZ] - b -> n[VZ]; 76} 77 78// Vector cross product 79void CMSEXPORT _cmsVEC3cross(cmsVEC3* r, const cmsVEC3* u, const cmsVEC3* v) 80{ 81 r ->n[VX] = u->n[VY] * v->n[VZ] - v->n[VY] * u->n[VZ]; 82 r ->n[VY] = u->n[VZ] * v->n[VX] - v->n[VZ] * u->n[VX]; 83 r ->n[VZ] = u->n[VX] * v->n[VY] - v->n[VX] * u->n[VY]; 84} 85 86// Vector dot product 87cmsFloat64Number CMSEXPORT _cmsVEC3dot(const cmsVEC3* u, const cmsVEC3* v) 88{ 89 return u->n[VX] * v->n[VX] + u->n[VY] * v->n[VY] + u->n[VZ] * v->n[VZ]; 90} 91 92// Euclidean length 93cmsFloat64Number CMSEXPORT _cmsVEC3length(const cmsVEC3* a) 94{ 95 return sqrt(a ->n[VX] * a ->n[VX] + 96 a ->n[VY] * a ->n[VY] + 97 a ->n[VZ] * a ->n[VZ]); 98} 99 100// Euclidean distance 101cmsFloat64Number CMSEXPORT _cmsVEC3distance(const cmsVEC3* a, const cmsVEC3* b) 102{ 103 cmsFloat64Number d1 = a ->n[VX] - b ->n[VX]; 104 cmsFloat64Number d2 = a ->n[VY] - b ->n[VY]; 105 cmsFloat64Number d3 = a ->n[VZ] - b ->n[VZ]; 106 107 return sqrt(d1*d1 + d2*d2 + d3*d3); 108} 109 110 111 112// 3x3 Identity 113void CMSEXPORT _cmsMAT3identity(cmsMAT3* a) 114{ 115 _cmsVEC3init(&a-> v[0], 1.0, 0.0, 0.0); 116 _cmsVEC3init(&a-> v[1], 0.0, 1.0, 0.0); 117 _cmsVEC3init(&a-> v[2], 0.0, 0.0, 1.0); 118} 119 120static 121cmsBool CloseEnough(cmsFloat64Number a, cmsFloat64Number b) 122{ 123 return fabs(b - a) < (1.0 / 65535.0); 124} 125 126 127cmsBool CMSEXPORT _cmsMAT3isIdentity(const cmsMAT3* a) 128{ 129 cmsMAT3 Identity; 130 int i, j; 131 132 _cmsMAT3identity(&Identity); 133 134 for (i=0; i < 3; i++) 135 for (j=0; j < 3; j++) 136 if (!CloseEnough(a ->v[i].n[j], Identity.v[i].n[j])) return FALSE; 137 138 return TRUE; 139} 140 141 142// Multiply two matrices 143void CMSEXPORT _cmsMAT3per(cmsMAT3* r, const cmsMAT3* a, const cmsMAT3* b) 144{ 145#define ROWCOL(i, j) \ 146 a->v[i].n[0]*b->v[0].n[j] + a->v[i].n[1]*b->v[1].n[j] + a->v[i].n[2]*b->v[2].n[j] 147 148 _cmsVEC3init(&r-> v[0], ROWCOL(0,0), ROWCOL(0,1), ROWCOL(0,2)); 149 _cmsVEC3init(&r-> v[1], ROWCOL(1,0), ROWCOL(1,1), ROWCOL(1,2)); 150 _cmsVEC3init(&r-> v[2], ROWCOL(2,0), ROWCOL(2,1), ROWCOL(2,2)); 151 152#undef ROWCOL //(i, j) 153} 154 155 156 157// Inverse of a matrix b = a^(-1) 158cmsBool CMSEXPORT _cmsMAT3inverse(const cmsMAT3* a, cmsMAT3* b) 159{ 160 cmsFloat64Number det, c0, c1, c2; 161 162 c0 = a -> v[1].n[1]*a -> v[2].n[2] - a -> v[1].n[2]*a -> v[2].n[1]; 163 c1 = -a -> v[1].n[0]*a -> v[2].n[2] + a -> v[1].n[2]*a -> v[2].n[0]; 164 c2 = a -> v[1].n[0]*a -> v[2].n[1] - a -> v[1].n[1]*a -> v[2].n[0]; 165 166 det = a -> v[0].n[0]*c0 + a -> v[0].n[1]*c1 + a -> v[0].n[2]*c2; 167 168 if (fabs(det) < MATRIX_DET_TOLERANCE) return FALSE; // singular matrix; can't invert 169 170 b -> v[0].n[0] = c0/det; 171 b -> v[0].n[1] = (a -> v[0].n[2]*a -> v[2].n[1] - a -> v[0].n[1]*a -> v[2].n[2])/det; 172 b -> v[0].n[2] = (a -> v[0].n[1]*a -> v[1].n[2] - a -> v[0].n[2]*a -> v[1].n[1])/det; 173 b -> v[1].n[0] = c1/det; 174 b -> v[1].n[1] = (a -> v[0].n[0]*a -> v[2].n[2] - a -> v[0].n[2]*a -> v[2].n[0])/det; 175 b -> v[1].n[2] = (a -> v[0].n[2]*a -> v[1].n[0] - a -> v[0].n[0]*a -> v[1].n[2])/det; 176 b -> v[2].n[0] = c2/det; 177 b -> v[2].n[1] = (a -> v[0].n[1]*a -> v[2].n[0] - a -> v[0].n[0]*a -> v[2].n[1])/det; 178 b -> v[2].n[2] = (a -> v[0].n[0]*a -> v[1].n[1] - a -> v[0].n[1]*a -> v[1].n[0])/det; 179 180 return TRUE; 181} 182 183 184// Solve a system in the form Ax = b 185cmsBool CMSEXPORT _cmsMAT3solve(cmsVEC3* x, cmsMAT3* a, cmsVEC3* b) 186{ 187 cmsMAT3 m, a_1; 188 189 memmove(&m, a, sizeof(cmsMAT3)); 190 191 if (!_cmsMAT3inverse(&m, &a_1)) return FALSE; // Singular matrix 192 193 _cmsMAT3eval(x, &a_1, b); 194 return TRUE; 195} 196 197// Evaluate a vector across a matrix 198void CMSEXPORT _cmsMAT3eval(cmsVEC3* r, const cmsMAT3* a, const cmsVEC3* v) 199{ 200 r->n[VX] = a->v[0].n[VX]*v->n[VX] + a->v[0].n[VY]*v->n[VY] + a->v[0].n[VZ]*v->n[VZ]; 201 r->n[VY] = a->v[1].n[VX]*v->n[VX] + a->v[1].n[VY]*v->n[VY] + a->v[1].n[VZ]*v->n[VZ]; 202 r->n[VZ] = a->v[2].n[VX]*v->n[VX] + a->v[2].n[VY]*v->n[VY] + a->v[2].n[VZ]*v->n[VZ]; 203} 204