Ayub Sherif
Ayub Sherif

Reputation: 1

Power Spectrum Analysis for a very large set of data

I have a voltage signal that I am trying to denoise. The signal comes in very large files (524288 cells). When I take the whole file and make n equal to length of data set, I get two extremely large peaks at frequency 0 and max.

sig = np.genfromtxt(directory + '/'+ file, skip_header=3, dtype=np.float64)
n = len(sig)
freq = np.arange(n)
fhat = np.fft.fft(sig, n)
PSD = fhat * np.conj(fhat) / n

plt.plot(freq,PSD)              
plt.show()

indices: bool = PSD > 100
fhat = indices * fhat
ffilt = np.fft.ifft(fhat)
plt.plot(ffilt)
plt.show()

Is there a way to analyse the whole file or I have to split it to smaller data sets?

Upvotes: 0

Views: 665

Answers (1)

Pedro
Pedro

Reputation: 330

You can analyze the whole file in one go. Those high peaks might a padding issue. Anyway, I would suggest you to use one of the already implemented methods in Python like the Welch method from Scipy. It would save you time and less headaches figuring out if your implementation is correct. Here is an example adapted from the Scipy lectures:

import numpy as np
from matplotlib import pyplot as plt    
from scipy import signal

# Seed the random number generator
np.random.seed(0)

time_step = .01
time_vec = np.arange(0, 70, time_step)

# A signal with a small frequency chirp
sig = np.sin(0.5 * np.pi * time_vec * (1 + .1 * time_vec))

plt.figure()
plt.plot(time_vec, sig)

# Compute the Power Spectral Density
freqs, psd = signal.welch(sig)

plt.figure()
plt.semilogx(freqs, psd)
plt.title('PSD: power spectral density')
plt.xlabel('Frequency')
plt.ylabel('Power')
plt.tight_layout()
plt.show()

And here are the results: enter image description here

Upvotes: 0

Related Questions