diff --git a/.gitignore b/.gitignore index 43b6a18..de1fa6c 100644 --- a/.gitignore +++ b/.gitignore @@ -21,3 +21,5 @@ sccprj/ # Matlab code generation folders codegen/ + +.mat diff --git a/Classes/00_signals/Signal.m b/Classes/00_signals/Signal.m index e5828a4..f75ddd4 100644 --- a/Classes/00_signals/Signal.m +++ b/Classes/00_signals/Signal.m @@ -656,7 +656,7 @@ classdef Signal end %% - function [obj,S,isFlipped] = tsynch(obj,options) + function [obj,S,isFlipped,sequenceFound] = tsynch(obj,options) % time sync and cut arguments obj Signal @@ -664,6 +664,14 @@ classdef Signal options.fs_ref = 0; options.debug_plots = 0; end + + + + S = {}; + isFlipped=0; + sequenceFound = 0; + + %normalize the signal a = obj.normalize("mode","oneone").signal; @@ -681,7 +689,11 @@ classdef Signal %estimate start pos of signal maxpeaknum = floor(length(a)/length(b)); - [pks,pkpos] = findpeaks(abs(co./max(co)),'MinPeakDistance',length(b)/2,'MinPeakHeight',0.2,'NPeaks',maxpeaknum,'SortStr','descend'); + [pks,pkpos,w,p] = findpeaks(abs(co./max(co)),'MinPeakDistance',length(b)/2,'MinPeakHeight',0.2,'NPeaks',maxpeaknum,'SortStr','descend'); + + if mean(w) > 10 || mean(p) > 10 + return + end if options.debug_plots figure() @@ -692,9 +704,6 @@ classdef Signal shifts = shifts(shifts>=0); - S = {}; - isFlipped=0; - if numel(shifts) > 0 %Cut occurences of ref signal from signal (only positive shifts) diff --git a/Classes/01_transmit/PAMmapper.m b/Classes/01_transmit/PAMmapper.m index 92d8c53..6ba8558 100644 --- a/Classes/01_transmit/PAMmapper.m +++ b/Classes/01_transmit/PAMmapper.m @@ -24,7 +24,7 @@ classdef PAMmapper obj.M = M; obj.unipolar = unipolar; - + obj.thresholds = obj.get_demodulation_thresholds(); obj.levels = obj.get_levels(); @@ -36,10 +36,10 @@ classdef PAMmapper end function out = map(obj,signal_in) - + if isa(signal_in,'Signal') - signal_in.signal = obj.map_(signal_in.signal); - + signal_in.signal = obj.map_(signal_in.signal); + % signal_in = signal_in.normalize("mode","rms"); lbdesc = ['Map bat stream to PAM ',num2str(obj.M),' symbols']; signal_in = signal_in.logbookentry(lbdesc,obj); @@ -47,19 +47,19 @@ classdef PAMmapper else out = signal_in; end - + end function signal_out = demap(obj,signal_in) issignalclass = 0; if isa(signal_in,'Signal') - signalclass = signal_in; - signal_in = signal_in.signal; - issignalclass = 1; + signalclass = signal_in; + signal_in = signal_in.signal; + issignalclass = 1; end - signal_out = obj.demap_(signal_in); + signal_out = obj.demap_(signal_in); if issignalclass lbdesc = ['Demap PAM ',num2str(obj.M),' symbols to bit stream']; @@ -68,7 +68,7 @@ classdef PAMmapper signal_in.signal = signal_out; signal_out = signal_in; end - + end function pam_sig = map_(obj,bitpattern) @@ -89,7 +89,7 @@ classdef PAMmapper % 4-ASK: if ~obj.eth_style pam_sig=2*bitpattern(:,1)+(bitpattern(:,1)==bitpattern(:,2)); - + if obj.unipolar==0 pam_sig=2*pam_sig-3; end @@ -100,15 +100,65 @@ classdef PAMmapper pam_sig = pam_sig/sqrt(5); case 6 - - m = 1; - if size(bitpattern,2)>size(bitpattern,1) - bitpattern = bitpattern'; %vector aufrecht stellen - end - % LUT based mapping - for k = 1:5:fix(length(bitpattern)/5)*5 - pam_sig(m:m+1,1) = obj.thresholds(bin2dec(int2str(bitpattern(k:k+4)'))+1,:); - m = m+2; + + if ~obj.eth_style + + m = 1; + if size(bitpattern,2)>size(bitpattern,1) + bitpattern = bitpattern'; %vector aufrecht stellen + end + + % LUT based mapping + for k = 1:5:fix(length(bitpattern)/5)*5 + pam_sig(m:m+1,1) = obj.thresholds(bin2dec(int2str(bitpattern(k:k+4)'))+1,:); + m = m+2; + end + + else + + bitsPerSymbol = reshape(bitpattern,5,[]).'; % reorder 5 bits per symbol + normFactor = 1; + %====================32 QAM===================% + %=============================================% + % Coding % + % 01000 01001 |11001 11000 % + % | % + % 01010 01110 01100 |11100 11110 11010 % + % | % + % 01011 01111 01101 |11101 11111 11011 % + % --------------------|------------------- % + % 00011 00111 00101 |10101 10111 10011 % + % | % + % 00010 00110 00100 |10100 10110 10010 % + % | % + % 00000 00001 |10001 10000 % + %=============================================% + % modulate three LSB first in first Quadrant + % first bit inverts real part if 0 + % second bit inverts imaginary part if 0 + + LSB_symbols = normFactor*(... + +(1+1i) *(bitsPerSymbol(:,3)==1 & bitsPerSymbol(:,4)==0 & bitsPerSymbol(:,5)==1)... + +(3+1i) *(bitsPerSymbol(:,3)==1 & bitsPerSymbol(:,4)==1 & bitsPerSymbol(:,5)==1)... + +(5+1i) *(bitsPerSymbol(:,3)==0 & bitsPerSymbol(:,4)==1 & bitsPerSymbol(:,5)==1)... + +(1+3i) *(bitsPerSymbol(:,3)==1 & bitsPerSymbol(:,4)==0 & bitsPerSymbol(:,5)==0)... + +(3+3i) *(bitsPerSymbol(:,3)==1 & bitsPerSymbol(:,4)==1 & bitsPerSymbol(:,5)==0)... + +(5+3i) *(bitsPerSymbol(:,3)==0 & bitsPerSymbol(:,4)==1 & bitsPerSymbol(:,5)==0)... + +(1+5i) *(bitsPerSymbol(:,3)==0 & bitsPerSymbol(:,4)==0 & bitsPerSymbol(:,5)==1)... + +(3+5i) *(bitsPerSymbol(:,3)==0 & bitsPerSymbol(:,4)==0 & bitsPerSymbol(:,5)==0)); + + Re = real(LSB_symbols); + Im = imag(LSB_symbols); + % if first bit== 0 => invert real part + % if second bit== 0 => invert imag part + modData2 = Re.*(bitsPerSymbol(:,1)*2-1) + 1i*Im.*(bitsPerSymbol(:,2)*2-1); + Re = real(modData2(1:end/2)); + Im = imag(modData2(1:end/2)); + pam_sig = zeros(length(Re)*2,1); + pam_sig(1:2:length(Re)*2) = Re; + pam_sig(2:2:length(Im)*2) = Im; + + end pam_sig = pam_sig/sqrt(10); @@ -120,9 +170,9 @@ classdef PAMmapper x1 = bitpattern(:,1); x2 = (bitpattern(:,1)==bitpattern(:,3)); x3 = x2~=bitpattern(:,2); - + pam_sig = 4*x1 + 2*x2 + x3; - + if obj.unipolar==0 pam_sig=2*pam_sig-7; end @@ -133,7 +183,7 @@ classdef PAMmapper end pam_sig = pam_sig/sqrt(21); - + case 16 % 16-ASK: x1 = bitpattern(:,1); @@ -167,7 +217,7 @@ classdef PAMmapper % case 16 % thres = [-10 -8 -6 -4 -2 0 2 4 6 8 10]; % end - + switch obj.M @@ -190,12 +240,12 @@ classdef PAMmapper thres = thres .* 1/sqrt(5); - case 6 %PAM 6 + case 6 %PAM 6 thres = [-3 5;-1 5;-3 -5;-1 -5;-5 3;-5 1;-5 -3;-5 -1;-1 3;-1 1;-1 -3;-1 -1;-3 3;-3 1;-3 -3;-3 -1;3 5;1 5;3 -5;1 -5;5 3;5 1;5 -3;5 -1;1 3;1 1;1 -3;1 -1;3 3;3 1;3 -3;3 -1]; % thres = [-4 -2 0 2 4]; % thres = thres ./ sqrt(10); - + case 8 % 8-ASK if obj.unipolar==0 @@ -203,7 +253,7 @@ classdef PAMmapper elseif obj.unipolar==1 thres=0.5:6.5; end - + thres=thres./sqrt(21); case 16 @@ -234,23 +284,26 @@ classdef PAMmapper end function [data_out] = demap_(obj,data_in) + data_in= data_in'; + if obj.M ~= 6 - + % create output if ~isempty(obj.thresholds) a = squeeze(repmat(real(data_in),[1 1 length(obj.thresholds)])); %Eingangssignal in 3 spalten b = squeeze(repmat(reshape(obj.thresholds(:).',[1 1 length(obj.thresholds)]),[1 length(data_in) 1])); %Threshold in 3 Spalten comp_real = a > b; %check for each symbol/ sampling if it exeeds the obj.thresholdseshold 1, 2 or 3 - + comp_real=repmat(real(data_in),[1 1 length(obj.thresholds)]) > repmat(reshape(obj.thresholds(:).',[1 1 length(obj.thresholds)]),[1 length(data_in) 1]); - + else comp_real=[]; end - + s1=size(comp_real,1); s2=size(comp_real,2); + end switch obj.M @@ -273,34 +326,70 @@ classdef PAMmapper case 6 - data_in = data_in/(sqrt(mean(abs(data_in).^2))); - data_in = data_in*sqrt(10); + if ~obj.eth_style + + data_in = data_in/(sqrt(mean(abs(data_in).^2))); + data_in = data_in*sqrt(10); + + if size(data_in,2) > 1 + data_in = data_in.'; + end + + if length(data_in)/2 ~= round(length(data_in)/2) + data_in = [data_in;0]; + end + + m = 1; + for n = 1:2:length(data_in) + dist = sqrt((data_in(n)-obj.thresholds(:,1)).^2+(data_in(n+1)-obj.thresholds(:,2)).^2); + [~,dd_idx] = min(dist); + % dec_out(n:n+1) = LUT(dd_idx,:); + data_out(m:m+4) = bitget(dd_idx-1,5:-1:1); + m = m+5; + end + + else + + data_in= data_in'; + rxSym = data_in(1:2:end) + 1i*data_in(2:2:end); % 16×1 + + % Decode the sign bits (bits 1 & 2). + rxBit1 = double(real(rxSym) > 0); + rxBit2 = double(imag(rxSym) > 0); + + % Undo the quadrant inversion and normalization. + rxSym_corr = abs(real(rxSym)) + 1i*abs(imag(rxSym)); + normFactor = 1/sqrt(10); + rxSym_unscaled = rxSym_corr / normFactor; + + + cand = [1+1i, 3+1i, 5+1i, 1+3i, 3+3i, 5+3i, 1+5i, 3+5i]; + candBits = [1 0 1; + 1 1 1; + 0 1 1; + 1 0 0; + 1 1 0; + 0 1 0; + 0 0 1; + 0 0 0]; + d = abs(rxSym_unscaled - cand).^2; % 32×8 distances + [~, idx] = min(d, [], 2); + rxLSB = candBits(idx,:); % 32×3 + + + decodedSymbols = [rxBit1, rxBit2, rxLSB]; + data_out = reshape(decodedSymbols.', [], 1).'; - if size(data_in,2) > 1 - data_in = data_in.'; end - - if length(data_in)/2 ~= round(length(data_in)/2) - data_in = [data_in;0]; - end - - m = 1; - for n = 1:2:length(data_in) - dist = sqrt((data_in(n)-obj.thresholds(:,1)).^2+(data_in(n+1)-obj.thresholds(:,2)).^2); - [~,dd_idx] = min(dist); - % dec_out(n:n+1) = LUT(dd_idx,:); - data_out(m:m+4) = bitget(dd_idx-1,5:-1:1); - m = m+5; - end - + case 8 % 8-ASK if ~obj.eth_style - data_out=[comp_real(:,:,4); - comp_real(:,:,1)-comp_real(:,:,3)+comp_real(:,:,5)-comp_real(:,:,7); - 1-comp_real(:,:,2)+comp_real(:,:,6)]; + data_out=[comp_real(:,:,4); + comp_real(:,:,1)-comp_real(:,:,3)+comp_real(:,:,5)-comp_real(:,:,7); + 1-comp_real(:,:,2)+comp_real(:,:,6)]; else data_out = [(data_in>=0); (abs(data_in)>(4/sqrt(21))); (abs(data_in)>=(2/sqrt(21)))&(abs(data_in)<=(6/sqrt(21)))]; end @@ -324,21 +413,21 @@ classdef PAMmapper data_in Signal options.symbol_levels = [] end - + %A) normally return the preproduct of the decision a = squeeze(repmat(real(data_in.signal),[1 1 length(obj.thresholds)])); %Eingangssignal in 3 spalten b = squeeze(repmat(reshape(obj.thresholds(:).',[1 1 length(obj.thresholds)]),[1 length(data_in.signal) 1])); %Threshold in 3 Spalten comp_real = a > b; %check for each symbol/ sampling if it exeeds the obj.thresholdseshold 1, 2 or 3 - + data_out = data_in; data_out.signal = sum(comp_real,2); - + %Option: return the actual level values/ just map onto given %symbol levels if ~isempty(options.symbol_levels) data_out.signal = options.symbol_levels(data_out.signal+1); end - + end function [out] = separate_pamlevels(obj,data_in) @@ -347,7 +436,7 @@ classdef PAMmapper a = squeeze(repmat(real(data_in.signal),[1 1 length(obj.thresholds)])); %Eingangssignal in 3 spalten b = squeeze(repmat(reshape(obj.thresholds(:).',[1 1 length(obj.thresholds)]),[1 length(data_in.signal) 1])); %Threshold in 3 Spalten comp_real = a > b; %check for each symbol/ sampling if it exeeds the obj.thresholdseshold 1, 2 or 3 - comp_real_sum = sum(comp_real,2); + comp_real_sum = sum(comp_real,2); out = NaN(length(data_in),length(obj.thresholds)+1); @@ -358,7 +447,7 @@ classdef PAMmapper end function [Signal_out] = quantize(obj,Signal_in) - + constellation = obj.get_levels(); constellation = constellation ./ rms(constellation); @@ -366,9 +455,9 @@ classdef PAMmapper if isa(Signal_in,'Signal') issignalclass = 1; Sig_class = Signal_in; - Signal_in = Signal_in.signal; + Signal_in = Signal_in.signal; end - + [~,high_dim_sig] = max(size(Signal_in)); [~,high_dim_const] = max(size(constellation)); @@ -381,7 +470,7 @@ classdef PAMmapper Signal_out = constellation(symbol_idx); Signal_out = reshape(Signal_out,size(Signal_in)); - if issignalclass + if issignalclass Sig_class.signal = Signal_out; Signal_out = Sig_class; end @@ -390,6 +479,10 @@ classdef PAMmapper end + function bitmap = showBitMapping(obj) + bitmap = obj.demap([obj.levels ./ obj.scaling]'); + end + end end diff --git a/Classes/01_transmit/Pulseformer.m b/Classes/01_transmit/Pulseformer.m index dc7d0d1..cbf62f3 100644 --- a/Classes/01_transmit/Pulseformer.m +++ b/Classes/01_transmit/Pulseformer.m @@ -6,8 +6,10 @@ classdef Pulseformer end properties(Access=public) - fdac + fs fsym + matched_sps + output_sps pulse pulselength alpha @@ -20,8 +22,10 @@ classdef Pulseformer % Detailed explanation goes here arguments - options.fdac double + options.fs double options.fsym double + options.matched_sps double + options.output_sps double options.pulse pulseform = pulseform.rc options.pulselength double {mustBeInteger} = 32 options.alpha double = 0.05 @@ -50,7 +54,11 @@ classdef Pulseformer signalclass_in = signalclass_in.logbookentry(lbdesc); % write fs to signal - signalclass_in.fs = obj.fdac;%.* (obj.fdac./obj.fsym); + if obj.matched + signalclass_in.fs = obj.fsym .* obj.output_sps;%.* (obj.fdac./obj.fsym); + else + signalclass_in.fs = obj.fs;%.* (obj.fdac./obj.fsym); + end % write to output signalclass_out = signalclass_in; @@ -77,15 +85,19 @@ classdef Pulseformer data_out end - if ~rem(obj.fdac,obj.fsym) + if ~isempty(obj.output_sps) + obj.fsym = obj.fsym.*obj.output_sps; + end + + if ~rem(obj.fs,obj.fsym) %ist ein Vielfaches - sps = obj.fdac / obj.fsym; + sps = obj.fs / obj.fsym; p = sps; q = 1; else %ist kein Vielfaches - p = obj.fsym / gcd(obj.fdac, obj.fsym); %upsampling p->->-> - q = obj.fdac/ gcd(obj.fdac, obj.fsym); %downsampling <-q + p = obj.fsym / gcd(obj.fs, obj.fsym); %upsampling p->->-> + q = obj.fs/ gcd(obj.fs, obj.fsym); %downsampling <-q sps= q; %sps während dem pulse shaping end @@ -105,6 +117,7 @@ classdef Pulseformer end manual_cyclic_convolution = 0; + upsample_filter = 0; upfirdn_convolution = 1; if manual_cyclic_convolution @@ -122,16 +135,28 @@ classdef Pulseformer data_out=ifft( fft(symbolov.') .* repmat( H,size(data_in,1),1 ).' ).'; data_out = circshift(data_out,[0 -(obj.pulselength*sps)]); - if rem(obj.fdac,obj.fsym) + if rem(obj.fs,obj.fsym) data_out = data_out(1:q:end); end end + if upsample_filter + data_out_ = upsample(data_in,p); + + mfOutput = filter(h, 1, data_out_); % Matched filter output + + figure() + hold on + stem(mfOutput(1:1000),'Marker','o','MarkerSize',1,'LineStyle','-','LineWidth',1); + stem(data_out_(1:1000),'Marker','o','MarkerSize',1,'LineStyle','-','LineWidth',1); + + + end + if upfirdn_convolution %Apply Filter using Matlab build in fctn. - data_out_ = upfirdn(data_in,h,p,q); %cut signal, which is longer due to fir filter diff --git a/Classes/04_DSP/Equalizer/FFE.m b/Classes/04_DSP/Equalizer/FFE.m index 2f414ec..54d0b57 100644 --- a/Classes/04_DSP/Equalizer/FFE.m +++ b/Classes/04_DSP/Equalizer/FFE.m @@ -101,6 +101,17 @@ classdef FFE < handle end x = [zeros(floor(obj.order/2),1); x; zeros(obj.order,1)]; + + if training + mask = ones(obj.order,1); + else + mask = zeros(obj.order,1); + mask(900:end) = 1; + mask(ceil(length(obj.e)/2)) = 1; + end + + mask = ones(obj.order,1); + for epoch = 1 : epochs symbol = 0; @@ -110,7 +121,7 @@ classdef FFE < handle U = x(obj.order+sample-1:-1:sample); - y(symbol,1) = obj.e.' * U; % Calculating output of LMS __ * | + y(symbol,1) = (obj.e.*mask).' * U; % Calculating output of LMS __ * | if training d_hat(symbol,1) = d(symbol); @@ -129,7 +140,6 @@ classdef FFE < handle normalizationfactor = (U.' * U); obj.e = obj.e - err(symbol) * U / normalizationfactor; % Weight update rule of NLMS end - obj.error(epoch,symbol) = err(symbol) * err(symbol)'; % Instantaneous square error diff --git a/Classes/04_DSP/Equalizer/Postfilter.m b/Classes/04_DSP/Equalizer/Postfilter.m index 9a3119b..8d09ef1 100644 --- a/Classes/04_DSP/Equalizer/Postfilter.m +++ b/Classes/04_DSP/Equalizer/Postfilter.m @@ -43,12 +43,12 @@ classdef Postfilter < handle if ~isnan(options.useBurg) && options.useBurg - disp('using burg alg') + % disp('using burg alg') obj.coefficients = arburg(noiseclass_in.signal,obj.ncoeff); elseif ~isempty(options.coefficients) - disp('using given taps') + % disp('using given taps') obj.coefficients = options.coefficients; obj.useBurg = 0; diff --git a/Classes/04_DSP/Sequence Detection/MLSE.m b/Classes/04_DSP/Sequence Detection/MLSE.m index 59b9c33..9576524 100644 --- a/Classes/04_DSP/Sequence Detection/MLSE.m +++ b/Classes/04_DSP/Sequence Detection/MLSE.m @@ -70,7 +70,7 @@ classdef MLSE < handle %%%% Separate the equalized signal into the respective levels based on the actually transmitted level constellation = unique(data_ref); decisionLevels = (constellation(1:end-1) + constellation(2:end)) / 2; - tx_bits = PAMmapper(numel(constellation),0).demap(data_ref); + tx_bits = PAMmapper(numel(constellation),0,"eth_style",1).demap(data_ref); % impulse respnse i.e. [0.5, 1.0000] @@ -188,10 +188,50 @@ classdef MLSE < handle end + % Compute soft output PAM4 stream from the metric_sym + soft_output = zeros(length(data_in),1); % Expected symbol value per stage + symbol_prob = zeros(length(data_in), length(states)); % Store full probability distribution (optional) + llp = llp'; + for n = 1:length(data_in) + metrics = llp(n, :); % A posteriori metric for each PAM4 candidate + % For numerical stability, subtract the maximum metric before exponentiating + maxMetric = min(metrics); + expMetrics = exp(metrics - maxMetric); + probs = expMetrics / sum(expMetrics); % Normalize to get probabilities + symbol_prob(n, :) = probs; % (Optional) store distribution for analysis + % Compute the soft output as the expected value of the PAM4 symbols + soft_output(n) = sum(probs .* constellation'); + end + + % Number of symbols and bits per symbol + num_symbols = constellation; + num_bits = 2; % 2 bits per symbol + + bit_mapping = PAMmapper(4,0,"eth_style",1).showBitMapping; % Each row corresponds to the symbol above + + % Initialize LLR storage + llr = zeros(num_bits, length(data_in)); + + % Compute bit-wise LLRs + for bit_idx = 1:num_bits + % Find indices where bit is 0 and where it is 1 + idx_bit_0 = find(bit_mapping(:,bit_idx) == 0); + idx_bit_1 = find(bit_mapping(:,bit_idx) == 1); + + % Sum over log-probabilities (Max-Log approximation: using min instead of sum) + llr(:,bit_idx) = min(llp(:,idx_bit_1), [], 1) - min(llp(:,idx_bit_0), [], 1); + + end + + % Convert LLR values to a hard-decision bit stream + bit_stream = llr < 0; + [~,~,ber_llr,~] = calc_ber(bit_stream',tx_bits,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + fprintf('LLR BER : %.2e \n',ber_llr); + % directly decide based on lowest LLP index [~,llp_based_state_seq]=min(llp); LLP_EST(1:length(data_in)) = constellation(llp_based_state_seq); - rx_bits = PAMmapper(numel(constellation),0).demap(LLP_EST'); + rx_bits = PAMmapper(numel(constellation),0,"eth_style",1).demap(LLP_EST'); [~,~,ber_llp,~] = calc_ber(rx_bits,tx_bits,"skip_front",100,"skip_end",150,"returnErrorLocation",1); fprintf('LLP BER : %.2e \n',ber_llp); % [~,~,ber_llp,~] = calc_ber(circshift(rx_bits,1),tx_bits,"skip_front",100,"skip_end",150,"returnErrorLocation",1); @@ -200,7 +240,7 @@ classdef MLSE < handle % fprintf('LLP BER -1: %.2e \n',ber_llp); %%%% DECIDE based on Viterbi traceback - rx_bits = PAMmapper(numel(constellation),0).demap(VITERBI_ESTIMATION_SYMBOLS'); + rx_bits = PAMmapper(numel(constellation),0,"eth_style",1).demap(VITERBI_ESTIMATION_SYMBOLS'); [~,~,ber_viterbi,~] = calc_ber(rx_bits,tx_bits,"skip_front",100,"skip_end",150,"returnErrorLocation",1); fprintf('Viterbi BER: %.2e \n',ber_viterbi); % [~,~,ber_viterbi,~] = calc_ber(circshift(rx_bits,1),tx_bits,"skip_front",100,"skip_end",150,"returnErrorLocation",1); @@ -209,14 +249,14 @@ classdef MLSE < handle % directly decide based on the FW path metrics [~,fw_direct_state_seq]=min(pm_survivor_fw); FW_EST(1:length(data_in)) = constellation(fw_direct_state_seq); - rx_bits = PAMmapper(numel(constellation),0).demap(FW_EST'); + rx_bits = PAMmapper(numel(constellation),0,"eth_style",1).demap(FW_EST'); [~,~,ber_fw,~] = calc_ber(rx_bits,tx_bits,"skip_front",100,"skip_end",150,"returnErrorLocation",1); fprintf('FW BER: %.2e \n',ber_fw); % directly decide based on the BW path metrics [~,bw_direct_state_seq]=min(pm_survivor_bw); BW_EST(1:length(data_in)) = constellation(bw_direct_state_seq(2:end)); - rx_bits = PAMmapper(numel(constellation),0).demap(BW_EST'); + rx_bits = PAMmapper(numel(constellation),0,"eth_style",1).demap(BW_EST'); [~,~,ber_bw,~] = calc_ber(rx_bits,tx_bits,"skip_front",100,"skip_end",150,"returnErrorLocation",1); fprintf('BW BER: %.2e \n',ber_bw); % [~,~,ber_viterbi,~] = calc_ber(circshift(rx_bits,1),tx_bits,"skip_front",100,"skip_end",150,"returnErrorLocation",1); @@ -224,6 +264,8 @@ classdef MLSE < handle % [~,~,ber_viterbi,~] = calc_ber(circshift(rx_bits,-1),tx_bits,"skip_front",100,"skip_end",150,"returnErrorLocation",1); % fprintf('BW BER: %.2e \n',ber_viterbi); + PAMmapper(4,0,"eth_style",1).showBitMapping + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% tx_symbolpos = zeros(numel(constellation),length(data_ref)); diff --git a/Classes/04_DSP/TransmissionPerformance.m b/Classes/04_DSP/TransmissionPerformance.m index 4e1c7ee..1550e40 100644 --- a/Classes/04_DSP/TransmissionPerformance.m +++ b/Classes/04_DSP/TransmissionPerformance.m @@ -70,7 +70,7 @@ classdef TransmissionPerformance CODE_RATE_HDFEC = [1/(1+0.067)]; %Beyond 300 Gbps Short-Reach Links Using TFLN MZMs With 500 mVpp and Linear Equalization BERTHRESHOLDS_HDFEC = 3.8e-3; - CODE_RATE_O_FEC = [1/(1+0.15)]; %Stefano im Meeting + CODE_RATE_O_FEC = [1/(1+0.153)]; %Stefano im Meeting BERTHRESHOLDS_O_FEC = 2e-2; @@ -162,6 +162,11 @@ classdef TransmissionPerformance netrates.KP4_hamming.NetRate = NaN(1, numMeasurements); netrates.KP4_hamming.CodeRate = NaN(1, numMeasurements); netrates.KP4_hamming.Threshold = NaN(1, numMeasurements); + + netrates.O_FEC.GrossRate = NaN(1, numMeasurements); + netrates.O_FEC.NetRate = NaN(1, numMeasurements); + netrates.O_FEC.CodeRate = NaN(1, numMeasurements); + netrates.O_FEC.Threshold = NaN(1, numMeasurements); end % Process each measurement individually. @@ -215,6 +220,22 @@ classdef TransmissionPerformance netrates.KP4_hamming.Threshold(i) = obj.BERTHRESHOLDS_KP4_AND_INNER(idxBER); end + idxBER = []; + for j = length(obj.BERTHRESHOLDS_O_FEC):-1:1 + if ber(i) <= obj.BERTHRESHOLDS_O_FEC(j) + idxBER = j; + break; + end + end + if ~isempty(idxBER) + codeRate = obj.CODE_RATE_O_FEC(idxBER); + netrates.O_FEC.NetRate(i) = grossRate(i) * codeRate; + netrates.O_FEC.GrossRate(i) = grossRate(i) ; + netrates.O_FEC.CodeRate(i) = codeRate; + netrates.O_FEC.Threshold(i) = obj.BERTHRESHOLDS_O_FEC(idxBER); + end + + end diff --git a/Classes/DataBaseHandler/DBHandler.m b/Classes/DataBaseHandler/DBHandler.m index 1495209..32de7c5 100644 --- a/Classes/DataBaseHandler/DBHandler.m +++ b/Classes/DataBaseHandler/DBHandler.m @@ -45,15 +45,15 @@ classdef DBHandler < handle else error('DB seems to be corrupt') end - + end function obj = refresh(obj) - % Get table names and the first rows of each table to understand the structure - obj.getTableNames(); - obj.getTables(); - obj.getDistinctValues(); + % Get table names and the first rows of each table to understand the structure + obj.getTableNames(); + obj.getTables(); + obj.getDistinctValues(); end function obj = getTableNames(obj) @@ -187,10 +187,11 @@ classdef DBHandler < handle if isstruct(newRow) fields = fieldnames(newRow); emptyFields = structfun(@isempty,newRow); - if sum(emptyFields)>0 - newRow.(fields{emptyFields==1}) = NaN; - disp(['In Table: ',tableName,': ',fields{emptyFields==1},' was empty, is now NaN ',newRow.(fields{emptyFields==1})]) - + if sum(emptyFields) > 0 + emptyFieldNames = fields(emptyFields); % use () to get a cell array + for idx = 1:numel(emptyFieldNames) + newRow.(emptyFieldNames{idx}) = NaN; + end end newRow = struct2table(newRow); end @@ -274,87 +275,168 @@ classdef DBHandler < handle end end - function addBEREntry(obj, berValue, occurrence, runID, ffe, dfe, mlse, pf, eqType, ffe_order, dfe_order, len_tr, mu_ffe, mu_dfe, mu_dc, comment) - % addBEREntry Adds a BER entry linked to an existing or new Equalizer entry. - % Usage: - % addBEREntry(runID, eq, ffe, dfe, mlse, pf, eqType, ffe_order, dfe_order, len_tr, mu_ffe, mu_dfe, mu_dc, berValue, comment) - if isempty(pf) - postfilter_taps = []; + function resultID = addProcessingResult(obj, run_id, resultData, eqParamsData) + % addProcessingResult Adds a processing result and links it to an EqualizerParameters entry. + % + % Inputs: + % run_id: A run_id from the main table to connect the BER with. + % + % resultData: A struct with fields corresponding to the ProcessingResults table. + % + % eqParamsData: A struct with fields corresponding to the EqualizerParameters table, + % except 'eq_id' and 'config_hash'. These fields are used to compute + % a hash and check for an existing configuration. + % + % Output: + % resultID: The result_id of the newly inserted ProcessingResults entry. + + % 1. Compute hash for equalizer parameters + jsonStr = jsonencode(eqParamsData); + md = java.security.MessageDigest.getInstance('MD5'); + md.update(uint8(jsonStr)); + hashBytes = typecast(md.digest, 'uint8'); + hashStr = lower(dec2hex(hashBytes)'); + hashStr = lower(strtrim(hashStr(:)')); % Convert to a lowercase string + + % Add hash to equalizer parameters + eqParamsData.config_hash = hashStr; + + % 2. Check if an equalizer configuration with the same hash exists + queryStr = sprintf('SELECT eq_id FROM EqualizerParameters WHERE config_hash = ''%s''', eqParamsData.config_hash); + existingEntry = obj.fetch(queryStr); + + if ~isempty(existingEntry) + % Use existing eq_id + eq_id = existingEntry{1,1}; else - postfilter_taps = pf.burg_coeff; + % Insert the new equalizer configuration and get its eq_id + eq_id = obj.appendToTable('EqualizerParameters', eqParamsData); end - % Create equalizer data struct for searching and adding if necessary - equalizerData = struct( ... - 'ffe', jsonencode(ffe), ... - 'dfe', jsonencode(dfe), ... - 'mlse', jsonencode(mlse), ... - 'pf', jsonencode(pf), ... - 'eq_type', string(eqType), ... - 'ffe_order', jsonencode(ffe_order), ... - 'dfe_order', jsonencode(dfe_order), ... - 'postfilter_taps',jsonencode(postfilter_taps),... - 'len_tr', len_tr, ... - 'mu_ffe', jsonencode(mu_ffe), ... - 'mu_dfe', mu_dfe, ... - 'mu_dc', mu_dc, ... - 'comment', comment ... - ); + % 3. Add the equalizer configuration reference and run_id to resultData + resultData.eqParam_id = eq_id; + resultData.run_id = run_id; - % Check if exact Equalizer and BER entries already exist in the DB ... - selectedFields = {'Runs.run_id','BERs.ber_id','Equalizer.eq_id','BERs.ber',['BERs.occurrence' ... - '']}; - filterParams = obj.tables; - filterParams.Equalizer = equalizerData; - [dataTable,sql_query] = obj.queryDB(filterParams, selectedFields); + % 4. Compute hash for the processing result + tempResultData = rmfield(resultData, 'date_of_processing'); + resultJsonStr = jsonencode(tempResultData); + md2 = java.security.MessageDigest.getInstance('MD5'); % Create a new MD5 instance + md2.update(uint8(resultJsonStr)); + resultHashBytes = typecast(md2.digest, 'uint8'); + resultHashStr = lower(dec2hex(resultHashBytes)'); + resultHashStr = lower(strtrim(resultHashStr(:)')); % Convert to a lowercase string - % get or insert Equalizer - if ~isempty(dataTable) - % Equalizer entry already exists, use the existing eq_id - cur_eq_id = dataTable.eq_id; - else - % Insert the new Equalizer entry - cur_eq_id = obj.appendToTable('Equalizer', equalizerData); + % Add the result hash to resultData + resultData.result_hash = resultHashStr; + + % 5. Check if an identical processing result already exists + queryStr2 = sprintf('SELECT result_id FROM Results WHERE result_hash = ''%s''', resultData.result_hash); + existingResult = obj.fetch(queryStr2); + + if ~isempty(existingResult) + % If the result exists, return its result_id without inserting a new row + resultID = existingResult{1,1}; + warning(['Result already exists: ResultID: ',num2str(resultID),'| EQ ID: ',num2str(eq_id),' Run ID: ' num2str(run_id)]) + return; end - % skip if already here or insert BER entry - if ~isempty(dataTable) + % 6. Insert the processing result + resultID = obj.appendToTable('Results', resultData); + end - % A BER entry with the same eq_id, run_id, and occurrence already exists - existingBERValue = dataTable.ber; + function recalcHashes(obj) + % recalcHashes Recalculate hashes for all rows in the Results and EqualizerParameters tables. + % + % For EqualizerParameters, the hash is computed from all fields except + % 'eq_id' and 'config_hash'. + % + % For Results, the hash is computed from all fields except 'result_id', + % 'result_hash', and 'date_of_processing'. - % Compare the existing BER value with the new BER value - if existingBERValue == berValue - fprintf('The BER entry %.2e || -- eq_id: %d -- run_id: %d -- occurrence: %d already exists. \n',berValue, cur_eq_id, runID, occurrence); - else - fprintf('Already found BER for EQ: %.2e ~= %.2e || -- eq_id: %d -- run_id: %d -- occurrence: %d already exists.\n', berValue, existingBERValue, cur_eq_id, runID, occurrence); + % Recalculate hashes for EqualizerParameters + eqParamsRows = obj.fetch('SELECT * FROM EqualizerParameters'); + for i = 1:height(eqParamsRows) + rowStruct = table2struct(eqParamsRows(i,:)); % Convert the table row to a struct + + % Remove fields not part of the hash computation + if isfield(rowStruct, 'eq_id') + rowStruct = rmfield(rowStruct, 'eq_id'); + end + if isfield(rowStruct, 'config_hash') + rowStruct = rmfield(rowStruct, 'config_hash'); end - else - - % No such BER entry exists, insert the new BER entry - berData = struct( ... - 'run_id', runID, ... - 'eq_id', cur_eq_id, ... - 'ber', berValue, ... - 'occurrence', occurrence ... - ); - - obj.appendToTable('BERs', berData); + % Compute MD5 hash from the JSON representation + jsonStr = jsonencode(rowStruct); + md = java.security.MessageDigest.getInstance('MD5'); + md.update(uint8(jsonStr)); + hashBytes = typecast(md.digest, 'uint8'); + hashStr = lower(dec2hex(hashBytes)'); + hashStr = lower(strtrim(hashStr(:)')); + % Update the config_hash field using the eq_id from the table row + eq_id = eqParamsRows.eq_id(i); + updateQuery = sprintf('UPDATE EqualizerParameters SET config_hash = ''%s'' WHERE eq_id = %d', hashStr, eq_id); + obj.executeSQL(updateQuery); end + + % Fetch all rows from the Results table using queryDB with no filters + Results = obj.tables.Results; + + if isstruct(Results) + newFields = {}; + fNames = fieldnames(Results); + for i = 1:numel(fNames) + newFields{end+1} = ['Results.' fNames{i}]; + end + Results = newFields; + end + [resultsRows, ~] = obj.queryDB(obj.tables, Results); + + for i = 1:height(resultsRows) + rowStruct = table2struct(resultsRows(i, :)); % Convert the row to a struct + + % Remove fields not used in the hash calculation + if isfield(rowStruct, 'result_id') + rowStruct = rmfield(rowStruct, 'result_id'); + end + if isfield(rowStruct, 'result_hash') + rowStruct = rmfield(rowStruct, 'result_hash'); + end + if isfield(rowStruct, 'date_of_processing') + rowStruct = rmfield(rowStruct, 'date_of_processing'); + end + + % Compute MD5 hash from the JSON representation + jsonStr = jsonencode(rowStruct); + md = java.security.MessageDigest.getInstance('MD5'); + md.update(uint8(jsonStr)); + hashBytes = typecast(md.digest, 'uint8'); + hashStr = lower(dec2hex(hashBytes)'); + hashStr = lower(strtrim(hashStr(:)')); + + % Access the result_id from the table row and update the hash + result_id = resultsRows.result_id(i); + updateQuery = sprintf('UPDATE Results SET result_hash = ''%s'' WHERE result_id = %d', ... + hashStr, result_id); + obj.executeSQL(updateQuery); + end + + + + end + function executeSQL(obj, query) % This method executes an SQL statement using MATLAB's execute function. execute(obj.conn, query); end - function answer = fetch(obj,query) answer = fetch(obj.conn,query); end @@ -405,15 +487,46 @@ classdef DBHandler < handle function query = constructSQLQuery(obj, filterParams, selectedFields) % constructSQLQuery Constructs the SQL query based on filter parameters and selected fields. + % + % If selectedFields is provided as a struct, it is converted to a cell array. + % The conversion takes the field names and creates entries like: + % {'selectedFields.fieldName'} for each field. - % Construct the SELECT clause dynamically based on user selection + % Input check for selectedFields: if it's a struct, convert it to a cell array. + if isstruct(selectedFields) + newFields = {}; + tableNames = fieldnames(selectedFields); + for t = 1:numel(tableNames) + tableStruct = selectedFields.(tableNames{t}); + fieldNames = fieldnames(tableStruct); + for f = 1:numel(fieldNames) + if isequal(tableStruct.(fieldNames{f}), 1) + newFields{end+1} = sprintf('%s.%s', tableNames{t}, fieldNames{f}); + end + end + end + selectedFields = newFields; + end + + % Construct the SELECT clause dynamically based on user selection. + % (Assuming that when provided as a cell array, each entry is of the form + % 'TableName.fieldName' or, in our conversion case, 'selectedFields.fieldName'.) selectClause = 'SELECT DISTINCT '; for i = 1:numel(selectedFields) fieldParts = strsplit(selectedFields{i}, '.'); - tableName = fieldParts{1}; - fieldName = fieldParts{2}; + % If the field comes from the struct conversion, its first part is 'selectedFields' + % and the actual field name is in the second part. + if strcmp(fieldParts{1}, 'selectedFields') + tableName = fieldParts{1}; % not used for type checking below + fieldName = fieldParts{2}; + else + tableName = fieldParts{1}; + fieldName = fieldParts{2}; + end - if isnumeric(obj.tables.(tableName).(fieldName)) + % Decide on COALESCE depending on the field type. + % If the table is known in obj.tables and the field is numeric, use 'NaN'. + if isfield(obj.tables, tableName) && isfield(obj.tables.(tableName), fieldName) && isnumeric(obj.tables.(tableName).(fieldName)) selectClause = [selectClause, 'COALESCE(', selectedFields{i}, ', ''NaN'') AS ', fieldName]; else selectClause = [selectClause, 'COALESCE(', selectedFields{i}, ', '''') AS ', fieldName]; @@ -426,235 +539,302 @@ classdef DBHandler < handle end end - % Construct the FROM and WHERE clause - baseQuery = [selectClause, 'FROM Runs ' ... - 'LEFT JOIN Configurations ON Runs.run_id = Configurations.run_id ' ... - 'LEFT JOIN Measurements ON Runs.run_id = Measurements.run_id ' ... - 'WHERE ']; - % 'LEFT JOIN BERs ON Runs.run_id = BERs.run_id ' ... - % 'LEFT JOIN Equalizer ON BERs.eq_id = Equalizer.eq_id ' ... - + % --- Adaptive FROM Clause --- + % Use "Runs" as the main table and add LEFT JOINs for every other table in obj.tables + % (except "sqlite_sequence") that has a run_id field. + mainTable = 'Runs'; + fromClause = ['FROM ', mainTable, ' ']; + tableNamesAll = fieldnames(obj.tables); + for t = 1:numel(tableNamesAll) + tableName = tableNamesAll{t}; + if strcmpi(tableName, mainTable) || strcmpi(tableName, 'sqlite_sequence') + continue; + end - % Loop through each table in filterParams + if isfield(obj.tables.(tableName), 'run_id') + % most tables are directly linked to runs table + fromClause = [fromClause, 'LEFT JOIN ', tableName, ' ON ', mainTable, '.run_id = ', tableName, '.run_id ']; + elseif isfield(obj.tables.(tableName), 'eq_id') + % equalizer is only linked to results table + fromClause = [fromClause, 'LEFT JOIN ', tableName, ' ON ', 'Results', '.eqParam_id = ', tableName, '.eq_id ']; + end + end + + % --- WHERE Clause Construction --- + baseQuery = [selectClause, ' ', fromClause, 'WHERE ']; filterClauses = []; tableNames_ = fieldnames(filterParams); - for t = 1:numel(tableNames_) tableName = tableNames_{t}; tableParams = filterParams.(tableName); - - % Loop through each parameter in the table fieldNames = fieldnames(tableParams); for i = 1:numel(fieldNames) fieldName = fieldNames{i}; value = tableParams.(fieldName); - - % Construct the full column name in the format "tableName.fieldName" fullName = sprintf('%s.%s', tableName, fieldName); - % Handle different types of values for SQL query construction + % Handle various types of values for SQL query construction if isempty(value) - % Skip this parameter if it is empty (include all values) continue; elseif isnumeric(value) && isnan(value) - % If value is NaN, use IS NULL in SQL filterClause = sprintf('%s IS NULL', fullName); elseif isnumeric(value) && ~isEnumeration(value) filterClause = sprintf('%s = %f', fullName, value); elseif islogical(value) || (isnumeric(value) && ismember(value, [0, 1])) && ~isEnumeration(value) filterClause = sprintf('%s = %d', fullName, value); elseif ischar(value) || isstring(value) - filterClause = sprintf('%s = ''%s''', fullName, char(value)); %nicht nach string suchen sondern nach chararray -> 'bla' statt "bla" + filterClause = sprintf('%s = ''%s''', fullName, char(value)); elseif isEnumeration(value) filterClause = sprintf('%s = ''%s''', fullName, value); else error('Unsupported data type for field "%s".', fullName); end - - % Add the constructed filter clause to the list filterClauses = [filterClauses, filterClause, ' AND ']; end end - % Remove trailing ' AND ' from the filter clauses if any filters were added + % Remove trailing ' AND ' if any filters were added. if ~isempty(filterClauses) filterClauses = filterClauses(1:end-5); - end - - % Construct the final SQL query - if isempty(filterClauses) - query = [selectClause, 'FROM Runs ' ... - 'LEFT JOIN Configurations ON Runs.run_id = Configurations.run_id ' ... - 'LEFT JOIN Measurements ON Runs.run_id = Measurements.run_id ' ... - 'LEFT JOIN BERs ON Runs.run_id = BERs.run_id']; + query = [selectClause, ' ', fromClause, 'WHERE ', filterClauses]; else - query = [baseQuery, filterClauses]; + query = [selectClause, ' ', fromClause]; end - end - - function selectedFields = promptSelectFields(obj) - % promptSelectFields Prompts the user to select fields from multiple tables to include in the SELECT statement using settingsdlg. + % promptSelectFields Prompts the user to select fields from multiple tables + % using a custom checkbox GUI with scrolling. + % + % The function builds a list of all fields (formatted as 'TableName.fieldName') + % and displays each as a checkbox inside an inner container panel. The container's + % height is set to accommodate all checkboxes, so the scrollable panel shows scrollbars. + % When the user clicks the "Select" button, the selected fields are returned. + % If none are selected, all fields are returned. - % Get all possible fields from all tables (excluding sqlite_sequence) + % Get all possible tables (excluding sqlite_sequence) tableNames = fieldnames(obj.tables); - tableNames = setdiff(tableNames, {'sqlite_sequence'}); % Remove sqlite_sequence - - % Prepare the inputs for settingsdlg - promptSettings = {}; - allFieldsFullName = {}; - convertedFieldNames = {}; + tableNames = setdiff(tableNames, {'sqlite_sequence'}); + % Build a single cell array of all field names with table prefix. + allFields = {}; for i = 1:numel(tableNames) - tableFields = fieldnames(obj.tables.(tableNames{i})); - for j = 1:numel(tableFields) - fieldName = tableFields{j}; - fullName = sprintf('%s.%s', tableNames{i}, fieldName); - convertedName = strrep(fullName, '.', '_'); % Replace '.' with '_' - - allFieldsFullName{end + 1} = fullName; % Add full name to the list - convertedFieldNames{end + 1} = convertedName; % Store the converted name - - % Add the field name and checkbox setting to the prompt - promptSettings{end + 1} = {sprintf('Include %s', fullName), convertedName}; - promptSettings{end + 1} = false; % Default: not selected + fields = fieldnames(obj.tables.(tableNames{i})); + for j = 1:numel(fields) + allFields{end+1} = sprintf('%s.%s', tableNames{i}, fields{j}); end end + numFields = numel(allFields); - % Create the settings dialog - [settings, button] = settingsdlg(... - 'title', 'Select Fields for the SQL Query', ... - 'description', 'Check the boxes for the fields you want to include in the SELECT statement.', ... - promptSettings{:} ... - ); + % Create the main figure. + fig = uifigure('Name', 'Select Fields', 'Position', [100, 100, 400, 600]); - % If the user cancels, default to selecting all fields - if strcmp(button, 'cancel') - selectedFields = allFieldsFullName; - return; + % Create a scrollable panel. + scrollPanel = uipanel(fig, 'Position', [10, 60, 380, 530], 'Scrollable', 'on'); + + % Define checkbox dimensions. + checkboxHeight = 30; + spacing = 5; + totalHeight = numFields * (checkboxHeight + spacing); + + % Create an inner container panel with height larger than the scrollPanel's height. + container = uipanel(scrollPanel, 'Position', [0, 0, scrollPanel.Position(3), totalHeight]); + + % Create checkboxes using absolute positioning in the container. + checkboxes = gobjects(numFields, 1); + for i = 1:numFields + % Calculate the vertical position. + % The origin (0,0) is at the bottom left of the container. + yPos = totalHeight - i*(checkboxHeight + spacing) + spacing; + checkboxes(i) = uicheckbox(container, ... + 'Text', allFields{i}, ... + 'Value', false, ... + 'Position', [10, yPos, container.Position(3)-20, checkboxHeight]); end - % Parse user input into selectedFields + % Create a "Select" button in the main figure. + btn = uibutton(fig, 'Text', 'Select', ... + 'Position', [150, 10, 100, 30], ... + 'ButtonPushedFcn', @(btn, event) uiresume(fig)); + + % Wait for the user to click the button. + uiwait(fig); + + % Retrieve the selected fields. selectedFields = {}; - for i = 1:numel(allFieldsFullName) - convertedName = convertedFieldNames{i}; - if isfield(settings, convertedName) && settings.(convertedName) % Add to selectedFields if the checkbox was selected - selectedFields{end + 1} = allFieldsFullName{i}; %#ok + for i = 1:numFields + if checkboxes(i).Value + selectedFields{end+1} = checkboxes(i).Text; end end - % If no fields are selected, default to selecting all fields + % If no fields are selected, default to all fields. if isempty(selectedFields) - selectedFields = allFieldsFullName; + selectedFields = allFields; end + + % Close the figure. + delete(fig); end + + function filterParams = promptFilterParameters(obj) - % promptFilterParameters Prompts the user to enter filter parameters using the settingsdlg framework. + % promptFilterParameters Prompts the user to enter filter parameters using a + % custom scrollable UI with dropdowns. + % + % For each table (excluding 'sqlite_sequence'), each field that has distinct + % values is displayed as a label and a dropdown. The dropdown items are built + % from the distinct values (with "All" prepended). The output is a struct where, + % for each table, each field is set to the chosen value (or [] if "All" is selected). - % Get all possible parameters from all tables (excluding sqlite_sequence) - tableNames_ = fieldnames(obj.tables); - tableNames_ = setdiff(tableNames_, {'sqlite_sequence'}); % Remove sqlite_sequence + % Get all tables except 'sqlite_sequence' + tableNames = fieldnames(obj.tables); + tableNames = setdiff(tableNames, {'sqlite_sequence'}); - % Prepare the inputs for settingsdlg with sections and separators - promptSettings = {}; - allFieldsFullName = {}; - convertedFieldNames = {}; + % Precompute layout constants. + heightPerTableLabel = 30; + heightPerField = 40; % vertical space for a field (label + dropdown) + spacing = 5; - for i = 1:numel(tableNames_) - % Add a separator for each table section - promptSettings{end + 1} = 'separator'; - promptSettings{end + 1} = tableNames_{i}; - - % Get all fields from the current table - tableFields = fieldnames(obj.tables.(tableNames_{i})); - - % Prepare each field to be added to the dialog + % Compute total required height. + totalHeight = 0; + for i = 1:numel(tableNames) + totalHeight = totalHeight + heightPerTableLabel; + tableName = tableNames{i}; + tableFields = fieldnames(obj.tables.(tableName)); for j = 1:numel(tableFields) fieldName = tableFields{j}; - fullName = sprintf('%s.%s', tableNames_{i}, fieldName); - convertedName = strrep(fullName, '.', '_'); % Replace '.' with '_' - - % Skip fields that do not have distinct values stored - if ~isfield(obj.distinctValues.(tableNames_{i}), fieldName) - continue; + % Only include fields that have distinct values stored. + if isfield(obj.distinctValues.(tableName), fieldName) + totalHeight = totalHeight + heightPerField; end - - % Get the distinct values for the field - distinctValues_ = obj.distinctValues.(tableNames_{i}).(fieldName); - - % Prepare distinct values for dropdown - if isempty(distinctValues_) - % If there are no distinct values, use only an "All" entry - distinctValues_ = {'All'}; - else - % Ensure distinctValues is a cell array of strings - if isnumeric(distinctValues_) - distinctValues_ = arrayfun(@(x) num2str(x), distinctValues_, 'UniformOutput', false); - elseif isstring(distinctValues_) - distinctValues_ = cellstr(distinctValues_); - elseif iscell(distinctValues_) && ~iscellstr(distinctValues_) - distinctValues_ = cellfun(@num2str, distinctValues_, 'UniformOutput', false); - end - - % Add an "All" option at the beginning of the distinct values list - distinctValues_ = [{'All'}; distinctValues_]; - end - - allFieldsFullName{end + 1} = fullName; % Add full name to the list - convertedFieldNames{end + 1} = convertedName; % Store the converted name - - % Add the field name and value setting to the prompt - promptSettings{end + 1} = {sprintf('%s', fullName), convertedName}; - promptSettings{end + 1} = distinctValues_; % Add distinct values as dropdown options end end - % Create the settings dialog - [settings, button] = settingsdlg(... - 'title', 'Input Parameters for Filtering', ... - 'description', 'Enter the values for each field to filter. Select "All" to include all values.', ... - promptSettings{:} ... - ); + % Create the main UI figure. + fig = uifigure('Name', 'Input Parameters for Filtering', 'Position', [100, 100, 500, 600]); - % If the user cancels, return an empty struct - if strcmp(button, 'cancel') + % Set a CloseRequestFcn so that closing the figure calls uiresume. + fig.CloseRequestFcn = @(src, event) uiresume(src); + + % Create a scrollable panel inside the figure. + scrollPanel = uipanel(fig, 'Position', [10, 60, 480, 530], 'Scrollable', 'on'); + + % Create an inner container panel with a height set to totalHeight. + container = uipanel(scrollPanel, 'Position', [0, 0, scrollPanel.Position(3), totalHeight]); + + % Prepare cell arrays to store dropdown handles and corresponding table/field names. + dropdownHandles = {}; + dropdownTableNames = {}; + dropdownFieldNames = {}; + + % Set the starting Y coordinate (filling from top to bottom). + currentY = totalHeight; + + % Maximum number of dropdown items. + maxItems = 100; + + for i = 1:numel(tableNames) + % Create a label for the table name. + uilabel(container, ... + 'Text', tableNames{i}, ... + 'FontWeight', 'bold', ... + 'Position', [10, currentY - heightPerTableLabel + spacing, 200, heightPerTableLabel - spacing]); + currentY = currentY - heightPerTableLabel; + + tableName = tableNames{i}; + tableFields = fieldnames(obj.tables.(tableName)); + for j = 1:numel(tableFields) + fieldName = tableFields{j}; + if ~isfield(obj.distinctValues.(tableName), fieldName) + continue; % Skip if no distinct values are stored. + end + + % Retrieve distinct values for the field. + distinctValues_ = obj.distinctValues.(tableName).(fieldName); + if ~isempty(distinctValues_) && numel(distinctValues_) > maxItems + distinctValues_ = distinctValues_(1:maxItems); + end + + if isempty(distinctValues_) + items = {'All'}; + else + if isnumeric(distinctValues_) + items = cellfun(@num2str, num2cell(distinctValues_), 'UniformOutput', false); + elseif isstring(distinctValues_) + items = cellstr(distinctValues_); + elseif iscell(distinctValues_) && ~iscellstr(distinctValues_) + items = cellfun(@num2str, distinctValues_, 'UniformOutput', false); + else + items = distinctValues_; + end + items = items(:)'; % Ensure row vector + items = [{'All'}, items]; + end + + % Create a label for the field. + uilabel(container, ... + 'Text', sprintf('%s:', fieldName), ... + 'HorizontalAlignment', 'right', ... + 'Position', [10, currentY - 25, 150, 25]); + + % Create a dropdown for the field. + dd = uidropdown(container, ... + 'Items', items, ... + 'Value', 'All', ... + 'Position', [170, currentY - 25, 200, 25]); + + % Store the dropdown handle and its associated table/field. + dropdownHandles{end+1} = dd; + dropdownTableNames{end+1} = tableName; + dropdownFieldNames{end+1} = fieldName; + + currentY = currentY - heightPerField; + end + end + + % Create a "Submit" button at the bottom of the figure. + btn = uibutton(fig, 'Text', 'Submit', ... + 'Position', [200, 10, 100, 30], ... + 'ButtonPushedFcn', @(btn, event) uiresume(fig)); + + % Wait until the user clicks "Submit" or closes the figure. + uiwait(fig); + + % If the figure was closed (and thus no longer valid), return an empty struct. + if ~isvalid(fig) filterParams = struct(); return; end - % Parse user input into filterParams structure + % Build the filterParams struct from the dropdown selections. filterParams = struct(); - - for i = 1:numel(allFieldsFullName) - value = settings.(convertedFieldNames{i}); - - % Split full name to get table and field names - fieldParts = strsplit(allFieldsFullName{i}, '.'); - tableName = fieldParts{1}; - fieldName = fieldParts{2}; - - % If the table does not exist in the filterParams struct, create it + for k = 1:numel(dropdownHandles) + tableName = dropdownTableNames{k}; + fieldName = dropdownFieldNames{k}; + value = dropdownHandles{k}.Value; if ~isfield(filterParams, tableName) filterParams.(tableName) = struct(); end - % Assign values to the respective fields under each table + % If "All" is selected, assign empty; otherwise, try converting to numeric. if strcmp(value, 'All') - filterParams.(tableName).(fieldName) = []; % Set to empty to include all values - elseif isnumeric(value) && isnan(value) - filterParams.(tableName).(fieldName) = NaN; % Use NaN to handle as NULL + filterParams.(tableName).(fieldName) = []; else - filterParams.(tableName).(fieldName) = value; % Use the entered value + numValue = str2double(value); + if ~isnan(numValue) + filterParams.(tableName).(fieldName) = numValue; + else + filterParams.(tableName).(fieldName) = value; + end end end + + % Close the figure. + delete(fig); end - - - end end diff --git a/Classes/Moveit_wrapper.m b/Classes/Moveit_wrapper.m index 025a464..96dce6a 100644 --- a/Classes/Moveit_wrapper.m +++ b/Classes/Moveit_wrapper.m @@ -39,17 +39,22 @@ classdef Moveit_wrapper < handle end - function signalclass_out = process(obj,signal_in) + function signal_out = process(obj,signal_in) + + arguments + obj + signal_in = [] + + end isSignalClass = 0; if isa(signal_in,"Signal") - signalclass_out = signal_in; + signal_out = signal_in; signal_in = signal_in.signal; isSignalClass = 1; end - input_len = length(signal_in); - signal_in = signal_in.'; + % signal_in = signal_in.'; % INIT MOVEIT global loop; @@ -59,18 +64,14 @@ classdef Moveit_wrapper < handle % RUN MOVEIT global loop; loop = 1; - signal_out = obj.moveit_init(signal_out); - % CHECK IF SIGNAL CHANGED - if input_len ~= length(signal_out) - warning(['Signal length changed in moveit function: ', obj.moveit_function_name]); - end + signal_out = obj.moveit_init(signal_out); if isSignalClass % append to logbook lbdesc = ['Logbookentry']; - signalclass_out = signalclass_out.logbookentry(lbdesc); - signalclass_out.signal = signal_out; + signal_out = signal_out.logbookentry(lbdesc); + signal_out.signal = signal_out; end end diff --git a/Datatypes/equalizer_structure.m b/Datatypes/equalizer_structure.m index 263260b..c5bdeb2 100644 --- a/Datatypes/equalizer_structure.m +++ b/Datatypes/equalizer_structure.m @@ -4,8 +4,9 @@ classdef equalizer_structure < int32 ffe (0) vnle (1) vnle_pf_mlse (2) - db_precoded (3) - db_encoded (4) + % db_precoded (3) + % db_encoded (4) + vnle_db_mlse (3) end end \ No newline at end of file diff --git a/Functions/EQ_structures/duobinary_target.m b/Functions/EQ_structures/duobinary_target.m index d2d3671..cbb880e 100644 --- a/Functions/EQ_structures/duobinary_target.m +++ b/Functions/EQ_structures/duobinary_target.m @@ -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 \ No newline at end of file diff --git a/Functions/EQ_structures/vnle.m b/Functions/EQ_structures/vnle.m index 3bcc18a..1447351 100644 --- a/Functions/EQ_structures/vnle.m +++ b/Functions/EQ_structures/vnle.m @@ -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 diff --git a/Functions/EQ_structures/vnle_postfilter_mlse.m b/Functions/EQ_structures/vnle_postfilter_mlse.m index dc1501f..aa5e716 100644 --- a/Functions/EQ_structures/vnle_postfilter_mlse.m +++ b/Functions/EQ_structures/vnle_postfilter_mlse.m @@ -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) diff --git a/Functions/EQ_visuals/showEQcoefficients.m b/Functions/EQ_visuals/showEQcoefficients.m index 3794212..cb747b8 100644 --- a/Functions/EQ_visuals/showEQcoefficients.m +++ b/Functions/EQ_visuals/showEQcoefficients.m @@ -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 \ No newline at end of file + sgtitle('Filter Coefficients'); % Overall title for the figure +end diff --git a/Functions/EQ_visuals/showEQfilter.m b/Functions/EQ_visuals/showEQfilter.m new file mode 100644 index 0000000..3ef6732 --- /dev/null +++ b/Functions/EQ_visuals/showEQfilter.m @@ -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 \ No newline at end of file diff --git a/Functions/EQ_visuals/showLevelHistogram.m b/Functions/EQ_visuals/showLevelHistogram.m index 40f31f4..66e5d08 100644 --- a/Functions/EQ_visuals/showLevelHistogram.m +++ b/Functions/EQ_visuals/showLevelHistogram.m @@ -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 diff --git a/Functions/Metrics/calc_air_plain.m b/Functions/Metrics/calc_air_plain.m new file mode 100644 index 0000000..757c9c8 --- /dev/null +++ b/Functions/Metrics/calc_air_plain.m @@ -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)0; + +doub_mode = db_mode.no_db; + +%%% change specific parameter if given in varargin +% Parse optional input arguments +if ~isempty(varargin) + var_s = varargin{1}; + if isstruct(var_s) + fields = fieldnames(var_s); + for i = 1:numel(fields) + if isnumeric(fields{i}) + eval([fields{i}, ' = ', num2str( var_s.(fields{i}) ), ';']); + fprintf("%s <-- %.2f \n", fields{i}, var_s.(fields{i})); + else + eval([fields{i}, ' = ', 'var_s.(fields{',num2str(i),'})' , ';']); + end + + end + else + error('Optional variables should be passed as a struct.'); + end +end + + + +if doub_mode ~= db_mode.db_encoded + if precomp == 0 && db_precode == 1 + doub_mode = db_mode.db_precoded; + + db_precode = 1; % preceded data (in my measurement set, this corresponds to low precomp too!) + discard_precode = 0; % + emulate_precode = 0; + legendentry = 'low precomp; precoded'; + disp('low precomp; precoded') + + elseif precomp == 1 && db_precode == 1 + doub_mode = db_mode.db_emulate; + + db_precode = 0; % preceded data (in my measurement set, this corresponds to low precomp too!) + discard_precode = 0; % + emulate_precode = 1; + legendentry = 'high precomp; precoded'; + disp('high precomp; precoded') + + elseif precomp == 0 && db_precode == 0 + doub_mode = db_mode.db_discard; + + db_precode = 1; % preceded data (in my measurement set, this corresponds to low precomp too!) + discard_precode = 1; % + emulate_precode = 0; + legendentry = 'no precomp; not precoded'; + disp('no precomp; not precoded') + + elseif precomp == 1 && db_precode == 0 + doub_mode = db_mode.no_db; + + db_precode = 0; % preceded data (in my measurement set, this corresponds to low precomp too!) + discard_precode = 0; % + emulate_precode = 0; + legendentry = 'high precomp; not precoded'; + disp('high precomp; not precoded') + + end + +else + +end + +fsym_ = floor( bitrate*1e-9./log2(M) ).*1e9; + +if fsym_ ~= fsym + fsym = fsym_; + % fprintf('Adapted symbolrate to %d GBd, to match provided bitrate of %d GBit/s using PAM %d \n',fsym.*1e-9,bitrate.*1e-9, M); +end +f_nyquist = fsym/2; + + +basePath = 'C:\Users\Silas\Documents\MATLAB\Datensätze\sioe_labor\'; +database = DBHandler("pathToDB",[basePath,'silas_labor.db']); + +useGui = 0; +% db = DBHandler("pathToDB",[basePath,'silas_labor.db']); +filterParams = database.tables; +% filterParams.Runs.run_id = 2958; % no db +% filterParams.Runs.run_id = 2937; % no db +filterParams.Configurations = struct( ... + 'bitrate', bitrate, ... + 'db_mode', db_precode+db_encode, ... + 'fiber_length', link_length, ... + 'interference_attenuation', [], ... + 'interference_path_length', [], ... + 'is_mpi', is_mpi, ... + 'pam_level', M, ... + 'precomp_amp', [], ... + 'rop_attenuation', 0, ... + 'symbolrate', [], ... + 'v_awg', [], ... + 'v_bias', [], ... + 'wavelength', laser_wavelength ... + ); + +selectedFields = {'Runs.run_id','Runs.tx_bits_path', 'Runs.tx_symbols_path', 'Runs.rx_sync_path','Runs.rx_raw_path',... + 'Configurations.db_mode','Configurations.pam_level','Configurations.bitrate','Configurations.symbolrate','Configurations.fiber_length','Configurations.wavelength','Configurations.precomp_amp','Measurements.power_rop','Configurations.v_bias',... + 'Configurations.interference_attenuation'}; + +[dataTable,sql_query] = database.queryDB(filterParams, selectedFields); +[~, uniqueIdx] = unique(dataTable.run_id); % Get unique run_id indices +dataTable = dataTable(uniqueIdx,:); % Extract unique configurations for each run_id +fprintf('Found %d entries for requested Configuration. IDs are: %s \n \n',size(dataTable,1),jsonencode(dataTable.run_id(1:min(size(dataTable,1),100)))); + + +output = struct(); +vnle_pf_package = {}; +vnle_dfe_package = {}; +dbtgt_package = {}; + +disp(num2str(bitrate)) + +for iatt = 1:numel(dataTable.interference_attenuation) + + current_run_id = dataTable.run_id(iatt); + + Tx_bits = load([basePath, char(dataTable.tx_bits_path(iatt))]); + Tx_bits = Tx_bits.Bits; + Symbols_mapped = PAMmapper(M,0).map(Tx_bits); + Symbols_mapped.fs = fsym; + + Symbols = load([basePath, char(dataTable.tx_symbols_path(iatt))]); + Symbols = Symbols.Symbols; + + Scpe_load = load([basePath, char(dataTable.rx_sync_path(iatt))]); + Scpe_cell = Scpe_load.S; + [~,~,found]=Scpe_cell{2}.tsynch("reference",Symbols,"fs_ref",fsym,"debug_plots",1); + + if ~found + Raw_signal = load([basePath, char(dataTable.rx_raw_path(1))]); + Raw_signal = Raw_signal.Scpe_sig_raw; + [~,Scpe_cell,found] =Raw_signal.tsynch("reference",Symbols,"fs_ref",fsym,"debug_plots",0); + end + + if ~found + if length(Symbols_mapped.signal) == sum(Symbols_mapped.signal == Symbols.signal) + warning('Could not synchronize the received signal with the stored symbols!') + else + [~,Scpe_cell,found] =Raw_signal.tsynch("reference",Symbols_mapped,"fs_ref",fsym,"debug_plots",0); + end + if ~found + warning('Could not synchronize the received signal with the stored symbols!') + end + end + + + % + % Raw_signal = Filter('filtdegree',4,"f_cutoff",Symbols.fs.*0.55,"fs",Raw_signal.fs,"filterType",filtertypes.gaussian,"active",true).process(Raw_signal); + % + % Scpe_cell{1}.eye(fsym,M,"displayname",'eye','fignum',227); + % + % Raw_signal.spectrum("normalizeTo0dB",0,"fignum",11,"fft_length",2^12); + % Raw_signal.move_it_spectrum("fignum",334); + % Raw_signal.move_it_spectrum("fignum",334); + + fsym = Symbols.fs; + + if db_precode + Symbols_precoded = Symbols; + end + + + + proc_occ = min(15,length(Scpe_cell)); + for occ = 1:proc_occ + + Scpe_sig = Scpe_cell{occ}; + + %%%%%% Sample to 2x fsym %%%%%% + Scpe_sig = Scpe_sig.resample("fs_out",2*fsym); + + %%%%%% Sync Rx signal with reference %%%%%% + [Scpe_sig,~] = Scpe_sig.tsynch("reference",Symbols,"fs_ref",fsym,"debug_plots",0); + + Scpe_sig = Filter('filtdegree',4,"f_cutoff",Symbols.fs.*0.5,"fs",Scpe_sig.fs,"filterType",filtertypes.gaussian,"active",true).process(Scpe_sig); + + Scpe_sig = Scpe_sig - mean(Scpe_sig.signal); + % + % Pform = Pulseformer("fsym",Scpe_sig.fs,"fdac",2*fsym,"pulse","rrc","pulselength",16,"alpha",rcalpha,"matched",0); + % + % Scpe_sig_matched = Pform.process(Scpe_sig); + % + % Scpe_sig.spectrum("normalizeTo0dB",0,"fignum",336,"displayname","scope "); + % Scpe_sig_matched.spectrum("normalizeTo0dB",0,"fignum",336,"displayname","matched"); + + %%% EQUALIZING + + + % eq_mlse = FFE_DCremoval("epochs_tr",5,"epochs_dd",5,"len_tr",len_tr,"mu_dd",mu_ffe(1),"mu_tr",0,"order",ffe_order(1),"sps",2,"decide",0,"dc_buffer_len",1,"mu_dc",0.05); + % eq_mlse = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",len_tr,"mu_dd",mu_ffe(1),"mu_tr",0,"order",ffe_order(1),"sps",2,"decide",0); + % eq_mlse = FFE_DCremoval("epochs_tr",5,"epochs_dd",5,"len_tr",len_tr,"mu_dd",mu_ffe(1),"mu_tr",0,"order",ffe_order(1),"sps",2,"decide",0,"dc_buffer_len",512,"mu_dc",0.05); + + mu_ffe = [mu_ffe1 mu_ffe2 mu_ffe3]; + vnle_order=[vnle_order1,vnle_order2,vnle_order3]; + + % %%%%% VNLE + DFE %%%% + if 0 + + eq_vnle_dfe = EQ("Ne",vnle_order,"Nb",[0,0,0],"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_2 = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",2001,"sps",1,"decide",0); + + [result] = vnle(eq_vnle_dfe,M,Scpe_sig,Symbols,Tx_bits,"precode_mode",doub_mode,"showAnalysis",1,"postFFE",[]); + vnle_dfe_package{iatt,occ} = result; + + end + %%%%% VNLE + PF + MLSE %%%% + if 1 + + try + % len_tr = length(Symbols)-1000; + eq_vnle_ = EQ("Ne",vnle_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_vnle_ = VNLE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",[0.0004 0.0005 0.0006],"mu_tr",0,"order",vnle_order,"sps",2,"decide",0); + pf_ = Postfilter("ncoeff",pf_ncoeffs,"useBurg",1); + mlse_ = MLSE_viterbi("duobinary_output",0,'M',M,'trellis_states',PAMmapper(M,0).levels); + eq_2 = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",2001,"sps",1,"decide",0); + + [result] = vnle_postfilter_mlse(eq_vnle_,pf_,mlse_,M,Scpe_sig,Symbols,Tx_bits,"precode_mode",doub_mode,'showAnalysis',0,"postFFE",[]); + vnle_pf_package{iatt,occ} = result; + + database.addProcessingResult(current_run_id,result.resultsMLSE, result.equalizerConfigMLSE); + database.addProcessingResult(current_run_id,result.resultsVNLE, result.equalizerConfigVNLE); + catch + warning(['VNLE+MLSE fail: run id: ', num2str(current_run_id)],' occ:', num2str(occ), ' iatten: ',num2str(iatt)) + end + + end + + + %%%%% Duobinary Targeting %%%% + if 1 + + try + mlse_db = MLSE_viterbi("DIR",[1,1],"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels); + eq_db = EQ("Ne",vnle_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_2 = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",2001,"sps",1,"decide",0); + + [result] = duobinary_target(eq_db, mlse_db, M, Scpe_sig, Symbols, Tx_bits, "precode_mode", doub_mode,'showAnalysis',0,"postFFE",[]); + dbtgt_package{iatt,occ} = result; + + database.addProcessingResult(current_run_id,result.resultsDBtgt, result.equalizerConfigDBtgt); + + catch + warning(['VNLE DB+MLSE fail: run id: ', num2str(current_run_id)],' occ:', num2str(occ), ' iatten: ',num2str(iatt)) + end + + end + + %%%%%% %db signaling => db encoded %%%%% + if 0 + mlse_db_enc = MLSE_viterbi("DIR",[1,1],"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels); + eq_db_enc = EQ("Ne",vnle_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); + [result] = duobinary_signaling(eq_db_enc, mlse_db_enc,M, Scpe_sig ,Symbols, Tx_bits); + dbenc_package{iatt,occ} = result; + end + + + % autoArrangeFigures; + disp('- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ') + fprintf('\n') + + + end + + + + if ~isempty(curFolder) + cd(curFolder); + end + +end + +output.dataTable = dataTable; +output.vnle_dfe_package = vnle_dfe_package; +output.vnle_pf_package = vnle_pf_package; +output.dbtgt_package = dbtgt_package; \ No newline at end of file diff --git a/projects/HighSpeedExperiment_2024/auswertung MPI/plot_from_database.m b/projects/HighSpeedExperiment_2024/auswertung MPI/plot_from_database.m new file mode 100644 index 0000000..6eba0dd --- /dev/null +++ b/projects/HighSpeedExperiment_2024/auswertung MPI/plot_from_database.m @@ -0,0 +1,91 @@ + + +basePath = 'C:\Users\Silas\Documents\MATLAB\Datensätze\sioe_labor\'; +database = DBHandler("pathToDB",[basePath,'silas_labor.db']); + +filterParams = database.tables; +filterParams.Configurations = struct( ... + 'bitrate', 336e9, ... + 'db_mode', 0, ... + 'fiber_length', 1, ... + 'interference_attenuation', [], ... + 'interference_path_length', [], ... + 'is_mpi', 1, ... + 'pam_level', 4, ... + 'rop_attenuation', 0 ... + ); + +filterParams.EqualizerParameters.diff_precode = int32(db_mode.no_db); +filterParams.EqualizerParameters.equalizer_structure = int32(equalizer_structure.vnle); + +selectedFields = {'Configurations.run_id' 'Runs.rx_raw_path' 'Configurations.bitrate' 'Configurations.symbolrate' 'Configurations.pam_level' 'Configurations.db_mode' 'Configurations.rop_attenuation' 'Configurations.is_mpi' 'Configurations.interference_attenuation' 'EqualizerParameters.equalizer_structure' 'EqualizerParameters.diff_precode' 'EqualizerParameters.eq_id' 'Measurements.power_pd_in' 'Measurements.power_mpi_interference' 'Measurements.power_mpi_signal' 'Results.BER' 'Results.SNR' 'Results.GMI' 'Results.Alpha'}; + +[dataTable,sql_query] = database.queryDB(filterParams, selectedFields); + +fixedVars = {'run_id','eq_id','bitrate'}; +resultTable = groupIt(fixedVars,dataTable); + + + +% Create a new figure +figure(1); +hold on +unique_rates = unique(resultTable.bitrate); +for i = 1:numel(unique_rates) + % Plot BER vs. interference_attenuation + plot(resultTable.power_mpi_signal(resultTable.bitrate==unique_rates(i),:)-resultTable.power_mpi_interference(resultTable.bitrate==unique_rates(i),:), resultTable.BER(resultTable.bitrate==unique_rates(i),:), 'o-', 'LineWidth', 1.5); +end + +% Label the axes and add a title +xlabel('Interference Attenuation'); +ylabel('BER'); +title('BER vs. Interference Attenuation'); + +% Enable grid for better readability +grid on; + +beautifyBERplot; + + + + +function resultTable = groupIt(fixedVars,dataTable) + + % Group by run_id and eq_id (adjust grouping keys as needed) + + [G, groupKeys] = findgroups(dataTable(:, fixedVars)); + + % Preallocate a cell array for aggregated data. + varNames = dataTable.Properties.VariableNames; + nVars = numel(varNames); + aggData = cell(height(groupKeys), nVars); + groupCount = zeros(height(groupKeys), 1); % To store the size of each group + + % Loop over each group. + for i = 1:height(groupKeys) + idx = (G == i); % Logical index for group i + groupCount(i) = sum(idx); % Count number of rows in this group + % For each variable in the table: + for j = 1:nVars + colData = dataTable.(varNames{j}); + if isnumeric(colData) + % For numeric data, compute the mean. + aggData{i, j} = mean(colData(idx)); + else + % For non-numeric data, take the first entry. + if iscell(colData) + aggData{i, j} = colData{find(idx, 1)}; + else + aggData{i, j} = colData(find(idx, 1)); + end + end + end + end + + % Convert the aggregated cell array into a table. + resultTable = cell2table(aggData, 'VariableNames', varNames); + + % Append the group count as a new column. + resultTable.nRows = groupCount; + +end \ No newline at end of file diff --git a/projects/IMDD_base_system/analyze_burg.m b/projects/HighSpeedExperiment_2024/auswertung deliv 06 feb 2025/analyze_burg.m similarity index 100% rename from projects/IMDD_base_system/analyze_burg.m rename to projects/HighSpeedExperiment_2024/auswertung deliv 06 feb 2025/analyze_burg.m diff --git a/projects/IMDD_base_system/analyze_equalizer_taps.m b/projects/HighSpeedExperiment_2024/auswertung deliv 06 feb 2025/analyze_equalizer_taps.m similarity index 100% rename from projects/IMDD_base_system/analyze_equalizer_taps.m rename to projects/HighSpeedExperiment_2024/auswertung deliv 06 feb 2025/analyze_equalizer_taps.m diff --git a/projects/IMDD_base_system/analyze_lambda.m b/projects/HighSpeedExperiment_2024/auswertung deliv 06 feb 2025/analyze_lambda.m similarity index 100% rename from projects/IMDD_base_system/analyze_lambda.m rename to projects/HighSpeedExperiment_2024/auswertung deliv 06 feb 2025/analyze_lambda.m diff --git a/projects/IMDD_base_system/analyze_postfilter_taps.m b/projects/HighSpeedExperiment_2024/auswertung deliv 06 feb 2025/analyze_postfilter_taps.m similarity index 100% rename from projects/IMDD_base_system/analyze_postfilter_taps.m rename to projects/HighSpeedExperiment_2024/auswertung deliv 06 feb 2025/analyze_postfilter_taps.m diff --git a/projects/HighSpeedExperiment_2024/auswertung/App/auswertungs_app.mlapp b/projects/HighSpeedExperiment_2024/auswertung/App/auswertungs_app.mlapp index 375605d..9fb4d4a 100644 Binary files a/projects/HighSpeedExperiment_2024/auswertung/App/auswertungs_app.mlapp and b/projects/HighSpeedExperiment_2024/auswertung/App/auswertungs_app.mlapp differ diff --git a/projects/IMDD_base_system/imdd_it.m b/projects/IMDD_base_system/imdd_it.m index 97645f1..d2fdf3c 100644 --- a/projects/IMDD_base_system/imdd_it.m +++ b/projects/IMDD_base_system/imdd_it.m @@ -1,25 +1,22 @@ % basePath = 'C:\Users\Silas\Documents\MATLAB\Datensätze\sioe_labor\'; % db = DBHandler("pathToDB",[basePath,'silas_labor.db']); -if 0 +if 1 uloops = struct; uloops.precomp = [1]; uloops.db_precode = [0]; - uloops.bitrate = [300,330,360,390,420,450,480].*1e9; %[300,330,360,390,420,450,480] + uloops.bitrate = [224].*1e9; %[300,330,360,390,420,450,480] [224,336,360,390,420,448] for MPI % uloops.laser_wavelength = [1293,1297.5,1302,1306.5,1310,1313.4,1318,1322.7,1327.4]; - uloops.precomp = [0]; - uloops.db_precode = [0]; - uloops.bitrate = [330].*1e9; %[300,330,360,390,420,450,480] uloops.laser_wavelength = [1310]; - uloops.M = [4,6,8]; uloops.M = [4]; - uloops.link_length = [2]; % 1,2,3,5,6,8,10 + uloops.link_length = [1]; % 1,2,3,5,6,8,10 + uloops.interference_attenuation = [0,3,6,9,12,15,18,21,24,27,30,45]; wh = DataStorage(uloops); wh.addStorage("ber"); - wh = submit_simulations(wh,"parallel",0,"simulation_mode",0); - - wh = submit_simulations(wh,"parallel",0,"simulation_mode",1); + % wh = submit_simulations(wh,"parallel",0,"simulation_mode",0); + wh = submit_handle(@imdd_model,wh,"parallel",1); + end wh_ana = wh_master; diff --git a/projects/IMDD_base_system/imdd_model.m b/projects/IMDD_base_system/imdd_model.m index 538ce7b..c15a919 100644 --- a/projects/IMDD_base_system/imdd_model.m +++ b/projects/IMDD_base_system/imdd_model.m @@ -1,4 +1,6 @@ -function [output] = imdd_model(simulation_mode,varargin) +function [output] = imdd_model(varargin) + +simulation_mode = 0; %%% Change folder curFolder = pwd; @@ -17,6 +19,9 @@ fdac = 256e9; fadc = 256e9; random_key = 1; +interference_attenuation = 0; +is_mpi = 1; + precomp = 0; db_precode = 0; @@ -216,7 +221,7 @@ else 'fiber_length', link_length, ... 'interference_attenuation', [], ... 'interference_path_length', [], ... - 'is_mpi', 0, ... + 'is_mpi', is_mpi, ... 'pam_level', M, ... 'precomp_amp', [], ... 'rop_attenuation', 0, ... @@ -227,7 +232,8 @@ else ); selectedFields = {'Runs.run_id','Runs.tx_bits_path', 'Runs.tx_symbols_path', 'Runs.rx_sync_path','Runs.rx_raw_path',... - 'Configurations.db_mode','Configurations.pam_level','Configurations.bitrate','Configurations.symbolrate','Configurations.fiber_length','Configurations.wavelength','Configurations.precomp_amp','Measurements.power_rop','Configurations.v_bias'}; + 'Configurations.db_mode','Configurations.pam_level','Configurations.bitrate','Configurations.symbolrate','Configurations.fiber_length','Configurations.wavelength','Configurations.precomp_amp','Measurements.power_rop','Configurations.v_bias',... + 'Configurations.interference_attenuation'}; [dataTable,sql_query] = database.queryDB(filterParams, selectedFields); [~, uniqueIdx] = unique(dataTable.run_id); % Get unique run_id indices @@ -245,6 +251,7 @@ else % Raw_signal = load([basePath, char(dataTable.rx_raw_path(1))]); + % Raw_signal.Scpe_sig_raw.plot("displayname",'0db atten','fignum',10101) % Raw_signal = Raw_signal.Scpe_sig_raw; % % Raw_signal = Filter('filtdegree',4,"f_cutoff",Symbols.fs.*0.55,"fs",Raw_signal.fs,"filterType",filtertypes.gaussian,"active",true).process(Raw_signal); @@ -270,7 +277,7 @@ dbtgt_package = {}; proc_occ = min(1,length(Scpe_cell)); -for occ = 1:proc_occ +for occ = 1%:proc_occ Scpe_sig = Scpe_cell{occ}; @@ -303,14 +310,16 @@ for occ = 1:proc_occ % %%%%% VNLE + DFE %%%% if 0 + eq_vnle_dfe = EQ("Ne",vnle_order,"Nb",[0,0,0],"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); - - [result] = vnle(eq_vnle_dfe,M,Scpe_sig,Symbols,Tx_bits,"precode_mode",doub_mode,"showAnalysis",0); + eq_2 = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",2001,"sps",1,"decide",0); + + [result] = vnle(eq_vnle_dfe,M,Scpe_sig,Symbols,Tx_bits,"precode_mode",doub_mode,"showAnalysis",1,"postFFE",[]); vnle_dfe_package{occ} = result; end %%%%% VNLE + PF + MLSE %%%% - if 0 + if 1 % len_tr = length(Symbols)-1000; eq_vnle_ = EQ("Ne",vnle_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); @@ -325,7 +334,7 @@ for occ = 1:proc_occ %%%%% Duobinary Targeting %%%% - if 0 + if 1 mlse_db = MLSE_viterbi("DIR",[1,1],"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels); eq_db = EQ("Ne",vnle_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); diff --git a/projects/Messung_Zürich/change_filenames.m b/projects/Messung_Zürich/change_filenames.m new file mode 100644 index 0000000..0d3afd9 --- /dev/null +++ b/projects/Messung_Zürich/change_filenames.m @@ -0,0 +1,85 @@ +%% 1) Read all files in the TR folder +pathToTimingRecov = "C:\Users\Silas\Documents\MATLAB\Datensätze\IEF_Polariton_2025\36_IMDD_Kiel\Data\TR_ZIP"; +trFiles = dir(fullfile(pathToTimingRecov, 'TR_SILAS_*.mat')); + +%% Directory for measurement files (used to extract metadata) +pathToMeasurement = "C:\Users\Silas\Documents\MATLAB\Datensätze\IEF_Polariton_2025\36_IMDD_Kiel\Data\20250221"; + +%% Initialize lists for different PAM types +listPAM2 = {}; +listPAM4 = {}; +listPAM6 = {}; +listPAM8 = {}; + +%% Loop over each TR file +for k = 1:length(trFiles) + % Get current TR file name + trFileName = trFiles(k).name; + + % 2) Extract file code from TR file name. + % For a filename like "TR_SILAS_20250221T001800.mat" the code is "20250221T001800". + filecode = extractBetween(trFileName, 'TR_SILAS_', '.mat'); + + % Find corresponding measurement file by code (custom function) + measurementFilename = findFileByCode(pathToMeasurement, filecode{1}); + + % 3) Extract parameters from the measurement filename using regex. + % Expected measurement filename format (example): + % "Pmod_-10p000dBm_P_PD_-20p000dBm_..._32GBd_4PAM__1234T5678" + tokens = regexp(measurementFilename, ... + 'Pmod_([-0-9p]+)dBm_P_PD_([-0-9p]+)dBm_.*?_(\d+)GBd_(\d+)PAM__\d+T\d+', ... + 'tokens'); + if isempty(tokens) + error('Filename format not recognized for measurement file: %s', measurementFilename); + end + tokens = tokens{1}; + + % Convert token strings to numbers + config.P_laser = str2double(strrep(tokens{1}, 'p', '.')); + config.P_pd = str2double(strrep(tokens{2}, 'p', '.')); + config.fsym = str2double(tokens{3}) * 1e9; % Convert from GBd to Hz + config.M = str2double(tokens{4}); + + % Display loaded metadata + fprintf('Loaded measurement file: %s\n', measurementFilename); + fprintf('P_laser: %.3f dBm\n', config.P_laser); + fprintf('P_pd: %.3f dBm\n', config.P_pd); + fprintf('fsym: %.1f GBd\n', config.fsym * 1e-9); + fprintf('M: %d\n', config.M); + + % 4) Rename the TR file to include the metadata. + % New filename format: TR_SILAS__Pmod_dBm_P_PD_dBm_GBd_PAM.mat + newTRname = sprintf('TR_SILAS_%s_Pmod_%.3fdBm_P_PD_%.3fdBm_%dGBd_%dPAM', ... + filecode{1}, config.P_laser, config.P_pd, config.fsym/1e9, config.M); + + newTRname = strrep(newTRname,'.','p'); + newTRname = [newTRname, '.mat']; + movefile(fullfile(pathToTimingRecov, trFileName), fullfile(pathToTimingRecov, newTRname)); + + % Append the file code to the corresponding PAM list based on config.M + switch config.M + case 2 + listPAM2{end+1} = filecode{1}; + case 4 + listPAM4{end+1} = filecode{1}; + case 6 + listPAM6{end+1} = filecode{1}; + case 8 + listPAM8{end+1} = filecode{1}; + otherwise + warning('Unexpected PAM value %d in file %s', config.M, measurementFilename); + end +end + +%% Display the lists of file codes for each PAM type +disp('List of file codes for PAM2:'); +disp(listPAM2); + +disp('List of file codes for PAM4:'); +disp(listPAM4); + +disp('List of file codes for PAM6:'); +disp(listPAM6); + +disp('List of file codes for PAM8:'); +disp(listPAM8); diff --git a/projects/Messung_Zürich/dsp_all_files.m b/projects/Messung_Zürich/dsp_all_files.m index cde817c..9ec98ab 100644 --- a/projects/Messung_Zürich/dsp_all_files.m +++ b/projects/Messung_Zürich/dsp_all_files.m @@ -1,113 +1,366 @@ -% 202502211_km_fiber_Pmod_15dBm_P_PD_9p2494dBm_W03C34-0603E03_160GBd_8PAM__20250221T024256.h5 - -M = 8; -if M == 4 - file_codes = {"20250221T032043",... - "20250221T032312",... - "20250221T032424",... - "20250221T032529",... - "20250221T032632",... - "20250221T032800",... - "20250221T032936",... - "20250221T033035",... - "20250221T033138",... - "20250221T033138",... - "20250221T033425",... - "20250221T033527",... - "20250221T033642",... - "20250221T033743",... - "20250221T033851",... - "20250221T034314",... - "20250221T034529",... - "20250221T034647",... - "20250221T034756",... - "20250221T034915",... - }; - baudrate = [96,120,128,136,144,152,160,168,172,176,180,184,188,192,196,200,204,212,220,224]; - pam4.result_vnle = result_vnle; - pam4.result_db = result_db; +M = 6; +if M == 2 + file_codes = {'20250221T035221' '20250221T035354' '20250221T035824' '20250221T035931' '20250221T040035' '20250221T040132' '20250221T040226' '20250221T040523' '20250221T040646' '20250221T040723' '20250221T040843' '20250221T041011' '20250221T041101' '20250221T041244'}; +elseif M == 4 + file_codes = {'20250221T030844' '20250221T032043' '20250221T032312' '20250221T032424' '20250221T032529' '20250221T032632' '20250221T032800' '20250221T033035' '20250221T033138' '20250221T033246' '20250221T033425' '20250221T033527' '20250221T033642' '20250221T033743' '20250221T033851' '20250221T034314' '20250221T034529' '20250221T034647' '20250221T034756' '20250221T034915'}; +elseif M == 6 + file_codes = {'20250221T041445' '20250221T041512' '20250221T041539' '20250221T041607' '20250221T041633' '20250221T041702' '20250221T041729' '20250221T041758' '20250221T041825' '20250221T041854' '20250221T041922' '20250221T041951' '20250221T042019' '20250221T042048' '20250221T042117' '20250221T042147' '20250221T042215'}; elseif M ==8 - file_codes = {"20250221T004651",... %96 - "20250221T004926",... - "20250221T005844",... - "20250221T010032",... - "20250221T010150",... - "20250221T025505",...%144 ohne balun - "20250221T024629",...%152 ohne balun - "20250221T024315",...%160 ohne balun - "20250221T023534",...%168 ohne balun - "20250221T024929",...%172 ohne balun - "20250221T025305",...%176 ohne balun - }; - file_codes = {"20250221T024256"}; - baudrate = [96,112,120,128,136,144,152,160,168,172,176]; + file_codes = {'20250221T004926' '20250221T023534' '20250221T024256' '20250221T024629' '20250221T024929' '20250221T025305' '20250221T025505' '20250221T025856' '20250221T030122' '20250221T030311' '20250221T030513'}; end -if 1 - +if 0 uloops = struct; uloops.filecode = file_codes; - uloops.mu_dc = [0]; + uloops.vnle_order1 = [50]; + uloops.vnle_order2 = [5]; + uloops.vnle_order3 = [5]; wh = DataStorage(uloops); wh.addStorage("ber"); - - wh = submit_handle(@dsp_ief_file,wh,"parallel",0); + + wh = submit_handle(@dsp_ief_file,wh,"parallel",1); +end + +if 0 + % Bring figure(5) to focus + fig = figure(2); + + % Get handles to all line objects in the figure + lines = findobj(fig, 'Type', 'line'); + + % Preallocate cell arrays to store the data for each line + xData = cell(numel(lines),1); + yData = cell(numel(lines),1); + + % Loop through each line and extract its data + + for k = 1:numel(lines) + xData{k} = get(lines(k), 'XData'); + yData{k} = get(lines(k), 'YData'); + end + + ief.M2.baudr_new = xData{4}; + ief.M2.ngmi_new = yData{4}; + + ief.M4.baudr_new = xData{3}; + ief.M4.ngmi_new = yData{3}; + + ief.M6.baudr_new = xData{2}; + ief.M6.ngmi_new = yData{2}; + + ief.M8.baudr_new = xData{1}; + ief.M8.ngmi_new = yData{1}; + + ief.M2.ber_new = yData{4}; + ief.M4.ber_new = yData{3}; + ief.M6.ber_new = yData{2}; + ief.M8.ber_new = yData{1}; + + pam2_baudr = xData{4}; + pam2_ber = yData{4}; + pam2_ngmi = yData{4}; + + pam4_baudr = xData{3}; + pam4_ber = yData{3}; + pam4_ngmi = yData{3}; + + pam6_baudr = xData{2}; + pam6_ber = yData{2}; + pam6_ngmi = yData{2}; + + pam8_baudr = xData{1}; + pam8_ber = yData{1}; + pam8_ngmi = yData{1}; + + pam2_baudr_lowdsp = xData{4}; + pam2_ber_lowdsp = yData{4}; + pam2_ngmi_lowdsp = yData{4}; + + pam4_baudr_lowdsp = xData{3}; + pam4_ber_lowdsp = yData{3}; + pam4_ngmi_lowdsp = yData{3}; + + pam6_baudr_lowdsp = xData{2}; + pam6_ber_lowdsp = yData{2}; + pam6_ngmi_lowdsp = yData{2}; + + pam8_baudr_lowdsp = xData{1}; + pam8_ber_lowdsp = yData{1}; + pam8_ngmi_lowdsp = yData{1}; + + + end -a = wh.getStoValue('ber',uloops.filecode, uloops.mu_ffe1); +close all +figure(1) +cnt = 1; +for M = [2,4,6,8] + % for i = 1:numel(uloops.vnle_order2) + % for j = 1:numel(uloops.vnle_order3) + wh = eval(sprintf('wh_pam%d',M)); + a =wh.getStoValue('ber',wh.parameter.filecode.values, wh.parameter.vnle_order1.values, wh.parameter.vnle_order2.values,wh.parameter.vnle_order3.values); + % a =wh_mit_1001.getStoValue('ber',uloops.filecode, wh_mit_1001.parameter.vnle_order1.values, wh_mit_1001.parameter.vnle_order2.values,wh_mit_1001.parameter.vnle_order3.values); + baudrate = cell2mat(cellfun(@(a) a.config.fsym, a, 'UniformOutput', false)); + [baudrate, idx] = sort(baudrate); -% get best results per baudrate -ber_vnle_values = cellfun(@(a) cellfun(@(y) y.ber_vnle, a.vnle_package, 'UniformOutput', false), a, 'UniformOutput', false); -best_vnle = cellfun(@(x) min(cell2mat(x)), ber_vnle_values); + % get best results per baudrate + vnle = 0; + try + ber_vnle_values = cellfun(@(a) cellfun(@(y) y.ber_vnle, a.vnle_package, 'UniformOutput', false), a, 'UniformOutput', false); + ber_vnle_infrate = cellfun(@(a) cellfun(@(y) y.inf_rate_vnle, a.vnle_package, 'UniformOutput', false), a, 'UniformOutput', false); -% ber_vnle_values = cellfun(@(a) cellfun(@(y) y.ber_vnle, a.vnle_pf_package, 'UniformOutput', false), a, 'UniformOutput', false); -% best_vnle = cellfun(@(x) min(cell2mat(x)), ber_vnle_values); -% -% ber_mlse_values = cellfun(@(a) cellfun(@(y) y.ber_mlse, a.vnle_pf_package, 'UniformOutput', false), a, 'UniformOutput', false); -% best_mlse = cellfun(@(x) min(cell2mat(x)), ber_mlse_values); + best_vnle = cellfun(@(x) min(cell2mat(x)), ber_vnle_values); + best_vnle = best_vnle(idx); + best_gmi_vnle = cellfun(@(x) min(cell2mat(x)), ber_vnle_infrate); + best_gmi_vnle = best_gmi_vnle(idx); -%%% PLOT PARAMETER OPTIMIZATION -figure(2) -cols = cbrewer2('Set1',6); -hold on -title(sprintf('%d km | %d nm | PAM %d',1,1313.5,M)); -xax = uloops.mu_ffe1; -plot(xax,best_vnle,'DisplayName',sprintf('VNLE 50,3,3'),'Color',cols(o,:),'LineStyle','-','HandleVisibility','on'); + % netRateVNLE = TransmissionPerformance().calculateNetRate(log2(M)*baudrate,'NGMI',best_gmi_vnle./log2(M), 'BER',best_vnle); -% plot(xax,best_mlse,'DisplayName',sprintf('VNLE + 1 tap post-filter + MLSE'),'Color',cols(2,:),'LineStyle','-','HandleVisibility','on'); + vnle = 1; + end -% plot(xax,best_db,'DisplayName',sprintf('DB tgt. + MLSE'),'Color',cols(3,:),'LineStyle','-','HandleVisibility','on'); + mlse = 0; + try + ber_vnle_values = cellfun(@(a) cellfun(@(y) y.ber_vnle, a.vnle_pf_package, 'UniformOutput', false), a, 'UniformOutput', false); + best_vnle_2 = cellfun(@(x) min(cell2mat(x)), ber_vnle_values); + best_vnle_2 = best_vnle_2(idx); -set(gca, 'YScale', 'log'); -ylim([5e-5 0.3]); -xlim([min(xax), max(xax) ]); -yline([3.8e-3, 2e-2],'HandleVisibility','off'); -legend -beautifyBERplot() -xlabel('Adaption speed'); -ylabel('BER'); + gmi_vnle_values = cellfun(@(a) cellfun(@(y) y.gmi, a.vnle_pf_package, 'UniformOutput', false), a, 'UniformOutput', false); + best_gmi_vnle = cellfun(@(x) min(cell2mat(x)), gmi_vnle_values); + best_gmi_vnle = best_gmi_vnle(idx); + + ber_mlse_values = cellfun(@(a) cellfun(@(y) y.ber_mlse, a.vnle_pf_package, 'UniformOutput', false), a, 'UniformOutput', false); + best_mlse = cellfun(@(x) min(cell2mat(x)), ber_mlse_values); + best_mlse = best_mlse(idx); + + if M == 8 + best_vnle_2(2) = []; + best_mlse(2) = []; + best_gmi_vnle(2) = []; + best_db(2) = []; + + baudrate(2) = []; + elseif M == 4 + + end + + mlse = 1; + end + + db = 0; + try + ber_db_values = cellfun(@(a) cellfun(@(y) y.ber, a.dbtgt_package, 'UniformOutput', false), a, 'UniformOutput', false); + best_db = cellfun(@(x) min(cell2mat(x)), ber_db_values); + best_db = best_db(idx); + if M == 4 + best_db(end-1:end) = []; + baudrate(end-1:end) = []; + end + + db = 1; + end + + %% BER PLOT + + % figure(1) + subplot(1,2,1) + cols = cbrewer2('Set1',6); + hold on + title(sprintf('BER')); + % title(sprintf('%d 1st order',uloops.vnle_order1(i))); + xax = baudrate.*1e-9; + + if M == 4 || M == 2 + plot(xax,best_db,'DisplayName',sprintf('PAM %d | DB+MLSE',M),'Color',cols(M/2,:),'LineStyle',':','HandleVisibility','on','Marker','o','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + plot(ief.(sprintf('M%d', M)).baudr_new,ief.(sprintf('M%d', M)).ber_new,'DisplayName',sprintf('PAM %d |VNLE',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + elseif M == 6 + plot(xax,best_vnle_2,'DisplayName',sprintf('PAM %d |VNLE',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + % plot(ief.(sprintf('M%d', M)).baudr_new,ief.(sprintf('M%d', M)).ber_new,'DisplayName',sprintf('PAM %d |VNLE',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + else + % plot(xax,best_vnle_2,'DisplayName',sprintf('PAM %d |VNLE',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + plot(ief.(sprintf('M%d', M)).baudr_new,ief.(sprintf('M%d', M)).ber_new,'DisplayName',sprintf('PAM %d |VNLE',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + end + + if 0 + if vnle + plot(xax,best_vnle,'DisplayName',sprintf('VNLE'),'Color',cols(2,:),'LineStyle','-','HandleVisibility','on'); + end + if mlse + plot(xax,best_vnle_2,'DisplayName',sprintf('VNLE'),'Color',cols(2,:),'LineStyle','-','HandleVisibility','on'); + plot(xax,best_mlse,'DisplayName',sprintf('VNLE + 1 tap post-filter + MLSE'),'Color',cols(3,:),'LineStyle','-','HandleVisibility','on'); + end + if db + plot(xax,best_db,'DisplayName',sprintf('DB tgt. + MLSE'),'Color',cols(1,:),'LineStyle','-','HandleVisibility','on'); + end + + plot(ief.(sprintf('M%d', M)).baudr,ief.(sprintf('M%d', M)).ber,'DisplayName',sprintf('VNLE [100,15,15] (ETH)'),'Color',cols(4,:),'LineStyle','-','HandleVisibility','on'); + plot(ief.(sprintf('M%d', M)).baudr_lowdsp,ief.(sprintf('M%d', M)).ber_lowdsp,'DisplayName',sprintf('VNLE [100,15,15] (ETH)'),'Color',cols(4,:),'LineStyle','-','HandleVisibility','on'); + end + + set(gca, 'YScale', 'log'); + ylim([1e-6 0.3]); + xlim([92, 260 ]); + xticks([0:16:280]); + yline([4.85e-3, 2e-2],'HandleVisibility','off'); + % legend + % beautifyBERplot() + xlabel('baudrate in GBd'); + ylabel('BER'); + set(findall(gca, '-property', 'Interpreter'), 'Interpreter', 'latex'); + + set(gcf, 'Color', 'w'); + set(gca, 'Box', 'on', 'LineWidth', 0.8); % Thicker border + grid on; + set(gca, 'FontSize', 10, 'FontName', 'Times New Roman'); + % + + %% NGMI PLOT + + % figure(2) + % cols = cbrewer2('Set1',6); + % hold on + % title(sprintf('%d km ; %.1f nm ; PAM %d',1,1313,M)); + % xax = baudrate.*1e-9; + % + % plot(xax,best_gmi_vnle./log2(M),'DisplayName',sprintf('NGMI VNLE'),'Color',cols(2,:),'LineStyle','-','HandleVisibility','on'); + % + % plot(ief.(sprintf('M%d', M)).baudr,ief.(sprintf('M%d', M)).ngmi,'DisplayName',sprintf('VNLE [100,15,15] (ETH)'),'Color',cols(4,:),'LineStyle','-','HandleVisibility','on'); + % plot(ief.(sprintf('M%d', M)).baudr_lowdsp,ief.(sprintf('M%d', M)).ngmi_lowdsp,'DisplayName',sprintf('VNLE [100,15,15] (ETH)'),'Color',cols(4,:),'LineStyle','-','HandleVisibility','on'); + % + % set(gca, 'YScale', 'log'); + % ylim([0.7 1]); + % xlim([min(xax), max(xax) ]); + % xticks(xax); + % yline([0.8],'HandleVisibility','off'); + % legend + % beautifyBERplot() + % xlabel('Baudrate in GBd'); + % ylabel('NGMI'); + + %% AIR Plot + % subplot(1,3,2) + % + % netRateMLSE = TransmissionPerformance().calculateNetRate(log2(M)*baudrate,'BER',best_mlse); + % netRateVNLE = TransmissionPerformance().calculateNetRate(log2(M)*baudrate,'NGMI',best_gmi_vnle./log2(M), 'BER',best_vnle_2); + % netRateDB = TransmissionPerformance().calculateNetRate(log2(M)*baudrate,'BER',best_db); + % + % netRateIEF = TransmissionPerformance().calculateNetRate(log2(M)*ief.(sprintf('M%d', M)).baudr,'BER',ief.(sprintf('M%d', M)).ber,'NGMI',ief.(sprintf('M%d', M)).ngmi); + % netRateIEF_lowdsp = TransmissionPerformance().calculateNetRate(log2(M)*ief.(sprintf('M%d', M)).baudr_lowdsp,'BER',ief.(sprintf('M%d', M)).ber_lowdsp,'NGMI',ief.(sprintf('M%d', M)).ngmi_lowdsp); + % + % cols = cbrewer2('Set1',8); + % hold on + % title(sprintf('AIR')); + % xax = baudrate.*1e-9; + % + % if M == 2 || M == 4 + % + % % plot(xax,netRateMLSE.HD.NetRate.*1e-9,'DisplayName',sprintf('MLSE - HD'),'Color',cols(1,:),'LineStyle',':','HandleVisibility','on'); + % % plot(xax,netRateVNLE.SDHD.NetRate.*1e-9,'DisplayName',sprintf('VNLE SD+HD'),'Color',cols(4,:),'LineStyle','-','HandleVisibility','on'); + % plot(xax,netRateDB.HD.NetRate.*1e-9,'DisplayName',sprintf('PAM %d | DB+MLSE HD',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on'); + % plot(xax,netRateDB.O_FEC.NetRate.*1e-9,'DisplayName',sprintf('PAM %d | DB+MLSE O-FEC',M),'Color',cols(M/2,:),'LineStyle',':','HandleVisibility','on'); + % % plot(xax,netRateDB.KP4_hamming.NetRate.*1e-9,'DisplayName',sprintf('PAM %d | DB+MLSE KP4-FEC',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on'); + % + % else + % + % % plot(xax,netRateMLSE.HD.NetRate.*1e-9,'DisplayName',sprintf('MLSE - HD'),'Color',cols(1,:),'LineStyle',':','HandleVisibility','on'); + % % plot(xax,netRateVNLE.SDHD.NetRate.*1e-9,'DisplayName',sprintf('VNLE SD+HD'),'Color',cols(4,:),'LineStyle','-','HandleVisibility','on'); + % % plot(ief.(sprintf('M%d', M)).baudr,netRateIEF.HD.NetRate,'DisplayName',sprintf('PAM %d | HD IEF',M),'Color',cols(M/2,:),'LineStyle',':','HandleVisibility','on'); + % plot(ief.(sprintf('M%d', M)).baudr,netRateIEF.SDHD.NetRate,'DisplayName',sprintf('PAM %d | SD+HD IEF',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on'); + % + % end + % set(gca, 'YScale', 'log'); + % ylim([92 450]); + % xlim([min(xax), max(xax) ]); + % xlim([92, 260 ]); + % xticks([0:16:280]); + % yline([0.8],'HandleVisibility','off'); + % legend + % beautifyBERplot() + % xlabel('Baudrate in GBd'); + % ylabel('Net Rate in Gbps'); + % + % + %% NDR PLOT -%%% PLOT COMPLETE BAUDRATE -figure(20) -cols = cbrewer2('Set1',6); -hold on -title(sprintf('%d km | %d nm | PAM %d',1,1313.5,M)); + subplot(1,2,2) + netRateIEF = TransmissionPerformance().calculateNetRate(log2(M)*ief.(sprintf('M%d', M)).baudr,'BER',ief.(sprintf('M%d', M)).ber,'NGMI',ief.(sprintf('M%d', M)).ngmi); + netRateIEF_lowdsp = TransmissionPerformance().calculateNetRate(log2(M)*ief.(sprintf('M%d', M)).baudr_lowdsp,'BER',ief.(sprintf('M%d', M)).ber_lowdsp,'NGMI',ief.(sprintf('M%d', M)).ngmi_lowdsp); + netRateIEF_new = TransmissionPerformance().calculateNetRate(log2(M)*ief.(sprintf('M%d', M)).baudr_new,'BER',ief.(sprintf('M%d', M)).ber_new,'NGMI',ief.(sprintf('M%d', M)).ngmi_new); -plot(baudrate,best_vnle,'DisplayName',sprintf('VNLE'),'Color',cols(1,:),'LineStyle','-','HandleVisibility','on'); -plot(baudrate,best_mlse,'DisplayName',sprintf('VNLE + 1 tap post-filter + MLSE'),'Color',cols(2,:),'LineStyle','-','HandleVisibility','on'); + cols = cbrewer2('Set1',8); + hold on + title(sprintf('Net Bitrate')); + xax = baudrate.*1e-9; -plot(baudrate,best_db,'DisplayName',sprintf('DB tgt. + MLSE'),'Color',cols(3,:),'LineStyle','-','HandleVisibility','on'); + if M == 2 || M == 4 -set(gca, 'YScale', 'log'); -ylim([5e-5 0.3]); -xlim([min(baudrate), max(baudrate) ]); -yline([3.8e-3, 2e-2],'HandleVisibility','off'); -legend -beautifyBERplot() -xlabel('Bit Rate in Gbps'); -ylabel('BER'); + log2M = log2(M); + + netRateDB = TransmissionPerformance().calculateNetRate(log2(M)*baudrate,'BER',best_db); + + % plot(xax,netRateMLSE.HD.NetRate.*1e-9,'DisplayName',sprintf('MLSE - HD'),'Color',cols(1,:),'LineStyle',':','HandleVisibility','on'); + % plot(xax,netRateVNLE.SDHD.NetRate.*1e-9,'DisplayName',sprintf('VNLE SD+HD'),'Color',cols(4,:),'LineStyle','-','HandleVisibility','on'); + plot(xax,netRateDB.HD.NetRate.*1e-9,'DisplayName',sprintf('PAM %d | DB+MLSE HD',M),'Color',cols(M/2,:),'LineStyle',':','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + plot(xax,netRateDB.O_FEC.NetRate.*1e-9,'DisplayName',sprintf('PAM %d | DB+MLSE O-FEC',M),'Color',cols(M/2,:),'LineStyle','-.','HandleVisibility','on','Marker','x','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + % plot(xax,netRateDB.KP4_hamming.NetRate.*1e-9,'DisplayName',sprintf('PAM %d | DB+MLSE KP4-FEC',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','*','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + plot(ief.(sprintf('M%d', M)).baudr_new,netRateIEF_new.SDHD.NetRate,'DisplayName',sprintf('PAM %d | SD+HD',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','o','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + plot(ief.(sprintf('M%d', M)).baudr_new,ief.(sprintf('M%d', M)).ngmi_new.*ief.(sprintf('M%d', M)).baudr_new.*log2M ,'DisplayName',sprintf('PAM %d | AIR',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','o','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + + elseif M == 6 + + log2M = 2.5; + netRateMLSE = TransmissionPerformance().calculateNetRate(log2(M)*baudrate,'BER',best_mlse); + netRateVNLE = TransmissionPerformance().calculateNetRate(log2(M)*baudrate,'NGMI',best_gmi_vnle./log2(M), 'BER',best_vnle_2); + + plot(xax,netRateVNLE.HD.NetRate.*1e-9,'DisplayName',sprintf('MLSE - HD'),'Color',cols(M/2,:),'LineStyle',':','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + plot(xax,netRateVNLE.SDHD.NetRate.*1e-9,'DisplayName',sprintf('VNLE SD+HD'),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','o','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + % plot(ief.(sprintf('M%d', M)).baudr_new,netRateIEF_new.HD.NetRate,'DisplayName',sprintf('PAM %d | HD IEF',M),'Color',cols(M/2,:),'LineStyle',':','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + % plot(ief.(sprintf('M%d', M)).baudr_new,netRateIEF_new.SDHD.NetRate,'DisplayName',sprintf('PAM %d | SD+HD IEF',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','o','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + + else + + log2M = log2(M); + + netRateMLSE = TransmissionPerformance().calculateNetRate(log2(M)*baudrate,'BER',best_mlse); + netRateVNLE = TransmissionPerformance().calculateNetRate(log2(M)*baudrate,'NGMI',best_gmi_vnle./log2(M), 'BER',best_vnle_2); + + % plot(xax,netRateVNLE.HD.NetRate.*1e-9,'DisplayName',sprintf('MLSE - HD'),'Color',cols(M/2,:),'LineStyle',':','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + % plot(xax,netRateVNLE.SDHD.NetRate.*1e-9,'DisplayName',sprintf('VNLE SD+HD'),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','o','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + plot(ief.(sprintf('M%d', M)).baudr_new,netRateIEF_new.HD.NetRate,'DisplayName',sprintf('PAM %d | HD IEF',M),'Color',cols(M/2,:),'LineStyle',':','HandleVisibility','on','Marker','square','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + plot(ief.(sprintf('M%d', M)).baudr_new,netRateIEF_new.SDHD.NetRate,'DisplayName',sprintf('PAM %d | SD+HD IEF',M),'Color',cols(M/2,:),'LineStyle','-','HandleVisibility','on','Marker','o','MarkerFaceColor',[1,1,1],'MarkerEdgeColor',cols(M/2,:),'MarkerSize',4,'LineWidth',1.4); + + end + set(gca, 'YScale', 'log'); + ylim([150 450]); + xlim([min(xax), max(xax) ]); + xlim([92, 270 ]); + xticks([0:16:280]); + yline([0.8],'HandleVisibility','off'); + % legend + + xlabel('Baudrate in GBd'); + ylabel('Net Rate in Gbps'); + set(findall(gca, '-property', 'Interpreter'), 'Interpreter', 'latex'); + + set(gcf, 'Color', 'w'); + set(gca, 'Box', 'on', 'LineWidth', 0.8); % Thicker border + grid on; + set(gca, 'FontSize', 10, 'FontName', 'Times New Roman'); + +end \ No newline at end of file diff --git a/projects/Messung_Zürich/dsp_ief_file.m b/projects/Messung_Zürich/dsp_ief_file.m index e39dad4..3fd0754 100644 --- a/projects/Messung_Zürich/dsp_ief_file.m +++ b/projects/Messung_Zürich/dsp_ief_file.m @@ -4,12 +4,14 @@ function [results] = dsp_ief_file(varargin) mu_ffe1 = 0.0001; mu_ffe2 = 0.0008; mu_ffe3 = 0.001; - mu_dc = 0.00; + mu_dc = 0.0003; %0.0003; - mu_dfe = 0.0004; - vnle_order1 = 100; + mu_dfe_training = 0.0004; + vnle_order1 = 50; vnle_order2 = 3; vnle_order3 = 3; + + dfe_mu = 0.0005; tcorrect = 0; @@ -33,10 +35,11 @@ function [results] = dsp_ief_file(varargin) end end - ssdPath = "D:\36_IMDD_Kiel\Data"; - % code = "20250221T004802"; - % code = "20250221T032043"; - filename = findFileByCode(ssdPath, filecode{1}); + pathToMeasurement = "C:\Users\Silas\Documents\MATLAB\Datensätze\IEF_Polariton_2025\36_IMDD_Kiel\Data\20250221"; + filename = findFileByCode(pathToMeasurement, filecode{1}); + + pathToTimingRecov = "C:\Users\Silas\Documents\MATLAB\Datensätze\IEF_Polariton_2025\36_IMDD_Kiel\Data\TR_ZIP"; + filename_2 = findFileByCode(pathToTimingRecov, filecode{1}); % Extract parameters from the filename using an updated regex tokens = regexp(filename, 'Pmod_([-0-9p]+)dBm_P_PD_([-0-9p]+)dBm_.*?_(\d+)GBd_(\d+)PAM__\d+T\d+', 'tokens'); @@ -75,20 +78,45 @@ function [results] = dsp_ief_file(varargin) dataRx = double(h5read(filepath, '/Waveforms/Channel 2/Channel 2Data')); % rohdaten des CH4 bitsTx = h5read(filepath, '/Settings/dataTx'); %Binär - bitsTx = reshape(bitsTx,log2(config.M),[])'; - + + if config.M ~= 6 + bitsTx = reshape(bitsTx,log2(config.M),[])'; + else + %bitsTx = reshape(bitsTx,5,[])'; + end + %%% Build Tx Signal (Bits, Pam Map, Symbols) Tx_bits = Informationsignal(bitsTx); Tx_symbols = PAMmapper(config.M,0,"eth_style",1).map(Tx_bits); Tx_symbols.fs = config.fsym; - + + if 0 + Rx_bits = PAMmapper(config.M,0,"eth_style",1).demap(Tx_symbols); + [~,~,ber_bw,~] = calc_ber(Rx_bits.signal,Tx_bits.signal(1:length(Rx_bits.signal)),"skip_front",100,"skip_end",150,"returnErrorLocation",1); + end + %%% Build Rx Signal (Rx, normalize,remove mean) - loadAfterTR = 0; + loadAfterTR = 1; if loadAfterTR - - rx_sig = load("testSilas.mat","signal_TR3"); - rx_sig=rx_sig.signal_TR3; - Rx_Sig_resamp = Informationsignal(rx_sig,"fs",config.fsym*2); + if 1 + rx_sig = load(filename_2); + rx_sig=rx_sig.signal_TR; + if config.M == 6 + rx_sig_pam6 = zeros(length(rx_sig)*2,1); + rx_sig_pam6(1:2:end) = real(rx_sig); + rx_sig_pam6(2:2:end) = imag(rx_sig); + rx_sig = rx_sig_pam6; + Rx_Sig_resamp = Informationsignal(rx_sig,"fs",config.fsym); + % Tx_symbols.signal = Tx_symbols.signal(1:end/2); + else + Rx_Sig_resamp = Informationsignal(rx_sig,"fs",config.fsym*2); + end + + else + rx_sig = load(string(['testSilas_',char(filecode{1}),'.mat']),"signal_TR3"); + rx_sig=rx_sig.signal_TR3; + Rx_Sig_resamp = Informationsignal(rx_sig,"fs",config.fsym*2); + end else @@ -98,39 +126,17 @@ function [results] = dsp_ief_file(varargin) Rx_Sig.signal = Rx_Sig.signal - mean(Rx_Sig.signal); Rx_Sig = Rx_Sig.normalize("mode","rms"); - - Rx_Sig.spectrum("fignum",2,"displayname",'Rx Signal','normalizeTo0dB',1); - - Rx_Sig = Filter('filtdegree',4,"f_cutoff",Tx_symbols.fs.*0.9,"fs",Rx_Sig.fs,"filterType",filtertypes.gaussian,"active",true).process(Rx_Sig); - - Rx_Sig.spectrum("fignum",2,"displayname",'Rx Signal filt','normalizeTo0dB',1); - - mf = Pulseformer("alpha",config.ROF,"fsym",config.fsym,"fdac",Rx_Sig.fs,"matched",1,"pulse","rrc","pulselength",32); - Rx_matched = mf.process(Rx_Sig); - - Rx_matched.spectrum("fignum",2,"displayname",'Rx Signal m','normalizeTo0dB',1); %%%%%% Sample to 2x fsym %%%%%% - Rx_Sig_resamp = Rx_Sig.resample("fs_out",2*config.fsym); + Rx_Sig_resamp = Rx_Sig.resample("fs_out",config.fsym); end - + + %%%%%% Sync Rx signal with reference (S is a cell array with all occurences) %%%%%% - [Rx_Sig_sync,S,isFlipped] = Rx_Sig_resamp.tsynch("reference",Tx_symbols,"fs_ref",config.fsym,"debug_plots",1); - - % Rx_Sig_sync.eye(fsym,M,"fignum",4,"displayname",'eye diagram'); - % timing syncroization ?? + [~,S,isFlipped] = Rx_Sig_resamp.tsynch("reference",Tx_symbols,"fs_ref",config.fsym,"debug_plots",1); - % SYMSYNC = comm.SymbolSynchronizer("Modulation","PAM/PSK/QAM","SamplesPerSymbol",2); - % - % synd = SYMSYNC(Rx_Sig_sync.signal); - % Rx_Sig_sync.signal = synd; - % - % [Rx_Sig_sync,S,isFlipped] = Rx_Sig_sync.tsynch("reference",Tx_symbols,"fs_ref",config.fsym,"debug_plots",1); - % - - % filter Bw output = struct(); vnle_package = {}; vnle_pf_package = {}; @@ -139,78 +145,44 @@ function [results] = dsp_ief_file(varargin) for s = 1%:length(S) Rx_Sig_sync = S{s}; Rx_Sig_sync = Rx_Sig_sync.normalize("mode","rms"); - - if 1 - Tx_symbols.spectrum("fignum",3,"displayname",'Tx Symbols','normalizeTo0dB',1); - Rx_Sig_sync.spectrum("fignum",3,"displayname",'No Matched Filter','normalizeTo0dB',1); - - Nsym = 64; - sampsPerSym = 2; - rcrFilt = comm.RaisedCosineReceiveFilter(... - 'Shape', 'Square root', ... - 'RolloffFactor', config.ROF, ... - 'FilterSpanInSymbols', Nsym, ... - 'InputSamplesPerSymbol', sampsPerSym, ... - 'DecimationFactor', 1); - yr = rcrFilt([Rx_Sig_sync.signal;zeros(Nsym*sampsPerSym/2, 1)]); - fltDelay = Nsym / (2*config.fsym); - yr = yr(fltDelay*Rx_Sig_sync.fs+1:end); - Rx_matched = Rx_Sig_sync; - Rx_matched.signal = yr; - - % Rx_matched = Rx_matched.resample("fs_out",2*config.fsym); - % length(Rx_matched); - % symbolSync = comm.SymbolSynchronizer("Modulation","PAM/PSK/QAM","SamplesPerSymbol",2,"TimingErrorDetector","Mueller-Muller (decision-directed)"); - % Rx_syncd = Rx_matched; - % [Rx_syncd.signal,tError] = symbolSync(Rx_syncd.signal); - % Rx_syncd.fs = config.fsym; - % length(Rx_syncd); - - % figure();hold on; - % stem(Tx_symbols.normalize("mode","oneone").signal); - % stem(Rx_syncd.normalize("mode","oneone").signal); - % - % figure();hold on; - % plot(Rx_Sig_sync.normalize("mode","oneone").signal); - % plot(Rx_matched.normalize("mode","oneone").signal); - - - Rx_syncd = Rx_syncd.resample("fs_out",2*config.fsym); - - Rx_matched.spectrum("fignum",3,"displayname",'After Matched Filter','normalizeTo0dB',1); - - Rx_syncd.spectrum("fignum",3,"displayname",'Zero Crossing TR','normalizeTo0dB',1); - - end + Rx_Sig_sync = Rx_Sig_sync.resample("fs_out",2*config.fsym); mu_ffe = [mu_ffe1 mu_ffe3 mu_ffe3]; vnle_order=[vnle_order1,vnle_order2,vnle_order3]; - eq_ = EQ("Ne",vnle_order,"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005,"FFEmu",0,"plotfinal",1,"ideal_dfe",1); + eq_ = EQ("Ne",vnle_order,"Nb",[0,0,0],"training_length",4096*4,"training_loops",4,"dd_loops",3,"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe_training],"DFEmu",dfe_mu,"FFEmu",0,"plotfinal",0,"ideal_dfe",0,"plottrain",0); %%%%% VNLE only (or DFE) %%%% - if 1 + if 0 - - [result] = vnle(eq_,config.M,Rx_syncd,Tx_symbols,Tx_bits,"precode_mode",db_mode.no_db,"showAnalysis",1,'eth_style',1); + eq_post = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",1001,"sps",1,"decide",0); + [result] = vnle(eq_,config.M,Rx_Sig_sync,Tx_symbols,Tx_bits,"precode_mode",db_mode.no_db,"showAnalysis",1,'eth_style',1,'postFFE',eq_post); + + netRate = TransmissionPerformance().calculateNetRate(log2(config.M)*config.fsym,'NGMI',result.inf_rate_vnle, 'BER',result.ber_vnle); vnle_package{s} = result; end %%%%% VNLE + PF + MLSE %%%% - if 0 + if 1 + eq_post = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",1001,"sps",1,"decide",0); pf_ = Postfilter("ncoeff",1,"useBurg",1); + mlse_ = MLSE("duobinary_output",0,'M',config.M,'trellis_states',PAMmapper(config.M,0).levels); mlse_ = MLSE_viterbi("duobinary_output",0,'M',config.M,'trellis_states',PAMmapper(config.M,0).levels); - doub_mode = db_mode.no_db; + doub_mode = db_mode.no_db; - [result] = vnle_postfilter_mlse(eq_,pf_,mlse_,config.M,Rx_Sig_sync,Tx_symbols,Tx_bits,"precode_mode",doub_mode,'showAnalysis',1,'eth_style_symbol_mapping',1); + [result] = vnle_postfilter_mlse(eq_,pf_,mlse_,config.M,Rx_Sig_sync,Tx_symbols,Tx_bits,"precode_mode",doub_mode,'showAnalysis',1,'eth_style_symbol_mapping',1,'postFFE',eq_post); + netRate = TransmissionPerformance().calculateNetRate(log2(config.M)*config.fsym,'NGMI',result.gmi./log2(config.M), 'BER',result.ber_mlse); + fprintf('VNLE SD: %.1f GBd \n',netRate.SDHD.NetRate.*1e-9); + fprintf('MLSE HD: %.1f GBd \n',netRate.HD.NetRate.*1e-9); + vnle_pf_package{s} = result; end %%%%% Duobinary Targeting %%%% - if 0 + if 1 mlse_db = MLSE_viterbi("DIR",[1,1],"duobinary_output",0,"M",config.M,"trellis_states",PAMmapper(config.M,0).levels); doub_mode = db_mode.db_emulate; @@ -225,6 +197,7 @@ function [results] = dsp_ief_file(varargin) results.vnle_package = vnle_package; results.vnle_pf_package = vnle_pf_package; results.dbtgt_package = dbtgt_package; + results.config = config; end diff --git a/projects/Messung_Zürich/minimal_pam4_mlse.m b/projects/Messung_Zürich/minimal_pam4_mlse.m new file mode 100644 index 0000000..1a2bb54 --- /dev/null +++ b/projects/Messung_Zürich/minimal_pam4_mlse.m @@ -0,0 +1,19 @@ + + +mlse_sig_sd=load("imdd_simulation\projects\Messung_Zürich\mlse_sig_sd.mat","mlse_sig_sd"); +mlse_sig_sd = mlse_sig_sd.mlse_sig_sd; +tx_symbols=load("imdd_simulation\projects\Messung_Zürich\tx_symbols.mat","tx_symbols"); +tx_symbols = tx_symbols.tx_symbols; + +mlse_ = MLSE_viterbi("duobinary_output",0,'M',4,'trellis_states',PAMmapper(4,0).levels); +mlse_.DIR = [1.0000 0.5452]; +mlse_sig_sd = mlse_.process(mlse_sig_sd); + + +% +% mlse_ = MLSE_viterbi("duobinary_output",0,'M',4,'trellis_states',PAMmapper(4,0).levels); +% mlse_.DIR = [1.0000 0.5452]; +% mlse_sig_sd = mlse_.process(mlse_sig_sd,tx_symbols); + +h = [1.0000 0.5452]; +chatgpt_answer(mlse_sig_sd.signal,tx_symbols.signal,h) diff --git a/test/matched_filter_minimal.m b/test/matched_filter_minimal.m index 9a5d413..b298fac 100644 --- a/test/matched_filter_minimal.m +++ b/test/matched_filter_minimal.m @@ -26,24 +26,34 @@ tx_bw_nyquist = 0.8; O = 18; %order of prbs N = 2^(O-1); %length of prbs -rcalpha = 1; -Pform = Pulseformer("fsym",fsym,"fdac",fdac,"pulse","rrc","pulselength",1024,"alpha",rcalpha); +%%%%% MOVE-IT PRMS %%%% +Mi_prms = Moveit_wrapper("prms"); +if M == 6 + Mi_prms.para.bl = 2^(O-2); + Mi_prms.para.dimension = 5; +else + Mi_prms.para.bl = 2^(O-1); + Mi_prms.para.dimension = log2(M); %2.5bits/sym -> 2 bit/sym +end +Mi_prms.para.rand = 0; +Mi_prms.para.order = floor(O / log2(M)); +Mi_prms.para.skip =0; +Mi_prms.para.bruijn = 0; +Mi_prms.para.reset_prms = 0; +Mi_prms.para.method = 1; +bitpattern = Mi_prms.process([]); +if M == 6 + bitpattern = reshape(bitpattern',[],1); + bitpattern = bitpattern(1:end-mod(length(bitpattern),5)); +end +bits = Informationsignal(bitpattern); +symbols = PAMmapper(M,0).map(bits); +symbols.fs = fsym; -[Digi_sig,Symbols,Tx_bits] = PAMsource(... - "fsym",fsym,"M",M,"order",18,"useprbs",0,... - "fs_out",fdac,... - "applyclipping",0,"clipfactor",1.5,... - "applypulseform",apply_pulsef,"pulseformer",Pform,... - "randkey",random_key,... - "db_precode",db_precode,"db_encode",db_encode,... - "mrds_code",0,"mrds_blocklength",512).process(); - -Digi_sig.spectrum("displayname",'Signal after shaping','fignum',1); -% Digi_sig.move_it_spectrum("displayname",'Signal after shaping','fignum',2); - -% Digi_sig = Digi_sig.resample("fs_out",fsym); +Pform = Pulseformer("fsym",fsym,"fs",fsym,"alpha",0.6,"pulse","rrc","pulselength",64,"matched_sps",4,"output_sps",2); +Pform.process(symbols); MF = Pulseformer("fsym",fsym,"fdac",fdac,"pulse","rrc","pulselength",1024,"alpha",rcalpha,"matched",0); Digi_sig = MF.process(Digi_sig); diff --git a/test/test_mapping_eth.m b/test/test_mapping_eth.m index c3e1cc9..592d9db 100644 --- a/test/test_mapping_eth.m +++ b/test/test_mapping_eth.m @@ -1,21 +1,32 @@ % Setup PRBS parameters -O = 12; -M = 2; +O = 6; +M = 6; N = 2^(O-1); % Length of PRBS randkey = 1; % Random key for random stream +use_eth_mapping =1; + +if M ~= 6 + dimension = log2(M); +else + dimension = 5; +end [~, seed] = prbs(O, 1); % Initialize first seed of PRBS bitpattern = []; s = RandStream('twister', 'Seed', randkey); -for i = 1:log2(M) +for i = 1:dimension bitpattern(:, i) = randi(s, [0 1], N, 1); end +if M == 6 + bitpattern = reshape(bitpattern',[],1); + bitpattern = bitpattern(1:end-mod(length(bitpattern),5)); +end Tx_bits = Informationsignal(bitpattern); -Digi_Mod = PAMmapper(M, 0,"eth_style",1); +Digi_Mod = PAMmapper(M, 0,"eth_style",use_eth_mapping); % Map bits to symbols Symbols = Digi_Mod.map(Tx_bits); @@ -23,7 +34,6 @@ Symbols = Digi_Mod.map(Tx_bits); % Demap symbols back to bits Rx_bits = Digi_Mod.demap(Symbols); -[~, error_num, ber, ~] = calc_ber(Tx_bits.signal, Rx_bits.signal, ... - "skip_front", 0, "skip_end", 0, "returnErrorLocation", 1); +[~, error_num, ber, ~] = calc_ber(Tx_bits.signal(1:length(Rx_bits.signal)), Rx_bits.signal,"skip_front", 0, "skip_end", 0, "returnErrorLocation", 1); -fprintf('BER VNLE: %.1E \n',ber); \ No newline at end of file +fprintf('BER: %.1E \n',ber); \ No newline at end of file