% thinwire_rcs.m% =========================================================================% MATLAB code for the calculation of%	Frequency Response of the Polarization Scattering Matrix % of thin wires at broadside, by converting 2-D (infinite bodies) % solutions of the wave equation in cylindrical coordinates to the 3-D % (finite bodies) cases. Simplified formulas are used. % Syntax:%		[thh,tvv,thv]=thinwire_rcs(a,l,tau,fmin,fmax,nf)% Input:%		a -- radius of the thin wire (cm)%		l -- length of the thin wire (cm)% 		tau-- angle tilted above the local horizontal (deg.)%		fmin -- minimum frequency (MHz)%		fmax -- maximum frequency (MHz)%		nf -- No. of frequency samples% Output:%		thh -- complex RCS for H-H polarization %		tvv -- complex RCS for V-V polarization%		thv -- complex RCS for H-V polarization (thv=tvh)% Conditions when holding accuracy: %		a<=wavelength/10.%		l>=5*wavelength% ==========================================================================function [thh,tvv,thv]=thinwire_rcs(a,l,tau,fmin,fmax,nf)l=l/100;    % (m)tau=tau*pi/180.;  % (rad.)j=sqrt(-1);gama=0.57721;c=30000.0;  % (cm/s)if nf~=1	fstep=(fmax-fmin)/(nf-1);else    fstep=0.0;endfor i=1:nf   freq(i)=fmin+(i-1)*fstep;   wavelength=c/freq(i);  % (cm)   k=2*pi/wavelength;      te=j*l*sqrt(pi)/(pi/2.+j*(gama+log(k*a/2.)));   th=-l*1.5*sqrt(pi)*(k*a)^2;	thh(i)=.5*(te-th)+.5*(te+th)*cos(2*tau);	tvv(i)=.5*(te-th)-.5*(te+th)*cos(2*tau);	thv(i)=.5*(te+th)*sin(2*tau);endgresult=thinwire_graph(thh,tvv,thv,freq);% Graphics% --------------------------------------------------function gresult=thinwire_graph(thh,tvv,thv,freq)   % --------------------------------------------------nf=max(size(freq));if nf~=1   figure      thha=20*log10(abs(thh)+eps);   tvva=20*log10(abs(tvv)+eps);   thva=20*log10(abs(thv)+eps);   amax=max([max(thha) max(tvva) max(thva)]);   amin=min([min(thha) min(tvva) min(thva)]);   if amin<-65      amin=-65;   end      thhp=atan2(imag(thh),real(thh))*180/pi;   tvvp=atan2(imag(tvv),real(tvv))*180/pi;   thvp=atan2(imag(thv),real(thv))*180/pi;      subplot(321)   plot(freq,thha)   axis([freq(1) freq(nf) amin-5 amax+5])   grid on   title('Magnitude of RCS (HH)')   %xlabel('Frequency (MHz)')   ylabel('RCS (dBsm)')      subplot(322)   plot(freq,thhp)   axis([freq(1) freq(nf) -180 180])   grid on   title('Phase of RCS (HH)')   %xlabel('Frequency (MHz)')   ylabel('Phase (deg.)')      subplot(323)   plot(freq,tvva)   axis([freq(1) freq(nf) amin-5 amax+5])   grid on   title('Magnitude of RCS (VV)')   %xlabel('Frequency (MHz)')   ylabel('RCS (dBsm)')      subplot(324)   plot(freq,tvvp)   axis([freq(1) freq(nf) -180 180])   grid on   title('Phase of RCS (VV)')   %xlabel('Frequency (MHz)')   ylabel('Phase (deg.)')      subplot(325)   plot(freq,thva)   axis([freq(1) freq(nf) amin-5 amax+5])   grid on   title('Magnitude of RCS (HV)')   xlabel('Frequency (MHz)')   ylabel('RCS (dBsm)')      subplot(326)   plot(freq,thvp)   axis([freq(1) freq(nf) -180 180])   grid on   title('Phase of RCS (HV)')   xlabel('Frequency (MHz)')   ylabel('Phase (deg.)')	gresult=1;   else   thha=20*log10(abs(thh)+eps)   thhp=atan2(imag(thh),real(thh))*180/pi   tvva=20*log10(abs(tvv)+eps)   tvvp=atan2(imag(tvv),real(thh))*180/pi   thva=20*log10(abs(thv)+eps)   thvp=atan2(imag(thv),real(thv))*180/pi   gresult=0;end