Issiue with implementation of 2D Discrete Cosine Transform in Python

Question

I'm trying to rewrite Zhao Koch steganography method from matlab into python and I am stuck right at the start.

The first two procedures as they are in matlab:

Step 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.

Step 2:

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

Python code analog:

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

However, when I try to execute the code I get the following error:

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)

Perhaps someone could kindly explain to me what am I doing wrong?

Additional Info: OS: Windows 10 Pro 64 bit Python: 2.7.12 scipy:0.18.1 numpy:1.11.2 pillow: 3.4.1

From numpy and scipy you have direct access to [`dct()`](http://stackoverflow.com/questions/34890585/in-scipy-why-doesnt-idctdcta-equal-to-a). Are you interested instead in rolling your own, slow approah? — Reti43, Oct 18 '16 at 10:09
Thank you, didn't know about it just tried it. Seems it works just fine. — Caiphas Kain, Oct 18 '16 at 10:19

Reti43 · Accepted Answer · 2016-10-18T13:26:56.990

Your code works fine, but it is designed to only accept a 2D array, just like dct2() in Matlab. Since your arr is a 3D array, you want to do

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

As mentioned in my comment, instead or reinventing the wheel, use the (fast) built-in dct() from the scipy package.

The code from the link in my comment effectively provides you this:

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')

But again, I must stress that you have to call this function for each colour plane individually. Scipy's dct() will happily accept any N-dimensional array and will apply the transform on the last axis. Since that's your colour planes and not your rows and columns of your pixels, you'll get the wrong result. Yes, there is a way to address this with the axis input parameter, but I won't unnecessarily overcomplicate this answer.

Regarding the various DCT implementations involved here, your version and scipy's implementation give the same result if you omit the norm='ortho' parameter from the snippet above. But with that parameter included, scipy's transform will agree with Matlab's.

Issiue with implementation of 2D Discrete Cosine Transform in Python

1 Answers1