Wave Form ( End DCF )  ( wav ) 
( T=7.5ps,  40Gb/s,  5000/50km )
 1995.10.6  by K.Shimoura

1998. 8. 21   last  revised
2001. 9. 15  slightly revised

(*-----Initial Conditions-----------------------------*)
(* LD parameters *)
tfw=7.5    ;Pm =10.0   ;GHz=40     ;
wl =1552   ;dnm=1.6    ;ch =+0.4   ;

(* Fiber parameters *)
s  =50.0   ;aa =0.20   ;Sl =0.07   ;Rt =5.0    ;
Dav=+0.03  ;dDp=0.0    ;kd =5      ;
Dc =+30    ;Nc =1      ;n2 =2.24 10^-20        ;

(* Amplification/Receiver parameters *)
ka =10     ;NF =+5.0   ;fw =6      ;
Rfil=1.4   ;Rthr=0.3   ;Rdec=0.6   ;

(* Calculation parameters *)
h  =5000   ;dt =1      ;NM =32     ;ym =1.0    ; 
kmax=1000  ;kw =100    ;kgt=20     ;geta=0     ;

(* Graphics Array *)
Array[   wav1,Round[kmax/kgt]];
Array[ snwav1,Round[kmax/kgt]];
Array[ pwwav1,Round[kmax/kgt]];
Print[ Date[ ]];

Lmax= kmax + geta;
dm  = Round[1000/GHz/dt];
M   = 2 Quotient[NM dm,2];

(*----Pulse shape definition--------------------------*)

(* Hyperbolic-secant 
kt =1.763 ; t0= tfw /kt ;  

EE0=Table[ Sech[ t dt/t0] ,{t,-M/2,M/2-1} ];
Chp=Table[ Exp[-((ch/2)(t dt/t0)^2 )I ],{t,-M/2,M/2-1}];
E0 =EE0 * Chp;
*)

(* Hyperbolic-secant Series 
kt =1.763 ; t0= tfw /kt ;  

E0= Table[ Sech[(t+4.5dm) dt/t0]  *0
		  +Sech[(t+3.5dm) dt/t0]  *1
          +Sech[(t+2.5dm) dt/t0]  *0
		  +Sech[(t+1.5dm) dt/t0]  *1
          +Sech[(t+0.5dm) dt/t0]  *1
          +Sech[(t-0.5dm) dt/t0]  *1
          +Sech[(t-1.5dm) dt/t0]  *0
          +Sech[(t-2.5dm) dt/t0]  *1
          +Sech[(t-3.5dm) dt/t0]  *1
          +Sech[(t-4.5dm) dt/t0]  *0
                                       ,{t,-M/2,M/2-1}];
*)

(* WDM of Hyperbolic-secant
kt =1.763 ; t0= tfw /kt ; 
dth=(2 Pi dt dnm)/8.04  ;
dn =dm/NM ; phd=0.2     ;

EE0=Table[ ( Sech[(t-1.5dm) dt/t0]*Exp[ phd Random[] I ] *1
            +Sech[(t-0.5dm) dt/t0]*Exp[ phd Random[] I ] *0
            +Sech[ t dt/t0]  *0
            +Sech[(t+0.5dm) dt/t0]*Exp[ phd Random[] I ] *1
            +Sech[(t+1.5dm) dt/t0]*Exp[ phd Random[] I ] *1)
                      
          *( Exp[-2t dth I ]  *1
            +Exp[-1t dth I ]  *1
            +      1          *0
            +Exp[+1t dth I ]  *1
            +Exp[+2t dth I ]  *1 )
                                     ,{t,-M/2,M/2-1} ];

SLOT=Table[ If[Abs[t dt]<=dm/2,1,1 ],{t,-M/2,M/2-1} ];
E0 =EE0;
*)


(* Chirped Gaussian *)  
kt= 1.665 ; t0= tfw /kt ; m= 1.0;  

