problem with GSAFT (backprojection) algorithm for nearfield image reconstruction on my simple planar array

Question

I have a planar array containing 16 transmitters and 88 receivers located on the same plane with dimension 30cm-by-30cm. I have an object of size 6cm-by-6cm. Object and antenna array plane are in parallel, and their distance is 15cm. GSAFT formula is the same as below:

This formula is extracted from the following article.

Its implementation in MATLAB is the same as below. Electrical fields are obtained by FEKO software:

clc;
clear all;
close all;
numT = 16; % number of transmitter
numR = 88; % number of receiver
numF = 11; % number of frequency

load Efield  % a (numT*numR , 6 , numF) array  
% Efield(:,1:3,:) = x,y,z dimension of receiver
% Efield(:,4,:) = complex electrical field Ex in each receiver position and frequency
% Efield(:,5,:) = complex electrical field Ey in each receiver position and frequency 
% Efield(:,6,:) = complex electrical field Ez in each receiver position and frequency


% defining transmitter coordinates

Tx(1)=29.250000e+00;
Ty(1)=5.250000e+00;
...
...
...
Tx(16)=-5.250000e+00;
Ty(16)=-29.250000e+00;

x = (-8:0.2:8); % x pixel in centimeter
y = (-8:0.2:8); % y pixel in centimeter
coor = zeros(numel(x)*numel(y),2); % coordinates of pixels of object in a row



for x1 = x,
    for y1 = y,
       coor(count,1) = x1;
       coor(count,2) = y1;
       count = count + 1;
    end
end

freq = 10e9:(10e9)/(numF-1):20e9; 
reflectivity = zeros(numel(x)*numel(y),3); % object pixel reflectivity

c = 3e8; % light speed
z = 15; % distance between antenna array plane and object 

for i = 1:numel(x)*numel(y),    
    for tx = 1:numT  ,
        rt = sqrt((Tx(tx)-coor(i,1))^2 + (Ty(tx)-coor(i,2))^2 + z^2); % difference between transmitter and object pixel
        for rx = 1:numR
            for f = 1:numF,
                rr = sqrt((Efield((tx-1)*numR + rx,1,f)-coor(i,1))^2 + (Efield((tx-1)*numR + rx,2,f)-coor(i,2))^2 + z^2); % distance between receiver and object pixel
                r = (rr + rt)/100; % converting total distance to meter

                
                reflectivity(i,1) = or(i,1) + Efield((tx-1)*numR + rx,4,f)*exp(1j*(2*pi*freq(f)/c)*r); % reflectivity resulted from Ex
                reflectivity(i,2) = or(i,2) + Efield((tx-1)*numR + rx,5,f)*exp(1j*(2*pi*freq(f)/c)*r); % reflectivity resulted from Ey
                reflectivity(i,3) = or(i,3) + Efield((tx-1)*numR + rx,6,f)*exp(1j*(2*pi*freq(f)/c)*r); % reflectivity resulted from Ez

            end
        end
    end    
end



image = zeros(numel(x)*numel(y),1);
image(:) = nthroot(abs(reflectivity(:,1)).^2 + abs(reflectivity(:,2)).^2 + abs(reflectivity(:,3)).^2 , 2);
image = image/max(image); % 
image = reshape(image,[numel(x),numel(y)]);

figure(1);
h = surf(x,y,20*log10(image));
colorbar;
set(h,'LineStyle','none');
view(2);
caxis([-15 0]);

When I reconstruct my image, it is like this. unfortunately, I do not why edges of my image are not reconstructed correctly:

Can you please help me on what is my problem with this method of millimeter wave image reconstruction?