1 //--------------------------------------------------------------------------------- 2 // 3 // Little Color Management System 4 // Copyright (c) 1998-2016 Marti Maria Saguer 5 // 6 // Permission is hereby granted, free of charge, to any person obtaining 7 // a copy of this software and associated documentation files (the "Software"), 8 // to deal in the Software without restriction, including without limitation 9 // the rights to use, copy, modify, merge, publish, distribute, sublicense, 10 // and/or sell copies of the Software, and to permit persons to whom the Software 11 // is furnished to do so, subject to the following conditions: 12 // 13 // The above copyright notice and this permission notice shall be included in 14 // all copies or substantial portions of the Software. 15 // 16 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, 17 // EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO 18 // THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND 19 // NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE 20 // LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION 21 // OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION 22 // WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. 23 // 24 //--------------------------------------------------------------------------------- 25 // 26 27 #include "lcms2_internal.h" 28 29 30 #define DSWAP(x, y) {cmsFloat64Number tmp = (x); (x)=(y); (y)=tmp;} 31 32 33 // Initiate a vector 34 void CMSEXPORT _cmsVEC3init(cmsVEC3* r, cmsFloat64Number x, cmsFloat64Number y, cmsFloat64Number z) 35 { 36 r -> n[VX] = x; 37 r -> n[VY] = y; 38 r -> n[VZ] = z; 39 } 40 41 // Vector subtraction 42 void CMSEXPORT _cmsVEC3minus(cmsVEC3* r, const cmsVEC3* a, const cmsVEC3* b) 43 { 44 r -> n[VX] = a -> n[VX] - b -> n[VX]; 45 r -> n[VY] = a -> n[VY] - b -> n[VY]; 46 r -> n[VZ] = a -> n[VZ] - b -> n[VZ]; 47 } 48 49 // Vector cross product 50 void CMSEXPORT _cmsVEC3cross(cmsVEC3* r, const cmsVEC3* u, const cmsVEC3* v) 51 { 52 r ->n[VX] = u->n[VY] * v->n[VZ] - v->n[VY] * u->n[VZ]; 53 r ->n[VY] = u->n[VZ] * v->n[VX] - v->n[VZ] * u->n[VX]; 54 r ->n[VZ] = u->n[VX] * v->n[VY] - v->n[VX] * u->n[VY]; 55 } 56 57 // Vector dot product 58 cmsFloat64Number CMSEXPORT _cmsVEC3dot(const cmsVEC3* u, const cmsVEC3* v) 59 { 60 return u->n[VX] * v->n[VX] + u->n[VY] * v->n[VY] + u->n[VZ] * v->n[VZ]; 61 } 62 63 // Euclidean length 64 cmsFloat64Number CMSEXPORT _cmsVEC3length(const cmsVEC3* a) 65 { 66 return sqrt(a ->n[VX] * a ->n[VX] + 67 a ->n[VY] * a ->n[VY] + 68 a ->n[VZ] * a ->n[VZ]); 69 } 70 71 // Euclidean distance 72 cmsFloat64Number CMSEXPORT _cmsVEC3distance(const cmsVEC3* a, const cmsVEC3* b) 73 { 74 cmsFloat64Number d1 = a ->n[VX] - b ->n[VX]; 75 cmsFloat64Number d2 = a ->n[VY] - b ->n[VY]; 76 cmsFloat64Number d3 = a ->n[VZ] - b ->n[VZ]; 77 78 return sqrt(d1*d1 + d2*d2 + d3*d3); 79 } 80 81 82 83 // 3x3 Identity 84 void CMSEXPORT _cmsMAT3identity(cmsMAT3* a) 85 { 86 _cmsVEC3init(&a-> v[0], 1.0, 0.0, 0.0); 87 _cmsVEC3init(&a-> v[1], 0.0, 1.0, 0.0); 88 _cmsVEC3init(&a-> v[2], 0.0, 0.0, 1.0); 89 } 90 91 static 92 cmsBool CloseEnough(cmsFloat64Number a, cmsFloat64Number b) 93 { 94 return fabs(b - a) < (1.0 / 65535.0); 95 } 96 97 98 cmsBool CMSEXPORT _cmsMAT3isIdentity(const cmsMAT3* a) 99 { 100 cmsMAT3 Identity; 101 int i, j; 102 103 _cmsMAT3identity(&Identity); 104 105 for (i=0; i < 3; i++) 106 for (j=0; j < 3; j++) 107 if (!CloseEnough(a ->v[i].n[j], Identity.v[i].n[j])) return FALSE; 108 109 return TRUE; 110 } 111 112 113 // Multiply two matrices 114 void CMSEXPORT _cmsMAT3per(cmsMAT3* r, const cmsMAT3* a, const cmsMAT3* b) 115 { 116 #define ROWCOL(i, j) \ 117 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] 118 119 _cmsVEC3init(&r-> v[0], ROWCOL(0,0), ROWCOL(0,1), ROWCOL(0,2)); 120 _cmsVEC3init(&r-> v[1], ROWCOL(1,0), ROWCOL(1,1), ROWCOL(1,2)); 121 _cmsVEC3init(&r-> v[2], ROWCOL(2,0), ROWCOL(2,1), ROWCOL(2,2)); 122 123 #undef ROWCOL //(i, j) 124 } 125 126 127 128 // Inverse of a matrix b = a^(-1) 129 cmsBool CMSEXPORT _cmsMAT3inverse(const cmsMAT3* a, cmsMAT3* b) 130 { 131 cmsFloat64Number det, c0, c1, c2; 132 133 c0 = a -> v[1].n[1]*a -> v[2].n[2] - a -> v[1].n[2]*a -> v[2].n[1]; 134 c1 = -a -> v[1].n[0]*a -> v[2].n[2] + a -> v[1].n[2]*a -> v[2].n[0]; 135 c2 = a -> v[1].n[0]*a -> v[2].n[1] - a -> v[1].n[1]*a -> v[2].n[0]; 136 137 det = a -> v[0].n[0]*c0 + a -> v[0].n[1]*c1 + a -> v[0].n[2]*c2; 138 139 if (fabs(det) < MATRIX_DET_TOLERANCE) return FALSE; // singular matrix; can't invert 140 141 b -> v[0].n[0] = c0/det; 142 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; 143 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; 144 b -> v[1].n[0] = c1/det; 145 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; 146 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; 147 b -> v[2].n[0] = c2/det; 148 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; 149 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; 150 151 return TRUE; 152 } 153 154 155 // Solve a system in the form Ax = b 156 cmsBool CMSEXPORT _cmsMAT3solve(cmsVEC3* x, cmsMAT3* a, cmsVEC3* b) 157 { 158 cmsMAT3 m, a_1; 159 160 memmove(&m, a, sizeof(cmsMAT3)); 161 162 if (!_cmsMAT3inverse(&m, &a_1)) return FALSE; // Singular matrix 163 164 _cmsMAT3eval(x, &a_1, b); 165 return TRUE; 166 } 167 168 // Evaluate a vector across a matrix 169 void CMSEXPORT _cmsMAT3eval(cmsVEC3* r, const cmsMAT3* a, const cmsVEC3* v) 170 { 171 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]; 172 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]; 173 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]; 174 } 175 176 177
