在 Python 中实现二维离散余弦变换的问题答案

【问题标题】：Issiue with implementation of 2D Discrete Cosine Transform in Python在 Python 中实现二维离散余弦变换的问题
【发布时间】：2016-10-18 09:24:28
【问题描述】：

我正在尝试将 Zhao Koch 隐写方法从 matlab 重写为 python，但我一开始就卡住了。

matlab 中的前两个程序：

第 1 步：

A = imread(casepath); # Reading stegonography case image and aquiring it's RGB values. In my case it's a 400x400 PNG image, so it gives a 400x400x3 array.

第 2 步：

D = dct2(A(:,:,3)); # Applying 2D DCT to blue values of the image

Python 代码模拟：

from scipy import misc
from numpy import empty,arange,exp,real,imag,pi
from numpy.fft import rfft,irfft

arr = misc.imread('casepath')# 400x480x3 array (Step 1)
arr[20, 30, 2] # Getting blue pixel value

def dct(y): #Basic DCT build from numpy
    N = len(y)
    y2 = empty(2*N,float)
    y2[:N] = y[:]
    y2[N:] = y[::-1]

    c = rfft(y2)
    phi = exp(-1j*pi*arange(N)/(2*N))
    return real(phi*c[:N])


def dct2(y): #2D DCT bulid from numpy and using prvious DCT function
    M = y.shape[0]
    N = y.shape[1]
    a = empty([M,N],float)
    b = empty([M,N],float)

    for i in range(M):
        a[i,:] = dct(y[i,:])
    for j in range(N):
        b[:,j] = dct(a[:,j])

    return b

 D = dct2(arr) # step 2 anlogue

但是，当我尝试执行代码时，出现以下错误：

Traceback (most recent call last):
File "path to .py file", line 31, in <module>
D = dct2(arr)
File "path to .py file", line 25, in dct2
a[i,:] = dct(y[i,:])
File "path to .py file", line 10, in dct
y2[:N] = y[:]
ValueError: could not broadcast input array from shape (400,3) into shape (400)

也许有人可以向我解释我做错了什么？

附加信息：操作系统：Windows 10 专业版 64 位蟒蛇：2.7.12 scipy：0.18.1 麻木：1.11.2 枕头：3.4.1

【问题讨论】：

从 numpy 和 scipy 你可以直接访问dct()。您是否有兴趣滚动自己的缓慢方法？
谢谢，不知道刚刚试了一下。看起来效果很好。

标签： python matlab dct

【解决方案1】：

您的代码工作正常，但它被设计为只接受二维数组，就像 Matlab 中的dct2()。由于您的 arr 是一个 3D 数组，您想要这样做

D = dct2(arr[...,2])

正如我在评论中提到的，使用 scipy 包中的（快速）内置 dct() 代替或重新发明轮子。

我评论中链接中的代码有效地为您提供了这个：

import numpy as np
from scipy.fftpack import dct, idct

def dct2(block):
    return dct(dct(block.T, norm='ortho').T, norm='ortho')

def idct2(block):
    return idct(idct(block.T, norm='ortho').T, norm='ortho')

但我必须再次强调，您必须为每个颜色平面单独调用此函数。 Scipy 的dct() 很乐意接受任何 N 维数组，并将在最后一个轴上应用变换。由于这是您的颜色平面，而不是像素的行和列，因此您将得到错误的结果。是的，有一种方法可以使用 axis 输入参数来解决这个问题，但我不会不必要地过度复杂化这个答案。

关于这里涉及的各种 DCT 实现，如果您从上面的 sn-p 中省略 norm='ortho' 参数，您的版本和 scipy 的实现会给出相同的结果。但包含该参数后，scipy 的转换将与 Matlab 的一致。

【讨论】：