calc_tail_jd

 This function generates the series of rays to represent the spectral tail

 INPUTS :

0001 function [Npar0_tail,dNpar0_tail,P0_2piRp_tail] = calc_tail_jd(Npar0,dNpar0,P0_2piRp,Nparmax_tail,n_tail,P_tail,opt_tail,display_mode)
0002 %
0003 % This function generates the series of rays to represent the spectral tail
0004 %
0005 % INPUTS :
0006 %
0007 %
0008 %
0009 if nargin < 8,
0010     display_mode = 0;
0011 end
0012 %
0013 if nargin < 7,
0014     opt_tail = 1;%linear
0015 end
0016 %
0017 if imag(Npar0) ~= 0,% case where nz is prescribed instead of Npar - only nz is tailed out
0018     %
0019     Ny = real(Npar0);
0020     Npar0 = imag(Npar0);% nz
0021     %
0022 else
0023     Ny = NaN;
0024 end
0025 %
0026 if Npar0 < 0,% working in positive spectrum
0027     s = -1;
0028     Npar0 = - Npar0;
0029 else
0030     s = 1;
0031 end
0032 %
0033 fac = exp(-1);%overlap factor
0034 %
0035 facln = sqrt(-log(fac));
0036 %
0037 if dNpar0 == 0 || imag(opt_tail) == 1,% uniform dNpar0 in the tail
0038     dNpar0_new = abs(Nparmax_tail - Npar0)/(1+2*n_tail)/facln;
0039 else
0040     dNpar0_new = dNpar0;
0041 end
0042 %
0043 N_ini = Npar0 + dNpar0_new*facln;
0044 %
0045 if Nparmax_tail <= N_ini,% no tail is required as single pass absorption is ensured
0046     n_tail = 0;
0047 end
0048 %
0049 if n_tail == 0,% no tail
0050     %
0051     Npar0_tail = s*Npar0;
0052     dNpar0_tail = dNpar0;
0053     P0_2piRp_tail = P0_2piRp;
0054     %
0055     if ~isnan(Ny),
0056         Npar0_tail = Ny + 1i*Npar0_tail;
0057     end
0058     %
0059     return
0060 end
0061 %
0062 N_sep = linspace(N_ini,Nparmax_tail,n_tail + 1);
0063 %
0064 Npar0_tail = [Npar0,(N_sep(2:end) + N_sep(1:end-1))/2];% lobes positions
0065 %
0066 dNpar0_tail = [dNpar0_new,diff(N_sep)/(2*facln)];% lobes width
0067 %
0068 h_tail = real(opt_tail);% final power fraction of the tail
0069 %
0070 if isnan(P_tail),% power in first lobe follows tail evolution
0071     %
0072     P_tail = n_tail*(n_tail + h_tail)/(n_tail+1)^2;
0073     %
0074 %     P0_2piRp_tail = P0_2piRp*((1-h_tail)*(2*n_tail+1:-2:1)/(n_tail+1)^2 ...
0075 %       + h_tail*ones(1,n_tail+1)/(n_tail+1));% lobes power (linear distribution)
0076 end
0077 %
0078 P0_2piRp_tail = P0_2piRp*[(1-P_tail),...
0079     (1-h_tail)*P_tail*(2*n_tail-1:-2:1)/n_tail^2 ...
0080   + h_tail*P_tail*ones(1,n_tail)/n_tail];% lobes power (linear distribution)
0081 %
0082 if display_mode == 2,
0083     %
0084     nx = 101;
0085     xlim = [Npar0-2*dNpar0_new,Nparmax_tail+dNpar0_new];
0086     %
0087     xN = linspace(xlim(1),xlim(2),nx);
0088     %
0089     xdPdN_0 = exp(-(xN - Npar0).^2/dNpar0.^2)/(sqrt(pi)*dNpar0);
0090     %
0091     for in = 1:n_tail + 1,
0092         xdPdN(in,:) = P0_2piRp_tail(in)/P0_2piRp*exp(-(xN - Npar0_tail(in)).^2/dNpar0_tail(in).^2)/(sqrt(pi)*dNpar0_tail(in));
0093     end
0094     %
0095     [sxdPdN] = feval(@fcalc_tail_yp,{Npar0_tail,dNpar0_tail,P0_2piRp_tail,P0_2piRp},xlim,xN);
0096     %
0097     leg = {'Orig.','Tail'};
0098     %
0099     figure(1),clf
0100     %
0101     graph1D_jd(xN,xdPdN_0,0,0,'','','',NaN,xlim,NaN,'-','none','b',2,20+14i);
0102     graph1D_jd(xN,sxdPdN,0,0,'N_{||}','dPdN_{||} (norm.)','',leg,xlim,NaN,'-','none','r',2,20+14i,gca,0.9,0.7,0.7);
0103     graph1D_jd(xN,xdPdN,0,0,'','','',NaN,xlim,NaN,'-','none','m',0.5,20+14i);
0104     %
0105 end
0106 %
0107 Npar0_tail = s*Npar0_tail;
0108 %
0109 if ~isnan(Ny),
0110     Npar0_tail = Ny + 1i*Npar0_tail;
0111 end
0112 %