Files
imdd_silas/Theory/Dissertation/dispersion_contour_bandwidth_lambda.m
2026-03-25 10:57:48 +01:00

184 lines
6.1 KiB
Matlab
Raw Permalink Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
%% ------------------------------------------------------------
% Contour plot: λ_null as function of bandwidth (f_target) and reach (L)
% ------------------------------------------------------------
% Parameters
lambda0 = 1310e-9; % [m]
S0 = 0.09; % [ps/(nm²·km)]
c = physconst('lightspeed');
% Sweep dimensions
f_targets = linspace(50e9, 130e9, 200); % [Hz] (x-axis)
L_values = linspace(0.5e3, 15e3, 200); % [m] (y-axis)
lambda_surface = zeros(numel(L_values), numel(f_targets));
Dacc_surface = zeros(numel(L_values), numel(f_targets));
% Outer loop over fiber length (since L must be scalar)
for iL = 1:numel(L_values)
L = L_values(iL);
[lambda_vec, Dacc_vec] = lambda_for_first_null_full(f_targets, L, lambda0, S0);
% Store the 1x absolute offset |lambda - lambda0|
lambda_surface(iL, :) = abs(lambda0 - lambda_vec);
Dacc_surface(iL, :) = Dacc_vec;
end
% Convert for plotting
lambda_surface_nm = lambda_surface * 1e9; % [nm]
L_km = L_values / 1000; % [km]
f_GHz = f_targets / 1e9; % [GHz]
%% Contour plot
figure('Color','w');
hold on;
% Define wavelength contour levels [nm]
% Focus on a clean range of 1x offset values
lambda_levels = unique([50:-10:30, 30:-5:5]);
% Contour plot
[C,h] = contourf(f_GHz, L_km, lambda_surface_nm, lambda_levels, ...
'LineWidth', 1.2, ...
'ShowText', 'off');
% Colormap: Modern Blue palette with light colors removed for visibility
cmap_full = cbrewer2('Blues', 40);
colormap(cmap_full(10:end, :));
clim([min(lambda_levels) max(lambda_levels)]);
cb = colorbar;
ylabel(cb, '$\Delta \lambda$ [nm]', 'Interpreter', 'latex');
% --- MANUAL TEXTBOX ANNOTATIONS ---
% Find placement along the first-null curve for each level
for i = 1:length(lambda_levels)
lvl = lambda_levels(i);
% Re-calculate the specific (f, L) curve for this delta-lambda
lambda_target = lambda0 - (lvl * 1e-9);
LHS = -( (S0*1e3) / 4 ) * (lambda_target - (lambda0^4)/(lambda_target^3)) * lambda_target^2;
const_val = (c*0.5) / LHS;
f_curve_GHz = linspace(min(f_GHz), max(f_GHz), 500);
L_curve_km = const_val ./ (f_curve_GHz * 1e9).^2 / 1000;
% Filter for points within the plot axes
in_bounds = find(L_curve_km >= min(L_km)*1.1 & L_curve_km <= max(L_km)*0.9 & ...
f_curve_GHz >= min(f_GHz)*1.1 & f_curve_GHz <= max(f_GHz)*0.9);
if ~isempty(in_bounds)
% Specific alternating pattern for weight to minimize overlapping
if i < 9
weight = 0.05;
else
weight = 0.05 + 0.1 * mod(i, 2);
end
idx = in_bounds(max(1, min(length(in_bounds), round(length(in_bounds) * weight))));
text(f_curve_GHz(idx), L_curve_km(idx), sprintf('%g nm', lvl), ...
'Color', 'k', 'BackgroundColor', 'w', 'Margin', 1.5, ...
'HorizontalAlignment', 'center', 'VerticalAlignment', 'middle', ...
'EdgeColor', 'k', 'FontSize', 9);
end
end
% Axis formatting
xlabel('Signal Bandwidth [GHz]', 'FontSize', 11);
ylabel('Fiber length [km]', 'FontSize', 11);
xticks(min(f_GHz):10:max(f_GHz));
yticks(min(L_km):2.5:max(L_km));
grid on; box on;
axis([min(f_GHz) max(f_GHz) min(L_km) max(L_km)]);
%% Optional: overlay accumulated-dispersion contours
if 0
hold on;
min_D = min(Dacc_surface(:), [], 'omitnan');
max_D = max(Dacc_surface(:), [], 'omitnan');
% Calculate 3 integer levels well within the data range
D_levels = unique(round(linspace(min_D*0.8, max_D*0.8, 3)));
[CS, h] = contour(f_GHz, L_km, Dacc_surface, D_levels, 'k--', 'LineWidth', 0.8);
clabel(CS, h, 'Color','k', 'FontSize',8);
end
%% Export
% Hier erzwingen wir die rote Colormap für pgfplots, damit mat2tikz es nicht blau exportiert!
% mat2tikz_improved("C:/Users/Silas/Documents/6971e0b65b380ca6d71c837f/02_IMDD_System/tikz/dispersion/dispersion_power_fading_contour2.tikz");
function [lambda_vec, Dacc_vec] = lambda_for_first_null_full(f_target, L, lambda0, S0)
% lambda_for_first_null_full (stable, single-branch + validity checks)
% --------------------------------------------------------------------
% Computes the wavelength(s) at which the first IM/DD fading null
% occurs at frequency/ies f_target using the full dispersion model:
%
% D(lambda) = (S0/4)*(lambda - lambda0^4 / lambda^3)
%
% Restricted to the NORMAL-dispersion branch (λ < λ0),
% and valid only in the O-band (12601360 nm).
%
% Inputs:
% f_target - scalar or vector of target null frequencies [Hz]
% L - fiber length [m]
% lambda0 - zero-dispersion wavelength (ZDW) [m]
% S0 - dispersion slope at ZDW [ps/(nm²·km)]
%
% Outputs:
% lambda_vec - wavelength(s) [m] where first null occurs (clamped to O-band)
% Dacc_vec - accumulated dispersion(s) [ps/nm] (NaN if out of valid range)
% --------------------------------------------------------------------
c = physconst('lightspeed');
S0_si = S0 * 1e3; % ps/(nm²·km) -> s/(m³)
% Define O-band boundaries (in meters)
lambda_min = 1255e-9;
lambda_max = 1361e-9;
% Force column vector
f_target = f_target(:);
N = numel(f_target);
lambda_vec = NaN(N,1);
Dacc_vec = NaN(N,1);
for k = 1:N
RHS = c * 0.5 / (f_target(k)^2 * L);
% Normal-dispersion branch (λ < λ0)
fun = @(lambda) -(S0_si/4).*(lambda - (lambda0^4)./(lambda.^3)).*lambda.^2 - RHS;
% Limit the search to [λ_min, λ0)
try
lambda_sol = fzero(fun, [lambda_min, lambda0 * 0.999]);
catch
% If the zero is not within bounds, skip this point
lambda_sol = NaN;
end
% Validate solution
if isnan(lambda_sol) || lambda_sol < lambda_min || lambda_sol > lambda_max
lambda_vec(k) = NaN;
Dacc_vec(k) = NaN;
continue
end
% Compute D(lambda) and accumulated dispersion
D_lambda = (S0_si/4) * (lambda_sol - (lambda0^4)/(lambda_sol^3)) / 1e-6; % ps/(nm·km)
Dacc_val = D_lambda * (L/1000); % ps/nm
% Sanity bound on dispersion (avoid unphysical > ±100 ps/nm)
if abs(Dacc_val) > 100
lambda_vec(k) = NaN;
Dacc_vec(k) = NaN;
else
lambda_vec(k) = lambda_sol;
Dacc_vec(k) = Dacc_val;
end
end
end