%
%   BESSPLT.M    Sample script file to plot some low-order Bessel Functions
%
%   This is a sample file to generate plots of the zero and first order Bessel 
%   functions - J0(x), J1(x)   and   Y0(x), Y1(x)
%             - I0(x), I1(x)   and   K0(x), K1(x)
%
%   File written by J. R. White, UMass-Lowell (Aug. 2003)
%

%
%   getting started
      clear all, close all,  nfig = 0;
%
%   setup independent variable, but don't evaluate at exactly zero since some
%   of the functions have a singular point at zero
      Nx1 = 201;   x1 = linspace(eps,10,Nx1);   % range for ordinary BF plots
      Nx2 = 201;   x2 = linspace(eps,4,Nx2);    % range for modified BF plots
%
%   evaluate ordinary Bessel functions
      J0 = besselj(0,x1);    Y0 = bessely(0,x1);
      J1 = besselj(1,x1);    Y1 = bessely(1,x1);
%
%   evaluate modified Bessel functions
      I0 = besseli(0,x2);    K0 = besselk(0,x2);
      I1 = besseli(1,x2);    K1 = besselk(1,x2);
%
%   now let's plot these curves
      nfig = nfig+1;   figure(nfig)
      subplot(2,2,1),plot(x1,[J0;Y0],'LineWidth',2),grid
      axis([0 10 -2 2]);
      gtext('J_0(x)'),gtext('Y_0(x)')
      gtext('J_0 and Y_0 Bessel Functions')
%
      subplot(2,2,3),plot(x1,[J1;Y1],'LineWidth',2),grid
      axis([0 10 -2 2]);
      gtext('J_1(x)'),gtext('Y_1(x)')
      gtext('J_1 and Y_1 Bessel Functions')
%
      subplot(1,2,2),plot(x2,[I0;I1;K0;K1],'LineWidth',2),grid
      axis([0 4 0 10]);
      gtext('I_0(x)'),gtext('I_1(x)'),gtext('K_0(x)'),gtext('K_1(x)')
      gtext('Modified Bessel Functions')
%
%   end of demo
