clear all
pi=3.14159;
%c=3*10^8;
j=sqrt(-1);
z=[0.1,0.06,0.04]; %光纤长度
gamma=[20,15,30]; %非线性系数
beta2=[-0.0065,0.001,0.008]; %二阶色散
beta4=[2.0*10^(-6),2*10^(-6),1*10^(-6)]; %四阶色散
p1=[1,0.5,0.2];
p2=[1,0.5,0.2]; %泵浦功率
ps=0.001; %信号功率
ppi=0.001; %闲频功率
f0=0.04; %零色散频率Hz
omega0=2*pi*f0; %零色散频率
% f1=c/(1.540*10^(-6))*10^(-12);%THz
% f2=c/(1.560*10^(-6))*10^(-12);%THz
% omega1=2*pi*f1;
% omega2=2*pi*f2; %两泵浦频率
% omegac=(omega1+omega2)/2;
fs=(-50:0.02:50); %信号频率THz,波长1428nm--1578nm
omegas=2*pi*fs;
n=length(fs);
bs0=0.03;
bi0=0.03; %信号与闲频的初始功率
% B0=1;
for m=1:3
for i=1:n
% lams(i)=c/fs(i)*10^(-12); %信号波长
dbeta(m,i)=beta2(m)*(omegas(i))^2+1/12*beta4(m)*(omegas(i))^4;
%dbeta(m,i)=beta2(m)*(omegas(i)-omega0)^2+1/12*beta4(m)*(omegas(i)-omega0)^4;
k(m,i)=dbeta(m,i)-gamma(m)*(p1(m)+p2(m))-1/2*beta2(m)*omegas(i)^2+1/24*beta4(m)*omegas(i)^4;
g(m,i)=sqrt(4*gamma(m)^2*p1(m)*p2(m)-1/4*k(m,i).^2);
bs(m,i)=(cosh(g(m,i)*z(m))+j*k(m,i)/(2*g(m,i))*sinh(g(m,i)*z(m)));
bi(m,i)=(j*2*sqrt(p1(m)*p2(m))/g(m,i))*sinh(g(m,i)*z(m));
end
end
B1=bs(1,:)*bs0+bi(1,:)*bi0;
G1=20*log10((abs(B1))/(bs0)); %一段的增益

B2=bs(1,:).*bs(2,:)*bs0+bi(1,:).*bi(2,:)*bi0;
G2=20*log10((abs(B2))/(bs0)); %两段级联的增益

B3=bs(1,:).*bs(2,:).*bs(3,:)*bs0+bi(1,:).*bi(2,:).*bi(3,:)*bi0;
G3=20*log10((abs(B3))/(bs0)); %三段级联的增益

figure(1)
plot(fs,G1,fs,G2,'-',fs,G3,'--')
legend('G1','G2','G3')

想得到g和光纤长度Z的关系曲线(z分别取3个值,gamma=20,beta=0.001,beta=2*10^(-6),p1=p2=1)
g和非线性系数gamma曲线(gamma取3个值,其他跟上面一样)
g和泵浦功率p的关系曲线(p1=p2取3个值,其他跟上面一样)
哪位大哥帮个忙,感激不尽~~~~~~~~~~~~