Merge branch 'main' of https://cau-git.rz.uni-kiel.de/nt/mitarbeiter/silas/imdd_simulation

Classes/00_signals/Signal.m

@@ -344,13 +344,14 @@ classdef Signal
 
             % spectrum_plot(obj.signal,options.fsamp,options.figurename,options.displayname);
 
-            N = 2^(nextpow2(length(obj.signal))-8);
+            N = 2^(nextpow2(length(obj.signal))-10);
 
             if options.normalizeToNyquist==0
                 [p_lin,w] = pwelch(obj.signal,hanning(N),N/2,N,obj.fs,"centered","power","mean");
                 w=w.*1e-9;
             else
                 [p_lin,w] = pwelch(obj.signal,hanning(N),N/2,N,"centered","power","mean");
+                p_lin = smooth(p_lin,0.05,'rloess');
             end
 
             

projects/Lab_2024/offline_dsp_analysis/load_n_dsp.m

@@ -1,10 +1,13 @@
 
 % load data points
-
+if 0
-foldername = 'C:\Users\Silas\Nextcloud4\Dokumente\02_Ablage_Office\Lab_Data_24\sir_sweep_pam4';
+    foldername = '/Users/silasoettinghaus/Nextcloud/Dokumente/02_Ablage_Office/Lab_Data_24/sir_sweep_sd40/';
-filename = 'PAM4_10km_';
+    filename = 'PAM4_56_';
-
+    filename = 'PAM4_56_v2';
-
+else
+    foldername = '/Users/silasoettinghaus/Nextcloud/Dokumente/02_Ablage_Office/Lab_Data_24/sir_sweep_sd40/DB_coded';
+    filename = 'PAM4_68_v2';
+end
 stuff = load([foldername,filesep,filename,'wh']);
 %stuff = load([foldername,filesep,'PAM4_v2_10km_wh']);
 wh = stuff.obj;
@@ -19,9 +22,9 @@ M = 4;
 % Tx Symbols
 % Rx Signal
 
-Bits = load([foldername, filesep, filename, 'bits'],'Bits');
+Bits = load([foldername, filesep, filename, 'bits_iatten_0'],'Bits');
 Bits = Bits.Bits;
-Symbols = load([foldername, filesep, filename, 'symbols'],'Symbols');
+Symbols = load([foldername, filesep, filename, 'symbols_iatten_0'],'Symbols');
 Symbols = Symbols.Symbols;
 
 % all_signals = load([foldername, filesep, filename, 'best_rx_signals']);
@@ -29,101 +32,125 @@ Symbols = Symbols.Symbols;
 
 Eq = EQ("Ne",[50,7,7],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.0,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
 
-cnt = 0;
+Eq = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",50,"sps",2,"decide",0);
-for atten = 40:-4:0
-    cnt = cnt+1;
 
+ber_mlse_ = zeros(11,2,3);
+ber_db_ = zeros(11,2,1);
+ber_vnle_ = zeros(11,2,1);
+
+cnt = 0;
+attens = 40:-4:0;
+for cnt = 10:numel(attens)
+    atten = attens(cnt);
+    ber_mlse = [];
+    ber_vnle = [];
+    ber_db = [];
     Signal_cell = load([foldername, filesep, filename, 'rx_signal_iatten_',num2str(atten)]);
     
     % Scpe_sig = all_signals(cnt);
 
-    ber_ffe_only_ = zeros(numel(Signal_cell.S),1);
+    for i = 1:2
-    ber_mlse_ = zeros(numel(Signal_cell.S),1);
+        cols = cbrewer2('Paired',8);
-    ber_db_ = zeros(numel(Signal_cell.S),1);
-    
-    parfor i = 1:numel(Signal_cell.S)
-        
         Scpe_sig = Signal_cell.S{i};
 
         sir_measured = wh.getStoValue('sir',v_bias,awg_vpp,eqmode,atten);
 
-        if eqmode == 2
+        if eqmode == 2 %postfilter
-            [EQ_sig] = Eq.process(Scpe_sig,Symbols);
 
-            EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1);
+            [EQ_sig, Noi] = Eq.process(Scpe_sig,Symbols);
-        
-            Noi = EQ_sig-Symbols;
 
             Rx_bits = PAMmapper(M,0).demap(EQ_sig);
