On the Impulse Response of Singular Discrete LTI Systems and Three Fourier Transform Pairs

: A basic tenet of linear invariant systems is that they are sufficiently described by either the impulse response function or the frequency transfer function. This implies that we can always obtain one from the other. However, when the transfer function contains uncanceled poles, the impulse function cannot be obtained by the standard inverse Fourier transform method. Specifically, when the input consists of a uniform train of pulses and the output sequence has a finite duration, the transfer function contains multiple poles on the unit cycle. We show how the impulse function can be obtained from the frequency transfer function for such marginally stable systems. We discuss three interesting discrete Fourier transform pairs that are used in demonstrating the equivalence of the impulse response and transfer functions for such systems. The Fourier transform pairs can be used to yield various trigonometric sums involving sin ( π k / N ) sin ( π Lk / N ) , where k is the integer summing variable and N is a multiple of integer L .


Introduction
A linear time-invariant (LTI) system is sufficiently described by the unit impulse response function (IRF).Given an LTI system with an unknown IRF, we can use any input function to obtain the IRF.A conceptually simple way to determine a system's IRF is to use a unit impulse and unit step function as input.However, as a unit impulse has an infinitely short duration and a unit step function has an infinitely long duration, they are not always practical to produce in reality.A more practical and convenient input test function is the rectangular function.For a discrete system, this corresponds to a sequence of uniform pulses.For example, let us have the input to a discrete LTI system as x 1 [n] = 1, n = 0, 1 0, otherwise and the corresponding output is y 1 [n] = 1, n = 0, 1, 2 0, otherwise .It can be seen that the IRF is which is an infinite sequence with the repeating pattern underlined.For an LTI system, the standard method to find the system's IRF is to inverse Fourier transform the frequency transfer function where X[k] and Y[k] are the discrete Fourier transform (DFT) of the input and output sequences, respectively.In this example, however, inverse discrete Fourier Transform (IDFT) cannot be applied since a singularity (pole) exists.In this article, we are largely concerned only with the singular frequency transfer function (SFTF) of marginally stable systems.In general, deconvolution via IDFT is not doable for SFTF and part of the information is considered lost [1].
To be more general, let the input sequence x[n] have L consecutive unit impulses, i.e., x recursive sum (e.g., [2,3]).Alternatively, the IRF for the above defined x[n] can be found through the convolution sum y Since y[n] = 0 for n ≥ L y , h L y + 1 = h L y − L + 1 ].h[n] becomes periodic with a period of L starting at max L y − L + 1, 0 or sooner, where max is the maximum operator.The number of independent points in h[n] is M = max L y + 1, L .h[n] can alternatively be expressed as where , and u[m] = 0 m < 0 1 m ≥ 0 is the unit step function.In general, such an LTI system has an infinite IRF.Fundamental to LTI systems is that they are sufficiently described by either h[n] or its Fourier transform, the frequency transfer function H[k] (e.g., [3,4]).If h[n] and H[k] are indeed equivalent in characterizing an LTI system, a logical conclusion is that we can always obtain one from the other.The fact that we cannot find h[n] from a singular H[k] via IDFT may not be because H[k] has incomplete information but because we have not found the right way to carry it out.Otherwise, we have to admit H[k] and h[n] are not always equivalent.The problem appears to be simple and yet fundamental to LTI systems.Although LTI theory is very mature, obtaining IRF from an SFTF has not received much attention.Despite extensive searches in the Web of Science database and the open internet, the author cannot find any references on the direct conversion from SFTF to IRF.There are possibly two reasons for this: (1) we cannot apply DFT or IDFT when there are singularities; (2) alternative ways, such as polynomial division, can be used to obtain h[n].General discussions on the zeros of the discrete-time system can be found in [5,6].The stability of fractional oscillators is analyzed by Li et al. [7].Schluter and Darup [8] showed that linear dynamic controllers with integer coefficients are usually unstable.Park et al. [9] discussed uniformly marginally stable zeros of linear sampled data systems and its application to stable inversion-based control.In the following, we demonstrate that we can recover h[n] from H[k] when X[k], the DFT of the input function, is a sinc function.
To describe the problem more precisely, let us consider the details of H[k].To perform the discrete Fourier transform (DFT), we zero-pad x[n] and y[n] to N elements.Let X[k] and Y[k] be the DFT of the zero-padded sequence of x[n] and y[n], respectively.The N-point DFT of x[n] is defined as In this article, we restrict our input to we have where π N = π N and W = e −j2π N .X[k] has zeros at k = l N L , where l ∈ {1, . . . ,L − 1}.We denote the set of the zeros of X[k] as k o = N L , . . .
The poles, which characterize the stability of a system, are at k o , which are all on the unit cycle in the corresponding Z-transform and H[0] = Y[0]/L.Such a system is considered to be marginally stable [10].IDFT cannot be directly performed on H[k] to obtain h[n] when there are singularities.
If we limit our objective to obtain h[n] from H[k], we may attempt to avoid poles in H[k].As we are concerned with discrete operations, we can sample H[k] in such a way that there are no poles.In fact, when N is not divisible by L, H[k] does not necessarily have a singularity.In this first approach, the inverse transform, however, does not yield h[n], as shown in Figure 1.In this example, N = 17, L = 3, y[0 : 4] = 0, and y[5 : N − 1] = 0.Alternatively, when N is divisible by L, we can assign a large value at k = k o to attempt for an approximation of h[n].We show an example of this second approach in Figure 2, where N = 16, L = 2, y[0 : 6] = 1, y[7 : N − 1] = 0, and H[k o ] = 10 10 .The approach does not work in two aspects.The h[n] computed from the IDFT method yields values proportional H[k o ] we artificially set and it is not even proportionally correct at the beginning indices.This invalidates the notion that a large enough real value is a reasonable numerical approximation when there is a zero in the denominator.
to be marginally stable [10].IDFT cannot be directly performed on [] to obtain ℎ[] when there are singularities.
If we limit our objective to obtain ℎ[] from [], we may attempt to avoid poles in [].As we are concerned with discrete operations, we can sample [] in such a way that there are no poles.In fact, when  is not divisible by , [] does not necessarily have a singularity.In this first approach, the inverse transform, however, does not yield ℎ[], as shown in we artificially set and it is not even proportionally correct at the beginning indices.This invalidates the notion that a large enough real value is a reasonable numerical approximation when there is a zero in the denominator.In the following, we show that when the input function [] is a rectangular function and the output [] is a finite length sequence, ℎ[] can indeed be recovered from [] and [] when  is divisible by  .Before we delve into the proofs, we outline our method of obtaining ℎ[] from [] in Section 2. The proof of our method is facilitated by three interesting discrete Fourier transform pairs, which are discussed in Section 3.Each of the three DFT pairs is presented as a lemma.Section 4 contains the proof of our method and we point out some potential questions to be explored further in Section 5.

