invsinc2.m

很多matlab的源代码

function b=invsinc2(x,tol)
% INVSINC2 Inverse sinc squared function.
%
%	C = INVSINC2(X,TOL) returns the inverse of the sinc squared function 
% 	TOL = required tolerance, TOL [Default: TOL = 1e-7]
%
%	NOTE: INVSINC2 requires |X| < 1 
%	NOTE: Value of C is returned over the mainlobe (0,1)
%
%	INVSINC2 (with no input arguments) invokes the following example:
%
%	% Let xa=sinc^2(0.25), Then invsinc2(xa) should equal 0.25.
%         >>xa = sinc2(0.25)
%         >>xb = invsinc2(xa)


% ADSP Toolbox: Version 2.0 
% For use with "Analog and Digital Signal Processing", 2nd Ed.
% Published by PWS Publishing Co.
%
% Ashok Ambardar, EE Dept. MTU, Houghton, MI 49931, USA
% http://www.ee.mtu/faculty/akambard.html
% e-mail: [e-mail]
% Copyright (c) 1998


if nargin==0,help invsinc2,disp('Strike a key to see results of the example')
pause,xa=sinc2(0.25),ya=invsinc2(xa),return,end

if nargin<2,tol=1e-7;end
if any(any(abs(x)>1)),error('invsinc2 requires |x| < 1'),break,end 

a=0*x;b=0.5+a;c=1+a;kmax=100;k=0;
while any(any(abs(x-sinc2(b))>tol)) & k<=kmax
k=k+1;
fb=sinc2(b);f=(x-fb);
c=b.*(f>0)+c.*(f<=0);
a=b.*(f<=0)+a.*(f>0);
b=0.5*(a+c);
end
if k>kmax,disp(['No convergence in ' num2str(kmax) ' iterations']),end
b=b.*(x>0 & x<1)+(x==0);