FFT de una señal coseno

Señal de entrada: x(t)= cos(ω0 t)
ω0 = (2π fs) en radianes, para el ejemplo fs=20 Hz

x(ω) =π[δ(ω-ω0)+δ(ω+ω0)] con eje en radianes.
La gráfica se muestra en Hz. fs0/2π
Se muestran las partes:
– real e imaginaria de la transformada y
– magnitud y fase

Script de python para usar la transformada rapida de Fourier.
En la seccion de función de entrada, se muestra para cuatro funciones: coseno, constante, impulso y sinc. Para usar una de ellas se elimina el simbolo de comentario (##)

# Transformadas Rapida de Fourier
# entrada funcion par
# Propuesta: edelros@espol.edu.ec
import numpy as np
import scipy.fftpack as fourier
import matplotlib.pyplot as plt

# Definir la funcion de Entrada par
def entradax(t,fs):
    # función matemática CAMBIAR AQUI
    x=np.cos(2*np.pi*fs*t)

##    x=1  #constante
    
##    x=0   # impulso
##    if abs(t)<(1/fs):
##        x=1

##    if t==0:  #para sinc aumentar k=8
##        t=(1/fs)/100000
##    x=(np.sin(2*np.pi*fs*t)/(np.pi*t))

    # función matemática CAMBIAR AQUI
    return(x)

# PROGRAMA
# INGRESO
fs=20       #float(input('frecuencia Señal Hz:'))
ancho=1     #ancho del intervalo
t0=0        #intervalo de t
tn=ancho*1/fs
resolucion=1    #resolucion de transformada

# PROCEDIMIENTO
n= 1024*resolucion 
dt=(tn-t0)/n   # intervalo de muestreo
# Analógica Referencia para mostrar que es par
t=np.arange(-tn,tn,dt)  # eje tiempo analógica
m=len(t)
xanalog=np.zeros(m, dtype=float)
for i in range(0,m):
    xanalog[i]=entradax(t[i],fs)

# FFT: Transformada Rapida de Fourier
# Analiza la parte t>=0 de xnalog muestras[n:2*n]
# Transformada de Fourier
xf=fourier.fft(xanalog[n:2*n])
xf=fourier.fftshift(xf)
# Rango de frecuencia para eje
frq=fourier.fftfreq(n, dt)
frq=fourier.fftshift(frq)

# x[w] real
xfreal=(1/n)*np.real(xf)
# x[w] angulo
xfimag=(1/n)*np.imag(xf)
# x[w] magnitud
xfabs=(1/n)*np.abs(xf)
# x[w] angulo en radianes
xfangle=np.angle(xf)# np.unwrap(np.angle(xf))


#SALIDA
plt.figure(1)   # define la grafica
plt.suptitle('Señal de Entrada')
plt.ylabel('xanalog[t]')
plt.xlabel('tiempo')
plt.plot(t,xanalog)
plt.margins(0,0.05)
plt.grid()

plt.figure(2)   # define la grafica
plt.suptitle('Transformada Rápida Fourier FFT')
ventana=ancho*0.01 # ventana de frecuencia alrededor f=0
ra=int((n/2)*(1-ventana)) #rango de indices
rb=int((n/2)*(1+ventana))+1
a=frq[ra] #rango de ejes uniformes entre graficos
b=frq[rb]
c=np.min([-0.05,1.1*np.min(xfreal)])
d=1.1*np.max(xfabs)
ejes=[a,b,c,d]

plt.subplot(221)    # grafica de 2x2, subgrafica 1
plt.ylabel('x[f] real')
plt.xlabel(' frecuencia (Hz)')
plt.axis(ejes)
plt.stem(frq[ra:rb],xfreal[ra:rb])
plt.grid()

plt.subplot(223)    # grafica de 2x2, subgrafica 2
plt.ylabel('x[f] imag')
plt.xlabel(' frecuencia (Hz)')
c=np.min([-0.05,1.1*np.min(xfimag)])
ejes=[a,b,c,d]
plt.axis(ejes)
plt.stem(frq[ra:rb],xfimag[ra:rb])
plt.grid()

plt.subplot(222)    # grafica de 2x2, subgrafica 3
plt.ylabel('x[f] magnitud')
plt.xlabel(' frecuencia (Hz)')
c=np.min([-0.05,1.1*np.min(xfabs)])
ejes=[a,b,c,d]
plt.axis(ejes)
plt.stem(frq[ra:rb],xfabs[ra:rb])

plt.grid()

plt.subplot(224)    # grafica de 2x2, subgrafica 4
plt.ylabel('x[f] fase radianes')
plt.xlabel(' frecuencia (Hz)')
c=1.1*np.min(xfangle)
d=1.1*np.max(xfangle)
ejes=[a,b,c,d]
plt.axis(ejes)
plt.stem(frq[ra:rb],xfangle[ra:rb])
plt.grid()

plt.show()

En el caso del seno(ω0 t), las graficas resultantes son:


con transformada:

continuar el ejercicio con las otras funciones.

Publicado por

Edison Del Rosario

edelros@espol.edu.ec / Profesor del FIEC/FCNM-ESPOL