当前位置:首页 >> 信息与通信 >>

Lesson02


Digital Signal Processing

? Dr. Fred J. Taylor, Professor

Introduction to Discrete-Time Systems Lesson Number: 02 [UF], not in IEEE

What’s it all about?
? ? ? What is a linear system? What are the theoretical underpinnings of linear systems? What’s with the impulse response thing?

This Lesson content is considered to be review material and would normally be studied in an introductory circuits or linear systems course. It is presented at this time for two purposes. First, to appreciate the analysis efficiency to be gained by using transform methods (Lesson 6) over the classic methods (Lesson 2). Second, to establish a framework in which discrete-time simulations can be used to study discrete-time and continuous-time systems.

Introduction
The function of a system is to process or create signals. As with signals, electronic system come in three basic flavors, continuous-time (analog), discrete-time (sampled-data), and digital (quantized). Digital systems give rise to the concept of digital signal processing or DSP. DSP has evolved to the point where it is consider an infrastructure technology, embedded into virtually every application that involves signals, images, or video. The DSP revolution is producing two classes of systems as developed below. ? As a facilitated technology: In a facilitating role, DSP systems have the outward behavior of the analog systems they replace. DSP systems are replacing significant parts of control, communications, and instrumentation systems that once were exclusively analog. The benefits gained are added features, lower cost, higher performance, and packaging advantages. As an enabling technology: DSP technology and methods are used to create new products that have defined areas that previously did not exist in entertainment, manmachine interface, software defined radio, instrumentation, biomedicine, and so forth. As a disruptive technology: DSP technology has created industries and products that where never envisioned (out of the blue), such as MP3 players, digital cameras, and so forth.

?

?

Linear and Non-linear Discrete-time Systems 线性和非线性系 统
Fundamentally, discrete-time systems can be partitioned into two broad classes called:

Lesson 02 - 1

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

? ?

linear and non-linear.

The differences between these two are profound. Linear systems are based on linear mathematical models and operations which are very familiar to the engineers and scientists. Linear system theory is supported with a myriad of analysis and synthesis tools. Non-linear systems, however, are often mathematically intractable and can frustrate even the most serious scholar. The difficulties of working with non-linear systems are often so severe that engineers and scientists will often overlook reality and assume the system being studied is linear. At other times, engineers simply give up and use simulations, neural nets, statistics, or intelligent guessing in an attempt to predict the performance of a non-linear system. It is therefore incumbent on the DSP practitioner to fully appreciate what a linear system is and what it is not.

Superposition Principle 线性叠加原理
Linear discrete-time systems exhibit a number of important properties. One of the most important is the superposition property. The superposition principle assumes that there exists a discrete-time system S having an input-output relationship xi [k ] ? S y i [k ] , where xi[k] is an ?→ input and yi[k] is an output caused by xi[k]. A system satisfying the superposition principle maps a linear combination of L inputs, xi[k],
x[k ] =

∑ α x [k ] .
i i k =1

L

1

into an output:
y [k ] =
k =1

∑ α i S (xi [k ]) = S (x[k ]),

L

2

This process is graphically interpreted in Figure 1.

Figure 1: Superposition Property If a system exhibits the superposition property then it is also said to be a linear. If a system is not linear, then it is said to be non-linear. These conditions are often casually interpreted, yielding to the temptation that “if it looks linear, it must be linear.” The next example offers a strong caveat to that presumption.

Example: Superposition
Consider the discrete-time S characterized by the linear algebraic equation

Lesson 02 - 2

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

y [k ] = ax [k ] + b , where x[k] is a causal input signal. The system response to inputs x1[k] and x2[k] are:
y 1[k ] = S( x1[k ]) = ax1[k ] + b y 2 [k ] = S( x2 [k ]) = ax 2 [k ] + b

The sum of the responses, namely y[k]=y1[k]+y2[k], is:and their sum is:
y [k ] = y1[k ] + y 2 [k ] = S ( x1[k ]) + S ( x2 [k ]) = (ax1[k ] + b) + (ax 2 [k ] + b) = ax1[k ] + ax 2 [k ] + 2b If the input is changed to x[k]=x1[k]+x2[k], the system response (labeled w[k]) is:
w [k ] = S( x1[k ] + x 2 [k ]) = ax1[k ] + ax 2 [k ] + b

If the system possesses the superposition property, then it must follow that y[k] must also equal w[k] which it does not for all b≠0. If b=0, however, superposition holds and the system is then classified as being linear.

End of example ---------------------This example points out an important observation that a system defined by a linear equation is only linear with respect to zero initial conditions. Many apparently linear systems that are missclassified as being linear simply because engineers and scientists have a natural tendency to assume that a system is at-rest (zero initial conditions and/or bias free). This assumption is sometimes based on experiences with analog systems which have a natural “tendency” to converge to a zero state in the absence of any input. It should be appreciated that digital signal processing systems are the antithesis of analog signal processors. Digital systems always remember! The results from a previous calculation may remain loitering in some digital register and assimilated unwittingly into a DSP algorithm during run-time. It is therefore strongly recommended that unless it can be absolutely guaranteed that the working registers of a DSP system are zero, they be explicitly initialized.

Time-Invariant and Time-Varying Systems 时不变和时变系统
It is generally considered to be a positive system attribute if a system behaves the same way every time it is used. To illustrate, suppose a causal signal x[k] is presented to a system S and produces a response y[k], where:

x [k ] → y [ k ]

S

3

A discrete-time system S is said to be shift-invariant or time-invariant if a time shift, or delay, in the input produces an identical time shift or delay in the output. This can be mathematically expressed as: x [k ? m ] → y [ k ? m ] ,
S

4

Lesson 02 - 3

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

