UFDC Home  Search all Groups  UF Institutional Repository  UF Institutional Repository  UF Theses & Dissertations  Vendor Digitized Files   Help 
Material Information
Subjects
Notes
Record Information

Full Text 
THE GENERALIZED MATRIX PRODUCT AND FAST FOURIER TRANSFORM FOR PERMUTOHEDRAL AGGREGATES By JAIME L. ZAPATA A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2000 To my sister, Luz Helena. In her presence I get to be everything that is good, powerful, and worthy. ACKNOWLEDGEMENTS First of all, I would like to thank Dr. Ritter, friend and mentor, for introducing me to image algebra, for giving me the opportunity to do independent research, and for serving as chair of my doctoral committee. Thanks to the other members of my committee: Dr. David Wilson, Dr. Joseph Wilson, Dr. James Keesling, and Dr. Murali Rao. I am indebted to them for their support and encouragement and for everything they taught me. Special thanks to Dr. Andrew Vince for carefully going over everything I had written and for helping me to simplify the presentation of ideas. Thanks to the faculty of the Department of Mathematics at the University of Florida. Nearly everyone was a teacher or friend at one point or another. Special thanks to Dr. Beverly Brechner for giving me a deep appreciation of topology; to Dr. Neil White for listening at a time when I was talking about the protein folding problem; to Dr. Gang Bao for his encouragement and mentorship; and to Mrs. TownsendI enjoyed working with you. More thanks to Sandy Gagnon, Sharon Easter, Sonia Simmons, and Julia Porchiazzo, all with the staff at the Department of Mathematics. Very special thanks to Gretchen Garrett. All of you helped me with the little things and the big things that make a graduate career possible. Thanks to my parents for creating an openended possibility that bears my name, for their unconditional love and support, and for being such powerful witnesses. Thanks to Lily, Rosana, and Luz Helena for witnessing me as a brother and a scholar, among other things. Finally, I would like to thank Haylee for coming into my life at the right time, for making me laugh, for her love and help during the difficult times, and for having the courage to say goodbye. TABLE OF CONTENTS ACKNOWLEDGMENTS ............................. iii LIST OF FIGURES ................................ vii ABSTRACT .................................... viii CHAPTERS 1 INTRODUCTION ............................. 1 1.1 Hexagonal and Truncated Octahedral Aggregates ........... 4 1.2 Overview ............................ ... 5 1.3 Summary of Results ......................... 6 2 THE DISCRETE FOURIER TRANSFORM .............. 7 2.1 Introduction ............... ... .. .. .. ..... 7 2.2 Vector Spaces Over Finite Abelian Groups .............. 7 2.3 The Discrete Fourier Transform and Inverse DFT ........ 11 2.4 The Convolution Theorem .......................... 17 2.5 The CooleyTukey Factorization .................... 20 3 QUOTIENT LATTICES ......................... 23 3.1 Introduction ............................. 23 3.2 The DFT and Inverse DFT for Quotient Lattices .......... 23 3.3 The CooleyTukey Factorization for Quotient Lattices ..... 28 4 PERMUTOHEDRAL AGGREGATES ................... 31 4.1 Introduction ............................. 31 4.2 Geometric Description ........... ..................... 32 4.3 Quotient Lattices Associated with Permutohedral Aggregates . 35 4.4 Aggregates as Sets of Coset Representatives . . . . .... ..41 4.5 The CooleyTukey Factorization for Permutohedral Aggregates 48 5 THE GENERALIZED MATRIX PRODUCT .... . . . . .. 50 5.1 Introduction ....................................50 5.2 Digital Indexing Sets . . . . . . . . . . ... .. 51 5.3 The Generalized Matrix Product for Permutohedral Aggregates 52 5.4 Digital Addressing for Permutohedral Aggregates . . . ... ..54 6 DIGITAL ALGORITHMS FOR PERMUTOHEDRAL AGGREGATES 57 6.1 Introduction . . . . . . . . . . . . . ... .. 57 6.2 The Discrete Fourier Transform and Inverse DFT . . . ... ..58 6.3 CooleyTukey Factorization Algorithms . . . . . . . 59 6.4 The Fast Fourier Transform and Inverse FFT Algorithms . 64 7 SUMMARY AND SUGGESTIONS FOR FURTHER RESEARCH . 67 REFERENCES . . . . . . . . . . . . . . . . ... .. 69 BIOGRAPHICAL SKETCH . . . . . . . . . . . . ... .. 71 LIST OF FIGURES 4.1 Level 1 hexagonal aggregate . . . . . . . . . . ... .. 33 4.2 Hierarchical structure of a level 3 hexagonal aggregate . . ... .. 34 4.3 Rhombic dodecahedron . . . . . . . . . . .... .. 35 4.4 Truncated octahedron . . . . . . . . . . . .... ..36 4.5 Level 1 truncated octahedral aggregate . . . . . . . ... .. 37 4.6 Level 1 hexagonal aggregate in the spatial domain . . . ... .. 45 4.7 Level 2 hexagonal aggregate in the spatial domain . . . ... .. 46 4.8 Level 2 hexagonal aggregate in the frequency domain . . . ... .. 47 5.1 Indexed level 3 hexagonal aggregate in the spatial domain ..... .. 56 Abstract of a Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy THE GENERALIZED MATRIX PRODUCT AND FAST FOURIER TRANSFORM FOR PERMUTOHEDRAL AGGREGATES By Jaime L. Zapata May 2000 Chairman: Dr. Gerhard X. Ritter Major Department: Mathematics Permutohedral aggregates are hierarchical arrangements of cells. Special cases include hexagonal aggregates and truncated octahedral aggregates in dimensions 2 and 3, respectively. These cells can be efficiently addressed using a scheme known as generalized balanced ternary, or GBT. The cells of a Ddimensional permutohedral aggregate are addressed, or indexed, using finite strings of a set of n = 2D+1 1 digits. The objects of interest in this work are functions whose domains are permutohedral aggregates. The generalized matrix product, also known as the pproduct, is a generaliza tion of standard matrix multiplication. The present work defines the discrete Fourier transform and the pproduct for permutohedral aggregates, and derives a radixn, decimationinspace, fast Fourier transform (FFT) for permutohedral aggregates in terms of the pproduct. Data reordering (also known as shuffle permutations) is generally associated with FFT algorithms. However, use of the pproduct makes data reordering unnecessary. CHAPTER 1 INTRODUCTION The goal of this dissertation is to describe algorithms that efficiently compute the discrete Fourier transform and the corresponding inverse transform for functions with discrete domains known as permutohedral aggregates. The Fourier transform is one of the most fundamental transforms in signal analysis and image processing. This transform extracts the frequency content of a function. Permutohedral aggregates are hierarchical arrangements of certain points in Ddimensional Euclidean space. These aggregates are of interest because they possess several useful properties: 1. they can be addressed efficiently using finite strings of digits; 2. they have a ring structure and addition and multiplication can be done entirely in terms of bit string operations; 3. in dimension two, aggregates approximate the region of lightsensing cells in the human visual system; 4. in dimensions two and three, aggregates represent optimal sampling strategies in the spatial domain. Algebraically, Alevel permutohedral aggregates are sets of coset representa tives of quotient lattices of the form L/NAL, where L is a Ddimensional lattice having a ring structure, N is an endomorphism of L, A is a positive integer, and N'L = NA(L). These quotient lattices consist of n' elements, where n = 2D+1 1. The Voronoi cell of a point x in a lattice I is the set {y C RD : y xa < ly zJ, for all z E I}. The Voronoi regions or cells associated with permutohe dral aggregates are known as permutohedra. From an image processing perspective, the Voronoi cells associated with aggregate elements play the role of picture elements, or pixels. Important special cases include hexagonal cells for D = 2, and truncated octahedral cells for D = 3. The cells of a hexagonal aggregate are centered on a hexagonal lattice and the cells of a truncated octahedral aggregate are centered on a bodycentered cubic lattice. The simplest (level 1) hexagonal aggregate can be visual ized as a hexagonal cell with six neighboring cells, while a level 2 hexagonal aggregate consists of 49 cells. A central cell in a truncated octahedral aggregate has fourteen neighbors. See the figures in Chapter 4. The cells in a Alevel permutohedral aggre gate can be addressed using strings of length A from the set of digits {0, 1,..., n 1}. This form of addressing, or indexing, is known as generalized balanced ternary (GBT), and was developed by Gibson and Lucas [8]. The majority of image processing techniques and algorithms can be broadly classified as either spatial domain or frequency domain techniques. Spatial domain techniques include morphological transforms, while frequency domain techniques in clude highpass and lowpass filtering (see Gonzalez and Woods [9], or Ritter and Wilson [18]). Frequency domain techniques depend on the discrete Fourier trans form (DFT), a linear transform that extracts the frequency content of images with discrete domains. This dissertation is concerned with the computation of the dis crete Fourier transform and the inverse DFT of functions, or images, whose domains are permutohedral aggregates. The geometric duals of permutohedral aggregates are function domains in the frequency space. The fast Fourier transform (FFT) is an efficient algorithm for the computation of the discrete Fourier transform. Gilbert Strang [19] calls the FFT "the most valuable numerical algorithm in our lifetime." The complexity of the FFT is O(M log M), where M is the cardinality of the discrete image domain. Mersereau and Speak [15], and Dudgeon and Mersereau [7] describe the FFT for periodic functions defined on multidimensional lattices, and Tolimieri, An, and Lu [20, 21], describe the FFT for functions defined over Zn, the ring of integers modulo n. The FFT is based on the CooleyTukey factorization of the discrete Fourier transform, which was published in 1965 [5]. However, some unpublished manuscripts by Gauss contained the essentials of the CooleyTukey factorization and the FFT algorithm around 1805 (see Heideman et al. [11]). Data reordering, also known as shuffle permutations, is generally associated with FFT algorithms. Zhu and Ritter [24, 26], derived algorithms for the computa tion of the FFT of images defined on rectangular arrays that did not require a data reordering stage. They accomplished this by using the generalized matrix product, a generalization of standard matrix multiplication. This product, also known as the pproduct, was first described by Ritter [17]. The present work defines the gener alized matrix product for arrays of data indexed by the GBT addressing scheme. Such a definition is necessary for the description of image processing algorithms for permutohedral aggregates and their duals. This dissertation exploits the hierarchical structure of permutohedral aggre gates and the GBT addressing scheme to derive pyramidal data processing algorithms that compute the FFT and the inverse FFT for aggregates. These algorithms have complexity O(M log M), where M = n' and n = 2D+l 1; they are expressed in terms of the generalized matrix product for permutohedral aggregates. 1.1 Hexagonal and Truncated Octahedral Aggregates We include some further remarks regarding hexagonal and truncated octahe dral aggregates since properties of these two special cases of permutohedral aggregates motivated the present study. Watson and Ahumada, Jr [23] derived a mathematical transform that models spacefrequency coding of images in the primary visual cortex of the brain. They noted that while retinal ganglion cells represent visual images with a spatial code, the cells of the primary visual cortex code spatial, frequency, and orientation information about an image. Furthermore, ganglion cell receptive fields form an approximate hexagonal lattice near the fovea, the region in the central portion of the retina that allows human beings to resolve fine details in an image. The input lattice for the image transform presented by Watson and Ahumada, Jr. is a hexagonal array that is equivalent to the hexagonal aggregates described here. Petersen and Middleton [16] provide further motivation for the study of multi dimensional functions sampled on hexagonal and bodycentered cubic lattices. They showed that the most efficient sampling lattices are not in general rectangular. Ef ficient sampling lattices attempt to use a minimum number of sampling points to exactly reconstruct bandlimited functions. In particular, sampling twodimensional isotropic functions on a hexagonal lattice and sampling threedimensional isotropic functions on a bodycentered cubic lattice is significantly more efficient than sam pling such functions on simple square and cubic lattices, respectively. The defining property of an isotropic function is that the Fourier transform of the function cuts off at the same magnitude in all directions. Thus, the support of the Fourier transform of an isotropic function is a hypersphere or ball. Petersen and Middleton specify minimum sampling lattices for isotropic functions in dimensions one through eight. 1.2 Overview Chapter 2 presents a number of standard definitions and results from finite abelian group harmonic analysis: the discrete Fourier transform, the inverse DFT, the inversion theorem, the convolution theorem, and the CooleyTukey factorization theorem. Chapter 3 covers the discrete Fourier transform, the inverse DFT, and the CooleyTukey factorizations for finite abelian groups that arise as quotient lattices and their geometric duals, respectively. Permutohedral aggregates and the quotient lattices AN associated with per mutohedral aggregates are described in Chapter 4. Permutohedral aggregates are particular sets of coset representatives of the quotients A\. Chapter 4 also discusses the subgroups B# in A\, which lead to quotients of the form Qa = A\/BO, and con tains statements of the CooleyTukey factorizations for the quotient lattices A,\ and the geometric duals A', respectively. Chapter 5 defines Findiceslexicographically ordered finite strings of digits and the generalized matrix product over Findices. The Findices are used to address the elements of permutohedral aggregates and their duals. Two theorems in Chapter 6 are formulations of the DFT and the inverse DFT for permutohedral aggregates in terms of the generalized matrix product over F indices, or digital indexing sets. Chapter 6 also states four main digital algorithms in terms of the generalized matrix product for permutohedral aggregates: the Cooley Tukey (CT) factorization, the inverse CT factorization, the fast Fourier transform, and the inverse FFT. Chapter 7 contains suggestions for further research. 1.3 Summary of Results The main results established in this work are as follows: 1. algebraic formulation of the discrete Fourier transform (DFT) and the inverse DFT for permutohedral aggregates; 2. definition of the generalized matrix product for permutohedral aggregates; 3. discrete Fourier transform and inverse DFT algorithms in terms of the gener alized matrix product for permutohedral aggregates; 4. CooleyTukey factorization and inverse CT factorization algorithms in terms of the generalized matrix product for permutohedral aggregates; 5. a decimationinspace, fast Fourier transform (FFT) algorithm for permutohe dral aggregates; 6. a decimationinfrequency, inverse FFT algorithm for permutohedral aggre gates. CHAPTER 2 THE DISCRETE FOURIER TRANSFORM 2.1 Introduction The purpose of this chapter is to present background material from finite abelian group harmonic analysis. Section 2.2 describes the normed vector space L'(A) of complexvalued functions on a finite abelian group A, the algebraic dual A, and the space L2(A). Section 2.3 presents some standard definitions and results for the discrete Fourier transform and the inverse DFT for finite abelian groups. We recall that the Fourier transform is an isometry and state the inversion theorem. Section 2.4 defines cyclic convolution for functions in L(A) and presents the convolution theorem for finite abelian groups. In Section 2.5 we prove the CooleyTukey factorization theorem. Denote the set of positive integers by Z+ and the field of complex numbers byC. 2.2 Vector Spaces Over Finite Abelian Groups Denote by card(S) the cardinality of a set S. Let A be a finite abelian group. The set L(A) of all complexvalued functions defined on A is a vector space over C with addition and multiplication given by (a + b)(r) = a(r) + b(r), for all a,b E L(A),r E A, (Ca)(r) = Ca(r), for all a E L(A), E C,r E A. For a,b E L(A), define ab E L(A) by (ab)(r) = a(r)b(r), for all r E A. Define an inner product on L(A) as follows. For a, b E L(A), let (a,b) = Z a(r)b(r), rEA where b(r) denotes the complex conjugate of b(r). For a E L(A), define the norm of a by / \ 1/2 ail = V(a,a) = ( Ja(r)2 12 Denote by L2(A) the inner product space L(A) together with the norm ]. We now describe a basis of L2(A). For r, s C A, define the Kronecker delta over indices in A by 3 1 if r =s ur= 0 if r s. For r E A, define er E L(A) by er(s) = Sr,, for all s E A. We will also have occasion to refer to the support of a function, defined as follows. If f is any function with domain A, R is any ring, and f : A 4 R, then the support of f, written supp(f), is defined by supp(f) = {d E A f(d) 0}. Thus, supp(er) = r, for all r E A. Proposition 1 The set C = {er r E A} is an orthonormal basis of L2(A). Proof. If a E L2(A), then a = Z a(r)er. rEA Hence, the set of vectors {er : r E A} spans the space L 2(A). If r and s are distinct elements of A, then (er, e.) = E er(t)e.(t) = e,(s) + e.(r) = J'. + S.r = 0. tEA Thus, is a set of (nonzero) orthogonal vectors, which implies that it is a linearly independent set. Since C also spans L2(A), it is a basis of L2(A). If r = s, then (er,e) = er,) = Zer() e r(() = 1. tEA Hence, is an orthonormal basis of L2(A). D Refer to as the standard basis of L2(A). If card(A) = n, then L2(A) is an ndimensional normed vector space, since card(E) = card(A). Section 2.3 describes a second basis of L2(A). Let U = {z E C: jzj = 1}, the unit circle, and for n E Z+, let U(n) denote the multiplicative group of all the nth roots of unity. Define the algebraic dual A of a finite abelian group A of cardinaliy n by A = Hom(A, U(n)), the set of all homomorphisms from A into U(n). Thus, 4 E A implies 0 : A + U(n) such that 4)(r + s) = 4)(r)}4(s), for all r, s E A. Define addition in A as follows: for ,), b E A, let (0 + 0)(r) = O(r)t(r), for all r E A. Note that kOc(r) = O(r)k, for all k E Z, r E A. This addition renders A an abelian group. Denote by 0o the zero element or trivial homomorphism in A. Hence, 00o(r) = 1, for all r E A. A standard result for finite abelian groups follows. Theorem 2 If A is a finite abelian group with card(A) = n and A = Hom(A, U(n)), then A A. In particular, card(A) = n. The elements of the dual group A enjoy the following properties. Proposition 3 If A is a finite abelian group and 4) E A, then 1. 0(O) = 1; 2. I1)(r)12 = O(r)O(r) = 1, for all r E A; 3. The additive inverse 4) of 4) satisfies O(r) = 4)(r)' = O)(r), for all r E A. Proof. Any homomorphism 4): A * U(n) maps the identity element 0 E A to 1, the identity element of U(n). Hence, 0(0) = 1,for all E A. which proves (1). If r E A, then ](r)l2 is the norm of O(r) E U(n). The elements of U(n) lie on the unit circle. Hence, IO(r)l2 = 1,for all r E A, which proves (2). If 0 is the additive inverse of c A, then (0 + (O))(r) = Oo(r) = 1,for all r E A. Since (0 + ())(r) = ((r))((r)) = 1, O(r) is the multi plicative inverse of O(r) E U(n). Hence, O(r) = 0(r)1 and 0(r)1 = O(r), for all r E A, since .O(r) E U, for all r E A, which proves (3). 0 Note that E A implies : A + U(n) C C, so that 4 L(A), for all ) E A. Thus, A C L(A). The set L(A) = {a: A + C} of all complexvalued functions on A is a vector space over C with addition and scalar multiplication given by (a+ b)(0) =a(0) +b(O), for all a,bE L(A),4 E A, ()(O) = Ca(), for all a E L(A), E C, E A. For a,b E L(A), define a. b E L (A) by (. b)(0) = a()b()), for all E A. Define an inner product and a norm on L(A) as follows. For a,b 6 E A, the inner product of a and b is given by (a,b) = a(4)b4), and the norm of a is given by 1/2 a=^2)= ja(0)12 In Section 2.3 we define a basis of L(A) with cardinality n. Denote by L2(A) the inner product space L(A) together with the norm I1 . Thus, L2(A) is a normed vector space over C of dimension n. In order to define the Fourier transform for finite abelian groups in Section 2.3 we need a bilinear form (.,.): A x A * U(n). Proposition 4 If A is a finite abelian group of cardinality n, then the map (.):A x A * U(n) defined by (r,) =(r), for all r E A and 0 E A, is a bilinear form, that is, 1. (r + s, b) = (r, 0)(s, ), for all r, s E A, E A; 2. (r, + 0) = (r, 4)(s, 4), for all r E A, EA; 3. (kr, ) = (r,))k = (r,ko), for all r E A,4) A. In particular, Proposition 4 implies (r,4) = (r,)' = (r, 4) = (r,), for all r E A, 0 E A. 2.3 The Discrete Fourier Transform and Inverse DFT Definition 1 (Discrete Fourier Transform) If A is a finite abelian group with card(A) = n, then the Fourier transform of a E L(A) is the function a E L(A) defined by (4) = a(r)(r), (2.1) rEA for all ) E A. We shall write TA for the map TA : L(A) + L(A) such that A (a) = i. The complex numbers a(O) are known as the Fourier coefficients of a E L(A). Computation of the DFT of a E L(A) via equation 2.1 requires no more than n2 complex multiplications (not counting multiplication by the normalization constant). The factors (r, 0) can be kept in a lookup table. Hence, their computation does not contribute to the complexity of equation 2.1. Example. For a E L(A) and s E A, the shift of a by s is the function a. E L(A) defined by a.(r) = a(r s), for all r E A. Applying the Fourier transform to a. yields a,(0) = (s, 4) a( () = 0(s) fa(O), for all q E A. One may easily verify linearity of the Fourier transform, which we state for mally. Proposition 5 If A is a finite abelian group with card(A) = n, then FA is a linear transform, that is, .A(a + b) = YA(a) + YA(b), for all a,b E L(A), and YA(a) = YA(a), for all a E L(A),C E C. For an application of the Fourier transform, we determine a basis of L2(A) by computing the Fourier transforms of the basis elements in E = {er : r E A}. Let er = A(er), for all r E A. Lemma 6 If r E A, then 1 (or al ( A.r) for all 0 E A. Proof. If r E A then = =Eer(s)(S,( ) sEA 1 = 1v for all E A. Noteworthy is the function eo, where 0 is the zero element of A. Since eo() = 0(0) = 1 for all 0 E A, o is constant on A. Lemma 7 If r E A, then Proof. If r = 0, then (r, from the cardinality of A. is, (r, 4) # 0. Thus, Sn, ifr=0 (r ) = 0, otherwise. 0) = (0, 0) = 0(0) = 1, for all 4 E A. The result follows If r # 0, then there exists 0 E A such that O(r) # 1, that ,EA OEA OEA which implies D(,) =0. 'kEA Proposition 8 The set 6 = {fr : r E A} is an orthonormal basis of L2(A). Proof. Suppose r, s E A. on A x A, we have Using the inner product in L2(A) and the bilinear form (Fr,es) = Fr(;b)?s() EA 1 E 4>(r)4O(s) (by Lemma 6) n  4'EA E !(rs0) (s n  _ EA n OA n O6A 1, ifr=s S 0, otherwise, the last equality by Lemma 7. Thus, (e,, ',) = Sr,, for all r, s E A, which implies & is an orthonormal basis of L2(A). 0 Since the cardinality of E is n, L(A) has dimension n. Furthermore, since TA maps an orthonormal basis to an orthonormal basis, we have the following result. Theorem 9 The Fourier transform FA : L2(A) + L2(A) is an isometry, that is, IIYA(a)II = \all, for all a E L2(A). Before defining the inverse Fourier transform, we present a dual to Lemma 7. Lemma 10 If q5 E A, then D r ) n, if b= 0 EA' f 0= otherwise, where O is the trivial homomorphism. 15 Proof. If 4 = o, then (r, ) = 0o(r) = 1, for all r E A. The result follows from the cardinality of A. If 0 # 0o, then O(s) 5 1, for some s E A, and rEA rEA rEA = (s,4))Z(r,4) reA implies D(r,4) = 0. rEA 0 The next example yields another basis of L2(A). Example. If 0, l E A, then 0, E L(A). We compute the inner product of and 4': ( 4 = S (r)O(r) rEA = D(r, 0)(r, V) rEA = 5(r,4)(r,4') rEA = 5(r,44) rEA n, if= = 0, otherwise, the last equality by Lemma 10. Normalization yields an orthonormal basis e = A =  fVnO 0 v}= n of L 2(A). Example. If 4 E A, then 4 E L(A). What are the Fourier coefficients of FA= = (O)? If E C A, then =E ) V rEA 1 rEA 1+ VnrEA Vfn, ifb= = 0, otherwise, the last equality by Lemma 10. Definition 2 (Inverse Discrete Fourier Transform) If A is a finite abelian group with card(A) = n, A is the algebraic dual of A, and a E L(A), then the inverse Fourier transform of a is the function a = )1(ad) E L(A) given by 1 a(r) = 4)(r,), (2.2) OEA for all r E A. Computation of the inverse DFT of a E L(A) via equation 2.2 requires at most n2 complex multiplications. The composition of FA and F;'1 is the identity transformation on L(A), as stated in the next result. Theorem 11 (Inversion Theorem) If a E L(A) and a = FA(a), then a = F'(a). Proof. Let a E L(A) and b = J^(ia). If s E A, then b(s) = ()s,) sEA V (rEA OEA r A OEA tEA = a(s) n (by Lemma 7) G~A rEA Z:a (r) E(r s n rEA ( a(s).n (by Lemma?7) = a(s), which implies b = a. E 2.4 The Convolution Theorem Definition 3 (Convolution) If a, b E L(A), then the (cyclic) convolution c E L(A) of a and b, written c = a b, is defined by c(t) = a(r)b(t r), rEA for all t E A. Example. We compute the convolution of a basis element e, E C with a function a E L(A) at t E A: (e, a)(t) = E e,(r)a(t r) = a(t s) = a,(t), rEA where a, is the shift of a by s A. Thus, e, a = a,. Example. If , 0 E A, then , 0 E L(A) and (0 )(0) = 0(r)7(0 r) rEA rEA rEA rEAi f n, if 0=0 0, otherwise, the last equality by Lemma 10. The next example shows that we may regard the Fourier transform as a con volution. Example. If a E L(A) and 4 E A, then the properties of the bilinear form on A x A imply (a (0))(0) = a(r)(0(0 r)) rEA E a(r)(r,*) rGA rEA = /n. (<). Thus, a40) = (a (k))(0). The Fourier transform relates convolution and pointwise multiplication. This relation is known as the convolution theorem. Theorem 12 (Convolution Theorem) If A is an abelian group with card(A) =n, a, b E L2(A) and c = a* b, then F= V/nii a. b. Proof. If a, b, c E L(A), with c = a b, and E A, then (a.6)() = b(1.( ( 1 a(r) (r ) ( Zbss)( i/n rA ) ( / E rEA SsEAA Za~r) (Zb(s)(sc~) (r,4, n rEA sEA n rEA sEA = 1TZ a(r)b(s)(s + r,,), n rEA sEA using bilinearity. Letting t = s + r and interchanging the summations, (a. b)() = a(r)b(t r)(t, 0) rEA tEA = (ca(r)b(t,r)) (t,) n tEA tEA 1  v"i V ,/ CIO). 0 Example. Since a, = e, a from a previous example, the Fourier transform of a, is a, = v/n a, by the convolution theorem. 2.5 The CooleyTukey Factorization The CooleyTukey (CT) factorization is an algorithm that computes the dis crete Fourier transform of a function a E L(A), where A is a finite abelian group. An immediate corollary computes the inverse DFT of a function a E L(A). For the spatial domain, consider the short exact sequence 0  B  A  AIB  U, where A is a finite abelian group and B is a subgroup of A. Note that B does not necessarily split in A, that is, it is not generally true that there exists a subgroup C of A such that A = B E C (internal direct sum). However, if R C A is a complete and nonredundant set of coset representatives of the quotient group A/B such that B n R = {0}, then A = B + R, and for every r E A there exist unique r3 E B and r, E R such that r = r# + r. We refer to such a set as a complete system of coset representatives (or residues) of A/B. Let B. = { E A : (r,>) = 1, for all r E B}. It is easy to see that B. is a subgroup of A. In the frequency domain consider the short exact sequence Oo  B.  A * A/B.  o where Oo is the trivial homomorphism in A. If R. C A is a complete system of coset representatives of A/B., then A = B. + R. and for every 0 E A there exist unique E5 C B. and 0# E R. such that 4 = (f + 0,. One can show that A/B. B and B. A/B. Let n = card(A), np = card(B), and n, = card(A/B). We have n = nnp = card(A), card(R) = n, = card(B.), and card(A/B.) = np = card(R.). An et al. [1], give a detailed discussion of the next result, although their approach depends on an isomorphism from A onto A. Theorem 13 (CooleyTukey Factorization) Let A be a finite abelian group, B a subgroup of A, and R C A a complete system of coset representatives of A/B. If a E L(A) and a = FA(a), then a()= ,),a(r)(r# O) *(r,,)) (r ) (2.3) 'rBER for all 4 = 0. + # c B. + R. = A, where R. C A is a complete system of coset representatives of A/B., and r = rp + r. E B + R = A. Proof. Let 4> = 0,, + o E A with 4O E B. and O E R.. If rp E B, then (r, 4o) = 1 by the definition of B.. Hence, (r,4) = (r +r,. + Op) = (r/3, o)(r/3, )/)(ro,(fa)(ro,(f>0) = (r/5, )(ro, )(ro, a), for all r = rp + ra E A and 4 = + Op E A. The Fourier transform a() = a(r)(r,'), (2.4) V rEA where n = nnn is the cardinality of A, can therefore be written as in the statement of the theorem (equation 2.3). 0 Corollary 14 Let A be a finite abelian group, B a subgroup of A, and R. C A a complete system of coset representatives of A/B.. If a E L(A) and a = ,'J(a), then a^ 1 =O 1 f, ) ).) v OPAvER 7 76 EB for all r = ro + r. E B + R = A, where R C A is a complete system of coset representatives of A/B, and 4 = 0, + 4p E B. + R. = A. The numbers (r,, 0) and (r,, ) in Theorem 13 and Corollary 14, respec tively, are known as twiddle factors. The CooleyTukey factorization is a threestage algorithm. The first stage is the inner summation in equation 2.3 which requires no more than npn,, = n complex multiplications for each Op, for a total of n np complex multiplications. The second stage is an entrybyentry multiplication by the twiddle factors consisting of na"n complex multiplications. The last stage is the outer sum mation in equation 2.3 which requires no more than n,"np complex multiplications for each 0,, for a total of n na complex multiplications. The total number of complex multiplications is bounded by n n + n n, a + n 6n. Consider the computational savings for the case n = 104, n = 100 = n,. Direct computation using equation 2.4 involves n2 = 108 complex multiplications. The CooleyTukey factorization algorithm requires no more that 2 106 + 104 complex multiplications. CHAPTER 3 QUOTIENT LATTICES 3.1 Introduction In this chapter we consider the discrete Fourier transform for functions defined on finite abelian groups that arise as quotients of Ddimensional lattices. Instead of working with algebraic duals in the frequency domain, we will take a more intuitive approach and work with the equivalent geometric duals. This approach differs from standard presentations of the discrete Fourier transform that emphasize the algebraic structure of the underlying indexing sets for function domains; our emphasis is on the geometric structure of the underlying indexing sets or coset representatives. Section 3.2 presents the discrete Fourier transform and the inverse DFT for quotient lattices. Since quotient lattices have the structure of finite abelian groups, all the results of Chapter 2 are immediately applicable. In particular, the inversion and convolution theorems apply to the quotients under consideration in this chapter. Section 3.3 states the CooleyTukey factorization for quotient lattices in the spatial and frequency domains. Denote the field of real numbers by R, Ddimensional Euclidean space by RD, and the set of integers by Z. Denote the inner product of x, y E RD by (x, y). 3.2 The DFT and Inverse DFT for Quotient Lattices A Ddimensional lattice is the set of all linear integer combinations of D linearly independent vectors in RD. Let L be a Ddimensional lattice. The geometric dual (or reciprocal lattice) L* of L is defined by L = {s E RD: (r,s) E Z, for all r E L}. Note that L* is also a Ddimensional lattice. Denote a sublattice relationship by L1 < Lo for two Ddimensional lattices L0 and L1. Note that "<" is also an abelian subgroup relationship. Proposition 15 (Lattice Duality Relationship) If Lo and L1 are Ddimensional latttices, then L, < Lo ==> L* < L\ Remark 1 From a signal processing perspective, we regard Lo in Proposition 15 as a spatial domain sampling lattice and L, as a spatial domain periodicity lattice. Similarly, we regard L\ as a frequency domain sampling lattice and LQ as a frequency domain periodicity lattice. The sublattice relationship also applies to sublattices of quotient lattices of the form Lo/Li. Quotient sublattices arise from sequences of the form L, < L2 < L0. Corollary 16 If Lo, L1, and L2 are Ddimensional latttices such that L, < L2 < Lo, then 1. L* < L < L\; 2. L2/L1 < Lo/L, and L2/L; < L*/Lo; 3. (Lo/L1)/(L2/Lj) Lo/L2 and (L/L*)/(L*/L ) LI/L;. Define quotient lattices for the spatial and frequency domains by A = Lo/LI and A' = L /L respectively. If n is the index of L1 in Lo, then A is a finite abelian group with cardinality n, and nr E L1,for all r E L0. By Proposition 18 (below), the index of L; in L\ is also n. Hence, ns E L, for all s E L\. The spatial and frequency domain short exact sequences corresponding to A and A' are 0 + L,  Lo  A 0 and 0 *+ L L A' 0, respectively. Recall that U(n) denotes the multiplicative group of all the nth roots of unity. The next result defines a bilinear form on A x A'. Proposition 17 Let Lo and L1 be Ddimensional lattices with L, < Lo. If A = Lo/LI1, A' = L*/Lo, and n = card(A), then the map (, A x A' + U(n) defined by (:) = exp[27ri(r,s)], where i = V/2 and F and represent the costs containing r and s, respectively, is a bilinear form. Proof. Since A = Lo/L1 and nr E L1, for all r E Lo, we have (nr, s) E Z, for all s L*, that is, (r,s) E IZ, for all r E L0 and s E L which implies (F,) = (r + L1,s + L ) = exp[27ri(r,s)] E U(n). To show that (r,s) is independent of the choice of coset representatives, let rl, r2 E L0 with ri = r2 mod L1. Then r, r2 = u, for some u E L1. Now, if s E L*, then exp[27ri(ri,s)] = exp[27ri(r2, s)] 4== exp[27ri(r1 r2,s)] = 1 4= exp[27ri(u,s)] = 1, which is true for all s E L\ since u E L1. A similar argument shows that (F,3) is independent of the choice of coset representative of the equivalence class X. Thus, (, .) is well defined. The identities (OF + 'F2, ) = (TI,5)(F2, ), for all F ,T2 E A, j E A', (F, 3 +s2) = (F, 1)(F,52), for all T E A, 31,32 E A', and (cj) = (r,3)c = ( ,cs), for all c E Z, r E A, s E A' are all easily verified. ] Note that (r,s) = 1 if and only if (r, s) E Z, where r is any coset representative of F E A and s is any coset representative of 5 E A'. The quotient lattice A' is the geometric analogue of the dual group A = Hom(A, U(n)). Two abelian group isomorphisms follow. Proposition 18 If Lo and L, are Ddimensional lattices such that L, < Lo, A = Lo/L1, A' = L*/L;, and n = card(A) (as in Proposition 17), then (1) A' A, and (2) A' ^ A. Proof. Part (1). For every = s + L; E A', define a map < : A + U(n), where n = card(A), by 0,(F) = (F, ), for all F = r + Li E A. The linearity of (,3) implies A,, that is, 4> is a homomorphism. Define 'P : A' 4 A by 'P(3) = ^, for all E A'. The linearity of (F, .) implies %F is a homomorphism. Since card(A) = n = card(A) (Theorem 2), we need only show T is injective. If i,S2 E A', then T(1) = 'I(52) == (T,;i) = (f,S2), for all F E A ==# (F(,31 52) = 1 === exp[27ri(r,si S2)] = 1, for all r E Lo, = S1 = 52 == 5I =S2 Hence, T is an isomorphism. Part (2). By Theorem 2, A A; by the first part of the proposition, A' A. Hence, A' A. 0 Denote by L(A) the ndimensional space of all complexvalued functions on A and by L(A') the ndimensional space of all complexvalued functions on A'. The geometric analogue of the Fourier transform (Definition 1) of a function a E L(A) is the function a' = .A(a), where FA : L(A) + L(A'), defined by a,(P) = = a(f)(T,) = = a(f)exp[27ri(r,s)], for all 5 6 A'. (3.1) v EA V 7EA The inverse Fourier transform of a function a' E L(A') (the geometric analogue of Definition 2) is the function a = Jc'(a') E L(A) defined by a(T) = a'(3)(T, ) = n a'(3) exp[27ri(r,s)], for all F E A. (3.2) E A' 3EA' In the sequel we will work with sets of coset representatives of A whenever A is a quotient lattice, that is, whenever A = Lo/L1 for Ddimensional lattices L0 and LI. Similarly for the dual quotient A' = LI/L In this context, let G C L0 be a complete system of residues for A and let G' C L\ be a complete system of residues for A'. Assume 0 E G and 0 E G', where 0 is the zero vector in RD. Define a map (, ) on Lo x L; by (r,s) = exp[27ri(r,s)], for all r E Lo, s E L*. (3.3) Since Lo and L1 are Zmodules, and from the proof of Proposition 17, we have (., .): Lo x L* U(n), where n = card(A), is a bilinear form. Now, define a: G 4 C by a(r) = a(f), for all r E G, and a': G' + C by a'(s) = a'(3), for all s E G'. Note that G x G' C Lo x L,. In view of Proposition 17, compute the DFT of a E L(A) as 1 1 a'(s) = E a(r)(r, s) = E a(r)exp[2ri(r, s)], for all s E G', (3.4) SrEG V rEG and the inverse DFT of a' E L(A') as a(r) = 1: > a'(s)(r,s) = = j a'(s)exp[27ri(r,s)], for all r E G. (3.5) BEG' sEG' 3.3 The CoolevTukey Factorization for Quotient Lattices The purpose of this section is to state the CooleyTukey (CT) factorization (Theorem 13) and the inverse CT factorization (Corollary 14) for the pair of sequences { B + A + A/B, B' + A' * A'/B'}, respectively, where A and A' arise as quotients of Ddimensional lattices, B is a subgroup of A, and B' is a subgroup of A'. If n = card(B) and n, = card(A/B), then n = nono = card(A). In the spatial domain, let L0, L1, and L2 be Ddimensional lattices such that L, < L2 < L0. Consider the short exact sequence 0 * B + A + A/B 0, where A = Lo/L1 and B = L2/Li < A. If G C Lo is a complete system of residues for A, Kp C L2 a complete system of residues for B, and Gc, C Lo such that G = KO+G, and KfnGl = {0}, then for every r E G there exist unique rj3 E KG and ro E G, such that r = ro + r,. Furthermore, card(G) = n, card(Kp) = n6, and card(Ga) = n,. Example. Let L0 = Z, L2 = 10Z, and L1 = 1000Z. Then A = Z/1000Z and B = 10Z/1000Z. The set G = {0, 1,..., 999} C Z is a complete system of residues for A. Furthermore, if K# = {0,10,20,... ,990} C 10Z and G, = {0,1,...,9}, then K# is a complete system of residues for B, Kp n G, = {0} and G = K# + Ga,. In the frequency domain, we have Lo < L2 < LI (Corollary 16). Consider the short exact sequence 0 + B'  A' + A'/B'  , where A' = LI/L* and B' = L'/L < A'. Note that B' # L/L*. (To understand the relationship between B and B', let F = r + L1 E B = L2/L1 and = s + L; E B' = L2/L;. Since r E L2 and s E L, we have (r,s) E Z. Hence, (r,s) = exp[27ri(r,s)] = 1. Thus, T E B and E B' implies (r,s) = 1.) If G' C LI is a complete system of residues for A', K" C L a complete system of residues for B', and G9 C L, such that K"f n G6 = {0} and G' = K" + G', then for every s C G' there exist unique sa E K' and sp E GO such that s = sa + sp. Corollary 16 and Proposition 18 imply B' A/B and A'/B' B. Since A' A, we have card(G') = n,card(K') = n, and card(GC) = n = card(A'/B'). Example. If Lo = Z, L2 = 10Z, and L1 = 1000Z, then L; = Z,L = 0.1Z, and L\ = 0.001Z. Furthermore, A' = 0.001Z/Z and B' = 0.1Z/Z. The set G' = {0,0.001,0.002,... ,0.999} C 0.001Z is a complete system of residues for A' and K' = {0,0.1,0.2,... ,0.9} C 0.1Z is a complete system of residues for B'. If G# = {0,0.001,0.002,..., 0.099} C 0.001Z, then K", Gn = {0} and G' = K" + G'p. If a E L(A) and a' = 17A(a) E L(A'), then the spatial domain CT factoriza tion for quotient lattices is determined by a'() = ( a a(r)(r, sf) (r.,so)) (r.,s.), (3.6) raEG .( rpEKe for all s = a+s + s E K + G = G', where r = r# +r, E K +G,, = G. If a' E L(A') and a = .F'(a') E L(A), then the inverse or frequency domain CT factorization for quotient lattices is determined by a(r) l 1 E a'(s)(r, sa)) (r.,so)) (rs), (3.7) for all r = ro + r,, E Kc + G,, = G, where s = so, + so E K + G' = G'. CHAPTER 4 PERMUTOHEDRAL AGGREGATES 4.1 Introduction The purpose of this chapter is to present a brief introduction to certain hier archical arrangements of points (or cells) known as permutohedral aggregates. These aggregates and their geometric duals are the objects of interest in this dissertation. This chapter applies the material from Chapter 3 to quotient lattices that are as sociated with permutohedral aggregates. Specifically, we describe the CooleyTukey factorization algorithms in the spatial and frequency domains, respectively, for these quotients. Section 4.2 gives a geometric description of permutohedral aggregates. In dimensions two and three they are hexagonal and truncated octahedral aggregates, respectively. See Conway and Sloane [4] for the necessary background. Section 4.3 discusses quotient lattices of the form A\ = L/N"L, where L is a Ddimensional lattice that has a ring structure, N is an endomorphism of L, A is a positive integer, and N'L = NA(L). A Alevel permutohedral aggregate is a particular set of coset representatives of the quotient AA. These representatives are discussed in greater detail in Section 4.4. Section 4.5 contains statements of the CooleyTukey (CT) factorization and inverse CT factorization in terms of coset representatives of A\ and the geometric duals A'. Although the quotient lattices A,\ have a ring structure, the factorization algorithms discussed here do not depend on the multiplicative structure of these quotients. 4.2 Geometric Description Permutohedral aggregates can be defined constructively starting with a poly tope, known as a permutohedron, associated with the classical Ddimensional root lattice LD = {(X0,X1,...,XD) E ZD+l1 :X o + *" +zrD = 0}, which lies in the hyperplane E xi = 0 in RD+l. A Ddimensional permutohedron is the Voronoi cell of the dual lattice LB; its vertices are the (D + 1)! points obtained by permuting the coordinates of (G(D), (D + 2), (D + 4),..., (D 2), tD). Permutohedra tile Ddimensional Euclidean space. In dimension one, L, and L\ can be represented by Z in R. The Voronoi cells are integer translates of the interval [1/2,1/2]. A level 1 aggregate consists of the points G = {1,0,1}. For A E Z+, a Alevel aggregate is the set of 3x integers G\ A ^=0 3kG. We also refer to the set of Voronoi cells associated with GA as a Alevel aggregate. In dimension two, L2 and L2 are hexagonal lattices. The hexagonal lattice is an optimal packing strategy for disks in the plane. The Voronoi cells of a hexago nal lattice are regular hexagons. A level 1 hexagonal aggregate consists of a central hexagon together with six neighboring hexagonal tiles as shown in Figure 4.1. A level 2 hexagonal aggregate consists of seven neighboring level 1 aggregates. In this manner, higherlevel aggregates are constructed from lowerlevel aggregates. The hi erarchical structure of a level 3 hexagonal aggregate consisting of 73 = 343 hexagonal cells is shown in Figure 4.2. In dimension three, L3 is a facecentered cubic (fcc) lattice and LZ is a body centered cubic (bcc) lattice. The fcc lattice is an optimal packing strategy for spheres Figure 4.1: Level 1 hexagonal aggregate. Figure 4.2: Hierarchical structure of a level 3 hexagonal aggregate. of equal radius in threespace. The Voronoi cells of an fcc lattice are rhombic do decahedra (Figure 4.3). The Voronoi cells of a bcc lattice are truncated octahedra Figure 4.3: Rhombic dodecahedron. (Figure 4.4). The bcc lattice can also be defined as the set of vectors (x,y,z) E R3 where the coordinates x, y, and z are all even or all odd integers. Figure 4.5 shows a level 1 truncated octahedral aggregate, which consists of a central truncated octahe dron and 14 neighboring cells. Incidentally, beehive chambers are hexagonal prisms capped by rhombic dodecahedral ends (see Kaproff [12]). In the next section we give an alternate description of permutohedral aggre gates. 4.3 Quotient Lattices Associated with Permutohedral Aggregates Let p represent the polynomial p(x) = 1 + x + ... + xD, where D E Z+, and let A = Z[x]/(p), a commutative ring with identity T = 1 + (p). Let = = Figure 4.4: Truncated octahedron. x + (p) and let r = 2 C. The quotient ring A is a free abelian group of rank D with ordered basis X = (1, ,... ,i"D), where 1 is the identity element. We have: D = (1 + + ... + CD1) and CD+1 = 1. The quotient ring A can be realized as a Ddimensional lattice in RD by mapping the basis elements to any set of D linearly independent vectors in RD. We will simply identify the basis elements in X with a chosen set of D linearly independent vectors in RD and regard the elements of X as labels or names for the assigned vectors. The element r determines an endomorphism N, on A defined by N, : p + rp, for all p E A. For example, if D = 2, then N,r(C) = r7 = (2 ) = 2_ C2 = 2C + 1 + = 1 + 3. Let NA = Nr(A). Note that N7A = rA = (r) is a principal ideal of A. Furthermore, rX = (r, T, ..., rTDl) is a Figure 4.5: Level 1 truncated octahedral aggregate. basis of rA. With respect to the basis X, S2 00 ... 01 1 2 0 ... 0 1 0 1 2 ... 0 1 Nr= : : (4.1) 0 0 0 .** 21 0 0 0 ... 1 3 Furthermore, det N, = 2D+1 1 (see Vince [22]). More generally, for A E Z+, NAA = NA(A) = TrA = (r7), is an ideal of A, and det N,' = n. Kitto, Vince, and Wilson [13] prove the following isomorphism of rings. Theorem 19 If D, A E Z+, then A/NAA A Z,, (ring isomorphism), where n = 2D+1 1. In particular, A/NAA is cyclic and card(A/NAA) = nA. Let {uj : j = 0,1,..., D 1} be a set of D linearly independent vectors in RD and denote by U the ordered basis (uo, U1,... u D1) of RD. Identify Ci with uj, for j = 0,1,..., D 1, and denote by L the Ddimensional lattice generated by U. Let N be the endomorphism of L such that the matrix of N with respect to the basis U is given by the matrix in equation 4.1. Thus, NL = N(L) is a Ddimensional sublattice of L and the index of NL in L is n = det N = 2D+1 1. More generally, if A E Z+ and NAL = NA(L), then NAL is a Ddimensional sublattice of L that is generated by the set of vectors {NA(uj) : j = 0,1,..., D 1}. Furthermore, the index of N'L in L is det NA = (det N) = n. We turn our attention to quotient lattices derived from L and the sublattices NAL. Let A = L/NL, that is, A is the cokernel of N E End(L), the set of endo morphisms of L. We have card(A) = n = 2D+1 1. More generally, for A E Z+, let A. = L/N'L, that is, A, is the cokernel of N' C End(L). We have card(A,) = nA. Remark 2 A Alevel permutohedral aggregate (in the spatial domain) is a particular set of coset representatives of A,\ as a quotient lattice. For permutohedral aggregates, the uj E U are chosen such that the set of D+ 1 points {u, :j = 0,1,...,D}, where UD = ]jDo1 u3, is the set of vertices of a regular Dsimplex with barycenter at the origin. See Lucas [14] for a particular choice of coordinates for these vectors. By the lattice duality relationship (Proposition 15), NL < L implies L* < (NL)*. Hence, the geometric counterpart of A in the frequency domain is the quotient lattice A' = (NL)*/L*. Similarly, for A E Z+, the geometric counterpart of AA in the frequency domain is A\ = (NAL)*/L*. We can regard the (injective) endomorphisms NA as Linvariant, invertible, linear operators on RD (and therefore bounded). Recall that every bounded linear operator on a Hilbert space has a unique adjoint. Theorem 20 If H is a Hilbert space and M is a bounded linear operator on H, then there exists a unique linear operator M*, known as the adjoint of M, satisfying (M(x),y) = (x, M*(y)), for all x,y E H. For the elementary results from linear operator theory in this section and for a proof of Theorem 20, which uses the Riesz representation theorem for a Hilbert space, see DeVito [6]. We shall use the following results for a linear operator M: if M is an invertible linear operator and M* is the adjoint of M, then M* is invertible and (M*)1 = (M1)*; if M is given by a realvalued square matrix with respect to a particular basis, then M* = Mt, the transpose matrix with respect to the same basis; if A E Z+, then (MA)* = (M*)A; finally, det M* = det M. The next proposition will allow us to express A' and A' in terms of N*, the unique adjoint of the Linvariant, linear operator N. In general, for I a Ddimensional lattice and M an Iinvariant linear operator on RD, let MI = M(I). Note that if M is invertible, then MI is a Ddimensional lattice and MI < I (sublattice and subgroup relationship). Proposition 21 If I is a Ddimensional lattice, I* is the dual lattice of I, M is an Iinvariant, invertible linear operator on RD, and M* is the unique adjoint of M, then (1) (MI)* = (M*)1I*, and (2) if I' = (MI)*, then M" is I'invariant. Proof. Part (1). We first show (MI)* C (M*)1I*. Let u E (MI). Then (MA(r),u) = (r, M*(u)) E Z, for all r e I. Hence, M*(u) E I*, that is, there exists s E I* such that M*(u) = s, or u = (M")(s). Thus, u E (M*)1I* and (MI)* C (M*)'I*. To show (M*)1I* C (MI)*, let s E I*. We must show (M*)'(s) E (MI)*, that is, (M(r), (AM)'(s)) E Z, for all r E I. Let r E I. Then M(r) E MI and (A4(r),(M*A)'(s)) = (M(r),(M)')(s)) = (M'(MA(r)),s) = (r, s) E Z, since r E I and s E 1*. Hence, (M )'1(s) E (MAI) and (M*)'I* C (MI)*. Part (2). The hypotheses imply MI < I. Let I' = (MI)* = (M,)I*, using the first part of the proposition. Then I* = M*I'. Now, MI < I implies I* < (MI)* (Proposition 15), that is, M*I' < I'. Hence, MA E End(F'). 0 For D, A E Z+ fixed, let L' = (NAL)*. Part (1) of the proposition implies L'\ = ((NA')*)L* = (N*)L*; part (2) implies (NA)* = (N*)A is L'invariant. Hence, L* = (N*)AL' and A' = L/1(N))'L'. Just as AA is the cokernel of N' E End(L), A' is the cokernel of (N*)' E End(LZ\). The corresponding short exact sequencesfor A\ and A', respectivelyare 0 NL 4 L +L/NL 0, in the spatial domain, and 0 4 (N)'L' + L L/(N)L', 40, in the frequency domain. Remark 3 Note that N* is L')invariant. To see this, consider the chain of lattices NAL < NIL < ... < NL < L. In particular, N'L < NL implies (NIL)" .< (NAL)*, or (N*)'L* < (N*)AL", by the first part of Proposition 21. Hence, if L' = (N*)AL*, then N'L'\ < L". Denote by L(A\) the nAdimensional space of complexvalued functions on A\, and denote by L(A') the n'dimensional space of complexvalued functions on A\. Furthermore, denote the discrete Fourier transform (equation 3.1) of a function a E L(AA) by a' = .Fx(a) E L(A'), and denote the inverse DFT (equation 3.2) of a function a' E L(A') by a = .T'(a') E L(AA,). 4.4 Aggregates as Sets of Coset Representatives Let I be a Ddimensional lattice in RD and M E End(I) such that I det M\ > 1. A set R C I is called a residue system for M if 0 E R and I = U{r + MI : r E R} (see Bandt [3]). Thus, R is a set of coset representatives of I/MI and card(R) = det M. A finite address of a lattice point r E I is a finite sequence rAl "'riro (A E Z+) such that rk E R, for all k, and r = E\= M k(rk). Vince [22] shows that the finite address of a lattice point, if it exists, is unique. A Alevel aggregate (not necessarily a permutohedral aggregate) is the set of lattice points having addresses of length no greater than A. We also refer to the corresponding set of Voronoi cells as a Alevel aggregate. The following lemma states that Alevel aggregates are particular sets of coset representatives of I/MAI. Lemma 22 Let D, A E Z+. If I is a Ddimensional lattice, M E End(I), I det MI > 1, and R is a residue system for M, then AI R, = Mk R) k=O is a set of coset representatives of IIM/I. Proof. Since R is a set of coset representatives of I/MI and MI fl R = {0}, we have I = MI + R. Now, I = .MI+R = M(M/I+R))+R = M21 +M (R) + R = MAI + MAl(R)+... +M(R) + R = M\,I+Rx. Hence, R\ is a set of coset representatives of IIMAI. 0 For D E Z+, let n = 2D+I 1, p(x) = 1 + x +... + xD, A = Z[x]/(p), and 7 = 2 , where = x + (p), as before. Kitto, Vince, and Wilson [13] prove the next result (we omit the use of overbars for elements of A/NAA; Nr is the endomorphism of A defined in Section 4.3). Lemma 23 Each of the following subsets of A is a set of coset representatives of the quotient ring A/NTA: 1. R={0,1,...,n1}; 2. S = {bDD +. + bj. + bo: bj E {O,1}, forj = O,1,...,D, not all bI = }. Note that elements in the second set in Lemma 23 can be written as strings (of length D+1) of binary digits, or bits. Thus, we can represent s = bDCD+.. +b1+bo0 E S by the string (bD ... b1b0). Conversely, if (bD .. b1bo) is the binary representation of r E R, then r = s mod NA, where s = bDD +... + bj + bo = (bD... bbio) E S. Example. If D = 1, then A = Z[x]/(1 + x) is isomorphic to Z, R = {0,1, 2} and S = {0, 1, }. Note that 2 = mod NAA. Example. If D = 2, then A = Z[x]/(1 + x + x2), R = {0,1,2,3,4,5,6} and S = {0,1, + 1, 2, 2 + }. The following equivalences (modulo NIA) hold between the elements of R and the elements of S: 2 = , 3 + +1, 4 = C2, 5 = 2+1, 6 +. Note that the binary representation of 6 E R is (110), which represents the equiva lence class 1 2 + 1 + 0. 1. Similarly for the other elements of R. We now consider sets of coset representatives of the quotients A, and A'. Let U = (uo, u1i,..., UD1) be an ordered set of D linearly independent vectors in RD and let UD = ,1D_0 uj. Under the identifications ++ uj, for j = 0,1,..., D, the set G D G= bju : bj EO, 1}, for = 0, L., D, not all b = 1} (4.2) j=O I is a residue system for N : L + L, where L is the Ddimensional lattice generated by U and N is the endomorphism of L specified by the matrix in equation 4.1 with respect to the basis U. The set G is equivalent to the second set in Lemma 23 under the identifications j ++* uj. Recall that det N = 2D+1 1. By Lemma 22, the set G, = Z Nk(G) k=O is a set of coset representatives of A\ = L/NAL. In the frequency domain, let 1, = (ul, u,... i u_) be an ordered set of D linearly independent vectors in RD that generate the Ddimensional lattice L' = (NAL)* = (N*)AL, and let u' = E. 'o u such that S1} Ga = bju1: bj E {0,1}, for j = 01,..., D, not all b6 = 1 (4.3) Ij=o I is a residue system for N* : L\ +* L'. Note that det N* = det N. By Lemma 22, the set AI G' = (N*)k(G') k=O is a set of coset representatives of A\ = L\1(N*)L'\. Example. Let D = 2 and let L be the 2dimensional lattice in the spatial domain generated byU = (uo, ul) and let u2 = (u0 +ul). We have G = {b2u2 +bu i+bouo: bj {0,1}, not all bj = 1}, N=( 1 3 with respect to the ordered basis U, and det N = 7 = card(G). If A = 2, then L' is the lattice (N2L)* = (N*)2L* in the frequency domain. Suppose L'2 is generated by U' = (uo,U ) and let u' = (u' + u'). Then G' = {b2z + blu + bo0uo : bj E {0,1}, not all bj = 1} is a set of coset representatives of A'2 = L'/(N*)2L'. If uo = (1,0)t and ul = (l/2,v/3/2)t, then L is a hexagonal lattice in the plane (congruent to the lattice L2 in Section 4.2); the associated Voronoi cells are hexagons. Figure 4.6 depicts a level 1 hexagonal aggregate (the set G) in the spatial domain. The figure shows the seven lattice points and the corresponding hexagonal tiles. Note that the ordered pair (uO, ul) defines a righthanded coordinate system in the plane. Level two hexagonal aggregates in the spatial and frequency domains are shown in Figures 4.7 and 4.8, respectively. Geometrically, the two domains differ in scale, orientation, and chirality (or handedness). Figure 4.6: Level 1 hexagonal aggregate in the spatial domain. Figure 4.7: Level 2 hexagonal aggregate in the spatial domain. 0.6 Figure 4.8: Level 2 hexagonal aggregate in the frequency domain. Example. If D = 3, then ( 2 0 1) N= 1 2 0 0 1 3 and detN = 15 = card(G). Choose u0 = (1,0,0)t,ui = (1/3,2v=/3,O),u = (1/3,v/2/3,2/v/6)t. These vectors generate a bodycentered cubic lattice (con gruent to the lattice L; from Section 4.2). The associated Voronoi cells are truncated octahedra. A level 1 truncated octahedral aggregate is shown in Figure 4.5. 4.5 The CooleyTukey Factorization for Permutohedral Aggregates Fix D, A E Z+ and let n = 2D+1 1. In the spatial domain, consider the short exact sequence O + B 4 A) + Q,  U, where AA = L/NAL, BO is the subgroup N'L/NAL of AA, and Qa, = AA,/BO. For 1 < 7 A, let Gy = E"= Nk(G), where G is a residue system for the endomorphism N : L + L. If Ko = No(G9), then K n G,, = {0}, and G\ = K# + G, is a complete set of coset representatives of A\. Hence, for every r E G\, there exist unique ro E K# and r, E G, such that r = r6 + ra,,. We regard KO as the set of highorder representatives in G\; similarly, we regard Ga as the set of loworder representatives in G\. In the frequency domain, let L'\ = (NAL)* and consider the short exact se quence 0 B A', Q 0, where A' = L'/(N*)'L', B, is the subgroup (N*)OL'/(N*)'L of A', and Q' = AI/B^. For 1 < 7y < A, let G'v = E'>=(N*)k(G!), where G' is a residue system for the endomorphism N* : L' * L'. If K', = (N*)O(G'), then K' n G = {0}, and G' = K' + Gn is a complete set of coset representatives of A'. Hence, for every s E G, there exist unique s. E K. and sp E G' such that s = S, + sp. We regard K' as the set of highorder representatives in GC, and G6 as the set of loworder representatives in G'. If a E L(A\), then the discrete Fourier transform a' = F\(a) E L(a) can be computed by the spatial domain CooleyTukey factorization (Theorem 13 and equation 3.6) as a' (s) o/1 = ;/ 1 a (r)(r, s,)) (r,,),) (r,S), raEGa rpEK / for all s = s, +sp E K +G, = G' where = r + r, E K +G = G\. If a' E L(A'), then the inverse DFT a = .T'(a') E L(AA) can be computed by the frequency domain CT factorization (Corollary 14 and equation 3.7) as aro) 1 z(( 1 ) s ,() = *PEG a EK a'(K)(r,,S) s. (,, (, ), for all r = r6 + ra, E K, + Ga = GA, where s = S, + so E K' + G = G'A. The ideals Bp do not split in A\. For example, in dimension two, A,\ Z/7AZ, which does not split since 7 is prime. Thus, the GoodThomas algorithm (as described by An et al. [1], for example) is not applicable to discrete functions defined on the quotients A\. The GoodThomas algorithm is related to the Chinese remainder theorem and is used to compute the DFT of functions with rectangular domains, that is, domains of the form Z,, x*..* x ZL. The algorithm exploits the linear separability of the Fourier transform by first computing the onedimensional transform over the rows of an image and then computing the onedimensional transform over the columns of the resulting intermediate image (similarly for higher dimensions). CHAPTER 5 THE GENERALIZED MATRIX PRODUCT 5.1 Introduction The goal of this chapter is to introduce the generalized matrix product over digital indexing sets. The generalized matrix product, formulated by Ritter [17], generalizes standard matrix multiplication. By using this product, data reordering stages in the FFT algorithm are eliminated. Data reordering is a feature of stan dard FFT algorithms. Zhu and Ritter [26] examined the computational advantages of using the generalized matrix product in signal processing applications. Digital indices address the elements of permutohedral aggregates and their geometric duals. This form of addressing is known as generalized balanced ternary (GBT) and was developed by Gibson and Lucas [8]. The generalized matrix product described by Ritter and Zhu is not easily adapted to data arrays indexed by strings of digits, as in the GBT scheme for permutohedral aggregates. While the generalized matrix prod uct and GBT addressing for permutohedral aggregates have been studied previously and independently, the product defined in this chapter is tailored to processing data arrays that are associated with permutohedral aggregates and their duals. Let F be a ring and denote by F,,xm the set of I x m matrices with entries in F. The definition of standard matrix multiplication, which is a map Fixm x F,,,mxr + Fxr, requires the number of columns in the left factor to equal the number of rows in the right factor. With p E Z+, the generalized matrix product of order p (or pproduct) is a map p : Fixmp x Fpqxr + Fq~mrn. Hence, the pproduct A p B requires that p divide the number of columns in the left factor and the number of rows in the right factor. Note the dimensions of the resulting matrix. We will also make use of the Hadamard product, a componentwise operation on matrices having the same dimensions, defined as follows. For A, B E Flxm, the Hadamard product C = (cij) of A = (a1j) and B = (bj), denoted by A o B, is given by cj = aijbj, for all indices i and j. Hence, the Hadamard product is a map o : Fixm,,, x Fixm + Fixm. Section 5.2 defines Findices, indexing sets consisting of finite strings of digits that address the elements of permutohedral aggregates and their geometric duals (that is, coset representatives of A\ and A', respectively). Addressing with Findices can be implemented in terms of strings of bits, or binary digits, and yields very efficient data processing algorithms (see Gibson and Lucas [8], and Kitto, Vince, and Wilson [13]). Section 5.3 defines the generalized matrix product, or pproduct, over rindices, and Section 5.4 covers digital addressing for permutohedral aggregates and their duals. 5.2 Digital Indexing Sets The elements of Ddimensional permutohedral aggregates and their duals can be addressed using finite strings of digits. We refer to the underlying digital indexing sets as Findices, defined as follows. Definition 4 (Gamma Indices) Let D E Z+ and n = 2D+1 1. 1. Define the set together with the order relationship 2. For A E Z+, define the set FA = {gAg1A2 ... gigo0 : gk E Fi, for k = 0,1,..., A 1} and extend the order relationship on F1 to lexicographic order on r\, where the rightmost digit go is the loworder, that is, least significant, digit; 3. For A,a,,3 E Z+ such that A = a + 0, define the set FrF = {gh :g E r,,,3 E Fr}. We refer to the sets F,\ as Findices. Remark 4 The following properties of Findices follow directly from Definition 4. Fix D E Z+ and let n = 2D+1 1. 1. The set F\ has nA elements; 2. If A, a, # E Z+ such that A = a +/3, then F F, = Fr = FrF,. 5.3 The Generalized Matrix Product for Permutohedral Aggregates The generalized matrix product, or pproduct, over Findices has a particularly simple form. Fix D E Z+ and let n = 2+1 1. Let F be any ring and a,,3 E Z+. Denote by Fr. xr, the set of all no x nO, Fvalued matrices with rows indexed by the set of indices F,,, and columns indexed by the set of indices Fp. In this notation, x E Fr.,, is a column vector of length n' indexed by Fr, and y E Flxr, is a row vector of length n3 indexed by Fr. If transposition is denoted by superscript t, then x E Fr.xI implies xt E Fixr0. Definition 5 (Generalized Matrix Product Over Gamma Indices) Let F be any ring, D E Z+, and n = 2D+1 1. Let A E Fr xrp, where a,3 E Z+, and B E Fr,xr,, where 7,5 E Z+. If A E Z+ such that A < {/,7}, then the pproduct of A and B, with p = nA, is the n'+" x nOA+ matrix C = A p B E Fr.r,_ xrg r, defined by C(gt,rm) = A(g,rs)B(st,m), sErF for all gt E rr.,,rm E 1A s. Remark 5 For full matrices A and B as in Definition 5, computation of the n\ product C = A E,, B requires n++"1+6' multiplications. The case A = = 7 corresponds to standard matrix multiplication. Example. If D = 1, then n = 3 and Fi = {0,1,2}. If x E Rixra and y E Rrxi, then x = (x(000), x(001), x(002), x(010),..., x(221),x(222)) and y = (y(00),y(O1),y(02),y(10),..., y(21),y(22))t. If A = 2 and z = x ,\ y = x 32 y, then z E Rixr, and z(r) = x(rs)y(s), sEF'2 for r= 0,1,2. Example. If D = 2, then n = 7 and F[ = {0,1,...,6}. Suppose A E Cr2, xr and B E Crxr1. If C= A&7B (A= 1), then CE Crsxr2 and C(gt,rm) = 1 A(g,rs)B(st,m), sEri for all (gt,rm) F2F3 x F2F1 = F5 x F3. 5.4 Digital Addressing for Permutohedral Aggregates Our goal in this section is to express the elements of a permutohedral aggregate (and of the corresponding geometric dual) entirely in terms of digits, or Findices. We will need to distinguish between digital addresses for elements of permutohedral ag gregates in the spatial domain and digital addresses for elements of the corresponding duals in the frequency domain. To this end, let Fr = Fr and Fr, = Fr (A E Z+). Fix D, A E Z+. Let L be the Ddimensional lattice generated by U = (uo, u1,..., UD1i) and let N be the endomorphism of L specified by the D x D matrix in equation 4.1 with respect to the basis U. Let U,\ = (u', u',..., un) be an ordered set of D linearly independent vectors that generate the lattice L' = (N*)AL. Let UD = D_ o uj and u'= _o u. Recall the sets of coset representatives G, G', GA, and G', described in Section 4.4. Note that any element of F1, that is, any digit, can be expressed in the form bD2D + ... + b12 + bo, where the bj E {0, 1}, for all j, and not all bj = 1. Furthermore, the bj are unique. Assign the digit bD2D +. +b12+b0 E F1 to the vector bDUD + +b6u +bouo E G. This assignment of digits establishes a onetoone correspondence between Fx and G. Now, assign the string of digits g\.. *.. gigo E FA to the vector r = N\ (r,\) +. + N(ri) +ro C G,\, where gk E F1 is the digit corresponding to the vector rk E G, for k = 0, 1,..., A 1. In this manner, we establish a onetoone correspondence between the digital in dexing set F\ and G\ (a permutohedral aggregate in the spatial domain). Write 9k = rk whenever gk E F1 is the digit corresponding to rk E G and refer to the string 9gA\ gi9o = N'A(gAI) + + N(gi) + go as the digital address of r. In the frequency domain, assign the digit bDo2 + ... + b12 + bo E to bDuoD + S+ blux + bouo E GC'. This assignment of digits defines a onetoone correspondence between Fr and G'. Assign the string of digits h,\\l h1ho E r, to the vector s = NAl(s\i) + + N(si) + so E G', where hk E F' is the digit assigned to sk E G', for k = 0,1,..., A 1. This assignment of digital strings defines a one toone correspondence between the digital indexing set Fr and G' (a permutohedral aggregate in the frequency domain). Write hk = Sk whenever hk E Fr is the digit assigned to Sk E G'; the string hA, ... h1ho = EiI(N*) k (hk) is the digital address ofs. Recall from Section 4.4 that the finite address of a lattice point, if it exists, is unique. Hence, the digital addressing scheme we have described assigns unique addresses to the elements of permutohedral aggregates in the spatial domain; similarly for their counterparts in the frequency domain. Also, note that the elements of Li or the corresponding vectors in Gcan be expressed as bit patterns of the form (bD ... bibo). Similarly for the elements of r' and the corresponding vectors in G'. Example. Let D = 2 and A = 3. The digital address of the vector N2(u2 + Uo) + N(u2)+uo = N2(1.u2+0.u+l1.uo)+N(1.u2+0.ui+O.uo)+(O.u2+0.ui+l1.uo) = N2(101) + N(100) + (001) = N2(5) + N(4) + 1 E G3 is 541 E r3. Note, however, that 541 E rF is the digital address of the vector (N')2(ul + Uo) + N*(u') + U'o E G'3 in the frequency domain. Example. Figure 5.4 shows an indexed level 3 hexagonal aggregate in the spatial domain. Figure 5.1: Indexed level 3 hexagonal aggregate in the spatial domain. CHAPTER 6 DIGITAL ALGORITHMS FOR PERMUTOHEDRAL AGGREGATES 6.1 Introduction This chapter brings together the concepts developed in earlier chapters to give statements of the discrete Fourier transform, the CooleyTukey factorization, and the fast Fourier transformtogether with the respective inverse transformsfor complexvalued functions defined on permutohedral aggregates and their geometric duals. The algorithms in this chapter are expressed in terms of the generalized matrix product over digital indexing sets. The two theorems in Section 6.2 define the discrete Fourier transform (DFT) and the inverse DFT for permutohedral aggregates. The theorems and lemmas in Section 6.3 are digital CooleyTukey (CT) factorization algorithms for permutohedral aggregates. Section 6.4 gives statements of the fast Fourier transform (FFT) and the inverse FFT algorithms for permutohedral aggregates addressed by digital indexing sets. Throughout, assume D, A E Z+ fixed and let n = 2D+1 1. Recall that n) is the cardinality of a Ddimensional, Alevel permutohedral aggregate. We will use the digital indexing sets FA and FP to address spatial and frequency domain data, respectively, as discussed in Section 5.4. As stated earlier, very efficient data processing results from the use of digital indexing sets to address permutohedral aggregates and their geometric duals. Define the map G, x G'  U(nA) by (r,s) = exp[27ri(r, s)], where G\ is a set of coset representatives of A\ = L/NAL, G' is a set of coset representatives of A' = L'\/(N*)A'L', and (r,s) is the inner product of r,s E RD (see Proposition 17); (.,.) is a bilinear form on L x L'. For a E L(A\), define a : GA 4 C by a(r) = a(f), for all r E GA; for a' E L(A'), define a': G'\ + C by a'(s) = a'(), for all s 6 G'.. If g E rA is the digital address of r E GA, and h E r, is the digital address of s E G', then define (g, h) = (r, s). Hence, in the context of digital addressing, (., .) is a map r\ x Fr 4 U(nA). 6.2 The Discrete Fourier Transform and Inverse DFT Let a,/# E Z+ such that A = a + P3. For the direct transform, define e0 E Cixr rr (a row vector of length n2") by e,4(vaq) = (N/(q'),v') E U(nA), q" E r,V r :,. Let e = el. Note that e\(h\gA) = (g\,hA) and e(hg) = (N1(g),h). For the inverse DFT, define e' E Cixr.r; (a row vector of length n2a) by e'(qvc') = (NA(q'),v')1 = (NIA(qa),v) U(nA), q rC r'. Let e' = e'. Note that e'(g\hA) = (gA, hA) and e'(gh) = (N1(g), h). Theorem 24 (Digital DFT for Permutohedral Aggregates) If D, A E Z+, n= 2D+1 1, a E L(AA), and a' = .F\(a), then 1 a' = ;h2 (e\ Gx a)t. Proof. The discrete Fourier transform of a E L(A)) is determined by a'(s) E a(r)(r,s), s E G. (6.1) rEG\ Let g" E r,\ be the digital address of r E G\ and hA E 1F the digital address of s E G'. Then equation 6.1 is equivalent to aG'(hA) = n/' a(g\)(g\, h\) g9'Er n/ 9A E rA n 1~I E e.\(h'g')a(gA), hA E rj\* In terms of the nAproduct, the last summation can be written a' = (e\ a)t (6.2) to yield a column vector of length n. 0 The nAproduct in equation 6.2 requires no more than n2A complex multipli cations. A similar proof yields the corresponding result for computation of the digital inverse DFT, which recovers a E L(A\) from a' E L(A'A). Theorem 25 (Digital Inverse DFT for Permutohedral Aggregates) If D, A E Z+,n = 2D+1 1, a' E L(A'J), and a = ."(a'), then a = 1 (e a,)t. 6.3 CooleyTukey Factorization Algorithms Let a, / E Z+ such that A = a + fl. For the CooleyTukey (CT) factorization of the direct transform, define the n" x no matrix Da,0 E Cr. xr, by DC,(qO',wf) = (qa,wf) E U(n'), q E rw E rf. For factorization of the inverse transform, define the n6 x n' matrix D' E Cr, xr by D',wq") = (qr,ww) E U(n'), wO 6 ,'( ,q' E F,. Lemma 26 Fix D,A E Z+ and let a, 3 E Z such that A = a +3. 1. If p E Fr and q E Fr,, then N0(p) + qa = pqa E r FA; 2. ifv E r' and w1 E r', then (N*)(v') + w, = vau E r,. Proof. Part (1). If p3 = pi pipo and qO = q,_ " qiqo, then paq = P0I . pipoqI ... q qo0 a+131  = Z k(pk.)+ k(qk) k=a k=O k=O k=0 '81 CI1 = NE Nk (pk)+ E Nk(qk) k=O k=O = NA(9)+ q. Part (2) is proved similarly. 0 Lemma 27 Fix D,A E Z+. If r E L and s E L'\, then (NA(r),s) = 1. Proof. We have NA(r) E NAL. Since s E L' = (N'L)*,(NA(r),s) E Z, which implies (N'(r),s) = exp[27ri(N(r),s)] = 1. 0 Theorem 28 (Digital CT Factorization) If D, A, a, 0 E Z+ with A = a + f3, n = 2n+l 1, and a E L(AA), then eA e,\ a = e, en (D0,, o (eO pn a)). Proof. The nAproduct e% E,. a is equivalent to the summation E e\(hAgA)a(gA), hA E r. 9" Er. Ifg\ = 9qo' E rFr, = FA and h\ = v=wo E 'FrF = rF, then e,(hVgA) = (gA, h\) = (pyqaIvaw) = (Na(pl) + q+ (N*)~ (va) + w) (by Lemma 26) = (Na (P),(N*)o(0'))* (N(p,),w'") (qQ,(qN*)fi(v0)) (q,w) = (NA(p9),v). (N"(po),w). (No(q'),v0) (qa, wo) = 1 eo(wlo) e,(vq') D.c,,(q', w), the last equality by Lemma 27. Hence, Z eA(hAg)a(g) = e(v"wpq)a(pq) (6.3) 9AEr, qEra p ErI, = E e,(vaq') [D'.,(q"',wI)Co(q`,w/o)], (6.4) qa'Era where Co(q,uw) = ) ep(wVp)a(p9 q'). (6.5) pOEFo Summation 6.4 is equivalent to e. E,,. (D ,, o Co) and equation 6.5 is equivalent to C# = e# (,o a. Hence, eA Gn a = @n (Da,, o Cp) = o e,1 (D,,, o (eo (O a)), a row vector of length n. 0 The inside noproduct in the CT factorization algorithm ea (no (Da,# o ,(eED a)) requires no more than n2" complex multiplications; the Hadamard product takes no more than nA complex multiplications. Hence, the total number of complex mul tiplications required by the CT factorization algorithm is bounded by n2" + n20 + n'. Remark 6 Theorems 24 and 28 imply a' = 1 (e,, ,,. (D(,0 o (ep ,, a)))t. Arguments similar to those in the proof of Theorem 28 establish the digital inverse CT factorization algorithm. Theorem 29 (Digital Inverse CT Factorization) If D, A, a, /3 E Z+ with A = a + /3, n = 2D+1 1, a' E L(A'), and a = .1(a'), then el eLA a' = e' n (DI, o (el ena a')). Remark 7 Theorems 25 and 29 imply 1 a= (ee ^ (DN,. o (e' (Bn, a')))'. Theorems 28 and 29 are digital versions of the CT factorization for permutohe dral aggregates. However, the next two lemmas, which use variations of the matrices D,,,3 and D,a,, respectively, yield digital FFT algorithms (direct and inverse) for permutohedral aggregates. Define T,,+i,i E Cr+, xr'__ by T.+,, 1(pMq',w1) = (N'(pM),w'01) E U(n'), pMq' E r ++,w' 1 e r, and define Tp++,_1 Cr.+ xr.i by T+1,Q1_(VMWqa1) = (q', (N*)P(vM))l E U(n'), VMW c F +,q'q 1 E I,1i. In the signal processing literature, the entries in matrices DQ,,3,D#,,,T +i,p, and T#+, 1_ are known as twiddle factors. Lemma 30 Let D,A E Z+, n = 2+1 1, and a E L(AA). If 3 E Z+ with 2 < t < A, then eo E)p a = e D (To+1,1 o (e_i .i a)), where a = A 3. Proof. The naproduct e nep a is equivalent to the summation E eo(w'p')a(pAq), q' E r,,w E r,. If P = P'PM E r,r, = ro and w3 = WMW31 E rr_ = rL, then ep(wlp) = (Ng(p),w)3) = (N'('pM),wMw'3^1) = (Na(N(p1) + PM), (N*)3(WM) + wf1) = (Ng+'(p') + Na(pM), (N*)I(wM) + wl1) = (Na+l(pl), (N*)/l(WM)) (N+'(p91), w1) (Na(pM),(N*)Il(wM)). (N2(pM), w1) = (NA(p1),WM). (N 1(p1),w1) (N''(pM),wM) (N'(pM),W 1) = 1 ei(w'p') e(wMpM) o+r ,_,(pm0, "), the last equality by Lemma E ew(w'p3)a(pAq') p Erp 27. Hence, = E E e(WMW ^'8'pM)a(^PlPMq') PMErI pAe1Ero_ = E [Ti + ,w)(PMq'w 1l)l(pMqc,w 1)] e(wpMprM), "e(WMPM), (6.6) 64 where C#(pMqaw8l)= S e#j(wlWpa)a(pPMq). (6.7) PO'Eroi Summation 6.6 is equivalent to e e, (Ta+l,pi o Cpi) and equation 6.7 is equivalent to C,I1 = ej1 n,1 a, a matrix of dimensions no x nO. Hence, ep Enf9a = e (n (Tr+i,p.i oC#_i) = e e, (TM+4,,61 o (ef_i (nPi a)). Similar arguments establish the corresponding inverse digital algorithm. Lemma 31 LetD, A E Z+,n = 2D+11, anda' E L(A'). If a E Z+ with < a < A, then e no a' = e' E (T_+1OI o (e'1_, E, a')), where 3 = A a. 6.4 The Fast Fourier Transform and Inverse FFT Algorithms The following algorithm is a radixn, decimationinspace, fast Fourier trans form for permutohedral aggregates addressed by digital indices. Theorem 32 (Digital FFT for Permutohedral Aggregates) Let D, A E Z+, n = 2D+1 l,a E L(AA), and a' = F(a). If C, = e (. a, (6.8) and C'o = e ,, (TA2+j,.1 o C0i), forf/ = 2,3,..., A, (6.9) 65 then a=Lc Proof. We claim Cp = e, &,, a, for 1 hypothesis. Arguing by mathematical induction, if C, = e eny a, for 1 < 7 < A 1, then Cf+l = e E) (TA_,o 0 CU,) = e n (Ty,y o (e ,, a)) = e.Y+1 .n+i a, (by Lemma 30) for 1 < 7y < A 1. Hence, CD = e3 nO a, for 1 < /3 < A. In particular, CA = e,\ ED a, which implies a' = nC/2 C,\, by Theorem 24. 0 A similar argument proves the corresponding digital inverse FFT for permu tohedral aggregates, a radixn, decimationinfrequency algorithm. Theorem 33 (Digital Inverse FFT for Permutohedral Aggregates) Let D, A E Z+,n = 2D+1 1,a' E L(A'), and a = .7' (a'). If C = e' ,n a', and '= e'. (T o C.), fora = 2,3,..., A, then ] t a = (C\) . Computation of the DFT of a function a E L(A\) consisting of n lattice points or cells requires A stages, as specified by equations 6.8 and 6.9. Computation of the nproduct in each of these equations requires no more than nA+1 complex multiplications. The Hadamard product in each of the A 1 stages specified by equation 6.9 requires at most n' complex multiplications. Hence, the total number of complex multiplications required to compute the DFT of a using the FFT algorithm is bounded by An+l +(A 1)nx. The algorithm, therefore, has complexity O(M log M), where M = n. Furthermore, since the FFT algorithm (Theorem 32) is written entirely in terms of the generalized matrix product and Hadamard multiplication, no explicit reordering of the data is required. Data reordering is a common feature of standard FFT algorithms. Similar remarks apply to the computation of the digital inverse FFT algorithm (Theorem 33). CHAPTER 7 SUMMARY AND SUGGESTIONS FOR FURTHER RESEARCH In this dissertation we showed how the additive structure of permutohedral aggregates can be exploited to derive fast algorithms for the computation of the discrete Fourier transform of functions whose domains are aggregates of any given level in any given dimension. In these algorithms, the data is addressed in terms of an efficient digital addressing scheme known as generalized balanced ternary (GBT) addressing. The algorithms described here are not based on the Chinese remainder theorem, which is used to derive the fast Fourier transform (FFT) for functions de fined on rectangular arrays. Furthermore, we defined the generalized matrix product for permutohedral aggregates and stated the FFT and the inverse FFT for permuto hedral aggregates entirely in terms of digital addressing and the generalized matrix product. The notation described in the text renders very succint algorithms that are amenable to computer implementation. Note, however, that we have not directly exploited the multiplicative structure of permutohedral aggregates. A number of algorithms exist that compute the discrete Fourier transform for domains that possess a multiplicative structure (see Tolimieri, An, and Lu [20]). Hence, we propose the following problems for further research: 1. characterization of ideals and quotients of permutohedral aggregates in terms of GBT addressing; 2. derivation of the discrete Fourier transform of periodic and decimated functions defined on permutohedral aggregates; 68 3. derivation of DFT algorithms that exploit the multiplicative structure of per mutohedral aggregates; 4. extension of the definition of the generalized matrix product for permutohedral aggregates to more general addressing schemes; 5. description of selfsimilar tiles derived from permutohedral aggregates and deriva tion of a discrete wavelet transform and a multiresolution analysis for such tiles. REFERENCES [1] M. An, I. Gertner, M. Rofheart, and R. Tolimieri, Discrete fast Fourier transform algorithms: a tutorial survey, Advances in Electronics and Electron Physics, Vol. 80, Academic Press, Inc., New York, NY, 1991, pp. 167. [2] L. Auslander and R. Tolimieri, Is computing with the Fourier transform pure or applied mathematics?, Bulletin (New Series) of the American Mathematical Society, Vol. 1, No. 6, November 1979, pp. 847897. [3] C. Bandt, Selfsimilar sets 5. Integer matrices and fractal tilings of R", Pro ceedings of the American Mathematical Society, Vol. 112, No. 2, June 1991, pp. 549562. [4] J. H. Conway and N. J. A. Sloane, Sphere packing, lattices, and groups, SpringerVerlag, New York, NY, 1988. [5] J. W. Cooley and J. W. Tukey, An algorithm for the machine calculation of complex Fourier series, Math. Comp., Vol. 19, No. 2, 1965, pp. 297301. [6] C. L. DeVito, Functional analysis and linear operator theory, AddisonWesley Publishing Company, New York, NY, 1990. [7] D. E. Dudgeon and R. M. Mersereau, Multidimensional digital signal processing, Prentice Hall, Englewood Cliffs, NJ, 1984. [8] L. Gibson and D. Lucas, Spatial data processing using generalized balanced ternary, Proceedings of the IEEE Conference on Computer Vision and Pat tern Recognition, June 1982, pp. 566571. [9] R. C. Gonzalez and R. E. Woods, Digital image processing, AddisonWesley, Reading, MA, 1992. [10] T. W. Hungerford, Algebra, SpringerVerlag, Berlin, 1980. [11] M. T. Heideman, D. H. Johnson, and C. S. Burrus, Gauss and the history of the fast Fourier transform, IEEE ASSP Magazine, October 1984. [12] J. Kappraf, Connections: the geometric bridge between art and science, McGraw Hill, Inc., New York, NY, 1991. [13] W. Z. Kitto, A. Vince, and D. C. Wilson, An isomorphism theorem between the padic integers and a ring associated with a tiling of Nspace by permutohedra, Discrete Applied Mathematics 52, 1994, pp. 3951. [14] D. Lucas, A multiplication in Nspace, Proceedings of the American Mathemat ical Society 74, 1979, pp. 18. [15] R. M. Mersereau and T. S. Speake, The processing of periodically sampled mul tidimensional signals, IEEE Transactions on Acoustics, Speech, and Signal Pro cessing, ASSP31, No. 1, 1983, pp. 188194. [16] D. P. Petersen and D. Middleton, Sampling and reconstruction of wavenumber limited functions in ndimensional euclidean spaces, Information and Control, Vol. 5, 1962, pp. 279323. [17] G. X. Ritter, Heterogeneous matrix products, in Image algebra and morpholog ical image processing II, Vol. 1568 of Proceedings of SPIE, San Diego, CA, July 1991, pp. 92100. [18] G. X. Ritter and J. N. Wilson, Handbook of computer vision algorithms in image algebra, CRC Press, Boca Raton, 1996. [19] G. Strang, Wavelet transforms versus Fourier transforms, Bulletin of the Amer ican Mathematical Society, Vol. 28, No. 2, April 1993, pp. 288305. [20] R. Tolimieri, M. An., and C. Lu, Algorithms for discrete Fourier transform and convolution, SpringerVerlag, New York, 1989. [21] R. Tolimieri, M. An., and C. Lu, Mathematics of multidimensional Fourier transform algorithms, 2nd ed., SpringerVerlag, New York, 1997. [22] A. Vince, Replicating tesselations, SIAM Journal of Discrete Mathematics, Vol. 6, No. 3, August 1993, pp. 501521. [23] A. B. Watson and A. J. Ahumada, Jr., A hexagonal orthogonaloriented pyramid as a model of image representation in the visual cortex, IEEE Transactions on Biomedical Engineering, Vol. 36, No. 1, 1989, pp. 97106. [24] H. Zhu, The generalized matrix product and its applications in signal processing, Ph.D. Thesis, 1993. [25] H. Zhu and G. X. Ritter, The generalized matrix product and the wavelet trans form, Journal of Mathematical Imaging and Vision, Vol. 3, No. 1, 1993, pp. 95104. [26] H. Zhu and G. X. Ritter, The pproduct and its applications in signal processing, SIAM Journal of Matrix Analysis and Applications, Vol. 16, No. 2, 1995, pp. 579601. BIOGRAPHICAL SKETCH Jaime L. Zapata was born in Bogota, Colombia. He received the Bachelor of Science degree in engineering with High Honors from the University of Florida in August 1985, and the Master of Science degree in biomedical engineering from the University of Miami in August 1986. From 1986 to 1990, he worked as a software engineer with Coulter Electronics, a major biomedical engineering firm in Hialeah, Florida. His research interests include image algebra, image processing, wavelets and multiresolution analysis, and bioinformatics and biocomputing. I certify that I have read this study and that in my opinion it cnforms to acceptable standards of scholarly presentation and is fully adequate/ i scope and quality, as a dissertation for the degree of Doctor of Philosophy. Gerwrd X. Ritter, Chairman P essor of Mathematics I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. 694se. L&P David C. Wilson Professor of Mathematics I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Murali Rao Professor of Mathematics I certify that I have read this study and that in my opini n it conforms to acceptable standards of scholarly presentation and is fully adequ te, in scope and quality, as a dissertation for the degree of actor of Philosophy. ] 'Ame sE. Keesling f /4 or of Mathematics I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. J ph N. Wilson Assistant Professor of Computer and Information Science and Engineering This dissertation was submitted to the Graduate Faculty of the Department of Mathematics in the College of Liberal Arts and Sciences and to the Graduate School and was accepted as partial fulfillment of the requirements for the degree of Doctor of Philosophy. May 2000 Dean, Graduate School UNIVERSITY OF FLORIDA 3 1262 08555 361W I .. .  I. / "' 