PR filter viz

Zürich data
Refactoring/ Organization

Theory/Optical/Dispersion/dispersion_power_fading.m

@@ -17,7 +17,7 @@
 %% Fiber and wavelength parameters
 lambda0  = 1310e-9;   % Zero-dispersion wavelength (ZDW) [m]
 S0       = 0.08;      % Dispersion slope at ZDW [ps/(nm^2·km)]
-L        = 10000;      % Fiber length [m]
+L        = 8000;      % Fiber length [m]
 alpha_dB = 0;         % Attenuation [dB/m] (ignored here)
 
 %% Target null frequency

Theory/Optical/Dispersion/power_fading.m

@@ -0,0 +1,50 @@
+%% ============================================================
+%  Minimal IM/DD Power Fading Plot
+% ============================================================
+
+
+%% Fiber and system parameters
+lambda0 = 1309e-9;       % zero-dispersion wavelength [m]
+lambda  = 1370e-9;       % operating wavelength [m]
+S0      = 0.092;          % dispersion slope [ps/(nm²·km)]
+L       = 30e3;          % fiber length [m]
+c       = physconst('lightspeed');
+
+%% Derived quantities
+S0_si = S0 * 1e3;                                    % → s/m³
+D_lambda = (S0/4) * (lambda*1e9 - (lambda0*1e9)^4/(lambda*1e9)^3); % ps/(nm·km)
+D_si = D_lambda * 1e-6;                              % → s/m²
+b2   = -D_si * lambda^2 / (2*pi*c);                  % s²/m
+
+Dacc = D_lambda * L;
+fprintf('Accumulated Dispersion: %.2f ps/nm \n', Dacc / 1e3);
+
+%% Frequency grid
+f_max = 200e9;
+f = linspace(0, f_max, 5000);                        % [Hz]
+
+%% IM/DD transfer function (power fading)
+phi = 2*pi^2 * b2 * f.^2 * L;
+H = abs(cos(phi));
+
+%% Plot
+figure('Color','w');
+plot(f/1e9, 10*log10(H), 'LineWidth', 1,'Color','black');
+grid on; box on;
+xlabel('Frequency [GHz]');
+ylabel('Magnitude [dB]');
+% title(sprintf('IM/DD Power Fading: 10 km; 1275nm', lambda*1e9, L/1000),"Interpreter","latex");
+ylim([-20 0]);
+
+%% Mark analytic null frequencies up to order 5
+max_order = 3;
+for n = 0:max_order
+    f_null_n = sqrt( c*(2*n + 1)/(2*abs(D_si)*lambda^2*L) );
+    xline(f_null_n/1e9, '--', 'LineWidth', 1.2, ...
+        'Color', [0.1216, 0.4706, 0.7059]);
+    text(f_null_n/1e9, -15, sprintf('$f_{\\mathrm{null}, %d}=%.1f$ GHz', n, f_null_n/1e9), ...
+        'BackgroundColor', 'w', 'EdgeColor', 'k', 'Interpreter', 'latex', ...
+        'HorizontalAlignment', 'center', 'VerticalAlignment', 'middle');
+end
+
+% mat2tikz_improved('C:\Users\Silas\Documents\6971e0b65b380ca6d71c837f\02_IMDD_System\tikz\dispersion\power_fading.tikz')

Theory/Signal_Processing/partial_response_filter.m

@@ -0,0 +1,100 @@
+
+% DUOBINARY CLASSES
+% Define the filter taps
+h_ = {1,[1 1],[1 2 1],[1 3 3 1]};
+desc = {'$1$','$(1+D)$','$(1+D)^2$','$(1+D)^3$'};
+
+for i = 1:length(h_)
+    h = h_{i};
+    
+    [H, w] = freqz(h, 1, 1024, 1);
+    
+    figure(1);
+    hold on
+    plot(w, 10*log10(abs(H)), 'LineWidth', 1, 'DisplayName',desc{i}); %todo
+    xlabel('Normalized Frequency');
+    ylabel('Amplitude in dB');
+    grid on;
+    ylim([-20,10]);
+    
+end
+
+legend
+beautifyBERplot("logscale",false,"setmarkers",false);
+% mat2tikz_improved("C:/Users/Silas/Documents/6971e0b65b380ca6d71c837f/03_DSP_Techniques/tikz/duobinary_response.tikz")
+
+%%
+
+% 1+ alpha D 
+% Define the filter taps
+h_ = {[1, 0],[1 0.2],[1 0.4],[1 0.6],[1 0.8],[1 1]};
+
+for i = 1:length(h_)
+    h = h_{i};
+    
+    [H, w] = freqz(h, 1, 1024, 1);
+    
+    figure(2);
+    hold on
+    plot(w, 10*log10(abs(H)), 'LineWidth', 2,'DisplayName',sprintf('$\\alpha = %.1f$', h(end))); %todo
+    xlabel('Normalized Frequency');
+    ylabel('Amplitude in dB');
+    grid on;
+    ylim([-15,5])
+end
+legend
+beautifyBERplot("logscale",false,"setmarkers",false);
+% mat2tikz_improved("C:/Users/Silas/Documents/6971e0b65b380ca6d71c837f/03_DSP_Techniques/tikz/partial_response.tikz")
+
+%%
+
+useprbs = 0;
+M = 4;
+randkey = 2;
+fsym = 112e9;
+
+%%%%% PRBS Generation in correct shape for Modulation Format %%%%%%
+O = 19; %order of prbs
+N = 2^(O-1); %length of prbs
+[~,seed] = prbs(O,1); %initialize first seed of prbs
+bitpattern=[];
+
+if useprbs
+    for i = 1:log2(M)
+        [bitpattern(:,i),seed] = prbs(O,N,seed);
+    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
+
+Tx_bits = Informationsignal(bitpattern);
+
+Digi_Mod = PAMmapper(M,0);
+Symbols_tx = Digi_Mod.map(Tx_bits);
+Symbols_tx.fs = fsym;
+
+cnt = 1;
+
+Symbols = Symbols_tx;
+
+Symbols = Duobinary().precode(Symbols);
+Symbols = Duobinary().encode(Symbols);
+
+figure()
+histogram()
+
+% Symbols = Duobinary().decode(Symbols);
+
+
+% coeff = [1,0.5];
+% 
+% Symbols.signal = filter(coeff, 1, Symbols.signal);
+% Symbols.spectrum("fignum",2,"displayname",['coeff:',num2str(coeff)],"normalizeTo0dB",1);

