A very important and fundamental operation in discrete-time signal processing is that of sampling. Discrete-time signals are often obtained from continuous-time signal by simple sampling. This is mathematically modeled as the evaluation of a function of a real variable at discrete values of time . Physically, it is a more complicated and varied process which might be modeled as convolution of the sampled signal by a narrow pulse or an inner product with a basis function or, perhaps, by some nonlinear process. The sampling of continuous-time signals is reviewed in the recent books by Marks which is a bit casual with mathematical details, but gives a good overview and list of references. He gives a more advanced treatment in . Some of these references are , , , , , , . These will discuss the usual sampling theorem but also interpretations and extensions such as sampling the value and one derivative at each point, or of non uniform sampling. Multirate discrete-time systems use sampling and sub sampling for a variety of reasons , . A very general definition of sampling might be any mapping of a signal into a sequence of numbers. It might be the process of calculating coefficients of an expansion using inner products. A powerful tool is the use of periodically time varying theory, particularly the bifrequency map, block formulation, commutators, filter banks, and multidimensional formulations. One current interest follows from the study of wavelet basis functions. What kind of sampling theory can be developed for signals described in terms of wavelets? Some of the literature can be found in , , , , . Another relatively new framework is the idea of tight frames , , . Here signals are expanded in terms of an over determined set of expansion functions or vectors. If these expansions are what is called a tight frame, the mathematics of calculating the expansion coefficients with inner products works just as if the expansion functions were an orthonormal basis set. The redundancy of tight frames offers interesting possibilities. One example of a tight frame is an over sampled band limited function expansion.

