http://www.developer.com/
2D Fourier Transforms using Java, Part 2August 9, 2005 Java Programming, Notes # 1491
PrefaceThis is the second part of a two-part lesson. The first part published earlier was entitled 2D Fourier transforms using Java. In this lesson, I will teach you how to perform two-dimensional (2D) Fourier transforms using Java. I will
Two separate programs I will present and explain two separate programs. One program consists of a single class named ImgMod30. The purpose of this class is to satisfy the computational requirements for forward and inverse 2D Fourier transforms. This class also provides a method for rearranging the spectral data into a more useful format for plotting. The second program named ImgMod31 will be used to test the 2D Fourier transform class, and also to illustrate the use of 2D Fourier transforms for some well known sample surfaces. A third class named ImgMod29 will be used to display various 3D surfaces resulting from the application of the 2D Fourier transform. I explained this class in an earlier lesson entitled Plotting 3D Surfaces using Java. Using the class named ImgMod30 The 2D Fourier transform class couldn't be easier to use. To perform a forward transform execute a statement similar to the following: ImgMod30.xform2D(spatialData,realSpect, imagSpect,amplitudeSpect); The first parameter in the above statement is a reference to an array object containing the data to be transformed. The other three parameters refer to array objects that will be populated with the results of the transform. To perform an inverse transform execute a statement similar to the following: ImgMod30.inverseXform2D(realSpect,imagSpect, recoveredSpatialData); The first two parameters in the above statement refer to array objects containing the complex spectral data to be transformed. The third parameter refers to an array that will be populated with the results of the inverse transform. To rearrange the spectral data for plotting, execute a statement similar to the following where the parameter refers to an array object containing the spectral data to be rearranged. double[][] shiftedRealSpect = ImgMod30.shiftOrigin(realSpect); Digital signal processing (DSP) This lesson will cover some technically difficult material in the general area of Digital Signal Processing, or DSP for short. As usual, the better prepared you are, the more likely you are to understand the material. For example, it would be well for you to already understand the one-dimensional Fourier transform before tackling the 2D Fourier transform. If you don't already have that knowledge, you can learn about one-dimensional Fourier transforms by studying the following lessons:
In addition, I strongly recommend that you study Part 1 of this lesson before embarking on this part. You might also enjoy studying my other lessons on DSP as well. Will use in subsequent lessons The 2D Fourier transform has many uses. I will use the 2D Fourier transform in several future lessons involving such diverse topics as:
Viewing tip You may find it useful to open another copy of this lesson in a separate browser window. That will make it easier for you to scroll back and forth among the different figures and listings while you are reading about them. Supplementary material I recommend that you also study the other lessons in my extensive collection of online Java tutorials. You will find those lessons published at Gamelan.com. However, as of the date of this writing, Gamelan doesn't maintain a consolidated index of my Java tutorial lessons, and sometimes they are difficult to locate there. You will find a consolidated index at www.DickBaldwin.com. General DiscussionThe space domain In Part 1 of this lesson, I extended the concept of the Fourier transform from the time domain into the space domain. I pointed out that while the time domain is one-dimensional, the space domain is thee-dimensional. However, in order to keep the complexity of this lesson in check, we will assume that space is only two-dimensional. This will serve us well later for such tasks as image processing.
A purely real space domain Although the space domain can be (and often is) complex, many interesting problems, (such as photographic image processing) can be dealt with under the assumption that the space domain is purely real. We will make that assumption in this lesson. This assumption will allow us to simplify our computations when performing the 2D Fourier transform to transform our data from the space domain into the wavenumber domain. PreviewI will present and explain two complete Java programs in this lesson. The first program is a single class named ImgMod30, which provides the capability to perform forward and inverse Fourier transforms on three-dimensional surfaces. In addition, the class provides a method that can be used to reformat the wavenumber spectrum to make it more suitable for display. The second program is named ImgMod31. This is an executable program whose purpose is to exercise and to test the ImgMod30 class using several examples for which the results should already be known. Sample ProgramsThe class named ImgMod30 This class provides 2D Fourier transform capability that can be used for image processing and other purposes. The class provides three static methods:
The class was tested using J2SE 5.0 and WinXP. The class uses the static import capability that was introduced in J2SE 5.0. Therefore, it should not compile using earlier versions of the compiler. The xform2D method The beginning of the class and the beginning of the static method named xform2D is shown in Listing 1. This method computes a forward 2D Fourier transform from the space domain into the wavenumber domain. The number of points produced for the wavenumber domain matches the number of points received for the space domain in both dimensions. Note that the input data must be purely real. In other words, the program assumes that there are no imaginary values in the space domain. Therefore, this is not a general purpose 2D complex-to-complex transform.
Parameters The first parameter is a reference to a 2D double array object containing the data to be transformed. The remaining three parameters are references to 2D double array objects of the same size that will be populated with the following transform results:
Listing 1 also determines and displays the dimensions of the incoming 2D array of data to be transformed. I won't bore you with the details as to how and why the 2D Fourier transform does what it does. Neither will I bore you with the details of the code that implements the 2D Fourier transform. If you understand the material that I have previously published on Fourier transforms in one dimension, this code and these concepts should be a straightforward extension from one dimension to two dimensions. The remainder of the xform2D method The remainder of the xform2D method is shown in Listing 2. Note that it was necessary to sacrifice indentation in order to force these very long equations to be compatible with this narrow publication format and still be somewhat readable.
The inverseXform2D method The inverseXform2d method is shown in its entirety in Listing 3. This method computes an inverse 2D Fourier transform from the wavenumber domain into the space domain. The number of points produced for the space domain matches the number of points received for the wavenumber domain in both dimensions. This method assumes that the inverse transform will produce purely real values in the space domain. Therefore, in the interest of computational efficiency, the method does not compute the imaginary output values. Therefore, this is not a general purpose 2D complex-to-complex transform. For correct results, the input complex data must match that obtained by performing a forward transform on purely real data in the space domain. Once again it was necessary to sacrifice indentation to force this very long equation to be compatible with this narrow publication format and still be readable.
Parameters This method requires three parameters. The first two parameters are references to 2D arrays containing the real and imaginary parts of the complex wavenumber spectrum that is to be transformed into a surface in the space domain. The third parameter is a reference to a 2D array of the same size that will be populated with the transform results. The shiftOrigin method This method deserves some explanation. The reason that this method is needed is illustrated by Figure 1. Figure 1 Two views of the same wavenumber spectrum Both of the images in Figure 1 represent the same wavenumber spectrum, but they are plotted against different coordinate systems. How the wavenumber spectrum is actually computed The top image shows how the wavenumber spectrum is actually computed. The wavenumber spectrum is computed covering an area of wavenumber space with the 0,0 origin in the upper left corner. The computation extends to twice the Nyquist folding wave number along each axis. Computationally sound but not visually pleasing While this format is computationally sound, it isn't what most of us are accustomed to seeing in plots wavenumber space. Rather, we are accustomed to seeing wavenumber spectra plotted with the 0,0 origin at the center. The wavenumber spectrum is periodic Knowing that the area of wavenumber space shown in the top image of Figure 1 covers one complete period of a periodic surface, the shiftOrigin method rearranges the results (for display purposes only) to that shown in the bottom image in Figure 1. The origin is at the center in the bottom image of Figure 1. The edges of the lower image in Figure 1 are the Nyquist folding wave numbers. The shiftOrigin method code The shiftOrigin method is shown in its entirety in Listing 4. Although this method is rather long, it is also completely straightforward. Therefore, it shouldn't require a further explanation. You may be able to develop a much shorter algorithm for accomplishing the same task.
End of the ImgMod30 class Listing 4 also signals the end of the class definition for the class named ImgMod30. The program named ImgMod31 The purpose of this program is to exercise and to test the 2D Fourier transform methods and the axis shifting method provided by the class named ImgMod30. Command line parameters The main method in this class reads two command line parameters and uses them to select:
Forward and inverse Fourier transforms The program performs a 2D Fourier transform on that surface followed by an inverse 2D Fourier transform. Six different plots are produced in this process showing different aspects of the transform and the inverse transform. Fourteen cases There are 14 different cases built into the program with case numbers ranging from 0 to 13 inclusive. Each of the cases is designed such that the results of the analysis should be known in advance by a person familiar with 2D Fourier analysis and the wavenumber domain. Thus, these cases can be used to confirm that the transform code was properly written. The cases are also designed to illustrate the impact of various space domain characteristics on the wavenumber spectrum. This information will be useful later when analyzing the results of performing 2D transforms on photographic images. A stack of output images Each time the program is run, it produces a stack of six output images in the upper left corner of the screen. A brief description of each of the output images is provided in the following list. The top-to-bottom order of the stack is:
To view the images near the bottom of the stack, you must physically move those on top to get them out of the way. Numeric output In addition, the program produces some numeric output on the command line screen that may be useful in confirming the validity of the forward and inverse transformation processes. Figure 2 shows an example of the numeric output.
The size of the surfaces The first two lines of numeric output in Figure 2 show the size of the spatial surface for the forward transform. The second two lines show the size of the wavenumber surface for the inverse transform. The quality of the transformation process The remaining nine lines indicate something about the quality of the forward and inverse transforms in terms of the ability of the inverse transform to replicate the original spatial surface. These lines also indicate something about the correctness of the overall scaling from original input to final output. Matching pairs of values Each of the last nine lines contains a pair of values. The first value is a sample from the original spatial surface. The second value is a sample from the corresponding location on the spatial surface produced by performing an inverse transform on the wavenumber spectrum. The two values in each pair of values should match. If they match, this indicates the probability of a valid result.
The match is very good in the example shown above. This example is from Case #12. How to use the program named ImgMod31 Usage information for the program is shown in Figure 3.
A description of each case is provided by the comments in this program. In addition, each case will be discussed in detail in this lesson. See ImgMod29 in the earlier lesson entitled Plotting 3D Surfaces using Java for a definition of DisplayType. You can terminate the program by clicking on the close button on any of the display frames produced by the program. Let's see some code The beginning of the class and the beginning of the main method is shown in Listing 5.
The code in Listing 5 gets the input parameters and uses them to set the case and the display format. A default case and a default display format are used if this information is not provided by the user. Create and save the test surface Listing 6 invokes the method named getSpatialData to create a test surface that matches the specified case. This surface will be used for testing the transformation process.
I will discuss the method named getSpatialData in detail later. For now, just assume that the 2D array object referred to by spatialData contains the test surface when this method returns. Display the test surface Listing 7 instantiates an object of the class named ImgMod29 to display the test surface in the display format indicated by the value of displayType.
The value of false in the third parameter indicates that the axes should not be displayed.
An example test surface plot The upper left image in Figure 4 is an example of the output produced by the code in Listing 7 for a displayType value of 0.
Figure 4 Figure 4 shows the results for a test surface switchCase value of 2. I will discuss the particulars of this case in detail later. Perform the forward Fourier transform Listing 8 performs the forward Fourier transform to transform the test surface into the wavenumber domain.
Prepare array objects to receive the transform results Listing 8 begins by preparing some array objects to receive the transform results. The forward transform receives an incoming surface array and returns the real and imaginary parts of the complex wavenumber spectrum along with the amplitude spectrum by populating three array objects passed as parameters to the method. Perform the forward transform Then Listing 8 invokes the static xform2D method of the ImgMod30 class to perform the forward transform, returning the results by way of the parameters to the method. Display unshifted amplitude spectrum The upper right image in Figure 4 is an example of the type of display produced by the code in Listing 9. This is a plot of the amplitude spectrum without the wavenumber origin being shifted to place it at the center.
This image also shows the result of passing true as the third parameter causing the red axes to be plotted on top of the spectral data.
Need to shift the origin for display The upper right image in Figure 4 is in a format that is not particularly good for viewing. In particular, the origin is at the upper left corner. The horizontal Nyquist folding wavenumber is near the horizontal center of the plot. The vertical Nyquist folding wave number is near the vertical center of the plot. It is much easier for most people to understand the plot when the wavenumber origin is shifted to the center of the plot with the Nyquist folding wave numbers at the edges of the plot. The method named shiftOrigin can be used to rearrange the data and shift the origin to the center of the plot. Shift the origin and display the results Listing 10 shifts the origin to the center of the plot and displays:
The axes are displayed in all three cases.
Example displays Examples of the displays produced by the code in Listing 10 are shown in Figure 4. The real part of the shifted wavenumber spectrum is shown in the image in the top center of Figure 4. The imaginary part of the shifted wavenumber spectrum is shown in the bottom center of Figure 4. The shifted amplitude spectrum is shown in the bottom right image in Figure 4. The origin has been shifted to the center in all three cases.
Perform an inverse transform Listing 11 performs an inverse Fourier transform to transform the complex wavenumber surface into a real surface in the space domain. Ideally, the result should exactly match the space domain surface that was transformed into the wavenumber domain in Listing 8. However, because of small arithmetic errors that accumulate in the forward and inverse transform computations, it is unusual for an exact match to be achieved.
Prepare an array object to store the results Listing 11 begins by preparing an array object to store the results of the inverse transformation process. Invoke the inverseXform2D method Then Listing 11 invokes the inverseXform2D method to transform the complex wavenumber spectrum into a real space function. The inverseXform2D method requires the real and imaginary parts of the complex wavenumber spectrum as input parameters.
The inverseXform2D method also receives an incoming array object in which to store the real result of the transformation process. Display the result of the inverse transform Finally, Listing 12 displays the result of the inverse transformation process as a surface in the space domain. This surface should compare favorably with the original surface that was transformed into the wavenumber domain in Listing 8
The output produced by Listing 12 is shown in the lower left image in Figure 4. Compare this with the input surface shown in the upper left image in Figure 4. As you can see, they do compare favorably. In fact, they appear to be identical in this grayscale plotting format. We will see later that when a more sensitive plotting format is used, small differences in the two may become apparent. Display some numeric results As discussed earlier, the code in Listing 13 samples and displays a few corresponding points on the original surface and the surface produced by the inverse transformation process. The results can be used to evaluate the overall quality of the process as well as the correctness of the overall scaling.
Each line of output text contains two values, and ideally the two values should be exactly the same. Realistically, because of small computational errors in the transform and inverse transform process, it is unlikely that the two values will be exactly the same except in computationally trivial cases. Unless the two values are very close, however, something probably went wrong in the transformation process and the results should not be trusted. Listing 13 also signals the end of the main method. What we know so far ... Now we know how to use the ImgMod30 class and the ImgMod29 class to:
It is time to for us to take a look at the method named getSpatialData that can be used to create any of fourteen standard surfaces in the space domain. The getSpatialData method This method constructs and returns a specific 3D surface in a 2D array of type double. The surface is identified by the value of an incoming parameter named switchCase. There are 14 possible cases. The allowable values for switchCase range from 0 through 13 inclusive. The other two input parameters specify the size of the surface that will be produced in units of rows and columns. The code for the getSpatialData method The getSpatialData method begins in Listing 14.
Listing 14 begins by creating a 2D array object of type double in which to store the surface. Then Listing 14 shows the beginning of a switch statement that will be used to select the code to create a surface that matches the value of the incoming parameter named switchCase. For switchCase = 0 Listing 15 shows the code that is executed for a value of switchCase equal to 0.
This case places a single non-zero point at the origin in the space domain. The origin is at the upper left corner. The surface produced by this case is shown in the leftmost image in Figure 5, and the non-zero value can be seen as the small white square in the upper left corner. In signal processing terminology, this point can be viewed as an impulse in space. It is well known that such an impulse produces a flat spectrum in wavenumber space. Figure 5 The output surface The rightmost image in Figure 5 shows the result of:
You can see the impulse as the small white square in the upper left corner of both images. The wavenumber spectrum is flat Because the wavenumber spectrum is flat, plots of the spectrum are completely featureless. Therefore, I did not include them in Figure 5. A very small error The numeric output shows that the final output surface matches the input surface to within an error that is less than about one part in ten to the fourteenth power. The program produces the expected results for this test case. If you were to go back to the equations in Listing 2 and Listing 3 and work this case out by hand, you would soon discover that the computational requirements are almost trivial. Most of the computation involves doing arithmetic using values of 1 and 0. Thus, there isn't a lot of opportunity for computational errors in this case. For switchCase = 1 Now we are going to a case that is more significant from a computational viewpoint. The input surface in this case will consist of a single impulse that is not located at the origin in the space domain. Rather, it is displaced from the origin. The wavenumber amplitude spectrum of a single impulse in the space domain should be flat regardless of the location of the impulse in the space domain. However, the real and imaginary parts of the wavenumber spectrum are flat only when the impulse is located at the origin in space. This case does not satisfy that requirement. Regardless of the fact that the real and imaginary parts are not flat, the square root of the sum of the squares of the real and imaginary parts (the amplitude) should be the same for every point in wavenumber space for this case. Thus, the real and imaginary parts are related in a very special way. The code for switchCase = 1 The code that is executed when switchCase equals 1 is shown in Listing 16.
This case places a single impulse close to but not at the origin in space. This produces a flat amplitude spectrum in wavenumber space just like in case 0. However, the real and imaginary parts of the spectrum are different from case 0. They are not flat. The computations are probably more subject to errors in this case than for case 0. The visual output Figure 6 shows the six output images produced by the program for a switchCase value of 1. Figure 6 The input and output surfaces match visually The input and output surfaces showing the single impulse are in the two leftmost images in Figure 6. From a visual viewpoint, the output at the bottom appears to be an exact match for the input at the top. The real and imaginary parts The real and imaginary parts of the wavenumber spectrum are shown in the two center images. The real part is at the top and the imaginary part is at the bottom. For a single impulse in the space domain, we would expect each of these surfaces to be a 3D sinusoidal wave (similar to a piece of corrugated sheet metal). That appears to be what we are seeing, with almost two full cycles of the sinusoidal wave between the origin and the bottom right corner of the image.
We know that the real part of a wavenumber spectrum resulting from the Fourier transform of a real space function is symmetric about the origin. We also know that the imaginary part is asymmetric about the origin. The symmetry/asymmetry requirements appear to be satisfied by this case. The color bands in the real part at the top are symmetric on either side of the origin. The imaginary part is asymmetric about the origin (the centers of the red/white and the blue/black bands appear to be equidistant from and on opposite sides of the origin). The amplitude spectrum is ugly The amplitude spectrum is shown in the two rightmost images in Figure 6. The unshifted amplitude spectrum is shown at the top. The amplitude spectrum with the origin shifted to the center is shown at the bottom. The ugliness of these two plots is an artifact of the 3D plotting scheme implemented by the class named ImgMod29. In order to maximize the use of the available dynamic range in the plot, each surface that is plotted is normalized such that:
This normalization is applied even when the distance between the highest and lowest elevation is very small. As a result of computational errors, the amplitude spectrum is not perfectly flat. Rather there are very small variations from one point to the next. As a result, the colors used to plot the surface switch among the full range of available colors even for tiny deviations from perfect flatness. A very small error The total error for this case is very small. The numeric output shows that the final output surface matches the input surface to within an error that is less than about one part in ten to the thirteenth power. The program produces the expected results for this test case. For switchCase = 2 Now we are going to take a look at another case for which we know in advance generally what the outcome should be. This will allow us to compare the outcome with our expectations to confirm proper operation of the program. A box on the diagonal in space This case places a box that is one unit tall on the diagonal near the origin in the space domain as shown in the upper left image in Figure 7. Figure 7 What do we know? On the basis of prior experience, we know that the amplitude spectrum of this surface along the horizontal and vertical axes of the wavenumber spectrum should have a rectified sin(x)/x shape (all negative values are converted to positive values). We know that the peak in this amplitude spectrum should appear at the origin in wavenumber space, and that the width of the peak should be inversely proportional to the size of the box. The code for switchCase = 2 The code that constructs the space domain surface for this case is shown in Listing 17.
This code is completely straightforward. It sets the value of each of nine adjacent points on the surface to a value of 1, while the values of all other points on the surface remain at zero. The arrangement of those nine points forms a square whose sides are parallel to the horizontal and vertical axes. The real and imaginary parts of the spectrum There isn't a lot that I can tell you about what to expect regarding the real and imaginary parts of this spectrum, other than that they should exhibit the same symmetry and asymmetry conditions that I described earlier for the real and imaginary parts in general. These requirements appear to be satisfied by the real part at the top center of Figure 7 and the imaginary part at the bottom center of Figure 7. Otherwise, the shape of the real and imaginary wavenumber spectra will depend on the location of the box in space and the size and orientation of the box. A different plotting color scheme Note that the plotting color scheme that I used for Figure 7 is different from any of the plots previously shown in this lesson. This scheme quantizes the range from the lowest to the highest elevation into 23 levels, coloring the lowest elevation black, the highest elevation white, and assigning very specific colors to the 21 levels in between. The colors and the levels that they represent are shown in the calibration scales under the plots in Figure 7. The lowest elevation is on the left end of the calibration scale. The highest elevation is on the right end of the calibration scale. The amplitude spectrum As before, the wavenumber amplitude spectrum with the origin in the upper left corner is shown in the upper right image in Figure 7. The amplitude spectrum with the origin shifted to the center is shown in the lower right image in Figure 7. If you were to use the calibration scale to convert the colors along the horizontal and vertical axes in the lower right image into numeric values, you would find that they approximate a rectified sin(x)/x shape as expected. The output surface The output surface produced by performing an inverse Fourier transform on the complex wavenumber spectrum is shown in the lower left image in Figure 7. This surface appears to match the input surface shown in the upper left image in Figure 7. The overall results The numeric output for this case isn't very useful because none of the samples for which numeric data is provided fall within the square. However, because the real and imaginary parts exhibit the correct symmetry, the shape of the amplitude spectrum is generally what we expect, and the output from the inverse Fourier transform appears to match the original input causes us to conclude that the program is working properly in this case. For switchCase = 3 This case places a raised box at the top near the origin in the space domain, but the box is not on the diagonal as it was in case 2.
Figure 8 Amplitude spectrum should not change As long as the size and the orientation of the box doesn't change, the wavenumber amplitude spectrum should be the same as case 2 regardless of the location of the box in space. Since the size and orientation of this box is the same as in case 2, the amplitude spectrum for this case should be the same as for case 2. The real and imaginary parts of the spectrum may change However, the real and imaginary parts (or the phase) change as the location of the box changes relative to the origin in space. A hypothetical example The purpose of this case is to illustrate a hypothetical example. If two different photographic images contain a picture of the same object in the same size and the same orientation in space, that object will contribute the same values to the amplitude spectrum of both images regardless of where the object is located in the different images. For example, assume that a photographic image includes a picture of a vase. Assume that the original image is cropped twice along two different borders producing two new images. Assume that both of the new images contain the picture of the vase, but in different locations. That vase will contribute the same values to the amplitude spectra of the two images regardless of the location of the vase in each of the images. This knowledge will be useful to us in future lessons when we begin using 2D Fourier transforms to process photographic images. Amplitude spectrum is the same If you compare Figure 8 with Figure 7, you will see that the amplitude spectrum is the same for both surfaces despite the fact that the box is in a different location in each of the two surfaces. However, the real and imaginary parts of the spectrum in Figure 8 are considerably different from the real and imaginary parts of the spectrum in Figure 7. The code that was used to create the surface for this case is straightforward. You can view that code in Listing 22 near the end of the lesson. For switchCase = 4 This case draws a short line containing eight points along the diagonal from upper left to lower right in the space domain. You can view this surface in the upper left image in Figure 9. You can view the code that generated this surface in Listing 22 near the end of the lesson. Figure 9 Another example of sin(x)/x On the basis of prior experience, we would expect the wavenumber amplitude spectrum, (when viewed along any line in wavenumber space parallel to the line in space), to have a rectified sin(x)/x shape. We would expect the peak of that shape to be centered on the origin in wavenumber space. We would expect the width of the peak in wavenumber space to be inversely proportional to the length of the line in space. We would expect the amplitude spectrum when viewed along any line in wavenumber space perpendicular to the line in space to have a constant value. Our expectations are borne out The shape of the amplitude spectrum shown in the lower right image in Figure 9 agrees with our expectations. Although not shown here, if we were to make the line of points longer, the width of the peak in the rectified sin(x)/x would become narrower. If we were to make the line of points shorter, the peak would become wider, as we will demonstrate in case 5. Our expectations regarding symmetry and asymmetry for the real and imaginary parts shown in the center images of Figure 9 are borne out. The real part is at the top center and the imaginary part is at the bottom center. The output from the inverse Fourier transform shown in the bottom left of Figure 9 matches the original space domain surface in the top left of Figure 9. For switchCase = 5 This case draws a short line consisting of only four points perpendicular to the diagonal from upper left to lower right. This line of points is perpendicular to the direction of the line of points in case 4. You can view the surface for this case in the upper left image of Figure 10. You can view the code that generated this surface in Listing 22 near the end of the lesson. Figure 10 Rotated by ninety degrees If you compare Figure 10 with Figure 9, you will see that the spectral result is rotated ninety degrees relative to that shown for case 4 where the line was along the diagonal. In other words, rotating the line of points by ninety degrees also rotated the structure in the wavenumber spectrum by ninety degrees. A wider peak In addition, the line of points for case 5 is shorter than the line of points for case 4 resulting in a wider peak in the rectified sin(x)/x shape for case 5. The real and imaginary parts While the real and imaginary parts of the spectrum shown in the center of Figure 10 are considerably different from anything that we have seen prior to this, they still satisfy the symmetry and asymmetry conditions that we expect for the real and imaginary parts. The final output matches the input The output from the inverse Fourier transform in the bottom left image in Figure 10 matches the input surface in the top left image in Figure 10. All of this matches our expectations for this case. For switchCase = 6 This case is considerably more complicated than the previous cases. You can view the surface for this case in the upper left image in Figure 11. You can view the code that generated this surface in Listing 22 near the end of the lesson. Figure 11 Many weighted lines of points This case draws horizontal lines, vertical lines, and lines on both diagonals. Each individual point on each line is given a value of either +1 or -1. The weights of the individual points are adjusted so that the sum of all the weights is 0. The weight at the point where the lines intersect is also 0. Black is -1, white is +1 The small black squares in the upper left image in Figure 11 represent points with a weight of -1. The small white squares represent points with a weight of +1. The green background color represents a value of 0. Symmetries on four different axes The wavenumber amplitude spectrum is shown in the bottom right image in Figure 11. As you can see from that image, performing a 2D Fourier transform on this surface produces a wavenumber amplitude spectrum that is symmetrical along lines drawn at 0, 45, 90, and 135 degrees to the horizontal. There is a line of symmetry in the amplitude spectrum for every line of points on the space domain surface. Must be zero at the wavenumber origin Because the sum of all the points is 0, the value of the wavenumber spectrum at the origin must also be zero. This is indicated by the black square at the origin in the lower right image. Peaks at the folding wave numbers This amplitude spectrum has major peaks at the folding wave number on each of the 45-degree axes. In addition, there are minor peaks at various other points in the spectrum. The real and imaginary parts As expected, the real and imaginary parts of the spectrum, shown in the center of Figure 11 exhibit the required symmetry and asymmetry that I discussed earlier. The final output The output produced by performing an inverse Fourier transform on the complex wavenumber spectrum is shown in the lower left image in Figure 11. This image matches the input surface shown in the top left image in Figure 11. For switchCase = 7 Now we are going to make a major change in direction. All of the surfaces from cases 0 through 6 consisted of a few individual points located in specific geometries in the space domain. All of the remaining points on the surface had a value of zero. This resulted in continuous (but sampled) surfaces in the wavenumber domain. Now we are going to generate continuous (but sampled) surfaces in the space domain. We will generate these surfaces as sinusoidal surfaces (similar to a sheet of corrugated sheet metal) or the sums of sinusoidal surfaces. Performing Fourier transforms on these surfaces will produce amplitude spectra consisting of a few non-zero points in wavenumber space with the remaining points in the spectrum having values near zero. Need to change the surface plotting scale In order to make these amplitude spectra easier to view, I have modified the program to cause the square representing each point in the amplitude spectrum to be five pixels on each side instead of three pixels on each side. To keep the overall size of the images under control, I reduced the width and the height of the surfaces from 41 points to 23 points. Display fewer results I suspect that you have seen all the real parts, imaginary parts, and unshifted amplitude spectra that you want to see. Therefore, at this point, I will begin displaying only the input surface, the amplitude spectrum, and the output surface that results from performing an inverse Fourier transform on the complex spectrum. A zero frequency sine wave The first example in this category is shown in Figure 12. The input surface for this example is a sinusoidal wave with a frequency of zero. This results in a perfectly flat surface in the space domain as shown in the leftmost image in Figure 12. This surface is perfectly flat and featureless. Figure 12 The code for this case The code that was used to generate this surface is shown in Listing 18. For the case of a sinusoidal wave with zero frequency, every point on the surface has a value of 1.0.
A single point at the origin As shown by the center image in Figure 12, the Fourier transform of this surface produces a single point at the origin in wavenumber space. This is exactly what we would expect. The inverse transform output is ugly The result of performing an inverse Fourier transform on the complex spectrum is shown in the rightmost image in Figure 12. As was the case earlier in Figure 6, the ugliness of this plot is an artifact of the 3D plotting scheme implemented by the class named ImgMod29. The explanation that I gave there applies here also. A very small error Once again, the total error is very small. The numeric output shows that the final output surface matches the input surface to within an error that is less than about one part in ten to the thirteenth power. Thus, the program produces the expected results for this test case. For switchCase = 8 This case draws a sinusoidal surface along the horizontal axis with one sample per cycle.
Thus, it is impossible to distinguish this surface from a surface consisting of a sinusoid with a frequency of zero. The code that was used to produce this surface is shown in Listing 19. This code is typical of the code that I will be using to produce the remaining surfaces in this lesson. This code is straightforward and shouldn't require further explanation.
The graphic output The Fourier transform of this surface produces a single peak at the origin in the wavenumber spectrum just like in Figure 12. I didn't provide a display of the graphic output for this case because it looks just like the graphic output shown for the zero frequency sinusoid in Figure 12. For switchCase = 9 This case draws a sinusoidal surface along the horizontal axis with two samples per cycle as shown in the leftmost image in Figure 13. This is the Nyquist folding wavenumber. The wavenumber spectrum The center image in Figure 13 shows the wavenumber amplitude spectrum for this surface. The wavenumber spectrum has white peak values at the positive and negative folding wave numbers. The colors in between these two peaks are green, blue, and gray indicating very low values. Figure 13 The inverse Fourier transform output The output from the inverse Fourier transform performed on the complex wavenumber spectrum for this case is shown in the rightmost image in Figure 13. The output is a good match for the input. You can view the code that was used to create this surface in Listing 22 near the end of the lesson. For switchCase = 10 This case draws a sinusoidal surface along the vertical axis with two samples per cycle. Again, this is the Nyquist folding wave number but the sinusoid appears along the vertical axis instead of appearing along the horizontal axis. If you run this case and view the results, you will see that it replicates the results from case 9 except that everything is rotated by ninety degrees in both the space domain and the wavenumber domain. For switchCase = 11 This case draws a sinusoidal surface along the horizontal axis with eight samples per cycle as shown in the leftmost image of Figure 14. Figure 14 The wavenumber spectrum Performing a forward Fourier transform on this surface produces symmetrical peaks on the horizontal axis on either side of the wavenumber origin. The two peaks are indicated by the small white and red squares on the horizontal axis in the center image in Figure 14.
For a sinusoidal surface with eight samples per cycle, we would expect the peaks to occur in the wavenumber spectrum about one-fourth of the distance from the origin to the folding wavenumber. Figure 14 meets that expectation. The peaks are surrounded on both sides by blue and aqua colors, indicating very low values. The inverse Fourier transform output The output from the inverse Fourier transformed performed on the complex spectrum is shown in the rightmost image in Figure 14. This output compares very favorably with the input surface shown in the leftmost image. The difference between the two is that the input has white vertical bands whereas the output has red vertical bands (with a single white spot). The above explanation of white versus red applies here also. You can view the code that created this surface in Listing 22 near the end of the lesson. For switchCase = 12 This case draws a sinusoidal surface on the horizontal axis with three samples per cycle plus a sinusoidal surface on the vertical axis with eight samples per cycle as shown by the leftmost image in Figure 15. Figure 15 The wavenumber spectrum Performing a forward Fourier transform produces symmetrical peaks on the horizontal and vertical axes on all four sides of the wave number origin. These peaks are indicated by the red and white squares in the center image in Figure 15.
The peaks on the vertical axis should be about one-fourth of the way between the origin and the folding wavenumber. This appears to be the case. The peaks on the horizontal axis should be about two-thirds of the way between the origin and the folding wavenumber, which they also appear to be. Inverse Fourier transform output The output produced by performing an inverse Fourier transform on the complex spectrum is shown in the rightmost image in Figure 15. Taking the red versus white issue into account, this output compares favorably with the input surface shown in the leftmost image in Figure 15. You can view the code that created this surface in Listing 22 near the end of the lesson. For switchCase = 13 This case draws a sinusoidal surface at an angle of approximately 45 degrees relative to the horizontal as shown in the leftmost image in Figure 16. This sinusoid has approximately eight samples per cycle. Figure 16 The wavenumber spectrum Performing a forward Fourier transform on this surface produces a pair of peaks in the wavenumber spectrum that are symmetrical about the origin at approximately 45 degrees relative to the horizontal axis. These peaks are indicated by the red and white squares in the center image in Figure 16. The inverse Fourier transform output The output produced by performing an inverse Fourier transform on the complex wavenumber spectrum is shown in the rightmost image in Figure 16. This output compares favorably with the input surface shown in the leftmost image in Figure 16. You can view the code that created this surface in Listing 22 near the end of the lesson. The end of the getSpatialData method Listing 20 shows the end of the method named getSpatialData and the end of the class named ImgMod31.
A default case is provided in the switch statement to deal with the possibility that the user may specify a case that is not included in the allowable limits of 0 through 13 inclusive. The method ends by returning a reference to the array object containing the 3D surface that was created by the selected case in the switch statement. Run the ProgramI encourage you to copy, compile, and run the programs that you will find in Listing 21 and Listing 22 near the end of the lesson.
Modify the programs and experiment with them in order to learn as much as you can about 2D Fourier transforms. Create some different test cases and work with them until you understand why they produce the results that they do. SummaryThen I introduced you to some practical examples showing how 2D Fourier transforms and wavenumber spectra can be useful in solving engineering problems involving antenna arrays. In Part 2, I provided and explained a class that can be used to perform forward and inverse 2D Fourier transforms, and can also be used to shift the wavenumber origin from the upper left to the center for a more pleasing plot of the wavenumber spectral data. Finally, I provided and explained a program that is used to:
What's Next?I will explain general purpose 2D convolution and will explain and demonstrate the relationships that exist between 2D Fourier transforms and 2D convolution in the next lesson in this series. Complete Program ListingsComplete listings of the classes presented in this lesson are provided in Listing 21 and Listing 22 below.
Copyright 2005, Richard G. Baldwin. Reproduction in whole or in part in any form or medium without express written permission from Richard Baldwin is prohibited. About the authorRichard Baldwin is a college professor (at Austin Community College in Austin, TX) and private consultant whose primary focus is a combination of Java, C#, and XML. In addition to the many platform and/or language independent benefits of Java and C# applications, he believes that a combination of Java, C#, and XML will become the primary driving force in the delivery of structured information on the Web.Richard has participated in numerous consulting projects and he frequently provides onsite training at the high-tech companies located in and around Austin, Texas. He is the author of Baldwin's Programming Tutorials, which have gained a worldwide following among experienced and aspiring programmers. He has also published articles in JavaPro magazine. In addition to his programming expertise, Richard has many years of practical experience in Digital Signal Processing (DSP). His first job after he earned his Bachelor's degree was doing DSP in the Seismic Research Department of Texas Instruments. (TI is still a world leader in DSP.) In the following years, he applied his programming and DSP expertise to other interesting areas including sonar and underwater acoustics. Richard holds an MSEE degree from Southern Methodist University and has many years of experience in the application of computer technology to real-world problems. |