MATLAB: How to plot Frequency domain of optical nonlinear Schrodinger equation

nonlinear optical fibernonlinear schrodingher equationsplit step method

Hello, This code solves the nonlinear in optical communication but, it only plots the time domain, I need to plot that frequency domain as well for both 3D and 2D plot. Thank you for help
if true
 clc; 
clear all; 
close all;
clf;
cputime=0;
tic;
p0=6;
alpha=0.00;%Fiber loss value in dB/km
gamma=0;%fiber non linearity in /W/m
t0=1e-10;%%total initial pulse width in second
C=0;%Chirp parameter
s=0;
b2=-20e-22;%second oder dispersion 
N=1;% soliton oder
m=0;
ln=1;
i=sqrt(-1);
pi=3.1415926535;
alph=alpha/(4.343); % Fiber loss equation, Ref page#55 eqn 2.5.3 (Agrawal)
Ld=(t0^2)/(abs(b2)); %%Dispersion
Ao=sqrt(p0); %power amplitude
tau =- 4096e-12:1e-12: 4095e-12;%  dt=t/to
dt=1e-12;
h1=1000;%%step size
for ii=0.1:0.1:1.0;
z=ii*Ld;
u1=Ao*sech(tau/t0).*exp(-i*C*(tau/t0).^2); % 
% u1=N*exp(-((1+i*(-C))/2)*(tau/t0).^2)
u2=Ao*sech(tau/t0);
%     u=N*sech(tau/to);%fundamental soliton pulse % first order soliton N=3
%     u=Ao*exp(-((1+i*(-C))/2)*(tau/to).^2) %%%N=1
%     u = Ao*sech(tau/to+1.5).*exp((-i*C/2).*(tau/to+1.5).^2)% first order with Chirp parameter
%%%%%%%%%%Plot Input Pulse%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
A11=u1;
A12=u2;
%    figure(1)
%    plot(tau/1e-10,abs(u1).^2,'r'); % Plot input pulse (with dispersion)
% %    axis([min(tau)/1e-9 max(tau)/1e-9,0 Ao])
grid on;
% figure(2)
%    plot(tau/1e-9,abs(u2).^2,'r'); % Plot input pulse (without pulse)
%    axis([min(tau)/1e-10 max(tau)/1e-10,0 Ao])
h=h1/Ld;%soliton conditions

Z=z/Ld;%soliton conditions
l=max(size(u1));  
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fwhm1=find(abs(u1)>abs(max(u1)/2));
fwhm1=length(fwhm1);
spectrum1=fft(fftshift(u1)); %Pulse spectrum

spectrum2=fft(fftshift(u2)); %Pulse spectrum
dw=(1/l)/dt*2*pi;
w=(-1*l/2:1:l/2-1)*dw;
 w=fftshift(w);
d=0;
for jj=h:h:Z
spectrum1=spectrum1.*exp(gamma*(h/2)+i*s/2*w.^2*(h/2)) ; 
f1=ifft(spectrum1);
f1=f1.*exp(i*(N^2)*((abs(f1)).^2)*(h));
% f=fftshift(f);

spectrum1=fft(f1);
spectrum1=spectrum1.*exp(gamma*(h/2)+i*s/2*w.^2*(h/2)) ; 
spectrum2=spectrum2.*exp(gamma*(h/2)+i*s/2*w.^2*(h/2)) ; 
f2=ifft(spectrum2);
f2=f2.*exp(i*(N^2)*((abs(f2)).^2)*(h));
% f=fftshift(f);
spectrum2=fft(f2);
spectrum2=spectrum2.*exp(gamma*(h/2)+i*s/2*w.^2*(h/2)) ; 
d=d+1;
end
f1=ifft(spectrum1);
f1=fftshift(f1);
op_pulse1(ln,:)=abs(f1);%saving output pulse at all intervals

fwhm=find(abs(f1)>abs(max(f1)/2));
fwhm=length(fwhm);
ratio=fwhm/fwhm1; %PBR at every value

pbratio1(ln)=ratio;%saving PBR at every step size

dd=atand((abs(imag(f1)))/(abs(real(f1))));
phadisp1(ln)=dd;%saving pulse phase

f2=ifft(spectrum2);
f2=fftshift(f2);
op_pulse2(ln,:)=abs(f2);%saving output pulse at all intervals
fwhm=find(abs(f2)>abs(max(f2)/2));
fwhm=length(fwhm);
ratio=fwhm/fwhm1; %PBR at every value
pbratio2(ln)=ratio;%saving PBR at every step size
dd=atand((abs(imag(f2)))/(abs(real(f2))));
phadisp2(ln)=dd;%saving pulse phase
ln=ln+1;
end
toc;
cputime=toc;
figure(3);
mesh(tau/1e-12,(1:1:ln-1)./12,op_pulse1(1:1:ln-1,:));
% axis([min(tau)/1e-12 max(tau)/1e-12,0 Ao])
axis([-500 500 0 1 0 inf]);
colormap([0,0,1]);
title('Pulse Evolution');
xlabel('Time (ps)'),ylabel('Distance(km)'); zlabel('Power (W)');
set(gca,'FontSize',24);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%frequency
% figure(110)
spectrum1=fft(fftshift(u1)).*(N*tau)/sqrt(2*pi); % spectrum
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Add optisystem data%%%%%%%%%%%
figure(10)
plot(tau/1e-12,abs(op_pulse2(ln-1,:)).^2,'r-',tau/1e-12,abs(A12).^2,':b','Linewidth',12);
axis([min(tau)/1e-12 max(tau)/1e-12,0 Ao]);
axis([-300,300 0 6]);
grid on;
title('Pulse Evolution');
xlabel('Time (ps) ');
ylabel('Power (W)');
set(gca,'FontSize',24);
legend ('Matlab IN','Matlab SPM-OUT');
legend boxon;
grid on;
  end

Best Answer

Unfortunately, I could not run and check your code.
Previously, I have implemented split-step Fourier method (SSFM) to solve the following nonlinear Schrodinger equation (NLSE) with 3rd order chromatic dispersion effect in optical fiber. Since I believe this code would be some help for your task, I would attach the code here.
If you run the SSFMtest.m, you can find basic soliton waveform before and after transmission in optical fiber where beta1, beta3 and alpha were set to 0.
By slightly modifying the code, you can plot spectrum of signal for any distance. But before going forward, let me clarify your definition of '3D plot'. Do you mean 'waterfall plot', like in the following link?
https://jp.mathworks.com/help/matlab/ref/waterfall.html
Related Question