Equation 4 is the test for time-invariance. If a system in not time-invariant then it is said to be time-varying.

Example: Time-Invariance
Suppose a discrete-time system is given by y[k]=(αk+β)x[k]. The output to a delayed input x[k?m] shall be denoted w[k]=(αk+β)x[k–m]. If the system is time-invariant, then it must follow that w[k]=y[k+m]. However, y[k–m]=(α(k–m)+β)x[k–m]≠(αk+β)x[k–m]=w[k] for all α≠0 and m≠0. The system is therefore time-varying. The system is, however, timeinvariant if α=0. If a system is both linear and time-invariant, then it is said to be linear time-invariant or (LTI). LTI systems are of fundamental important to digital filtering and are generally far easier to analyze and design than time-varying and/or non-linear systems. However, it is important to remember that a system can be linear and not be time-invariant (i.e., linear systems with time-varying coefficients).

Causal and Non-Causal Systems 因果系统和非因果系统
A system is said to be causal, or non-anticipatory, if the output does not precede the input (i.e., nothing happens until the system is turned-on). An at-rest (zero initial condition) LTI system is causal if and only its response to a causal input is zero for all sample values preceding the application of the input (i.e., k<0). If a system is not causal, it is said to be non-causal, or anticipatory. A non-causal system is prophetic.

Difference Equations,差分方程 ,
One of the common means of modeling a discrete-time signal is as a difference equation (DE). A homogenous DE represents a system having no external input, driven only by internal initial conditions. Specifically, an Nth order unforced difference equation has the general form:
y [k ] ?



N

an y [k ? n ] =

n =1

n =0

∑ b x [k ? n ] ,
n

N

5

where the N initial conditions are:
y [0] = y 0 ; y [ ?1] = y ?1; ... ; y [-(N ? 1)] = y ? ( N ?1) , 6

The shape or envelope of the solution y[k] is determined by the system’s input x[k], coefficients {ai, bi}, and initial conditions y[-i].

Lesson 02 - 4

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

Example: Difference Equation Production
Difference equations can often be extracted from continuous-time system models. For example, the forced response behavior of an at-rest continuous time signal is defined in terms of the system’s impulse response h(t)冲激响应 and the convolution integral:卷积
y (t ) =


?∞

∫ x (τ )h(t ? τ )dτ

卷积

Suppose h(t) = λe-λt u(t) = 105.361 e-105.361t u(t), then from elementary continuous-time linear system theory, the step response of the system is given by: x(t)=u(t)
y (t ) =


?∞

λ ∫ x (τ )h(t ? τ )dτ = (1? e )u (t )
? t

阶跃响应

The response to x(t)=u(t)-u(t-0.1), a 1/10th of a second pulse, is:

? 1? e ? λt for t ∈ [0,0.1) y (t ) = ? ? λ (t ? 0.1) ? λ (0.1) e for t ≥ 0.1 ?e

(

(

)

)

In order to emulate the response of the continuous time-system with a discrete-time model, a sample rate is chosen which is sufficiently high to insure that a continuous-time signal can be approximated with a dense cover of samples. For a chosen sample rate: h[k] ? Ts h(kTs)

对连续信号采样 采样率 Ts

For fs=1k Sa/s, for example, Ts=1 ms (10-3). At t=kTs, h(t) = λe-λt = λe-λkTs = λe(-λTs)k = 105.361 e0.10536k = 105.361 (0.9)k. Therefore the discrete-time model of the 1st order continuous-time system’s impulse response is: h[k] ? Ts h(kTs) = 0.105361(0.9)k. 乘了 Ts For computational purposes, it is assumed that system’s impulse response is extinguished after 50 samples ((0.9)50 ~ 0.0) as shown in Figure 2 (first panel). The simulation of the discrete-time systems to a unit pulse of duration 0.1 second by a linear system having an impulse response h[k] = 0.105361(0.9)k u[k], is simulated using the MATLAB using the code shown below. The response (second panel) exhibits a transient “build-up,” followed by at steady state period, followed by a “build-down” to the zero state. MatLab code ? Ts=0.001; lamda=105.361;%设置采样周期 Ts 和λ ? n=0:50; scale=Ts*lamda; %生成时间序列,计算刻度比例 ? h=scale*(0.9).^n; %计算冲激响应序列 ? subplot(1,2,1); %选择子图 ? stem(h); %用 stem(h);绘出离散图 ? x=ones(1,100); %输入方波 ? y=conv(h,x); %求输出

Lesson 02 - 5

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

? subplot(1,2,2); %选择子图 ? stem(y) %绘出输出离散图
0.12 0.1 0.08 0.06 0.6 0.04 0.02 0 0 10 20 30 40 50 60 0.4 0.2 0 0 50 100 150 1.4 1.2 1 0.8

y[k], k∈[0,150] Steady State Transient

h[k]=(0.105361)0.9k,

Figure 2: Simulated response of a system h(t)=e-105.361t to x(t)=u(t) and x(t)=u(t)-u(t-0.1) using a discrete-time model sampled at 1 kHz with Ts=0.001 seconds.

Example: Difference Equation
Assume that the solution to the homogeneous, or unforced difference equation y [k + 1] ? ay [k ] = 0 , for y[0]=y0, is y [k ] = y 0a k u[k ] , where u[k] is a unit step. In this case u[k] denotes signal and system causality. The candidate solution can be verified using substitution. First observe that y[0]=y0a0u[0]=y0, as required. Inductively, note that y[k+1]= y0a(k+1)u[k+1]= a(y0ak)=ay[k].

Homogeneous Solution
The solution of a homogeneous linear difference equation can be approached in the manner historically used to solve an Nth order homogeneous ordinary differential equations (ODE). This process begins with the definition of an ODE form,

∑a
i =0

N

i

