A Comparative Study of Stability Testing Approaches of Two-Dimensional Recursive Digital Filters

: There are many problems in science and engineering whose solution is applied in the design of Multi-Dimensional (MD) digital filters. Digital filtering finds an important position in the field of digital signal and image processing. (cid:1) Recently there had been a great deal of interest in the design and stability analysis of Two-Dimensional (2-D) recursive digital filters. The design techniques for stable One Dimensional (1-D) digital filters are relatively well developed; but their extension to 2-D had been tried with difficulties. The stability problem of 2D filters had been the subject of intensive research during the recent decades. The trend is towards finding sufficient conditions, because necessary and sufficient conditions are difficult to find or to test. In this study, we make a detailed comparison of the various approaches available for testing the stability of 2-D recursive digital filters. Though researchers had done a lot on the stability testing methods, a comparative study by listing the computational complexities of various methods had not been done so far. This study would help the researchers to choose proper stability testing methods based on the computational complexities.


INTRODUCTION
The growing interest for the design of Two-Dimensional (2-D) filters is due to a variety of applications in fields such as digital image processing, medical data processing, artificial vision, radar and sonar data processing, remote sensing, pattern recognition, numerical stereoscopy, astronomy, biomedical engineering, biochemistry and robotics [1] . For a given response characteristic, recursive filters have less hardware requirements and so, wherever linear phase is not a requirement, recursive filters are preferred over the nonrecursive filters. Recursive filters have lower sidelobes in the stopband than nonrecursive filters with the same number of parameters. Despite these merits, the problem of stability is associated with the design of recursive digital filters. In the design of 2-D filters, the normal practice is to use first and/or second order sections and to cascade these realizations to obtain any higher order filter. The advantages of a cascade realization are reduction of quantization noise and sensitivity to coefficient truncation [2] . Even while designing a second order filter section from a given magnitude or phase specification using lp design technique [3] , one has to check at the end of each iteration whether the resulting filter is stable or not.
A system is considered stable in the Bounded-Input Bounded-Output (BIBO) sense if and only if a bounded input always leads to a bounded output [4,10] . Stability is often a desirable constraint to impose, since an unstable system can generate an unbounded output, which can cause system overload or other difficulties. If a filter is unstable, any noise, including round-off errors to computation will propagate to the output and be amplified. In the spatial domain, the necessary and sufficient condition for a 2-D Linear Shift Invariant (LSI) system (or equivalently, the filter) to be stable is that its impulse response h ij is absolutely summable [4,10,14] . Mathematically, the stability of a 2-D causal recursive filter is: Alternatively, in the frequency domain, the stability of 2-D recursive filter is determined by the coefficients of the denominator polynomial B(z 1 , z 2 ) of the transfer function H(z 1 , z 2 ) described as follows [5] : The stability of the designed digital filters is therefore essential for their practical implementation. However, most of the existing filter design algorithms may result in an unstable filter [1,3] . In 1-D recursive filters, determination of stability is quite straightforward. It is only necessary to locate a finite set of roots in the z-plane. As for 1-D, any higher-order polynomial can be factored as a product of lower-order polynomials, which is not feasible in the case of 2-D. This set back that occurs in 2-D has a major impact on many results in signal processing [4] . Therefore testing of 2-D recursive filters for stability is extremely cumbersome. However, testing of recursive digital filters for stability is a very important at the design stage. Numerous approaches can be found in the literature for testing the 2-D recursive digital filters for stability [5,[6][7][8][9]11,12,14] . In this study, we make a comparison of the various approaches for testing the stability of 2-D recursive digital filters which are devoid of Nonessential Singularities of the Second Kind (NSSK) [4] .
The organization of this study is as follows. We first present the existing Theorems on stability and approaches for testing the stability of 2-D recursive digital filters. Then we briefly review the various approaches for testing the stability conditions of 2-D recursive filters. A summarized table with the computational complexities of the various stability testing methods is prepared.
Brief review of stability theorems and stability testing approaches on 2-D recursive digital filters: In this research, we review the various stability testing Theorems on 2-D recursive digital filter and approaches for testing its stability.
A 2-D causal recursive digital filter is characterized by the two-dimensional z-transform function as in (2). In this, a ij and b ij are real constants. Without loss of generality, we can assume that b 00 >0. The variables z 1 and z 2 are defined as , where s 1 and s2 are horizontal and vertical spatial complex frequency variables and T 1 and T 2 are the constants representing the sampling period. In (2), we assume A(z 1 , z 2 ) and B(z 1 , z 2 ) are mutually prime and also the transfer function H(z 1 , z 2 ) has no nonessential singularities of the second kind. A variety of Theorems relate the stability of a recursive filter to the zero region of the polynomial B(z 1 , z 2 ). Most of these Theorems have two parts in general [2] . The first step usually consists of checking the zero distribution of a onedimensional polynomial. The second step consists of determining whether B(z 1 , z 2 ) is zero free on some higher dimensional surface. In certain Theorems, this surface has two Euclidean dimensions, three in some and four in other Theorems [2,4] . The popular 2-D recursive filter stability Theorems available in the literature are mainly due to Shanks and Huang [4,10] .
Shanks' theorem is historically one of the first Theorems developed, which is conceptually very simple and has led to other stability Theorems. According to this theorem: This means that if there are any values (real or complex) of z 1 and z 2 for which B(z 1 , z 2 ) is zero and for which z 1 and z 2 are simultaneously less than or equal to one in magnitude, then the filter described by the transfer function H(z 1 , z 2 ) will be unstable, or else, the filter will be stable. The condition in (3) suggests a procedure in which B(z 1 , z 2 ) is evaluated in the 4-D space. This stability condition implies that there are no pole surfaces for any 1 z 1 The search in the 4-D space is of course a tremendous work and hence this procedure cannot be used in practice. In this, one way to test the stability is to map the unit disk d 1 (z 1 , |z 1 | 1) in the z 1 plane into the z 2 plane by the implicit mapping relation B(z 1 , z 2 ) = 0. The filter is stable if and only if the image of d 1 in the z 2 plane does not overlap the unit disk d 2 (z 2 , |z 2 | 1). Shanks' Theorem is simplified considerably by [12] . According to Huang's Theorem [11] , the system (filter) described by the transfer function H(z 1 , z 2 ) as in (2) is stable if and only if: (Condition 1): (Condition 2) [1,2,4,5,7,17] : Condition-1 in (4) is relatively easy to check using any 1-D stability test [5] . Condition-2 in (5) is more difficult since it includes two variables.
In general, the 2-D recursive filter stability conditions can be tested with algebraic and mapping methods [4,[5][6][7][8][9]11,12] . Theoretically, only algebraic methods can exactly determine the stability conditions in (5) in a finite number of steps, while mapping methods lead to approximate results and need an infinite number of steps to acquire that reliability. However, there are some difficulties which limit applications of algebraic methods. In algebraic methods, the stability tests rely on certain algebraic properties of complex polynomials. Huang developed the first algebraic procedure [12] to test the stability conditions based on double bilinear transformation and Ansell method. In a while, a number of other algebraic procedures were developed. Maria and Fahmy used the table form to test the stability of two-dimensional filters. The positivity of self-inverse polynomial is checked for the stability. Two schemes employed for checking the positivity of such polynomials are Sturm's and Cohn's methods. Anderson and Jury came up with a new algebraic procedure [11] for testing the stability conditions. This method does not involve the use of bilinear transformations. In this method, the 2-D polynomial is considered as a 1-D polynomial and then Schur-Cohn method is employed. Siljak [16] has modified the method given by [11] proving that only one self inverse polynomial need to be tested for positivity in the course of testing condition instead of 'q' such polynomials as required in Anderson and Jury method.
In [5] , an alternative necessary and sufficient condition for the stability of a class of 2-D recursive digital filters is discussed. This method makes use of Schur-Cohn Hermitian matrix formulation and bilinear transformation in a single variable. A new necessary stability conditions for 2-D recursive digital filter are given in [6] . It has been shown here that the inversion of the necessary conditions gives sufficient conditions. This approach uses the minimal delay property of 2-D polynomials. In [15] , an improved stability test algorithm for two-dimensional digital filters is discussed. It has been shown the tabular algorithm given by Jury for obtaining inner determinants with numerical entries is also valid for the case of polynomial entries. A parametric algorithm is suggested here for the construction of stability array for 2-D systems that avoids the exponentially growth in the order of polynomial entries.
A new form for the Jury's tabular test was proposed by [13] . It is shown that it is possible to telescope the last polynomial of the table by interpolation and circumvent the construction of the 2-D table. In [8] , a computationally simpler algebraic procedure is proposed. The stability of a two or multidimensional digital filter can be analyzed using the simple and straightforward approach of Routh's stability test. This method involves testing one necessary condition and the sufficient condition involves testing four polynomials using Routh test for 2-D recursive filters.
Most of the algebraic methods need to deal with polynomial entries in their procedures. The orders of those polynomials increase very rapidly as the order of the system under test increases. In [18] , a modified polynomial array is developed for testing the stability conditions, which eliminates redundant information introduced in the ordinary table technique and leads to a greater reduction in the order of polynomial entries of the array.