Fourier Techniques We first start with the most basic sampling ideas based on various forms of Fourier transforms , , . The Spectrum of a Continuous-Time Signal and the Fourier Transform Although in many cases digital signal processing views the signal as simple sequence of numbers, here we are going to pose the problem as originating with a function of continuous time. The fundamental tool is the classical Fourier transform defined by F ( ω ) = ∫ f ( t ) e - j ω t d t and its inverse f ( t ) = 1 2 π ∫ F ( ω ) e j ω t d ω . where j=-1. The Fourier transform of a signal is called its spectrum and it is complex valued with a magnitude and phase. If the signal is periodic with period f(t)=f(t+P), the Fourier transform does not exist as a function (it may as a distribution) therefore the spectrum is defined as the set of Fourier series coefficients C ( k ) = 1 P ∫ 0 P f ( t ) e - j 2 π k t / P d t with the expansion having the form f ( t ) = ∑ k C ( k ) e j 2 π k t / P . The functions gk(t)=ej2πkt/P form an orthogonal basis for periodic functions and () is the inner product C(k)=〈f(t),gk(t)〉. For the non-periodic case in () the spectrum is a function of continuous frequency and for the periodic case in (), the spectrum is a number sequence (a function of discrete frequency). The Spectrum of a Sampled Signal and the DTFT The discrete-time Fourier transform (DTFT) as defined in terms samples of a continuous function is F d ( ω ) = ∑ n f ( T n ) e - j ω T n and its inverse f ( T n ) = T 2 π ∫ - π / T π / T F d ( ω ) e j ω T n d ω can be derived by noting that Fd(ω) is periodic with period P=2π/T and, therefore, it can be expanded in a Fourier series with () resulting from calculating the series coefficients using (). The spectrum of a discrete-time signal is defined as the DTFT of the samples of a continuous-time signal given in (). Samples of the signal are given by the inverse DTFT in () but they can also be obtained by directly sampling f(t) in () giving f ( T n ) = 1 2 π ∫ - ∞ ∞ F ( ω ) e j ω T n d ω which can be rewritten as an infinite sum of finite integrals in the form f ( T n ) = 1 2 π ∑ ℓ ∫ 0 2 π / T F ( ω + 2 π ℓ / T ) e j ( ω + 2 π ℓ / T ) T n d ω = 1 2 π ∫ 0 2 π / T ∑ ℓ F ( ω + 2 π ℓ / T ) e j ( ω + 2 π ℓ / T ) T n d ω where Fp(ω) is a periodic function made up of shifted versions of F(ω) (aliased) defined in () Because () and () are equal for all Tn and because the limits can be shifted by π/T without changing the equality, the integrands are equal and we have F d ( ω ) = 1 T ∑ ℓ F ( ω + 2 π ℓ / T ) = 1 T F p ( ω ) . where Fp(ω) is a periodic function made up of shifted versions of F(ω) as in (). The spectrum of the samples of f(t) is an aliased version of the spectrum of f(t) itself. The closer together the samples are taken, the further apart the centers of the aliased spectra are. This result is very important in determining the frequency domain effects of sampling. It shows what the sampling rate should be and it is the basis for deriving the sampling theorem. Samples of the Spectrum of a Sampled Signal and the DFT Samples of the spectrum can be calculated from a finite number of samples of the original continuous-time signal using the DFT. If we let the length of the DFT be N and separation of the samples in the frequency domain be Δ and define the periodic functions F p ( ω ) = ∑ ℓ F ( ω + N Δ ℓ ) and f p ( t ) = ∑ m f ( t + N T m ) then from () and () samples of the DTFT of f(Tn) are F p ( Δ k ) = T ∑ n f ( T n ) e - j T Δ n k = T ∑ m ∑ n = 0 N - 1 f ( T n + T N m ) e - j Δ ( T n + T N m ) k = T ∑ n = 0 N - 1 ∑ m f ( T n + T N m ) e - j Δ ( T n + T N m ) k , therefore, F p ( Δ k ) = DFT { f p ( T n ) } if ΔTN=2π. This formula gives a method for approximately calculating values of the Fourier transform of a function by taking the DFT (usually with the FFT) of samples of the function. This formula can easily be verified by forming the Riemann sum to approximate the integrals in () or (). Samples of the DTFT of a Sequence If the signal is discrete in origin and is not a sampled function of a continuous variable, the DTFT is defined with T=1 as H ( ω ) = ∑ n h ( n ) e - j ω n with an inverse h ( n ) = 1 2 π ∫ - π π H ( ω ) e j ω n d ω . If we want to calculate H(ω), we must sample it and that is written as H ( Δ k ) = ∑ n h ( n ) e - j Δ k n which after breaking the sum into an infinite sum of length-N sums as was done in () becomes H ( Δ k ) = ∑ m ∑ n = 0 N - 1 h ( n + N m ) e - j Δ k n if Δ=2π/N. This allows us to calculate samples of the DTFT by taking the DFT of samples of a periodized h(n). H ( Δ k ) = DFT { h p ( n ) } . This a combination of the results in () and in (). Fourier Series Coefficients from the DFT If the signal to be analyzed is periodic, the Fourier integral in () does not converge to a function (it may to a distribution). This function is usually expanded in a Fourier series to define its spectrum or a frequency description. We will sample this function and show how to approximately calculate the Fourier series coefficients using the DFT of the samples. Consider a periodic signal f˜(t)=f˜(t+P) with N samples taken every T seconds to give Tn˜(t) for integer n such that NT=P. The Fourier series expansion of f˜(t) is f ˜ ( t ) = ∑ k = - ∞ ∞ C ( k ) e 2 π k t / P with the coefficients given in (). Samples of this are f ˜ ( T n ) = ∑ k = - ∞ ∞ C ( k ) e 2 π k T n / P = ∑ k = - ∞ ∞ C ( k ) e 2 π k n / N which is broken into a sum of sums as f ˜ ( T n ) = ∑ ℓ - ∞ ∞ ∑ k = 0 N - 1 C ( k + N ℓ ) e 2 π ( k + N ℓ ) n / N = ∑ k = 0 N - 1 ∑ ℓ - ∞ ∞ C ( k + N ℓ ) e 2 π k n / N . But the inverse DFT is of the form f ˜ ( T n ) = 1 N ∑ k = 0 N - 1 F ( k ) e j 2 π n k / N therefore, DFT { f ˜ ( T n ) } = N ∑ ℓ C ( k + N ℓ ) = N C p ( k ) . and we have our result of the relation of the Fourier coefficients to the DFT of a sampled periodic signal. Once again aliasing is a result of sampling. Shannon's Sampling Theorem Given a signal modeled as a real (sometimes complex) valued function of a real variable (usually time here), we define a bandlimited function as any function whose Fourier transform or spectrum is zero outside of some finite domain | F ( ω ) | = 0 for | ω | > W for some W<∞. The sampling theorem states that if f(t) is sampled f s ( n ) = f ( T n ) such that T<2π/W, then f(t) can be exactly reconstructed (interpolated) from its samples fs(n) using f ( t ) = ∑ n = - ∞ ∞ f s ( n ) sin ( π t / T - π n ) π t / T - π n . This is more compactly written by defining the sinc function as sinc ( x ) = sin ( x ) x which gives the sampling formula () the form f ( t ) = ∑ n f s ( n ) sinc ( π t / T - π n ) . The derivation of () or () can be done a number of ways. One of the quickest uses infinite sequences of delta functions and will be developed later in these notes. We will use a more direct method now to better see the assumptions and restrictions. We first note that if f(t) is bandlimited and if T<2π/W then there is no overlap or aliasing in Fp(ω). In other words, we can write () as f ( t ) = 1 2 π ∫ - ∞ ∞ F ( ω ) e j ω t d ω = 1 2 π ∫ - π / T π / T F p ( ω ) e j ω t d ω but F p ( ω ) = ∑ ℓ F ( ω + 2 π ℓ / T ) = T ∑ n f ( T n ) e - j ω T n therefore, f ( t ) = 1 2 π ∫ - π / T π / T T ∑ n f ( T n ) e - j ω T n e j ω t d ω = T 2 π ∑ n f ( T n ) ∫ - π / T π / T e j ( t - T n ) ω d ω = ∑ n f ( T n ) sin ( π T t - π n ) π T t - π n which is the sampling theorem. An alternate derivation uses a rectangle function and its Fourier transform, the sinc function, together with convolution and multiplication. A still shorter derivation uses strings of delta function with convolutions and multiplications. This is discussed later in these notes. There are several things to notice about this very important result. First, note that although f(t) is defined for all t from only its samples, it does require an infinite number of them to exactly calculate f(t). Also note that this sum can be thought of as an expansion of f(t) in terms of an orthogonal set of basis function which are the sinc functions. One can show that the coefficients in this expansion of f(t) calculated by an inner product are simply samples of f(t). In other words, the sinc functions span the space of bandlimited functions with a very simple calculation of the expansion coefficients. One can ask the question of what happens if a signal is “under sampled". What happens if the reconstruction formula in () is used when there is aliasing and () is not true. We will not pursue that just now. In any case, there are many variations and generalizations of this result that are quite interesting and useful.

