Skip to content
Snippets Groups Projects
Commit cc2a8409 authored by Eric Kooistra's avatar Eric Kooistra
Browse files

Added run_pfb_complex.m to provide complex input stimuli for HDL PFB.

parent 9bcec3fc
Branches
No related tags found
No related merge requests found
Show whitespace changes
Inline Side-by-side
applications/apertif/matlab/apertif_matlab_readme.txt
+ 9
− 7
View file @ cc2a8409
......@@ -36,15 +36,13 @@ This directory contains Matlab code for modelling DSP aspects of Apertif:
pfir_coeff.m Function Compute polyphase filterbank coefficients (same purpose as FIR_LPF_ApertifCF.m but with more options)
run_pfir_coeff.m Script Run pfir_coeff.m to show filter response in time and frequency domain
run_pfir.m Script Run the PFIR with real input to create HDL reference data (based on one_pfb.m)
4) FFT
sinusoids_and_fft_frequency_bins.m Script Code on DFT and DTFT from Matlab blog on internet
try_fft.m Script Try DFT and DTFT
5) Data path HDL functions
5) Data path functions
wg.m Function Waveform generator per block of data
dt.m Function Delay tracking per block of data
bsn_source.m Function Block Sequence Number source per block of data
......@@ -54,10 +52,6 @@ This directory contains Matlab code for modelling DSP aspects of Apertif:
reorder_serial.m Function Select data from block of data
corner_turn.m Function Collect M blocks of in_data(1:K) and then transpose [M] and [K] for corner turn
run_pfft.m Script Run the PFFT with real input to create HDL reference data (based on one_pfb.m)
run_pfft_complex.m Script Run the PFFT with complex input to create HDL reference data (based on one_pfb.m)
6) Polyphase filterbank
one_pfb.m Script Run data path model with subband polyphase filterbank (real input)
......@@ -66,3 +60,11 @@ This directory contains Matlab code for modelling DSP aspects of Apertif:
two_pfb.m Script Run data path model with subband polyphase filterbank followed by a channel polyphase filterbank
7) HDL reference data
run_pfir.m Script Run the PFIR with real input to create HDL reference data (based on one_pfb.m)
run_pfft.m Script Run the PFFT with real input to create HDL reference data (based on one_pfb.m)
run_pfft_complex.m Script Run the PFFT with complex input to create HDL reference data (based on one_pfb.m)
run_pfb_complex.m Script Run the PFB with complex input to create HDL reference data (based on one_pfb.m)
applications/apertif/matlab/run_pfb_complex.m 0 → 100644
+ 488
− 0
View file @ cc2a8409
%-----------------------------------------------------------------------------
%
% Copyright (C) 2016
% ASTRON (Netherlands Institute for Radio Astronomy) <http://www.astron.nl/>
% P.O.Box 2, 7990 AA Dwingeloo, The Netherlands
%
% This program is free software: 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
% MERCHANTABILITY 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. If not, see <http://www.gnu.org/licenses/>.
%
%-----------------------------------------------------------------------------
% Author: Eric Kooistra, 2016 (for Apertif)
%
% Purpose : Run the PFB = PFIR + PFFT with complex input:
% To create input and output reference data for HDL implementation in
% data/, e.g.:
%
% . PFIR coeffcients: run_pfb_complex_m_pfir_coeff_fircls1_16taps_128points_16b.dat
% . PFIR coeffcients, WG, PFIR, PFFT data: run_pfb_complex_m_chirp_8b_16taps_128points_16b_16b.dat
%
% Description :
%
% * See one_pfb.m for general description.
% This run_pfb_complex.m is almost the same as run_pfft_complex.m,
% except that it includes the PFIR.
clear all;
close all;
rng('default');
fig=0;
%% Selectable test bench settings
% Input signal
tb.model_signal = 'sinusoid'; % Use sinusoid to check the frequency response of the FFT with frequency set via tb.channel_wg
tb.model_signal = 'phasor'; % Use phasor to check the frequency response of the FFT with frequency set via tb.channel_wg
%tb.model_signal = 'noise'; % Use noise to check all aspects in one signal
%tb.model_signal = 'square'; % Use square to create a square wave with period set via tb.channel_wg
%tb.model_signal = 'impulse'; % Use impulse to check the impulse response for different ctrl_wg.sop
tb.model_quantization = 'floating point';
tb.model_quantization = 'fixed point';
tb.nof_channels = 32;
% Carrier frequency
tb.channel_wg = 1; % channel range -tb.nof_channels/2:tb.nof_channels/2-1, can be fraction to have any sinusoid frequency
%tb.channel_wg = 1.55;
%tb.channel_wg = 12;
tb.phase = 0; % initial 'sinusoid' phase offset normalized to 2pi
tb.sop = 1; % initial 'impulse' start index in range ctrl_wg.block_size = tb.channel_fft_size
% Model a frequency sweep of the 'sinusoid'
tb.chirp = 0; % 0 = use fixed tb.channel_wg frequency or pulse period equal to block_size
tb.chirp = 1; % else increment WG frequency every block to have chirp frequency sweep or slide the pulse
if strcmp(tb.model_signal, 'noise')
tb.nof_tchan = 10 ; % number of channel periods to simulate
elseif tb.chirp
tb.nof_tchan = 200; % number of channel periods to simulate
else
tb.nof_tchan = 5; % number of channel periods to simulate
end
tb.plot_per_block = 0; % 1 = plot spectrum for each block in time, else skip this plot to save time
%tb.plot_per_block = 1;
%% Derived test bench settings
tb.nof_complex = 2;
tb.channel_i = floor(tb.channel_wg + tb.nof_channels/2); % natural channel index in range 0:tb.nof_channels-1
tb.channel_I = tb.channel_i+1; % equivalent Matlab indexing for range 1:tb.nof_channels
if tb.channel_i<tb.nof_channels-1 % determine index of next neighbour channel
tb.channel_iplus1 = tb.channel_i+1;
else
tb.channel_iplus1 = 0; % wrap channel index tb.nof_channels to channel index 0
end
tb.channel_Iplus1 = tb.channel_iplus1+1; % equivalent Matlab index for next neighbour channel in range 1:tb.nof_channels
tb.channel_fft_size = tb.nof_channels; % channel complex FFT
fs = 1; % normalized sample frequency
fchan = fs/tb.channel_fft_size; % channel frequency relative to fs
disp(sprintf (['Test bench settings:\n', ...
'. Model coefficients and data = %s\n', ...
'. Number of channels = %d\n', ...
'. Number of channel periods to simulate = %d\n', ...
'. Channel filterbank complex FFT size = %d\n', ...
'. WG input signal = %s\n', ...
'. WG channel frequency = %f\n'], ...
tb.model_quantization, tb.nof_channels, tb.nof_tchan, tb.channel_fft_size, tb.model_signal, tb.channel_wg))
%% Waveform generator
% Common WG settings
ctrl_wg.complex = false;
if strcmp(tb.model_signal, 'phasor') || strcmp(tb.model_signal, 'noise')
ctrl_wg.complex = true;
end
ctrl_wg.block_nr = 0;
ctrl_wg.block_size = tb.channel_fft_size;
if strcmp(tb.model_quantization, 'floating point')
ctrl_wg.data_w = 0;
lsb = 1/2^8; % Assume 8 bit ADC to have reference for AGWN level in floating point
else
ctrl_wg.data_w = 8;
lsb = 1/2^ctrl_wg.data_w;
end
ctrl_wg.q_full_scale = 2^(ctrl_wg.data_w-1);
ctrl_wg.agwn_sigma = 2.5*lsb;
ctrl_wg.agwn_sigma = 0; % AGWN sigma
ctrl_wg.ampl = 1; % full scale is 1, ampl > 1 causes analogue clipping
if ctrl_wg.agwn_sigma>0
ctrl_wg.ampl = 1-5*ctrl_wg.agwn_sigma;
ctrl_wg.ampl = 0.9;
%ctrl_wg.ampl = 1*lsb;
end
ctrl_wg.offset = 0; % DC offset
%ctrl_wg.offset = 0.1;
% Signal specific WG settings
if strcmp(tb.model_signal, 'sinusoid') || strcmp(tb.model_signal, 'phasor')
ctrl_wg.signal = 'sinusoid';
ctrl_wg.psk = 0; % no PSK modulation
if tb.chirp
ctrl_wg.df = 0.49*fchan; % increment freq by df per block to create chirp
%ctrl_wg.df = 0.5*fchan;
else
ctrl_wg.df = 0;
end
if ctrl_wg.df == 0
ctrl_wg.freq = tb.channel_wg*fchan; % start WG at tb.channel_wg
else
ctrl_wg.freq = (tb.channel_i-1)*fchan; % start WG about one channel before tb.channel_wg to account for impulse response of several tchan
end
ctrl_wg.phase = tb.phase; % Phase offset normalized to 2pi
if ctrl_wg.freq == 0 && ctrl_wg.df == 0
ctrl_wg.offset = 1; % Force DC offset for zero Hz
end
elseif strcmp(tb.model_signal, 'noise')
ctrl_wg.signal = 'noise';
ctrl_wg.ampl = 0;
ctrl_wg.agwn_sigma = 1/3; % AGWN sigma is 1/3 full scale
else
ctrl_wg.signal = 'pulse';
ctrl_wg.sop = tb.sop; % Initial start of pulse index normalized to ctrl_wg.block_size = tb.channel_fft_size
ctrl_wg.period = ctrl_wg.block_size; % Pulse period
if tb.chirp
ctrl_wg.period = ctrl_wg.block_size+1; % Pulse period for pulse moving in block like a sop chirp
end
if strcmp(tb.model_signal, 'impulse')
ctrl_wg.width = 1; % Pulse width for FFT impulse response
else % default to 'square')
ctrl_wg.ampl = 2*ctrl_wg.ampl;
ctrl_wg.offset = -ctrl_wg.ampl/2;
ctrl_wg.period = ctrl_wg.block_size/tb.channel_wg; % Pulse period equal to channel_wg period
ctrl_wg.width = ctrl_wg.period/2; % Pulse width for square wave
end
ctrl_wg.width_remaining = 0;
end
%% Channel PFIR filter parameters
ctrl_pfir_channel.bypass = 0;
ctrl_pfir_channel.nof_polyphases = tb.channel_fft_size;
if strcmp(tb.model_quantization, 'floating point')
ctrl_pfir_channel.data_w = 0; % no quantization
ctrl_pfir_channel.coeff_w = 0; % no quantization
else
ctrl_pfir_channel.data_w = 16;
ctrl_pfir_channel.backoff_w = 1; % attenuate PFIR data output by 2^backoff_w before quantization to account for PFIR overshoot
ctrl_pfir_channel.coeff_w = 16; % nof bits per coefficient including 1 sign bit
ctrl_pfir_channel.coeff_q_max = 2^(ctrl_pfir_channel.coeff_w-1)-1;
ctrl_pfir_channel.coeff_q_full_scale = 2^(ctrl_pfir_channel.coeff_w-1);
ctrl_pfir_channel.data_q_full_scale = 2^(ctrl_pfir_channel.data_w-1-ctrl_pfir_channel.backoff_w);
end
% Calculate free choice PFIR coefficients (can be floating point or fixed point)
ctrl_pfir_channel.nof_taps = 16; % Number of taps
ctrl_pfir_channel.nof_coefficients = ctrl_pfir_channel.nof_polyphases*ctrl_pfir_channel.nof_taps; % Number of filter coefficients (taps)
ctrl_pfir_channel.config.design = 'fir1';
ctrl_pfir_channel.config.design = 'fircls1'; % 'fir1', 'fircls1'
ctrl_pfir_channel.config.design_flag = ''; % only for fircls1: 'trace', ''
ctrl_pfir_channel.config.interpolate = 'interpft'; % only for fircls1: 'resample', 'fourier', 'interpft'
ctrl_pfir_channel.r_pass = 1e-3; % only for fircls1
ctrl_pfir_channel.r_stop = 1e-4; % only for fircls1
ctrl_pfir_channel.hp_factor = 1.050; % Adjust channel half power bandwidth
ctrl_pfir_channel.hp_factor = 1;
ctrl_pfir_channel.BWchan = ctrl_pfir_channel.hp_factor / tb.channel_fft_size; % Channel bandwidth
hfir_channel_coeff = pfir_coeff(ctrl_pfir_channel.nof_polyphases, ...
ctrl_pfir_channel.nof_taps, ...
ctrl_pfir_channel.BWchan, ...
ctrl_pfir_channel.r_pass, ...
ctrl_pfir_channel.r_stop, ...
ctrl_pfir_channel.coeff_w, ...
ctrl_pfir_channel.config);
ctrl_pfir_channel.coeff = reshape(hfir_channel_coeff, ctrl_pfir_channel.nof_polyphases, ctrl_pfir_channel.nof_taps);
ctrl_pfir_channel.coeff = flipud(ctrl_pfir_channel.coeff);
ctrl_pfir_channel.Zdelays = zeros(ctrl_pfir_channel.nof_polyphases, ctrl_pfir_channel.nof_taps-1);
ctrl_pfir_channel.gain = 1; % no gain adjustment in PFIR, just apply the coefficients to the input data
disp(sprintf(['Channel PFIR filter settings:\n', ...
'. Bypass = %d\n', ...
'. Number of polyphases = %d\n', ...
'. Number of taps = %d\n'], ...
ctrl_pfir_channel.bypass, ...
ctrl_pfir_channel.nof_polyphases, ...
ctrl_pfir_channel.nof_taps));
if strcmp(tb.model_quantization, 'fixed point')
disp(sprintf(['. Coefficient width = %d\n', ...
'. Coefficient q_max = %d\n', ...
'. Data width = %d\n', ...
'. Data backoff width = %d\n', ...
', Data q_full_scale = %d\n'], ...
ctrl_pfir_channel.coeff_w, ...
ctrl_pfir_channel.coeff_q_max, ...
ctrl_pfir_channel.data_w, ...
ctrl_pfir_channel.backoff_w, ...
ctrl_pfir_channel.data_q_full_scale));
end
%% Channel FFT parameters
ctrl_pfft_channel.fft_size = tb.channel_fft_size;
ctrl_pfft_channel.complex = true; % complex input FFT
ctrl_pfft_channel.gain = ctrl_pfft_channel.fft_size;
if strcmp(tb.model_quantization, 'floating point')
ctrl_pfft_channel.data_w = 0;
else
ctrl_pfft_channel.data_w = 16;
end
ctrl_pfft_channel.data_q_full_scale = 2^(ctrl_pfft_channel.data_w-1);
%% Run the data path processing
t_start = cputime;
data_wg = zeros(tb.nof_tchan, tb.channel_fft_size);
data_pfir_channel = zeros(tb.nof_tchan, tb.channel_fft_size);
data_pfft_channel = zeros(tb.nof_tchan, tb.nof_channels);
for bi = 0:tb.nof_tchan-1
bI = bi+1;
% Data path (DP)
[ctrl_wg, block_wg] = wg( ctrl_wg);
[ctrl_pfir_channel, block_pfir_channel] = pfir( ctrl_pfir_channel, block_wg);
[ block_pfft_channel] = pfft( ctrl_pfft_channel, block_pfir_channel);
% Capture data at each DP interface
data_wg(bI, :) = block_wg;
data_pfir_channel(bI, :) = block_pfir_channel;
data_pfft_channel(bI, :) = block_pfft_channel;
end
t_stop = cputime;
disp(sprintf('Total processing time: %f seconds', t_stop-t_start));
%% Save coeffcients and data reference files for stimuli and verification of HDL implementation
if strcmp(tb.model_quantization, 'fixed point')
% Quantize channel PFIR filter coefficients
int_hfir = round(ctrl_pfir_channel.coeff_q_full_scale * hfir_channel_coeff); % map largest coefficient to q_max
% Quantize WG output
wg_q_data = round(ctrl_wg.q_full_scale * data_wg);
% Quantize PFIR channel output
pfir_channel_q_data = round(ctrl_pfir_channel.data_q_full_scale * data_pfir_channel);
% Quantize PFFT channel output real and imaginary
pfft_channel_q_data = round(ctrl_pfft_channel.data_q_full_scale * data_pfft_channel);
% Save the quantized integer channel PFIR filter coefficients in a separate file
L = ctrl_pfir_channel.nof_taps;
N = ctrl_pfir_channel.nof_polyphases; % = tb.channel_fft_size;
file_name_prefix = 'data/run_pfb_complex_m_';
file_name = sprintf([file_name_prefix, 'pfir_coeff_', ctrl_pfir_channel.config.design, '_%dtaps_%dpoints_%db.dat'], L, N, ctrl_pfir_channel.coeff_w);
fid = fopen(file_name, 'w');
fprintf(fid,'%d\n', int_hfir);
fclose(fid);
% Save PFIR filter coefficients again together with the data path signals for WG, PFIR and PFFT
tbegin = 1;
tend = tb.nof_tchan; % <= tb.nof_tchan
file_name = sprintf([file_name_prefix, '%s_%db_%dtaps_%dpoints_%db_%db.dat'], ctrl_wg.name, ctrl_wg.data_w, L, N, ctrl_pfir_channel.data_w, ctrl_pfft_channel.data_w);
fid = fopen(file_name, 'w');
fprintf(fid,'Nof lines PFIR coefficients: %d\n', length(int_hfir)); % save PFIR filter coefficients again also in the data file to have a complete set
fprintf(fid,'Nof lines WG output: %d\n', length(wg_q_data(:)));
fprintf(fid,'Nof lines PFIR output: %d\n', length(pfir_channel_q_data(:)));
fprintf(fid,'Nof lines PFFT output: %d\n', length(pfft_channel_q_data(:)));
% Channel PFIR coefficients
fprintf(fid,'%d\n', int_hfir);
% WG output
for bI = tbegin:tend
for bJ = 1:N
fprintf(fid,'%8d%8d\n', real(wg_q_data(bI, bJ)), imag(wg_q_data(bI, bJ)));
end
end;
% PFIR channel output
for bI = tbegin:tend
for bJ = 1:N
fprintf(fid,'%12d%12d\n', real(pfir_channel_q_data(bI, bJ)), imag(pfir_channel_q_data(bI, bJ)));
end
end
% PFFT channel output real and imaginary
for bI = tbegin:tend
for bJ = 1:tb.nof_channels
fprintf(fid,'%10d%10d\n', real(pfft_channel_q_data(bI, bJ)), imag(pfft_channel_q_data(bI, bJ)));
end
end
fclose(fid);
end
%% Plot data path results
xfig = 300;
yfig = 200;
xfigw = 1000;
yfigw = 800;
dfig = 20;
ts = (0:tb.channel_fft_size*tb.nof_tchan-1)/tb.channel_fft_size; % time of ADC / WG samples in channel periods
tchan_all = (0: tb.nof_channels*tb.nof_tchan-1)/tb.nof_channels; % time in channel periods for block of channels
tchan_one = (0: tb.nof_tchan-1); % time in channel periods for one channel
chan_I = tb.channel_I + [0: tb.nof_channels: tb.nof_channels*tb.nof_tchan-1]; % get Matlab indices of all channel_I
chan_Iplus1 = tb.channel_Iplus1 + [0: tb.nof_channels: tb.nof_channels*tb.nof_tchan-1]; % get Matlab indices of all next neighbour channel_I+1
%% Plot WG output
fig=fig+1;
figure('position', [xfig+fig*dfig yfig-fig*dfig xfigw yfigw]);
figure(fig);
data = data_wg.';
plot(ts, real(data(:)), 'r', ts, imag(data(:)), 'b')
ylim([-1.3 1.3]);
title(sprintf('WG output data (WG channel %6.3f)', tb.channel_wg));
xlabel(sprintf('Time 0:%d [Tchan]', tb.nof_tchan-1));
ylabel('Voltage');
grid on;
%% Plot FIR-filter coefficients
fig=fig+1;
figure('position', [xfig+fig*dfig yfig-fig*dfig xfigw yfigw]);
figure(fig);
h = hfir_channel_coeff;
NL = ctrl_pfir_channel.nof_coefficients;
N = ctrl_pfir_channel.nof_polyphases;
L = ctrl_pfir_channel.nof_taps;
hi = 1:NL; % coefficients index
hx = hi / N; % coefficients axis in tap units
if ctrl_pfir_channel.coeff_w == 0
plot(hx, h);
title(['FIR filter coefficients for ', ctrl_pfir_channel.config.design, ' (floating point)']);
else
plot(hx, h);
title(['FIR filter coefficients for ', ctrl_pfir_channel.config.design, ' (', num2str(ctrl_pfir_channel.coeff_w), ' bits)']);
end
ylo = min(h)*1.2;
yhi = max(h)*1.2;
ylim([ylo yhi]);
grid on;
%% Plot channel PFIR transfer function
h = hfir_channel_coeff;
NL = ctrl_pfir_channel.nof_coefficients;
L = ctrl_pfir_channel.nof_taps;
hf_abs = abs(fftshift(fft(h / sum(h), NL)));
hi = 1:NL; % coefficients index
fx = (hi - NL/2-1) / L; % frequency axis in channel units
fig=fig+1;
figure('position', [xfig+fig*dfig yfig-fig*dfig xfigw yfigw]);
figure(fig);
plot(fx, db(hf_abs)); % db() = 20*log10 for voltage
%xlim([-3 3]);
grid on;
title(['Channel PFIR filter transfer function for ', ctrl_pfir_channel.config.design, ' (nof taps = ', int2str(ctrl_pfir_channel.nof_taps), ')']);
xlabel('Frequency [channels]');
ylabel('Magnitude [dB]');
%% Plot channel PFIR output
fig=fig+1;
figure('position', [xfig+fig*dfig yfig-fig*dfig xfigw yfigw]);
figure(fig);
data = data_pfir_channel.';
plot(ts, real(data(:)), 'r', ts, imag(data(:)), 'b')
title(sprintf('Channel PFIR filter output - FFT input data (WG channel %6.3f)', tb.channel_wg));
ylim([-2 2]); % Delay tracking step when tb.channel_wg is .5 causes double range
xlabel(sprintf('Time 0:%d [Tchan]', tb.nof_tchan-1));
ylabel('Voltage');
grid on;
%% Plot PFFT channels spectrum and phase for all tb.nof_tchan in one plot
chan_ampl = abs(data_pfft_channel);
chan_ampl_max = max(chan_ampl(:));
chan_phase = angle(data_pfft_channel)*180/pi;
x = chan_ampl < 0.1*chan_ampl_max;
chan_phase(x) = 0; % force phase of too small signals to 0
fig=fig+1;
figure('position', [xfig+fig*dfig yfig-fig*dfig xfigw yfigw]);
figure(fig);
subplot(2,1,1);
data = chan_ampl.';
data = data(:);
plot(tchan_all, data, 'k', tchan_all(chan_I), data(chan_I), 'ko', tchan_all(chan_Iplus1), data(chan_Iplus1), 'kx');
title(sprintf('Channel data - amplitude (o,x = channel %d,%d for WG channel = %6.3f)', tb.channel_i, tb.channel_iplus1, tb.channel_wg));
xlabel(sprintf('Channels 0:%d at time 0:%d [Tchan]', tb.nof_channels-1, tb.nof_tchan-1));
ylabel('Voltage');
grid on;
subplot(2,1,2);
data = chan_phase.';
data = data(:);
plot(tchan_all, data, 'k', tchan_all(chan_I), data(chan_I), 'ko', tchan_all(chan_Iplus1), data(chan_Iplus1), 'kx');
ylim([-180 180])
title(sprintf('Channel data - amplitude (o,x = channel %d,%d for WG channel = %6.3f)', tb.channel_i, tb.channel_iplus1, tb.channel_wg));
xlabel(sprintf('Channels 0:%d at time 0:%d [Tchan]', tb.nof_channels-1, tb.nof_tchan-1));
ylabel('Phase [degrees]');
grid on;
%% Plot phase for the channel that is set in the WG
wg_chan_phase = angle(data_pfft_channel(:, tb.channel_I))*180/pi;
fig=fig+1;
figure('position', [xfig+fig*dfig yfig-fig*dfig xfigw yfigw]);
figure(fig);
plot(tchan_one, wg_chan_phase, '-o')
ylim([-180 180])
title(sprintf('Channel phase for channel %d in range 0:%d', tb.channel_i, tb.nof_channels-1));
xlabel(sprintf('Time 0:%d [Tchan]', tb.nof_tchan-1));
ylabel('Phase [degrees]');
grid on;
%% Plot PFFT channels spectrum and phase for all tb.nof_tchan in separate plots
chan_ampl = abs(data_pfft_channel);
chan_ampl_max = max(chan_ampl(:));
chan_phase = angle(data_pfft_channel)*180/pi;
x = chan_ampl < 0.1*chan_ampl_max;
chan_phase(x) = 0; % force phase of too small signals to 0
fig=fig+1;
figure('position', [xfig+fig*dfig yfig-fig*dfig xfigw yfigw]);
figure(fig);
if tb.plot_per_block
for bi = 0:tb.nof_tchan-1
bI = bi+1;
subplot(2,1,1);
data = chan_ampl(bI, :);
plot(0:tb.nof_channels-1, data, 'k', tb.channel_i, data(tb.channel_I), 'ko', tb.channel_iplus1, data(tb.channel_Iplus1), 'kx');
ylim([0 chan_ampl_max])
title(sprintf('Channel data - amplitude (o,x = channel %d,%d for WG channel = %6.3f)', tb.channel_i, tb.channel_iplus1, tb.channel_wg));
xlabel(sprintf('Channels 0:%d at time %d [Tchan]', tb.nof_channels-1, bi));
ylabel('Voltage');
grid on;
subplot(2,1,2);
data = chan_phase(bI, :);
plot(0:tb.nof_channels-1, data, 'k', tb.channel_i, data(tb.channel_I), 'ko', tb.channel_iplus1, data(tb.channel_Iplus1), 'kx');
ylim([-180 180])
title(sprintf('Channel data - amplitude (o,x = channel index %d,%d for WG channel = %6.3f)', tb.channel_i, tb.channel_iplus1, tb.channel_wg));
xlabel(sprintf('Channels 0:%d at time %d [Tchan]', tb.nof_channels-1, bi));
ylabel('Phase [degrees]');
grid on;
drawnow;
%pause(0.3)
end
end
%% Plot channel spectrogram for all tb.nof_tchan in one plot
data = db(chan_ampl); % no need to scale data, range is already normalized
if strcmp(tb.model_quantization, 'floating point')
db_low = -150;
else
db_low = -20 - floor(6.02 * ctrl_pfft_channel.data_w);
end
%db_low = floor(min(data(data~=-Inf)))
data(data<db_low) = db_low;
mymap = jet(-db_low);
%mymap(1,:) = [0 0 0]; % force black for dB(0)
colormap(mymap);
imagesc(data',[db_low 0]);
colorbar;
title(sprintf('Channel spectogram (max value = %f = %.2f dB)', chan_ampl_max, db(chan_ampl_max)));
xlabel(sprintf('Time 0:%d [Tchan]', tb.nof_tchan-1));
ylabel(sprintf('Channels 0:%d', tb.nof_channels-1));
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment