Last commit after 400G Lab Measurement

Classes/04_DSP/Coding/Duobinary.m

@@ -125,6 +125,13 @@ classdef Duobinary
             %make bipolar
             data = (data-round(mean(data),1));
 
+            [unique_points, ~, idx] = unique(data);
+            counts = accumarray(idx, 1);
+            total_samples = numel(data);
+            probabilities = counts / total_samples;
+            mean_power = sum((unique_points .^ 2) .* probabilities);
+            scaling_factor = sqrt(mean_power);
+
             if M == 4
                 data = data ./ sqrt(2.5); % 7-level constellation weighted with probability after DB code i.e. mean([-3 3 -2 -2 2 2 -1 -1 -1 1 1 1 0 0 0 0].^2) = 2.5 --> sqrt(2.5) == rms(constellation)
             elseif M == 6

Classes/04_DSP/Equalizer/EQ.m

@@ -1,4 +1,4 @@
-classdef EQ
+classdef EQ %< handle
     %EQ Summary of this class goes here
     %   Detailed explanation goes here
 

Classes/Warehouse_class/classes/DataStorage.m

@@ -129,7 +129,7 @@ classdef DataStorage < handle
                         tmp = obj.sto.(storageVarName){lin_idx(i)};
                         if ~isempty(tmp)
 
-                            if isa(tmp,'Signal') || isa(tmp,'struct') ||  isa(tmp,'Exfo_laser')
+                            if isa(tmp,'Signal') || isa(tmp,'struct') ||  isa(tmp,'Exfo_laser') ||  isa(tmp,'DC_supply')
                                 if i == 1
                                     value = {};
                                 end

Functions/showCurrentMeasurement.m

@@ -41,7 +41,7 @@ function showCurrentMeasurement(varargin)
     for i = 1:length(values)
         varNameLower = names{i};  % Convert variable name to lowercase for case-insensitive comparison
         if isnumeric(values{i})
-            if any(contains(varNameLower, {'ber','mlse','ffe','db'}))
+            if any(contains(varNameLower, {'ber','mlse','ffe','db'},"IgnoreCase",true))
                 % Format 'ber' values in exponential notation with two decimal places
                 values{i} = sprintf('%.2e', values{i});
             else

Tests/04_DSP/Coding/Duobinary_test.m

@@ -12,9 +12,9 @@ classdef Duobinary_test < matlab.unittest.TestCase
     properties (MethodSetupParameter)
         % Define method-level parameters for PRBS and bit pattern
         % i.e.: useprbs = {0,1};
-        useprbs = struct('true', 1,'false', 0);  % variations: {0, 1}
+        useprbs = struct('true', 1);  % variations: {0, 1}
         M = struct('M2', 2, 'M4', 4, 'M6', 6, 'M8', 8);  % variations: {2, 4, 6, 8}
-        O = struct('O10', 10,'O11', 11,'O12', 12,'O13', 13, 'O15', 15, 'O17', 17);  % variations: {10, 15, 17}
+        O = struct('O10', 10,'O11', 11,'O12', 12,'O13', 13, 'O15', 15, 'O17', 17, 'O18', 18, 'O19', 19);  % variations: {10, 15, 17}
     end
 
     properties (TestParameter)
@@ -31,14 +31,42 @@ classdef Duobinary_test < matlab.unittest.TestCase
             fsym = round(datarate / log2(M)) ;
 
             %%%%% PRBS Generation in correct shape for Modulation Format %%%%%%
-            N = 2^O; %length of prbs
+            N = 2^(O-1); %length of prbs
             [~,seed] = prbs(O,1); %initialize first seed of prbs
             bitpattern=[];
 
             if useprbs
-                for i = 1:log2(M)
+                 %%%%% MOVE-IT PRMS %%%%
-                    [bitpattern(:,i),seed] = prbs(O,N,seed);
+
+                state = struct();
+
+                para = struct();
+
+                if M == 6
+                    para.bl = 2^(O-2);
+                    para.dimension = 5;
+                else
+                    para.bl = 2^(O-1);
+                    para.dimension = log2(M); %2.5bits/sym -> 2 bit/sym
                 end
+
+                para.rand = 0;
+
+                para.order = floor(O / log2(M));
+                para.skip =0;
+                para.bruijn = 0;
+                para.reset_prms = 0;
+                para.method = 1;
+
+                data_in = [];
+                global loop;
+                loop = 0;
+                [data_out,state_] = prms(data_in, state, para);
+                loop = 1;
+                [data_out,state_out] = prms(data_in, state_, para);
+                bitpattern = data_out';
+
+                %%%%% END MOVE-IT %%%%%%%
             else
                 s = RandStream('twister','Seed',randkey);
                 for i = 1:log2(M)
@@ -46,8 +74,9 @@ classdef Duobinary_test < matlab.unittest.TestCase
                 end
             end
 
+
             if M == 6
