www.pudn.com > array.zip > circ_array.m, change:2003-07-07,size:5959b

function [pattern,amn] = circ_array(N,dolxr,dolyr,theta0,phi0,winid,win,nbits);
%%%%%%%%%%%%%%%%%%%% ************************ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This function computes the 3-D directive gain patterns for a circular plannar array
% This function uses the fft2 to compute its output. It assumes that there are the same number
% of elements along the major x- and y-axes
%%%%%%%%%%%%%%%%%%%% ************  INPUTS ************ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% N ==> number of elements along x-aixs or y-axis
% dolxr ==> element spacing in x-direction; dolyr ==> element spacing in y-direction. Both are in lambda units
% theta0 ==> elevation steering angle in degrees, phi0 ==> azimuth steering angle in degrees
% This function uses the function (rec_to_circ) which computes the circular array from a square
% array (of size NXN) using the notation developed by ALLEN,J.L.,"The Theory of Array Antennas
% (with Emphasis on Radar Application)" MIT-LL Technical Report No. 323,July, 25 1965.
% winid ==> window identifier; winid negative ==> no window ; winid positive ==> use window given by win
% win ==> input window function (2-D window) MUST be of size (Nxr X Nyr)
% nbits is the number of nbits used in phase quantization; nbits negative ==> NO quantization
%%%%%%%%%%%%%%%%%%%% *********** OUTPUTS ************* %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% amn ==> array of ones and zeros; ones indicates true element location on the grid
% zeros mean no elements at that location; pattern ==> directive gain pattern
%%%%%%%%%%%%%%%%%%%% ************************ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

eps = 0.0001;
nx = 0:N-1;
ny = 0:N-1;
i = sqrt(-1);

% check that window size is the same as the array size
[nw,mw] = size(win);
if winid >0
if mw ~= N
fprintf('STOP == Window size must be the same as the array')
return
end
if nw ~= N
fprintf('STOP == Window size must be the same as the array')
return
end
end

%if dol is > 0.5 then; choose dol = 0.5 and compute new N
if(dolxr <=0.5)
ratiox = 1 ;
dolx = dolxr ;
Nx = N ;
else
ratiox = ceil(dolxr/.5) ;
Nx = (N-1) * ratiox + 1 ;
dolx = 0.5 ;
end
if(dolyr <=0.5)
ratioy = 1 ;
doly = dolyr ;
Ny = N ;
else
ratioy = ceil(dolyr/.5);
Ny = (N-1)*ratioy + 1 ;
doly = 0.5 ;
end

% choose proper size fft, for minimum value choose 256X256
Nrx = 10 * Nx;
Nry = 10 * Ny;
nfftx = 2^(ceil(log(Nrx)/log(2)));
nffty = 2^(ceil(log(Nry)/log(2)));
if nfftx < 256
nfftx = 256;
end
if nffty < 256
nffty = 256;
end

% generate array of elements with or without window
if winid < 0
array = ones(N,N);
else
array = win;
end

% convert steering angles (theta0, phi0) to radians
theta0 = theta0 * pi / 180;
phi0 = phi0 * pi / 180;
% convert steering angles (theta0, phi0) to U-V sine-space
u0 = sin(theta0) * cos(phi0);
v0 = sin(theta0) * sin(phi0);

% Use formula thetal = (2*pi*n*dol) * sin(theta0) divided into 2^m levels
% and rounded to the nearest qunatization level
if nbits < 0
phasem = exp(i*2*pi*dolx*u0 .* nx * ratiox);
phasen = exp(i*2*pi*doly*v0 .* ny * ratioy);
else
levels = 2^nbits;
qlevels = 2.0*pi / levels; % compute quantization levels
sinthetaq = round(dolx .* nx * u0 * levels * ratiox) .* qlevels; % vector of possible angles
sinphiq = round(doly .* ny * v0 * levels *ratioy) .* qlevels; % vector of possible angles
phasem = exp(i*sinthetaq);
phasen = exp(i*sinphiq);
end

% add the phase shift terms
array = array .* (transpose(phasem) * phasen) ;

% determine if interpolation is needed (i.e N > Nr)
if (Nx > N )| (Ny > N)
for xloop = 1 : N
temprow = array(xloop, :) ;
w( (xloop-1)*ratiox+1, 1:ratioy:Ny) =  temprow ;
end
array = w;
else
w(1:Nx, :) = array(1:N,:);
end

% Convert rectangular array into circular using function rec_to_circ
[m,n] = size(w) ;
NC = max(m,n);  % Use Allens algorithm
if Nx == Ny
temp_array = w;
else
midpoint = (NC-1)/2 +1 ;
midwm = (m-1)/2 ;
midwn = (n-1)/2 ;
temp_array = zeros(NC,NC);
temp_array(midpoint-midwm:midpoint+midwm, midpoint-midwn:midpoint+midwn) = w ;
end
amn = rec_to_circ(NC);  % must be rectangular array (Nx=Ny)
amn = temp_array .* amn ;

% Compute array pattern
arrayfft = abs(fftshift(fft2(amn,nfftx,nffty))).^2 ;

%compute [su,sv] matrix
U = [-nfftx/2:(nfftx/2)-1] ./(dolx*nfftx);
indexx = find(abs(U) <= 1);
U = U(indexx);
V = [-nffty/2:(nffty/2)-1] ./(doly*nffty);
indexy = find(abs(V) <= 1);
V = V(indexy);
[SU,SV] = meshgrid(V,U);
indx = find((SU.^2 + SV.^2) >1);
arrayfft(indx) = eps/10;

%Normalize to generate gain patern
rbar=sum(sum(arrayfft(indexx,indexy))) / dolx/doly/4./nfftx/nffty;
arrayfft = arrayfft(indexx,indexy) ./rbar;
[SU,SV] = meshgrid(V,U);
indx = find((SU.^2 + SV.^2) >1);
arrayfft(indx) = eps/10;
pattern = 10*log10(arrayfft +eps);

figure(1)
mesh(V,U,pattern);
xlabel('V')
ylabel('U');
zlabel('Gain pattern - dB')

figure(2)
contour(V,U,pattern)
axis image
grid
xlabel('V')
ylabel('U');
axis([-1 1 -1 1])

figure(3)

x0 = (NC+1)/2 ;
y0 = (NC+1)/2 ;
radiusx = dolx*((NC-1)/2 + 0.05/dolx) ;
radiusy = doly*((NC-1)/2 + 0.05/dolx) ;
theta = 5  ;
[xxx, yyy]=find(abs(amn)>0);
xxx = xxx-x0 ;
yyy = yyy-y0 ;
plot(yyy*doly, xxx*dolx,'rx')
axis equal
hold on
grid
title('antenna spacing pattern');
xlabel('y - \lambda units')
ylabel('x - \lambda units')

plot(y, x) ;
return