%This is a file to look at the antenna response of a patch and a linear
%dipole antenna. The dipole is in the plane of the patch, which means we
%have to rotate the typical toroidal response by 90 degrees, but keep the
%same theta and phi as used to describe the patch. On top of this the
%outputs of the dipole antenna can have an arbitrary phase and amplitude
%shift before it is added to the patch. These antennas are close together,
%but there is still a small array factor to include.
clear;
N = 40;%even number
%Here we rotate the toroidal response of a dipole but we have to do this
%and keep the same phi and theta used for the other antenna responses. We
%do this by taking the phi and theta we want to use in the rotated version,
%and then rotate it 90 degrees to be the figure we can easily generate in
%spherical coordinates, with symmetry in theta.
% A more general rotation function could be generated by deriving the
% conversion functions from spherical to cart then back, with an arbitrary
% rotation angle, but I only needed 90 degrees so I used a shortcut.
%generate the theta and phi for the rotated response, this corresponds to
%the theta and phi of the other antenna that is different and not rotated
%swapping x and z results in a rotation around the y axis
[THETA,PHI] = meshgrid(0:((2*pi)/N):2*pi,0:((2*pi)/N):2*pi);
[m,n] = size(PHI);
R1=ones(m,n);%generate an R matrix for the rotation, this does not matter
%because R cancels out in the conversion back to spherical, this is a handy
%response to have for troubleshooting
[x,y,z]=sph2cart(THETA,PHI,R1);%Convert to cartestian
[THETA1,PHI1,R2]=cart2sph(z,y,x);%flip x and z and convert back - R2 is a
%throw-away set of variables you can overwrite it, or not.
% use the new values of theta and phi to generate the R matrix, in this
% case we only need the new phi. To use the result we use the new R values
% with the original theta and phi
%R3 = (cos((pi/2)*cos(pi/2 - PHI1.*1)));
%*********************************************************************** R3 AMPLITUDE
R3 =.5*abs(cos(PHI1));%simplified toroidal response with a gain value
% Now we generate the response of the other antenna element, this is a
% circular patch type response
%***********************************************************************PATCH AMPLITUDE RATIO TERMS
R4(1:N/2,1:N+1) = 2*abs(sin((pi/2)*cos(pi/2 - PHI(1:N/2,1:N+1))));
R4(N/2+1:N+1,1:N+1) = .3*abs(sin((pi/2)*cos(pi/2 - PHI(N/2+1:N+1,1:N+1))));
%The following section shows how to make an array factor for the 2 antennas,
%but what we really want to do in this case is to get an array phase term
%which we can add to the phase shift we use in the combining process
%Let us generate the array curve. If we have a spacing of "spacing" then this is
%expressed in wavelengths lambda as spacing/lambda. Orthogonal to the array will
%have a gain of 1, and in the plane of the array the reduction as a
%function of the angle P will be
% 1/sqrt(2)*sqrt(1+cos[2*pi*spacing/lambda*sin(P)])
%Once again this is easy to describe in spherical coordinates but at 90
%degrees from the way we want. We can use the same trick to generate the
%array factor that we used to rotate the dipole response, but this time we
%will do it around the X axis.
%*************************************************************************** D
D=.1;%spacing in terms of wavelengths
%rotate around the X axis, similar to above
[x,y,z]=sph2cart(THETA,PHI,R1);%Convert to cartestian
[THETA1,PHI1,R2]=cart2sph(x,z,y);
%below is the array amplitude factor, which we don't really need
%ArrayFactor=1/sqrt(2)*sqrt(1+cos(2*pi*D*sin(PHI1)));
%similar to above, let us generate an array phase term
PhaseArray=2*pi*D*sin(PHI1);
PhaseAdjust=pi;%the phase shift for the combiner***************************COMBINER PHASE
AmplitudeAdjust=1;%Amplitude adjustment for the combiner*******************COMBINER AMPLITUDE
PhaseTotal=PhaseAdjust+PhaseArray;
R=sqrt((R4+AmplitudeAdjust*R3.*cos(PhaseTotal)).^2+(AmplitudeAdjust*R3.*sin(PhaseTotal)).^2);
%use the unity sphere to debug the process
%R=sqrt((R1+AmplitudeAdjust*R3.*cos(PhaseTotal)).^2+(AmplitudeAdjust*R3.*sin(PhaseTotal)).^2);
%convert to dB, sort of. We search for values of R plus 20 dB that are less
%than 0, then replace the negative values with 0. The plot maximums will be
%20 dB high
R=20+20*log10(R);
[r,c,v]=find(R<0);
R(r,c)=0;
% This next part does the coordinate conversion.
[x,y,z] = sph2cart(THETA,PHI,R);
surf(x,y,z);
xlabel('x');
ylabel('y');
zlabel('z');
title('F(theta)');
zlim([-15,25]);
ylim([-25,25]);
xlim([-25,25]);