E0= Table[(Exp[ -.5(1 + ch I)((t+15.5dm) dt/t0) ^(2m)] *0
          +Exp[ -.5(1 + ch I)((t+14.5dm) dt/t0) ^(2m)] *0
	      +Exp[ -.5(1 + ch I)((t+13.5dm) dt/t0) ^(2m)] *1
	      +Exp[ -.5(1 + ch I)((t+12.5dm) dt/t0) ^(2m)] *0
	      +Exp[ -.5(1 + ch I)((t+11.5dm) dt/t0) ^(2m)] *1
	      +Exp[ -.5(1 + ch I)((t+10.5dm) dt/t0) ^(2m)] *1
	      +Exp[ -.5(1 + ch I)((t+9.5dm) dt/t0)  ^(2m)] *0
	      +Exp[ -.5(1 + ch I)((t+8.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t+7.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t+6.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t+5.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t+4.5dm) dt/t0) ^(2m)]  *0
	      +Exp[ -.5(1 + ch I)((t+3.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t+2.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t+1.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t+0.5dm) dt/t0) ^(2m)]  *0
	      
	      +Exp[ -.5(1 + ch I)((t-0.5dm) dt/t0) ^(2m)]  *0
	      +Exp[ -.5(1 + ch I)((t-1.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t-2.5dm) dt/t0) ^(2m)]  *0
	      +Exp[ -.5(1 + ch I)((t-3.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t-4.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t-5.5dm) dt/t0) ^(2m)]  *0
	      +Exp[ -.5(1 + ch I)((t-6.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t-7.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t-8.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t-9.5dm) dt/t0) ^(2m)]  *1
	      +Exp[ -.5(1 + ch I)((t-10.5dm) dt/t0) ^(2m)] *0
	      +Exp[ -.5(1 + ch I)((t-11.5dm) dt/t0) ^(2m)] *1
	      +Exp[ -.5(1 + ch I)((t-12.5dm) dt/t0) ^(2m)] *1
	      +Exp[ -.5(1 + ch I)((t-13.5dm) dt/t0) ^(2m)] *1
	      +Exp[ -.5(1 + ch I)((t-14.5dm) dt/t0) ^(2m)] *0
	      +Exp[ -.5(1 + ch I)((t-15.5dm) dt/t0) ^(2m)] *0)
	      
	      *(+Exp[ +2t dth I ]  *0
            +Exp[ +1t dth I ]  *0
            +      1           *1
            +Exp[ -1t dth I ]  *0
            +Exp[ -2t dth I ]  *0)	                                              
                                        ,{t,-M/2,M/2-1} ];
	                                              
(* Chirped Gaussian 2  
kt= 1.665 ; t0= tfw /kt ; m= 1.0;  
E0= Table[ Exp[ -.5(1 + ch I)( t dt/t0) ^(2m) ]		
                                        ,{t,-M/2,M/2-1} ];
*)

(* Dark soliton 
kt =1.763 ; t0= tfw /kt ; 
EE0=Table[ Tanh[t dt/t0] * 
                Exp[-(t dt/(4t0) )^4] ,{t,-M/2,M/2-1}];
Chp=Table[ Exp[-((ch/2)(t dt/t0)^2 )I],{t,-M/2,M/2-1}];
E0 =EE0 * Chp; 
*)

(* Dispersion length -----------------------*)
Lt = kmax h/1000;
La = ka h/1000;
att= 10^( -aa h /40000);
Dp = Dav - ( 1000 Dc )/( h ka Nc );

Ld = t0^2 /Abs[ 5.31 10^-7 wl^2 Dp ];
Lq = t0^3 /Abs[ 2.82 10^-13 wl^4 Sl ];
Ls = 1.571 * Ld;
Print["Dp, dDp      = ",Dp," +- ",dDp,"( ps/nm/km )"];
Print["Ld, Lq       = ",Ld,"( km )  ",Lq,"( km ) "];
Print["Lt, Lt/Ls    = ",Lt,"( km )  ",Lt/Ls];
Print["La, La/Ld    = ",La,"( km )  ",La/Ld];
      
(* Initial width & power -------------------*)
uu = Table[ t ,{t,-M/2,M/2-1} ];
EET= Interpolation[ Re[ E0 ] ];
    Plot[ EET[t],{t,1,M},PlotRange-> {-ym,ym},
                  AxesLabel ->{"time","E-field"}];
AE0= (Abs[E0])^2;
EEU= Interpolation[ AE0 ];
TEU= Interpolation[ uu * AE0 ];
TEV= Interpolation[ uu^2 * AE0 ];

IEU= NIntegrate[ EEU[t],{t,1,M},AccuracyGoal->3];
ITU= NIntegrate[ TEU[t],{t,1,M},AccuracyGoal->3];
ITV= NIntegrate[ TEV[t],{t,1,M},AccuracyGoal->3];
Pow= IEU Pm dt GHz /(1000 NM);
trms= Sqrt[ (ITV/IEU) - (ITU/IEU)^2 ] * dt;

