El teorema de muestreo y sus aplicaciones

FACULTAD DE MATEM´ ATICAS DEPARTAMENTO DE AN´ ALISIS MATEM´ ATICO EL TEOREMA DE MUESTREO Y SUS APLICACIONES Memoria realizada por Jos´e Garc´ıa Fern´andez Dirigido por: Dr. D. Renato ´ Alvarez Nodarse Resumen En este trabajo nuestro objeto de estudio ser´a la teor´ıa de se˜nales, el objetivo es llegar al teorema de muestreo de Shannon-Whittaker-Kotelnikov, con su correspondiente demostraci´on, y aplicarlo a distintas situaciones, en particular al problema de recuperaci´on de se˜nales. Para ello necesitaremos trabajar con herramientas de an´alisis Fourier, as´ı que realizaremos un estudio previo sobre la serie de Fourier, con una bre- ve introducci´on hist´orica, adem´as de ejemplos, propiedades y convergencia. Seguiremos el estudio con otra de las herramientas del an´alisis de Fourier, la transformada de Fourier, que ser´a un instrumento fundamental en este trabajo, y al que por supuesto dedicaremos un estudio previo con ejemplos y propiedades fundamentales. Por ´ultimo realizaremos una breve descripci´on de lo que consideraremos se˜nales, veremos algunos ejemplos y tipos de se˜nales que nos podemos encontrar. Para finalizar, tras realizar el estudio te´orico, efectuaremos un estudio num´erico, contrastando nuestro estudio te´orico con ejemplos pr´acticos. Para esto ´ultimo nos apoyaremos en herramientas de software libre como son Maxima u Octave. Todos los programas que se han realizado para esta labor est´an a disposici´on de todos los lectores en: http://euler.us.es/~renato/clases/tfg/jgf/. 4 Abstract The aim of this work is to present an introduction to the Mathematical Signal Theory. We will prove the Shannon-Whitaker-Kotelnikov Sampling Theorem and we will apply it to different situations, in particular to the problem of signal recovering. In order to do that, we will use some tools from the Fourier Analysis. So we will carry out a previous study of the Fourier Series, with a brief historical introduction, examples, properties and convergence Theorems. The study will continue with the Fourier Transform, which will be a key instrument along this work, and, therefore, we will include a previous study including some examples and fundamental properties. Finally, we will briefly describe what we consider to be signals and we present some examples of them. To conclude our analysis, we will a numerical study, contrasting our provide theoretical results with examples. In order to do this, we will use free software tools such as Maxima and Octave. All the programs used along this research are available to all readers and can be found on the web site of the advisor of this work: http://euler.us.es/~renato/clases/tfg/jgf/. ´ Indice general 1. La serie trigonom´etrica de Fourier 9 1.1. Introducci´on ............................ 9 1.2. Definiciones principales y ejemplos ............... 11 1.3. Convergencia puntual de la serie de Fourier ........... 15 1.4. Convergencia uniforme de la serie de Fourier .......... 16 1.5. La serie exponencial de Fourier ................. 17 1.6. La serie de Fourier en un intervalo arbitrario .......... 18 1.7. El f´enomeno de Gibbs ...................... 19 2. La transformada de Fourier 23 2.1. Definiciones principales y ejemplos ............... 23 2.2. Propiedades fundamentales .................... 26 2.3. Identidad de Plancherel-Parseval ................ 27 2.4. Una normalizaci´on diferente ................... 28 3. Se˜nales 31 3.1. Definiciones principales ...................... 31 3.2. Primeros ejemplos ......................... 32 3.3. Tipos de se˜nales .......................... 33 3.4. Teorema de muestreo ....................... 34 3.5. Convergencia ........................... 36 3.6. Cota de error cometido ...................... 37 3.7. El fen´omeno del aliasing ..................... 39 3.8. Muestreo Irregular ........................ 39 4. Programas 41 4.1. Serie de Fourier, paquete fourie ................ 41 4.2. Serie de Fourier, la funci´on aprfou ............... 43 7 8´ INDICE GENERAL 4.3. Transformada de Fourier con Octave .............. 46 4.4. Espectro de una se˜nal con Octave ................ 48 4.5. Espectro de una se˜nal con ruido ................. 49 4.6. Reconstrucci´on de se˜nales con Maxima ............. 50 4.7. Reconstrucci´on de se˜nales con Octave .............. 56 4.8. Aliasing .............................. 58 4.9. Lista de programas ........................ 61 Bibliograf´ıa 63 Cap´ıtulo 1 La serie trigonom´etrica de Fourier En este apartado vamos a estudiar con alg´un detalle la serie trigonom´etrica de Fourier. En particular, veremos algunos teoremas de convergencia puntual y uniforme de las mismas. En adelante, asumiremos que la funci´on fes de cuadrado integrable en [−π, π], es decir, f∈L2[−π, π], aunque la mayor´ıa de los resultados son ciertos para fabsolutamente integrable. N´otese que el sistema de funciones {1}∪{sin (nx),cos (nx)}∞ n=1 es un sistema ortogonal dos a dos respecto al producto escalar: hf, gi=Zπ −π f(x)g(x)dx. 1.1. Introducci´on El 21 de diciembre de 1807, Joseph Fourier, present´o al Institut de France una memoria titulada M´emoire sur la propagation de la chaleur dans les corps solides. Cuatro miembros, uno m´as de lo que era habitual, Lagrange, Laplace, Lacroix y Monge, fueron designados para emitir un informe, que nunca se lleg´o a escribir a pesar de las insistencias de Fourier, que deseaba un juicio sobre su trabajo. A cambio, el Instituto lo propuso como tema para el premio que se deb´ıa otorgar en 1812 (la teor´ıa matem´atica de las leyes de propagaci´on del calor y comparar los resultados de esta teor´ıa con la de los experimentos). A finales de 1811 Fourier entreg´o una nueva memoria como 9 16 CAP´ ITULO 1. LA SERIE TRIGONOM ´ ETRICA DE FOURIER 1.4. Convergencia uniforme de la serie de Fourier Teorema 1.4.1 Sea f:R→Runa funci´on continua y 2π-peri´odica y supongamos que fes casi-continuamente derivable en [−π, π]. Entonces: 1. La serie de Fourier f0se obtiene derivando t´ermino a t´ermino Sf(x). 2. P∞ k=1 |ak|+|bk|<+∞. 3. Sf(x) converge uniformemente a fen [−π, π]. Demostraci´on: Calculemos los coeficientes an(f0) y bn(f0) de f0en funci´on de los coeficientes an(f) y bn(f) de f. ak(f) = 1 πZπ −π f(x) cos(kx)dx =−1 πk Zπ −π f0(x) cos(kx)dx =−1 kbk(f0), bk(f) = 1 πZπ −π f(x) sin(kx)dx =1 πk Zπ −π f0(x) sin(kx)dx =−1 kak(f0), por tanto la serie de Fourier de f0es: Sf0(x) = ∞ X k=1 −kak(f) sin(kx) + kbk(f) cos(kx), de donde se deduce 1. Para probar 2 usamos que: |ak(f)|=1 kbk(f0)≤1 21 k2+b2 k(f0), |bk(f)|=1 kak(f0)≤1 21 k2+a2 k(f0). En la desigualdad anterior hemos usado que (|a|−|b|)2≥0. Como f0es casi-derivable, tenemos por la desigualdad de Bessel [3]: 1.5. LA SERIE EXPONENCIAL DE FOURIER 17 ∞ X k=1 ak(f0)2+bk(f0)2≤1 πZπ −πf0(x)2dx < +∞,luego ∞ X k=1 |ak(f)|+|bk(f)| ≤ ∞ X k=1 1 k2+1 2ak(f0)2+bk(f0)2<+∞. De la desigualdad anterior aplicando el Criterio de Weierstrass para series de funciones vemos la convergencia uniforme de Sf(x).  Corolario 1.4.2 Sea f:R→Runa funci´on continua y 2π-peri´odica y supongamos que fes casi-continuamente derivable en [−π, π]. Entonces Sf(x) converge absolutamente e uniformemente a fen [−π, π]. 1.5. La serie exponencial de Fourier Sea f:R→C, y sean las funciones: en(x) = einx, n ∈Z. N´otese que: Zπ −π einxe−imxdx =Zπ −π ei(n−m)xdx =f(x) =      2π, si n=m, 2 sin(n−m)π n−m,si n6=m. Es decir, el conjunto (en)nes un conjunto ortogonal respecto al producto escalar: hf, gi=Zπ −π f(x)g(x)dx. Definici´on 1.5.1 Dada una funci´on f:R→C, de cuadrado integrable en [-π,π] y 2π-peri´odica definiremos la serie exponencial de Fourier por: Sf(x) = ∞ X n=−∞ cneinx,(1.5.1) donde los coeficientes vienen dados por la siguiente expresi´on: 18 CAP´ ITULO 1. LA SERIE TRIGONOM ´ ETRICA DE FOURIER cn=1 2πZπ −π f(x)e−inxdx. Obviamente esta serie es equivalente a la serie de (1.2.1). Para ello basta usar la f´ormula de euler einx = cos(nx) + sin i(nx) lo que nos resulta: a0= 2c0, an=cn+c−n, bn=i(cn+c−n). cn=a0 2, cn=an−ibn 2, c−n=an+ibn 2. Adem´as, de las relaciones anteriores podemos deducir que cn=c−n. Ejemplo 1.5.2 Encontrar la serie exponencial de Fourier de la funci´on: f(x) =    Asi 0 ≤x≤π, 0 si −π≤x≤0. Vemos claramente que c0es A/2. Para n6= 0 tenemos: cn=A 2πZπ 0 e−inxdx =A(1 −e−inπ) 2πin =     0 si npar, A πni si nimpar. 1.6. La serie de Fourier en un intervalo arbitrario Si en vez de trabajar con funciones 2π-peri´odicas, trabajamos con periodos arbitrarios, por ejemplo 2ltendremos que usar los siguientes sistemas, que son ortogonales y completos en [−l, l] : 1,sin kπx l,cos kπx ln∈N , O en la forma exponencial: nekπxi lo.n∈N 1.7. EL F ´ ENOMENO DE GIBBS 19 As´ı, tenemos las series equivalentes: Sf(x) = ∞ X n=−∞ cnenπix l,donde cn=1 2πZπ −π f(x)e−nπxi ldx, n ∈Z. Sf(x) = a0 2+ n X k=1 akcos nπx l+bksin nπx l, an=1 πZπ −π f(x) cos nπx ldx, n = 0,1,2, ..., bn=1 πZπ −π f(x) sin nπx ldx, n = 0,1,2, .... 1.7. El f´enomeno de Gibbs Vamos a ver un ejemplo interesante. Consideramos la funci´on signo(x) en [−π, π] y extendamos peri´odicamente la funci´on a todo R. Su serie de Fourier es: signo(x) = 4 π ∞ X k=0 sin(2k+ 1)x 2k+ 1 . Recordamos que: signo(x) =            −1 si x < 0, 0 si x= 0, 1 si x > 1. Vemos que adem´as converge puntualmente en (−1/2,1/2) y lo hace uniformemente en cualquier compacto [δ, 1/2−δ]⊂(0,1/2). Obviamente la convergencia no puede ser uniforme en (−1/2,1/2) pues signo(x) no es continua en dicho intervalo. Si dibujamos las gr´afica de signo de xy de algunas de sus sumas parciales vemos que cerca del cero aparece un peque˜no pico que no parece disminuir incluso aumentando el n´umero de sumados en las sumas parciales. 20 CAP´ ITULO 1. LA SERIE TRIGONOM ´ ETRICA DE FOURIER Figura 1.4: Serie de Fourier con 20 t´erminos de signo(x). Figura 1.5: Serie de Fourier con 200 t´erminos de signo(x). De hecho, al aumentar nvemos en un entorno de x= 0 que el pico no disminuye, sino que mantiene su altura. 1.7. EL F ´ ENOMENO DE GIBBS 21 Para aclarar lo que ocurre calculemos el m´aximo de la sumas S2n−1f: S2n−1f(x) = 4 π n−1 X k=0 sin 2π(2k+ 1)x 2k+ 1 , que se obtienen al tomar la derivada en x: d dxS2n−1f(x) = 8 n−1 X k=0 cos 2π(2k+ 1)x=4 sin 4nπx sin 2πx .(1.7.1) La derivada se anula en xk=k 4n, k = 1, ..., 2n−1,y podemos comprobar que: d dxS2n−1fk 4n= 16πcsc k 2ncos(kπ), Por tanto, en k= 1 tenemos un m´aximo mientras que en k= 2 un m´ınimo, y as´ı, sucesivamente. Adem´as de (1.7.1) tenemos que, al ser S2n−1f(0) = 0, S2n−1f(x) = 4 Zx 0 sin 4nπz sin 2πz dz. A partir de esta expresi´on, podemos comprobar que precisamente el primer m´aximo es el m´aximo global de la suma parcial, es decir, m´ax x∈[0,δ]S2n−1f(x) = S2n−1fk 4n. Ahora bien, S2n−1fk 4n=4 π n−1 X k=0 sin π(2k+ 1) k 2n 2k+ 1 = n X k=0 2 π sin πξk ξk 1 n, ξk=2k+ 1 2n. N´otese que la ´ultima suma, constituye una suma de Riemann de la integral Z1 0 2 π sin πx xdx, en la partici´on xk=k/n, k = 0,1, ..., n, al evaluar la funci´on integrando en ξk=xk+xk+1 2. Luego : 22 CAP´ ITULO 1. LA SERIE TRIGONOM ´ ETRICA DE FOURIER l´ımn→∞S2n−1fk 4n= l´ımn→∞ n X k=0 2 π sin πξk ξk 1 n=Z1 0 2 π sin πx xdx = 1,17897... As´ı, el m´aximo de las sumas parciales de la serie signo(x) no tiende a 1 sino a 1,178979744..., es decir: m´ax signo(x)−m´ax l´ım n→∞Snf(x)= 1,17897... Lo anterior es generalizable a cualquier funci´on fcon una discontinuidad en alg´un punto. Supongamos que en cierto intervalo I= [α, β] la funci´on tiene un punto de discontinuidad en x0y sea µ=f(x+ 0)−f(x− 0) 2y supongamos que fes tal que la serie converge uniformemente a ambos lados de x0. Entonces la funci´on g(x) = f(x)−µsigno(x) es continua en x0y su serie converge uniformemente en un entorno de x0, por lo tanto el comportamiento de fen un entorno de x0es el mismo que el de µsigno(x), es decir m´ax f(x)−m´ax l´ım n→∞Snf(x)=µ1,17897... El efecto es v´alido para toda funci´on discontinua se conoce como efecto de Gibbs quien estudi´o el problema en 1899. Cap´ıtulo 2 La transformada de Fourier En este cap´ıtulo, introduciremos la transformada de Fourier, que va a tener una gran importancia a la hora de demostrar el teorema de muestreo como ya veremos m´as adelante. Intuitivamente y adelant´andonos un poco, la transformada de Fourier transforma una se˜nal que depende del tiempo en una se˜nal que depende de la frecuencia. 2.1. Definiciones principales y ejemplos En adelante, asumiremos que f∈L1(R). Definici´on 2.1.1 Dada una funci´on f:R→C, definiremos su transformada de Fourier que denotaremos por b f(λ) o F[f] a la funci´on: b f(λ) = F[f](λ) := 1 2πZ∞ −∞ f(x)e−iλxdx. (2.1.1) Por comodidad definiremos tambi´en la transformaci´on: F[f](λ) := 1 2πZ∞ −∞ f(x)eiλxdx. (2.1.2) Bajo ciertas condiciones sobre f, 2πF[f](λ) es la inversa de F[f](λ), y se verifica: f(x) = 2πF[b f](x), F−1[f] = 2πF[b f]. 23 24 CAP´ ITULO 2. LA TRANSFORMADA DE FOURIER Definici´on 2.1.2 A la funci´on F(λ)=b f(λ), se le denomina espectro de la funci´on f. Ejemplos Calcular las transformadas de Fourier de las siguientes funciones. 1.f:R→R,f(x) =    1 si |x|<1, 0 si |x| ≥ 1. F[f](λ) = 1 2πZ∞ −∞ e−iλxdx =1 2πZ1 −1 (cos(λx)−isin(λx))dx =sin(λ) πλ . Figura 2.1: Pulso rectangular (arriba) y su transformada (abajo). 2.1. DEFINICIONES PRINCIPALES Y EJEMPLOS 25 2.h:R→R,h(x) = e−x2/2. b h(λ) = F[h](λ) = 1 2πZ∞ −∞ e−x2/2e−iλxdx =Z∞ −∞ e−x2/2cos(λx)dx. (2.1.3) Entonces, ∂b h(λ) ∂λ =1 2πZ∞ −∞ e−x2/2(−x) sin(λx)dx =−λ 2πZ∞ −∞ e−x2/2cos(λx)dx =−λb h(λ). En la segunda igualdad hemos integrado por partes, luego: b h(λ) = Ce−λ2/2⇒C=b h(0) = 1 2πZ∞ −∞ e−x2/2dx =1 √2π⇒F[h](λ) = 1 √2πe−λ2/2. Figura 2.2: Gaussiana (arriba) y su transformada (abajo). 32 CAP´ ITULO 3. SE ˜ NALES 3.2. Primeros ejemplos Escal´on de Heaviside: Heaviside (t) =    1 si t≥0, 0 si t < 0. Pulso rectangular: Rect(t) =    1 si |t| ≤ 1, 0 si |t|>1. Figura 3.1: Escal´on de Heaviside (izquierda) y pulso Rectangular Rect(t) (derecha). Se˜nal sinusoidal: f(t) = Asin(wt +ϕ), w > 0. Se˜nal con ruido: f(t) = Asin(wt +ϕ) + a0random(t), w > 0. 3.3. TIPOS DE SE ˜ NALES 33 Figura 3.2: Se˜nal sinusoidal f(t) = cos(20πt) + sin(40πt) (izquierda) y con ruido f(t) = cos(20πt) + sin(40πt) + 2 random(t) (derecha). 3.3. Tipos de se˜nales En adelante, denotaremos por Xal espacio de las se˜nales. Uno de los objetivos del an´alisis de se˜nales, es poder distinguir las se˜nales, trasmitirlas, codificarlas, etc. Las funciones anteriores pueden representar el voltaje o la intensidad de corriente en un circuito, intensidad de la m´usica en un reproductor de CD, potencia de un altavoz, etc. Una se˜nal que sea casi-continua, se dice que es anal´ogica. Si la se˜nal es un conjunto de valores dados s1, s2, ..., sN, ... se dice que la se˜nal es digital. As´ı: f:R→R, f ∈C(R) =⇒se˜nal anal´ogica. Si en vez de ftenemos sus muestras tomadas con una frecuencia 1/a, diremos que la se˜nal es digital: (f(xk))k:xk=ka, k ∈Z=⇒se˜nal digital. En la naturaleza las se˜nales suelen ser en general, anal´ogicas. Uno de los principales problemas consiste en convertir una se˜nal anal´ogica en digital (discreta) para, por ejemplo, transmitirla (fotos digitales), grabarla (CD), y luego poder recuperarla sin perder informaci´on. El proceso de conversi´on se denomina “muestreo”de la se˜nal. Nos interesar´a el problema de la recuperaci´on de una se˜nal a partir de sus muestras. 34 CAP´ ITULO 3. SE ˜ NALES 3.4. Teorema de muestreo Vamos a enunciar y posteriormente probar, el resultado fundamental de este trabajo: el teorema de muestreo de Whittaker-Kotelnikov-Shannon, que como veremos posteriormente tiene una gran utilidad en la teor´ıa de se˜nales. Teorema 3.4.1 Sea una se˜nal x(t) de energ´ıa finita y de banda limitada a B, entonces: x(t) = 1 π ∞ X n=−∞ xn 2Bsin (π(2Bt −n)) 2Bt −n.(3.4.1) Es decir, para recuperar una se˜nal banda limitada a Bbasta conocer sus muestras tomadas con una frecuencia igual a 2B. Dicha frecuencia se denomina frecuencia de Nyquist. Demostraci´on Primero vamos a considerar la funci´on f(u) = e2πiut en el intervalo (−B, B). Vamos a extenderla peri´odicamente y a calcular su serie exponencial de Fourier en la base e2πi(n/2B)un≥0,que es una base ortogonal de L2(−B, B). Nuestra funci´on es continua y derivable luego su serie de Fourier ser´a uniformemente convergente como hemos visto en el apartado de convergencia. f(u) = ∞ X n=−∞ cne2πi(n/2B)u,(3.4.2) donde: cn=1 2BZB −B e2πiut(e−2πi(n/2B)u)du =1 2BZB −B e2πi(t−n/2B)udu = 1 2BZB −B cos 2πt−n 2Budu +i1 2BZB −B sin 2πt−n 2Budu = =1 BZB 0 cos 2πt−n 2Budu =1 B sin (2πt−n BB 2πt−n B  B 0 = 3.4. TEOREMA DE MUESTREO 35 =sin (π(2tB −n) π(2Bt −n). Vemos que: i1 2BZB −B sin (2πt−n 2Bu)du = 0, ya que es una integral en la cual el integrando es una funci´on impar en un dominio de integraci´on sim´etrico. En adelante usaremos la notaci´on sinc(x) = sin (πx) πx ,que denota a la funci´on seno cardinal. Recordamos que definimos la trasformada de Fourier y su inversa como: X(u) = Z∞ −∞ x(t)e−2πiutdt, x(t) = ZB −B X(u)e2πiutdu. Sustituimos e2πiut por su serie de Fourier: x(t) = ZB −B ∞ X n=−∞ X(u)eiπnu/B sinc(2Bt −n)du = = ∞ X n=−∞ ZB −B X(u)eiπnu/B sinc(2Bt −n)du = = ∞ X n=−∞ sinc(2Bt −n)ZB −B X(u)eiπnu/Bdu = ∞ X n=−∞ sinc(2Bt −n)xn 2B= =1 π ∞ X −∞ xn 2Bsin (π(2Bt −n)) 2Bt −n.  Nota: La integraci´on t´ermino a t´ermino, la podemos realizar debido a que nuestra serie es uniformemente convergente. La igualdad se entiende en el sentido de L2. 36 CAP´ ITULO 3. SE ˜ NALES Nota: En teor´ıa de se˜nales se define la serie cardinal de la se˜nal xcomo: 1 π ∞ X n=−∞ xn 2Bsin (π(2Bt −n)) 2Bt −n= (3.4.3) =1 π ∞ X n=−∞ x(na) sin πt a−n t a−n =1 π n=∞ X n=−∞ x(na) sin π a(t−na) π a(t−na). Para trabajar con (3.4.3) de una manera m´as c´omoda, hemos realizado un cambio de notaci´on llamando a=1 2B. 3.5. Convergencia La serie cardinal converge en L2. Vamos ahora a probar que nuestra convergencia no solo es en L2, sino que tenemos convergencia uniforme. Para ello, probaremos un lema auxiliar y posteriormente que la serie cardinal converge uniformemente. Lema 3.5.1   X N≤|n|≤M  sin π a(t−na) π a(t−a) 2  = 1. Demostraci´on Consideremos la funci´on k(t) = e2πiλtχ[−1 2a,1 2a], λ ∈R, y extend´amosla peri´odicamente a todo t∈R.Su serie de Fourier es: k(t) = X n∈R cne2πinat, cn=aZ1 2a −1 2a e2πiλte−2πinatdt = sin π a(λ−na) π a(λ−na), adem´as (cn)n∈L2.Si ahora usamos la identidad de Parseval obtenemos: 3.6. COTA DE ERROR COMETIDO 37 X n∈R|cn|2=aZ1 2a −1 2ae2πiλt2dt = 1. Teorema 3.5.2 La serie cardinal (3.4.1) converge uniformemente. Demostraci´on Sean M, N ∈N, N ≤M: X |n|≤M x(na) sin π a(t−na) π a(t−na)−X |n|≤N x(na) sin π a(t−na) π a(t−na)≤ ≤X N≤|n|≤M|x(na)| sin π a(t−na) π a(t−na)≤ ≤ X N≤|n|≤M|x(na)|2  1/2  X N≤|n|≤M  sin π a(t−na) π a(t−a) 2   1/2 . Por el lema anterior, la serie   X N≤|n|≤M  sin π a(t−na) π a(t−a) 2  est´a acotada, luego: ≤ X N≤|n|≤M|x(na)|2  1/2  X N≤|n|≤M  sin π a(t−na) π a(t−a) 2   1/2 < ε, ya que si MyNson suficientemente grandes, la serie Pn∈Z|x(na)|2<∞. Por tanto, por el criterio de Cauchy para la convergencia uniforme, la serie converge uniformemente en R. 3.6. Cota de error cometido Es importante recordar que cuando hagamos los correspondientes programas no podremos tomar la serie infinita, sino una serie truncada, con lo cual nos ser´a de inter´es una cota del error que cometemos. Para ello, enunciamos el siguiente teorema. 38 CAP´ ITULO 3. SE ˜ NALES Teorema 3.6.1 Sea h(t) una se˜nal con energ´ıa finita y de banda limitada a B, se verifica que: |h(t)| ≤ rBE π.(3.6.1) Demostraci´on Sea H(w) la transformada de Fourier de la funci´on h(t). Ahora consideramos la inversa de la transformada de Fourier de h: h(t) = 1 2πZB −B H(w)eiwtdw , y aplicamos la desigualdad de Cauchy- Schwarz: |h(t)|2≤1 4π2ZB −B|H(w)|2dw ZB −Beiwt2dw = =B 2π2ZB −B|H(w)|2dw =B π1 2πZB −B|H(w)|2dw= =B πZB −B|h(w)|2dw =B πE⇒ |h(t)| ≤ rBE π.  Nota En la ´ultima igualdad, hemos usado la f´ormula de Parseval. Corolario 3.6.2 Tomando h(t) = x(t)−xN(t),es inmediato de (3.6.1) que obtenemos la siguiente cota del error: |x(t)−xN(t)| ≤ rBE π. Vemos como se cumplen las condiciones del teorema anterior, dado que hes una funci´on de L2al ser una resta de dos funciones de L2, y adem´as hes una se˜nal de banda limitada al ser resta de funciones de banda limitada (esta propiedad se obtiene gracias a la linealidad de la integral). 3.7. EL FEN ´ OMENO DEL ALIASING 39 3.7. El fen´omeno del aliasing El aliasing es el efecto que causa que se˜nales continuas distintas no sean distinguibles cuando realizamos un muestreo digital. En este caso, no podemos recuperar la se˜nal original de forma ´unica a partir de la se˜nal digital. En otras palabras, tenemos se˜nales distintas, pero que coinciden en las muestras tomadas, por lo que la recuperaci´on de nuestra se˜nal no es posible. Nos encontramos con este problema cuando nuestras muestras no son las “indicadas”, es decir, cuando no verifican que νs= 1/a ≥2B, donde νs representa nuestra frecuencia de muestreo. Veamos un ejemplo de aliasing. Consideremos x(t) = sin(2πt): t Muestras Reconstrucción Original -3 -2 -1 0 1 2 3 0 0.5 1 1.5 2 t Muestras Reconstrucción Original -3 -2 -1 0 1 2 3 0 0.5 1 1.5 2 Figura 3.3: Recuperaci´on de nuestra se˜nal tomando como frecuencia de muestreo 4 (izquierda) y 1 (derecha). 3.8. Muestreo Irregular Es interesante resaltar que hemos estado trabajando en nuestra teor´ıa de se˜nales con lo que se conoce como muestreo regular, es decir, que nuestras muestras se toman de forma equidistante. Hay diferentes tipos de muestreo, entre los cuales uno de los m´as destacamos es el muestreo irregular. ¿Cu´ales son las ventajas y desventajas con respecto al muestreo regular? 40 CAP´ ITULO 3. SE ˜ NALES Como desventaja, es muy importante el hecho de que el Teorema de Muestreo no es aplicable, y para la reconstrucci´on tenemos que usar m´etodos que son computacionalmente intensivos, como por ejemplo la transformada r´apida de Fourier, que en el caso del muestreo irregular no es ni siquiera invertible. Como ventaja, en el muestreo irregular es menos probable que se produzca el fen´omeno del aliasing, que es uno de los inconvenientes m´as importantes del muestreo regular. Para m´as informaci´on acerca del muestreo irregular, consultar [5] y [8] en la bibliograf´ıa. Cap´ıtulo 4 Programas En este cap´ıtulo mostraremos algunos de los programas que a lo largo del trabajo hemos usado y expondremos algunos ejemplos nuevos de inter´es. Comenzaremos mostrando varias formas de construir las series de Fourier. En primer lugar, mostraremos un ejemplo donde podemos realizar los c´alculos de manera anal´ıtica, usando un programa simb´olico (en este caso de Maxima), tambi´en veremos alg´un programa que realice los c´alculos num´ericamente. Adem´as de la transformada de Fourier, veremos programas de recuperaci´on de se˜nales y el problema del aliasing. Para consultar los comandos utilizados tanto en Maxima como en Octave, ver referencias [1] y [10]. 4.1. Serie de Fourier, paquete fourie Vamos a usar Maxima para el caso de f(x) = xen [−π, π], extendida peri´odicamente a todo R. Ante todo debemos cargar el paquete fourie, que nos permite trabajar simb´olicamente con series de Fourier. (%i1) kill(all)$ (%i1) load(fourie)$ Definimos la funci´on que vamos a desarrollar en serie, en este caso f(x) = x. (%i2) f1(x) := x; (%o2) f1(x):= x 41 48 CAP´ ITULO 4. PROGRAMAS 4.4. Espectro de una se˜nal con Octave Este programa nos ser´a muy ´util a la hora de reconstruir se˜nales, ya que podemos ver el espectro de una se˜nal mediante su transformada. El programa es el siguiente: T=1; % periodo de muestreo fs=256; % frecuencia de muestreo t=(0:1/fs:T); % puntos de muestreo N=length(t); f=sin(5*2*pi*t); % se~nal original figure, plot(t,f), ylim=[-1.5,1.5], xlabel(’Tiempo’), title(’Se~nal original’) F=fft(f)/sqrt(N); % trasformada de Fourier de f % Para dibujarlo, despreciamos la % mitad del dominio debido a la simetr´ıa omega=0.5*fs*linspace(0,1,floor(N/2)+1); % vector de % frecuencias discretas range=(1:floor(N/2)+1); % rango del espectro de potencia P=abs(F(range)).^2; % espectro de potencia de la se~nal f figure, plot(omega,P),xlim([0,10]), xlabel(’Frecuencia’), title(’Espectro de la se~nal f’) 4.5. ESPECTRO DE UNA SE ˜ NAL CON RUIDO 49 Figura 4.5: Se˜nal original f(t) = sin(10πt) (izquierda). Vemos que la frecuencia de Nyquist se alcanza en 6 (derecha). 4.5. Espectro de una se˜nal con ruido Vamos a realizar un programa similar al anterior, pero con una funci´on con ruido: T=1; % periodo de muestreo fs=2^9; % frecuencia de muestreo t=(0:1/fs:T); % puntos de muestreo N=length(t); f=cos(10*2*pi*t)+cos(20*2*pi*t)+4*rand(size(t)); % se~nal original figure, plot(t,f), ylim([-4,8.5]), xlabel(’Tiempo’), title(’Se~nal original’) F=fft(f)/sqrt(N); % trasformada de Fourier de f % Para dibujarlo, despreciamos % la mitad del dominio debido % a la simetr´ıa omega=0.5*fs*linspace(0,1,floor(N/2)+1); % vector de frecuencias discretas range=(1:floor(N/2)+1); % rango del espectro de potencia P=abs(F(range)).^2; % espectro de potencia de la se~nal f 50 CAP´ ITULO 4. PROGRAMAS figure, plot(omega,P), xlabel(’Frecuencia’), title(’Espectro de la se~nal f’) Figura 4.6: Se˜nal original f(t) = 2 cos(20πt) + 3 sin(20πt) + 4 random(t) (izquierda) y espectro de nuestra se˜nal (derecha). 4.6. Reconstrucci´on de se˜nales con Maxima Una de las utilidades m´as destacadas del teorema de muestreo es la reconstrucci´on de se˜nales. La idea es que una se˜nal anal´ogica puede ser reconstruirla a trav´es de sus muestras. Como consecuencia del teorema de muestreo 3.4.1, para poder recuperar nuestra se˜nal banda limitada a Bse debe verificar que nuestra frecuencia de muestreo νs,sea mayor o igual que el doble de la frecuencia de Nyquist, es decir: νs=1 a≥2B. Veamos un ejemplo con la se˜nal x(t) = sin(2πt). Si aplicamos el apartado anterior a nuestra se˜nal para ver su frecuencia obtenemos: 4.6. RECONSTRUCCI ´ ON DE SE ˜ NALES CON MAXIMA 51 T=1; % periodo de muestreo fs=256; % frecuencia de muestreo t=(0:1/fs:T); % puntos de muestreo N=length(t); f=sin(2*pi*t); % se~nal original figure, plot(t,f), ylim([-1.5,1.5]), xlabel(’Tiempo’), title(’Se~nal original’) F=fft(f)/sqrt(N); % trasformada de Fourier de f % Para dibujarlo, despreciamos la % mitad del dominio debido a la simetr´ıa omega=0.5*fs*linspace(0,1,floor(N/2)+1); % vector de frecuencias discretas range=(1:floor(N/2)+1); % rango del espectro de potencia P=abs(F(range)).^2; % espectro de potencia de la se~nal f figure, plot(omega,P), xlim([-0.5 2.5]), xlabel(’Frecuencia’), title(’Espectro de la se~nal f’) 52 CAP´ ITULO 4. PROGRAMAS Figura 4.7: x(t) = sin(2πt) (izquierda). Espectro de x(t) = sin(2πt) (derecha). Vemos como la frecuencia de Nyquist es 2, as´ı que tomamos como frecuencia de muestreo 4, es decir, a= 1/4. Vamos a realizar un programa en Maxima que tome unas muestras de nuestra se˜nal y, utilizando el teorema de muestreo, la reconstruya. El programa es el siguiente: define(f(t),sin(2*%pi*t)); NN:10$ a:1/4$ tt:makelist( (k)*a , k,0,NN-1); x:makelist(f((k)*a),k,0,NN-1); dis:makelist( [tt[k],x[k]], k,1,NN); ll:length(x); sum(x[n]*sin((%pi/a)*(t-(n-1)*a))/((%pi/a)*(t-(n-1)*a )),n,1,ll); define(fr(t),%); plot2d([[discrete,dis],fr(t),f(t)],[t,0,2], [style,points,lines, lines],[y,-3,3], [legend,"Muestras","Reconstrucci´on","Original"], [gnuplot_pdf_term_command, "set term pdfcairo color solid lw 2 size 15 cm, 15 cm font \",24\""], [pdf_file, "C:/Users/PC/AppData/Local/Temp/recuperacion.pdf"], [style,points,[lines,3.5],[lines,3.5]]); 4.6. RECONSTRUCCI ´ ON DE SE ˜ NALES CON MAXIMA 53 t Muestras Reconstrucción Original -3 -2 -1 0 1 2 3 0 0.5 1 1.5 2 Figura 4.8: Recuperaci´on de nuestra se˜nal a partir de sus muestras. Ahora realizamos una reconstrucci´on con dos arm´onicos. Sea x(t) = sin(2πt) + cos(3πt). Aplicando los programas anteriores obtenemos: T=1; % periodo de muestreo fs=256; % frecuencia de muestreo t=(0:1/fs:T); % puntos de muestreo N=length(t); f=sin(2*pi*t)+cos(3*pi*t); % se~nal original figure, plot(t,f), ylim([-1.5,1.5]), xlabel(’Tiempo’), title(’Se~nal original’) 54 CAP´ ITULO 4. PROGRAMAS F=fft(f)/sqrt(N); % trasformada de Fourier de f % Para dibujarlo, despreciamos la % mitad del dominio debido a la simetr´ıa omega=0.5*fs*linspace(0,1,floor(N/2)+1); % vector de frecuencias discretas range=(1:floor(N/2)+1); % rango del espectro de potencia P=abs(F(range)).^2; % espectro de potencia de la se~nal f figure, plot(omega,P), xlim([-0.5 12.5]),,ylim([0,38]), xlabel(’Frecuencia’), title(’Espectro de la se~nal f’) Figura 4.9: x(t) = sin(2πt)+cos(3πt) (izquierda). Espectro de x(t) (derecha). La frecuencia de Nyquist est´a en torno a 8, as´ı que podemos tomar por ejemplo, a= 1/16. Aplicando el programa anterior obtenemos: define(f(t),sin(2*%pi*t)+cos(3*%pi*t)); NN:30$ a:1/16$ 4.6. RECONSTRUCCI ´ ON DE SE ˜ NALES CON MAXIMA 55 tt:makelist( (k)*a , k,0,NN-1); x:makelist(f((k)*a),k,0,NN-1); dis:makelist( [tt[k],x[k]], k,1,NN); ll:length(x); sum(x[n]*sin((%pi/a)*(t-(n-1)*a))/((%pi/a)*(t-(n-1)*a )),n,1,ll); define(fr(t),%); plot2d([[discrete,dis],fr(t),f(t)],[t,0,2], [style,points,lines, lines],[y,-3,3], [legend,"Muestras","Reconstrucci´on","Original"], [gnuplot_pdf_term_command, "set term pdfcairo color solid lw 2 size 15 cm, 15 cm font \",24\""], [pdf_file, "C:/Users/PC/AppData/Local/Temp/recuperacion.pdf"], [style,points,[lines,3.5],[lines,3.5]]); t Muestras Reconstrucción Original -3 -2 -1 0 1 2 3 0 0.2 0.4 0.6 0.8 1 Figura 4.10: Recuperaci´on de nuestra se˜nal a partir de sus muestras. 56 CAP´ ITULO 4. PROGRAMAS 4.7. Reconstrucci´on de se˜nales con Octave Realizamos un programa similar al anterior, pero esta vez usando Octave: Vamos a reconstruir x(t) = 2 cos(10πt) : T=0.1; f=5; no=0:0.0001:1; xo=cos(2*pi*f*no); n=(0:T:1)’; xs=cos(2*pi*f*n); t=linspace(0,1,1/0.0001)’; ya=sinc((1/T)*t(:,ones(size(n)))-(1/T)*n(:,ones(size(t)))’)*xs; plot(n,xs,’or’,t,ya,’r’,no,xo,’b’);grid; legend(’Muestras’,’Reconstruida’,’Original’); xlabel(’Tiempo’); axis([0 1 -1.9 1.9]); Figura 4.11: Recuperaci´on de nuestra se˜nal a partir de sus muestras. 4.7. RECONSTRUCCI ´ ON DE SE ˜ NALES CON OCTAVE 57 Por ´ultimo, vamos a realizar un ejemplo m´as complejo con nuestro programa en Octave. Para ello, sea x(t) = cos(2πt) + sin(68πt)−cos(5πt). T=0.01; f=1; no=0:0.0001:10; xo=cos(2*pi*f*no)+sin(34*pi*2*f*no)-cos(5*pi*f*no); n=(0:T:1)’; xs=cos(2*pi*f*n) +sin(34*pi*2*f*n)-cos(5*pi*f*n); t=linspace(0,1,1/0.0001)’; ya=sinc((1/T)*t(:,ones(size(n)))-(1/T)*n(:,ones(size(t)))’)*xs; plot(n,xs,’or’,t,ya,’r’,no,xo,’b’);grid; legend(’Muestras’,’Reconstruida’,’Original’); xlabel(’Tiempo’); axis([0 1 -3.3 3.3]); Figura 4.12: Recuperaci´on de nuestra se˜nal a partir de sus muestras. 64 BIBLIOGRAF´ IA [9] Athanasios Papoulis. Signal Analysis. McGraw-Hill Companies City, 1977. [10] Jos´e Mar´ıa Valiente Cifuentes. Manual de iniciaci´on a GNU Octave. E.U. Polit´ecnica de Teruel, 2006. http://softlibre.unizar.es/manuales/aplicaciones/octave/ manual_octave.pdf.