Summary of the Solution
Let  be divisible by  throughout the following discussion, and we define to be marginally stable [10].IDFT cannot be directly performed on [] to obtain ℎ[] when there are singularities.
If we limit our objective to obtain ℎ[] from [], we may attempt to avoid poles in [].As we are concerned with discrete operations, we can sample [] in such a way that there are no poles.In fact, when  is not divisible by , [] does not necessarily have a singularity.In this first approach, the inverse transform, however, does not yield ℎ[], as shown in we artificially set and it is not even proportionally correct at the beginning indices.This invalidates the notion that a large enough real value is a reasonable numerical approximation when there is a zero in the denominator.In the following, we show that when the input function [] is a rectangular function and the output [] is a finite length sequence, ℎ[] can indeed be recovered from [] and [] when  is divisible by  .Before we delve into the proofs, we outline our method of obtaining ℎ[] from [] in Section 2. The proof of our method is facilitated by three interesting discrete Fourier transform pairs, which are discussed in Section 3.Each of the three DFT pairs is presented as a lemma.Section 4 contains the proof of our method and we point out some potential questions to be explored further in Section 5.

Summary of the Solution
Let  be divisible by  throughout the following discussion, and we define In the following, we show that when the input function x[n] is a rectangular function and the output y[n] is a finite length sequence, h[n] can indeed be recovered from X[k] and Y[k] when N is divisible by L. Before we delve into the proofs, we outline our method of obtaining h[n] from H[k] in Section 2. The proof of our method is facilitated by three interesting discrete Fourier transform pairs, which are discussed in Section 3.Each of the three DFT pairs is presented as a lemma.Section 4 contains the proof of our method and we point out some potential questions to be explored further in Section 5.

