Matlab - how to do a near field / far field transformation

I would like to do a near field - far field transformation with MATLAB. I am using the software tool Feko to model my patch antenna and to get the near field data. With this data, I can display a plane that shows the near E-field of my microstrip antenna (at z = 5mm).

At the moment, I just put Re(Ex)+i.Im(Ex), Re(Ey)+i.Im(Ey) and Re(Ez)+i.Im(Ez) in the same matrix.

I know I should use fft or something like "plane wave spectrum" but I don't know how.

In the end, I would like to display a plane that shows the far E-field (z = 500mm for example).

Any idea to help me would be welcome !

Thank you

Here is the code with what I have done, and comments to explain it.

clear all

%% GET NEAR-FIELD DATAS FROM FEKO AND DELETE COMMENTS

[ficname,pathname] = uigetfile('*.efe','fichier ''.efe'' a convertir ?');
nomfic = [pathname ficname];
i0 = find(ficname=='.');
system(['sed -e "/^#/d;/^*/d" ',' "',nomfic,'"| tr -s " " " " > result.txt']);  
A = load('result.txt');

%% DISPLAY THE E NEAR-FIELD  

X = A(:,1);
Y = A(:,2);
Z = A(:,3);
aux = A(:,4:9);
Ex = sqrt((aux(:,1).^2 + aux(:,2).^2));
Ey = sqrt((aux(:,3).^2 + aux(:,4).^2));
%%Ez = sqrt((aux(:,5).^2 + aux(:,6).^2));
Etot = sqrt(Ex.^2+Ey.^2);

%% Grid interpolation (meshgrid)
patchd = 28 .* 0.001;
patchw = 43 .* 0.001;
pas = 0.5 .* 0.001;

xmin = -patchw/2;
xmax = patchw/2;
pasx = pas;

ymin = -patchd/2;
ymax = patchd/2;
pasy = pas;

[xq,yq] = meshgrid(xmin:pasx:xmax,ymin:pasy:ymax);
EInterp = griddata(X,Y,Etot,xq,yq);

figure(1);
surf(xq,yq,EInterp);
caxis([1 90]);
view(0,90);

%% TRANSFORMATION NEAR-FIELD FAR-FIELD


%% This is where I need to do an FFT on my near-field datas

Ex_im = aux(:,1)+i*aux(:,2);
Ey_im = aux(:,3)+i*aux(:,4);
Ez_im = aux(:,5)+i*aux(:,6);

N = 1024;
E_im = [Ex_im Ey_im];


%%C = fftshift(abs(fft2(E_im,N)));
%%C = fftshift(fft(fftshift(E_im,N)));


%% DISPLAY THE E FAR-FIELD

%% I will see later

%%surf(20*log10(C));

Solution

Only two field components are needed when implementing plane wave spectrum method.

I suppose you are going to use Ex and Ey from frequency domain at a height of z=5mm to calculate Ex, Ey and Ez at z'=500mm.

Calculate kx,ky and kz:

kx range from -pi/X to pi/X and ky from -pi/Y to pi/Y. X and Y are lengths of the plane. kz can be calculated use following equation.

Equation of kz
Calculate FEx and FEy at z=5mm:

This can be easily done using Matlab.
```
 FEx=ifftshift(ifft2(Ex));

 FEy=ifftshift(ifft2(Ey));
```
Calculate FEz at z=5mm:

Equation of FEz

Calculate FEx, FEy, FEz at z'=500mm:

 FEx_500=FEx.*exp(-1i*kz*(z'-z));

 FEy_500=FEy.*exp(-1i*kz*(z'-z));

 FEz_500=FEz.*exp(-1i*kz*(z'-z));

Calculate Ex, Ey, Ez at z'=500mm:

 Ex_500=fft2(fftshift(FEx_500));

 Ey_500=fft2(fftshift(FEy_500));

 Ez_500=fft2(fftshift(FEz_500));

There are some references you may feel interested to take a look:

Investigation on Planar Near-to-Far-Field Transformations for EMC Applications
Plane wave spectrum theory applied to nearfield measurements for electromagnetic compatibility investigations