Multichannel Interpolation for Periodic Signals via FFT, Error Analysis and Image Scaling

1 Introduction

In practical applications, one often has a number of data points, obtained by sampling or experimentation, which represent the values of a continuous signal. Interpolation by simple functions (e.g., trigonometric functions or rational functions), as one of the important mathematical techniques used in physics and engineering sciences, deals with the problem of reconstructing or approximating the continuous signals from a series of discrete points. An ideal interpolation usually tells us that a continuous signal can be exactly and uniquely reconstructed from the discrete points if it satisfies some suitable conditions (e.g., bandlimited). Such an interpolation is also referred to as a sampling theorem in signal processing [1]

. In general, sampling formulas interpolate the given data even if the specific conditions for perfect reconstruction cannot be met. The error analysis for interpolation or sampling formulas is of great importance because the continuous signals that need to be recovered do not satisfy the conditions for ideal interpolation in most circumstances. Due to the wide range applications of interpolation, finding new interpolation or sampling formulas with error estimations as well as developing their fast algorithms for implementation have received considerable attentions in recently years (see e.g.,

[2, 3, 4, 5, 6, 7]).

Most of classical sampling formulas have centered around reconstructing a signal from its own samples [1, 6]. In fact, reconstructing a signal from data other than the samples of the original signal is possible. Papoulis [8] first proposed generalized sampling expansion (GSE) for bandlimited signals defined on real line $R$ . The GSE indicates that a $σ$ -bandlimited signal $f$ can be reconstructed from the samples of output signals of $M$ linear systems. That is, one can reconstruct $f$ from the samples $g_{1} (n T), \dots, g_{M} (n T)$ of the output signals

g_{m} (t) = \frac{1}{\sqrt{2 π}} \int_{- σ}^{σ} F (ω) H_{m} (ω) e^{i t ω} d ω, m = 1, \dots, M

where $F$ is Fourier transform of $f$ and $H_{1}, \dots, H_{M}$ are the system functions and $T = M π / σ$ . The idea of GSE has been extended in different directions. Cheung [9] introduced GSE for real-valued multidimensional signals associated with Fourier transform, while Wei, Ran and Li [10, 11] presented the GSE with generalized integral transformation, such as fractional Fourier transform and linear canonical transform. In [12], the authors studied GSE for quaternion-valued signals associated with quaternion Fourier transform. Importantly, the applications based on GSE as well as its generalizations have been conducted by researchers [3, 13].

A finite duration signal is commonly assumed to be periodic in practice. Accordingly, there are certain studies concerning the interpolation theory of periodic signals. A sampling theorem for periodic signals was first introduced by Goldman [14]. The author in [15] proposed pseudo-orthogonal bases for reconstructing periodic signals. In a series of papers [16, 17, 18], researchers extensively discussed the sinc interpolation of discrete periodic signals and they also derived several equivalent interpolation formulas in distinct forms. As an extension of periodic sinc interpolation, decomposing a periodic signal in a basis of shifted and scaled versions of a generating function was studied in [19]. Moreover, they further presented an error analysis for their approximation method. Recently, the non-uniform sampling theorems for periodic signals were also presented [20, 21].

In this paper, a multichannel interpolation for periodic signals is studied. We derive a general interpolation formula that depends on the transfer functions of the linear systems. The formula bears a resemblance to the classical GSE defined on real line. Nevertheless, the recipe of derivation is different from the traditional one, and moreover, the proposed formula is given by a finite summation, as opposed to a infinite summation in the traditional case. The contributions of this paper are summarized as follows:

We propose a multichannel interpolation for periodic signals (MIP). The proposed MIP is capable of generating various useful interpolation formulas by selecting suitable parameters according to the types and the amount of collected data. In addition, it would restore the original signal $f$ and some integral transformations (such as Hilbert transform) of $f$ .
A fast algorithm based on FFT is presented which brings high computational efficiency and reliability of the proposed algorithm of MIP.
To the authors’ knowledge, there seems no error analysis for classical GSE in approximating general finite energy signal on real line in the literature. By contrast, we analyze the errors that arise in approximating non-bandlimited periodic signals by the proposed MIP.
The proposed MIP is applied to image scaling problem. Its main advantage makes good use of multifaceted information of image, so that the scaled image retains lots of information of the original unscaled image. The experimental simulations show the superiority of our method over the other three popularly used conventional methods: nearest neighbour, bilinear and bicubic.

