Le théorème de Parseval ne tient pas pour FFT d'une sinusoïde + bruit?

Question

Merci d'avance pour toute aide sur ce sujet. J'ai récemment essayé d'élaborer le théorème de Parseval pour les transformées de Fourier discrètes quand le bruit est inclus. J'ai basé mon code de this code. Ce que je m'attendais à voir est que (comme quand aucun bruit n'est inclus) la puissance totale dans le domaine fréquentiel est la moitié de la puissance totale dans le domaine temporel, comme j'ai coupé les fréquences négatives. Toutefois, plus le bruit est ajouté au signal dans le domaine temporel, plus la puissance totale de la transformée de Fourier du signal + bruit devient très inférieure à la moitié de la puissance totale du signal + bruit.

Mon code est le suivant:

import numpy as np 
import numpy.fft as nf 
import matplotlib.pyplot as plt 

def findingdifference(randomvalues): 

    n    = int(1e7) #number of points 
    tmax   = 40e-3 #measurement time 
    f1    = 30e6 #beat frequency 

    t    = np.linspace(-tmax,tmax,num=n) #define time axis 
    dt    = t[1]-t[0] #time spacing 

    gt    = np.sin(2*np.pi*f1*t)+randomvalues #make a sin + noise 

    fftfreq   = nf.fftfreq(n,dt) #defining frequency (x) axis 
    hkk    = nf.fft(gt) # fourier transform of sinusoid + noise 
    hkn    = nf.fft(randomvalues) #fourier transform of just noise 

    fftfreq   = fftfreq[fftfreq>0] #only taking positive frequencies 
    hkk    = hkk[fftfreq>0] 
    hkn    = hkn[fftfreq>0] 

    timedomain_p = sum(abs(gt)**2.0)*dt #parseval's theorem for time 
    freqdomain_p = sum(abs(hkk)**2.0)*dt/n # parseval's therom for frequency 
    difference  = (timedomain_p-freqdomain_p)/timedomain_p*100 #percentage diff 

    tdomain_pn = sum(abs(randomvalues)**2.0)*dt #parseval's for time 
    fdomain_pn = sum(abs(hkn)**2.0)*dt/n # parseval's for frequency 
    difference_n = (tdomain_pn-fdomain_pn)/tdomain_pn*100 #percent diff 

    return difference,difference_n 

def definingvalues(max_amp,length): 

    noise_amplitude  = np.linspace(0,max_amp,length) #defining noise amplitude 
    difference   = np.zeros((2,len(noise_amplitude))) 
    randomvals   = np.random.random(int(1e7)) #defining noise 

    for i in range(len(noise_amplitude)): 
     difference[:,i] = (findingdifference(noise_amplitude[i]*randomvals)) 

    return noise_amplitude,difference 

def figure(max_amp,length): 

    noise_amplitude,difference = definingvalues(max_amp,length) 

    plt.figure() 
    plt.plot(noise_amplitude,difference[0,:],color='red') 
    plt.plot(noise_amplitude,difference[1,:],color='blue') 
    plt.xlabel('Noise_Variable') 
    plt.ylabel(r'Difference in $\%$') 
    plt.show() 

    return 
figure(max_amp=3,length=21) 

Mon dernier graphique ressemble à ce figure. Est-ce que je fais quelque chose de mal en travaillant cela? Y at-il une raison physique que cette tendance se produit avec un bruit supplémentaire? S'agit-il de faire une transformation de Fourier sur un signal pas parfaitement sinusoïdal? La raison pour laquelle je fais cela est de comprendre un signal sinusoïdal très bruyant pour lequel j'ai des données réelles.

Answer 1

Le théorème de Parseval est généralement utilisé si vous utilisez le spectre entier (fréquences positives et) pour calculer la puissance.

La raison de l'écart est le composant DC (f = 0), qui est traité quelque peu spécial.

D'abord, d'où vient le composant DC? Vous utilisez np.random.random pour générer des valeurs aléatoires entre 0 et 1. Donc, en moyenne, vous augmentez le signal de 0.5 * noise_amplitude, ce qui implique beaucoup de puissance. Ce pouvoir est correctement calculé dans le domaine temporel. Cependant, dans le domaine fréquentiel, il n'y a qu'un seul bac FFT qui correspond à f = 0. La puissance de toutes les autres fréquences est répartie sur deux cases, seule la puissance DC est contenue dans une seule case.

En augmentant le bruit, vous ajoutez une alimentation en courant continu. En supprimant les fréquences négatives, vous supprimez la moitié de la puissance du signal, mais la majeure partie de la puissance du bruit se trouve dans le composant CC qui est utilisé complètement.

Vous avez plusieurs options:

1. Utilisez toutes les fréquences pour calculer la puissance.
2. Utiliser le bruit sans composante continue: randomvals = np.random.random(int(1e7)) - 0.5
3. « Fix » le calcul de la puissance en supprimant la moitié de la puissance DC: hkk[fftfreq==0] /= np.sqrt(2)

je partirais avec l'option 1. La seconde pourrait être OK et je ne recommande pas vraiment 3.

Enfin, il y a un problème mineur avec le code:

fftfreq   = fftfreq[fftfreq>0] #only taking positive frequencies 
hkk    = hkk[fftfreq>0] 
hkn    = hkn[fftfreq>0] 

Cela ne fait pas vraiment de sens. Mieux le changer à

hkk    = hkk[fftfreq>=0] 
hkn    = hkn[fftfreq>=0] 

ou l'enlever complètement pour l'option 1.