Let g(x) = |g(x)| exp(if(x)) be the amplitude function of the electric field in the object plane of the first lens of an optical processor, E(x) µ g(x) exp[i(kx-wt).
Let the processor have lenses with a focal length of 500 mm. Let the image in the object plane of the first lens be a "grating" with 3 slits of width 100 micrometers each, separated by 200 micrometers. Let g(x) be real, g(x) = 1 over the slits, g(x) = 0 everywhere else.
We can find E(kx'') and I(kx'') in the frequency plane of the processor using using Microsoft Excel's FFT function.
The Fourier Analysis tool is a part of the Analysis ToolPak. It can be used to analyzes periodic data by using the Fast Fourier Transform (FFT) method to transform the data. This tool also supports inverse transformations, in which the inverse of transformed data returns the original data. The tool can be used to analyze one-dimensional data.
To access the tool click "Tools, Data Analysis, Fourier Analysis" on Excel's toolbar.
 |
Excel works with discrete data. One can define a Fourier transform for a discrete series of points called the discrete Fourier transform (DFT). Excel has a build in Fast Fourier Transform (FFT) algorithm. The FFT is simply a DFT that is faster to calculate on a computer. Assume that we have a periodic function with period L which we are sampling at N discrete points.
f(xn) = S-¥+¥C(km)exp(ikmxn)
~ Sm=0N-1C(km)exp(ikmxn) If we choose our units so that L = 1 and xn = n(L/N) = n/N and km = m2p, then we can write
f(n) = Sm=0N-1C(m)exp(i2pmn/N) 0 £ n £ N-1, 0 £ m £ N-1 For N = 2n data, Excel's FFT and inverse FFT implement the relationships
f(n) = (1/N)Sm=0N-1C'(m)exp(i2pmn/N) If f(n) is real, we can rewrite the above equations in terms of a summation of sine and cosine functions with real coefficients.
f(n) = (1/N)Sm=0N-1[a(m)cos(2pmn/N) + b(m)sin(2pmn/N)] The C(m) are complex. C(0) is the constant contribution and C(N/2) corresponds to the Nyquist frequency. The Nyquist criterion states when sampling discrete data points at some sampling frequency, fs, one can obtain reliable frequency information only for frequencies less than fs/2.
[Assume you sample an object over a period L and your
sample contains feature with period L/2. The frequency of this feature
is f2 = 2/L or k2 = 2*2p/L.
With only 2 evenly spaced data points, it is impossible to recognize features that have a frequency f2.
With 4 evenly spaced data points, it is possible to recognize features that have a frequency f2, except for certain phase differences between sample and feature.] The C(k>N/2) correspond to negative frequencies and for a real sequence x, they are complex conjugates of the C(k<N/2). The magnitude plot of C(k) versus k is perfectly symmetrical about the Nyquist frequency. The useful information in the signal is found in the range 0 to fs/2. |
The input range for the Fourier Analysis tool can be range of real or complex data. Complex data must be in either x+yi or x+yj format. If x is a negative number, precede it with an apostrophe ( ' ). The number of input range values must be an even power of 2. The maximum number of values is 4096.
For the analysis of periodic data the input range should hold one period of the data. For the analysis of non-periodic data, the data to be analyzed should be restricted to a fraction of the input range and the rest of the input range should be zero. In this approximation the input range represents all space, (i.e. space is restricted to a box), and the data occupy a fraction of space.
Let the input range contain N values labeled by x(0) to x(N-1). After calculating the Fourier transform of the input range the output range will contain N values. The absolute values of the entries x(1) to x((N/2)-1) will be mirrored by the absolute values of the entries x((N/2)+1)) to x(N-1). But the entries are complex numbers and, in general, the phases will not be mirrored.
If we use this output range as the input range for an inverse FFT, the output range of this "inverse transform" will contain the original data.
To simulate our optical processor let us use an input range with N = 1024 values. Let the maximum object size (i.e. the maximum size of the image to be analyzed) be 10.24 mm. Our sampling "frequency" therefore is 1/(0.01) mm = 100 mm-1.
If our range of x-values is L = 10.24 mm, then the range of values kx'' Excel will include is kx'' = (2p/L)m = (2p/(10.24 mm))m, where m is an integer ranging from 0 to 1023.
Since kx'' @ kq and x' @ fq, the value x'' corresponding to kx'' in the frequency plane is x'' = (f/k)kx'' = (fl/2p)(m2p/L) = mfl/L for the first N/2 values of kx''. For the second N/2 values of kx'', whose absolute values mirror those of the first N/2 values, we choose x'' = -(N-m)fl/L, i.e. we let x'' be negative to produce a plot the absolute value of the Fourier transform that is symmetric about x'' = 0.
The linked Excel spreadsheet contains columns for x, kx'', and x'', so that it is easy to plot intensities in the object plane, frequency plane, and image plane versus distance and to observe the simulated patterns.
The spreadsheet contains columns for the object, the transform, and the inverse transform. The intensity in the object plane is represented by "object^2", the intensity in the frequency plane is represented by "trans^2", and the intensity in the image plane is represented by "inv. trans^2". The spreadsheet also contains a column called "mask". In this column we can manipulate 1 mm in the middle of the 10.24 mm object without having to page down to the middle of the object column. The changes we make in the "mask" column are automatically copied to the middle of the object column.
The graphs above shows the intensities in the object, frequency, and image planes of our optical processor, as simulated using Excel's FFT using the "grating" described above as the object.
We can now manipulate the transform in frequency plane and observe how this affects the synthesized image in the image plane for different objects in the object plane.
To change the pattern in the object plane, manipulate the mask. To calculate the Fourier transform of the new pattern, click "Tools, Data Analysis, Fourier Analysis" on Excel's toolbar. Choose the "object" column (F3:F1026) for the input range and the "transform" column (H3) for the output range. To manipulate the the synthesized image, apply filters and masks in the transform column. For example, to apply low-pass filter, set everything except the central peak in the transform column equal to zero. (Click the LP Filter button.) To apply high-pass filter, set the central region in the transform column equal to zero. (Click the HP Filter button.) To calculate the inverse Fourier transform, click "Tools, Data Analysis, Fourier Analysis" on Excel's toolbar. Choose the "transform" column (H3:H1026) for the input range and the "inv. trans" column (J3) for the output range.
[For the spreadsheet to work the "Add-Ins" Analysis ToolPak and Analysis ToolPak-VBA have to be installed and checked. The spreadsheet contains macros for the FFT and the inverse FFT. If macros are enabled, you can just click the appropriate button to calculate the FFT and the inverse FFT. If the macros do not run on your computer, click Tools, Macro, Visual Basic Editor. In the editor check that funcres.xla, procdb.xla, and atpvbaen.xla are loaded. Click Tools, Preferences, and browse for the missing files on your hard disk. You should find them in a directory with a name similar C:\Program Files\Microsoft Office\Office10\Library\Analysis.]
A low- pass filter applied to our grating will yield the results shown below.
Spatial filtering is the manipulation of an image with masks in the frequency plane. Low-pass and high-pass filters are common spatial filters. A simple low-pass filter consists of a pinhole located on the optical axis in the frequency plane. It blocks all the high-frequency (short wavelength) components in the Fourier transform of the object. A simple high-pass filter consists of a small, opaque spot located on the optical axis in the frequency plane. It blocks all the low-frequency (long wavelength) components in the Fourier transform of the object.
Why does anyone build an optical image processor if it is so easy to manipulate digital images?
In ordinary microscopy, for example, most of the detail of living cells is undetectable because there is too little contrast between structures with similar transparency and there is insufficient natural pigmentation. Edge enhancement of the digital images does not work, because the images are very smooth. However the various organelles show wide variation in index of refraction. These index of refraction variations can be transformed into brightness variations by means of phase contrast microscopy.
Phase contrast microscopy uses phase plates, which when placed in the frequency plane change the phase of the low-frequency components by p/2, while leaving the phase of the high frequency components unchanged.
Assume a sample is transparent and only changes the phase of the transmitted wave by a small amount because of small index of refraction changes.
Aexp(i(kx-wt)) --> Aexp(-id(x))exp(i(kx-wt))
If d(x) is small (d(x)
<< 1) the transmitted amplitude may be written as A(1- id(x)).
The Fourier transform of such a function contains a large central peak around k
= 0 representing the constant term and much smaller peaks at larger values of k
representing the small variations d(x). The
synthesized image amplitude is A(1- id(x')) and the
image intensity is
|A(1- id(x'))|2 ~ A2.
Assume we use a quarter wave plate in the frequency plane (the focal plane of a lens) to change the phase of the central spot by exp(-ip/2) = -i. The synthesized image amplitude then is approximately A(-i - id(x')) = -iA(1 + d(x')) and the image intensity is
|iA(1+ d(x'))|2 ~ A2(1+2d(x'))
The phase contrast has been converted into an amplitude contrast and the object becomes visible.
Link: The excel spreadsheet from above with a phase contrast button.
