It depends on the accuracy needed and how important it is to you. You can use short-time Fourier transforms (STFT) to get rough estimates of the frequency and decay time, or you can put effort into getting high-quality estimates. I'll outline how I'd do the latter.
There are 2N+2 parameters that need to be estimated, where N is the number of excitation impulses:
- Oscillation frequency f
- Decay time τ
- Number of excitation impulses N
- Time of each impulse, tn
- Amplitude of each impulse, An
Any real system will have noise present, so the accuracy you can achieve will be determined by the SNR of the smallest impulse. Good estimates require SNR much greater than one. I would use a Bayesian estimator with these steps.
- Form the autocorrelation function. There will be N2 correlation peaks, assuming that the intervals between pulses are unique. This gives you N.
- I would produce K-1 frequency estimates from the phase differences between successive pairs of samples, where K is the number of samples in your data record. The resulting distribution of values can be used as the Bayesian prior distribution for f.
- The Fourier transform (FT) of an exponentially damped sinusoid is a Lorentzian line centered at f. (A Lorentzian line has a real absorption peak and an imaginary dispersive shape that is bimodal.) The FT of a sum of decaying exponentials with random time offsets is a sum of Lorentzian lines at f but with random phase rolls. The random phases will cause the sum to tend towards zero, but if N is not too large, you should be able to see the basic line shape. Estimate τ from the width of the real absorption peak. For the prior probability distribution, use a Gaussian centered on τ with a standard deviation (width) that represents the state of uncertainty of this fairly crude initial guess.
- Use Bretthorst's Bayesian methodology to refine f and τ, https://bayes.wustl.edu/glb/book.pdf , using his multiple sinusoid model. The pulse times and amplitudes are "nuisance parameters" that can be integrated away so they do not need to be estimated.
- Repeat the entire analysis for N-1 and N+1 pulses (and maybe more) to be sure you've gotten it right. This is particularly important if SNR is low or if the time between pulses is short compared to the decay time.
If you know the noise power level, you can compute the Cramer Rao lower bound on the variance of the frequency estimate to get an idea of the best that can be possibly accomplished, and then use that to determine how many frequencies you should test for. Let's say you decide to test 100 frequencies and also try 100 decay time values, then you will compute 10,000 posterior probabilities for each N value. The parameter set [f, τ, N] with the highest probability (the maximum
a posteriori probability or MAP) is your best estimate. A nice thing about Bayesian estimation is it provides information about the confidence in the estimate. If the MAP is 95%, you can be confident that you have fit the data well, but if it's 22% (due to low SNR, e.g.), you know not to trust it too much.
If that's too much, just play around with the FT and use the estimates in #2 and #3 as your estimates instead of as prior probability information for the Bayesian inference.