d i y (t ) = 0 , N阶常系数齐次微分方程 dt i

7

Classically, the solution to Equation 7 is based on the intelligent use of eigen functions. Continuous-time eigen functions are assumed to have the form φn (t ) = es n t , n∈[1,N], where sn is an eigenvalue of the ODE’s characteristic equation. The general solution is given by: y (t ) =

∑α iφi (t ) =∑α i esi t , 通解
i =1 i =1

N

N

8

where the coefficients αi, and therefore the solution, can be obtained using the method of indeterminate coefficients (solving N equations in N unknowns).

Lesson 02 - 6

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

A similar approach can be applied to the problem of determining the response of a system having a homogeneous difference equations of the form:

∑ a y [k-i ] = 0 .
i i =0

N

N阶常系数差分方程

9

The eigenfunctions 特征函数 of an Nth order difference equation takes on the form Φ n [k ] = λk , n n∈[0,N], where λn is a root (eigenvalue) 特征值 of the difference equations characteristic equation
i =0

∑ ai λin = 0 .

N

10

The difference equation’s homogeneous solution can be expressed in terms of a linear combination of eigenfunctions Φ n [k ] = λk , namely n

y [k ] =

∑α i Φ i [k ] = ∑α i λki , 通解
i =1 i =1

N

N

11

where the coefficients αi are computed using the method of indeterminate coefficients.

Example: Eigenfunctions
Consider the unforced 2nd order difference equation: y[k]=(3/4)y[k–1] –(1/8)y[k–2], y[0]=1, y[–1]=0. The characteristic equation is given by λ2–(3/4)λ+1/8=(λ–1/2)(λ –1/4) and the eigenvalues are λ1=1/2, and λ2=1/4. From Equation 11, it follows that 解出特征值 y[k]=(α1λ1k+α2λ2k)=(α1(1/2)k+α2(1/4)k).

写出通解

Substituting the known initial conditions into the above equation results in y[0]=(α1+α2)=1

套用初始条件

y[–1]= (α1(1/2)–1+α2(1/4)–1)=(2α1+4α2)=0, which can be expressed in matrix-vector form as:
? 1 1 ?? α1 ? ? 1 ? ? ? 2 4 ?? α ? = ? 0 ? . ?? ? ? ? ? ?? 2 ? ? ?

解线性方程

Solving the matrix equation results in α1=2 and α2= –1. Therefore, y[k]=2(1/2)k–1(1/4)k which satisfies the conditions y[0]=1 and y[–1]=0 as required. Upon substituting y[k] into the given DE, one obtains the required: 得到解

Lesson 02 - 7

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

y [k ] =

3 1 y [k ? 1] ? y [k ? 2] 4 8 k ?1 k ?1 k ?2 k ?2 ? ? 3 ? ? 1? ? 1? ? ? 2? ? ? 1? 1 ? ? ? 1 ? 2? 1 ? = ? 1? ? ? ? ? ? ? 4? ?2? ?4? ? 8? ?2? ?4? ? ? ? ? ?
k k

? 1? ? 1? = 2? ? ? 1? ? , 2? ? ?4?

LTI Discrete-Time Systems 线性时不变离散时间系统
Most contemporary digital filters are based on an Nth order LTI system model having the form:

m =0



N

am y [k ? m] =

m=0

∑b

M

m x [k

? m ] , 线性时不变非齐次差分方程

12

with the appropriate initial conditions specified. The system is said to be proper if N≥M, and strictly proper if N>M. If a0=1, then the system is said to be monic. 归一化 Equation 12 can also be expressed as
y [k ] =
M 1? N ? ? ? ∑ am y [k ? m ] + ∑ bm x [k ? m ] ? . a0 ? m =1 m =0 ?

13

The solution to Equation 13 can be modeled as the sum of two solutions called the homogeneous (unforced)非强迫(齐次) and inhomogeneous (forced 强迫)非齐次 solutions. The homogeneous equation satisfies:
m =0

∑ am y h [k ? m] = 0 .

N

齐次,只有初始条件的影响

14

where yh[k] denotes the homogeneous solution which is “powered” only by the system’s initial conditions. The homogeneous solution can be defined in terms of roots of the system’s characteristic equation:
m =0

∑ am r N ? m = 0 .

N

特征方程

15

These roots were previously are called eigenvalues 特征值. There may be up to N distinct solutions to this Equation 15. If a root ri is unique (i.e., non-repeated), then the root is said to be distinct or have multiplicity one. If a root ri occurs ni times, it is said to be repeated and has a multiplicity ni. The solution to the homogeneous equation, namely Equation 14, is in general: y h [k ] = ∑ γ i k ni ?1ri k ,
i =1 N

16

Lesson 02 - 8

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

where the coefficients γi, for i∈{1, … ,N}, are constants whose value is dependent upon knowing N initial conditions.

Example: Homogeneous Solution 齐次解
Consider the homogeneous LTI system equation given by
y [k ] ? 1 1 y [k ? 1] ? y [k ? 2] = 0; y [ ?1] = 1; y [ ?2] = 1 4 8

Beginning with the characteristic equation, observe that: r2 ? 1 1 r ? = 0 特征方程 4 8

or
1 ?? 1? ? ? r ? ?? r + ? = 0 解特征方程 2 ?? 4? ?

The roots of the characteristic equations are therefore seen to be distinct and are given by r1 = 1/2 and r2 = –1/4 (i.e., n1=n2=1). Therefore, the homogeneous solution has the general form:
? 1? ? 1? y [k ] = γ 1? ? + γ 2 ? ? ? , 齐次解 ?2? ? 4? The task therefore becomes on of determining the values of γ1 and γ2. To solve for γ1 and γ2 , two independent equations in two unknowns are needed. Two data points are known in this case and they are the given initial conditions y[–1] and y[–2]. This information can be used to specify ? 1? ? 1? y [ ?1] = γ 1? ? + γ 2 ? ? ? ? 2? ? 4? ? 1? y [ ?2] = γ 1? ? ?2?
?2 ?1 ?1 k k

