Toolbox for the Computation of 2D Fourier Transforms in Polar Coordinates via Maple

multi-dimensions and curvilinear coordinates. A toolbox of functions for the computation of two dimensional Fourier transforms in polar coordinates with symbolic computer algebra (Maple) was developed. The implementation of the 2D Fourier transform in polar coordinates within the toolbox is a combination of two significantly simpler transforms. A modular approach is used along with the idea of lookup tables to help avoid the issue of indeterminate results when attempting to directly evaluate the transform. This concept helps prevent unnecessary computation of already known transforms thereby saving memory and processing time.


Introduction
Symbolic computation allows a wider range of expression for mathematical formulae and their various transformation rules, while computer algebra admits greater algorithmic precision as it constructs algorithms for computing algebraic quantities in various arithmetic domains. The simplest form of this Fourier transform is the one dimensional case and symbolic computer algebra has successfully been applied to this case.
In Cartesian coordinates, the 2D case is simply two one-dimensional cases (one in each Cartesian direction). However, there are occasions when a system is best expressed in polar coordinates, prompting the need for a Fourier transform in polar coordinates. Recently, the development of the polar coordinate version of a 2D Fourier transform was documented along with the corresponding primary rules [1].
In this paper, the development of a Symbolic Computer Algebra (Maple) toolbox for the computation of algebraic, closed-form versions of the two-dimensional Fourier transform in polar coordinates is outlined.

Brief outline of the theory of 2D Fourier transforms in polar coordinates
The 2D Fourier transform of a function f (x, y) is defined as [1]: .
The inverse Fourier transform is given by where the shorthand notation of ( , ), ( , ) x y r x y w w w = =   has been used. Polar coordinates can be introduced as x = r cos θ, y = r sin θ and similarly in the spatial frequency domain as ω x = ρ cos ψ ω y = ρ sin ψ otherwise written as, r 2 = x 2 + y 2 , θ = arctan (y/x) and ρ 2 = ω 2 x + ω 2 y , ψ = arctan (ω y /ω x ). It then follows that the two-dimensional Fourier transform can be written as In terms of polar coordinates, the Fourier transform operation transforms the spatial position radius and angle (r, θ) to the frequency radius and angle (ρ, ψ). The corresponding 2D inverse Fourier transform is written as

Toolbox for the Computation of 2D Fourier Transforms in Polar Coordinates via Maple
The Fourier transform is one of the most useful tools in science and engineering and can be expanded to multi-dimensions and curvilinear coordinates. A toolbox of functions for the computation of two dimensional Fourier transforms in polar coordinates with symbolic computer algebra (Maple) was developed. The implementation of the 2D Fourier transform in polar coordinates within the toolbox is a combination of two significantly simpler transforms. A modular approach is used along with the idea of lookup tables to help avoid the issue of indeterminate results when attempting to directly evaluate the transform. This concept helps prevent unnecessary computation of already known transforms thereby saving memory and processing time.  Therefore, the nth term in the Fourier series for the original function will Hankel transform into the nth term of the Fourier series of the Fourier transform function. However, it is an nth order Hankel transform for the nth term, so that all the terms are not equivalently transformed. For reference, the nth order Hankel transform is defined by the integral [2] where J n (z) is the nth order Bessel function with the overhat indicating a Hankel transform as shown in equation (11). Here, n may be an arbitrary real or complex number. The Hankel transform is self-reciprocating and the inversion formula is the same as that given by the forward formula.
It is noted that this definition of the Hankel transform is not the same as that used by the built-in Maple function for defining Hankel transforms, therefore another function was written in order to match the definition given in equation (11).

Implementation of 2D Fourier transform in polar coordinates
Most importantly, it can be seen that the operation of finding the 2D Fourier transform F (ρ, ψ) of a function F (r , θ) is equivalent to 1) First finding its Fourier series coefficients in the angular variable f n (r), given by equation (6). 2) Finding the Fourier series coefficient of the Fourier transform, F n (ρ) via F n (ρ) = 2π i -n H n {f n (r)}. That is, by finding the nth order Hankel transform (of the spatial radial variable to the spatial frequency radial variable) of the nth coefficient in the Fourier series and appropriately scaling the result. 3) Finally, taking the inverse Fourier series transform (summing the series) with respect to the frequency angular variable, given by equation (7).

Implementation and architecture
The toolbox created for the evaluation of 2D polar Fourier transforms is called SCAToolbox. It is a symbolic computer algebra toolbox that provides a comprehensive collection of interactive tools suitable for computing mainly the two dimensional Fourier transform of expressions in polar coordinates. SCAToolbox consists of several procedures (this is the term used for "functions" in Maple) and operations as well as tables. One of the packages in SCAToolbox is named the IntegralTrans package. This package contains the procedures, tables and operations necessary for computing some integral transforms.

Set-up of the SCAToolbox
The SCAToolbox is broken into sub packages allowing it to be easy to use and maintain. The code structure is made up of four main sections: 1) creating the toolbox, 2) supporting functions, 3) integral transforms and 4) testing and verification. Creating the toolbox depends on the CAS software being used, in this case Maple.

