%script template for ECE123 Assignment 3 part 2
%created 02/01/2008 by WAC for array of perpendicular dipoles
%array axis is z and dipole axis is x
%
%Need N=number of elements and d=element spacing in wavelengths
%and off-axis steering angle in deg = steer
thetas=90-steer;
cts=cosd(thetas);

figure(1)
nsamp=256;              %can increase if a narrower beam is being analyzed
[x,y,z]=sphere(nsamp);  %samples sphere in mesh equally spaced in theta and phi
dtheta=pi/nsamp;dphi=2*pi/nsamp;    %get increments in theta and phi
%gain=(1-x.^2).*(sin(pi*N*d*(z-cts+eps))./sin(pi*d*(z-cts+eps))).^2;
gain=(sin(pi*N*d*(z-cts+eps))./sin(pi*d*(z-cts+eps))).^2;
integrand=gain.*sqrt(1-z.^2);   %multiply by sin(thetaz) and make an array
int=sum(integrand(:))*dtheta*dphi;  %sum(a(:)) sums over entire multidimensional array
gain=(4*pi/int)*gain;   %normalize gain
10*log10(max(gain(:)))  %convert max gain to dBi
xx=x.*gain; yy=y.*gain; zz=z.*gain; %radial scaling by gain (could be done in loggain)
mesh(xx,yy,zz)
axis equal
%note theta=0 is vertical on the mesh plot