silas/imdd_silas
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Changes from mwork PC.

Browse Source 
PDP 2025

MPI analysis

new focus on database and SQL
This commit is contained in:
Silas Oettinghaus
2025-03-21 08:11:40 +01:00
parent 402e491506
commit 74066d0669
36 changed files with 2234 additions and 620 deletions
Download Patch File Download Diff File Expand all files Collapse all files

141
Functions/EQ_structures/duobinary_target.m
View File

@@ -1,36 +1,42 @@
function [eq_package] = duobinary_target(eq_, mlse_,M, rx_signal, tx_symbols, tx_bits, options)
arguments
eq_
mlse_
M
rx_signal
tx_symbols
tx_bits
eq_
mlse_
M
rx_signal
tx_symbols
tx_bits
options.precode_mode db_mode
options.showAnalysis = 0;
options.eth_style_symbol_mapping = 0;
options.postFFE = [];
end
%Duobinary Targeting
[eq_signal, eq_noise] = eq_.process(rx_signal,Duobinary().encode(tx_symbols));
%Duobinary Targeting
db_ref_sequence = Duobinary().encode(tx_symbols);
db_ref_constellation = unique(db_ref_sequence.signal);
[eq_signal, eq_noise] = eq_.process(rx_signal,db_ref_sequence);
% dir = [1,1];
mlse_sig_sd = mlse_.process(eq_signal);
if ~isempty(options.postFFE)
[eq_signal,eq_noise] = options.postFFE.process(eq_signal,db_ref_sequence);
end
mlse_sig_hd = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).quantize(mlse_sig_sd);
% precoding to mitigate error propagation, most prominently used in
% combination with duobinary signaling to avoid catastrophic error
% behavior (see J.W.M. Bergmans, Digital Baseband Transmission and Recording -> partial response signaling)
% dir = [1,1];
mlse_sig_sd = mlse_.process(eq_signal);
% takes:
% -> eq_signal_hd: hard decision signal after eq
% -> tx_symbols: that where used as reference for eq
mlse_sig_hd = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).quantize(mlse_sig_sd);
switch options.precode_mode
case db_mode.db_emulate
% precoding to mitigate error propagation, most prominently used in
% combination with duobinary signaling to avoid catastrophic error
% behavior (see J.W.M. Bergmans, Digital Baseband Transmission and Recording -> partial response signaling)
% takes:
% -> eq_signal_hd: hard decision signal after eq
% -> tx_symbols: that where used as reference for eq
switch options.precode_mode
case db_mode.db_emulate
mlse_sig_hd = Duobinary().encode(mlse_sig_hd,"M",M);
mlse_sig_hd = Duobinary().decode(mlse_sig_hd,"M",M);
@@ -40,43 +46,102 @@ end
tx_bits = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(tx_symbols_precoded);
case db_mode.db_discard
case db_mode.db_discard
% normal dsp for precoded sequence == discard/omit/ignore precode
tx_bits = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(tx_symbols);
case db_mode.db_encoded
case db_mode.db_encoded
% normal DB encoded data (only for 10KM)
case db_mode.db_precoded
case db_mode.db_precoded
mlse_sig_hd = Duobinary().encode(mlse_sig_hd,"M",M);
mlse_sig_hd = Duobinary().decode(mlse_sig_hd,"M",M);
end
end
% M = numel(unique(tx_symbols.signal));
rx_bits = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(mlse_sig_hd);
% M = numel(unique(tx_symbols.signal));
rx_bits = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(mlse_sig_hd);
[~,numErrors,ber,~] = calc_ber(rx_bits.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
[bits_db,errors_db,ber_db,errorIndice_db] = calc_ber(rx_bits.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
eq_package.ber = ber_db;
eq_package.ber = ber;
resultsDBtgt = struct( ...
'result_id', NaN, ... %
'run_id', NaN, ... % Beispielhafte Run-ID
'eqParam_id', NaN, ... % Beispielhafter Fremdschlüssel zur EqualizerParameters-Tabelle
'date_of_processing', datetime('now'), ... % Aktuelles Datum und Uhrzeit
'numBits', bits_db, ... % Beispiel: 1.000.000 Bits
'numBitErr', errors_db, ... % Beispiel: 120 Bitfehler
'BER', ber_db, ... % BER = 120 / 1.000.000
'SNR', [], ... % Beispielhafte SNR
'SNR_level', jsonencode([]), ... % SNR-Level als JSON-codiertes Array
'GMI', [], ... % Beispielhafter GMI-Wert
'AIR', [], ... % Beispielhafter AIR-Wert
'EVM', [], ... % Beispielhafte EVM
'EVM_level', jsonencode([]), ... % EVM-Level als JSON-codiertes Array
'Alpha', [] ... % Beispielhafter Alpha-Wert
);
if options.showAnalysis
eq_noise = eq_noise - mean(eq_noise.signal);
if ~isempty(options.postFFE)
npostFFE = options.postFFE.order;
else
npostFFE = 0;
end
rx_signal.spectrum("normalizeTo0dB",1,"fignum",250,"displayname","Rx Spectrum");
equalizerConfigDBtgt = struct( ...
'eq_id', NaN, ... % Auto-Inkrement, wird in der DB gesetzt
'equalizer_structure', int32(equalizer_structure.vnle_db_mlse), ... % Beispiel: 1 (z.B. für vnle)
'M', M, ... % Ordnung der PAM-Konstellation
'target_constellation', jsonencode(round(db_ref_constellation,5)), ... % Beispielhafter Target-String
'db_target', 1, ... % 0 oder 1
'diff_precode', int32(options.precode_mode), ... % 0 oder 1
'postFFE', ~isempty(options.postFFE), ... % Beispielwert
'NpostFFE', npostFFE, ... % Beispielwert
'Ne1', eq_.Ne(1), ... % Feedforward Koeffizienten 1. Ordnung
'Ne2', eq_.Ne(2), ... % Feedforward Koeffizienten 2. Ordnung
'Ne3', eq_.Ne(3), ... % Feedforward Koeffizienten 3. Ordnung
'Nb1', eq_.Nb(1), ... % Decision Feedback Koeffizienten 1. Ordnung
'Nb2', eq_.Nb(2), ... % Decision Feedback Koeffizienten 2. Ordnung
'Nb3', eq_.Nb(3), ... % Decision Feedback Koeffizienten 3. Ordnung
'K', eq_.K, ... % Samples pro Symbol
'DCmu', eq_.DCmu, ... % Anpassungsrate für DC-Tap
'ideal_dfe', eq_.ideal_dfe, ... % Flag für ideal DFE (0 oder 1)
'training_length', eq_.training_length, ... % Anzahl Trainingssymbole
'training_loops', eq_.training_loops, ... % Anzahl Trainingsdurchläufe
'TRmu1', eq_.FFEmu, ... % mu für DD-Modus (1. Ordnung)
'TRmu2', eq_.FFEmu, ... % mu für DD-Modus (2. Ordnung)
'TRmu3', eq_.FFEmu, ... % mu für DD-Modus (3. Ordnung)
'TRmuDFE', eq_.DFEmu, ... % mu für DFE-Modus im DD
'dd_loops', 5, ... % Anzahl Durchläufe im DD-Modus
'DDmu1', eq_.DDmu(1), ... % mu für DD-Modus (1. Ordnung)
'DDmu2', eq_.DDmu(2), ... % mu für DD-Modus (2. Ordnung)
'DDmu3', eq_.DDmu(3), ... % mu für DD-Modus (3. Ordnung)
'DDmuDFE', eq_.DDmu(4), ... % mu für DFE-Modus im DD
'MLSE_mode', 'viterbi', ... % Beispiel: MLSE-Modus als String
'MLSE_trellis_states', jsonencode(mlse_.trellis_states), ... % Trellis-States, z.B. als JSON-String oder kommasepariert
'comment', 'function: duobinary_target.m', ... % Zusätzliche Kommentare
'config_hash', NaN ...
);
Duobinary().encode(tx_symbols).spectrum("normalizeTo0dB",1,"fignum",250,"displayname","DB encoded reference");
eq_package.resultsDBtgt = resultsDBtgt;
eq_package.equalizerConfigDBtgt = equalizerConfigDBtgt;
showEQNoisePSD(eq_noise,"fignum",250,"displayname",'Duobinary Target Noise after Equalization');
if options.showAnalysis
eq_noise = eq_noise - mean(eq_noise.signal);
fprintf('DB tgt BER: %.2e \n',ber);
rx_signal.spectrum("normalizeTo0dB",1,"fignum",250,"displayname","Rx Spectrum");
end
Duobinary().encode(tx_symbols).spectrum("normalizeTo0dB",1,"fignum",250,"displayname","DB encoded reference");
showEQNoisePSD(eq_noise,"fignum",250,"displayname",'Duobinary Target Noise after Equalization');
fprintf('DB tgt BER: %.2e \n',ber);
end
end

30
Functions/EQ_structures/vnle.m
View File

@@ -23,6 +23,7 @@ function [eq_package] = vnle(eq_,M,rx_signal,tx_symbols,tx_bits,options)
options.precode_mode db_mode
options.showAnalysis = 0
options.eth_style = 0;
options.postFFE = [];
end
%FFE or VNLE
@@ -32,6 +33,10 @@ function [eq_package] = vnle(eq_,M,rx_signal,tx_symbols,tx_bits,options)
end
[eq_signal_sd,eq_noise] = eq_.process(rx_signal,tx_symbols);
if ~isempty(options.postFFE)
[eq_signal_sd,eq_noise] = options.postFFE.process(eq_signal_sd,tx_symbols);
end
eq_signal_hd = PAMmapper(M,0).quantize(eq_signal_sd);
% precoding to mitigate error propagation, most prominently used in
@@ -89,15 +94,34 @@ function [eq_package] = vnle(eq_,M,rx_signal,tx_symbols,tx_bits,options)
fprintf('SNR: %d dB \n',snr(eq_signal_sd.signal,eq_noise.signal));
if M == 6
logm = 2.5;
else
logm = log2(M);
end
fprintf('NGMI: %.4f \n', inf_rate/logm);
fprintf(['VNLE EVM lvl: ',repmat('%.3f ',1,numel(evm_lvl)),' \n'],evm_lvl);
fprintf('VNLE BER: %.2e \n',ber);
showEQcoefficients(eq_.e,eq_.e2,eq_.e3,"displayname",'Coefficients');
disp("%%%%%%%%%%%%%%%%%%%%%")
showEQNoiseSNR(tx_symbols,eq_noise,"displayname",'SNR after VNLE','fignum',301);
% showEQcoefficients('n1',eq_.e,'n2',eq_.e2,'n3',eq_.e3,"displayname",'Coefficients');
%
% if ~isempty(options.postFFE)
% showEQcoefficients('n1',options.postFFE.e,"displayname",'Coefficients');
% end
%
% showEQNoisePSD(eq_noise);
%
% showEQfilter(eq_.e,eq_signal_sd.fs.*2)
% noiselessness(tx_symbols,eq_noise,"displayname",'SNR after VNLE','fignum',301);
% showLevelHistogram(eq_signal_sd,tx_symbols,"fignum",302);
showLevelHistogram(eq_signal_sd,tx_symbols,"fignum",302);
end

144
Functions/EQ_structures/vnle_postfilter_mlse.m
View File

@@ -11,11 +11,17 @@ arguments
options.precode_mode db_mode
options.showAnalysis = 0;
options.eth_style_symbol_mapping = 0;
options.postFFE = [];
options.database = [];
end
%FFE or VNLE
[eq_signal_sd,eq_noise] = eq_.process(rx_signal,tx_symbols);
if ~isempty(options.postFFE)
[eq_signal_sd,eq_noise] = options.postFFE.process(eq_signal_sd,tx_symbols);
end
eq_signal_hd = PAMmapper(M,0).quantize(eq_signal_sd);
mlse_sig_sd = pf_.process(eq_signal_sd,eq_noise);
@@ -70,28 +76,152 @@ end
% METRICS OF VNLE %
rx_bits_vnle = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(eq_signal_hd);
[~,~,ber_vnle,~] = calc_ber(rx_bits_vnle.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
[bits_vnle,errors_vnle,ber_vnle,errorIndice_vnle] = calc_ber(rx_bits_vnle.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
[snr_vnle,snr_vnle_lvl] = calc_snr(tx_symbols.signal,eq_noise.signal);
% correct TUM implementation of AIR
[inf_rate_vnle] = calc_air(eq_signal_sd,tx_symbols,"skip_front",10000,"skip_end",10000);
[gmi_vnle] = calc_air(eq_signal_sd,tx_symbols,"skip_front",10000,"skip_end",10000);
gmi_vnle = min(max(abs(gmi_vnle),0.1),0); %set to zero if no convergence of gmi below 0.1 to avoid negative or any other "dumb" value
air_vnle = tx_symbols.fs .* floor(log2(8)*10)/10 .* gmi_vnle;
[evm_vnle_total,evm_vnle_lvl] = calc_evm(eq_signal_sd,tx_symbols);
% METRICS OF MLSE (HD-VITERBI)
rx_bits_mlse = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(mlse_sig_hd);
[~,~,ber_mlse,~] = calc_ber(rx_bits_mlse.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
[bits_mlse,errors_mlse,ber_mlse,errorIndice_mlse]= calc_ber(rx_bits_mlse.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
pf_.ncoeff = 1;
pf_.process(eq_signal_sd,eq_noise);
alpha = pf_.coefficients(2);
eq_package.ber_mlse = ber_mlse;
eq_package.ber_vnle = ber_vnle;
eq_package.evm_vnle_total = evm_vnle_total;
eq_package.evm_vnle_lvl = evm_vnle_lvl;
eq_package.air = inf_rate_vnle;
eq_package.gmi = gmi_vnle;
eq_package.eq = eq_;
eq_package.pf = pf_;
eq_package.mlse = mlse_;
resultsVNLE = struct( ...
'result_id', NaN, ... %
'run_id', NaN, ... % Beispielhafte Run-ID
'eqParam_id', NaN, ... % Beispielhafter Fremdschlüssel zur EqualizerParameters-Tabelle
'date_of_processing', datetime('now'), ... % Aktuelles Datum und Uhrzeit
'numBits', bits_vnle, ... % Beispiel: 1.000.000 Bits
'numBitErr', errors_vnle, ... % Beispiel: 120 Bitfehler
'BER', ber_vnle, ... % BER = 120 / 1.000.000
'SNR', snr_vnle, ... % Beispielhafte SNR
'SNR_level', jsonencode(snr_vnle_lvl), ... % SNR-Level als JSON-codiertes Array
'GMI', gmi_vnle, ... % Beispielhafter GMI-Wert
'AIR', air_vnle, ... % Beispielhafter AIR-Wert
'EVM', evm_vnle_total, ... % Beispielhafte EVM
'EVM_level', jsonencode(evm_vnle_lvl), ... % EVM-Level als JSON-codiertes Array
'Alpha', [] ... % Beispielhafter Alpha-Wert
);
if ~isempty(options.postFFE)
npostFFE = options.postFFE.order;
else
npostFFE = 0;
end
equalizerConfigVNLE = struct( ...
'eq_id', NaN, ... % Auto-Inkrement, wird in der DB gesetzt
'equalizer_structure', int32(equalizer_structure.vnle), ... % Beispiel: 1 (z.B. für vnle)
'M', M, ... % Ordnung der PAM-Konstellation
'target_constellation', jsonencode(round(unique(tx_symbols.signal),5)), ... % Beispielhafter Target-String
'db_target', 0, ... % 0 oder 1
'diff_precode', int32(options.precode_mode), ... % 0 oder 1
'postFFE', ~isempty(options.postFFE), ... % Beispielwert
'NpostFFE', npostFFE, ... % Beispielwert
'Ne1', eq_.Ne(1), ... % Feedforward Koeffizienten 1. Ordnung
'Ne2', eq_.Ne(2), ... % Feedforward Koeffizienten 2. Ordnung
'Ne3', eq_.Ne(3), ... % Feedforward Koeffizienten 3. Ordnung
'Nb1', eq_.Nb(1), ... % Decision Feedback Koeffizienten 1. Ordnung
'Nb2', eq_.Nb(2), ... % Decision Feedback Koeffizienten 2. Ordnung
'Nb3', eq_.Nb(3), ... % Decision Feedback Koeffizienten 3. Ordnung
'K', eq_.K, ... % Samples pro Symbol
'DCmu', eq_.DCmu, ... % Anpassungsrate für DC-Tap
'ideal_dfe', eq_.ideal_dfe, ... % Flag für ideal DFE (0 oder 1)
'training_length', eq_.training_length, ... % Anzahl Trainingssymbole
'training_loops', eq_.training_loops, ... % Anzahl Trainingsdurchläufe
'TRmu1', eq_.FFEmu, ... % mu für DD-Modus (1. Ordnung)
'TRmu2', eq_.FFEmu, ... % mu für DD-Modus (2. Ordnung)
'TRmu3', eq_.FFEmu, ... % mu für DD-Modus (3. Ordnung)
'TRmuDFE', eq_.DFEmu, ... % mu für DFE-Modus im DD
'dd_loops', 5, ... % Anzahl Durchläufe im DD-Modus
'DDmu1', eq_.DDmu(1), ... % mu für DD-Modus (1. Ordnung)
'DDmu2', eq_.DDmu(2), ... % mu für DD-Modus (2. Ordnung)
'DDmu3', eq_.DDmu(3), ... % mu für DD-Modus (3. Ordnung)
'DDmuDFE', eq_.DDmu(4), ... % mu für DFE-Modus im DD
'comment', 'function: vnle_postfilter_mlse', ... % Zusätzliche Kommentare
'config_hash', NaN ...
);
resultsMLSE = struct( ...
'result_id', NaN, ... %
'run_id', NaN, ... % Beispielhafte Run-ID
'eqParam_id', NaN, ... % Beispielhafter Fremdschlüssel zur EqualizerParameters-Tabelle
'date_of_processing', datetime('now'), ... % Aktuelles Datum und Uhrzeit
'numBits', bits_mlse, ... % Beispiel: 1.000.000 Bits
'numBitErr', errors_mlse, ... % Beispiel: 120 Bitfehler
'BER', ber_mlse, ... % BER = 120 / 1.000.000
'SNR', [], ... % Beispielhafte SNR
'SNR_level', jsonencode([]), ... % SNR-Level als JSON-codiertes Array
'GMI', [], ... % Beispielhafter GMI-Wert
'AIR', [], ... % Beispielhafter AIR-Wert
'EVM', [], ... % Beispielhafte EVM
'EVM_level', jsonencode([]), ... % EVM-Level als JSON-codiertes Array
'Alpha', alpha, ... % Beispielhafter Alpha-Wert
'MLSE_dir', jsonencode([mlse_.DIR])...
);
equalizerConfigMLSE = struct( ...
'eq_id', NaN, ... % Auto-Inkrement, wird in der DB gesetzt
'equalizer_structure', int32(equalizer_structure.vnle_pf_mlse), ... % Beispiel: 1 (z.B. für vnle)
'M', M, ... % Ordnung der PAM-Konstellation
'target_constellation', jsonencode(round(unique(tx_symbols.signal),5)), ... % Beispielhafter Target-String
'db_target', 0, ... % 0 oder 1
'diff_precode', int32(options.precode_mode), ... % 0 oder 1
'postFFE', ~isempty(options.postFFE), ... % Beispielwert
'NpostFFE', npostFFE, ... % Beispielwert
'Ne1', eq_.Ne(1), ... % Feedforward Koeffizienten 1. Ordnung
'Ne2', eq_.Ne(2), ... % Feedforward Koeffizienten 2. Ordnung
'Ne3', eq_.Ne(3), ... % Feedforward Koeffizienten 3. Ordnung
'Nb1', eq_.Nb(1), ... % Decision Feedback Koeffizienten 1. Ordnung
'Nb2', eq_.Nb(2), ... % Decision Feedback Koeffizienten 2. Ordnung
'Nb3', eq_.Nb(3), ... % Decision Feedback Koeffizienten 3. Ordnung
'K', eq_.K, ... % Samples pro Symbol
'DCmu', eq_.DCmu, ... % Anpassungsrate für DC-Tap
'ideal_dfe', eq_.ideal_dfe, ... % Flag für ideal DFE (0 oder 1)
'training_length', eq_.training_length, ... % Anzahl Trainingssymbole
'training_loops', eq_.training_loops, ... % Anzahl Trainingsdurchläufe
'TRmu1', eq_.FFEmu, ... % mu für DD-Modus (1. Ordnung)
'TRmu2', eq_.FFEmu, ... % mu für DD-Modus (2. Ordnung)
'TRmu3', eq_.FFEmu, ... % mu für DD-Modus (3. Ordnung)
'TRmuDFE', eq_.DFEmu, ... % mu für DFE-Modus im DD
'dd_loops', 5, ... % Anzahl Durchläufe im DD-Modus
'DDmu1', eq_.DDmu(1), ... % mu für DD-Modus (1. Ordnung)
'DDmu2', eq_.DDmu(2), ... % mu für DD-Modus (2. Ordnung)
'DDmu3', eq_.DDmu(3), ... % mu für DD-Modus (3. Ordnung)
'DDmuDFE', eq_.DDmu(4), ... % mu für DFE-Modus im DD
'MLSE_mode', 'viterbi', ... % Beispiel: MLSE-Modus als String
'MLSE_trellis_states', jsonencode(mlse_.trellis_states), ... % Trellis-States, z.B. als JSON-String oder kommasepariert
'comment', 'function: vnle_postfilter_mlse', ... % Zusätzliche Kommentare
'config_hash', NaN ...
);
eq_package.resultsVNLE = resultsVNLE;
eq_package.resultsMLSE = resultsMLSE;
eq_package.equalizerConfigVNLE = equalizerConfigVNLE;
eq_package.equalizerConfigMLSE = equalizerConfigMLSE;
% eq_package.vnle_out = eq_signal_sd;
if options.showAnalysis
% fprintf(['VNLE EVM lvl: ',repmat('%.3f ',1,numel(evm_lvl)),' \n'],evm_lvl);
@@ -110,13 +240,13 @@ end
showLevelHistogram(eq_signal_sd,tx_symbols)
% showLevelHistogram(mlse_sig_sd,tx_symbols)
showEQcoefficients(eq_.e,eq_.e2,eq_.e3,"displayname",'Coefficients');
showEQcoefficients('n1',eq_.e,'n2',eq_.e2,'n3',eq_.e3,"displayname",'Coefficients');
showEQNoiseSNR(tx_symbols,eq_noise,"displayname",'vnle snr','fignum',101);
%%% EQ SNR Spectrum %230
%snr
snr_vnle = snr(tx_symbols.signal,eq_noise.signal);
% showErrorBurstCount(eq_signal_sd,tx_symbols)

77
Functions/EQ_visuals/showEQcoefficients.m
View File

@@ -1,58 +1,67 @@
function showEQcoefficients(n1, n2, n3, options)
% Show filter coefficients as stem plot
% n1, n2, and n3 in different subplots
% Scale all y-axis to -1 and 1
function showEQcoefficients(options)
% Show filter coefficients as stem plots.
% Only the provided coefficient arrays (n1, n2, n3) are shown,
% each in its own subplot. The y-axis is scaled to [-1, 1].
arguments
n1
n2
n3
options.fignum (1,1) double = NaN % Default to NaN if not provided
options.displayname (1,:) char = '' % Default to an empty string if not provided
options.n1 = [];
options.n2 = [];
options.n3 = [];
options.fignum (1,1) double = NaN; % Default: create new figure if NaN
options.displayname (1,:) char = ''; % Default: empty string
options.color = [0.2157, 0.4941, 0.7216];
options.clf = 0; % Clear figure before plotting new
options.clf = 0; % Clear figure before plotting if set to 1
end
% Determine the figure number to use or create a new figure
% Determine the figure number to use or create a new one.
if isnan(options.fignum)
fig = figure; % Create a new figure and get its handle
fig = figure;
else
fig = figure(options.fignum); % Use the specified figure number
fig = figure(options.fignum);
end
if options.clf
clf(fig); % Clear the figure if requested
clf(fig);
end
hold on
ax = gca;
N = numel(ax.Children);
% Set up a colormap for consistent coloring
% Set up a colormap for consistent coloring.
cmap = linspecer(8);
options.color = cmap(mod(N, size(cmap, 1)) + 1, :);
% Create subplots for n1, n2, n3
for i = 1:3
subplot(3, 1, i);
switch i
case 1
stem(n1, 'Color', options.color, 'LineWidth', 1,'Marker','.','MarkerSize',10);
title(sprintf('1st order Filter Coefficients: %d',numel(n1)));
case 2
stem(n2, 'Color', options.color, 'LineWidth', 1,'Marker','.','MarkerSize',10);
title(sprintf('2nd order Filter Coefficients: %d',numel(n2)));
case 3
stem(n3, 'Color', options.color, 'LineWidth', 1,'Marker','.','MarkerSize',10);
title(sprintf('3rd order Filter Coefficients: %d',numel(n3)));
end
ylim([-1, 1]); % Scale y-axis to -1 and 1
% Build cell arrays for coefficients and their corresponding titles.
coeffs = {};
titles = {};
if ~isempty(options.n1)
coeffs{end+1} = options.n1;
titles{end+1} = sprintf('1st order Filter Coefficients: %d', numel(options.n1));
end
if ~isempty(options.n2)
coeffs{end+1} = options.n2;
titles{end+1} = sprintf('2nd order Filter Coefficients: %d', numel(options.n2));
end
if ~isempty(options.n3)
coeffs{end+1} = options.n3;
titles{end+1} = sprintf('3rd order Filter Coefficients: %d', numel(options.n3));
end
numSubplots = numel(coeffs);
for i = 1:numSubplots
subplot(1, numSubplots, i);
stem(coeffs{i}, 'Color', options.color, 'LineWidth', 1, ...
'Marker', '.', 'MarkerSize', 10);
title(titles{i});
ylim([-1, 1]); % Set y-axis limits to [-1, 1]
grid on;
grid minor
grid minor;
xlabel('Coefficient Index');
ylabel('Amplitude');
end
% Ensure the layout is tight for better visibility
sgtitle('Filter Coefficients'); % Overall title
end
sgtitle('Filter Coefficients'); % Overall title for the figure
end

36
Functions/EQ_visuals/showEQfilter.m Normal file
View File

@@ -0,0 +1,36 @@
function showEQfilter(coefficients,fs)
% Assuming that obj.e contains the final FFE filter coefficients.
% Set the number of frequency points and sampling frequency.
nfft = 1024; % Number of frequency points
% Compute the frequency response of the FFE filter.
[H, f] = freqz(coefficients, 1, nfft, fs);
% Keep only the first half of the frequency response (up to the Nyquist frequency).
half_nfft = floor(nfft/2) + 1;
f = f(1:half_nfft);
H = H(1:half_nfft);
% Plot the magnitude and phase responses.
figure;
% Magnitude response (in dB)
subplot(2,1,1);
hold on
plot(f.*1e-9, 20*log10(abs(1./H)));
title('(Inverted) Magnitude Response of FFE Filter');
xlabel('Frequency (Hz)');
ylabel('Magnitude (dB)');
grid on;
% Phase response
subplot(2,1,2);
plot(f.*1e-9, unwrap(angle(H)));
title('Phase Response of FFE Filter');
xlabel('Frequency (Hz)');
ylabel('Phase');
grid on;
end

2
Functions/EQ_visuals/showLevelHistogram.m
View File

@@ -39,7 +39,9 @@ end
intermediate = received_sd(lvl,:);
cnt(lvl) = round(numel(intermediate(~isnan(intermediate)))./length(eq_signal),3).*100;
hold on
warning off
histogram(received_sd(lvl,:),1000,"EdgeAlpha",0,'DisplayName',['Lvl ',num2str(lvl),' | ',num2str(cnt(lvl)),' %'],'FaceColor',lvlcol(lvl,:),'Normalization','pdf');
warning on
end
legend
grid on

40
Functions/Metrics/calc_air_plain.m Normal file
View File

@@ -0,0 +1,40 @@
function [airs] = calc_air_plain(noisy_signal,reference_signal,options)
% Calculation of AIR acc. to J. Kozesnik, „Numerically Computing Achievable Rates of Memoryless Channels“, Francisco Javier Garcıa-Gomez, doi: 10.1007/978-94-009-9857-5.
% Implementation is not accessible, I mailed TUM to get the code...
arguments(Input)
noisy_signal;
reference_signal;
options.skip_front = 0;
options.skip_end = 0;
options.returnErrorLocation = 0;
end
options.skip_end = abs(options.skip_end);
options.skip_front = abs(options.skip_front);
assert((options.skip_end+options.skip_front)<length(noisy_signal),"You can not skip more bits than overall length of data! Set skip_front or skip_end to lower value or check data_in");
% TRIM
[noisy_signal,reference_signal]=trimseq(noisy_signal,reference_signal,options.skip_front,options.skip_end);
% CALC EVM
%%% new implementation of AIR
constellation = unique(reference_signal);
reference_idx = arrayfun(@(x) find(constellation == x, 1), reference_signal);
air = air_garcia_implementation(constellation',noisy_signal',reference_idx');
function [data_,reference_]=trimseq(data,reference,skipstart,skip_end)
data_ = data(skipstart+1:end-skip_end,:);
delta_bits = length(reference) - length(data);
skip_end = delta_bits + skip_end;
reference_ = reference(skipstart+1:end-skip_end,:);
end
end

36
Functions/Metrics/calc_snr.m Normal file
View File

@@ -0,0 +1,36 @@
function [snr_all, snr_per_level] = calc_snr(tx_signal, eq_noise)
% CALC_SNR Calculates overall SNR and level-wise SNR for a PAM-M constellation.
%
% [snr_all, snr_per_level] = calc_snr(tx_signal, eq_noise)
%
% Inputs:
% tx_signal - Vector of transmitted signal values.
% eq_noise - Vector of corresponding noise samples.
%
% Outputs:
% snr_all - Overall SNR computed using all signal values.
% snr_per_level - A vector where each element is the SNR computed
% for a unique amplitude level in tx_signal.
%
% The function first computes the overall SNR using the full signal vectors.
% Then it uses the unique levels in tx_signal to calculate the SNR for
% the symbols corresponding to each level separately.
% Calculate overall SNR using the complete signals
snr_all = snr(tx_signal, eq_noise);
% Get the unique amplitude levels in the transmitted signal
levels = unique(tx_signal);
% Preallocate an array to store the SNR for each unique level
snr_per_level = zeros(size(levels));
% Loop over each unique level to compute the SNR for that level
for i = 1:length(levels)
% Find indices where tx_signal equals the current level
idx = (tx_signal == levels(i));
% Compute the SNR for these indices
snr_per_level(i) = snr(tx_signal(idx), eq_noise(idx));
end
end

2
Functions/beautifyBERplot.m
View File

@@ -25,7 +25,7 @@ function beautifyBERplot()
% Set logarithmic scale for y-axis, but only if it makes sense.
% If this is not always desired, you could condition this on the presence of lines or data.
% set(gca, 'YScale', 'log');
set(gca, 'YScale', 'log');
% Customize grid and box appearance
set(gca, 'Box', 'on', 'LineWidth', 0.8); % Thicker border