Permalink
Switch branches/tags
Nothing to show
Find file Copy path
Fetching contributors…
Cannot retrieve contributors at this time
102 lines (91 sloc) 3.37 KB
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
%% Filename: showspectrogram.m
%%
%% Project: A Wishbone Controlled PWM (audio) controller
%%
%% Purpose: To generate a spectrogram image which can then be used to
%% evaluate the wavfp.dbl file produced by the pdmdemo executable
%% in this directory.
%%
%% This script has been tested with Octave, and it is known to work with
%% Octave. While it may work within Matlab as well, no representation is
%% being made to that effect.
%%
%% Creator: Dan Gisselquist, Ph.D.
%% Gisselquist Technology, LLC
%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
%% Copyright (C) 2017, Gisselquist Technology, LLC
%%
%% This program is free software (firmware): you can redistribute it and/or
%% modify it under the terms of the GNU General Public License as published
%% by the Free Software Foundation, either version 3 of the License, or (at
%% your option) any later version.
%%
%% This program is distributed in the hope that it will be useful, but WITHOUT
%% ANY WARRANTY; without even the implied warranty of MERCHANTIBILITY or
%% FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
%% for more details.
%%
%% You should have received a copy of the GNU General Public License along
%% with this program. (It's in the $(ROOT)/doc directory. Run make with no
%% target there if the PDF file isn't present.) If not, see
%% <http://www.gnu.org/licenses/> for a copy.
%%
%% License: GPL, v3, as defined and found on www.gnu.org,
%% http://www.gnu.org/licenses/gpl.html
%%
%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%
%%
function [img] = showspectrogram(fname)
% Read the waveform into memory
fid = fopen(fname,'r');
dat = fread(fid, inf, 'double');
fclose(fid);
% You can plot it at this point if you would like
% plot(dat);
sample_rate = 44.1e3;
window_length = 512; % Size of the FFT
Nc = 256; % Number of colors
h = hanning(window_length); % About 1 ms
rows = (length(dat)-length(h)) / length(h) * 2
img = zeros(length(h)/2, rows);
size(img)
for n=1:rows
cut = dat((1:length(h)) + (n-1)*length(h)/2);
cutf = abs(fft(cut.*h)).^2;
img(:,n) = cutf(1:(length(h)/2));
end
mintime = 0;
maxtime = rows / sample_rate * window_length / 2;
minfreq = 0;
maxfreq = 0.5 * sample_rate;
% Normalize the image so that the maximum value is one.
img = img ./ max(max(img));
%
% Turn out spectrogram image into dB
%
% Adding 1e-8 is useful for artificially forcing the floor of the
% image to be at a particular value, as well as for keeping the log
% from reporting values too negative to plot successfully.
img = 10 * log(img + 1e-8)/log(10);
img = (img + 80)/80;
colormap(mymap(256));
image([mintime, maxtime], [minfreq, maxfreq/1e3], img * Nc);
ylabel('Frequency (kHz)');
xlabel('Time (s)');
% Trim this image output to the relevant portion of the output
axis([0, 6.2, 0, maxfreq/1e3]);
% For the plot, comment out the img ./ max(max(img)) line, as well as the
% img = (img+80)/80 line. Then you can call this as:
% img = showspectrogram('pdm8k-weak.dbl');
% imgpwm = showspectrogram('pwm8k-weak.dbl');
%
% freq = (0:(window_length-1))/(window_length)*(sample_rate/2/1e3)
% plot(freq,img(:,1400),'b;PDM;',freq,imgpwm(:,1400),'r;PWM;');
% grid on; xlabel('Frequency (kHz)'); ylabel('Volume (dB)');
% axis([0,22,-70,50]);