Home > atmlab > mie > mie.m

mie

PURPOSE ^

Computation of Mie Efficiencies for given

SYNOPSIS ^

function result = mie(m, x)

DESCRIPTION ^

 Computation of Mie Efficiencies for given 
 complex refractive-index ratio m=m'+im" 
 and size parameter x=k0*a, where k0= wave number in ambient 
 medium, a=sphere radius, using complex Mie Coefficients
 an and bn for n=1 to nmax,
 s. Bohren and Huffman (1983) BEWI:TDD122, p. 103,119-122,477.
 Result: m', m", x, efficiencies for extinction (qext), 
 scattering (qsca), absorption (qabs), backscattering (qb), 
 asymmetry parameter (asy=<costeta>) and (qratio=qb/qsca).
 Uses the function "Mie_ab" for an and bn, for n=1 to nmax.
 C. M�zler, May 2002, revised July 2002.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

DOWNLOAD ^

mie.m

SOURCE CODE ^

0001 function result = mie(m, x)
0002 % Computation of Mie Efficiencies for given
0003 % complex refractive-index ratio m=m'+im"
0004 % and size parameter x=k0*a, where k0= wave number in ambient
0005 % medium, a=sphere radius, using complex Mie Coefficients
0006 % an and bn for n=1 to nmax,
0007 % s. Bohren and Huffman (1983) BEWI:TDD122, p. 103,119-122,477.
0008 % Result: m', m", x, efficiencies for extinction (qext),
0009 % scattering (qsca), absorption (qabs), backscattering (qb),
0010 % asymmetry parameter (asy=<costeta>) and (qratio=qb/qsca).
0011 % Uses the function "Mie_ab" for an and bn, for n=1 to nmax.
0012 % C. M�zler, May 2002, revised July 2002.
0013 
0014 if x==0                 % To avoid a singularity at x=0
0015     result=[0 0 0 0 0 1.5];
0016 elseif x>0              % This is the normal situation
0017     nmax=round(2+x+4*x.^(1/3));
0018     n1=nmax-1;
0019     n=(1:nmax);cn=2*n+1; c1n=n.*(n+2)./(n+1); c2n=cn./n./(n+1);
0020     x2=x.*x;
0021     f=mie_ab(m,x);
0022     anp=(real(f(1,:))); anpp=(imag(f(1,:)));
0023     bnp=(real(f(2,:))); bnpp=(imag(f(2,:)));
0024     g1(1:4,nmax)=[0; 0; 0; 0]; % displaced numbers used for
0025     g1(1,1:n1)=anp(2:nmax);    % asymmetry parameter, p. 120
0026     g1(2,1:n1)=anpp(2:nmax);
0027     g1(3,1:n1)=bnp(2:nmax);
0028     g1(4,1:n1)=bnpp(2:nmax);   
0029     dn=cn.*(anp+bnp);
0030     q=sum(dn);
0031     qext=2*q/x2;
0032     en=cn.*(anp.*anp+anpp.*anpp+bnp.*bnp+bnpp.*bnpp);
0033     q=sum(en);
0034     qsca=2*q/x2;
0035     qabs=qext-qsca;
0036     fn=(f(1,:)-f(2,:)).*cn;
0037     gn=(-1).^n;
0038     f(3,:)=fn.*gn;
0039     q=sum(f(3,:));
0040     qb=q*q'/x2;
0041     asy1=c1n.*(anp.*g1(1,:)+anpp.*g1(2,:)+bnp.*g1(3,:)+bnpp.*g1(4,:));
0042     asy2=c2n.*(anp.*bnp+anpp.*bnpp);
0043     asy=4/x2*sum(asy1+asy2)/qsca;
0044     qratio=qb/qsca;
0045     result=[qext qsca qabs qb asy qratio];
0046 end;

Generated on Mon 15-Sep-2014 13:31:28 by m2html © 2005