1/* 2 * Copyright (c) 2002-2007, Communications and Remote Sensing Laboratory, Universite catholique de Louvain (UCL), Belgium 3 * Copyright (c) 2002-2007, Professor Benoit Macq 4 * Copyright (c) 2001-2003, David Janssens 5 * Copyright (c) 2002-2003, Yannick Verschueren 6 * Copyright (c) 2003-2007, Francois-Olivier Devaux and Antonin Descampe 7 * Copyright (c) 2005, Herve Drolon, FreeImage Team 8 * All rights reserved. 9 * 10 * Redistribution and use in source and binary forms, with or without 11 * modification, are permitted provided that the following conditions 12 * are met: 13 * 1. Redistributions of source code must retain the above copyright 14 * notice, this list of conditions and the following disclaimer. 15 * 2. Redistributions in binary form must reproduce the above copyright 16 * notice, this list of conditions and the following disclaimer in the 17 * documentation and/or other materials provided with the distribution. 18 * 19 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS `AS IS' 20 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 21 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 22 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE 23 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 24 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 25 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS 26 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN 27 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 28 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 29 * POSSIBILITY OF SUCH DAMAGE. 30 */ 31 32#ifdef __SSE__ 33#include <xmmintrin.h> 34#endif 35 36#include "opj_includes.h" 37 38/* <summary> */ 39/* This table contains the norms of the basis function of the reversible MCT. */ 40/* </summary> */ 41static const double mct_norms[3] = { 1.732, .8292, .8292 }; 42 43/* <summary> */ 44/* This table contains the norms of the basis function of the irreversible MCT. */ 45/* </summary> */ 46static const double mct_norms_real[3] = { 1.732, 1.805, 1.573 }; 47 48/* <summary> */ 49/* Foward reversible MCT. */ 50/* </summary> */ 51void mct_encode( 52 int* restrict c0, 53 int* restrict c1, 54 int* restrict c2, 55 int n) 56{ 57 int i; 58 for(i = 0; i < n; ++i) { 59 int r = c0[i]; 60 int g = c1[i]; 61 int b = c2[i]; 62 int y = (r + (g * 2) + b) >> 2; 63 int u = b - g; 64 int v = r - g; 65 c0[i] = y; 66 c1[i] = u; 67 c2[i] = v; 68 } 69} 70 71/* <summary> */ 72/* Inverse reversible MCT. */ 73/* </summary> */ 74void mct_decode( 75 int* restrict c0, 76 int* restrict c1, 77 int* restrict c2, 78 int n) 79{ 80 int i; 81 for (i = 0; i < n; ++i) { 82 int y = c0[i]; 83 int u = c1[i]; 84 int v = c2[i]; 85 int g = y - ((u + v) >> 2); 86 int r = v + g; 87 int b = u + g; 88 c0[i] = r; 89 c1[i] = g; 90 c2[i] = b; 91 } 92} 93 94/* <summary> */ 95/* Get norm of basis function of reversible MCT. */ 96/* </summary> */ 97double mct_getnorm(int compno) { 98 return mct_norms[compno]; 99} 100 101/* <summary> */ 102/* Foward irreversible MCT. */ 103/* </summary> */ 104void mct_encode_real( 105 int* restrict c0, 106 int* restrict c1, 107 int* restrict c2, 108 int n) 109{ 110 int i; 111 for(i = 0; i < n; ++i) { 112 int r = c0[i]; 113 int g = c1[i]; 114 int b = c2[i]; 115 int y = fix_mul(r, 2449) + fix_mul(g, 4809) + fix_mul(b, 934); 116 int u = -fix_mul(r, 1382) - fix_mul(g, 2714) + fix_mul(b, 4096); 117 int v = fix_mul(r, 4096) - fix_mul(g, 3430) - fix_mul(b, 666); 118 c0[i] = y; 119 c1[i] = u; 120 c2[i] = v; 121 } 122} 123 124/* <summary> */ 125/* Inverse irreversible MCT. */ 126/* </summary> */ 127void mct_decode_real( 128 float* restrict c0, 129 float* restrict c1, 130 float* restrict c2, 131 int n) 132{ 133 int i; 134#ifdef __SSE__ 135 __m128 vrv, vgu, vgv, vbu; 136 vrv = _mm_set1_ps(1.402f); 137 vgu = _mm_set1_ps(0.34413f); 138 vgv = _mm_set1_ps(0.71414f); 139 vbu = _mm_set1_ps(1.772f); 140 for (i = 0; i < (n >> 3); ++i) { 141 __m128 vy, vu, vv; 142 __m128 vr, vg, vb; 143 144 vy = _mm_load_ps(c0); 145 vu = _mm_load_ps(c1); 146 vv = _mm_load_ps(c2); 147 vr = _mm_add_ps(vy, _mm_mul_ps(vv, vrv)); 148 vg = _mm_sub_ps(_mm_sub_ps(vy, _mm_mul_ps(vu, vgu)), _mm_mul_ps(vv, vgv)); 149 vb = _mm_add_ps(vy, _mm_mul_ps(vu, vbu)); 150 _mm_store_ps(c0, vr); 151 _mm_store_ps(c1, vg); 152 _mm_store_ps(c2, vb); 153 c0 += 4; 154 c1 += 4; 155 c2 += 4; 156 157 vy = _mm_load_ps(c0); 158 vu = _mm_load_ps(c1); 159 vv = _mm_load_ps(c2); 160 vr = _mm_add_ps(vy, _mm_mul_ps(vv, vrv)); 161 vg = _mm_sub_ps(_mm_sub_ps(vy, _mm_mul_ps(vu, vgu)), _mm_mul_ps(vv, vgv)); 162 vb = _mm_add_ps(vy, _mm_mul_ps(vu, vbu)); 163 _mm_store_ps(c0, vr); 164 _mm_store_ps(c1, vg); 165 _mm_store_ps(c2, vb); 166 c0 += 4; 167 c1 += 4; 168 c2 += 4; 169 } 170 n &= 7; 171#endif 172 for(i = 0; i < n; ++i) { 173 float y = c0[i]; 174 float u = c1[i]; 175 float v = c2[i]; 176 float r = y + (v * 1.402f); 177 float g = y - (u * 0.34413f) - (v * (0.71414f)); 178 float b = y + (u * 1.772f); 179 c0[i] = r; 180 c1[i] = g; 181 c2[i] = b; 182 } 183} 184 185/* <summary> */ 186/* Get norm of basis function of irreversible MCT. */ 187/* </summary> */ 188double mct_getnorm_real(int compno) { 189 return mct_norms_real[compno]; 190} 191