개시하다

배경.

재밌어 보여서 이산 부립엽 변환의 변환 행렬을 만들고 싶어요.

목표

이산 부립엽 변환의 변환 행렬을 만들고 실제 변환과 반변환이 가능한지 확인합니다.

이산 부립엽 변환

이산 부립엽 변환에 대한 설명
또한wiki의 공식에 따라 변환된 코드는 이전wiki에 실시되었으니 참고하시기 바랍니다.

변환 행렬

우선, 이산 부립엽 변환의 방정식은 다음과 같다.
${\displaystyle F(t)=\sum_{x=0}^{N-1} f(x) e^{-i\frac{2\pi t x}{N}}\quad\quad }$
t와 x를 T, X로 바꿔 분위기에서의 느낌(X와 T 모두 0~N-1의 1차원 배열로 변화)
${\displaystyle F(T)= f(X) e^{-i\frac{2\pi T X}{N}}=f(X)M}$
이제 변환 매트릭스 함수 만들기↓
dft.py

import numpy as np
import cmath

def dft_matrix(N):
    A = np.arange(N)
    T = A.reshape(1, -1)
    X = A.reshape(-1, 1)
    M = cmath.e**(-1j * 2 * cmath.pi * T * X  / N)
    return M

N = 3의 예에서 각각 무엇을 나타낸다.

N = 3
A = np.arange(N)
T = A.reshape(1, -1)
X = A.reshape(-1, 1)
M = cmath.e**(-1j * 2 * cmath.pi * T * X  / N)
print(A)
print('------------------')
print(T * X)
print('------------------')
print(M)

실행하다 보니까 이런 느낌이에요.

[0 1 2]
------------------
[[0 0 0]
 [0 1 2]
 [0 2 4]]
------------------
[[ 1. +0.j         1. +0.j         1. +0.j       ]
 [ 1. +0.j        -0.5-0.8660254j -0.5+0.8660254j]
 [ 1. +0.j        -0.5+0.8660254j -0.5-0.8660254j]]

T * X N×N의 배열, (i, j)의 요소(T * X)[i][j]는 ${2^i×2^j}달러입니다.
그리고 ${M=e^{-i\racc{2\piTX}{N}의 나머지 부분M을 계산합니다.↓ 여기 부분

M = cmath.e**(-1j * 2 * cmath.pi * T * X  / N)

실험

변환원이 되는 파형y이 필요하기 때문에 적당히 인의 합성파로 제작한다.

import numpy as np
import matplotlib.pyplot as plt

import dft

x = np.linspace(-1, 1, 200)
y  = np.sin(2 * np.pi * x) + np.sin(2 * np.pi * 2.5 * x)

print('元波形')
plt.plot(x, y)
plt.show()

이런 느낌.
이산 부립엽 변환의 실현
샘플 수dft_matrix를 다음 생성된 함수에 전달하여 변환 행렬N을 얻습니다.
변환 매트릭스M를 사용하는 매트릭스 계산에서 변환M.행렬 연산은numbyy이기 때문에 사용합니다.

# 変換行列
N = len(y)
M = dft.dft_matrix(N)
# 波形 ⇒　M　⇒ 周波数表現
fy = np.dot(y, M)
print('変換後:実部')
plt.plot(fy.real)
plt.show()
print('変換後:虚部')
plt.plot(fy.imag)
plt.show()

변환 결과

이어서 반변환을 진행하는데 이때 M의 반행렬을 구해야 한다.
역 행렬 변환도 numbynumpy.dot이기 때문에 사용합니다.
그리고 변환된 numpy.linalg.inv를 이전과 같은 요령에 따라 반변환한다.결과를 원래의 파형과 비교해 보자.

inv_M =  np.linalg.inv(M)
# 周波数表現 ⇒　-M　⇒ 波形
yd = np.dot(fy, inv_M)
print('元波形')
plt.plot(x, y)
plt.show()
print('逆変換後')
plt.plot(yd.real)
plt.show()

너무 완벽해(·)ω・)b

최후

dft를 묻고 바로 이동하는 코드를 거기에 두세요.(numby와 matplotlib이 없을 때는 pip로 살짝 넣으세요)

import cmath

import numpy as np
import matplotlib.pyplot as plt


def dft_matrix(N):
    A = np.arange(N)
    T = A.reshape(1, -1)
    X = A.reshape(-1, 1)
    M = cmath.e**(-1j * 2 * cmath.pi * T * X  / N)
    return M

# 元波形生成部
x = np.linspace(-1, 1, 200)
y  = np.sin(2 * np.pi * x) + np.sin(2 * np.pi * 2.5 * x)
print('元波形')
plt.plot(x, y)
plt.show()

# 変換行列の作成と変換
M = dft_matrix(len(y))
fy = np.dot(y, M)
print('変換後:実部')
plt.plot(fy.real)
plt.show()
print('変換後:虚部')
plt.plot(fy.imag)
plt.show()

# 逆変換行列の作成と逆変換
inv_M =  np.linalg.inv(M)
yd = np.dot(fy, inv_M)
print('逆変換後')
plt.plot(yd.real)
plt.show()

총결산

이산 부립엽 변환의 변환 행렬을 만들어 실제 변환과 반변환이 가능하다는 것을 증명했다.
결과:재미있다ρ(・ω・)