From bf94e3dc2f6bf9f6276a2a35b5d96091652eaffc Mon Sep 17 00:00:00 2001 From: Silas Labor Zizou Date: Tue, 15 Oct 2024 08:43:27 +0200 Subject: [PATCH] Work on lab PC - many new automations - scripts to record and save MPI and bias optimizations... --- Classes/04_DSP/Coding/Duobinary.m | 2 +- Classes/05_Lab/Exfo_laser.m | 74 ++++ Classes/Warehouse_class/classes/DataStorage.m | 16 +- Functions/showCurrentMeasurement.m | 136 +++++++ Functions/updateWaitbar.m | 45 +++ projects/Lab_2024/bias_sweep_evaluation.m | 145 +++++++ projects/Lab_2024/bias_sweep_evaluation_2.m | 124 ++++++ projects/Lab_2024/lab_baudrate_sweep.m | 274 +++++++++++++ projects/Lab_2024/lab_bias_sweep.m | 269 +++++++++++++ projects/Lab_2024/lab_bias_sweep_gigantisch.m | 320 +++++++++++++++ projects/Lab_2024/lab_db_precode.m | 197 +++++----- projects/Lab_2024/lab_modulator_tf_sweep.m | 68 ++++ projects/Lab_2024/lab_precompensation_sweep.m | 265 +++++++++++++ projects/Lab_2024/lab_rop_sweep.m | 261 +++++++++++++ projects/Lab_2024/lab_sir_sweep.m | 368 ++++++++++++++++++ projects/Lab_2024/sweep_laser_vs_power.m | 88 +++++ 16 files changed, 2548 insertions(+), 104 deletions(-) create mode 100644 Classes/05_Lab/Exfo_laser.m create mode 100644 Functions/showCurrentMeasurement.m create mode 100644 Functions/updateWaitbar.m create mode 100644 projects/Lab_2024/bias_sweep_evaluation.m create mode 100644 projects/Lab_2024/bias_sweep_evaluation_2.m create mode 100644 projects/Lab_2024/lab_baudrate_sweep.m create mode 100644 projects/Lab_2024/lab_bias_sweep.m create mode 100644 projects/Lab_2024/lab_bias_sweep_gigantisch.m create mode 100644 projects/Lab_2024/lab_modulator_tf_sweep.m create mode 100644 projects/Lab_2024/lab_precompensation_sweep.m create mode 100644 projects/Lab_2024/lab_rop_sweep.m create mode 100644 projects/Lab_2024/lab_sir_sweep.m create mode 100644 projects/Lab_2024/sweep_laser_vs_power.m diff --git a/Classes/04_DSP/Coding/Duobinary.m b/Classes/04_DSP/Coding/Duobinary.m index e9888b9..259e04d 100644 --- a/Classes/04_DSP/Coding/Duobinary.m +++ b/Classes/04_DSP/Coding/Duobinary.m @@ -108,7 +108,7 @@ classdef Duobinary data = data - b; data = data ./ 2; - assert(isequal((0:M-1)',unique(data)),'Check Duobinary Precoding'); %seems the signal is not unipolar + % assert(isequal((0:M-1)',unique(data)),'Check Duobinary Precoding'); %seems the signal is not unipolar % duobinary coding (1+D) % coeff = [1,1]; diff --git a/Classes/05_Lab/Exfo_laser.m b/Classes/05_Lab/Exfo_laser.m new file mode 100644 index 0000000..1fa6acf --- /dev/null +++ b/Classes/05_Lab/Exfo_laser.m @@ -0,0 +1,74 @@ +classdef Exfo_laser + + properties(Access=public) + wavelength + power + end + + methods (Access=public) + + function obj = Exfo_laser(options) + + + arguments + options.wavelength = 1310; %dbm + options.power = -10; %dbm + end + + % + fn = fieldnames(options); + for n = 1:numel(fn) + try + obj.(fn{n}) = options.(fn{n}); + end + end + + end + + function success = set(obj,options) + + + arguments + obj + options.wavelength = obj.wavelength; %dbm + options.power = obj.power; %dbm + end + + % Connect to the laser + o = serialport("COM8", 9600); + configureTerminator(o, "CR"); % Set the terminator to carriage return (CR) + writeline(o, "*IDN?"); + pause(1); + if o.NumBytesAvailable ~= 0 + disp(['Laser Mainframe: ', readline(o)]); + else + error('No connection to the mainframe'); + clear o; + end + + + + end + + % Function to set the wavelength of the laser + function setLaserWavelength(~,serialObj, channel, wavelength) + command = ['CH', num2str(channel), ':L=', num2str(wavelength)]; + writeline(serialObj, command); + pause(0.5); % Allow time for the wavelength to change + writeline(serialObj, ['CH', num2str(channel), ':L?']); % Query current wavelength + current_wavelen = readline(serialObj); + disp(['Current Wavelength: ', current_wavelen]); + end + + % Function to set the laser power + function setLaserPower(~,serialObj, channel, power_dBm) + command = ['CH', num2str(channel), ':P=', num2str(power_dBm)]; + writeline(serialObj, command); + pause(0.2); % Allow time for power to adjust + writeline(serialObj, ['CH', num2str(channel), ':P?']); % Query current power + current_power = readline(serialObj); + disp(['Current Power: ', current_power, ' dBm']); + end + + end +end diff --git a/Classes/Warehouse_class/classes/DataStorage.m b/Classes/Warehouse_class/classes/DataStorage.m index 6fd1e97..dd9eb4b 100644 --- a/Classes/Warehouse_class/classes/DataStorage.m +++ b/Classes/Warehouse_class/classes/DataStorage.m @@ -128,7 +128,21 @@ classdef DataStorage < handle try tmp = obj.sto.(storageVarName){lin_idx(i)}; if ~isempty(tmp) - value(i,:) = tmp ; + if isa(tmp,'Signal') + 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} ; + end + else + value(i,:) = tmp ; + end else errcnt = errcnt+1; diff --git a/Functions/showCurrentMeasurement.m b/Functions/showCurrentMeasurement.m new file mode 100644 index 0000000..b9710b8 --- /dev/null +++ b/Functions/showCurrentMeasurement.m @@ -0,0 +1,136 @@ +function showCurrentMeasurement(varargin) + % showCurrentMeasurement displays measurement data in a figure with variable names + % as column headers and values listed below. Calling the function multiple times + % with the same variable names but different values adds more data points. + % + % Usage: + % showCurrentMeasurement('VariableName1', VariableValue1, 'VariableName2', VariableValue2, ...) + % + % Example: + % % First measurement + % voltage = 5.12; + % ber = 3.86e-5; + % power = voltage * 0.85; + % showCurrentMeasurement('Voltage', voltage, 'ber', ber, 'Power', power); + % + % % Second measurement + % voltage = 5.15; + % ber = 2.54e-5; + % power = voltage * 0.90; + % showCurrentMeasurement('Voltage', voltage, 'ber', ber, 'Power', power); + + % Validate that inputs are in name-value pairs + if mod(nargin, 2) ~= 0 + error('Inputs must be provided as name-value pairs.'); + end + + % Extract variable names and values + numPairs = nargin / 2; + names = varargin(1:2:end); + values = varargin(2:2:end); + + % Ensure variable names are strings + for i = 1:length(names) + if ~ischar(names{i}) && ~isstring(names{i}) + error('Variable names must be strings.'); + end + names{i} = char(names{i}); + end + + % Convert values to strings for display + for i = 1:length(values) + varNameLower = lower(names{i}); % Convert variable name to lowercase for case-insensitive comparison + if isnumeric(values{i}) + if strcmp(varNameLower, 'ber') + % Format 'ber' values in exponential notation with two decimal places + values{i} = sprintf('%.2e', values{i}); + else + values{i} = num2str(values{i}); + end + else + values{i} = char(values{i}); + end + end + + % Define a unique tag for the figure to locate it later + figTag = 'CurrentMeasurementsFigure'; + + % Try to find an existing figure with the specified tag + hFig = findobj('Type', 'figure', 'Tag', figTag); + + if isempty(hFig) + % Create a new figure and table + hFig = figure('Name', 'Current Measurements', 'NumberTitle', 'off', ... + 'MenuBar', 'none', 'ToolBar', 'none', 'Resize', 'on', ... + 'Tag', figTag); + + % Initialize data and column names + data = values; + columnNames = names; + + % Create the uitable + hTable = uitable('Parent', hFig, 'Data', data, ... + 'ColumnName', columnNames, ... + 'FontSize', 14, ... + 'RowName', [], ... + 'Units', 'normalized', ... + 'Position', [0, 0, 1, 1]); + % Adjust column widths + setColumnWidths(hTable); + + % Store the table handle for future use + setappdata(hFig, 'DataTable', hTable); + else + % Retrieve the existing table handle + hTable = getappdata(hFig, 'DataTable'); + + % Get current data and column names + currentData = get(hTable, 'Data'); + columnNames = get(hTable, 'ColumnName'); + + % Ensure that variable names are consistent + if ~isequal(columnNames, names') + error('Variable names must be consistent with previous calls.'); + end + + % Append new data to the existing data + updatedData = [currentData; values]; + + % Update the table data + set(hTable, 'Data', updatedData); + + % Adjust column widths + setColumnWidths(hTable); + end + + % Adjust the figure size to fit the table content without changing its position + drawnow; + tableExtent = get(hTable, 'Extent'); + % Get the current figure position + figPosition = get(hFig, 'Position'); + % Update the figure size while preserving the position + figPosition(3) = max(figPosition(3), tableExtent(3) + 20); % Width + figPosition(4) = max(figPosition(4), tableExtent(4) + 20); % Height + set(hFig, 'Position', figPosition); +end + +function setColumnWidths(hTable) + % Helper function to adjust column widths based on content + data = get(hTable, 'Data'); + columnNames = get(hTable, 'ColumnName'); + numColumns = length(columnNames); + columnWidths = cell(1, numColumns); + + % Calculate the maximum width needed for each column + for col = 1:numColumns + maxContentLength = max(cellfun(@length, data(:, col))); + headerLength = length(columnNames{col}); + maxLength = max(maxContentLength, headerLength); + + % Estimate pixel width (approximate, adjust as needed) + pixelWidth = maxLength * 14; % 8 pixels per character as an estimate + columnWidths{col} = pixelWidth; + end + + set(hTable, 'ColumnWidth', columnWidths); +end diff --git a/Functions/updateWaitbar.m b/Functions/updateWaitbar.m new file mode 100644 index 0000000..94accd9 --- /dev/null +++ b/Functions/updateWaitbar.m @@ -0,0 +1,45 @@ +function updateWaitbar(currentIteration, totalIterations) + +if currentIteration == 1 + + % Check if the waitbar already exists using its unique Tag + hWaitbar = findobj('Tag', 'MyUniqueWaitbar'); + + if isempty(hWaitbar) || ~ishandle(hWaitbar) + % Create a waitbar with a unique Tag if it doesn't exist + hWaitbar = waitbar(0, 'Starting process...', 'Name', 'Processing Progress', 'Tag', 'MyUniqueWaitbar'); + else + % Waitbar exists, reset the progress bar + waitbar(0, hWaitbar, 'Resuming process...'); + end + +elseif currentIteration == totalIterations+1 + % Check if the waitbar already exists using its unique Tag + hWaitbar = findobj('Tag', 'MyUniqueWaitbar'); + + % Close the waitbar after the loop is completed + if ishandle(hWaitbar) + close(hWaitbar); + end + +else + + % Check if the waitbar already exists using its unique Tag + hWaitbar = findobj('Tag', 'MyUniqueWaitbar'); + + % Calculate the progress fraction + progressFraction = currentIteration / totalIterations; + + % Update the waitbar's progress and message + if ishandle(hWaitbar) + waitbar(progressFraction, hWaitbar, ... + sprintf('Progress: %d/%d', currentIteration, totalIterations)); + else + % % If the waitbar was closed, recreate it + % hWaitbar = waitbar(progressFraction, 'Resuming process...', 'Name', 'Processing Progress', 'Tag', 'MyUniqueWaitbar'); + end + +end + + + diff --git a/projects/Lab_2024/bias_sweep_evaluation.m b/projects/Lab_2024/bias_sweep_evaluation.m new file mode 100644 index 0000000..08ffcac --- /dev/null +++ b/projects/Lab_2024/bias_sweep_evaluation.m @@ -0,0 +1,145 @@ + + +filename = "C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\bias_sweep\PAM6_10km_ffe__wh.mat"; +a = load(filename); +wh2 = a.obj; + +m = wh2.getStoValue('m',wh2.parameter.vbias.values(1),wh2.parameter.awg_vpp.values(1)); + +v_bias_vals = wh2.parameter.vbias.values; +awg_vpp_vals = wh2.parameter.awg_vpp.values; + +bers = []; +rop_measured = []; +cnt = 0; +for awg_vpp_cur = awg_vpp_vals + cnt = cnt+1; + bers(cnt,:) = wh2.getStoValue('ber',v_bias_vals,awg_vpp_cur); + rop_measured(cnt,:) = wh2.getStoValue('rop',v_bias_vals,awg_vpp_cur); +end + +[bestber,bestindex] = min(bers,[],'all'); +[awg_pos,v_bias_pos]=ind2sub(size(bers),bestindex); +bestawgvpp=awg_vpp_vals(awg_pos); +bestvbias=v_bias_vals(v_bias_pos); + +disp(['Best Vpp: ',num2str(bestvbias),' V; Best Vpp AWG: ',num2str(bestawgvpp),' V' ]) + +figure(); +sgtitle(['PAM ', num2str(m)]) +subplot1 = subplot(1,2,1); + +% Compute the logarithm of BER data +% Adding a small epsilon to avoid log(0) +epsilon = 1e-12; +log_bers = log10(bers + epsilon); + +% Set limits for z-data scaling in log scale +zmin = log10(1e-4 + epsilon); +zmax = log10(0.5 + epsilon); + +% Plot the filled contour plot with log-scaled z-data +contourf_handle = contourf(v_bias_vals, awg_vpp_vals, log_bers, 'Parent', subplot1, "ShowText",true,"LabelFormat", @mylabelfun); + +% Set x and y labels with subscripts for clarity +xlabel('V_{bias}'); +ylabel('V_{pp} AWG'); +title('BER (Mind used Equalizer!)'); + +% Adjust the grid to display white lines +grid on; +set(subplot1, 'GridColor', [1 1 1]); % Set grid color to white + +% Set limits for z-data scaling +clim([zmin zmax]); + +% Adjust the colormap +colormap(flipud(cbrewer2('RdBu',64))); + +% Add a colorbar and adjust its ticks to represent actual BER values +c = colorbar; +% Set colorbar ticks at log-spaced intervals +tick_values = [1e-4 1e-3 1e-2 1e-1 0.5]; +tick_positions = log10(tick_values + epsilon); +set(c, 'Ticks', tick_positions, 'TickLabels', arrayfun(@num2str, tick_values, 'UniformOutput', false)); + +% Store variables in the figure's application data for use in the data tip function +setappdata(gcf, 'v_bias_vals', v_bias_vals); +setappdata(gcf, 'awg_vpp_vals', awg_vpp_vals); +setappdata(gcf, 'bers', bers); +setappdata(gcf, 'power', rop_measured); % Store the Power data + +% Set up the data cursor mode to display custom data tips +dcm_obj = datacursormode(gcf); +set(dcm_obj, 'UpdateFcn', @customDataTip); +hold on +scatter(bestvbias,bestawgvpp,100,"red",'Marker','x','LineWidth',2); + + + +subplot2 = subplot(1,2,2); + +% Plot the filled contour plot +contourf_handle = contourf(v_bias_vals, awg_vpp_vals, rop_measured, 'Parent', subplot2); + +% Set x and y labels +xlabel('V_{bias}'); +ylabel('V_{pp} AWG'); +title('Power at Rx'); + +% Adjust the grid to display white lines +grid on; +set(subplot2, 'GridColor', [1 1 1]); % Set grid color to white + +% Set limits for z-data scaling +clim([-5 -3]); + +% Adjust the colormap +colormap(flipud(cbrewer2('RdBu',64))); + +% Add a colorbar +colorbar; + +% Store variables in the figure's application data for use in the data tip function +setappdata(gcf, 'v_bias_vals', v_bias_vals); +setappdata(gcf, 'awg_vpp_vals', awg_vpp_vals); +setappdata(gcf, 'bers', bers); +setappdata(gcf, 'power', rop_measured); % Store the Power data + +% Set up the data cursor mode to display custom data tips +dcm_obj = datacursormode(gcf); +set(dcm_obj, 'UpdateFcn', @customDataTip); + +function labels = mylabelfun(vals) + lab = 10.^vals; + labels = arrayfun(@(x) num2str(x, '%.1e'), lab, 'UniformOutput', false); +end + + +% Define the custom data tip function +function txt = customDataTip(~, event_obj) + % Retrieve stored variables + v_bias_vals = getappdata(gcf, 'v_bias_vals'); + awg_vpp_vals = getappdata(gcf, 'awg_vpp_vals'); + bers = getappdata(gcf, 'bers'); + power = getappdata(gcf, 'power'); % Retrieve the Power data + + % Get the position of the data cursor + pos = event_obj.Position; + xdata = pos(1); + ydata = pos(2); + + % Find the nearest indices in the data arrays + [~, xInd] = min(abs(v_bias_vals - xdata)); + [~, yInd] = min(abs(awg_vpp_vals - ydata)); + + % Get the corresponding BER and Power values + berValue = bers(yInd, xInd); + powerValue = power(yInd, xInd); % Get the Power value + + % Format the text for the data tip + txt = {['V_{bias} = ', num2str(xdata)], ... + ['V_{pp} AWG = ', num2str(ydata)], ... + ['BER = ', num2str(berValue, '%.1e')], ... + ['Power = ', num2str(powerValue)]}; +end \ No newline at end of file diff --git a/projects/Lab_2024/bias_sweep_evaluation_2.m b/projects/Lab_2024/bias_sweep_evaluation_2.m new file mode 100644 index 0000000..c310afd --- /dev/null +++ b/projects/Lab_2024/bias_sweep_evaluation_2.m @@ -0,0 +1,124 @@ + + +filename = "C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\bias_sweep_gigantisch\wh_pam4.mat"; +a = load(filename); +wh2 = a.wh; + + +v_bias_vals = wh2.parameter.vbias.values; +awg_vpp_vals = wh2.parameter.awg_vpp.values; +eq_mode_vals = wh2.parameter.eq_mode.values; +eq_mode_show = eq_mode_vals(2); +eq_modes = ["FFE","FFE+MLSE","DB precoded","DB encoded"]; + +precomp_amp_max_vals = wh2.parameter.precomp_amp_max.values; +precomp_amp_max_show = precomp_amp_max_vals(2); + +m = wh2.getStoValue('m',wh2.parameter.vbias.values(1),wh2.parameter.awg_vpp.values(1),wh2.parameter.eq_mode.values(1),wh2.parameter.precomp_amp_max.values(1)); + +figure(); +sgtitle(['PAM ', num2str(m),' | EQ: ', char(eq_modes(eq_mode_show))]) + + +for p = 1:numel(precomp_amp_max_vals) + precomp_amp_max_show = precomp_amp_max_vals(p); + subplot1 = subplot(2,3,p); + + bers = []; + rop_measured = []; + cnt = 0; + for awg_vpp_cur = awg_vpp_vals + cnt = cnt+1; + bers(cnt,:) = wh2.getStoValue('ber',v_bias_vals,awg_vpp_cur,eq_mode_show,precomp_amp_max_show); + rop_measured(cnt,:) = wh2.getStoValue('rop',v_bias_vals,awg_vpp_cur,eq_mode_show,precomp_amp_max_show); + end + + [bestber,bestindex] = min(bers,[],'all'); + [awg_pos,v_bias_pos]=ind2sub(size(bers),bestindex); + bestawgvpp=awg_vpp_vals(awg_pos); + bestvbias=v_bias_vals(v_bias_pos); + + disp(['Best Vpp: ',num2str(bestvbias),' V; Best Vpp AWG: ',num2str(bestawgvpp),' V' ]) + + % Compute the logarithm of BER data + % Adding a small epsilon to avoid log(0) + epsilon = 1e-12; + log_bers = log10(bers + epsilon); + + % Set limits for z-data scaling in log scale + zmin = log10(1e-4 + epsilon); + zmax = log10(0.5 + epsilon); + + % Plot the filled contour plot with log-scaled z-data + contourf_handle = contourf(v_bias_vals, awg_vpp_vals, log_bers, 'Parent', subplot1, "ShowText",true,"LabelFormat", @mylabelfun); + + % Set x and y labels with subscripts for clarity + xlabel('V_{bias}'); + ylabel('V_{pp} AWG'); + title(['Prec. Ampl.: ',num2str(precomp_amp_max_show), 'dB']); + + % Adjust the grid to display white lines + grid on; + set(subplot1, 'GridColor', [1 1 1]); % Set grid color to white + + % Set limits for z-data scaling + clim([zmin zmax]); + + % Adjust the colormap + colormap(flipud(cbrewer2('RdBu',64))); + + % Add a colorbar and adjust its ticks to represent actual BER values + c = colorbar; + % Set colorbar ticks at log-spaced intervals + tick_values = [1e-4 1e-3 1e-2 1e-1 0.5]; + tick_positions = log10(tick_values + epsilon); + set(c, 'Ticks', tick_positions, 'TickLabels', arrayfun(@num2str, tick_values, 'UniformOutput', false)); + + % Store variables in the figure's application data for use in the data tip function + setappdata(gcf, 'v_bias_vals', v_bias_vals); + setappdata(gcf, 'awg_vpp_vals', awg_vpp_vals); + setappdata(gcf, 'bers', bers); + setappdata(gcf, 'power', rop_measured); % Store the Power data + + % Set up the data cursor mode to display custom data tips + dcm_obj = datacursormode(gcf); + set(dcm_obj, 'UpdateFcn', @customDataTip); + hold on + scatter(bestvbias,bestawgvpp,100,"red",'Marker','x','LineWidth',2); + +end + + +function labels = mylabelfun(vals) + lab = 10.^vals; + labels = arrayfun(@(x) num2str(x, '%.1e'), lab, 'UniformOutput', false); +end + + +% Define the custom data tip function +function txt = customDataTip(~, event_obj) + % Retrieve stored variables + v_bias_vals = getappdata(gcf, 'v_bias_vals'); + awg_vpp_vals = getappdata(gcf, 'awg_vpp_vals'); + bers = getappdata(gcf, 'bers'); + power = getappdata(gcf, 'power'); % Retrieve the Power data + + % Get the position of the data cursor + pos = event_obj.Position; + xdata = pos(1); + ydata = pos(2); + + % Find the nearest indices in the data arrays + [~, xInd] = min(abs(v_bias_vals - xdata)); + [~, yInd] = min(abs(awg_vpp_vals - ydata)); + + % Get the corresponding BER and Power values + berValue = bers(yInd, xInd); + powerValue = power(yInd, xInd); % Get the Power value + + % Format the text for the data tip + txt = {['V_{bias} = ', num2str(xdata)], ... + ['V_{pp} AWG = ', num2str(ydata)], ... + ['BER = ', num2str(berValue, '%.1e')], ... + ['Power = ', num2str(powerValue)]}; +end \ No newline at end of file diff --git a/projects/Lab_2024/lab_baudrate_sweep.m b/projects/Lab_2024/lab_baudrate_sweep.m new file mode 100644 index 0000000..43b88a2 --- /dev/null +++ b/projects/Lab_2024/lab_baudrate_sweep.m @@ -0,0 +1,274 @@ + +folderpath = 'C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\baudrate_sweep\'; +experiment_name = 'PAM4_10km_ffe_'; + +ffe_only = 0; +postfilter_approach = 0; +db_channel_approach = 1; +db_coding_approach = 0; + +db_precode = db_coding_approach || db_channel_approach; + +%%% SIR Sweep for MPI Experiment %%% +params = struct; +params.fsym = [56,68,80,92].*1e9; +params.fsym = [92].*1e9; + +wh = DataStorage(params); + +wh.addStorage("ber"); +wh.addStorage("pd_in"); +wh.addStorage("rop"); +wh.addStorage("m"); +wh.addStorage("signals"); + +precomp_path = "C:\Users\sioe\Documents\MATLAB\imdd_simulation\projects\standard_system\"; +precomp_fn = "lab_mpi_setup_2"; +precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active +precomp_amp_max = 2; + +M = 4; +pn_key = 2; +usemrds = 0; +fdac = 92e9; +awg_vpp = 0.15; +fadc = 160e9; +rrcalpha = 0.05; +v_bias = 2.25; +pd_in_set = 4; +rop_atten = 0; + +disp(['Start Measurement of ',num2str(prod(wh.dim)),' loops...']) + +for fsym = wh.parameter.fsym.values + + loop_name = ['_fsym_',num2str(fsym)]; + + %%%%% SET Voltages %%%%%% + dcs = DC_supply("active",[1,1],"voltage",[v_bias, 9]); + dcs.set("voltage",[v_bias, 9]); + + %%%%% 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 %%%%%% + SCP = ScopeKeysight("model","DSAZ634A",'autoscale',1,"fadc",'GSa_160',"channel",[1,0],"recordLen",2000000,"removeDC",1); + AWG = AwgKeysight("model","M8196A","fdac",fdac,"scaletodac",[1,1,1,1],"skews",[0,0,0,0],"voltages",[0,0,0,awg_vpp]); + A2S = Awg2Scope(AWG,SCP,[0,0,0,1]); + + + %%%%% Symbol Generation %%%%%% + Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rrcalpha); + + [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",18,"useprbs",1,... + "fs_out",fdac,"applyclipping",0,"clipfactor",1.7,... + "applypulseform",0,"pulseformer",Pform,"randkey",pn_key,... + "db_precode",db_precode,... + "mrds_code",usemrds,"mrds_blocklength",512,"db_encode",db_coding_approach).process(); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 % measure channel + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.buildOFDM(); + elseif precomp_mode == 2 % apply precomp + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.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",10); + + save([folderpath,[experiment_name,'bits'],loop_name],"Bits"); + save([folderpath,[experiment_name,'symbols'],loop_name],"Symbols"); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %%%%% AWG --> Scope %%%%%% + [~,~,~,Scpe_sig] = A2S.process("signal4",Digi_sig); + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + Scpe_sig.spectrum("displayname","Scope PSD","fignum",20); + + Scpe_sig.plot("displayname","Scope raw signal","fignum",25); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + + %%%%%% Sample to 2x fsym %%%%%% + Scpe_sig = Scpe_sig.resample("fs_in",160e9,"fs_out",2*fsym); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 + freqresp.estimate(Scpe_sig,"save",true,"savePath",precomp_path,"fileName",precomp_fn); + freqresp.plot(); + end + + %%%%%% Sync Rx signal with reference %%%%%% + [Scpe_sig,S] = Scpe_sig.tsynch("reference",Symbols,"fs_ref",fsym); + + %%%%%% SNR CHEAT - Avg. the measured signal occurences %%%%%% + 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; + end + scope_mean = scope_mean ./ n; + Scpe_sig.signal = scope_mean; + end + + %%%%% Plot and Save Routine 2 %%%%%%%%%%%%%%%%%%%%%%%%% + save([folderpath,experiment_name,'rx_signal',loop_name],"S"); + Scpe_sig.eye(fsym,M,"fignum",40,"displayname",' after Scope'); + + voa.readvals(); + rop = voa.power_state(1); + pd_in = voa.power_state(2); + + + %%%%% 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); + + if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ'); + + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber),'| ROP: ',num2str(rop),' dB | PD_in: ',num2str(pd_in),' dBm']); + + if 1 + figure(53); + constellation = unique(Symbols.signal); + received = NaN(numel(constellation),length(Symbols)); + for lvl = 1:numel(constellation) + received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig.signal(Symbols.signal==constellation(lvl)); + hold on + histogram(received(lvl,:),1000,"EdgeAlpha",0); + end + end + + elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + + Noi = EQ_sig-Symbols; + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber_ffe_only,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); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber_ffe_only),' -> PF -> MLSE: ',sprintf('%.1E',ber),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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); + [~,errors_bm,ber,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,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,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),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + + end + + wh.addValueToStorage(ber,'ber',fsym); + wh.addValueToStorage(rop,'rop',fsym); + wh.addValueToStorage(pd_in,'pd_in',fsym); + wh.addValueToStorage(Rx_bits,'signals',fsym); + wh.addValueToStorage(M,'m',fsym); + + showCurrentMeasurement('BER', ber, 'ROP', rop, 'PD in', pd_in, 'PAM',M, 'Vbias', v_bias, 'AWG Vpp', Awg_vpp); + autoArrangeFigures(3,3,2); +end + +wh.save([folderpath,experiment_name,'_wh']); + +cols = linspecer(8); + +fsym_vals = wh.parameter.fsym.values; + +bers = wh.getStoValue('ber',fsym_vals); +rop_measured = wh.getStoValue('rop',fsym_vals); +pd_in_measured = wh.getStoValue('pd_in',fsym_vals); + +figure(90); +hold on; % Retain the plot so new points can be added without complete redraw + +% Plot the data and get the line handle +hLine = plot(fsym_vals.*1e-9, bers, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", experiment_name); + +% Store pd_in_measured in the ZData property +hLine.ZData = pd_in_measured; + +% Customize the data tips +% Set labels for existing data tip rows +hLine.DataTipTemplate.DataTipRows(1).Label = 'Fsym'; +hLine.DataTipTemplate.DataTipRows(2).Label = 'BER'; +hLine.DataTipTemplate.DataTipRows(2).Format = '%.2e'; % Format BER as "3e-4" + + +% Add a new data tip row for PDin +pdinRow = dataTipTextRow('PDin', 'ZData'); +hLine.DataTipTemplate.DataTipRows(3) = pdinRow; + +% Continue with the rest of your plot settings +yline(3.8e-3, 'DisplayName', 'HD-FEC', 'LineStyle', '--', 'HandleVisibility', 'off'); +xlabel('Symbol Rate'); +ylabel('Bit Error Rate (BER)'); +title('Bit Error Rate vs. ROP'); +set(gca, 'yscale', 'log'); +set(gca, 'Box', 'on'); +grid on; +grid minor; +legend('Interpreter', 'none'); + +autoArrangeFigures(3,3,2) + diff --git a/projects/Lab_2024/lab_bias_sweep.m b/projects/Lab_2024/lab_bias_sweep.m new file mode 100644 index 0000000..b22f415 --- /dev/null +++ b/projects/Lab_2024/lab_bias_sweep.m @@ -0,0 +1,269 @@ + +folderpath = 'C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\bias_sweep\'; +experiment_name = 'PAM4_DB_precoded_10km_ffe_'; + +% a = load([folderpath,experiment_name,'_wh']); +% wh2 = a.obj; + + +ffe_only = 0; +postfilter_approach = 0; +db_channel_approach = 1; +db_coding_approach = 0; + +db_precode = db_coding_approach || db_channel_approach; + +%%% SIR Sweep for MPI Experiment %%% +params = struct; + +params.vbias = [2.1:0.05:2.4]; +params.awg_vpp = [0.15:0.05:0.6]; + +wh = DataStorage(params); + +wh.addStorage("ber"); +wh.addStorage("pd_in"); +wh.addStorage("rop"); +wh.addStorage("m"); +wh.addStorage("signals"); + +precomp_path = "C:\Users\sioe\Documents\MATLAB\imdd_simulation\projects\standard_system\"; +precomp_fn = "lab_mpi_setup_2"; +precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active +precomp_amp_max = -2; + +M = 4; +pn_key = 2; +usemrds = 0; +fsym = 92e9; +fdac = 92e9; +Awg_vpp = 0.35; +fadc = 160e9; +rrcalpha = 0.05; +v_bias = 2.25; +pd_in_set = 6; +rop_atten = 0; + +looptatal = prod(wh.dim); +disp(['Start Measurement of ',num2str(looptatal),' loops...']) +hWaitbar = waitbar(0, 'Starting measurement...', 'Name', 'Processing Progress'); +loopcnt = 0; + +for v_bias = wh.parameter.vbias.values + for awg_vpp = wh.parameter.awg_vpp.values + + loopcnt = loopcnt+1; + progressFraction = loopcnt / looptatal; + waitbar(progressFraction, hWaitbar, ... + sprintf('Progress: %d/%d', loopcnt, looptatal)); + + try + a + % ber = wh2.getStoValue('ber', v_bias,awg_vpp); + % rop = wh2.getStoValue('rop', v_bias,awg_vpp); + % pd_in = wh2.getStoValue('pd_in', v_bias,awg_vpp); + % Rx_bits = wh2.getStoValue('signals', v_bias,awg_vpp); + % Rx_bits = Rx_bits{1}; + % M = wh2.getStoValue('m', v_bias,awg_vpp); + + + catch + + loop_name = ['_fsym_',num2str(fsym)]; + + %%%%% SET Voltages %%%%%% + dcs = DC_supply("active",[1,1],"voltage",[v_bias, 9]); + dcs.set("voltage",[v_bias, 9]); + + %%%%% 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 %%%%%% + SCP = ScopeKeysight("model","DSAZ634A",'autoscale',1,"fadc",'GSa_160',"channel",[1,0],"recordLen",1000000,"removeDC",1); + AWG = AwgKeysight("model","M8196A","fdac",fdac,"scaletodac",[1,1,1,1],"skews",[0,0,0,0],"voltages",[0,0,0,awg_vpp]); + A2S = Awg2Scope(AWG,SCP,[0,0,0,1]); + + + %%%%% Symbol Generation %%%%%% + Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rrcalpha); + + [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",18,"useprbs",1,... + "fs_out",fdac,"applyclipping",0,"clipfactor",1.7,... + "applypulseform",0,"pulseformer",Pform,"randkey",pn_key,... + "db_precode",db_precode,... + "mrds_code",usemrds,"mrds_blocklength",512,"db_encode",db_coding_approach).process(); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 % measure channel + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.buildOFDM(); + elseif precomp_mode == 2 % apply precomp + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.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",10); + + % save([folderpath,[experiment_name,'bits'],loop_name],"Bits"); + % save([folderpath,[experiment_name,'symbols'],loop_name],"Symbols"); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %%%%% AWG --> Scope %%%%%% + [~,~,~,Scpe_sig] = A2S.process("signal4",Digi_sig); + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + % Scpe_sig.spectrum("displayname","Scope PSD","fignum",20); + + % Scpe_sig.plot("displayname","Scope raw signal","fignum",25); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + + %%%%%% Sample to 2x fsym %%%%%% + Scpe_sig = Scpe_sig.resample("fs_in",160e9,"fs_out",2*fsym); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 + freqresp.estimate(Scpe_sig,"save",true,"savePath",precomp_path,"fileName",precomp_fn); + freqresp.plot(); + end + + %%%%%% Sync Rx signal with reference %%%%%% + [Scpe_sig,S] = Scpe_sig.tsynch("reference",Symbols,"fs_ref",fsym); + + %%%%%% SNR CHEAT - Avg. the measured signal occurences %%%%%% + 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; + end + scope_mean = scope_mean ./ n; + Scpe_sig.signal = scope_mean; + end + + %%%%% Plot and Save Routine 2 %%%%%%%%%%%%%%%%%%%%%%%%% + save([folderpath,experiment_name,'rx_signal',loop_name],"S"); + Scpe_sig.eye(fsym,M,"fignum",40,"displayname",' after Scope'); + + voa.readvals(); + rop = voa.power_state(1); + pd_in = voa.power_state(2); + + + %%%%% 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); + + if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ'); + + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber),'| ROP: ',num2str(rop),' dB | PD_in: ',num2str(pd_in),' dBm']); + + if 1 + figure(53); + constellation = unique(Symbols.signal); + received = NaN(numel(constellation),length(Symbols)); + for lvl = 1:numel(constellation) + received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig.signal(Symbols.signal==constellation(lvl)); + hold on + histogram(received(lvl,:),1000,"EdgeAlpha",0); + end + end + + elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + + Noi = EQ_sig-Symbols; + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber_ffe_only,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); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber_ffe_only),' -> PF -> MLSE: ',sprintf('%.1E',ber),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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); + [~,errors_bm,ber,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),' | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_coding_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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,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),' | PD_in: ',num2str(pd_in),' dBm']); + + end + + end + + wh.addValueToStorage(ber,'ber',v_bias,awg_vpp); + wh.addValueToStorage(rop,'rop',v_bias,awg_vpp); + wh.addValueToStorage(pd_in,'pd_in',v_bias,awg_vpp); + wh.addValueToStorage(Rx_bits,'signals',v_bias,awg_vpp); + wh.addValueToStorage(M,'m',v_bias,awg_vpp); + + showCurrentMeasurement('BER', ber, 'ROP', rop, 'PD in', pd_in, 'PAM',M, 'Vbias', v_bias, 'AWG Vpp', awg_vpp, 'Precomp MaxAmp',precomp_amp_max); + autoArrangeFigures(3,3,2); + + + end +end + +close(hWaitbar); + +wh.save([folderpath,experiment_name,'_wh']); + +autoArrangeFigures(3,3,2) + diff --git a/projects/Lab_2024/lab_bias_sweep_gigantisch.m b/projects/Lab_2024/lab_bias_sweep_gigantisch.m new file mode 100644 index 0000000..d00c0ca --- /dev/null +++ b/projects/Lab_2024/lab_bias_sweep_gigantisch.m @@ -0,0 +1,320 @@ + +folderpath = 'C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\bias_sweep_4db_pdin\'; +experiment_name = 'PAM6_alles_10km_ffe_'; + +wh2 = obj; + +ffe_only = 0; +postfilter_approach = 0; +db_channel_approach = 1; +db_coding_approach = 0; + +db_precode = db_coding_approach || db_channel_approach; + +%%% SIR Sweep for MPI Experiment %%% +params = struct; + +params.vbias = [2.2:0.05:2.5]; +params.awg_vpp = [0.15:0.05:0.6]; +params.eq_mode = [1]; +params.precomp_amp_max = [0:2:5]; + + +wh = DataStorage(params); + +wh.addStorage("ber"); +wh.addStorage("pd_in"); +wh.addStorage("rop"); +wh.addStorage("m"); +wh.addStorage("signals"); + +precomp_path = "C:\Users\sioe\Documents\MATLAB\imdd_simulation\projects\standard_system\"; +precomp_fn = "lab_mpi_setup_2"; +precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active +precomp_amp_max = -2; + +M = 6; +pn_key = 2; +usemrds = 0; +fsym = 70e9; +fdac = 92e9; +awg_vpp = 0.35; +fadc = 160e9; +rrcalpha = 0.05; +v_bias = 2.25; +pd_in_set = 4; +rop_atten = 0; + +looptatal = prod(wh.dim); +iterationTimes = zeros(looptatal, 1); % Preallocate for speed + +disp(['Start Measurement of ',num2str(looptatal),' loops...']) +hWaitbar = waitbar(0, 'Starting measurement...', 'Name', 'Processing Progress'); +loopcnt = 0; + +estimatedTimeRemaining = 0; +estimatedTotalTime = 0; +for eq_mode = wh.parameter.eq_mode.values + for precomp_amp_max = wh.parameter.precomp_amp_max.values + for v_bias = wh.parameter.vbias.values + for awg_vpp = wh.parameter.awg_vpp.values + + iterationStartTime = tic; + experiment_name = ['PAM6_alles_10km_eq',num2str(eq_mode),'maxamp',num2str(precomp_amp_max),'vbias',num2str(v_bias),'awgvpp',num2str(awg_vpp)]; + experiment_name = strrep(experiment_name,'.','_'); + + loopcnt = loopcnt+1; + progressFraction = loopcnt / looptatal; + waitbar(progressFraction, hWaitbar, ... + sprintf('Progress: %d/%d\nEstimated time remaining: %.2f hours\nEstimated time remaining: %.2f hours', ... + loopcnt, looptatal, estimatedTimeRemaining/60/60, estimatedTotalTime/60/60)); + + + loop_name = ['_fsym_',num2str(fsym)]; + + try + + a + ber = wh2.getStoValue('ber',v_bias,awg_vpp,eq_mode,precomp_amp_max); + assert(~isempty(ber),'err') + rop = wh2.getStoValue('rop',v_bias,awg_vpp,eq_mode,precomp_amp_max); + assert(~isempty(rop)) + pd_in = wh2.getStoValue('pd_in',v_bias,awg_vpp,eq_mode,precomp_amp_max); + assert(~isempty(pd_in)) + M = wh2.getStoValue('m',v_bias,awg_vpp,eq_mode,precomp_amp_max); + assert(~isempty(M)) + + catch + + + + + switch eq_mode + + case 2 + ffe_only = 0; + postfilter_approach = 1; + db_channel_approach = 0; + db_coding_approach = 0; + + db_precode = db_coding_approach || db_channel_approach; + case 3 + ffe_only = 0; + postfilter_approach = 0; + db_channel_approach = 1; + db_coding_approach = 0; + + db_precode = db_coding_approach || db_channel_approach; + case 4 + ffe_only = 0; + postfilter_approach = 0; + db_channel_approach = 0; + db_coding_approach = 1; + + db_precode = db_coding_approach || db_channel_approach; + end + + + + %%%%% SET Voltages %%%%%% + dcs = DC_supply("active",[1,1],"voltage",[v_bias, 9]); + dcs.set("voltage",[v_bias, 9]); + + %%%%% 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 %%%%%% + SCP = ScopeKeysight("model","DSAZ634A",'autoscale',1,"fadc",'GSa_160',"channel",[1,0],"recordLen",1000000,"removeDC",1); + AWG = AwgKeysight("model","M8196A","fdac",fdac,"scaletodac",[1,1,1,1],"skews",[0,0,0,0],"voltages",[0,0,0,awg_vpp]); + A2S = Awg2Scope(AWG,SCP,[0,0,0,1]); + + %%%%% Symbol Generation %%%%%% + Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rrcalpha); + + [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",18,"useprbs",1,... + "fs_out",fdac,"applyclipping",0,"clipfactor",1.7,... + "applypulseform",0,"pulseformer",Pform,"randkey",pn_key,... + "db_precode",db_precode,... + "mrds_code",usemrds,"mrds_blocklength",512,"db_encode",db_coding_approach).process(); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 % measure channel + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.buildOFDM(); + elseif precomp_mode == 2 % apply precomp + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.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",10); + + % save([folderpath,[experiment_name,'bits'],loop_name],"Bits"); + % save([folderpath,[experiment_name,'symbols'],loop_name],"Symbols"); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %%%%% AWG --> Scope %%%%%% + [~,~,~,Scpe_sig] = A2S.process("signal4",Digi_sig); + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + % Scpe_sig.spectrum("displayname","Scope PSD","fignum",20); + + % Scpe_sig.plot("displayname","Scope raw signal","fignum",25); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + + %%%%%% Sample to 2x fsym %%%%%% + Scpe_sig = Scpe_sig.resample("fs_in",160e9,"fs_out",2*fsym); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 + freqresp.estimate(Scpe_sig,"save",true,"savePath",precomp_path,"fileName",precomp_fn); + freqresp.plot(); + end + + %%%%%% Sync Rx signal with reference %%%%%% + [Scpe_sig,S] = Scpe_sig.tsynch("reference",Symbols,"fs_ref",fsym); + + %%%%%% SNR CHEAT - Avg. the measured signal occurences %%%%%% + 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; + end + scope_mean = scope_mean ./ n; + Scpe_sig.signal = scope_mean; + end + + %%%%% Plot and Save Routine 2 %%%%%%%%%%%%%%%%%%%%%%%%% + save([folderpath,experiment_name,'rx_signal',loop_name],"S"); + % Scpe_sig.eye(fsym,M,"fignum",40,"displayname",' after Scope'); + + voa.readvals(); + rop = voa.power_state(1); + pd_in = voa.power_state(2); + + + %%%%% 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); + + if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + % EQ_sig.plot("fignum",50,"displayname",'After EQ'); + + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber),'| ROP: ',num2str(rop),' dB | PD_in: ',num2str(pd_in),' dBm']); + + if 0 + figure(53); + constellation = unique(Symbols.signal); + received = NaN(numel(constellation),length(Symbols)); + for lvl = 1:numel(constellation) + received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig.signal(Symbols.signal==constellation(lvl)); + hold on + histogram(received(lvl,:),1000,"EdgeAlpha",0); + end + end + + elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + % EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + + Noi = EQ_sig-Symbols; + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber_ffe_only,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); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber_ffe_only),' -> PF -> MLSE: ',sprintf('%.1E',ber),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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); + [~,errors_bm,ber,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),' | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_coding_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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,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),' | PD_in: ',num2str(pd_in),' dBm']); + + end + + end + + wh.addValueToStorage(ber,'ber',v_bias,awg_vpp,eq_mode,precomp_amp_max); + wh.addValueToStorage(rop,'rop',v_bias,awg_vpp,eq_mode,precomp_amp_max); + wh.addValueToStorage(pd_in,'pd_in',v_bias,awg_vpp,eq_mode,precomp_amp_max); + %wh.addValueToStorage(Rx_bits,'signals',v_bias,awg_vpp); + wh.addValueToStorage(M,'m',v_bias,awg_vpp,eq_mode,precomp_amp_max); + + showCurrentMeasurement('BER', ber, 'ROP', rop, 'PD in', pd_in, 'PAM',M, 'Vbias', v_bias, 'AWG Vpp', awg_vpp, 'Precomp MaxAmp',precomp_amp_max); + + iterationTimes(loopcnt) = toc(iterationStartTime); + averageTimePerIteration = mean(iterationTimes(1:loopcnt)); + estimatedTotalTime = averageTimePerIteration * looptatal; + estimatedTimeRemaining = estimatedTotalTime - sum(iterationTimes(1:loopcnt)); + %autoArrangeFigures(3,3,2); + + end + end + end +end + +close(hWaitbar); + +wh.save([folderpath,experiment_name,'wh']); + +autoArrangeFigures(3,3,2) + +disp("measurement done") \ No newline at end of file diff --git a/projects/Lab_2024/lab_db_precode.m b/projects/Lab_2024/lab_db_precode.m index 1cdcd1a..91900e7 100644 --- a/projects/Lab_2024/lab_db_precode.m +++ b/projects/Lab_2024/lab_db_precode.m @@ -1,11 +1,11 @@ -folderpath = 'C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\mpi_ofc_2024\'; -experiment_name = '10km_db_transmit_no_mpi_'; +folderpath = 'C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\no_mpi_2024\'; +experiment_name = '10km_ffe_no_mpi_'; only_dsp = 0; -ffe_only = 0; -postfilter_approach = 1; +ffe_only = 1; +postfilter_approach = 0; db_channel_approach = 0; db_coding_approach = 0; @@ -25,34 +25,35 @@ wh.addStorage("signals"); precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active precomp_amp_max = 3; -M = 4; +M = 6; pn_key = 2; usemrds = 0; fdac = 92e9; -fsym = 92e9; +fsym = 56e9; fadc = 160e9; rrcalpha = 0.05; v_bias = 2.25; i_atten = params.i_atten(1); -pd_in_desired = 7; +pd_in_desired = 6; -if ~only_dsp - disp(['Start Measurement of ',num2str(prod(wh.dim)),' loops...']) - +disp(['Start Measurement of ',num2str(prod(wh.dim)),' loops...']) + +for i = 1 + %%%%% SET Volatges %%%%%% dcs = DC_supply("active",[1,1],"voltage",[v_bias, 9]); dcs.set("voltage",[v_bias, 9]); - + voa = OptAtten("active",[1,2,1,1],"value",[0,pd_in_desired,0,i_atten],"wavelength",[1310,1310,1310,1310]); - % voa.set('active',[1,2,1,1],'value',[0,pd_in_desired,0,i_atten]); + %voa.set('active',[1,2,1,1],'value',[0,pd_in_desired,0,i_atten]); % voa.readvals(); - + SCP = ScopeKeysight("model","DSAZ634A",'autoscale',1,"fadc",'GSa_160',"channel",[1,0],"recordLen",2000000,"removeDC",1); AWG = AwgKeysight("model","M8196A","fdac",fdac,"scaletodac",[1,1,1,1],"skews",[0,0,0,0],"voltages",[0,0,0,0.62]); A2S = Awg2Scope(AWG,SCP,[0,0,0,1]); - - + + if 1 precomp_path = "C:\Users\sioe\Documents\MATLAB\imdd_simulation\projects\standard_system\"; precomp_fn = "lab_mpi_setup_2"; @@ -60,17 +61,17 @@ if ~only_dsp precomp_path = "C:\Users\sioe\Documents\MATLAB\model-collection\sioe_models\Labor_2024\Lab_PAM4\"; precomp_fn = "precomp_bla__loop1_1"; end - - + + %%%%% Symbol Generation %%%%%% Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rrcalpha); - + [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",18,"useprbs",1,... "fs_out",fdac,"applyclipping",0,"clipfactor",1.7,... "applypulseform",0,"pulseformer",Pform,"randkey",pn_key,... "db_precode",db_precode,... "mrds_code",usemrds,"mrds_blocklength",512,"db_encode",db_coding_approach).process(); - + if precomp_mode == 1 % measure channel freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); Digi_sig = freqresp.buildOFDM(); @@ -78,50 +79,42 @@ if ~only_dsp freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); Digi_sig = freqresp.precomp(Digi_sig,'maxampdb',precomp_amp_max,'loadPath',precomp_path,'fileName',precomp_fn); end - + Digi_sig = Digi_sig.resample("fs_out",AWG.fdac); Digi_sig.spectrum("displayname","Normal Tx","fignum",10); - % Digi_sig = Filter('filtdegree',1,"f_cutoff",45e9,"fs",Digi_sig.fs,"filterType",filtertypes.butterworth,"active",true).process(Digi_sig); - % Digi_sig.spectrum("displayname","Lowpass Tx","fignum",999); - % Digi_sig.eye(fsym,M); - - save([folderpath,[experiment_name,'bits']],"Bits"); save([folderpath,[experiment_name,'symbols']],"Symbols"); -end -for i_atten = wh.parameter.i_atten.values - - if ~only_dsp + for i_atten = wh.parameter.i_atten.values %%%%% SET ATTENUATOR %%%%%% voa.set('active',[1,2,1,1],'value',[0,7,0,i_atten]); - + %%%%% AWG --> Scope %%%%%% [~,~,~,Scpe_sig] = A2S.process("signal4",Digi_sig); Scpe_sig.spectrum("displayname","Scope PSD","fignum",20); - - + + % Scpe_sig.spectrum("displayname",'Rx Signal','fignum',10); - + %%%%%% Sample to 2x fsym %%%%%% Scpe_sig = Scpe_sig.resample("fs_in",160e9,"fs_out",2*fsym); - + if precomp_mode == 1 freqresp.estimate(Scpe_sig,"save",true,"savePath",precomp_path,"fileName",precomp_fn); freqresp.plot(); end - + %%%%%% Sync Rx signal with reference %%%%%% [Scpe_sig,S] = Scpe_sig.tsynch("reference",Symbols,"fs_ref",fsym); save([folderpath,experiment_name,'rx_signal_iatten_',num2str(i_atten),''],"S"); - + average_signals = 0; if average_signals scope_mean = zeros(size(S{1}.signal)); @@ -135,96 +128,96 @@ for i_atten = wh.parameter.i_atten.values Scpe_sig.plot("displayname","Scope PSD","fignum",30); Scpe_sig.eye(fsym,M,"fignum",40,"displayname",' after Scope'); - sir = voa.power_state(3)-voa.power_state(4); pd_in = voa.power_state(2); - end - %%%%% 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); + %%%%% 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.0,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1); - if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% + if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% - [EQ_sig] = Eq.process(Scpe_sig,Symbols); + [EQ_sig] = Eq.process(Scpe_sig,Symbols); - EQ_sig.plot("fignum",50,"displayname",'After EQ'); + EQ_sig.plot("fignum",50,"displayname",'After EQ'); - Rx_bits = PAMmapper(M,0).demap(EQ_sig); - [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - disp(['FFE: ',sprintf('%.1E',ber),'| SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + disp(['FFE: ',sprintf('%.1E',ber),'| SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); - elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% - [EQ_sig] = Eq.process(Scpe_sig,Symbols); + [EQ_sig] = Eq.process(Scpe_sig,Symbols); - EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); - Noi = EQ_sig-Symbols; + Noi = EQ_sig-Symbols; - Rx_bits = PAMmapper(M,0).demap(EQ_sig); - [~,errors_bm,ber_ffe_only,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber_ffe_only,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - nc = 2; - burg_coeff = arburg(Noi.signal,nc); + nc = 2; + burg_coeff = arburg(Noi.signal,nc); - EQ_sig = EQ_sig.filter(burg_coeff,1); + 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); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber_ffe_only),' -> PF -> MLSE: ',sprintf('%.1E',ber),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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); + [~,errors_bm,ber,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),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_coding_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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,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),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); - 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); - [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - disp(['FFE: ',sprintf('%.1E',ber_ffe_only),' -> PF -> MLSE: ',sprintf('%.1E',ber),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); - - elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% - - [EQ_sig, Noi] = Eq.process(Scpe_sig,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); - [~,errors_bm,ber,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),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); - - elseif db_coding_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% - - [EQ_sig, Noi] = Eq.process(Scpe_sig,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,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),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + wh.addValueToStorage(ber,'ber',i_atten); + wh.addValueToStorage(sir,'sir',i_atten); + wh.addValueToStorage(pd_in,'pd_in',i_atten); + wh.addValueToStorage(Rx_bits,'signals',i_atten); end - - - - wh.addValueToStorage(ber,'ber',i_atten); - wh.addValueToStorage(sir,'sir',i_atten); - wh.addValueToStorage(pd_in,'pd_in',i_atten); - wh.addValueToStorage(Rx_bits,'signals',i_atten); - end wh.save([folderpath,experiment_name,'_wh']); diff --git a/projects/Lab_2024/lab_modulator_tf_sweep.m b/projects/Lab_2024/lab_modulator_tf_sweep.m new file mode 100644 index 0000000..6c41e4a --- /dev/null +++ b/projects/Lab_2024/lab_modulator_tf_sweep.m @@ -0,0 +1,68 @@ +% Initialize a structure to hold parameters +params = struct; + +% Define the bias voltage range from 1.2V to 2.8V with 0.01V increments +params.v_bias = 1.2:0.01:2.8; + +% Create a DataStorage object with the defined parameters +wh = DataStorage(params); + +% Add a storage field for output power measurements +wh.addStorage("p_out"); + +% Display the total number of measurement loops to be executed +disp(['Start Measurement of ', num2str(prod(wh.dim)), ' loops...']); + +% Initialize the Optical Attenuator (VOA) with specified settings +voa = OptAtten(... + "active", [1, 0, 0, 0], ... % Activate only the first channel + "value", [0, 0, 0, 0], ... % Set attenuation values to 0 dB + "wavelength", [1310, 1310, 1310, 1310]); % Set the wavelength for each channel + +% Initialize the DC Power Supply with specified settings +dcs = DC_supply(... + "active", [1, 1], ... % Activate the first two channels + "voltage", [v_bias, 9]); % Set initial voltages for channels + +% Set the VOA active channels and attenuation values +voa.set('active', [1, 0, 0, 0], 'value', [0, 0, 0, 0]); + +% Loop over each bias voltage value to perform measurements +for v_bias = wh.parameter.v_bias.values + try + % Try to retrieve existing output power measurement to avoid repetition + p_out = wh2.getStoValue('p_out', v_bias); + wh.addValueToStorage(p_out, 'p_out', v_bias); + catch + % If no existing measurement, proceed with the measurement + + %%%%% SET Voltages %%%%%% + % Update the DC supply voltage for the current bias voltage + dcs.set("voltage", [v_bias, 9]); + + %%%%% Measure VOA %%%%%% + % Read the current values from the VOA + voa.readvals(); + + % Store the measured output power in the DataStorage object + wh.addValueToStorage(voa.power_state(1), 'p_out', v_bias); + end + + % Retrieve all stored output power measurements up to the current point + p_out = wh.getStoValue('p_out', wh.parameter.v_bias.values); + + % Plot the measured output power in dBm versus the negative bias voltage + figure(90); + plot(-wh.parameter.v_bias.values(1:numel(p_out)), p_out, 'DisplayName', 'Measured Output Power in dBm'); + xlim([-max(params.v_bias), -min(params.v_bias)]); % Set x-axis limits + grid on; % Enable grid for better readability +end + +% After completing the measurements, retrieve all output power data +p_out = wh.getStoValue('p_out', wh.parameter.v_bias.values); + +% Plot the output power converted from dBm to linear scale (Watts) +figure(91); +plot(-wh.parameter.v_bias.values(1:numel(p_out)), db2pow(p_out), 'DisplayName', 'Measured Output Power in Watts'); +xlim([-max(params.v_bias), -min(params.v_bias)]); % Set x-axis limits +grid on; % Enable grid for better readability diff --git a/projects/Lab_2024/lab_precompensation_sweep.m b/projects/Lab_2024/lab_precompensation_sweep.m new file mode 100644 index 0000000..1786ebc --- /dev/null +++ b/projects/Lab_2024/lab_precompensation_sweep.m @@ -0,0 +1,265 @@ + +folderpath = 'C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\precompensation_sweep\'; +experiment_name = 'PAM4_DBencode_10km_'; + +% a = load([folderpath,experiment_name,'_wh']); +% wh2 = a.obj; + + +ffe_only = 0; +postfilter_approach = 0; +db_channel_approach = 0; +db_coding_approach = 1; + +db_precode = db_coding_approach || db_channel_approach; + +%%% SIR Sweep for MPI Experiment %%% +params = struct; + +params.precomp_amp_max = [-6:6]; + +wh = DataStorage(params); + +wh.addStorage("ber"); +wh.addStorage("pd_in"); +wh.addStorage("rop"); +wh.addStorage("m"); +wh.addStorage("signals"); + +precomp_path = "C:\Users\sioe\Documents\MATLAB\imdd_simulation\projects\standard_system\"; +precomp_fn = "lab_mpi_setup_2"; +precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active +% precomp_amp_max = -2; + +M = 4; +pn_key = 2; +usemrds = 0; +fsym = 92e9; +fdac = 92e9; +awg_vpp = 0.2; +fadc = 160e9; +rrcalpha = 0.05; +v_bias = 2.35; +pd_in_set = 6; +rop_atten = 0; + +looptatal = prod(wh.dim); +disp(['Start Measurement of ',num2str(looptatal),' loops...']) +hWaitbar = waitbar(0, 'Starting measurement...', 'Name', 'Processing Progress'); +loopcnt = 0; + +for precomp_amp_max = wh.parameter.precomp_amp_max.values + + loopcnt = loopcnt+1; + progressFraction = loopcnt / looptatal; + waitbar(progressFraction, hWaitbar, ... + sprintf('Progress: %d/%d', loopcnt, looptatal)); + + + loop_name = ['_fsym_',num2str(fsym)]; + + %%%%% SET Voltages %%%%%% + dcs = DC_supply("active",[1,1],"voltage",[v_bias, 9]); + dcs.set("voltage",[v_bias, 9]); + + %%%%% 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 %%%%%% + SCP = ScopeKeysight("model","DSAZ634A",'autoscale',1,"fadc",'GSa_160',"channel",[1,0],"recordLen",1000000,"removeDC",1); + AWG = AwgKeysight("model","M8196A","fdac",fdac,"scaletodac",[1,1,1,1],"skews",[0,0,0,0],"voltages",[0,0,0,awg_vpp]); + A2S = Awg2Scope(AWG,SCP,[0,0,0,1]); + + + %%%%% Symbol Generation %%%%%% + Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rrcalpha); + + [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",18,"useprbs",1,... + "fs_out",fdac,"applyclipping",0,"clipfactor",1.7,... + "applypulseform",0,"pulseformer",Pform,"randkey",pn_key,... + "db_precode",db_precode,... + "mrds_code",usemrds,"mrds_blocklength",512,"db_encode",db_coding_approach).process(); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 % measure channel + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.buildOFDM(); + elseif precomp_mode == 2 % apply precomp + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.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",10); + + % save([folderpath,[experiment_name,'bits'],loop_name],"Bits"); + % save([folderpath,[experiment_name,'symbols'],loop_name],"Symbols"); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %%%%% AWG --> Scope %%%%%% + [~,~,~,Scpe_sig] = A2S.process("signal4",Digi_sig); + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + % Scpe_sig.spectrum("displayname","Scope PSD","fignum",20); + + % Scpe_sig.plot("displayname","Scope raw signal","fignum",25); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + + %%%%%% Sample to 2x fsym %%%%%% + Scpe_sig = Scpe_sig.resample("fs_in",160e9,"fs_out",2*fsym); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 + freqresp.estimate(Scpe_sig,"save",true,"savePath",precomp_path,"fileName",precomp_fn); + freqresp.plot(); + end + + %%%%%% Sync Rx signal with reference %%%%%% + [Scpe_sig,S] = Scpe_sig.tsynch("reference",Symbols,"fs_ref",fsym); + + %%%%%% SNR CHEAT - Avg. the measured signal occurences %%%%%% + 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; + end + scope_mean = scope_mean ./ n; + Scpe_sig.signal = scope_mean; + end + + %%%%% Plot and Save Routine 2 %%%%%%%%%%%%%%%%%%%%%%%%% + save([folderpath,experiment_name,'rx_signal',loop_name],"S"); + Scpe_sig.eye(fsym,M,"fignum",40,"displayname",' after Scope'); + + voa.readvals(); + rop = voa.power_state(1); + pd_in = voa.power_state(2); + + + %%%%% 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); + + if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ'); + + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber),'| ROP: ',num2str(rop),' dB | PD_in: ',num2str(pd_in),' dBm']); + + if 1 + figure(53); + constellation = unique(Symbols.signal); + received = NaN(numel(constellation),length(Symbols)); + for lvl = 1:numel(constellation) + received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig.signal(Symbols.signal==constellation(lvl)); + hold on + histogram(received(lvl,:),1000,"EdgeAlpha",0); + end + end + + elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + + Noi = EQ_sig-Symbols; + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber_ffe_only,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); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber_ffe_only),' -> PF -> MLSE: ',sprintf('%.1E',ber),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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); + [~,errors_bm,ber,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),' | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_coding_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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,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),' | PD_in: ',num2str(pd_in),' dBm']); + + end + + + wh.addValueToStorage(ber,'ber',precomp_amp_max); + wh.addValueToStorage(rop,'rop',precomp_amp_max); + wh.addValueToStorage(pd_in,'pd_in',precomp_amp_max); + wh.addValueToStorage(Rx_bits,'signals',precomp_amp_max); + wh.addValueToStorage(M,'m',precomp_amp_max); + + showCurrentMeasurement('BER', ber, 'ROP', rop, 'PD in', pd_in, 'PAM',M, 'Vbias', v_bias, 'AWG Vpp', awg_vpp, 'Precomp MaxAmp',precomp_amp_max); + autoArrangeFigures(3,3,2); + + +end + +close(hWaitbar); + +wh.save([folderpath,experiment_name,'_wh']); + +autoArrangeFigures(3,3,2) + +figure(90) +amp_vals = wh.parameter.precomp_amp_max.values; +ber = wh.getStoValue('ber',wh.parameter.precomp_amp_max.values); +plot(amp_vals,ber,'DisplayName',['PAM ',num2str(M)]); +xlabel('Precompensation Max Amp'); +ylabel('BER'); +grid on +grid minor +title('Bit Error Rate vs. ROP'); +set(gca, 'yscale', 'log'); +legend diff --git a/projects/Lab_2024/lab_rop_sweep.m b/projects/Lab_2024/lab_rop_sweep.m new file mode 100644 index 0000000..0e52e1a --- /dev/null +++ b/projects/Lab_2024/lab_rop_sweep.m @@ -0,0 +1,261 @@ + +folderpath = 'C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\no_mpi_2024\'; +experiment_name = '10km_ffe_no_mpi_'; + +ffe_only = 1; +postfilter_approach = 0; +db_channel_approach = 0; +db_coding_approach = 0; + +db_precode = db_coding_approach || db_channel_approach; + +%%% SIR Sweep for MPI Experiment %%% +params = struct; +params.rop_atten = [7:-1:0]; % high atten to low atten to make sure there is no sudden opening of VOA +params.pd_in_set = [6]; % desired P_out (outp. power mode=2) + +wh = DataStorage(params); + +wh.addStorage("ber"); +wh.addStorage("pd_in"); +wh.addStorage("rop"); +wh.addStorage("signals"); + +precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active +precomp_amp_max = 3; + +M = 4; +pn_key = 2; +usemrds = 0; +fdac = 92e9; +fsym = 92e9; +fadc = 160e9; +rrcalpha = 0.05; +v_bias = 2.25; + +disp(['Start Measurement of ',num2str(prod(wh.dim)),' loops...']) + +for rop_atten = wh.parameter.rop_atten.values + for pd_in_set = wh.parameter.pd_in_set.values + + loop_name = ['_ropatten_',num2str(rop_atten),'_pdin_',num2str(pd_in_set)]; + + %%%%% SET Voltages %%%%%% + dcs = DC_supply("active",[1,1],"voltage",[v_bias, 9]); + dcs.set("voltage",[v_bias, 9]); + + 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(); + + SCP = ScopeKeysight("model","DSAZ634A",'autoscale',1,"fadc",'GSa_160',"channel",[1,0],"recordLen",1000000,"removeDC",1); + AWG = AwgKeysight("model","M8196A","fdac",fdac,"scaletodac",[1,1,1,1],"skews",[0,0,0,0],"voltages",[0,0,0,0.62]); + A2S = Awg2Scope(AWG,SCP,[0,0,0,1]); + + if 1 + precomp_path = "C:\Users\sioe\Documents\MATLAB\imdd_simulation\projects\standard_system\"; + precomp_fn = "lab_mpi_setup_2"; + else + precomp_path = "C:\Users\sioe\Documents\MATLAB\model-collection\sioe_models\Labor_2024\Lab_PAM4\"; + precomp_fn = "precomp_bla__loop1_1"; + end + + + %%%%% Symbol Generation %%%%%% + Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rrcalpha); + + [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",18,"useprbs",1,... + "fs_out",fdac,"applyclipping",0,"clipfactor",1.7,... + "applypulseform",0,"pulseformer",Pform,"randkey",pn_key,... + "db_precode",db_precode,... + "mrds_code",usemrds,"mrds_blocklength",512,"db_encode",db_coding_approach).process(); + + if precomp_mode == 1 % measure channel + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.buildOFDM(); + elseif precomp_mode == 2 % apply precomp + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.precomp(Digi_sig,'maxampdb',precomp_amp_max,'loadPath',precomp_path,'fileName',precomp_fn); + end + + Digi_sig = Digi_sig.resample("fs_out",AWG.fdac); + + Digi_sig.spectrum("displayname","Normal Tx","fignum",10); + + + save([folderpath,[experiment_name,'bits'],loop_name],"Bits"); + save([folderpath,[experiment_name,'symbols'],loop_name],"Symbols"); + + %%%%% AWG --> Scope %%%%%% + [~,~,~,Scpe_sig] = A2S.process("signal4",Digi_sig); + + Scpe_sig.spectrum("displayname","Scope PSD","fignum",20); + + + + % Scpe_sig.spectrum("displayname",'Rx Signal','fignum',10); + + %%%%%% Sample to 2x fsym %%%%%% + Scpe_sig = Scpe_sig.resample("fs_in",160e9,"fs_out",2*fsym); + + if precomp_mode == 1 + freqresp.estimate(Scpe_sig,"save",true,"savePath",precomp_path,"fileName",precomp_fn); + freqresp.plot(); + end + + %%%%%% Sync Rx signal with reference %%%%%% + [Scpe_sig,S] = Scpe_sig.tsynch("reference",Symbols,"fs_ref",fsym); + save([folderpath,experiment_name,'rx_signal',loop_name],"S"); + + 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; + end + scope_mean = scope_mean ./ n; + Scpe_sig.signal = scope_mean; + end + + Scpe_sig.plot("displayname","Scope PSD","fignum",30); + Scpe_sig.eye(fsym,M,"fignum",40,"displayname",' after Scope'); + + voa.readvals(); + rop = voa.power_state(1); + pd_in = voa.power_state(2); + + + %%%%% 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.0,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1); + + if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ'); + + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber),'| ROP: ',num2str(rop),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + + Noi = EQ_sig-Symbols; + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber_ffe_only,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); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber_ffe_only),' -> PF -> MLSE: ',sprintf('%.1E',ber),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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); + [~,errors_bm,ber,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),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_coding_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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,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),' | SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + + end + + wh.addValueToStorage(ber,'ber',rop_atten,pd_in_set); + wh.addValueToStorage(rop,'rop',rop_atten,pd_in_set); + wh.addValueToStorage(pd_in,'pd_in',rop_atten,pd_in_set); + wh.addValueToStorage(Rx_bits,'signals',rop_atten,pd_in_set); + + showCurrentMeasurement('BER', ber, 'ROP', rop, 'PD in', pd_in); + autoArrangeFigures(3,3,2); + end +end + +wh.save([folderpath,experiment_name,'_wh']); + +cols = linspecer(8); + +rop_vals = wh.parameter.rop_atten.values; +pd_in_set = wh.parameter.pd_in_set.values(1); + +bers = wh.getStoValue('ber',rop_vals,pd_in_set); +rop_measured = wh.getStoValue('rop',rop_vals,pd_in_set); +pd_in_measured = wh.getStoValue('pd_in',rop_vals,pd_in_set); + +s = wh.getStoValue('signals',rop_vals,pd_in_set); + +figure(90); +hold on; % Retain the plot so new points can be added without complete redraw + +% Plot the data and get the line handle +hLine = plot(rop_measured, bers, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", experiment_name); + +% Store pd_in_measured in the ZData property +hLine.ZData = pd_in_measured; + +% Customize the data tips +% Set labels for existing data tip rows +hLine.DataTipTemplate.DataTipRows(1).Label = 'ROP'; +hLine.DataTipTemplate.DataTipRows(2).Label = 'BER'; +hLine.DataTipTemplate.DataTipRows(2).Format = '%.2e'; % Format BER as "3e-4" + + +% Add a new data tip row for PDin +pdinRow = dataTipTextRow('PDin', 'ZData'); +hLine.DataTipTemplate.DataTipRows(3) = pdinRow; + +% Continue with the rest of your plot settings +yline(3.8e-3, 'DisplayName', 'HD-FEC', 'LineStyle', '--', 'HandleVisibility', 'off'); +xlabel('Received Optical Power (dBm)'); +ylabel('Bit Error Rate (BER)'); +title('Bit Error Rate vs. ROP'); +set(gca, 'yscale', 'log'); +set(gca, 'Box', 'on'); +grid on; +grid minor; +legend('Interpreter', 'none'); + +autoArrangeFigures(3,3,2) + diff --git a/projects/Lab_2024/lab_sir_sweep.m b/projects/Lab_2024/lab_sir_sweep.m new file mode 100644 index 0000000..2bb3efc --- /dev/null +++ b/projects/Lab_2024/lab_sir_sweep.m @@ -0,0 +1,368 @@ + +folderpath = 'C:\Users\sioe\Nextcloud\Dokumente\02_Ablage_Office\Lab_Data_24\sir_sweep_pam4\'; +experiment_name = 'PAM4_DB_encoded_10km_'; + +%%% SIR Sweep for MPI Experiment %%% +params = struct; + +params.vbias = [2.45]; +params.awg_vpp = [0.25]; +params.eq_mode = [4]; +params.i_atten = [0:4:40]; + +wh = DataStorage(params); + +wh.addStorage("ber"); +wh.addStorage("pd_in"); +wh.addStorage("m"); +wh.addStorage("sir"); +wh.addStorage("s_pow"); +wh.addStorage("i_pow"); +wh.addStorage("signals"); + +precomp_path = "C:\Users\sioe\Documents\MATLAB\imdd_simulation\projects\standard_system\"; +precomp_fn = "lab_mpi_setup_2"; +precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active +precomp_amp_max = 5; + +M = 4; +pn_key = 2; +usemrds = 0; +fsym = 92e9; +fdac = 92e9; +awg_vpp = 0.35; +fadc = 160e9; +rrcalpha = 0.05; +v_bias = 2.25; +pd_in_set = 6; + +looptotal = prod(wh.dim); +iterationTimes = zeros(looptotal, 1); % Preallocate for speed + +disp(['Start Measurement of ',num2str(looptotal),' loops...']) +hWaitbar = waitbar(0, 'Starting measurement...', 'Name', 'Processing Progress'); +loopcnt = 0; + +estimatedTimeRemaining = 0; +estimatedTotalTime = 0; +for eq_mode = wh.parameter.eq_mode.values + for i_atten = wh.parameter.i_atten.values + for v_bias = wh.parameter.vbias.values + for awg_vpp = wh.parameter.awg_vpp.values + + iterationStartTime = tic; + + loop_name = ['_iatten_',num2str(i_atten)]; + + loopcnt = loopcnt+1; + progressFraction = loopcnt / looptotal; + waitbar(progressFraction, hWaitbar, ... + sprintf('Progress: %d/%d\nEstimated time remaining: %.2f hours\nEstimated time remaining: %.2f hours', ... + loopcnt, looptotal, estimatedTimeRemaining/60/60, estimatedTotalTime/60/60)); + + try + + z + ber = wh2.getStoValue('ber',v_bias,awg_vpp,eq_mode,precomp_amp_max); + assert(~isempty(ber),'err') + rop = wh2.getStoValue('rop',v_bias,awg_vpp,eq_mode,precomp_amp_max); + assert(~isempty(rop)) + pd_in = wh2.getStoValue('pd_in',v_bias,awg_vpp,eq_mode,precomp_amp_max); + assert(~isempty(pd_in)) + M = wh2.getStoValue('m',v_bias,awg_vpp,eq_mode,precomp_amp_max); + assert(~isempty(M)) + + catch + + + + + switch eq_mode + case 1 + ffe_only = 1; + postfilter_approach = 0; + db_channel_approach = 0; + db_coding_approach = 0; + + db_precode = db_coding_approach || db_channel_approach; + + case 2 + ffe_only = 0; + postfilter_approach = 1; + db_channel_approach = 0; + db_coding_approach = 0; + + db_precode = db_coding_approach || db_channel_approach; + case 3 + ffe_only = 0; + postfilter_approach = 0; + db_channel_approach = 1; + db_coding_approach = 0; + + db_precode = db_coding_approach || db_channel_approach; + case 4 + ffe_only = 0; + postfilter_approach = 0; + db_channel_approach = 0; + db_coding_approach = 1; + + db_precode = db_coding_approach || db_channel_approach; + end + + + + %%%%% SET Voltages %%%%%% + dcs = DC_supply("active",[1,1],"voltage",[v_bias, 9]); + dcs.set("voltage",[v_bias, 9]); + + %%%%% SET Attenuator %%%%%% + voa = OptAtten("active",[1,2,1,1],"value",[0,pd_in_set,0,i_atten],"wavelength",[1310,1310,1310,1310]); + voa.set('active',[1,2,1,1],'value',[0,pd_in_set,0,i_atten]); + % voa.readvals(); + + %%%%% Construct AWG and Scope Modules %%%%%% + SCP = ScopeKeysight("model","DSAZ634A",'autoscale',1,"fadc",'GSa_160',"channel",[1,0],"recordLen",2000000,"removeDC",1); + AWG = AwgKeysight("model","M8196A","fdac",fdac,"scaletodac",[1,1,1,1],"skews",[0,0,0,0],"voltages",[0,0,0,awg_vpp]); + A2S = Awg2Scope(AWG,SCP,[0,0,0,1]); + + %%%%% Symbol Generation %%%%%% + Pform = Pulseformer("fsym",fsym,"fdac",4*fsym,"pulse","rrc","pulselength",16,"rrcalpha",rrcalpha); + + [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",19,"useprbs",1,... + "fs_out",fdac,"applyclipping",0,"clipfactor",1.7,... + "applypulseform",0,"pulseformer",Pform,"randkey",pn_key,... + "db_precode",db_precode,... + "mrds_code",usemrds,"mrds_blocklength",512,"db_encode",db_coding_approach).process(); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 % measure channel + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.buildOFDM(); + elseif precomp_mode == 2 % apply precomp + freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); + Digi_sig = freqresp.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",10); + + if loopcnt == 1 + save([folderpath,[experiment_name,'bits'],loop_name],"Bits"); + save([folderpath,[experiment_name,'symbols'],loop_name],"Symbols"); + end + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + %%%%% AWG --> Scope %%%%%% + [~,~,~,Scpe_sig] = A2S.process("signal4",Digi_sig); + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + Scpe_sig.spectrum("displayname","Scope PSD","fignum",20); + + % Scpe_sig.plot("displayname","Scope raw signal","fignum",25); + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + + + %%%%%% Sample to 2x fsym %%%%%% + Scpe_sig = Scpe_sig.resample("fs_in",160e9,"fs_out",2*fsym); + + %%%%% Precompensation Routine %%%%%% + if precomp_mode == 1 + freqresp.estimate(Scpe_sig,"save",true,"savePath",precomp_path,"fileName",precomp_fn); + freqresp.plot(); + end + + %%%%%% Sync Rx signal with reference %%%%%% + [Scpe_sig,S] = Scpe_sig.tsynch("reference",Symbols,"fs_ref",fsym); + + %%%%%% SNR CHEAT - Avg. the measured signal occurences %%%%%% + 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; + end + scope_mean = scope_mean ./ n; + Scpe_sig.signal = scope_mean; + end + + %%%%% Plot and Save Routine 2 %%%%%%%%%%%%%%%%%%%%%%%%% + save([folderpath,experiment_name,'rx_signal',loop_name],"S"); + % Scpe_sig.eye(fsym,M,"fignum",40,"displayname",' after Scope'); + + voa.readvals(); + + pd_in = voa.power_state(2); + s_pow = voa.power_state(3); + i_pow = voa.power_state(4); + + sir = s_pow- i_pow; + + %%%%% 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.0,"DDmu",[0.0004 0.0004 0.0004 0.0004 ],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1); + + if ffe_only %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ'); + + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber),'| SIR: ',num2str(sir),' dB | PD_in: ',num2str(pd_in),' dBm']); + + if 0 + figure(53); + constellation = unique(Symbols.signal); + received = NaN(numel(constellation),length(Symbols)); + for lvl = 1:numel(constellation) + received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig.signal(Symbols.signal==constellation(lvl)); + hold on + histogram(received(lvl,:),1000,"EdgeAlpha",0); + end + end + + elseif postfilter_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + EQ_sig.plot("fignum",50,"displayname",'After EQ','clear',1); + + Noi = EQ_sig-Symbols; + + Rx_bits = PAMmapper(M,0).demap(EQ_sig); + [~,errors_bm,ber_ffe_only,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 1 + 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 + + if 0 + figure(53); + constellation = unique(Symbols.signal); + received = NaN(numel(constellation),length(Symbols)); + for lvl = 1:numel(constellation) + received(lvl,Symbols.signal==constellation(lvl)) = EQ_sig.signal(Symbols.signal==constellation(lvl)); + hold on + histogram(received(lvl,:),1000,"EdgeAlpha",0); + end + end + + EQ_sig = MLSE("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); + [~,errors_bm,ber,errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); + + disp(['FFE: ',sprintf('%.1E',ber_ffe_only),' -> PF -> MLSE: ',sprintf('%.1E',ber),' dB | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_channel_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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); + [~,errors_bm,ber,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),' | PD_in: ',num2str(pd_in),' dBm']); + + elseif db_coding_approach %%%%%%%%%%%%%%%%%%%%%%%%%%% + + [EQ_sig, Noi] = Eq.process(Scpe_sig,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,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),' | PD_in: ',num2str(pd_in),' dBm']); + + end + + end + + wh.addValueToStorage(ber,'ber',v_bias,awg_vpp,eq_mode,i_atten); + wh.addValueToStorage(pd_in,'pd_in',v_bias,awg_vpp,eq_mode,i_atten); + wh.addValueToStorage(Rx_bits,'signals',v_bias,awg_vpp,eq_mode,i_atten); + wh.addValueToStorage(sir,'sir',v_bias,awg_vpp,eq_mode,i_atten); + wh.addValueToStorage(s_pow,'s_pow',v_bias,awg_vpp,eq_mode,i_atten); + wh.addValueToStorage(i_pow,'i_pow',v_bias,awg_vpp,eq_mode,i_atten); + wh.addValueToStorage(M,'m',v_bias,awg_vpp,eq_mode,i_atten); + + showCurrentMeasurement('BER', ber, 'PD in', pd_in, 'PAM',M, 'Vbias', v_bias, 'AWG Vpp', awg_vpp, 'Precomp MaxAmp',precomp_amp_max); + + iterationTimes(loopcnt) = toc(iterationStartTime); + averageTimePerIteration = mean(iterationTimes(1:loopcnt)); + estimatedTotalTime = averageTimePerIteration * looptotal; + estimatedTimeRemaining = estimatedTotalTime - sum(iterationTimes(1:loopcnt)); + %autoArrangeFigures(3,3,2); + + end + end + end +end + +close(hWaitbar); + +wh.save([folderpath,experiment_name,'wh']); + + +cols = linspecer(8); + +i_atten_vals = wh.parameter.i_atten.values; +v_bias = wh.parameter.vbias.values(1); +awg_vpp = wh.parameter.awg_vpp.values(1); +eq_mode = wh.parameter.eq_mode.values(1); + +bers = wh.getStoValue('ber',v_bias,awg_vpp,eq_mode,i_atten_vals); + +figure(90); +hold on; % Retain the plot so new points can be added without complete redraw + +% Plot the data and get the line handle +hLine = plot(i_atten_vals, bers, "LineWidth", 0.5, "LineStyle", "-", "Marker", ".", "MarkerSize", 15, "DisplayName", experiment_name); + +% Customize the data tips +% Set labels for existing data tip rows +hLine.DataTipTemplate.DataTipRows(1).Label = 'Fsym'; +hLine.DataTipTemplate.DataTipRows(2).Label = 'BER'; +hLine.DataTipTemplate.DataTipRows(2).Format = '%.2e'; % Format BER as "3e-4" + + +% Continue with the rest of your plot settings +yline(3.8e-3, 'DisplayName', 'HD-FEC', 'LineStyle', '--', 'HandleVisibility', 'off'); +xlabel('Signal to Interference Ratio in dB'); +ylabel('Bit Error Rate (BER)'); +title('Bit Error Rate vs. SIR'); +set(gca, 'yscale', 'log'); +set(gca, 'Box', 'on'); +grid on; +grid minor; +legend('Interpreter', 'none'); + + +autoArrangeFigures(3,3,2) + +disp("measurement done") \ No newline at end of file diff --git a/projects/Lab_2024/sweep_laser_vs_power.m b/projects/Lab_2024/sweep_laser_vs_power.m new file mode 100644 index 0000000..60df6e1 --- /dev/null +++ b/projects/Lab_2024/sweep_laser_vs_power.m @@ -0,0 +1,88 @@ + +% 1) Establish connection to laser +o = serialport("COM9",9600); %per USB angeschlossen +configureTerminator(o,"CR") +writeline(o,"*IDN?"); +wait(1) +if o.NumBytesAvailable ~= 0 + disp(['Laser Mainframe: ',readline(o)]); +else + error('Keine Verbindung zum Mainframe mglich?') + clear o + +end + +% 2) connect to device +v = visa('keysight', 'TCPIP0::134.245.243.248::inst0::INSTR'); +fopen(v); +fprintf(v, '*IDN?;'); +disp(['Powermeter: ' fscanf(v)]); + +% Define channels +laser_channel = 7; +powermeter_slot = 1; +l = 1300:0.1:1320; + +% turn on the laser +writeline(o,['CH',num2str(laser_channel),':ENABLE']); +wait(0.2) +current_wavelen = readline(o); + +clear power +clear lambda + +for n = 1:length(l) + + % 2 change wavelength in laser slot + command = string(['CH',num2str(laser_channel),':L=',num2str(l(n))]); + writeline(o,command); + wait(0.5); + readline(o); + + % query and check wavelength + writeline(o,['CH',num2str(laser_channel),':L?']); + wait(0.2) + current_wavelen = readline(o); + current_wavelen = str2double(strrep(regexp(current_wavelen,'([CH7:L=])+([\d]*)+([.])+([\d]*)','match'),'CH7:L=','')); + + if l(n) ~= current_wavelen + clear o + fclose(v); + delete(v); + clear v + error('Wellenlnge wurde nicht bernommen'); + end + + wait(0.75); + + % get current power in slot + slot = 1; + fprintf(v, [':READ' num2str(slot) ':POW?']); + power(n) = sscanf(fscanf(v),'%f'); + lambda(n) = current_wavelen; + + if mod(n,10)==1 + disp(['Measured ',num2str(power(n)),' dBm at ',num2str(lambda(n)),' nm']) + end +end + +figure(2) +hold on +plot(lambda,power,'Marker','*'); +xlabel('Wavelngth in nm') +ylabel('Power in dBm') +grid minor + +figure(211) +hold on +plot(lambda,10.^(power/10),'Marker','*'); +xlabel('Wavelngth in nm') +ylabel('Power in mW') +grid minor + +%close the serial connection + +clear o +fclose(v); +delete(v); +clear v