% This program simulates a single-channel fiber transmission link
% using the symmetrized split-step Fourier algorithm.
%
% written by Jong-Hyung Lee
clear all %=============================================
% Define Time Window and Frequency Window
%============================================= taum = 2000; dtau = 2*taum/2^11; tunit= 1e-12; % make time unit in psec
tau = (-taum:dtau:(taum-dtau))*tunit; fs = 1/(dtau*tunit); tl = length(tau)/2; w = 2*pi*fs*(-tl:(tl-1))/length(tau); % w=angular freq.
wst = w(2)-w(1); %============================================= % Define Physical Parameters
%============================================= c = 3e5; %[km/sec] speed of light
ram0 = 1.55e-9; %[km] center wavelength
k0 = 2*pi/ram0; n2 = 6e-13 ; %[1/mW]
gamm = k0*n2 ; %[1/(km*mW)]
alphaDB = 0.2 ; % [dB/km] Power Loss
alpha = alphaDB/(10*log10(exp(1))); %[1/km] Power Loss in linear scale
% Dispersion parameters (beta3 term ignored)
Dp = -2; % [ps/nm.km]
beta2 = -(ram0)^2*Dp/(2*pi*c); % [sec^2/km]
%============================================= % Define Input Signal
%============================================= % A single Gaussian pulse is assumed.
Po = 2; % [mW] initial peak power of signal source
C = 0; % Chirping Parameter
m = 1; % Super Gaussian parameter (m=1 ==> Gaussian)
t0 = 50e-12; %[sec] initial pulse width
at = sqrt(Po)*exp(-0.5*(1+1i*C)*(tau./t0).^(2*m)); % Input field in the
time domain a0 = fft(at(1,:)); af = fftshift(a0); % Input field in the frequency domain
%============================================= % Define Simulation Distance and Step Size
%============================================= zfinal = 100; %[km] propagation distance
pha_max = 0.01; %[rad] maximum allowable phase shift due to the
nonlinear operator % pha_max = h*gamma*Po (h = simulation step length)
h = fix(pha_max/(gamm*Po)); % [km] simulation step length
M = zfinal/h; % Partition Number
% Define Dispersion Exp. operator
% Dh = exp((h/2)*D^), D^=-(1/2)*i*sgnb2*P, P=>(-i*w)^2
Dh = exp((h/2)*(-alpha/2+(1i/2)*beta2*w.^2)); % %================================================%
% Propagation Through Fiber %
%================================================% % Call the subroutine, sym_ssf.m for the symmetrized split-step Fourier
method [bt,bf] = sym_ssf(M,h,gamm,Dh,af); % Preamplifier at the receiver
% Optical amplifier is assumed ideal (flat frequency response and no noise)
GdB = 20; % [dB] optical amplifier power gain
gainA = sqrt(10^(GdB/10)); % field gain in linear scale
rt = gainA*bt; % plot the received power signal
figure(1) plot(tau,abs(rt).^2,'r') function [to,fo] =symssf(M,h,gamma,Dh,uf0) % Symmetrized Split-Step Fourier Algorithm
% % ==Inputs==
% M = Simulation step number ( M*h = simulation distance )
% h = Simulation step
% gamma = Nonlinearity coefficient
% Dh = Dispersion operator in frequency domain
% uf0 = Input field in the frequency domain
% % ==Outputs==
% to = Output field in the time domain
% fo = Output field in the frequency domain
% % written by Jong-Hyung Lee for k = 1:M %=============================================================
% Propagation in the first half dispersion region, z to z+h/2
%============================================================= Hf = Dh.*uf0; %==========================================================
% Initial estimate of the nonlinear phase shift at z+(h/2)
%========================================================== % Initial estimate value
ht = ifft(Hf); % time signal after h/2 dispersion region
pq = ht.*conj(ht); % intensity in time
u2e = ht.*exp(h*1i*gamma*pq); %Time signal
%============================================================= % Propagation in the second Dispersion Region, z+(h/2) to z+h
%============================================================= u2ef = fft(u2e); u3ef = u2ef.*Dh; u3e = ifft(u3ef); u3ei = u3e.*conj(u3e); %========================================================
% Iteration for the nonlinear phase shift(two iterations)
%======================================================== u2 = ht.*exp((h/2)*1i*gamma*(pq+u3ei)); u2f = fft(u2) ; u3f = u2f.* Dh; u4 = ifft(u3f); u4i = u4.*conj(u4); u5 = ht.*exp((h/2)*1i*gamma*(pq+u4i)); u5f = fft(u5); uf0 = u5f.*Dh; u6 = ifft(uf0); u6i = u6.*conj(u6); %============================================================= % Maximum allowable tolerance after the two iterations
etol = 1e-5; if abs(max(abs(u6i))-max(abs(u4i)))/max(abs(u6i)) > etol disp('Peak value is not converging! Reduce Step Size'), break end %============================================================= end to = u6; fo = uf0;
MATLAB: Fix the mistake line 91 function
error
Best Answer