I am working with a modular spectral processing STFT toolkit I have written in C. I often need to carry spectral peak detection and analisys (frequency estimation, peak displacement, etc) but I am into troubles when there are two or more peaks merging in a same apparent peak, i.e many overlapping peaks resulting in a single peak detectable as a single local maximum. I am aware that this problem has been addressed copiously in literature, in particular in chromatography or mass spectrometry. A commonly suggested approach is analysing the 2nd derivative looking for valleys, but this trick is not viable when the spectral resolution is insufficient, as when you have perhaps two or three peaks sharing few bins. Also I found many cases where a very low frequency sawtooth spectrum results in a smooth and large single peak, for which even the 2nd derivative of magnitude looks monotonous. So what to do ? I suspect that proper analisys of phases other than of magnitudes alone could give more clues of the actual underlieing peak structure but I am not a maths guru despite skilled enough. I would like to be able, once detected a whole peak "body", to decompose in its constituting overlapping peaks (if more than one) and possibly isolate them, or at least having a good estimate. Is there a way to do that, which is not too computationally intensive and which can be done in real time ? Thanks for any relevant answers!
so, adding to both Damien's a hotpaw's answer, i recommend using a Gaussian window and making your window as wide as practical. the Fourier Transform of a Gaussian function is a Gaussian function. so the wider your Gaussian window in the time domain, the narrower the Gaussian peaks will be in the frequency domain. you can use quadratic interpolation of the log of the peak magnitude to locate the true peak frequency (between FFT bins) more precisely.
the Gaussian function does not have sidelobe bumps. but, because a true Gaussian goes on forever, you actually have to truncate the Gaussian to a finite length, which does give it tiny sidelobes in the frequency domain.
if you want to estimate both frequency and frequency sweep rate, this paper of mine from 2001 discusses how to do that with a Gaussian window.
Update: I found a trick which has the potential to work. It has still to be perfectioned though. I am explaining it here in the hope I can turn out useful or interesting for somebody else in future.
-Start from bin 0 of spectrum thru last bin;
-For every local maximum (peak tip) found, carry estimation of frequency (in fractionary bin units), phase and magnitude (such estimations can easily be made by parabola fitting or more advanced sinc-fitting methods or ok even by gaussian fitting but this is another topic);
-Basing on said estimates, resynthesize the peak by a suitable sinc function and subtract iteratively the whole synthesized peak (even just its right side) from the spectrum (for debugging purposes, the synthesized peaks can even be added to a separate spectral buffer);
-This way all hidden peaks will eventually come out.
This method can't by any means be defined mathematically sound and perfect of course, but is promising. A demonstration of that is that converting to audio all stft spectral frames obtained in the debug buffer said above, the result has minimal artifacts and is 99% close to the original (by ear). I also make the algorithm marking all detected peaks for debug purposes and I would say it does its job (it has still to be perfectioned though as I said but it basically works). It is even capable of separating peaks with very close frequency like a sinewave with chorus.
Try a different windowing function - with windowing you can either have good spectral resolution or low spectral leakage.
So if you look at frequency spectrum of each window, it has a mainlobe and many sidelobes. You can window with functions that have a large main lobe and small side lobes, thus reducing spectral leakage but at the cost of lower frequency resolution. Or if you desire greater frequency resolution you look for a window with a thinner main lobe at the cost of more spectral leakage.
The reason for this is when you truncate a signal for fourier analysis you are convoluting your signal with a window of that width - applying different windows will therefore affect your resultant output.
Astronomers are talking about the Rayleigh distance. The (angular) distance between two point sources (stars) where they can be resolved with some robustness using diffraction limited optics: https://courses.lumenlearning.com/physics/chapter/27-6-limits-of-resolution-the-rayleigh-criterion/
This article by Julius Orion Smith discuses spectral interpolation of (well resolved) sinoidals: https://ccrma.stanford.edu/~jos/sasp/Quadratic_Interpolation_Spectral_Peaks.html
If you want more frequency resolution and time resolution (as in: you cannot trust the input to be stationary for any longer), then you are in a tough situation. Perhaps fitting a smoothly sweeping sparse set of sinoids where you restrict amplitude or phase relationship?