• 理解DCT与DST【一】：离散傅里叶变换
• 理解DCT与DST【二】：离散余弦变换

## 3 离散正弦变换（Discrete Sine Transform, DST）

### 3.2 DST-I

DST-I 延拓的公式化表示如下，

### 3.4 DST-III

DST-III 延拓的公式化表示如下，

### 3.6 DST-V

DST-V 延拓的公式化表示如下，

DST-I

DST-II

DST-III

DST-IV

DST-V

DST-VI

DST-VII

DST-VIII

## DST 基图像生成脚本

``````# -*- coding: utf-8 -*-
from __future__ import division
import numpy as np
import matplotlib.pyplot as plt

def dst1_kern(N):
x = np.zeros([N, N])
for k in range(N):
x[k, :] = np.sin(np.pi * (k + 1) * (np.arange(N) + 1) / (N + 1))
x *= np.sqrt(2. / (N + 1))
return x

def dst2_kern(N):
x = np.zeros([N, N])
for k in range(N):
x[k, :] = np.sin(np.pi * (k + 1) * (2 * np.arange(N) + 1) / (2 * N))
x[-1, :] /= np.sqrt(2)
x *= np.sqrt(2. / N)
return x

def dst3_kern(N):
x = np.zeros([N, N])
for k in range(N):
x[k] = np.sin(np.pi * (2 * k + 1) * (np.arange(N) + 1) / (2 * N))
x[:, -1] /= np.sqrt(2)
x *= np.sqrt(2. / N)
return x

def dst4_kern(N):
x = np.zeros([N, N])
for k in range(N):
x[k] = np.sin(np.pi * (2 * k + 1) * (2 * np.arange(N) + 1) / (4 * N))
x *= np.sqrt(2. / N)
return x

def dst5_kern(N):
x = np.zeros([N, N])
for k in range(N):
x[k] = np.sin(2 * np.pi * (k + 1) * (np.arange(N) + 1) / (2 * N + 1))
x *= np.sqrt(4. / (2 * N + 1))
return x

def dst6_kern(N):
x = np.zeros([N, N])
for k in range(N):
x[k] = np.sin(np.pi * (k + 1) * (2 * np.arange(N) + 1) / (2 * N + 1))
x *= np.sqrt(4. / (2 * N + 1))
return x

def dst7_kern(N):
x = np.zeros([N, N])
for k in range(N):
x[k] = np.sin(np.pi * (2 * k + 1) * (np.arange(N) + 1) / (2 * N + 1))
x *= np.sqrt(4. / (2 * N + 1))
return x

def dst8_kern(N):
x = np.zeros([N, N])
for k in range(N):
x[k] = np.sin(np.pi * (2 * k + 1) * (2 * np.arange(N) + 1) / (4 * N - 2))
x[:, -1] /= np.sqrt(2)
x[-1, :] /= np.sqrt(2)
x *= np.sqrt(4. / (2 * N - 1))
return x

if __name__ == '__main__':
path = 'dst_img/basis/'
num = ['I', 'II', 'III', 'IV', 'V', 'VI', 'VII', 'VIII']

for n in range(1, 9):
kern = globals().get('dst{}_kern'.format(n))(8)
'''
把变换矩阵按行分块，记第k行的转置为ak，那么 h=[a1.T, a2.T, ..., a8.T]，
用矩阵分块的知识得到 v.dot(h) 的每一个子块是 ai.dot(aj.T)，
这样一次性就可以得到64个基图像拼起来的64x64的总图像了。
'''
h = kern.reshape([1, -1])
v = h.T
img = v.dot(h)

plt.figure(figsize=(8, 8))
plt.imshow(img, 'gray')
plt.savefig(path+'/dst8x8_{}_all.png'.format(n), dpi='figure')

imgs = []
for i in range(8):
for j in range(8):
imgs.append(img[8*i:8*i+8, 8*j:8*j+8])
#imgs.append(kern[[i], :].T.dot(kern[[j], :])) # 这种方法就是基图像一个一个地求

imgs = np.stack(imgs)
vmin, vmax = np.min(imgs), np.max(imgs)

plt.figure(figsize=(32, 33))
for i in range(8):
for j in range(8):
plt.subplot(8, 8, i*8+j+1)
plt.imshow(imgs[i*8+j], 'gray', vmin=vmin, vmax=vmax)
plt.axis('off')
plt.tight_layout()
plt.subplots_adjust(0.02, 0.02, 0.98, 0.95, 0.1, 0.1)
plt.suptitle('DST-'+num[n-1]+' basis', fontsize=48)
plt.savefig(path+'/dst8x8_{}.png'.format(n), dpi=32)``````