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

Simulation preps for high speed.

Browse Source
This commit is contained in:
Silas Oettinghaus
2025-11-14 11:45:13 +01:00
parent 0080cb2264
commit 888cbbd23e
9 changed files with 384 additions and 535 deletions
Download Patch File Download Diff File Expand all files Collapse all files

4
Classes/00_signals/Signal.m
View File

@@ -965,8 +965,8 @@ classdef Signal
maxA = max(sig(100:end-100))*1.3; maxA = max(sig(100:end-100))*1.3;
minA = min(sig(100:end-100))*1.3; minA = min(sig(100:end-100))*1.3;
maxA = 0.0025; % maxA = 0.0025;
minA = 0; % minA = 0;
difference= maxA-minA; difference= maxA-minA;

493
Classes/04_DSP/Equalizer/ML_MLSE.m
View File

@@ -1,44 +1,16 @@
classdef ML_MLSE < handle classdef ML_MLSE < handle
% ALGORITHM DESCRIBED IN: % ---------------------------------------------------------------------
% W. Lanneer and Y. Lefevre, “Machine Learning-Based Pre-Equalizers for % W. Lanneer and Y. Lefevre,
% Maximum Likelihood Sequence Estimation in High-Speed PONs,” % “Machine Learning-Based Pre-Equalizers for Maximum Likelihood
% in 2023 31st European Signal Processing Conference % Sequence Estimation in High-Speed PONs,” EUSIPCO 2023
% ---------------------------------------------------------------------
% Further ML Refs: % This implementation reproduces the closed-loop ML-based
% https://machinelearningmastery.com/cross-entropy-for-machine-learning/ % pre-equalizer training for MLSE, supporting both training and
% https://docs.pytorch.org/docs/stable/generated/torch.nn.CrossEntropyLoss.html % detection (decision-directed) modes.
% ---------------------------------------------------------------------
% The central idea is to overcome the (white-) noise assumption within the previously described
% Viterbi algorithm, more precisely a closed-loop optimization is proposed that finds a suitable
% filter-set to directly compute the branch metrics c_k (s,s^' ). These can directly be used to
% carry out the conventional Viterbi algorithm. The system consists of S^L S=F linear FIR filters,
% combined with one bias coefficient respectively. These filters take the received input samples to
% compute the branch metrics estimates (c_k ) ̂(s,s^' ) according toThe central idea is to overcome
% the (white-) noise assumption within the previously described Viterbi algorithm, more precisely
% a closed-loop optimization is proposed that finds a suitable filter-set to directly compute the
% branch metrics c_k (s,s^' ). These can directly be used to carry out the conventional Viterbi
% algorithm. The system consists of S^L S=F linear FIR filters, combined with one bias coefficient
% respectively. These filters take the received input samples to compute the branch metrics
% estimates. Finally, the usual Viterbi is carried out...
% Recommended Settings and some findings:
% Requires many training epochs. According to ML people, 100,200 or
% even up to 1000 epochs are normal for ML-convergence
% The mu parameter _can_ be adaptive - using the cross entropy and when
% analyzing the isolated training it looks very promisig. However, is
% later use I found this is not as stable as a fixed learning rate.
% mu = 0.1 worked good for me
% Longer orders/ filter length are not always better. For me order=11
% was good.
% Delay factor (delta) is good when the order is also increased. With
% order = 11, a delta of =4 shows good results
properties properties
sps % usually 2 sps
order order
e e
e_tr e_tr
@@ -48,27 +20,23 @@ classdef ML_MLSE < handle
mu_tr mu_tr
epochs_tr epochs_tr
dd_mode % 1 or 0 to set DD-mode on or off dd_mode
mu_dd %weight update in dd mode mu_dd
epochs_dd epochs_dd
adaptive_mu adaptive_mu
constellation constellation
L
L %viterbi memory length
alpha alpha
DIR DIR
DIR_flip DIR_flip
trellis_states trellis_states
traceback_depth traceback_depth
delta
% --- Added internal class variables used later --- % Internal variables
S S
Nf Nf
delta
nStates nStates
nFeasible nFeasible
combs combs
@@ -79,38 +47,32 @@ classdef ML_MLSE < handle
valid_from_idx valid_from_idx
w w
% --- New: fast state lookup --- % Fast lookup
nSym
key_table
trans_index
true_to_state_idx true_to_state_idx
state_dict % containers.Map: key(sequence)->state index
key_fmt = '%.8g_'; % key format for sequence strings
nSym % |constellation|
% Debug metrics
ber = [] ber = []
ce = ones(1,1); ce = ones(1,1)
end end
methods methods
function obj = ML_MLSE(options) function obj = ML_MLSE(options)
arguments(Input) arguments(Input)
options.sps = 2; options.sps = 2;
options.order = 15; options.order = 15;
options.len_tr = 4096; options.len_tr = 4096;
options.mu_tr = 0; options.mu_tr = 0.001;
options.epochs_tr = 5; options.epochs_tr = 5;
options.dd_mode = 1; options.dd_mode = 1;
options.mu_dd = 1e-5; options.mu_dd = 1e-5;
options.epochs_dd = 5; options.epochs_dd = 5;
options.adaptive_mu = 1; options.adaptive_mu = 1;
options.delta = 0; options.delta = 0;
options.traceback_depth = 1024; options.traceback_depth = 1024;
options.L = 1;
options.L = 1
end end
fn = fieldnames(options); fn = fieldnames(options);
@@ -122,13 +84,12 @@ classdef ML_MLSE < handle
obj.error = 0; obj.error = 0;
end end
% ==============================================================
% PROCESS
% ==============================================================
function [X,X_viterbi] = process(obj, X, D) function [X,X_viterbi] = process(obj, X, D)
% Normalize input RMS
% actual processing of the signal (steps 1. - 3.)
% 1 normalize RMS
X = X.normalize("mode","rms"); X = X.normalize("mode","rms");
% Use sorted constellation for deterministic mapping
obj.constellation = sort(unique(D.signal),'ascend'); obj.constellation = sort(unique(D.signal),'ascend');
obj.nSym = numel(obj.constellation); obj.nSym = numel(obj.constellation);
@@ -136,22 +97,17 @@ classdef ML_MLSE < handle
warning('Signal length does not fit to reference!'); warning('Signal length does not fit to reference!');
end end
% ==============================================================
% INITIALIZATION (only before final epoch and detection mode)
% ==============================================================
% --- Parameters % --- Parameters
obj.S = numel(obj.constellation); % alphabet size obj.S = obj.nSym;
obj.Nf = obj.order*obj.sps; % filter length obj.Nf = obj.order * obj.sps;
% obj.delta = 3;%ceil(obj.Nf/2); % delay parameter
obj.nStates = obj.S^obj.L; obj.nStates = obj.S^obj.L;
obj.nFeasible = obj.nStates*obj.S; obj.nFeasible = obj.nStates * obj.S;
% --- Trellis mapping % --- Trellis mapping
obj.trellis_states = reshape(obj.constellation,1,[]); obj.trellis_states = reshape(obj.constellation,1,[]);
pre_comb_mat = repmat(obj.trellis_states, obj.L, 1); pre_comb_mat = repmat(obj.trellis_states, obj.L, 1);
pre_comb_cell = mat2cell(pre_comb_mat, ones(1,obj.L), size(pre_comb_mat,2)); pre_comb_cell = mat2cell(pre_comb_mat, ones(1,obj.L), size(pre_comb_mat,2));
obj.combs = fliplr(combvec(pre_comb_cell{:}).'); % rows: states, columns: [x_k, x_{k-1}, ...] obj.combs = fliplr(combvec(pre_comb_cell{:}).');
obj.first_sym = obj.combs(:,1); obj.first_sym = obj.combs(:,1);
obj.last_sym = obj.combs(:,end); obj.last_sym = obj.combs(:,end);
obj.nStates = size(obj.combs,1); obj.nStates = size(obj.combs,1);
@@ -165,328 +121,235 @@ classdef ML_MLSE < handle
end end
end end
end end
[obj.valid_to_idx, obj.valid_from_idx] = find(obj.valid); [obj.valid_to_idx,obj.valid_from_idx] = find(obj.valid);
% --- Allocate vectors and weights % --- Initialize weights
% !! IF SHAPE FIT, then we already have smth there an we want
% to start with the existing fitler-set
if isempty(obj.w) || any(size(obj.w) ~= [obj.Nf+1,obj.nFeasible]) if isempty(obj.w) || any(size(obj.w) ~= [obj.Nf+1,obj.nFeasible])
obj.w = zeros(obj.Nf+1,obj.nFeasible); % filter weights per transition + bias tap
obj.w = randn(obj.Nf+1,obj.nFeasible); obj.w = randn(obj.Nf+1,obj.nFeasible);
end end
% --- Precompute dictionary for fast state lookup (sequence -> state) % --- Fast lookup tables
keys = cell(obj.nStates,1); [~, sym_idx_mat] = ismember(obj.combs, obj.constellation);
for i = 1:obj.nStates key_vals = 1 + sum((sym_idx_mat - 1) .* (obj.nSym .^ (0:obj.L-1)), 2);
keys{i} = obj.seq_key(obj.combs(i,:)); % combs row is already [x_k, x_{k-1}, ...] max_key = obj.nSym^obj.L;
obj.key_table = zeros(max_key,1,'uint32');
obj.key_table(key_vals) = 1:obj.nStates;
obj.trans_index = sparse(obj.nStates,obj.nStates);
for i = 1:length(obj.valid_from_idx)
f = obj.valid_from_idx(i);
t = obj.valid_to_idx(i);
obj.trans_index(t,f) = i;
end end
obj.state_dict = containers.Map(keys, 1:obj.nStates);
% ============================================================== % ==============================================================
% TRAINING % TRAINING
% ============================================================== % ==============================================================
fprintf('\n--- Training mode ---\n');
% Training Mode obj.equalize(X.signal, D.signal, obj.mu_tr, obj.epochs_tr, obj.len_tr, true);
n = obj.len_tr;
training = 1;
obj.equalize(X.signal, D.signal,obj.mu_tr,obj.epochs_tr,n,training);
obj.e_tr = obj.e; obj.e_tr = obj.e;
% ============================================================== % ==============================================================
% DD-Mode / Fixed Mode % DECISION-DIRECTED / TESTING
% ============================================================== % ==============================================================
fprintf('--- Decision-directed / detection mode ---\n');
% Decision Directed Mode [y, y_vit] = obj.equalize(X.signal, D.signal, obj.mu_dd, obj.epochs_dd, X.length, false);
n = X.length;
training = 0;
[y,y_vit]=obj.equalize(X.signal, D.signal,obj.mu_dd,obj.epochs_dd,n,training);
X_viterbi = X; X_viterbi = X;
X.signal = y; X.signal = y;
X.fs = D.fs; %change sampling frequency of outgoing signal from fdac e.g. 2 sps to symbol spaced = fsym
lbdesc = [num2str(obj.order),' tap FFE'];
X = X.logbookentry(lbdesc); % append to logbook
X_viterbi.signal = y_vit; X_viterbi.signal = y_vit;
X_viterbi.fs = D.fs; %change sampling frequency of outgoing signal from fdac e.g. 2 sps to symbol spaced = fsym
lbdesc = [num2str(obj.order),'order FFE + PF + Viterbi'];
X_viterbi = X_viterbi.logbookentry(lbdesc); % append to logbook
end end
% ==============================================================
% EQUALIZE
% ==============================================================
function [y,y_ref] = equalize(obj,x,d,mu,epochs,N,training) function [y,y_ref] = equalize(obj,x,d,mu,epochs,N,training)
% ==============================================================
% FFE + Whitening + ML-Based Branch Metric Estimation + Viterbi
% ==============================================================
debug = 1; debug = 1;
showPlots = 1; showPlots = 1;
% --- Input padding and preallocation
y = zeros(N,1); y = zeros(N,1);
% number of symbol steps in this block
nSymbols = ceil(N/obj.sps); nSymbols = ceil(N/obj.sps);
for epoch = 1:epochs for epoch = 1:epochs
pm = zeros(obj.nStates,1);
% state metrics (log-domain costs): keep as column [nStates×1] pred = zeros(nSymbols,obj.nStates,'uint32');
pm = zeros(obj.nStates,1); % v_{k-1}(s′) pm_sto = nan(obj.nStates,nSymbols,'like',pm);
c_hat = zeros(1,obj.nFeasible);
v_tilde = zeros(1,obj.nFeasible);
pred = zeros(nSymbols, obj.nStates, 'uint32');
pm_sto = nan(obj.nStates, nSymbols,'like',pm);
CE_accum = 0; CE_accum = 0;
start_sample = 1;
end_sample = N;
start_symbol = 1 + floor((start_sample - 1)/obj.sps);
%%% START IDX % --- initialize true state
if training
max_start = length(x) - ( (ceil(N/obj.sps)-1)*obj.sps + 1 );
max_start = max(1, max_start); % safety
start_sample = randi([1, max_start], 1); %rnd training; not really good
start_sample = 1;
end_sample = start_sample + (ceil(N/obj.sps)-1)*obj.sps;
else
start_sample = 1;%obj.len_tr;
end_sample = N;
end
start_symbol = 1 + floor((start_sample - 1)/obj.sps); % ABSOLUTE symbol index
if numel(d) >= obj.L && start_symbol >= obj.L if numel(d) >= obj.L && start_symbol >= obj.L
init_seq = d(start_symbol-obj.L+1 : start_symbol); % [d_k-L+1 ... d_k] init_seq = d(start_symbol-obj.L+1:start_symbol);
true_to_state_idx = obj.state_dict(obj.seq_key(flip(init_seq))); % [d_k ... d_k-L+1] key_init = obj.seq2key(init_seq);
true_to_state_idx = obj.key_table(key_init);
if true_to_state_idx==0, true_to_state_idx=1; end
else else
% Not enough history – fall back to state 1
true_to_state_idx = uint32(1); true_to_state_idx = uint32(1);
end end
symbol = 0;
for sample = start_sample:obj.sps:end_sample for sample = start_sample:obj.sps:end_sample
symbol = symbol + 1; symbol = (sample - start_sample)/obj.sps + 1;
k = symbol;
sym_idx = start_symbol + (symbol - 1); sym_idx = start_symbol + (symbol - 1);
% --- Build Δ-delayed observation window y_k % --- Observation window (with delta)
i1 = sample - obj.Nf + 1 + obj.delta; i1 = sample - obj.Nf + 1 + obj.delta;
i2 = sample + obj.delta; i2 = sample + obj.delta;
buf = x(max(1,i1):min(length(x),i2)); buf = x(max(1,i1):min(length(x),i2));
padL = max(0,1 - i1); padL = max(0,1 - i1);
padR = max(0,i2 - length(x)); padR = max(0,i2 - length(x));
yk = [zeros(padL,1); buf(:); zeros(padR,1)]; % Nf×1 yk = [zeros(padL,1); buf(:); zeros(padR,1)];
yk = [yk;1]; yk = [yk;1];
% --- Predict branch metrics for all feasible transitions: c_hat % --- Branch metrics
c_hat = (yk.' * obj.w); % [1×nFeasible] c_hat = (yk.' * obj.w).';
c_hat = c_hat.'; % [nFeasible×1]
% --- Extended path metrics: v_tilde = pm(from) + c_hat
% normalize pm to avoid growth (invariant to additive const)
pm = pm - min(pm); pm = pm - min(pm);
v_tilde = pm(obj.valid_from_idx) + c_hat; % [nFeasible×1] v_tilde = pm(obj.valid_from_idx) + c_hat;
% ===== Gradient update (Algorithm 1) ===== % --- allocate once
if epoch==1 && symbol==1
obj.true_to_state_idx = ones(ceil(N/obj.sps),1,'uint32');
end
if 1 %training % --- previous "to" becomes "from"
% --- allocate storage once if symbol>1
if epoch == 1 && symbol == 1 true_from_state_idx = obj.true_to_state_idx(symbol-1);
obj.true_to_state_idx = ones(ceil(N/obj.sps),1,'uint32'); else
end true_from_state_idx = 1;
end
% --- previous "to" becomes current "from" % --- compute or reuse "to" state
if symbol > 1 if epoch==1
true_from_state_idx = obj.true_to_state_idx(symbol-1); if sym_idx>=obj.L
else key_to = obj.seq2key(d(sym_idx-obj.L+1:sym_idx));
true_from_state_idx = 1; state_idx = obj.key_table(key_to);
end if state_idx==0
state_idx = true_from_state_idx;
% --- compute or reuse "to" state
if epoch == 1
% only compute in first epoch
if sym_idx >= obj.L
key_to = obj.seq_key(flip(d(sym_idx-obj.L+1 : sym_idx)));
if isKey(obj.state_dict, key_to)
obj.true_to_state_idx(symbol) = obj.state_dict(key_to);
else
obj.true_to_state_idx(symbol) = true_from_state_idx;
end
else
obj.true_to_state_idx(symbol) = true_from_state_idx;
end end
end obj.true_to_state_idx(symbol) = state_idx;
% --- reuse cached state from second epoch onward
true_to_state_idx = obj.true_to_state_idx(symbol);
% --- ensure valid (from,to)
dirac = zeros(obj.nFeasible,1);
mask = obj.valid_from_idx==true_from_state_idx & ...
obj.valid_to_idx ==true_to_state_idx;
if any(mask)
dirac(mask) = 1;
else else
idx = find(obj.valid_from_idx==true_from_state_idx,1,'first'); obj.true_to_state_idx(symbol) = true_from_state_idx;
dirac(idx) = 1;
obj.true_to_state_idx(symbol) = obj.valid_to_idx(idx);
end end
end
true_to_state_idx = obj.true_to_state_idx(symbol);
% --- fast Dirac creation
dirac = zeros(obj.nFeasible,1);
trans_idx = obj.trans_index(true_to_state_idx,true_from_state_idx);
if trans_idx~=0
dirac(trans_idx)=1;
end
% ===================================================================
% TRAINING MODE (weight update)
% softmax over -v_tilde (numerically safe shift) % ===================================================================
v_shift = -(v_tilde - min(v_tilde)); % shift to small positive numbers if training
v_shift = min(v_shift, 100); % clamp exponent argument (≈ exp(50)=3e21) % --- Softmax and CE
v_shift = -(v_tilde - min(v_tilde));
v_shift = min(v_shift,100);
expv = exp(v_shift); expv = exp(v_shift);
p = expv ./ (sum(expv) + eps); p = expv./(sum(expv)+eps);
CE_symbol(symbol) = -log(p(dirac==1)+eps);
% for logging only: % --- CE smoothing and adaptive μ
CE_symbol(symbol) = -log(p(dirac==1) + eps); if sym_idx>obj.L
CE_smooth(symbol)=0.01*CE_symbol(symbol)+0.99*CE_symbol(symbol-1);
if sym_idx > obj.L
CE_smooth(symbol) = 0.01*CE_symbol(symbol) + 0.99*CE_smooth(symbol-1);
else else
if epoch > 1 CE_smooth(symbol)=CE_symbol(symbol);
CE_smooth(symbol) = obj.ce(end); %use ce from last epoch or =1 for very first round?!
else
CE_smooth(symbol) = CE_symbol(symbol);
end
end end
CE_accum=CE_accum+CE_symbol(symbol);
CE_accum = CE_symbol(symbol) + CE_accum; % --- Gradient update
dmp=(dirac-p)';
dL_Dw=(yk).*dmp;
% gradient term (t - p) if sym_idx>=obj.L
dmp = (dirac - p)'; % 1×nFeasible
% Per-feature gradient; implicit expansion gives (Nf+1)×nFeasible
dL_Dw = (yk) .* dmp;
% Start updates only when the ABSOLUTE symbol index has ≥ L history
if sym_idx >= obj.L
if obj.adaptive_mu if obj.adaptive_mu
mu_eff = CE_smooth(sym_idx); mu_eff=CE_smooth(symbol);
mu_eff = max(min(mu_eff, 0.2), 1e-4); mu_eff=max(min(mu_eff,0.2),1e-4);
else else
mu_eff = mu; mu_eff=mu;
end end
obj.w=obj.w - mu_eff.*dL_Dw;
obj.w = obj.w - mu_eff .* dL_Dw; % (Nf+1)×nFeasible
end end
% if debug && epoch > 2
% figure(100);
% subplot(4,1,1);
% heatmap(p');
% title('Probs')
% subplot(4,1,2);
% heatmap(dmp);
% title('Update')
% subplot(4,1,3);
% heatmap(dL_Dw);
% title('Update')
% subplot(4,1,4);
% heatmap(bj.w);
% title('Update')
%
% end
end end
% ===================================================================
% DECODING MODE (Viterbi only)
% --- Compare-Select (matrix form, min of costs) % ===================================================================
v_tilde_mat = inf(obj.nStates, obj.nStates); % Compare-Select (always executed)
v_tilde_mat(obj.valid) = v_tilde; vmat=inf(obj.nStates,obj.nStates);
[pm_next, pred(k,:)] = min(v_tilde_mat, [], 2); vmat(obj.valid)=v_tilde;
[pm_next,pred(symbol,:)]=min(vmat,[],2);
% re-center to keep metrics bounded (decision-invariant) pm_next=pm_next-min(pm_next);
pm_next = pm_next - min(pm_next); pm=pm_next;
pm_sto(:,symbol)=pm;
pm = pm_next;
pm_sto(:,symbol) = pm;
end end
% --- Traceback (full; you can window with traceback_depth if desired) % --- Traceback
[~, s_end] = min(pm); [~,s_end]=min(pm);
viterbi_path = zeros(symbol,1,'uint32'); vpath=zeros(symbol,1,'uint32');
viterbi_path(symbol) = s_end; vpath(symbol)=s_end;
for n = symbol:-1:2 for n=symbol:-1:2
viterbi_path(n-1) = pred(n, viterbi_path(n)); vpath(n-1)=pred(n,vpath(n));
end end
y_ref = d(start_symbol:end); y_ref=d(start_symbol:end);
y = obj.first_sym(viterbi_path); y=obj.first_sym(vpath);
if debug && training
sym_start = start_symbol;
sym_end = start_symbol + symbol - 1;
ref_slice = d(sym_start : sym_end);
err = sum(y ~= ref_slice(1:numel(y)));
% --- BER/CE reporting and plots
if training
err=sum(y~=y_ref(1:length(y)));
ser=err/length(y);
try try
ref_bits = PAMmapper(obj.S,0).demap(ref_slice); ref_bits=PAMmapper(obj.S,0).demap(y_ref(1:length(y)));
eq_bits = PAMmapper(obj.S,0).demap(y); eq_bits=PAMmapper(obj.S,0).demap(y);
[~, ~, ber, ~] = calc_ber(ref_bits, eq_bits, "skip_front", 10, "skip_end", 10, "returnErrorLocation", 1); [~,~,ber,~]=calc_ber(ref_bits,eq_bits,"skip_front",10,"skip_end",10,"returnErrorLocation",1);
fprintf('Epoch: %d - BER: %.1e \n',epoch, ber); fprintf('Epoch %d - BER: %.2e\n',epoch,ber);
obj.ber(epoch) = ber; obj.ber(epoch)=ber;
catch catch
ser = err./length(y); fprintf('Epoch %d - SER: %.2e\n',epoch,ser);
fprintf('Epoch: %d - SER: %.1e \n',epoch, ser); obj.ber(epoch)=ser;
end end
obj.ce(epoch)=CE_accum/symbol;
obj.ce(epoch) = CE_accum./symbol; if debug && mod(epoch,10)==1 && showPlots
if showPlots
figure(10);clf figure(10);clf
subplot(3,2,1:2); subplot(3,2,1:2);
heatmap(obj.w); imagesc(obj.w);axis xy;colorbar;title('Filter W');
title('Filter')
subplot(3,2,3); subplot(3,2,3);
v_tildemat = NaN(obj.nStates, obj.nStates); vtilde_mat=NaN(obj.nStates,obj.nStates);
v_tildemat(obj.valid) = v_tilde; % log-domain scores vtilde_mat(obj.valid)=v_tilde;
heatmap(v_tildemat); imagesc(vtilde_mat);axis xy;colorbar;title('Path Metrics (v\_tilde)');
title('Path Metrics (v_tilde)')
subplot(3,2,4); subplot(3,2,4);
scatter(1:symbol,pm_sto,1,'.') plot(1:symbol,pm_sto);title('Path Metric Evolution');
title('Path Metric Winners') subplot(3,2,5);hold on;
subplot(3,2,5);hold on
scatter(1:symbol,CE_symbol,1,'.'); scatter(1:symbol,CE_symbol,1,'.');
scatter(1:symbol,CE_smooth,1,'.') scatter(1:symbol,CE_smooth,1,'.');
title('Cross Entropy') title('Cross Entropy');
subplot(3,2,6);hold on;
subplot(3,2,6); hold on
% Left y-axis: Cross Entropy (linear)
yyaxis left yyaxis left
scatter(1:length(obj.ce), obj.ce, 10, 's', 'filled') scatter(1:length(obj.ce),obj.ce,10,'s','filled');
ylabel('Cross Entropy') ylabel('Cross Entropy');
% Right y-axis: BER (logarithmic)
yyaxis right yyaxis right
scatter(1:length(obj.ber), obj.ber, 10, 'd', 'filled') scatter(1:length(obj.ber),obj.ber,10,'d','filled');
set(gca, 'YScale', 'log') set(gca,'YScale','log');
ylabel('BER (log scale)') ylabel('BER (log)');
xlabel('Epoch');grid on;
xlim([1, epochs]) title('Convergence');
xlabel('Epoch') drawnow;
title('Cross Entropy // BER')
grid on
drawnow
end end
end end
end end
end end
end
methods (Access=private) % ==============================================================
function k = seq_key(obj, seq) % Helper: Sequence → key (always scalar)
% Build a stable key string for a sequence row vector in the *same order as combs rows* ([x_k, x_{k-1}, ...]) % ==============================================================
% Use rounding via sprintf to avoid floating-point issues. function key = seq2key(obj, seq)
% seq must be a row vector. [~, idx] = ismember(flip(seq), obj.constellation);
k = sprintf(obj.key_fmt, seq); pow = (obj.nSym .^ (0:obj.L-1)).';
key = 1 + sum((idx(:) - 1) .* pow);
end end
end end
end end

11
Functions/EQ_structures/dsp_runid.m
View File

@@ -97,10 +97,10 @@ try
use_ffe = 0; use_ffe = 0;
use_dfe = 0; use_dfe = 0;
use_vnle_mlse = 0; use_vnle_mlse = 1;
use_dbtgt = 0; use_dbtgt = 0;
use_dbenc = 0; use_dbenc = 0;
use_ml_mlse = 1; use_ml_mlse = 0;
addProcessingResultToDatabase = 0; addProcessingResultToDatabase = 0;
@@ -140,7 +140,7 @@ try
% Preprocess signal % Preprocess signal
Scpe_sig = preprocessSignal(Scpe_cell{r}, Symbols, fsym); Scpe_sig = preprocessSignal(Scpe_cell{r}, Symbols, fsym);
% Scpe_sig.spectrum("fignum",2223,"normalizeTo0dB",1,"displayname",'Rx'); Scpe_sig.spectrum("fignum",2223,"normalizeTo0dB",1,"displayname",'Rx');
% Scpe_sig.spectrum("fignum",22233,"normalizeTo0dB",0,"displayname",'Rx'); % Scpe_sig.spectrum("fignum",22233,"normalizeTo0dB",0,"displayname",'Rx');
% Scpe_sig.eye(fsym,M,"fignum",1024); % Scpe_sig.eye(fsym,M,"fignum",1024);
@@ -194,7 +194,8 @@ try
pf_ncoeffs = 1; pf_ncoeffs = 1;
ffe_order = [50, 5, 5]; ffe_order = [50, 5, 5];
eq_ = EQ("Ne",ffe_order,"Nb",dfe_order,"training_length",len_tr,"training_loops",5,"dd_loops",5,"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1); eq_ = EQ("Ne",ffe_order,"Nb",dfe_order,"training_length",len_tr,"training_loops",5,"dd_loops",5,"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",0);
% eq_ = VNLE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",[0.0004 0.0005 0.0006],"mu_tr",0.0004,"order",[50,5,5],"sps",2,"decide",0);
pf_ = Postfilter("ncoeff",pf_ncoeffs,"useBurg",1); pf_ = Postfilter("ncoeff",pf_ncoeffs,"useBurg",1);
useviterbi = 0; useviterbi = 0;
@@ -288,6 +289,8 @@ try
if duob_mode == db_mode.db_encoded if duob_mode == db_mode.db_encoded
mlse_db_enc = MLSE_viterbi("duobinary_output",0,'M',M,'trellis_states',PAMmapper(M,0).levels); mlse_db_enc = MLSE_viterbi("duobinary_output",0,'M',M,'trellis_states',PAMmapper(M,0).levels);
mlse_db_enc = MLSE("DIR", [1,1], "duobinary_output", 0, "M", M, "trellis_states", PAMmapper(M,0).levels);
db_results = duobinary_signaling(eq_db_enc, mlse_db_enc, M, Scpe_sig, Symbols, Tx_bits, "precode_mode",duob_mode, "showAnalysis",0,"postFFE",[]); db_results = duobinary_signaling(eq_db_enc, mlse_db_enc, M, Scpe_sig, Symbols, Tx_bits, "precode_mode",duob_mode, "showAnalysis",0,"postFFE",[]);
output.dbenc_package{r} = db_results; output.dbenc_package{r} = db_results;
if options.append_to_db if options.append_to_db

22
Functions/EQ_structures/duobinary_signaling.m
View File

@@ -35,8 +35,26 @@ if ~isempty(options.postFFE)
[eq_signal, eq_noise] = options.postFFE.process(eq_signal, tx_symbols); [eq_signal, eq_noise] = options.postFFE.process(eq_signal, tx_symbols);
end end
% Process through MLSE if isa(mlse_,'MLSE_viterbi')
[mlse_signal] = mlse_.process(eq_signal); [mlse_signal] = mlse_.process(eq_signal);
else
% Aufpassen mit welcher Sequenz man hier vergleicht für LLR stuff...
% gespeichtere "Symbols" sind schon DB codiert, das wollen wir hier
% nicht! Sondern die precoded aber nicht db-encoded müssen als ref in
% die LLR berechnung gehen!
ref_sym = PAMmapper(M,0).map(tx_bits); %ist klar
ref_sym_dpc = Duobinary().precode(ref_sym); % precoded
% ref_sym_dbenc = Duobinary().encode(ref_sym_dpc); %encoded - das wurde gesendet!
% ref_sym_dec = Duobinary().decode(ref_sym_dbenc); %ref_sym wieder zurück!
mlse_.trellis_states = PAMmapper(M,0).levels;
mlse_.trellis_state_mode = 1;
[mlse_signal,LLR,GMI_MLSE] = mlse_.process(eq_signal,ref_sym_dpc);
end
% tx_symbols_ = Duobinary().decode(tx_symbols); % tx_symbols_ = Duobinary().decode(tx_symbols);
% [mlse_signal,~,GMI_MLSE] = mlse_.process(eq_signal,tx_symbols); % [mlse_signal,~,GMI_MLSE] = mlse_.process(eq_signal,tx_symbols);

55
Functions/EQ_structures/ml_mlse.m
View File

@@ -44,31 +44,35 @@ try
end end
ml_mlse_results.config = Equalizerstruct(); ml_mlse_results.config = Equalizerstruct();
eq_small = strip_eq(eq_, 10);
json_str = jsonencode(eq_small);
ml_mlse_results.config.eq = jsonencode(eq_); ml_mlse_results.config.eq = jsonencode(eq_);
ml_mlse_results.config.equalizer_structure = int32(equalizer_structure.ml_mlse); ml_mlse_results.config.equalizer_structure = int32(equalizer_structure.ml_mlse);
ml_mlse_results.config.comment = 'function: ML-based MLSE'; ml_mlse_results.config.comment = 'function: ML-based MLSE';
ml_mlse_results.metrics = Metricstruct; ml_mlse_results.metrics = Metricstruct;
ml_mlse_results.metrics.result_id = NaN; % ml_mlse_results.metrics.result_id = NaN;
ml_mlse_results.metrics.run_id = NaN; % ml_mlse_results.metrics.run_id = NaN;
ml_mlse_results.metrics.eqParam_id = NaN; % ml_mlse_results.metrics.eqParam_id = NaN;
ml_mlse_results.metrics.date_of_processing = datetime('now'); ml_mlse_results.metrics.date_of_processing = datetime('now');
ml_mlse_results.metrics.BER = ber; ml_mlse_results.metrics.BER = ber;
ml_mlse_results.metrics.numBits = bits; ml_mlse_results.metrics.numBits = bits;
ml_mlse_results.metrics.numBitErr = errors; ml_mlse_results.metrics.numBitErr = errors;
ml_mlse_results.metrics.BER_precoded = ber_precoded; ml_mlse_results.metrics.BER_precoded = ber_precoded;
ml_mlse_results.metrics.numBitErr_precoded = errors_precoded; ml_mlse_results.metrics.numBitErr_precoded = errors_precoded;
ml_mlse_results.metrics.SNR = NaN; % ml_mlse_results.metrics.SNR = NaN;
ml_mlse_results.metrics.SNR_level = NaN; % ml_mlse_results.metrics.SNR_level = NaN;
ml_mlse_results.metrics.STD = NaN; % ml_mlse_results.metrics.STD = NaN;
ml_mlse_results.metrics.STD_level = NaN; % ml_mlse_results.metrics.STD_level = NaN;
ml_mlse_results.metrics.STDrx = NaN; % ml_mlse_results.metrics.STDrx = NaN;
ml_mlse_results.metrics.STDrx_level = NaN; % ml_mlse_results.metrics.STDrx_level = NaN;
ml_mlse_results.metrics.GMI = NaN; % ml_mlse_results.metrics.GMI = NaN;
ml_mlse_results.metrics.AIR = NaN; % ml_mlse_results.metrics.AIR = NaN;
ml_mlse_results.metrics.EVM = NaN; % ml_mlse_results.metrics.EVM = NaN;
ml_mlse_results.metrics.EVM_level = NaN; % ml_mlse_results.metrics.EVM_level = NaN;
ml_mlse_results.metrics.Alpha = NaN; % ml_mlse_results.metrics.Alpha = NaN;
end end
@@ -111,3 +115,26 @@ switch precode_mode
[bits, errors, ber, error_pos] = calc_ber(rx_bits.signal, tx_bits_demapped.signal, "skip_front", 30000, "skip_end", 150, "returnErrorLocation", 1); [bits, errors, ber, error_pos] = calc_ber(rx_bits.signal, tx_bits_demapped.signal, "skip_front", 30000, "skip_end", 150, "returnErrorLocation", 1);
end end
end end
function eq_out = strip_eq(eq_, max_elems)
% strip_eq removes all large fields from the ML_MLSE object
% eq_out = strip_eq(eq_, max_elems)
% max_elems ... maximum number of elements to keep (default = 10)
if nargin < 2
max_elems = 10; % default threshold
end
props = properties(eq_);
for i = 1:numel(props)
val = eq_.(props{i});
if ~isempty(val)
% Count total number of elements
if numel(val) > max_elems
eq_.(props{i}) = [];
end
end
end
eq_out = eq_;
end

23
Functions/EQ_structures/vnle_postfilter_mlse.m
View File

@@ -45,7 +45,26 @@ end
eq_signal_hd = PAMmapper(M, 0).quantize(eq_signal_sd); eq_signal_hd = PAMmapper(M, 0).quantize(eq_signal_sd);
% Process through postfilter and MLSE % Process through postfilter and MLSE
[mlse_sig_sd,whitened_noise] = pf_.process(eq_signal_sd, eq_noise); [mlse_sig_sd,whitened_noise] = pf_.process(eq_signal_sd, eq_noise);
if 0 %tx_symbols.fs > 190e9
if pf_.ncoeff == 1
if pf_.coefficients(2) < 0
% coeff is negative for too high/ bad VNLE convergence
pf_.coefficients(2) = 0.9;
end
else
%long memory / pf respinse - not sure what to set here in a worst
%case :-)
end
%do it again:
pf_.useBurg = 0;
[mlse_sig_sd,whitened_noise] = pf_.process(eq_signal_sd, eq_noise);
end
mlse_.DIR = pf_.coefficients; mlse_.DIR = pf_.coefficients;
GMI_MLSE = NaN; GMI_MLSE = NaN;
@@ -249,7 +268,8 @@ figure(336);
hold on; hold on;
eq_signal_sd.spectrum("displayname",'Equalized Signal','fignum',336,'normalizeTo0dB',0); eq_signal_sd.spectrum("displayname",'Equalized Signal','fignum',336,'normalizeTo0dB',0);
eq_noise.spectrum("displayname",'Equalized Signal','fignum',336,'normalizeTo0dB',0); eq_noise.spectrum("displayname",'Equalized Signal','fignum',336,'normalizeTo0dB',0);
showEQNoisePSD(eq_noise, "fignum", 336, "displayname", 'Residual Noise after VNLE', 'postfilter_taps', pf_.coefficients);
showEQNoisePSD(eq_noise, "fignum", 338, "displayname", 'Residual Noise after VNLE', 'postfilter_taps', pf_.coefficients);
for t = 1:4 for t = 1:4
pf_.ncoeff = t; pf_.ncoeff = t;
@@ -263,6 +283,7 @@ if ~isempty(postFFE)
showEQcoefficients('n1', postFFE.e, "displayname", 'Coefficients', 'fignum', 338); showEQcoefficients('n1', postFFE.e, "displayname", 'Coefficients', 'fignum', 338);
end end
showEQcoefficients('n1', eq_.e,'n2', eq_.e2,'n3', eq_.e3, "displayname", 'Coefficients', 'fignum', 339);
showEQfilter(eq_.e, eq_signal_sd.fs.*2); showEQfilter(eq_.e, eq_signal_sd.fs.*2);
figure(340); clf; figure(340); clf;

2
Functions/Job_Processing/preprocessSignal.m
View File

@@ -16,7 +16,7 @@ Scpe_sig = Scpe_sig.resample("fs_out", 2*fsym);
[Scpe_sig, ~] = Scpe_sig.tsynch("reference", Symbols, "fs_ref", fsym, "debug_plots", 0); [Scpe_sig, ~] = Scpe_sig.tsynch("reference", Symbols, "fs_ref", fsym, "debug_plots", 0);
% Apply Gaussian filter % Apply Gaussian filter
Scpe_sig = Filter('filtdegree', 4, "f_cutoff", Symbols.fs.*0.6, ... Scpe_sig = Filter('filtdegree', 8, "f_cutoff", Symbols.fs.*0.52, ...
"fs", Scpe_sig.fs, "filterType", filtertypes.gaussian, ... "fs", Scpe_sig.fs, "filterType", filtertypes.gaussian, ...
"active", true).process(Scpe_sig); "active", true).process(Scpe_sig);

287
projects/HighSpeedExperiment_2024/Auswertung_JLT/run_dsp_from_db.m
View File

@@ -1,6 +1,6 @@
% === SETTINGS === % === SETTINGS ===
dsp_options.append_to_db = 1; dsp_options.append_to_db = 0;
dsp_options.max_occurences = 2; dsp_options.max_occurences = 15;
experiment = "highspeed_2024"; experiment = "highspeed_2024";
dsp_options.mode = "load_run_id"; % 'simulate' & 'load_files' dsp_options.mode = "load_run_id"; % 'simulate' & 'load_files'
@@ -39,19 +39,19 @@ end
% === Get Run ID's === % === Get Run ID's ===
fp = QueryFilter(); fp = QueryFilter();
% fp.where('Runs', 'run_id','EQUALS', 987); % fp.where('Runs', 'run_id','EQUALS', 2776);
M = 4; M = 4;
fp.where('Runs', 'pam_level','EQUALS', M); fp.where('Runs', 'pam_level','EQUALS', M);
fp.where('Runs', 'bitrate','EQUALS', 300e9);%360,390 % fp.where('Runs', 'bitrate','EQUALS', 300e9);%360,390
% fp.where('Runs', 'symbolrate','EQUALS', 195e9); % fp.where('Runs', 'symbolrate','EQUALS', 195e9);
% fp.where('Runs', 'fiber_length','EQUALS', 1); fp.where('Runs', 'fiber_length','EQUALS', 2);
fp.where('Runs', 'is_mpi','EQUALS', 0); % fp.where('Runs', 'is_mpi','EQUALS', 0);
% fp.where('Runs', 'interference_path_length','EQUALS', 1000); % fp.where('Runs', 'interference_path_length','EQUALS', 1000);
% fp.where('Runs', 'loop_id','GREATER_THAN', 11); % fp.where('Runs', 'loop_id','GREATER_THAN', 11);
% fp.where('Runs', 'sir','EQUALS',18); % fp.where('Runs', 'sir','EQUALS',18);
fp.where('Runs', 'wavelength','EQUAL', 1310); fp.where('Runs', 'wavelength','EQUAL', 1310);
fp.where('Runs', 'db_mode','EQUALS', 1); fp.where('Runs', 'db_mode','EQUALS', 1);
% fp.where('Runs', 'rop_attenuation','EQUAL', 0); fp.where('Runs', 'rop_attenuation','EQUAL', 0);
% fp.where('Runs', 'power_pd_in','GREATER_THAN', 7); % fp.where('Runs', 'power_pd_in','GREATER_THAN', 7);
[dataTable,~] = db.queryDB(fp, db.getTableFieldNames('Runs')); [dataTable,~] = db.queryDB(fp, db.getTableFieldNames('Runs'));
@@ -72,203 +72,115 @@ wh.addStorage("mlmlse_package");
[results,wh] = submitJobs(dataTable.run_id(:), dsp_options, "serial", 'wh', wh, 'waitbar', true); [results,wh] = submitJobs(dataTable.run_id(:), dsp_options, "serial", 'wh', wh, 'waitbar', true);
%%
results_db = results(dataTable.db_mode==1);
results_nodb = results(dataTable.db_mode==0);
%BER for ML based MLSE
for i = 1:length(results_db)
ber_m = cellfun(@(c) c.metrics.BER_precoded, results_db{1,i}.mlmlse_package);
[BER_MLMLSE_pre_emph(i), ~] = min(ber_m);
ber_m = cellfun(@(c) c.metrics.BER_precoded, results_nodb{1,i}.mlmlse_package);
[BER_MLMLSE(i), ~] = min(ber_m);
ber_m = cellfun(@(c) c.metrics.BER_precoded, results_nodb{1,i}.mlmlse_package);
[BER_MLMLSE(i), ~] = min(ber_m);
baudrate(i) = dataTable.symbolrate(i);
rop_db(i) = dataTable.power_rop((2*i)-1); %% =========================================================================
% LOAD METADATA
% =========================================================================
[dataTable, ~] = db.queryDB(fp, db.getTableFieldNames('Runs'));
results = results(:).'; % ensure row vector
N = numel(results);
% =========================================================================
% PREALLOCATE METRIC ARRAYS
% =========================================================================
BER_VNLE = nan(1,N);
BER_MLSE = nan(1,N);
BER_DB = nan(1,N);
BER_DB_PREC = nan(1,N);
BER_MLMLSE = nan(1,N);
BER_MLMLSE_PREC = nan(1,N);
% =========================================================================
% EXTRACT METRICS (ONE LOOP, ROBUST)
% =========================================================================
for i = 1:N
r = results{i};
% ---- VNLE (no-DB mode) ----
if isfield(r, 'vnle_package') && ~isempty(r.vnle_package)
pkg = r.vnle_package;
BER_VNLE(i) = min(cellfun(@(c) c.metrics.BER, pkg));
end
% ---- Classical MLSE (DB mode) ----
if isfield(r, 'mlse_package') && ~isempty(r.mlse_package)
pkg = r.mlse_package;
BER_MLSE(i) = min(cellfun(@(c) c.metrics.BER, pkg));
BER_MLSE_PREC(i) = min(cellfun(@(c) c.metrics.BER_precoded, pkg));
end
% ---- DB Target (DB mode) ----
if isfield(r, 'dbtgt_package') && ~isempty(r.dbtgt_package)
pkg = r.dbtgt_package;
BER_DB(i) = min(cellfun(@(c) c.metrics.BER, pkg));
BER_DB_PREC(i) = min(cellfun(@(c) c.metrics.BER_precoded, pkg));
end
% ---- ML-based MLSE (both modes) ----
if isfield(r, 'mlmlse_package') && ~isempty(r.mlmlse_package)
pkg = r.mlmlse_package;
% raw BER
BER_MLMLSE(i) = min(cellfun(@(c) c.metrics.BER, pkg));
% precoded BER
if isfield(pkg{1}.metrics, 'BER_precoded')
BER_MLMLSE_PREC(i) = min(cellfun(@(c) c.metrics.BER_precoded, pkg));
end
end
end end
figure(11);hold on %% =========================================================================
plot(sort(rop_db),sort(BER_MLMLSE_pre_emph)) % METADATA (ALWAYS INDEX-ALIGNED WITH RESULTS)
plot(sort(rop_db),sort(BER_MLMLSE)) % =========================================================================
beautifyBERplot bitrate = dataTable.bitrate(:).';
baudrate = dataTable.symbolrate(:).';
%% rop_atten = dataTable.rop_attenuation(2:2:end).';
results_db = results(dataTable.db_mode==1); rop_pre = dataTable.power_rop(1:2:end).';
results_nodb = results(dataTable.db_mode==0); rop = dataTable.power_rop(2:2:end).';
for i = 1:numel(results_nodb) % =========================================================================
% PLOT STYLE
% =========================================================================
STYLE_BASE = 2;
MARKER_SIZE = STYLE_BASE;
LINE_WIDTH = max(2, STYLE_BASE/3);
% VNLE (from results_nodb) cols = cbrewer2('Paired', 8);
gmi_v = cellfun(@(c) c.metrics.GMI, results_nodb{1,i}.vnle_package);
ber_v = cellfun(@(c) c.metrics.BER, results_nodb{1,i}.vnle_package);
air_v = cellfun(@(c) c.metrics.AIR, results_nodb{1,i}.vnle_package);
snr_v = cellfun(@(c) c.metrics.SNR, results_nodb{1,i}.vnle_package);
[BER_VNLE(i), idx_ber] = min(ber_v);
GMI_VNLE(i) = gmi_v(idx_ber);
AIR_VNLE(i) = air_v(idx_ber);
SNR_VNLE(i) = max(snr_v);
idx_gmi_min_vnle(i) = find(gmi_v == min(gmi_v), 1);
idx_air_max_vnle(i) = find(air_v == max(air_v), 1);
% MLSE (from results_db) cm.VNLE = cols(1,:);
gmi_m = cellfun(@(c) c.metrics.GMI, results_db{1,i}.mlse_package); cm.MLSE = cols(2,:);
ber_m = cellfun(@(c) c.metrics.BER, results_nodb{1,i}.mlse_package); cm.DB_PREC = cols(3,:);
air_m = cellfun(@(c) c.metrics.AIR, results_db{1,i}.mlse_package); cm.DB = cols(4,:);
[BER_MLSE(i), idx_ber] = min(ber_m); cm.ML_MLSE = cols(6,:);
GMI_MLSE(i) = gmi_m(idx_ber);
AIR_MLSE(i) = air_m(idx_ber);
idx_gmi_min_mlse(i) = find(gmi_m == min(gmi_m), 1);
idx_air_max_mlse(i) = find(air_m == max(air_m), 1);
% DB (from results_db, BER_precoded) mk = @(col,shape) {'Marker',shape,'MarkerFaceColor',col,'MarkerEdgeColor',col,'MarkerSize',MARKER_SIZE};
gmi_db = cellfun(@(c) c.metrics.GMI, results_db{1,i}.dbtgt_package);
ber_db = cellfun(@(c) c.metrics.BER, results_db{1,i}.dbtgt_package);
ber_db_prec = cellfun(@(c) c.metrics.BER_precoded, results_db{1,i}.dbtgt_package);
air_db = cellfun(@(c) c.metrics.AIR, results_db{1,i}.dbtgt_package);
[BER_DB(i), idx_ber] = min(ber_db);
[BER_DB_PREC(i), idx_ber] = min(ber_db_prec);
GMI_DB(i) = gmi_db(idx_ber); % =========================================================================
AIR_DB(i) = air_db(idx_ber); % FIGURE 1: BER vs BAUDRATE
idx_gmi_min_db(i) = find(gmi_db == min(gmi_db), 1); % =========================================================================
idx_air_max_db(i) = find(air_db == max(air_db), 1);
%BER for ML based MLSE
ber_m = cellfun(@(c) c.metrics.BER, results_nodb{1,i}.mlmlse_package);
[BER_MLMLSE(i), idx_ber] = min(ber_m);
ber_m = cellfun(@(c) c.metrics.BER_precoded, results_db{1,i}.mlmlse_package);
[BER_MLMLSE_PREC(i), idx_ber] = min(ber_m);
% metadata
bitrate(i) = dataTable.bitrate(i);
baudrate(i) = dataTable.symbolrate(i);
rop_atten(i) = dataTable.rop_attenuation(2*i);
rop_pre(i) = dataTable.power_rop((2*i)-1);
rop(i) = dataTable.power_rop((2*i));
end
%%
STYLE_BASE = 2; % adjust this single number to scale markers & lines
MARKER_SIZE = STYLE_BASE; % marker size (MATLAB MarkerSize)
LINE_WIDTH = max(2, STYLE_BASE/3); % line width (keeps lines reasonable when STYLE_BASE large)
% --- color map / method -> color assignment (keeps colors consistent) ---
cols = cbrewer2('Paired',8);
cols = linspecer(6);
d = 0;
cm.VNLE = cols(1 + d, :);
cm.MLSE = cols(2 + d, :);
cm.DB_precode = cols(3 + d, :);
cm.DB = cols(4 + d, :); % duobinary
cm.ML_MLSE = cols(6 + d, :);
% prepare x values in GBd
xGHz = baudrate .* 1e-9;
xticks_vals = xGHz;
xtick_labels = arrayfun(@(v) sprintf('%d', round(v)), xticks_vals, 'UniformOutput', false);
% common marker settings (filled, same face+edge color)
mk.VNLE = {'Marker','o','MarkerFaceColor',cm.VNLE,'MarkerEdgeColor',cm.VNLE,'MarkerSize',MARKER_SIZE};
mk.MLSE = {'Marker','*','MarkerFaceColor',cm.MLSE,'MarkerEdgeColor',cm.MLSE,'MarkerSize',MARKER_SIZE};
mk.DB_precode = {'Marker','^','MarkerFaceColor',cm.DB_precode,'MarkerEdgeColor',cm.DB_precode,'MarkerSize',MARKER_SIZE};
mk.DB = {'Marker','d','MarkerFaceColor',cm.DB,'MarkerEdgeColor',cm.DB,'MarkerSize',MARKER_SIZE};
mk.ML_MLSE = {'Marker','^','MarkerFaceColor',cm.ML_MLSE,'MarkerEdgeColor',cm.ML_MLSE,'MarkerSize',MARKER_SIZE};
% ---------------- FIGURE : BER ----------------
figure(112+M); clf; hold on; figure(112+M); clf; hold on;
plot(xGHz, BER_VNLE, ... xGHz = baudrate * 1e-9;
'DisplayName','VNLE', ...
mk.VNLE{:}, 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.VNLE); plot(xGHz, BER_VNLE, 'DisplayName','VNLE', 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.VNLE);
plot(xGHz, BER_MLSE, ... plot(xGHz, BER_MLSE, 'DisplayName','MLSE', 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.MLSE);
'DisplayName','MLSE', ... % plot(xGHz, BER_MLSE_PREC, 'DisplayName','MLSE', 'LineStyle',':','LineWidth',LINE_WIDTH,'Color',cm.MLSE);
mk.MLSE{:}, 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.MLSE); plot(xGHz, BER_DB_PREC, 'DisplayName','Diff. Precode + DB', 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.DB_PREC);
plot(xGHz, BER_DB_PREC, ... plot(xGHz, BER_MLMLSE_PREC, 'DisplayName','ML-based MLSE', 'LineStyle',':','LineWidth',LINE_WIDTH,'Color',cm.ML_MLSE);
'DisplayName','Diff. Precode + DB tgt.', ...
mk.DB{:}, 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.DB);
plot(xGHz, BER_MLMLSE_PREC, ...
'DisplayName','ML-based MLSE (L=2)', ...
mk.ML_MLSE{:}, 'LineStyle',':','LineWidth',LINE_WIDTH,'Color',cm.ML_MLSE);
yline(2e-2,'LineWidth',1,'HandleVisibility','off'); yline(2e-2,'LineWidth',1,'HandleVisibility','off');
yline(4.85e-3,'LineWidth',1,'HandleVisibility','off'); yline(4.85e-3,'LineWidth',1,'HandleVisibility','off');
yline(2.2e-4,'LineWidth',1,'HandleVisibility','off'); yline(2.2e-4,'LineWidth',1,'HandleVisibility','off');
xlabel('Baudrate in GBd'); xlabel('Baudrate in GBd');
ylabel('BER'); ylabel('BER');
set(gca, 'yscale', 'log'); set(gca, 'YScale', 'log'); grid on; legend('Location','best');
set(gca, 'XTick', xticks_vals(1:end), 'XTickLabel', xtick_labels(1:end)); % beautifyBERplot;
grid on;
legend('Location','best');
beautifyBERplot
%%
%%
STYLE_BASE = 2; % adjust this single number to scale markers & lines
MARKER_SIZE = STYLE_BASE; % marker size (MATLAB MarkerSize)
LINE_WIDTH = max(2, STYLE_BASE/3); % line width (keeps lines reasonable when STYLE_BASE large)
% --- color map / method -> color assignment (keeps colors consistent) ---
cols = cbrewer2('Paired',8);
cols = linspecer(6);
d = 0;
cm.VNLE = cols(1 + d, :);
cm.MLSE = cols(2 + d, :);
cm.DB_precode = cols(3 + d, :);
cm.DB = cols(4 + d, :); % duobinary
cm.ML_MLSE = cols(6 + d, :);
% prepare x values in GBd
xdbm = flip(sort(rop));
xticks_vals = xdbm;
xtick_labels = arrayfun(@(v) sprintf('%d', round(v)), xticks_vals, 'UniformOutput', false);
% common marker settings (filled, same face+edge color)
mk.VNLE = {'Marker','o','MarkerFaceColor',cm.VNLE,'MarkerEdgeColor',cm.VNLE,'MarkerSize',MARKER_SIZE};
mk.MLSE = {'Marker','*','MarkerFaceColor',cm.MLSE,'MarkerEdgeColor',cm.MLSE,'MarkerSize',MARKER_SIZE};
mk.DB_precode = {'Marker','^','MarkerFaceColor',cm.DB_precode,'MarkerEdgeColor',cm.DB_precode,'MarkerSize',MARKER_SIZE};
mk.DB = {'Marker','d','MarkerFaceColor',cm.DB,'MarkerEdgeColor',cm.DB,'MarkerSize',MARKER_SIZE};
mk.ML_MLSE = {'Marker','^','MarkerFaceColor',cm.ML_MLSE,'MarkerEdgeColor',cm.ML_MLSE,'MarkerSize',MARKER_SIZE};
% ---------------- FIGURE : BER ----------------
figure(112+M); clf; hold on;
plot(flip(sort(rop)), sort(BER_VNLE), ...
'DisplayName','VNLE', ...
mk.VNLE{:}, 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.VNLE);
plot(flip(sort(rop)), sort(BER_MLSE), ...
'DisplayName','MLSE', ...
mk.MLSE{:}, 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.MLSE);
plot(flip(sort(rop)), sort(BER_MLSE), ...
'DisplayName','MLSE', ...
mk.MLSE{:}, 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.MLSE);
plot(flip(sort(rop_pre)), sort(BER_DB_PREC), ...
'DisplayName','Diff. Precode + DB tgt.', ...
mk.DB{:}, 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.DB);
plot(flip(sort(rop_pre)), sort(BER_MLMLSE_PREC), ...
'DisplayName','ML-based MLSE Diff. Prec. (L=2)', ...
mk.ML_MLSE{:}, 'LineStyle','-','LineWidth',LINE_WIDTH,'Color',cm.ML_MLSE);
plot(flip(sort(rop_pre)), sort(BER_MLMLSE), ...
'DisplayName','ML-based MLSE (L=2)', ...
mk.ML_MLSE{:}, 'LineStyle',':','LineWidth',LINE_WIDTH,'Color',cm.DB_precode);
yline(2e-2,'LineWidth',1,'HandleVisibility','off');
yline(4.85e-3,'LineWidth',1,'HandleVisibility','off');
yline(2.2e-4,'LineWidth',1,'HandleVisibility','off');
xlabel('ROP in dBm');
ylabel('BER');
set(gca, 'yscale', 'log');
% set(gca, 'XTick', xticks_vals(1:end), 'XTickLabel', xtick_labels(1:end));
grid on;
legend('Location','best');
beautifyBERplot
%%
@@ -283,7 +195,10 @@ beautifyBERplot
% ---------------- FIGURE 15 : GMI ----------------
%% ---------------- FIGURE 15 : GMI ----------------
figure(113+M); clf; hold on; figure(113+M); clf; hold on;
plot(xGHz, GMI_VNLE, ... plot(xGHz, GMI_VNLE, ...
'DisplayName','VNLE', ... 'DisplayName','VNLE', ...

18
projects/HighSpeedExperiment_2024/Auswertung_JLT/run_plot_measurements.m
View File

@@ -5,16 +5,16 @@ db = DBHandler("dataBase", [dataBase], "type", database_type);
M = 4; M = 4;
fp = QueryFilter(); fp = QueryFilter();
% fp.where('Runs', 'run_id','EQUALS', 987); % fp.where('Runs', 'run_id','EQUALS', 987);
% fp.where('Runs', 'pam_level','EQUALS', M); fp.where('Runs', 'pam_level','EQUALS', M);
% fp.where('Runs', 'symbolrate','EQUALS', 165e9); %150, 165, 180, 195, 210, 225, 240 % fp.where('Runs', 'symbolrate','EQUALS', 165e9); %150, 165, 180, 195, 210, 225, 240
% fp.where('Runs', 'fiber_length','EQUALS', 10); fp.where('Runs', 'fiber_length','EQUALS', 10);
% fp.where('Runs', 'is_mpi','EQUALS', 0); % fp.where('Runs', 'is_mpi','EQUALS', 0);
% fp.where('Runs', 'power_pd_in','GREATER_THAN', 7); % fp.where('Runs', 'power_pd_in','GREATER_THAN', 7);
% fp.where('Runs', 'interference_path_length','EQUALS', 1000); % fp.where('Runs', 'interference_path_length','EQUALS', 1000);
% fp.where('Runs', 'loop_id','GREATER_THAN', 11); % fp.where('Runs', 'loop_id','GREATER_THAN', 11);
% fp.where('Runs', 'sir','EQUALS',18); % fp.where('Runs', 'sir','EQUALS',18);
% fp.where('Runs', 'wavelength','EQUALS', 1310); fp.where('Runs', 'wavelength','EQUALS', 1310);
fp.where('Runs', 'db_mode','EQUALS', 2); % 0 == high preemphasis // 1 == low preemphasis fp.where('Runs', 'db_mode','EQUALS', 1); % 0 == high preemphasis // 1 == low preemphasis
% fp.where('Runs', 'rop_attenuation','EQUALS', 0); % fp.where('Runs', 'rop_attenuation','EQUALS', 0);
fields = db.getTableFieldNames('power_state_info'); fields = db.getTableFieldNames('power_state_info');
@@ -25,15 +25,15 @@ fields = [fields; db.getTableFieldNames('dashboard_ungrouped_new')];
cfg = struct; cfg = struct;
cfg.x_axis = 'symbolrate'; % 'symbolrate' | 'bitrate' | 'wavelength' cfg.x_axis = 'symbolrate'; % 'symbolrate' | 'bitrate' | 'wavelength'
cfg.y_axis = 'BER'; % 'BER' | 'GMI' | 'AIR' | ... cfg.y_axis = 'BER'; % 'BER' | 'GMI' | 'AIR' | ...
cfg.group_by = {'wavelength'}; cfg.group_by = {'wavelength','equalizer_structure','pre_emph'};
cfg.filters = struct('is_mpi',0,'pam_level',M,'equalizer_structure',[equalizer_structure.vnle_db_mlse]);%,equalizer_structure.ffe,equalizer_structure.vnle_pf_mlse,equalizer_structure.vnle_db_mlse,equalizer_structure.dfe]); cfg.filters = struct('is_mpi',0,'pam_level',M,'equalizer_structure',[equalizer_structure.ml_mlse,equalizer_structure.vnle,equalizer_structure.vnle_pf_mlse,equalizer_structure.vnle_db_mlse,equalizer_structure.dfe]);
cfg.y_scale = 'auto'; % auto -> log for BER*, linear otherwise cfg.y_scale = 'auto'; % auto -> log for BER*, linear otherwise
cfg.outlier = 'mad'; % simple, robust; 'none' or 'pctl' also available cfg.outlier = 'mad'; % simple, robust; 'none' or 'pctl' also available
cfg.show_raw = true; cfg.show_raw = true;
cfg.show_spread = 'iqr'; % 'none' or 'iqr' cfg.show_spread = 'none'; % 'none' or 'iqr'
cfg.agg = 'mean'; % or 'median' cfg.agg = 'mean'; % or 'median'
cfg.show_precoded = 0; cfg.show_precoded = 1;
cfg.fec_lines = [2.2e-4 4.85e-3 2e-2]; % optional cfg.fec_lines = [2.2e-4 4.85e-3 2e-2]; % optional
% New styling knobs % New styling knobs
@@ -47,3 +47,5 @@ cfg.plot.legendLocation = 'best';
cfg.plot.fecLineWidth = 2.4; % thicker FEC limits cfg.plot.fecLineWidth = 2.4; % thicker FEC limits
plot_measurements_gpt(dataTable, cfg); plot_measurements_gpt(dataTable, cfg);
% beautifyBERplot()