Print["W(FWHM),chirp = ",tfw,"( ps )  ",ch ];
Print["Peak Power    = ",Pm,"( mW ) ",
                     N[ 10 Log[10,Pm]],"( dBm )"];
Print["Average Power = ",Pow,"( mW ) ",
                        10 Log[10,Pow],"( dBm )"];

(* Dynamic soliton parameters --------------*)
P0 = 2.63 10^-28 wl^3 s Dp/(n2 tfw^2);
Ga = 10^(aa ka h/10000);
gam= Ga Log[Ga]/(Ga-1);
Pam= P0 * gam;
Ne = Sqrt[ Pm /P0 ];
Nt = Sqrt[ gam ];
Print["P0, Pa, gam   = ",P0,"( mW )  ",Pam,"( mW )  "
      ,gam];
Print["N-exp,N-th    = ",Ne,"   ",Nt ];
Print[" "];

(* Spectrum parameters ---------------------*)
F0 = InverseFourier[ E0 ] ;
tt = Table[ If[ t <= M/2,(t+M/2),(t-M/2)],{t,1,M}];
AF0= Part[ Abs[F0],tt ]^2 ; 
FFU= Interpolation[ Reverse[ 10 Log[10,AF0] ] ];
ff = Table[ If[ t <= M/2,(t-1),(t-M-1) ],{t,1,M} ];
fr = ff /(M dt);

(* Dispersion slope ------------------------*)
dl = 3.34 10^-6 ( wl^2 /(M dt) );
Dpl= Table[ Dp +Sl(t dl),{t,-M/2,M/2-1} ];
ltf= Table[ -t+M+1,{t,1,M} ];    
Ddf= Part[ Dpl,ltf ];
DDT= Interpolation[ Reverse[ Ddf ] ];

(* Dispersion compensation parameters-------*)
Print["Dc, Nc      = ",Dc," ( ps/nm ) ",Nc];
Print["Dav, Dp     = ",Dav,Dp," ( ps/nm/km )"];

(* Initial shape drawing--------------------*)
Plot[ DDT[t],{t,1,M},PlotRange-> {-2,+2},
               AxesLabel ->{"W.L.","D( ps/nm/km )"} ];
Plot[ EEU[t],{t,1,M},PlotRange-> {0,ym},
               AxesLabel ->{"time","Intensity"}];
Plot[ FFU[t],{t,1,M},PlotRange-> {-30,20},
               Frame -> True,GridLines -> Automatic];
               

(*---- Pulse propagation by S.S.F.------------------*)

