2

This a much discussed topic, but I happen to have an issue that has not been answered yet. My problem is not the method it self but rather it's applicability: My image's f(x,y) represent physical values that can be negative or positive. When I mask the peaks corresponding with, say the median, i get, after application of the inverse FFT, an image which is complex.

This seams logical as image != ifft(fft(image)) if image != image, thus it can very well be complex result?

I thus take the absolute value of my image array, and get a nicely cleaned image. But by taking the abs of the image i have lost the negative values!

My code is complex and uses multiple images to find the correct positions where to mask so I will break down to the essentials:

def everything(fft,fftImage,sizeOfField,shapeOfFFT):
max_x = []
max_y = []
median = np.median(fft)

threshold = 500
#correctLocalMax() holds several subfunctions that look for the propper max_x and max_y. This works fine and returns 2 lists max_x,max_Y that contain the coordiantes of the max's
max_x,max_y = correctLocalMax(iStart = 0,iStop = 30, jStart =0 , jStop = shapeOfFFT[1],threshold=threshold, max_x = max_x, max_y = max_y)

for i in range(len(max_x)):
    for k in range(sizeOfField):
        for l in range(sizeOfField):
            fftImage[max_x[i]+k][max_y[i]+l] = median

return(fftImage)

image, coverage, stdev = pickleOpener(dataDir,i)
field = getROI(image,area,i0,j0)

fftImage = np.fft.fft2(image)
fftImage = np.fft.fftshift(fftImage)

fft = np.fft.fft2(coverage)
fft = np.fft.fftshift(fft)

fftMod = everything(fft, fftImage, sizeOfField, shapeOfFFT)
imageBack = np.fft.ifft2(fftMod)
imageBack = np.abs(imageBack)
field = getROI(imageBack,area,i0,j0)

The images I have and get after processing look like this: enter image description here The stripe pattern are the ones I wish to remove

enter image description here These are the masks applied to the FFT

enter image description here The stripe pattern is mostly removed, however now the image is purely positive!

You can find the proper solution to the problem in the comments!

Community
  • 1
  • 1
Sebastiano1991
  • 867
  • 1
  • 10
  • 26
  • This is a really great question, and there’s a good answer for it. This blog post explains what’s happening and how to fix it in very generic terms: http://blogs.mathworks.com/steve/2010/07/16/complex-surprises-from-fft/ but it’s 1D and doesn’t address the problem you’re trying to solve. Can you post the original image data somewhere (with the stripes)? I can then show you how to adjust your mask to get real-only outputs out of the IFFT. – Ahmed Fasih Oct 18 '16 at 13:49
  • Never mind about asking the data, I reread the post about how `correctLocalMax` is a complicated function. The basic problem is that when you mask the spectrum to zeros, you need to ensure you preserve the spectrum’s conjugate-symmetry. – Ahmed Fasih Oct 18 '16 at 13:55
  • 1
    What’s confusing me is: your spectrum (the output of `fft2`, your `fftMod`) should already be conjugate-symmetric: peaks in it (which here correspond to your stripes) should be symmetric left-to-right and anti-symmetric up-to-down. But your `correctLocalMax` seems to return masks that aren't conjugate-symmetric? If you fix that, then your code should work. – Ahmed Fasih Oct 18 '16 at 15:17
  • Amazing answer as well! this indeed made the deal, I just had to squeeze `correctLocalMax` to return symmetrical patches and it worked just fine! – Sebastiano1991 Oct 19 '16 at 10:22

1 Answers1

1

You could try two different approaches: Either you scale your image between the original values first, and rescale later, somewhat like this:

max_val = max(max(A))
min_val = min(min(A))
% normalize to [0,1]
image = norm(image)
% do your stuff here
% then rescale to original values
image = min_val + (max_val - min_val).*image / (max_val - min_val)

An alternative would be to save which values where negative in first place. Although I'd recommend to check whether they got changed during your function call to avoid reinstating your noise

dennlinger
  • 9,890
  • 1
  • 42
  • 63
  • Hello, concerning the first idea: I was thinking about the same, however does this overcome the general issue - As far as I understand if I want to generate an image from a complex array generated from a FFT the absolute value is the way to go? The second one seams interesting as well, I will look into that a bit. Thanks! I saw in other posts (i.e. http://stackoverflow.com/questions/34027840/removing-periodic-noise-from-an-image-using-the-fourier-transform) that people used np.real to generate the output image. This however discards all the complex values which method seams more reasonable? – Sebastiano1991 Oct 18 '16 at 11:56
  • I am beginning to think that one should prefer np.real. My heuristic mathematical argument goes like that: If an the fft(image) is back transformed to the original image there should not be any imaginary parts. The changing of the fft introduces imaginary parts, while at the same time removing the pattern: The back transform closest to the original image would have no imaginary part --> thus drop it? Is this too much experimental physicists math to be true? – Sebastiano1991 Oct 18 '16 at 12:22
  • np.real could resolve your issue. Although I'm no expert with the numpy, that's what i believe is told in the MATLAB documentary regarding this topic (the syntax is pretty similar to begin with, so maybe you can look for answers there, too). Tbh, i can't really help you much more than I already did. Maybe one final thing: I noticed you don't use `ifftshift`, so you might look into that as well – dennlinger Oct 19 '16 at 12:14