离散余弦变换一维Matlab
Discrete Cosine Transformation 1D Matlab
我正在尝试使用矩阵向量乘法在 MATLAB 中实现 DCT。具体来说,我想创建系数的 DCT 矩阵,然后用它与一维信号相乘来计算一维 DCT。
这是我生成 DCT 矩阵的代码:
function D=dct1d(n)
for u=0:n-1
if u==0
au=sqrt(1/n);
else
au=sqrt(2/n);
end
for x=0:n-1
D(u+1,x+1)=au*cos(((2*x+1)*u*pi)/2*n);
end
end
在此之后,我尝试使用 x = [1 2 3 4 5 6 7 8]
:
的测试信号执行 DCT
x=[1 2 3 4 5 6 7 8];
y=dct1(8)*x';
它给出的答案是:
12.7279
18.0000
18.0000
18.0000
18.0000
18.0000
18.0000
18.0000
然而,正确答案是:
12.7279
-6.4423
-0.0000
-0.6735
0
-0.2009
-0.0000
-0.0507
您的代码中的错误非常轻微但很重要。您计算系数的行:
D(u+1,x+1)=au*cos(((2*x+1)*u*pi)/2*n);
看看这行的最后一部分:
D(u+1,x+1)=au*cos(((2*x+1)*u*pi)/2*n);
%// ^^^
Because multiplication and division are equal in precedence,这和做完全一样:
D(u+1,x+1)=au*cos((((2*x+1)*u*pi)/2)*n);
因此,您没有除以 2n
。你除以 2 然后 乘以 n
这是不正确的。您只需要用方括号将 2*n
操作括起来:
D(u+1,x+1)=au*cos(((2*x+1)*u*pi)/(2*n));
一旦你这样做了,我们就得到了正确的 DCT 矩阵。顺便说一句,检查您是否有正确答案的一种方法是使用 dctmtx
function, which calculates the N x N
DCT coefficient matrix that you're seeking. However, this function is part of the Signal Processing Toolbox, so if you don't have that then you unfortunately can't use this function, but if I can suggest an alternative answer instead of using for
loops, I would build a 2D grid of coordinates with meshgrid
,然后计算元素乘积。
这样的东西会起作用:
function D = dct1d(n)
[x,u] = meshgrid(0:n-1);
D = sqrt(2/n)*cos(((2*x+1).*u*pi)/(2*n));
D(1,:) = D(1,:) / sqrt(2);
end
不是使用 if
语句来确定每行需要应用的权重,我们可以只使用 sqrt(2/n)
然后除以第一行的 sqrt(2)
这样你而是除以 sqrt(1/n)
。此代码应产生与更正代码相同的结果1.
无论如何,一旦我做了这些更正,我就比较了你的代码给出的答案和 dctmtx
给出的答案,它是正确的:
>> dct1d(8)
ans =
0.3536 0.3536 0.3536 0.3536 0.3536 0.3536 0.3536 0.3536
0.4904 0.4157 0.2778 0.0975 -0.0975 -0.2778 -0.4157 -0.4904
0.4619 0.1913 -0.1913 -0.4619 -0.4619 -0.1913 0.1913 0.4619
0.4157 -0.0975 -0.4904 -0.2778 0.2778 0.4904 0.0975 -0.4157
0.3536 -0.3536 -0.3536 0.3536 0.3536 -0.3536 -0.3536 0.3536
0.2778 -0.4904 0.0975 0.4157 -0.4157 -0.0975 0.4904 -0.2778
0.1913 -0.4619 0.4619 -0.1913 -0.1913 0.4619 -0.4619 0.1913
0.0975 -0.2778 0.4157 -0.4904 0.4904 -0.4157 0.2778 -0.0975
>> dctmtx(8)
ans =
0.3536 0.3536 0.3536 0.3536 0.3536 0.3536 0.3536 0.3536
0.4904 0.4157 0.2778 0.0975 -0.0975 -0.2778 -0.4157 -0.4904
0.4619 0.1913 -0.1913 -0.4619 -0.4619 -0.1913 0.1913 0.4619
0.4157 -0.0975 -0.4904 -0.2778 0.2778 0.4904 0.0975 -0.4157
0.3536 -0.3536 -0.3536 0.3536 0.3536 -0.3536 -0.3536 0.3536
0.2778 -0.4904 0.0975 0.4157 -0.4157 -0.0975 0.4904 -0.2778
0.1913 -0.4619 0.4619 -0.1913 -0.1913 0.4619 -0.4619 0.1913
0.0975 -0.2778 0.4157 -0.4904 0.4904 -0.4157 0.2778 -0.0975
一旦我们得到校正后的 DCT 矩阵,我们可以通过与您使用的 1:8
的测试向量进行矩阵乘法来检查实际的 DCT 计算:
>> dct1d(8)*((1:8).')
ans =
12.7279
-6.4423
-0.0000
-0.6735
0
-0.2009
-0.0000
-0.0507
1.这段代码实际上是在 dctmtx
中所做的,一旦你去掉所有的错误检查和输入一致性检查。
我正在尝试使用矩阵向量乘法在 MATLAB 中实现 DCT。具体来说,我想创建系数的 DCT 矩阵,然后用它与一维信号相乘来计算一维 DCT。
这是我生成 DCT 矩阵的代码:
function D=dct1d(n)
for u=0:n-1
if u==0
au=sqrt(1/n);
else
au=sqrt(2/n);
end
for x=0:n-1
D(u+1,x+1)=au*cos(((2*x+1)*u*pi)/2*n);
end
end
在此之后,我尝试使用 x = [1 2 3 4 5 6 7 8]
:
x=[1 2 3 4 5 6 7 8];
y=dct1(8)*x';
它给出的答案是:
12.7279
18.0000
18.0000
18.0000
18.0000
18.0000
18.0000
18.0000
然而,正确答案是:
12.7279
-6.4423
-0.0000
-0.6735
0
-0.2009
-0.0000
-0.0507
您的代码中的错误非常轻微但很重要。您计算系数的行:
D(u+1,x+1)=au*cos(((2*x+1)*u*pi)/2*n);
看看这行的最后一部分:
D(u+1,x+1)=au*cos(((2*x+1)*u*pi)/2*n);
%// ^^^
Because multiplication and division are equal in precedence,这和做完全一样:
D(u+1,x+1)=au*cos((((2*x+1)*u*pi)/2)*n);
因此,您没有除以 2n
。你除以 2 然后 乘以 n
这是不正确的。您只需要用方括号将 2*n
操作括起来:
D(u+1,x+1)=au*cos(((2*x+1)*u*pi)/(2*n));
一旦你这样做了,我们就得到了正确的 DCT 矩阵。顺便说一句,检查您是否有正确答案的一种方法是使用 dctmtx
function, which calculates the N x N
DCT coefficient matrix that you're seeking. However, this function is part of the Signal Processing Toolbox, so if you don't have that then you unfortunately can't use this function, but if I can suggest an alternative answer instead of using for
loops, I would build a 2D grid of coordinates with meshgrid
,然后计算元素乘积。
这样的东西会起作用:
function D = dct1d(n)
[x,u] = meshgrid(0:n-1);
D = sqrt(2/n)*cos(((2*x+1).*u*pi)/(2*n));
D(1,:) = D(1,:) / sqrt(2);
end
不是使用 if
语句来确定每行需要应用的权重,我们可以只使用 sqrt(2/n)
然后除以第一行的 sqrt(2)
这样你而是除以 sqrt(1/n)
。此代码应产生与更正代码相同的结果1.
无论如何,一旦我做了这些更正,我就比较了你的代码给出的答案和 dctmtx
给出的答案,它是正确的:
>> dct1d(8)
ans =
0.3536 0.3536 0.3536 0.3536 0.3536 0.3536 0.3536 0.3536
0.4904 0.4157 0.2778 0.0975 -0.0975 -0.2778 -0.4157 -0.4904
0.4619 0.1913 -0.1913 -0.4619 -0.4619 -0.1913 0.1913 0.4619
0.4157 -0.0975 -0.4904 -0.2778 0.2778 0.4904 0.0975 -0.4157
0.3536 -0.3536 -0.3536 0.3536 0.3536 -0.3536 -0.3536 0.3536
0.2778 -0.4904 0.0975 0.4157 -0.4157 -0.0975 0.4904 -0.2778
0.1913 -0.4619 0.4619 -0.1913 -0.1913 0.4619 -0.4619 0.1913
0.0975 -0.2778 0.4157 -0.4904 0.4904 -0.4157 0.2778 -0.0975
>> dctmtx(8)
ans =
0.3536 0.3536 0.3536 0.3536 0.3536 0.3536 0.3536 0.3536
0.4904 0.4157 0.2778 0.0975 -0.0975 -0.2778 -0.4157 -0.4904
0.4619 0.1913 -0.1913 -0.4619 -0.4619 -0.1913 0.1913 0.4619
0.4157 -0.0975 -0.4904 -0.2778 0.2778 0.4904 0.0975 -0.4157
0.3536 -0.3536 -0.3536 0.3536 0.3536 -0.3536 -0.3536 0.3536
0.2778 -0.4904 0.0975 0.4157 -0.4157 -0.0975 0.4904 -0.2778
0.1913 -0.4619 0.4619 -0.1913 -0.1913 0.4619 -0.4619 0.1913
0.0975 -0.2778 0.4157 -0.4904 0.4904 -0.4157 0.2778 -0.0975
一旦我们得到校正后的 DCT 矩阵,我们可以通过与您使用的 1:8
的测试向量进行矩阵乘法来检查实际的 DCT 计算:
>> dct1d(8)*((1:8).')
ans =
12.7279
-6.4423
-0.0000
-0.6735
0
-0.2009
-0.0000
-0.0507
1.这段代码实际上是在 dctmtx
中所做的,一旦你去掉所有的错误检查和输入一致性检查。