diff --git a/applications/apertif/matlab/apertif_matlab_readme.txt b/applications/apertif/matlab/apertif_matlab_readme.txt
index 4a9ace53c93919a69e5a878ce95c110fa1be7bf4..1e4164c3cfd841826325107281e36d22af5ab057 100644
--- a/applications/apertif/matlab/apertif_matlab_readme.txt
+++ b/applications/apertif/matlab/apertif_matlab_readme.txt
@@ -40,7 +40,7 @@ This directory contains Matlab code for modelling DSP aspects of Apertif:
try_fft.m Script Try DFT and DTFT
-5) Data path functions
+5) Data path HDL 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
diff --git a/applications/apertif/matlab/one_pfb.m b/applications/apertif/matlab/one_pfb.m
index cc560b30f8559fb7527096c819954fa91ed1e79f..702f6258e483587050f66b099a22aaa612b1e4ae 100644
--- a/applications/apertif/matlab/one_pfb.m
+++ b/applications/apertif/matlab/one_pfb.m
@@ -21,7 +21,7 @@
% Author: W. Lubberhuizen, 2015 (top.m for LOFAR Station, uses classes)
% Eric Kooistra, 2016 (for Apertif)
%
-% Purpose : Model one PFB with real input:
+% Purpose : Model one PFB = PFIR + PFFT with real input:
% 1) to show the transfer function of a polyphase filterbank with real
% input like the subband filterbank
% 2) to create reference data for HDL implementation in data/:
@@ -95,7 +95,7 @@ if strcmp(tb.model_filterbank, 'LOFAR')
tb.nof_subbands = 512;
else
tb.model_quantization = 'floating point';
- tb.model_quantization = 'fixed point';
+ %tb.model_quantization = 'fixed point';
tb.nof_subbands = 64;
end
@@ -106,7 +106,7 @@ tb.subband_wg = 47; % subband range 0:tb.nof_subbands-1, can be
% Model a frequency sweep of the 'sinusoid'
tb.chirp = 0; % 0 = use fixed tb.subband_wg frequency, else increment WG frequency every block to have chirp frequency sweep
-tb.chirp = 1;
+%tb.chirp = 1;
if tb.chirp
tb.nof_tsub = 200; % number of subband periods to simulate
else
@@ -139,7 +139,7 @@ disp(sprintf (['Test bench settings:\n', ...
'. Number of subband periods to simulate = %d\n', ...
'. Subband filterbank real FFT size = %d\n', ...
'. WG subband frequency = %f\n'], ...
- tb.model_filterbank, tb.model_quantization, tb.nof_subbands, tb.nof_tsub, tb.subband_fft_size, tb.subband_wg))
+ tb.model_filterbank, tb.model_quantization, tb.nof_subbands, tb.nof_tsub, tb.subband_fft_size, tb.subband_wg));
%% Waveform generator
@@ -156,9 +156,10 @@ 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.6*lsb;
ctrl_wg.agwn_sigma = 0; % AGWN sigma
-ctrl_wg.ampl = 1; % full scale is 1, ampl > 1 causes analogue clipping
+ctrl_wg.ampl = 1; % full scale maps to 1, so 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;
@@ -181,28 +182,32 @@ else
end
%% Subband PFIR filter parameters
-ctrl_pfir_subband.nof_polyphases= tb.subband_fft_size;
+ctrl_pfir_subband.bypass = 0;
+ctrl_pfir_subband.nof_polyphases = tb.subband_fft_size;
+
+if strcmp(tb.model_quantization, 'floating point')
+ ctrl_pfir_subband.data_w = 0; % no quantization
+ ctrl_pfir_subband.coeff_w = 0; % no quantization
+else
+ ctrl_pfir_subband.data_w = 16;
+ ctrl_pfir_subband.backoff_w = 1; % attenuate PFIR data output by 2^backoff_w before quantization to account for PFIR overshoot
+ ctrl_pfir_subband.coeff_w = 16; % nof bits per coefficient including 1 sign bit
+ ctrl_pfir_subband.coeff_q_max = 2^(ctrl_pfir_subband.coeff_w-1)-1;
+ ctrl_pfir_subband.coeff_q_full_scale = 2^(ctrl_pfir_subband.coeff_w-1);
+ ctrl_pfir_subband.data_q_full_scale = 2^(ctrl_pfir_subband.data_w-1-ctrl_pfir_subband.backoff_w);
+end
+
if strcmp(tb.model_filterbank, 'LOFAR')
- % Load the quantized PFIR coefficients of LOFAR station
+ % Load the quantized PFIR coefficients of LOFAR station (fixed point)
ctrl_pfir_subband.nof_taps = 16; % Number of taps
ctrl_pfir_subband.nof_coefficients = ctrl_pfir_subband.nof_polyphases*ctrl_pfir_subband.nof_taps; % Number of filter coefficients (taps)
- ctrl_pfir_subband.backoff_w = 1;
ctrl_pfir_subband.data_w = 16;
ctrl_pfir_subband.config.design = 'lofar file';
hfir_subband_coeff = load('data/Coeffs16384Kaiser-quant.dat');
else
- % Calculate free choice PFIR coefficients
+ % Calculate free choice PFIR coefficients (can be floating point or fixed point)
ctrl_pfir_subband.nof_taps = 16; % Number of taps
ctrl_pfir_subband.nof_coefficients = ctrl_pfir_subband.nof_polyphases*ctrl_pfir_subband.nof_taps; % Number of filter coefficients (taps)
- if strcmp(tb.model_quantization, 'floating point')
- ctrl_pfir_subband.coeff_w = 0;
- ctrl_pfir_subband.backoff_w = 1;
- ctrl_pfir_subband.data_w = 0;
- else
- ctrl_pfir_subband.coeff_w = 16; % bits/coefficient = 1 sign + (w-1) mantissa
- ctrl_pfir_subband.backoff_w = 1; % normalize to 2^backoff_w before quantization to account for PFIR overshoot
- ctrl_pfir_subband.data_w = 16;
- end
ctrl_pfir_subband.config.design = 'fir1';
ctrl_pfir_subband.config.design = 'fircls1'; % 'fir1', 'fircls1'
ctrl_pfir_subband.config.design_flag = ''; % only for fircls1: 'trace', ''
@@ -220,28 +225,41 @@ else
ctrl_pfir_subband.coeff_w, ...
ctrl_pfir_subband.config);
end
-
ctrl_pfir_subband.coeff = reshape(hfir_subband_coeff, ctrl_pfir_subband.nof_polyphases, ctrl_pfir_subband.nof_taps);
ctrl_pfir_subband.coeff = flipud(ctrl_pfir_subband.coeff);
ctrl_pfir_subband.Zdelays = zeros(ctrl_pfir_subband.nof_polyphases, ctrl_pfir_subband.nof_taps-1);
-ctrl_pfir_subband.gain = sum(ctrl_pfir_subband.coeff(:)) / ctrl_pfir_subband.nof_polyphases;
-ctrl_pfir_subband.bypass = 0;
+ctrl_pfir_subband.gain = 1; % no gain adjustment in PFIR, just apply the coefficients to the input data
disp(sprintf(['Subband PFIR filter settings:\n', ...
'. Bypass = %d\n', ...
'. Number of polyphases = %d\n', ...
'. Number of taps = %d\n'], ...
- ctrl_pfir_subband.bypass, ctrl_pfir_subband.nof_polyphases, ctrl_pfir_subband.nof_taps))
-
+ ctrl_pfir_subband.bypass, ...
+ ctrl_pfir_subband.nof_polyphases, ...
+ ctrl_pfir_subband.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_subband.coeff_w, ...
+ ctrl_pfir_subband.coeff_q_max, ...
+ ctrl_pfir_subband.data_w, ...
+ ctrl_pfir_subband.backoff_w, ...
+ ctrl_pfir_subband.data_q_full_scale));
+end
+
%% Subband FFT parameters
ctrl_pfft_subband.fft_size = tb.subband_fft_size;
ctrl_pfft_subband.complex = false; % real input PFB
-ctrl_pfft_subband.gain = ctrl_pfft_subband.fft_size;
+ctrl_pfft_subband.gain = ctrl_pfft_subband.fft_size; % compensate default gain of FFT in PFFT
if strcmp(tb.model_quantization, 'floating point')
ctrl_pfft_subband.data_w = 0;
else
ctrl_pfft_subband.data_w = 16;
end
+ctrl_pfft_subband.data_q_full_scale = 2^(ctrl_pfft_subband.data_w-1);
%% Run the data path processing
@@ -263,42 +281,39 @@ for bi = 0:tb.nof_tsub-1
data_pfft_subband(bI, :) = block_pfft_subband;
end
t_stop = cputime;
-disp(sprintf('Total processing time: %f seconds', t_stop-t_start));
+disp(sprintf('Total processing time: %f seconds\n', t_stop-t_start));
+
%% Save coeffcients and data reference files for stimuli and verification of HDL implementation
if strcmp(tb.model_quantization, 'fixed point')
file_name_prefix = 'data/one_pfb_m_';
L = ctrl_pfir_subband.nof_taps;
- N = tb.subband_fft_size;
- coeff_w = ctrl_pfir_subband.coeff_w;
+ N = ctrl_pfir_subband.nof_polyphases; % = tb.subband_fft_size;
+
+ % Quantize subband PFIR filter coefficients
+ int_hfir = round(ctrl_pfir_subband.coeff_q_full_scale * hfir_subband_coeff); % map largest coefficient to q_max
+
+ % Quantize WG output
+ wg_q_data = round(ctrl_wg.q_full_scale * data_wg);
- % Save the quantized integer subband PFIR filter coefficients
- hfir = hfir_subband_coeff;
- hfir_norm = 1/ctrl_pfir_subband.gain;
- hfir_q_fule_scale = 2^(coeff_w-1)-1;
- hfir_q_norm = 1/max(hfir);
- int_hfir = round(hfir_q_fule_scale * hfir * hfir_q_norm); % use full scale hfir_q_fule_scale for largest coefficient normalized by hfir_q_norm
- file_name = sprintf([file_name_prefix, 'pfir_coeff_', ctrl_pfir_subband.config.design, '_%dtaps_%dpoints_%db.dat'], L, N, coeff_w);
+ % Quantize PFIR subband output
+ pfir_subband_q_data = round(ctrl_pfir_subband.data_q_full_scale * data_pfir_subband);
+
+ % Quantize PFFT subband output real and imaginary
+ pfft_subband_q_data = round(ctrl_pfft_subband.data_q_full_scale * data_pfft_subband);
+
+ % Save the quantized integer subband PFIR filter coefficients in a separate file
+ file_name = sprintf([file_name_prefix, 'pfir_coeff_', ctrl_pfir_subband.config.design, '_%dtaps_%dpoints_%db.dat'], L, N, ctrl_pfir_subband.coeff_w);
fid = fopen(file_name, 'w');
fprintf(fid,'%d\n', int_hfir);
fclose(fid);
- % Save the data path signals
+ % Save PFIR filter coefficients again together with the data path signals for WG, PFIR and PFFT
tbegin = 1;
tend = tb.nof_tsub; % <= tb.nof_tsub
- % WG output
- wg_q_full_scale = 2^(ctrl_wg.data_w-1);
- wg_q_data = round(wg_q_full_scale * data_wg); % use full scale q_max for data_wg value 1.0
- % PFIR subband output
- pfir_subband_q_full_scale = 2^(ctrl_pfir_subband.data_w-1-ctrl_pfir_subband.backoff_w);
- pfir_subband_q_norm = hfir_q_norm / hfir_norm;
- pfir_subband_q_data = round(pfir_subband_q_full_scale * data_pfir_subband * pfir_subband_q_norm);
- % PFFT subband output real and imaginary
- pfft_subband_q_full_scale = 2^(ctrl_pfft_subband.data_w-1);
- pfft_subband_q_data = round(pfft_subband_q_full_scale * data_pfft_subband);
file_name = sprintf([file_name_prefix, '%s_%db_%dtaps_%dpoints_%db_%db.dat'], ctrl_wg.name, ctrl_wg.data_w, L, N, ctrl_pfir_subband.data_w, ctrl_pfft_subband.data_w);
fid = fopen(file_name, 'w');
- fprintf(fid,'Nof lines PFIR coefficients: %d\n', length(int_hfir)); % save again also in the data file to have a complete set
+ 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_subband_q_data(:)));
fprintf(fid,'Nof lines PFFT output: %d\n', length(pfft_subband_q_data(:)));
@@ -352,6 +367,28 @@ xlabel(sprintf('Time 0:%d [Tsub]', tb.nof_tsub-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_subband_coeff;
+NL = ctrl_pfir_subband.nof_coefficients;
+N = ctrl_pfir_subband.nof_polyphases;
+L = ctrl_pfir_subband.nof_taps;
+hi = 1:NL; % coefficients index
+hx = hi / N; % coefficients axis in tap units
+if ctrl_pfir_subband.coeff_w == 0
+ plot(hx, h);
+ title(['FIR filter coefficients for ', ctrl_pfir_subband.config.design, ' (floating point)']);
+else
+ plot(hx, h);
+ title(['FIR filter coefficients for ', ctrl_pfir_subband.config.design, ' (', num2str(ctrl_pfir_subband.coeff_w), ' bits)']);
+end
+ylo = min(h)*1.2;
+yhi = max(h)*1.2;
+ylim([ylo yhi]);
+grid on;
+
%% Plot subband PFIR transfer function
h = hfir_subband_coeff;
NL = ctrl_pfir_subband.nof_coefficients;
diff --git a/applications/apertif/matlab/run_pfft.m b/applications/apertif/matlab/run_pfft.m
index 54fe1ae103d376eb17492448cdc3f9f9a3439ae0..f1fe458fde1cd7edd08ae1a175651b577196afa1 100644
--- a/applications/apertif/matlab/run_pfft.m
+++ b/applications/apertif/matlab/run_pfft.m
@@ -18,8 +18,7 @@
% along with this program. If not, see <http://www.gnu.org/licenses/>.
%
%-----------------------------------------------------------------------------
-% Author: W. Lubberhuizen, 2006 (top.m for LOFAR Station, uses classes)
-% Eric Kooistra, 2016 (for Apertif)
+% Author: Eric Kooistra, 2016 (for Apertif)
%
% Purpose : Run the PFFT with real input:
% To create input and output reference data for HDL implementation in
diff --git a/applications/apertif/matlab/run_pfir_coeff.m b/applications/apertif/matlab/run_pfir_coeff.m
index 4ca5e1547e17874ee4a71419ba401e9554fd027a..094b8ad4613f5dfd54f2ecca4f912495fcec5f82 100644
--- a/applications/apertif/matlab/run_pfir_coeff.m
+++ b/applications/apertif/matlab/run_pfir_coeff.m
@@ -132,16 +132,16 @@ hx = hi / N; % coefficients axis in tap units
fx = (hi - NL/2-1) / L; % frequency axis in channel units
% Print FIR-filter parameters
-sprintf([ ...
- 'Filter design : %s\n', ...
- 'FFT size : N = %d\n', ...
- 'Number of FIR taps : L = %d\n', ...
- 'Number of FIR coefficients : N*L = %d\n', ...
- 'Sample rate : fs = %6.2f [Hz]\n' , ...
- 'Channel frequency : f_chan = %6.2f [Hz]\n' , ...
- 'Channel bandwidth : BWchan = %6.2f [Hz]\n' , ...
- 'Channel half power factor : hp_factor = %6.3f\n'] , ...
- config.design, N, L, NL, fs, f_chan, BWchan*fs, hp_factor)
+disp(sprintf([ ...
+ 'Filter design : %s\n', ...
+ 'Number of polyphases (FFT size) : N = %d\n', ...
+ 'Number of FIR taps : L = %d\n', ...
+ 'Number of FIR coefficients : N*L = %d\n', ...
+ 'Sample rate : fs = %6.2f [Hz]\n' , ...
+ 'Channel frequency : f_chan = %6.2f [Hz]\n' , ...
+ 'Channel bandwidth : BWchan = %6.2f [Hz]\n' , ...
+ 'Channel half power factor : hp_factor = %6.3f\n'] , ...
+ config.design, N, L, NL, fs, f_chan, BWchan*fs, hp_factor));
if strcmp(config.design, 'fircls1')
sprintf([ ...
'Pass band deviation : r_pass = %9.7f = %6.2f dB ripple\n' , ...
@@ -158,7 +158,7 @@ if coeff_w>0
int_coeff = coeff;
else
q_max = 2^(coeff_w-1)-1;
- int_coeff = round(coeff/max(coeff) * q_max);
+ int_coeff = round(q_max * coeff/max(coeff));
end
file_name_prefix = 'data/run_pfir_coeff_m_';
file_name = sprintf([file_name_prefix, config.design, '_%dtaps_%dpoints_%db.dat'], L, N, coeff_w);