From cf4e0f2b12e55a2fb1610af91bbdea93afabfb3b Mon Sep 17 00:00:00 2001 From: Silas Labor Zizou Date: Tue, 29 Oct 2024 14:38:22 +0100 Subject: [PATCH] measurement state --- Classes/00_signals/Signal.m | 4 +- Classes/01_transmit/ChannelFreqResp.m | 3 +- Classes/04_DSP/Coding/Duobinary.m | 4 +- Classes/05_Lab/Exfo_laser.m | 2 +- Classes/05_Lab/OptAtten.m | 31 +- Classes/05_Lab/ScopeKeysight.m | 7 +- Classes/Warehouse_class/classes/DataStorage.m | 36 +- Functions/holdAndShowValue.m | 48 ++ Functions/showCurrentMeasurement.m | 2 +- Functions/waitUntilClick.m | 20 + .../bias_evaluation.m | 147 +++- .../bias_optimization.m | 633 ++++++++++-------- .../master_evaluation.m | 433 ++++++++++++ .../mpi_measurement.m | 418 ++++++++++++ 14 files changed, 1469 insertions(+), 319 deletions(-) create mode 100644 Functions/holdAndShowValue.m create mode 100644 Functions/waitUntilClick.m create mode 100644 projects/HighSpeedExperiment_2024/master_evaluation.m create mode 100644 projects/HighSpeedExperiment_2024/mpi_measurement.m diff --git a/Classes/00_signals/Signal.m b/Classes/00_signals/Signal.m index a3a5007..500eb76 100644 --- a/Classes/00_signals/Signal.m +++ b/Classes/00_signals/Signal.m @@ -27,7 +27,7 @@ classdef Signal [~,obj.gitSHA] = system('git rev-parse HEAD'); [~,obj.gitStatus] = system('git status --porcelain'); - [~,obj.gitPatch] = system('git diff'); + % [~,obj.gitPatch] = system('git diff'); %%% Stuff for Logbook %%% SignalType = []; @@ -344,7 +344,7 @@ classdef Signal % spectrum_plot(obj.signal,options.fsamp,options.figurename,options.displayname); - N = 2^(nextpow2(length(obj.signal))-8); + N = 2^(nextpow2(length(obj.signal))-10); if options.normalizeToNyquist==0 [p_lin,w] = pwelch(obj.signal,hanning(N),N/2,N,obj.fs,"centered","power","mean"); diff --git a/Classes/01_transmit/ChannelFreqResp.m b/Classes/01_transmit/ChannelFreqResp.m index 5bf2ede..d003650 100644 --- a/Classes/01_transmit/ChannelFreqResp.m +++ b/Classes/01_transmit/ChannelFreqResp.m @@ -241,10 +241,11 @@ classdef ChannelFreqResp < handle xlim([0.2 .5*max(obj.faxis)*1e-9]); grid on; %%% plot for publication - figure(20);hold on;box on;title('Magnitude Freq. Response'); + figure(30);hold on;box on;title('Magnitude Freq. Response'); % xlim([0 max(obj.faxis)*1e-9]); % ylim([-20, 10]); fax = obj.faxis - obj.f_ref/2; + Havg = Havg ./ max(abs(Havg)); plot(fax/1e9, 20*log10(abs(fftshift(Havg)))+7,'LineWidth',2); grid on; diff --git a/Classes/04_DSP/Coding/Duobinary.m b/Classes/04_DSP/Coding/Duobinary.m index 259e04d..191766a 100644 --- a/Classes/04_DSP/Coding/Duobinary.m +++ b/Classes/04_DSP/Coding/Duobinary.m @@ -164,11 +164,11 @@ classdef Duobinary elseif I == 11 %todo data = data .* sqrt(5.8); - warning('Check if PAM16 implementation, mapping and scaling is correct!') + warning('Check db decode implementation, mapping and scaling is correct!') elseif I == 15 data = data .* sqrt(10.5); elseif I == 16 - warning('Check if PAM16 implementation, mapping and scaling is correct!') + warning('Check db decode implementation, mapping and scaling is correct!') end data = round(data); diff --git a/Classes/05_Lab/Exfo_laser.m b/Classes/05_Lab/Exfo_laser.m index c6ae8bc..ce8f1d3 100644 --- a/Classes/05_Lab/Exfo_laser.m +++ b/Classes/05_Lab/Exfo_laser.m @@ -52,7 +52,7 @@ classdef Exfo_laser < handle if obj.safety_mode warning("Safety_mode ON: Display all information. Ask when switching laser. Turn safety_mode off in class instance or during initialization Exfo_laser(...,'safety_mode',0)"); else - warning("safety_mode OFF: EVERYTHING IS EXECUTED WITHOUT ASKING :-) BE SHURE WHAT YOU DO"); + warning("safety_mode OFF: EVERYTHING IS EXECUTED WITHOUT ASKING :-) BE SURE WHAT YOU DO"); end end diff --git a/Classes/05_Lab/OptAtten.m b/Classes/05_Lab/OptAtten.m index 93deccc..8d65e8a 100644 --- a/Classes/05_Lab/OptAtten.m +++ b/Classes/05_Lab/OptAtten.m @@ -53,8 +53,11 @@ classdef OptAtten < handle try %connect to device + + v = visadev('TCPIP::134.245.243.248::INSTR'); + debug = 0; if debug disp(['Connected to Instrument: ',char(v.Vendor),' ',char(v.Model),' SerNo:',char(v.SerialNumber)]); @@ -213,7 +216,7 @@ classdef OptAtten < handle ch_info_txt = ['OptAtten: ',num2str(i),' -> Attenuated by: ',num2str(options.value(i)),' dB; Cur Output: ',num2str(state.outputpower(i)),'dBm']; disp(ch_info_txt); end - + end end @@ -258,11 +261,6 @@ classdef OptAtten < handle answer = sscanf(line,'%f'); obj.power_state(cnt) = answer; - writeline(v, [':READ' num2str(s) ':POW?']); - line = readline(v); - answer = sscanf(line,'%f'); - obj.power_state(cnt) = answer; - writeline(v, [':INP' num2str(s) ':ATT?']); line = readline(v); answer = sscanf(line,'%f'); @@ -273,7 +271,7 @@ classdef OptAtten < handle answer = sscanf(line,'%f'); obj.wavelength_state(cnt) = answer; - writeline(v, [':INP' num2str(s) ':WAV?']); + writeline(v, [':INP' num2str(s) ':ATT:SPE?']); line = readline(v); answer = sscanf(line,'%f'); obj.speed_state(cnt) = answer; @@ -294,6 +292,25 @@ classdef OptAtten < handle end end + function [p1,p2,p3,p4]=readPower(obj) + + v = visadev('TCPIP::134.245.243.248::INSTR'); + cnt = 1; + for s = 1:2:7 + writeline(v, [':READ' num2str(s) ':POW?']); + line = readline(v); + answer = sscanf(line,'%f'); + obj.power_state(cnt) = answer; + + cnt = cnt+1; + end + p1 = obj.power_state(1); + p2 = obj.power_state(2); + p3 = obj.power_state(3); + p4 = obj.power_state(4); + + end + end end diff --git a/Classes/05_Lab/ScopeKeysight.m b/Classes/05_Lab/ScopeKeysight.m index 8eee6ba..fc4b336 100644 --- a/Classes/05_Lab/ScopeKeysight.m +++ b/Classes/05_Lab/ScopeKeysight.m @@ -131,7 +131,7 @@ classdef ScopeKeysight %use this to finetune autoscaling - VERY helpful % higher value leads to higher scaling - fintunefactor = 1.4; % within [1,...,2] + fintunefactor = 1.5; % within [1,...,2] obj.writeNcheck(v,sprintf(':CHANnel%u:RANGe %.3f',n,range/fintunefactor)); end else @@ -141,6 +141,8 @@ classdef ScopeKeysight % After Autoscale, let the user adjust the scope scaling... % "options.waitUntilClick" was a shit name but now this is it :-) if options.waitUntilClick + + if 0 % Create a dialog box with the desired text d = dialog('Name', 'Scale Scope then press continue', 'Position', [300, 300, 300, 180]); @@ -160,6 +162,9 @@ classdef ScopeKeysight % Pause execution until the dialog box is closed uiwait(d); + else + holdAndShowValue + end end if obj.extRef diff --git a/Classes/Warehouse_class/classes/DataStorage.m b/Classes/Warehouse_class/classes/DataStorage.m index dd9eb4b..6f0d2ed 100644 --- a/Classes/Warehouse_class/classes/DataStorage.m +++ b/Classes/Warehouse_class/classes/DataStorage.m @@ -125,23 +125,49 @@ classdef DataStorage < handle lin_idx = obj.getIndicesByPhys(varargin); errcnt = 0; for i=1:numel(lin_idx) - try + tmp = obj.sto.(storageVarName){lin_idx(i)}; if ~isempty(tmp) - if isa(tmp,'Signal') + + if isa(tmp,'Signal') || isa(tmp,'struct') || isa(tmp,'Exfo_laser') if i == 1 value = {}; end value{i} = tmp ; + elseif isa(tmp,'cell') if isa(tmp{1},'Signal') if i == 1 value = {}; end value{i} = tmp{1} ; + else + value{i} = tmp{1} ; end + else - value(i,:) = tmp ; + + try + value(i,:) = tmp ; + catch + % value(i,:) = tmp(1:size(value,2)) ; + + if size(value,2) < size(tmp,2) + + diff = size(tmp,2) - size(value,2); + value(:,end+1:end+diff) = NaN(size(value,1),diff); + value(i,:) = tmp ; + + elseif size(value,2) > size(tmp,2) + + diff = size(value,2) - size(tmp,2); + tmp(:,end+1:end+diff) = NaN(1,diff); + value(i,:) = tmp ; + + end + + end + end else errcnt = errcnt+1; @@ -166,9 +192,7 @@ classdef DataStorage < handle % warning([num2str(errcnt),' requested datapoint(s) not in warehouse.']); end - catch - error('Error in Datastorage: Something happened while looking up in warehouse.') - end + end else error('Wrong Request using ExampleWarehouse.getStoValue(*parameter set*). Give me all the Parameters! Please!') diff --git a/Functions/holdAndShowValue.m b/Functions/holdAndShowValue.m new file mode 100644 index 0000000..e45c588 --- /dev/null +++ b/Functions/holdAndShowValue.m @@ -0,0 +1,48 @@ +function holdAndShowValue() + % Create the UI figure + fig = uifigure('Name', 'Multi-Channel Power Monitor', 'Position', [100 100 600 150]); + + % Create labels to display the power values for each channel in a horizontal layout + powerLabels = gobjects(4, 1); + for i = 1:4 + powerLabels(i) = uilabel(fig, 'Position', [50 + (i-1)*130, 70, 120, 30], ... + 'FontSize', 18, 'HorizontalAlignment', 'center'); + powerLabels(i).Text = sprintf('CH %d: Fetching...', i); + end + + % Create the "OK" button + okButton = uibutton(fig, 'push', 'Text', 'OK', 'Position', [250 20 100 40], ... + 'ButtonPushedFcn', @(src, event) closeWindow()); + + % Initialize the timer + updateTimer = timer('ExecutionMode', 'fixedRate', 'Period', 1, ... + 'TimerFcn', @(src, event) updatePowerValues()); + + % Start the timer + start(updateTimer); + + % Function to close the window and stop the timer + function closeWindow() + stop(updateTimer); + delete(updateTimer); + delete(fig); + end + + % Function to fetch and update power values for all channels + function updatePowerValues() + powerValues = fetchPowerValues(); % Replace with your data-fetching function + for i = 1:4 + powerLabels(i).Text = sprintf('CH %d: %.2f dBm', i, powerValues(i)); + end + end + + % Wait until the figure is closed to resume execution + waitfor(fig); +end + +% Example function to simulate fetching power values for four channels +function powerValues = fetchPowerValues() + % Replace this with actual code to get the power values from the power meter + [p1,p2,p3,p4] = OptAtten().readPower(); + powerValues = [p1,p2,p3,p4]; +end diff --git a/Functions/showCurrentMeasurement.m b/Functions/showCurrentMeasurement.m index 7205f6f..d97db44 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(strcmpi(varNameLower, {'ber','mlse','ffe','db'})) + if any(contains(varNameLower, {'ber','mlse','ffe','db'})) % Format 'ber' values in exponential notation with two decimal places values{i} = sprintf('%.2e', values{i}); else diff --git a/Functions/waitUntilClick.m b/Functions/waitUntilClick.m new file mode 100644 index 0000000..6e06737 --- /dev/null +++ b/Functions/waitUntilClick.m @@ -0,0 +1,20 @@ +function waitUntilClick() + % Create the UI figure + fig = uifigure('Name', 'Pause Execution', 'Position', [100 100 300 150]); + + % Create a label to inform the user + uilabel(fig, 'Position', [50 80 200 40], 'Text', 'Click "Continue" to proceed', ... + 'FontSize', 14, 'HorizontalAlignment', 'center'); + + % Create the "Continue" button + continueButton = uibutton(fig, 'push', 'Text', 'Continue', 'Position', [100 20 100 40], ... + 'ButtonPushedFcn', @(src, event) closeWindow()); + + % Function to close the window when "Continue" is clicked + function closeWindow() + delete(fig); % Close the figure window + end + + % Wait until the figure is closed to resume execution + waitfor(fig); +end \ No newline at end of file diff --git a/projects/HighSpeedExperiment_2024/bias_evaluation.m b/projects/HighSpeedExperiment_2024/bias_evaluation.m index 2c27189..babd881 100644 --- a/projects/HighSpeedExperiment_2024/bias_evaluation.m +++ b/projects/HighSpeedExperiment_2024/bias_evaluation.m @@ -1,41 +1,125 @@ -% filename = "C:\Users\sioe\Documents\High_Speed_Measurement_2024\bias_testing\PAM4_b2b_bias_sweep_20241023_130342_wh.mat"; -% -% a = load(filename); -% wh = a.obj; +filename = "C:\Users\sioe\Documents\High_Speed_Measurement_2024\bias_5km\PAMX_5km_20241025_204334_wh.mat"; + +a = load(filename); +wh = a.obj; v_bias_vals = wh.parameter.vbias.values; awg_vpp_vals = wh.parameter.awg_vpp.values; precomp_amp_max_vals = wh.parameter.precomp_amp_max.values; rop_atten_vals = wh.parameter.rop_atten.values; -PAM = wh.parameter.M.values(1); +lambda_vals = wh.parameter.lambda.values; +M_vals = wh.parameter.M.values; +fsym_vals = [168e9, 144e9, 120e9]; -m = wh.getStoValue('m',v_bias_vals(1),awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),PAM); - - - -ber_ffe= []; -ber_mlse= []; -rop_measured= []; -pd_in_measured= []; +ber_ffe = []; +ber_mlse = []; +rop_measured = []; +pd_in_measured = []; rop_measured = []; cnt = 0; -for PAM = wh.parameter.M.values - cnt = cnt+1; - ber_ffe(cnt,:) = wh.getStoValue('ber_ffe',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),PAM)'; - ber_mlse(cnt,:) = wh.getStoValue('ber_mlse',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),PAM); - rop_measured(cnt,:) = wh.getStoValue('rop',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),PAM); - pd_in_measured(cnt,:) = wh.getStoValue('pd_in',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),PAM); + +figure(252) +clf +hold on +cols = cbrewer2('Set1',3); +for l = 1:numel(lambda_vals) + for m = 1:numel(M_vals) + + ber_ffe = wh.getStoValue('ber_ffe',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),M_vals(m),lambda_vals(l)); + ber = wh.getStoValue('ber_collect',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),M_vals(m),lambda_vals(l)); + exfo = wh.getStoValue('exfo',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),M_vals(m),lambda_vals(l)); + + for e = 1:numel(exfo) + laser_pow(e) = exfo{e}.cur_power; + end + + rop_measured = wh.getStoValue('rop',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),M_vals(m),lambda_vals(l)); + pd_in_measured(l,m,:) = wh.getStoValue('pd_in',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),M_vals(m),lambda_vals(l)); + + rx_logbook = wh.getStoValue('rx_logbook',v_bias_vals(1),awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),M_vals(1),lambda_vals(1)); + + + subplot(1,3,l) + hold on + a = scatter(v_bias_vals,min(ber,[],2),40,'LineWidth',2,'Marker','.','DisplayName',['PAM ',num2str(M_vals(m))],'MarkerEdgeColor',cols(m,:)); + title([num2str(lambda_vals(l)),'nm']) + + a.DataTipTemplate.DataTipRows(1).Label = 'Vbias'; + + a.DataTipTemplate.DataTipRows(2).Label = 'BER'; + a.DataTipTemplate.DataTipRows(2).Format ='%.1e'; + + a.DataTipTemplate.DataTipRows(3).Label = 'P_{out}'; + a.DataTipTemplate.DataTipRows(3).Value = rop_measured; + a.DataTipTemplate.DataTipRows(3).Format = ['auto']; + + a.DataTipTemplate.DataTipRows(4).Label = 'Baudr'; + a.DataTipTemplate.DataTipRows(4).Value = repmat(fsym_vals(m).*1e-9,size(ber_ffe)); + a.DataTipTemplate.DataTipRows(4).Format = ['%d',' GBd']; + + a.DataTipTemplate.DataTipRows(5).Label = 'L_{out}'; + a.DataTipTemplate.DataTipRows(5).Value = laser_pow; + a.DataTipTemplate.DataTipRows(5).Format = ['auto']; + + % Polynomial fit (e.g., second-order polynomial) + [woutliers,n] = rmoutliers( min(ber,[],2) ); + p = polyfit( v_bias_vals(~n), log10(woutliers), 4); % Adjust order as needed + BER_fit = polyval(p, v_bias_vals); + + + % Plot the fitted curve + plot(v_bias_vals, 10.^(BER_fit), '-r', 'LineWidth', 1.5,'Color',cols(m,:),'HandleVisibility','off'); + + % Continue with the rest of your plot settings + yline(3.8e-3, 'DisplayName', 'HD-FEC', 'LineStyle', '--', 'HandleVisibility', 'off'); + xlabel('Bias Voltage'); + ylabel('Bit Error Rate (BER)'); + sgtitle('Bit Error Rate vs. ROP'); + set(gca, 'yscale', 'log'); + set(gca, 'Box', 'on'); + grid on; + grid minor; + legend('Interpreter', 'none'); + + ylim([1e-3,0.5]); + xlim([-16, -2]); + + ylim([1e-3,0.5]); + xlim([min(v_bias_vals) max(v_bias_vals)]); + + end end +filename = "C:\Users\sioe\Documents\High_Speed_Measurement_2024\bias_testing_and_b2b\PAM4_b2b_bias_sweep_20241023_191202_wh_BB_BIAS_FINAL.mat"; + +a = load(filename); +wh = a.obj; + +v_bias_vals = wh.parameter.vbias.values; +awg_vpp_vals = wh.parameter.awg_vpp.values; +precomp_amp_max_vals = wh.parameter.precomp_amp_max.values; +rop_atten_vals = wh.parameter.rop_atten.values; +M_vals = wh.parameter.M.values; + + +ber_ffe = []; +ber_mlse = []; +rop_measured = []; +pd_in_measured = []; figure(2024) for i = 1:3 - + + ber_ffe(i,:) = wh.getStoValue('ber_ffe',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),M_vals(i)); + + rop_measured(i,:) = wh.getStoValue('rop',v_bias_vals,awg_vpp_vals(1),precomp_amp_max_vals(1),rop_atten_vals(1),M_vals(i)); + + [bestber,bestindex] = min(ber_ffe(i,:),[],'all'); [awg_pos,v_bias_pos]=ind2sub(size(ber_ffe(i,:)),bestindex); bestawgvpp=awg_vpp_vals(awg_pos); @@ -43,24 +127,35 @@ for i = 1:3 disp(['Best Vpp: ',num2str(bestvbias),' V; Best Vpp AWG: ',num2str(bestawgvpp),' V' ]); + % Polynomial fit (e.g., second-order polynomial) + [woutliers,n] = rmoutliers( ber_ffe(i,:) ); + p = polyfit( v_bias_vals(~n), log10(woutliers), 8); % Adjust order as needed + BER_fit = polyval(p, v_bias_vals); + + % Plot the fitted curve + plot(v_bias_vals, 10.^(BER_fit), '-r', 'LineWidth', 1.5,'Color',cols(i,:),'HandleVisibility','off'); hold on - a = scatter(rop_measured(i,:),ber_ffe(i,:),'Marker','+','DisplayName',['PAM ',num2str(wh.parameter.M.values(i))]); - a.DataTipTemplate.DataTipRows(1).Label = 'P_{out}'; + a = scatter(v_bias_vals,ber_ffe(i,:),'Marker','+','DisplayName',['PAM ',num2str(wh.parameter.M.values(i))],'MarkerEdgeColor',cols(i,:)); + a.DataTipTemplate.DataTipRows(1).Label = 'Vbias'; a.DataTipTemplate.DataTipRows(2).Label = 'BER'; a.DataTipTemplate.DataTipRows(2).Format ='%.1e'; - a.DataTipTemplate.DataTipRows(3).Label = 'Vbias'; - a.DataTipTemplate.DataTipRows(3).Value = v_bias_vals; + a.DataTipTemplate.DataTipRows(3).Label = 'P_{out}'; + a.DataTipTemplate.DataTipRows(3).Value = rop_measured(i,:); a.DataTipTemplate.DataTipRows(3).Format = 'auto'; end + % Continue with the rest of your plot settings yline(3.8e-3, 'DisplayName', 'HD-FEC', 'LineStyle', '--', 'HandleVisibility', 'off'); -xlabel('Measured MZM Output Power (dBm)'); +xlabel('Bias Voltage'); ylabel('Bit Error Rate (BER)'); -title('Bit Error Rate vs. ROP'); +title('Bit Error Rate vs. ROP | MI->DO | B2B'); set(gca, 'yscale', 'log'); set(gca, 'Box', 'on'); grid on; grid minor; legend('Interpreter', 'none'); +ylim([1e-4,0.5]); +xlim([1.6 3.2]); + diff --git a/projects/HighSpeedExperiment_2024/bias_optimization.m b/projects/HighSpeedExperiment_2024/bias_optimization.m index 431373d..88985c0 100644 --- a/projects/HighSpeedExperiment_2024/bias_optimization.m +++ b/projects/HighSpeedExperiment_2024/bias_optimization.m @@ -5,9 +5,9 @@ currentTime = datetime('now', 'Format', 'yyyyMMdd_HHmmss'); timeStr = char(currentTime); experiment_name = [experiment_name, timeStr]; -ffe_only = 1; +ffe_only = 0; postfilter_approach = 0; -db_channel_approach = 0; +db_channel_approach = 1; db_coding_approach = 0; db_precode = db_coding_approach || db_channel_approach; @@ -15,19 +15,20 @@ db_precode = db_coding_approach || db_channel_approach; %%% SIR Sweep for MPI Experiment %%% params = struct; -params.vbias = [1.7:0.02:3.2]; % PAM6=2.3V %PAM8=2.68V -params.awg_vpp = [2.7]; -params.precomp_amp_max = [5]; -params.rop_atten = [0]; -params.M = [4,6,8]; -params.lambda = [1292,1310,1327]; %calcWavelengthPlan(16, 400e9 , 1310); - -% params.vbias = [2.23]; % PAM6=2.3V %PAM8=2.68V +% params.vbias = [1.7:0.02:3.2]; % PAM6=2.3V %PAM8=2.68V % params.awg_vpp = [2.7]; % params.precomp_amp_max = [5]; % params.rop_atten = [0]; -% params.M = [6]; -% params.lambda = [1310]; %calcWavelengthPlan(16, 400e9 , 1310); +% params.M = [4,6,8]; +% params.lambda = [1293,1310,1327.4]; %calcWavelengthPlan(16, 400e9 , 1310); + +params.vbias = [2.3]; %PAM4=2.3 V PAM6=2.3V %PAM8=2.6V +params.awg_vpp = [2.7]; +params.precomp_amp_max = [-50]; +params.rop_atten = [0]; +params.M = [4]; +params.lambda = [1310]; %calcWavelengthPlan(16, 400e9 , 1310); +params.rcalpha = [0.05]; wh = DataStorage(params); @@ -39,14 +40,17 @@ wh.addStorage("pd_in"); wh.addStorage("rop"); wh.addStorage("m"); wh.addStorage("rx_logbook"); +wh.addStorage("dcs"); +wh.addStorage("pdfa"); +wh.addStorage("exfo"); precomp_path = "C:\Users\sioe\Documents\High_Speed_Measurement_2024\precomp\"; precomp_fn = "lab_high_speed"; precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active -precomp_amp_max = 4; +precomp_amp_max = -34; random_key = 2; -pd_in_set = 7; +pd_in_set = 8; looptotal = prod(wh.dim); @@ -62,301 +66,385 @@ loopcnt = 0; estimatedTimeRemaining = 0; estimatedTotalTime = 0; -for lambda = wh.parameter.lambda.values +for rcalpha = wh.parameter.rcalpha.values + 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; + 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) + if ~(exfo.cur_wavelength == lambda) - % 1) - pdfa.disablePDFA; + % 1) + pdfa.disablePDFA; - % 2) - exfo.setWavelength(lambda); + % 2) + exfo.setWavelength(lambda); - end + % 3) + pdfa.enablePDFA(); - % 3) - pdfa.enablePDFA(); + % 4) + pdfa.setPumpLevel(100); - % 4) - pdfa.setPumpLevel(100); - - % 5) SET to first vbias and wait 30 minutes - v_bias_first = wh.parameter.vbias.values(1); - dcs = DC_supply("active",[1,0],"voltage",[v_bias_first, 0]); - dcs.set("voltage",[v_bias_first, 0]); - % pause(30*60); %wait 30 minutes for stable bias - - for v_bias = wh.parameter.vbias.values - for rop_atten = wh.parameter.rop_atten.values - for precomp_amp_max = wh.parameter.precomp_amp_max.values - for M = wh.parameter.M.values - for awg_vpp = wh.parameter.awg_vpp.values - - iterationStartTime = tic; - loopcnt = loopcnt+1; - - if M == 4 - fsym = 168e9; - pulsef = 1; - elseif M == 6 - fsym = 144e9; - pulsef = 0; - elseif M == 8 - fsym = 120e9; - pulsef = 0; - end - - %%%%% Loop Preps - %fsym = round(targetrate/log2(M)); - loop_name = ['_fsym_',num2str(fsym)]; - - %%%%% SET Voltages %%%%%% - dcs = DC_supply("active",[1,0],"voltage",[v_bias, 0]); - dcs.set("voltage",[v_bias, 0]); - - %%%%% SET Attenuator %%%%%% - voa = OptAtten("active",[1,2,1,1],"value",[rop_atten,pd_in_set,0,0],"wavelength",[1310,1310,1310,1310]); - voa.set('active',[1,2,1,1],'value',[rop_atten,pd_in_set,0,0]); - % voa.readvals(); - - %%%%% Construct AWG and Scope Modules %%%%%% - fdac = 256e9; - fadc = 256e9; - SCP = ScopeKeysight("model","UXR1104B",'autoscale',1,"fadc","GSa_256","channel",[0,1,0,0],"recordLen",2000000,"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",0); % - - %%%%% Symbol Generation %%%%%% - Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",0.05); - - [Digi_sig,Symbols,Bits] = PAMsource(... - "fsym",fsym,"M",M,"order",19,"useprbs",1,... - "fs_out",fdac,... - "applyclipping",0,"clipfactor",1.5,... - "applypulseform",pulsef,"pulseformer",Pform,... - "randkey",random_key,... - "db_precode",db_precode,"db_encode",db_coding_approach,... - "mrds_code",0,"mrds_blocklength",512).process(); - - Digi_sig.spectrum("displayname","Normal Tx","fignum",20,"normalizeToNyquist",0,"normalizeTo0dB",0); - - %%%%% Precompensation Routine %%%%%% - if precomp_mode == 1 % measure channel - precomp_est = ChannelFreqResp("Nacq",1024,"Navg",64,"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); - 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); + end - %%%%% Plot and Save Routine 1 %%%%%%%%%%%%%%%%%%%%%%%%% - Digi_sig.spectrum("displayname","Normal Tx","fignum",10,"normalizeToNyquist",0,"normalizeTo0dB",0); + % 5) SET to first vbias and wait 30 minutes + v_bias_first = wh.parameter.vbias.values(1); + dcs = DC_supply("active",[1,0],"voltage",[v_bias_first, 0]); + dcs.set("voltage",[v_bias_first, 0]); + dcs.readVals(); - save([folderpath,[experiment_name,'_bits'],loop_name],"Bits"); - save([folderpath,[experiment_name,'_symbols'],loop_name],"Symbols"); - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + % pause(30*60); %wait 30 minutes for stable bias - %%%%% AWG --> Scope %%%%%% - [~,Scpe_sig_raw,~,D] = A2S.process("signal2",Digi_sig,"waitUntilClick",0); + for v_bias = wh.parameter.vbias.values + for rop_atten = wh.parameter.rop_atten.values + for precomp_amp_max = wh.parameter.precomp_amp_max.values + for M = wh.parameter.M.values + for awg_vpp = wh.parameter.awg_vpp.values - % Scpe_sig_raw.spectrum("displayname","Scope PSD","fignum",20); - % Scpe_sig_raw.plot("displayname","Scope raw signal","fignum",25,"clear",1); - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + iterationStartTime = tic; + loopcnt = loopcnt+1; - %%%%%% Sample to 2x fsym %%%%%% - Scpe_sig_resampled = Scpe_sig_raw.resample("fs_in",fadc,"fs_out",2*fsym); - - %%%%% Precompensation Routine %%%%%% - if precomp_mode == 1 - 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); - pd_in = voa.power_state(2); - disp(['ROP: ',num2str(rop),' dBm || PD in: ',num2str(pd_in), ' dBm']); - - %%%%%% 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); - - %%%%%% SNR CHEAT - Avges the measured signal occurences found after correlation in "tsynch" %%%%%% - average_signals = 0; - if average_signals - scope_mean = zeros(size(S{1}.signal)); - for n=1:numel(S) - scope_mean = scope_mean + S{n}.signal; + if M == 4 + fsym = 220e9; + pulsef = 0; + elseif M == 6 + fsym = 180e9; + pulsef = 0; + elseif M == 8 + fsym = 160e9; + pulsef = 0; end - scope_mean = scope_mean ./ n; - Scpe_sig_syncd.signal = scope_mean; - end - %%%%% Plot and Save Routines: SAVE RECEIVED SIGNALS %%%%%%%%%%%%%%%%%%%%%%%%% - save([folderpath,experiment_name,'_rx_signal',loop_name],"S"); + %%%%% Loop Preps + %fsym = round(targetrate/log2(M)); + loop_name = ['_fsym_',num2str(fsym)]; - Scpe_sig_syncd.eye(fsym,M,"fignum",40,"displayname",' after Scope'); + %%%%% SET Voltages %%%%%% + dcs = DC_supply("active",[1,0],"voltage",[v_bias, 0]); + dcs.set("voltage",[v_bias, 0]); - %%%%% EQUALIZE %%%%%% - Eq = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",50,"sps",2,"decide",0); - % Eq = VNLE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",[0.0004 0.0005 0.0006],"mu_tr",0,"order",[50,7,7],"sps",2,"decide",1); - Eq = 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); + %%%%% SET Attenuator %%%%%% + voa = OptAtten("active",[1,2,1,1],"value",[rop_atten,pd_in_set,0,0],"wavelength",[1310,1310,1310,1310],"speed",[1000,100,1000,1000]); + voa.set('active',[1,2,1,1],'value',[rop_atten,pd_in_set,0,0]); + % voa.readvals(); - % set to minus one not zero not avoid confusion if BER is acutally zero - ber_ffe = -1; - ber_mlse = -1; - ber_db = -1; + %%%%% Construct AWG and Scope Modules %%%%%% + fdac = 256e9; + fadc = 256e9; + SCP = ScopeKeysight("model","UXR1104B",'autoscale',1,"fadc","GSa_256","channel",[0,1,0,0],"recordLen",3000000,"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",0); % - if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% - ber = []; + %%%%% Symbol Generation %%%%%% + Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rcalpha); - parfor i = 1:numel(S) + [Digi_sig,Symbols,Bits] = PAMsource(... + "fsym",fsym,"M",M,"order",19,"useprbs",1,... + "fs_out",fdac,... + "applyclipping",0,"clipfactor",1.5,... + "applypulseform",pulsef,"pulseformer",Pform,... + "randkey",random_key,... + "db_precode",db_precode,"db_encode",db_coding_approach,... + "mrds_code",0,"mrds_blocklength",512).process(); - [EQ_sig] = Eq.process(S{i},Symbols); + Digi_sig.spectrum("displayname","Normal Tx","fignum",10,"normalizeToNyquist",0,"normalizeTo0dB",0); - % EQ_sig.plot("fignum",50,"displayname",'After EQ'); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 % measure channel + precomp_est = ChannelFreqResp("Nacq",1024,"Navg",64,"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); + 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); + + + %%%%% Plot and Save Routine 1 %%%%%%%%%%%%%%%%%%%%%%%%% + Digi_sig.spectrum("displayname","Normal Tx","fignum",14,"normalizeToNyquist",0,"normalizeTo0dB",0); + + % save([folderpath,[experiment_name,'_bits'],loop_name],"Bits"); + % save([folderpath,[experiment_name,'_symbols'],loop_name],"Symbols"); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %%%%% AWG --> Scope %%%%%% + [~,Scpe_sig_raw,~,D] = A2S.process("signal2",Digi_sig,"waitUntilClick",0); + + Scpe_sig_raw.spectrum("displayname","Scope PSD","fignum",20,"normalizeTo0dB",1); + Scpe_sig_raw.plot("displayname","Scope raw signal","fignum",29,"clear",1); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %%%%%% Sample to 2x fsym %%%%%% + Scpe_sig_resampled = Scpe_sig_raw.resample("fs_in",fadc,"fs_out",2*fsym); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 + 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); + pd_in = voa.power_state(2); + disp(['ROP: ',num2str(rop),' dBm || PD in: ',num2str(pd_in), ' dBm']); + + %%%%%% 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); + + %%%%%% SNR CHEAT - Avges the measured signal occurences found after correlation in "tsynch" %%%%%% + average_signals = 0; + if average_signals + Scpe_sig_avg = Scpe_sig_syncd; + scope_mean = zeros(size(S{1}.signal)); + for n=1:numel(S) + scope_mean = scope_mean + S{n}.signal; + end + scope_mean = scope_mean ./ n; + Scpe_sig_avg.signal = scope_mean; + + Scpe_sig_avg.spectrum("displayname","Scope PSD","fignum",20,"normalizeTo0dB",1); + Scpe_sig_avg.plot("displayname","Scope raw signal","fignum",27,"clear",1); + Scpe_sig_avg.eye(fsym,M,"fignum",41,"displayname",' Eye of AVG Signal'); + end + + % Optfilter = Filter('filtdegree',6,"f_cutoff",100e9,"fs",Scpe_sig_avg.fs,"filterType",filtertypes.gaussian,"active",true); + % Scpe_sig_syncd = Optfilter.process(Scpe_sig_syncd); + % Scpe_sig_syncd.spectrum("displayname","Scope PSD","fignum",20,"normalizeTo0dB",1); + + %%%%% Plot and Save Routines: SAVE RECEIVED SIGNALS %%%%%%%%%%%%%%%%%%%%%%%%% + % save([folderpath,experiment_name,'_rx_signal',loop_name],"S"); + + Scpe_sig_syncd.eye(fsym,M,"fignum",40,"displayname",' after Scope'); + + %%%%% EQUALIZE %%%%%% + Eq = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",50,"sps",2,"decide",0); + % Eq = VNLE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",[0.0004 0.0005 0.0006],"mu_tr",0,"order",[50,7,7],"sps",2,"decide",1); + Eq = 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); + + % set to minus one not zero not avoid confusion if BER is acutally zero + ber_ffe = -1; + ber_mlse = -1; + ber_db = -1; + + if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% + ber = []; + + parfor i = 1:numel(S) + + [EQ_sig] = Eq.process(S{i},Symbols); + + % EQ_sig.plot("fignum",50,"displayname",'After EQ'); + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + + [~,errors_bm,ber(i),errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + % disp(['FFE: ',sprintf('%.1E',ber(i)),'| ROP: ',num2str(rop),' dB | PD_in: ',num2str(pd_in),' dBm']); + end + + disp(['FFE EQ: BEST BER: ',sprintf('%.1E',min(ber)),' AVG BER: ',sprintf('%.1E',mean(ber)),' WORST:',sprintf('%.1E',max(ber)),'. Out of ',num2str(numel(ber))]); + + try + ber_ffe = mean(rmoutliers(ber)); + catch + ber_ffe = min(ber); + end + + if 0 + + EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + + + figure(56); + clf + title(sprintf('PAM %d ; BER: %1.2e',M, ber_ffe)); + 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.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 + + elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + ber_vnle = []; + ber_vnle_mlse = []; + ber_ffe_mlse =[]; + ber_ffe = []; + + parfor s = 1:numel(S) + + if 1 + %FFE LINEAR + Eq = 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); + Scpe_sig_syncd = S{s}; + [EQ_ffe] = Eq.process(Scpe_sig_syncd,Symbols); + Noi = EQ_ffe-Symbols; + Rx_bits = PAMmapper(M,0).demap(EQ_ffe); + [~,num_errors,ber_ffe(s),pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + end + + if 1 + %FFE + MLSE + nc = 2; + burg_coeff = arburg(Noi.signal,nc); + EQ_ffe = EQ_ffe.filter(burg_coeff,1); + + EQ_mlse = MLSE("DIR",burg_coeff,"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_ffe); + Rx_bits = PAMmapper(M,0).demap(EQ_mlse); + [~,num_errors,ber_ffe_mlse(s),pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + end + + %VNLE + Eq = 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); + Scpe_sig_syncd = S{s}; + [EQ_vnle] = Eq.process(Scpe_sig_syncd,Symbols); + Noi = EQ_vnle-Symbols; + Rx_bits = PAMmapper(M,0).demap(EQ_vnle); + [~,num_errors,ber_vnle(s),pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + + %VNLE + MLSE + if 0 + nc = 2; + burg_coeff = arburg(Noi.signal,nc); + EQ_mlse = EQ_vnle.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); + [~,num_errors,ber_vnle_mlse(s),pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + end + 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))]); + + + if 0 + window = 100; + Noi_ = Noi; + Noi.signal = Noi.signal - movmean(Noi.signal,[floor(window/2),ceil(window/2)]); + + EQ_vnle.spectrum('displayname','EQ out PSD','fignum',123); + Noi.spectrum('displayname','Noise PSD','fignum',123); + + nc = 1; + burg_coeff = arburg(Noi.signal,nc); + [h,w] = freqz(1,burg_coeff,length(Noi),"whole",Noi.fs); + h = h/max(abs(h)); + hold on + w_ = (w - Noi.fs/2); + plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']); + end + if 0 + figure(57); + clf + title(sprintf('PAM %d ; BER: %1.2e',M, ber_ffe)); + 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.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: ',sprintf('%.1E',ber_ffe),' -> PF -> MLSE: ',sprintf('%.1E',ber_mlse),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + parfor s = 1:numel(S) + Scpe_sig_syncd = S{s}; + [EQ_sig, Noi] = Eq.process(Scpe_sig_syncd,Duobinary().encode(Symbols)); + + % EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + + 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); + + %disp([' DB Precode -> Channel -> FFE -> Decode/ Mod ',sprintf('%.1E',ber_db(s)),' | PD_in: ',num2str(pd_in),' dBm']); + 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))]); + ber = min(ber_db); + + elseif db_coding_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig_syncd,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); + [~,errors_bm,ber_db,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - [~,errors_bm,ber(i),errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - disp(['FFE: ',sprintf('%.1E',ber(i)),'| ROP: ',num2str(rop),' dB | PD_in: ',num2str(pd_in),' dBm']); + EQ_sig.plot("fignum",50,"displayname",'After EQ'); + disp([' DB Precode -> DB Code -> Channel -> FFE -> Decode/ Mod ',sprintf('%.1E',ber_db),' | PD_in: ',num2str(pd_in),' dBm']); end - - try - ber_ffe = mean(rmoutliers(ber)); - catch - ber_ffe = min(ber); - end - - if 0 - - EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); - figure(56); - clf - title(sprintf('PAM %d ; BER: %1.2e',M, ber_ffe)); - 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.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 + %%%%% Store measurement into measurement "warehouse" %%%%%% + wh.addValueToStorage(ber,'ber_collect',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); + wh.addValueToStorage(ber_ffe,'ber_ffe',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); + wh.addValueToStorage(ber_mlse,'ber_mlse',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); + wh.addValueToStorage(ber_db,'ber_db',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); - elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + wh.addValueToStorage(rop,'rop',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); + wh.addValueToStorage(pd_in,'pd_in',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); + % Rx_bits.logbook.SignalCopy = []; + % wh.addValueToStorage(Rx_bits,'rx_logbook',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); + wh.addValueToStorage(M,'m',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); - [EQ_sig] = Eq.process(Scpe_sig_syncd,Symbols); + wh.addValueToStorage(dcs,'dcs',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); + wh.addValueToStorage(pdfa,'pdfa',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); + wh.addValueToStorage(exfo,'exfo',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda,rcalpha); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + %%%%% Plot stuff into Table (feel free to add own values in -> 'name',value <- notation. Must be closed when table is changed)%%%%%%%%%%%%%%%%%%%%%% + % showCurrentMeasurement('BER', min(ber_vnle),'BER',mean(ber_vnle),'Alpha',rcalpha, 'Fsym',fsym.*1e-9, 'ROP', rop,'pulsef',pulsef,'rrcalpha',rrcalpha, 'PAM',M, 'Vbias', v_bias, 'AWG Vpp', awg_vpp, 'Precomp MaxAmp',precomp_amp_max); - Noi = EQ_sig-Symbols; + %%%%% Arrange Figures %%%%%%%%%%%%%%%%%%%%%% + autoArrangeFigures(3,3,2); - Rx_bits = PAMmapper(M,0).demap(EQ_sig); + iterationTimes(loopcnt) = toc(iterationStartTime); + averageTimePerIteration = mean(iterationTimes(1:loopcnt)); + estimatedTotalTime = averageTimePerIteration * looptotal; + estimatedTimeRemaining = estimatedTotalTime - sum(iterationTimes(1:loopcnt)); + progressFraction = loopcnt / looptotal; + waitbar(progressFraction, hWaitbar, ... + sprintf('Loop: %d of %d \n Runtime: %.1f min | %.1f sec per Loop |Time to go: %.1f min ', ... + loopcnt, looptotal, sum(iterationTimes(1:loopcnt))/60, averageTimePerIteration, estimatedTimeRemaining/60 )); - [~,num_errors,ber_ffe,pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - - nc = 2; - burg_coeff = arburg(Noi.signal,nc); - - EQ_sig = EQ_sig.filter(burg_coeff,1); - - if 0 - Noi.spectrum('displayname','Noise PSD','fignum',123) - [h,w] = freqz(1,burg_coeff,length(Noi),"whole",Noi.fs); - h = h/max(abs(h)); - hold on - w_ = (w - Noi.fs/2); - plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']); - end - - EQ_sig = MLSE("DIR",burg_coeff,"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig); - - Rx_bits = PAMmapper(M,0).demap(EQ_sig); - - [~,num_errors,ber_mlse,pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - disp(['FFE: ',sprintf('%.1E',ber_ffe),' -> PF -> MLSE: ',sprintf('%.1E',ber_mlse),' dB | PD_in: ',num2str(pd_in),' dBm']); - - elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% - - [EQ_sig, Noi] = Eq.process(Scpe_sig_syncd,Duobinary().encode(Symbols)); - - EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); - - 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,pos_errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - disp([' DB Precode -> Channel -> FFE -> Decode/ Mod ',sprintf('%.1E',ber_db),' | PD_in: ',num2str(pd_in),' dBm']); - - elseif db_coding_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% - - [EQ_sig, Noi] = Eq.process(Scpe_sig_syncd,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); - [~,errors_bm,ber_db,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - EQ_sig.plot("fignum",50,"displayname",'After EQ'); - - disp([' DB Precode -> DB Code -> Channel -> FFE -> Decode/ Mod ',sprintf('%.1E',ber_db),' | PD_in: ',num2str(pd_in),' dBm']); + wh.save([folderpath,experiment_name,'_wh']); end - - - %%%%% Store measurement into measurement "warehouse" %%%%%% - wh.addValueToStorage(ber,'ber_collect',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda); - wh.addValueToStorage(ber_ffe,'ber_ffe',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda); - wh.addValueToStorage(ber_mlse,'ber_mlse',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda); - wh.addValueToStorage(ber_db,'ber_db',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda); - - wh.addValueToStorage(rop,'rop',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda); - wh.addValueToStorage(pd_in,'pd_in',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda); - Rx_bits.logbook.SignalCopy = []; - wh.addValueToStorage(Rx_bits,'rx_logbook',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda); - wh.addValueToStorage(M,'m',v_bias,awg_vpp,precomp_amp_max,rop_atten,M,lambda); - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - - %%%%% Plot stuff into Table (feel free to add own values in -> 'name',value <- notation. Must be closed when table is changed)%%%%%%%%%%%%%%%%%%%%%% - showCurrentMeasurement('FFE', ber_ffe,'MLSE',ber_mlse, 'Fsym',fsym.*1e-9, 'ROP', rop, 'PD in', pd_in, 'PAM',M, 'Vbias', v_bias, 'AWG Vpp', awg_vpp, 'Precomp MaxAmp',precomp_amp_max); - - %%%%% Arrange Figures %%%%%%%%%%%%%%%%%%%%%% - autoArrangeFigures(3,3,2); - - iterationTimes(loopcnt) = toc(iterationStartTime); - averageTimePerIteration = mean(iterationTimes(1:loopcnt)); - estimatedTotalTime = averageTimePerIteration * looptotal; - estimatedTimeRemaining = estimatedTotalTime - sum(iterationTimes(1:loopcnt)); - progressFraction = loopcnt / looptotal; - waitbar(progressFraction, hWaitbar, ... - sprintf('Loop: %d of %d \n Runtime: %.1f min | %.1f sec per Loop |Time to go: %.1f min ', ... - loopcnt, looptotal, sum(iterationTimes(1:loopcnt))/60, averageTimePerIteration, estimatedTimeRemaining/60 )); - - wh.save([folderpath,experiment_name,'_wh']); - end end end @@ -364,11 +452,12 @@ for lambda = wh.parameter.lambda.values end end + close(hWaitbar); wh.save([folderpath,experiment_name,'_wh']); -autoArrangeFigures(3,3,2) +% autoArrangeFigures(3,3,2) %%% LAMBDA PLOT diff --git a/projects/HighSpeedExperiment_2024/master_evaluation.m b/projects/HighSpeedExperiment_2024/master_evaluation.m new file mode 100644 index 0000000..cf09489 --- /dev/null +++ b/projects/HighSpeedExperiment_2024/master_evaluation.m @@ -0,0 +1,433 @@ + +folderpath = 'C:\Users\sioe\Documents\High_Speed_Measurement_2024\8km_bitrate_rop_master\'; +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.rop_atten = [0]; +end + +if 0 + %%% ROP SWEEP + awg_vpp = 2.7; + pd_in_set = 6; + random_key = 0; + params = struct; + params.M = [8,6,4]; + params.lambda = [1310]; %calcWavelengthPlan(16, 400e9 , 1310); + params.bitrate = [300:30:480].*1e9; + params.duobinary = [1,0]; + params.rop_atten = [0:1.5:7.5]; +end + +wh = DataStorage(params); + +wh.addStorage("ber_ffe"); +wh.addStorage("ber_ffe_mlse"); +wh.addStorage("ber_vnle"); +wh.addStorage("ber_vnle_mlse"); +wh.addStorage("ber_db"); + +wh.addStorage("FFE"); +wh.addStorage("VNLE"); + +wh.addStorage("pd_in"); +wh.addStorage("rop"); +wh.addStorage("m"); + +wh.addStorage("dcs"); +wh.addStorage("pdfa"); +wh.addStorage("exfo"); +wh.addStorage("voa"); + +precomp_path = "C:\Users\sioe\Documents\High_Speed_Measurement_2024\precomp\"; +precomp_fn = "lab_high_speed"; +precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active + +looptotal = prod(wh.dim); + +disp(['Start Measurement of ',num2str(looptotal),' loops...']) +iterationTimes = zeros(looptotal, 1); % Preallocate for speed +if ~exist('hWaitbar', 'var') || ~isvalid(hWaitbar) + hWaitbar = waitbar(0, sprintf('Starting %d measurements',looptotal), 'Name', 'Processing Progress'); +else + waitbar(0, hWaitbar, sprintf('Starting %d measurements',looptotal)); +end + +loopcnt = 0; +estimatedTimeRemaining = 0; +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]); + 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); + + end + + for bitrate = wh.parameter.bitrate.values + + fsym = floor( bitrate*1e-9./log2(M) ).*1e9; + + for db = wh.parameter.duobinary.values + + 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 + + + + %%%%% Construct AWG and Scope Modules %%%%%% + fdac = 256e9; + fadc = 256e9; + SCP = ScopeKeysight("model","UXR1104B",'autoscale',1,"fadc","GSa_256","channel",[0,1,0,0],"recordLen",4000000,"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",0); % + + %%%%% Symbol Generation %%%%%% + rcalpha = 0.05; + Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rcalpha); + + Pamsource = PAMsource(... + "fsym",fsym,"M",M,"order",19,"useprbs",1,... + "fs_out",fdac,... + "applyclipping",0,"clipfactor",1.5,... + "applypulseform",pulsef,"pulseformer",Pform,... + "randkey",random_key,... + "db_precode",db_precode,"db_encode",db_coding_approach,... + "mrds_code",0,"mrds_blocklength",512); + + [Digi_sig,Symbols,Bits] = Pamsource.process(); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 % measure channel + precomp_est = ChannelFreqResp("Nacq",1024,"Navg",64,"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); + 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); + + 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,'.','_'); + + %%%%% READ Voltages %%%%%% + dcs.readVals(); + + %%%%% SET Attenuator %%%%%% + voa = OptAtten("active",[1,2,1,1],"value",[rop_atten,pd_in_set,0,0],"wavelength",[1310,1310,1310,1310],"speed",[1000,100,1000,1000]); + voa.set('active',[1,2,1,1],'value',[rop_atten,pd_in_set,0,0]); + voa.readvals(); + + %%%% SIGNAL USUALLY HERE, NOW ABOVE ROP_ATTEN %%% + + %%%%% Plot and Save Routine 1 - same for all rops, thus save only once %%%%%%%%%%%%%%%%%%%%%%%%% + if rop_atten == 0 + save([folderpath,experiment_name,loop_name,'_bits'],"Bits"); + save([folderpath,experiment_name,loop_name,'_symbols'],"Symbols"); + end + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %%%%% AWG --> Scope %%%%%% + [~,Scpe_sig_raw,~,D] = A2S.process("signal2",Digi_sig,"waitUntilClick",0); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + + % 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.spectrum("displayname","Scope PSD","fignum",30,"normalizeTo0dB",1); + % Scpe_sig_raw.eye(fsym,M,"displayname",'eye','fignum',200); + + %%%%%% Sample to 2x fsym %%%%%% + Scpe_sig_resampled = Scpe_sig_raw.resample("fs_in",fadc,"fs_out",2*fsym); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 + 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); + pd_in = voa.power_state(2); + + %%%%%% 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 + ber_vnle = [-1]; + ber_vnle_mlse = [-1]; + ber_ffe_mlse =[-1]; + ber_ffe = [-1]; + ber_db = [-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,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 + 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 + 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 + 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 + + + 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))); + hold on + histogram(received(lvl,:),1000,"EdgeAlpha",0,'DisplayName',['Lvl ',num2str(lvl),' | ',num2str(cnt(lvl)),' entries']); + end + legend + 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); + + if 0 + parfor s = 1:numel(S) + + Scpe_sig_syncd = S{s}; + + [EQ_sig, Noi] = 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 = 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 + + 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') + + 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 + + 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 + + + %%%%% 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); + wh.addValueToStorage({ber_vnle_mlse},'ber_vnle_mlse',M,lambda,bitrate,db,rop_atten); + wh.addValueToStorage({ber_db},'ber_db',M,lambda,bitrate,db,rop_atten); + + wh.addValueToStorage(rop,'rop',M,lambda,bitrate,db,rop_atten); + wh.addValueToStorage(pd_in,'pd_in',M,lambda,bitrate,db,rop_atten); + wh.addValueToStorage(M,'m',M,lambda,bitrate,db,rop_atten); + + wh.addValueToStorage(dcs,'dcs',M,lambda,bitrate,db,rop_atten); + wh.addValueToStorage(pdfa,'pdfa',M,lambda,bitrate,db,rop_atten); + wh.addValueToStorage(exfo,'exfo',M,lambda,bitrate,db,rop_atten); + wh.addValueToStorage(voa,'voa',M,lambda,bitrate,db,rop_atten); + + wh.addValueToStorage(ffe,'FFE',M,lambda,bitrate,db,rop_atten); + wh.addValueToStorage(vnle,'VNLE',M,lambda,bitrate,db,rop_atten); + + + iterationTimes(loopcnt) = toc(iterationStartTime); + averageTimePerIteration = mean(iterationTimes(1:loopcnt)); + estimatedTotalTime = averageTimePerIteration * looptotal; + estimatedTimeRemaining = estimatedTotalTime - sum(iterationTimes(1:loopcnt)); + progressFraction = loopcnt / looptotal; + waitbar(progressFraction, hWaitbar, ... + sprintf('Loop: %d of %d \n Runtime: %.1f min | %.1f sec per Loop |Time to go: %.1f min ', ... + loopcnt, looptotal, sum(iterationTimes(1:loopcnt))/60, averageTimePerIteration, estimatedTimeRemaining/60 )); + + wh.save([folderpath,experiment_name,'_wh']); + + end + end + end + end +end + +close(hWaitbar); + +wh.save([folderpath,experiment_name,'_wh']); + + diff --git a/projects/HighSpeedExperiment_2024/mpi_measurement.m b/projects/HighSpeedExperiment_2024/mpi_measurement.m new file mode 100644 index 0000000..37e1b04 --- /dev/null +++ b/projects/HighSpeedExperiment_2024/mpi_measurement.m @@ -0,0 +1,418 @@ + +folderpath = 'C:\Users\sioe\Documents\High_Speed_Measurement_2024\mpi_measurement\'; +experiment_name = 'testen'; +currentTime = datetime('now', 'Format', 'yyyyMMdd_HHmmss'); +timeStr = char(currentTime); +experiment_name = [experiment_name, timeStr]; + + +%%% BITRATE Sweep for MPI Experiment %%% +awg_vpp = 2.7; +rop_atten = 0; %VOA 1 -> %nicht angeschlossen +pd_in_set = 8; %VOA 2 -> PD in +%VOA 3 -> Signal path +%VOA 4 -> Interference path + +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]; + +wh = DataStorage(params); + +wh.addStorage("ber_ffe"); +wh.addStorage("ber_ffe_mlse"); +wh.addStorage("ber_vnle"); +wh.addStorage("ber_vnle_mlse"); +wh.addStorage("ber_db"); + +wh.addStorage("FFE"); +wh.addStorage("VNLE"); + +wh.addStorage("pd_in"); +wh.addStorage("rop"); +wh.addStorage("s_power"); +wh.addStorage("i_power"); +wh.addStorage("sir"); + + +wh.addStorage("m"); + +wh.addStorage("dcs"); +wh.addStorage("pdfa"); +wh.addStorage("exfo"); +wh.addStorage("voa"); + +precomp_path = "C:\Users\sioe\Documents\High_Speed_Measurement_2024\precomp\"; +precomp_fn = "lab_high_speed"; +precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active + +looptotal = prod(wh.dim); + +disp(['Start Measurement of ',num2str(looptotal),' loops...']) +iterationTimes = zeros(looptotal, 1); % Preallocate for speed +if ~exist('hWaitbar', 'var') || ~isvalid(hWaitbar) + hWaitbar = waitbar(0, sprintf('Starting %d measurements',looptotal), 'Name', 'Processing Progress'); +else + waitbar(0, hWaitbar, sprintf('Starting %d measurements',looptotal)); +end + +loopcnt = 0; +estimatedTimeRemaining = 0; +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 + + 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 + + %%%%% Construct AWG and Scope Modules %%%%%% + fdac = 256e9; + fadc = 256e9; + + + %%%%% Symbol Generation %%%%%% + rcalpha = 0.05; + Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rcalpha); + + Pamsource = PAMsource(... + "fsym",fsym,"M",M,"order",19,"useprbs",1,... + "fs_out",fdac,... + "applyclipping",0,"clipfactor",1.5,... + "applypulseform",pulsef,"pulseformer",Pform,... + "randkey",random_key,... + "db_precode",db_precode,"db_encode",db_coding_approach,... + "mrds_code",0,"mrds_blocklength",512); + + [Digi_sig,Symbols,Bits] = Pamsource.process(); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 % measure channel + precomp_est = ChannelFreqResp("Nacq",1024,"Navg",64,"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); + 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 = 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; + + 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); + 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 + + %%%%% Loop Preps + iterationStartTime = tic; + loopcnt = loopcnt+1; + loop_name = ['_PAM_',num2str(M),'_R_',num2str(bitrate),'_DB_',num2str(db),'_I_atten_',num2str(interference_atten)]; + loop_name = strrep(loop_name,'.','_'); + + %%%%% READ Voltages %%%%%% + dcs.readVals(); + + %%%%% SET Attenuator %%%%%% + voa = OptAtten("active",[1,2,1,1],"value",[rop_atten,pd_in_set,0,interference_atten],"wavelength",[1310,1310,1310,1310],"speed",[1000,100,1000,1000]); + voa.set('active',[1,2,1,1],'value',[rop_atten,pd_in_set,0,interference_atten]); + + %%%% SIGNAL USUALLY HERE, NOW ABOVE ROP_ATTEN %%% + + %%%%% Plot and Save Routine 1 - same for all rops, thus save only once %%%%%%%%%%%%%%%%%%%%%%%%% + if interference_atten == 0 + save([folderpath,experiment_name,loop_name,'_bits'],"Bits"); + save([folderpath,experiment_name,loop_name,'_symbols'],"Symbols"); + end + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %%%%% 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); + + 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.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']); + + %%%%%% 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 + ber_vnle = [-1]; + ber_vnle_mlse = [-1]; + ber_ffe_mlse =[-1]; + ber_ffe = [-1]; + ber_db = [-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,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 %%%%%%%%%%%%%%%%%%%%%%%%%%% + + + 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 + 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 + 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 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); + 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; + end + + 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))]); + + 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); + + if 0 + + 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{s}, Noi{s}] = ffe.process(Scpe_sig_syncd,Duobinary().encode(Symbols)); + + EQ_sig{s}.signal = EQ_sig{s}.signal-mean(EQ_sig{s}.signal); + + EQ_sig_mlse = MLSE("DIR",[1,1],"duobinary_output",1,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig{s}); + + 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 + Noi{1}.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); + 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') + + end + + end + + + %%%%% Store measurement into measurement "warehouse" %%%%%% + wh.addValueToStorage({ber_ffe},'ber_ffe',M,bitrate,db,interference_atten); + wh.addValueToStorage({ber_ffe_mlse},'ber_ffe_mlse',M,bitrate,db,interference_atten); + wh.addValueToStorage({ber_vnle},'ber_vnle',M,bitrate,db,interference_atten); + wh.addValueToStorage({ber_vnle_mlse},'ber_vnle_mlse',M,bitrate,db,interference_atten); + wh.addValueToStorage({ber_db},'ber_db',M,bitrate,db,interference_atten); + + wh.addValueToStorage(rop,'rop',M,bitrate,db,interference_atten); + wh.addValueToStorage(pd_in,'pd_in',M,bitrate,db,interference_atten); + wh.addValueToStorage(s_power,'s_power',M,bitrate,db,interference_atten); + wh.addValueToStorage(i_power,'i_power',M,bitrate,db,interference_atten); + wh.addValueToStorage(sir,'sir',M,bitrate,db,interference_atten); + + wh.addValueToStorage(M,'m',M,bitrate,db,interference_atten); + + wh.addValueToStorage(dcs,'dcs',M,bitrate,db,interference_atten); + + exfo = Exfo_laser("serialport_number",'COM8','mainframe_channel',1,'safety_mode',0); + exfo.getLaserInfo; + + wh.addValueToStorage(exfo,'exfo',M,bitrate,db,interference_atten); + wh.addValueToStorage(voa,'voa',M,bitrate,db,interference_atten); + + wh.addValueToStorage(ffe,'FFE',M,bitrate,db,interference_atten); + wh.addValueToStorage(vnle,'VNLE',M,bitrate,db,interference_atten); + + iterationTimes(loopcnt) = toc(iterationStartTime); + averageTimePerIteration = mean(iterationTimes(1:loopcnt)); + estimatedTotalTime = averageTimePerIteration * looptotal; + estimatedTimeRemaining = estimatedTotalTime - sum(iterationTimes(1:loopcnt)); + progressFraction = loopcnt / looptotal; + waitbar(progressFraction, hWaitbar, ... + sprintf('Loop: %d of %d \n Runtime: %.1f min | %.1f sec per Loop |Time to go: %.1f min ', ... + 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); + +wh.save([folderpath,experiment_name,'_wh']); + +disp('Measurement complete')