= 2γ 1 ? 4γ 2 = 1
?2

? 1? + γ 2? ? ? ? 4?

= 4γ 1 + 16γ 2 = 1 带入初始条件

which, in turn, can be expressed in matrix-vector form as

? 2 ? 4? ? γ 1 ? ?1? ? ? ? ? = ? ? . 求出系数 ? 4 16 ? ?γ 2 ? ?1? The solution is γ1=5/12 and γ2=–1/24. For k≥0, it follows that

y [k ] =

10 ? 1 ? 1 ? 1? ? ? ? ? ? ? . 得到解 24 ? 2 ? 24 ? 4 ?
Lesson 02 - 9

k

k

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

The homogeneous solution, by itself, has limited value. What is of far greater interest, in general, is knowledge of the forced or inhomogeneous response or response to an arbitrary input x[k]. It can be shown that the inhomogeneous response of a discrete-time LTI system is defined in terms of the system’s impulse response which is closely related to the homogeneous system response.

Impulse Response 冲激响应
The impulse response of a causal LTI system is, as the name implies, the output trajectory of an at-rest (i.e., y[k]=0 for all k<0) system to an impulse input (i.e., x[k]=δ[k]). The system’s impulse response, denoted h[k], of the system defined in Equation 12, is: 冲激响应就是在 K 时刻有一个脉冲输入,其余时刻输入都为 0,形成的输出数字序列 h[k]

m =0

∑a

N

m y [k

? m] =

m =0

∑b

M

m x[k

? m]

12

令 x[k]=δ[k] 则 y[k]=h[k] 得到下面方程:

m =0

∑a

N

m h[ k

? m] =

m =0

∑b

M

m

δ [ k ? m] .

17

The production of the impulse response can be patterned after the process used to produce the homogeneous system response, which was based on the knowledge of the system’s initial conditions. The inhomogeneous solution, however, assumes the system’s initial conditions are zero. Therefore the procedure must be “tricked” into finding replacements for the missing initial conditions. Consider using a set of pseudo-initial conditions, denoted {h[0],h[1], … ,h[M]}, which represents the first M+1 samples of output y[k] driven by a unit impulse δK[k]. The pseudo-initial conditions can be obtained by solving the matrix-vector equation

脉冲响应的产生可以用处理齐次方程的方法解决问题。基于系统的初始条件。非齐次解则假定初 始条件为 0。因而,找到初始条件是解决问题的诀窍。考虑一组假定的初始条件,表示 为;{h[0],h[1], … ,h[M]},代表在单位脉冲δK[k]的激励下,y[k]的第 M+1 个采样。这个假定的初始条 件可以通过解矩阵方程获得:
? a0 ? a Ah = ? 1 ? M ? ?aM 0 a0 M aM -1 0 ? ? h[0] ? ? b0 ? ? L 0 ? ? h[1] ? ? b1 ? ? ? = ? ? =b. O M ?? M ? ? M ? ?? ? ? ? L a0 ? ?h[M ]? ?bM ? L

18

The solution to Equation 18 can be solved for h by

h = A ?1 b .

19

produces the an array of data (i.e., h) that can serve the role of the pseudo initial conditions used to solve the homogeneous problem. Once the pseudo initial-conditions are determined, h[k] can be computed for all time as: 求得 h 后就可以用齐次方程的方法,解决问题。

Lesson 02 - 10

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

h[k ] = γ 0δ [k ] + ∑ γ i k ni ?1ri k , 齐次方程解和特解求得h[k]
i =1

N

20

which satisfies the system equation (Equation 17) for all sample index k>M, namely:
m =0

∑ am h[k ? m] = 0 .

N

满足方程17

21

Example: Impulse Response
An at-rest 3rd order LTI is given by

y [k ] ?

1 1 y [k ? 1] ? y [k ? 2] = 2 x [k ] ? x [k ? 1] + x [k ? 2] . 3阶线性时不变差分方程 4 8

If the input x[k] is an impulse (i.e., x[k]=δ[k]), the impulse response, for k∈{0,1,2} can be defined as the solution to Ah=b, namely:

? ? 1 ? 1 Ah = ?? ? 4 ? 1 ?? 8 ?

? 0? ? ?h[0]? ? 2 ? 1 0? ? h[1] ? = ? ? 1? = b . 求假定的初始条件h[0]、h[1]、h[2] ? ? ? ?? 1 ? ?h[2]? ? 1 ? ? ? ? ? ? 1? 4 ? 0

Alternately, the pseudo initial conditions can be produced using procedure motivated below. k 0 1 2 x[k] 1 0 0 x[k-1] 0 1 0 x[k-2] 0 0 1 y[k] 2 -1/2 9/8 y[k-1] 0 2 -1/2 y[k-2] 0 0 2

Solving the matrix equation results in a set of pseudo initial conditions h[0]=2, h[1]=–1/2, and h[2]=9/8. The assumed form of the impulse response is:

用初等变换求得 h[0]=2, h[1]=–1/2、h[2]=9/8,然后假定冲激响应表达式如下:

h[ k ] = γ 0δ [ k ] + γ 1r1k + γ 2 r2k , 带入h[0]h[1]h[3]
where r1 and r2 are the roots of the characteristic equation (eignevalues), namel

再求得方程 y [k ] ?

1 1 y [k ? 1] ? y [k ? 2] = 2 x[k ] ? x[k ? 1] + x [k ? 2] 的特征根 4 8 1 1 ; r2 = ? . 2 4