Calculation of the Fourier Transform and Fourier Series using the FFT Most theoretical and mathematical analysis of signals and systems use the Fourier series, Fourier transform, Laplace transform, discrete-time Fourier transform (DTFT), or the z-transform, however, when we want to actually evaluate transforms, we calculate values at sample frequencies. In other words, we use the discrete Fourier transform (DFT) and, for efficiency, usually evaluate it with the FFT algorithm. An important question is how can we calculate or approximately calculate these symbolic formula-based transforms with our practical finite numerical tool. It would certainly seem that if we wanted the Fourier transform of a signal or function, we could sample the function, take its DFT with the FFT, and have some approximation to samples of the desired Fourier transform. We saw in the previous section that it is, in fact, possible provided some care is taken. Summary For the signal that is a function of a continuous variable we have FT: f ( t ) → F ( ω ) DTFT: f ( T n ) → 1 T F p ( ω ) = 1 T ∑ ℓ F ( ω + 2 π ℓ / T ) DFT: f p ( T n ) → 1TFp(Δk) for ΔTN=2π

For the signal that is a function of a discrete variable we have DTFT: h ( n ) → H ( ω ) DFT: h p ( n ) → H(Δk) for ΔN=2π

For the periodic signal of a continuous variable we have FS: g ˜ ( t ) → C ( k ) DFT: g ˜ ( T n ) → NCp(k) for TN=P