The rest of the paper is organized as follows: Section 2 recalls some preliminaries of Fourier series. Section 3 formulates MIP and presents some examples to illustrate how to use MIP flexibly. The error analysis is drawn in Section 4. In Section 5, the effectiveness of the proposed MIP for approximating signals is demonstrated by several numerical examples and the application of MIP to image scaling is also addressed. Finally, conclusions are made in Section 6.

2 Preliminaries

This part recalls some preparatory knowledge of Fourier series (see e.g. [22]). Without loss of generality, we restrict attention to the $2 π$ -periodic signals, or alternatively, we will consider the signals defined on unit circle $T$ .

Let $L^{p} (T)$ be the totality of functions such that

{∥ f ∥}_{p} := {(\frac{1}{2 π} \int_{T} | f (t) |^{p} d t)}_{T}^{\frac{1}{p}} < \infty .

For $f \in L^{2} (T) \subset L^{1} (T)$ , it can be written as

f (t) = \sum n \in Z a (n) e^{i n t}

(1)

with $\sum_{n \in Z} {| a (n) |}^{2} < \infty$ , where the Fourier series is convergent to $f$ in $L^{2}$ norm. It is well known that $L^{2} (T)$ is a Hilbert space with the inner product

(f, h) := \frac{1}{2 π} \int_{T} f (t) ¯ ¯¯¯¯¯¯¯ ¯ h (t) d t, \forall f, h \in L^{2} (T) .

If $f, h \in L^{2} (T)$ and their Fourier coefficients are respectively given by

f \sim a (n), h \sim b (n) .

By Hölder inequality, we have ${a (n) b (n)}, {a (n) ¯ ¯¯¯¯¯¯¯¯ ¯ b (n)} \in l^{1}$ and the general version of Parseval’s identity is of the form:

(f, h) = \sum n \in Z a (n) ¯ ¯¯¯¯¯¯¯¯ ¯ b (n) .

Moreover, the convolution theorem gives

(f * h) (t) := \frac{1}{2 π} \int_{T} f (s) h (t - s) d s = \sum n \in Z a (n) b (n) e^{i n t} .

(2)

The analytic signal defined by the linear combination of original signal and its (circular) Hilbert transform is regarded as an useful representation from which the phase, energy and frequency may be estimated [23, 24]. The circular Hilbert transform [25, 26] for periodic signal $f$ is given by

\begin{matrix} H f (t) & = \sum n \in Z (- i sgn (n)) a (n) e^{i n t} = \frac{1}{2 π} p . v . \int_{T} f (s) cot (\frac{}{t - s} 2) d s, \end{matrix}

where sgn is the signum function taking values $1$ , $- 1$ or $0$ for $n > 0$ , $n < 0$ or $n = 0$ respectively. By definition, we have

H^{2} f (t) = - f (t) + a (0) .

(3)

We will first concentrate on reconstruction problem of finite order trigonometric polynomials. To maintain consistent terminology with the classical case, in what follows, a finite order trigonometric polynomial is called a bandlimited signal. Specifically, a periodic signal $f (t)$ is said to be bandlimited if its sequence of Fourier coefficients possesses finite nonzero elements. Let $N = (N_{1}, N_{2}) \in Z^{2}$ , in the sequel we denote by $B_{N}$ the totality of bandlimited signals with the following form:

f (t) = \sum n \in I^{N} a (n) e^{i n t}, I^{N} = {n : N_{1} \leq n \leq N_{2}} .

(4)

3 Formulation of MIP

Throughout the paper, let $N_{1}, N_{2} \in Z$ , $M \in N^{+}$ and assume that $N_{2} - N_{1} + 1$ is divisible by $M$ , namely, $\frac{N_{2} - N_{1} + 1}{M} = L \in N^{+}$ . We cut the set of integers into pieces for convenience. Let us set

I_{k} = {n : N_{1} + (k - 1) L \leq n \leq N_{1} + k L - 1}, J_{k} = M + k ⋃ l = k + 1 I_{l} .

Then we have $I^{N} = ⋃_{k = 1}^{M} I_{k} = J_{0}$ and $Z = ⋃_{k} I_{k}$ .

