0
$\begingroup$

I am trying to calculate the instantaneous power of a 1-d signal within a particular frequency band in an on-line, real-time, streaming application. I have tried using sliding window FFTs but these require computing an average power over a window rather than the instantaneous power at each incoming sample. Wavelet transforms also appear to not be suitable for on-line or real-time instantaneous power because of edge effects in the analysis window. I have considered using Hilbert-Huang Transforms, but these seem to have the same issue as wavelet transforms with edge effects and compute time. RMS analysis is fastest, but still returns average power over a window of samples rather than the power at a particular sample.

Problem statement: Given a data stream data_stream = x_0, x_1, x_2, x_3..., x_n where x_i is a single, real-valued number, how to compute the the instantaneous power of the signal at x_i for a particular frequency band?

Attempts at FFT and RMS solutions (pseudocode):

  1. Sliding window FFT:
window_length = 50 # number of samples required to capture 2 cycles of lowest frequency in frequency band of interest
low_f = 0.25 # 0.25 Hz lowest frequency of interest
high_f = 0.45 # 0.45 Hz highest frequency of interest
band = [low_f, high_f] # band of interest includes all components of signal from 0.25-0.45 Hz
start_index = 0
end_index = start_index + window_length
while end_index < len(data_stream):
    window = data_stream[start_index: end_index]
    fft = fft(window) # compute fft on window
    power = fft[band]
    start_index += 1
    end_index += 1

Here, I keep the window length constant at 50 samples, and move the window one sample forward each time I receive a new sample. The band power is the average power of the band from all 50 samples in the window data_stream[start_index: end_index] rather than the instantaneous power of the signal at data_stream[end_index - 1].

  1. Sliding window RMS:
filter = bandpass_filter(low_f, high_f)
window_length = 50 # number of samples required to capture 2 cycles of lowest frequency in frequency band of interest
start_index = 0
end_index = start_index + window_length
while end_index < len(data_stream):
    window = data_stream[start_index: end_index]
    filtered_window = filter(window) # reduce signal to only frequencies of interest
    power = RMS(window)
    start_index += 1
    end_index += 1

Here, I bandpass filter each window of the data_stream to keep only the frequency band of interest. I keep the filtered window length constant at 50 samples, and move the filtered window one sample forward each time I receive a new sample. The band power is the average power of the band from all 50 samples in the filtered window filtered(data_stream[start_index: end_index]) rather than the instantaneous power of the signal at data_stream[end_index - 1]. This sliding window RMS approach should give me approximately the same results as the sliding window FFT as FFT power is ~ RMS.

  1. "Instantaneous Power":
online_filter = online_bandpass_filter(low_f, high_f) # see https://stackoverflow.com/questions/40483518/how-to-real-time-filter-with-scipy-and-lfilter
current_index = 0
while current_index < len(data_stream):
    current_sample = data_stream(current_index) # get new sample
    filtered_sample = online_filter(current_sample)
    power = filtered_sample ** 2
    current_index += 1

Here I use an on-line bandpass filter (as discussed in this post)on each sample as it is received. Assuming that filtered_sample represents only the frequency band of interest, I then square the value of that sample to get the "instantaneous power".

As far as I understand, the way real-time oscilloscopes / spectroscopes work is computing the FFTs on a sliding window of incoming samples, and then presenting the FFT PSD for each window as the "real-time" PSD.

For example, if a window of 2.0s is required to capture all frequencies of interest, and the sampling rate is 2 Hz, then the oscilloscope /spectroscope computes the FFT PSD on the samples from [t = -2.0s, t = 0.0s] and presents that PSD as the "instant, real-time" PSD at time t = 0.0s. Another sample is obtained, and the window shifts to [t = -1.5s, t = 0.5s]. FFT PSD is computed on the new window, and that PSD is reported as the "instant, real-time" PSD at time t = 0.5s.

In reality, the PSD computed from [t = -2.0s, t = 0.0s] is the PSD over the interval rather than the PSD at t = 0.0s, and likewise for the PSD computed at the next sample acquisition from [t = -1.5s, t = 0.5s].

Am I understanding this correctly? Is there a better way to calculate the power of a signal at the latest available data point while new data arrives in a on-line / streaming manner? My attempt #3 ("Instantaneous Power") seems to be the closest, but I'm not sure if it's theoretically sound. The output is super high variance, noisy, and when averaged, not equivalent to the FFT output.

Thanks for your thoughts!

EDIT: Edited wrt to comments and previous answer.

$\endgroup$
2
  • $\begingroup$ You need to clarify your requirements. Your problem statement doesn't include the word frequency and the (only) answer as written would be simply $x_i^2$. 'Instantaneous Power of a Frequency Band' is non-sensical if you interpret Instantaneous to mean "only dependent on the current time domain sample". $\endgroup$
    – Hilmar
    Commented Mar 30, 2022 at 11:37
  • $\begingroup$ Thanks @Hilmar. Clarified above. $\endgroup$ Commented Mar 30, 2022 at 18:04

1 Answer 1

1
$\begingroup$

If you want to calculate the power in real time and only for the current sample, you can do so but only for whole frequency range, not for single bands. Also, such a measure would not yield any useful information, as it would just show a squared version of the waveform itself if you plotted the values in line.

Any method of calculating power in a certain frequency band needs a certain amount of samples to function correctly. This is kind of intuitive, since you cannot assess the spectral contents of a signal by just looking at one sample, you will always need a history, which in this case means some kind of sliding window, or, from a different viewpoint, some kind of buffer. Put in other words: frequency is the reciprocal of time. If the time range you consider is virtually zero (one sample), there cannot be any frequency information in this sample by itself.

The PSD is a well known and established method for looking at the spectral contents of a signal. The term "instantaneous" was coined for a very short buffer/sliding window and is misleading in this regard. "short time power" would probably be more apt.

As a sidenote: this is true also for analog signals and filterbanks. While there is no buffer/sliding window per se, the capacitors in those circuits act as kind of an analog memory, when the observed time range is in the same order of magnitude as the time constants of the circuit.

$\endgroup$
3
  • $\begingroup$ Thanks @Max. Clarified above. I am using an on-line bandpass filter so theoretically only the frequency band I care about is left in the signal after filtering? In which case, then I care only about power, and not frequency information so time range being zero should not be an issue? $\endgroup$ Commented Mar 30, 2022 at 18:06
  • $\begingroup$ Yes indeed. But this filter also has a history of samples and is not instantaneous in this sense. That the output is high variance and noisy is possible. Try feeding a pure sine signal through it, then you should see a constant value. With, say, music signals, high variance is usual. $\endgroup$
    – Max
    Commented Mar 30, 2022 at 18:52
  • $\begingroup$ Got it. Squaring each filtered sample as it comes through sounds like the best bet. $\endgroup$ Commented Mar 30, 2022 at 20:07

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.