<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 FourierCats.ipynb
‹ FourierCats.ipynb › Copyright (C) ‹ 2017 › ‹ Anna Scaife - anna.scaife@manchester.ac.uk ›
This program is free software
Step1:  We need to import some libraries
Step2:  Open up the input image
Step3:  ...and let's take a look at it
Step4:  Next we'll dive straight in and Fourier Transform our cat. Note that there are a couple of steps to this, the first is the fft itself and the second is an fft "shift". This is necessary because of the way that the frequency space is ordered in the Fourier Transform.
Step5:  Let's see how that looks
Step6:  Now we can set things up for filtering our Fourier Cat. We need to know (1) the dimenions of the image and (2) where the centre is.
Step7:  To start with, let's make a filter function that separates the inner most 40 x 40 pixels from everything else.
Step8:  We can then use this filter to, firstly, mask out the inner most 40 x 40 pixels in Fourier space. This removes our small Fourier frequencies, i.e. the large scale information in our image.
Step9:  We can then Fourier transform this back into image space.
Step10:  ...and, see how it looks.
Step11:  Now let's filter out the large Fourier frequencies
Step12:  and Fourier transform that back into image space
Step13:  This looks like
Step14:  We can also take the Fourier Transform of the filter function to see how the PSF looks
Step15:  Let's use that to make the plot from the lecture
Step16:  What about... instead of filtering out a contiguous range of Fourier frequencies, we rnadomly selected the Fourier components?
<ASSISTANT_TASK:>
Python Code:
%matplotlib inline
Explanation: FourierCats.ipynb
‹ FourierCats.ipynb › Copyright (C) ‹ 2017 › ‹ Anna Scaife - anna.scaife@manchester.ac.uk ›
This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with this program. If not, see http://www.gnu.org/licenses/.
[AMS - 170905] Notebook created for TIARA Astrostatistics Summer School, Taipei, September 2017.
This notebook presents some simple properties of Fourier space and accompanies the lecture: "Imaging the Invisible"
It uses a number of Python libraries, which are all installable using pip.
We want plotting inline in the notebook:
End of explanation
import numpy as np   # for array manipulation and the fft
import pylab as pl   # for plotting
import cv2           # for image file handling
Explanation: We need to import some libraries:
End of explanation
cat = cv2.imread('./FIGURES/cat1.jpg',0)
Explanation: Open up the input image:
End of explanation
pl.imshow(cat,cmap='gray')
pl.show()
Explanation: ...and let's take a look at it:
End of explanation
cat_squiggle = np.fft.fft2(cat)
cat_squiggle_shifted = np.fft.fftshift(cat_squiggle)
cat_spectrum = 20*np.log(np.abs(cat_squiggle_shifted))
Explanation: Next we'll dive straight in and Fourier Transform our cat. Note that there are a couple of steps to this, the first is the fft itself and the second is an fft "shift". This is necessary because of the way that the frequency space is ordered in the Fourier Transform.
End of explanation
pl.subplot(121),pl.imshow(cat, cmap = 'gray')
pl.title('Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(122),pl.imshow(cat_spectrum, cmap = 'gray')
pl.title('Fourier Cat'), pl.xticks([]), pl.yticks([])
pl.show()
Explanation: Let's see how that looks:
End of explanation
rows, cols = cat.shape
crow,ccol = rows/2 , cols/2
Explanation: Now we can set things up for filtering our Fourier Cat. We need to know (1) the dimenions of the image and (2) where the centre is.
End of explanation
filter_fnc = np.zeros(cat_squiggle_shifted.shape)
filter_fnc[crow-20:crow+20, ccol-20:ccol+20] = 1.0
pl.imshow(filter_fnc)
pl.show()
Explanation: To start with, let's make a filter function that separates the inner most 40 x 40 pixels from everything else.
End of explanation
cat_squiggle_hpf = np.copy(cat_squiggle_shifted)
cat_squiggle_hpf[np.where(filter_fnc==1.)] = 0.0+0*1j
Explanation: We can then use this filter to, firstly, mask out the inner most 40 x 40 pixels in Fourier space. This removes our small Fourier frequencies, i.e. the large scale information in our image.
End of explanation
cat_filtered = np.fft.ifftshift(cat_squiggle_hpf)
cat_filtered_hpf = np.fft.ifft2(cat_filtered)
cat_filtered_hpf = np.abs(cat_filtered_hpf)
Explanation: We can then Fourier transform this back into image space.
End of explanation
#pl.subplot(121),pl.imshow(cat, cmap = 'gray')
#pl.title('Image'), pl.xticks([]), pl.yticks([])
pl.subplot(121),pl.imshow(20*np.log(np.abs(cat_squiggle_hpf)), cmap = 'gray')
pl.title('Filtered Fourier Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(122),pl.imshow(cat_filtered_hpf)
pl.title('HPF Cat'), pl.xticks([]), pl.yticks([])
pl.show()
Explanation: ...and, see how it looks.
End of explanation
cat_squiggle_lpf = np.copy(cat_squiggle_shifted)
cat_squiggle_lpf[np.where(filter_fnc==0.)] = 0.+0.*1j
Explanation: Now let's filter out the large Fourier frequencies:
End of explanation
cat_filtered = np.fft.ifftshift(cat_squiggle_lpf)
cat_filtered_lpf = np.fft.ifft2(cat_filtered)
cat_filtered_lpf = np.abs(cat_filtered_lpf)
Explanation: and Fourier transform that back into image space:
End of explanation
#pl.subplot(121),pl.imshow(cat, cmap = 'gray')
#pl.title('Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(121),pl.imshow(20*np.log(np.abs(cat_squiggle_lpf)), cmap = 'gray')
pl.title('Filtered Fourier Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(122),pl.imshow(cat_filtered_lpf)
pl.title('LPF Cat'), pl.xticks([]), pl.yticks([])
pl.show()
Explanation: This looks like:
End of explanation
psf = np.fft.ifft2(filter_fnc)
psf = np.fft.ifftshift(psf)
psf = np.abs(psf)
Explanation: We can also take the Fourier Transform of the filter function to see how the PSF looks:
End of explanation
pl.subplot(231),pl.imshow(20*np.log(np.abs(cat_squiggle_lpf)), cmap = 'gray')
pl.title('Filtered Fourier Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(232),pl.imshow(20*np.log(np.abs(cat_squiggle_shifted)), cmap = 'gray')
pl.title('Fourier Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(233),pl.imshow(filter_fnc, cmap = 'gray')
pl.title('Filter'), pl.xticks([]), pl.yticks([])
pl.subplot(234),pl.imshow(cat_filtered_lpf, cmap = 'gray')
pl.xlabel('LPF Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(235),pl.imshow(cat, cmap = 'gray')
pl.xlabel('Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(236),pl.imshow(psf, cmap = 'gray')
pl.xlabel('PSF'), pl.xticks([]), pl.yticks([])
pl.show()
Explanation: Let's use that to make the plot from the lecture:
End of explanation
filter_mask = np.random.randint(2, size=cat_squiggle_shifted.shape)
pl.imshow(filter_mask)
pl.show()
cat_squiggle_msk = np.copy(cat_squiggle_shifted)
cat_squiggle_msk[np.where(filter_mask==0.)] = 0.+0.*1j
cat_filtered = np.fft.ifftshift(cat_squiggle_msk)
cat_filtered_msk = np.fft.ifft2(cat_filtered)
cat_filtered_msk = np.abs(cat_filtered_msk)
pl.subplot(121),pl.imshow(cat, cmap = 'gray')
pl.title('Cat'), pl.xticks([]), pl.yticks([])
#pl.subplot(132),pl.imshow(20*np.log(np.abs(cat_squiggle_msk)), cmap = 'gray')
#pl.title('Filtered Fourier Cat'), pl.xticks([]), pl.yticks([])
pl.subplot(122),pl.imshow(cat_filtered_msk)
pl.title('Masked Cat'), pl.xticks([]), pl.yticks([])
pl.show()
Explanation: What about... instead of filtering out a contiguous range of Fourier frequencies, we rnadomly selected the Fourier components?
End of explanation
<END_TASK>