Ei = E0;
Do[

If[ Mod[k,kgt] == 1 && k > geta,
   wav1[(k-geta-1)/kgt +1]= Ei ];
      
   Dpx= Dp + If[ EvenQ[ Quotient[k,kd]],+dDp,-dDp];
   
   gm = -0.0252 Dpx fr^2 +0.0673 Sl fr^3; 
   A  = att Exp[ -(gm h/2)I ];
   B  = -( 6.283*10^18 (n2/wl) (Pm/s) )I;

   F0 = InverseFourier[Ei];
   E1 = Fourier[A * F0];

(* Raman *)   
   AE1= (Abs[E1])^2;
   DE1= Append[ Rest[AE1],First[AE1] ];
   RAM= ( DE1 - AE1 ) * 0.001 Rt/dt ;
   B1 = Exp[ h B ( AE1 - RAM ) ];
   E2 = B1 * E1;
   F2 = InverseFourier[E2];
   
(* B1 = Exp[ h B (Abs[E1])^2 ];
   E2 = B1 * E1;
   F2 = InverseFourier[E2];  *)

   F3 = A * F2;
   E3 = Fourier[F3];


(*------Dispersion Compensation before Amplicication------*)
If[ Mod[k,ka*Nc] == 0,
   A3 = Exp[ 25.2 fr^2 Dc I ];
   E3 = Fourier[A3 * F3];
   F3 = InverseFourier[E3]     
   ];
   

(*------------Amplification and Filtering----------*)
If[ Mod[k,ka] == 0,
  
(* Filter *)
   ktf= 0.6931;
   fm = 1     ;
   fc = dnm/ dl;
   fwt= fw / dl;
   
   FIL= Table[ 
               +Exp[-ktf Abs[2(f+8fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f+7fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f+6fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f+5fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f+4fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f+3fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f+2fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f+ fc) /fwt]^(2fm)] *0
               
               +Exp[-ktf Abs[2f/fwt]^(2fm)] *1
               
               +Exp[-ktf Abs[2(f- fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f-2fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f-3fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f-4fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f-5fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f-6fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f-7fc) /fwt]^(2fm)] *0
               +Exp[-ktf Abs[2(f-8fc) /fwt]^(2fm)] *0
               
                                      ,{f,-M/2,M/2-1}];
   FIL= FIL / Max[FIL];
                                 
(*
(* Filter Shapes *) 
   FI = Interpolation[ FIL ];
   Plot[ FI[t],{t,1,M},PlotRange->{0,1},
             Frame-> True, GridLines-> Automatic ];
   FIP= Interpolation[ 10 Log[10,FIL] ];
   Plot[ FIP[t],{t,1,M},PlotRange->{-50,0},
             Frame-> True, GridLines-> Automatic ];
   dlr= dl * 50;  Print[ dlr,"  (nm/div)" ];
*)      
 
   F4 = F3 * Sqrt[ Part[FIL,tt] ];
   E4 = Fourier[ F4 ];

(* Amplifier *)
   dAp= 0.0;
   GE4= Interpolation[ (Abs[E4])^2 ];
   IE4= NIntegrate[ GE4[t],{t,1,M},AccuracyGoal->3];
   dmp= If[ EvenQ[ Quotient[k,ka] ]
                         ,10^(-dAp/10),10^(dAp/10)];
   AMP= IEU/IE4 * dmp;                        
   F5 = F4 * Sqrt[ AMP ];
   E5 = Fourier[ F5 ];

(* Noise *)
   ASE= 6.41 10^-8 * 10^(NF/10) * AMP * (1.0 GHz); 
   WNE= Table[ ( ASE/Pm ) * Exp[6.28 Random[] I],{M}];
   E3 = E5 + Sqrt[ WNE ];
   F3 = InverseFourier[ E3 ];  	     
   ];   
   

(*----Pulse Width----------------*)
If[ Mod[k,kgt] == 0 && k > geta,      
   
   AE3 = (Abs[E3])^2;
   Pmax= Max[ AE3 ];
   GE3 = Interpolation[ AE3 ];
      
(* FWHM *)   
Do[
   x1 = x;
   S1 = GE3[x1]  -Pmax/2;
   S2 = GE3[x1+1]-Pmax/2;
   If[ S1*S2<0 && S1<0, Break[] ]
,{x,1,M-1}];
   x1 = x1 -S1/(-S1+S2);
Do[
   x2 = M-x;
   S1 = GE3[x2]  -Pmax/2;
   S2 = GE3[x2+1]-Pmax/2;
   If[ S1*S2<0 && S1>0, Break[] ]
,{x,1,M-1}];
   x2 = x2 +S1/( S1-S2);
      
   Wid1= If[x2>x1, (x2-x1)*dt, 999];
   Wid2= (M+1-x1+x2) *dt;
   Wfwh= Min[ Wid1,Wid2 ];
   Wfwh= If[ Wfwh>dt, Wfwh, dt ];

 (*Print[ "Wfwh   =  ", Wfwh ]; *)      
   pwwav1[ (k-geta)/kgt ]= Wfwh /tfw;
   ];


(*---- Q-factor Calculation ----------------------------*)
If[ Mod[k,kgt] == 0 && k > geta,      
   
(* E-filter *)
   ktf= 0.6931;
   fwR= Rfil * GHz ( wl^2 10^-8 /2.998 );
   fwQ= fwR / dl; 
   FIE= Table[ Exp[ -ktf Abs[ 2f/fwQ ]^(2fm) ]
                                      ,{f,-M/2,M/2-1}];
   F4 = F3 * Sqrt[ Part[ FIE,tt ] ];
   E4 = Fourier[ F4 ];   
   
   AE4 = (Abs[ E4 ])^2;
   Ymax= Max[ AE4 ];
   GE4 = Interpolation[ AE4 ];
   
(* Auto Search ------------*)   
   ww  = AE4;   
   AEF = 0;
   Do[
      AEF = AEF + Take[ ww, dm ];
      ww  = Drop[ ww, dm ]
   ,{t,1,NM}];
   
   GEF = Interpolation[ AEF ];
   IX  = NIntegrate[ GEF[t],{t,1,dm},AccuracyGoal->3];
   AMAX= Max[ AEF ];
 (*Plot[ GEF[t]/Pmax,{t,1,dm},PlotRange->{0,NM*ym},
                    AxesLabel ->{"time","Intensity"}];*)
   Do[
      CEN = t;
      If[ Part[ AEF,t ] > 0.99 AMAX, Break[ ]]
   ,{t,1,dm} ];   
      
   ww  = AE4;   
   If[ CEN > dm/2, ww = Drop[ ww, Round[CEN -dm/2] ]
                 , ww = Drop[ ww, Round[CEN +dm/2] ]];
 
 (* eye-pattern -----------*)
   nn  = Table[ t,{t,1,M} ];
   cn1 = Table[ Mod[(t +dm -CEN), 2dm],{t,1,M} ];
   yey1= { Part[cn1,nn], Part[AE4,nn] };
   eye1= Transpose[ yey1 ];
   
   cn2 = Table[ Mod[(t -CEN), 2dm],{t,1,M} ];
   yey2= { Part[cn2,nn], Part[AE4,nn] };
   eye2= Transpose[ yey2 ];
   
   eye = Join[ eye1,eye2 ];
 (*ListPlot[ eye, Prolog->AbsolutePointSize[3]
                , PlotRange->{0,Ymax} 
                (*PlotJoined -> True*) ]; *)
 
 (* Q-factor --------------*)     
   iww= Interpolation[ ww ];
   
 (* energy  
   WD = Rdec * dm;
   PWD= Table[ NIntegrate[ iww[t],{t, (n-.5)dm -WD/2,
                       (n-.5)dm +WD/2}] ,{n,1,NM-1}];*)
                       
 (* level *)
   PWD= Table[ iww[(n-.5)dm ] ,{n,1,NM-1}];  
   PWE = Sort[ PWD ];
 (*Print[ PWE ];*)
   thr = Rthr * Max[ PWE ];
  
   For[ i=1, PWE[[i]] <thr, i++ ];
   zer = Take[ PWE, i-1];
   one = Drop[ PWE, i-1];
   n0  = Part[ Dimensions[ zer ],1] +0.001;
   n1  = Part[ Dimensions[ one ],1] +0.001;

   zrm= Apply[Plus,zer] /n0;
   zrd= (zer-zrm)^2;
   zrs= Sqrt[ Apply[Plus,zrd] /n0 ];
   
   oem= Apply[Plus,one] /n1;
   oed= (one-oem)^2;
   oes= Sqrt[ Apply[Plus,oed] /n1 ];
   
   soz = zrs + oes + 0.0001;
   Qfac= Abs[ (oem-zrm)/soz ];
   
 (*Print[ n0,"  ",zrm,"  ",zrs ];
   Print[ n1,"  ",oem,"  ",oes ];
   Print[ "Qfac=  ", Qfac ];
   Print[ "Rthr,thr =  ", Rthr,"  ",thr ];*)
 
   snwav1[ (k-geta)/kgt ]= Qfac;
   ];

(*
(*---Dispersion Compensation after Amplicication---------*)
If[ Mod[k,ka*Nc] == 0,
   E3 = Fourier[ A3 * F3 ]
   ];
*)

(*---- Pulse Shape Drawing ------------------------------*)
If[ Mod[k,kw] == 0 && k > geta,
   
   L   = k h /1000;
   AE3 = (Abs[ E3 ])^2;
   Pmax= Max[ AE3 ];
   GE3 = Interpolation[ AE3 ];
   IE3 = NIntegrate[ GE3[t],{t,1,M},AccuracyGoal->3];
   Po3 = IE3 Pm dt GHz /(1000 NM);

   AF3 = Part[ Abs[F3],tt ]^2;
   GF3 = Interpolation[ AF3 ];
   LF3 = Interpolation[ Reverse[ 10 Log[10,AF3] ] ];
   
   Print[" "];
   Print["L, L/Ls       = ",L,"( km )  ",L/Ls];
 (*Print["Peak Power    = ",Pmax * Pm,"( mW )  ",
                     N[ 10 Log[10,Pmax*Pm]],"( dBm )"];
   Print["Average Power = ",Po3,"( mW )  ",
                              10 Log[10,Po3]"( dBm )"];*)
      
   Plot[ GE3[t]/Pmax,{t,1,M},PlotRange->{0,ym},
                    AxesLabel ->{"time","Intensity"}];
   Plot[ LF3[t],{t,1,M},PlotRange->{-30,20},
                 Frame -> True,GridLines -> Automatic];
   ListPlot[ eye, Prolog->AbsolutePointSize[3]
                , PlotRange->{0,Ymax} 
                (*PlotJoined -> True*) ];              
  ];  

 Ei= E3
 
,{k,1,Lmax} ];

