Matlab谱图函数的手动实现
Manual implementation of Matlab spectogram function
我正在尝试实现我自己的函数,它提供与 Matlab 频谱图函数相同的结果。
到目前为止,我已经完成了这样的功能:
function out = manulaSpectogram(x, win, noverlap, nfft)
x = x(:);
n = length(x);
wlen = length(win);
nUnique = ceil((1+nfft)/2); % number of uniqure points
L = fix((n-noverlap)/(wlen-noverlap)); % number of signal frames
out = zeros(L, nUnique);
index = 1:wlen;
for i = 0:L-1
xw = win.*x(index);
X = fft(xw, nfft);
out(i+1, :) = X(1:nUnique);
index = index + (wlen - noverlap);
end
end
在我的测试中,当参数 nfft 大于或等于 window.
的长度时,它可以完美运行并给出与频谱图函数相同的结果
% first test (nnft = window length):
A = [1,2,3,4,5,6];
window = 6;
overlap = 2;
nfft = 6;
s = spectrogram(A, hamming(window), overlap, nfft)'
s2 = manulaSpectogram(A, hamming(window), overlap, nfft)
% results:
s =
9.7300 + 0.0000i -5.2936 + 0.9205i 0.7279 - 0.3737i -0.1186 + 0.0000i
s2 =
9.7300 + 0.0000i -5.2936 - 0.9205i 0.7279 + 0.3737i -0.1186 + 0.0000i
% second test (nfft > window length):
A = [1,2,3,4,5,6];
window = 3;
overlap = 2;
nfft = 6;
s = spectrogram(A, hamming(window), overlap, nfft)'
s2 = manulaSpectogram(A, hamming(window), overlap, nfft)
% results:
s =
2.3200 + 0.0000i 0.9600 + 1.9399i -1.0400 + 1.5242i -1.6800 + 0.0000i
3.4800 + 0.0000i 1.5000 + 2.8752i -1.5000 + 2.3209i -2.5200 + 0.0000i
4.6400 + 0.0000i 2.0400 + 3.8105i -1.9600 + 3.1177i -3.3600 + 0.0000i
5.8000 + 0.0000i 2.5800 + 4.7458i -2.4200 + 3.9144i -4.2000 + 0.0000i
s2 =
2.3200 + 0.0000i 0.9600 - 1.9399i -1.0400 - 1.5242i -1.6800 + 0.0000i
3.4800 + 0.0000i 1.5000 - 2.8752i -1.5000 - 2.3209i -2.5200 + 0.0000i
4.6400 + 0.0000i 2.0400 - 3.8105i -1.9600 - 3.1177i -3.3600 + 0.0000i
5.8000 + 0.0000i 2.5800 - 4.7458i -2.4200 - 3.9144i -4.2000 + 0.0000i
如果window的长度小于nfft,则结果不同。
% third test (nfft < window length):
A = [1,2,3,4,5,6];
window = 6;
overlap = 2;
nfft = 3;
s = spectrogram(A, hamming(window), overlap, nfft)'
s2 = manulaSpectogram(A, hamming(window), overlap, nfft)
% results:
s =
9.7300 + 0.0000i 0.7279 - 0.3737i
s2 =
3.6121 + 0.0000i -1.6861 + 1.6807i
那么,即使在 nnft 小于 window 长度的情况下,我如何改进我的函数以接收相同的结果? Matlab的spectogram是如何计算这种情况的?
我正在尝试实现我自己的函数,因为频谱图函数是我需要从 Matlab 到 C# 语言实现的大型算法的一部分,所以我想知道频谱图“黑匣子”的作用..
我注意到当 window 大小大于 nfft 标量数时,必须以某种方式转换数据。最后我找到了一个可能在原始 spectogram
Matlab 函数中调用的内部 Matlab 函数。它被命名为 datawrap
并包装输入数据模 nfft.
所以在我的函数中,我必须在调用 fft 之前转换数据段(与 datawrap 函数的处理方式相同)。
改进功能:
function out = manulaSpectogram(x, win, noverlap, nfft)
x = x(:);
n = length(x);
wlen = length(win);
nUnique = ceil((1+nfft)/2); % number of uniqure points
L = fix((n-noverlap)/(wlen-noverlap)); % number of signal frames
out = zeros(L, nUnique);
index = 1:wlen;
for i = 0:L-1
xw = win.*x(index);
% added transformation
if length(xw) > nfft
xw = sum(buffer(xw, nfft), 2);
end
% end of added transformation
X = fft(xw, nfft);
out(i+1, :) = X(1:nUnique);
index = index + (wlen - noverlap);
end
end
我相信它工作正常,因为它给出了与 Matlab spectogram
函数相同的结果。
使用 window 计算块的末端为零。这需要作为强制性要求,window 和 FFT-Block 的矢量长度相同。 Matlab 允许使用默认参数不相等的值:nsc 不等于 nff(请参阅默认的 Matlab docu 频谱图)。如果 window 和 nff 的块长度不同,则需要一条错误消息。我认为默认参数的 Matlab 频谱图是错误的。与 LabView、Dasylab、Hewlett-Packard 和 http://www.schmid-werren.ch/hanspeter/publications/2012fftnoise.pdf
比较
我正在尝试实现我自己的函数,它提供与 Matlab 频谱图函数相同的结果。 到目前为止,我已经完成了这样的功能:
function out = manulaSpectogram(x, win, noverlap, nfft)
x = x(:);
n = length(x);
wlen = length(win);
nUnique = ceil((1+nfft)/2); % number of uniqure points
L = fix((n-noverlap)/(wlen-noverlap)); % number of signal frames
out = zeros(L, nUnique);
index = 1:wlen;
for i = 0:L-1
xw = win.*x(index);
X = fft(xw, nfft);
out(i+1, :) = X(1:nUnique);
index = index + (wlen - noverlap);
end
end
在我的测试中,当参数 nfft 大于或等于 window.
的长度时,它可以完美运行并给出与频谱图函数相同的结果% first test (nnft = window length):
A = [1,2,3,4,5,6];
window = 6;
overlap = 2;
nfft = 6;
s = spectrogram(A, hamming(window), overlap, nfft)'
s2 = manulaSpectogram(A, hamming(window), overlap, nfft)
% results:
s =
9.7300 + 0.0000i -5.2936 + 0.9205i 0.7279 - 0.3737i -0.1186 + 0.0000i
s2 =
9.7300 + 0.0000i -5.2936 - 0.9205i 0.7279 + 0.3737i -0.1186 + 0.0000i
% second test (nfft > window length):
A = [1,2,3,4,5,6];
window = 3;
overlap = 2;
nfft = 6;
s = spectrogram(A, hamming(window), overlap, nfft)'
s2 = manulaSpectogram(A, hamming(window), overlap, nfft)
% results:
s =
2.3200 + 0.0000i 0.9600 + 1.9399i -1.0400 + 1.5242i -1.6800 + 0.0000i
3.4800 + 0.0000i 1.5000 + 2.8752i -1.5000 + 2.3209i -2.5200 + 0.0000i
4.6400 + 0.0000i 2.0400 + 3.8105i -1.9600 + 3.1177i -3.3600 + 0.0000i
5.8000 + 0.0000i 2.5800 + 4.7458i -2.4200 + 3.9144i -4.2000 + 0.0000i
s2 =
2.3200 + 0.0000i 0.9600 - 1.9399i -1.0400 - 1.5242i -1.6800 + 0.0000i
3.4800 + 0.0000i 1.5000 - 2.8752i -1.5000 - 2.3209i -2.5200 + 0.0000i
4.6400 + 0.0000i 2.0400 - 3.8105i -1.9600 - 3.1177i -3.3600 + 0.0000i
5.8000 + 0.0000i 2.5800 - 4.7458i -2.4200 - 3.9144i -4.2000 + 0.0000i
如果window的长度小于nfft,则结果不同。
% third test (nfft < window length):
A = [1,2,3,4,5,6];
window = 6;
overlap = 2;
nfft = 3;
s = spectrogram(A, hamming(window), overlap, nfft)'
s2 = manulaSpectogram(A, hamming(window), overlap, nfft)
% results:
s =
9.7300 + 0.0000i 0.7279 - 0.3737i
s2 =
3.6121 + 0.0000i -1.6861 + 1.6807i
那么,即使在 nnft 小于 window 长度的情况下,我如何改进我的函数以接收相同的结果? Matlab的spectogram是如何计算这种情况的?
我正在尝试实现我自己的函数,因为频谱图函数是我需要从 Matlab 到 C# 语言实现的大型算法的一部分,所以我想知道频谱图“黑匣子”的作用..
我注意到当 window 大小大于 nfft 标量数时,必须以某种方式转换数据。最后我找到了一个可能在原始 spectogram
Matlab 函数中调用的内部 Matlab 函数。它被命名为 datawrap
并包装输入数据模 nfft.
所以在我的函数中,我必须在调用 fft 之前转换数据段(与 datawrap 函数的处理方式相同)。 改进功能:
function out = manulaSpectogram(x, win, noverlap, nfft)
x = x(:);
n = length(x);
wlen = length(win);
nUnique = ceil((1+nfft)/2); % number of uniqure points
L = fix((n-noverlap)/(wlen-noverlap)); % number of signal frames
out = zeros(L, nUnique);
index = 1:wlen;
for i = 0:L-1
xw = win.*x(index);
% added transformation
if length(xw) > nfft
xw = sum(buffer(xw, nfft), 2);
end
% end of added transformation
X = fft(xw, nfft);
out(i+1, :) = X(1:nUnique);
index = index + (wlen - noverlap);
end
end
我相信它工作正常,因为它给出了与 Matlab spectogram
函数相同的结果。
使用 window 计算块的末端为零。这需要作为强制性要求,window 和 FFT-Block 的矢量长度相同。 Matlab 允许使用默认参数不相等的值:nsc 不等于 nff(请参阅默认的 Matlab docu 频谱图)。如果 window 和 nff 的块长度不同,则需要一条错误消息。我认为默认参数的 Matlab 频谱图是错误的。与 LabView、Dasylab、Hewlett-Packard 和 http://www.schmid-werren.ch/hanspeter/publications/2012fftnoise.pdf
比较