Summary of the Solution
Let N be divisible by L throughout the following discussion, and we define and to zeros but otherwise is the same as From where L L, with int being the integer operator.We construct h 3 [k] in the following manner: where m is a positive integer, L h = int L 2 , and N h = int N 2 .Alternatively, h 3 [n] can be computed via the following operation: where and make the following modifications: where In the next section, we present three Fourier transform pairs that are used for the proofs.

Three Fourier Transform Pairs and Their Properties
where i.e., q 0 [m] and Z 0 ) form a DFT pair.In (7a), the index of q 0 can be negative.Noting that W −mN = 1 and where * is for complex conjugate, q 0 [m] can be alternatively expressed as follows: As q 0 [m] and y[n] are real, h 1 [n] is thus real as well.Since q 0 [m] is an IDFT, all the DFT/IDFT properties apply, including the periodicity property, q 0 where is the comb function with a period of L, and I is the set of all integers.
Proof.Equation (11) can be written as where ∑ N−1 i=1:N/L indicates that the summation starts at i with an incremental of N/L and ends at a value up to but not larger than N − 1.Since W iL = W (N/L+i)L for i ∈ ⋖, the denominator in each summation does not vary with the summation index.The above equation can be simplified to If m is divisible by L and for m ∈ N, we have Using the fact that As q 0 [m] = q 0 [N + m], the recursive Equations ( 12) and (13) lead to The expression can be alternatively written as Equation (11b).□ Expression (11b) consists of two linear functions in one period for m = [0, N − 1].Two example plots of q 0 [m] expressed as the sum of (11) and the closed form of (11a) are shown in Figure 3 for N = 20, L = 2 and N = 4 5, L = 5.Lemma (1) shows that the sum expressed by Equation (11a), which happens to be the Fourier transform of Z 0 [k], has a simple expression.For instance, letting m = 0, we have Figure 3. Example plots to show that the summation expressed by (11a) is the same as the closed form expressed in (11b).Note that q 0 consists of two linear functions.
Other than being real and periodic, for i, k, m ∈ I, q 0 has the following properties: Proof.From Lemma 1, (7a) can be expressed as Letting n 1 = kL + l and n 2 = n 1 + L, we have The condition for h 1 [kL + l] to be a linear function of k is that the above expression is independent of k.This condition is satisfied when Combin- ing the cases for n 1 ≥ m and n 1 < m proves the first part of the theorem.
When N − L + L y is even, we have For finite L y and bounded y[m], all of the terms in the summation are bounded.As When N + L y − L or N is odd, the same argument holds when N is replaced by N − 1.This proves the second part of the theorem.□ is typically smaller than around n = N/2.As the periodic linearity starts at L y − L + 1, the middle point of the linear region is at , which is the converging point of

Inverse Fourier Transform of Z
Proof.
The N-point The denominators of the four geometric sums are the same.The sums of the geometric sequences can be reduced to As q 1 [m] = q 1 [N + m], a general expression for q 1 [m] is Equation (20a), which can be alternatively expressed as where

□
Lemma 2 can also be proven by following the approach used in proving Lemma 1.Like q 0 [m], q 1 [m] is real and periodic.Additional properties of ) .Lemma 2 shows that this sum of trigonometric functions has only three possible values.For example, ∑ Two example plots for q 0 [m] are shown in Figure 5 for N = 20, L = 2 and N = 30, L = 5.  20) is the same as the closed form expressed in (20a).q 1 has only three possible values.

Inverse Fourier Transform of 1 − e jπk cos(π N Lk) Z[k]
We define i.e., q 2 [m] and 1 − e jπk cos(π N Lk) Z[k] form a DFT pair. where Proof.As in the proof for Lemma 2, we take the N-point DFT of q 2 [m] as expressed by Equation (23a) for m ∈ N The sums of the geometric sequences can be reduced to As Z[k o ] = 0, the above equation can be simply written as The IDFT of Q 2 [k] leads back to Equation (23a) for m ∈ N, i.e., a more general expression for any m is Equation (23a).Alternatively, we can express Equation (23a) as where m] zero at these two points.Analogous to Equation (20a), Equation (23a) can also be written as where b Other than the fact that q 2 is real and periodic, additional properties include the following: , where N 1 and N 2 are used to generate q 1 [m] and q 2 [m], respectively, and min is the minimum operator.Like q 1 [m], q 2 [m] has only three values.The lengthiness in Expression (23a) is to account for the positions of the three values.The number of positive values in q 1 per period is N/L, while there are N/L − 1 positive values in q 2 per period.Similar to Lemmas 1 and 2, Lemma 3 can also be used to find sums containing sin(πk/N) sin(πLk/N) .Two example plots of q 2 [m] are shown in Figure 6 for N = 18, L = 2 and N = 25, L = 5.

