从功率谱生成热图,Matlab
Generating heatmap from power spectrum, Matlab
我有一个非常大的局部场电位(原始电压)数据集,我已经对其进行了预处理以去除噪声和异常值。我对数据进行了排列,使每一行代表 30 秒的样本。我生成的功率谱如下:
Fs = 1024
LFP = 1075x30720 double
pxx = 1075x4097 double
for k = 1:1075;
pxx(kk,:) = pwelch(LFP(k,:));
end
目标:生成热图,使 pxx 的每一行都是生成的热图上的一列,所以我应该在 x 轴上有 1075 个 bin,我希望将 Y 轴的频率限制在 0 - 120 之间赫兹。我试过使用 imagesc 但遇到困难,谢谢。
要绘制结果,您需要做几件事:
- select 列对应最高 120Hz 的频率;
- 转置
pxx
使 pxx
的行显示为生成图像的列;
- 如果你想让最高频率出现在顶部,用
imagesc
; 用flipud
把数据颠倒过来
- 可选择转换为对数分贝标度;
- 绘图使用
imagesc
或 pcolor
;
- 使用
caxis
选择要显示的值的动态范围,以便您在颜色图中获得适当的值分布;
- select 热图样式的颜色图,例如
colormap(hot)
。
这可以通过以下方式完成:
% 1) Compute maximum frequency index
Fmax = 120; % Hz
M = 1 + round(Fmax/(0.5*Fs/(size(pxx,2)-1)));
% select displayed section
pxx_select = pxx(:,1:M);
% 2) transpose matrix
pxx_reshape = transpose(pxx_select);
% 3) flip data upside down for imagesc
pxx_reshape = flipud(pxx_reshape);
% 4) convert to decibel scale
pxx_dB = 10*log10(pxx_reshape);
% 5) generate plot
figure(1);
imagesc(pxx_dB);
% 6) choose dynamic range
% assign e.g. 80dB below peak power to the first value in the colormap
% and the peak power to the last value in the colormap
caxis([max(max(pxx_dB))-80 max(max(pxx_dB))]);
% 7) select colormap
colormap(hot);
或者,如果您想控制显示的轴:
% 1) Compute maximum frequency index
Fmax = 120; % Hz
M = 1 + round(Fmax/(0.5*Fs/(size(pxx,2)-1)));
% select displayed section
pxx_select = pxx(:,1:M);
% 2) transpose matrix
pxx_reshape = transpose(pxx_select);
% 3) flipud not needed with pcolor, instead set t & f axis:
t = (size(LPF,2)/Fs)*[0:size(LPF,1)];
f = [0:M]*Fmax/(M-1);
% 4) convert to decibel scale
pxx_dB = 10*log10(pxx_reshape);
% 5) generate plot
figure(2);
% Note: extend by one row & column since pcolor does not show the last row/col
P2 = [pxx_dB pxx_dB(:,end); pxx_dB(end,:) pxx_dB(end,end)];
pcolor(t,f,P2); shading flat;
% 6) choose dynamic range
% assign e.g. 80dB below peak power to the first value in the colormap
% and the peak power to the last value in the colormap
caxis([max(max(pxx_dB))-80 max(max(pxx_dB))]);
% 7) select colormap
colormap(hot);
xlabel("time (s)");
ylabel("frequency (Hz)");
作为示例,您会得到类似于
的图表
对于通过以下方式生成的简单的缓慢频率变化的音调:
T = size(LPF,1)-1;
phi = 0;
n = [0:size(LPF,2)-1];
for k=1:size(LPF,1)
f = 0.5*(fmin+fmax) + 0.5*(fmax-fmin)*sin(2*pi*k/T);
LPF(k,:) = sin(2*pi*f*n/Fs + phi);
phi = mod(phi + 2*pi*f*size(LPF,2)/Fs, 2*pi);
end
我有一个非常大的局部场电位(原始电压)数据集,我已经对其进行了预处理以去除噪声和异常值。我对数据进行了排列,使每一行代表 30 秒的样本。我生成的功率谱如下:
Fs = 1024
LFP = 1075x30720 double
pxx = 1075x4097 double
for k = 1:1075;
pxx(kk,:) = pwelch(LFP(k,:));
end
目标:生成热图,使 pxx 的每一行都是生成的热图上的一列,所以我应该在 x 轴上有 1075 个 bin,我希望将 Y 轴的频率限制在 0 - 120 之间赫兹。我试过使用 imagesc 但遇到困难,谢谢。
要绘制结果,您需要做几件事:
- select 列对应最高 120Hz 的频率;
- 转置
pxx
使pxx
的行显示为生成图像的列; - 如果你想让最高频率出现在顶部,用
imagesc
; 用 - 可选择转换为对数分贝标度;
- 绘图使用
imagesc
或pcolor
; - 使用
caxis
选择要显示的值的动态范围,以便您在颜色图中获得适当的值分布; - select 热图样式的颜色图,例如
colormap(hot)
。
flipud
把数据颠倒过来
这可以通过以下方式完成:
% 1) Compute maximum frequency index
Fmax = 120; % Hz
M = 1 + round(Fmax/(0.5*Fs/(size(pxx,2)-1)));
% select displayed section
pxx_select = pxx(:,1:M);
% 2) transpose matrix
pxx_reshape = transpose(pxx_select);
% 3) flip data upside down for imagesc
pxx_reshape = flipud(pxx_reshape);
% 4) convert to decibel scale
pxx_dB = 10*log10(pxx_reshape);
% 5) generate plot
figure(1);
imagesc(pxx_dB);
% 6) choose dynamic range
% assign e.g. 80dB below peak power to the first value in the colormap
% and the peak power to the last value in the colormap
caxis([max(max(pxx_dB))-80 max(max(pxx_dB))]);
% 7) select colormap
colormap(hot);
或者,如果您想控制显示的轴:
% 1) Compute maximum frequency index
Fmax = 120; % Hz
M = 1 + round(Fmax/(0.5*Fs/(size(pxx,2)-1)));
% select displayed section
pxx_select = pxx(:,1:M);
% 2) transpose matrix
pxx_reshape = transpose(pxx_select);
% 3) flipud not needed with pcolor, instead set t & f axis:
t = (size(LPF,2)/Fs)*[0:size(LPF,1)];
f = [0:M]*Fmax/(M-1);
% 4) convert to decibel scale
pxx_dB = 10*log10(pxx_reshape);
% 5) generate plot
figure(2);
% Note: extend by one row & column since pcolor does not show the last row/col
P2 = [pxx_dB pxx_dB(:,end); pxx_dB(end,:) pxx_dB(end,end)];
pcolor(t,f,P2); shading flat;
% 6) choose dynamic range
% assign e.g. 80dB below peak power to the first value in the colormap
% and the peak power to the last value in the colormap
caxis([max(max(pxx_dB))-80 max(max(pxx_dB))]);
% 7) select colormap
colormap(hot);
xlabel("time (s)");
ylabel("frequency (Hz)");
作为示例,您会得到类似于
的图表对于通过以下方式生成的简单的缓慢频率变化的音调:
T = size(LPF,1)-1;
phi = 0;
n = [0:size(LPF,2)-1];
for k=1:size(LPF,1)
f = 0.5*(fmin+fmax) + 0.5*(fmax-fmin)*sin(2*pi*k/T);
LPF(k,:) = sin(2*pi*f*n/Fs + phi);
phi = mod(phi + 2*pi*f*size(LPF,2)/Fs, 2*pi);
end