-                bitpattern = reshape(bitpattern,[],1);
+                bitpattern = reshape(bitpattern',[],1);
                 bitpattern = bitpattern(1:end-mod(length(bitpattern),5));
             end
 

projects/HighSpeedExperiment_2024/master_auswertung_10km.m

@@ -0,0 +1,80 @@
+
+wh = load('C:\Users\sioe\Documents\High_Speed_Measurement_2024\10km_bitrate_complete\20241030_170224_wh.mat');
+wh = wh.obj;
+
+M_vals = wh.parameter.M.values;
+M_choose = M_vals(1);
+lambda_vals = wh.parameter.lambda.values;
+bitrate_vals = wh.parameter.bitrate.values;
+duobinary_vals = wh.parameter.duobinary.values;
+rop_atten_vals = wh.parameter.rop_atten.values;
+
+
+figure(177)
+
+
+for M_choose = [8]
+    sgtitle(['PAM',num2str(M_choose)])
+    for l = 1:numel(lambda_vals)
+        for b = 1:numel(bitrate_vals)
+
+            cel = wh.getStoValue('ber_vnle',M_choose(1),lambda_vals(l),bitrate_vals(b),duobinary_vals(1),rop_atten_vals(1));
+            ber_vnle(b)=min(cel{1});
+
+            cel = wh.getStoValue('ber_vnle_mlse',M_choose(1),lambda_vals(l),bitrate_vals(b),duobinary_vals(1),rop_atten_vals(1));
+            ber_vnle_mlse(b)=min(cel{1});
+
+            cel = wh.getStoValue('ber_db',M_choose(1),lambda_vals(l),bitrate_vals(b),duobinary_vals(2),rop_atten_vals(1));
+            ber_db(b)=min(cel{1});
+
+            cel = wh.getStoValue('ber_db',M_choose(1),lambda_vals(l),bitrate_vals(b),duobinary_vals(3),rop_atten_vals(1));
+            ber_db_enc(b)=min(cel{1});
+
+            dcs_ = wh.getStoValue('dcs',M_choose(1),lambda_vals(l),bitrate_vals(b),duobinary_vals(2),rop_atten_vals(1));
+
+
+        end
+
+        cols = linspecer(4);
+        subplot(3,3,l)
+
+        if M_choose == 4
+            lst = '-';
+            mkr = 'o';
+            hv = 'on';
+        elseif M_choose == 6
+            lst = '-';
+            mkr = 'x';
+            hv = 'on';
+        elseif M_choose == 8
+            lst = '-';
+            mkr = 'diamond';
+            hv = 'on';
+        end
+
+
+        fsym_vals = floor( bitrate_vals*1e-9./log2(M_choose) );
+        hold on
+
+        plot(bitrate_vals*1e-9,ber_db,'Color',cols(1,:),'Marker',mkr,'MarkerFaceColor','auto','DisplayName','DB pre','LineStyle',lst,'HandleVisibility',hv);
+        plot(bitrate_vals*1e-9,ber_db_enc,'Color',cols(2,:)','Marker',mkr,'MarkerFaceColor','auto','DisplayName','DB enc','LineStyle',lst,'HandleVisibility',hv);
+        plot(bitrate_vals*1e-9,ber_vnle,'Color',cols(3,:),'Marker',mkr,'MarkerFaceColor','auto','DisplayName','VNLE','LineStyle',lst,'HandleVisibility',hv);
+        plot(bitrate_vals*1e-9,ber_vnle_mlse,'Color',cols(4,:),'Marker',mkr,'MarkerFaceColor','auto','DisplayName','VNLE+PF+MLSE','LineStyle',lst,'HandleVisibility',hv);
+
+        % Continue with the rest of your plot settings
+
+        yline(3.8e-3, 'DisplayName', 'HD-FEC', 'LineStyle', '--', 'HandleVisibility', 'off');
+        yline(2e-2, 'DisplayName', '20%', 'LineStyle', '--','LineWidth',1, 'HandleVisibility', 'off');
+        xlabel('Bitrate');
+        ylabel('Bit Error Rate (BER)');
+        title([num2str(lambda_vals(l)),' nm']);
+        set(gca, 'yscale', 'log');
+        set(gca, 'Box', 'on');
+        grid on;
+        grid minor;
+        legend('Interpreter', 'none','Location','southwest');
+        ylim([1e-4,1e-1]);
+        xlim([bitrate_vals(1)*1e-9,bitrate_vals(end)*1e-9])
+
+    end
+end

projects/HighSpeedExperiment_2024/master_evaluation.m

@@ -1,22 +1,24 @@
 
-folderpath = 'C:\Users\sioe\Documents\High_Speed_Measurement_2024\8km_bitrate_rop_master\';
+folderpath = 'C:\Users\sioe\Documents\High_Speed_Measurement_2024\10km_bitrate_complete\';
 experiment_name = '';
 currentTime = datetime('now', 'Format', 'yyyyMMdd_HHmmss');
 timeStr = char(currentTime);
 experiment_name = [experiment_name, timeStr];
 
 if 1
+
     %%% BITRATE Sweep for MPI Experiment %%%
     awg_vpp = 2.7;
     pd_in_set  = 8;
     random_key = 0;
 
     params = struct;
-    params.M = [4,6,8];
+    params.M = [8];
-    params.lambda = flip([1293, 1302, 1310, 1318, 1327.4]); %calcWavelengthPlan(16, 400e9 , 1310);
+    params.lambda = [1293, 1297.5, 1302, 1306.5, 1310, 1313.4, 1318, 1322.7, 1327.4]; %calcWavelengthPlan(16, 400e9 , 1310);
-    params.bitrate = [300:30:480].*1e9;
+    params.bitrate = [270,300,330,360,390,400,410,420,430,440,450,460,470,480].*1e9;
-    params.duobinary = [0,1];
+    params.duobinary = [0,1,2];
     params.rop_atten = [0];
+
 end
 
 if 0
@@ -51,6 +53,8 @@ wh.addStorage("dcs");
 wh.addStorage("pdfa");
 wh.addStorage("exfo");
 wh.addStorage("voa");
+wh.addStorage("filename");
+
 
 precomp_path = "C:\Users\sioe\Documents\High_Speed_Measurement_2024\precomp\";
 precomp_fn = "lab_high_speed";
@@ -72,94 +76,118 @@ estimatedTotalTime = 0;
 
 for M = wh.parameter.M.values
     
-
     dcs = DC_supply("active",[1,0],"voltage",[2.3, 0]);
 
-    if M == 4
+    for db = wh.parameter.duobinary.values
 
-        v_bias_for_pam = 2.3;
+        if db == 1
-        dcs.set("voltage",[v_bias_for_pam, 0]);
+            ffe_only            = 0;
-        pulsef = 1;
+            postfilter_approach = 0;
+            db_channel_approach = 1;
+            db_coding_approach  = 0;
+            db_precode = db_coding_approach || db_channel_approach;
+            if M == 4
+                pulsef=1;
+                precomp_amp_max = -50;
+                v_bias_for_pam = 2.3;
+                dcs.set("voltage",[v_bias_for_pam, 0]);
+                pulsef = 1;
+            elseif M == 6
+                pulsef=0;
+                precomp_amp_max = -50;
+                v_bias_for_pam = 2.3;
+                dcs.set("voltage",[v_bias_for_pam, 0]);
+                pulsef = 1;
+            elseif M == 8
+                pulsef=0;
+                precomp_amp_max = -50;
+                v_bias_for_pam=2.6;
+                dcs.set("voltage",[v_bias_for_pam, 0]);
+                pause(7*60); %wait 30 minutes for stable bias
+                pulsef = 0;
+            end
 
-    elseif M == 6
+        elseif db == 2
-        %pause(7*60); %wait 30 minutes for stable bias
-        v_bias_for_pam=2.3;
-        dcs.set("voltage",[v_bias_for_pam, 0]);
-        pulsef = 0;
 
-    elseif M == 8
+            ffe_only            = 0;
+            postfilter_approach = 0;
+            db_channel_approach = 0;
+            db_coding_approach  = 1;
+            db_precode = db_coding_approach || db_channel_approach;
+            if M == 4
+                pulsef=1;
+                precomp_amp_max = -38;
+                v_bias_for_pam = 2.8;
+                dcs.set("voltage",[v_bias_for_pam, 0]);
+                pulsef = 1;
+            elseif M == 6
+                pulsef=0;
+                precomp_amp_max = -38;
+                v_bias_for_pam = 2.8;
+                dcs.set("voltage",[v_bias_for_pam, 0]);
+                pulsef = 1;
+            elseif M == 8
+                pulsef=0;
+                precomp_amp_max = -38;
+                v_bias_for_pam = 2.8;
+                dcs.set("voltage",[v_bias_for_pam, 0]);
+                pulsef = 1;
+            end
 
-        v_bias_for_pam=2.6;
+        elseif db == 0
-        dcs.set("voltage",[v_bias_for_pam, 0]);
-        pause(7*60); %wait 30 minutes for stable bias
-        pulsef = 0;
-
-    end
-
-    for lambda = wh.parameter.lambda.values
-
-        exfo = Exfo_laser("serialport_number",'COM8','mainframe_channel',1,'safety_mode',0);
-        pdfa = Thor_PDFA("safety_mode",0);
-        exfo.getLaserInfo;
-
-        if ~(exfo.cur_wavelength == lambda)
-
-            % 1)
-            pdfa.disablePDFA;
-
-            % 2)
-            exfo.setWavelength(lambda);
-
-            % 3)
-            pdfa.enablePDFA();
-
-            % 4)
-            pdfa.setPumpLevel(100);
 
+            ffe_only            = 0;
+            postfilter_approach = 1;
+            db_channel_approach = 0;
+            db_coding_approach  = 0;
+            db_precode = db_coding_approach || db_channel_approach;
+            if M == 4
+                pulsef=1;
+                precomp_amp_max = -37;
+                v_bias_for_pam = 2.3;
+                dcs.set("voltage",[v_bias_for_pam, 0]);
+                pulsef = 1;
+            elseif M == 6
+                pulsef=0;
+                precomp_amp_max = -34;
+                v_bias_for_pam = 2.3;
+                dcs.set("voltage",[v_bias_for_pam, 0]);
+                pulsef = 1;
+            elseif M == 8
+                pulsef=0;
+                precomp_amp_max = -34;
+                v_bias_for_pam=2.6;
+                dcs.set("voltage",[v_bias_for_pam, 0]);
+                pause(7*60); %wait 30 minutes for stable bias
+                pulsef = 0;
+            end
         end
 
-        for bitrate = wh.parameter.bitrate.values
+        for lambda = wh.parameter.lambda.values
 
-            fsym = floor( bitrate*1e-9./log2(M) ).*1e9;
+            exfo = Exfo_laser("serialport_number",'COM8','mainframe_channel',1,'safety_mode',0);
+            pdfa = Thor_PDFA("safety_mode",0);
+            exfo.getLaserInfo;
 
