1 pièce(s) jointe(s)
[Hann Window se déplaçant]
Bonjour,
Je souhaite réaliser le spectre d'un jeu de données afin de déterminer les caractéristiques du signal (Tp, Amplitude caractéristique etc.)
J'ai déjà fait un premier code mais les résultats sont bofs car il faut que j'ajoute une fenêtre de Hann.
J'ai donc l'algorithme suivant pour réaliser l'analyse de fourier :
- Temps d'acquisition : 30 min
- Fréquence acquisition : 2.56 Hz
L'algorithme demande d'étudier Fourier sur 200 secondes et faire 17 demi recouvrements:
- Première boucle : analyse de fourier sur les 200 premières secondes (de t=0 à t=200)
- Deuxième boucle : analyse de Fourier de t=100 à t=300
...
- Dernière boucle : t=1600 à t=1800
Pour chaque analyse, il faut faire l'algorithme suivant : (cf pièce jointe)Pièce jointe 278020
J'ai réalisé le code suivant où je stock les données dans une liste mais après je n'arrive pas à faire le recouvrement.
Code:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
| import pandas as pd
import math
import cmath
import matplotlib.pyplot as plt
import matplotlib #only needed to determine Matplotlib version number
import numpy as np
file=r'C:\Users\NathanS\Dropbox\Weptos Nathan\DWR4_74026{displ}2015-03-05.csv'
df = pd.read_csv(file, sep='\t', header=None, usecols=(0, 2))
df.columns = ['Time', 'Elevation']
#df=pd.concat([df.loc[df.index>=4607+4608+4608],pd.DataFrame()],ignore_index=True) #To start at t=60 min (to see if there is a change)
t=0
T=200
fs=2.56
ts=1/fs
N=int(T/ts)
PSD_moy=[]
m=0
nstart=0
df_temp=pd.concat([df.loc[(df.index>=m*4608) & (df.index<(1+m)*4608)],pd.DataFrame()],ignore_index=True) #To take only 30 minutes
while nstart!=1700:
H=[]
somme = 0
for i in range(0,N):
w=[]
h=[]
for k in range(0,N):
w.append((1/2)*(1-math.cos(2*math.pi*(k)/(N-1)))) #Hann Window
h.append(df_temp.Elevation[k+t*T/(2*ts)])
w2=0
for n in range(0,len(w)):
w2=w2+w[n]**2
wnorm=np.sqrt(fs*w2)
for k in range(0,len(w)):
somme=somme+(w[k]/wnorm)*h[k]*cmath.exp(2*math.pi*k*i*1j/N) #Create each H
H.append(somme)
freq=[]
for i in range(int(N/2)):
freq.append(i*fs/N) #Creation of a vector with frequency from 0 to fs/2 with resoltuion of 0.005
PSD=[]
# Creation of Power Spectral Density
PSD.append(abs(H[0])**2)
for i in range(1,len(freq)-1):
PSD.append(abs(H[i])**2+abs(H[N-i])**2)
PSD.append(abs(H[len(freq)])**2)
PSD=np.asarray(PSD)
PSD_moy.append(PSD)
nstart=nstart+100
t=t+1 |
Si quelqu'un a une idée pour arriver à écrire cet algorithme je vous en serai très reconnaissant !