Array[ wav1,kmax/kgt]   >> C:\soliton\wav1.m;
Array[ snwav1,kmax/kgt] >> C:\soliton\snwav1.m;
Array[ pwwav1,kmax/kgt] >> C:\soliton\pwwav1.m;
Print[ Date[ ]]
(*----------------------------------------------------*)

(* T = 7.5  ps, 5000/50  km, fw = 6  nm, NF = +5 dB,
   Dav = +0.03, Pm = 10  dBm, Dc = +30, Nc = 2 *)

\(wav1\  = << C:\soliton\wav1.m;\)
Power3D  = Abs[ wav1 ^2];
\(ListPlot3D[\ 
          Power3D, \ \ PlotRange -> {0, 
              1.1}, 
	\ \ \ \ \ \ \ \ \ \ \ BoxRatios -> {1, 2, 0.4}, 
          ViewPoint -> {1, \(-2\),  .4}\ ];\)
\(ListPlot3D[\ 
        Log[10, Power3D], \ \ PlotRange -> {\(-5\), 
            1}, 
	\ \ \ \ \ \ \ \ \ \ \ BoxRatios -> {1, 2, 0.4}, 
        ViewPoint -> {1, \(-2\),  .4}\ ];\)

\(
\(pwwav1\  = << C:\soliton\pwwav1.m;\)
\)
\(\(ListPlot[\ pwwav1, \ 
        PlotRange -> {0.8, 1.5}, 
	\ \ \ \ \ \ \ \ \ PlotJoined\  -> 
          True, \ AxesLabel\  -> {"Lt", "t(FWHM)"}\ ];\)\)