-            for db = wh.parameter.duobinary.values
+            if ~(exfo.cur_wavelength == lambda)
 
-                if db == 1
+                % 1)
-                    ffe_only            = 0;
+                pdfa.disablePDFA;
-                    postfilter_approach = 0;
-                    db_channel_approach = 1;
-                    db_coding_approach  = 0;
-                    db_precode = db_coding_approach || db_channel_approach;
-                    if M == 4
-                        pulsef=1;
-                        precomp_amp_max = -50;
-                    elseif M == 6
-                        pulsef=0;
-                        precomp_amp_max = -50;
-                    elseif M == 8
-                        pulsef=0;
-                        precomp_amp_max = -50;
-                    end
-                elseif db == 0
-                    ffe_only            = 0;
-                    postfilter_approach = 1;
-                    db_channel_approach = 0;
-                    db_coding_approach  = 0;
-                    db_precode = db_coding_approach || db_channel_approach;
-                    if M == 4
-                        pulsef=1;
-                        precomp_amp_max = -38;
-                    elseif M == 6
-                        pulsef=0;
-                        precomp_amp_max = -34;
-                    elseif M == 8
-                        pulsef=0;
-                        precomp_amp_max = -34;
-                    end
-                end
 
+                % 2)
+                exfo.setWavelength(lambda);
 
+                % 3)
+                pdfa.enablePDFA();
 
+                % 4)
+                pdfa.setPumpLevel(86);
+
+            end
+
+            for bitrate = wh.parameter.bitrate.values
+
+                fsym = floor( bitrate*1e-9./log2(M) ).*1e9;
                 %%%%% Construct AWG and Scope Modules %%%%%%
                 fdac        = 256e9;
                 fadc        = 256e9;
@@ -174,7 +202,7 @@ for M = wh.parameter.M.values
                 Pamsource = PAMsource(...
                     "fsym",fsym,"M",M,"order",19,"useprbs",1,...
                     "fs_out",fdac,...
-                    "applyclipping",0,"clipfactor",1.5,...
+                    "applyclipping",0,"clipfactor",1.2,...
                     "applypulseform",pulsef,"pulseformer",Pform,...
                     "randkey",random_key,...
                     "db_precode",db_precode,"db_encode",db_coding_approach,...
@@ -182,19 +210,39 @@ for M = wh.parameter.M.values
 
                 [Digi_sig,Symbols,Bits] = Pamsource.process();
 
+                % Digi_sig.plot("displayname","Digi_sig clipped","fignum",21,"clear",1);
+
                 %%%%% Precompensation Routine %%%%%%
                 if precomp_mode == 1 % measure channel
-                    precomp_est = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',fdac);
+                    precomp_est = ChannelFreqResp("Nacq",2048,"Navg",100,"Ncp",63,'f_ref',fdac);
                     Digi_sig = precomp_est.buildOFDM();
                 elseif precomp_mode == 2 % apply precomp
-                    precomp_est = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs);
+                    precomp_est = ChannelFreqResp("Nacq",2048,"Navg",100,"Ncp",63,'f_ref',Digi_sig.fs);
                     Digi_sig = precomp_est.precomp(Digi_sig,'maxampdb',precomp_amp_max,'loadPath',precomp_path,'fileName',precomp_fn);
                 end
 
                 %%%%% Resample to DAC rate %%%%%%
                 Digi_sig = Digi_sig.resample("fs_out",AWG.fdac);
 
-                % Digi_sig.spectrum("displayname","TX After precomp","fignum",30,"normalizeToNyquist",0,"normalizeTo0dB",1);
+                % Digi_sig.spectrum("displayname","TX After precomp","fignum",222,"normalizeToNyquist",0,"normalizeTo0dB",1);
+
+                if 0 %negative performance...
+                    % X: design FIR filter for sinc precomp
+                    % https://www.dsprelated.com/showarticle/1191.php
+                    ntaps = 13;
+                    npts = 32;
+                    % least-squares FIR design
+                    fmax = AWG.fdac*0.5;
+                    ff = linspace(0,fmax,npts);
+                    hsinc = sin(pi*ff/AWG.fdac)./(pi*ff/AWG.fdac + eps); % transfer function of sample and hold DAC
+                    hsinc(1) = 1;
+                    h_goal= 1./hsinc; % goal function
+                    f = 2.*ff./AWG.fdac; %vector between 0 and 1, where 1 is nyquist is fsamp/2
+                    b       = firls(ntaps-1,f,h_goal);
+                    Digi_sig.signal = conv(Digi_sig.signal,b,"same");
+                end
+
+                % Digi_sig.spectrum("displayname","TX After SINC precomp","fignum",222,"normalizeToNyquist",0,"normalizeTo0dB",1);
 
                 for rop_atten = wh.parameter.rop_atten.values
 
@@ -224,12 +272,12 @@ for M = wh.parameter.M.values
                     %%%%% AWG --> Scope %%%%%%
                     [~,Scpe_sig_raw,~,D] = A2S.process("signal2",Digi_sig,"waitUntilClick",0);
                     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
+                    save([folderpath,experiment_name,loop_name,'_raw_signal'],"Scpe_sig_raw");
 
                     % Scpe_sig_raw = Filter('filtdegree',5,"f_cutoff",0.55.*fsym,"fs",fadc,"filterType",filtertypes.gaussian,"active",true).process(Scpe_sig_raw);
 
                     %
-                    % Scpe_sig_raw.plot("displayname","Scope raw signal","fignum",20,"clear",1);
+                    Scpe_sig_raw.plot("displayname","Scope raw signal","fignum",20,"clear",1);
                     % Scpe_sig_raw.spectrum("displayname","Scope PSD","fignum",30,"normalizeTo0dB",1);
                     % Scpe_sig_raw.eye(fsym,M,"displayname",'eye','fignum',200);
 
@@ -252,7 +300,7 @@ for M = wh.parameter.M.values
 
                     %%%%% Plot and Save Routines: SAVE RECEIVED SIGNALS %%%%%%%%%%%%%%%%%%%%%%%%%
                     save([folderpath,experiment_name,loop_name,'_rx_signal'],"S");
-                    save([folderpath,experiment_name,loop_name,'_raw_signal'],"Scpe_sig_raw");
+
 
                     %%%%% EQUALIZE %%%%%%
                     % set to minus one not zero not avoid confusion if BER is acutally zero