For the sampled bandlimited signal we have Sinc: f ( t ) → f ( T n ) f ( t ) = ∑ n f ( T n ) sinc ( 2 π t / T - π n ) if F(ω)=0 for |ω|>2π/T

These formulas summarize much of the relations of the Fourier transforms of sampled signals and how they might be approximately calculate with the FFT. We next turn to the use of distributions and strings of delta functions as tool to study sampling.

Sampling Functions — the Shah Function Th preceding discussions used traditional Fourier techniques to develop sampling tools. If distributions or delta functions are allowed, the Fourier transform will exist for a much larger class of signals. One should take care when using distributions as if they were functions but it is a very powerful extension. There are several functions which have equally spaced sequences of impulses that can be used as tools in deriving a sampling formula. These are called “pitch fork" functions, picket fence functions, comb functions and shah functions. We start first with a finite length sequence to be used with the DFT. We define ⨿ M ( n ) = ∑ m = 0 L - 1 δ ( n - M m ) where N=LM. D F T { ⨿ M ( n ) } = ∑ n = 0 N - 1 ∑ m = 0 L - 1 δ ( n - M m ) e - j 2 π n k / N = ∑ m = 0 L - 1 ∑ n = 0 N - 1 δ ( n - M m ) e - j 2 π n k / N = ∑ m = 0 L - 1 e - j 2 π M m k / N = ∑ m = 0 L - 1 e - j 2 π m k / L = { 0 otherwise L k L = 0 = L ∑ l = 0 M - 1 δ ( k - L l ) = L ? L ( k ) For the DTFT we have a similar derivation: D T F T { ⨿ M ( n ) } = ∑ n = - ∞ ∞ ∑ m = 0 L - 1 δ ( n - M m ) e - j ω n = ∑ m = 0 L - 1 ∑ n = - ∞ ∞ δ ( n - M m ) e - j ω n = ∑ m = 0 L - 1 e - j ω M m = { L ω = k 2 π / M 0 otherwise = ∑ l = 0 M - 1 δ ( ω - 2 π l / M l ) = K ? 2π/M ( ω ) where K is constant.An alternate derivation for the DTFT uses the inverse DTFT. I D T F T { ⨿ 2 π / M ( ω ) } = 1 2 π ∫ - π π ⨿ 2 π / M ( ω ) e j ω n d ω = 1 2 π ∫ - π π ∑ l δ ( ω - 2 π l / M ) e j ω n d ω = 1 2 π ∑ l ∫ - π π δ ( ω - 2 π l / M ) e j ω n d ω = 1 2 π ∑ l = 0 M - 1 e 2 π l n / M = { 0 otherwise M / 2 π n = M = ( M 2π ) ? 2π/M ( ω ) Therefore, ? M ( n ) → ( 2π M ) ? 2π/T ( ω ) For regular Fourier transform, we have a string of impulse functions in both the time and frequency. This we see from: F T { ⨿ T ( t ) } = ∫ - ∞ ∞ ∑ n δ ( t - n T ) e - j ω t d t = ∑ n ∫ δ ( t - n T ) e - j ω t d t = ∑ n e - j ω n T = { 0 otherwise ∞ ω = 2 π / T = 2π T ? 2π/T ( ω ) The multiplicative constant is found from knowing the result for a single delta function.These “shah functions" will be useful in sampling signals in both the continuous time and discrete time cases.