-            [~,errors_bm,ber_ffe_only_(i),errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+            [~,~,ber_vnle(i),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
 
-            nc = 2;
-            burg_coeff = arburg(Noi.signal,nc);
 
-            EQ_sig = EQ_sig.filter(burg_coeff,1);
 
-            if 1
+            for nc = 1:3
-                Noi.spectrum('displayname',['SIR: ', num2str(sir_measured), ' dB '],'fignum',123)
-                [h,w] = freqz(1,burg_coeff,length(Noi),"whole",Noi.fs);
-                h = h/max(abs(h));
-                hold on
-                w_ = (w - Noi.fs/2);
-                plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']);
-            end
 
-            if 0
+                burg_coeff = arburg(Noi.signal,nc);
-                figure(53);
+
-                clf
+                EQ_sig_filt = EQ_sig.filter(burg_coeff,1);
-                constellation = unique(Symbols.signal);
+
-                received = NaN(numel(constellation),length(Symbols));
+                % EQ_sig.spectrum("displayname","Signal Spectrum after Postfilter","fignum",1234);
-                for lvl = 1:numel(constellation)
+
-                    received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig.signal(Symbols.signal==constellation(lvl));
+                EQ_sig_mlse = MLSE("DIR",burg_coeff,"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig_filt);
+
+                % EQ_sig.spectrum("displayname","Signal Spectrum after MLSE","fignum",1234);
+
+                if nc == 3
+                    c = 6;
+
+                    EQ_sig.spectrum('displayname',['FFE out',num2str(sir_measured)],'fignum',11,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(c,:));
+                    
+                    Noi.spectrum('displayname',['Noise',num2str(sir_measured)],'fignum',12,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(c,:));
+                    
+                    EQ_sig_filt.spectrum('displayname',['PF ',num2str(sir_measured)],'fignum',13,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(c,:));
+
+                    EQ_sig_mlse.spectrum('displayname',['MLSE out ',num2str(sir_measured)],'fignum',14,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(c,:));
+
+                    [h,w] = freqz(1,burg_coeff,length(Noi),"whole");
+                    h = 1./h;
+                    h = h/max(abs(h));
                     hold on
-                    histogram(received(lvl,:),1000,"EdgeAlpha",0);
+                    w_ = (w - pi);
+                    figure(13)
+                    plot(w_,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']);
+
+                    Scpe_sig.plot("displayname","FFE out","fignum",12);
+
                 end
+
+                Rx_bits = PAMmapper(M,0).demap(EQ_sig_mlse);
+
+                [~,errors_bm,ber_mlse(i,nc),errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+               
+
             end
 
-            EQ_sig = MLSE("DIR",burg_coeff,"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig);
+            % disp(['FFE: ',sprintf('%.1E',ber_vnle(cnt,i)),' -> PF -> MLSE: ',sprintf('%.1E',ber_mlse(i,nc)),'  dB']);
 
-            Rx_bits = PAMmapper(M,0).demap(EQ_sig);
-            [~,errors_bm,ber_mlse_(i),errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
-        
-            disp(['FFE: ',sprintf('%.1E',ber_ffe_only_(i)),' -> PF -> MLSE: ',sprintf('%.1E',ber_mlse_(i)),'  dB ']);
         elseif eqmode == 3
 
 
             [EQ_sig, Noi] = Eq.process(Scpe_sig,Duobinary().encode(Symbols));
 
-            % EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1);
+
+
+            c = 6;
+
+            EQ_sig.spectrum('displayname',['DB out',num2str(sir_measured)],'fignum',11,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(c,:));
+
+            Noi.spectrum('displayname',['DB Noise',num2str(sir_measured)],'fignum',12,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(c,:));
+
+
 
             EQ_sig = MLSE("DIR",[1,1],"duobinary_output",1,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig);
             EQ_sig = Duobinary().decode(EQ_sig);
 
             Rx_bits = PAMmapper(M,0).demap(EQ_sig);
-            [~,errors_bm,ber_db_(i),errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+            [~,errors_bm,ber_db(i),errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
 
-            disp([' DB Precode -> Channel -> FFE -> Decode/ Mod  ',sprintf('%.1E',ber_db_(i)),' ']);
+            disp([' DB Precode -> Channel -> FFE -> Decode/ Mod  ',sprintf('%.1E',ber_db(i)),' ']);
 
 
 
         end
 
+
     end
 
     if eqmode == 2
-        [ber_ffe_2(cnt),idx] = min(ber_ffe_only_);
+        ber_mlse_(cnt,:,:) = ber_mlse(1:size(ber_mlse_,2),:);
-        [ber_mlse_2(cnt),idx] = min(ber_mlse_);
+        ber_vnle_(cnt,:,:) = ber_vnle(1:size(ber_mlse_,2));
-
     elseif eqmode == 3
-        [ber_db(cnt),idx] = min(ber_db_);
+        ber_db_(cnt,:,:) = ber_db;
     end
     
-    all_signals(:,cnt) = Signal_cell.S{idx};
-
 end
 
 
-
+ber_vnle_min = flip(squeeze(min(ber_vnle_,[],2)));
+ber_mlse_min = flip(squeeze(min(ber_mlse_,[],2)));
+ber_db_min = flip(min(ber_db_,[],2));
 
 cols = linspecer(8);
 
@@ -140,13 +167,14 @@ hold on; % Retain the plot so new points can be added without complete redraw
 
 % Plot the data and get the line handle
 if eqmode == 2
-    hLine1 = plot(sirs, ber_ffe_2, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'VNLE 2');
+    plot(sirs, ber_vnle_min, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'VNLE');
-    hLine2 = plot(sirs, ber_mlse_2, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'VNLE+MLSE 2');
+    plot(sirs, ber_mlse_min(:,1), "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'MLSE 1');
- elseif eqmode == 3
+    plot(sirs, ber_mlse_min(:,2), "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'MLSE 2');
-    hLine2 = plot(sirs, ber_db, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'DB ');
+    plot(sirs, ber_mlse_min(:,3), "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'MLSE 3');
+elseif eqmode == 3
+    plot(sirs, ber_db_min, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'DB ');
 end
 
-hLine3 = plot(sirs, bers, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'from measurement');
 
 
 % Continue with the rest of your plot settings

projects/Lab_2024/offline_dsp_analysis/process_vnle_mlse.m

@@ -0,0 +1,115 @@
+
+
+foldername = '/Users/silasoettinghaus/Nextcloud/Dokumente/02_Ablage_Office/Lab_Data_24/sir_sweep_sd40/';
+filename = 'PAM4_56_';
+
+M = 4;
+
+Bits = load([foldername, filesep, filename, 'bits_iatten_0'],'Bits');
+Bits = Bits.Bits;
+
+Symbols = load([foldername, filesep, filename, 'symbols_iatten_0'],'Symbols');
+Symbols = Symbols.Symbols;
+
+% all_signals = load([foldername, filesep, filename, 'best_rx_signals']);
+% all_signals = all_signals.all_signals;
+
+stuff = load([foldername,filesep,filename,'wh']);
+wh = stuff.obj;
+i_atten_vals = wh.parameter.i_atten.values;
+v_bias = wh.parameter.vbias.values(1);
+awg_vpp = wh.parameter.awg_vpp.values(1);
+eqmode = wh.parameter.eq_mode.values(1);
+
+
+ber_vnle = zeros(1,numel(all_signals));
+ber_mlse = zeros(numel(all_signals),3);
+
+Eq = EQ("Ne",[25,0,0],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.0,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.00,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
+
+Eq  = FFE_adaptive_decision("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",25,"sps",2,"decide",1,"buffer_length",80);
+
+Eq     = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",25,"sps",2,"decide",1);
+
+col = linspecer(11);
+atten = 40:-4:0;
+for i = 1:numel(atten)
+
+    sir_measured = wh.getStoValue('sir',v_bias,awg_vpp,eqmode,atten(i));
+    
+    Scpe_sig = all_signals(i);
+
+    [EQ_sig, Noi] = Eq.process(Scpe_sig,Symbols);
+
+    Noi.spectrum('displayname',['SIR: ', num2str(sir_measured), ' dB '],'fignum',12345,'color',col(i,:))   
+
+    % EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1);
+
+    Rx_bits = PAMmapper(M,0).demap(EQ_sig);
+
+    [~,~,ber_vnle(i),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+
+    % for nc = 1:3
+    % 
+    %     burg_coeff = arburg(Noi.signal,nc);
+    % 
+    %     EQ_sig_mlse = EQ_sig.filter(burg_coeff,1);
+    % 
+    % 
+    %     if 1
+    %         Noi.spectrum('displayname',['SIR: ', num2str(sir_measured), ' dB '],'fignum',123)
+    %         [h,w] = freqz(1,burg_coeff,length(Noi),"whole",Noi.fs);
+    %         h = h/max(abs(h));
+    %         hold on
+    %         w_ = (w - Noi.fs/2);
+    %         plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']);
+    %     end
+    % 
+    %     if 0
+    %         figure(53);
+    %         clf
+    %         constellation = unique(Symbols.signal);
+    %         received = NaN(numel(constellation),length(Symbols));
+    %         for lvl = 1:numel(constellation)
+    %             received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig.signal(Symbols.signal==constellation(lvl));
+    %             hold on
+    %             histogram(received(lvl,:),1000,"EdgeAlpha",0);
+    %         end
+    %     end
+    % 
+    %     EQ_sig_mlse = MLSE("DIR",burg_coeff,"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig_mlse);
+    % 
+    %     Rx_bits = PAMmapper(M,0).demap(EQ_sig_mlse);
+    %     [~,~,ber_mlse(i,nc),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+    % 
+    %     disp(['FFE: ',sprintf('%.1E',ber_vnle(i)),' -> PF -> MLSE: ',sprintf('%.1E',ber_mlse(i,nc)),'  dB ']);
+    % end
+
+
+end
+
+
+
+bers = wh.getStoValue('ber',v_bias,awg_vpp,eqmode,i_atten_vals);
+sirs = wh.getStoValue('sir',v_bias,awg_vpp,eqmode,i_atten_vals);
+
+figure(90);
+hold on; % Retain the plot so new points can be added without complete redraw
+
+% Plot the data and get the line handle
+hLine1 = plot(sirs, (ber_vnle), "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'VNLE 2');
+hLine2 = plot(sirs, (ber_mlse(:,1)), "LineWidth", 0.5, "LineStyle", "-", "Marker", "o", "MarkerSize", 15, "DisplayName", ['VNLE + PF ',num2str(1),' taps + MLSE']);
+hLine2 = plot(sirs, (ber_mlse(:,2)), "LineWidth", 0.5, "LineStyle", "-", "Marker", "x", "MarkerSize", 15, "DisplayName", ['VNLE + PF ',num2str(2),' taps + MLSE']);
+hLine2 = plot(sirs, (ber_mlse(:,3)), "LineWidth", 0.5, "LineStyle", "-", "Marker", "square", "MarkerSize", 15, "DisplayName", ['VNLE + PF ',num2str(3),' taps + MLSE']);
+% hLine3 = plot(sirs, bers, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", 'from measurement');
+
+% Continue with the rest of your plot settings
+yline(3.8e-3, 'DisplayName', 'HD-FEC', 'LineStyle', '--', 'HandleVisibility', 'off');
+xlabel('Signal to Interference Ratio in dB');
+ylabel('Bit Error Rate (BER)');
+title('Bit Error Rate vs. SIR');
+set(gca, 'yscale', 'log');
+set(gca, 'Box', 'on');
+grid on;
+grid minor;
+legend('Interpreter', 'none');