Supporting functions
Many operations and procedures have been designed and implemented as part of this toolbox. The supporting functions are those operations that do not constitute the core of the toolbox but are vital for the complete and effective functionality of the system. These supporting functions help manage other structures within the toolbox. There are two types of supporting functions within the toolbox. The first type consists of procedures that perform operations on functions. All the procedures in this group are convolutions (one of the major reasons to implement a Fourier transform). The convolutions implemented include one dimensional and two dimensional convolutions in Cartesian coordinates, angular/circular convolution, radial convolution, two dimensional convolution in polar coordinates and series convolution. As their names suggest, these convolutions convolve different function types and so have different rules that apply to their evaluation. A "Bracket" convolution is implemented to make these convolutions easily accessible and usable -that is, one convolution procedure is defined and the type of convolution desired is passed as an option to the procedure. This makes for less error when trying to evaluate a convolution since the syntax becomes uniform for six different types of convolution. Convolutions can be stand-alone procedures and as such are not placed in the IntegralTrans package of the SCAToolbox. Using these operations therefore does not require loading the IntegralTrans package; downloading and loading the SCAToolbox is sufficient.
Procedures that help manage/manipulate other structures in the toolbox form the other type of supporting functions. These are important to the inner workings of the integral transforms and are stored in the IntegralTrans package within the SCAToolbox. Among these procedures are the takeAlook procedure, addToTable procedure, Ekronecker procedure and EdDirac procedure. The former two procedures help manage the toolbox by manipulating tables.
The takeAlook procedure manipulates tables by comparing the entered function with the list of functions in the first column of the lookup table of interest and returning the mapping result in the second column of the table. A process named pattern matching is used to compare the functions. The entered function has particular properties (arithmetic and variable types, etc.) that are matched to a fitting pattern in the table. With the concept of a lookup table comes the need to add to the table. This is accomplished by implementing the addToTable procedure. This procedure adds an expression and its corresponding transform to a given transform table and makes extending an existing table quite simple.
The Ekronecker procedure provides the Kronecker delta function given the appropriate variables. The need for this operation arose when the built-in kronecker function did not work as expected. In Maple, the integral of a shifted Dirac function is halved at the end points. However, in our theoretical work, the desired result at the endpoints is not halved. The EdDirac procedure therefore redefines the built-in Dirac function so that the result matches the theoretical development.

Integral transforms
The integral transforms implemented in this toolbox include the forward nth order Hankel transform (equation (11)) and its inverse, the forward and inverse one dimensional Fourier series transform (equation (5) and (6)) and the forward and inverse two dimensional Fourier transform in polar coordinates (equations (3) and (4)). When lookup and remember tables fail to return a result, the transforms are evaluated directly by integration. A separate function that calls the direct 2D Fourier transform (which attempts to implement the transform by direct integration, equation (3)) is also included, although we found this approach to be ineffective. If evaluation is unsuccessful, the outputs of the transform procedures are written as they are entered except in the case of the Fourier Series transform and its inverse where the actual integral and sum respectively are returned.
The Hankel transform and 1D Fourier series are implemented by using lookup tables. The 2D Fourier transform in polar coordinates is implemented via the two preceding transforms. Here, it draws on the modular nature of the code to implement the transform. This method is part of the core foundation of this work and involves breaking the 2D polar Fourier transform into three steps.
A summary of the contents of the toolbox is shown in Table 1 where f and g are functions of x, ff and gg are functions of x and y, p and q are functions of r, u and v is a function of ϑ, and pp and qq are functions of r and ϑ. Also, s and w are series in n. Additionally, expr is an arbitrary expression, tablename is the name of the table to consider, procname is the name of the procedure to consider, and ans is the result that maps to expr.

Quality control
To ensure that these operations and procedures work well, they have been tested and verified with known examples/data.
Below are examples of some simple functions and their transforms. The results shown in this section involve Maple outputs. 'MapleIn' is used to indicate the queries and 'MapleOut' gives the output of the toolbox. 'MapleCommand' represents executable commands.

Example 1
First to be tested is the 1D convolution in Cartesian coordinates. Theory: A Dirac function convolved (single dimension Cartesian convolution) with any other function returns the function i.e. f (x)*δ(x) = f (x). This is confirmed via the procedure OneDCartConv. The equivalent statement using the "Bracket" version of the convolution is OneDCartConv The above computation shows that the "Bracket" convolution works and that the operation does indeed commute.

Example 2
A list of functions is entered for which the Hankel transform of order 0 of the contents is sought. The output list has the correct order of the transforms within it.   Figure 1 gives a summary of a few functions (including the example above) that were tested and their transforms.

Operating system
Windows XP and higher.

Programming language
Maple version 12 and higher.

Additional system requirements
If using Maple 12, minimum system requirements are 512MB of RAM and 1GB of hard disk space. Higher versions of Maple will required additional memory and disk space. The persistent identifier provided above is a link to a complete fileset that includes the toolbox, instructions for downloading and using the toolbox and sample Maple code for using the toolbox.

License
This software is released under the MIT license.

Support
Please send email to nbaddour@uottawa.ca or senaedem@gmail.com.