Equivalence to the Unit Impulse Response Function
Proof.Case 1: N L is an even integer.From the definition of h 1 [n] and h 2 [n], we have From Equation ( 9) and Lemma 2, we have The minimum value for as L l is zero or a multiple of L. Thus, This is the same as Equation (2).Therefore, if In order to produce M points (corresponding to index n = 0 : M − 1) accurately in Since N/L needs to be an even integer, N should be 18 to authentically produce the first five elements.
Case 2: N L is an odd integer Similar to Case 1 and applying Equation (10) for odd N/L and Lemma 3, h 3 [n] can be written as Using the same process as in case 1, h 3 [n] can be shown to be identical to h Proof.Let n = kM 1 +n 1 , where M 1 = N/2 − L l , k = int(n/M 1 ), and n 1 = mod(n, M 1 ).
The recursive Equation (10) for even N/L can be written as as demonstrated in Theorem 2. We just need to consider k ≥ 1 and N 2 − L l ≤ n < N 2 in the following.Note that mod(n 1 + (k − i)M 1 , L) = mod(n 1 , L) for all integer k and n 1 .Using property (iv) of q 1 [n] and Equation (27), the above summation results in When The above is the same equation as (25a) except the applicable range for n is improved to n < N/2.□ When N L is an odd integer and n ≤ N h− − 1, the following recursive equation appears to be true: However, the author is unable to obtain a proof at this point.

Examples
In our first example (N/L is even), input x 7A, while h 3 [n] and h 4 [n] are plotted along with the true impulse response function, h[n], in Figure 7B.As seen in Figure 7A, h 1 [n] starts to be periodically linear at index 3 i.e., L y − L +

Conclusions
The frequency transfer function, H[k], of a linear time-invariant system is generally considered to be equivalent to the unit impulse response function h[n].The normal process of obtaining the impulse function from the frequency transfer function via the DFT, however, breaks down when the latter has singularities for discrete systems.This happens when the input time-domain function is a uniform train of pulses and the transfer function contains zeros on the unit cycle.We show that the frequency transfer function contains complete information to determine the IRF in such a case.A broader fundamental question to be answered is whether the IRF can always be determined from a frequency transfer function containing singularities.If the answer is positive, the ensuring question is how to find h[n] from H[k] with any kind of pole in general.
Our method to obtain h[n] from H[k] involves performing the inverse Fourier transform of the transfer function by excluding the poles to obtain a time-domain function, h 1 [n].We show that a linear combination of h 1 [n] and its shifts can yield h[n].It is important to note that the poles in the transfer function discussed here are all on the unit circle in the z-domain.Such systems are marginally stable and h[n] is periodic after max L y − L + 1, 0 points.The number of independent points in h[n], M, is the larger of L y + 1 or L. Our method requires the length of the discrete transfer function, N, to be at least 2M to correctly capture the first M points in h[n].It is of interest to explore whether it is possible to obtain the first M points of h[n] using N < 2M.
In demonstrating that h[n] can be obtained using the DFT method, we have explored three interesting Fourier transform pairs.These DFT pairs are used not only for obtaining h[n] from H[k] but can be used to yield various trigonometric sums involving sin(πk/N) sin(πLk/N) , where N is a multiple of integer L and k is the summing integer variable.
Funding: The study is partially supported by the U.S. National Science Foundation grant AGS-2152109.
[n] = 1, 0 ≤ n ≤ L − 1 0, n ≥ L .The corresponding output of an LTI system is a sequence of real values, y[n].The last non-zero value in y[n] occurs at index L y − 1, i.e., y[n] = 0, ∀n ≥ L y and y[n] ̸ = 0 for n = L y − 1.In this study, we are confined to finite L and L y .To find the IRF, h[n], we can carry out a time-domain deconvolution via polynomial division or a

