|
2 | 2 | # The DCT converts spatial image data into frequency domain for compression. |
3 | 3 | module CrImage::JPEG |
4 | 4 | # Forward Discrete Cosine Transform (DCT) implementation |
5 | | - # Uses direct formula for correctness |
6 | | - # Operates on 8x8 blocks |
| 5 | + # Uses a separable 8x8 transform with precomputed coefficients. |
7 | 6 | module DCT |
| 7 | + SCALE = [ |
| 8 | + 1.0 / ::Math.sqrt(8.0), |
| 9 | + ::Math.sqrt(2.0 / 8.0), |
| 10 | + ::Math.sqrt(2.0 / 8.0), |
| 11 | + ::Math.sqrt(2.0 / 8.0), |
| 12 | + ::Math.sqrt(2.0 / 8.0), |
| 13 | + ::Math.sqrt(2.0 / 8.0), |
| 14 | + ::Math.sqrt(2.0 / 8.0), |
| 15 | + ::Math.sqrt(2.0 / 8.0), |
| 16 | + ] |
| 17 | + |
| 18 | + BASIS = [ |
| 19 | + [0.3535533905932738, 0.3535533905932738, 0.3535533905932738, 0.3535533905932738, 0.3535533905932738, 0.3535533905932738, 0.3535533905932738, 0.3535533905932738], |
| 20 | + [0.4903926402016152, 0.4157348061512726, 0.27778511650980114, 0.09754516100806417, -0.0975451610080641, -0.277785116509801, -0.4157348061512727, -0.4903926402016152], |
| 21 | + [0.46193976625564337, 0.19134171618254492, -0.19134171618254486, -0.46193976625564337, -0.4619397662556434, -0.19134171618254517, 0.191341716182545, 0.46193976625564326], |
| 22 | + [0.4157348061512726, -0.0975451610080641, -0.4903926402016152, -0.2777851165098011, 0.2777851165098009, 0.4903926402016153, 0.0975451610080644, -0.41573480615127256], |
| 23 | + [0.35355339059327384, -0.35355339059327373, -0.35355339059327384, 0.3535533905932737, 0.35355339059327384, -0.35355339059327334, -0.3535533905932733, 0.35355339059327323], |
| 24 | + [0.27778511650980114, -0.4903926402016152, 0.09754516100806415, 0.41573480615127273, -0.41573480615127256, -0.09754516100806429, 0.49039264020161516, -0.27778511650980076], |
| 25 | + [0.19134171618254495, -0.4619397662556434, 0.46193976625564326, -0.19134171618254528, -0.19134171618254495, 0.46193976625564315, -0.4619397662556437, 0.19134171618254314], |
| 26 | + [0.09754516100806417, -0.2777851165098011, 0.41573480615127273, -0.4903926402016153, 0.4903926402016152, -0.415734806151272, 0.27778511650980076, -0.09754516100806251], |
| 27 | + ] |
| 28 | + |
8 | 29 | # Perform 2D DCT on an 8x8 block |
9 | 30 | # Input: 64-element array of pixel values (0-255 range) |
10 | 31 | # Output: 64-element array of DCT coefficients |
11 | 32 | def self.transform(block : Array(Int32)) : Array(Int32) |
| 33 | + result = Array(Int32).new(64, 0) |
| 34 | + transform_into(block, result) |
| 35 | + result |
| 36 | + end |
| 37 | + |
| 38 | + def self.transform_into(block : Indexable(Int32), result : Array(Int32)) : Nil |
12 | 39 | raise ArgumentError.new("Block must have 64 elements") unless block.size == 64 |
| 40 | + raise ArgumentError.new("Result must have 64 elements") unless result.size == 64 |
13 | 41 |
|
14 | | - # Create a working copy and shift values from [0, 255] to [-128, 127] |
15 | | - temp = block.map { |val| (val - 128).to_f64 } |
| 42 | + shifted = StaticArray(Float64, 64).new(0.0) |
| 43 | + temp = StaticArray(Float64, 64).new(0.0) |
16 | 44 |
|
17 | | - # Apply 2D DCT |
18 | | - result = Array(Float64).new(64, 0.0) |
| 45 | + 64.times do |i| |
| 46 | + shifted[i] = (block[i] - 128).to_f64 |
| 47 | + end |
19 | 48 |
|
20 | | - 8.times do |ver| |
21 | | - 8.times do |hor| |
| 49 | + 8.times do |y| |
| 50 | + row_base = y * 8 |
| 51 | + 8.times do |u| |
22 | 52 | sum = 0.0 |
23 | | - 8.times do |ypos| |
24 | | - 8.times do |xpos| |
25 | | - sum += temp[ypos * 8 + xpos] * |
26 | | - ::Math.cos((2 * xpos + 1) * hor * ::Math::PI / 16.0) * |
27 | | - ::Math.cos((2 * ypos + 1) * ver * ::Math::PI / 16.0) |
28 | | - end |
| 53 | + basis = BASIS[u] |
| 54 | + 8.times do |x| |
| 55 | + sum += shifted[row_base + x] * basis[x] |
29 | 56 | end |
| 57 | + temp[row_base + u] = sum |
| 58 | + end |
| 59 | + end |
30 | 60 |
|
31 | | - # Apply normalization: (2/N) * C(u) * C(v) where N=8 |
32 | | - cu = hor == 0 ? 1.0 / ::Math.sqrt(2.0) : 1.0 |
33 | | - cv = ver == 0 ? 1.0 / ::Math.sqrt(2.0) : 1.0 |
34 | | - result[ver * 8 + hor] = sum * 0.25 * cu * cv # 2/8 = 0.25 |
| 61 | + 8.times do |v| |
| 62 | + 8.times do |u| |
| 63 | + sum = 0.0 |
| 64 | + 8.times do |y| |
| 65 | + sum += BASIS[v][y] * temp[y * 8 + u] |
| 66 | + end |
| 67 | + result[v * 8 + u] = sum.round.to_i32 |
35 | 68 | end |
36 | 69 | end |
37 | 70 |
|
38 | | - result.map(&.round.to_i32) |
| 71 | + nil |
39 | 72 | end |
40 | 73 | end |
41 | 74 | end |
0 commit comments