For $1 \leq m \leq M$ , let

	$h_{m} (t)$	$= \sum n \in Z b_{m} (n) e^{i n t},$		(5)
	$g_{m} (t)$	$= (f * h_{m}) (t) = \frac{1}{2 π} \int_{T} f (s) h_{m} (t - s) d s .$

It follows from (2) that

g_{m} (t) = \sum n \in Z c_{m} (n) e^{i n t},

where $c_{m} (n) = a (n) b_{m} (n)$ . We particularly mention that the series (5) may not be convergent in general. Nevertheless, $g_{m} (t)$ is well defined when ${c_{m} (n)} \in l^{1}$ . In this case, $h_{m} (t)$ may be regarded as a distribution.

The proposed MIP is to reconstruct $f (t)$ from the samples of $g_{1} (t), g_{2} (t), \dots, g_{M} (t)$ . To achieve this, it is naturally to expect that the simultaneous sampling of $M$ signals will reduces sampling rate by $1 / M$ . we shall note that the signal expressed as Eq. (4

) may be represented by a shorter length of summation as long as we introduce the following vectors.

For $n \in I_{1}$ , let

	$A_{n} = [a (n), a (n + L), a (n + 2 L), \dots, a (n + (M - 1) L)],$
	$E_{n} (t) = {[e^{i n t}, e^{i (n + L) t}, e^{i (n + 2 L) t}, \dots, e^{i (n + (M - 1) L) t}]}^{T} .$

It follows that

\begin{matrix} f (t) & = \sum n \in I^{N} a (n) e^{i n t} = \sum n = I_{1} M - 1 \sum k = 0 a (n + k L) e^{i (n + k L) t} = \sum n = I_{1} A_{n} E_{n} (t) . \end{matrix}

(6)

Similar considerations applying to $g_{m} (t)$ , we have

g_{m} (t) = \sum n \in I_{1} C_{m, n} E_{n} (t),

where

C_{m, n} = [c (n), c (n + L), c (n + 2 L), \dots, c (n + (M - 1) L)] .

Owing to the periodicity of $E_{n} (t)$ , we easily obtain

E_{n} (\frac{2 π p}{L}) = e^{i n \frac{2 π p}{L}} [1, 1, \dots, 1]

for every $0 \leq p \leq L - 1$ . This leads to a simple expression for samples of $g_{m} (t)$ , that is

\begin{matrix} g_{m} (\frac{2 π p}{L}) & = \sum n \in I_{1} C_{m, n} E_{n} (\frac{2 π p}{L}) = \sum n \in I_{1} e^{i n \frac{2 π p}{L}} M - 1 \sum k = 0 c_{m} (n + k L) = \sum n \in I_{1} d_{m} (n) e^{i n \frac{2 π p}{L}}, \end{matrix}

(7)

where $d_{m} (n) = \sum_{k = 0}^{M - 1} c_{m} (n + k L)$ . The factor $e^{i n \frac{2 π p}{L}}$ gives us a hint that we may be able to express the samples of $g_{m} (t)$ in terms of the DFT matrix or its inverse matrix.

Lemma 1

Let $t_{p} = \frac{2 π p}{L}$ and ${˜ D}_{m} = [d_{m} (N_{1}), d_{m} (N_{1} + 1), \dots, d_{m} (N_{1} + L - 1)] .$ There is a matrix representation for samples of $g_{m} (t)$ in terms of the inverse DFT matrix. That is,

\frac{1}{L} [g_{m} (t_{0}), g_{m} (t_{1}), \dots, g_{m} (t_{L - 1})] = {˜ D}_{m} F_{L}^{- 1} U_{L}^{- 1},

(8)

where $F_{L}$ is the $L$ -th order DFT matrix

F_{L} = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ \begin{matrix} ω^{0} & ω^{0} & ω^{0} & \dots & ω^{0} ω^{0} & ω^{1} & ω^{2} & \dots & ω^{L - 1} ω^{0} & ω^{2} & ω^{4} & \dots & ω^{2 (L - 1)} ⋮ & ⋮ & ⋮ & ⋱ & ⋮ ω^{0} & ω^{L - 1} & ω^{2 (L - 1)} & \dots & ω^{(L - 1)^{2}} \end{matrix} ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

(9)

with $ω = e^{- 2 π i / L}$ and $U_{L}$ is a diagonal matrix