@@ -264,39 +312,20 @@ for M = wh.parameter.M.values
                     ffe = EQ("Ne",[50,0,0],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.05,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
                     vnle = EQ("Ne",[50,7,7],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.05,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
 
-
                     if postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%%
 
 
-                        eq_values = min(numel(S),8);
+                        if 1
-                        Noi = cell(eq_values,1);
+
-                        EQ_vnle= cell(eq_values,1);
+                            eq_values = min(numel(S),8);
-                        EQ_ffe= cell(eq_values,1);
+                            Noi = cell(eq_values,1);
-                        if 0
+                            EQ_vnle= cell(eq_values,1);
+                            EQ_ffe= cell(eq_values,1);
+
                             parfor s = 1:eq_values
-
-                                if 0
-                                    %FFE LINEAR
-                                    Scpe_sig_syncd = S{s};
-                                    [EQ_ffe{s}] = ffe.process(Scpe_sig_syncd,Symbols);
-                                    Noi{s} = EQ_ffe{s}-Symbols;
-                                    Rx_bits = PAMmapper(M,0).demap(EQ_ffe{s});
-                                    [~,~,ber_ffe(s),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
-                                end
-
-                                if 0
-                                    %FFE + MLSE
-                                    nc = 2;
-                                    burg_coeff = arburg(Noi{s}.signal,nc);
-                                    EQ_ffe{s} = EQ_ffe{s}.filter(burg_coeff,1);
-
-                                    EQ_mlse = MLSE("DIR",burg_coeff,"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_ffe{s});
-                                    Rx_bits = PAMmapper(M,0).demap(EQ_mlse);
-                                    [~,~,ber_ffe_mlse(s),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
-                                end
-
                                 if 1
                                     %VNLE
+                                    vnle = 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);
                                     Scpe_sig_syncd = S{s};
                                     [EQ_vnle{s}] = vnle.process(Scpe_sig_syncd,Symbols);
                                     Noi{s} = EQ_vnle{s}-Symbols;
@@ -306,6 +335,7 @@ for M = wh.parameter.M.values
 
                                 %VNLE + MLSE
                                 if 1
+                                    Noi{s}.signal = Noi{s}.signal - mean(Noi{s}.signal);
                                     nc = 2;
                                     burg_coeff = arburg(Noi{s}.signal,nc);
                                     EQ_mlse = EQ_vnle{s}.filter(burg_coeff,1);
@@ -317,76 +347,149 @@ for M = wh.parameter.M.values
                             end
 
 
-                            disp(['FFE EQ: BEST BER: ',sprintf('%.1E',min(ber_ffe)),' AVG BER: ',sprintf('%.1E',mean(ber_ffe)),' WORST:',sprintf('%.1E',max(ber_ffe)),'. Out of ',num2str(numel(ber_ffe))]);
+                            if 0
-                            disp(['FFE + MLSE EQ: BEST BER: ',sprintf('%.1E',min(ber_ffe_mlse)),' AVG BER: ',sprintf('%.1E',mean(ber_ffe_mlse)),' WORST:',sprintf('%.1E',max(ber_ffe_mlse)),'. Out of ',num2str(numel(ber_ffe_mlse))]);
+                                [~,s]=min(ber_vnle_mlse);
-                            disp(['VNLE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle)),' WORST:',sprintf('%.1E',max(ber_vnle)),'. Out of ',num2str(numel(ber_vnle))]);
+                                nc = 2;
-                            disp(['VNLE+MLSE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle_mlse)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle_mlse)),' WORST:',sprintf('%.1E',max(ber_vnle_mlse)),'. Out of ',num2str(numel(ber_ffe))]);
+                                burg_coeff = arburg(Noi{s}.signal,nc);
-
+                                Noi{s}.spectrum('displayname','Noise PSD','fignum',123);
-                            [~,i] = min(ber_vnle);
+                                [h,w] = freqz(1,burg_coeff,length(Noi{s}),"whole",Noi{s}.fs);
-                            figure(56);
+                                h = h/max(abs(h));
-                            clf
-                            title(sprintf('PAM %d ; BER: %1.2e',M, ber_vnle(i)));
-                            constellation = unique(Symbols.signal);
-                            received = NaN(numel(constellation),length(Symbols));
-                            for lvl = 1:numel(constellation)
-                                %Separate the equalized signal into the
-                                %respective levels based on the actually
-                                %transmitted level!
-                                received(lvl,Symbols.signal==constellation(lvl)) = EQ_vnle{i}.signal(Symbols.signal==constellation(lvl));
-                                intermediate = received(lvl,:);
-                                cnt(lvl) = numel(intermediate(~isnan(intermediate)));
                                 hold on
-                                histogram(received(lvl,:),1000,"EdgeAlpha",0,'DisplayName',['Lvl ',num2str(lvl),' | ',num2str(cnt(lvl)),' entries']);
+                                w_ = (w - Noi{s}.fs/2);
+                                plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']);
                             end
-                            legend
+
+                            if 0
+                                figure(55);
+                                clf
+                                title(sprintf('PAM %d ; BER: %1.2e',M, ber_vnle_mlse(s) ));
+                                constellation = unique(Symbols.signal);
+                                received = NaN(numel(constellation),length(Symbols));
+                                for lvl = 1:numel(constellation)
+                                    %Separate the equalized signal into the
+                                    %respective levels based on the actually
+                                    %transmitted level!
+                                    received(lvl,Symbols.signal==constellation(lvl)) = EQ_vnle{s}.signal(Symbols.signal==constellation(lvl));
+                                    intermediate = received(lvl,:);
+                                    cnt(lvl) = numel(intermediate(~isnan(intermediate)));
+                                    hold on
+                                    histogram(received(lvl,:),1000,"EdgeAlpha",0,'DisplayName',['Lvl ',num2str(lvl),' | ',num2str(cnt(lvl)),' entries']);
+                                end
+                                legend
+                            end
+
+
+
+                            % disp(['FFE EQ: BEST BER: ',sprintf('%.1E',min(ber_ffe)),' AVG BER: ',sprintf('%.1E',mean(ber_ffe)),' WORST:',sprintf('%.1E',max(ber_ffe)),'. Out of ',num2str(numel(ber_ffe))]);
+                            % disp(['FFE + MLSE EQ: BEST BER: ',sprintf('%.1E',min(ber_ffe_mlse)),' AVG BER: ',sprintf('%.1E',mean(ber_ffe_mlse)),' WORST:',sprintf('%.1E',max(ber_ffe_mlse)),'. Out of ',num2str(numel(ber_ffe_mlse))]);
+                            disp(['VNLE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle)),' WORST:',sprintf('%.1E',max(ber_vnle)),'. Out of ',num2str(numel(ber_vnle))]);
+                            % disp(['VNLE+MLSE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle_mlse)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle_mlse)),' WORST:',sprintf('%.1E',max(ber_vnle_mlse)),'. Out of ',num2str(numel(ber_vnle_mlse))]);
                         end
 
                     elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%%
 
                         ffe = EQ("Ne",[50,0,0],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.05,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
+                        ffe = 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);
 
-                        if 0
+                        if 1
-                            parfor s = 1:numel(S)
+
+                            eq_values = min(numel(S),8);
+                            Noi = cell(eq_values,1);
+                            EQ_sig = cell(eq_values,1);
+
+                            parfor s = 1:eq_values
 
                                 Scpe_sig_syncd = S{s};
 
-                                [EQ_sig, Noi] = ffe.process(Scpe_sig_syncd,Duobinary().encode(Symbols));
+                                [EQ_sig{s}, Noi{s}] = ffe.process(Scpe_sig_syncd,Duobinary().encode(Symbols));
 
-                                EQ_sig = MLSE("DIR",[1,1],"duobinary_output",1,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig);
+                                EQ_sig{s}.signal = EQ_sig{s}.signal-mean(EQ_sig{s}.signal);
 
-                                EQ_sig = Duobinary().decode(EQ_sig);
+                                EQ_sig_mlse = MLSE("DIR",[1,1],"duobinary_output",1,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig{s});
 
-                                Rx_bits = PAMmapper(M,0).demap(EQ_sig);
+                                EQ_sig_mlse = Duobinary().decode(EQ_sig_mlse);
+
+                                Rx_bits = PAMmapper(M,0).demap(EQ_sig_mlse);
 
                                 [~,num_errors,ber_db(s),pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
 
                             end
 
+                            if 0
+                                [~,s]=min(ber_vnle_mlse);
+
+                                Noi{s}.spectrum('displayname',['Noise; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1);
+
+                                Duobinary().encode(Symbols).spectrum('displayname',['DB coded symbols; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1);
+
+                                EQ_sig{s}.spectrum('displayname',['EQ; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1);
+
+                                EQ_sig{s}.signal = EQ_sig{1}.signal-mean(EQ_sig{s}.signal);
+                            end
+
                             disp(['DB EQ: BEST BER: ',sprintf('%.1E',min(ber_db)),' AVG BER: ',sprintf('%.1E',mean(ber_db)),' WORST:',sprintf('%.1E',max(ber_db)),'. Out of',num2str(numel(ber_db))]);
 
                         else
 
-                            disp('Disabled MLSE for DB in all cases, due to time in measurement loop')
+                            % disp('Disabled MLSE for DB in all cases, due to time in measurement loop')
 
                         end
 
                     elseif db_coding_approach
 
-                        parfor s = 1:numel(S)
+                        ffe = 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);
-                            S{s}=S{s}.normalize("mode","rms");
+
-                            ffe = EQ("Ne",[50,7,7],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.05,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
+                        if 1
-                            [EQ_sig, Noi] = ffe.process(S{s},Symbols);
+
-                            EQ_sig.plot("displayname",'After EQ','fignum',112);
+                            eq_values = min(numel(S),8);
-                            EQ_sig = MLSE("DIR",[1,1],"duobinary_output",1,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig);
+                            Noi = cell(eq_values,1);
-                            EQ_sig = Duobinary().decode(EQ_sig);
+                            EQ_sig = cell(eq_values,1);
-                            Rx_bits = PAMmapper(M,0).demap(EQ_sig);
+                            EQ_sig_mlse = cell(eq_values,1);
-                            [~,num_errors,ber_db(s),pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+
+                            parfor s = 1:eq_values
+                                [EQ_sig{s}, Noi{s}] = ffe.process(S{s},Symbols);
+                                EQ_sig_mlse{s} = MLSE("DIR",[1,1],"duobinary_output",1,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig{s});
+                                EQ_sig_mlse{s} = Duobinary().decode(EQ_sig_mlse{s});
+                                Rx_bits = PAMmapper(M,0).demap(EQ_sig_mlse{s});
+                                [~,num_errors,ber_db(s),pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+                            end
+
+                            disp(['DB EQ: BEST BER: ',sprintf('%.1E',min(ber_db)),' AVG BER: ',sprintf('%.1E',mean(ber_db)),' WORST:',sprintf('%.1E',max(ber_db)),'. Out of',num2str(numel(ber_db))]);
+
+                            if 1
+                                [~,s]=min(ber_db);
+                                Noi{s}.spectrum('displayname',['Noise '],'fignum',50,'normalizeTo0dB',1);
+                                EQ_sig{s}.spectrum('displayname',['EQLZD '],'fignum',50,'normalizeTo0dB',1);
+                                Symbols.spectrum('displayname',['Symbols '],'fignum',222,'normalizeTo0dB',1);
+                            end
+
+                            if 1
+                                figure(51);
+                                clf
+                                title(sprintf('DB coded PAM after EQ ; BER: %1.2e',M, ber_db(s) ));
+                                constellation = unique(Symbols.signal);
+                                received = NaN(numel(constellation),length(Symbols));
+                                for lvl = 1:numel(constellation)
+                                    %Separate the equalized signal into the
+                                    %respective levels based on the actually
+                                    %transmitted level!
+                                    received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig{s}.signal(Symbols.signal==constellation(lvl));
+                                    intermediate = received(lvl,:);
+                                    cnt(lvl) = numel(intermediate(~isnan(intermediate)));
+                                    hold on
+                                    histogram(received(lvl,:),1000,"EdgeAlpha",0,'DisplayName',['Lvl ',num2str(lvl),' | ',num2str(cnt(lvl)),' entries']);
+                                end
+                                legend
+                            end
+
                         end
 
-                          disp(['DB EQ: BEST BER: ',sprintf('%.1E',min(ber_db)),' AVG BER: ',sprintf('%.1E',mean(ber_db)),' WORST:',sprintf('%.1E',max(ber_db)),'. Out of',num2str(numel(ber_db))]);
-
                     end
 
+                    % showCurrentMeasurement('Vbias', v_bias_for_pam,'MIN BER', min(ber_db),'MEAN BER',mean(ber_db),'MAX BER',max(ber_db), 'Fsym',fsym.*1e-9, 'ROP', rop, 'Precomp MaxAmp',precomp_amp_max);
+
+
+
 
                     %%%%% Store measurement into measurement "warehouse" %%%%%%
 
@@ -408,6 +511,8 @@ for M = wh.parameter.M.values
                     wh.addValueToStorage(ffe,'FFE',M,lambda,bitrate,db,rop_atten);
                     wh.addValueToStorage(vnle,'VNLE',M,lambda,bitrate,db,rop_atten);
 
+                    wh.addValueToStorage(string([experiment_name,loop_name]),'filename',M,lambda,bitrate,db,rop_atten);
+
 
                     iterationTimes(loopcnt) = toc(iterationStartTime);
                     averageTimePerIteration = mean(iterationTimes(1:loopcnt));

projects/HighSpeedExperiment_2024/mpi_measurement.m

@@ -1,6 +1,6 @@
 
-folderpath = 'C:\Users\sioe\Documents\High_Speed_Measurement_2024\mpi_measurement\';
+folderpath = 'C:\Users\sioe\Documents\High_Speed_Measurement_2024\MPI_duobinary_encoded\';
-experiment_name = 'testen';
+experiment_name = '';
 currentTime = datetime('now', 'Format', 'yyyyMMdd_HHmmss');
 timeStr = char(currentTime);
 experiment_name = [experiment_name, timeStr];
@@ -16,10 +16,10 @@ pd_in_set       = 8; %VOA 2 -> PD in
 random_key      = 0;
 
 params          = struct;
-params.M        = [8];
+params.M        = [4];
-params.bitrate  = [336].*1e9;%[90:60:480].*1e9;%[300:30:480].*1e9;
+params.bitrate  = [336,360,390,420,448].*1e9;%[90:60:480].*1e9;%[300:30:480].*1e9;
-params.duobinary = [0];
+params.duobinary = [2];
-params.interference_atten = [45];
+params.interference_atten = [0:1:20,45];
 
 wh = DataStorage(params);
 
@@ -38,7 +38,7 @@ wh.addStorage("s_power");
 wh.addStorage("i_power");
 wh.addStorage("sir");
 
-
+wh.addStorage("filename");
 wh.addStorage("m");
 
 wh.addStorage("dcs");
@@ -66,37 +66,13 @@ estimatedTotalTime = 0;
 
 for M = wh.parameter.M.values
 
-
-    dcs = DC_supply("active",[1,0],"voltage",[2.3, 0]);
-
-    if M == 4
-
-        v_bias_for_pam = 2.3;
-        dcs.set("voltage",[v_bias_for_pam, 0]);
-        pulsef = 1;
-
-    elseif M == 6
-        %pause(7*60); %wait 30 minutes for stable bias
-        v_bias_for_pam=2.3;
-        dcs.set("voltage",[v_bias_for_pam, 0]);
-        pulsef = 0;
-
-    elseif M == 8
-
-        v_bias_for_pam=2.6;
-        dcs.set("voltage",[v_bias_for_pam, 0]);
-        disp('waiting...');
-        %pause(5*60); %wait 5 minutes for stable bias
-        pulsef = 0;
-
-    end
-
     for bitrate = wh.parameter.bitrate.values
 
         fsym = floor( bitrate*1e-9./log2(M) ).*1e9;
 
         for db = wh.parameter.duobinary.values
 
+            dcs = DC_supply("active",[1,0],"voltage",[2.3, 0]);
             if db == 1
                 ffe_only            = 0;
                 postfilter_approach = 0;
@@ -106,13 +82,50 @@ for M = wh.parameter.M.values
                 if M == 4
                     pulsef=1;
                     precomp_amp_max = -50;
+                    v_bias_for_pam = 2.3;
+                    dcs.set("voltage",[v_bias_for_pam, 0]);
+                    pulsef = 1;
                 elseif M == 6
                     pulsef=0;
                     precomp_amp_max = -50;
+                    v_bias_for_pam = 2.3;
+                    dcs.set("voltage",[v_bias_for_pam, 0]);
+                    pulsef = 1;
                 elseif M == 8
                     pulsef=0;
                     precomp_amp_max = -50;
+                    v_bias_for_pam=2.6;
+                    dcs.set("voltage",[v_bias_for_pam, 0]);
+                    pause(7*60); %wait 30 minutes for stable bias
+                    pulsef = 0;
                 end
+
+            elseif db == 2
+                ffe_only            = 0;
+                postfilter_approach = 0;
+                db_channel_approach = 0;
+                db_coding_approach  = 1;
+                db_precode = db_coding_approach || db_channel_approach;
+                if M == 4
+                    pulsef=1;
+                    precomp_amp_max = -38;
+                    v_bias_for_pam = 2.8;
+                    dcs.set("voltage",[v_bias_for_pam, 0]);
+                    pulsef = 1;
+                elseif M == 6
+                    pulsef=0;
+                    precomp_amp_max = -38;
+                    v_bias_for_pam = 2.8;
+                    dcs.set("voltage",[v_bias_for_pam, 0]);
+                    pulsef = 1;
+                elseif M == 8
+                    pulsef=0;
+                    precomp_amp_max = -38;
+                    v_bias_for_pam = 2.8;
+                    dcs.set("voltage",[v_bias_for_pam, 0]);
+                    pulsef = 1;
+                end
+
             elseif db == 0
                 ffe_only            = 0;
                 postfilter_approach = 1;
@@ -121,13 +134,23 @@ for M = wh.parameter.M.values
                 db_precode = db_coding_approach || db_channel_approach;
                 if M == 4
                     pulsef=1;
-                    precomp_amp_max = -38;
+                    precomp_amp_max = -37;
+                    v_bias_for_pam = 2.3;
+                    dcs.set("voltage",[v_bias_for_pam, 0]);
+                    pulsef = 1;
                 elseif M == 6
                     pulsef=0;
                     precomp_amp_max = -34;
+                    v_bias_for_pam = 2.3;
+                    dcs.set("voltage",[v_bias_for_pam, 0]);
+                    pulsef = 1;
                 elseif M == 8
                     pulsef=0;
                     precomp_amp_max = -34;
+                    v_bias_for_pam=2.6;
+                    dcs.set("voltage",[v_bias_for_pam, 0]);
+                    pause(7*60); %wait 30 minutes for stable bias
+                    pulsef = 0;
                 end
             end
 
@@ -161,24 +184,23 @@ for M = wh.parameter.M.values
             end
 
             %%%%% Resample to DAC rate %%%%%%
-            Digi_sig = Digi_sig.resample("fs_out",AWG.fdac);
+            Digi_sig = Digi_sig.resample("fs_out",fdac);
 
             % Digi_sig = Filter('filtdegree',5,"f_cutoff",0.75*fsym,"fs",fadc,"filterType",filtertypes.gaussian,"active",true).process(Digi_sig);
 
             % Digi_sig.spectrum("displayname","TX After precomp","fignum",10,"normalizeToNyquist",0,"normalizeTo0dB",0);
 
-
+            % holdAndShowValue;
-            holdAndShowValue;
 
             scopeAutoScale = 1;
 
             for interference_atten = wh.parameter.interference_atten.values
 
-                SCP = ScopeKeysight("model","UXR1104B",'autoscale',scopeAutoScale,"fadc","GSa_256","channel",[0,1,0,0],"recordLen",10000000,"removeDC",1);
+                SCP = ScopeKeysight("model","UXR1104B",'autoscale',scopeAutoScale,"fadc","GSa_256","channel",[0,1,0,0],"recordLen",6000000,"removeDC",1);
                 AWG = AwgKeysight("model","M8199B","fdac",fdac,"scaletodac",[1,1],"skews",[0,0],"voltages",[0,awg_vpp]);
                 A2S = Awg2Scope(AWG,SCP,[0,2,0,0],"waitUntilClick",1); %
 
-                scopeAutoScale = 0; %until is set to 1 in next db change and then bitrate
+                % scopeAutoScale = 0; %until is set to 1 in next db change and then bitrate
 
                 %%%%% Loop Preps
                 iterationStartTime = tic;
@@ -204,35 +226,34 @@ for M = wh.parameter.M.values
 
                 %%%%% AWG --> Scope %%%%%%
                 [~,Scpe_sig_raw,~,D] = A2S.process("signal2",Digi_sig,"waitUntilClick",0);
-                %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+                save([folderpath,experiment_name,loop_name,'_raw_signal'],"Scpe_sig_raw");
-
-                % Scpe_sig_raw.spectrum("displayname","Scope PSD before filter","fignum",30,"normalizeTo0dB",1);
-
-                Scpe_sig_raw = Filter('filtdegree',5,"f_cutoff",0.65.*fsym,"fs",fadc,"filterType",filtertypes.gaussian,"active",true).process(Scpe_sig_raw);
-
-                % Scpe_sig_raw.spectrum("displayname","Scope PSD after filter","fignum",30,"normalizeTo0dB",1);
-
-                %%%%%% Sample to 2x fsym %%%%%%
-                Scpe_sig_resampled = Scpe_sig_raw.resample("fs_in",fadc,"fs_out",2*fsym);
-
                 voa.readvals();
                 rop = voa.power_state(1);
                 pd_in = voa.power_state(2);
                 i_power = voa.power_state(4);
                 s_power = voa.power_state(3);
                 sir = s_power-i_power;
+                %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+                % Scpe_sig_raw.spectrum("displayname","Scope PSD before filter","fignum",30,"normalizeTo0dB",1);
+
+                % Scpe_sig_raw = Filter('filtdegree',5,"f_cutoff",0.65.*fsym,"fs",fadc,"filterType",filtertypes.gaussian,"active",true).process(Scpe_sig_raw);
+
+                % Scpe_sig_raw.spectrum("displayname","Scope PSD after filter","fignum",30,"normalizeTo0dB",1);
+
+                %%%%%% Sample to 2x fsym %%%%%%
+                Scpe_sig_resampled = Scpe_sig_raw.resample("fs_in",fadc,"fs_out",2*fsym);
 
                 Scpe_sig_raw.plot("displayname",['SIR: ',sprintf('%.2f',sir),' dB'],"fignum",20,"clear",1);
 
-
+                % disp(['PDin: ',sprintf('%.2f',pd_in),' dB | S: ',sprintf('%.2f',s_power),' dB | I: ',sprintf('%.2f',i_power),' dB | -> SIR: ',sprintf('%.2f',sir),' dB']);
-                disp(['PDin: ',sprintf('%.2f',pd_in),' dB | S: ',sprintf('%.2f',s_power),' dB | I: ',sprintf('%.2f',i_power),' dB | -> SIR: ',sprintf('%.2f',sir),' dB']);
 
                 %%%%%% Sync Rx signal with reference (S is a cell array with all occurences) %%%%%%
                 [Scpe_sig_syncd,S,isFlipped] = Scpe_sig_resampled.tsynch("reference",Symbols,"fs_ref",fsym);
 
                 %%%%% Plot and Save Routines: SAVE RECEIVED SIGNALS %%%%%%%%%%%%%%%%%%%%%%%%%
                 save([folderpath,experiment_name,loop_name,'_rx_signal'],"S");
-                save([folderpath,experiment_name,loop_name,'_raw_signal'],"Scpe_sig_raw");
+                
 
                 %%%%% EQUALIZE %%%%%%
                 % set to minus one not zero not avoid confusion if BER is acutally zero
@@ -244,7 +265,6 @@ for M = wh.parameter.M.values
                 ffe = EQ("Ne",[50,0,0],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.05,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
                 vnle = EQ("Ne",[50,7,7],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.05,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
 
-
                 if postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%%
 
 
@@ -256,29 +276,9 @@ for M = wh.parameter.M.values
                         EQ_ffe= cell(eq_values,1);
 
                         parfor s = 1:eq_values
-
-                            if 0
-                                %FFE LINEAR
-                                Scpe_sig_syncd = S{s};
-                                [EQ_ffe{s}] = ffe.process(Scpe_sig_syncd,Symbols);
-                                Noi{s} = EQ_ffe{s}-Symbols;
-                                Rx_bits = PAMmapper(M,0).demap(EQ_ffe{s});
-                                [~,~,ber_ffe(s),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
-                            end
-
-                            if 0
-                                %FFE + MLSE
-                                nc = 2;
-                                burg_coeff = arburg(Noi{s}.signal,nc);
-                                EQ_ffe{s} = EQ_ffe{s}.filter(burg_coeff,1);
-
-                                EQ_mlse = MLSE("DIR",burg_coeff,"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_ffe{s});
-                                Rx_bits = PAMmapper(M,0).demap(EQ_mlse);
-                                [~,~,ber_ffe_mlse(s),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
-                            end
-
                             if 1
                                 %VNLE
+                                vnle = 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);
                                 Scpe_sig_syncd = S{s};
                                 [EQ_vnle{s}] = vnle.process(Scpe_sig_syncd,Symbols);
                                 Noi{s} = EQ_vnle{s}-Symbols;
@@ -288,6 +288,7 @@ for M = wh.parameter.M.values
 
                             %VNLE + MLSE
                             if 1
+                                Noi{s}.signal = Noi{s}.signal - mean(Noi{s}.signal);
                                 nc = 2;
                                 burg_coeff = arburg(Noi{s}.signal,nc);
                                 EQ_mlse = EQ_vnle{s}.filter(burg_coeff,1);
@@ -298,29 +299,50 @@ for M = wh.parameter.M.values
                             end
                         end
 
-                        if 0
+
-                            nc=1;
+                        if 1
-                            Noi{1}.spectrum('displayname',['Noise; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1);
+                            [~,s]=min(ber_vnle_mlse);
-                            burg_coeff = arburg(Noi{1}.signal,nc);
+                            nc = 2;
-                            [h,w] = freqz(1,burg_coeff,length(Noi{1}),"whole",Noi{1}.fs);
+                            burg_coeff = arburg(Noi{s}.signal,nc);
+                            Noi{s}.spectrum('displayname','Noise PSD','fignum',123);
+                            [h,w] = freqz(1,burg_coeff,length(Noi{s}),"whole",Noi{s}.fs);
                             h = h/max(abs(h));
                             hold on
-                            w_ = (w - Noi{1}.fs/2);
+                            w_ = (w - Noi{s}.fs/2);
-                            plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',[num2str(nc), ' burg; SIR:',sprintf('%.2f',sir)]);
+                            plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']);
-                            drawnow;
                         end
 
-                    end
+                        if 1
+                            figure(55);
+                            clf
+                            title(sprintf('PAM %d ; BER: %1.2e',M, ber_vnle_mlse(s) ));
+                            constellation = unique(Symbols.signal);
+                            received = NaN(numel(constellation),length(Symbols));
+                            for lvl = 1:numel(constellation)
+                                %Separate the equalized signal into the
+                                %respective levels based on the actually
+                                %transmitted level!
+                                received(lvl,Symbols.signal==constellation(lvl)) = EQ_vnle{s}.signal(Symbols.signal==constellation(lvl));
+                                intermediate = received(lvl,:);
+                                cnt(lvl) = numel(intermediate(~isnan(intermediate)));
+                                hold on
+                                histogram(received(lvl,:),1000,"EdgeAlpha",0,'DisplayName',['Lvl ',num2str(lvl),' | ',num2str(cnt(lvl)),' entries']);
+                            end
+                            legend
+                        end
 
-                    % disp(['FFE EQ: BEST BER: ',sprintf('%.1E',min(ber_ffe)),' AVG BER: ',sprintf('%.1E',mean(ber_ffe)),' WORST:',sprintf('%.1E',max(ber_ffe)),'. Out of ',num2str(numel(ber_ffe))]);
+
-                    % disp(['FFE + MLSE EQ: BEST BER: ',sprintf('%.1E',min(ber_ffe_mlse)),' AVG BER: ',sprintf('%.1E',mean(ber_ffe_mlse)),' WORST:',sprintf('%.1E',max(ber_ffe_mlse)),'. Out of ',num2str(numel(ber_ffe_mlse))]);
+
-                    disp(['VNLE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle)),' WORST:',sprintf('%.1E',max(ber_vnle)),'. Out of ',num2str(numel(ber_vnle))]);
+                        % disp(['FFE EQ: BEST BER: ',sprintf('%.1E',min(ber_ffe)),' AVG BER: ',sprintf('%.1E',mean(ber_ffe)),' WORST:',sprintf('%.1E',max(ber_ffe)),'. Out of ',num2str(numel(ber_ffe))]);
-                    % disp(['VNLE+MLSE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle_mlse)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle_mlse)),' WORST:',sprintf('%.1E',max(ber_vnle_mlse)),'. Out of ',num2str(numel(ber_vnle_mlse))]);
+                        % disp(['FFE + MLSE EQ: BEST BER: ',sprintf('%.1E',min(ber_ffe_mlse)),' AVG BER: ',sprintf('%.1E',mean(ber_ffe_mlse)),' WORST:',sprintf('%.1E',max(ber_ffe_mlse)),'. Out of ',num2str(numel(ber_ffe_mlse))]);
+                        disp(['VNLE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle)),' WORST:',sprintf('%.1E',max(ber_vnle)),'. Out of ',num2str(numel(ber_vnle))]);
+                        % disp(['VNLE+MLSE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle_mlse)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle_mlse)),' WORST:',sprintf('%.1E',max(ber_vnle_mlse)),'. Out of ',num2str(numel(ber_vnle_mlse))]);
+                    end
 
                 elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%%
 
                     ffe = EQ("Ne",[50,0,0],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.05,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
-                    ffe = EQ("Ne",[50,7,7],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.05,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
+                    ffe = 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);
 
                     if 0
 
@@ -346,24 +368,88 @@ for M = wh.parameter.M.values
 
                         end
 
-                        if 0
+                        if 1
-                            Noi{1}.spectrum('displayname',['Noise; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1);
+                            [~,s]=min(ber_vnle_mlse);
+
+                            Noi{s}.spectrum('displayname',['Noise; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1);
 
                             Duobinary().encode(Symbols).spectrum('displayname',['DB coded symbols; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1);
-                            EQ_sig{1}.spectrum('displayname',['EQ; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1);
+
-                            EQ_sig{1}.signal = EQ_sig{1}.signal-mean(EQ_sig{1}.signal);
+                            EQ_sig{s}.spectrum('displayname',['EQ; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1);
+
+                            EQ_sig{s}.signal = EQ_sig{1}.signal-mean(EQ_sig{s}.signal);
                         end
 
                         disp(['DB EQ: BEST BER: ',sprintf('%.1E',min(ber_db)),' AVG BER: ',sprintf('%.1E',mean(ber_db)),' WORST:',sprintf('%.1E',max(ber_db)),'. Out of',num2str(numel(ber_db))]);
 
                     else
 
-                        disp('Disabled MLSE for DB in all cases, due to time in measurement loop')
+                        % disp('Disabled MLSE for DB in all cases, due to time in measurement loop')
+
+                    end
+
+                elseif db_coding_approach
+
+                    ffe = 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);
+
+                    if 1
+
+                        eq_values = min(numel(S),8);
+                        Noi = cell(eq_values,1);
+                        EQ_sig = cell(eq_values,1);
+                        EQ_sig_mlse = cell(eq_values,1);
+
+                        parfor s = 1:eq_values
+                            [EQ_sig{s}, Noi{s}] = ffe.process(S{s},Symbols);
+                            EQ_sig_mlse{s} = MLSE("DIR",[1,1],"duobinary_output",1,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig{s});
+                            EQ_sig_mlse{s} = Duobinary().decode(EQ_sig_mlse{s});
+                            Rx_bits = PAMmapper(M,0).demap(EQ_sig_mlse{s});
+                            [~,num_errors,ber_db(s),pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+                        end
+
+                        disp(['DB EQ: BEST BER: ',sprintf('%.1E',min(ber_db)),' AVG BER: ',sprintf('%.1E',mean(ber_db)),' WORST:',sprintf('%.1E',max(ber_db)),'. Out of',num2str(numel(ber_db))]);
+
+                        if 1
+                            [~,s]=min(ber_db);
+                            Noi{s}.spectrum('displayname',['Noise '],'fignum',50,'normalizeTo0dB',1);
+                            EQ_sig{s}.spectrum('displayname',['EQLZD '],'fignum',50,'normalizeTo0dB',1);
+                            Symbols.spectrum('displayname',['Symbols '],'fignum',50,'normalizeTo0dB',1);
+
+                            nc = 5;
+                            burg_coeff = arburg(Noi{s}.signal,nc);
+                            [h,w] = freqz(1,burg_coeff,length(Noi{s}),"whole",Noi{s}.fs);
+                            h = h/max(abs(h));
+                            hold on
+                            w_ = (w - Noi{s}.fs/2);
+                            plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']);
+
+                        end
+
+                        if 1
+                            figure(51);
+                            clf
+                            title(sprintf('DB coded PAM after EQ ; BER: %1.2e',M, ber_db(s) ));
+                            constellation = unique(Symbols.signal);
+                            received = NaN(numel(constellation),length(Symbols));
+                            for lvl = 1:numel(constellation)
+                                %Separate the equalized signal into the
+                                %respective levels based on the actually
+                                %transmitted level!
+                                received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig{s}.signal(Symbols.signal==constellation(lvl));
+                                intermediate = received(lvl,:);
+                                cnt(lvl) = numel(intermediate(~isnan(intermediate)));
+                                hold on
+                                histogram(received(lvl,:),1000,"EdgeAlpha",0,'DisplayName',['Lvl ',num2str(lvl),' | ',num2str(cnt(lvl)),' entries']);
+                            end
+                            legend
+                        end
 
                     end
 
                 end
 
+                 showCurrentMeasurement('Att.',interference_atten,'SIR',sir,'MIN BER', min(ber_db),'MEAN BER',mean(ber_db),'MAX BER',max(ber_db), 'Fsym',fsym.*1e-9, 'ROP', rop, 'Precomp MaxAmp',precomp_amp_max);
+
 
                 %%%%% Store measurement into measurement "warehouse" %%%%%%
                 wh.addValueToStorage({ber_ffe},'ber_ffe',M,bitrate,db,interference_atten);
@@ -378,6 +464,8 @@ for M = wh.parameter.M.values
                 wh.addValueToStorage(i_power,'i_power',M,bitrate,db,interference_atten);
                 wh.addValueToStorage(sir,'sir',M,bitrate,db,interference_atten);
 
+                wh.addValueToStorage(string([experiment_name,loop_name]),'filename',M,bitrate,db,interference_atten);
+
                 wh.addValueToStorage(M,'m',M,bitrate,db,interference_atten);
 
                 wh.addValueToStorage(dcs,'dcs',M,bitrate,db,interference_atten);
@@ -402,13 +490,10 @@ for M = wh.parameter.M.values
 
                 wh.save([folderpath,experiment_name,'_wh']);
               
-                showCurrentMeasurement('SIR', sir, 'I Att', interference_atten,'BER V', min(ber_vnle),'BER M',min(ber_vnle_mlse),'BER DB',min(ber_db),'Fsym',fsym.*1e-9, 'ROP', rop, 'PAM',M);
-
 
             end
         end
     end
-
 end
 
 close(hWaitbar);

projects/HighSpeedExperiment_2024/offline_dsp/mpi_dsp.m

@@ -0,0 +1,97 @@
+
+% Set your folder path here
+folderPath = 'C:\Users\sioe\Documents\High_Speed_Measurement_2024\mpi_measurement\';
+
+fileTimeStmp = 'testen20241028_163511';
+fileBody = [fileTimeStmp,'_PAM_4_R_336000000000_DB_1_I_atten_'];
+
+fileBits = [folderPath,fileBody,'0_bits'];
+fileSymbols = [folderPath,fileBody,'0_symbols'];
+fileWareHouse = [folderPath,fileTimeStmp,'_wh'];
+
+Bits = load(fileBits,"Bits");Bits = Bits.Bits;
+Symbols = load(fileSymbols);Symbols = Symbols.Symbols;
+wh = load(fileWareHouse);wh = wh.obj;
+M = wh.parameter.M.values(1);
+
+cnt = 1;
+for i_atten = 24%flip([0:3:30,45])
+
+    fileRx = [folderPath,fileBody,num2str(i_atten),'_rx_signal'];
+    S = load(fileRx);S = S.S;
+
+    Scpe_sig_raw = load([folderPath,fileBody,num2str(i_atten),'_raw_signal']);
+    Scpe_sig_raw = Scpe_sig_raw.Scpe_sig_raw;
+
+    Scpe_sig_raw.spectrum("displayname",'Raw signal','fignum',80,'normalizeTo0dB',1);
+
+    eq_values = min(numel(S),8);
+    Noi = cell(eq_values,1);
+    EQ_vnle= cell(eq_values,1);
+    EQ_ffe= cell(eq_values,1);
+
+    %ffe = EQ("Ne",[50,0,0],"Nb",[0,0,0],"training_length",4096*2,"training_loops",5,"dd_loops",5,"K",2,"DCmu",0.05,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
+    vnle = EQ("Ne",[50,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.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1);
+
+    parfor s = 1:eq_values
+        if 1
+            %VNLE
+            Scpe_sig_syncd = S{s};
+            [EQ_vnle{s}] = vnle.process(Scpe_sig_syncd,Symbols);
+            Noi{s} = EQ_vnle{s}-Symbols;
+            Rx_bits = PAMmapper(M,0).demap(EQ_vnle{s});
+            [~,~,ber_vnle(s),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+        end
+
+        %VNLE + MLSE
+        if 1
+            Noi{s}.signal = Noi{s}.signal - mean(Noi{s}.signal);
+            nc = 2;
+            burg_coeff = arburg(Noi{s}.signal,nc);
+            EQ_mlse = EQ_vnle{s}.filter(burg_coeff,1);
+
+            EQ_mlse = MLSE("DIR",burg_coeff,"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_mlse);
+            Rx_bits = PAMmapper(M,0).demap(EQ_mlse);
+            [~,~,ber_vnle_mlse(s),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1);
+        end
+    end
+
+    if 1
+        [~,s]=min(ber_vnle_mlse);
+        nc = 2;
+        burg_coeff = arburg(Noi{s}.signal,nc);
+        Noi{s}.spectrum('displayname','Noise PSD','fignum',123);
+        [h,w] = freqz(1,burg_coeff,length(Noi{s}),"whole",Noi{s}.fs);
+        h = h/max(abs(h));
+        hold on
+        w_ = (w - Noi{s}.fs/2);
+        plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']);
+    end
+
+    disp(['VNLE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle)),' WORST:',sprintf('%.1E',max(ber_vnle)),'. Out of ',num2str(numel(ber_vnle))]);
+    disp(['VNLE+MLSE EQ: BEST BER: ',sprintf('%.1E',min(ber_vnle_mlse)),' AVG BER: ',sprintf('%.1E',mean(ber_vnle_mlse)),' WORST:',sprintf('%.1E',max(ber_vnle_mlse)),'. Out of ',num2str(numel(ber_vnle_mlse))]);
+
+    vnle_result(cnt) = min(ber_vnle);
+    mlse_result(cnt) = min(ber_vnle_mlse);
+    cnt=cnt+1;
+end
+
+
+
+i_atten = flip([0:3:30,45]);
+figure(90)
+plot(i_atten,mlse_result);
+% Continue with the rest of your plot settings
+title('MPI')
+yline(3.8e-3, 'DisplayName', 'HD-FEC', 'LineStyle', '--', 'HandleVisibility', 'off');
+xlabel('Bit Rate in GBps');
+ylabel('Attenuation of I Branch');
+set(gca, 'yscale', 'log');
+set(gca, 'Box', 'on');
+grid on;
+grid minor;
+legend('Interpreter', 'none');
+
+
+
+