r1 =

Lesson 02 - 11

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

Evaluating, with respect to the pseudo-initial conditions, results in
? 1 ?1 ? ? 1 ?0 2 ? ? ? 1 ?2 ?0 ? ? ? ?2? ? produces: 就可以求得 ? ? ? ? ? ?γ ? ? 2 ? 0 1 ?? ? ? 1? ? ? γ 1 = ?? ? , 4 ?? ? ? 2? 2 ?γ 2 ? ? ? ? 9 ? ? 1? ? ? 8 ? ? ? ? ? ? ? ? 4? ? ? 1

γ 0 = ?8; γ 1 = ; γ 2 =

8 3

22 . 3

It therefore follows that the impulse response is given by
k k ? ? ? ? 8δ [k ] + 8 ? 1 ? + 22 ? ? 1 ? ?u[k ] .得到冲激响应 h[k ] = ? ? ? ? ? 3?2? 3 ? 4? ? ? ?

Notice that h[k] agrees with the pseudo-initial conditions at k={0,1,2}. Unfortunately, it can also be seen that this method is computationally intensive and may be unsuitable for use with highorder systems. Fortunately more computational efficient paradigm exists in the form of transforms and, in particular, the z-transform. When introduced (Lesson 6), it will be used to produce homogeneous, inhomogeneous, and impulse responses in a relatively straightforward manner. Nevertheless, it can be seen that the impulse response is the gateway to the study of discrete-time systems.

Example: Classical Discrete-Time System Analysis
Consider the 2nd order difference equation given by:2 阶差分方程 y[k+2] -5 y[k+1] + 6 y[k] = x[k+1] – 5 x[k] with an input 输入序列 x[k] = (3k+5)u[k] and known responses (initial conditions): y[0]=4, y[1]=13 The classical analysis process begins with the determination of the system’s eigenvalues. (γ2-5γ+6) = 0, γ1=2, γ2=3 (eigenvalues) 首先解出特征根 From the eigenvalues, the general form of the solution can be stipulated and it consists of two parts called the natural and particular solutions. For the system under study:

初始条件

Lesson 02 - 12

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

y[k] = {b1(γ1)k+ b2(γ2)k}natural solution + {c1k+c2}particular solution. 通解 特解

The natural solution is associated with the eigenvalue solution and the particular solution takes on the form of the system’s forcing function x[k]. The assumed particular solution form is then substituted into the defining difference equation y[k+2] -5 y[k+1] + 6 y[k] = x[k+1] – 5 x[k], x[k] = (3k+5)u[k], resulting in:用假定的特解{c1k+c2}带入方程 (c1k+2c1+c2) -5(c1k+c1+c2)+6(c1k+c2)=(3k+8)-5(3k+5) Solving for the unknown coefficients c1 and c2, one obtains: (2c1)=(-12); (-3c1+c2) = -(17) ? c1= -6; c2= -35/2 Therefore, the particular solution is y[k]particular= -6k -35/2. The complete solution, which now has the form y[k] = {b1(γ1)k+ b2(γ2)k}natural solution + {-6k-35/2}particular solution can now be computed by invoking the initial conditions: 利用初始条件,带入特解 k=0,k=1 y[0] = {b1+ b2}natural solution + {-35/2}particular solution.= 4 y[1] = {b1(2)+ b2(3)}natural solution + {-6 -35/2}particular solution = 13 Solving, one obtains b1 = 28; b2= -13/2. Therefore 求得 b1,b2 得到差分方程解 y[k] = {28(2)k – (13/2) (3)k}natural solution + {-6k-35/2}particular solution.

Difference Equation Diagrams 差分方程图
A common means of studying LTI system is through the use of a computer simulation. Simulations can be created from the defining difference equation. Associated with each difference equation is discrete-time block diagram. The basic building blocks are shift register delays, and adders, and parametric multipliers. For example, y[k]=0.5x[k]+0.4x[k-1]-0.2x[k-2], has a block diagram description shown in Figure 3.

x[k]

0.5 Delay 0.4 Delay -0.2

Σ

y[k]

Σ

Lesson 02 - 13

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

Figure 3: Block diagram representation of y[k]=0.5x[k]+0.4x[k-1]-0.2x[k-2]. Observe that the system shown in Figure 3 possess only feedforward paths. As a result, the system is sometimes referred to as a non-recursive filter or finite impulse response filter. The propagation of an impulse through the system, and the attendant impulse response can be directly inferred from the block diagram shown in Figure 3. The impulse response is initially h[0]=0.5, h[1]=0.4, and h[2]=-0.2 which sets the filter’s coefficients (called tap weight coefficients).

Example: MATLAB IIR Simulation 无限冲激响应滤波器仿真
MATLAB provides a mechanism to compute the response of an FIR system to an arbitrary input. The MATLAB function is called conv (convolve) and is summarized below. Program conv is used to compute the impulse response of an LTI shown in Figure 3. MATLAB 提供了一套计算任 意输入的滤波器响应的方法。函数名 conv,卷积 ? help conv 解决具有这样表达式的情况 y[k]=0.5x[k]+0.4x[k-1]-0.2x[k-2]

CONV Convolution and polynomial multiplication. 卷积和多项式乘法 C = CONV(A, B) convolves vectors A and B. The resulting vector is length LENGTH(A)+LENGTH(B)-1. 结果矢量的长度 If A and B are vectors of polynomial coefficients, convolving them is equivalent to multiplying the two polynomials. 如果 A 和 B 是多项式系数矢量 See also XCORR, DECONV, CONV2, FILTER, and CONVMTX in the Signal Processing Toolbox.MATLAB Simulation The following script is used to create two coefficient arrays and an impulse which are combined to produce the system’s impulse response. The following script is used to create two coefficient arrays and an impulse which are combined to produce the system’s impulse response. ? a=[.5, .4, -.2]; {filer coefficients} ? b=[1,0,0,0,0]; ? conv(a,b) ans = 0.5000 0.4000 -0.2000 0 0 0 0 {simulated impulse}