UL=⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣ω0ωN1{\huge 0}ω2N1{\huge 0}⋱ω(L−1)N1⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦.

(10)

Proof. From Eq. (7), we can rewrite $g_{m} (t_{p})$ as

{⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ \begin{matrix} d_{m} (N_{1}) d_{m} (N_{1} + 1) ⋮ d_{m} (N_{1} + L - 1) \end{matrix} ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦}^{T} ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ \begin{matrix} e^{i N_{1} \frac{2 π p}{L}} e^{i (1 + N_{1}) \frac{2 π p}{L}} ⋮ e^{i (L - 1 + N_{1}) \frac{2 π p}{L}} \end{matrix} ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ .

Note that $ω = e^{- 2 π i / L}$ , $U_{L}^{- 1} =_{L}$ and $F_{L}^{- 1} = \frac{1}{L}_{L}$ , it follows that

\begin{matrix} \frac{1}{L} {[e^{i N_{1} \frac{2 π p}{L}}, e^{i (1 + N_{1}) \frac{2 π p}{L}}, \dots, e^{i (L - 1 + N_{1}) \frac{2 π p}{L}}]}^{T} = & Product of F_{L}^{- 1} % and (p + 1) -th column of U_{L}^{- 1} . \end{matrix}

Hence we immediately obtain Eq. (8) which completes the proof. $□$

Following the definition of $A_{n}$ , for $n \in I_{1}$ , we further set $M$ by $M$ matrix

H_{n} = {[b_{k} (n + j L - L)]}_{j k},

(11)

which means the $j k$ element of $H_{n}$ is $b_{k} (n + j L - L)$ . Suppose that $H_{n}$ is invertible for every $n \in I_{1}$ and we denote the inverse matrix as

H_{n}^{- 1} = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ \begin{matrix} q_{11} (n) & q_{12} (n) & \dots & q_{1 M} (n) q_{21} (n) & q_{22} (n) & \dots & q_{2 M} (n) ⋮ & ⋮ & ⋮ q_{M 1} (n) & q_{M 2} (n) & \dots & q_{M M} (n) \end{matrix} ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ .

Next we use the elements of $H_{n}^{- 1}$ to construct the interpolation functions. Let

