Here is the code for plotting a 5x5 3D planar array 
% Defining the size of the array
M = 5;
N = 5;
theta = linspace(0, pi, 180); % Taking sample values for theta
phi = linspace(0, 2*pi, 360); % Taking sample values for phi
% Creating a meshgrid for theta and phi
[Theta, Phi] = meshgrid(theta, phi);
% Defining the wave number and element spacing
k = (2 * pi) / 30; % Assuming wavelength lambda = 30
d_x = 7.5; % Given lambda/4
d_y = 7.5; % Given d_x = d_y
% Calculating the phase shifts
psi_x = k * d_x * sin(Theta) .* cos(Phi);
psi_y = k * d_y * sin(Theta) .* sin(Phi);
% Calculating the array factor
AF_x = (1/M) * (sin(M/2 * psi_x) ./ sin(psi_x/2));
AF_y = (1/N) * (sin(N/2 * psi_y) ./ sin(psi_y/2));
AF = abs(AF_x .* AF_y);
figure;
surf(psi_x, psi_y, AF, 'EdgeColor', 'none');
xlabel('$\psi_x$','Interpreter', 'latex');
ylabel('$\psi_y$', 'Interpreter', 'latex');
zlabel('Array Factor');
title('3D Planar Array Factor');
colorbar;
