Asked by abed
on 3 Dec 2013

% 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;

Answer by Image Analyst
on 3 Dec 2013

Accepted Answer

You're calling sym_ssf() but the function is actually called symssf() with no underline.

Opportunities for recent engineering grads.

Apply Today
## 3 Comments

## Jos (10584) (view profile)

## Direct link to this comment

https://www.mathworks.com/matlabcentral/answers/108424-fix-the-mistake-line-91-function#comment_183213

## Ghislaine Flore Kabadiang ngon (view profile)

## Direct link to this comment

https://www.mathworks.com/matlabcentral/answers/108424-fix-the-mistake-line-91-function#comment_738399

## Image Analyst (view profile)

## Direct link to this comment

https://www.mathworks.com/matlabcentral/answers/108424-fix-the-mistake-line-91-function#comment_738577

Sign in to comment.