【離散コサイン変換】DFTを使ってDCTしてみる

背景

【離散コサイン変換】DCT-IIの式をぼんやり理解するの続きです。

使用したライブラリ

【画像処理100本ノックに挑戦】独自の離散フーリエ変換クラスを作る

DFTを使ってDCT-IIしてみる

$N$点のデータ$x_n$があったとして、$4N$点の$y_n$を

\begin{eqnarray}
y_{2n+1} &=& y_{4N-2n-1}=x_n\;(0\leq n < N)\\
\end{eqnarray}

で定義してやってDFTすればOKです。上の点以外はゼロとします。

std::vector<double> dct(std::vector<double> in)
{
    int N = in.size();
    std::vector<double> out(N);
    for (int k = 0; k < N; k++)
    {
        out[k] = 0;
        for (int n = 0; n < N; n++)
        {
            out[k] += in[n] * cos(pi / N * (n + 0.5) * k);
        }
    }
    return out;
}

int main()
{
    int N = 128;
    std::vector<double> data(N), out(N);
    for (int i = 0; i < N; i++)
    {
        data[i] = cos(2 * pi * i / 10);
    }

    out = dct(data);

    std::vector<std::complex<double>> data2(4 * N);
    for (int i = 0; i < 4 * N; i++) data2[i] = 0;
    for (int i = 0; i < N; i++) data2[2 * i + 1] = data2[4 * N - 2 * i - 1] = data[i] / 2.;

    DFT dft;
    std::vector<std::complex<double>> out2 = dft.forward(data2);

    for (int i = 0; i < N; i++)
    {
        std::cout << out[i] << " " << out2[i].real() << std::endl;
    }

    return 0;
}

定義どおりにDCTした場合と結果が一致することが確認できました。

DFTを使ってDCT-IIIしてみる

$N$点のデータ$x_n$があったとして、$4N$点の$y_n$を

\begin{eqnarray}
0\leq n<N&\to& y_n = x_n\\
n=N&\to& y_n = 0\\
N<n\leq 2N &\to& y_n = -x_{2N-n}\\
2N<n<3N &\to& y_n = -x_{n-2N}\\
n=3N &\to& y_n = 0\\
3N<n<4N &\to& y_n = x_{4N-n}
\end{eqnarray}

で定義してやってDFTすればOKです。

std::vector<double> dct3(std::vector<double> in)
{
    int N = in.size();
    std::vector<double> out(N);
    for (int k = 0; k < N; k++)
    {
        out[k] = in[0] / 2.;
        for (int n = 1; n < N; n++)
        {
            out[k] += in[n] * cos(pi / N * (k + 0.5) * n);
        }
    }
    return out;
}

int main()
{
    int N = 128;
    std::vector<double> data(N), out(N);
    for (int i = 0; i < N; i++)
    {
        data[i] = cos(2 * pi * i / 10);
    }

    out = dct3(data);

    std::vector<std::complex<double>> data2(4 * N);
    for (int i = 0; i < 4 * N; i++) data2[i] = 0;
    for (int n = 0; n < N; n++) data2[n] = data[n];
    for (int n = N + 1; n <= 2 * N; n++) data2[n] = -data[2 * N - n];
    for (int n = 2 * N + 1; n < 3 * N; n++) data2[n] = -data[n - 2 * N];
    for (int n = 3 * N + 1; n < 4 * N; n++) data2[n] = data[4 * N - n];

    DFT dft;
    std::vector<std::complex<double>> out2 = dft.forward(data2);

    for (int i = 0; i < N; i++)
    {
        std::cout << out[i] << " " << out2[2 * i + 1].real() / 4. << std::endl;
    }

    return 0;
}

定義どおりにDCTした場合と結果が一致することが確認できました。

次回は二次元DCTについて少し考えてみるつもりです。