The system’s impulse response is therefore y[k]={0.5, 0.4, -0.2, 0, 0, … }.

Lesson 02 - 14

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

Example: MATLAB IIR Simulation 递归的情况
The difference equation y[k] + 0.5y[k-1] + 0.5y[k-2] =0.5x[k]+0.4x[k-1]-0.2x[k-2], has a block diagram description shown in Figure 4.

x[k]

0.5

Σ

Σ -

y[k]

Delay 0.4 Delay -0.2 0.5 Σ Σ 0.5

Delay

Delay

Figure 4: Block diagram representation of y[k]=0.5x[k]+0.4x[k-1]-0.2x[k-2]. Observe that in this case the system possess both feedforward 前馈 and feedback paths 反馈. As a result, the system is sometimes referred to as a recursive 递归 filter. The presence of feedback paths give rise to the possibility of an infinitely long impulse response and the infinite impulse response (IIR) filter. The propagation of an impulse through the system, and the attendant impulse response are now far more complicated than in the non-recursive case. This mathematical challenge can be overcome with a well designed simulation.

Example: MATLAB IIR Simulation
MATLAB provides a mechanism to compute the response of an IIR system to an arbitrary input. The MATLAB function is called filter and is summarized below. Program filter is used to compute the impulse response presented in Figure 3. ? help filter FILTER One-dimensional digital filter. Y = FILTER(B,A,X) filters the data in vector X with the filter described by vectors A and B to create the filtered data Y. The filter is a "Direct Form II Transposed" implementation of the standard difference equation:解决如下表达形式的解 a(1)*y(n) =b(1)*x(n)+b(2)*x(n-1)+ ... +b(nb+1)*x(n-nb) -a(2)*y(n-1)- ... -a(na+1)*y(n-na) If a(1) is not equal to 1, FILTER normalizes 归一化 the filter coefficients by a(1). When X is a matrix, FILTER operates on the columns of X. When X is an N-D array, FILTER operates along the first non-singleton dimension.

Lesson 02 - 15

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

[Y,Zf] = FILTER(B,A,X,Zi) gives access to initial and final conditions, Zi and Zf, of the delays. Zi is a vector of length MAX(LENGTH(A),LENGTH(B))-1 or an array of such vectors, one for each column of X. 带出条件的解 FILTER(B,A,X,[],DIM) or FILTER(B,A,X,Zi,DIM) operates along the dimension DIM. See also FILTER2, FILTFILT (in the Signal Processing Toolbox). The following script is used to create two coefficient arrays and an impulse which are combined to produce the system’s impulse response. ? b=[.5, .4, -.2]; a=[1,.5, .5]; x=[1,0,0,0,0,0]; {simulated impulse} ? y=filter(b,a,x) y = 0.5000 0.1500 -0.5250 0.1875 0.1688 -0.1781 The system’s impulse response is therefore y[k]={0.5, 0.15, -0.525, 0.1875, … }. One of the primary uses of systems described by difference equations is as a filter (a.k.a., digital filter). These systems were previously modeled as Equation 13 which is repeated below of convenience.
y [k ] =
M 1? N ? ? ? ∑ am y [k ? m] + ∑ bm x [k ? m] ? . a0 ? m =1 m =0 ?

22

Equation 22 describes a general filter class called a recursive filter. If a0≠0, and all the other ai=0, the filter simplifies to the aforementioned non-recursive or FIR filter. Such a system would satisfy the difference equation:
y [k ] = 1 a0 ? M ? ? bm x [k ? m ] ? .非递归滤波器y输出只和输入x有关 ? ? ? m =0 ?



23

Two simple but important member of this class are called the Nth order moving average (MA) and the comb filter. They are mathematically represented by the difference equations: y [k ] = ? 1 ? N ?1 ? x[k ? m] ? .(Moving Average Filter) 移动平均滤波器 ? N ? m =0 ? ?



23 22

y [k ] = x [k ] ± x [k ? (N ? 1)] .梳状滤波器

and interpreted in block diagram form in Figure 5.移动平均滤波器和梳状滤波器结构图

Σ

Σ

Lesson 02 - 16

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

Figure 5: Moving Average and Comb FIRs It should be apparent that the MA FIR computes a "running average" of the input signal contained in a moving N-sample uniform window while the comb FIR computes the sum or difference of samples of the input that are separated by N-samples.

