Intégration de spectre et écart-type
Bonjour,
Je suis doctorant en aérodynamique et je cherche à exploiter des données obtenus en soufflerie. J'essaie d'étudier le taux de turbulence que l'on trouve étonnant. Je possède un fichier avec la tension d'une sonde (250 000 échantillons à 25 kHz, avec un filtre passe-bas pour éviter l'aliasing à 11.3 kHz), cette tension peut être convertie en vitesse avec une simple loi, j'ai donc de manière équivalente 250 000 échantillons de vitesse.
Je cherche à étudier l'influence des très basses fréquences sur la turbulence que l'on définit ainsi: Tu = Ecart-Type(U) / U_moyen, U étant la vitesse. Pour faire ça j'essaie de faire la chose suivante:
-Effectuer la PSD de la partie fluctuante de la vitesse (pour les 250 000 échantillons, on appelle cette composante u[i] = U[i] - U_moyen)
-Intégrer cette PSD sur toute la gamme de fréquence sauf la plus basse (environ 6Hz seront coupés), pour retomber sur l'écart-type(U) sans les basses fréquences (sauf erreur, intégrer la PSD de u[i] = U[i] - U_moyen doit me faire retomber sur la valeur de l'écart-type)
Dans une premier temps, j'essaie d'intégrer l'écart-type en prenant toute la gamme de fréquence pour vérifier qu'on obtient bien la même valeur que par la commande np.std(U), mais ceci ne marche pas. Sauriez-vous d'où vient l'erreur? Je n'arrive pas à savoir si elle vient d'un souci de raisonnement de ma part ou simplement d'implémentation.
Je mets ci-dessous ma boucle principale simplifiée (j'ai 36 points avec 250 000 échantillons de vitesse)
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
|
for i in range (nbdepoints):
#Initialisation des variables temporaires
tension_FC=[]
vitesse_FC=[]
#Récupération de la position et des valeurs de tension
position=np.genfromtxt(filchaud,dtype=str,usecols=1, unpack=True, skip_header=3, max_rows=1)
tension_FC =np.genfromtxt(filchaud,dtype=str,usecols=(1), unpack=True, skip_header=2,max_rows=nb_echantillons)
position=str(position)
#Conversion virgules en point et tension en vitesses
position=convertisseur1(position)
tension_FC=convertisseur(tension_FC)
vitesse_FC=convertisseurvitesse(tension_FC)
for j in range(nb_echantillons):
U_inst[i][j]=vitesse_FC[j]
#alimentation des données à garder
U[i]=np.mean(vitesse_FC)
u_rms[i]=np.std(vitesse_FC)
position_longi[i]=position+position_depart
Tu[i]=u_rms[i]/U[i]
if i==36: ######## c'est une version simplifiée où je n'étude que le dernier point pour aller plus vite
f1,PSD1=signal.welch(vitesse_FC,frequence_acqui, nperseg=4096)
#On calcule u' adimensionné par U ######## vitesse fluctuante appelée u ci-dessus
u_prime_ad=[]
for j in range(len(vitesse_FC)):
u_prime_ad.append((vitesse_FC[j] - U[i]))
f2,PSD2=signal.welch(u_prime_ad,frequence_acqui, nperseg=4096)
#Intégraton par la méthode des trapèzes
a=np.sqrt(np.trapz(PSD2,f2)) #Sur toutes les fréquences ##### Logiquement je pensais obtenir ici la même valeur que par np.stp(Vitesse)
a1=np.sqrt(np.trapz(PSD2[1:],f2[1:])) #On enlève la première bande de fréquences
a2=np.sqrt(np.trapz(PSD2[2:],f2[2:])) #On enlève les deux premières |
Merci d'avance