From 492f889dce1828684f56f9ab10a6d8dc271068c6 Mon Sep 17 00:00:00 2001 From: Silas Labor Zizou Date: Mon, 4 Nov 2024 08:57:53 +0100 Subject: [PATCH 1/2] Last commit after 400G Lab Measurement --- Classes/04_DSP/Coding/Duobinary.m | 7 + Classes/04_DSP/Equalizer/EQ.m | 2 +- Classes/Warehouse_class/classes/DataStorage.m | 2 +- Functions/showCurrentMeasurement.m | 2 +- Tests/04_DSP/Coding/Duobinary_test.m | 41 +- .../master_auswertung_10km.m | 80 ++++ .../master_evaluation.m | 423 +++++++++++------- .../mpi_measurement.m | 281 ++++++++---- .../offline_dsp/mpi_dsp.m | 97 ++++ 9 files changed, 669 insertions(+), 266 deletions(-) create mode 100644 projects/HighSpeedExperiment_2024/master_auswertung_10km.m create mode 100644 projects/HighSpeedExperiment_2024/offline_dsp/mpi_dsp.m diff --git a/Classes/04_DSP/Coding/Duobinary.m b/Classes/04_DSP/Coding/Duobinary.m index 191766a..5801163 100644 --- a/Classes/04_DSP/Coding/Duobinary.m +++ b/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 diff --git a/Classes/04_DSP/Equalizer/EQ.m b/Classes/04_DSP/Equalizer/EQ.m index a92f827..9f528ab 100644 --- a/Classes/04_DSP/Equalizer/EQ.m +++ b/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 diff --git a/Classes/Warehouse_class/classes/DataStorage.m b/Classes/Warehouse_class/classes/DataStorage.m index 6f0d2ed..06b77b9 100644 --- a/Classes/Warehouse_class/classes/DataStorage.m +++ b/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 diff --git a/Functions/showCurrentMeasurement.m b/Functions/showCurrentMeasurement.m index d97db44..dc316f6 100644 --- a/Functions/showCurrentMeasurement.m +++ b/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 diff --git a/Tests/04_DSP/Coding/Duobinary_test.m b/Tests/04_DSP/Coding/Duobinary_test.m index 33de45e..d6afa42 100644 --- a/Tests/04_DSP/Coding/Duobinary_test.m +++ b/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) - [bitpattern(:,i),seed] = prbs(O,N,seed); + %%%%% MOVE-IT PRMS %%%% + + 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 diff --git a/projects/HighSpeedExperiment_2024/master_auswertung_10km.m b/projects/HighSpeedExperiment_2024/master_auswertung_10km.m new file mode 100644 index 0000000..16b5039 --- /dev/null +++ b/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 diff --git a/projects/HighSpeedExperiment_2024/master_evaluation.m b/projects/HighSpeedExperiment_2024/master_evaluation.m index cf09489..7df54fe 100644 --- a/projects/HighSpeedExperiment_2024/master_evaluation.m +++ b/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.lambda = flip([1293, 1302, 1310, 1318, 1327.4]); %calcWavelengthPlan(16, 400e9 , 1310); - params.bitrate = [300:30:480].*1e9; - params.duobinary = [0,1]; + params.M = [8]; + params.lambda = [1293, 1297.5, 1302, 1306.5, 1310, 1313.4, 1318, 1322.7, 1327.4]; %calcWavelengthPlan(16, 400e9 , 1310); + params.bitrate = [270,300,330,360,390,400,410,420,430,440,450,460,470,480].*1e9; + 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"; @@ -71,95 +75,119 @@ estimatedTimeRemaining = 0; 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; - dcs.set("voltage",[v_bias_for_pam, 0]); - pulsef = 1; + if db == 1 + ffe_only = 0; + 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 - %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 db == 2 - 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; - 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); + 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 = -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 - ffe_only = 0; - 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 + % 1) + pdfa.disablePDFA; - + % 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,27 +210,47 @@ 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 %%%%% Loop Preps iterationStartTime = tic; loopcnt = loopcnt+1; loop_name = ['_PAM_',num2str(M),'_L_',num2str(lambda),'_R_',num2str(bitrate),'_DB_',num2str(db),'_ROP_',num2str(rop_atten)]; - loop_name = strrep(loop_name,'.','_'); + loop_name = strrep(loop_name,'.','_'); %%%%% READ Voltages %%%%%% dcs.readVals(); @@ -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); @@ -241,7 +289,7 @@ for M = wh.parameter.M.values precomp_est.estimate(Scpe_sig_resampled,"save",true,"savePath",precomp_path,"fileName",precomp_fn); precomp_est.plot(); end - + voa.readvals(); rop = voa.power_state(1); @@ -252,8 +300,8 @@ 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 ber_vnle = [-1]; @@ -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); - Noi = cell(eq_values,1); - EQ_vnle= cell(eq_values,1); - EQ_ffe= cell(eq_values,1); - if 0 + + if 1 + + eq_values = min(numel(S),8); + Noi = cell(eq_values,1); + 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,79 +347,152 @@ 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))]); - 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_ffe))]); - - [~,i] = min(ber_vnle); - figure(56); - 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))); + if 0 + [~,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 - 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 - parfor s = 1:numel(S) + if 1 + + 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) - 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); - [EQ_sig, Noi] = ffe.process(S{s},Symbols); - EQ_sig.plot("displayname",'After EQ','fignum',112); - 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); - [~,num_errors,ber_db(s),pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - end + 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); - 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 + + 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',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 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" %%%%%% - + wh.addValueToStorage({ber_ffe},'ber_ffe',M,lambda,bitrate,db,rop_atten); wh.addValueToStorage({ber_ffe_mlse},'ber_ffe_mlse',M,lambda,bitrate,db,rop_atten); wh.addValueToStorage({ber_vnle},'ber_vnle',M,lambda,bitrate,db,rop_atten); @@ -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)); diff --git a/projects/HighSpeedExperiment_2024/mpi_measurement.m b/projects/HighSpeedExperiment_2024/mpi_measurement.m index 37e1b04..754c99f 100644 --- a/projects/HighSpeedExperiment_2024/mpi_measurement.m +++ b/projects/HighSpeedExperiment_2024/mpi_measurement.m @@ -1,6 +1,6 @@ -folderpath = 'C:\Users\sioe\Documents\High_Speed_Measurement_2024\mpi_measurement\'; -experiment_name = 'testen'; +folderpath = 'C:\Users\sioe\Documents\High_Speed_Measurement_2024\MPI_duobinary_encoded\'; +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.bitrate = [336].*1e9;%[90:60:480].*1e9;%[300:30:480].*1e9; -params.duobinary = [0]; -params.interference_atten = [45]; +params.M = [4]; +params.bitrate = [336,360,390,420,448].*1e9;%[90:60:480].*1e9;%[300:30:480].*1e9; +params.duobinary = [2]; +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); - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - - % 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); - + save([folderpath,experiment_name,loop_name,'_raw_signal'],"Scpe_sig_raw"); 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; - Noi{1}.spectrum('displayname',['Noise; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1); - burg_coeff = arburg(Noi{1}.signal,nc); - [h,w] = freqz(1,burg_coeff,length(Noi{1}),"whole",Noi{1}.fs); + + 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{1}.fs/2); - plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',[num2str(nc), ' burg; SIR:',sprintf('%.2f',sir)]); - drawnow; + w_ = (w - Noi{s}.fs/2); + plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']); 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(['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 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.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 - Noi{1}.spectrum('displayname',['Noise; SIR:',sprintf('%.2f',sir)],'fignum',40,'normalizeTo0dB',1); - + if 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); @@ -401,14 +489,11 @@ for M = wh.parameter.M.values loopcnt, looptotal, sum(iterationTimes(1:loopcnt))/60, averageTimePerIteration, estimatedTimeRemaining/60 )); 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); diff --git a/projects/HighSpeedExperiment_2024/offline_dsp/mpi_dsp.m b/projects/HighSpeedExperiment_2024/offline_dsp/mpi_dsp.m new file mode 100644 index 0000000..cbce0f1 --- /dev/null +++ b/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'); + + + + From 9ed763d094c08ea7cbc799610e86fe1ccb93d7ab Mon Sep 17 00:00:00 2001 From: sioe Date: Mon, 4 Nov 2024 09:00:51 +0100 Subject: [PATCH 2/2] CHnages from silas PC --- Classes/00_signals/Signal.m | 8 +- Classes/01_transmit/ChannelFreqResp.m | 12 +- projects/400G_FTN_setups/imdd_mpi_dsp.m | 132 ++++++++---------- .../offline_dsp_analysis/load_n_dsp.m | 5 +- .../load_n_plot_transfer_characteristic.m | 8 ++ 5 files changed, 80 insertions(+), 85 deletions(-) create mode 100644 projects/Lab_2024/offline_dsp_analysis/load_n_plot_transfer_characteristic.m diff --git a/Classes/00_signals/Signal.m b/Classes/00_signals/Signal.m index a2f70ae..c2a7244 100644 --- a/Classes/00_signals/Signal.m +++ b/Classes/00_signals/Signal.m @@ -344,14 +344,14 @@ classdef Signal % spectrum_plot(obj.signal,options.fsamp,options.figurename,options.displayname); - N = 2^(nextpow2(length(obj.signal))-10); + N = 2^(nextpow2(length(obj.signal))-2); 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'); + % p_lin = smooth(p_lin,0.05,'rloess'); end @@ -376,10 +376,10 @@ classdef Signal if options.normalizeToNyquist==0 xlabel("Frequency in GHz"); - %xlim([-obj.fs/2 obj.fs/2].*1e-9) + % xlim([-obj.fs/2 obj.fs/2].*1e-9) edgetick = 2^(nextpow2(obj.fs*1e-9)); xticks([-edgetick:16:edgetick]); - xlim([100*round( min(w.*1e-9)/100,1)-10,100*round( max(w.*1e-9)/100,1)+10]) + xlim([100*round( min(w)/100,1)-10,100*round( max(w)/100,1)+10]) else xlabel("Normalized Frequency"); xlim([-pi, pi]); diff --git a/Classes/01_transmit/ChannelFreqResp.m b/Classes/01_transmit/ChannelFreqResp.m index e77cf43..f204112 100644 --- a/Classes/01_transmit/ChannelFreqResp.m +++ b/Classes/01_transmit/ChannelFreqResp.m @@ -241,10 +241,14 @@ classdef ChannelFreqResp < handle xlim([0.2 .5*max(obj.faxis)*1e-9]); grid on; %%% plot for publication - figure(98989);hold all;box on;title('Magnitude Freq. Response'); - xlim([0.2 .5*max(obj.faxis)*1e-9]); - ylim([-20, 2]); - plot(obj.faxis/1e9, 20*log10(abs(Havg)),'LineWidth',2); + figure(1234);hold all;box on;title('Magnitude Freq. Response'); + %xlim([0.2 .5*max(obj.faxis)*1e-9]); + ylim([-40, 2]); + Havg_smooth = smooth(Havg,50); + symaxis = (obj.faxis-(obj.f_ref/2))/1e9; + Havg = fftshift(Havg); + Havg = smooth(Havg); + plot(symaxis, 20*log10(abs(Havg)),'LineWidth',0.5); grid on; diff --git a/projects/400G_FTN_setups/imdd_mpi_dsp.m b/projects/400G_FTN_setups/imdd_mpi_dsp.m index d9464c2..2016158 100644 --- a/projects/400G_FTN_setups/imdd_mpi_dsp.m +++ b/projects/400G_FTN_setups/imdd_mpi_dsp.m @@ -5,17 +5,17 @@ params = struct; params.M = [4]; params.datarate = [448]; params.rop = [0]; -params.sir = 40;%15:1:40; +params.sir = 45; params.random_key_laser_phase = 10:20; precomp_mode = 0; %0=do nothing ; 1= measure; 2=precomp active -postfilter = 0; % noise whiten. approach -> Postfilter + MLSE +postfilter = 1; % noise whiten. approach -> Postfilter + MLSE -db_precode = 1; +db_precode = 0; db_encode = 0; -db_channelapproach = 1; +db_channelapproach = 0; -laser_linewidth = 50e5; +laser_linewidth = 5e5; random_key_sequence = 15; random_key_laser_phase = 66; sir = 20; @@ -61,14 +61,14 @@ for M = wh.parameter.M.values % MAIN SIGNAL %%%%% Symbol Generation MAIN %%%%%% - [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",18,"useprbs",1,... + [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",19,"useprbs",1,... "fs_out",M8199.fdac,"applyclipping",0,"clipfactor",1.5,... "applypulseform",0,"pulseformer",Pform,"randkey",random_key_sequence,... "db_precode",db_precode,"db_encode",db_encode,... "mrds_code",usemrds,"mrds_blocklength",512).process(); %%%%% Symbol Generation INTERFERENCE %%%%%% - [Digi_sig_I,Symbols_I,Bits_I] = PAMsource("fsym",fsym,"M",M,"order",18,"useprbs",0,... + [Digi_sig_I,Symbols_I,Bits_I] = PAMsource("fsym",fsym,"M",M,"order",19,"useprbs",0,... "fs_out",M8199.fdac,"applyclipping",0,"clipfactor",1.5,... "applypulseform",0,"pulseformer",Pform,"randkey",random_key_sequence+1,... "db_precode",db_precode,"db_encode",db_encode,... @@ -185,8 +185,8 @@ for M = wh.parameter.M.values end % Scpe_sig_normalized = Scpe_sig.normalize("mode","rms"); - - % Scpe_sig.normalize("mode","rms").spectrum("displayname",'After Scope','fignum',10); + % + % Scpe_sig_normalized.normalize("mode","rms").spectrum("displayname",'After Scope','fignum',23); %%%%%% Sample to 2x fsym %%%%%% Scpe_sig = Scpe_sig.resample("fs_in",fadc,"fs_out",2*fsym); @@ -213,11 +213,17 @@ for M = wh.parameter.M.values if db_channelapproach % ref symbols and transm. sequence are precoded [EQ_sig, Noi] = Eq.process(Scpe_sig,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 = Duobinary().decode(EQ_sig); Rx_bits = PAMmapper(M,0).demap(EQ_sig); [~,~,ber_vnle(i,j),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - + + + EQ_sig.spectrum('displayname','EQ DB Out','fignum',12345,'normalizeTo0dB',0,'normalizeToNyquist',0,'color',cols(3,:)); + Noi.spectrum('displayname','Noise PSD optimal','fignum',1234,'normalizeTo0dB',0,'normalizeToNyquist',0,'color',cols(4,:)); + + elseif db_encode [EQ_sig, Noi] = Eq.process(Scpe_sig,Symbols); @@ -227,21 +233,39 @@ for M = wh.parameter.M.values elseif postfilter [EQ_sig, Noi] = Eq.process(Scpe_sig,Symbols); - - % EQ_sig.plot("displayname",'After VNLE','fignum',90,'clear',1); - % Quantization is too far from orig. symbols -> - % error psd is quite different - % Sym_ = PAMmapper(M,0).quantize(EQ_sig); - % Noi_ = Sym_-EQ_sig; - % Noi_.normalize('mode','rms').spectrum('displayname','Noise PSD','fignum',1234,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(nc+1,:)); + %%% REMOVE DC peak from Noi PSD + S = Noi.signal; + N1 = 1001; + + % recursion + % Initialize the moving sum for the first window + half_window = (N1 - 1) / 2; + moving_sum = sum(S(1:N1)); + + % Calculate the first element of R1 + S_(1:half_window+1) = S(1:half_window+1) - (moving_sum / N1); + + % Loop over the signal and apply the recursive moving average subtraction + for n = (half_window+2):(length(S)-half_window) + % Update the moving sum by subtracting the oldest value and adding the new one + moving_sum = moving_sum - S(n-half_window-1) + S(n+half_window); + + % Calculate the new value of R1 + S_(n) = S(n) - (moving_sum / N1); + end + + S_(n+1:length(S)) = S(n+1:end) - (moving_sum / N1); + Noi.signal = S_; + %%% END REMOVE DC PEAK %%% Rx_bits = PAMmapper(M,0).demap(EQ_sig); [~,~,ber_vnle(i,j),~] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - EQ_sig.normalize('mode','rms').spectrum('displayname','EQ Out','fignum',1234,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(nc,:)); - Noi.normalize('mode','rms').spectrum('displayname','Noise PSD optimal','fignum',1234,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(nc+1,:)); + EQ_sig.spectrum('displayname','EQ Out','fignum',12345,'normalizeTo0dB',0,'normalizeToNyquist',0,'color',cols(1,:)); + + + Noi.spectrum('displayname','Noise PSD optimal','fignum',22,'normalizeTo0dB',1,'normalizeToNyquist',0,'color',cols(2,:)); for nc = 1:3 @@ -258,24 +282,16 @@ for M = wh.parameter.M.values if 1 cols = linspecer(12); + % EQ_sig_filt.normalize('mode','rms').spectrum('displayname','Noise PSD','fignum',1234,'normalizeTo0dB',1,'normalizeToNyquist',0,'color',cols(nc+2,:)); - - EQ_sig_filt.normalize('mode','rms').spectrum('displayname','Noise PSD','fignum',1234,'normalizeTo0dB',1,'normalizeToNyquist',1,'color',cols(nc+2,:)); - - % [h,w] = freqz(1,burg_coeff,length(Noi),"whole",Noi.fs); - % h = h/max(abs(h)); - % hold on - % w_ = (w - Noi.fs/2); - % figure(123) - % plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']); - - [h,w] = freqz(1,burg_coeff,length(Noi),"whole"); + [h,w] = freqz(1,burg_coeff,length(Noi),"whole",Noi.fs); + % h = 1./h; h = h/max(abs(h)); hold on - w_ = (w - pi); - plot(w_,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']); + w_ = (w - Noi.fs/2); + figure(22) + plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']); - end Rx_bits = PAMmapper(M,0).demap(EQ_sig_mlse); @@ -286,41 +302,7 @@ for M = wh.parameter.M.values else - % S = Scpe_sig.signal; - % N1 = 101; - % - % % Initialize the running sum with the first window's sum - % running_sum = mean( S(1:N1) ); - % - % % Calculate the first output value - % S_(1) = S(1) - running_sum; - % - % % Recursive running sum filter - % for n = 2 : length(S) - N1 - % % Update running sum by removing the oldest sample and adding the newest - % avg_win(n) = mean( S(n:n+N1) ); - % S_(n) = S(n) - avg_win(n); - % end - % - % % movmean - % S__ = S - movmean(S,[floor(N1/2),ceil(N1/2)]); - % - % % recursion - % % Initialize the moving sum for the first window - % half_window = (N1 - 1) / 2; - % moving_sum = sum(S(1:N1)); - % - % % Calculate the first element of R1 - % S___(half_window+1) = S(half_window+1) - (moving_sum / N1); - % - % % Loop over the signal and apply the recursive moving average subtraction - % for n = (half_window+2):(length(S)-half_window) - % % Update the moving sum by subtracting the oldest value and adding the new one - % moving_sum = moving_sum - S(n-half_window-1) + S(n+half_window); - % - % % Calculate the new value of R1 - % S___(n) = S(n) - (moving_sum / N1); - % end + [EQ_sig, Noi] = Eq.process(Scpe_sig,Symbols); @@ -379,16 +361,16 @@ ber_mlse=squeeze(mean(ber_mlse,1)); ber_vnle = mean(ber_vnle,1); % Create the initial plot -figure(44); +figure(466); a = gca; hold on; % Retain the plot so new points can be added without complete redraw dispname = ['Lw: ',num2str(laser_linewidth.*1e-6),' MHz']; - -plot(wh.parameter.sir.values,ber_vnle,"LineWidth",0.5,"LineStyle","-","Marker",".","MarkerSize",15,"DisplayName",['VNLE ',dispname]); -plot(wh.parameter.sir.values,ber_mlse(:,1),"LineWidth",0.5,"LineStyle","-","Marker",".","MarkerSize",15,"DisplayName",['MLSE 1 ',dispname]); -plot(wh.parameter.sir.values,ber_mlse(:,2),"LineWidth",0.5,"LineStyle","-","Marker",".","MarkerSize",15,"DisplayName",['MLSE 2',dispname]); -plot(wh.parameter.sir.values,ber_mlse(:,3),"LineWidth",0.5,"LineStyle","-","Marker",".","MarkerSize",15,"DisplayName",['MLSE 3',dispname]); +cols = linspecer(6); +plot(wh.parameter.sir.values,ber_vnle,"LineWidth",0.5,"LineStyle","--","Marker",".","MarkerSize",15,"DisplayName",['PAM',num2str(M),' VNLE ',dispname],'Color',cols(1,:)); +plot(wh.parameter.sir.values,ber_mlse(:,1),"LineWidth",0.5,"LineStyle","--","Marker",".","MarkerSize",15,"DisplayName",['PAM',num2str(M),'MLSE 1 ',dispname],'Color',cols(2,:)); +plot(wh.parameter.sir.values,ber_mlse(:,2),"LineWidth",0.5,"LineStyle","--","Marker",".","MarkerSize",15,"DisplayName",['PAM',num2str(M),'MLSE 2',dispname],'Color',cols(3,:)); +plot(wh.parameter.sir.values,ber_mlse(:,3),"LineWidth",0.5,"LineStyle","--","Marker",".","MarkerSize",15,"DisplayName",['PAM',num2str(M),'MLSE 3',dispname],'Color',cols(4,:)); yline(3.8e-3,'DisplayName','HD-FEC','LineStyle','--','HandleVisibility','off'); xlabel('Received Optical Power (dBm)'); diff --git a/projects/Lab_2024/offline_dsp_analysis/load_n_dsp.m b/projects/Lab_2024/offline_dsp_analysis/load_n_dsp.m index fe8434f..784127b 100644 --- a/projects/Lab_2024/offline_dsp_analysis/load_n_dsp.m +++ b/projects/Lab_2024/offline_dsp_analysis/load_n_dsp.m @@ -1,13 +1,14 @@ % load data points if 0 - foldername = '/Users/silasoettinghaus/Nextcloud/Dokumente/02_Ablage_Office/Lab_Data_24/sir_sweep_sd40/'; + foldername = 'C:\Users\Silas\Nextcloud4\Dokumente\02_Ablage_Office\Lab_Data_24\sir_sweep_sd40'; filename = 'PAM4_56_'; filename = 'PAM4_56_v2'; else - foldername = '/Users/silasoettinghaus/Nextcloud/Dokumente/02_Ablage_Office/Lab_Data_24/sir_sweep_sd40/DB_coded'; + foldername = 'C:\Users\Silas\Nextcloud4\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; diff --git a/projects/Lab_2024/offline_dsp_analysis/load_n_plot_transfer_characteristic.m b/projects/Lab_2024/offline_dsp_analysis/load_n_plot_transfer_characteristic.m new file mode 100644 index 0000000..a23d25a --- /dev/null +++ b/projects/Lab_2024/offline_dsp_analysis/load_n_plot_transfer_characteristic.m @@ -0,0 +1,8 @@ + + +precomp_path = "C:\Users\Silas\Documents\MATLAB\imdd_simulation\projects\standard_system"; +precomp_filename = "lab_mpi_setup_2"; + +freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',92e9); +freqresp.load('loadPath',precomp_path,'fileName',precomp_filename); +freqresp.plot();