import math
import numpy as np
from matplotlib.pyplot import *
from numpy.fft import fft 

def signal(t):
    return 0.1+np.cos(2*np.pi*t)+0.5*np.cos(2*2*np.pi*t) \
    +0.25*np.sin(3*2*np.pi*t)+0.01*np.sin(4*2*np.pi*t)
            

T=100.345821
fe = 10
t = np.arange(start=0.0,stop=T,step=1.0/fe)
echantillons = signal(t)
            

tfd = fft(echantillons)
N=len(echantillons)
spectre = np.absolute(tfd)*2/N
            

freq=np.arange(N)*1.0/T
            

figure(figsize=(10,4))
plot(freq,spectre,'r')
xlabel('f')
ylabel('A')
axis([-0.1,fe/2,0,spectre.max()])
grid()
            

spectre_db = 20*np.log10(spectre/spectre.max())
figure(figsize=(10,4))
plot(freq,spectre_db,'r')
xlabel('f')
ylabel('AdB')
axis([-0.1,fe/2,spectre_db.min(),spectre_db.max()])
grid()
            

zeros=np.zeros(N)
echantillons_z = np.concatenate((zeros,echantillons,zeros))
spectre_z = np.absolute(fft(echantillons_z))*2/N
N_z=len(echantillons_z)
freq_z=np.arange(N_z)*fe/N_z
figure(figsize=(10,4))
plot(freq_z,spectre_z,'r')
xlabel('f')
ylabel('A')
axis([-0.1,fe/2,0,spectre_z.max()])
grid()
             

axis([0.8,1.2,0,spectre_z.max()])
             

def fenetre(t):
    if t < T:
        return 1.0
    else:
        return 0.0
Tf = T*10
t = np.arange(start=0.0,stop=Tf,step=1.0/fe)
Nf = t.size
echantillons_fenetre = np.zeros(Nf)
for k in range(Nf):
    echantillons_fenetre[k] = fenetre(t[k])
spectre_fenetre = np.absolute(fft(echantillons_fenetre))
spectre_fenetre_db = 20*np.log10(spectre_fenetre/spectre_fenetre.max())
Nf = spectre_fenetre_db.size
freq_fenetre = np.zeros(Nf)
for k in range(Nf):
    freq_fenetre[k] = 1.0/Tf*k
figure(figsize=(8,4))
plot(freq_fenetre,spectre_fenetre_db,'r')
xlabel('f')
ylabel('AdB')
axis([0,0.2,-80,spectre_fenetre_db.max()])
grid()
        

from scipy import signal 
hamming = signal.hamming(N)
figure(figsize=(8,4))
plot(hamming)
xlabel('k')
ylabel('A')
ylim(0,1)
grid()
        

spectre_fenetre = np.absolute(fft(hamming,Nf))
spectre_fenetre_db = 20*np.log10(spectre_fenetre/spectre_fenetre.max())
figure(figsize=(10,4))
plot(freq_fenetre,spectre_fenetre_db,'r')
xlabel('f')
ylabel('AdB')
axis([0,0.2,-100,spectre_fenetre_db.max()])
grid()
        

echantillons = echantillons*hamming
echantillons_z = np.concatenate((zeros,echantillons,zeros))
spectre_z = np.absolute(fft(echantillons_z))*2/N/0.54
N_z=len(echantillons_z)
freq_z=np.arange(N_z)*fe/N_z
figure(figsize=(10,4))
plot(freq_z,spectre_z,'r')
xlabel('f')
ylabel('A')
axis([-0.1,fe/2,0,1])
grid()
        

def calculerSpectre(echantillons,fe,fenetre,nz=2,db=False):
    N = echantillons.size
    p = int(np.log(nz*N)/np.log(2))
    NN = 2**p
    zeros=np.zeros(NN-N)
    echantillons = np.concatenate((echantillons*signal.get_window(fenetre,N),zeros))
    spectre = np.absolute(fft(echantillons))
    spectre = spectre / spectre.max()
    if db:
        spectre = 20*np.log10(spectre)
    freq = np.arange(NN)*fe/NN
    return [freq,spectre]
            

[freq,spectre] = calculerSpectre(echantillons,fe,"hamming",nz=2,db=True)
figure(figsize=(10,4))
plot(freq,spectre,'r')
xlabel('f')
ylabel('AdB')
axis([-0.1,fe/2,-160,0])
grid()
            

[freq,spectre] = calculerSpectre(echantillons,fe,"hann",nz=2,db=True)
figure(figsize=(10,4))
plot(freq,spectre,'r')
xlabel('f')
ylabel('AdB')
axis([-0.1,fe/2,-160,0])
grid()