r_{m} (n) = {\begin{matrix} q_{m k} (n + L - k L), & if n \in I_{k}, k = 1, 2, \dots, M, 0 & if n \notin I^{N}, \end{matrix}

(12)

and define

y_{m} (t) = \sum n \in I^{N} r_{m} (n) e^{i n t} .

(13)

As with most conventional bandlimited interpolation methods, the interpolation functions of MIP are generated from some fixed functions by translations. The following result provides a matrix representation of shifted functions generated by $y_{m} (t)$ .

Lemma 2

Let $t_{p} = \frac{2 π p}{L}$ and $y_{m} (t)$ be defined by (13). The shifted functions of $y_{m} (t)$ can be expressed by

[y_{m} (t - t_{0}), y_{m} (t - t_{1}), \dots, y_{m} (t - t_{L - 1})] = {˜ V}_{m} (t) F_{L} U_{L}

(14)

where $F_{L}, U_{L}$ are given respectively by (9), (10) and ${˜ V}_{m} (t)$ is defined as

{˜ V}_{m} (t) = [v_{m, N_{1}} (t), v_{m, N_{1} + 1} (t), \dots, v_{m, L + N_{1} - 1} (t)]

with $v_{m, n} (t) = \sum_{k = 1}^{M} q_{m k} (n) e^{i (n + k L - L) t}$ , for $n \in I_{1}$ .

Proof. By the definition of $y_{m}$ and rearranging the summation terms of Eq. (13), we get

\begin{matrix} y_{m} (t - t_{p}) = & \sum n \in I^{N} r_{m} (n) e^{i n (t - \frac{2 π p}{L})} = & M \sum k = 1 \sum n \in I_{1} r_{m} (n + k L - L) e^{i (n + k L - L) (t - \frac{2 π p}{L})} = & \sum n \in I_{1} (M \sum k = 1 r_{m} (n + k L - L) e^{i (n + k L - L) t}) e^{- i n \frac{2 π p}{L}} . \end{matrix}

From the piecewise defined $r_{m} (n)$ in Eq. (12), it is easy to verify that $r_{m} (n + k L - L) = q_{m k} (n)$ for all $n \in I_{1}$ and $1 \leq k \leq M$ . It follows that $y_{m} (t - t_{p})$ can be rewritten as

y_{m} (t - t_{p}) = \sum n \in I_{1} v_{m, n} (t) e^{- i n \frac{2 π p}{L}} .

(15)

In view of the resemblance of Eq. (15) and Eq. (7) and by similar arguments to Lemma 1, we immediately have Eq. (14) which completes the proof. $□$

Until now, in addition to representing the sampled values of $g_{m} (t)$ by ${˜ D}_{m}$ , the interpolation functions $y_{m} (t - t_{p})$ ( $1 \leq m \leq M$ , $0 \leq p \leq L - 1$ ) have been constructed from the elements of $H_{n}^{- 1}$ . we shall now make a connection between $f (t)$ and $y_{m} (t - t_{p})$ . To achieve this, let

	$D_{n}$	$= [d_{1} (n), d_{2} (n), d_{3} (n), \dots, d_{M} (n)],$
	$V_{n} (t)$	$= [v_{1, n} (t), v_{2, n} (t), v_{3, n} (t), \dots, v_{M, n} (t)] .$

By straightforward computations, we get

	$A_{n} H_{n}$	$= D_{n},$		(16)
	$H_{n} {[V_{n} (t)]}_{n}^{T}$	$= E_{n} (t) .$		(17)

Having introduced above notations, we are in position to show our main result which is referred to multichannel interpolation for periodic signals (MIP).

Theorem 1

Let $f \in B_{N}$ and $g_{m}$ , $H_{n}$ be given above. If $H_{n}$ is invertible for all $n \in I_{1}$ , then

f (t) = \frac{1}{L} M \sum m = 1 L - 1 \sum p = 0 g_{m} (t_{p}) y_{m} (t - t_{p})

(18)

where $t_{p} = \frac{2 π p}{L}$ and $y_{m}$ is given by Eq. (13).

Proof. Multiplying Eq. (17) by $A_{n}$ and using Eq. (16), the expression (4) can be rewritten as

$f (t)$	$= \sum n \in I_{1} D_{n} {[V_{n} (t)]}_{n}^{T} = \sum n \in I_{1} M \sum m = 1 d_{m} (n) v_{m, n} (t)$
	$= M \sum m = 1 \sum n \in I_{1} d_{m} (n) v_{m, n} (t) = M \sum m = 1 {˜ D}_{m} [{˜ V}_{m} (t)]^{T}$
	$= M \sum m = 1 {˜ D}_{m} F_{L}^{- 1} U_{L}^{- 1} U_{L} F_{L} [{˜ V}_{m} (t)]^{T}$	(19)
	$= \frac{1}{L} M \sum m = 1 L - 1 \sum p = 0 g_{m} (t_{p}) y_{m} (t - t_{p})$

The last equality is a consequence of Lemma 1 and 2. The proof is complete. $□$

Remark 1

In general, $H_{n}$ has to be invertible for every $n \in I_{1}$ because the definition of $v_{m, n} (t)$ depends on $q_{m k} (n)$ ( $1 \leq k \leq M$ ) which are the elements of $H_{n}^{- 1}$ . However, if there is an index set $Λ \subset I_{1}$ such that $A_{n} = 0$ for every $n \in Λ$ , then $H_{n}$ is not necessary to be invertible for $n \in Λ$ . This is due to the fact that the zero terms of (6) may be removed from the summation. In this case, we just need to set $q_{m k} (n) = 0$ ( $1 \leq m, k \leq M$ ) for $n \in Λ$ .

Now we use some examples to show how Theorem 1 can be flexibly used to derive various sampling formulas for periodic signals. For simplicity we shall restrict to the case where $N_{2} = - N_{1} = N$ .

Example 1

We first present the most basic example for our results. Let $M = 1$ , $t_{p} = \frac{2 π p}{2 N + 1}$ , $b (n) = 1$ . By Theorem 1, we conclude that $f \in B_{N}$ can be recovered from its samples. That is,

f (t) = \frac{1}{2 N + 1} 2 N \sum p = 0 f (t_{p}) D_{N} (t - t_{p})

(20)

where

D_{N} (t) = ⎧ ⎪ ⎨ ⎪ ⎩ \begin{matrix} \frac{sin (\frac{1}{2} + N) t}{sin \frac{1}{2} t}, & t \neq 2 k π, 2 N + 1 & t = 2 k π; \end{matrix}

is the $N$ -th order Dirichlet kernel. This formula is referred to trigonometric interpolation for data $f (t_{p})$ ( $p = 0, 1, \dots, 2 N$ ). It should be stressed that the interpolation formulas given in [16, 17, 18] are mathematically equivalent to Eq. (20). Note that the expression (20) is not numerically stable at which the function takes $0 / 0$ form. Based on the DFT, the author in [17] obtained a numerically stable formulation which is somehow equivalent to interpolation by the FFT [27]. In fact, the DFT or FFT based expressions for trigonometric interpolation can be subsumed in Eq. (19).

Example 2

Theorem 1 permits us to express $f \in B_{N}$ from its samples and samples of its derivatives. Note that

f^{'} (t) = \sum n \in I^{N} i n a (n) e^{i n t}, f^{''} (t) = \sum n \in I^{N} - n^{2} a (n) e^{i n t} .

Let $g_{1} (t) = f (t), g_{2} (t) = f^{'} (t), g_{3} (t) = f^{''} (t)$ , then $H_{n}$ should be defined by

H_{n} = ⎡ ⎢ ⎣ \begin{matrix} 1 & i n & - n^{2} 1 & i (n + L) & - (n + L)^{2} 1 & i (n + 2 L) & - (n + 2 L)^{2} \end{matrix} ⎤ ⎥ ⎦ .

By straightforward computations,

H_{n}^{- 1} = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ \begin{matrix} \frac{2 L^{2} + 3 L n + n^{2}}{2 L^{2}} & - \frac{n (2 L + n)}{L^{2}} & \frac{n (L + n)}{2 L^{2}} \frac{i (3 L + 2 n)}{2 L^{2}} & - \frac{2 i (L + n)}{L^{2}} & \frac{i (L + 2 n)}{2 L^{2}} - \frac{1}{2 L^{2}} & \frac{1}{L^{2}} & - \frac{1}{2 L^{2}} \end{matrix} ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ .

In order to satisfy the assumption that $2 N + 1$ is divisible by $M = 3$ , we have to assume that $N = 1 (m o d 3)$ . Suppose that $N = 3 N_{0} + 1$ , then $L = \frac{2 N + 1}{3} = 2 N_{0} + 1$ . Consequently, we obtain $v_{n, m} (t)$ as follows:

	$v_{1, n} (t) =$	$\frac{2 L^{2} + 3 L n + n^{2}}{2 L^{2}} e^{i n t} - \frac{n (2 L + n)}{L^{2}} e^{i (n + L) t}$
		$+ \frac{n (L + n)}{2 L^{2}} e^{i (n + 2 L) t},$
	$v_{2, n} (t) =$	$\frac{i (3 L + 2 n)}{2 L^{2}} e^{i n t} - \frac{2 i (L + n)}{L^{2}} e^{i (n + L) t}$
		$+ \frac{i (L + 2 n)}{2 L^{2}} e^{i (n + 2 L) t},$
	$v_{3, n} (t) =$	$\frac{- 1}{2 L^{2}} e^{i n t} + \frac{1}{L^{2}} e^{i (n + L) t} - \frac{1}{2 L^{2}} e^{i (n + 2 L) t} .$

By using $L = 2 N_{0} + 1$ and letting $p = 0$ in Eq. (15), we have that $y_{1} (t)$ equals

\frac{{sin}^{3} ((N_{0} + \frac{1}{2}) t) ({N_{0}}^{2} + N_{0} + 1 - (N_{0} + 1) N_{0} cos t)}{(2 N_{0} + 1)^{2} {sin}^{3} (\frac{}{t} 2)}

and

	$y_{2} (t)$	$= \sum n \in I_{1} v_{2, n} (t) = \frac{sin (t) {sin}^{3} ((N_{0} + \frac{1}{2}) t)}{(2 N_{0} + 1)^{2} {sin}^{3} (\frac{t}{2})},$
	$y_{3} (t)$	$= \sum n \in I_{1} v_{3, n} (t) = \frac{2 {sin}^{3} ((N_{0} + \frac{1}{2}) t)}{(2 N_{0} + 1)^{2} sin (\frac{t}{2})} .$

Note that for $b_{m} (n) = (i n)^{m - 1}$ ( $1 \leq m \leq M$ ), the determinant of $H_{n}$ defined by (11) is a Vandermonde determinant. Therefore $H_{n}$ is invertible for all $n \in I_{1}$ . This would imply that we can further derive a sampling formula for recovering $f$ from its samples as well as the samples of its first $M - 1$ derivatives.

Having considered interpolation for the samples of original signal and its derivatives. We now continue to show how the samples of Hilbert transform of $f$ can be used to reconstruct $f$ itself. Interestingly, based on anti-involution property (3), MIP can be applied to compute Hilbert transform as well.

Example 3

Let $M = 1$ and $b (n) = - i sgn (n)$ . Note that $H_{n} = b (n)$ does not satisfy the conditions of Theorem 1 because it is not invertible at $n = 0$ . Nevertheless, By Remark 1, this inconsistence can be avoided if $a (0) = 0$ . Since $b (n)^{- 1} = i sgn (n)$ for $n \neq 0$ , then $y (t)$ equals

\sum 0 < | n | \leq N i sgn (n) e^{i n t} = - 2 csc (\frac{t}{2}) sin (\frac{N t}{2}) sin (\frac{1}{2} (1 + N) t)

and the sampling formula

\begin{matrix} f (t) = & \frac{- 2}{2 N + 1} 2 N \sum p = 0 H f (t_{p}) csc (\frac{t - t_{p}}{2}) \times sin (\frac{N (t - t_{p})}{2}) sin (\frac{1 + N}{2} (t - t_{p})) \end{matrix}

(21)

holds if $f \in B_{N}$ and $a (0) = 0$ .

Note that $H f (t) = \sum_{n \neq 0} (- i sgn (n)) a_{n} e^{i n t}$ satisfies the preconditions for establishment of Eq. (21) inherently. By substituting $f$ with $H f$ in (21) and using (3), we conclude that the Hilbert transform of $f \in B_{N}$ can be computed by

\begin{matrix} H f (t) = & \frac{2}{2 N + 1} 2 N \sum p = 0 (f (t_{p}) - a (0)) csc (\frac{t - t_{p}}{2}) \times sin (\frac{N (t - t_{p})}{2}) sin (\frac{1 + N}{2} (t - t_{p})) . \end{matrix}

(22)

Unlike (21), the formula (22) is valid for the case of $a (0) \neq 0$ . In fact, if $f \in B_{N}$ , we have

a (0) = \frac{1}{2 π} \int_{T} f (t) d t = \frac{1}{2 N + 1} 2 N \sum p = 0 f (t_{p})

(23)

by direct computation.

Remark 2

Given a real-valued $(2 N + 1)$ -point discrete signal ${x_{r} (p)}$ , $p = 0, 1, \dots, 2 N$ . By replacing $f (t_{p})$ with $x_{r} (p)$ and letting $t = t_{j} = \frac{2 π j}{2 N + 1}, j = 0, 1, \dots, 2 N$ in the right hand side of (22), we get a discrete signal $x_{i} (j) = H f (t_{j})$ . The DFT for ${z (j) = x_{r} (j) + i x_{i} (j)}$ is

Z (k) = ⎧ ⎨ ⎩ \begin{matrix} X (0), & k = 0, 2 X (k), & 1 \leq k \leq N, 0, & N + 1 \leq k \leq 2 N, \end{matrix}

where $X (\cdot)$ is the DFT of $x_{r} (\cdot)$ . That means that $x_{i}$ is the discrete Hilbert transform of $x_{r}$ [28]. This fact reveals the relationship between discrete Hilbert transform and continuous circular Hilbert transform. For the discrete signal of even number, similar arguments can be made, we omit the details.

Example 4

The analytic signal associated with $f$ is defined by $f_{A} (t) = f (t) + i H f (t)$ , it can be rewritten as

f_{A} (t) = N \sum n = 0_{n} e^{i n t}

where $_{n} = a_{n}$ for $n = 0$ and $_{n} = 2 a_{n}$ for $n \neq 0$ . Applying Theorem 1 to $f_{A} (t)$ , we have the following. If $f \in B_{N}$ , then

f (t) + i H f (t) = \frac{1}{1 + N} N \sum p = 0 (f (t_{p}) + i H f (t_{p})) y (t - t_{p}),

where $t_{p} = \frac{2 π p}{N + 1}$ and $y (t) = \sum_{n = 0}^{N} e^{i n t} = y_{r} (t) + i y_{i} (t)$ with

y_{r} (t) = csc (\frac{t}{2}) cos (\frac{N t}{2}) sin (\frac{1}{2} (1 + N) t)

and

y_{i} (t) = csc (\frac{t}{2}) sin (\frac{N t}{2}) sin (\frac{1}{2} (1 + N) t) .

Moreover, if $f$ is real-valued, then

	$f (t)$	$= \frac{1}{1 + N} N \sum p = 0 f (t_{p}) y_{r} (t - t_{p}) - H f (t_{p}) y_{i} (t - t_{p}),$		(24)
	$H f (t)$	$= \frac{1}{1 + N} N \sum p = 0 f (t_{p}) y_{i} (t - t_{p}) + H f (t_{p}) y_{r} (t - t_{p}) .$		(25)

Example 5

In Example 4, we express $f$ (or $H f$ ) by the samples of $f$ and $H f$ . However, the formulas (24) and (25) are only valid for the real-valued signals. By Theorem 1, it is natural to expect that we can deduce a new expression of $f$ (resp. $H f$ ) in terms of the samples of $f$ and $H f$ by letting $M = 2, b_{1} (n) = 1, b_{2} (n) = - i sgn (n)$ . Rewrite $f \in B_{N}$ as $f (t) = \sum_{n = - N}^{N + 1} a_{n} e^{i n t}$ with $a_{N + 1} = 0$ . Let $M = 2, L = N + 1$ and $b_{1} (n) = 1, b_{2} (n) = - i sgn (n)$ . Then $g_{1} (t) = f (t)$ , $g_{2} (t) = H f (t)$ and

H_{n} = [\begin{matrix} 1 & - i sgn (n) 1 & - i sgn (n + L) \end{matrix}] .

It is clear that

H_{n}^{- 1} = [\begin{matrix} \frac{1}{2} & \frac{1}{2} - \frac{i}{2} & \frac{i}{2} \end{matrix}] for - N \leq n \leq - 1, H_{0}^{- 1} = [\begin{matrix} 1 & 0 \end{matrix}

Multichannel Interpolation for Periodic Signals via FFT, Error Analysis and Image Scaling

Experimental Demonstration of Data Transmission Based on the Exact Inverse Periodic Nonlinear Fourier Transform

Fast approximation by periodic kernel-based lattice-point interpolation with application in uncertainty quantification

A Fast Template Periodogram for Detecting Non-sinusoidal Fixed-shape Signals in Irregularly Sampled Time Series

Processing of large sets of stochastic signals: filtering based on piecewise interpolation technique

Time Encoding Sampling of Bandpass Signals

Speaking-Rate-Controllable HiFi-GAN Using Feature Interpolation

ZeLiC and ZeChipC: Time Series Interpolation Methods for Lebesgue or Event-based Sampling

1 Introduction

2 Preliminaries

3 Formulation of MIP

Lemma 1

Lemma 2

Theorem 1

Remark 1

Example 1

Example 2

Example 3

Remark 2

Example 4

Example 5

Multichannel Interpolation for Periodic Signals via FFT, Error Analysis and Image Scaling

Related Research

Experimental Demonstration of Data Transmission Based on the Exact Inverse Periodic Nonlinear Fourier Transform

Fast approximation by periodic kernel-based lattice-point interpolation with application in uncertainty quantification

A Fast Template Periodogram for Detecting Non-sinusoidal Fixed-shape Signals in Irregularly Sampled Time Series

Processing of large sets of stochastic signals: filtering based on piecewise interpolation technique

Time Encoding Sampling of Bandpass Signals

Speaking-Rate-Controllable HiFi-GAN Using Feature Interpolation

ZeLiC and ZeChipC: Time Series Interpolation Methods for Lebesgue or Event-based Sampling

1 Introduction

2 Preliminaries

3 Formulation of MIP

Lemma 1

Lemma 2

Theorem 1

Remark 1

Example 1

Example 2

Example 3

Remark 2

Example 4

Example 5