projects/Messung_Zürich/check_ber_to_ndr_plot.m

@@ -0,0 +1,91 @@
+% load data points from BER plot C:/Users/Silas/Documents/MATLAB/imdd_simulation/Projects/Messung_Zürich/BER_Silas.fig
+% load data points from NDR plot C:/Users/Silas/Documents/MATLAB/imdd_simulation/Projects/Messung_Zürich/NDR_Silas.fig
+
+% use Functions/Metrics/TransmissionPerformance.m to calculate the NDR from the BER values
+
+% 09. April 2026 - während plasmonic JLT review prozess
+
+ber = get_xy_from_fig("C:/Users/Silas/Documents/MATLAB/imdd_simulation/Projects/Messung_Zürich/BER_Silas.fig");
+tp = TransmissionPerformance;
+
+% 1 = PAM2 VNLE ;3=PAM2 DBt
+% 2 = PAM4 VNLE ;4 PAM4 DBt
+
+for n = 1:4
+
+    bps = 2 - mod(n,2);
+
+    figure(3);hold on
+    plot(ber(n).X,ber(n).Y,'DisplayName','ber')
+    set(gca,'YScale','log');
+    xlim([96,256]);
+    ylim([1e-4, 0.5]);
+    beautifyBERplot;
+    
+    ndr = tp.calculateNetRate(ber(n).X.*bps,'BER',ber(n).Y);
+    
+    figure(4); hold on
+    plot(ber(n).X,ndr.STAIR.NetRate,'DisplayName','ber')
+    xlim([96,256]);
+    ylim([100,330]);
+    beautifyBERplot;
+
+    figure(5); hold on
+    plot(ber(n).X,ndr.KP4_hamming.NetRate,'DisplayName','ber')
+    xlim([96,256]);
+    ylim([100,330]);
+    beautifyBERplot;
+
+    figure(6); hold on
+    plot(ber(n).X,ndr.O_FEC.NetRate,'DisplayName','ber')
+    xlim([96,256]);
+    ylim([100,330]);
+    beautifyBERplot;
+end
+
+function data = get_xy_from_fig(figPath)
+%GET_XY_FROM_FIG Silently open a .fig file and extract X/Y data
+%
+% Usage:
+%   data = get_xy_from_fig("C:\path\to\myfigure.fig");
+%
+% Output:
+%   data(k).Type
+%   data(k).DisplayName
+%   data(k).X
+%   data(k).Y
+
+    if ~isfile(figPath)
+        error('File not found: %s', figPath);
+    end
+
+    % Make sure figures stay invisible
+    oldVisible = get(groot, 'DefaultFigureVisible');
+    cleanupVisible = onCleanup(@() set(groot, 'DefaultFigureVisible', oldVisible));
+    set(groot, 'DefaultFigureVisible', 'off');
+
+    % Open figure invisibly
+    fig = openfig(figPath, 'invisible');
+    cleanupFig = onCleanup(@() close(fig));
+
+    % Find all objects that have XData and YData
+    objs = findall(fig, '-property', 'XData', '-and', '-property', 'YData');
+
+    % Optional: reverse order so it is closer to plotting order
+    objs = flipud(objs);
+
+    data = struct('Type', {}, 'DisplayName', {}, 'X', {}, 'Y', {});
+
+    for k = 1:numel(objs)
+        data(k).Type = get(objs(k), 'Type');
+
+        if isprop(objs(k), 'DisplayName')
+            data(k).DisplayName = get(objs(k), 'DisplayName');
+        else
+            data(k).DisplayName = '';
+        end
+
+        data(k).X = get(objs(k), 'XData');
+        data(k).Y = get(objs(k), 'YData');
+    end
+end