\(
\(snwav1\  = << C:\soliton\snwav1.m;\)
\)
Qfactor  =  Sqrt[ snwav1 ];
\(\(ListPlot[\ Qfactor, \ 
        PlotRange -> {0, 
            20}, 
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ PlotJoined\  -> 
          True, \ AxesLabel\  -> {"Lt", "Q-fac"}\ ];\)
 (*ListPlot[\ 
        Log[10,  .5  Erfc[Qfactor/1.414]\ ], \ 
        PlotRange -> {0, \(-30\)}, 
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ PlotJoined\  -> True, \ AxesLabel\  -> {"Lt", "BER"}\ ]; 

      ListPlot[\ 20\ Log[10, Qfactor]\ , \ 
        PlotRange -> {10, 30}, 
	\ \ \ \ \ \ \ \ PlotJoined\  -> True, \ 
        AxesLabel\  -> {"Lt", "20 LogQ"}\ ];*) 
\)

\(\(wav1\  = << C:\soliton\wav1.m;\)
\)
Power3D  = Abs[ wav1 ^2];
\(\(ListContourPlot[\ \(-Power3D\), \ PlotRange -> {0, \(-1.4\)}, 
	
        Contours -> {0, \(-0.2\), \(-0.4\), \(-0.6\), \(-0.8\), \(-1\), \(-1.2\), \(-1.4\)}\ ]\ ;\)\)
\(\(ListContourPlot[\ \(-Log[10, Power3D]\), \ 
        PlotRange -> {\(-1\), \(+5\)}, 
	
        Contours -> {5, 4, 3, 2, 1, 
            0, \(-1\)}\ ];\)\ \ \ \ \ \ \ \)

\(
 (*\ 
      Filter\ Shapes\ *) \ 
\(FI\  = \ Interpolation[\ FIL\ ];\)\)
\(\(Plot[\ FI[t], {t, 1, M}, 
        PlotRange -> {0, 1}, 
\ \ \ \ \ \ \ \ \ \ \ \ \ Frame -> \ True, \ 
        GridLines -> \ Automatic\ ];\)\)
FIP =  Interpolation[ 10 Log[10, FIL] ];
\(\(Plot[\ FIP[t], {t, 1, M}, 
        PlotRange -> {\(-50\), 0}, 
\ \ \ \ \ \ \ \ \ \ \ \ \ Frame -> \ 
          True, \ GridLines -> \ Automatic\ ];\)\)
\(\(\(dlr = \ dl\ *\ 200\)\(;\)\(\ \ \)\(Print[\ 
      dlr, "  (nm/div)"\ ]\)\(
\)\)\)