diff --git a/Classes/00_signals/Signal.m b/Classes/00_signals/Signal.m index dd4cf9e..eb41f39 100644 --- a/Classes/00_signals/Signal.m +++ b/Classes/00_signals/Signal.m @@ -693,6 +693,8 @@ classdef Signal return end + pkpos = sort(pkpos); + if mean(w) > 10 || mean(p) > 10 return else @@ -709,7 +711,7 @@ classdef Signal shifts = lags(pkpos); sequenceStarts = shifts; - % shifts = shifts(shifts>=0); + shifts = shifts(shifts>=0); if numel(shifts) > 0 @@ -720,17 +722,17 @@ classdef Signal for c = shifts sig = obj.delay(-c,'mode','samples'); - sig.signal = sig.signal(1:length(b)) .* -inverted; + sig.signal = sig.signal(1:length(b));% .* -inverted; S{end+1,1} = sig; end - %return/keep the sinal with the highest correlation (only within positive shifts) - [~,idx]=max(pks); - obj.signal = S{idx}.signal; - %put signal with highest corr. to first index in S array - swap = S{1}; - S{1} = S{idx}; - S{idx} = swap; + % %return/keep the sinal with the highest correlation (only within positive shifts) + % [~,idx]=max(pks); + % obj.signal = S{idx}.signal; + % %put signal with highest corr. to first index in S array + % swap = S{1}; + % S{1} = S{idx}; + % S{idx} = swap; for c = 1:numel(shifts) S{c}.logbook = []; @@ -962,7 +964,7 @@ classdef Signal % add information - if 1 + if 0 pwr_dbm = round(obj.power,3); pwr_lin = obj.power("unit",power_notation.W); diff --git a/Classes/DataBaseHandler/DBHandler.m b/Classes/DataBaseHandler/DBHandler.m index 6a77aec..926d016 100644 --- a/Classes/DataBaseHandler/DBHandler.m +++ b/Classes/DataBaseHandler/DBHandler.m @@ -83,7 +83,7 @@ classdef DBHandler < handle result = fetch(obj.conn, 'SELECT name FROM sqlite_master WHERE type="table"'); obj.tableNames = result.name; end - + catch e error('Failed to retrieve table names: %s', e.message); end @@ -169,14 +169,6 @@ classdef DBHandler < handle function healthyDB = dbIsHealthy(obj) healthyDB = false; - num_runs = obj.fetch('SELECT COUNT(*) AS total_runs FROM Runs'); - - num_configs = obj.fetch('SELECT COUNT(*) AS total_configurations FROM Configurations'); - - num_meas = obj.fetch('SELECT COUNT(*) AS total_measurements FROM Measurements'); - - assert((num_runs{1,1}==num_configs{1,1})&&(num_configs{1,1}==num_meas{1,1}),'Different num of entries per table'); - %Check for any duplicate paths duplictae_raw = obj.fetch("SELECT COALESCE(Runs.rx_raw_path,'NaN') AS rx_raw_path, COUNT(*) AS occurrences FROM Runs GROUP BY rx_raw_path HAVING COUNT(*) > 1"); duplictae_sync = obj.fetch("SELECT COALESCE(Runs.rx_sync_path,'NaN') AS rx_sync_path, COUNT(*) AS occurrences FROM Runs GROUP BY rx_sync_path HAVING COUNT(*) > 1"); @@ -190,7 +182,6 @@ classdef DBHandler < handle end - function lastID = appendToTable(obj, tableName, newRow) % appendToTable Appends a new row to the specified table % @@ -206,48 +197,52 @@ classdef DBHandler < handle error('Table %s does not exist in the database or has not been fetched.', tableName); end - % Convert newRow to a table if it is a struct + % Handle struct preprocessing before table conversion if isstruct(newRow) fields = fieldnames(newRow); - emptyFields = structfun(@isempty,newRow); + + % Handle empty fields + emptyFields = structfun(@isempty, newRow); if sum(emptyFields) > 0 - emptyFieldNames = fields(emptyFields); % use () to get a cell array + emptyFieldNames = fields(emptyFields); for idx = 1:numel(emptyFieldNames) newRow.(emptyFieldNames{idx}) = NaN; end end + + % Convert arrays to JSON strings + for i = 1:length(fields) + fieldValue = newRow.(fields{i}); + if isnumeric(fieldValue) && length(fieldValue) > 1 + newRow.(fields{i}) = jsonencode(fieldValue); + end + end + + % Convert to table newRow = struct2table(newRow); end % Ensure the new row matches the structure of the existing table existingTableStructure = obj.tables.(tableName); - % Perform data type checks and conversions + % Perform remaining data type checks and conversions for colName = newRow.Properties.VariableNames - % Extract the value and its intended column type value = newRow.(colName{1}); - existingValue = existingTableStructure.(colName{1}); - % If the value is a class object, convert it to JSON format - if isobject(value) && ~isdatetime(value) && ~isa(value,"string") + if iscell(value) && ~isempty(value) + if ischar(value{1}) || isstring(value{1}) + newRow.(colName{1}) = value{1}; + end + elseif isobject(value) && ~isdatetime(value) && ~isa(value, "string") newRow.(colName{1}) = string(jsonencode(value)); - - % If the value is a character array, convert it to a string elseif ischar(value) newRow.(colName{1}) = string(value); - - elseif isdatetime(value) - - % newRow.(colName{1}) = string(value); - end - end % Parameters for retry logic maxRetries = 50; basePause = 0.05; % seconds - attempt = 0; success = false; @@ -257,36 +252,40 @@ classdef DBHandler < handle newRow_ = obj.convertTableToCellStrings(newRow); end sqlwrite(obj.conn, tableName, newRow); - success = true; % Write successful + success = true; catch e if contains(e.message, 'database is locked') || contains(e.message, 'cannot rollback transaction') attempt = attempt + 1; - pauseTime = basePause * (1 + rand()); % Add random jitter to reduce collision chance - fprintf('Database locked. Retry %d/%d after %.3f seconds...\n', attempt, maxRetries, pauseTime); + pauseTime = basePause * (1 + rand()); + fprintf('Database locked. Retry %d/%d after %.3f seconds...\n', ... + attempt, maxRetries, pauseTime); pause(pauseTime); else + fprintf('Error details:\n'); + fprintf('Column values:\n'); + disp(newRow); error('Failed to append to the table %s: %s', tableName, e.message); end end end if ~success - error('Failed to append to table %s after %d retries due to database lock.', tableName, maxRetries); + error('Failed to append to table %s after %d retries due to database lock.', ... + tableName, maxRetries); end - % Retrieve the measurement_id of the newly inserted row for linking other tables - + % Retrieve the ID of the newly inserted row if obj.type == "mysql" - queue = "SELECT LAST_INSERT_ID()"; + query = "SELECT LAST_INSERT_ID()"; elseif obj.type == "sqlite" - queue = "SELECT last_insert_rowid()"; + query = "SELECT last_insert_rowid()"; end - result = fetch(obj.conn, queue); - lastID = result{1, 1}; % Access the value directly from the table + result = fetch(obj.conn, query); + lastID = result{1, 1}; end - function exists = checkIfRunExists(obj, table2check, column2check, value2check) + function [exists, count] = checkIfRunExists(obj, table2check, column2check, value2check) % checkIfRunExists Checks if a specific value exists in a specified column of a table % % Usage: @@ -311,7 +310,11 @@ classdef DBHandler < handle end % Construct the query to check for the value in the specified column - query = sprintf('SELECT COUNT(*) FROM %s WHERE %s = "%s"', table2check, column2check, value2check); + if isnumeric(value2check) + query = sprintf('SELECT COUNT(*) FROM %s WHERE %s = %d', table2check, column2check, value2check); + else + query = sprintf('SELECT COUNT(*) FROM %s WHERE %s = "%s"', table2check, column2check, value2check); + end % Execute the query and pass the value2check to avoid SQL injection issues try @@ -325,12 +328,21 @@ classdef DBHandler < handle exists = count > 0; if exists - disp(['The value "', value2check, '" already exists in the column "', column2check, '" of the table "', table2check, '".']); + % disp(['The value "', num2str(value2check), '" already exists in the column "', column2check, '" of the table "', table2check, '".']); else % disp(['The value "', value2check, '" does not exist in the column "', column2check, '" of the table "', table2check, '".']); end end + function hashStr = calcHash(~,object) + jsonStr = jsonencode(object); + md = java.security.MessageDigest.getInstance('MD5'); + md.update(uint8(jsonStr)); + hashBytes = typecast(md.digest, 'uint8'); + hashStr = lower(dec2hex(hashBytes)'); + hashStr = lower(strtrim(hashStr(:)')); % Convert to a lowercase string + end + function resultID = addProcessingResult(obj, run_id, resultData, eqParamsData) % addProcessingResult Adds a processing result and links it to an EqualizerParameters entry. @@ -348,18 +360,19 @@ classdef DBHandler < handle % resultID: The result_id of the newly inserted ProcessingResults entry. % 1. Compute hash for equalizer parameters - jsonStr = jsonencode(eqParamsData); - md = java.security.MessageDigest.getInstance('MD5'); - md.update(uint8(jsonStr)); - hashBytes = typecast(md.digest, 'uint8'); - hashStr = lower(dec2hex(hashBytes)'); - hashStr = lower(strtrim(hashStr(:)')); % Convert to a lowercase string + % jsonStr = jsonencode(eqParamsData); + % md = java.security.MessageDigest.getInstance('MD5'); + % md.update(uint8(jsonStr)); + % hashBytes = typecast(md.digest, 'uint8'); + % hashStr = lower(dec2hex(hashBytes)'); + % hashStr = lower(strtrim(hashStr(:)')); % Convert to a lowercase string + hashStr = obj.calcHash(eqParamsData); % Add hash to equalizer parameters - eqParamsData.config_hash = hashStr; + eqParamsData.hash = hashStr; % 2. Check if an equalizer configuration with the same hash exists - queryStr = sprintf('SELECT eq_id FROM EqualizerParameters WHERE config_hash = ''%s''', eqParamsData.config_hash); + queryStr = sprintf('SELECT eq_id FROM Equalizer WHERE hash = ''%s''', eqParamsData.hash); existingEntry = obj.fetch(queryStr); if ~isempty(existingEntry) @@ -367,27 +380,30 @@ classdef DBHandler < handle eq_id = existingEntry{1,1}; else % Insert the new equalizer configuration and get its eq_id - eq_id = obj.appendToTable('EqualizerParameters', eqParamsData); + eqParamsData = eqParamsData.toStruct; + eq_id = obj.appendToTable('Equalizer', eqParamsData); end % 3. Add the equalizer configuration reference and run_id to resultData - resultData.eqParam_id = eq_id; + resultData = resultData.toStruct; + resultData.eq_id = eq_id; resultData.run_id = run_id; % 4. Compute hash for the processing result tempResultData = rmfield(resultData, 'date_of_processing'); - resultJsonStr = jsonencode(tempResultData); - md2 = java.security.MessageDigest.getInstance('MD5'); % Create a new MD5 instance - md2.update(uint8(resultJsonStr)); - resultHashBytes = typecast(md2.digest, 'uint8'); - resultHashStr = lower(dec2hex(resultHashBytes)'); - resultHashStr = lower(strtrim(resultHashStr(:)')); % Convert to a lowercase string + % resultJsonStr = jsonencode(tempResultData); + % md2 = java.security.MessageDigest.getInstance('MD5'); % Create a new MD5 instance + % md2.update(uint8(resultJsonStr)); + % resultHashBytes = typecast(md2.digest, 'uint8'); + % resultHashStr = lower(dec2hex(resultHashBytes)'); + % resultHashStr = lower(strtrim(resultHashStr(:)')); % Convert to a lowercase string + resultHashStr = obj.calcHash(tempResultData); % Add the result hash to resultData - resultData.result_hash = resultHashStr; + resultData.hash = resultHashStr; % 5. Check if an identical processing result already exists - queryStr2 = sprintf('SELECT result_id FROM Results WHERE result_hash = ''%s''', resultData.result_hash); + queryStr2 = sprintf('SELECT result_id FROM Results WHERE hash = ''%s''', resultData.hash); existingResult = obj.fetch(queryStr2); if ~isempty(existingResult) @@ -397,7 +413,66 @@ classdef DBHandler < handle return; end - % 6. Insert the processing result + % 6. check if obj.tables.Results matches Metricstruct + % Fields to exclude from comparison + excludeFields = {'result_id', 'run_id', 'eq_id'}; + + % 7. Chack for new fields in Metric struct and append to + % database if necessary + ms=Metricstruct; + metricFields = setdiff(fieldnames(ms), excludeFields); + tableFields = setdiff(fieldnames(obj.tables.Results), excludeFields); + + % Check matches and find missing fields + matches = all(ismember(metricFields, tableFields)); + + if ~matches + missingFields = setdiff(metricFields, tableFields); + + + % If there are missing fields, add them to the SQL table + if ~isempty(missingFields) + for i = 1:length(missingFields) + fieldName = missingFields{i}; + + % Determine SQL data type based on MATLAB class + fieldValue = ms.(fieldName); + if isnumeric(fieldValue) + if isinteger(fieldValue) + sqlType = 'INTEGER'; + else + sqlType = 'REAL'; + end + elseif ischar(fieldValue) || isstring(fieldValue) + sqlType = 'TEXT'; + elseif isdatetime(fieldValue) + sqlType = 'DATETIME'; + elseif iscell(fieldValue) || isstruct(fieldValue) || islogical(fieldValue) + sqlType = 'TEXT'; % Store as JSON + else + sqlType = 'TEXT'; % Default to TEXT for unknown types + end + + % Create ALTER TABLE query + queryStr = sprintf('ALTER TABLE Results ADD COLUMN %s %s', fieldName, sqlType); + + try + % Execute the query + obj.fetch(queryStr); + fprintf('Added field "%s" of type %s to Results table\n', fieldName, sqlType); + catch ME + fprintf('Error adding field "%s": %s\n', fieldName, ME.message); + end + end + else + fprintf('No missing fields to add.\n'); + end + + end + + + + % 8. Insert the processing result resultID = obj.appendToTable('Results', resultData); end @@ -486,8 +561,6 @@ classdef DBHandler < handle end - - function executeSQL(obj, query) % This method executes an SQL statement using MATLAB's execute function. execute(obj.conn, query); @@ -518,22 +591,6 @@ classdef DBHandler < handle selectedFields = []; end - % Step 1: Prompt the user to input filter parameters if not provided - if isempty(filterParams) - filterParams = obj.promptFilterParameters(); - end - - % Step 2: Prompt the user to select fields to include in the SELECT statement - if isempty(selectedFields) - selectedFields = obj.promptSelectFields(); - else - if iscell(selectedFields) - - elseif isstruct(selectedFields) - - end - end - % Step 3: Construct the SQL query based on the inputs query = obj.constructSQLQuery(filterParams, selectedFields); @@ -619,16 +676,20 @@ classdef DBHandler < handle end end - function query = constructSQLQuery(obj, filterParams, selectedFields) % constructSQLQuery Constructs the SQL query based on filter parameters and selected fields. - % - % If selectedFields is provided as a struct, it is converted to a cell array. - % The conversion takes the field names and creates entries like: - % {'selectedFields.fieldName'} for each field. - % Input check for selectedFields: if it's a struct, convert it to a cell array. - if isstruct(selectedFields) + arguments + obj + filterParams + selectedFields + end + + % Step 1: Handle selectedFields conversion + if isempty(selectedFields) || (ischar(selectedFields) && strcmpi(selectedFields, 'all')) + selectedFields = obj.getTableFieldNames('Runs'); % Default to Runs table + elseif isstruct(selectedFields) + % Convert struct to cell array of 'Table.field' format newFields = {}; tableNames = fieldnames(selectedFields); for t = 1:numel(tableNames) @@ -643,39 +704,10 @@ classdef DBHandler < handle selectedFields = newFields; end - % Construct the SELECT clause dynamically based on user selection. - % (Assuming that when provided as a cell array, each entry is of the form - % 'TableName.fieldName' or, in our conversion case, 'selectedFields.fieldName'.) - selectClause = 'SELECT DISTINCT '; - for i = 1:numel(selectedFields) - fieldParts = strsplit(selectedFields{i}, '.'); - % If the field comes from the struct conversion, its first part is 'selectedFields' - % and the actual field name is in the second part. - if strcmp(fieldParts{1}, 'selectedFields') - tableName = fieldParts{1}; % not used for type checking below - fieldName = fieldParts{2}; - else - tableName = fieldParts{1}; - fieldName = fieldParts{2}; - end + % Step 2: Generate COALESCE string for SELECT clause + selectClause = obj.generateCoalesceString(selectedFields); - % Decide on COALESCE depending on the field type. - % If the table is known in obj.tables and the field is numeric, use 'NaN'. - if isfield(obj.tables, tableName) && isfield(obj.tables.(tableName), fieldName) && isnumeric(obj.tables.(tableName).(fieldName)) - selectClause = [selectClause, 'COALESCE(', selectedFields{i}, ', ''NaN'') AS ', fieldName]; - else - selectClause = [selectClause, 'COALESCE(', selectedFields{i}, ', '''') AS ', fieldName]; - end - - if i < numel(selectedFields) - selectClause = [selectClause, ', ']; - else - selectClause = [selectClause, ' ']; - end - end - - - % --- Adaptive FROM Clause --- + % Step 3: Generate FROM clause with joins mainTable = 'Runs'; fromClause = ['FROM ', mainTable, ' ']; @@ -693,58 +725,90 @@ classdef DBHandler < handle if isfield(obj.tables.(tableName), 'run_id') % Direct join to Runs - normalJoins = [normalJoins, 'LEFT JOIN ', tableName, ' ON ', mainTable, '.run_id = ', tableName, '.run_id ']; + normalJoins = [normalJoins, 'LEFT JOIN ', tableName, ... + ' ON ', mainTable, '.run_id = ', tableName, '.run_id ']; elseif isfield(obj.tables.(tableName), 'eq_id') % Equalizer depends on Results, collect this join separately - equalizerJoin = ['LEFT JOIN ', tableName, ' ON Results.eqParam_id = ', tableName, '.eq_id ']; + equalizerJoin = ['LEFT JOIN ', tableName, ... + ' ON Results.eq_id = ', tableName, '.eq_id ']; end end % Combine joins: normal joins first, Equalizer last fromClause = [fromClause, normalJoins, equalizerJoin]; + % Step 4: Generate WHERE clause if filter parameters exist + if ~isempty(filterParams) + whereClause = obj.generateWhereClause(filterParams); + if ~isempty(whereClause) + query = [selectClause, ' ', fromClause, 'WHERE ', whereClause]; + else + query = [selectClause, ' ', fromClause]; + end + else + query = [selectClause, ' ', fromClause]; + end + end - - - % --- WHERE Clause Construction --- - baseQuery = [selectClause, ' ', fromClause, 'WHERE ']; + function whereClause = generateWhereClause(obj, filterParams) filterClauses = []; - tableNames_ = fieldnames(filterParams); - for t = 1:numel(tableNames_) - tableName = tableNames_{t}; + + % If filterParams is a DbFilterParameter object, get its internal structure + if isa(filterParams, 'QueryFilter') + filterParams = filterParams.toStruct(); + end + + % Now proceed with the structure + tableNames = fieldnames(filterParams); + + for t = 1:numel(tableNames) + tableName = tableNames{t}; tableParams = filterParams.(tableName); fieldNames = fieldnames(tableParams); + for i = 1:numel(fieldNames) fieldName = fieldNames{i}; value = tableParams.(fieldName); fullName = sprintf('%s.%s', tableName, fieldName); - % Handle various types of values for SQL query construction - if isempty(value) + % Skip empty values + if isempty(value) || (isa(value, 'QueryFilter') && isempty(value.value)) continue; - elseif isnumeric(value) && isnan(value) - filterClause = sprintf('%s IS NULL', fullName); - elseif isnumeric(value) && ~isEnumeration(value) - filterClause = sprintf('%s = %f', fullName, value); - elseif islogical(value) || (isnumeric(value) && ismember(value, [0, 1])) && ~isEnumeration(value) - filterClause = sprintf('%s = %d', fullName, value); - elseif ischar(value) || isstring(value) - filterClause = sprintf('%s = ''%s''', fullName, char(value)); - elseif isEnumeration(value) - filterClause = sprintf('%s = ''%s''', fullName, value); - else - error('Unsupported data type for field "%s".', fullName); end - filterClauses = [filterClauses, filterClause, ' AND ']; + + % Handle Filter class + if isa(value, 'SqlFilter') + if isnumeric(value.value) + filterClause = sprintf('%s %s %f', ... + fullName, value.operator, value.value); + elseif ischar(value.value) || isstring(value.value) + filterClause = sprintf('%s %s ''%s''', ... + fullName, value.operator, char(value.value)); + else + continue; % Skip unsupported types + end + filterClauses = [filterClauses, filterClause, ' AND ']; + else + % Handle direct values (legacy support) + if isnumeric(value) && isnan(value) + filterClause = sprintf('%s IS NULL', fullName); + elseif isnumeric(value) + filterClause = sprintf('%s = %f', fullName, value); + elseif ischar(value) || isstring(value) + filterClause = sprintf('%s = ''%s''', fullName, char(value)); + else + continue; % Skip unsupported types + end + filterClauses = [filterClauses, filterClause, ' AND ']; + end end end - % Remove trailing ' AND ' if any filters were added. + % Remove trailing ' AND ' if any filters were added if ~isempty(filterClauses) - filterClauses = filterClauses(1:end-5); - query = [selectClause, ' ', fromClause, 'WHERE ', filterClauses]; + whereClause = filterClauses(1:end-5); else - query = [selectClause, ' ', fromClause]; + whereClause = ''; end end @@ -980,5 +1044,59 @@ classdef DBHandler < handle end + function fieldNames = getTableFieldNames(obj, tableName) + % Returns all field names for a given table as a cell array in the format {'TableName.fieldName'} + if isfield(obj.tables, tableName) + % Get raw field names + rawFields = fieldnames(obj.tables.(tableName)); + + % Create cell array with table name prefix + fieldNames = cellfun(@(x) [tableName, '.', x], ... + rawFields, ... + 'UniformOutput', false); + else + error('Table "%s" not found in obj.tables.', tableName); + end + end + + + function coalesceStr = generateCoalesceString(obj, selectedFields) + % Generates a COALESCE string for selected fields + % Input: + % selectedFields: cell array of strings in format {'Table.field'} + % e.g., {'Runs.run_id', 'Runs.bitrate'} + % Output: + % coalesceStr: string with COALESCE statements + + arguments + obj + selectedFields cell + end + + % Initialize cell array to store each COALESCE statement + coalesceStatements = cell(length(selectedFields), 1); + + % Generate COALESCE statement for each field + for i = 1:length(selectedFields) + % Split table and field name + parts = strsplit(selectedFields{i}, '.'); + if length(parts) ~= 2 + error('Field name must be in format "Table.field": %s', selectedFields{i}); + end + tableName = parts{1}; + fieldName = parts{2}; + + % Generate COALESCE statement + coalesceStatements{i} = sprintf('COALESCE(%s.%s, ''NaN'') AS %s ', ... + tableName, fieldName, fieldName); + end + + % Join with comma, newline and MATLAB string continuation + coalesceStr = strjoin(coalesceStatements, [', ' sprintf('\n ')]); + % Add initial newline and spacing for formatting + coalesceStr = [sprintf('SELECT DISTINCT \n ') coalesceStr]; + end + + end end diff --git a/Classes/DataBaseHandler/Equalizerstruct.m b/Classes/DataBaseHandler/Equalizerstruct.m new file mode 100644 index 0000000..3a72307 --- /dev/null +++ b/Classes/DataBaseHandler/Equalizerstruct.m @@ -0,0 +1,68 @@ +classdef Equalizerstruct + % Equalizerstruct - Class to store and manage equalizer structure data + + properties + eq_id (1,1) double {mustBeNumeric} = NaN + equalizer_structure equalizer_structure = equalizer_structure.ffe + eq + mlse + comment char = string.empty() % Changed to string with proper empty initialization + hash char = string.empty() % Changed to string with proper empty initialization + end + + methods + function obj = Equalizerstruct(varargin) + % Constructor method for Equalizerstruct + % Can be called empty or with name-value pairs + + if nargin > 0 + for i = 1:2:nargin + if isprop(obj, varargin{i}) + obj.(varargin{i}) = varargin{i+1}; + else + error('Property %s does not exist in Equalizerstruct', varargin{i}); + end + end + end + end + + function s = toStruct(obj) + % Convert the object to a struct + s = struct(); + props = properties(obj); + for i = 1:length(props) + s.(props{i}) = obj.(props{i}); + end + end + + function str = toString(obj) + % Convert the object to a formatted string + s = obj.toStruct(); + str = sprintf('Equalizerstruct:\n'); + fields = fieldnames(s); + for i = 1:length(fields) + val = s.(fields{i}); + if isempty(val) + str = sprintf('%s%s: []\n', str, fields{i}); + elseif isnumeric(val) && length(val) > 1 + str = sprintf('%s%s: [%s]\n', str, fields{i}, num2str(val')); + else + str = sprintf('%s%s: %s\n', str, fields{i}, string(val)); + end + end + end + end + + methods (Static) + function obj = fromStruct(s) + % Create an Equalizerstruct object from a struct + obj = Equalizerstruct(); + fields = fieldnames(s); + for i = 1:length(fields) + if isprop(obj, fields{i}) + obj.(fields{i}) = s.(fields{i}); + end + end + end + end +end \ No newline at end of file diff --git a/Classes/DataBaseHandler/Metricstruct.m b/Classes/DataBaseHandler/Metricstruct.m new file mode 100644 index 0000000..9161e47 --- /dev/null +++ b/Classes/DataBaseHandler/Metricstruct.m @@ -0,0 +1,87 @@ +classdef Metricstruct + % ResultData - Class to store and manage metric data from signal processing results + + properties + result_id (1,1) double {mustBeNumeric} = NaN + run_id (1,1) double {mustBeNumeric} = NaN + eqParam_id (1,1) double {mustBeNumeric} = NaN + date_of_processing (1,1) datetime = datetime('now') + + numBits (1,1) double {mustBeInteger, mustBeNonnegative} = 0 + BER (1,1) double {mustBeNumeric, mustBeNonnegative, mustBeLessThanOrEqual(BER,1)} = 0 + numBitErr (1,1) double {mustBeInteger, mustBeNonnegative} = 0 + BER_precoded (1,1) double {mustBeNumeric, mustBeNonnegative, mustBeLessThanOrEqual(BER_precoded,1)} = 0 + numBitErr_precoded (1,1) double {mustBeInteger, mustBeNonnegative} = 0 + + SNR (1,1) double {mustBeNumeric} = NaN + SNR_level (:,1) double {mustBeNumeric} = [] + STD (1,1) double {mustBeNumeric} = NaN + STD_level (:,1) double {mustBeNumeric, mustBeNonnegative} = [] + STDrx (1,1) double {mustBeNumeric} = NaN + STDrx_level (:,1) double {mustBeNumeric, mustBeNonnegative} = [] + EVM (1,1) double {mustBeNumeric} = NaN + EVM_level (:,1) double {mustBeNumeric} = [] + + GMI (1,1) double {mustBeNumeric} = NaN + AIR (1,1) double {mustBeNumeric} = NaN + Alpha (:,1) double {mustBeNumeric} = [] + MLSE_dir (:,1) double {mustBeNumeric} = [] + end + + methods + function obj = Metricstruct(varargin) + % Constructor method for ResultData + % Can be called empty or with name-value pairs + + % Process name-value pairs if provided + if nargin > 0 + for i = 1:2:nargin + if isprop(obj, varargin{i}) + obj.(varargin{i}) = varargin{i+1}; + else + error('Property %s does not exist in ResultData', varargin{i}); + end + end + end + end + + function s = toStruct(obj) + % Convert the object to a struct + s = struct(); + props = properties(obj); + for i = 1:length(props) + s.(props{i}) = obj.(props{i}); + end + end + + function str = toString(obj) + % Convert the object to a formatted string + s = obj.toStruct(); + str = sprintf('ResultData:\n'); + fields = fieldnames(s); + for i = 1:length(fields) + val = s.(fields{i}); + if isempty(val) + str = sprintf('%s%s: []\n', str, fields{i}); + elseif isnumeric(val) && length(val) > 1 + str = sprintf('%s%s: [%s]\n', str, fields{i}, num2str(val')); + else + str = sprintf('%s%s: %s\n', str, fields{i}, string(val)); + end + end + end + end + + methods (Static) + function obj = fromStruct(s) + % Create a ResultData object from a struct + obj = Metricstruct(); + fields = fieldnames(s); + for i = 1:length(fields) + if isprop(obj, fields{i}) + obj.(fields{i}) = s.(fields{i}); + end + end + end + end +end \ No newline at end of file diff --git a/Classes/DataBaseHandler/QueryFilter.m b/Classes/DataBaseHandler/QueryFilter.m new file mode 100644 index 0000000..86dfe01 --- /dev/null +++ b/Classes/DataBaseHandler/QueryFilter.m @@ -0,0 +1,169 @@ +% QueryFilter - Class for building SQL WHERE conditions for database queries. +% +% Usage: +% qf = QueryFilter(); +% qf.where('Runs', 'fiber_length', SqlOperator.GREATER_THAN, 5); +% qf.where('Runs', 'bitrate', SqlOperator.IN, [224, 336, 448]); +% qf.where('Runs', 'is_mpi', SqlOperator.EQUALS, 0); +% +% % Convert to struct for use in DBHandler or other query functions: +% filterStruct = qf.toStruct(); +% +% % Display current filters: +% disp(qf); + +classdef QueryFilter < handle + + + properties (Access = private) + filters struct = struct() + end + + methods + function obj = QueryFilter(oldFormat) + % Constructor - optionally convert from old filter format + if nargin > 0 && isstruct(oldFormat) + tableNames = fieldnames(oldFormat); + for t = 1:length(tableNames) + tableName = tableNames{t}; + if isstruct(oldFormat.(tableName)) + fields = fieldnames(oldFormat.(tableName)); + for f = 1:length(fields) + fieldName = fields{f}; + value = oldFormat.(tableName).(fieldName); + if ~isempty(value) + obj.where(tableName, fieldName, value); + end + end + end + end + end + end + + function where(obj, tableName, fieldName, operator, value) + % Add a WHERE condition to the query + % + % Inputs: + % tableName - Name of the database table (char) + % fieldName - Name of the field/column (char) + % operator - SqlOperator enum or value if using EQUALS + % value - Value to filter by + % + % Example: + % filter.where('Runs', 'fiber_length', SqlOperator.GREATER_THAN, 5) + arguments + obj + tableName char + fieldName char + operator SqlOperator + value = [] + end + + % Handle case where operator is the value (using default EQUALS) + if nargin < 5 + value = operator; + operator = SqlOperator.EQUALS; + end + + % Validate operator type + if ~isa(operator, 'SqlOperator') + error('Operator must be a SqlOperator enumeration'); + end + + % Validate value based on operator + obj.validateValue(operator, value); + + % Create filter + filter = SqlFilter(value, operator.toSqlString()); + + % Add to filters + if ~isfield(obj.filters, tableName) + obj.filters.(tableName) = struct(); + end + obj.filters.(tableName).(fieldName) = filter; + end + + function s = toStruct(obj) + % Convert filters to struct format for database query + s = obj.filters; + end + + function display(obj) + % Custom display of filter conditions + fprintf('QueryFilter with conditions:\n'); + if isempty(fieldnames(obj.filters)) + fprintf(' No filters set\n'); + return; + end + + tables = fieldnames(obj.filters); + for t = 1:length(tables) + tableName = tables{t}; + if ~isempty(fieldnames(obj.filters.(tableName))) + fprintf('\nTable: %s\n', tableName); + fields = fieldnames(obj.filters.(tableName)); + for f = 1:length(fields) + fieldName = fields{f}; + filter = obj.filters.(tableName).(fieldName); + if isnumeric(filter.value) + if length(filter.value) > 1 + valueStr = ['[', num2str(filter.value), ']']; + else + valueStr = num2str(filter.value); + end + elseif isempty(filter.value) + valueStr = 'empty'; + else + valueStr = char(filter.value); + end + fprintf(' %s %s %s\n', fieldName, filter.operator, valueStr); + end + end + end + end + + function clear(obj, tableName) + % Clear all filters or filters for a specific table + if nargin < 2 + obj.filters = struct(); + else + if isfield(obj.filters, tableName) + obj.filters = rmfield(obj.filters, tableName); + end + end + end + + function remove(obj, tableName, fieldName) + % Remove a specific filter + if isfield(obj.filters, tableName) && ... + isfield(obj.filters.(tableName), fieldName) + obj.filters.(tableName) = rmfield(obj.filters.(tableName), fieldName); + end + end + end + + methods (Access = private) + function validateValue(~, operator, value) + % Validate value based on operator type + switch operator + case SqlOperator.BETWEEN + if ~isnumeric(value) || length(value) ~= 2 + error('BETWEEN operator requires array of 2 numbers'); + end + case SqlOperator.IN + if ~isnumeric(value) || isempty(value) + error('IN operator requires non-empty array'); + end + case SqlOperator.LIKE + if ~ischar(value) && ~isstring(value) + error('LIKE operator requires string value'); + end + otherwise + % For other operators, just ensure value is not empty + if isempty(value) + error('Value cannot be empty'); + end + end + end + end +end \ No newline at end of file diff --git a/Classes/DataBaseHandler/SqlFilter.m b/Classes/DataBaseHandler/SqlFilter.m new file mode 100644 index 0000000..3dbc07f --- /dev/null +++ b/Classes/DataBaseHandler/SqlFilter.m @@ -0,0 +1,22 @@ +% SqlFilter - Internal class for holding a value and SQL operator for a filter condition. +% +% Usage: +% f = SqlFilter(10, SqlOperator.GREATER_THAN.toSqlString()); +% % f.value == 10, f.operator == '>' + + +classdef SqlFilter + properties + value + operator char = '=' + end + + methods + function obj = SqlFilter(value, operator) + obj.value = value; + if nargin > 1 + obj.operator = operator; + end + end + end +end \ No newline at end of file diff --git a/Classes/DataBaseHandler/SqlOperator.m b/Classes/DataBaseHandler/SqlOperator.m new file mode 100644 index 0000000..4e4beea --- /dev/null +++ b/Classes/DataBaseHandler/SqlOperator.m @@ -0,0 +1,48 @@ +% SqlOperator - Enumeration of supported SQL operators for query building. +% +% Usage: +% op = SqlOperator.GREATER_THAN; +% opStr = op.toSqlString(); % returns '>' +% +% Supported operators: +% EQUALS, GREATER_THAN, LESS_THAN, GREATER_EQUAL, LESS_EQUAL, +% NOT_EQUAL, IN, BETWEEN, LIKE + +classdef SqlOperator < uint32 + enumeration + EQUALS (1) % == + GREATER_THAN (2) % > + LESS_THAN (3) % < + GREATER_EQUAL (4) % >= + LESS_EQUAL (5) % <= + NOT_EQUAL (6) % != + IN (7) % IN + BETWEEN (8) % BETWEEN + LIKE (9) % LIKE + end + + methods + function op = toSqlString(obj) + switch obj + case SqlOperator.EQUALS + op = '='; + case SqlOperator.GREATER_THAN + op = '>'; + case SqlOperator.LESS_THAN + op = '<'; + case SqlOperator.GREATER_EQUAL + op = '>='; + case SqlOperator.LESS_EQUAL + op = '<='; + case SqlOperator.NOT_EQUAL + op = '!='; + case SqlOperator.IN + op = 'IN'; + case SqlOperator.BETWEEN + op = 'BETWEEN'; + case SqlOperator.LIKE + op = 'LIKE'; + end + end + end +end \ No newline at end of file diff --git a/Classes/DataBaseHandler/copyStylingFrom.m b/Classes/DataBaseHandler/copyStylingFrom.m deleted file mode 100644 index aaf110f..0000000 --- a/Classes/DataBaseHandler/copyStylingFrom.m +++ /dev/null @@ -1,76 +0,0 @@ -function copyStylingFrom(figNumSource, figNumTgt) - % Get handles to the source and target figures - sourceFig = figure(figNumSource); - targetFig = figure(figNumTgt); - - % Get axes of source and target figures - sourceAxes = findall(sourceFig, 'type', 'axes'); - targetAxes = findall(targetFig, 'type', 'axes'); - - % Ensure the number of axes match - if length(sourceAxes) ~= length(targetAxes) - error('Number of axes in source and target figures must be the same.'); - end - - % Loop through each pair of axes and copy styling properties - for i = 1:length(sourceAxes) - copyAxesProperties(sourceAxes(i), targetAxes(i)); - end - - % Apply general figure properties if desired - targetFig.Color = sourceFig.Color; % Background color -end - -function copyAxesProperties(sourceAx, targetAx) - % List of properties to copy from source to target axes - propsToCopy = {'XColor', 'YColor', 'ZColor', 'FontSize', 'FontName', ... - 'GridColor', 'GridLineStyle', 'MinorGridColor', 'Box', ... - 'XGrid', 'YGrid', 'ZGrid', 'XMinorGrid', 'YMinorGrid', 'ZMinorGrid', ... - 'LineWidth', 'TitleFontSizeMultiplier', 'LabelFontSizeMultiplier'}; - - % Copy properties from source to target - for i = 1:length(propsToCopy) - try - targetAx.(propsToCopy{i}) = sourceAx.(propsToCopy{i}); - catch - % Skip property if it doesn't exist or can't be copied - end - end - - % Copy axis labels and titles - targetAx.Title.String = sourceAx.Title.String; - targetAx.XLabel.String = sourceAx.XLabel.String; - targetAx.YLabel.String = sourceAx.YLabel.String; - targetAx.ZLabel.String = sourceAx.ZLabel.String; - - % Copy children elements like lines, patches, etc. - sourceChildren = allchild(sourceAx); - targetChildren = allchild(targetAx); - - % Ensure the number of children elements match - if length(sourceChildren) ~= length(targetChildren) - warning('Number of elements in source and target axes differ. Styling may not be applied completely.'); - end - - % Copy properties of children (like lines, patches, etc.), except colors and legends - for i = 1:min(length(sourceChildren), length(targetChildren)) - copyObjectProperties(sourceChildren(i), targetChildren(i)); - end -end - -function copyObjectProperties(sourceObj, targetObj) - % List of common properties to copy for plot elements (lines, patches, etc.) - propsToCopy = {'LineStyle', 'LineWidth', 'Marker', 'MarkerSize', ... - 'MarkerEdgeColor', 'MarkerFaceColor', 'DisplayName'}; - - % Copy properties from source to target, excluding colors - for i = 1:length(propsToCopy) - try - if ~contains(propsToCopy{i}, 'Color') % Skip color properties - targetObj.(propsToCopy{i}) = sourceObj.(propsToCopy{i}); - end - catch - % Skip property if it doesn't exist or can't be copied - end - end -end diff --git a/Classes/DataBaseHandler/related functions/cleanUpTable.m b/Classes/DataBaseHandler/related functions/cleanUpTable.m new file mode 100644 index 0000000..239ca30 --- /dev/null +++ b/Classes/DataBaseHandler/related functions/cleanUpTable.m @@ -0,0 +1,32 @@ +function cleanedTable = cleanUpTable(inputTable) +% Converts string numbers to numeric, 'NaN' to NaN, and tries to convert date strings to datetime. + +cleanedTable = inputTable; +varNames = cleanedTable.Properties.VariableNames; + +for i = 1:numel(varNames) + col = cleanedTable.(varNames{i}); + if iscell(col) + numericCol = str2double(col); + if all(isnan(numericCol) == strcmpi(col, 'NaN') | cellfun(@isempty, col)) + cleanedTable.(varNames{i}) = numericCol; + else + try + cleanedTable.(varNames{i}) = datetime(col, 'InputFormat', 'yyyy-MM-dd HH:mm:ss.SSSSSS'); + catch + cleanedTable.(varNames{i}) = string(col); + end + end + elseif isstring(col) + numericCol = str2double(col); + if all(isnan(numericCol) == strcmpi(col, "NaN")) + cleanedTable.(varNames{i}) = numericCol; + else + try + cleanedTable.(varNames{i}) = datetime(col, 'InputFormat', 'yyyy-MM-dd HH:mm:ss'); + catch + end + end + end +end +end \ No newline at end of file diff --git a/Classes/DataBaseHandler/related functions/groupIt.m b/Classes/DataBaseHandler/related functions/groupIt.m new file mode 100644 index 0000000..85a1851 --- /dev/null +++ b/Classes/DataBaseHandler/related functions/groupIt.m @@ -0,0 +1,53 @@ +function resultTable = groupIt(fixedVars, dataTable, aggregationFunction) +% groupIt Groups data in a table based on fixedVars and applies aggregationFunction to numeric data. +% +% resultTable = groupIt(fixedVars, dataTable, aggregationFunction) +% +% Inputs: +% fixedVars - Cell array of variable names to group by +% dataTable - Input MATLAB table +% aggregationFunction - Function handle (e.g., @mean, @min, @max) +% +% Output: +% resultTable - Grouped and aggregated table + +[G, groupKeys] = findgroups(dataTable(:, fixedVars)); +varNames = dataTable.Properties.VariableNames; +nVars = numel(varNames); +aggData = cell(height(groupKeys), nVars); +groupCount = zeros(height(groupKeys), 1); + +for i = 1:height(groupKeys) + idx = (G == i); + groupCount(i) = sum(idx); + for j = 1:nVars + colData = dataTable.(varNames{j}); + if isnumeric(colData) + if any(idx) + aggData{i, j} = aggregationFunction(colData(idx)); + else + aggData{i, j} = NaN; + end + else + if iscell(colData) + nonEmptyIdx = find(idx & ~cellfun(@isempty, colData), 1); + if ~isempty(nonEmptyIdx) + aggData{i, j} = colData{nonEmptyIdx}; + else + aggData{i, j} = []; + end + else + nonEmptyIdx = find(idx, 1); + if ~isempty(nonEmptyIdx) + aggData{i, j} = colData(nonEmptyIdx); + else + aggData{i, j} = []; + end + end + end + end +end + +resultTable = cell2table(aggData, 'VariableNames', varNames); +resultTable.nRows = groupCount; +end \ No newline at end of file diff --git a/Classes/DataBaseHandler/related functions/removeGroupOutliers.m b/Classes/DataBaseHandler/related functions/removeGroupOutliers.m new file mode 100644 index 0000000..cd287da --- /dev/null +++ b/Classes/DataBaseHandler/related functions/removeGroupOutliers.m @@ -0,0 +1,61 @@ +function [cleanedTable, outliersTable] = removeGroupOutliers(dataTable, fixedVars, y_var) +% removeGroupOutliers removes outliers in y_var within each group defined by fixedVars. +% +% [cleanedTable, outliersTable] = removeGroupOutliers(dataTable, fixedVars, y_var) +% +% Inputs: +% dataTable - Input MATLAB table +% fixedVars - Cell array of variable names to group by +% y_var - Name of the variable to check for outliers (string or char) +% +% Outputs: +% cleanedTable - Table with outliers removed +% outliersTable - Table of removed outlier rows + +[G, groupKeys] = findgroups(dataTable(:, fixedVars)); +keepIdx = true(height(dataTable), 1); +outlierRecords = []; + +for groupIdx = 1:height(groupKeys) + groupRows = (G == groupIdx); + y_values = dataTable.(y_var)(groupRows); + + % Skip groups with fewer than 3 points + if sum(groupRows) < 3 + continue; + end + + % Detect outliers in log10 space (robust for BER, etc.) + y_log = log10(y_values); + outlierMask = isoutlier(y_log, 'quartiles', 1); + + if any(outlierMask) + groupData = dataTable(groupRows, :); + outlierGroupTable = groupData(outlierMask, :); + + % Optionally, add group key values for traceability + for k = 1:numel(fixedVars) + outlierGroupTable.(['Group_', fixedVars{k}]) = repmat(groupKeys{groupIdx, k}, height(outlierGroupTable), 1); + end + + outlierRecords = [outlierRecords; outlierGroupTable]; %#ok + end + + % Mark outliers for removal + groupRowIdx = find(groupRows); + keepIdx(groupRowIdx(outlierMask)) = false; +end + +cleanedTable = dataTable(keepIdx, :); + +if isempty(outlierRecords) + outliersTable = table(); +else + outliersTable = outlierRecords; +end + +nRemoved = sum(~keepIdx); +nTotalOriginal = height(dataTable); +fprintf('Removed %d outliers from the data table (%.2f%% of total %d entries).\n', ... + nRemoved, 100*nRemoved/nTotalOriginal, nTotalOriginal); +end \ No newline at end of file diff --git a/Datatypes/processingMode.m b/Datatypes/processingMode.m new file mode 100644 index 0000000..c47dc7a --- /dev/null +++ b/Datatypes/processingMode.m @@ -0,0 +1,8 @@ + +classdef processingMode < int32 + enumeration + serial (1) + parallel (2) + end + +end \ No newline at end of file diff --git a/Functions/EQ_structures/dsp_runid.m b/Functions/EQ_structures/dsp_runid.m new file mode 100644 index 0000000..c3d6a63 --- /dev/null +++ b/Functions/EQ_structures/dsp_runid.m @@ -0,0 +1,137 @@ +function [output] = dsp_runid(run_id, options) + +arguments + run_id + options.append_to_db = 0; + options.max_occurences = 4; + options.parameters = struct(); + options.database_path + options.database_name + options.storage_path +end + +% Initialize output structures +output.ffe_package = {}; +output.mlse_package = {}; +output.vnle_package = {}; +output.dbtgt_package = {}; +output.dbenc_package = {}; + +% Initialize database connection +database = DBHandler("pathToDB", [options.database_path, options.database_name], "type", "sqlite"); +dataTable = queryRunid(run_id, database); +fsym = dataTable.symbolrate; +M = double(dataTable.pam_level); +duob_mode = db_mode(dataTable.db_mode); + +if database.checkIfRunExists('Results','run_id',run_id) + disp(['Already got at least one reulst for run id: ',num2str(run_id),' ']) + return +end + + +% Handle Settings and argument replacement + +len_tr = 4096*2; + +ffe_order = [50, 5, 5]; +dfe_order = [0, 0, 0]; +pf_ncoeffs = 1; +mu_ffe = [0.0001, 0.0008, 0.001]; +mu_dfe = 0.0004; +mu_dc = 0.00; + +use_ffe = 1; +use_vnle_mlse = 1; +use_dbtgt = 1; +use_dbenc = 0; + +addProcessingResultToDatabase = 0; + +% Overwrite default parameters if given in options.parameters +paramStruct = options.parameters; +if ~isempty(paramStruct) + paramNames = fieldnames(paramStruct); + for i = 1:numel(paramNames) + thisName = paramNames{i}; + thisValue = paramStruct.(thisName); + eval([thisName ' = thisValue;']); + end +end + +% Configure equalizers +eq_lin = EQ("Ne",[50,0,0],"Nb",dfe_order,"training_length",len_tr,"training_loops",5,"dd_loops",5,"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1); + +eq_ = EQ("Ne",ffe_order,"Nb",dfe_order,"training_length",len_tr,"training_loops",5,"dd_loops",5,"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1); +pf_ = Postfilter("ncoeff",pf_ncoeffs,"useBurg",1); +mlse_ = MLSE_viterbi("duobinary_output",0,'M',M,'trellis_states',PAMmapper(M,0).levels); +mlse_db_ = MLSE_viterbi("DIR",[1,1],"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels); +eq_post = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",2001,"sps",1,"decide",0); +% Duobinary signaling (db encoded) +mlse_db_enc = MLSE_viterbi("DIR", [1,1], "duobinary_output", 0, "M", M, "trellis_states", PAMmapper(M,0).levels); +eq_db_enc = EQ("Ne", ffe_order, "Nb", dfe_order, "training_length", len_tr, ... + "training_loops", 5, "dd_loops", 5, "K", 2, "DCmu", mu_dc, ... + "DDmu", [mu_ffe mu_dfe], "DFEmu", 0.005, "FFEmu", 0, "plotfinal", 0, "ideal_dfe", 1); + + +% Load signal data +[Tx_bits, Symbols, Scpe_cell, ~] = loadSignalData(dataTable, options); + + +for r = 1:numel(Scpe_cell) + + % Preprocess signal + Scpe_sig = preprocessSignal(Scpe_cell{r}, Symbols, fsym); + + if duob_mode ~= db_mode.db_encoded + if use_ffe + ffe_results = ffe(eq_lin,M,Scpe_sig,Symbols,Tx_bits,... + "precode_mode",duob_mode,... + 'showAnalysis',1,... + "postFFE",[],... + "eth_style_symbol_mapping",0); + output.ffe_package{r} = ffe_results; + if options.append_to_db + database.addProcessingResult(run_id, ffe_results.metrics, ffe_results.config); + end + end + + if use_vnle_mlse + [ffe_results, mlse_results] = vnle_postfilter_mlse(eq_, pf_, mlse_, M, Scpe_sig, Symbols, Tx_bits, ... + "precode_mode", duob_mode,... + 'showAnalysis', 1, ... + "postFFE", [],... + "eth_style_symbol_mapping", 0); + output.mlse_package{r} = mlse_results; + output.vnle_package{r} = ffe_results; + if options.append_to_db + database.addProcessingResult(run_id, mlse_results.metrics, mlse_results.config); + database.addProcessingResult(run_id, ffe_results.metrics, ffe_results.config); + end + end + + if use_dbtgt + dbt_results = duobinary_target(eq_, mlse_db_, M, Scpe_sig, Symbols, Tx_bits, ... + "precode_mode", duob_mode, ... + 'showAnalysis', 0,... + "postFFE", []); + output.dbtgt_package{r} = dbt_results; + if options.append_to_db + database.addProcessingResult(run_id, dbt_results.metrics, dbt_results.config); + end + end + end + + if duob_mode == db_mode.db_encoded + db_results = duobinary_signaling(eq_db_enc, mlse_db_enc, M, Scpe_sig, Symbols, Tx_bits, "precode_mode",duob_mode, "showAnalysis",0,"postFFE",[]); + output.dbenc_package{r} = db_results; + if options.append_to_db + database.addProcessingResult(run_id, db_results.metrics, db_results.config); + end + end + + +end + + +end \ No newline at end of file diff --git a/Functions/EQ_structures/duobinary_signaling.m b/Functions/EQ_structures/duobinary_signaling.m index 69c0449..013e184 100644 --- a/Functions/EQ_structures/duobinary_signaling.m +++ b/Functions/EQ_structures/duobinary_signaling.m @@ -1,96 +1,129 @@ -function [eq_package] = duobinary_signaling(eq_, mlse_,M ,rx_signal, tx_symbols, tx_bits,options) -%Duobinary Signaling +function [db_results] = duobinary_signaling(eq_, mlse_, M, rx_signal, tx_symbols, tx_bits, options) +% DUOBINARY_SIGNALING Processes signals through duobinary signaling +% +% Inputs: +% eq_ - Equalizer object +% mlse_ - MLSE object +% M - Modulation order +% rx_signal - Received signal +% tx_symbols - Transmitted symbols +% tx_bits - Transmitted bits +% options - Optional parameters +% +% Outputs: +% db_results - Results from duobinary signaling processing + arguments - eq_ - mlse_ - M - rx_signal - tx_symbols - tx_bits - options.postFFE = []; + eq_ + mlse_ + M + rx_signal + tx_symbols + tx_bits + options.precode_mode db_mode + options.showAnalysis = 0 + options.eth_style_symbol_mapping = 0 + options.postFFE = [] + options.database = [] end +%% Process signals through equalizer +[eq_signal, eq_noise] = eq_.process(rx_signal, tx_symbols); -[eq_signal, eq_noise] = eq_.process(rx_signal,tx_symbols); - +% Apply post-FFE if provided if ~isempty(options.postFFE) - [eq_signal,eq_noise] = options.postFFE.process(eq_signal,tx_symbols); + [eq_signal, eq_noise] = options.postFFE.process(eq_signal, tx_symbols); end - +% Process through MLSE eq_signal = mlse_.process(eq_signal); +% Apply duobinary encoding and decoding eq_signal = Duobinary().encode(eq_signal); eq_signal = Duobinary().decode(eq_signal); -% M = numel(unique(eq_signal.signal)); -rx_bits = PAMmapper(M,0).demap(eq_signal); +% Demap symbols to bits +rx_bits = PAMmapper(M, 0, "eth_style", options.eth_style_symbol_mapping).demap(eq_signal); -[bits_db,errors_db,ber_db,~] = calc_ber(rx_bits.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); +%% Calculate BER and metrics +[bits_db, errors_db, ber_db, error_pos] = calc_ber(rx_bits.signal, tx_bits.signal, "skip_front", 100, "skip_end", 150, "returnErrorLocation", 1); -eq_package.ber = ber_db; - -resultsDBsignaling = struct( ... - 'result_id', NaN, ... % - 'run_id', NaN, ... % Beispielhafte Run-ID - 'eqParam_id', NaN, ... % Beispielhafter Fremdschlüssel zur EqualizerParameters-Tabelle - 'date_of_processing', datetime('now'), ... % Aktuelles Datum und Uhrzeit - 'numBits', bits_db, ... % Beispiel: 1.000.000 Bits - 'numBitErr', errors_db, ... % Beispiel: 120 Bitfehler - 'BER_precoded', [], ... % BER = 120 / 1.000.000 - 'numBitErr_precoded', [], ... % Beispiel: 120 Bitfehler - 'BER', ber_db, ... % BER = 120 / 1.000.000 - 'SNR', [], ... % Beispielhafte SNR - 'SNR_level', jsonencode([]), ... % SNR-Level als JSON-codiertes Array - 'GMI', [], ... % Beispielhafter GMI-Wert - 'AIR', [], ... % Beispielhafter AIR-Wert - 'EVM', [], ... % Beispielhafte EVM - 'EVM_level', jsonencode([]), ... % EVM-Level als JSON-codiertes Array - 'Alpha', [] ... % Beispielhafter Alpha-Wert - ); +% Calculate performance metrics +[snr, snr_lvl] = calc_snr(tx_symbols.signal, eq_noise.signal); +[gmi] = calc_air(eq_signal, tx_symbols, "skip_front", 10000, "skip_end", 10000); +air = tx_symbols.fs .* floor(log2(double(M))*10)/10 .* gmi ./ log2(double(M)); +[evm_total, evm_lvl] = calc_evm(eq_signal, tx_symbols); +[std_total, std_lvl] = calc_std(eq_signal, tx_symbols); +[std_rxraw_total, std_rxraw_lvl] = calc_std(rx_signal.resample("fs_out", tx_symbols.fs), tx_symbols); +%% Prepare output structure +% Determine postFFE order if ~isempty(options.postFFE) npostFFE = options.postFFE.order; else npostFFE = 0; end -equalizerConfigDBsignaling = struct( ... - 'eq_id', NaN, ... % Auto-Inkrement, wird in der DB gesetzt - 'equalizer_structure', int32(equalizer_structure.db_encoded), ... % Beispiel: 1 (z.B. für vnle) - 'M', M, ... % Ordnung der PAM-Konstellation - 'target_constellation', jsonencode(round(unique(tx_symbols.signal),5)), ... % Beispielhafter Target-String - 'db_target', 1, ... % 0 oder 1 - 'diff_precode', 1, ... % 0 oder 1 - 'postFFE', ~isempty(options.postFFE), ... % Beispielwert - 'NpostFFE', npostFFE, ... % Beispielwert - 'Ne1', eq_.Ne(1), ... % Feedforward Koeffizienten 1. Ordnung - 'Ne2', eq_.Ne(2), ... % Feedforward Koeffizienten 2. Ordnung - 'Ne3', eq_.Ne(3), ... % Feedforward Koeffizienten 3. Ordnung - 'Nb1', eq_.Nb(1), ... % Decision Feedback Koeffizienten 1. Ordnung - 'Nb2', eq_.Nb(2), ... % Decision Feedback Koeffizienten 2. Ordnung - 'Nb3', eq_.Nb(3), ... % Decision Feedback Koeffizienten 3. Ordnung - 'K', eq_.K, ... % Samples pro Symbol - 'DCmu', eq_.DCmu, ... % Anpassungsrate für DC-Tap - 'ideal_dfe', eq_.ideal_dfe, ... % Flag für ideal DFE (0 oder 1) - 'training_length', eq_.training_length, ... % Anzahl Trainingssymbole - 'training_loops', eq_.training_loops, ... % Anzahl Trainingsdurchläufe - 'TRmu1', eq_.FFEmu, ... % mu für DD-Modus (1. Ordnung) - 'TRmu2', eq_.FFEmu, ... % mu für DD-Modus (2. Ordnung) - 'TRmu3', eq_.FFEmu, ... % mu für DD-Modus (3. Ordnung) - 'TRmuDFE', eq_.DFEmu, ... % mu für DFE-Modus im DD - 'dd_loops', 5, ... % Anzahl Durchläufe im DD-Modus - 'DDmu1', eq_.DDmu(1), ... % mu für DD-Modus (1. Ordnung) - 'DDmu2', eq_.DDmu(2), ... % mu für DD-Modus (2. Ordnung) - 'DDmu3', eq_.DDmu(3), ... % mu für DD-Modus (3. Ordnung) - 'DDmuDFE', eq_.DDmu(4), ... % mu für DFE-Modus im DD - 'MLSE_mode', 'viterbi', ... % Beispiel: MLSE-Modus als String - 'MLSE_trellis_states', jsonencode(mlse_.trellis_states), ... % Trellis-States, z.B. als JSON-String oder kommasepariert - 'comment', 'function: duobinary_target.m', ... % Zusätzliche Kommentare - 'config_hash', NaN ... - ); +% Create results structure +db_results = struct(); +db_results.metrics = Metricstruct; +db_results.metrics.result_id = NaN; +db_results.metrics.run_id = NaN; +db_results.metrics.eqParam_id = NaN; +db_results.metrics.date_of_processing = datetime('now'); +db_results.metrics.BER = ber_db; +db_results.metrics.numBits = bits_db; +db_results.metrics.numBitErr = errors_db; +db_results.metrics.SNR = snr; +db_results.metrics.SNR_level = snr_lvl; +db_results.metrics.STD = std_total; +db_results.metrics.STD_level = std_lvl; +db_results.metrics.STDrx = std_rxraw_total; +db_results.metrics.STDrx_level = std_rxraw_lvl; +db_results.metrics.GMI = gmi; +db_results.metrics.AIR = air; +db_results.metrics.EVM = evm_total; +db_results.metrics.EVM_level = evm_lvl; +db_results.metrics.MLSE_dir = mlse_.DIR; -eq_package.resultsDBsignaling = resultsDBsignaling; -eq_package.equalizerConfigDBsignaling = equalizerConfigDBsignaling; +% Create configuration structure +eq_.e = []; +eq_.e2 = []; +eq_.e3 = []; +db_results.config = Equalizerstruct(); +db_results.config.eq = jsonencode(eq_); +mlse_.DIR = []; +db_results.config.mlse = jsonencode(mlse_); +db_results.config.equalizer_structure = int32(equalizer_structure.db_encoded); +db_results.config.comment = 'function: duobinary_signaling'; +%% Display analysis if requested +if options.showAnalysis + displayAnalysis(eq_noise, eq_signal, rx_signal, eq_, tx_symbols, M, options.postFFE); +end + +end + +%% Helper Function +function displayAnalysis(eq_noise, eq_signal, rx_signal, eq_, tx_symbols, M, postFFE) +% Display analysis plots and metrics +figure(336); +showEQNoisePSD(eq_noise, "fignum", 336, "displayname", 'Residual Noise after Duobinary'); + +if ~isempty(postFFE) + showEQcoefficients('n1', postFFE.e, "displayname", 'Coefficients', 'fignum', 338); +end + +showEQfilter(eq_.e, eq_signal.fs.*2); + +figure(341); clf; +showLevelHistogram(eq_signal, tx_symbols, "fignum", 341); + +warning off +figure(400); clf; +showLevelScatter(eq_signal, tx_symbols, "fignum", 400); + +figure(401); clf; +showLevelScatter(rx_signal.resample("fs_out", tx_symbols.fs), tx_symbols, "fignum", 401); +warning on end \ No newline at end of file diff --git a/Functions/EQ_structures/duobinary_target.m b/Functions/EQ_structures/duobinary_target.m index 8d66efb..b5ff07c 100644 --- a/Functions/EQ_structures/duobinary_target.m +++ b/Functions/EQ_structures/duobinary_target.m @@ -1,4 +1,4 @@ -function [eq_package] = duobinary_target(eq_, mlse_,M, rx_signal, tx_symbols, tx_bits, options) +function [db_results] = duobinary_target(eq_, mlse_,M, rx_signal, tx_symbols, tx_bits, options) arguments eq_ @@ -22,7 +22,7 @@ if ~isempty(options.postFFE) [eq_signal,eq_noise] = options.postFFE.process(eq_signal,db_ref_sequence); end -% dir = [1,1]; +mlse_.DIR = [1,1]; mlse_sig_sd = mlse_.process(eq_signal); mlse_sig_hd = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).quantize(mlse_sig_sd); @@ -73,70 +73,36 @@ end rx_bits = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(mlse_sig_hd); [bits_db,errors_db,ber_db,errorIndice_db] = calc_ber(rx_bits.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); -eq_package.ber = ber_db; -resultsDBtgt = struct( ... - 'result_id', NaN, ... % - 'run_id', NaN, ... % Beispielhafte Run-ID - 'eqParam_id', NaN, ... % Beispielhafter Fremdschlüssel zur EqualizerParameters-Tabelle - 'date_of_processing', datetime('now'), ... % Aktuelles Datum und Uhrzeit - 'numBits', bits_db, ... % Beispiel: 1.000.000 Bits - 'numBitErr', errors_db, ... % Beispiel: 120 Bitfehler - 'BER_precoded', ber_db_diff_precoded, ... % BER = 120 / 1.000.000 - 'numBitErr_precoded', errors_db_diff_precoded, ... % Beispiel: 120 Bitfehler - 'BER', ber_db, ... % BER = 120 / 1.000.000 - 'SNR', [], ... % Beispielhafte SNR - 'SNR_level', jsonencode([]), ... % SNR-Level als JSON-codiertes Array - 'GMI', [], ... % Beispielhafter GMI-Wert - 'AIR', [], ... % Beispielhafter AIR-Wert - 'EVM', [], ... % Beispielhafte EVM - 'EVM_level', jsonencode([]), ... % EVM-Level als JSON-codiertes Array - 'Alpha', [] ... % Beispielhafter Alpha-Wert - ); +db_results = struct(); +db_results.metrics = Metricstruct; +db_results.metrics.result_id = NaN; +db_results.metrics.run_id = NaN; +db_results.metrics.eqParam_id = NaN; +db_results.metrics.date_of_processing = datetime('now'); +db_results.metrics.BER = ber_db; +db_results.metrics.numBits = bits_db; +db_results.metrics.numBitErr = errors_db; +db_results.metrics.BER_precoded = ber_db_diff_precoded; +db_results.metrics.numBitErr_precoded = errors_db_diff_precoded; +db_results.metrics.SNR = NaN; +db_results.metrics.SNR_level = NaN; +db_results.metrics.GMI = NaN; +db_results.metrics.AIR = NaN; +db_results.metrics.EVM = NaN; +db_results.metrics.EVM_level = NaN; +db_results.metrics.MLSE_dir = mlse_.DIR; -if ~isempty(options.postFFE) - npostFFE = options.postFFE.order; -else - npostFFE = 0; -end - -equalizerConfigDBtgt = struct( ... - 'eq_id', NaN, ... % Auto-Inkrement, wird in der DB gesetzt - 'equalizer_structure', int32(equalizer_structure.vnle_db_mlse), ... % Beispiel: 1 (z.B. für vnle) - 'M', M, ... % Ordnung der PAM-Konstellation - 'target_constellation', jsonencode(round(db_ref_constellation,5)), ... % Beispielhafter Target-String - 'db_target', 1, ... % 0 oder 1 - 'diff_precode', int32(options.precode_mode), ... % 0 oder 1 - 'postFFE', ~isempty(options.postFFE), ... % Beispielwert - 'NpostFFE', npostFFE, ... % Beispielwert - 'Ne1', eq_.Ne(1), ... % Feedforward Koeffizienten 1. Ordnung - 'Ne2', eq_.Ne(2), ... % Feedforward Koeffizienten 2. Ordnung - 'Ne3', eq_.Ne(3), ... % Feedforward Koeffizienten 3. Ordnung - 'Nb1', eq_.Nb(1), ... % Decision Feedback Koeffizienten 1. Ordnung - 'Nb2', eq_.Nb(2), ... % Decision Feedback Koeffizienten 2. Ordnung - 'Nb3', eq_.Nb(3), ... % Decision Feedback Koeffizienten 3. Ordnung - 'K', eq_.K, ... % Samples pro Symbol - 'DCmu', eq_.DCmu, ... % Anpassungsrate für DC-Tap - 'ideal_dfe', eq_.ideal_dfe, ... % Flag für ideal DFE (0 oder 1) - 'training_length', eq_.training_length, ... % Anzahl Trainingssymbole - 'training_loops', eq_.training_loops, ... % Anzahl Trainingsdurchläufe - 'TRmu1', eq_.FFEmu, ... % mu für DD-Modus (1. Ordnung) - 'TRmu2', eq_.FFEmu, ... % mu für DD-Modus (2. Ordnung) - 'TRmu3', eq_.FFEmu, ... % mu für DD-Modus (3. Ordnung) - 'TRmuDFE', eq_.DFEmu, ... % mu für DFE-Modus im DD - 'dd_loops', 5, ... % Anzahl Durchläufe im DD-Modus - 'DDmu1', eq_.DDmu(1), ... % mu für DD-Modus (1. Ordnung) - 'DDmu2', eq_.DDmu(2), ... % mu für DD-Modus (2. Ordnung) - 'DDmu3', eq_.DDmu(3), ... % mu für DD-Modus (3. Ordnung) - 'DDmuDFE', eq_.DDmu(4), ... % mu für DFE-Modus im DD - 'MLSE_mode', 'viterbi', ... % Beispiel: MLSE-Modus als String - 'MLSE_trellis_states', jsonencode(mlse_.trellis_states), ... % Trellis-States, z.B. als JSON-String oder kommasepariert - 'comment', 'function: duobinary_target.m', ... % Zusätzliche Kommentare - 'config_hash', NaN ... - ); - -eq_package.resultsDBtgt = resultsDBtgt; -eq_package.equalizerConfigDBtgt = equalizerConfigDBtgt; +% Create DB results structure +db_results.config = Equalizerstruct(); +eq_.e = []; +eq_.e2 = []; +eq_.e3 = []; +db_results.config.eq = jsonencode(eq_); +mlse_.DIR = []; +db_results.config.mlse = jsonencode(mlse_); +db_results.config.equalizer_structure = int32(equalizer_structure.vnle_db_mlse); +db_results.config.comment = 'function: Duobinary tgt. (VNLE -> MLSE)'; if options.showAnalysis eq_noise = eq_noise - mean(eq_noise.signal); diff --git a/Functions/EQ_structures/ffe.m b/Functions/EQ_structures/ffe.m new file mode 100644 index 0000000..f9ff57f --- /dev/null +++ b/Functions/EQ_structures/ffe.m @@ -0,0 +1,163 @@ +function [ffe_results] = ffe(eq_, M, rx_signal, tx_symbols, tx_bits, options) +% FFE Processes signals through FFE equalizer +% +% Inputs: +% eq_ - Equalizer object +% M - Modulation order +% rx_signal - Received signal +% tx_symbols - Transmitted symbols +% tx_bits - Transmitted bits +% options - Optional parameters +% +% Outputs: +% ffe_results - Results from FFE processing + +arguments + eq_ + M + rx_signal + tx_symbols + tx_bits + options.precode_mode db_mode + options.showAnalysis = 0; + options.eth_style_symbol_mapping = 0; + options.postFFE = []; + options.database = []; +end + +%% Process signals through equalizer +% FFE or VNLE +[eq_signal_sd, eq_noise] = eq_.process(rx_signal, tx_symbols); + +% Apply post-FFE if provided +if ~isempty(options.postFFE) + tic + [eq_signal_sd, eq_noise] = options.postFFE.process(eq_signal_sd, tx_symbols); + toc +end + +% Hard decision on FFE output +eq_signal_hd = PAMmapper(M, 0).quantize(eq_signal_sd); + +%% Calculate BER based on precoding mode +[bits, errors, ber, error_pos, errors_precoded, ber_precoded] = calculateBER(eq_signal_hd, tx_symbols, tx_bits, options.precode_mode, M, options.eth_style_symbol_mapping); + +%% Calculate performance metrics +[snr, snr_lvl] = calc_snr(tx_symbols.signal, eq_noise.signal); +[gmi] = calc_air(eq_signal_sd, tx_symbols, "skip_front", 10000, "skip_end", 10000); +air = tx_symbols.fs .* floor(log2(double(M))*10)/10 .* gmi ./ log2(double(M)); +[evm_total, evm_lvl] = calc_evm(eq_signal_sd, tx_symbols); +[std_total, std_lvl] = calc_std(eq_signal_sd, tx_symbols); +[std_rxraw_total, std_rxraw_lvl] = calc_std(rx_signal.resample("fs_out", tx_symbols.fs), tx_symbols); + +%% Display analysis if requested +if options.showAnalysis + displayAnalysis(eq_noise, eq_signal_sd, rx_signal, eq_, tx_symbols, M, options.postFFE); +end + + +%% Prepare output structure +% Determine postFFE order +if ~isempty(options.postFFE) + npostFFE = options.postFFE.order; +else + npostFFE = 0; +end + +% Create FFE results structure +ffe_results = struct(); + +eq_.e = []; +eq_.e2 = []; +eq_.e3 = []; +ffe_results.config = Equalizerstruct(); +ffe_results.config.eq = jsonencode(eq_); +ffe_results.config.equalizer_structure = int32(equalizer_structure.ffe); +ffe_results.config.comment = 'function: ffe'; + +ffe_results.metrics = Metricstruct; +ffe_results.metrics.result_id = NaN; +ffe_results.metrics.run_id = NaN; +ffe_results.metrics.eqParam_id = NaN; +ffe_results.metrics.date_of_processing = datetime('now'); +ffe_results.metrics.BER = ber; +ffe_results.metrics.numBits = bits; +ffe_results.metrics.numBitErr = errors; +ffe_results.metrics.BER_precoded = ber_precoded; +ffe_results.metrics.numBitErr_precoded = errors_precoded; +ffe_results.metrics.SNR = snr; +ffe_results.metrics.SNR_level = snr_lvl; +ffe_results.metrics.STD = std_total; +ffe_results.metrics.STD_level = std_lvl; +ffe_results.metrics.STDrx = std_rxraw_total; +ffe_results.metrics.STDrx_level = std_rxraw_lvl; +ffe_results.metrics.GMI = gmi; +ffe_results.metrics.AIR = air; +ffe_results.metrics.EVM = evm_total; +ffe_results.metrics.EVM_level = evm_lvl; + + + +end + +%% Helper Functions +function [bits, errors, ber, error_pos, errors_precoded, ber_precoded] = calculateBER(eq_signal_hd, tx_symbols, tx_bits, precode_mode, M, eth_style) +% Calculate BER based on precoding mode +mapper = PAMmapper(M, 0, "eth_style", eth_style); + +switch precode_mode + case db_mode.no_db + % TX Data is not precoded + % A) Emulate diff precoding + eq_signal_hd_precoded = Duobinary().encode(eq_signal_hd, "M", M); + eq_signal_hd_precoded = Duobinary().decode(eq_signal_hd_precoded, "M", M); + + tx_symbols_precoded = Duobinary().encode(tx_symbols); + tx_symbols_precoded = Duobinary().decode(tx_symbols_precoded); + + tx_bits_precoded = mapper.demap(tx_symbols_precoded); + + rx_bits = mapper.demap(eq_signal_hd_precoded); + [~, errors_precoded, ber_precoded, ~] = calc_ber(rx_bits.signal, tx_bits_precoded.signal, "skip_front", 30000, "skip_end", 150, "returnErrorLocation", 1); + + % B) Just determine BER + rx_bits = mapper.demap(eq_signal_hd); + [bits, errors, ber, error_pos] = calc_ber(rx_bits.signal, tx_bits.signal, "skip_front", 30000, "skip_end", 150, "returnErrorLocation", 1); + + case db_mode.db_precoded + % Data is precoded on TX side + % A) Decode at Rx if no DB targeting was applied + eq_signal_hd_decoded = Duobinary().encode(eq_signal_hd, "M", M); + eq_signal_hd_decoded = Duobinary().decode(eq_signal_hd_decoded, "M", M); + rx_bits_decoded = mapper.demap(eq_signal_hd_decoded); + [~, errors_precoded, ber_precoded, ~] = calc_ber(rx_bits_decoded.signal, tx_bits.signal, "skip_front", 30000, "skip_end", 150, "returnErrorLocation", 1); + + % B) Omit the Coding by comparing with demapped TX symbol sequence + tx_bits_demapped = mapper.demap(tx_symbols); + rx_bits = mapper.demap(eq_signal_hd); + [bits, errors, ber, error_pos] = calc_ber(rx_bits.signal, tx_bits_demapped.signal, "skip_front", 30000, "skip_end", 150, "returnErrorLocation", 1); +end +end + +function displayAnalysis(eq_noise, eq_signal_sd, rx_signal, eq_, tx_symbols, M, postFFE) +% Display analysis plots and metrics +figure(336); +showEQNoisePSD(eq_noise, "fignum", 336, "displayname", 'Residual Noise after FFE'); + +if ~isempty(postFFE) + showEQcoefficients('n1', postFFE.e, "displayname", 'Coefficients', 'fignum', 338); +end + +showEQfilter(eq_.e, eq_signal_sd.fs.*2); + +figure(341); clf; +showLevelHistogram(eq_signal_sd, tx_symbols, "fignum", 341); + +warning off +figure(400); clf; +showLevelScatter(eq_signal_sd, tx_symbols, "fignum", 400); + +figure(401); clf; +showLevelScatter(rx_signal.resample("fs_out", tx_symbols.fs), tx_symbols, "fignum", 401); +warning on +end \ No newline at end of file diff --git a/Functions/EQ_structures/vnle.m b/Functions/EQ_structures/vnle.m deleted file mode 100644 index f143452..0000000 --- a/Functions/EQ_structures/vnle.m +++ /dev/null @@ -1,192 +0,0 @@ -function [eq_package] = vnle(eq_,M,rx_signal,tx_symbols,tx_bits,options) - -arguments - eq_ - M - rx_signal - tx_symbols - tx_bits - options.precode_mode db_mode - options.showAnalysis = 0; - options.eth_style_symbol_mapping = 0; - options.postFFE = []; - options.database = []; -end - -mudc_given = eq_.mu_dc; - -%FFE or VNLE -[eq_signal_sd,eq_noise] = eq_.process(rx_signal,tx_symbols); - -if ~isempty(options.postFFE) - tic - [eq_signal_sd,eq_noise] = options.postFFE.process(eq_signal_sd,tx_symbols); - toc -end - -eq_signal_hd = PAMmapper(M,0).quantize(eq_signal_sd); - -% precoding to mitigate error propagation, most prominently used in -% combination with duobinary signaling to avoid catastrophic error -% behavior (see J.W.M. Bergmans, Digital Baseband Transmission and Recording -> partial response signaling) - -switch options.precode_mode - - case db_mode.no_db - % TX Data is not precoded: - - % A) Emulate diff precoding - eq_signal_hd_precoded = Duobinary().encode(eq_signal_hd,"M",M); - eq_signal_hd_precoded = Duobinary().decode(eq_signal_hd_precoded,"M",M); - - tx_symbols_precoded = Duobinary().encode(tx_symbols); - tx_symbols_precoded = Duobinary().decode(tx_symbols_precoded); - - tx_bits_precoded = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(tx_symbols_precoded); - - rx_bits_vnle = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(eq_signal_hd_precoded); - [~,errors_vnle_diff_precoded,ber_vnle_diff_precoded,~] = calc_ber(rx_bits_vnle.signal,tx_bits_precoded.signal,"skip_front",30000,"skip_end",150,"returnErrorLocation",1); - - %B) Just determine BER - rx_bits_vnle = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(eq_signal_hd); - [bits_vnle,errors_vnle,ber_vnle,error_pos_vnle] = calc_ber(rx_bits_vnle.signal,tx_bits.signal,"skip_front",30000,"skip_end",150,"returnErrorLocation",1); - - max_burst_length = 10; - burst_count = count_error_bursts(error_pos_vnle, max_burst_length); - - case db_mode.db_precoded - - % Daten SIND TATSÄCHLICH precoded auf TX Seite: - - % A) Decode at Rx if no DB targeting was applied (we are in VNLE or MLSE EQ structure here! - eq_signal_hd_decoded = Duobinary().encode(eq_signal_hd,"M",M); - eq_signal_hd_decoded = Duobinary().decode(eq_signal_hd_decoded,"M",M); - rx_bits_vnle_decoded = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(eq_signal_hd_decoded); - [~,errors_vnle_diff_precoded,ber_vnle_diff_precoded,~] = calc_ber(rx_bits_vnle_decoded.signal,tx_bits.signal,"skip_front",30000,"skip_end",150,"returnErrorLocation",1); - - % B) Omit the Coding by comparing with demapped TX symbol sequence - - tx_bits = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(tx_symbols); - rx_bits_vnle = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(eq_signal_hd); - [bits_vnle,errors_vnle,ber_vnle,~] = calc_ber(rx_bits_vnle.signal,tx_bits.signal,"skip_front",30000,"skip_end",150,"returnErrorLocation",1); - - -end - -% METRICS OF VNLE SD Signal: -[snr_vnle,snr_vnle_lvl] = calc_snr(tx_symbols.signal,eq_noise.signal); -[gmi_vnle] = calc_air(eq_signal_sd,tx_symbols,"skip_front",10000,"skip_end",10000); -air_vnle = tx_symbols.fs .* floor(log2(double(M))*10)/10 .* gmi_vnle ./ log2(double(M)); -[evm_vnle_total,evm_vnle_lvl] = calc_evm(eq_signal_sd,tx_symbols); -[std_vnle_total,std_vnle_lvl] = calc_std(eq_signal_sd,tx_symbols); -[std_rxraw_total,std_rxraw_lvl] = calc_std(rx_signal.resample("fs_out",tx_symbols.fs),tx_symbols); - -eq_package.ber_vnle = ber_vnle; -eq_package.evm_vnle_total = evm_vnle_total; -eq_package.evm_vnle_lvl = evm_vnle_lvl; -eq_package.gmi = gmi_vnle; - -eq_package.eq = eq_; - -resultsVNLE = struct( ... - 'result_id', NaN, ... % - 'run_id', NaN, ... % Beispielhafte Run-ID - 'eqParam_id', NaN, ... % Beispielhafter Fremdschlüssel zur EqualizerParameters-Tabelle - 'date_of_processing', datetime('now'), ... % Aktuelles Datum und Uhrzeit - 'BER', ber_vnle, ... % BER = 120 / 1.000.000 - 'numBits', bits_vnle, ... % Beispiel: 1.000.000 Bits - 'numBitErr', errors_vnle, ... % Beispiel: 120 Bitfehler - 'BER_precoded', ber_vnle_diff_precoded, ... % BER = 120 / 1.000.000 - 'numBitErr_precoded', errors_vnle_diff_precoded, ... % Beispiel: 120 Bitfehler - 'SNR', snr_vnle, ... % Beispielhafte SNR - 'SNR_level', jsonencode(snr_vnle_lvl), ... % SNR-Level als JSON-codiertes Array - 'STD', std_vnle_total, ... - 'STD_level', jsonencode(std_vnle_lvl),... - 'STDrx' , std_rxraw_total, ... - 'STDrx_level', jsonencode(std_rxraw_lvl),... - 'GMI', gmi_vnle, ... % Beispielhafter GMI-Wert - 'AIR', air_vnle, ... % Beispielhafter AIR-Wert - 'EVM', evm_vnle_total, ... % Beispielhafte EVM - 'EVM_level', jsonencode(evm_vnle_lvl), ... % EVM-Level als JSON-codiertes Array - 'Alpha', [] ... % Beispielhafter Alpha-Wert - ); - -if ~isempty(options.postFFE) - npostFFE = options.postFFE.order; -else - npostFFE = 0; -end - -equalizerConfigVNLE = struct( ... - 'eq_id', NaN, ... % Auto-Inkrement, wird in der DB gesetzt - 'equalizer_structure', int32(equalizer_structure.vnle), ... % Beispiel: 1 (z.B. für vnle) - 'M', M, ... % Ordnung der PAM-Konstellation - 'target_constellation', jsonencode(round(unique(tx_symbols.signal),5)), ... % Beispielhafter Target-String - 'db_target', 0, ... % 0 oder 1 - 'diff_precode', int32(options.precode_mode), ... % 0 oder 1 - 'postFFE', int32(~isempty(options.postFFE)), ... % Beispielwert - 'NpostFFE', npostFFE, ... % Beispielwert - 'Ne1', eq_.order, ... % Feedforward Koeffizienten 1. Ordnung - 'K', eq_.sps, ... % Samples pro Symbol - 'DCmu', mudc_given, ... % Anpassungsrate für DC-Tap - 'training_length', eq_.len_tr, ... % Anzahl Trainingssymbole - 'training_loops', eq_.epochs_tr, ... % Anzahl Trainingsdurchläufe - 'TRmu1', eq_.mu_tr, ... % mu für DD-Modus (1. Ordnung) - 'dd_loops', eq_.epochs_dd, ... % Anzahl Durchläufe im DD-Modus - 'DDmu1', eq_.mu_dd, ... % mu für DD-Modus (1. Ordnung) - 'comment', 'function: ffe dc removal', ... % Zusätzliche Kommentare - 'ffe_buffer_len', eq_.ffe_buffer_len, ... - 'smoothing_buffer_len', eq_.smoothing_buffer_length, ... - 'smoothing_buffer_update', eq_.smoothing_buffer_update, ... - 'dc_buffer_len', eq_.dc_buffer_len, ... - 'config_hash', NaN ... - ); - -eq_package.resultsVNLE = resultsVNLE; -eq_package.equalizerConfigVNLE = equalizerConfigVNLE; - - -% eq_package.vnle_out = eq_signal_sd; -if options.showAnalysis - - - % fprintf(['VNLE EVM lvl: ',repmat('%.3f ',1,numel(evm_lvl)),' \n'],evm_lvl); - - % fprintf('VNLE BER: %.2e \n',ber_vnle); - % - % fprintf('MLSE BER: %.2e \n',ber_mlse); - - figure(336); - showEQNoisePSD(eq_noise,"fignum",336,"displayname",'Residual Noise after VNLE'); - - % figure(337);clf; - % rx_signal.spectrum("normalizeTo0dB",1,"fignum",337,"displayname",'Rx Signal'); - - % figure(338);clf; - % showEQcoefficients('n1',eq_.e,'n2',eq_.e2,'n3',eq_.e3,"displayname",'Coefficients','fignum',338); - - if ~isempty(options.postFFE) - showEQcoefficients('n1',options.postFFE.e,"displayname",'Coefficients','fignum',338); - end - - showEQfilter(eq_.e,eq_signal_sd.fs.*2); - - % figure(340);clf; - % eq_signal_sd.eye(eq_signal_sd.fs,M,"fignum",340); - - figure(341);clf; - showLevelHistogram(eq_signal_sd,tx_symbols,"fignum",341); - - figure(400);clf; - warning off - showLevelScatter(eq_signal_sd,tx_symbols,"fignum",400); - - showLevelScatter(rx_signal.resample("fs_out",tx_symbols.fs),tx_symbols,"fignum",400); - warning on - % autoArrangeFigures(3,3,2) - - -end - -end - diff --git a/Functions/EQ_structures/vnle_postfilter_mlse.m b/Functions/EQ_structures/vnle_postfilter_mlse.m index 85a17e1..90691d8 100644 --- a/Functions/EQ_structures/vnle_postfilter_mlse.m +++ b/Functions/EQ_structures/vnle_postfilter_mlse.m @@ -1,4 +1,19 @@ -function [eq_package] = vnle_postfilter_mlse(eq_,pf_,mlse_,M,rx_signal,tx_symbols,tx_bits,options) +function [ffe_results, mlse_results] = vnle_postfilter_mlse(eq_, pf_, mlse_, M, rx_signal, tx_symbols, tx_bits, options) +% VNLE_POSTFILTER_MLSE Processes signals through VNLE, postfilter, and MLSE +% +% Inputs: +% eq_ - Equalizer object +% pf_ - Postfilter object +% mlse_ - MLSE object +% M - Modulation order +% rx_signal - Received signal +% tx_symbols - Transmitted symbols +% tx_bits - Transmitted bits +% options - Optional parameters +% +% Outputs: +% ffe_results - Results from FFE/VNLE processing +% mlse_results - Results from MLSE processing arguments eq_ @@ -15,274 +30,202 @@ arguments options.database = []; end -%FFE or VNLE -[eq_signal_sd,eq_noise] = eq_.process(rx_signal,tx_symbols); +%% Process signals through equalizers +% FFE or VNLE +[eq_signal_sd, eq_noise] = eq_.process(rx_signal, tx_symbols); +% Apply post-FFE if provided if ~isempty(options.postFFE) tic - [eq_signal_sd,eq_noise] = options.postFFE.process(eq_signal_sd,tx_symbols); + [eq_signal_sd, eq_noise] = options.postFFE.process(eq_signal_sd, tx_symbols); toc end -eq_signal_hd = PAMmapper(M,0).quantize(eq_signal_sd); - -mlse_sig_sd = pf_.process(eq_signal_sd,eq_noise); +% Hard decision on VNLE output +eq_signal_hd = PAMmapper(M, 0).quantize(eq_signal_sd); +% Process through postfilter and MLSE +mlse_sig_sd = pf_.process(eq_signal_sd, eq_noise); mlse_.DIR = pf_.coefficients; -% [mlse_sig_hd,mlse_sig_sd] = mlse_.process(mlse_sig_sd,tx_symbols); mlse_sig_sd = mlse_.process(mlse_sig_sd); +mlse_sig_hd = PAMmapper(M, 0, "eth_style", options.eth_style_symbol_mapping).quantize(mlse_sig_sd); -mlse_sig_hd = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).quantize(mlse_sig_sd); +%% Calculate BER based on precoding mode +[numbits, errors, bers, ~] = calculateBER(eq_signal_hd, mlse_sig_hd, tx_symbols, tx_bits, options.precode_mode, M, options.eth_style_symbol_mapping); -% precoding to mitigate error propagation, most prominently used in -% combination with duobinary signaling to avoid catastrophic error -% behavior (see J.W.M. Bergmans, Digital Baseband Transmission and Recording -> partial response signaling) +%% Calculate performance metrics +% VNLE metrics +[snr_vnle, snr_vnle_lvl] = calc_snr(tx_symbols.signal, eq_noise.signal); +[gmi_vnle] = calc_air(eq_signal_sd, tx_symbols, "skip_front", 10000, "skip_end", 10000); +air_vnle = tx_symbols.fs .* floor(log2(double(M))*10)/10 .* gmi_vnle ./ log2(double(M)); +[evm_vnle_total, evm_vnle_lvl] = calc_evm(eq_signal_sd, tx_symbols); +[std_vnle_total, std_vnle_lvl] = calc_std(eq_signal_sd, tx_symbols); +[std_rxraw_total, std_rxraw_lvl] = calc_std(rx_signal.resample("fs_out", tx_symbols.fs), tx_symbols); -switch options.precode_mode - - case db_mode.no_db - % TX Data is not precoded: - - % A) Emulate diff precoding - eq_signal_hd_precoded = Duobinary().encode(eq_signal_hd,"M",M); - eq_signal_hd_precoded = Duobinary().decode(eq_signal_hd_precoded,"M",M); - - mlse_sig_hd_precoded = Duobinary().encode(mlse_sig_hd,"M",M); - mlse_sig_hd_precoded = Duobinary().decode(mlse_sig_hd_precoded,"M",M); - - tx_symbols_precoded = Duobinary().encode(tx_symbols); - tx_symbols_precoded = Duobinary().decode(tx_symbols_precoded); - - tx_bits_precoded = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(tx_symbols_precoded); - - rx_bits_vnle = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(eq_signal_hd_precoded); - [~,errors_vnle_diff_precoded,ber_vnle_diff_precoded,~] = calc_ber(rx_bits_vnle.signal,tx_bits_precoded.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - rx_bits_mlse = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(mlse_sig_hd_precoded); - [~,errors_mlse_diff_precoded,ber_mlse_diff_precoded,~] = calc_ber(rx_bits_mlse.signal,tx_bits_precoded.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - %B) Just determine BER - rx_bits_vnle = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(eq_signal_hd); - [bits_vnle,errors_vnle,ber_vnle,~] = calc_ber(rx_bits_vnle.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - rx_bits_mlse = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(mlse_sig_hd); - [bits_mlse,errors_mlse,ber_mlse,~] = calc_ber(rx_bits_mlse.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - case db_mode.db_precoded - - % Daten SIND TATSÄCHLICH precoded auf TX Seite: - - % A) Decode at Rx if no DB targeting was applied (we are in VNLE or MLSE EQ structure here! - eq_signal_hd_decoded = Duobinary().encode(eq_signal_hd,"M",M); - eq_signal_hd_decoded = Duobinary().decode(eq_signal_hd_decoded,"M",M); - rx_bits_vnle_decoded = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(eq_signal_hd_decoded); - [~,errors_vnle_diff_precoded,ber_vnle_diff_precoded,~] = calc_ber(rx_bits_vnle_decoded.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - mlse_sig_hd_decoded = Duobinary().encode(mlse_sig_hd,"M",M); - mlse_sig_hd_decoded = Duobinary().decode(mlse_sig_hd_decoded,"M",M); - rx_bits_mlse_decoded = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(mlse_sig_hd_decoded); - [~,errors_mlse_diff_precoded,ber_mlse_diff_precoded,~] = calc_ber(rx_bits_mlse_decoded.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - % B) Omit the Coding by comparing with demapped TX symbol sequence - - tx_bits = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(tx_symbols); - rx_bits_vnle = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(eq_signal_hd); - [bits_vnle,errors_vnle,ber_vnle,~] = calc_ber(rx_bits_vnle.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - - rx_bits_mlse = PAMmapper(M,0,"eth_style",options.eth_style_symbol_mapping).demap(mlse_sig_hd); - [bits_mlse,errors_mlse,ber_mlse,~] = calc_ber(rx_bits_mlse.signal,tx_bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); - -end - -% METRICS OF VNLE SD Signal: -[snr_vnle,snr_vnle_lvl] = calc_snr(tx_symbols.signal,eq_noise.signal); -[gmi_vnle] = calc_air(eq_signal_sd,tx_symbols,"skip_front",10000,"skip_end",10000); -air_vnle = tx_symbols.fs .* floor(log2(double(M))*10)/10 .* gmi_vnle ./ log2(double(M)); -[evm_vnle_total,evm_vnle_lvl] = calc_evm(eq_signal_sd,tx_symbols); -[std_vnle_total,std_vnle_lvl] = calc_std(eq_signal_sd,tx_symbols); -[std_rxraw_total,std_rxraw_lvl] = calc_std(rx_signal.resample("fs_out",tx_symbols.fs),tx_symbols); - -% METRICS OF MLSE (HD-VITERBI) -pf_.ncoeff = 1; -pf_.process(eq_signal_sd,eq_noise); +% MLSE metrics alpha = pf_.coefficients(2); -eq_package.ber_mlse = ber_mlse; -eq_package.ber_vnle = ber_vnle; -eq_package.evm_vnle_total = evm_vnle_total; -eq_package.evm_vnle_lvl = evm_vnle_lvl; -eq_package.gmi = gmi_vnle; - -eq_package.eq = eq_; -eq_package.pf = pf_; -eq_package.mlse = mlse_; - -resultsVNLE = struct( ... - 'result_id', NaN, ... % - 'run_id', NaN, ... % Beispielhafte Run-ID - 'eqParam_id', NaN, ... % Beispielhafter Fremdschlüssel zur EqualizerParameters-Tabelle - 'date_of_processing', datetime('now'), ... % Aktuelles Datum und Uhrzeit - 'BER', ber_vnle, ... % BER = 120 / 1.000.000 - 'numBits', bits_vnle, ... % Beispiel: 1.000.000 Bits - 'numBitErr', errors_vnle, ... % Beispiel: 120 Bitfehler - 'BER_precoded', ber_vnle_diff_precoded, ... % BER = 120 / 1.000.000 - 'numBitErr_precoded', errors_vnle_diff_precoded, ... % Beispiel: 120 Bitfehler - 'SNR', snr_vnle, ... % Beispielhafte SNR - 'SNR_level', jsonencode(snr_vnle_lvl), ... % SNR-Level als JSON-codiertes Array - 'STD', std_vnle_total, ... - 'STD_level', jsonencode(std_vnle_lvl),... - 'STDrx' , std_rxraw_total, ... - 'STDrx_level', jsonencode(std_rxraw_lvl),... - 'GMI', gmi_vnle, ... % Beispielhafter GMI-Wert - 'AIR', air_vnle, ... % Beispielhafter AIR-Wert - 'EVM', evm_vnle_total, ... % Beispielhafte EVM - 'EVM_level', jsonencode(evm_vnle_lvl), ... % EVM-Level als JSON-codiertes Array - 'Alpha', [] ... % Beispielhafter Alpha-Wert - ); - +%% Prepare output structures +% Determine postFFE order if ~isempty(options.postFFE) npostFFE = options.postFFE.order; else npostFFE = 0; end -equalizerConfigVNLE = struct( ... - 'eq_id', NaN, ... % Auto-Inkrement, wird in der DB gesetzt - 'equalizer_structure', int32(equalizer_structure.vnle), ... % Beispiel: 1 (z.B. für vnle) - 'M', M, ... % Ordnung der PAM-Konstellation - 'target_constellation', jsonencode(round(unique(tx_symbols.signal),5)), ... % Beispielhafter Target-String - 'db_target', 0, ... % 0 oder 1 - 'diff_precode', int32(options.precode_mode), ... % 0 oder 1 - 'postFFE', ~isempty(options.postFFE), ... % Beispielwert - 'NpostFFE', npostFFE, ... % Beispielwert - 'Ne1', eq_.Ne(1), ... % Feedforward Koeffizienten 1. Ordnung - 'Ne2', eq_.Ne(2), ... % Feedforward Koeffizienten 2. Ordnung - 'Ne3', eq_.Ne(3), ... % Feedforward Koeffizienten 3. Ordnung - 'Nb1', eq_.Nb(1), ... % Decision Feedback Koeffizienten 1. Ordnung - 'Nb2', eq_.Nb(2), ... % Decision Feedback Koeffizienten 2. Ordnung - 'Nb3', eq_.Nb(3), ... % Decision Feedback Koeffizienten 3. Ordnung - 'K', eq_.K, ... % Samples pro Symbol - 'DCmu', eq_.DCmu, ... % Anpassungsrate für DC-Tap - 'ideal_dfe', eq_.ideal_dfe, ... % Flag für ideal DFE (0 oder 1) - 'training_length', eq_.training_length, ... % Anzahl Trainingssymbole - 'training_loops', eq_.training_loops, ... % Anzahl Trainingsdurchläufe - 'TRmu1', eq_.FFEmu, ... % mu für DD-Modus (1. Ordnung) - 'TRmu2', eq_.FFEmu, ... % mu für DD-Modus (2. Ordnung) - 'TRmu3', eq_.FFEmu, ... % mu für DD-Modus (3. Ordnung) - 'TRmuDFE', eq_.DFEmu, ... % mu für DFE-Modus im DD - 'dd_loops', 5, ... % Anzahl Durchläufe im DD-Modus - 'DDmu1', eq_.DDmu(1), ... % mu für DD-Modus (1. Ordnung) - 'DDmu2', eq_.DDmu(2), ... % mu für DD-Modus (2. Ordnung) - 'DDmu3', eq_.DDmu(3), ... % mu für DD-Modus (3. Ordnung) - 'DDmuDFE', eq_.DDmu(4), ... % mu für DFE-Modus im DD - 'comment', 'function: vnle_postfilter_mlse', ... % Zusätzliche Kommentare - 'config_hash', NaN ... - ); +ffe_results = struct(); +ffe_results.metrics = Metricstruct; +ffe_results.metrics.result_id = NaN; +ffe_results.metrics.run_id = NaN; +ffe_results.metrics.eqParam_id = NaN; +ffe_results.metrics.date_of_processing = datetime('now'); +ffe_results.metrics.BER = bers.vnle; +ffe_results.metrics.numBits = numbits.vnle; +ffe_results.metrics.numBitErr = errors.vnle; +ffe_results.metrics.BER_precoded = bers.vnle_precoded; +ffe_results.metrics.numBitErr_precoded = errors.vnle_precoded; +ffe_results.metrics.SNR = snr_vnle; +ffe_results.metrics.SNR_level = snr_vnle_lvl; +ffe_results.metrics.STD = std_vnle_total; +ffe_results.metrics.STD_level = std_vnle_lvl; +ffe_results.metrics.STDrx = std_rxraw_total; +ffe_results.metrics.STDrx_level = std_rxraw_lvl; +ffe_results.metrics.GMI = gmi_vnle; +ffe_results.metrics.AIR = air_vnle; +ffe_results.metrics.EVM = evm_vnle_total; +ffe_results.metrics.EVM_level = evm_vnle_lvl; +ffe_results.metrics.Alpha = []; -resultsMLSE = struct( ... - 'result_id', NaN, ... % - 'run_id', NaN, ... % Beispielhafte Run-ID - 'eqParam_id', NaN, ... % Beispielhafter Fremdschlüssel zur EqualizerParameters-Tabelle - 'date_of_processing', datetime('now'), ... % Aktuelles Datum und Uhrzeit - 'BER', ber_mlse, ... % BER = 120 / 1.000.000 - 'numBits', bits_mlse, ... % Beispiel: 1.000.000 Bits - 'numBitErr', errors_mlse, ... % Beispiel: 120 Bitfehler - 'BER_precoded', ber_mlse_diff_precoded, ... % BER = 120 / 1.000.000 - 'numBitErr_precoded', errors_mlse_diff_precoded, ... % Beispiel: 120 Bitfehler - 'SNR', [], ... % Beispielhafte SNR - 'SNR_level', jsonencode([]), ... % SNR-Level als JSON-codiertes Array - 'GMI', [], ... % Beispielhafter GMI-Wert - 'AIR', [], ... % Beispielhafter AIR-Wert - 'EVM', [], ... % Beispielhafte EVM - 'EVM_level', jsonencode([]), ... % EVM-Level als JSON-codiertes Array - 'Alpha', alpha, ... % Beispielhafter Alpha-Wert - 'MLSE_dir', jsonencode([mlse_.DIR])... - ); +% Create FFE results structure +eq_.e = []; +eq_.e2 = []; +eq_.e3 = []; +ffe_results.config = Equalizerstruct(); +ffe_results.config.eq = jsonencode(eq_); +ffe_results.config.equalizer_structure = int32(equalizer_structure.vnle); +ffe_results.config.comment = 'function: vnle_postfilter_mlse - FFE part'; -equalizerConfigMLSE = struct( ... - 'eq_id', NaN, ... % Auto-Inkrement, wird in der DB gesetzt - 'equalizer_structure', int32(equalizer_structure.vnle_pf_mlse), ... % Beispiel: 1 (z.B. für vnle) - 'M', M, ... % Ordnung der PAM-Konstellation - 'target_constellation', jsonencode(round(unique(tx_symbols.signal),5)), ... % Beispielhafter Target-String - 'db_target', 0, ... % 0 oder 1 - 'diff_precode', int32(options.precode_mode), ... % 0 oder 1 - 'postFFE', ~isempty(options.postFFE), ... % Beispielwert - 'NpostFFE', npostFFE, ... % Beispielwert - 'Ne1', eq_.Ne(1), ... % Feedforward Koeffizienten 1. Ordnung - 'Ne2', eq_.Ne(2), ... % Feedforward Koeffizienten 2. Ordnung - 'Ne3', eq_.Ne(3), ... % Feedforward Koeffizienten 3. Ordnung - 'Nb1', eq_.Nb(1), ... % Decision Feedback Koeffizienten 1. Ordnung - 'Nb2', eq_.Nb(2), ... % Decision Feedback Koeffizienten 2. Ordnung - 'Nb3', eq_.Nb(3), ... % Decision Feedback Koeffizienten 3. Ordnung - 'K', eq_.K, ... % Samples pro Symbol - 'DCmu', eq_.DCmu, ... % Anpassungsrate für DC-Tap - 'ideal_dfe', eq_.ideal_dfe, ... % Flag für ideal DFE (0 oder 1) - 'training_length', eq_.training_length, ... % Anzahl Trainingssymbole - 'training_loops', eq_.training_loops, ... % Anzahl Trainingsdurchläufe - 'TRmu1', eq_.FFEmu, ... % mu für DD-Modus (1. Ordnung) - 'TRmu2', eq_.FFEmu, ... % mu für DD-Modus (2. Ordnung) - 'TRmu3', eq_.FFEmu, ... % mu für DD-Modus (3. Ordnung) - 'TRmuDFE', eq_.DFEmu, ... % mu für DFE-Modus im DD - 'dd_loops', 5, ... % Anzahl Durchläufe im DD-Modus - 'DDmu1', eq_.DDmu(1), ... % mu für DD-Modus (1. Ordnung) - 'DDmu2', eq_.DDmu(2), ... % mu für DD-Modus (2. Ordnung) - 'DDmu3', eq_.DDmu(3), ... % mu für DD-Modus (3. Ordnung) - 'DDmuDFE', eq_.DDmu(4), ... % mu für DFE-Modus im DD - 'MLSE_mode', 'viterbi', ... % Beispiel: MLSE-Modus als String - 'MLSE_trellis_states', jsonencode(mlse_.trellis_states), ... % Trellis-States, z.B. als JSON-String oder kommasepariert - 'comment', 'function: vnle_postfilter_mlse', ... % Zusätzliche Kommentare - 'config_hash', NaN ... - ); -eq_package.resultsVNLE = resultsVNLE; -eq_package.resultsMLSE = resultsMLSE; -eq_package.equalizerConfigVNLE = equalizerConfigVNLE; -eq_package.equalizerConfigMLSE = equalizerConfigMLSE; +mlse_results = struct(); +mlse_results.metrics = Metricstruct; +mlse_results.metrics.result_id = NaN; +mlse_results.metrics.run_id = NaN; +mlse_results.metrics.eqParam_id = NaN; +mlse_results.metrics.date_of_processing = datetime('now'); +mlse_results.metrics.BER = bers.mlse; +mlse_results.metrics.numBits = numbits.mlse; +mlse_results.metrics.numBitErr = errors.mlse; +mlse_results.metrics.BER_precoded = bers.mlse_precoded; +mlse_results.metrics.numBitErr_precoded = errors.mlse_precoded; +mlse_results.metrics.SNR = NaN; +mlse_results.metrics.SNR_level = NaN; +mlse_results.metrics.GMI = NaN; +mlse_results.metrics.AIR = NaN; +mlse_results.metrics.EVM = NaN; +mlse_results.metrics.EVM_level = NaN; +mlse_results.metrics.Alpha = alpha; +mlse_results.metrics.MLSE_dir = mlse_.DIR; +% Create MLSE results structure +mlse_results.config = Equalizerstruct(); +mlse_results.config.eq = jsonencode(eq_); +mlse_.DIR = []; +mlse_results.config.mlse = jsonencode(mlse_); +mlse_results.config.equalizer_structure = int32(equalizer_structure.vnle_pf_mlse); +mlse_results.config.comment = 'function: vnle_postfilter_mlse - MLSE part'; -% eq_package.vnle_out = eq_signal_sd; +%% Display analysis if requested if options.showAnalysis - - - % fprintf(['VNLE EVM lvl: ',repmat('%.3f ',1,numel(evm_lvl)),' \n'],evm_lvl); - - % fprintf('VNLE BER: %.2e \n',ber_vnle); - % - % fprintf('MLSE BER: %.2e \n',ber_mlse); - - figure(336); - showEQNoisePSD(eq_noise,"fignum",336,"displayname",'Residual Noise after VNLE','postfilter_taps',pf_.coefficients); - - % figure(337);clf; - % rx_signal.spectrum("normalizeTo0dB",1,"fignum",337,"displayname",'Rx Signal'); - - % figure(338);clf; - % showEQcoefficients('n1',eq_.e,'n2',eq_.e2,'n3',eq_.e3,"displayname",'Coefficients','fignum',338); - - if ~isempty(options.postFFE) - showEQcoefficients('n1',options.postFFE.e,"displayname",'Coefficients','fignum',338); - end - - showEQfilter(eq_.e,eq_signal_sd.fs.*2); - - figure(340);clf; - eq_signal_sd.eye(eq_signal_sd.fs,M,"fignum",340); - - figure(341);clf; - showLevelHistogram(eq_signal_sd,tx_symbols,"fignum",341); - - warning off - showLevelScatter(eq_signal_sd,tx_symbols,"fignum",400); - - showLevelScatter(rx_signal.resample("fs_out",tx_symbols.fs),tx_symbols,"fignum",401); - warning on - drawnow; - - % autoArrangeFigures(3,3,2) - + displayAnalysis(eq_noise, eq_signal_sd, rx_signal, eq_, pf_, mlse_, tx_symbols, M, options.postFFE); +end end +%% Helper Functions +function [numbits, errors, ber, error_locations] = calculateBER(eq_signal_hd, mlse_sig_hd, tx_symbols, tx_bits, precode_mode, M, eth_style) +% Initialize output structure +numbits = struct('vnle', 0, 'mlse', 0); +errors = struct('vnle', 0, 'mlse', 0, 'vnle_precoded', 0, 'mlse_precoded', 0); +ber = struct('vnle', 0, 'mlse', 0, 'vnle_precoded', 0, 'mlse_precoded', 0); +error_locations = struct('vnle', [], 'mlse', []); + +% PAM mapper for demapping +mapper = PAMmapper(M, 0, "eth_style", eth_style); + +switch precode_mode + case db_mode.no_db + % TX Data is not precoded + % A) Emulate diff precoding + eq_signal_hd_precoded = Duobinary().encode(eq_signal_hd, "M", M); + eq_signal_hd_precoded = Duobinary().decode(eq_signal_hd_precoded, "M", M); + + mlse_sig_hd_precoded = Duobinary().encode(mlse_sig_hd, "M", M); + mlse_sig_hd_precoded = Duobinary().decode(mlse_sig_hd_precoded, "M", M); + + tx_symbols_precoded = Duobinary().encode(tx_symbols); + tx_symbols_precoded = Duobinary().decode(tx_symbols_precoded); + + tx_bits_precoded = mapper.demap(tx_symbols_precoded); + + rx_bits_vnle = mapper.demap(eq_signal_hd_precoded); + [~, errors.vnle_precoded, ber.vnle_precoded, ~] = calc_ber(rx_bits_vnle.signal, tx_bits_precoded.signal, "skip_front", 100, "skip_end", 150, "returnErrorLocation", 1); + + rx_bits_mlse = mapper.demap(mlse_sig_hd_precoded); + [~, errors.mlse_precoded, ber.mlse_precoded, ~] = calc_ber(rx_bits_mlse.signal, tx_bits_precoded.signal, "skip_front", 100, "skip_end", 150, "returnErrorLocation", 1); + + % B) Just determine BER + rx_bits_vnle = mapper.demap(eq_signal_hd); + [numbits.vnle, errors.vnle, ber.vnle, error_locations.vnle] = calc_ber(rx_bits_vnle.signal, tx_bits.signal, "skip_front", 100, "skip_end", 150, "returnErrorLocation", 1); + + rx_bits_mlse = mapper.demap(mlse_sig_hd); + [numbits.mlse, errors.mlse, ber.mlse, error_locations.mlse] = calc_ber(rx_bits_mlse.signal, tx_bits.signal, "skip_front", 100, "skip_end", 150, "returnErrorLocation", 1); + + case db_mode.db_precoded + % Data is precoded on TX side + % A) Decode at Rx if no DB targeting was applied + eq_signal_hd_decoded = Duobinary().encode(eq_signal_hd, "M", M); + eq_signal_hd_decoded = Duobinary().decode(eq_signal_hd_decoded, "M", M); + rx_bits_vnle_decoded = mapper.demap(eq_signal_hd_decoded); + [~, errors.vnle_precoded, ber.vnle_precoded, ~] = calc_ber(rx_bits_vnle_decoded.signal, tx_bits.signal, "skip_front", 100, "skip_end", 150, "returnErrorLocation", 1); + + mlse_sig_hd_decoded = Duobinary().encode(mlse_sig_hd, "M", M); + mlse_sig_hd_decoded = Duobinary().decode(mlse_sig_hd_decoded, "M", M); + rx_bits_mlse_decoded = mapper.demap(mlse_sig_hd_decoded); + [~, errors.mlse_precoded, ber.mlse_precoded, ~] = calc_ber(rx_bits_mlse_decoded.signal, tx_bits.signal, "skip_front", 100, "skip_end", 150, "returnErrorLocation", 1); + + % B) Omit the Coding by comparing with demapped TX symbol sequence + tx_bits_demapped = mapper.demap(tx_symbols); + rx_bits_vnle = mapper.demap(eq_signal_hd); + [numbits.vnle, errors.vnle, ber.vnle, error_locations.vnle] = calc_ber(rx_bits_vnle.signal, tx_bits_demapped.signal, "skip_front", 100, "skip_end", 150, "returnErrorLocation", 1); + + rx_bits_mlse = mapper.demap(mlse_sig_hd); + [numbits.mlse, errors.mlse, ber.mlse, error_locations.mlse] = calc_ber(rx_bits_mlse.signal, tx_bits_demapped.signal, "skip_front", 100, "skip_end", 150, "returnErrorLocation", 1); +end +end + + +function displayAnalysis(eq_noise, eq_signal_sd, rx_signal, eq_, pf_, mlse_, tx_symbols, M, postFFE) +% Display analysis plots and metrics +figure(336); +showEQNoisePSD(eq_noise, "fignum", 336, "displayname", 'Residual Noise after VNLE', 'postfilter_taps', pf_.coefficients); + +if ~isempty(postFFE) + showEQcoefficients('n1', postFFE.e, "displayname", 'Coefficients', 'fignum', 338); +end + +showEQfilter(eq_.e, eq_signal_sd.fs.*2); + +figure(340); clf; +eq_signal_sd.eye(eq_signal_sd.fs, M, "fignum", 340); + +figure(341); clf; +showLevelHistogram(eq_signal_sd, tx_symbols, "fignum", 341); + +warning off +showLevelScatter(eq_signal_sd, tx_symbols, "fignum", 400); +showLevelScatter(rx_signal.resample("fs_out", tx_symbols.fs), tx_symbols, "fignum", 401); +drawnow; +warning on end \ No newline at end of file diff --git a/Functions/EQ_visuals/showLevelHistogram.m b/Functions/EQ_visuals/showLevelHistogram.m index 66e5d08..437b65c 100644 --- a/Functions/EQ_visuals/showLevelHistogram.m +++ b/Functions/EQ_visuals/showLevelHistogram.m @@ -40,7 +40,7 @@ end cnt(lvl) = round(numel(intermediate(~isnan(intermediate)))./length(eq_signal),3).*100; hold on warning off - histogram(received_sd(lvl,:),1000,"EdgeAlpha",0,'DisplayName',['Lvl ',num2str(lvl),' | ',num2str(cnt(lvl)),' %'],'FaceColor',lvlcol(lvl,:),'Normalization','pdf'); + histogram(received_sd(lvl,:),1000,"EdgeAlpha",0,'DisplayName',['Lvl ',num2str(lvl),' ; ',num2str(cnt(lvl)),' '],'FaceColor',lvlcol(lvl,:),'Normalization','pdf'); warning on end legend diff --git a/Functions/EQ_visuals/showLevelScatter.m b/Functions/EQ_visuals/showLevelScatter.m index 1bc7e33..a5d34ef 100644 --- a/Functions/EQ_visuals/showLevelScatter.m +++ b/Functions/EQ_visuals/showLevelScatter.m @@ -1,4 +1,4 @@ -function showLevelScatter(eq_signal,ref_symbols,options) +function [symbols_for_lvl,avg_for_lvl] = showLevelScatter(eq_signal,ref_symbols,options) arguments eq_signal ref_symbols @@ -7,22 +7,27 @@ arguments options.f_sym = []; end +plot_shit = 0; + if isa(eq_signal,'Signal') options.f_sym = eq_signal.fs; eq_signal = eq_signal.signal; + assert(~isempty(options.f_sym),'No fsym given'); end if isa(ref_symbols,'Signal') ref_symbols = ref_symbols.signal; end -% Determine the figure number to use or create a new figure -if isnan(options.fignum) - fig = figure; % Create a new figure and get its handle -else - fig = figure(options.fignum); % Use the specified figure number + +if plot_shit + % Determine the figure number to use or create a new figure + if isnan(options.fignum) + fig = figure; % Create a new figure and get its handle + else + fig = figure(options.fignum); % Use the specified figure number + end end -assert(~isempty(options.f_sym),'No fsym given'); rx_symbols = eq_signal ./ rms(eq_signal); correct_symbols = ref_symbols; @@ -32,57 +37,54 @@ col = cbrewer2('Paired',numel(unique(correct_symbols))*2); ccnt = -1; levels = unique(correct_symbols); - +symbols_for_lvl = NaN(numel(levels),length(correct_symbols)); start = 1; ende = length(correct_symbols); -% start = 30000; -% ende = 40000; - for l = 1:numel(levels) ccnt = ccnt+2; level_amplitude = levels(l); - symbols_for_lvl = NaN(1,length(correct_symbols)); - symbols_for_lvl(correct_symbols==level_amplitude) = rx_symbols(correct_symbols==level_amplitude); - std_lvl(l) = std(symbols_for_lvl,'omitnan'); + + symbols_for_lvl(l,correct_symbols==level_amplitude) = rx_symbols(correct_symbols==level_amplitude); + std_lvl(l) = std(symbols_for_lvl(l,:),'omitnan'); xax_in_sec = ((1:length(correct_symbols)) / f_sym) * 1e6; - % xax_in_sec = 1:length(correct_symbols); -if 0 - scatter(xax_in_sec(start:ende),symbols_for_lvl(start:ende),10,'.','MarkerFaceAlpha',0.5,'MarkerEdgeAlpha',0.5,'MarkerEdgeColor',col(ccnt,:)); - hold on; -end + + if plot_shit + scatter(xax_in_sec(start:ende),symbols_for_lvl(l,start:ende),10,'.','MarkerFaceAlpha',0.5,'MarkerEdgeAlpha',0.5,'MarkerEdgeColor',col(ccnt,:)); + hold on; + end + end std_lvl = round(std_lvl,2); ccnt = 0; - +avg_for_lvl = NaN(numel(levels),length(correct_symbols)); % Add the windowed/ smoothed curves for l = 1:numel(levels) ccnt = ccnt+2; level_amplitude = levels(l); - symbols_for_lvl = NaN(1,length(correct_symbols)); - movmean = 1/250 .* movsum(rx_symbols(correct_symbols==level_amplitude),[250/2,250/2], 'Endpoints', 'fill'); - symbols_for_lvl(correct_symbols==level_amplitude) = movmean; + avg_for_lvl(l,correct_symbols==level_amplitude) = movmean; - nanx = isnan(symbols_for_lvl); - t = 1:numel(symbols_for_lvl); - symbols_for_lvl(nanx) = interp1(t(~nanx), symbols_for_lvl(~nanx), t(nanx)); + nanx = isnan(avg_for_lvl(l,:)); + t = 1:numel(avg_for_lvl(l,:)); + avg_for_lvl(l,nanx) = interp1(t(~nanx), avg_for_lvl(l,~nanx), t(nanx)); xax_in_sec = ((1:length(correct_symbols)) / f_sym) * 1e6; % xax_in_sec = 1:length(correct_symbols); - plot(xax_in_sec(start:ende),symbols_for_lvl(start:ende),'Color',col(ccnt,:)); + if plot_shit + plot(xax_in_sec(start:ende),avg_for_lvl(l,start:ende),'Color',col(ccnt,:)); + end hold on end -%yline(max(rx_symbols(correct_symbols==levels(2)))) -yline(levels); + if 0 annotation(fig,'textbox',... @@ -121,9 +123,9 @@ if 0 'FitBoxToText','off'); end -% xlim([0, 2.6]) -% ylim([-2 2]) +if plot_shit +yline(levels); xlabel('Time in $\mu$s'); ylabel('Normalized Amplitude'); - +end end diff --git a/Functions/Job_Processing/configureEqualizers_remove.m b/Functions/Job_Processing/configureEqualizers_remove.m new file mode 100644 index 0000000..f207088 --- /dev/null +++ b/Functions/Job_Processing/configureEqualizers_remove.m @@ -0,0 +1,40 @@ +function [eq_, pf_, mlse_, mlse_db_, eq_post] = configureEqualizers(M, len_tr, vnle_order, dfe_order, mu_dc, mu_ffe, mu_dfe, pf_ncoeffs) + % CONFIGUREEQUALIZERS Creates and configures equalizer objects + % + % Inputs: + % M - PAM level + % len_tr - Training length + % vnle_order - Array with orders for VNLE [order1, order2, order3] + % dfe_order - Array with orders for DFE + % mu_dc - DC adaptation rate + % mu_ffe - Array with adaptation rates for FFE [mu1, mu2, mu3] + % mu_dfe - Adaptation rate for DFE + % pf_ncoeffs - Number of coefficients for postfilter + % + % Outputs: + % eq_ - Configured EQ object + % pf_ - Configured Postfilter object + % mlse_ - Configured MLSE_viterbi object + % mlse_db_ - Configured MLSE_viterbi object for duobinary + % eq_post - Configured FFE object for post-processing + + % Configure main equalizer + eq_ = EQ("Ne", vnle_order, "Nb", dfe_order, ... + "training_length", len_tr, "training_loops", 5, "dd_loops", 5, ... + "K", 2, "DCmu", mu_dc, "DDmu", [mu_ffe mu_dfe], ... + "DFEmu", 0.005, "FFEmu", 0, "plotfinal", 0, "ideal_dfe", 1); + + % Configure postfilter + pf_ = Postfilter("ncoeff", pf_ncoeffs, "useBurg", 1); + + % Configure MLSE objects + mlse_ = MLSE_viterbi("duobinary_output", 0, 'M', M, ... + 'trellis_states', PAMmapper(M,0).levels); + mlse_db_ = MLSE_viterbi("DIR", [1,1], "duobinary_output", 0, ... + "M", M, "trellis_states", PAMmapper(M,0).levels); + + % Configure post-FFE + eq_post = FFE("epochs_tr", 5, "epochs_dd", 5, "len_tr", 4096*2, ... + "mu_dd", 1e-4, "mu_tr", 0, "order", 2001, ... + "sps", 1, "decide", 0); + end \ No newline at end of file diff --git a/Functions/Job_Processing/loadSignalData.m b/Functions/Job_Processing/loadSignalData.m new file mode 100644 index 0000000..cff4ff8 --- /dev/null +++ b/Functions/Job_Processing/loadSignalData.m @@ -0,0 +1,106 @@ +function [Bits, Symbols, Scpe_cell, found_sync] = loadSignalData(dataTable, options) +% LOADSIGNALDATA Loads and synchronizes signal data from storage +% +% Inputs: +% dataTable - Table with file paths and configuration +% options - Struct with storage_path and max_occurences +% +% Outputs: +% Symbols_mapped - Mapped symbols from bits +% Symbols - Original symbols +% Scpe_cell - Cell array of synchronized signals +% found_sync - Boolean indicating if synchronization was successful + +found_sync = 0; +tempLocalStorage = 1; + +% Part A: Check and load from local storage if available +if tempLocalStorage == 1 + local_filename = sprintf('sync_data_run_%s.mat', num2str(dataTable.run_id)); + if exist(local_filename, 'file') + % Load from local storage and return + load(local_filename, 'Bits', 'Symbols', 'Scpe_cell'); + found_sync = 1; + end +end + +% If not locally saved, load from storage +if ~found_sync + % Load transmitted bits + Bits = load([options.storage_path, char(dataTable.tx_bits_path)]); + Bits = Bits.Bits; + + % Map bits to symbols + M = double(dataTable.pam_level); + fsym = dataTable.symbolrate; + Symbols_mapped = PAMmapper(M,0).map(Bits); + Symbols_mapped.fs = fsym; + + % Load original symbols + Symbols = load([options.storage_path, char(dataTable.tx_symbols_path)]); + Symbols = Symbols.Symbols; + + found_sync = 0; + Scpe_cell = {}; + + % Try to load pre-synchronized data + try + Scpe_load = load([options.storage_path, char(dataTable.rx_sync_path)]); + Scpe_cell = Scpe_load.S; + [~,~,~,found_sync] = Scpe_cell{2}.tsynch("reference", Symbols, "fs_ref", fsym, "debug_plots", 1); + catch + % Continue to next method if this fails + end +end + +% If not found, try with raw data +if ~found_sync + try + Scpe_sig_raw = load([options.storage_path, char(dataTable.rx_raw_path(1))]); + Scpe_sig_raw = Scpe_sig_raw.Scpe_sig_raw; + Scpe_sig_resampled = Scpe_sig_raw.resample("fs_in", Scpe_sig_raw.fs, "fs_out", 2*fsym); + [~, Scpe_cell, ~, found_sync] = Scpe_sig_resampled.tsynch("reference", Symbols, "fs_ref", fsym, "debug_plots", 1); + catch + % Continue to next method if this fails + end +end + +% Last attempt with mapped symbols +if ~found_sync && exist('Scpe_sig_raw', 'var') + if length(Symbols_mapped.signal) ~= sum(Symbols_mapped.signal == Symbols.signal) + [~, Scpe_cell, ~, found_sync] = Scpe_sig_raw.tsynch("reference", Symbols_mapped, "fs_ref", fsym, "debug_plots", 0); + end +end + +% Part B: Save to local storage if data was loaded and synced +if tempLocalStorage == 1 && found_sync + local_filename = sprintf('sync_data_run_%s.mat', num2str(dataTable.run_id)); + + % Create directory if it doesn't exist + if ~exist('temp_sync_data', 'dir') + mkdir('temp_sync_data'); + end + + % List existing files and remove oldest if more than N + max_local_files = 10; % Store up to N files + files = dir('sync_data_run_*.mat'); + if length(files) >= max_local_files + % Sort by date + [~, idx] = sort([files.datenum]); + % Delete oldest file + delete(files(idx(1)).name); + end + + % Save current data + save(local_filename, 'Bits', 'Symbols', 'Scpe_cell'); +end + +% Limit number of occurrences +if found_sync + record_realizations = min(options.max_occurences, length(Scpe_cell)); + Scpe_cell = Scpe_cell(1:record_realizations); +else + warning('Could not synchronize the received signal with the stored symbols!'); +end + +end \ No newline at end of file diff --git a/Functions/Job_Processing/preprocessSignal.m b/Functions/Job_Processing/preprocessSignal.m new file mode 100644 index 0000000..4828bba --- /dev/null +++ b/Functions/Job_Processing/preprocessSignal.m @@ -0,0 +1,26 @@ +function Scpe_sig = preprocessSignal(Scpe_sig, Symbols, fsym) +% PREPROCESSSIGNAL Performs standard preprocessing on a signal +% +% Inputs: +% Scpe_sig - Input signal +% Symbols - Reference symbols for synchronization +% fsym - Symbol frequency +% +% Outputs: +% Scpe_sig - Preprocessed signal + +% Resample to 2x symbol rate +Scpe_sig = Scpe_sig.resample("fs_out", 2*fsym); + +% Synchronize with reference +[Scpe_sig, ~] = Scpe_sig.tsynch("reference", Symbols, "fs_ref", fsym, "debug_plots", 0); + +% Apply Gaussian filter +Scpe_sig = Filter('filtdegree', 4, "f_cutoff", Symbols.fs.*0.6, ... + "fs", Scpe_sig.fs, "filterType", filtertypes.gaussian, ... + "active", true).process(Scpe_sig); + +% Remove DC offset +Scpe_sig = Scpe_sig - mean(Scpe_sig.signal); + +end \ No newline at end of file diff --git a/Functions/Job_Processing/printResults.m b/Functions/Job_Processing/printResults.m new file mode 100644 index 0000000..4cd1ddc --- /dev/null +++ b/Functions/Job_Processing/printResults.m @@ -0,0 +1,47 @@ +function printResults(run_id, M, Symbols, vnle_pf_package, dbtgt_package, occ) + % PRINTRESULTS Prints formatted results from equalization + % + % Inputs: + % run_id - Run identifier + % M - PAM level + % Symbols - Symbol data with fs property + % vnle_pf_package - VNLE+PF results package + % dbtgt_package - DB target results package (optional) + % occ - Occurrence index + + % Print header + fprintf("==== EQUALIZATION RUN-ID %d | PAM-%d | %.2f GBd ====\n\n", run_id, M, Symbols.fs.*1e-9); + + % VNLE Results + if ~isempty(vnle_pf_package) + vnle_result = vnle_pf_package{occ}.resultsVNLE; + mlse_result = vnle_pf_package{occ}.resultsMLSE; + + fprintf(">> VNLE Results:\n"); + fprintf(" BER: %.2e\n", vnle_result.BER); + fprintf(" BER (pre-code): %.2e\n", vnle_result.BER_precoded); + fprintf(" SNR: %.2f dB\n", vnle_result.SNR); + fprintf(" GMI: %.4f\n", vnle_result.GMI); + fprintf(" Linerate: %.2f Gbps\n", Symbols.fs .* floor(log2(M)*10)/10 .*1e-9); + fprintf(" AIR: %.2f Gbps\n", vnle_result.AIR.*1e-9); + fprintf("\n"); + + % MLSE Results + fprintf(">> MLSE Results:\n"); + fprintf(" BER: %.2e\n", mlse_result.BER); + fprintf(" BER (pre-code): %.2e\n", mlse_result.BER_precoded); + fprintf(" Channel Alpha: %.2f\n", mlse_result.Alpha); + fprintf("\n"); + end + + % DB Target Results + if ~isempty(dbtgt_package) + dbtgt = dbtgt_package{occ}.resultsDBtgt; + fprintf(">> DB Target Results:\n"); + fprintf(" BER: %.2e\n", dbtgt.BER); + fprintf(" BER (pre-code): %.2e\n", dbtgt.BER_precoded); + fprintf("\n"); + end + + fprintf("- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - \n\n"); + end \ No newline at end of file diff --git a/Functions/Job_Processing/queryRunid.m b/Functions/Job_Processing/queryRunid.m new file mode 100644 index 0000000..7a0e048 --- /dev/null +++ b/Functions/Job_Processing/queryRunid.m @@ -0,0 +1,11 @@ +function dataTable = queryRunid(run_id, database) + + % Query database for run configuration + filterParams = database.tables; + filterParams.Runs = struct('run_id', run_id); + + [dataTable, ~] = database.queryDB(filterParams, database.getTableFieldNames('Runs')); + [~, uniqueIdx] = unique(dataTable.run_id); % Get unique run_id indices + dataTable = dataTable(uniqueIdx, :); % Extract unique configurations for each run_id + +end \ No newline at end of file diff --git a/Functions/Job_Processing/submitJobs.m b/Functions/Job_Processing/submitJobs.m new file mode 100644 index 0000000..a5a45a4 --- /dev/null +++ b/Functions/Job_Processing/submitJobs.m @@ -0,0 +1,213 @@ +function [results,wh] = submitJobs(run_ids, dsp_options, submit_mode, submit_options) +% SUBMITJOBS Submits dsp_runid jobs for processing +% +% Inputs: +% run_ids - Single run ID or array of Run IDs for processing +% dsp_options - Options for dsp_runid function +% submit_mode - Execution mode: 'parallel' or 'linear' +% options - Optional parameters +% +% Outputs: +% results - Cell array of results from each job +% wh - Updated DataStorage object + +% USAGE: +% % === SETTINGS === +% +% dsp_options.append_to_db = 1; +% dsp_options.max_occurences = 15; +% dsp_options.database_path = 'C:\Users\Silas\Documents\MATLAB\Datensätze\sioe_labor\'; +% dsp_options.database_name = 'silas_labor_newdsp_newstructure.db'; +% dsp_options.storage_path = 'Z:\2024\sioe_labor\'; +% dsp_options.parameters = struct(); +% dsp_options.parameters.mu_dc = [0.005]; +% +% % === Get Run ID's === +% db = DBHandler("pathToDB", [dsp_options.database_path, dsp_options.database_name], "type", "sqlite"); +% fp = QueryFilter(); +% % fp.where('Runs', 'run_id','EQUALS', 5108); +% fp.where('Runs', 'is_mpi','EQUALS', 0); +% fp.where('Runs', 'fiber_length','EQUALS', 1); +% fp.where('Runs', 'wavelength','EQUALS', 1310); +% fp.where('Runs', 'db_mode','EQUALS', 0); +% fp.where('Runs', 'rop_attenuation','EQUALS', 0); +% fp.where('Runs', 'pam_level','EQUALS', 4); +% fp.where('Runs', 'bitrate','EQUALS', 360e9); +% fp.where('Runs', 'power_pd_in','GREATER_THAN', 7); +% [dataTable,~] = db.queryDB(fp, db.getTableFieldNames('Runs')); + +% % === Initialize DataStorage === +% wh = DataStorage(dsp_options.parameters); +% wh.addStorage("ffe_package"); +% wh.addStorage("mlse_package"); +% wh.addStorage("vnle_package"); +% wh.addStorage("dbtgt_package"); +% wh.addStorage("dbenc_package"); +% +% % === RUN IT === +% [results,wh] = submitJobs(dataTable.run_id(:), dsp_options, "parallel", 'wh', wh, 'waitbar', true); +% wh.getStoValue('ffe_package',0.005); +% wh.getStoValue('mlse_package',0.005); + + +arguments + run_ids int32 = 0 + dsp_options struct = struct(); + submit_mode processingMode = processingMode.serial + submit_options.waitbar (1,1) logical = true + submit_options.wh = DataStorage(struct()); % Optional DataStorage object +end + +% Ensure run_ids is a row vector +run_ids = run_ids(:)'; + +% Get number of jobs per run_id +nJobsPerRunId = submit_options.wh.getLastLinIndice; +nRunIds = length(run_ids); +totalJobs = nJobsPerRunId * nRunIds; + +% Initialize results array for all jobs +results = cell(nJobsPerRunId, nRunIds); +futures = cell(totalJobs, 1); + +% Initialize waitbar if requested +if submit_options.waitbar + h = waitbar(0, 'Processing Jobs...'); + cleanupObj = onCleanup(@() delete(h)); +end + +% Process jobs based on mode +if submit_mode == processingMode.parallel + setupParallelPool(); + + % Submit jobs to parallel pool + jobCounter = 0; + for r = 1:nRunIds + for k = 1:nJobsPerRunId + jobCounter = jobCounter + 1; + optionalVars = buildOptionalVars(k, submit_options.wh); + + % Submit job with proper arguments + futures{jobCounter} = parfeval(@dsp_runid, 1, run_ids(r), ... + "database_name", dsp_options.database_name, ... + "append_to_db", dsp_options.append_to_db, ... + "database_path", dsp_options.database_path, ... + "max_occurences", dsp_options.max_occurences, ... + "storage_path", dsp_options.storage_path, ... + "parameters", optionalVars); + + fprintf('[RunID %d, Job %d] Submitted to pool.\n', run_ids(r), k); + end + end + + if submit_options.waitbar + setupParallelWaitbar(futures, totalJobs, h); + end + + % Fetch results + jobCounter = 0; + for r = 1:nRunIds + for k = 1:nJobsPerRunId + jobCounter = jobCounter + 1; + try + results{k,r} = fetchOutputs(futures{jobCounter}); + fprintf('[RunID %d, Job %d] Completed successfully.\n', run_ids(r), k); + storeResult(results{k,r}, k, submit_options.wh); + catch ME + handleError(ME, k, run_ids(r)); + results{k,r} = ME; + end + end + end + +elseif submit_mode == processingMode.serial + % Process jobs sequentially + jobCounter = 0; + for r = 1:nRunIds + for k = 1:nJobsPerRunId + jobCounter = jobCounter + 1; + optionalVars = buildOptionalVars(k, submit_options.wh); + + try + fprintf('[RunID %d, Job %d] Running in linear mode...\n', run_ids(r), k); + + results{k,r} = dsp_runid(run_ids(r),... + "database_name", dsp_options.database_name,... + "append_to_db", dsp_options.append_to_db,... + "database_path", dsp_options.database_path,... + "max_occurences", dsp_options.max_occurences,... + "storage_path", dsp_options.storage_path,... + "parameters", optionalVars); + + fprintf('[RunID %d, Job %d] Completed successfully.\n', run_ids(r), k); + storeResult(results{k,r}, k, submit_options.wh); + catch ME + handleError(ME, k, run_ids(r)); + results{k,r} = ME; + end + + % Update waitbar if requested + if submit_options.waitbar + waitbar(jobCounter/totalJobs, h, sprintf('Completed %d/%d jobs', jobCounter, totalJobs)); + end + end + end +else + error('') +end + +wh = submit_options.wh; + +% Helper functions remain the same except for handleError +% Modified handleError function to include run_id + function handleError(ME, jobIndex, run_id) + fprintf('[RunID %d, Job %d] ERROR [%s]: %s\n', run_id, jobIndex, ME.identifier, ME.message); + for st = ME.stack' + fprintf(' %s:%d (%s)\n', st.file, st.line, st.name); + end + fprintf('Full report:\n%s\n', getReport(ME,'extended')); + end + +% Other helper functions remain unchanged + function setupParallelPool() + curpool = gcp('nocreate'); + if isempty(curpool) + parpool; + else + if ~isempty(curpool.FevalQueue.QueuedFutures) || ~isempty(curpool.FevalQueue.RunningFutures) + oldq = length(curpool.FevalQueue.QueuedFutures) + length(curpool.FevalQueue.RunningFutures); + curpool.FevalQueue.cancelAll; + fprintf('Canceled %d unfetched jobs from old queue.\n', oldq); + end + end + end + + function optionalVars = buildOptionalVars(jobIndex, wh) + optionalVars = struct(); + if ~isempty(wh.getDimension) + [parametervalues, parameternames] = wh.getPhysIndicesByLinIndex(jobIndex); + for pidx = 1:numel(parameternames) + optionalVars.(parameternames{pidx}) = parametervalues{pidx}; + end + end + end + + function setupParallelWaitbar(futures, totalJobs, h) + updateWaitbar = @(~) waitbar(sum(cellfun(@(f) strcmp(f.State, 'finished'), futures))/totalJobs, h, ... + sprintf('Completed %d/%d jobs', sum(cellfun(@(f) strcmp(f.State, 'finished'), futures)), totalJobs)); + + futureArray = [futures{:}]; + afterEach(futureArray, updateWaitbar, 0); + end + + function storeResult(result, jobIndex, wh) + if ~isempty(wh) + wh.addValueToStorageByLinIdx(result.ffe_package, 'ffe_package', jobIndex); + wh.addValueToStorageByLinIdx(result.mlse_package, 'mlse_package', jobIndex); + wh.addValueToStorageByLinIdx(result.vnle_package, 'vnle_package', jobIndex); + wh.addValueToStorageByLinIdx(result.dbtgt_package, 'dbtgt_package', jobIndex); + wh.addValueToStorageByLinIdx(result.dbenc_package, 'dbenc_package', jobIndex); + end + end + +end \ No newline at end of file diff --git a/Functions/Metrics/calc_evm.m b/Functions/Metrics/calc_evm.m index c597732..a8e4b7e 100644 --- a/Functions/Metrics/calc_evm.m +++ b/Functions/Metrics/calc_evm.m @@ -28,39 +28,39 @@ end [evm_total,evm_lvl] = calc_evm_(test_signal,reference_signal); -function [evm_total,evm_lvl] = calc_evm_(test_signal,reference_signal) + function [evm_total, evm_lvl] = calc_evm_(test_signal, reference_signal) + % Validate input + assert(length(test_signal) == length(reference_signal), "Sequence length does not match"); - assert(length(test_signal) == length(reference_signal),"Sequence length does not match"); + % Calculate error vector + error_vector = test_signal - reference_signal; - error_vector = (test_signal-reference_signal); + % EVM (RMS) as percentage, per MathWorks definition + evm_total = sqrt(sum(abs(error_vector).^2) / sum(abs(reference_signal).^2)) * 100; - %%% Overall EVM - evm_total = rms(error_vector); - - try - %%% Per Level EVM + % Per-level EVM k = unique(reference_signal); + evm_lvl = NaN(1, length(k)); for lvl = 1:length(k) - lvl_errors = error_vector(reference_signal==k(lvl)); - evm_lvl(lvl) = rms(lvl_errors); + idx = reference_signal == k(lvl); + if any(idx) + lvl_errors = error_vector(idx); + lvl_refs = reference_signal(idx); + evm_lvl(lvl) = sqrt(sum(abs(lvl_errors).^2) / sum(abs(lvl_refs).^2)) * 100; + end end - catch - evm_lvl = NaN; - warning('No EVM per level calculated') end -end + function [data_,reference_]=trimseq(data,reference,skipstart,skip_end) -function [data_,reference_]=trimseq(data,reference,skipstart,skip_end) + data_ = data(skipstart+1:end-skip_end,:); - data_ = data(skipstart+1:end-skip_end,:); + delta_bits = length(reference) - length(data); - delta_bits = length(reference) - length(data); + skip_end = delta_bits + skip_end; - skip_end = delta_bits + skip_end; + reference_ = reference(skipstart+1:end-skip_end,:); - reference_ = reference(skipstart+1:end-skip_end,:); - -end + end end \ No newline at end of file diff --git a/Functions/Metrics/calc_std.m b/Functions/Metrics/calc_std.m index 5f3a2e5..927b18d 100644 --- a/Functions/Metrics/calc_std.m +++ b/Functions/Metrics/calc_std.m @@ -25,34 +25,46 @@ end [test_signal,reference_signal]=trimseq(test_signal,reference_signal,options.skip_front,options.skip_end); % CALC EVM -[std_total,std_lvl] = calc_std_(test_signal,reference_signal); +[std_total, std_lvl, nsd_lvl, d_min]= calc_std_(test_signal,reference_signal); +std_lvl = nsd_lvl; + function [std_total, std_lvl, nsd_lvl, d_min] = calc_std_(test_signal, reference_signal) + + % NSD (Normalized Standard Deviation) expresses the noise spread at each symbol level + % relative to the minimum distance between levels. A low NSD means low error risk, + % while NSD approaching 0.5 indicates a high chance of symbol errors due to noise. -function [std_total,std_lvl] = calc_std_(test_signal,reference_signal) - - assert(length(test_signal) == length(reference_signal),"Sequence length does not match"); - - error_vector = (test_signal-reference_signal); - - %%% Overall EVM - std_total = std(test_signal); - - test_signal = test_signal ./ rms(test_signal); - - try - %%% Per Level EVM + % Ensure input is column vector + test_signal = test_signal(:); + reference_signal = reference_signal(:); + + assert(length(test_signal) == length(reference_signal), "Sequence length does not match"); + + % Find unique levels and their minimum spacing k = unique(reference_signal); + d_min = min(diff(k)); % Minimum distance between adjacent levels + + % Normalize test signal to RMS=1 (if not already) + test_signal = test_signal / rms(test_signal); + + % Overall standard deviation (not normalized) + std_total = std(test_signal); + + % Per-level standard deviation and normalized std + std_lvl = zeros(1, length(k)); + nsd_lvl = zeros(1, length(k)); for lvl = 1:length(k) - % lvl_errors = error_vector(reference_signal==k(lvl)); - std_lvl(lvl) = std(test_signal(reference_signal==k(lvl))); + idx = reference_signal == k(lvl); + if any(idx) + std_lvl(lvl) = std(test_signal(idx)); + nsd_lvl(lvl) = std_lvl(lvl) / d_min; + else + std_lvl(lvl) = NaN; + nsd_lvl(lvl) = NaN; + end end - catch - std_lvl = NaN; - warning('No EVM per level calculated') end -end - function [data_,reference_] = trimseq(data,reference,skipstart,skip_end) data_ = data(skipstart+1:end-skip_end,:); diff --git a/Functions/Theory/propagation_time.m b/Functions/Theory/propagation_time.m new file mode 100644 index 0000000..9c64897 --- /dev/null +++ b/Functions/Theory/propagation_time.m @@ -0,0 +1,11 @@ +function prop_time = propagation_time(fiber_len_m) + + c = physconst('lightspeed'); % Speed of light in vacuum (m/s) + n = 1.46; % Refractive index of the fiber + + % Calculate the propagation speed in the fiber + propagation_speed = c / n; + + % Calculate the propagation time + prop_time = fiber_len_m ./ propagation_speed; +end \ No newline at end of file diff --git a/Functions/beautifyBERplot.m b/Functions/beautifyBERplot.m index e7a02ef..0d0a665 100644 --- a/Functions/beautifyBERplot.m +++ b/Functions/beautifyBERplot.m @@ -9,9 +9,9 @@ function beautifyBERplot() for i = 1:length(lines) lines(i).LineWidth = 1.3; % Thicker line width %lines(i).LineStyle = '-'; % Solid lines for simplicity - if string(lines(i).Marker) == "none" - lines(i).Marker = markers{mod(i-1, num_markers) + 1}; % Assign markers cyclically - end + % if string(lines(i).Marker) == "none" + % lines(i).Marker = markers{mod(i-1, num_markers) + 1}; % Assign markers cyclically + % end lines(i).MarkerSize = 4; % Marker size lines(i).MarkerFaceColor = 'auto'; % Use line color for marker face end diff --git a/projects/ECOC_2025/dsp_run_id.m b/projects/ECOC_2025/dsp_run_id.m index 5ef32d3..d1941a9 100644 --- a/projects/ECOC_2025/dsp_run_id.m +++ b/projects/ECOC_2025/dsp_run_id.m @@ -12,13 +12,13 @@ arguments options.storage_path end -% database = DBHandler("pathToDB",[options.database_path,options.database_name],"type","sqlite"); -database = DBHandler("type","mysql"); +database = DBHandler("pathToDB",[options.database_path,options.database_name],"type","sqlite"); +% database = DBHandler("type","mysql"); filterParams = database.tables; filterParams.Configurations = struct('run_id', run_id); -selectedFields = {'Runs.run_id','Runs.tx_bits_path','Runs.tx_signal_path','Runs.tx_symbols_path','Runs.rx_sync_path','Runs.rx_raw_path',... +selectedFields = {'Runs.run_id','Runs.tx_bits_path','Runs.tx_symbols_path','Runs.rx_sync_path','Runs.rx_raw_path',... 'Configurations.db_mode','Configurations.pam_level','Configurations.bitrate','Configurations.symbolrate','Configurations.fiber_length','Configurations.wavelength','Configurations.precomp_amp','Measurements.power_rop','Configurations.v_bias',... 'Configurations.interference_attenuation'}; @@ -28,8 +28,11 @@ dataTable = dataTable(uniqueIdx,:); % Extract unique configurations for each ru fsym = dataTable.symbolrate; M = double(dataTable.pam_level); +try duob_mode = db_mode.(strrep(char(dataTable.db_mode),'"','')); - +catch +duob_mode = db_mode(dataTable.db_mode); +end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% len_tr = 4096*2; @@ -82,8 +85,6 @@ dbtgt_package = {}; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -Tx_signal = load([options.storage_path, char(dataTable.tx_signal_path)]); -Tx_signal = Tx_signal.Digi_sig; Tx_bits = load([options.storage_path, char(dataTable.tx_bits_path)]); Tx_bits = Tx_bits.Bits; @@ -98,21 +99,21 @@ found_sync = 0; try Scpe_load = load([options.storage_path, char(dataTable.rx_sync_path)]); Scpe_cell = Scpe_load.S; - [~,~,found_sync] = Scpe_cell{2}.tsynch("reference",Symbols,"fs_ref",dataTable.symbolrate,"debug_plots",0); + [~,~,~,found_sync] = Scpe_cell{2}.tsynch("reference",Symbols,"fs_ref",dataTable.symbolrate,"debug_plots",1); end if ~found_sync Scpe_sig_raw = load([options.storage_path, char(dataTable.rx_raw_path(1))]); Scpe_sig_raw = Scpe_sig_raw.Scpe_sig_raw; Scpe_sig_resampled = Scpe_sig_raw.resample("fs_in",Scpe_sig_raw.fs,"fs_out",2*fsym); - [~,Scpe_cell,found_sync] =Scpe_sig_resampled.tsynch("reference",Symbols,"fs_ref",dataTable.symbolrate,"debug_plots",1); + [~,Scpe_cell,~,found_sync] =Scpe_sig_resampled.tsynch("reference",Symbols,"fs_ref",dataTable.symbolrate,"debug_plots",1); end if ~found_sync if length(Symbols_mapped.signal) == sum(Symbols_mapped.signal == Symbols.signal) warning('Could not synchronize the received signal with the stored symbols!') else - [~,Scpe_cell,found_sync] =Scpe_sig_raw.tsynch("reference",Symbols_mapped,"fs_ref",dataTable.symbolrate,"debug_plots",0); + [~,Scpe_cell,~,found_sync] =Scpe_sig_raw.tsynch("reference",Symbols_mapped,"fs_ref",dataTable.symbolrate,"debug_plots",0); end if ~found_sync warning('Could not synchronize the received signal with the stored symbols!') @@ -120,7 +121,7 @@ if ~found_sync end record_realizations = min(options.max_occurences,length(Scpe_cell)); -for occ = 4:record_realizations +for occ = 1:record_realizations Scpe_sig = Scpe_cell{occ}; @@ -136,10 +137,10 @@ for occ = 4:record_realizations if duob_mode ~= db_mode.db_encoded - vnle_pf = 0; + vnle_pf = 1; dbtgt = 0; % %%%%% VNLE + DFE %%%% - if 1 + if 0 eq_ = FFE_DCremoval_adaptive_mu("epochs_tr",5,"epochs_dd",3,"len_tr",4096*2,"mu_dd",... 0.0002,"mu_tr",0,"order",25,"sps",2,"decide",0,... @@ -149,6 +150,8 @@ for occ = 4:record_realizations "smoothing_buffer_length",smoothing_buffer_length,... "smoothing_buffer_update",smoothing_buffer_update); + eq_ = EQ("Ne",vnle_order,"Nb",dfe_order,"training_length",len_tr,"training_loops",5,"dd_loops",5,"K",2,"DCmu",mu_dc,"DDmu",[mu_ffe mu_dfe],"DFEmu",0.005,"FFEmu",0,"plotfinal",0,"ideal_dfe",1); + result = vnle(eq_,M,Scpe_sig,Symbols,Tx_bits,"precode_mode",duob_mode,'showAnalysis',1,"postFFE",[],"eth_style_symbol_mapping",0); vnle_package{occ} = result; fprintf("FFE Results: %.2e\n", result.ber_vnle); @@ -196,7 +199,7 @@ for occ = 4:record_realizations % Print header fprintf("==== EQUALIZATION RUN-ID %d | PAM-%d | %.2f GBd ====\n\n", run_id, M, Symbols.fs.*1e-9); -s + % VNLE Results fprintf(">> VNLE Results:\n"); fprintf(" BER %.2e\n", vnle_result.BER); diff --git a/projects/ECOC_2025/dsp_loop_id.m b/projects/ECOC_2025/dsp_standalone_2.m similarity index 98% rename from projects/ECOC_2025/dsp_loop_id.m rename to projects/ECOC_2025/dsp_standalone_2.m index 2d9ab8f..97007ae 100644 --- a/projects/ECOC_2025/dsp_loop_id.m +++ b/projects/ECOC_2025/dsp_standalone_2.m @@ -14,7 +14,7 @@ filterParams.Configurations = struct( ... 'fiber_length', 0, ... 'db_mode', '"no_db"', ... 'interference_attenuation', [], ... - 'interference_path_length', 10, ... + 'interference_path_length', 0, ... 'is_mpi', 1, ... 'pam_level', 4, ... 'wavelength', 1310, ... diff --git a/projects/ECOC_2025/load_signal_standalone.m b/projects/ECOC_2025/load_signal_standalone.m index cb27584..f58c77f 100644 --- a/projects/ECOC_2025/load_signal_standalone.m +++ b/projects/ECOC_2025/load_signal_standalone.m @@ -1,6 +1,5 @@ -run_id = 3966; savePath = 'Z:\2025\ECOC Silas\ecoc_2025\'; databasePath = 'C:\Users\Silas\Documents\MATLAB\imdd_simulation\projects\ECOC_2025\'; @@ -11,37 +10,38 @@ db = DBHandler("type","mysql"); filterParams = db.tables; % filterParams.Configurations = struct('run_id', run_id); filterParams.Configurations = struct( ... - 'symbolrate', [], ... %[224,336,360,390,420,448] + 'symbolrate', 112e9, ... %[224,336,360,390,420,448] 'fiber_length', 0, ... 'db_mode', '"no_db"', ... 'interference_attenuation', 4, ... - 'interference_path_length', 10, ... + 'interference_path_length', 300, ... 'is_mpi', 1, ... - 'pam_level', [], ... + 'pam_level', 4, ... 'wavelength', 1310, ... 'precomp_amp', [], ... 'signal_attenuation', [], ... - 'v_awg', [], ... - 'v_bias', 2.65 ... + 'v_awg', [], ... + 'v_bias', [] ... ); - selectedFields = {'Runs.run_id','Runs.tx_bits_path','Runs.tx_signal_path','Runs.tx_symbols_path','Runs.rx_sync_path','Runs.rx_raw_path',... 'Configurations.db_mode','Configurations.pam_level','Configurations.bitrate','Configurations.symbolrate','Configurations.fiber_length','Configurations.wavelength','Configurations.precomp_amp','Measurements.power_rop','Configurations.v_bias',... - 'Configurations.interference_attenuation'}; + 'Configurations.interference_attenuation', 'Measurements.power_mpi_interference'}; [dataTable,sql_query] = db.queryDB(filterParams, selectedFields); [~, uniqueIdx] = unique(dataTable.run_id); % Get unique run_id indices +dataTable.SIR = -7 - round(dataTable.power_mpi_interference); -for i = 1:4%size(dataTable,1) +%% +for i = 1:size(dataTable,1) dataTable_ = dataTable(i,:); % Extract unique configurations for each run_id - + fsym = dataTable_.symbolrate; M = double(dataTable_.pam_level); duob_mode = db_mode.(strrep(char(dataTable_.db_mode),'"','')); - + Tx_signal = load([savePath, char(dataTable_.tx_signal_path)]); Tx_signal = Tx_signal.Digi_sig; @@ -52,17 +52,52 @@ for i = 1:4%size(dataTable,1) Symbols = load([savePath, char(dataTable_.tx_symbols_path)]); Symbols = Symbols.Symbols; - + Scpe_sig_raw = load([savePath, char(dataTable_.rx_raw_path(1))]); Scpe_sig_raw = Scpe_sig_raw.Scpe_sig_raw; - Scpe_sig_raw.plot("displayname",['Scope Signal (Run ID: ',num2str(dataTable_.run_id)],"fignum",dataTable_.run_id,"clear",0); + % Scpe_sig_raw.plot("displayname",['Scope Signal (Run ID: ',num2str(dataTable_.run_id)],"fignum",dataTable_.run_id,"clear",0); - Scpe_sig_resampled = Scpe_sig_raw.resample("fs_in",Scpe_sig_raw.fs,"fs_out",2*fsym); - [~,Scpe_cell,found_sync,~,shifts] = Scpe_sig_resampled.tsynch("reference",Symbols,"fs_ref",dataTable_.symbolrate,"debug_plots",1); + Scpe_sig_resampled = Scpe_sig_raw.resample("fs_in",Scpe_sig_raw.fs,"fs_out",fsym); + [~,Scpe_cell,found_sync,test,shifts] = Scpe_sig_resampled.tsynch("reference",Symbols,"fs_ref",dataTable_.symbolrate,"debug_plots",0); shifts_mus = shifts./Scpe_sig_resampled.fs .*1e6; - Scpe_sig_raw.plot("displayname",['Scope Signal (Run ID: ',num2str(dataTable_.run_id)],"fignum",dataTable_.run_id,"clear",0); + Scpe_sig_resampled = Scpe_sig_resampled.normalize("mode","rms"); + % Scpe_sig_resampled.plot("displayname",['Scope Signal (Run ID: ',num2str(dataTable_.run_id)],"fignum",dataTable_.run_id,"clear",0); hold on; - xline(shifts_mus,'HandleVisibility','off'); + % xline(shifts_mus,'HandleVisibility','off'); + + shifts = shifts-shifts(1)+1; + sep_sig = NaN(M,length(Scpe_sig_resampled)); + avg_sig = NaN(M,length(Scpe_sig_resampled)); + for j = 1:size(Scpe_cell,1) + + [sep_sig_,avg_sig_]=showLevelScatter(Scpe_cell{j}.resample("fs_out",Symbols.fs),Symbols,"fignum",400); + s = shifts(j); + sep_sig(:,s+1:s+length(sep_sig_)) = sep_sig_; + avg_sig(:,s+1:s+length(avg_sig_)) = avg_sig_; + + end + + disp(num2str(filterParams.Configurations.interference_path_length)); + var(sep_sig,0,2,'omitnan') +%% + xax_in_sec = ((1:length(avg_sig)) / fsym) * 1e6; + figure();hold on; + cols = cbrewer2('Paired',8); + % for p = 1:size(avg_sig,1) + % sc=scatter(xax_in_sec,sep_sig(p,:),1,'.','MarkerEdgeColor',cols((2*p)-1,:),'MarkerEdgeAlpha',0.1); + % end + for p = 1:size(avg_sig,1) + sc=plot(xax_in_sec,avg_sig(p,:),'LineWidth',1,'Color',cols((2*p),:)); + end + + yline(unique(Symbols.signal),'HandleVisibility','off'); + % xline(shifts./ fsym .*1e6,'HandleVisibility','off'); + xlabel('time in $\mu$s'); + ylabel('Normalized Amplitude'); + xlim([0 25]); + ylim([-2 2]); + + drawnow; end \ No newline at end of file diff --git a/projects/ECOC_2025/plots_from_database/plot_mpi_trial.m b/projects/ECOC_2025/plots_from_database/plot_mpi_trial.m index eadb53c..3b2696c 100644 --- a/projects/ECOC_2025/plots_from_database/plot_mpi_trial.m +++ b/projects/ECOC_2025/plots_from_database/plot_mpi_trial.m @@ -16,7 +16,7 @@ filterParams.Configurations = struct( ... 'wavelength', 1310, ... 'precomp_amp', [], ... 'signal_attenuation', [], ... - 'v_awg', [], ... + 'v_awg', [], ... 'v_bias', [] ... ); @@ -25,7 +25,7 @@ filterParams.Configurations = struct( ... % % filterParams.EqualizerParameters.dc_buffer_len = 1; % filterParams.EqualizerParameters.ffe_buffer_len = 1; % filterParams.EqualizerParameters.smoothing_buffer_len = 4096; -% filterParams.EqualizerParameters.smoothing_buffer_update = 224; +% filterParams.EqualizerParameters.smoothing_buffer_update = 224; % filterParams.EqualizerParameters.DCmu = 0; % end @@ -34,220 +34,291 @@ selectedFields = {'Configurations.run_id' 'Runs.loop_id' 'Runs.date_of_run' 'Run 'EqualizerParameters.equalizer_structure' 'EqualizerParameters.diff_precode' 'EqualizerParameters.eq_id' 'EqualizerParameters.dc_buffer_len' 'EqualizerParameters.ffe_buffer_len' 'EqualizerParameters.smoothing_buffer_len' 'EqualizerParameters.smoothing_buffer_update' 'EqualizerParameters.DCmu' 'Measurements.power_pd_in' ... 'Measurements.power_mpi_interference' 'Measurements.power_mpi_signal' 'Results.BER' 'Results.BER_precoded' 'Results.EVM' 'Results.SNR' 'Results.GMI' 'Results.Alpha' 'Results.date_of_processing'}; -% [dataTable_raw,sql_query] = database.queryDB(filterParams, selectedFields); +[dataTable_raw,sql_query] = database.queryDB(filterParams, selectedFields); + +%% + +dataTable_clean = dataTable_raw; +dataTable_clean.SIR = -7 - round(dataTable_clean.power_mpi_interference); +dataTable_clean.NGMI = dataTable_clean.GMI ./ log2(dataTable_clean.pam_level); +dataTable_clean = cleanUpTable(dataTable_clean); + +dataTable_clean(dataTable_clean.BER>0.2,:) = []; %% - - -dataTable_clean = dataTable_raw; -dataTable_clean.SIR = round(-7.5 - dataTable_clean.power_mpi_interference); -dataTable_clean.NGMI = dataTable_clean.GMI ./ log2(dataTable_clean.pam_level); -dataTable_clean = cleanUpTable(dataTable_clean); - -%% - -dataTable = dataTable_clean; - -figure(26); -plotBoundaries = 1; -plotRealizations = 0; cols = linspecer(8); % Ensure color count matches +figure() +tiledlayout(1, 4, 'TileSpacing', 'compact', 'Padding', 'compact'); +y_here = 0; +figcnt = 0; +for int_len = [0,50,300,1000] + figcnt = figcnt+1; + % figure(int_len+1); + nexttile; + hold on + mode = 4; -% ideal DC tracking -if 0 - dataTable = dataTable(dataTable.dc_buffer_len == 1, :); - dataTable = dataTable(dataTable.ffe_buffer_len == 1, :); - dataTable = dataTable(dataTable.smoothing_buffer_len == 0, :); - dataTable = dataTable(dataTable.smoothing_buffer_update == 0, :); - dataTable = dataTable(dataTable.DCmu == 0.005, :); - cols = cols(1:1+1,:); - method = 'ideal dc tracking'; - % slow DC tracking -elseif 1 - dataTable = dataTable(dataTable.dc_buffer_len == 224, :); - dataTable = dataTable(dataTable.ffe_buffer_len == 1, :); - dataTable = dataTable(dataTable.smoothing_buffer_len == 0, :); - dataTable = dataTable(dataTable.smoothing_buffer_update == 0, :); - dataTable = dataTable(dataTable.DCmu == 0.005, :); - cols = cols(2:2+1,:); - method = 'parallelized dc tracking'; - % slow DC smoothing -elseif 0 + for mode = [1,4] + + hold on; + dataTable = dataTable_clean; - dataTable = dataTable(dataTable.dc_buffer_len == 1, :); - dataTable = dataTable(dataTable.ffe_buffer_len == 1, :); - dataTable = dataTable(dataTable.smoothing_buffer_len == 4096, :); - dataTable = dataTable(dataTable.smoothing_buffer_update == 224, :); - dataTable = dataTable(dataTable.DCmu == 0, :); - cols = cols(3:3+1,:); - method = 'dc smoothing'; -elseif 0 - % No compensation method - dataTable = dataTable(dataTable.dc_buffer_len == 1, :); - dataTable = dataTable(dataTable.ffe_buffer_len == 1, :); - dataTable = dataTable(dataTable.smoothing_buffer_len == 0, :); - dataTable = dataTable(dataTable.smoothing_buffer_update == 0, :); - dataTable = dataTable(dataTable.DCmu == 0, :); - cols = cols(4:4+1,:); - method = 'ffe only'; -end -dataTable(dataTable.eq_id==0,:) = []; -dataTable(dataTable.equalizer_structure~=1,:) = []; + plotBoundaries = 1; + plotRealizations = 1; + cols = linspecer(8); % Ensure color count matches -% dataTable.interference_path_length(dataTable.interference_path_length < 20, :) = 20; + if mode == 1 + % No compensation method + dataTable = dataTable(dataTable.dc_buffer_len == 1, :); + dataTable = dataTable(dataTable.ffe_buffer_len == 1, :); + dataTable = dataTable(dataTable.smoothing_buffer_len == 0, :); + dataTable = dataTable(dataTable.smoothing_buffer_update == 0, :); + dataTable = dataTable(dataTable.DCmu == 0, :); + cols = cols(1:1+1,:); + method = 'ffe only'; -dataTable = dataTable(dataTable.interference_path_length == 1000, :); + % slow DC smoothing + elseif mode == 2 -% dataTable(dataTable.interference_path_length ~= 50, :) = []; -% dataTable = dataTable(dataTable.interference_path_length < 51, :); + dataTable = dataTable(dataTable.dc_buffer_len == 1, :); + dataTable = dataTable(dataTable.ffe_buffer_len == 1, :); + dataTable = dataTable(dataTable.smoothing_buffer_len == 4096, :); + dataTable = dataTable(dataTable.smoothing_buffer_update == 224, :); + dataTable = dataTable(dataTable.DCmu == 0, :); + + cols = cols(4:4+1,:); + method = 'dc smoothing'; -% dataTable(dataTable.loop_id<200,:) = []; -% Filter by time -filter_by_time = 0; -if filter_by_time - startTime = datetime('2025-04-20 18:00:00', 'InputFormat', 'yyyy-MM-dd HH:mm:ss'); - stopTime = datetime('2025-04-30 19:30:00', 'InputFormat', 'yyyy-MM-dd HH:mm:ss'); - dataTable.date_of_run = datetime(dataTable.date_of_run, 'InputFormat', 'yyyy-MM-dd HH:mm:ss.SSSSSS'); - dataTable.date_of_processing = datetime(dataTable.date_of_processing, 'InputFormat', 'yyyy-MM-dd HH:mm:ss.SSSSSS'); - dataTable = dataTable(dataTable.date_of_processing > startTime, :); - dataTable = dataTable(dataTable.date_of_processing < stopTime, :); -end + % slow DC tracking + elseif mode == 3 -% Group by smth -y_var = 'BER'; -x_var = 'SIR'; + dataTable = dataTable(dataTable.dc_buffer_len == 224, :); + dataTable = dataTable(dataTable.ffe_buffer_len == 1, :); + dataTable = dataTable(dataTable.smoothing_buffer_len == 0, :); + dataTable = dataTable(dataTable.smoothing_buffer_update == 0, :); + dataTable = dataTable(dataTable.DCmu == 0.005, :); + + cols = cols(3:3+1,:); + method = 'parallelized dc tracking'; -loop_var = 'interference_path_length'; -fixedVars = {'equalizer_structure','interference_path_length',x_var}; + elseif mode == 4 -[dataTable, outliersTable] = removeGroupOutliers(dataTable, fixedVars, y_var); + % ideal DC tracking + dataTable = dataTable(dataTable.dc_buffer_len == 1, :); + dataTable = dataTable(dataTable.ffe_buffer_len == 1, :); + dataTable = dataTable(dataTable.smoothing_buffer_len == 0, :); + dataTable = dataTable(dataTable.smoothing_buffer_update == 0, :); + dataTable = dataTable(dataTable.DCmu == 0.005, :); + + cols = cols(2:2+1,:); + method = 'ideal dc tracking'; + end -dataTableGrpd_mean = groupIt(fixedVars, dataTable, @mean); -dataTableGrpd_min = groupIt(fixedVars, dataTable, @min); -dataTableGrpd_max = groupIt(fixedVars, dataTable, @max); + dataTable(dataTable.eq_id==0,:) = []; + dataTable(dataTable.equalizer_structure~=1,:) = []; -% Create a new figure -mkr = '.'; -hold on + % Modify values in 'interference_path_length' where the condition is met + dataTable.interference_path_length(dataTable.interference_path_length < 101 & dataTable.interference_path_length > 1) = 50; -unique_loop_var = unique(dataTable.(loop_var)); + dataTable.interference_path_length(dataTable.interference_path_length == 1) = 0; + dataTable = dataTable(dataTable.interference_path_length == int_len, :); -for i = 1:numel(unique_loop_var) - % Prepare filtered data for this loop variable - loopValue = unique_loop_var(i); + dataTable(dataTable.run_id == 3866, :) = []; + dataTable(dataTable.run_id == 3865, :) = []; + dataTable(dataTable.run_id == 3796, :) = []; + dataTable(dataTable.run_id == 3797, :) = []; + dataTable(dataTable.run_id == 3798, :) = []; + dataTable(dataTable.run_id == 4002, :) = []; + dataTable(dataTable.run_id == 4200, :) = []; + dataTable(dataTable.run_id == 4199, :) = []; - loopFiltGrpd = dataTableGrpd_mean.(loop_var) == loopValue; - if ~any(loopFiltGrpd) - continue; % Skip if no data for this loop var - end + % 0 + % 1 + % 10 + % 15 + % 20 + % 50 + % 100 + % 300 + % 1000 - % Extract values - x_values = dataTableGrpd_mean.(x_var)(loopFiltGrpd, :); - y_mean = dataTableGrpd_mean.(y_var)(loopFiltGrpd, :); - y_min = dataTableGrpd_min.(y_var)(loopFiltGrpd, :); - y_max = dataTableGrpd_max.(y_var)(loopFiltGrpd, :); + % dataTable(dataTable.interference_path_length ~= 50, :) = []; + % dataTable = dataTable(dataTable.interference_path_length < 51, :); - % Compute bounds: distance from mean - y_lower = y_mean - y_min; - y_upper = y_max - y_mean; - y_bounds = [y_lower, y_upper]; + % dataTable(dataTable.loop_id<200,:) = []; - % Display name (optional) - try - idx = find(dataTable.(loop_var) == loopValue, 1, 'first'); - % dispname = char(equalizer_structure(dataTable.equalizer_structure(idx))); - dispname = [method]; - dispname = [dispname, '/ ',num2str(unique_loop_var(i)) ,' m']; - % dispname = [dispname,'; ',num2str(unique(dataTable.interference_path_length)),' m']; - dispname = [dispname, '/ PAM ', num2str(filterParams.Configurations.pam_level)]; - dispname = [dispname, '/ ', num2str(filterParams.Configurations.symbolrate.*1e-9),' GBd']; - end + % Filter by time + filter_by_time = 0; + if filter_by_time + startTime = datetime('2025-04-20 18:00:00', 'InputFormat', 'yyyy-MM-dd HH:mm:ss'); + stopTime = datetime('2025-04-30 19:30:00', 'InputFormat', 'yyyy-MM-dd HH:mm:ss'); + dataTable.date_of_run = datetime(dataTable.date_of_run, 'InputFormat', 'yyyy-MM-dd HH:mm:ss.SSSSSS'); + dataTable.date_of_processing = datetime(dataTable.date_of_processing, 'InputFormat', 'yyyy-MM-dd HH:mm:ss.SSSSSS'); + dataTable = dataTable(dataTable.date_of_processing > startTime, :); + dataTable = dataTable(dataTable.date_of_processing < stopTime, :); + end - if plotBoundaries - % Plot bounded line - [hl, hp] = boundedline(x_values, y_mean, y_bounds, ... - 'alpha', 'transparency', 0.1, ... - 'cmap', cols(i,:), ... - 'nan', 'fill', ... - 'orientation', 'vert'); - - % % Style the main line: thinnest, dotted, no marker - set(hl, 'LineWidth', 1, 'LineStyle', '-', 'Marker', 'none', ... - 'Color', cols(i,:), 'DisplayName', string(dispname)); + % Group by smth + y_var = 'BER'; + x_var = 'SIR'; - plt = errorbar(x_values,y_mean,y_lower,y_upper,'LineWidth', 0.1, 'LineStyle', 'none', 'Marker', 'none','Color', cols(i,:), 'DisplayName', string(dispname),'HandleVisibility','off'); - - % Hide patch (shaded area) from legend - set(hp, 'HandleVisibility', 'off','LineStyle',':','LineWidth',0.5,'Marker','none'); - - % % Add invisible scatter for DataTips - % plt = scatter(x_values, y_mean, ... - % 'Marker', 'o', 'MarkerEdgeColor', 'none', 'MarkerFaceColor', 'none', ... - % 'HandleVisibility', 'off', 'PickableParts', 'all'); - else + loop_var = 'interference_path_length'; + fixedVars = {'equalizer_structure','interference_path_length',x_var}; - plt= plot(x_values,y_mean,'LineWidth', 1, 'LineStyle', '-', 'Marker', 'none', ... - 'Color', cols(i,:), 'DisplayName', string(dispname)); - % plt= errorbar(x_values,y_mean,y_lower,y_upper,'LineWidth', 1, 'LineStyle', '-', 'Marker', 'none','Color', cols(i,:), 'DisplayName', string(dispname)); + [dataTable, outliersTable] = removeGroupOutliers(dataTable, fixedVars, y_var); + + dataTableGrpd_mean = groupIt(fixedVars, dataTable, @mean); + dataTableGrpd_min = groupIt(fixedVars, dataTable, @min); + dataTableGrpd_max = groupIt(fixedVars, dataTable, @max); + + % Create a new figure + + hold on + + unique_loop_var = unique(dataTable.(loop_var)); + + for i = 1:numel(unique_loop_var) + + % Prepare filtered data for this loop variable + loopValue = unique_loop_var(i); + + loopFiltGrpd = dataTableGrpd_mean.(loop_var) == loopValue; + if ~any(loopFiltGrpd) + continue; % Skip if no data for this loop var + end + + % Extract values + x_values = dataTableGrpd_mean.(x_var)(loopFiltGrpd, :); + y_mean = dataTableGrpd_mean.(y_var)(loopFiltGrpd, :); + y_min = dataTableGrpd_min.(y_var)(loopFiltGrpd, :); + y_max = dataTableGrpd_max.(y_var)(loopFiltGrpd, :); + + % Compute bounds: distance from mean + y_lower = y_mean - y_min; + y_upper = y_max - y_mean; + y_bounds = [y_lower, y_upper]; + + % Display name (optional) + try + idx = find(dataTable.(loop_var) == loopValue, 1, 'first'); + % dispname = char(equalizer_structure(dataTable.equalizer_structure(idx))); + dispname = [method]; + dispname = [dispname, '/ ',num2str(unique_loop_var(i)) ,' m']; + % dispname = [dispname,'; ',num2str(unique(dataTable.interference_path_length)),' m']; + % dispname = [dispname, '/ PAM ', num2str(filterParams.Configurations.pam_level)]; + % dispname = [dispname, '/ ', num2str(filterParams.Configurations.symbolrate.*1e-9),' GBd']; + end + + if plotBoundaries + % Plot bounded line + [hl, hp] = boundedline(x_values, y_mean, y_bounds, ... + 'alpha', 'transparency', 0.1, ... + 'cmap', cols(i,:), ... + 'nan', 'fill', ... + 'orientation', 'vert'); + + % % Style the main line: thinnest, dotted, no marker + set(hl, 'LineWidth', 0.5, 'LineStyle', ':', 'Marker', 'none', ... + 'Color', cols(i,:), 'DisplayName', string(dispname)); + + plt = errorbar(x_values,y_mean,y_lower,y_upper,'LineWidth', 0.9, 'LineStyle', 'none', 'Marker', 'none','Color', cols(i,:), 'DisplayName', string(dispname),'HandleVisibility','off'); + + % Hide patch (shaded area) from legend + set(hp, 'HandleVisibility', 'off','LineStyle',':','LineWidth',0.5,'Marker','none'); + + % Fit a 4th-order polynomial to log10(BER) + p = polyfit(x_values, log10(y_mean), 3); % 4 is fitting order, adjust as needed + + % Evaluate the fitted polynomial + x_fit = linspace(min(x_values), max(x_values), 300); % Fine points + y_fit_log = polyval(p, x_fit); % Still in log10 domain + y_fit = 10.^y_fit_log; % Back to BER domain + + + + plot(x_fit,y_fit,'LineWidth', 1, 'LineStyle', '-', 'Marker', 'none', ... + 'Color', cols(i,:), 'DisplayName', string(dispname),'HandleVisibility','off'); + + % % Add invisible scatter for DataTips + % plt = scatter(x_values, y_mean, ... + % 'Marker', 'o', 'MarkerEdgeColor', 'none', 'MarkerFaceColor', 'none', ... + % 'HandleVisibility', 'off', 'PickableParts', 'all'); + else + + plt= plot(x_values,y_mean,'LineWidth', 1, 'LineStyle', '-', 'Marker', 'none', ... + 'Color', cols(i,:), 'DisplayName', string(dispname)); + % plt= errorbar(x_values,y_mean,y_lower,y_upper,'LineWidth', 1, 'LineStyle', '-', 'Marker', 'none','Color', cols(i,:), 'DisplayName', string(dispname)); + + end + % Add data tips to the invisible scatter + pair_one = {'Run ID', dataTableGrpd_mean.run_id(loopFiltGrpd, :)}; + pair_two = {'Rate', dataTableGrpd_mean.bitrate(loopFiltGrpd, :) * 1e-9}; + pair_three = {'PD in', round(dataTableGrpd_mean.power_pd_in(loopFiltGrpd, :), 2)}; + addDatatips(plt, pair_one, pair_two, pair_three); + + + + % Optionally: outline bounds for better visibility (optional) + % hnew = outlinebounds(hl, hp); + + % Tick marks (x-axis) + + + % Optional: scatter realizations + if plotRealizations + loopFiltSingle = dataTable.(loop_var) == loopValue; + x_single = double(dataTable.(x_var)(loopFiltSingle, :)); + y_single = double(dataTable.(y_var)(loopFiltSingle, :)); + + mkr = '.'; + sc = scatter(x_single+(mode*0.1)-0.2, y_single,15, 'Marker', mkr, 'MarkerEdgeColor', cols(i, :), ... + 'LineWidth', 0.5, 'HandleVisibility', 'off', 'DisplayName', string(dispname)); + + pair_one = {'Run ID', dataTable.run_id(loopFiltSingle, :)}; + pair_two = {'Rate', dataTable.bitrate(loopFiltSingle, :) * 1e-9}; + pair_three = {'PD in', round(dataTable.power_pd_in(loopFiltSingle, :), 2)}; + addDatatips(sc, pair_one, pair_two, pair_three); + end + end + + % Label axes and title + legend('Interpreter', 'latex'); + xlabel(x_var); + if ~y_here + ylabel(y_var); + yticklabels = []; + + y_here = 1; + end + if int_len ~= 0 + set(gca, 'YTick', []); + end + % title([x_var, ' vs. ', y_var]); + + if string(y_var) == "BER" + yline(2.2e-4, 'LineWidth', 1, 'LineStyle', '--', 'HandleVisibility', 'off'); + yline(3.8e-3, 'LineWidth', 1, 'LineStyle', '--', 'HandleVisibility', 'off'); + yline(2e-2, 'LineWidth', 1, 'LineStyle', '--', 'HandleVisibility', 'off'); + ylim([1e-5, 0.1]); + end + + xlim([15,35]); + ylim([9e-5 0.1 ]); + xticks([13:2:35]); + + % Enable grid and beautify + grid on; + beautifyBERplot; - end - % Add data tips to the invisible scatter - pair_one = {'Run ID', dataTableGrpd_mean.run_id(loopFiltGrpd, :)}; - pair_two = {'Rate', dataTableGrpd_mean.bitrate(loopFiltGrpd, :) * 1e-9}; - pair_three = {'PD in', round(dataTableGrpd_mean.power_pd_in(loopFiltGrpd, :), 2)}; - addDatatips(plt, pair_one, pair_two, pair_three); - - - - % Optionally: outline bounds for better visibility (optional) - % hnew = outlinebounds(hl, hp); - - % Tick marks (x-axis) - xticks(round(unique(x_values),2)); - - % Optional: scatter realizations - if plotRealizations - loopFiltSingle = dataTable.(loop_var) == loopValue; - x_single = double(dataTable.(x_var)(loopFiltSingle, :)); - y_single = double(dataTable.(y_var)(loopFiltSingle, :)); - - sc = scatter(x_single, y_single, 'Marker', mkr, 'MarkerEdgeColor', cols(i, :), ... - 'LineWidth', 0.5, 'HandleVisibility', 'off', 'DisplayName', string(dispname)); - - pair_one = {'Run ID', dataTable.run_id(loopFiltSingle, :)}; - pair_two = {'Rate', dataTable.bitrate(loopFiltSingle, :) * 1e-9}; - pair_three = {'PD in', round(dataTable.power_pd_in(loopFiltSingle, :), 2)}; - addDatatips(sc, pair_one, pair_two, pair_three); end end -% Label axes and title -legend('Interpreter', 'latex'); -xlabel(x_var); -ylabel(y_var); -title([x_var, ' vs. ', y_var]); - -if string(y_var) == "BER" - yline(4e-4, 'LineWidth', 1, 'LineStyle', '--', 'HandleVisibility', 'off'); - yline(3.8e-3, 'LineWidth', 1, 'LineStyle', '--', 'HandleVisibility', 'off'); - yline(2e-2, 'LineWidth', 1, 'LineStyle', '--', 'HandleVisibility', 'off'); - ylim([1e-5, 0.1]); -end - -xlim([13,35]); - -% Enable grid and beautify -grid on; -beautifyBERplot; - - - function resultTable = groupIt(fixedVars, dataTable, aggregationFunction) % groupIt Groups data in a table based on fixedVars and applies aggregationFunction to numeric data. @@ -463,9 +534,9 @@ for groupIdx = 1:height(groupKeys) continue; end - % Detect outliers in log space + % Detect outliers in log space y_log = log10(y_values); - outlierMask = isoutlier(y_log, 'quartiles'); % or 'median', 'grubbs', etc. + outlierMask = isoutlier(y_log, 'quartiles',1); % or 'median', 'grubbs', etc. % If any outliers found, collect their data if any(outlierMask) diff --git a/projects/HighSpeedExperiment_2024/auswertung MPI/plot_from_database.m b/projects/HighSpeedExperiment_2024/auswertung MPI/plot_from_database.m index 3093427..0aa8010 100644 --- a/projects/HighSpeedExperiment_2024/auswertung MPI/plot_from_database.m +++ b/projects/HighSpeedExperiment_2024/auswertung MPI/plot_from_database.m @@ -1,26 +1,27 @@ basePath = 'C:\Users\Silas\Documents\MATLAB\Datensätze\sioe_labor\'; -database = DBHandler("pathToDB",[basePath,'silas_labor.db']); +database = DBHandler("pathToDB",[basePath,'silas_labor_plain.db'],"type",'sqlite'); filterParams = database.tables; filterParams.Configurations = struct( ... 'bitrate', [], ... %[224,336,360,390,420,448] - 'db_mode', int32(db_mode.no_db), ... + 'db_mode', [], ... 'fiber_length', 1, ... 'interference_attenuation', [], ... 'interference_path_length', [], ... - 'is_mpi', 1, ... + 'is_mpi', 0, ... 'pam_level', 4, ... 'rop_attenuation', 0, ... 'wavelength', 1310 ... ); % filterParams.EqualizerParameters.diff_precode = int32(db_mode.no_db); - filterParams.EqualizerParameters.equalizer_structure = int32(equalizer_structure.vnle); + % filterParams.EqualizerParameters.equalizer_structure = int32(equalizer_structure.vnle); % filterParams.EqualizerParameters.DCmu = 0.005; selectedFields = {'Configurations.run_id' 'Runs.rx_raw_path' 'Configurations.bitrate' 'Configurations.symbolrate' 'Configurations.pam_level' 'Configurations.db_mode' 'Configurations.rop_attenuation' 'Configurations.is_mpi' 'Configurations.interference_attenuation' 'EqualizerParameters.equalizer_structure' 'EqualizerParameters.diff_precode' 'EqualizerParameters.eq_id' 'Measurements.power_pd_in' 'Measurements.power_mpi_interference' 'Measurements.power_mpi_signal' 'Results.BER' 'Results.SNR' 'Results.GMI' 'Results.Alpha'}; +% selectedFields = {'Configurations.run_id'}; [dataTable,sql_query] = database.queryDB(filterParams, selectedFields); diff --git a/projects/HighSpeedExperiment_2024/db_auswertung/rate_vs_ber.m b/projects/HighSpeedExperiment_2024/db_auswertung/rate_vs_ber.m index e53d7c6..a42bad8 100644 --- a/projects/HighSpeedExperiment_2024/db_auswertung/rate_vs_ber.m +++ b/projects/HighSpeedExperiment_2024/db_auswertung/rate_vs_ber.m @@ -1,13 +1,13 @@ basePath = 'C:\Users\Silas\Documents\MATLAB\Datensätze\sioe_labor\'; -database = DBHandler("pathToDB",[basePath,'silas_labor.db']); +database = DBHandler("pathToDB",[basePath,'silas_labor.db'],"type",'sqlite'); filterParams = database.tables; filterParams.Configurations = struct( ... 'bitrate', [], ... %[224,336,360,390,420,448] - 'db_mode', int32(db_mode.db_encoded), ... - 'fiber_length', 10, ... + 'db_mode', [], ... + 'fiber_length', 1, ... 'interference_attenuation', [], ... 'interference_path_length', [], ... 'is_mpi', 0, ... @@ -18,11 +18,11 @@ filterParams.Configurations = struct( ... % filterParams.EqualizerParameters.diff_precode = int32(db_mode.db_encoded); % filterParams.EqualizerParameters.equalizer_structure = int32(equalizer_structure.vnle); -filterParams.EqualizerParameters.DCmu = 0.00; +% filterParams.EqualizerParameters.DCmu = 0.00; selectedFields = {'Configurations.run_id' 'Runs.rx_raw_path' 'Configurations.bitrate' 'Configurations.symbolrate' 'Configurations.pam_level'... 'Configurations.db_mode' 'Configurations.rop_attenuation' 'Configurations.is_mpi' 'Configurations.interference_attenuation' ... - 'EqualizerParameters.equalizer_structure' 'EqualizerParameters.diff_precode' 'EqualizerParameters.eq_id' 'Measurements.power_pd_in' ... + 'EqualizerParameters.equalizer_structure' 'EqualizerParameters.diff_precode' 'EqualizerParameters.eq_id' 'EqualizerParameters.DCmu' 'Measurements.power_pd_in' ... 'Measurements.power_mpi_interference' 'Measurements.power_mpi_signal' 'Results.BER' 'Results.BER_precoded' 'Results.SNR' 'Results.GMI' 'Results.Alpha' 'Results.date_of_processing'}; [dataTable,sql_query] = database.queryDB(filterParams, selectedFields); diff --git a/projects/HighSpeedExperiment_2024/db_auswertung/realization_vs_ber.m b/projects/HighSpeedExperiment_2024/db_auswertung/realization_vs_ber.m index b941680..786605e 100644 --- a/projects/HighSpeedExperiment_2024/db_auswertung/realization_vs_ber.m +++ b/projects/HighSpeedExperiment_2024/db_auswertung/realization_vs_ber.m @@ -1,11 +1,11 @@ basePath = 'C:\Users\Silas\Documents\MATLAB\Datensätze\sioe_labor\'; -database = DBHandler("pathToDB",[basePath,'silas_labor.db']); +database = DBHandler("pathToDB",[basePath,'silas_labor.db'],"type",'sqlite'); filterParams = database.tables; filterParams.Configurations = struct( ... - 'bitrate', 450e9, ... %[224,336,360,390,420,448] + 'bitrate', 420e9, ... %[224,336,360,390,420,448] 'db_mode', int32(db_mode.no_db), ... 'fiber_length', 10, ... 'interference_attenuation', [], ... diff --git a/projects/Job_Processing/run_offline_dsp.m b/projects/Job_Processing/run_offline_dsp.m new file mode 100644 index 0000000..2544680 --- /dev/null +++ b/projects/Job_Processing/run_offline_dsp.m @@ -0,0 +1,117 @@ +% === SETTINGS === + +dsp_options.append_to_db = 1; +dsp_options.max_occurences = 15; +dsp_options.database_path = 'C:\Users\Silas\Documents\MATLAB\Datensätze\sioe_labor\'; +dsp_options.database_name = 'silas_labor_newdsp_newstructure.db'; +dsp_options.storage_path = 'Z:\2024\sioe_labor\'; +dsp_options.parameters = struct(); +dsp_options.parameters.mu_dc = [0.005]; + +% === Get Run ID's === +db = DBHandler("pathToDB", [dsp_options.database_path, dsp_options.database_name], "type", "sqlite"); +fp = QueryFilter(); +% fp.where('Runs', 'run_id','EQUALS', 5108); +fp.where('Runs', 'is_mpi','EQUALS', 0); +% fp.where('Runs', 'fiber_length','EQUALS', 10); +fp.where('Runs', 'wavelength','EQUALS', 1310); +fp.where('Runs', 'db_mode','EQUALS', 1); +fp.where('Runs', 'rop_attenuation','EQUALS', 0); +fp.where('Runs', 'pam_level','EQUALS', 4); +% fp.where('Runs', 'bitrate','EQUALS', 360e9); +% fp.where('Runs', 'power_pd_in','GREATER_THAN', 7); +[dataTable,~] = db.queryDB(fp, db.getTableFieldNames('Runs')); + +% === Initialize DataStorage === +wh = DataStorage(dsp_options.parameters); +wh.addStorage("ffe_package"); +wh.addStorage("mlse_package"); +wh.addStorage("vnle_package"); +wh.addStorage("dbtgt_package"); +wh.addStorage("dbenc_package"); + +% === RUN IT === + +% [results,wh] = submitJobs(dataTable.run_id(:), dsp_options, "parallel", 'wh', wh, 'waitbar', true); +% wh.getStoValue('ffe_package',0.005); +% wh.getStoValue('mlse_package',0.005); + +[dataTable,~] = db.queryDB(fp, [db.getTableFieldNames('Runs');db.getTableFieldNames('Results');db.getTableFieldNames('Equalizer')]); + +dataTable = cleanUpTable(dataTable); + +% === Look at it === +y_var = 'BER_precoded'; +x_var = 'bitrate'; +fixedVars = {'equalizer_structure', x_var}; + +[dataTableClean, outliersTable] = removeGroupOutliers(dataTable, fixedVars, y_var); + +% --- Group and aggregate --- +dataTableGrpd_mean = groupIt(fixedVars, dataTableClean, @mean); +dataTableGrpd_min = groupIt(fixedVars, dataTableClean, @min); +dataTableGrpd_max = groupIt(fixedVars, dataTableClean, @max); + +% Choose a color map +cols = linspecer(numel(unique(dataTableGrpd_mean.equalizer_structure))); + +figure; +hold on; + +% Get unique equalizer structures for grouping +unique_eq = unique(dataTableGrpd_mean.equalizer_structure); + +for i = 1:numel(unique_eq) + eq_val = unique_eq(i); + + % Filter grouped data for this equalizer structure + filt = dataTableGrpd_mean.equalizer_structure == eq_val; + + x = dataTableGrpd_mean.(x_var)(filt); + y_mean = dataTableGrpd_mean.(y_var)(filt); + y_min = dataTableGrpd_min.(y_var)(filt); + y_max = dataTableGrpd_max.(y_var)(filt); + + % Bounds for boundedline (distance from mean) + y_lower = y_mean - y_min; + y_upper = y_max - y_mean; + y_bounds = [y_lower, y_upper]; + + % --- Bounded line (mean ± min/max) --- + if exist('boundedline', 'file') + [hl, hp] = boundedline(x, y_mean, y_bounds, ... + 'alpha', 'transparency', 0.1, ... + 'cmap', cols(i,:), ... + 'nan', 'fill', ... + 'orientation', 'vert'); + set(hl, 'LineWidth', 1.2, 'DisplayName', sprintf('Eq %s', eq_val)); + set(hp, 'HandleVisibility', 'off'); + else + % If boundedline is not available, use errorbar + errorbar(x, y_mean, y_lower, y_upper, ... + 'o-', 'Color', cols(i,:), 'LineWidth', 1.2, ... + 'DisplayName', sprintf('Eq %d', eq_val),'HandleVisibility', 'off'); + end + + % --- Normal line (mean only) --- + plot(x, y_mean, '-', 'Color', cols(i,:), 'LineWidth', 1.5, ... + 'DisplayName', sprintf('Mean Eq %s', eq_val),'HandleVisibility', 'off'); + + % --- Scatter plot for individual points (from original data) --- + % Filter original data for this group + orig_filt = dataTableClean.equalizer_structure == eq_val; + x_scatter = dataTableClean.(x_var)(orig_filt); + y_scatter = dataTableClean.(y_var)(orig_filt); + + scatter(x_scatter, y_scatter, 10,cols(i,:), 'filled', ... + 'MarkerFaceAlpha', 0.5, 'DisplayName', sprintf('Scatter Eq %s', eq_val),'HandleVisibility', 'off'); +end + +yline([2.2e-4,4.85e-3,2e-2],'HandleVisibility', 'off','LineWidth',1,'LineStyle','--'); +set(gca, 'YScale', 'log'); % BER is usually plotted log-scale +xlabel(x_var, 'Interpreter', 'none'); +ylabel(y_var, 'Interpreter', 'none'); +legend('show', 'Location', 'best'); +grid on; +title(sprintf('%s vs. %s', y_var, x_var), 'Interpreter', 'none'); +hold off; \ No newline at end of file