Asked by Mona Mahboob Kanafi
on 20 May 2015

I'm trying to generate a randomly rough surface(isotropic surface). Suppose I have the absolute values of fourier components of the surface in "H" and H has symmetry with respect to center of the array. Now, I need to define some random phase for the fourier series which are between -pi to pi (due to fftshift) as in:

phi = -pi + (pi+pi)*rand(n); % random phase

Now, I don't know what to do to make the output of my surface real..

Here is the gist of the code I used so far:

[a,b] = pol2cart(phi,H);

H_complex = complex(a,b); % the complex fourier transform composed of H & phi

z = ifft2(ifftshift((H_complex)));

If I need to apply conjugate symmetry to the fourier components, how should I do it? The output is complex and I need a real surface.

Please help! Thanks a lot.

Answer by Mona Mahboob Kanafi
on 21 May 2015

Edited by Mona Mahboob Kanafi
on 21 May 2015

Accepted Answer

I found my answer, so maybe useful for someone else as well. I need to apply conjugate symmetry condition for 2D fourier transform, both in magnitude values and phase components. I applied it as below for a n*m matrix:

"H" is the magnitude:

H(1,1) = 0;

H(1,m/2+1) = 0;

H(n/2+1,m/2+1) = 0;

H(n/2+1,1) = 0;

H(2:end,2:m/2) = rot90(H(2:end,m/2+2:end),2);

H(1,2:m/2) = rot90(H(1,m/2+2:end),2);

H(n/2+2:end,1) = rot90(H(2:n/2,1),2);

H(n/2+2:end,m/2+1) = rot90(H(2:n/2,m/2+1),2);

Similar operation but with negative sign applies to "phi" (phase):

phi(1,1) = 0;

phi(1,m/2+1) = 0;

phi(n/2+1,m/2+1) = 0;

phi(n/2+1,1) = 0;

phi(2:end,2:m/2) = -rot90(phi(2:end,m/2+2:end),2);

phi(1,2:m/2) = -rot90(phi(1,m/2+2:end),2);

phi(n/2+2:end,1) = -rot90(phi(2:n/2,1),2);

phi(n/2+2:end,m/2+1) = -rot90(phi(2:n/2,m/2+1),2);

You can refer to this link for a good example of 2D fourier transform on a matrix and how the results should look like: http://fourier.eng.hmc.edu/e101/lectures/image_processing/node6.html

Ralf Mackenbach
on 10 Apr 2019

Small correction to this. The condition that

H(1,1) = 0;

H(1,m/2+1) = 0;

H(n/2+1,m/2+1) = 0;

H(n/2+1,1) = 0;

Does not need to hold necessarily. So a better (examplary) version would be

% Choose even m and n

m = 100 ;

n = 100 ;

% Choose average value of surface

z0 = 1 ;

% Fill random examplary matrices

H = rand(n,m) ;

phi = 2*pi*rand(n,m) ;

% Apply symmetry to matrices

H(1,1) = z0*m*n ;

H(2:end,2:m/2) = rot90(H(2:end,m/2+2:end),2);

H(1,2:m/2) = rot90(H(1,m/2+2:end),2);

H(n/2+2:end,1) = rot90(H(2:n/2,1),2);

H(n/2+2:end,m/2+1) = rot90(H(2:n/2,m/2+1),2);

phi(1,1) = 0;

phi(1,m/2+1) = 0;

phi(n/2+1,m/2+1) = 0;

phi(n/2+1,1) = 0;

phi(2:end,2:m/2) = -rot90(phi(2:end,m/2+2:end),2);

phi(1,2:m/2) = -rot90(phi(1,m/2+2:end),2);

phi(n/2+2:end,1) = -rot90(phi(2:n/2,1),2);

phi(n/2+2:end,m/2+1)= -rot90(phi(2:n/2,m/2+1),2);

% Convert to polar

[a,b] = pol2cart(phi,H);

H_complex = complex(a,b);

% Inverse and plot

S = ifft2(H_complex, 'symmetric') ;

surf(S)

Sign in to comment.

Answer by Thomas Koelen
on 20 May 2015

use abs() or real()? depending on if you want only the real number or the absolute value.

Mona Mahboob Kanafi
on 20 May 2015

Sign in to comment.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.