Up–Sampling, Signal Stretching, and Interpolation In several situations we would like to increase the data rate of a signal or, to increase its length if it has finite length. This may be part of a multi rate system or part of an interpolation process. Consider the process of inserting M-1 zeros between each sample of a discrete time signal. y ( n ) = { 0 otherwise x ( n / M ) n M = 0 ( or n = kM ) For the finite length sequence case we calculate the DFT of the stretched or up–sampled sequence by C s ( k ) = ∑ n = 0 M N - 1 y ( n ) W M N n k C s ( k ) = ∑ n = 0 M N - 1 x ( n / M ) ⨿ M ( n ) W M N n k where the length is now NM and k=0,1,⋯,NM-1. Changing the index variable n=Mm gives: C s ( k ) = ∑ m = 0 N - 1 x ( m ) W N m k = C ( k ) . which says the DFT of the stretched sequence is exactly the same as the DFT of the original sequence but over M periods, each of length N. For up–sampling an infinitely long sequence, we calculate the DTFT of the modified sequence in () as C s ( ω ) = ∑ n = - ∞ ∞ x ( n / M ) ⨿ M ( n ) e - j ω n = ∑ m x ( m ) e - j ω M m = C ( M ω ) where C(ω) is the DTFT of x(n). Here again the transforms of the up–sampled signal is the same as the original signal except over M periods. This shows up here as Cs(ω) being a compressed version of M periods of C(ω). The z-transform of an up–sampled sequence is simply derived by: Y ( z ) = ∑ n = - ∞ ∞ y ( n ) z - n = ∑ n x ( n / M ) ⨿ M ( n ) z - n = ∑ m x ( m ) z - M m = X ( z M ) which is consistent with a complex version of the DTFT in (). Notice that in all of these cases, there is no loss of information or invertibility. In other words, there is no aliasing.

Down–Sampling, Subsampling, or Decimation In this section we consider the sampling problem where, unless there is sufficient redundancy, there will be a loss of information caused by removing data in the time domain and aliasing in the frequency domain. The sampling process or the down sampling process creates a new shorter or compressed signal by keeping every Mth sample of the original sequence. This process is best seen as done in two steps. The first is to mask off the terms to be removed by setting M-1 terms to zero in each length-M block (multiply x(n) by ⨿M(n)), then that sequence is compressed or shortened by removing the M-1 zeroed terms. We will now calculate the length L=N/M DFT of a sequence that was obtained by sampling every M terms of an original length-N sequence x(n). We will use the orthogonal properties of the basis vectors of the DFT which says: ∑ n = 0 M - 1 e - j 2 π n l / M = { 0 otherwise. M if n is an integer multiple of M We now calculate the DFT of the down-sampled signal. C d ( k ) = ∑ m = 0 L - 1 x ( M m ) W L m k where N = L M and k = 0,1,...,L-1 . This is done by masking x ( n ) . C d ( k ) = ∑ n = 0 N - 1 x ( n ) x M ( n ) W L n k = ∑ n = 0 N - 1 x ( n ) [ 1 M ∑ l = 0 M - 1 e - j 2 π n l / M ] e - j 2 π n k / N = 1 M ∑ l = 0 M - 1 ∑ n = 0 N - 1 x ( n ) e j 2 π ( k + L l ) n / N = 1 M ∑ l = 0 M - 1 C ( k + L l ) The compression or removal of the masked terms is achieved in the frequency domain by using k = 0,1,...,L-1 This is a length-L DFT of the samples of x ( n ) . Unless C(k) is sufficiently bandlimited, this causes aliasing and x(n) is not unrecoverable. It is instructive to consider an alternative derivation of the above result. In this case we use the IDFT given by x ( n ) = 1 N ∑ k = 0 N - 1 C ( k ) W N - n k . The sampled signal gives y ( n ) = x ( M n ) = 1 N ∑ k = 0 N - 1 C ( k ) W N - M n k . for n=0,1,⋯,L-1. This sum can be broken down by y ( n ) = 1 N ∑ k = 0 L - 1 ∑ l = 0 M - 1 C ( k + L l ) W N - M n ( k + L l ) . = 1 N ∑ k = 0 L - 1 ∑ l = 0 M - 1 C ( k + L