Skip to content

Instantly share code, notes, and snippets.

@vbkaisetsu

vbkaisetsu/dct-steganography.pyx

Last active June 1, 2017 13:58
Show Gist options
  • Save vbkaisetsu/10909e07ecdcafd35dd9f56ff46d3eae to your computer and use it in GitHub Desktop.
Save vbkaisetsu/10909e07ecdcafd35dd9f56ff46d3eae to your computer and use it in GitHub Desktop.
Download ZIP
DCT Steganography Implemented in Cython
Raw
dct-steganography.pyx
#
# Copyright 2017 Koichi Akabe
# License: GPLv3+
#
import numpy as np
cimport numpy as np
#DTYPE = np.int32
#ctypedef np.int32_t DTYPE_t
DTYPE = np.double
ctypedef np.double_t DTYPE_t
cdef np.ndarray[DTYPE_t, ndim=5] dct(np.ndarray[DTYPE_t, ndim=5] data):
cdef int size_x, size_y
cdef int x, y, i, s, ch
cdef DTYPE_t in0, in1, in2, in3, in4, in5, in6, in7
cdef DTYPE_t t0, t1, t2, t3, t4, t5, t6, t7
cdef DTYPE_t p0, p3
cdef DTYPE_t c0, c1, c2, c3
#cdef DTYPE_t r1 = 1.38703984532215 #np.cos(0.0625 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r2 = 1.30656296487638 #np.cos(0.1250 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r3 = 1.17587560241936 #np.cos(0.1875 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r5 = 0.785694958387102 #np.cos(0.3125 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r6 = 0.541196100146197 #np.cos(0.3750 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r7 = 0.275899379282943 #np.cos(0.4375 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t invsqrt2 = 0.707106781186547 #1.0 / np.sqrt(2.0)
cdef np.ndarray[DTYPE_t, ndim=5] result
cdef DTYPE_t *result_p
cdef DTYPE_t *data_p
cdef DTYPE_t *result_subst
cdef DTYPE_t *data_subst
size_y = data.shape[1]
size_x = data.shape[2]
result = DTYPE(np.zeros((3, size_y, size_x, 8, 8)))
result_p = <DTYPE_t*>result.data
data_p = <DTYPE_t*>data.data
for ch in range(3):
for y in range(size_y):
for x in range(size_x):
s = ((size_y * ch + y) * size_x + x) * 64
for i in range(8):
result_subst = &result_p[s + 8 * i]
data_subst = &data_p[s + 8 * i]
in0 = data_subst[0]; in1 = data_subst[1]; in2 = data_subst[2]; in3 = data_subst[3]
in4 = data_subst[4]; in5 = data_subst[5]; in6 = data_subst[6]; in7 = data_subst[7]
t0 = in0 + in7
t1 = in1 + in6
t2 = in2 + in5
t3 = in3 + in4
t7 = in0 - in7
t6 = in1 - in6
t5 = in2 - in5
t4 = in3 - in4
c0 = t0 + t3
c1 = t1 + t2
c3 = t0 - t3
c2 = t1 - t2
#result_subst[2] = c2 * r6 + c3 * r2
#result_subst[6] = c3 * r6 - c2 * r2
result_subst[2] = c2 * 0.541196100146197 + c3 * 1.30656296487638
result_subst[6] = c3 * 0.541196100146197 - c2 * 1.30656296487638
result_subst[0] = c0 + c1
result_subst[4] = c0 - c1
#c3 = t4 * r3 + t7 * r5
#c2 = t5 * r1 + t6 * r7
#c1 = t6 * r1 - t5 * r7
#c0 = t7 * r3 - t4 * r5
c3 = t4 * 1.17587560241936 + t7 * 0.785694958387102
c2 = t5 * 1.38703984532215 + t6 * 0.275899379282943
c1 = t6 * 1.38703984532215 - t5 * 0.275899379282943
c0 = t7 * 1.17587560241936 - t4 * 0.785694958387102
result_subst[5] = c3 - c1
result_subst[3] = c0 - c2
#p0 = (c0 + c2) * invsqrt2
#p3 = (c3 + c1) * invsqrt2
p0 = (c0 + c2) * 0.707106781186547
p3 = (c3 + c1) * 0.707106781186547
result_subst[1] = p0 + p3
result_subst[7] = p0 - p3
for i in range(8):
result_subst = &result_p[s + i]
in0 = result_subst[8 * 0]; in1 = result_subst[8 * 1]; in2 = result_subst[8 * 2]; in3 = result_subst[8 * 3]
in4 = result_subst[8 * 4]; in5 = result_subst[8 * 5]; in6 = result_subst[8 * 6]; in7 = result_subst[8 * 7]
t0 = in0 + in7
t1 = in1 + in6
t2 = in2 + in5
t3 = in3 + in4
t7 = in0 - in7
t6 = in1 - in6
t5 = in2 - in5
t4 = in3 - in4
c0 = t0 + t3
c1 = t1 + t2
c3 = t0 - t3
c2 = t1 - t2
#result_subst[8 * 2] = c2 * r6 + c3 * r2
#result_subst[8 * 6] = c3 * r6 - c2 * r2
result_subst[8 * 2] = c2 * 0.541196100146197 + c3 * 1.30656296487638
result_subst[8 * 6] = c3 * 0.541196100146197 - c2 * 1.30656296487638
result_subst[8 * 0] = c0 + c1
result_subst[8 * 4] = c0 - c1
#c3 = t4 * r3 + t7 * r5
#c2 = t5 * r1 + t6 * r7
#c1 = t6 * r1 - t5 * r7
#c0 = t7 * r3 - t4 * r5
c3 = t4 * 1.17587560241936 + t7 * 0.785694958387102
c2 = t5 * 1.38703984532215 + t6 * 0.275899379282943
c1 = t6 * 1.38703984532215 - t5 * 0.275899379282943
c0 = t7 * 1.17587560241936 - t4 * 0.785694958387102
result_subst[8 * 5] = c3 - c1
result_subst[8 * 3] = c0 - c2
#p0 = (c0 + c2) * invsqrt2
#p3 = (c3 + c1) * invsqrt2
p0 = (c0 + c2) * 0.707106781186547
p3 = (c3 + c1) * 0.707106781186547
result_subst[8 * 1] = p0 + p3
result_subst[8 * 7] = p0 - p3
result_subst = &result_p[s]
for i in range(64):
result_subst[i] *= 0.125
return result
cdef np.ndarray[DTYPE_t, ndim=5] idct(np.ndarray[DTYPE_t, ndim=5] data):
cdef int w, h, size_x, size_y
cdef int x, y, i, s, ch
cdef DTYPE_t in0, in1, in2, in3, in4, in5, in6, in7
cdef DTYPE_t t0, t1, t2, t3, t4, t5, t6, t7
cdef DTYPE_t p0, p3
cdef DTYPE_t c0, c1, c2, c3
#cdef DTYPE_t r1 = 1.38703984532215 #np.cos(0.0625 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r2 = 1.30656296487638 #np.cos(0.1250 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r3 = 1.17587560241936 #np.cos(0.1875 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r5 = 0.785694958387102 #np.cos(0.3125 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r6 = 0.541196100146197 #np.cos(0.3750 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t r7 = 0.275899379282943 #np.cos(0.4375 * np.pi) * np.sqrt(2.0)
#cdef DTYPE_t invsqrt2 = 0.707106781186547 #1.0 / np.sqrt(2.0)
cdef np.ndarray[DTYPE_t, ndim=5] result
cdef DTYPE_t *result_p
cdef DTYPE_t *data_p
cdef DTYPE_t *result_subst
cdef DTYPE_t *data_subst
size_y = data.shape[1]
size_x = data.shape[2]
result = DTYPE(np.zeros((3, size_y, size_x, 8, 8)))
result_p = <DTYPE_t*>result.data
data_p = <DTYPE_t*>data.data
for ch in range(3):
for y in range(size_y):
for x in range(size_x):
s = ((size_y * ch + y) * size_x + x) * 64
for i in range(8):
result_subst = &result_p[s + 8 * i]
data_subst = &data_p[s + 8 * i]
in0 = data_subst[0]; in1 = data_subst[1]; in2 = data_subst[2]; in3 = data_subst[3]
in4 = data_subst[4]; in5 = data_subst[5]; in6 = data_subst[6]; in7 = data_subst[7]
#p0 = (in1 + in7) * invsqrt2
#p3 = (in1 - in7) * invsqrt2
p0 = (in1 + in7) * 0.707106781186547
p3 = (in1 - in7) * 0.707106781186547
c3 = p3 + in5
c0 = p0 + in3
c1 = p3 - in5
c2 = p0 - in3
#t7 = c3 * r5 + c0 * r3
#t6 = c2 * r7 + c1 * r1
#t4 = c3 * r3 - c0 * r5
#t5 = c2 * r1 - c1 * r7
t7 = c3 * 0.785694958387102 + c0 * 1.17587560241936
t6 = c2 * 0.275899379282943 + c1 * 1.38703984532215
t4 = c3 * 1.17587560241936 - c0 * 0.785694958387102
t5 = c2 * 1.38703984532215 - c1 * 0.275899379282943
#c3 = in2 * r2 + in6 * r6
#c2 = in2 * r6 - in6 * r2
c3 = in2 * 1.30656296487638 + in6 * 0.541196100146197
c2 = in2 * 0.541196100146197 - in6 * 1.30656296487638
c0 = in0 + in4
c1 = in0 - in4
t0 = c0 + c3
t1 = c1 + c2
t3 = c0 - c3
t2 = c1 - c2
result_subst[0] = t0 + t7
result_subst[1] = t1 + t6
result_subst[2] = t2 + t5
result_subst[3] = t3 + t4
result_subst[7] = t0 - t7
result_subst[6] = t1 - t6
result_subst[5] = t2 - t5
result_subst[4] = t3 - t4
for i in range(8):
result_subst = &result_p[s + i]
in0 = result_subst[8 * 0]; in1 = result_subst[8 * 1]; in2 = result_subst[8 * 2]; in3 = result_subst[8 * 3]
in4 = result_subst[8 * 4]; in5 = result_subst[8 * 5]; in6 = result_subst[8 * 6]; in7 = result_subst[8 * 7]
#p0 = (in1 + in7) * invsqrt2
#p3 = (in1 - in7) * invsqrt2
p0 = (in1 + in7) * 0.707106781186547
p3 = (in1 - in7) * 0.707106781186547
c3 = p3 + in5
c0 = p0 + in3
c1 = p3 - in5
c2 = p0 - in3
#t7 = c3 * r5 + c0 * r3
#t6 = c2 * r7 + c1 * r1
#t4 = c3 * r3 - c0 * r5
#t5 = c2 * r1 - c1 * r7
t7 = c3 * 0.785694958387102 + c0 * 1.17587560241936
t6 = c2 * 0.275899379282943 + c1 * 1.38703984532215
t4 = c3 * 1.17587560241936 - c0 * 0.785694958387102
t5 = c2 * 1.38703984532215 - c1 * 0.275899379282943
#c3 = in2 * r2 + in6 * r6
#c2 = in2 * r6 - in6 * r2
c3 = in2 * 1.30656296487638 + in6 * 0.541196100146197
c2 = in2 * 0.541196100146197 - in6 * 1.30656296487638
c0 = in0 + in4
c1 = in0 - in4
t0 = c0 + c3
t1 = c1 + c2
t3 = c0 - c3
t2 = c1 - c2
result_subst[8 * 0] = t0 + t7
result_subst[8 * 1] = t1 + t6
result_subst[8 * 2] = t2 + t5
result_subst[8 * 3] = t3 + t4
result_subst[8 * 7] = t0 - t7
result_subst[8 * 6] = t1 - t6
result_subst[8 * 5] = t2 - t5
result_subst[8 * 4] = t3 - t4
result_subst = &result_p[s]
for i in range(64):
result_subst[i] *= 0.125
return result
def encode(np.ndarray[np.uint8_t, ndim=3] image, bytes hidden_data, int intensity):
cdef int i, x, y, w, h, ch, datasize
cdef DTYPE_t a, b, tmp
cdef np.ndarray[DTYPE_t, ndim=3] image_arr = DTYPE(image)
cdef np.ndarray[DTYPE_t, ndim=5] image_arr_trans
cdef np.ndarray[DTYPE_t, ndim=5] image_dct
cdef DTYPE_t *image_dct_p
cdef unsigned char *data = hidden_data
h = image.shape[0]
w = image.shape[1]
if h % 8 != 0 or w % 8 != 0:
return None
datasize = len(hidden_data) * 8
image_arr_trans = image_arr.reshape(h // 8, 8, -1, 8, 3).transpose(4, 0, 2, 1, 3).copy()
image_dct = dct(image_arr_trans)
image_dct_p = <DTYPE_t*>image_dct.data
for ch in range(3):
i = 0
for y in range(h / 8):
for x in range(w / 8):
a = image_dct[ch, y, x, 6, 4]
b = image_dct[ch, y, x, 4, 2]
#print(data[i / 8] >> i % 8 & 1)
if data[i / 8] >> i % 8 & 1 == 0:
if a > b:
tmp = a
a = b
b = tmp
tmp = max(0, (intensity - (b - a)) / 2)
a -= tmp
b += tmp
else:
if a < b:
tmp = a
a = b
b = tmp
tmp = max(0, (intensity - (a - b)) / 2)
a += tmp
b -= tmp
#print(intensity)
#print(a, b)
#import sys
#sys.exit()
image_dct[ch, y, x, 6, 4] = a
image_dct[ch, y, x, 4, 2] = b
i += 1
if i >= datasize:
i = 0
image_arr_trans = idct(image_dct)
image_arr = image_arr_trans.transpose(1, 3, 2, 4, 0).reshape(h, w, 3).copy()
image_arr[image_arr < 0] = 0
image_arr[image_arr > 255] = 255
#print(np.max(image_arr))
#print(np.min(image_arr))
#import sys
#sys.exit()
return np.uint8(image_arr)
def decode(np.ndarray[np.uint8_t, ndim=3] image, int string_size=0):
cdef int i, x, y, w, h, ch
cdef DTYPE_t a, b
cdef np.ndarray[DTYPE_t, ndim=3] image_arr = DTYPE(image)
cdef np.ndarray[DTYPE_t, ndim=5] image_arr_trans
cdef np.ndarray[DTYPE_t, ndim=5] image_dct
cdef np.ndarray[np.int32_t, ndim=1] result_tmp
cdef np.ndarray[np.uint8_t, ndim=1] result
h = image.shape[0]
w = image.shape[1]
if h % 8 != 0 or w % 8 != 0:
return None
if string_size > 0:
result_tmp = np.int32(np.zeros(string_size * 8))
else:
result_tmp = np.int32(np.zeros(h * w / 64))
image_arr_trans = image_arr.reshape(h // 8, 8, -1, 8, 3).transpose(4, 0, 2, 1, 3).copy()
image_dct = dct(image_arr_trans)
for ch in range(3):
i = 0
for y in range(h / 8):
for x in range(w / 8):
a = image_dct[ch, y, x, 6, 4]
b = image_dct[ch, y, x, 4, 2]
if a > b:
result_tmp[i] += 1
elif a < b:
result_tmp[i] -= 1
i += 1
if i >= len(result_tmp):
i = 0
result = np.uint8(np.zeros(len(result_tmp) / 8))
for i in range(len(result_tmp)):
if result_tmp[i] > 0:
result[i / 8] |= 1 << i % 8
return result.tobytes()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment