diff --git a/Classes/00_signals/Signal.m b/Classes/00_signals/Signal.m index fa4c44a..7e1c3c6 100644 --- a/Classes/00_signals/Signal.m +++ b/Classes/00_signals/Signal.m @@ -192,6 +192,19 @@ classdef Signal end + %% Add signals from one signal to another, the first object will sustain + function Difference = minus(X,y) + + if isa(y,'Signal') + Difference = X; + Difference.signal = X.signal - y.signal; + elseif isnumeric(y) + Difference = X; + Difference.signal = X.signal - y; + end + + end + function Product = times(X,y) if isa(y,'Signal') @@ -271,14 +284,22 @@ classdef Signal N = 2^(nextpow2(length(obj.signal))-8); [p_lin,w] = pwelch(obj.signal,hanning(N),N/2,N,obj.fs,"centered","power","mean"); - p_dbm = 10*log10(p_lin)+30; %dB to dBm in case of "power" + normalize = 1; + if normalize + p_lin = p_lin./ max(p_lin); + p_dbm = 10*log10(p_lin); %dB to dBm in case of "power" + ylab = "normalized to 0 dB"; + else + p_dbm = 10*log10(p_lin)+30; %dB to dBm in case of "power" + ylab = "Power (dBm)"; + end figure(options.fignum); % If figure does not exist, create new figure hold on plot(w.*1e-9,p_dbm,'DisplayName',options.displayname,'LineWidth',1); xlabel("Frequency in GHz"); %ylabel("Power/frequency (dB/Hz)"); - ylabel("Power (dBm)"); + ylabel(ylab); xlim([-obj.fs/2 obj.fs/2].*1e-9) edgetick = 2^(nextpow2(obj.fs*1e-9)); % xticks([-edgetick:16:edgetick]); @@ -362,7 +383,7 @@ classdef Signal end - %% + %% Normalize function obj = normalize(obj,options) arguments @@ -381,7 +402,7 @@ classdef Signal end - %% + %% Delay function [obj] = delay(obj,delay,options) arguments @@ -402,7 +423,7 @@ classdef Signal end - %% + %% Synchronize with reference function [obj,D,cuts] = tsynch(obj,options) % time sync and cut arguments @@ -446,6 +467,10 @@ classdef Signal end + function obj = filter(obj,a,b) + obj.signal = filter(a,b,obj.signal); + end + function er = extinctionratio(obj,fsym,M) histpoints = 1024; %% verticale resolution histpoints = floor(histpoints/2)*2+1; %% to have the eye digram centered around one point make the vertical resolution uneven diff --git a/Classes/01_transmit/ChannelFreqResp.m b/Classes/01_transmit/ChannelFreqResp.m index 592e932..7c3da7f 100644 --- a/Classes/01_transmit/ChannelFreqResp.m +++ b/Classes/01_transmit/ChannelFreqResp.m @@ -184,12 +184,12 @@ classdef ChannelFreqResp < handle % it could be helpful to clip at linear 1 (=to keep comp from attenuating the signal) % iH(abs(iH)<1) = 1.*exp(1j*angle(iH(abs(iH)<1))); - + if 0 figure(7);hold on;plot(fnew,20*log10(abs(iH))) end - % iH(1) is DC ---> iH(end) is High Freq. + % iH(1) is DC ---> iH(end) is High Freq. H_inv = [iH(1) iH fliplr(conj(iH)) conj(iH(1))]; H_inv = [iH(1) iH iH(end) fliplr(conj(iH)) conj(iH(1))]; diff --git a/Classes/01_transmit/PAMsource.m b/Classes/01_transmit/PAMsource.m index e72033c..1566ee0 100644 --- a/Classes/01_transmit/PAMsource.m +++ b/Classes/01_transmit/PAMsource.m @@ -10,6 +10,7 @@ classdef PAMsource randkey db_precode + db_encode mrds_code mrds_blocklength @@ -36,6 +37,7 @@ classdef PAMsource options.randkey = 0; options.db_precode = 0; + options.db_encode = 0; options.mrds_code = 0; options.mrds_blocklength = 512; @@ -95,10 +97,16 @@ classdef PAMsource symbols = PAMmapper(obj.M,0).map(bits); symbols.fs = obj.fsym; + %%%%%% Duobinary %%%%%%%%%%% if obj.db_precode symbols = Duobinary().precode(symbols); + end + + if obj.db_encode symbols = Duobinary().encode(symbols); end + + % figure(12);hold on;histogram(symbols.signal,'Normalization','probability'); if obj.mrds_code symbols = MRDS_coding("blocklength",obj.mrds_blocklength).encode(symbols); diff --git a/Classes/04_DSP/Coding/Duobinary.m b/Classes/04_DSP/Coding/Duobinary.m index a9e76c1..e9888b9 100644 --- a/Classes/04_DSP/Coding/Duobinary.m +++ b/Classes/04_DSP/Coding/Duobinary.m @@ -47,10 +47,14 @@ classdef Duobinary % Pre coding bk = zeros(numel(data),1); - for k = 1:numel(data)-1 - bk(k+1) =mod(data(k)-bk(k),M); - end + % for k = 1:numel(data)-1 + % bk(k+1) =mod(data(k)-bk(k),M); + % end + for k = 2:numel(data) + bk(k) =mod(data(k)-bk(k-1),M); + end + %make bipolar bk = bk .* 2; bk = bk + b; @@ -107,9 +111,16 @@ classdef Duobinary assert(isequal((0:M-1)',unique(data)),'Check Duobinary Precoding'); %seems the signal is not unipolar % duobinary coding (1+D) + % coeff = [1,1]; + % + % data = conv(data,coeff,"same"); + coeff = [1,1]; - data = conv(data,coeff,"same"); + data = filter(coeff, 1, data); + %data = data ./ rms(data); + + % data = ifft(fft(real(data)).*fft(fliplr(coeff'),length(data))) + 1i*ifft(fft(imag(data)).*fft(fliplr(coeff'),length(data))); %make bipolar data = (data-round(mean(data),1)); diff --git a/Classes/04_DSP/Equalizer/FFE.m b/Classes/04_DSP/Equalizer/FFE.m index 111de61..fdd82e5 100644 --- a/Classes/04_DSP/Equalizer/FFE.m +++ b/Classes/04_DSP/Equalizer/FFE.m @@ -51,7 +51,7 @@ classdef FFE < handle end - function [X] = process(obj, X, D) + function [X,Noi] = process(obj, X, D) % actual processing of the signal (steps 1. - 3.) % 1 normalize RMS @@ -65,10 +65,10 @@ classdef FFE < handle obj.equalize(X.signal, D.signal,obj.mu_tr,obj.epochs_tr,obj.len_tr,training,showviz); % Decision Directed Mode - N = X.length; + n = X.length; training = 0; showviz = 0; - [signal,decision]=obj.equalize(X.signal, D.signal,obj.mu_dd,obj.epochs_dd,N,training,showviz); + [signal,decision]=obj.equalize(X.signal, D.signal,obj.mu_dd,obj.epochs_dd,n,training,showviz); % Output Signal if obj.decide @@ -76,9 +76,13 @@ classdef FFE < handle else X.signal = signal; end + X.fs = D.fs; %change sampling frequency of outgoing signal from fdac e.g. 2 sps to symbol spaced = fsym lbdesc = [num2str(obj.order),' tap FFE']; X = X.logbookentry(lbdesc); % append to logbook + + Noi = X; + Noi = X - D; end diff --git a/Classes/04_DSP/Sequence Detection/MLSE.m b/Classes/04_DSP/Sequence Detection/MLSE.m index 29d1f8e..aad4134 100644 --- a/Classes/04_DSP/Sequence Detection/MLSE.m +++ b/Classes/04_DSP/Sequence Detection/MLSE.m @@ -54,10 +54,12 @@ classdef MLSE obj.DIR(1:DIR_nonzero(1)-1) = []; end - if length(obj.DIR) == 1 + if isscalar(obj.DIR) obj.DIR = [0 obj.DIR]; end + obj.DIR = flip(obj.DIR); + % RMS normalization of input data data_in = data_in ./ rms(data_in); @@ -95,10 +97,65 @@ classdef MLSE % i.e. match the rms values of data_in to noise_free_received if isreal(data_in) if obj.M == round(obj.M) - data_in = data_in * rms(noise_free_received,'all','omitnan'); + data_in = data_in * rms(noise_free_received(noise_free_received ~= inf),'all','omitnan'); end end + + + % + % %% Optimized + % % Preallocate and initialize variables + % data_out = NaN(size(data_in)); % output vector + % sum_path_metrics_res = zeros(length(states), length(data_in)); + % path_idx = zeros(length(states), length(data_in)); + % + % % Precompute repmat size + % num_states = length(states); + % num_signals = numel(noise_free_received); + % + % % First trellis path (initialize) + % sum_path_metrics = zeros(num_states, num_states); + % path_metrics = (abs(data_in(1) - noise_free_received)).^2; % Use broadcasting instead of repmat + % sum_path_metrics = sum_path_metrics + path_metrics; % Compute initial path metrics + % + % [sum_path_metrics_res(:,1), path_idx(:,1)] = min(sum_path_metrics, [], 2); % Best path for first step + % + % % Preallocate path_metrics and sum_path_metrics to avoid reallocating in each loop + % path_metrics = zeros(num_states, num_signals); + % + % % Loop over remaining trellis paths + % for n = 2:length(data_in) + % % Avoid reallocation of sum_path_metrics, reuse the same matrix and update + % previous_sum_path_metrics = sum_path_metrics_res(:,n-1).'; % Transpose once for broadcasting + % sum_path_metrics = repmat(previous_sum_path_metrics, num_states, 1); % Avoid dynamic resizing + % + % % Calculate path metrics using broadcasting + % path_metrics = (abs(data_in(n) - noise_free_received)).^2; % Avoid repmat + % + % % Update sum path metrics + % sum_path_metrics = sum_path_metrics + path_metrics; + % + % % Find the best path for each state and store results + % [sum_path_metrics_res(:,n), path_idx(:,n)] = min(sum_path_metrics, [], 2); + % end + % + % %% Traceback + % ideal_path = NaN(1, length(data_in)+1); % Preallocate ideal path + % [~, ideal_path(length(data_in)+1)] = min(sum_path_metrics_res(:,length(data_in))); % Start from final state + % + % % Trace back through trellis + % for h = length(data_in):-1:1 + % ideal_path(h) = path_idx(ideal_path(h+1),h); % Follow the best path back + % end + % + % % Extract the output indices + % idx_out = ideal_path(2:end); + + + + + % OLD % initilaize the output vector data_out = NaN(size(data_in)); @@ -132,12 +189,17 @@ classdef MLSE ideal_path(h) = path_idx(ideal_path(h+1),h); end - idx_out = ideal_path(1:length(data_in)); + idx_out = ideal_path(2:length(data_in)+1); + + + + if obj.duobinary_output %use duobinary encoder, output is already scaled inside this %one data_out = Duobinary().encode(first_sym(idx_out)); + else % data_out(1:length(data_in)) = first_sym(idx_out); diff --git a/Tests/01_transmit/PAMmapper_test.m b/Tests/01_transmit/PAMmapper_test.m index 70b37d8..05b0b72 100644 --- a/Tests/01_transmit/PAMmapper_test.m +++ b/Tests/01_transmit/PAMmapper_test.m @@ -1,61 +1,68 @@ classdef PAMmapper_test < matlab.unittest.TestCase - + properties - Digi_Mod - Tx_bits - Rx_bits - % M = 8; % Modulation order - fsym = 112e9; % Symbol rate - + Tx_bits + PAMmapper + end + + properties (MethodSetupParameter) + % Define method-level parameters for PRBS and bit pattern + useprbs = struct('false', 0, 'true', 1); % variations: {0, 1} + M = struct('M2', 2, 'M4', 4, 'M6', 6, 'M8', 8); % variations: {2, 4, 6, 8} + O = struct('O10', 10, 'O15', 15, 'O17', 17); % variations: {10, 15, 17} end properties (TestParameter) - M = {4,6,8}; + % Define test-level parameters for M and O + end - + methods (TestMethodSetup) - function setupModulation(testCase) - % Setup PRBS and bit pattern for the test case - O = 10; % Order of PRBS + + function setupModulation(testCase, useprbs, M, O) + % Setup PRBS and bit pattern for the test case using parameters + + % Setup PRBS parameters N = 2^(O-1); % Length of PRBS - useprbs = 1; % Flag to use PRBS randkey = 1; % Random key for random stream [~, seed] = prbs(O, 1); % Initialize first seed of PRBS bitpattern = []; if useprbs - for i = 1:log2(testCase.M) + for i = 1:log2(M) [bitpattern(:, i), seed] = prbs(O, N, seed); end else s = RandStream('twister', 'Seed', randkey); - for i = 1:log2(testCase.M) + for i = 1:log2(M) bitpattern(:, i) = randi(s, [0 1], N, 1); end end - if testCase.M == 6 + if M == 6 bitpattern = reshape(bitpattern, [], 1); bitpattern = bitpattern(1:end-mod(length(bitpattern), 5)); end - + testCase.Tx_bits = Informationsignal(bitpattern); - testCase.Digi_Mod = PAMmapper(testCase.M, 0); + testCase.PAMmapper = PAMmapper(M, 0); end + end - - methods (Test) - function testBackToBackMapping(testCase) + + methods (Test, ParameterCombination = 'sequential') + % Test with sequential combination of parameters + function testBackToBackMapping(testCase, M, useprbs, O) % Test Bits -> Symbols -> Bits (Back-to-Back Mapping) % Map bits to symbols - Symbols = testCase.Digi_Mod.map(testCase.Tx_bits); - + Symbols = testCase.PAMmapper.map(testCase.Tx_bits); + % Demap symbols back to bits - testCase.Rx_bits = testCase.Digi_Mod.demap(Symbols); - + Rx_bits = testCase.PAMmapper.demap(Symbols); + % Validate BER is zero - [~, error_num, ber, ~] = calc_ber(testCase.Tx_bits.signal, testCase.Rx_bits.signal, ... + [~, error_num, ber, ~] = calc_ber(testCase.Tx_bits.signal, Rx_bits.signal, ... "skip_front", 0, "skip_end", 0, "returnErrorLocation", 1); % Assert that BER is zero @@ -63,5 +70,5 @@ classdef PAMmapper_test < matlab.unittest.TestCase testCase.verifyEqual(error_num, 0, 'No errors should occur in the mapping process'); end end - + end diff --git a/Tests/04_DSP/Coding/Duobinary_test.m b/Tests/04_DSP/Coding/Duobinary_test.m index c6be6b4..33de45e 100644 --- a/Tests/04_DSP/Coding/Duobinary_test.m +++ b/Tests/04_DSP/Coding/Duobinary_test.m @@ -1,11 +1,108 @@ classdef Duobinary_test < matlab.unittest.TestCase - % Test class for Duobinary + + properties + %genereal properties of the test + ber + errors + errorIndice + numbits + numsymbols + end + + properties (MethodSetupParameter) + % Define method-level parameters for PRBS and bit pattern + % i.e.: useprbs = {0,1}; + useprbs = struct('true', 1,'false', 0); % variations: {0, 1} + M = struct('M2', 2, 'M4', 4, 'M6', 6, 'M8', 8); % variations: {2, 4, 6, 8} + O = struct('O10', 10,'O11', 11,'O12', 12,'O13', 13, 'O15', 15, 'O17', 17); % variations: {10, 15, 17} + end + + properties (TestParameter) + % Define test-level parameters for M and O + + end + + methods (TestMethodSetup, ParameterCombination = 'pairwise') + + function setupTest(testCase, useprbs, M, O) + + randkey = 1; + datarate = 448e9; + fsym = round(datarate / log2(M)) ; + + %%%%% PRBS Generation in correct shape for Modulation Format %%%%%% + N = 2^O; %length of prbs + [~,seed] = prbs(O,1); %initialize first seed of prbs + bitpattern=[]; + + if useprbs + for i = 1:log2(M) + [bitpattern(:,i),seed] = prbs(O,N,seed); + end + else + s = RandStream('twister','Seed',randkey); + for i = 1:log2(M) + bitpattern(:,i) = randi(s,[0 1], N, 1); + end + end + + if M == 6 + bitpattern = reshape(bitpattern,[],1); + bitpattern = bitpattern(1:end-mod(length(bitpattern),5)); + end + + testCase.numbits = length(bitpattern); + + Tx_bits = Informationsignal(bitpattern); + + Symbols_tx = PAMmapper(M,0).map(Tx_bits); + Symbols_tx.fs = fsym; + + testCase.numsymbols = length(Symbols_tx); + + Symbols_tx = PAMmapper(M,0).map(Tx_bits); + + Symbols_tx.fs = fsym; + + Symbols = Duobinary().precode(Symbols_tx); + + Symbols = Duobinary().encode(Symbols); + + Symbols = Duobinary().decode(Symbols); + + Rx_bits = PAMmapper(M,0).demap(Symbols); + + [bits,testCase.errors,testCase.ber,testCase.errorIndice] = calc_ber(Tx_bits.signal,Rx_bits.signal,"skip_front",0,"skip_end",0,"returnErrorLocation",1); + + end + + end methods (Test) - function testExample(testCase) - % Example test case for Duobinary - % Add your test code here - testCase.verifyTrue(true); + % Test with sequential combination of parameters + + function checkBerZero(testCase, useprbs, M, O) + % Assert that BER is zero + testCase.verifyEqual(testCase.ber, 0, ['Numer of Errors is not zero, found ',num2str(testCase.errors), ' errors at position ', num2str(testCase.errorIndice),' out of ', num2str(2^O)]); + end + + function checkBerCloseToZero(testCase, useprbs, M, O) + + testCase.verifyLessThan(testCase.errors, 5, ['Found more than 5 errors. Found ',num2str(testCase.errors),' errors.']); + + % if ~isempty(testCase.errorIndice) + % testCase.verifyEqual(testCase.errorIndice, 2^O , ['Error is NOT at the End but at: ', num2str(testCase.errorIndice),' of ', num2str(2^O)]); + % testCase.verifyEqual(testCase.errorIndice,1, ['Error is NOT at the Beginning but at: ', num2str(testCase.errorIndice),' of ', num2str(2^O)]); + % end + + end + + function checkSystemFailure(testCase, useprbs, M, O) + testCase.verifyLessThan(testCase.ber, 0.1, ['BER is higher than 10%. Soemthing is completely wrong! (BER is ',num2str(testCase.ber),')']); + end + + end + end diff --git a/Tests/Testclass_template.m b/Tests/Testclass_template.m new file mode 100644 index 0000000..64058e2 --- /dev/null +++ b/Tests/Testclass_template.m @@ -0,0 +1,50 @@ +classdef Testclass_template < matlab.unittest.TestCase + + properties + %genereal properties of the test + end + + properties (MethodSetupParameter) + % Define method-level parameters for PRBS and bit pattern + % i.e.: useprbs = {0,1}; + useprbs = {0,1}; + M = {2,4,6}; + end + + properties (TestParameter) + % Define test-level parameters for M and O + + end + + methods (TestMethodSetup) + + function setupTest(testCase, useprbs, M) + % Setup everything that is somehow reguired for the test, + % include the MethodSetupParameters here as arguments + end + + end + + methods (Test) + % Test with sequential combination of parameters + function myTest(testCase, M, useprbs) + % Write the actual test, + % include the MethodSetupParameters here as arguments + + % Assert that BER is zero + testCase.verifyEqual(M, M, 'error message'); + testCase.verifyEqual(useprbs, useprbs, 'error message'); + end + + % Test with sequential combination of parameters + function anotherTest(testCase, M) + % Write the actual test, + % include the MethodSetupParameters here as arguments + + % Assert that BER is zero + testCase.verifyEqual(M, M, 'error message'); + testCase.verifyClass(M,'double'); + end + end + +end diff --git a/projects/400G_FTN_setups/duobinary_example.m b/projects/400G_FTN_setups/duobinary_example.m new file mode 100644 index 0000000..e69de29 diff --git a/projects/400G_FTN_setups/duobinary_example.mlx b/projects/400G_FTN_setups/duobinary_example.mlx deleted file mode 100644 index 1e310e7..0000000 Binary files a/projects/400G_FTN_setups/duobinary_example.mlx and /dev/null differ diff --git a/projects/AWG_characterization/output_power.m b/projects/AWG_characterization/output_power.m index 2fb812c..9fca42b 100644 --- a/projects/AWG_characterization/output_power.m +++ b/projects/AWG_characterization/output_power.m @@ -16,6 +16,7 @@ for i = 1:log2(M) [bitpattern(:,i),seed] = prbs(O,N,seed); end + if M == 6 bitpattern = reshape(bitpattern,[],1); bitpattern = bitpattern(1:end-mod(length(bitpattern),5)); diff --git a/projects/standard_system/400G_simulative_setup.mat b/projects/standard_system/400G_simulative_setup.mat index d60da59..f4a118f 100644 Binary files a/projects/standard_system/400G_simulative_setup.mat and b/projects/standard_system/400G_simulative_setup.mat differ diff --git a/projects/standard_system/imdd_minimal_2.m b/projects/standard_system/imdd_minimal_2.m index 1a5bf4d..2c6cddd 100644 --- a/projects/standard_system/imdd_minimal_2.m +++ b/projects/standard_system/imdd_minimal_2.m @@ -2,11 +2,17 @@ %% Parameter to simulate and save params = struct; -params.M = [4]; -params.datarate = [300]; -params.rop = [0]; +params.M = [8]; +params.datarate = [448]; +params.rop = [-12:0]; + +precomp_mode = 0; %0=do nothing ; 1= measure; 2=precomp active +postfilter = 0; % noise whiten. approach -> Postfilter + MLSE + +db_precode = 1; +db_encode = 1; +db_channelapproach = 0; -precomp_mode = 2; %0=do nothing ; 1= measure; 2=precomp active if ismac precomp_path = "/Users/silasoettinghaus/Documents/MATLAB/imdd_simulation/projects/standard_system"; @@ -24,7 +30,7 @@ wh = DataStorage(params); wh.addStorage("ber_ffe"); %% Init Params -link_length = 0; %meter +link_length = 1000; %meter pn_key = 2; laser_linewidth = 0; @@ -34,7 +40,6 @@ cnt=0; disp(['Start Simulation of ',num2str(endcnt),' loops...']) tic - for M = wh.parameter.M.values for datarate = wh.parameter.datarate.values @@ -50,32 +55,34 @@ for M = wh.parameter.M.values %%%%% Symbol Generation %%%%%% [Digi_sig,Symbols,Bits] = PAMsource("fsym",fsym,"M",M,"order",18,"useprbs",0,... - "fs_out",M8199.fdac,"applyclipping",1,"clipfactor",1.7,... + "fs_out",M8199.fdac,"applyclipping",1,"clipfactor",1.5,... "applypulseform",0,"pulseformer",Pform,"randkey",pn_key,... - "db_precode",0,... + "db_precode",db_precode,"db_encode",db_encode,... "mrds_code",usemrds,"mrds_blocklength",512).process(); - Digi_sig.eye(fsym,M); - Digi_sig.spectrum("fignum",11,"displayname",'before precomp'); + % Digi_sig.eye(fsym,M); + Digi_sig.spectrum("fignum",123434,"displayname",'Digital Tx Signal'); if precomp_mode == 1 freqresp = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs); Digi_sig = freqresp.buildOFDM(); elseif precomp_mode == 2 - Digi_sig = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs).precomp(Digi_sig,'maxampdb',3,'loadPath',precomp_path,'fileName',precomp_fn); + Digi_sig = ChannelFreqResp("Nacq",1024,"Navg",64,"Ncp",63,'f_ref',Digi_sig.fs).precomp(Digi_sig,'maxampdb',1,'loadPath',precomp_path,'fileName',precomp_fn); Digi_sig.spectrum("fignum",11,"displayname",'after precomp'); end %%%%% AWG %%%%%% El_sig = M8199.process(Digi_sig); + % El_sig.spectrum("displayname",'el','fignum',123434); + % El_sig.signal = awgn(El_sig.signal,-3,'measured',pn_key); %%%%% Lowpass el. components %%%%%% El_sig = Filter('filtdegree',2,"f_cutoff",60e9,"fs",fdac*kover,"filterType",filtertypes.butterworth,"active",true).process(El_sig); %%%%% Electrical Driver Amplifier %%%%%% - El_sig = Amplifier("amp_mode","ideal_no_noise","gain_mode","gain","amplification_db",9).process(El_sig); + El_sig = Amplifier("amp_mode","ideal_no_noise","gain_mode","gain","amplification_db",3).process(El_sig); % El_sig = El_sig.setPower(6,"dBm"); fprintf('Driver output power: %s dBm\n', num2str(El_sig.power)); @@ -90,9 +97,8 @@ for M = wh.parameter.M.values [Opt_sig] = EML("mode",eml_mode.im_cosinus,"power",3,"fsimu",El_sig.fs,"lambda",1290,"bias",vbias,"u_pi",u_pi,"linewidth",laser_linewidth,"randomkey",pn_key).process(El_sig); - - % Opt_sig.eye(fsym,M); - % % + % Opt_sig.eye(fsym,7); + % % figure(10) % hold on % scatter(El_sig.signal(1:100000)+vbias,(abs(Opt_sig.signal(1:100000)).^2)*1e3,0.1,'.','DisplayName','Modulator TF') @@ -104,6 +110,8 @@ for M = wh.parameter.M.values Optfilter = Filter('filtdegree',6,"f_cutoff",fsym.*0.7,"fs",fdac*kover,"filterType",filtertypes.gaussian,"active",true); Opt_sig = Optfilter.process(Opt_sig); + % Opt_sig.spectrum("fignum",122,"displayname",['Tx SPectrum; PAM ',num2str(M)]); + Opt_sig = Amplifier("amp_mode","ideal_no_noise","gain_mode","output_power","amplification_db",0).process(Opt_sig); i_ = wh.parameter.rop.length; @@ -143,19 +151,54 @@ for M = wh.parameter.M.values freqresp.plot(); end + Scpe_sig.spectrum("displayname",'After Scope','fignum',123434); + %%%%%% Sample to 2x fsym %%%%%% Scpe_sig = Scpe_sig.resample("fs_in",fadc,"fs_out",2*fsym); - %%%%%% Sync Rx signal with reference %%%%%% [Scpe_sig,D,cuts] = Scpe_sig.tsynch("reference",Symbols,"fs_ref",fsym); - %%%%% EQUALIZE %%%%%% - % - Eq = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",25,"sps",2,"decide",1); + Eq = FFE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",1e-4,"mu_tr",0,"order",25,"sps",2,"decide",0); %Eq = VNLE("epochs_tr",5,"epochs_dd",5,"len_tr",4096*2,"mu_dd",[0.0004 0.0005 0.0006],"mu_tr",0,"order",[50,7,7],"sps",2,"decide",1); - [EQ_sig] = Eq.process(Scpe_sig,Symbols); + + if db_channelapproach + % ref symbols and transm. sequence are precoded + [EQ_sig, Noi] = Eq.process(Scpe_sig,Duobinary().encode(Symbols)); + else + [EQ_sig, Noi] = Eq.process(Scpe_sig,Symbols); + end + + if db_encode || db_channelapproach + EQ_sig = MLSE("DIR",[1,1],"duobinary_output",1,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig); + EQ_sig = Duobinary().decode(EQ_sig); + end + + if postfilter + % Noi.spectrum("displayname",'Noise Spectrum','fignum',1234); + % EQ_sig.spectrum("displayname","Signal Spectrum","fignum",1234); + + nc = 2; + burg_coeff = arburg(Noi.signal,nc); + + EQ_sig = EQ_sig.filter(burg_coeff,1); + + % EQ_sig.spectrum("displayname","Signal Spectrum after Postfilter","fignum",1234); + tic + EQ_sig = MLSE("DIR",burg_coeff,"duobinary_output",0,"M",M,"trellis_states",PAMmapper(M,0).levels).process(EQ_sig); + toc + % EQ_sig.spectrum("displayname","Signal Spectrum after MLSE","fignum",1234); + + if 0 + [h,w] = freqz(1,burg_coeff,length(Noi),"whole",Noi.fs); + h = h/max(abs(h)); + hold on + w_ = (w - Noi.fs/2); + plot(w_.*1e-9,20*log10(fftshift(h)),'DisplayName',['', num2str(nc), ' coefficients for burg alg.']); + end + end + Rx_bits = PAMmapper(M,0).demap(EQ_sig); [~,errors_bm,ber_ffe(i),errors] = calc_ber(Rx_bits.signal,Bits.signal,"skip_front",100,"skip_end",150,"returnErrorLocation",1); disp(['BER: ',sprintf('%.1E',ber_ffe(i)),' - - ROP: ',num2str(patten(i)),'dBm - - PAM-',num2str(M),' - - ',num2str(fsym*1e-9),' GBd']); @@ -170,19 +213,21 @@ for M = wh.parameter.M.values end toc - disp(['Simulated: ',num2str(cnt/endcnt*100),' %']); + % wh.save('C:\Users\Silas\Documents\MATLAB\imdd_simulation\projects\MPI_August\auswertung\') end end +disp('Simulation Done!') + + cols = linspecer(8); %cnt = cnt+1; ber_ffe = wh.getStoValue('ber_ffe',M,datarate,wh.parameter.rop.values); - % Create the initial plot figure(44); diff --git a/projects/standard_system/sensitivity_analysis.m b/projects/standard_system/sensitivity_analysis.m index ed730a1..e88e059 100644 --- a/projects/standard_system/sensitivity_analysis.m +++ b/projects/standard_system/sensitivity_analysis.m @@ -2,6 +2,9 @@ M = wh.parameter.M.values(1); datarates = wh.parameter.datarate.values; rop = wh.parameter.rop.values; +clear("sens_ffe") + +col = flip(cbrewer2("RdYlGn",numel(datarates))); cnt = 1; for dr = datarates @@ -10,9 +13,9 @@ for dr = datarates sens_ffe(cnt) = getIntersection(ber_ffe,rop); - figure(112) + figure(116) hold on - plot(rop,ber_ffe,"LineWidth",0.5,"LineStyle","-","Marker",".","MarkerSize",15,"DisplayName",['FFE only; ',num2str(dr), ' Gbps']); + plot(rop,ber_ffe,"LineWidth",0.5,"LineStyle","-","Marker",".","MarkerSize",15,"DisplayName",['Postfilter + MLSE; ',num2str(dr), ' Gbps'],'Color',col(cnt,:)); yline(3.8e-3,'DisplayName','HD-FEC','LineStyle','--','HandleVisibility','off'); xlabel('Received Optical Power (dBm)'); ylabel('Bit Error Rate (BER)'); @@ -28,7 +31,7 @@ for dr = datarates end -figure(223); +figure(22); hold on; % Retain the plot so new points can be added without complete redraw plot(datarates,sens_ffe,"LineWidth",1,"LineStyle","-","Marker",".","MarkerSize",10,"DisplayName",['Linewidth: ',num2str(laser_linewidth.*1e-6),' MHz']); xlabel('Signal to Interference Ratio (dB)'); diff --git a/test/TestRand_examplecode.m b/test/TestRand_examplecode.m new file mode 100644 index 0000000..51ecb07 --- /dev/null +++ b/test/TestRand_examplecode.m @@ -0,0 +1,54 @@ +classdef TestRand_examplecode < matlab.unittest.TestCase + properties (ClassSetupParameter) + generator = {'twister','combRecursive','multFibonacci'}; + end + + properties (MethodSetupParameter) + seed = {0,123,4294967295}; + end + + properties (TestParameter) + dim1 = struct('small',1,'medium',2,'large',3); + dim2 = struct('small',2,'medium',3,'large',4); + dim3 = struct('small',3,'medium',4,'large',5); + type = {'single','double'}; + end + + methods (TestClassSetup) + function classSetup(testCase,generator) + orig = rng; + testCase.addTeardown(@rng,orig) + rng(0,generator) + end + end + + methods (TestMethodSetup) + function methodSetup(testCase,seed) + orig = rng; + testCase.addTeardown(@rng,orig) + rng(seed) + end + end + + methods (Test, ParameterCombination = 'sequential') + function testSize(testCase,dim1,dim2,dim3) + testCase.verifySize(rand(dim1,dim2,dim3),[dim1 dim2 dim3]) + end + end + + methods (Test, ParameterCombination = 'pairwise') + function testRepeatable(testCase,dim1,dim2,dim3) + state = rng; + firstRun = rand(dim1,dim2,dim3); + rng(state) + secondRun = rand(dim1,dim2,dim3); + testCase.verifyEqual(firstRun,secondRun) + end + end + + methods (Test) + function testClass(testCase,dim1,dim2,type) + testCase.verifyClass(rand(dim1,dim2,type),type) + end + end +end \ No newline at end of file diff --git a/test/duobinary_minimal_example.m b/test/duobinary_minimal_example.m index 0b7c407..bcd2de7 100644 --- a/test/duobinary_minimal_example.m +++ b/test/duobinary_minimal_example.m @@ -1,11 +1,12 @@ useprbs = 1; -M = 8; +M = 6; randkey = 1; -fsym = 112e9; +datarate = 448e9; +fsym = round(datarate / log2(M)) ; %%%%% PRBS Generation in correct shape for Modulation Format %%%%%% -O = 18; %order of prbs -N = 2^(O-1); %length of prbs +O = 17; %order of prbs +N = 2^(O); %length of prbs [~,seed] = prbs(O,1); %initialize first seed of prbs bitpattern=[]; @@ -34,12 +35,6 @@ Symbols = Duobinary().precode(Symbols_tx); Symbols = Duobinary().encode(Symbols); -% Symbols.spectrum("displayname","db","fignum",12) - -% Symbols = Duobinary().feedbackdetector(Symbols); - -% Symbols.eye(fsym,M); - Symbols = Duobinary().decode(Symbols); Rx_bits = PAMmapper(M,0).demap(Symbols); @@ -48,15 +43,35 @@ Rx_bits = PAMmapper(M,0).demap(Symbols); disp(['BER: ',sprintf('%.1E',ber),' - - PAM-',num2str(M)]); % -figure;hold on -subplot(2,1,1) -title('Bits Compare') -stairs(Rx_bits.signal(1:100,1),'LineWidth',2,'DisplayName','Rx Bits') -stairs(Tx_bits.signal(1:100,1),'LineStyle',':','LineWidth',2,'DisplayName','Tx Bits'); -legend -subplot(2,1,2) + + +figure(3); +clf +%sgtitle(['BER: ',num2str(ber),' // Error is at position: ',num2str(error_pos),'']) +subplot(2,2,1) hold on -title('Symbols Compare') +title('First Bits') +stairs(Tx_bits.signal(1:100,1),'LineStyle','-','LineWidth',2,'DisplayName','Tx Bits'); +stairs(Rx_bits.signal(1:100,1),'LineWidth',2,'DisplayName','Rx Bits','LineStyle',':') +legend + +subplot(2,2,2) +title('Last Bits') +hold on +stairs(Tx_bits.signal(end-50:end,1),'LineStyle','-','LineWidth',2,'DisplayName','Tx Bits'); +stairs(Rx_bits.signal(end-50:end,1),'LineWidth',2,'DisplayName','Rx Bits','LineStyle',':') +legend + +subplot(2,2,3) +hold on +title('First Symbols Compare') +stairs(Symbols_tx.signal(1:100,1),'LineWidth',2,'DisplayName','Tx Symbols','LineStyle','-') stairs(Symbols.signal(1:100,1),'LineStyle',':','LineWidth',2,'DisplayName','Rx Symbols'); -stairs(Symbols_tx.signal(1:100,1),'LineWidth',2,'DisplayName','Tx Symbols') +legend + +subplot(2,2,4) +hold on +title('Last Symbols Compare') +stairs(Symbols_tx.signal(end-50:end,1),'LineWidth',2,'DisplayName','Tx Symbols','LineStyle','-') +stairs(Symbols.signal(end-50:end,1),'LineStyle',':','LineWidth',2,'DisplayName','Rx Symbols'); legend \ No newline at end of file diff --git a/test/mapping_minimal_example.m b/test/mapping_minimal_example.m index b4f4f76..fdc86ac 100644 --- a/test/mapping_minimal_example.m +++ b/test/mapping_minimal_example.m @@ -30,8 +30,6 @@ end Tx_bits = Informationsignal(bitpattern); - - %%%%% ACTUAL TEST: Back to Back Mapping: Bits -> Symbols and Symbols -> Bits %%%%%% Digi_Mod = PAMmapper(M,0); @@ -41,14 +39,11 @@ Symbols = Digi_Mod.map(Tx_bits); Rx_bits = PAMmapper(M,0).demap(Symbols); - %%%%% VALIDATION: BER is required to be zero %%%%%% [~,error_num,ber,error_pos] = calc_ber(Tx_bits.signal,Rx_bits.signal,"skip_front",0,"skip_end",0,"returnErrorLocation",1); - - %%%%% For User: Show Debug Info and Results %%%%%% if viewresults diff --git a/test/mlse_minimal_example.m b/test/mlse_minimal_example.m index 506d94a..04e74db 100644 --- a/test/mlse_minimal_example.m +++ b/test/mlse_minimal_example.m @@ -1,6 +1,6 @@ -useprbs = 1; -M = 8; -randkey = 1; +useprbs = 0; +M = 4; +randkey = 2; fsym = 112e9; %%%%% PRBS Generation in correct shape for Modulation Format %%%%%% @@ -31,29 +31,50 @@ Digi_Mod = PAMmapper(M,0); Symbols_tx = Digi_Mod.map(Tx_bits); Symbols_tx.fs = fsym; -Symbols = Duobinary().precode(Symbols_tx); +if 0 -Symbols = Duobinary().encode(Symbols); + Symbols = Duobinary().precode(Symbols_tx); -Symbols = MLSE("DIR",[1,1],"duobinary_output",1,"trellis_states",Digi_Mod.levels,"M",M).process(Symbols); + Symbols = Duobinary().encode(Symbols); -Symbols = Duobinary().decode(Symbols); + Symbols = MLSE("DIR",[1 1],"duobinary_output",1,"trellis_states",Digi_Mod.levels,"M",M).process(Symbols); -Rx_bits = PAMmapper(M,0).demap(Symbols); + Symbols = Duobinary().decode(Symbols); -[~,error_num,ber,error_pos] = calc_ber(Tx_bits.signal,Rx_bits.signal,"skip_front",0,"skip_end",0,"returnErrorLocation",1); +else + cnt = 1; -disp(['BER: ',sprintf('%.1E',ber),' - - PAM-',num2str(M)]); -% -figure;hold on + Symbols = Symbols_tx; + + coeff = [1,0.5,0.2,0.1]; + + Symbols.signal = filter(coeff, 1, Symbols.signal); + + Symbols.spectrum("fignum",129,"displayname",['coeff:',num2str(coeff)]); + + Symbols = MLSE("DIR",coeff,"duobinary_output",0,"trellis_states",Digi_Mod.levels,"M",M).process(Symbols); + + Rx_bits = PAMmapper(M,0).demap(Symbols); + + [~,error_num,ber,error_pos] = calc_ber(Tx_bits.signal,Rx_bits.signal,"skip_front",0,"skip_end",0,"returnErrorLocation",1); + + disp(['BER: ',sprintf('%.1E',ber),' - - PAM-',num2str(M)]); + +end + + +% +figure(494) +clf subplot(2,1,1) title('Bits Compare') +hold on stairs(Rx_bits.signal(1:100,1),'LineWidth',2,'DisplayName','Rx Bits') stairs(Tx_bits.signal(1:100,1),'LineStyle',':','LineWidth',2,'DisplayName','Tx Bits'); legend subplot(2,1,2) hold on title('Symbols Compare') -stairs(Symbols.signal(1:100,1),'LineStyle',':','LineWidth',2,'DisplayName','Rx Symbols'); -stairs(Symbols_tx.signal(1:100,1),'LineWidth',2,'DisplayName','Tx Symbols') +stairs(Symbols.signal(1:100,1),'LineStyle','-','LineWidth',2,'DisplayName','Rx Symbols'); +stairs(Symbols_tx.signal(1:100,1),'LineWidth',2,'DisplayName','Tx Symbols','LineStyle',':') legend \ No newline at end of file diff --git a/test/test_modulation.m b/test/test_modulation.m index e080981..3c39bba 100644 --- a/test/test_modulation.m +++ b/test/test_modulation.m @@ -1,36 +1,74 @@ +classdef test_modulation < matlab.unittest.TestCase + properties + Tx_bits + Digi_Mod + end + properties (MethodSetupParameter) + % Define method-level parameters for PRBS and bit pattern + useprbs = struct('false', 0, 'true', 1); % variations: {0, 1} + M = struct('M2', 2, 'M4', 4, 'M6', 6, 'M8', 8); % variations: {2, 4, 6, 8} + O = struct('O10', 10, 'O15', 15, 'O17', 17); % variations: {10, 15, 17} + end + properties (TestParameter) + % Define test-level parameters for M and O -O = 10; %order of prbs -N = 2^(O-1); %length of prbs -[~,seed] = prbs(O,1); %initialize first seed of prbs -bitpattern=[]; -M=6; + end + + methods (TestMethodSetup) + + function setupModulation(testCase, useprbs, M, O) + % Setup PRBS and bit pattern for the test case using parameters + + % Setup PRBS parameters + N = 2^(O-1); % Length of PRBS + randkey = 1; % Random key for random stream + + [~, seed] = prbs(O, 1); % Initialize first seed of PRBS + bitpattern = []; + if useprbs + for i = 1:log2(M) + [bitpattern(:, i), seed] = prbs(O, N, seed); + end + else + s = RandStream('twister', 'Seed', randkey); + for i = 1:log2(M) + bitpattern(:, i) = randi(s, [0 1], N, 1); + end + end + + if M == 6 + bitpattern = reshape(bitpattern, [], 1); + bitpattern = bitpattern(1:end-mod(length(bitpattern), 5)); + end + + testCase.Tx_bits = Informationsignal(bitpattern); + testCase.Digi_Mod = PAMmapper(M, 0); + end + + end + + methods (Test, ParameterCombination = 'sequential') + % Test with sequential combination of parameters + function testBackToBackMapping(testCase, M, useprbs, O) + % Test Bits -> Symbols -> Bits (Back-to-Back Mapping) + + % Map bits to symbols + Symbols = testCase.Digi_Mod.map(testCase.Tx_bits); + + % Demap symbols back to bits + Rx_bits = testCase.Digi_Mod.demap(Symbols); + + % Validate BER is zero + [~, error_num, ber, ~] = calc_ber(testCase.Tx_bits.signal, Rx_bits.signal, ... + "skip_front", 0, "skip_end", 0, "returnErrorLocation", 1); + + % Assert that BER is zero + testCase.verifyEqual(ber, 0, 'BER should be zero'); + testCase.verifyEqual(error_num, 0, 'No errors should occur in the mapping process'); + end + end -for i = 1:log2(M) - [bitpattern(:,i),seed] = prbs(O,N,seed); end - -if M == 6 - bitpattern = reshape(bitpattern,[],1); -end - -% [bitpattern,seed] = prbs(O,N,seed); -% % -% bitpattern = bitpattern'; - -% 2 ) Build Inf. signal class -bits = Informationsignal(bitpattern); - - -% 3) Digi modulation -> PAM-M signal -digimod_out = PAMmapper(M,0).map(bits); - -Rx_Bits = PAMmapper(M,0).demap(digimod_out); - -[~,errors_bm,BER,errors] = calc_ber(Rx_Bits.signal,bitpattern,"skip_front",0,"skip_end",0,"returnErrorLocation",1); - -formatted_ber = sprintf('%.1e', BER); -disp(['BER: ',formatted_ber]); - \ No newline at end of file