GEOMETRIC DESIGN TECHNIQUE FOR TWO-DIMENSIONAL ANALOG AND DIGITAL FILTERS Cosk&n Zumrutkaya Thesis Z8 NAVAL POSTGRADUATE SCHOOL Monterey, California THESIS TWO- GEOMETRIC DESIGN TECHNIQUE FOR DIMENSIONAL ANALOG AND DIGITAL FILTERS by Coskun Zumrutkaya June 19 7 9 Thesis Advisor: S. R. Parker Approved for public release; distribution unlimited UNCLASSIFIED SECURITY CLASSIFICATION OF THIS PAGE (Whan Data Entered) REPORT DOCUMENTATION PAGE I a.£»OAT SUMII* READ INSTRUCTIONS BEFORE COMPLETING FORM 2. GOVT ACCESSION NO. 1. RECIPIENT'S CATALOG NUMBER 4. T|TlE fend Subtttl*) Geometric Design Technique for Two- Dimensional Analog and Digital Filters 5- TYPE OF REPORT a PERIOD COVERED Engineer's Thesis; June 197 9 « PERFORMING ORG. REPORT NUMBER 7. AUTHORS Coskun Zumrutkaya •. CONTRACT OR GRANT NUMBERS.) >. PERFORMING ORGANIZATION NAME AND ADDRESS Naval Postgraduate School Monterey, California 93940 10. PROGRAM ELEMENT. PROJECT, TASK AREA « WORK UNIT NUMBERS II. CONTROLLING OFFICE NAME ANO AOORESS Naval Postgraduate School Monterey, California 12. REPORT DATE June 19 7 9 IS. NUMBER OF PAGES 153 U. MONITORING AGENCY NAME ft AOORESSf different from Controlling Ottlee) IS. SECURITY CLASS, (of thio report) Unclassified <Sa. OCCLASSIFI CATION/ DOWNGRADING SCHEDULE 16. DISTRIBUTION STATEMENT for thit Report) Approved for public release; distribution unlimited. 17. DISTRIBUTION STATEMENT (o: thm omotrmct entered In Block 30, II different from Report) IS. SUPPLEMENTARY NOTES IS. KEY WORDS 'Continue on rereree tide it neceeeerr and Identity or block nimtber) Analog Filter Digital Filter Extended Bode approach 20. ABSTRACT (Continue an referee eide It I •ary and Identity ey bleek neer) A technique for the use of semi-infinite planes to approximate the log magnitude versus log frequency characteristics of two-dimensional filters is presented (Extended Bode approach) . This technique is applied to quarter plane filters and works well for separable transfer functions. Other non-separable canonic (basic) transfer functions are also studied in terms of their planar approximation. The general technique is shown DO , r ™ n 1473 (Page 1) EDITION OF 1 MOV «» IS OBSOLETE S/N 0103-014-6601 | UNCLASSIFIED SECURITY CLASSIFICATION OF THIS PAGE fWhan Dote Entered) UNCLASSIFIED f*eu*wTv claudication q» Tmt mtQtrw*,^ n»i« i.<, K< (20. ABSTRACT, Continued) to be useful for the insight it provides as well as being a simple approach to design. The double bilinear z-transform is studied for use with the two-dimensional analog transfer functions to convert them to the digital domain. DD Form 1473 „, 1 Jan 73 S/N 0102-014-6601 1 Jan 73 UNCLASSIFIED N 0102-014-6601 2 iieuAW classification o* ▼»••• natr**™ o«« *«#•*•*) Approved for public release; distribution unlimited. Geometric Design Technique for Two-Dimensional Analog and Digital Filters by Coskun Zumrutkaya Lieutenant /''Turkish Navy B.S.E.E., Naval Postgraduate School, 1977 M.S.E.E., Naval Postgraduate School, 1979 Submitted in partial fulfillment of the requirements for the degree of ELECTRICAL ENGINEER from the Naval Postgraduate School June 1979 ABSTRACT A technique for the use of semi-infinite planes to approximate the log magnitude versus log frequency charac- teristics of two-dimensional filters is presented (Extended Bode approach) . This technique is applied to quarter plane filters and works well for separable transfer functions. Other non-separable canonic (basic) transfer functions are also studied in terms of their planar approximation. The general technique is shown to be useful for the insight it provides as well as being a simple approach to design. The double bilinear z-transform is studied for use with the two-dimensional analog transfer functions to convert them to the digital domain. TABLE OF CONTENTS I. INTRODUCTION 8 A. DIGITAL SIGNAL PROCESSING 8 B. AREA OF INVESTIGATION AND OBJECTIVES 18 C. PREVIEW OF RESULTS 20 II. FUNDAMENTALS OF TWO-DIMENSIONAL RECURSIVE FILTERING 21 A. INTRODUCTION 21 B. DEFINITIONS : 22 1. Two-Dimensional Signals 22 2. Two-Dimensional Systems 28 3. Causality, Separability, Stability 31 4. Two-Dimensional Difference Equations - 36 5. Two-Dimensional z-Tranform 39 6. Two-Dimensional DFT 45 C. STRUCTURES OF TWO-DIMENSIONAL FILTERS 4 7 1. One-Dimensional Digital Filter Realization 47 2 . Two-Dimensional Digital Filter Realization 53 D. STABILITY 61 III. MATHEMATICAL BACKGROUND FROM ANALYTIC PLANAR GEOMETRY 69 IV. CANONIC FACTORS IN CONTINUOUS FREQUENCY DOMAIN 87 A. SEPARABLE FACTORS 87 B. NON-SEPARABLE FACTORS 10 8 C. CONCLUSIONS 118 V. TWO-DIMENSIONAL DIGITAL FILTER DESIGN TECHNIQUES 122 A. INTRODUCTION 122 B. TWO-DIMENSIONAL DIGITAL FILTER DESIGN PROCEDURE BY DOUBLE Z-TRANSFORM 125 1. Double z-Transform 125 2. Design Procedure 130 C. DIGITAL CANONIC FACTORS AND STABILITY CONSIDERATIONS 132 1. Introduction 132 2. Separable Factors 132 3. Non-Separable Factors 133 4. Conclusions 145 VI. CONCLUSIONS 148 LIST OF REFERENCES 150 INITIAL DISTRIBUTION LIST 153 ACKNOWLEDGMENT I am fortunate to have this opportunity to thank the many people who have helped in the preparation of this thesis. First of all, I am most indebted to my advisor. Dr. Sidney R. Parker, for his encouragement, patience, guidance and instruction. I want to thank the Turkish Navy for the generous opportunity to pursue this research. I want to thank Mrs. Rosemary Lande for her excellent typing. Last, but not least, I would like to thank my wife, Sumru and my son Devrim, for their patience and understanding I. INTRODUCTION A. DIGITAL SIGNAL PROCESSING Digital signal processing has its roots in 16th century mathematics, especially in the fields of astronomy and the compilation of mathematical tables. Today it has become a powerful tool in a multitude of diverse fields of science and technology. The applications of digital signal pro- cessing varies from low-frequency spectrum seismology through spectral analysis of speech and sonar into the video spectrum of radar systems [1] . In their book, Theory and Applications of Digital Signal Processing [2] , Rabiner and Gold have combined digi- tal signal processing theory, with a variety of applications ranging from sonar, radar, communication, music, seismic and medical signal processing, with digital component technology. This technology is the main driving force for progress in this field as well as the general area of com- puter design. It was also observed by them, that, although the formulation of engineering problems is often as vague as those of the "softer" sciences (such as anthropology, psychology, etc.), the execution of these problems appears to depend on greater and greater accuracy and reproducibility The capability of digital systems to achieve a guaranteed accuracy and essentially perfect reproducibility is very appealing to engineers. A digital filter is defined in [3] to be a computa- tional process in which a sequence of numbers acting as input, is transformed into a second sequence of numbers, or output digital signal; where the term "digital" implies that both time (the independent variable) and amplitude are quantized. The field of one-dimensional digital filtering, which is outlined in figure 1.1, encompasses recursive, non- recursive and Fast-Fourier Transform processing. The terms recursive and nonrecursive, instead of IIR (Infinite Impulse Response) and FIR (Finite Impulse Response) are preferred in this thesis. It is shown that recursive processing is much more efficient than nonrecursive processing. Stockham's method [4] to perform fast convolution, which later became known as the FFT method, improved the efficiency of non- recursive techniques, so that comparisons in one-dimension are no longer strongly biased toward recursive methods [2] . One-dimensional recursive digital filter design depends strongly on the effective and well developed continuous filter design theory. Moreover, the stability analysis of higher order one-dimensional recursive realizations can be solved by investigation of the root distribution of its factored form representation. It is important to realize that continuous domain design techniques and factorization property (the fundamental theorem of algebra) exists only in one-dimension. a CO Eh >H fa en Q i— i Q I 8 w Q B g £ us I H fa fa o „ z fa Eh CO fa fa , M O fa as II fa CD 9 fa H fa O fa S3 H Z en I H fa 388 S l_T 2 O H H Q u S R ^ las§ O H fa < i< H O H O M £ .3 2 H fa U Q 03 2 fa Eh 02 >H cn fa Eh H U H Q fa Z O H Z fa H a i m z o fa o fa Hi > fa fa > o fa fa D 10 Areas of application of one-dimensional digital signal processing are listed in figure 1.1. There are important military and specifically naval applications in addition to those listed, such as missile and torpedo guidance, launcher control system, fire control, combat information collection, dissemination and display. Digital filtering in several dimensions, which is over- viewed in figure 1.2, has gained considerable importance, especially for the two-dimensional case. Figure 1.3 presents an example of two-dimensional low-pass and high-pass filtering Until 196 6, two-dimensional digital filtering was implemented by nonrecursive or convolutional techniques. In this method the output is the weighted sum of unit sample responses of all past input values. The serious disadvantage of the convolutional method is the requirements of a very large number of arithmetic operations. The development of the Fast Fourier Transform (FFT) in 1966, reduced the number of arithmetic operations consider- ably and is used extensively today. Filtering via FFT is accomplished by computing the transform of the input func- tion, multiplying by the frequency response of the filter, and inverse transforming the result. The recursive algorithm, which has in general an infinite impulse response [3] , constitutes another technique for the realization of two-dimensional digital filters. Hall [5] compares the amount of computation required for the three filtering techniques. He shows that, in general, the FFT and recursive algorithms are preferable to the 11 02 T *—* co o S 1-1 s CO En 3 en W >H 2 co H Q J 1 < Z Eh H s O H M Q •r 1 K fa O O fa CO >H 2 05 O £ fa Eh K 1 Eh CS3 i 1 fa OS D Eh U D « Eh CO CO H CO >• < z < fa Eh < Eh CO 05 a ij H & b §£ H CO a z; o z eg fa Eh CO >H CO fa < Eh H U fa <: z o H CO fa S H Q I Z fa o fa HI > fa > o CM fa fa CJ H fa 12 Fig. 1.3. a. An example of two-dimensional filtering (An original photograph) [2] 13 Fig. 1.3.b. An example of two-dimensional filtering (Low pass filtered version of Fig. 1.3a) [2] 14 .'.-^■■'/r.-a"' ,\'^y' : ^*-"' '■' ■ \* ■ * S ' Fig. 1 . 3c. An example of two-dimensional filtering (High pass filtered version of Fig. 1.3a) [2] 15 nonrecursive one in number of computations and storage requirements. Hall demonstrates that the recursive filter algorithm constitutes the best method for large data, i.e., it is the fastest and cheapest. When processing two-dimensional data by recursive filters a fundamental problem exists due to inherent feedback, namely, the problem of numerical stability. Since, in general, two variable polynomials cannot be factored into a product of first and second order real coefficient polynomials in each of the variables, it is difficult to solve the stability problem for two-dimensional recursive filtering. Conse- quently, the majority of papers published in two-dimensional digital filtering deal with the design of nonrecursive filters which are inherently stable [6-10]. There are several papers discussing two-dimensional recursive digital filters. For example [11] and [12] formulate two-dimensional recursive filters by the z- transform and linear difference equations, and although they investigate problem areas related to stability and realization, the majority of problems remain to be solved. The most important applications of two-dimensional digital filtering of digital data are in picture processing and geophysical data analysis. Picture processing can be categorized into: a) Digital image restoration and enhancement, b) Computer pictorial pattern recognition. 16 The methods of image restoration and enhancement are applied to invert degradations, such as aberrations, atmospheric effect, scanning, motion, and to manipulate images to improve viewing phenomena experienced by the human eye. Image restoration and enhancement has been used with great success in biomedicine, i.e., in extraction of quan- titative information from x-ray films, chromosome counting, measuring the extent of arteriosclerosis from arteriograms, in forensic sciences application, i.e., fingerprint image enhancement for automatic classification, and in astronomy, i.e., removal of turbulence from astronomical photography [13]. Pictoral pattern recognition by digital computer has gained great importance in satellite surveillance and mili- tary reconnaissance. The application of two-dimensional digital filters permits the separation of different horizontal scans on magnetic and gravity maps in geophysics and can be used equally well for structural and topographic maps or for any other type of factor which is available in the format of a planar grid [14] . Although digital filtering in several dimensions is gaining importance in medicine, i.e., heart volume measurements, in fire control problems, and to analyze complex electronic circuitry, there exists no general theory concerning structures, analysis and design. These considerations can be summarized as follows. Since one-dimensional recursive digital analysis and design 17 techniques depend strongly on the fundamental theorem of algebra and on extensive utilization of continuous design theory, it represents a special, i.e., a non-generalizable field in the area of N-dimensional recursive digital filtering. Multi-dimensional digital filtering is well developed in two- dimensions but is based, due to the absence of recursive theory analysis and design tools, predominantly on nonrecursive filtering schemes. Because of their inherent advantages, much of the current research is directed towards recursive filtering algorithms and techniques. B. AREA OF INVESTIGATION AND OBJECTIVES One area of investigation of this thesis is to see how the semi-infinite straight line approximation technique of one- dimensional filters (based upon the log modulus/log frequency or decibel vs. log- frequency plots) can be extended to two (or higher) dimensional recursive filters. The advantages of log-modulus approach is that a transfer function in factored form becomes a linear sum of terms when the logarithm of its magnitude is taken, and each of these terms can be approxi- mated by semi-infinite straight line segments in the one- dimensional case, and semi-infinite planes in the multi- dimensional case, when expressed in terms of log-frequency. The following are specific objectives: 1. To postulate typical (canonic) factors for two-dimensional filters and to investigate their properties in the log-modulus/log frequency domain. 18 2. To develop an approximation technique for the design of two-dimensional recursive filters to meet given specifications in the frequency domain. The frequency response of two dimensional filters is inherently a volume in the space of magnitude, log cj, and log CO2 • This volume is the given or desired specification of the filter, which in the separable case can be approxi- mated by the combination of semi-infinite planes. In this thesis we investigate how we can approximate this volume by semi-infinite planes as an extension of the approximation used in one dimension, i.e., by semi- infinite lines. After determining the planar equations which perform the given specification of the filter, we investigate ways to get a transfer function; first, in the continuous frequency domain and second for two-dimensional recursive equations in the discrete variable domain. The results of this study will classify those charac- teristics that can be achieved by cascading simple canonic factors, and indicate limitations that are imposed by trying to achieve a particular frequency characteristic (such as the fan filter) in as simple form as possible. The approach taken here is to investigate the problem in the continuous domain first, and then study the trans- formation from continuous frequency (differential equation) domain to discrete time (sampled data) domain. 19 C. PREVIEW OF RESULTS We have proposed a geometric method to design two- dimensional recursive digital filters. This is a new method which is relatively simple and easy with respect to other design methods. We have investigated properties of some typical (canonic) factors and considered their stability properties. In Chapter II we go through the basic properties of two-dimensional filtering. In Chapter III we provide background in terms of planar geometry for the discussion in subsequent chapters. In Chapter IV we investigate the properties of typical factors and their stability in continu- ous variable domain. The new design technique, including stability considerations, starting with the continuous frequency domain and extending to two-dimensional recursive filters in the discrete domain, is presented in Chapter V. The results are summarized and some comments for further study are made in Chapter VI. It should be noted that the results presented apply to quarter plane filters, as defined in Chapter V. 20 II. FUNDAMENTALS OF TWO-DIMENSIONAL RECURSIVE FILTERING A. INTRODUCTION There are many signals that are inherently two-dimensional for which two-dimensional signal processing techniques are required. The most common examples are, of course, images in which the two variables are the spatial coordinates. Image processing [15], [16] plays an important role in many areas of scientific and technical research. Some examples are satellite imagery radiography, radar, and biomedical images such as acoustical holograms, medical x-rays, and electron micrographs [2]. Not all two-dimensional data come from an image. In prospecting for gas and oil, a common procedure is to monitor the seismic signals produced by detonating an explosive charge [17] . An array of geo- phones situated along a line radiating in a vertical or horizontal direction from the point of explosion is used to record these signals. In this case, the two variables are distance (from the point of explosion measured along the line) and time (the time at which the seismic signal reaches a given geophone or, equivalently , a given distance from the explosion) . Although two-dimensional signals may be processed by one-dimensional systems, it is generally preferable to con- sider using two-dimensional systems. Many of the basic ideas of one-dimensional signal processing may be readily 21 extended to the two-dimensional case. There are some very important concepts of one-dimensional systems, however, that are not directly extentable to two dimensions. It is the goal of this chapter to discuss the basic ideas and techniques of two-dimensional signal processing and to illustrate them into the context of two-dimensional filter design. A more detailed discussion of this theory plus an overview of recent work in two-dimensional filtering is contained in the excellent survey paper by Mersereau and Dudgeon [18] , and H. Chang and J.K. Aggarwal [19] . Other useful references are the first few chapters of [1] and chapter 7 of [2] . B. DEFINITIONS 1. Two-Dimensional Signals One-dimensional signals are functions of a single variable and two-dimensional signals are functions of two integer variables. Consider the two-dimensional sequence x(m,n) where m and n are integer variables, shown in figure 2.1. As in the one-dimensional case, the notation x(m,n) is often short hand for a sampled version of a continuous two-dimensional signal x(s,t): i.e., x(m,n) = x(mT,,nT 2 ) = x(s,t) (2.D s=mT 1 ,t=nT 2 One interpretation of x is a function which assigns a (generally complex) number to each integer ordered pair 22 - 7 y- it/a M * Fig. 2.1. Pictorial representation for two- dimensional discrete array X P £$ flit \-$l> J 23 (m,n) . The function x is not defined when either or both of its arguments (m,n) is not an integer. Such a signal will be interchangeably referred to as either a two-dimensional array, a two-dimensional sequence or simply as a two-dimensional signal . Some useful two-dimensional sequences are defined below and shown in figures 2.2 to 2.4. These include: 1. Digital Impulse or Unit Sample 1 m = n = u (m,n) = \ (2.2) elsewhere 2. Digital Step ( I 1 m,n >_ u ,(m,n) = \ (2.3) elsewhere 3. Exponential m n . n a, a~ m,n x(m,n) = X (2.4) otherwise 4. Sinusoid (Complex) j ((jo^m+oj^n) x(m,n) = e -°° < m,n < +=° (2.5) As seen above, the two-dimensional step is related to the two-dimensional impulse by 24 / ■+ + +-H-4 -t ■f- f- -h-h ■/ -h-h -h +- /-■/- +-/-++- -h- -h ■/- f- -H f 4 A -h A 4 +■ f +■ 4 + f 4 4 4 -h + + 4- 4- 4- -hif- 4 4- 4 -h-f- -h -4 -h -i-4- f- -f- + -h-h-h 4 i- -h -h 4- -h +- -f-f- 4 1-+- j-4- + j-4 4- f -h- 4- ■/- 4 -h -h i--l-4--/-4- -/-4-4--h+- -h4—h 4-/- 4- 4 +■ 4 4 I / / > rn. Fig. 2.2. Digital Impulse 25 - I * * -f 4 f +■ + -t-_zf- X XXX X -h-h + + f XXXXXX -h + + 4-4- XX xxxx -h-h-h h-h.XX XXXX, -h-t- 4- 4- -4 .XXX XXXX +.+ + + i .XiXXxxx / / / / / / y ) } y } ) > " 4--h4-4-4-4-4-4-i--h4- 4- -h -h 4- 4- f -4 4 h^4 4 4--A 4- 4 4 4 4 4 + 4--h -h 4 -h-j-h- 4- 4- -h 4- -h + + -h -/ 4 -L- -h 4 4 4 Fig. 2.3. Digital Step 26 -t- -h H- -h -j- -J- -f- ■f- f- -f- -*- + + ++• f--*- ff f -/- -f- ■+ 4-+ f 4- 4-f- *■ *■ -f- + + -h -h f 1 * f- *-* ■t- -t + -h HfA- J- j- -f- f- + m 4 4 -h -h 4 4 4- + -f- -h -h +■ 4- + 4 + + 4- 4- + -h +- -h -t- +-4-4 -h 4 4- 4 4- 4 4 f- -4 4- 7-/ -h -h A -h +~ 4- Fig. 2.4. Exponential 27 m n u_ 1 (m,n) =11 u o (m 1 ,n 1 ) (2.5a) m, =-<» m, =-oo 2 . Two- Dimensional Systems A two-dimensional discrete system is defined as the unique transformation T that maps an input sequence (x(m,n)} and the initial condition sequence {s(m,n)} into an output sequence {y(m,n)} y(m,n) = T[{x(m,n) },{s(m,n) }] (2.6) In signal processing applications, we assume zero initial conditions. With this assumption, the above relation may be shortened to y(m,n) = T[{x(m,n)}] (2.7) The system characterized by T is (a) linear if T[{ax 1 (m,n) + bx 2 (m,n) }] = aT [ {x 1 (m,n) } ] + bT [ {x 2 (m,n) }] (2.8) for arbitrary constants a and b. (b) shift invariant if y(m-m 1 ,n-n 1 ) = T [ {x (m-ir^ , n-n 1 ) } ] (2.9) for all m, and n, . 28 A hardware and software system which is equivalent to the transformation T of (2.6) is called a digital filter. The class of complex-valued sequences which are square summable will play the role of input and output signals of the digital filter. Linearity and shif t-invariance are indepen- dent properties of a system in that neither property implies the other. Thus the system T[x(m,n)] = h(m,n) x(m,n) is linear but not shift invariant and the system 2 T[x(m,n)] = x (m,n) is shift invariant but not linear For two-dimensional discrete linear systems, the impulse response is defined as the response of the system at (m,n) to a two-dimensional digital impulse input at (m, /n,) with zero initial conditions, i.e., h(m,n;m,,n, ) = T [ {u q (m-m 1 ,n-n 1 ) } ] (2.10) where u (m,n) denotes the digital impulse. In particular, the impulse response of a two-dimensional discrete LSI (Linear Shift Invariant) system is 29 h (m,n;m, ,n, ) = h (m-m, ,n-n, ) or h(m,n) = T[{u o (m,n) }] (2.11) The response of a two-dimensional LSI* system with zero initial conditions is completely characterized by its impulse response. In fact using the identity x(m,n) = I x(m, ,n , ) u (m-m, ,n-n 1 ) l'"l'"o % V "1' m 1 ,n 1 =-« I x(m-m 1 ,n-n 1 )u (m 1 ,n 1 ) m^n^- 00 we get y(m,n) = T[{ x (m, ,n 1 ) u (m-m, ,n-n, ) } ] ra r n r l x(m 1 ,n 1 )T[{u (m-m 1 ,n-n 1 ) }] m 1 /n 1 =- c x (m, ,n, ) h (m-m, ,n-n^) m l' n l = ~°° ' *LSI = Linear Shift Invariance 30 y(m,n) = £ x (m-m^n-r^) h (m, ,n,) m 1/ n 1 =-< = x(m,n) * h(m,n) (2.12) where "*" denotes two-dimensional convolution. In other words, for LSI systems the basic convolution theorem is valid. 3. Causality, Separability, Stability In discussing various features of two-dimensional discrete systems, it is convenient to have spectral nota- tion [19] to represent certain regions of the spatial domain. 1. Let S(a,3) denote a sector with the angular interval (a, 3). In polar form S(a,3) may be described by S(a,3) = ( (r,6) |r > 0,a < 9 < 3) (2.13) where (a, 3) are measured in radians. In particular, we give the following notation for the sectors frequently encountered: S++ = S[0,tt/2]; S-+ = S[Tr/2,ir] S— = S[-tt,-tt/2] ; S+- = S[-tt/2,0] (2.14) S*+ = S[0,tt) ; S*- = S[-tt,0) where "[" and "]" imply including the boundary and ")" implies X not including the boundary. 31 Notice that S++, S-+, S — , and S+- represent first, second, third, and fourth quadrants, respectively, and S*+ and S*- represent half-planes which are symmetric to each other with respect to the origin. (See fig. 2.5.) 2. For a two-dimensional sequence (x(m,n)}, the set E v = { (m,n) |x(m,n) ? 0} (2.15) is called the support of (x(m,n)}. 3. Let {h(m,n)} be the impulse response of a two- dimensional discrete LSI system, and let E, be the support of (h(m,n) }. Then, a two-dimensional LSI system is said to be causal if E, is the set S++, and semi-causal if E, is the set S* + . They are shown in Fig. 2.5. A two-dimensional discrete LSI system is said to be separable if its impulse response can be factored into a product of one-dimensional responses (Fig. 2.6); i.e., h(m,n) = h x (m) h 2 (n) (2.16) If the equation (2.16) is not satisfied, the filter is said to be nonseparable. The advantages of separable filters is that two dimensional convolution (2.12) may be carried out as a sequence of one-dimensional convolutions. This can be seen by rewriting (2.16) as follows. 32 fcm Fig. 2.5. Support Area of Impulse Response (a) Causal (b) Semi-causal S. (c) Semi-causal S*_ 33 X(«r«.) H,(-0 H 2 («.) ) 1' 2 H(z) = H 1 (z 1 )H 2 (z 2 ) Fig. 2.6. Separable Digital Systems 34 y(m,n) = j» h 1 (m^) h 2 (r^) x (m-m-^n-n-^ (2.17a) m, ,n =-<» I h 1 (m 1 ) [ I h 2 (n 1 )x(m-m 1/ n-n 1 ) ] (2.17b) m =-oo n,=-oo I h 1 (m 1 ) a(m-m.,n) (2.17c) m l =_c where a (m-m, , n) is a sequence of one-dimensional convolutions (obtained by evaluating the terms inside the bracket of (2.17b) for each fixed value of m, ) . Equation (2.17c) shows that y(m,n) may be obtained by a second sequence one-dimensional convolution. If both the input sequence x(m,n) and the filter impulse response h(m,n) are separable, then it is readily seen that the output sequence y(m,n) is also separable. In this case we obtain the result [using (2.12) and (2.16)] 00 00 y(m,n) = [ I h ± Cm^b (m-m 1 ) J I I h 2 (n 1 ) c (n-i^) ] (2.18a) m =-« n =-oo = a(m)3(n) (2.18b) where x(m,n) = b(m) c(n) 35 V A two-dimensional LSI system is said to be (BIBO) stable if every bounded input sequence produces a bounded output sequence. This condition is satisfied if and only if the two-dimensional impulse response of the system is absolutely summable, i.e., 00 I |h(m,n) | < c° (2.19) m, n=-°° As in the one-dimensional case, (2.19) can be shown to be a <e- necessary and sufficient condition for stability [20]. One problem with (2.19) is that it can be quite difficult to evaluate for an arbitrary h(m,n). 4 . Two-dimensional Difference Equations As in the one-dimensional case, two-dimensional LSI systems can often be described by a two-dimensional linear constant-coefficient difference equation relating the output y(m,n) of the filter to its input x(m,n). The most general form for such a difference equation is shown by I a(k 1 ,k 2 )x(m-k 1 ,n-k 2 ) - I b ( ^^ y (m-Jl^n-j^) = (k 1 ,k 2 )eR a (i 1 ,l 2 )eR h (2.20) wh ere R and R are finite sets of spatial grid points, called the input and output mask, respectively, and {a(k 1 ,k 2 )} and (bUwJU)} are the set of constant coefficients that 36 characterize the particular filter. Generally, a difference equation has a family of solutions as is the case with a differential equation. The difference equation is said to be recursive if there exists a sequence of computations which yields the output sequence serially from the input sequence and that portion of the output sequence which is already computed/ including the initial condition. An implicit ordering of grid points in the spatial domain is assumed for the serial computation of the output sequence. There are two cases of (2.20): 1. R , K in the S++ (first quadrant), b(0,0) ^ and y(m,n) = for outside of S++ (zero initial conditions), then (2.20) represents a causal system. 2. If R , R. in S*+, b(0,0) ? and y(m,n) = for a o outside of S*+ (zero initial conditions), then (2.20) represents a semi-causal system. In both cases we may assume b(0,0) = 1 without loss of generality and write (2.20) as y(m,n) = £ a (k-^k^ x (m-k-^n-k^ (k r k 2 )£R a I b(l 1 ,!L 2 )Y(m-l 1 ,n-l 2 ) (2.21) (l lf l 2 ) eR b -(0,0) It is clear that (2.21) is a recursive equation, and the implicit ordering, for example may be given as 37 {(0,0); (1,0), (0,1), (1,1); (2,0) ,(2,1) , (0,2) ,(1,2); (2,2);...} for a causal system, and {(0,0), (0,1) , (0,2) ,...; (-1,1) , (0,1) , (1,1) ;...;... ... (-1,2) , (0,2) , (1,2) ,.. .; . ..} for a semi-causal system. If {b(£ 1 ,£ 2 )} = {u Q (£ 1 ,£ 2 )} the equation (2.21) becomes y(m,n) = I a(k 1 ,k 2 )x(m-k 1 ,n-k 2 ) (2.22) (k 1 ,k 2 ) £ R a so that the output sequence is a weighted moving average of the input sequence. The impulse response of (2.22) is unique and is obviously {a(k,,k 2 )}. Such filters are called finite extent impulse filters (FIR) because the impulse response has only finitely many nonzero terms. FIR filters obviously are always BIBO stable because a(k,,k 2 ) is always absolutely summable. (It has only finitely many nonzero terms . ) If, on the other hand, the sequence (b(£,,£ 2 )} has several nonzero terms, then the impulse response will generally have infinitely many terms. We will say that such filters are infinite-extent impulse response filters (IIR) . 38 Such filters may be BIBO unstable, and so it is necessary to devise a stability test for such filters. FIR filters have the advantages that they are easier to design, that controlling the phase response is easier (in particular, it is easy to produce filters which have a linear phase response) , and that since the input- output relationship is a finite-extent convolution, Fast Fourier Transform (FFT) techniques can be used to speed the computation of the output sequence. On the other hand, recursive (IIR) filters generally require less high speed storage and usually require less computation time than FIR filters. Comparison of FIR and IIR filters is contained in [5]. 5 . Two-dimensional z-tranform A useful mathematical tool for representation of a sequence x(m,n) is its two-dimensional z-transform. Given a sequence x(m,n), the z-transform of this sequence is defined to be: < CO oo Z[{x(m,n)}] 4 X(z 1 ,z 2 ) = I £ x (m,n) z^ (2.23) m=-°° n=-°° where z, and z 2 are independent complex variables. The most useful property of the z-transform is the fact that the z- transform of the convolution of two sequences is the product of their z-transforms , i.e.: 39 Z[{x(m,n) }*{h(m,n) }] = X(z 1 ,z 2 )H(z 1 ,z 2 ) (2.24) Therefore, the input-output relationship of an LSI filter may be expressed in the z-transform domain as: Y(z 1 ,z 2 ) = H(z 1 ,z 2 )X(z 1 ,z 2 ) (2.25) where H(z-,,z 2 ), the z-transform of the impulse response is called the transfer function. A two-dimensional transfer function evaluated on the unit circle is called the frequency response of the system, i.e., jw, jco 9 ju>, jo> 9 jui, joo 9 Y(e ,e *) = H(e x ,e )X(e x ,e z ) (2.26) when a two-dimensional LSI system is describable by the difference equation of (2.20), and we take a z-transform >' of this, it follows that -k -k 2 [ I a(k 1 ,k 2 )z 1 z 2 ]X(.z 1 ,z 2 ) (k x ,k 2 ) eR a - [ I b(£ 1 ,£ 2 )z 1 z 2 Z ]Y(z L ,z 2 ) = (i lf l 2 ) el^ 40 or H(z 1 ,z 2 ) = Y(z 1 ,z 2 )/X(z 1 ,z 2 ) ' k l ~ k 2 I a(k 1 ,k 2 )z 1 z 2 (k l' k 2 )eR a -A, -l I b(£ 1 ,£ 2 )z 1 x z 2 z ( *■ -i / & ^ ) s^j-. = A(z 1 ,z 2 )/B(z 1 ,z 2 ) (2.27) where A and B are the z-transf orms of (a(k,,k 2 )} and {b (£■•* & 2 ) } respectively. In particular when (2.27) repre- sents a causal or semi-causal system, i.e., (2.27) is the z-transform of (2.21), we can give explicit descriptions for R and R, in each case. For the causal case we define a d (Fig. 2.5a) R a = * ( k i' k 2> I ° i k l i M l' °i k 2 - M 2 } (2.28) Rj^ = ( (l lf i 2 ) | <_l 1 <_N lf 0£^ 2 - N 2^ For the semi-causal case we define (Fig. 2.5b) R a = { (k i' k 2 } 1° i k l- r V k 2 = ° ; " M ct i k li M 6'° < k 2 - M 2 } (2.29) R b = {(£ 1 ,£ 2 )|0<^ ll N 3 ,£ 2 = 0;-N al £ 1 <N a ,Q<£ 2l N 2 } 41 HI A +) •H +j 3 a e o u en ■H a; o c <d cr w 4-1 3 a 4-> a M O 0) 4J (U H 4-1 <44 5 4J C 144 13 O <T5 3 s cr ft ^4 4J n3 ^4 •H -H Q 4-J (T3 CM en •H fa 42 C+ 2 — i -4 1 jC ■P •H T3 Q) -U a, g o o CO •H <D O c CD a 1 CO -u o a 4-> M o CD O 4-1 <-< w 3 e to (C U U I 03 g •H CD a co IT) CN CT> •H fa 43 and in both cases b(0,0) = 1 is assumed. The inverse z-transform is defined as x(m,n) = j j / X(z,,z ) z™ 1 z^ dz,dz„ (2.30) (2njr n „ 12 12 12 C l c 2 where c^ and c 2 are suitable closed contours in the z, and z 2 planes. The two-dimensional Fourier transform of a sequence is defined as DO)., 3u ? X(e ,e *) = F[x(m,n)] X(z 1 ,z 2 ) 2 (z 1 ,z 2 ) eF - j ( co-j m+cu^n) £ x(m,n)e x * (2.32) m r n=-°° For any given sequence the set of (z,,z 2 ) for which the z-transform converges is called the region of convergence. Uniform convergence of the z-transform requires that the sequence be absolutely summable, i.e. oo I |x(m,n) | iz 1 |" m |z 2 |" n < =o (2.32a) m, n=-°° 44 Similarly, the Fourier transform X(e ,e ) converges uniformly to a continuous function of co, and w_ if I |x(m,n) | < oo . (2.32) m,n=-°° A power series of the form of (2.23) is known as a two- dimensional Laurent series in the theory of complex functions in two variables . 6 . Two-dimensional DFT Let {x(m,n)} be a finite support sequence of the size MxN. Then the two-dimensional discrete Fourier transform (DFT) of (x(m,n)} is defined as M ~ 1N ~ X -^ ( R k i + £ k 2> X(k 1# k 2 ) =11 x(m,n)e M L N z (2.33) m=0 n=0 for <_ k. < M-l and <_ k p < N-l. Similarly, the inverse discrete Fourier transform (IDFT) is given by 45 X(m,n) M-l - I MN k,=0 N-l I k 2 =o x(k lf k 2 )e k, k 2 +:2TT( ir m + ir n) (2.34) Note that the DFT corresponds to sampling the Fourier trans- form of MxN points, i.e., X(k x ,k 2 ) = x (e Ua>i 300. /e k l k 2 a) l =27T lT' a) 2 =27T -N- (2.35) ]W. 30). More generally, if X(e ,e ) is the Fourier transform of (x(m,n)} whose size is greater than MxN, then the IDFT of {X(k 1 ,k )} will result in aliasing. In other words, the IDFT of (X(k,,k 2 )} will yield one period of (x(m,n)} where x(m,n) = m 1 ,n 1 =-<» x (m+m,M,n+n,N) (2.36) 46 C. STRUCTURES OF TWO-DIMENSIONAL FILTERS 1. One-dimensional Digital Filter Realization One-dimensional digital filters, in general, can be shown in the z-transform domain as follows: y 1=0 Y(z i ) H < z l> = ~ = xuJt (2 ' 37) 1 - I a m Z l m=l or , in the time domain M L y(n) I a y(n-m) + J b x(n-l) (2.38) m m=l 1=0 The realization implementing equation (2.38) directly, is known as the direct form 1 and shown in Figure 2.7. Equation (2.37) can be rewritten in a slightly different form by introducing a new variable, W(z 1 ) , such that W(z,) Y(z,) H / z ) = L_ . . x s (2.39) H(Z 1 J X(z x ) W z 1 The corresponding set of difference equations consists of: 47 I N J t— 1 2 t 3 I H as kl M r e o t^ +j o cu u •H Q CN 0^ •H fa o 5< 48 M w(n) = x(.n) + I a m w(n-m) (2.40) m=l and L y(n) I b £ w(n-£) (2.41) 1=0 Equations (2.40) and (2.41) are realized in the direct form 2 as shown in Figure 2.8 which can be visualized as a cascade realization of two digital filters, realizing the denominator and numerator polynomials, respectively. Figure 2.8 can be redrawn as shown in Figure 2.9. The resulting realization is called direct form 3 or canonic form, since it represents a structure with the number of delays equal to the order. At this point, an additional realization will be introduced, the direct form 4, which combines the charac- teristics of the direct forms 1 and 2 having only two entries in each summer, which corresponds directly to hardware implementation and by realizing the numerator and denominator polynomials separately in cascaded form (Figure 2. 10) . There are several other methods to realize a one- dimensional M order digital filter. For example, the cascade form of first and second order sections H(z,), as shown in Figure 2.11, where 49 H(z X ) 2 b £ Z l £=0 -£ M 1 " I a z. L ml m=l -m X(n) W(n) I «— 2 H «, r H rf 1 ^ 3 V 7-1 Fig. 2.8. Direct Form 2 50 H(z 1 ) £=0 * L M i - y a z" L ml m=0 -m Y(-) Fig. 2.9. Direct Form 3 (Canonic) 51 I X t~o II J ^ N WH l^H- u-jK JO cv^~~ ■ ^ a CjW- 3 — 4 _ — -± KT c X £ O fa -P o CD S-i •H Q o H •H fa 52 H(z 1 ) = a Q n H i (z 1 ), * = ^ (2.42) i=l Each component H^(z^) can be realized in one of the above outlined forms. Partial fraction expansion methods applied to equation (2.37) lead to a parallel form realization of first and second order sections (Figure 2.12), i.e., k H(z 1 ) - c+ J H i (z 1 ), k= [^+1] (2.43) i=l Other forms of realization include hybrid structures, i.e., parallel-cascade forms, the transpose configuration, which can be obtained for all previously outlined structures by reversing the direction of signal flow and by interchanging all branch and summing modes, wave-digital and continued fraction expansion filters. The direct forms discussed above and the series/parallel arrangement of lower order sections are by far the most often used ones in computer simulation and in digital hardware. 2. Two-dimensional Digital Filter Realization Two-dimensional Quarter plane digital filters, in general form, can be shown as 53 N ii n" x s u fa 0) o en u CM •H fa 54 HW = C+2H } U,) Fig. 2.12. Parallel Form 55 £ =0 £ =0 X Z H(z,,z ) = — - (2.44) M l M 2 -m.. -m 1 ~ I I a mm z, z L *< m.iru 1 2 m l =0 m 2 = ° The numerator of the two-dimensional transfer function (2.44) can be written as a one-variable polynomial in z~, weighted by coefficients which are functions in z,, i.e., L 2 L l I ( I «,,=0 J..=0 L 2 ' S V £ 1 _il 2 z 2 (2.45) " £ 2 :„ (2.46) The same reasoning applied to the denominator polynomials leads to M M -m, -itu D ( z r z 2 } = 1 - ^ ( l ^2 Zl )Z2 m,=0 m.=0 (2.47) M 2 = 1-1 (V. < z l»» Z 2 m (2 - 48) m 2 =0 56 To obtain the direct form 3 representation, a new function, W(z 1 z 2 ) is introduced: W(z,,z 9 ) Y(z,,z 9 ) H(z i z 2> - xi^Tzf) ■ w(^fr (2 - 49) where each factor can be written using equation (2.46) as W(z x ,z 2 ) = *2 £ 2 -0 / N I 1 £ 2 (Z 1 ))Z 2 *<V*2 J (2.50) and equation (2.48) as M 2 -m 2 Y(z 1 ,z 2 ) = W(z 1 ,z 2 ) + I (D M m (z 1 ))z 2 Z Y(z 1 z 2 ) m 2 =0 ^ 2 (2.51) The realization of equations (2.50) and (2.51) is shown in Figure 2.13 using compact notation, and Figure 2.14 where N L l ( z l^ anc ^ D M m ^ z l^ are ^ m Pl emented f° r each (L-,^ 2 ) and (nunu) by one-dimensional direct form 4 realizations. The general form of two-dimensional causal digital filter (2.44) can be rewritten by introducing a new variable, W(z-,z 2 ) or in the previous function, then, W ( z , z « ) Y ( z , , z 9 ) H <V Z 2> " X( Zl ,z 2 ) • W(z',z 2 ) (2 " 52 ' 57 e o •H -p +J 2 -u o (d a e o u a e o +J u aj u ■H Q m H •H 58 Q CO e S-l Pn -P U <U u •H Q ■H fa 59 Each factor can be written as : L l L 2 -I -I Y(l l'»2 J = I I W Z l 1 z 2 2 WCz r z 2 ) £ 1 =0 £ 2 =0 (2.53) and M M -m -m« W(z r z 2 ) =- XC V z 2 ) + I I a z x z 2 W(z r z 2 ) m 1 =0 m 2 =0 - 1 - (m 1+ m 2 )^0) (2<54) The equations (2.53) and (2.54) can be written as a differ- ence equation as follows, L l L 2 y(m,n) = I I b z z W(m-2. 1 ,n-£ 2 ) (2.55) l l =0 £ 2 =0 1 2 and M l M 2 w(m,n) = x(m,n) + I I a m m w ^ m " m 2 ' n ~ m 2 * (2.56) m 1 =0 m 2 =0 m 1 ,m 2 ^ 60 Without loss of generality, L. is chosen equal to M. for all i. The structure realizing equation (2.55) and (2.56) is shown in Figure 2.15. It is noted that the number of delays in Figure 2.15 equals the order of H(z..,z ). The structure is, therefore, by definition, canonic. The direct form 4 (Fig. 2.10) realization of H(z.,z 2 ) is shown in Figures 2.16 and 2.17, using compact notation and one-dimensional direct form 4 (Fig. 2.10) implementation of each N « and D , respectively. L 1^2 M l m 2 D. STABILITY Two major problems in the design of a recursive filter are approximation and stability. Since the output values are used through feedback by the recursive algorithm, it is possible for the output values to become arbitrarily large independent of the input. A filter of this kind is said to be unstable, which is an undesirable condition. Thus we need to know what constraints to put on the recursive filter coefficients, so 'that the filtering operation will be stable. This problem is similar to the stability problem for the one- dimensional case, except that the added dimension inherently increases the complexity of the analysis. The major diffi- culty is that the concept of zeros and poles does not hold in multi-dimensional systems so that simple algebraic systems extrapolation of one-dimensional results is not obvious. • Stability Theorems : Theorem 1 (Shank's theorem): A causal recursive filter with the z-transform H(ZwZ 2 ) = A(z 1 , z 2 ) /B (z 1 , z 2 ) , where A 61 (Q*— I — hQ m (hh— I — ►© M Qk 4 © t-j TCNJ c c U Q CN CO e o -U u CD U •H Q U") CN •H fa 62 rc\j Tsl og — @* — M rsj + i x c o ■H +J id -u c +J o td cu 2 o Q u o 4J o 0) S-l •H Q ^o CN Cn ■H 63 rsj H % J5' Q CM 6 S-i o fa u 0) CM t7> ■iH fa and B are polynomials in z, and z 2 , is stable if and only if there are no values of z, and z such that B(z ,z„) = for all z. _> 1 and z« >^ 1. (In general we assume B(Z 1' Z 2 } = b 00 + ^O 2 !" 1 + b 01 Z 2" 1 + ^l 2 !* 1 ^" 1 ••" that is an expression in increasing powers of z. and z 9 .) In other words, the theorem says that if there are any values (real or complex) of z, and z 2 for which B(z,,z ? ) is zero for z. and z 2 simultaneously greater than or equal to unity in magnitude, then the filter will be unstable. It is much more difficult to determine the stability of two-dimensional filters than of one-dimensional filters. In the one-dimensional case, it is only necessary to locate finite set of roots in the z-plane. In order to mechanize the theorem, we must, in general, find the values of z. and z_ for which B(z,,z~) = 0. There is no technique for locating the zeros of a general two-dimensional polynomial in (z,,z_). Shanks, Tretial and Justice, in their paper [11] proposed a technique to apply theorem 1 to determine the stability of a two-dimensional recursive filter, H(z,,z 2 ), by which the unit disk in the z, plane must be mapped into the z^ plane by solving the implicit two-variable denominator polynomial of H(z,,z 2 ). A necessary and sufficient condition for stability is determined if the map of the unit disk d, = (z-,z_ _> 1) does not overlap the unit disk d 2 E (z ? ,z > 1) in the z plane. This method requires an infinite number of mappings and, therefore, can not be applied exactly. 65 Huang [20] simplified Shank's stability theorem considerably by stating the stability theorem as follows Theorem 2: A causal filter with a z-transform A(z, ,z.J H(z w z n ) = 'l'"2' B(z 1 ,z 2 ) where A and B are two-variable polynomials in z,,z , is stable if and only if: 1) the map of the unit circle 6d = (z , | z, | =1 ) , according to B(z ,z 2 ) = 0, lies outside of the unit disk d_ = (z ,|z |_>1) in the z -plane, and 2) no point in the unit disk d.. = ( z ^ , | z , | _> 1 ) maps into the point z = by the relation B(z..,z 2 ) = 0. In order to apply this theorem as a stability test, we have to map the unit disk 5d, = (z , |z 1 |=1) into the z ? plane according to B(z, ,z 2 ) = or z 2 = f (z,) | i i-, and to see whether the contour in the z^ plane lies inside the unit disk d 2 = (z- , | z 2 |<_1) . Also, it is necessary to solve B(z,,0) = to see whether there are any roots with magnitude greater than 1. Although theorem 2 is much simpler than using Shank's original theorem, the procedure is still infinite in the sense that, in principle, we have to map the unit circle 5d, into the z ? -plane. Ansell [21] , whose main contribution is to couple the use of a Hermite test for checking stability with a series of 66 Sturm tests checking positivity, has reduced theorem 2. Although Ans ell's results enable us to test stabilty in a finite number of steps, it still is, unfortunately, very tedious. Make the change of variables: and, 1 - Z 2 1 + z 2 and let E(p, /P?) H(ZwZ ) = l'~2 fcp 1 ,p 2 ; where E and F are polynomials in p, and p~. We can restate Theorem 2 as follows. Theorem 3 (Ansell's theorem): The causal recursive filter H(z,,z 2 ) is stable if and only if: 1) for all real finite w , the complex polynomial in p„, F(joj,p 2 ) has no zeros in ReP 2 > 0, and 2) the real polynomial in p-,, F(p,,l) has no zeros in Rep, 0. This theorem is essentially the same as Theorem 2, but an advantage of Theorem 3 is that condition 1) can be tested using standard techniques of circuit theory. 67 Anderson and Jury [22] modified Huang's method by out- lining a procedure which replaces the bilinear transforma- tions by the construction of a Schur-Cohn matrix and checking for positivity of a set of self-inversive polynomials. It also replaces the Hermite test by a Schur-Cohn matrix test and requires a series of Sturm tests or equivalently , a system of tests establishing the roots distribution of a polynomial. These modifications represent a substantial reduction in computations as compared to Huang's method. We will discuss stability considerations which are related to our design technique in Chapter IV and Chapter V in detail. In this chapter, basic properties of two-dimensional digital filters have been reviewed and summarized including semi-causality to establish a fundamental theory of two- dimensional digital filters. 68 III. MATHEMATICAL BACKGROUND FROM ANALYTIC PLANAR GEOMETRY In this study, we propose using semi-infinite planes to approximate a given frequency response of a two-dimensional analog filters by extending one-dimensional semi-infinite straight line approximation to two-dimensional cases. For example, as we will see, the frequency response (dB vs. log- frequency) of the separable canonic factor H(ja) 1 ,j(u 2 ) = cl+jUiTi ) (l+ja) 2 T 2 ) can be approximated by a collection of planes. This is discussed in Chapter IV in detail. To prepare for the discussions in the next chapters, it is worthwhile to review and develop useful equations from analytic planar geometry which are going to be used throughout this study with some examples of how to use them. The properties of the analytic planar geometry can be summarized as follows: 1. The general equations for a plane are given by Ax+By+Cz+D=0 (3.1) where A, B and C are the direction numbers of the normal to the plane. The equation of a plane passing through the point (xwy^z,) and perpendicular to the line which has direction numbers of A, B and C is given by: 69 A(x-x 1 ) + B(y-y 1 ) + CU-z^ = (3.2) 2. The direction cosines of the normal to the plane are given by- cos a = — ^ = T-TTo (3.3a) (A Z +B Z +C ) 1/Z cos 6 = — ? ^ — TT7T (3.3b) (A z +B +C ) ±/Z cos y = {A 2 + B 2 + C 2 ) l/2 where a, 8/ and y are the angles that the normal makes with the x, y, z axes respectively. The relationship between direction cosines is given by 2 2 Cos a + Cos 3 + Cos y = 1 (3.3d) 3. For a line joining points (x, ,y^/Z-^) and (x 2 ,y2> z 2) the equation is given by Cos a _ Cos 3 _ Cos y ,^ ^, x 2 -x x y 2 -y! z 2 -z 1 4. The equation of a straight line with the direction numbers a, b, c and passing through the point (x^y^z^ is given by x- Xl y-y 1 z- Zl a ~ b = n — (3.5) 5. The angle 6, between two intersecting lines with direction numbers a, b, c and a', b f , c' respectively is given by: Cos = Cos a Cos a 1 + Cos 6 Cos 3' + Cos y Cos y 1 (3.6a) or Cos 9 = aa ' + bb ' + cc ' - (3.6b) Va 2 +b 2 + c 2 Va ,2 +b ,2 + c'2 6. For parallel lines, the relationship between direction numbers is given by a' b^ c* U * /J 7. For perpendicular lines: aa' + bb' + cc ' = (3.8) 8. The direction numbers of the line of the intersection between two given intersecting planes, with direction numbers A, B, C and A' , B', C respectively, is given by 71 a = b = B B* C A C A' (3.9a) (3.9b) c = B A' B' (3.9c) The angle 6 between these two planes is given by Cos AA' + BB 1 + CC* V A 2 +B' + C V A ' (3.10) 2 2 2 9. The equation for a plane in terms of its slope can be calculated as follows: As seen from Fig. 3.1, the slope, m, of the normal line is given by m = tan h g ^hFW (3.11) m = slope of the plane M = tan (90° - 9) = - 2. = . 1 h m (3.12) 72 normal to olane normal to nlane plane Fig. 3.1. A General Plane Configuration 73 The equation of the plane is given by one of the following: A(x-x 1 ) +By + Cz = (3.13a) Ax + B(y-y,) + Cz = (3.13b) and Ax + By + C(z-z 1 ) = (3.13c) A*l = By ± = Cz 1 = -D (3.13d) The intersection of the plane with the z-x plane is given by A(x -x 1 ) + Cz = (3.14) and the slope, nu, of this line is »i = a§ = -§ (3 - 14a) Similarly, the intersection of the plane with the z-y plane is given by B(y -y, ) + Cz = (3.15) 74 and the slope of this line, nu , is m 2 ~ " dy~ " " C (3.15a, Substituting (3.15a) and (3.14a) into (3.13) yields the following alternate equations for the plane. -m(x-x-j) - m 2 y + z = (3.16a) -m,x - m 2 (y-y 2 ) + z = (3.16b) -itux - m 2 y + Cz - Zj) = (3.16c) The direction cosines of the normal in terms of m-, / m 2 are -m, Cos a = (3.17a) Vm^ 2 2 -m~ Cos 6 = (3.17b) v^7 m 2 2 + l Cos Y = Vm-L 2 + m 2 2 +1 (3.17c) Substituting (3.14a) and (3.15a) into (3.11) and (3.12) yields m = i- ( 3 - 18a) m 2 75 and M = -V m i +m 2 (3.18b) 2 The angle, 8, between two intersecting planes from (3.10), (3.14a) and (3.15a), is given by m,m, * + irunu ' + 1 ,- ,„, Cos 9 = LJ: Li (3.19) ymj 2 + m 2 2 + 1 Vmp+mp + 1 The direction numbers of the line of intersection, from (3.9), (3.14a) and (3.15a) are given by a = (m£ - m 2 ) CC (3.20a) b = (m£ - m 1 ) CC* (3.20b) c = (m|m 2 - m^ir^) CC * (3.20c) The direction cosines of this intersection line are given by m 2 ~ m 2 (3.21a) Cos a = V2 2 (ml - nu) + (m| - m, ) + (m-jnu - mlm^) 2 m' - m, Cos 6 = (3.21b) V2 2 (mi -nu) + ^ m i ~ m i) + (m^m 2 -m^ir^) 76 m l m 2 " m 2 m l ,, ,. , Cos y = (3.21c) V^ 7 2 2 2 m 2 ) + (m* -m,) + (m|m - m'm, ) 10. When two planes are added in one dimension to form a new plane, the following results. Given two planes -m,x - nuy + z - z, = -m,'x - mly + z 1 - z' = Their sum in the z-direction with z" = z* + z is given by -(m 1 +m')x - (m 2 +m 2 ')y + z" - ^ + z£) = The slopes of the intercepts with the zx and zy planes m 1 " = -(m 1 + m|) m 2 " = -(m 2 + m 2 ') The slope of the normal is given by 2 , ,_ , .,2,-1/2 m " = [ (m, +m[) + (m 2 +m 2 ) ] (3.22 77 To aid in comprehension, two examples follow: Example (1) : Given m, (slope in z,x plane), m, (slope in x,y plane) , and x, (the point at which x-axis goes through the plane). Find the equation of the plane. The general plane equation (3.1) is given by Ax + By + Cz + D = When z = B D x A Y " A m 3 = If = " I (given) and When y = and x, = - j (given) A D " " C x " C ■l - E - "I (given) D Z l C 78 If we substitute these into the general equation, the result is or A Ax - Am-,y - — z - Ax, = 3 J m, 1 x - nuy - — z - x, = 3 m, 1 or z = nux - nunuy - x i m i (3.23) As shown here, we can derive the plane equation easily from the given specification. After finding this equation, we can find the canonic factor by simply letting x = log go. y = log a)- i i 2 z = 10 log |T| Alternately, we can go from the canonic factor to equations. Thus equation (3.23) can be written as 2 10 log |T(joo 1 , joo 2 ) | = m 1 log oj 1 - m^ log oj 2 - 10 log K (3.24) where 10 log K = ^i x i 79 This result assumes that i i2 w l |T(ju>,,ju> 2 ) | = (3.24a; Taking the logarithm of (3.24a) yields 2 10 log |T| = 10 p log u^ - 10 q log oj 2 - 10 log K (3.25) By equating (3.25) and (3.24), we get 10 p = m 1 (3.26a) 10 q = m 1 m 3 (3.26b] It follows from (3.26a) and (3.26b) that m 3 = q/p . Thus the transfer function of (3.24a) has a logarithmic characteristic as sketched in Fig. 3.2. The previous result can be extended to the following transfer function (1 + j Wl ) p/2 T(j<D, / ja> 9 ) = —t-t? 7TT? (3.27) 1 Z K ±/Z (l + jw 2 ) q/ By substituting (3.26a) and (3.26c) into (3.27), then 80 Z=io 1^1T| ■s loovJ *=: looi*/| Fig. 3.2. The logarithmic characteristic of P T(jgo 1 , jo) 2 ) jj- Kcu. 20 log |t| = m, log co, - m, m log aj_ - 10 log K 81 nu/20 (1 + jco,) L T( >l' jw 2 ) ■ — Hu/20 (3 * 28) K 1/Z (l + ja) 2 ) Z When w, << 1 and oj 2 << 1 , region I T(jo Jl /jco 2 ) = K~ 1/2 and 20 log |T| = - 10 log K When co, << 1 and oj- >> 1* region II T(j Ul ,ja> 2 ) = — ^-720 K 1/2 (ja J2 ) A and 20 log |T| = - 10 log K - nu log cj. When oj 2 << 1 and oj, >> 1, region III 111,720 (j«2.) T(ju) 1 ,jo) 2 ) = -j/2 K and 20 log |T| = - 10 log K + m, log 0). 82 When u>, >> 1 and oo- >> 1, region IV nu/20 j OJ- T(ju), ,jw~) = , 1 2. . nu/20 K X / 2 H 1 2 ' and 20 log |T| = - 10 logK + m 1 log oj, - m, log to,. The resulting planar approximation for the logarithmic transfer function is shown in Figure 3.3 Example 2 ; In this example, we derive planar equations from the given transfer function. Given a separable two-dimensional transfer function: T(jw, , ju 9 ) = ■= (3.29) (1 + JU) 1 T 1 ) P (1+ JU) 2 T 2 ) H The magnitude square of the transfer function (3.29) is given by K 2 |T(ja) 1 ,ja) 2 )| = 2 2~p" 2 2~q < 3 - 3 °) 1 A (1 + o Ji Z t 1 Z ) P (1 + ^ 2 T 2 } Take the logarithm of both sides of (3.30), 10 log |T| 2 = -10 p log(l+co 1 2 x 1 2 ) - 10 q log(l+u> 2 x 2 ) + 20 log K (3.31) when w,!, << 1 and co 2 x 2 << 1, the equation (3.31) becomes 83 A fcrtoUjITl t 7^ta>W IE Br 3- l°j^ > TUlaj^l Fig. 3.3. The logarithmic characteristic of (1 + jo)) Mi/20 X Z K ±/Z (l + jco 2 ) il 2/ zu 20 log |t I = -10 log K + m, log u>, - m- log go- 84 2 10 log |T| = 20 log K = C where C is some constant. When w^t^ >> 1 and oj 2 t 2 >> 1, or 10 log |T| = -20p log u^ - 20q log a, T + 20 log K i i 2 10 log |T| = -20p log o^ - 20q log a> 2 - 20p log t. Let -20q log x 2 + 20 log K 2 z = 10 log |t| x = log to. y = log a). and we get the planar equation, Z = - 20px - 20qy - 20p log t, - 20q log t 2 + 20 log K or 85 p q T T (-20p)x + (-20q)y - z = 20 log (-± ? ) K or n^x + m 2 y - z - v where m^ = slope of the line (x,z) plane = -20 p db/decade. m 2 = slope of the line (y,z) plane = -2 0q db/decade. In a similar manner the regions (co, t, << 1; CJ 2 T 2 >> "^ an< ^ ^l 1 ! >> ^' (J °2 T 2 <<c ^ can ^ e i nvest: '- < ? ate( i' The planar approximation for the logarithmic transfer function is sketched in Fig. 4.1 where q = p = 1. In this chapter we have developed some useful planar geometry formulas and shown how they can be used to determine a planar approximation for the two-dimensional logarithmic transfer function for separable cases. Some non-separable factors are discussed in the next chapter in which experimental verification of this approximation technique is also presented. 86 IV. CANONIC FACTORS IN CONTINUOUS FREQUENCY DOMAIN In Chapter III, we discussed some equations of planar geo- metry and showed how to use them for approximating separable two- dimensional transfer functions. In order to improve the scope of the proposed design technique, we now consider the behavior of other separable and non-separable factors insofar as their planar geometry approximation is concerned. This will provide insight into the frequency response characteristics of two-dimensional systems and indicate what can be expected by cascading these canonic factors, knowing their individual behavior. The logarithmic characteristic of cascaded factor is simply the sum of their individual logarithmic character- istics. This approach parallels the semi-infinite line approximation (dB vs. log w ) technique which has been used so successfully in one-dimensional design. A. SEPARABLE FACTORS The following several separable canonic factors are now considered. (1) (1 + jaj 1 T 1 ) P (l + jt 2 u3 2 ) q (4. a) i (2) (1 + jaJ 1 T 1 ) P (l + ju) 2 T 2 ) q (l + joo 2 T 3 ) q (4.b) (3) [1 + (ja Jl T 1 ) P ] a • (4.c) where p, q and a are integers. 87 (1) First, as an example of the technique, consider the following analog separable transfer function, the first separable canonic factor with p = q = -1; H(ju> ,ju> ) = ,-, 7 ■ w-i ■ • r (4.1) 1 2 (1 + joj 1 t 1 ) (1 + joo 2 t 2 ) The logarithm of magnitude-square of this function is given by *> *j _ 10 log |H(co 1 ,oo 2 ) | = -10 logd+u^ t^ 2 ) - 10 log(l+co 2 t 2 ) (4.2) Let us define the following for simplicity of notation. z = 10 log |H(oo 1 ,o) 2 ) I = 20 log |H(o) 1 ,a) 2 ) | (dB) (4.3a) x = 20 log oo-l (4.3b) y = 20 log w 2 (4.3c) Consider the regions in the (gj, ,u) 2 ) plane (as seen in Figure 4.2) as follows: Region I, when oo,x, < 1 and ^ 2 x 2 < 1, then z = db (4.4a) Region II, when co^t -^ > 1 and co 2 t 2 < 1 z = -20 log (w 1 t 1 ) db (4.4b) 88 Region III, when ^^l < 1 and W 2 T 2 > 1 z = -20 log(co 2 T 2 ) db (4.4c) and Region IV, when ou^ > 1 and u t~ > 1 z = = -20 log (oj 1 t 1 ) - 20 log (oj 2 t 2 ) (4.4d) The above equations (4.4a), (4.4b), (4.4c) and (4.4d) can be expressed as Region I, z = (4.5a) Region II, z = -x - 20 log t-, (4.5b) Region III, z = -y - 20 log t 2 (4.5c) Region IV, z = -x - y - 20 log t 2 - 20 log x, (4.5d) These equations describe four intersecting planes in x,y,z space as shown in Figure 4.1. Plane I is the horizontal plane at zero db. Plane II and Plane III have slopes of -20 db/decade passing through the corner lines as indicated in Figure 4.1. These lines are analogous to the break point frequencies of one-dimensional 89 Z = 20 I o g |H( w r w 2 ) I og w slope —-20 d B - 20dB Fig. 4.1. Planar Approximation for 1 H (co, ,u>~) 1+J0J 1 T 1 ) (1 + J00 2 T 2 ) 2 2 20 log I H (oi 1 ,co 2 ) | = -10 log (1+w, T, -10 log (1+co 2 2 t 2 2 ) 90 *2* 1 /^ Y = I O g w IT W I X = °g w. Fig. 4.2. The regions of H(j<VJ0) 2 ) " (1 + J03 1 T 1 ) (1 + J0) 2 T 2 ) 2 2 20 log |H(o) 1 ,o) 2 ) | = -10 log (1+u^ Tj ) 2 2 T 2 2 2 -10 log (1+oj„ t 9 ) in (log u), , log co 9 ) plane 91 log-modulus plots. Several properties of this plane are summarized. Pertinent equations from analytic geometry are given in Chapter III. The actual frequency response of this separable canonic factor with t, = t, = 1 is shown in Figure 4.3a (dB vs. log w,, log oj 2 ) which fits the planar approximation , In Fig. 4.3b/ the contour plot of the frequency response is given, which verifies the regions in log uj-, / log u>- plane shown in Fig. 4.2. The logarithmic plots versus linear frequency are given in Fig. 4.3c and Fig. 4.3d. It is apparent that the approximation regions are well defined when logarithmic frequency is used. The frequency response of this separable canonic factor with t, = 1, t- = 2 is given in Figures 4.4. These results again confirm the logarithmic approximation. (2) Next consider the canonic form H(j Ml ,ja> 2 ) = (1 + j^tjj (1 + ju) 2 x 2 ) (1 + ja> 2 T 3 ) where we have two separable factors of oj 2 . The frequency response of this factor with t, * 1, x 2 = 2 and x^ = 4 is shown in Figures 4.5. As seen from this figure, we have -40 dB/decade slope in the direction of u 2 , as expected from the one-dimensional case. The regions are shown in the contour plot with dotted line (Fig. 4.5c) . 92 In a similar manner we can analyze these regions as follows: The logarithm of magnitude square of this function is given by 10 log |H(u) 1 ,u) 2 ) | 2 = - 10 log(l + oj 1 2 t 1 2 ) - 10 log (1 + w 2 2 x 2 2 ) - 10 log(l + w 2 2 t 2 ) Consider the following regions in (to, , w 2 ) plane (as seen in Fig. 4.5a) as follows: Region I, when u-,t-, << 1 and co 2 x 2 < ^ 2 t 3 << 1 then z = dB. Region II, when u,t, >> 1 and o) 2 t 2 < ^t-d <k 1 z = - 20 log (oj- l t 1 ) dB Region II, when co-,t-, << 1 and oj 2 t 2 < ^2 T 3 <<; 1 z = -20 log (co 2 x 2 ) - 20 log(oj 2 x 3 ) dB Region IV, when m,T, >> 1, ^2 T 3 > C0 2 T 2 >:> ^ z = -20 log u,t, - 20 log ^ 2 T 2 " 20 l0g W 2 T 3 93 Region V, when co^ << 1, ^^^ >> 1, ^ T >>o) 2 t 2 >> 1 Z = -20 logU 2 T 3 ) - 20 log(o, 2 T 2 ) Region VI, when ^^l >> ^' W 2 T 3 > T 2 (J °2 >:> 1 Z = - 20 log a) 1 T t - 20 log oo 2 t 2 - 20 log oj 3 t 3 where Z = 20 log |H(u) 1 ,t 2 ) | . (3) Consider the separable canonic factor now, H(a) 1 ,o) 2 ) = [1 + (jw 1 T 1 ) P ] a We can show that the following regions occur by using the notation of (4.3a) and (4.3b) for simplicity: 94 Fig. 4.3a. The frequency response of . 1 H(jo) 1 ,:o) 2 ) - (i +jWi ) (l+jo) 2 ) 2 z = 2Q log|H(o) 1/ o) 2 ) | = -10 logd+u^ ) (dB vs. log bu,log o> 2 ) -10 log (l+w 2 ) 95 X-SCflLE=l .OOE+00 UNITS INCH T-SCRLE=I . OOE+00 UNITS INCH •ig. 4.3b. The contour plot of H^,^) = (1+j ) ( i +j(0 ) z = 20 log |H| = -10 log (l+o^ 2 ) - 10 log (l+io 2 2 ) 96 Fig. 4.3c. The frequency response of H(joJ 1 ,jco 2 ) (l+jai 1 ) (l+jo) 2 ) z = 20 log|H| = -10 log(l+oj 1 2 ) -10 logd+o^ 2 ) (dB vs. 03, ,0)-) (linear frequency) 97 10. >", a-, 3. I T Fig. 4.3d. The contour plot of H(oo..,go_) = . . , . , . = , . r r 1 2 (1+]U) 1 ) (l+:co 2 ) z = 20 log|H| = -10 logd+u^ 2 ) - 10 log(l+w 2 2 ) (linear frequency) 98 Fig. 4.4a. The frequency response of H(do) 1 ,do) 2 ) (l+joa 1 ) (l+j2co 2 ) 2 2 = 20 log|H| = -10 log(l+oj 1 ) - 10 log(l+4co 2 ) (dB vs. log co, , log co 2 ) 99 X-SCRLE-1 .OOE+00 UNITS INCH. T-SCflLE=l . OOE+00 UNITS INCH. Fig. 4.4b. The contour plot of the frequency response of H(o3 1 w 2 ) = (l+ja) 1 ) (l+j2u) 2 ) z = 20 log|H| = -10 log(l+co 1 2 ) - 10 log(l+4to 2 2 ) in (log to, , log to 2 ) plane 100 Fig. 4.4c. The frequency response of 1 H ( j co , , j co -, ) = ,-, , • - J 1 2' (l+jco,) (1 + j2co 2 ) z = 10 log|H| = -10 logd+co, 2 ) - 10 log(l+4co 2 2 (dB vs. co, ,co_) (linear frequency ccale) 101 Fig. 4.4d. The contour plot of the frequency response of H(jo) 1 ,joj 2 ) = (1+w, ) (l+4oo 2 ) z = 20 log|H| = -10 log(l+oo 1 2 ) - 10 log(l+4aj 2 ) (dB vs. o) 1 ,o) 2 ) (linear frequency scale) 102 ^ ' o g w 2 «H "2 YE nr 12 0)i*Tx i I UH-i -► !og w 1 Fig. 4.5a. The region of H(j(D 1 ,ju) 2 ) = (l+jo) 1 T 1 ) (1+joo t ) (l+jco 2 x 3 ) in (log a), /log co„) plane 103 Fig. 4.5b. The frequency response of 1 H(jo) 1 ,jo) 2 ) - (i+j 0Jl ) (l+j2co 2 ) (l+34w 2 ; z = 20 log|H| = -10 log(l+03 1 2 ) - 10 log(l+4w 2 ) 2 -10 log(l+16u> 2 ) (dB vs. log o) lf log w 2 ) 104 X-SCRLE=1 . OOE+00 UNITS INCH. Y-SCRLE=I . OOE+00 UNITS INCH. Fig. 4.5c. The contour plot of the frequency response of H(ju lf ju 2 ) = ( i+j Ul ) (l+j2oo 2 ) (l+34co 2 ; 2 2 z = 20 log|H| = -10 log(l+a) 1 ) - 10 log (1+4^ ) 2 -10 log(l+16oj 2 ) (log oo-L vs. log co 2 ) 105 Fig. 4.5d. z = The frequency response of H(joJ 1 ,jOi 2 ) = (1+j(Ui) (i +j 2co 2 ) (l+j4oo 2 ) 20 log|H| = -10 log(l+ca 1 2 ) - 10 log(l + 2aj 2 2 ) -10 log(l+16oo 2 2 ) (dB vs. (ji-.,uiy) (linear frequency scale) 106 10, 0. - X- /«u Fig. 4.5e. The contour plot of the frequency response of . 1 H(D(o 1 ,30) 2 J - (1+jc0i ) (l+j2oa 2 ) (l + j4co 2 ) ? 2 z = 20 log|H| = -10 log (1+w^) - 10 log(l + 4aj 1 ) -10 log (l+16oo 2 2 ) (linear frequency scale) 107 Region I, when go, > 1/t, , z = apx + 20ap log t-, dB Region II, when to, < 1/t-, z = dB These separable canonic factors are identical to the canonic factors of a one-dimensional transfer function. They are just an extension to the two-dimensional case. Therefore, they need not be discussed further. B. NON-SEPARABLE FACTORS The following non-separable canonic factors are now considered. (1) [1+ (jaj 1 x 1 ) P (joj 2 T 2 ) q ] a (4.d) (2) [1 + (ju 1 T 1 + joo 2 T 2 ) P ] a (4.e) where p, q and a are integers. (1) First, consider the non-separable canonic factor of the form: H(oJ 1 ,aj 2 ) = [1 + (jaj- L T 1 ) P (joJ 2 T 2 ) ] 108 Region I, when P q 1 1 2 T P T q z = dB Region II, when P q 1 2 T p , q T l T 2 z = apx + aqy + 20ap log t, + 20 ay log t- dB The cornerline between the two regions is given by the equation or p q p q u l w 2 T l T 2 px + qy + 20p log t, + 20q log t 2 = The cornerline is sketched in Figure 4.6. The slope of the plane in region II is given by -4 2 'p + q . Direction numbers of the normal are given by -ap, -aq and -1. 109 »,= « 2 = I (*)(*!* ^ x - lo g I T M±)H^ y -lo g Fig. 4.6. The regions of non-separable case in log to, , log co^ plane 110 Depending upon the sign of p and q and the magnitude of t, and x 2 several different orientations of the corner line are possible. Example (1) Discuss this canonic form in the special case with p = q = 1, a=+l, x, =0.5, and x 2 = 1. The transfer function gets the following form: H(o),,a) 2 ) = [1 + (.5jw,) ( ja) 2 ) ] Take the logarithm of the magnitude square. i 2 | i+2 z = +20 log|l - .5(jj,co 2 | In region I, w,t, <<: ^ an< ^ aj ? T 2 <<: "*" z = dB. In region II, id,t, >> 1 and w 2 t 2 >> 1 z=x + y+20 log (0.5) where x = 20 log co 1 , y = 20 log cu 2 and z = 20 log |H(oj 1 ,u 2 ) The actual frequency response of this factor is shown in Fig. 4.7 including a contour plot. Ill Fig. 4.7a. The frequency response of H( joo 1 , joo 2 ) = 1 + ( ja>,) ( jO. 5o3 2 ) z = 20 log|H| = 20 log (1 - 0.5oj,oj 2 ) (dB vs. log go,, log oo^) 112 X-SCRLE=l;OOE+00 UNITS INCH. T-SCRLE=1..0£iE + 00 units inch. Fig. 4.7b. The contour plot of the frequency response of H(jco L ,ja> 2 ) = 1 + ( jco 1 ) (j0.5oj 2 ) z = 20 log|H| = 20 log (1 - O.Su^u^) in (log co, , log ai 2 ) plane 113 When a = -1 we have a stability problem in this example The denominator of this transfer function becomes zero when or 1 — . 5o)-j 0)2 = o)-| cop — 2 This is the singularity of the transfer function in (w,,^) plane- It can be called a singular line of the system. It is a curve in the (ui,,to 2 ) plane to be compared with a pole in a one-dimensional transfer function. Example #2 ; In this example, we discuss the same canonic form with different coefficients, i.e., p = l, q = -1 , a = +l, x, = 0.5, x 2 = 1 Then the transfer gets the form of * 5uj l +1 H( Ul ,« 2 ) = (1 + ~— > 114 Taking the logarithm of its magnitude squared |H(a), ,u 9 ) = ll + ±\ Z . 5(jd, z = +20 log ll + w 2 . 5a).. In Region I when - << 1 a) 2 z = dB. • 5ai, In Region II, when >> 1 w 9 z = 20x - 20y + 20 log (.5) where x = log u^, y = log u 2 and z = 20 log |h(u>,/u> 2 ) |. The boundary between these regions is given by the equation . 5oi-| = oi~ We also have a stability problem for this canonic factor when u>2 s 0. This is also a singularity in (u), ,w 2 ) plane in the form of a straight line. This type of factor is unstable. (2) Now, consider the non-separable canonic factor H(o),,a> 2 ) = [1 + (ju),w, + ju^oO ] We have two regions again, as follows: Region I, when w,t, + w 2 t 2 < 1 z = dB 115 Region II, when u),x, + w-x- > 1 z = ap log (oa 1 x 1 + aj 2 T 2 ) • The separation line between the two regions is given by the equation: U 1 T 1 + 0) 2 T 2 = 1 or log (10, x, + o) 2 t 2 ) = . (4.6) The locus of points of (4.6) is shown in Figure 4.8, For region II, we have to think of two conditions, as follows: a) When w,x, > oo 2 x 2 z = apx + ap20 log x, b) When oj,t, < ^ 2 x 2 z = aqy + ap20 log x 2 . Thus, region II is split into two areas with the separation line given by W 1 T 1 = W 2 T 2 ' 116 Fig. 4.8. Locus points of logC^T^ + w 2 x 2 ) - 117 Frequency, response of this non-separable transfer function is shown in Fig. 4.9 with p - 1, a = -1, x = x = l. C. CONCLUSIONS In this chapter, we discussed the properties of the separable and non-separable two-dimensional properties of canonic factors in the analog domain (co,,^ plane). We can conclude that these canonic factors can be approximated by planes in the log magnitude versus log co, and log oj^ space. As a design technique we can obtain any given frequency specification by proper combinations of planar approximation of canonic factors. It is obvious that only certain special cases and symmetries can be obtained using the factors dis- cussed in this thesis. For example, a wedge-shaped character- istic in frequency domain could not be synthesized by means of the factors we have considered. The other result of this study is the stability problem which arises in certain transfer functions, that is, the singularity curve or line. For a separable case, we don't have a stability problem; but in non-separable cases we may have stability problems. 118 Fig. 4.9a. The frequency response of J 1 J 2 1 + joo, + jog. z = 20 log|H| = 10 log 1 1+ (oj, +oj 2 ) (dB vs. log to, /log co 2 ) -, 2 119 X-SCRLE=1 . OOE + 00 UNITS INCH T-SCflLE-1 . OOE+00 UNITS INCH Fig. 4.9b. The contour plot of the frequency response of H(ju>, ,jw ) = . —-. z = 20 log|H| = 10 log |l + (u^+u^) 2 ! in (log co, , log co-) plane 120 Fig. 4.9c. The frequency response of 1 H(3oa 1 ,:a) 2 ) = 1 + ju). JO). z = 20 log|H| = 10 log |l + (u^+a^) | (dB vs. a) 1 ,oo 2 ) (linear frequency scale) 121 v - TWO-DIMENSIONAL DIGITAL FILTER DESIGN TECHNIQUE A. INTRODUCTION The design of two and multi-dimensional recursive digital filters is difficult due to the absence of the Fundamental Theorem of Algebra in two or more dimensions [15]. That is, polynomials in two variables may not, in general, be factored into lower order polynomials. For this reason, several design procedures for two-dimensional recursive filters have been proposed in the literature. We then consider the use of the double bilinear z-transfer with the semi-infinite plane approximation in the previous chapters. These are briefly discussed in the following. Shank et al suggested three design techniques in their paper [11] . 1. The use of a series of stability-preserving transforma- tions to produce a stable two-dimensional recursive discrete filters from a stable one-dimensional continuous filter. 2. A space domain design technique in which the unit sample response of the filter being designed is made to approximate a given desired unit sample response. 3. A stabilization technique in which the denominator polynomial is replaced by its double planar least-squares inverse. Planar least-squares can be defined as follows: We have some given array C. We would like to find an array P such that C convolved with P approximates the unit pulse array U. That is, C * P « U 122 where the symbol * denotes two-dimensional convolution. In general it will not be possible to make C * P exactly equal to U. In actuality, C * P will be equal to some other array, G. If we choose P such that the sum of the squares of the elements of U-G is minimized, P is called a planar least square inverse of C. Planar least square inverse of P is called double planar least square inverse [11] . This new polynomial is conjuctured to have no zeros on -2 A , , , , u - {z^z^: | ZjJ ^_, |z 2 | >_ 1} and to have a frequency response which closely approximates that of the original transfer function. The first technique is known as a filter design with rotating axes and these filters are called rotated filters. This design procedure consists in considering a transfer function in the analog frequency domain, rotating the axes, and applying the bilinear transformation to obtain a two- dimensional recursive filter. In general, the starting point is a separable filter and stability of the resulting filter is not guaranteed even if the prototype filter is stable. This technique was used by Costa and Venetsanopoulus [23] to design low-pass filters having circular symmetry, and also they proved that rotated filters are stable if the angle of rotation is between 270° and 360°. A technique similar in character but rotated axes in digital frequency domain ((z,,z 2 ) domain) has been given by Chang and Aggarwal [24]. Hirona and Aggarwal [25] have accomplished the rota- tion through the introduction of rational powers of z^ and 123 z 2 , and the interpolation of input and output sequences. They have proved that, in this technique, the resulting filter is stable if the prototype filter is stable. Further, the distortion, a consequence of the bilinear transformation, is absent in the resulting filter when compared with Shank's et al. 's filter. The authors of [11] were unable to prove that the new denominator polynomial produced by the third technique, which is an extension of the one-dimensional case, was stable in the two-dimensional case. Nonetheless, several hundred filters have been designed using this technique without a counter example appearing [18]; furthermore, Jury, Kolavennu and Anderson [26] have shown it to be valid for certain low order cases. However, Genin and Kamp [27] have recently shown that counter-examples do, in fact, exist. Spectral factorization is another technique which can be used to produce a stable filter from an unstable one. By spectral factorization we mean a complex map that carries a stable rational transfer function into another stable transfer function exhibiting a different frequency response, at the same time maintaining some desirable characteristic [29] . Pendergrass, Mitra and Jury [28] used stability-preserving spectral transformations to modify the frequency response of quarter-plane filters. Mitra and Chakrabarti [29] extended this work. Read and Treitel [30] suggested a two-dimensional Hilbert trans- form to perform the spectral factorization, while Pister 124 [31] and Ekstrom and Woods [32] made use of the two- dimensional cepstrum. Ekstrom and Woods also suggested ways that the factors could be truncated and smoothed to preserve stability. They also showed that the usual quarter-plane filter is not the most general two-dimensional recursive filter. It is instead possible to have a recursive filter which in the limit has a unit sample response with support on an asymmetric half-plane, and also they developed stability conditions for this new class of filters. B. TWO-DIMENSIONAL DIGITAL FILTER DESIGN PROCEDURE BY DOUBLE Z-TRANSFORM 1. Double Z-transform In the one-dimensional recursive digital filter design, the bilinear z-transform is often applied to an analog filter transfer function. We proposed that this one- dimensional design approach can be extended to the two- dimensional case. Before going into the design procedure, two-dimensional recursive filter design via the double bilinear z-transform, it is better to prove that we can use double z-transform for two-dimensional systems. Each of the canonic factors, proposed in Chapter IV, either separable or non-separable, can be converted into discrete equations by the double z-transform, as follows: 2 1 - Z I' 125 2 X " z 2 1 + z 2 The single z-transform is equivalent to trapezoidal integration which can be shown as follows: X(s) = S Y(s) (5.2) .<*> = &J§i. or n Y n ~ y n _ x = / x(t) dt n-1 ?(x + x .) (5.4) 2 n n-1 where T is the sampling period. We write (5.4) in the z transform domain as (1 - z 1 )Y(z) = |d + z X )X(z) Hence X(z) = §( 1 " \[ ) Y(z) (5.5) 1 + z 126 Comparing (5.5) with (5.2), it is seen that 2 1 - z" 1 S = J(i ^- T ) (5.6) 1 + 2 Equation (5.6) is known as a trapezoidal bilinear z -trans form, We can drive the double z-transform in a similar way, as follows: or Then ,2 W(x,y) = 33^ Z(x,y) (5.7a) W(s,,s 2 ) = s i s 2 z ( s ]_' s 2^ (5.7b) / /d Z(x,y) = / /W(.x,y) dx, dy (5.8) by a discrete approximation, the left-hand side of (5.8) becomes (see Fig. 5.1) / /d 2 Z(x,y) = [Z(x n ,y n " z ( x n' y n-l )] " ^V-l'V " Z(x n-l' Y n-l )] -1 -1 . -1 -1, (1 - z x - z 2 + z 1 z 2 )Z(.z 1 ,z 2 ) (5.9) 127 and the right-hand side of (5.8) is approximated as follows: //W(x,y)dxdy = A^£ [w(Xn ,y n ) + W(x n ,y n-1 ) + "(Vl'V + W(x n-l'*n-l )] AXAVr, -1 , -1 , "I -1, TT/ V = — -^- [1 + z, + z 2 + z i z 2 ]W(z 1 ,z 2 ) (5.10) By equation (5.10) and (5.9) we have W(z 1 ,z 2 ) 1 - z" 1 - z" 1 + z~ z 2 Z(z,,z ) . -1 -1 -1 -1 AxAy 1 2 1 + z, + z 2 + z, z 2 2 (5.11) or -1 i - 1 2 1 " Z l 2 X ~ Z 2 = ^ ( T—*^ ) ^ ( TT7^ ) (5 * 12) 1 + z, J 1 + z 2 By comparing (5.1a) and (5.2b) it is seen that 2 1 " z l 1 " Z 2 Q.E.D. Thus, this is called a double z-transform. The double z-transform seems to be a reasonable substitution 128 y(t) x (O — n-1 n ► t (a) One-dimensional case Z(*,y) w(x,y) >W(x,y) (b) Two-dimensional case Fig. 5-1- Two-dimensional extension of the Bilinear Transformation 129 for each canonic factor. However, each factor must be checked for stability, which is discussed in detail later in this chapter. If each transformed factor is BIBO stable the overall transfer function will be BIBO stable. But sometimes the double z-transform may in certain cases lead to unstable solutions [33] . Finally, it is noted that the double z-transform results in a quarter-plane symmetry which is a severe limitation. 2 . Design Procedure The two-dimensional recursive digital filter design investigated in this thesis is the double bilinear z-transform. Proceed as follows: 1. Determine a two-dimensional analog transfer function, A(s 1 ,s 9 ) hi«i'-2> = niTTifr (5 - 12> where A(s,,s 2 ) and B(s,,s 2 ) are mutually prime two-dimensional polynomials and B(s 1 ,s 2 ) satisfies B(s 1 ,s 2 ) f (5.13) V(s 1 ,s 2 ) z r = {( Sl ,s 2 ): Re (s 1 ) >_0,Re(s 2 ) >_ } in such a way that, the frequency response of G(s 1 ,s 2 ) meets the given specifications (for example as discussed in Chapters III and IV) . 130 2. Apply the double bilinear z-transform 2 1 " z l 1 5 1 = A { ^T } (5.14a) 1 + z, 2 1 " z l 5 2 = -(- ij) (5.14b) 1 + z 2 where A is the sampling period, which is assumed the same in both directions. This yields the two-dimensional digital transfer function 2 1 -z, 1-2, C(z,,z ) H( f ( TI ) 'f ( =T )} = T ^V Z 2 ] ~ D(z z ) (5 ' 15) A 1 + z 2 X A 1 + Z2 1 X 2 D( Z;l ,z 2 j where C(z,,z 2 ) and D(z,,z 2 ) are mutually prime two-dimensional polynomials in z, and z 2 . (5.15) will have the desired frequency response with proper sampling time without aliasing effect. Since , 1 - z? 1 si = !< *r 1 A 1 + z. maps the closed unit disc onto the closed right half plane, one would expect (5.15) to satisfy the following sufficient stability conditions [11], [20]. D(z 1 ,z 2 ) ? 131 \AW £ u2 - {(.z v z 2 ): | Zl | >_ 1, z 2 | >_ 1} (5.16) Unfortunately, however, the conditions (5.13) and (5.16) as we will show later in this chapter, are not always satisfied. In summary, we have two conditions to be satisfied for stability. I. B(s lf s 2 ) ? o \/(s 1 ,s 2 ) E r = {(s 1 ,s 2 ): R e s 1 >_ 0,R e s 2 > 0} (5.17) II. D(2 1 ,z 2 ) ? Y(z 1 ,z 2 ) £ u = {(z 1 ,z 2 ): Iz-J > l/|z 2 | >_ 1} (5.18) where V stands for "for each", and u is the closed unit bidisc C. DIGITAL CANONIC FACTORS AND STABILITY CONSIDERATIONS 1 . Introduction Certain properties of two variable polynomials and rational functions which will be needed later will now be discussed. It is well known that a two variable polynomial is not, in general, factorable into first-order polynomials; but, a two variable polynomial can be factored into irreduci- ble factors which are themselves two variable polynomials but which cannot be further factored. Of course a given 132 polynomial may itself be irreducible. These irreducible polynomials are unique up to multiplicative constants [33] . Two polynomials which have no irreducible factors in common are said to be mutually prime. Consider two variable rational functions, in general, P (z, ,z 9 G(z, ,z ) = l'"2' Q(z 1 ,z 2 ) where P(z^,z 2 ) and Q(z,,z 2 ) are mutually prime. Using terminology of [34], a point (zj,zl) making Q(z|,z 2 ) = but P(z',z') ^ is called a pole or a nonessential singu- larity of the first kind (such a point is analogous to a pole in the one variable case). A point (z£,z£) making P(zV,z~) = Q(zV,z 2 ) = is called a nonessential singu- larity of the second kind (such a point has no analogue in the one variable case). Clearly, if (zi,zl) is a pole, G(z,',z') = °°, but if (zV,z") is a nonessential singularity of the second kind, then G(zV,z") is undefined. 2. Separable Factors The most common canonic separable factor in Chapter IV is H(<jj,,u) ) — "7t — —j, \ 1 1 . a r (5.19a) 1 2 (1 + Joj, x,) (1 + 3oj 2 t 2 ) assuming x = x_ = 1. Substituting jw-L - S 1 , 133 ju> 2 = S 2 We get H(s i' s 2» " (TTT[flTTip < 5 - 19b > Substituting equations (5.1a and 5.1b) into (5.19b), performs the double z-transformation. Thus, 8(1 + z" 1 ) (1 + z" 1 ) T(z r z 2 ) = i- £_ (5.20) (a + z 1 ) (a + z 2 ) where A 2 (A -2) 2 a = a + 2 a - 2 A = sampling period This separable canonic form obviously satisfies the condition of (5.17). We should be careful about the second condition of (5.18) . When z, = z 2 = -1, we have a nonessential singularity of the second kind if a = 1, but a = 1 never occurs. Therefore, there is no nonessential singularity of the second kind. Thus, this separable case is stable by satisfying 134 two essential necessities (5.17) and (5.18). The frequency response of this separable case is plotted in Fig. 5.2, with sampling time A = 0.08 which fits that of the analog frequency response. The double z-transform for equation (5.19a) is A 2 (l + z 1 1 ) (l + z^ 1 ) T(Z 1 /Z 2 ) = Ta-2 T;l ) (A-2t 2 ) A+2 T;l ~ aT2~^ ~ ( A^2T7 +Z 1 )( A^2TJ +Z 2 > (5.21) and this general case for the two-dimensional separable factor is stable by satisfying the stability necessities of (5.17) and (5.18) The second analog canonic factor in Chapter IV is H(j Ul rJ0. 2 ) (1 + jo) 1 T 1 ) (1 + JW 2 T 2 ) (1 + j w 2 T 3 ) ' When we perform the double bilinear transform for this, we get 3 -1 -12 A J (1 + z/) (1 + z 9 x r T(z n ,z„) = 1' 2' A+2t, , A+2T- , A+2t-. , (A-2x 1 )(A-2x 2 )(A-2x 3 )[^ + z- i ][ I=2 ^ + z 2 H I=n i+«; ] (5.22) The frequency response of (5.22) for t, = 1, t 2 = 2, t^ = 4 is shown in Fig. 5.3 which almost fits the analog case. 135 Fig. 5.2a. The frequency response of T(z 1 ,z 2 ) = where 3 = e(l+z 1 " 1 ) (l+z 2 1 ) (a+z^ ) (a+z~ ) A 2 A+2 ~ , a = (A-2)" (dB vs. log a), , log u}~) A-2 = 0.08 136 Fig. 5.2b. The frequency response of T(z 1 ,z 2 ) = 6(l+z 1 ~ 1 ) (l+z 2 1 ) (a+z 1 " 1 ) (a+z 2 -1 ) 3 = A' (A-2) 2 ' a = |i| , A = 0.08 (dB vs. a). ,0)-) (linear frequency scale) 137 IT). 0. -, *4 10.-, Fig. 5.2c. The contour plot of the frequency response of 3(l+z 1 " 1 ) (l+z 2 _1 ) T(z x ,z 2 ) 6 = 2 ' (a+z 1 ~ ) (a+z 2 ) a = f±§, A= O.OB where (A-2) (linear frequency scale) 138 We can summarize the m xn analog separable factor as M N H(ju> ,ju, ) = [ IT , .-*• ] [ n , , / ] (5.23) After performing double z-transform, (5.23) becomes M -1 N -1 M a (1 + z, ) N A (1 + z,) T(z,,z ) = [ n — =i ][ n — | ] (A-2t.) (-— ii+zT 1 ) 3 L (A-2x.) (_-^J-+z: 1 ) 1 A-2x -1 3 A-2t ■ 2. 15.24) Equations (5.23) and (5.24) satisfy the conditions of (5.17) and (5.18). As a result, a separable filter results in separable filters in digital domain and they are stable. 3. Non-separable Factors The first non-separable factor in Chapter IV is 4.d with a = -1, p = q = t, = t 2 = 1 as given below 1 H ( joj, t j oj ) - j— - ■ - 1 1 2 1 + 3oj-i](jo~ (5.25a) or «<V S 2> ■ TTT^ (5 - 25b: 139 Fig. 5.3a. The frequency response A 3 (l+z 1 " 1 ) d+z^ 1 ) T(z 1/ z 2 ) = (A-2) (A-4) (A-S, (f±| .z,- 1 ) (f±| .zf 1 ) <£±f +z 2 -h where A = 0.0 8 (dB vs. log co, , log u^) 140 Fig. 5.3b. The frequency response of T( 2l ,z 2 ) = A 3 (l+z 1 " 1 ) U+z^ 1 ) (A-2) (A-4) (A-8) (|±| +z 1 " 1 ) (|±i +z x X ) (f±f + z 2 1 ) where A = 0.08 (dB vs. u.,ii)J 141 !S.u Fig. 5.3c. The contour plot of the frequency response of T(z 1 ,z 2 ) = A 3 (l + z 1 " 1 ) (l+z^ 1 ) (A-2) (A-4) (A-8) (|±| +z 1 - 1 ) (|±J +a 1 " 1 ) <£±§ -hz^ 1 ) where A = 0.0 8 (linear frequency scale) 142 and by the double 2 -trans form, we get A 2 (l +Z, 1 ) (1 + Z' 1 ) T(Z 1' Z 2 ) = — 2 =T TT T Tt-ZT (5 ' 26) (A +4) (A -4) (z L + z l ) + (A +4)z 1 z X In the analog domain (5.25b) does not satisfy the condition (5.17). We have a singularity at Ul u) 2 = 1 For the second condition, we have a nonessential singularity of the second kind at z, = -1, z_ = -1. This does not tell us whether it is stable or not, When we check the frequency response, it is seen that this factor is unstable. The second non-separable canonic factor can be derived from (4.d) with a = q = -1 and T-,=T 2 = las follows. H(ju)i ,jui ) = 3 (5.27a) 1 + 2 tu. or H < s i< s 2> - — 4r = sr^r < 5 - 27b > 1 + IT s 2 143 Then, (5.27b) becomes after performing the double z-transform, as follows: (l + zT 1 ) (1 -z' 1 ) T(z ,z ) = _ 2 ( 5.28) 2U-« 1 r « 2 I ) This canonic factor does not satisfy either stability condition. Hence it is unstable. The other series of canonic factors can be derived from (4.e) . For example, with a = -1, t 1 = t- « p » 1, we get HtjoiwjwJ = -, , - 1 , ■ (5.29a; -L ^ 1 + 3 co, +3 u^ or H(S 1' S 2» = l + s^ + s 2 ( 5 - 29b> or in the digital domain A(l + z" 1 ) (1 + z" 1 ) T(z,,z ) = — —, _ (5.30) (A+4) +A(z x x + z 2 ) + (A-4)z 1 z 2 This canonic factor satisfies the first condition but it has a nonessential singularity of the second kind, at z, = z~ = -1 and so we should check the frequency 144 response. We see that there is a singularity at (u) 1 ,o) 2 ) = (it, it). Thus this non-separable canonic factor is unstable. 4 . Conclusions All of the canonic factors discussed in this chapter are collected in Table I. By comparing the results shown in this table, we can say that separable stable canonic factors result in stable separable filters. As seen from Table I, the separable transfer functions 1, 2 and 3 which start with stable separable ones in the two-dimensional continuous frequency domain end up as stable separable trans- fer functions in the two-dimensional discrete domain. This is also true for the mxn separable canonic factor 4. Non-separable canonic factors have very different characteristics. The unstable non-separable case ends up with a non-separable case in the discrete domain, which is also unstable. Although we could not prove this result generally, we also could not find any counter example to exist. The transfer functions of 5, and 6 start in the unstable case, and end up with the unstable case in the discrete domain. The transfer functions of 7 starts with the stable case in the continuous frequency domain, ends up with the unstable case in the discrete domain. Thus we have to be very careful with the double z-transform since it does not yield stable results every time. It may yield unstable results even if it starts from the stable case. 145 UJ _J CO < >- CO < co C\J M CM <f) J) LU GQ < ea H LU CD < CO I m w I NT I 4- CM LU co < uo J-flL If << V I i X I 00 + rO LU CD < CO r r ^-^ \ V ♦ *" s t-r 1 < ~ • *~<+ «- ^ .**• C«i *< II rJ WO *=* (/T + i ! i «■ 3iaVd3d 3S 146 LU CD < LU CQ < UJ —I cn < ui i « 4- i — m +• -< fW J" «i 1 H (N *— «. \ 1 m — hi -.jj • 1 ^ ^ r« i . i fV — + _ N '-" / I .4 > * "7. » • -. I* ru + + *W iT) CD + v/T 4» 1^ 319VyVd3SNn 147 VI. CONCLUSIONS A technique for the use of semi-infinite planes to approximate the log magnitude versus log frequency charac- teristics of two-dimensional filters is presented. This technique is applied to quarter plane filters and works well for separable transfer functions. Other non- separable canonic (basic) transfer functions are also studied in terms of their planar approximation. The following separable and nonseparable canonic factors have been studied. (1) (1 + jo) 1 x 1 ) p (l + jw 2 T 2 ) q (2) (1 + oo 1 x 1 ) P (l + jco 2 T 2 ) q (l + joj 2 T 3 ) q ' (3) [1 + (jco 1 t 1 ) P ] a q,a (4) [1 + (joo 1 T 1 )^(l + ja> 2 T 2 )*] p,a (5) [1 + (ja) 1 T 1 + jco 2 T 2 )^] The properties and stability of these canonic factors are investigated first in the continuous variable domain and then after the double bilinear z-transform is applied to extend them to two-dimensional recursive digital filters These are all summarized in Table I in Chapter V. An important result is that although the double bilinear z-transform works well for separable factors, it does not 148 always lead to a stable digital realization starting with a stable analog transfer function. The following are suggested for further research: 1. Other types of canonic factors need to be considered. For example, the general form: H(j» lf jw 2 ) = [1 + (ja3 1 T 1 ) P +(ja J2 T 2 ) q +(jco 1 T 3 ) r (jw 2 T 4 ) t ] a 2. The work should be extended to non-causal (half- plane) filters. 3. Transfer functions which have circular symmetry in the discrete domain need to be investigated. 4. Another approach is to consider half-plane recursive algorithms, by working backwards from the discrete domain (z, ,zj to continuous domain (s,,s 2 ). 5. Stability of canonic factors needs to be investigated in detail . 6. The use of predistortion should also be considered. 149 LIST OF REFERENCES 1. Oppenheim, A.V. and Schafer, R.W., Digital Signal Processing , Prentice Hall, Englewood Cliffs, New Jersey, 1975. 2. Rabiner, L. and Gold, B. , Theory and Application of Digital Signal Processing , Prentice Hall, Englewood Cliffs, New Jersey, 1975. 3. Rabiner, L.R., et al., "Terminology in Digital Signal Processing," IEEE Trans. Audio Electroacoust . , vol. AU-20, p. 322-334, December 1972. 4. Stockham, T.G., "High Speed Convolution and Correlation," AFIPS Conference Proceedings, 28, pp. 229-233, 1966. 5. Hall, E.L., "A Comparison of Computations for Spatial Frequency Filtering," Proc. IEEE (Special Issue on Digital Picture Processing), Vol. 60, pp. 887-891, July 1972. 6. Huang, T.S., "Two-Dimensional Windows," IEEE Trans. Audio Electroacoust. (Corresp.), Vol. AU-20, pp. 88-89, March 1972. 7. Hu, J.V. and Rabiner, L.R., "Design Techniques for Two- Dimensional Digital Filters," IEEE Trans. Audio Electro- acoust. (Special Issue on Digital Filtering), Vol. AU-20, pp. 249-257, October 1972. 8. Arazi, E., "Two-Dimensional Digital Processing of One- Dimensional Signal," IEEE Trans. Acoustics, Speech, and Signal Processing, Vol. ASSP-22, pp. 81-86, April 1974. 9. Kamp, Y. and Thiran, J. P., "Maximally Flat Nonrecursive Two-dimensional Digital Filters," IEEE Trans. Circuits and Systems , Vol. CAS-21, pp. 437-449, May 1974. 10. Rabiner, L.R., McClellan, J.H., and Parks, T.W. , "FIR Digital Filter Design Techniques Using Weighted Chebyshev Approximation," Proc. IEEE , Vol. 63, pp. 595-609, April 1975. 11. Shanks, J.L., Treitel, S. and Justice, J.H., "Stability and Synthesis of Two-Dimensional Recursive Filters," IEEE Trans. Audio Electroacoust. , Vol AU-20, pp. 115-128, June 1972. 150 12. Mitra, S.K., Sagar, A.D., and Pendergrass, N.A., "Realizations of Two- Dimensional Recursive Digital Filters," IEEE Trans. Circuits and Systems , Vol. CAS-22, pp. 177-184, March 1975. 13. Green, B.W., and O'Handley, A.D., "Recent Developments in Digital Image Processing at the Image Processing Laboratory at the Jet Propulsion Laboratory," Proc. IEEE , Vol. 60, pp. 821-827, July 1972. 14. Treitel, S., Shanks, J.L. , and Frasier, C.W., "Some Aspects of Far Filtering," Geophysics , Vol. 32, pp. 789-800, October 1967. 15. Huang, T.S., Schreiber, F.W., and Tretiak, O.J., "Image Processing," Proc. IEEE , Vol. 59, pp. 1586-1609, November 1971. 16. Hunt, B.R., "Digital Image Processing," Proc. IEEE , Vol. 63, pp. 693-708, April 1975. 17. Wood, L.C. and Treitel, S., "Seismic Signal Processing," Proc. IEEE , Vol. 63, pp. 649-661, April 1975. 18. Mersereau, Russell M. and Dudgeon, D.E., "Two- Dimensional Digital Filtering," Proc. IEEE , Vol. 63, pp. 610-623, April 1975. 19. Chong, H. and Aggarwal , J.K., "Fundamentals of Two- Dimensional Recursive Filtering," Digital Signal Processing , (Edited by Aggarwal, J.K.), pp. 211-260, Pt. Lobos Press, CA. , 1979. 20. Huang, T.S., "Stability of Two-Dimensional Recursive Filters," IEEE Trans. Audio Electroacoust . , Vol. AU-20, pp. 158-163, June 1972. 21. Ansell, H.G., "On Certain Two-Variable Generalizations of Circuit Theory, with Applications to Network of Transmission Lines and Lumped Reactances," IEEE Trans. Circuit Theory , VoL CT-11, pp. 214-223, June 1964. 22. Anderson, B.D.O. and Jury, E.I., "Stability Test for Two-Dimensional Recursive Filters," IEEE Trans. Audio Electroacoust . , Vol. AU-21, pp. 366-372, August 1973. 23. Costa, J.M., and Venetsanopoulos , A., "Design of Circularly Symmetric Two-Dimensional Recursive Filters, 1 IEEE Trans. Acoustics, Speech and Signal Processing , Vol. ASSP-22, pp. 432-443, December 1974. 151 24. Chang, H. and Aggarwal, J.K., "Design of Two-Dimensional Recursive Filters by Interpolation/' IEEE Trans. Circuits and Systems , Vol. CAS-24, pp. 281-291, June 1977. 25. Hirono, K. and Aggarwal, J.K., "Design of Two- Dimensional Recursive Digital Filters," IEEE Trans. Circuits and Systems , Vol. CAS-25, pp. 1066-1076, December 1978. 26. Jury, I., Kolavennu, V.R., and Anderson, B.D.O., "Stabilization of Certain Two-Dimensional Recursive Digital Filters," Proc. IEEE , Vol. 65, pp. 887-891, June 1977. 27. Genin, Y.V. and Kamp , Y.G., "Two-Dimensional Stability and Orthogonal Polynomials on the Hypercircle, " Proc. IEEE , Vol. 65, pp. 873-881, June 1977. 28. Pendergrass, N.A., Mitra, S.K., and Jury, E.I., "Spectral Transformations for Two-Dimensional Digital Filters", IEEE Trans. Circuits and Systems , Vol. CAS-23, pp. 26- 35, January 197 6. 29. Chakrabarti, S. and Mitra, S.K., "Design of Two-Dimensional Digital Filters via Spectral Transformations", Proc. IEEE , Vol. 65, pp. 905-914, June 1977. 30. Read, Randol. R. and Treitel , S., "The Stabilization of Two-Dimensional Recursive Filters via the Discrete Hilbert Transform" , IEEE Trans. Geoscience Elect . Vol. GE-11, pp. 153-160. July 1973. 31. Pister, P., "Stability Criterion for Recursive Filters", IBM Journal of Research and Development , Vol. 18, pp. 59-71, January 1974. 32. Ekstrom, P.M. and Woods, J.W., "Two-Dimensional Spectral Factorization with Applications in Recursive Digital Filtering", IEEE Trans. Acoustics, Speech and Signal Processing , Vol. ASSP-24, Pp. 115-128, April 1976. 33. Goodman, D., "Some Difficulties with the Double Bilinear Transformation in Two-Dimensional Recursive Filter Design", Proc. IEEE , Vol. 66, pp. 796-797, July 1978. 34. Goodman, D., "Some Stability Properties of Two-Dimensional Linear Swift-Invariant Digital Filters", IEEE Trans . Circuits and Systems , Vol. CAS-24, pp. 201-208, April 1977. 152 INITIAL DISTRIBUTION LIST No. Copies 1. Defense Documentation Center 2 Cameron Station Alexandria, Virginia 22314 2. Library, Code 014 2 2 Naval Postgraduate School Monterey, California 93940 3. Department Chairman, Code 62 1 Department of Electrical Engineering Naval Postgraduate School Monterey, California 93940 4. Professor S. R. Parker, Code 62Px 10 Department of Electrical Engineering Naval Postgraduate School Monterey, California 93940 5. Deniz Kuvvetleri K-ligi 1 Personel Sube Bsk-ligi ANKARA/TURKEY 6. Deniz Harb Okulu K-ligi 1 Heybeliada ISTANBUL/TURKEY 7. Coskun Zumrutkaya 3 Yuzbasiler cad. No=20 Degirmendere KOCAELI/TURKEY 153 I<2 AU^Qo 25787. Thesis Z8 Zumrutkaya ] 8381 6 Geometric design tech- nique for two-dimensional analog and digital filters. c.l l,c A^jdO - 25 7 8 7^ Thesis Z8 Zumrutkaya c.l Geometric design tech- nique for two-dimensional analog and digital filters 183816 thesZ8 Geometric design technique for two-dimen 3 2768 001 07029 5 DUDLEY KNOX LIBRARY