Source file src/pkg/image/jpeg/idct.go
1 // Copyright 2009 The Go Authors. All rights reserved. 2 // Use of this source code is governed by a BSD-style 3 // license that can be found in the LICENSE file. 4 5 package jpeg 6 7 // This is a Go translation of idct.c from 8 // 9 // http://standards.iso.org/ittf/PubliclyAvailableStandards/ISO_IEC_13818-4_2004_Conformance_Testing/Video/verifier/mpeg2decode_960109.tar.gz 10 // 11 // which carries the following notice: 12 13 /* Copyright (C) 1996, MPEG Software Simulation Group. All Rights Reserved. */ 14 15 /* 16 * Disclaimer of Warranty 17 * 18 * These software programs are available to the user without any license fee or 19 * royalty on an "as is" basis. The MPEG Software Simulation Group disclaims 20 * any and all warranties, whether express, implied, or statuary, including any 21 * implied warranties or merchantability or of fitness for a particular 22 * purpose. In no event shall the copyright-holder be liable for any 23 * incidental, punitive, or consequential damages of any kind whatsoever 24 * arising from the use of these programs. 25 * 26 * This disclaimer of warranty extends to the user of these programs and user's 27 * customers, employees, agents, transferees, successors, and assigns. 28 * 29 * The MPEG Software Simulation Group does not represent or warrant that the 30 * programs furnished hereunder are free of infringement of any third-party 31 * patents. 32 * 33 * Commercial implementations of MPEG-1 and MPEG-2 video, including shareware, 34 * are subject to royalty fees to patent holders. Many of these patents are 35 * general enough such that they are unavoidable regardless of implementation 36 * design. 37 * 38 */ 39 40 const ( 41 w1 = 2841 // 2048*sqrt(2)*cos(1*pi/16) 42 w2 = 2676 // 2048*sqrt(2)*cos(2*pi/16) 43 w3 = 2408 // 2048*sqrt(2)*cos(3*pi/16) 44 w5 = 1609 // 2048*sqrt(2)*cos(5*pi/16) 45 w6 = 1108 // 2048*sqrt(2)*cos(6*pi/16) 46 w7 = 565 // 2048*sqrt(2)*cos(7*pi/16) 47 48 w1pw7 = w1 + w7 49 w1mw7 = w1 - w7 50 w2pw6 = w2 + w6 51 w2mw6 = w2 - w6 52 w3pw5 = w3 + w5 53 w3mw5 = w3 - w5 54 55 r2 = 181 // 256/sqrt(2) 56 ) 57 58 // idct performs a 2-D Inverse Discrete Cosine Transformation, followed by a 59 // +128 level shift and a clip to [0, 255], writing the results to dst. 60 // stride is the number of elements between successive rows of dst. 61 // 62 // The input coefficients should already have been multiplied by the 63 // appropriate quantization table. We use fixed-point computation, with the 64 // number of bits for the fractional component varying over the intermediate 65 // stages. 66 // 67 // For more on the actual algorithm, see Z. Wang, "Fast algorithms for the 68 // discrete W transform and for the discrete Fourier transform", IEEE Trans. on 69 // ASSP, Vol. ASSP- 32, pp. 803-816, Aug. 1984. 70 func idct(dst []byte, stride int, src *block) { 71 // Horizontal 1-D IDCT. 72 for y := 0; y < 8; y++ { 73 // If all the AC components are zero, then the IDCT is trivial. 74 if src[y*8+1] == 0 && src[y*8+2] == 0 && src[y*8+3] == 0 && 75 src[y*8+4] == 0 && src[y*8+5] == 0 && src[y*8+6] == 0 && src[y*8+7] == 0 { 76 dc := src[y*8+0] << 3 77 src[y*8+0] = dc 78 src[y*8+1] = dc 79 src[y*8+2] = dc 80 src[y*8+3] = dc 81 src[y*8+4] = dc 82 src[y*8+5] = dc 83 src[y*8+6] = dc 84 src[y*8+7] = dc 85 continue 86 } 87 88 // Prescale. 89 x0 := (src[y*8+0] << 11) + 128 90 x1 := src[y*8+4] << 11 91 x2 := src[y*8+6] 92 x3 := src[y*8+2] 93 x4 := src[y*8+1] 94 x5 := src[y*8+7] 95 x6 := src[y*8+5] 96 x7 := src[y*8+3] 97 98 // Stage 1. 99 x8 := w7 * (x4 + x5) 100 x4 = x8 + w1mw7*x4 101 x5 = x8 - w1pw7*x5 102 x8 = w3 * (x6 + x7) 103 x6 = x8 - w3mw5*x6 104 x7 = x8 - w3pw5*x7 105 106 // Stage 2. 107 x8 = x0 + x1 108 x0 -= x1 109 x1 = w6 * (x3 + x2) 110 x2 = x1 - w2pw6*x2 111 x3 = x1 + w2mw6*x3 112 x1 = x4 + x6 113 x4 -= x6 114 x6 = x5 + x7 115 x5 -= x7 116 117 // Stage 3. 118 x7 = x8 + x3 119 x8 -= x3 120 x3 = x0 + x2 121 x0 -= x2 122 x2 = (r2*(x4+x5) + 128) >> 8 123 x4 = (r2*(x4-x5) + 128) >> 8 124 125 // Stage 4. 126 src[8*y+0] = (x7 + x1) >> 8 127 src[8*y+1] = (x3 + x2) >> 8 128 src[8*y+2] = (x0 + x4) >> 8 129 src[8*y+3] = (x8 + x6) >> 8 130 src[8*y+4] = (x8 - x6) >> 8 131 src[8*y+5] = (x0 - x4) >> 8 132 src[8*y+6] = (x3 - x2) >> 8 133 src[8*y+7] = (x7 - x1) >> 8 134 } 135 136 // Vertical 1-D IDCT. 137 for x := 0; x < 8; x++ { 138 // Similar to the horizontal 1-D IDCT case, if all the AC components are zero, then the IDCT is trivial. 139 // However, after performing the horizontal 1-D IDCT, there are typically non-zero AC components, so 140 // we do not bother to check for the all-zero case. 141 142 // Prescale. 143 y0 := (src[8*0+x] << 8) + 8192 144 y1 := src[8*4+x] << 8 145 y2 := src[8*6+x] 146 y3 := src[8*2+x] 147 y4 := src[8*1+x] 148 y5 := src[8*7+x] 149 y6 := src[8*5+x] 150 y7 := src[8*3+x] 151 152 // Stage 1. 153 y8 := w7*(y4+y5) + 4 154 y4 = (y8 + w1mw7*y4) >> 3 155 y5 = (y8 - w1pw7*y5) >> 3 156 y8 = w3*(y6+y7) + 4 157 y6 = (y8 - w3mw5*y6) >> 3 158 y7 = (y8 - w3pw5*y7) >> 3 159 160 // Stage 2. 161 y8 = y0 + y1 162 y0 -= y1 163 y1 = w6*(y3+y2) + 4 164 y2 = (y1 - w2pw6*y2) >> 3 165 y3 = (y1 + w2mw6*y3) >> 3 166 y1 = y4 + y6 167 y4 -= y6 168 y6 = y5 + y7 169 y5 -= y7 170 171 // Stage 3. 172 y7 = y8 + y3 173 y8 -= y3 174 y3 = y0 + y2 175 y0 -= y2 176 y2 = (r2*(y4+y5) + 128) >> 8 177 y4 = (r2*(y4-y5) + 128) >> 8 178 179 // Stage 4. 180 src[8*0+x] = (y7 + y1) >> 14 181 src[8*1+x] = (y3 + y2) >> 14 182 src[8*2+x] = (y0 + y4) >> 14 183 src[8*3+x] = (y8 + y6) >> 14 184 src[8*4+x] = (y8 - y6) >> 14 185 src[8*5+x] = (y0 - y4) >> 14 186 src[8*6+x] = (y3 - y2) >> 14 187 src[8*7+x] = (y7 - y1) >> 14 188 } 189 190 // Level shift by +128, clip to [0, 255], and write to dst. 191 for y := 0; y < 8; y++ { 192 for x := 0; x < 8; x++ { 193 c := src[y*8+x] 194 if c < -128 { 195 c = 0 196 } else if c > 127 { 197 c = 255 198 } else { 199 c += 128 200 } 201 dst[y*stride+x] = uint8(c) 202 } 203 } 204 }