Figure 1 .
Figure 1.An example showing that ℎ[] obtained by IDFT is not correct when N/L is not an integer.Red stems are the results from IDFT, while blue circles represent the true ℎ[].

Figure 2 .
Figure 2.An example showing that ℎ[] obtained by IDFT (red) is not correct when [  ] is set to 10 10 .Note that, in this case, the IDFT results, as indicated by the red symbols, are scaled down by 10 −9 .Blue circles represent the true ℎ[] values.

Figure 1 .
Figure 1.An example showing that h[n] obtained by IDFT is not correct when N/L is not an integer.Red stems are the results from IDFT, while blue circles represent the true h[n].

Figure 1 .
Figure 1.An example showing that ℎ[] obtained by IDFT is not correct when N/L is not an integer.Red stems are the results from IDFT, while blue circles represent the true ℎ[].

Figure 2 .
Figure 2.An example showing that ℎ[] obtained by IDFT (red) is not correct when [  ] is set to 10 10 .Note that, in this case, the IDFT results, as indicated by the red symbols, are scaled down by 10 −9 .Blue circles represent the true ℎ[] values.

Figure 2 .
Figure 2.An example showing that h[n] obtained by IDFT (red) is not correct when H[k o ] is set to 10 10 .Note that, in this case, the IDFT results, as indicated by the red symbols, are scaled down by 10 −9 .Blue circles represent the true h[n] values.

Figure 4 .
Two additional examples are provided in Section 4.

Figure 4 .
Figure 4.An example of h 1 [n].In this example, N = 42, L = 3 and y is a train of 11 contiguous unit impulses.

Figure 5 .
Figure 5. Example plots to show that the summation expressed by (20) is the same as the closed form expressed in (20a).q 1 has only three possible values.

Figure 6 .
Figure 6.Example plots to show that the summation expressed by (20) is the same as the closed form expressed in (20a).

Figure 9 , 2 −
we show an example of relatively large L (80) and L y (149) with N = 4L.The output sequence, y[n], is an integer array consisting of randomized numbers from 1 to 149.h 4 [n] is the same as h[n] for the first N 2 = 160 points, while h 4 [n] equals h[n] for the first N L l = 80 points.

Figure 7 .
Figure 7.The upper panel is for h 1 [n] and h 2 [n] while the lower panel is for h 3 [n] and h 4 [n] and the true impulse response, h[n].In this example, N = 36, L = 3 and L y = 5.Note that h[n] and h 3 [n] are the same for n = 0:14, while h[n] and h 4 [n] are the same for n = 0:17.

Figure 8 .
Figure 8. Same as Figure 7 except that, in this example, N = 36, L = 4 and L y = 6.h[n] and h 3 [n] are the same for n = 0:11, while h[n] and h 4 [n] are the same for n = 0:15.

Figure 9 .
Figure 9.An example plot for h[n], h 3 [n] and h 4 [n] with L = 80, L y = 149 and N = 320.Each element in y[n] is a random integer from 1 to 149.
for finite n andL l .Thus, h 2 [n] → h[n] for finite n as N → ∞ .The corollary states that h 2 [n] can be a good approximation of h[n] if N is sufficiently large.□ If N ≥ 2 L y + 1 and N/L is even, h 4 1 .h3[n]isthe same as h[n] = [1, 0, 0, 1, 0, −1...] for the first N/2 − int(Ly/L)L = 15 points and h 4 [n] is the same as h[n] for the first N/2 = 18 points in agreement with the above discussion.In ≤ n ≤ N , with N = 36.The impulse response function for this system ish[n] = [−1, 1, 1, 1, 0, 2, 2,2, 1, −5, . ..], with the periodic oscillations starting after six data points.The periodic linearity for h[n] starts at index 6, as seen from Figure 8A.As indicated in Figure 8A, h 3 [n] is the same as h[n] for the first (N − L)/2 − int(Ly/L)L = 8 data points.In this case, N is not large enough for h 3 [n] to obtain all 10 independent points in h[n].h 4 [n], which is computed from Equation (29), is the same as h[n] for the first N−L 2 − 1 = 15 points.In