-1

Let's say we have a time-transient with a chirped signal, $$V(t) = \cos\left( 2 \pi f(t) t \right)$$ where $f(t)$ is the frequency of the signal and is time dependant. My question is how would one describe the final frequency, as produced by an FFT spectrum? My intuition is that it would be the time average such that $$f^{\rm{FFT}}_{\rm{final}} = \frac{1}{T}\int^{T}_{0} f(t) \ {\rm{d}}t$$ where $T$ is the length of the time-transient signal used in the FFT. Strictly speaking I think would use a Gabor transform for non-stationary signals, but seeing as the spectrum analyser I am working with uses an FFT algorithm.

user27119
  • 199
  • 2
  • 17
  • 2
    Your first equation is wrong. Frequency is the derivative of the phase. Your first equation will not yield the frequency you expect unless frequency is constant.

    $V(t) = cos(2\pi \theta (t)) $ AND

    $ f(t) = \frac{d}{dt} \theta (t) * \frac{1}{2\pi}$

     – Ben Apr 09 '21 at 16:25
  • You're asking for "the" frequency, but that doesn't exist: the resulting spectrum is continuous, not discrete. – Marcus Müller Apr 09 '21 at 16:31
  • @Ben I think you made an important point, can you expand a bit more on how one would correct my expression? – user27119 Apr 09 '21 at 16:35
  • 2
    Related: this answer – Matt L. Apr 09 '21 at 16:39
  • @MarcusMüller of course you are right, the frequency is by definition changing throughout the entire acquisition time. But the final spectrum will have a peak signature which occurs at a frequency which must be some aggregate or average of all frequencies over the acquisition time. – user27119 Apr 09 '21 at 16:44
  • @MattL. Extremely helpful. Thanks. – user27119 Apr 09 '21 at 16:48
  • @Q.P. the final spectrum will have a peak signature no, that is not true, for most choices of $f(t)$, and especially not for things that I'd call "chirp". – Marcus Müller Apr 09 '21 at 17:52
  • 1
    @Q.P. You'd need to have a reproducible example in e.g. Octave and include both the time- and frequency-domain graphs. You claim that something happens according to the intuition - unnecessarily. It's easy to check if it actually happens. See for yourself first! :) – Kuba hasn't forgotten Monica Apr 09 '21 at 18:20
  • @MarcusMüller Maybe there is a distinction between what people really define as a chirp. In my case I have a very small frequency chip, $f(t) = f_{0} + \delta f e^{-t/\tau}$. Where $\delta f = 0.2 \ {\rm{Hz}}$ and $\tau = 200 \ {\rm{s}}$. Now for $\cos(2 \pi f(t) t)$ the FFT will produce a singled peak spectrum. – user27119 Apr 09 '21 at 18:38
  • 1
    no it really won't; an exponential has a very wideband spectrum, so you mathematically will get something that really resembles more higher-order Bessel functions than a peak. You still claim things on intuition that really aren't the case. Of course, I can define anything to be a "peak", but if you do that, you don't get to call it a chirp anymore, because you've "zoomed" out mentally far enough that the nature of the thing not being a constant-frequency tone got lost. – Marcus Müller Apr 09 '21 at 18:50
  • 2
    Does this answer your question? Simulation of a Frequency ramp – Dan Boschen Apr 10 '21 at 11:35

1 Answers1

0

Depending on the chirp, if you seek a closed form solution, there may not be one. For empirical estimation, one can use non-stationary methods.

Take exponential with frequency ranging from 0.005 to .125 of sampling frequency. Applying synchrosqueezed CWT and plotting energy norm per row vs cumulative sum of rows' associated frequency:

enter image description here

The red line marks where energy contributions grow insignificant (and thus max freq is likely reached), which can be found with a flatline detection algorithm. It's heuristics and engineering.

Code below, can be improved.

import numpy as np
from ssqueezepy import ssq_cwt, TestSignals
from ssqueezepy.visuals import plot, imshow

def flatline_after_max(x, count=10, rel_dist=.1, max_frac=.99):
    # pick slightly less than max
    max_idx = np.where(x > max_frac*x.max())[0][0]
    ctr = 0
    for i in range(max_idx, len(x)):
        if abs(x[i] - x[i - 1]) / x[i] < rel_dist:
            ctr += 1
        else:
            ctr = 0
        if ctr == count:
            return i - (count - 1)


define test signal


N = 4096  # sampling rate (or number of samples)
fmin = .005
fmax = .125
x = TestSignals().echirp(N=4096, fmin=fminN, fmax=fmaxN)[0]


synchrosqueeze & visualize


Tx, , ssq_freqs, * = ssq_cwt(x, ('gmw', {'beta': 4}))
imshow(Tx, abs=1, yticks=ssq_freqs[::-1], ylabel="Fraction of fs",
       title="abs(SSQ_CWT)")


compute row energies and find 'stopping point'


Tx_energies = np.sum(np.abs(Tx[::-1])**2, axis=1)
Tx_energies_cs = np.cumsum(Tx_energies)
max_idx = flatline_after_max(Tx_energies_cs)
est_freq = ssq_freqs[max_idx]


plot


plot(ssq_freqs, Tx_energies_cs, xlabel="Fraction of fs",
     title="cumsum(energy(Tx, axis=1)) | est_freq=%.3f" % est_freq,
     vlines=(est_freq, {'color': 'tab:red'}))
OverLordGoldDragon
  • 8,912
  • 5
  • 23
  • 74