Search code examples
matlabfftnumerical-methodsscilab

Algorithm behaves differently in matlab and scilab

I'm implementing the Gerchberg-Saxton-Algorithm for phase retrieval in Matlab as part of a report. When I didn't have access to a PC with a Matlab license, I thought I could just continue working with Scilab. However, the results differ noticably.

Matlab code:

function [ realphase ] = GerchSax( source, target, iterations )

[n, m] = size(target);

fourierspace = fft2(source);
for iter=1:iterations
    fourierphase = angle(fourierspace);
    fourierspace = sqrt(abs(target)) .* exp(1i * fourierphase);

    realspace = ifft2(fourierspace);
    realphase = angle(realspace);
    realspace = sqrt(abs(source)) .* exp(1i * realphase);

    fourierspace = fft2(realspace);
end

Scilab code: 

function [ res ] = angle(A)
    res = atan(imag(A), real(A))
endfunction


function [ realphase ] = GerchSax( source, target, iterations )

fourierspace = fft(source);
for iter=1:iterations
    fourierphase = angle(fourierspace);
    fourierspace = sqrt(abs(target)) .* exp(%i * fourierphase);

    realspace = ifft(fourierspace);
    realphase = angle(realspace);
    realspace = sqrt(abs(source)) .* exp(%i * realphase);

    fourierspace = fft(realspace);
end

The implementation of the "angle" function comes directly from Scilab documentation. As I understand it, fft in matlab performs the Fourier transformation on each column vector independently if given a matrix (thus the use of fft2 for actual 2D Fourier transformation), while scilab already does a 'real' 2D transformation with fft.

I used source = [1 1 1; 1 2 1; 1 1 1], target = [0 1 1; 0 1 1; 0 0 0] and 1000 iterations in both programs. Matlab returns the following phase distribution:

   -3.0589   -1.0472    0.9645
    3.1416         0    3.1416
    3.0589    1.0472   -0.9645

While I get this phase distribution from Scilab: 

    2.0943951  - 1.0471976    2.0943951  
    3.1415927  - 1.977D-16    3.1415927  
  - 2.0943951    1.0471976  - 2.0943951

Did I make any mistakes converting the Matlab code into its correct Scilab equivalent? Is the FFT implemented differently in these environments?

Solution

  • fft is different between Scilab and Matlab the link below will show you how to implement fft to have an equivalent result in either matlab or scilab.

    fft (Matlab function)