0
$\begingroup$

After reading the paper Phase retrieval algorithms: a comparison by J.R.Fienup. I implemented the Error-Reduction Algorithm in MATLAB. The initial phase input was randomly generated, so the final recovered image in each execution is different. The only information given is the Fourier Modulus of the original image.
Sometimes I can get a pretty well result, but most of the time the result is a fragmented version of the original image as the figures below shown:
The Recovered Image The Original Image:$\pi$
As you can see, the image of $\pi$ is "almost" recovered except that the head and the body of $\pi$ had been cut and swapped.
My questions are:
1. Is this a common problem using phase retrieval algorithm? What causes this?
2. Is there any method to avoid this?
3. Is there any other great methods in phase retrieval?

$\endgroup$
1
  • $\begingroup$ There is a $\pi /2$ phase shift in Matlab fft/fft2 $\endgroup$
    – lennon310
    Mar 18, 2015 at 13:26

1 Answer 1

0
$\begingroup$

I think the effect you're seeing can be explained in terms of the shift theorem of Fourier transforms and the periodicity of the discrete Fourier transform.

It looks like the recovered object has been shifted upwards, and because the DFT is periodic, the top part of the object that's left the field of view has wrapped around to re-appear at the bottom of the field. The reason the object is shifted is hard to see without knowing more about how your implementation of phase retrieval. I can say that, due to the shift theorem, we know a shift in the object domain causes a linear phase in the Fourier domain, but does not affect the Fourier modulus. So, the Fourier modulus alone can't help you find the correct shift.
I suspect you're leaving out some key information when you say "The only information given is the Fourier Modulus", because if that were strictly true you wouldn't be able to do phase retrieval at all. As discussed in the paper, you need to have some other piece of information, typically the modulus of the object or the object support (the object support is a closed boundary that you know contains the object). If you know the modulus of the object or at least its support, then reinforcing that constraint on every iteration of the error-reduction algorithm would prevent the object estimate from shifting.

As to your third question, I would point out that error-reduction is only the simplest algorithm discussed in the paper. Gradient searches or the input-output methods described in the later sections are usually faster and/or more reliable. In more recent years, successful methods that leverage nonlinear optimization algorithms such as BFGS have also become common. If you're interested, you might check out Prof. Fienup's list of publications.

$\endgroup$
4
  • $\begingroup$ Thank you @user966000 for your detailed answer. This experiment is done by following these step: $\endgroup$
    – Mr. Who
    Mar 30, 2015 at 1:32
  • $\begingroup$ 1) Draw the white color pattern $\pi$ in a 100by100 black image; 2) Use imread to load the .bmp image as a matrix p; 3) Calculated its Fourier modulus matrix as absP; 4) Use random phase theta to form a Fourier domain guess G_pr. 5) Inverse transform G_pr to g_pr; 6) As the phase of p is a zero matrix, we let the some elements of g_pr, whose phase values are too big, be zero; 7) Keep the real part of g_pr ; 8)FFT g_pr into G and only keep its phase as theta then back to step 4. $\endgroup$
    – Mr. Who
    Mar 30, 2015 at 1:54
  • $\begingroup$ Here is the Matlab code: <br/> format long p=imread('pi.bmp'); P=fftshift(fft2(p)); absP=abs(P); k=0; theta=2*pi.*rand(100)-pi; %[-pi,pi] figure(1) figure(2) while k<8000 G_pr=absP.*exp(1i.*theta); g_pr=ifft2(ifftshift(G_pr)); absPhase=abs(angle(g_pr)); maxPh=max(max(absPhase)); minPh=min(min(absPhase)); g_pr(absPhase>=(minPh+0.2*(maxPh-minPh)))=0; g_pr=real(g_pr); gg=255*g_pr/(max(max(g_pr))); figure(1),imshow(uint8(gg)); G=fftshift(fft2(g_pr)); G=G./abs(G); theta=angle(G); figure(2),imshow(theta),title(num2str(k)) k=k+1; end <br/> $\endgroup$
    – Mr. Who
    Mar 30, 2015 at 1:54
  • $\begingroup$ @Mr.Who Thanks for code and outline. From those I can see that, besides knowing the Fourier modulus, you also assume that the object is real-valued and non-negative. So, that's enough information that we can reasonably expect the results you're getting. However, you'll need some more data and/or assumptions if you want to resolve the ambiguity in the shift I described earlier. That could be as simple as including a support constraint, though. $\endgroup$
    – user966000
    Apr 1, 2015 at 19:14

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.