Example: Moving Average Filter 移动平均滤波器
A moving average filter (developed later in as a type of filter) is a simple but elegant means of smoothing data. Moving average filters of length N is intrinsically a lowpass filter having an impulse response: 1/N h[k]=(1/N)×[1, 1, … , 1] {scaled string of M 1’s} 1/N x This filter sums the N most current samples of an input signal x[k] and then scales the sum by (1/N) to obtain the average values of the N samples. MATLAB often employs a different indexing scheme that is contrary to common convention. In “MATLAB speak” the Nth order moving average filter is called an Mth order filter where M=N-1. An 12th order filter (M=11 in MATLAB speak), is simulated below for an input consisting of 2 sinusoidal inputs, one at a low frequency (5% of the sample rate) and one at a much higher frequency (45% of the sample rate). The two signals are added together and filtered (convolved using the MATLAB filter command). The output is shown to be essentially the low frequency signal where the high frequency signal has be averaged to a value near zero.
% Simulation of an M-point Moving Average Filter % Generate the input signal n = 0:100; %%取100个样本 s1 = cos(2*pi*0.05*n); % A low-frequency sinusoid 生成f=0.05的低频率离散正弦波序列 s2 = cos(2*pi*0.47*n); % A high frequency sinusoid 生成f=0.47的高频率离散正弦波序列 x = s1+s2; %%合成输入序列x % Implementation of the moving average filter 用移动平均滤波器 M = input('Desired length of the filter = '); %%输入滤波器长度,就是M个数平均 num = ones(1,M); %%生成滤波器 M个数平均,y(k)=x(k)+x(k-1)+…+x(k-M+1) y = filter(num,1,x)/M; %%用filter计算输出,B=num,A=1,X=x % Display the input and output signals clf; subplot(2,2,1); %%选择子图 plot(n, s1); %画出f=0.05的低频率离散正弦波序列 axis([0, 100, -2, 2]); %%设置坐标轴 xlabel('Time index n'); ylabel('Amplitude');%%显示坐标轴标签 title('Signal #1'); %%设置标题 subplot(2,2,2); plot(n, s2); %画出f=0.47的高频率离散正弦波序列 采样频率和正弦波频率不符合奈奎斯特要求的结果 axis([0, 100, -2, 2]); xlabel('Time index n'); ylabel('Amplitude'); title('Signal #2'); subplot(2,2,3); plot(n, x); %画出合成的离散正弦波序列 axis([0, 100, -2, 2]); xlabel('Time index n'); ylabel('Amplitude'); title('Input Signal'); subplot(2,2,4);

Lesson 02 - 17

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

plot(n, y); %%%画出滤波后离散正弦波序列 滤出了x1,抑制了高频率的x2 axis([0, 100, -2, 2]); xlabel('Time index n'); ylabel('Amplitude'); title('Output Signal'); axis;

1 0.5 0 -0.5 -1 2 Amplitude Amplitude 1 0 -1 -2 0 20 40 60 Time index n 80 100 0 20 40 60 Input Signal 80 100

1

x1[n]

0.5 0 -0.5

x2[n] -1
2 1 0 -1 -2

0

20

40 60 Output Signal

80

100

0

20

40 60 Time index n

80

100

Example: Moving Average Filter (again) 例二移动平均滤波器
Suppose an 50th order Moving Average (MA) FIR is presented with a sampled signal given by x[k]=1 + cos(0.2πk) + cos(0.4πk). The 50th order moving average FIR should intuitively behave as an averaging (i.e., smoothing or lowpass filter), attenuating anything having a high frequency (relative to the sample rate). Since the only low-frequency signal contained in the test signal is DC, it sands to reason that at steady state the DC will signal will be passed and all other signal components will be attenuated. This hypothesis can be tested using simulation. Shown in Figure 5 is the result of processing a 150-sample signal record x[k] with a 50th order MA filter. The filtering is performed using the MATLAB filter function. Observe that the output has transient and stead-state regimes.
FILTER One-dimensional digital filter. Y = FILTER(B,A,X) filters the data in vector X with the filter described by vectors A and B to create the filtered data Y. The filter is a "Direct Form II Transposed" implementation of the standard difference equation: a(1)*y(n) = b(1)*x(n) + b(2)*x(n-1) + ... + b(nb+1)*x(n-nb) - a(2)*y(n-1) - ... - a(na+1)*y(n-na) If a(1) is not equal to 1, FILTER normalizes the filter coefficients by a(1). When X is a matrix, FILTER operates on the columns of X. When X is an N-D array, FILTER operates along the first non-singleton dimension. [Y,Zf] = FILTER(B,A,X,Zi) gives access to initial and final conditions, Zi and Zf, of the delays. Zi is a vector of length MAX(LENGTH(A),LENGTH(B))-1 or an array of such vectors, one for

Lesson 02 - 18

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

each column of X. FILTER(B,A,X,[],DIM) or FILTER(B,A,X,Zi,DIM) operates along the dimension DIM. See also FILTER2, FILTFILT (in the Signal Processing Toolbox).

The MATLAB script is shown below.
n=0:150; //取150个样本时间, F=.1; //信号频率F=0.1赫兹 x1=1+cos(2*pi*n*F)+cos(2*pi*2.5*F*n); //生成输入序列=直流+F赫兹余弦+2.5*F余弦 x=[zeros(1,50), x1, zeros(1,100)];//前后扩展时间轴 Nx=length(x); //计算输入序列长度 nx=0:Nx-1; //生成时间轴 subplot(1,2,1) //选择子图 plot(nx,x) //绘出输入序列x,见图6左边 h=(1/50)*ones(1,50); //构建50阶移动平均滤波器,y(k)=(1/50)*(x(k)+x(k-1)+...+x(k=50)) y=filter(h,1,x); 用filter函数计算输出,A=1,B=h,X=x Ny=length(y); ny=0:Ny-1; //按长度生成时间轴序列 subplot(1,2,2) //选择子图 plot(ny,y) //画出输出Y的图,见图7保留了直流信号,衰减了交流(看稳态区)

and the outcome reported in Figure 6.
3 2.5 2 1.5 1 0.6 0.5 0 -0.5 -1 0.4 0.2 0 1.4 1.2 1

y[k] Steady State = 1.0

x[k]

0.8

Transient Response
0 50 100 150 200 250 300

0

50

100

150

200

250

300

Figure 6: Moving Average FIR experiment.