Comparison of stability theorems and stability testing approaches on 2-D recursive digital filters:
In this research, we compare the stability Theorems and various approaches to test the stability conditions of 2-D recursive filters.
According to Shanks' Theorem, stability implies that the denominator polynomial of the transfer function,  [4] . From the search perspective, this Theorem is considerably simpler than Huang's Theorem which requires a 3-D search. The Huang's Theorem which has been restated by Anderson and Jury, says that the system is stable if and only if conditions-1 and 2 are satisfied. All these tests rely on certain algebraic properties of 1-D complex polynomial and therefore they will be called algebraic methods. In 2-D algebraic testing methods, the 2-D denominator polynomial B(z 1 , z 2 ) is viewed as a 1-D polynomial with respect to one variable say, z 2 and the complex coefficients of the polynomial are themselves 1-D real polynomials with respect to the other variable z 1 . Then applying 1-D algebraic tests to this case result in a sequence of real polynomials that contain the stability information. This sequence of real polynomials must be checked as a function of z 1 to determine the system stability.
Huang developed the first algebraic procedure for testing stability conditions. Soon after this, a number of other algebraic procedures were developed. Testing the stability of a two-dimensional recursive filter using Huang Theorem is much simpler than using Shanks' original Theorem. However, this procedure too is still infinite. Huang Theorem can be reduced to a result due to Ansell [12] . Ansell's result enables us to test the stability in a finite number of steps, which can still be very tedious. The procedure of Huang with Ansell's result though finite, requires the application of two bilinear transformations. In essence, Ansell's main contribution is to couple the use of a Hermite test for checking stability with a series of Sturm tests. This test involves the construction of a Schur-Cohn matrix and checking for positivity on the unit circle of a set of selfinverse polynomials, but the computation of coefficients of the determinant of a polynomial matrix is very difficult. On the other hand, the 1-D s-plane procedure of Routh can be modified so that they give the same information as the Hermite test in a much more programmable manner. In this procedure a series of polynomials are iteratively computed and tested for positivity. The disadvantage of this approach being the order of the polynomials to be checked for positivity doubles for each row so the last polynomial has order approximately 2 M M, where M is the order of the polynomial. Hence, even a tenth order filter will present a substantial computational burden. In order to efficiently implement the bilinear transformation, Huang method requires O(M 3 ) calculations. It is worthwhile to note that methods such as Maria and Fahmy anderson and Jury and Siljak and Bose avoid the bilinear transformation and they seem to be more straightforward way of attacking the problem.
Maria and Fahmy method tests the stability conditions by employing the Marden-Jury table. This table performs the same functions for determining a complex polynomial's root distribution with respect to the unit circle as Routh array does for the left-half plane. In fact, an M th order complex polynomial has all its zeros outside the unit circle if and only if the first column of the table is positive. In this method, two schemes can be employed for checking the positivity of self inversive polynomials. They are Sturm's and Cohn's methods and they both require approximately O(N 2 ) operations where N is the order of the selfinverse polynomials. In order to carry out Maria and Fahmy test completely, the positivity of M polynomials of order M to 2 M-1 M must be checked. Therefore, it takes less than O(M 2 4 M ) real calculations to guarantee that a given M th order filter is stable. This means a 16×16 2-D causal recursive digital filter stability test requires over a trillion calculations to guarantee stability. It is slightly easier to implement as no bilinear transform is needed.
The procedure presented by [11] is much simpler than by [12,14] . The procedure by Anderson and Jury, like Huang procedure is finite. It requires no bilinear transformation and replaces the Hermite test component of the main part of Ansell's procedure by Schur-Cohn matrix. Then it allows either a series of Sturm tests, or, what turns out to be equivalent, a series of tests for establishing the root distribution polynomial. The checking of stability conditions falls into two distinct phases. First, we apply a Schur-Cohn test, after a fashion. Secondly, we check the positivity of a number of polynomials on |z 1 |. Anderson and Jury method is really to be compared against that of Huang. In general, these methods involve the same sort of calculations, save that we avoid the bilinear transformation component of Huang method [12] . Undoubtedly, this represents a substantial computational load for any but the simplest two-variable polynomials. On the other hand, as soon as the two variable polynomials under test become at all complex, presumably powerful computers will be used to do the checking and it might well prove the case that programming considerations determine which is the better method.
In the method by Reddy et al. [5] , the computational complexity is further reduced. It is because this method makes use of both Schur-Cohn Hermitian matrix formulation and bilinear transformation in a single variable. Compared to Huang Anderson and Jury and Siljak methods this method is found to be computationally superior in all cases except a special case when no Sturm test is performed. The minimal delay property that gives the useful necessary conditions for 1-D polynomial, when extended to 2-D gives the necessary conditions for testing the stability of 2-D polynomial [6] . In a method proposed by [7] , the results of Barnett have been used to test the stability conditions. It is shown here that for the computation of the polynomial matrix, C, we find that the total cost involved is O(3n 2 ) multiplications and additions. On the other hand, the method of Schur-Cohn includes O(n 3 ) multiplications and O(n 3 /2) additions. From this it is evident that the method given in [7] is better than the classical Schur-Cohn method with respect to computational complexity.
In the algorithm proposed by [15] , there is a greater reduction in the polynomial entries of the array from 1 n 1 2 n 2 − to n 1 n 2 . The proposed algorithm is used to construct a polynomial array, whose degrees are comparatively lower than the established methods and hence lower in complexity. The telepolation of the tabular test by [13] replaces the table by n 1 n 2 +1 stability tests of 1-D polynomial of degree n 1 and n 2 of certain form. The resulting new 2-D stability testing reduces the count of operations to O(N 4 ). A new algebraic test procedure reported in [8] relies on Routh criterion for testing the stability conditions. It has been shown that the number of necessary condition to be tested for the stability of 2-D filter is only one and the sufficient condition involves 4 different cases.
Most of the stability testing methods need to deal with polynomial entries in their procedures. The orders of those polynomials increase very rapidly as the order of the system under test increases. For example, to test a 10th order filter, the table technique involves computing and positivity testing of a sequence of polynomials up to 5120th order. Many proposed methods, however, are too complicated to program. Furthermore, both round-off error and coefficient quantization effects may rise to affect the accuracy of the result because of such a large number of calculations. In [18] , a modified polynomial array is developed which leads to a great reduction in the order of the polynomial entries from 1 n 1 2 n 2 − to n 1 n 2 . Bose test requires approximately O(6M 4 ) operations.
Anderson and Jury stability test with Siljak simplification is more efficient algebraic stability test for 2-D causal recursive digital filters. However, it is stressed that this and any other algebraic methods are infeasible for all but low-order filters. The algorithms that are needed for this test are complicated and difficult to program. Although algebraic methods can, in theory, determine the stability of any recursive filter, there exist a number of practical difficulties in their implementation. Since even a low-order filter requires a large number of calculations, a powerful digital computer must be utilized to implement the algorithm. Therefore, the computer's finite work length restriction means that both round-off error and coefficient quantization effects must be considered. Mapping algorithms are superior to the algebraic techniques as they are easy to program and are able to handle virtually any order filter without any modifications.
Computational complexities: Based on the above review, a shown in Table 1 with the computational complexities of the various stability testing methods is prepared which, could be used by researchers as a quick reference. Table 1 shows the various approaches available for testing the stability of 2-D recursive digital filters, the techniques/procedures used and the computations involved.  [12] Double bilinear O(n 2 4 n ) 1 transformation with Ansell's result Maria and Sturm procedure O(n 2 4 n ) Fahmy [3] Schur-cohn procedure Anderson and Schur-cohn procedure O(n 2 4 n ) Jury [11] (Bilinear transformation eliminated) Reddy et al. [5] Bilinear transformation O(n 2 4 n ) in a single variable with Schur-cohn hermitian procedure Mastorakis [6] Barnett results O(3n 2 ) (for computation of polynomial matrix alone) Yang and Parametric algorithm O(n 6 ) Hwang [15] developed based on Jury's modified table Bistritz [13] Reformulated Jury's Testing last polynomial 2-D tabular test alone requires O(n 4 ) computations while overall complexity is O(n 6 ) Telepolation technique Overall complexity is O(n 4 ) Sivanandam and Routh Criterion O(n 6 ) Sivakumar [8] Hu and Jury [18] 2-D Tabular test O(n 6 ) Bose [9] Local positivity testing O(6n 4 ) 1 : It is assumed that the degrees of two variable polynomial in both the directions are n and n

CONCLUSION
In this study we have presented in detail the various stability theorems available in the literature on 2-D recursive digital filters and the approaches available for testing 2-D filters for stability. The algebraic method available in the literature is not efficient particularly when the order of the filter polynomial is greaten than four. This is because of the complexity of computations involved and also the inaccuracies that creep in a big way because of the finite precision involved in the computers. It has been concluded that for lower filters, algebraic methods offer advantages, whereas for higher order filters approximation methods are efficient.