Example: Comb FIR 梳状滤波器
Suppose an 20th order Comb FIR (using the minus option in Figure 5) is presented with a sampled signal x[k]=1 + cos(2πk/20) + cos(2πk/40). That is, the input consists of a DC term plus two sinusoids of period 20 and 40 respectively. The configured 20th order Comb FIR will subtract signal samples located 20 samples apart. The comb FIR is presented with a sampled signal x[k]=1 + cos(2πk/20) + cos(2πk/40). The difference between the samples (illustrated in Lesson 02 - 19

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

Figure 7), separated 20 sample instances, would be zero for both DC and the period 20 sinusoid. The samples are, however, additive for a period 40 sinusoid. As a result, in an ideal setting, the DC and period 20 signals would be blocked ("nulled") while the period 40 signal would be passed with a doubling of amplitude.

Period 20

Period 40

20 samples

20 samples Figure 7: Comb FIR filter analysis.

The transient and steady-state regimes hold the same position as shown in Figure 7. It can be seen that at steady-state the comb filter output contains no DC or period 20 signal information. The system response was simulated using knowledge of the comb filter’s difference equation and MATLAB. The MATLAB script file is shown below.
n=0:200; //取200个样本, x1=1+cos(2*pi*n/20)+cos(2*pi*n/40); //生成输入离散数字序列 x=[zeros(1,50), x1, zeros(1,100)]; // 前扩展50个时间单位,后扩展100个 Nx=length(x); //计算输入序列长度 nx=0:Nx-1; // 生成前后扩展后的样本时间序列 subplot(1,2,1) //指定子图 plot(nx,x) //绘出输入序列图 见图8左边 h=[1, zeros(1,18), -1]; y=filter(h,1,x); Ny=length(y); ny=0:Ny-1; subplot(1,2,2) plot(ny,y) //够建梳状滤波器,延时20相减y[k]=x[k]-x[k-20] //用filter函数计算输出,A=1,B=h,X=x //计算序列长度 //按长度生成时间轴序列 //选择子图2 //画出输出Y的图,见图7滤出了直流和周期20的信号

The outcome is shown in Figure 8.

Lesson 02 - 20

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

3 2.5 2

3

Transient Response
2

1 1.5 1 -1 0.5 0 -0.5 0 50 100 150 200 250 300 350 -2

x[k]
0

Steady State Average=0

-3

Transient Response
0 50 100 150 200 250 300 350

Figure 8: Comb FIR (with minus option) experiment. Notice that the amplitude of the period 40 signal envelope is doubled.

Challenge 挑战
1st-Order System
A simple first-order RC circuit is show in Figure 1.

Figure 1: RC Circuit The relationship between the input forcing function v(t) and voltage developed across the capacitor, denoted vo(t), is defined by the ordinary differential equation
dv 0 (t ) 1 1 + v 0 (t ) = v (t ) 用电流相等的原理即可得到此式 dt RC RC From basic circuit theory, the transfer function associated with the RC circuit model satisfies:
? 1 ? ? ? ? RC ? 传递函数和冲激响应是一对积分变换,拉普拉斯变换 H (s ) = 1 s+ RC which has a DC gain of H(0)=1. It immediately follows that the system’s impulse response is given by:

Lesson 02 - 21

Digital Signal Processing

? Dr. Fred J. Taylor, Professor

h(t ) =

1 - t/RC e u (t ) RC

What is the equivalent discrete-time model of the circuit’s impulse response based on a sample rate of fs=1 k Sa/s, and RC=10-2? h(k)=Ts*e(k*Ts)

h(k ) = 10?3 ? 102 e ? (100 / 1000) k = 0.1e ?0.1k = 0.1 * 0.9048k
取 N=50,
h(k-m)=

k=0,1……N;

Lesson 02 - 22


赞助商链接
相关文章:
Leadership_Lesson02- Leadership
Leadership_Lesson02- Leadership - Store of Learning 培训店 第二课 领导艺术阅读练习 “How can your company gro...
Lesson02_CreatingSystemConnections_chi
Lesson02_CreatingSystemConnections_chi - 2. 钢结构系统组件 Tekla Structures 13.0 基本培训 2007-4-16 Copyrig...
test_lesson02
test_lesson02 - 测试 --1. 查询工资大于 1200 的员工姓名和工资 --2. 查询员工号为 7654 的员工的姓名和部门号 --3. 选择工资不在 1500 到 20...
Lesson02_CreatingSystemConnections_chi
lesson02 12页 1财富值 lesson02 35页 2财富值 LESSON02 8页 免费 Leadership Lesson02-Lead... 9页 5财富值 Lesson02--属性与方法 暂无评价 26页 免费 Le...
test_lesson02过滤和排序数据
test_lesson02过滤和排序数据 - 测试 1. 查询工资大于 12000 的员工姓名和工资 a) select last_name,salary b) from employees ...
Lesson 02 There be与have has
Lesson 02 There be与have has_英语考试_外语学习_教育专区。英语语法系列练习题 There be 句型与 have, has 的区别一、There be 句型 二、have,has 1.表示:...
小学英语冀教版三年级上册Unit1《lesson02 Boy,Girl an...
小学英语冀教版三年级上册Unit1《lesson02 Boy,Girl and Teacher》优质课公开课教案教师资格证面试试讲教案_三年级英语_英语_小学教育_教育专区。小学英语冀教版三...
新概念二册笔记Lesson 02 Breakfast or lunch
新概念二册笔记Lesson 02 Breakfast or lunch_英语学习_外语学习_教育专区。Lesson 2 Breakfast or lunch? 【New words and expressions】 ★until prep.直到 ...
新概念英语名师精讲笔记第一册lesson01-02
新概念英语名师精讲笔记第一册lesson01-02 - 重点和难点 一般疑问句 词汇 pen n.钢笔 【词组】pen name 笔名 /pen pal 笔友/ pen pusher 文员 pe...
lesson 02 Look at the tree
lesson 02 Look at the tree_英语_小学教育_教育专区。国家地理儿童百科入门级第二课 Look at the tree 课文 lesson 02 Look at the tree Look at the tree...
更多相关标签: