# Estimating the Most Likely Harmonic Signal in a Spectrogram

Experimental data description: below image is from spectrogram of doppler radar return as I walk toward sensor. Bright sinusoid shape modulated with frequency F is from motion of small retroreflector swinging in my right hand. Lower level ripple at freq. 2*F from my overall walking motion. Other noisy signals at F and 2F with other phase shifts may be from multipath reflections. Uncorrelated noise may be from trees blowing in wind, RFI, thermal receiver noise, and other sources. This is just one example, but I am interested in the general problem of tracking motion in noisy conditions. [edited for clarity]

When I look at this spectrogram image by eye, I immediately see a sinewave-like signal amid noise, containing about 12 cycles, and I also notice it is changing slightly in amplitude, frequency, and harmonic content over time. I think I could trace by hand a single line overlaying the sine-like signal I see, and with a lot of care I believe I could get a repeatability near one pixel doing this, given several images containing the same signal but different background noise of comparable amplitude to that shown here.

Is there a signal processing algorithm that could do this in a robust way? I imagine I could easily find a best-fit to one fixed frequency, single-amplitude sine, but I want to do better. If I do the simplest thing I could think of; independently select the brightest pixel from each vertical column, I mostly get the peaks of the sine-like signal, but randomly capture other noise toward the middle of each cycle. If I first do a 2D smoothing, it's not much better, and with a 1D vertical smoothing, still not nearly as good as my intuition suggests should be possible. Is there some defined algorithm analogous to what my visual system does to find and trace out this signal?

My limited understanding of Kalman filters is that you need to start with a model for your signal. However I saw the sinewave in this image without necessarily expecting one, and I'm wondering if there's a way to formalize whatever process that is.

Simply selecting brightest pixels at each timestep gives this. I would say the sine signal is only sort of visible here.

Here is another version of the first image showing only the brightest parts. I think I can still visually sketch out the same sine-like signal from this image, starting at the right-hand edge and then working backwards into fainter areas.

Ideally, I am imagining some algorithm that could produce the line shown in white (which here was traced by hand),

and could also work with other continuous-trace signals with similar signal-to-noise conditions, which are not so periodic. For example, this trace

where the faint segment at 8..12 seconds would be more dubious without the connecting and more visible portions before and after.

To try this yourself, here is the signal as mp3 file with audio used to create the sinusoidal-signal spectrogram image. To generate a plot containing the image segment, convert the .mp3 file to .wav format, and then run this plotIQ script, as shown with these two commands on a Linux system:

sox 2023-01-03_Doppler_Segment.mp3 2023-01-03_Doppler_Segment.wav
python plotIQ.py 2023-01-03_Doppler_Segment.wav

• My guess is visually, I start with the brighter parts of the image that form recognizable curves (the tops and bottoms of the sine wave). Then, I draw a line connecting these together in a smooth way like spline fitting, but by this point I've recognized a sine shape, and so I do a piecewise fitting to form one connected curve like a sine, but smoothly deformed to pass through image local maxima that are "near enough" to the expected path of the sine curve. Is something like this already coded somewhere? Commented Jan 3, 2023 at 21:56
• Kalman filters and many other approaches require you to think you know something about how the signal was generated. Do you know what you're expecting to see? Some background information here (not about the Kalman filter) might help you think about the problem.
– Peter K.
Commented Jan 4, 2023 at 0:10
• what can you do to AVOID all that noise? what causes this noise? edit your question to include the context you give in the comments below. Commented Jan 5, 2023 at 8:48
• @user3732008, Any chance you share the data and the script to create the image so we can try work on it?
– Royi
Commented Jan 11, 2023 at 6:55
• Edited post to include link to data and script. Commented Jan 11, 2023 at 22:08

Given this description:

When I look at this spectrogram [...] I immediately see a sinewave-like signal amid the noise, containing about 12 cycles, and I also notice it is changing slightly in amplitude, frequency, and harmonic content over time.

what you are trying to do is track a frequency modulated signal.

A typical frequency modulation process looks like this:

$$\begin{matrix} x(t) = u_{mod} \times sin(2 \pi f_{mod} t) \\ y(t) = u_{car} \times cos(2 \pi (f_{car} + x(t)) t) \end{matrix}$$

Where:

• $$u_{mod}, u_{car}$$ is the amplitude of $$x(t), y(t)$$ respectively, but also notice here that the amplitude of $$x(t)$$ determines "how far" off $$f_{car}$$ will be offset in $$y(t)$$.
• $$f_{mod}$$ is the modulation frequency (and $$x(t)$$ is the modulation signal)
• $$f_{car}$$ is the carrier frequency (and $$y(t)$$ is the carrier signal)

"Carrier" and "modulation" usually refer to a communications context where we are trying to get a baseband signal (e.g. audio, typically between 20Hz - 20kHz) to modulate a carrier signal (e.g. RF at the range of MHz) in order to better couple it over a medium (e.g. RF, water, optical fibres and other media).

The correspondence with your descrption here is this:

...I immediately see a sinewave-like signal amid the noise...

That is your carrier signal. A high frequency signal that creates the "...brightest pixel from each vertical column..." in your time-frequency (spectrogram) representation.

That is your modulation signal. Its $$f_{mod} \approx 12$$ Hz ... and;

... I also notice it is changing slightly in amplitude, frequency, and harmonic content over time

...this part refers to the dynamics of the characteristics of $$x(t)$$, or, in other words turning amplitude and frequency to functions of time ($$u_{mod}(t)$$, $$f_{mod}(t)$$).

The easiest way to track a frequency modulated signal is using a Phase Locked Loop.

A PLL is a closed loop automatic control system that can track phase (frequency) changes of a given signal at its input.

The "error signal" of the PLL (at the output of the low pass filter) will return the $$f_{car} + x(t)$$ component which is the dominant "shape" that you see in the spectrogram.

The PLL can track a single frequency within its design specification. Looking at your spectrogram, it seems that your modulation produces harmonics, which might be a product of the dynamics of $$u_{mod}(t)$$ or creeping in to your signal from other sources.

Personally, what I see here is a dominant sinusoid but at the same time there is a lower amplitude, approximately same frequency sinusoid that is like the first derivative of the dominant one. Its peaks coincide with the slopes of the dominant sunsoid (and with a negative coefficient, when the first derivative hits positive maximum the dominant sinusoid is decreasing its frequency). Towards the end it is clearer that there are instances where you have two relatively bright pixels in the same column. If this was audible it would sound like an ambulance siren with an additional faint component at the back...sometimes you get this kind of distortion from very cheap, small buzzers trying to reproduce anything other than a single beep.

The PLL (on its own) is not going to help you with explaining all of this variation and some knowledge of what this signal is and some domain knowledge about it would definitely help, as it has already been pointed out in the comments, to home in to additional methods of tracking it.

Hope this helps.

EDIT:

What is hinted at in the additions to the original question and comments to this answer is a way of optimising a model from the spectrogram by viewing the spectrogram as a grayscale image that contains a form to be approximated.

This method of dropping "...a springy, flexible metal strip on top, shake it around and find the average position. I imagine the lowest-energy shape of that strip traces out the curve I'm looking for." is called an Active Contour Model.

An active contour model deforms a spline iteratively until an equilibrium is reached. The equilibrium is between forces that represent the dynamics of the spline itself (e.g. stiffness) and forces that deform the spline by pulling it towards features of interest (using some form of gradient).

Although the word "model" is in the name of this technique, it does not directly relate the form of the spline to characteristics that were described in the original question (e.g. Depth of modulation, modulation frequency).

Active contours would present a number of challenges for this application:

1. Selection of appropriate parameters for the spline itself (too soft and it will be attracted easily by small perturbations, too stiff and it will fit in a more...general sense.)
2. The initial state affects the final solution. For the same image, different initial states of the spline can lead to it being attracted to different image features.
3. Convergence can be lengthy

As far as the "look-ahead" is concerned, it will have to be taken into account as a set of additional constraints over the way the active contour is optimised. Otherwise, the original model does not imply some kind of "look-ahead". In other words, the fact that one part of the image would seemingly fit a sinusoid-like form does not mean that the model would "figure out on its own" that it has to extend that form through areas of less good fitness in a particular way.

Finally, a spline over a spectrogram implies a single frequency being associated with any given point along the "trajectory" of the spline. Therefore, it does not add something to the simpler PLL method which can also track a single frequency.

... is the PLL optimal in some sense?

I do not think I can answer this question because we do not have enough details about the dynamics of the possible targets that might show up within the field of view of the radar. But, at the very least, the PLL matches the operation / output of the Doppler radar module.

• Thanks for comments. Data is a doppler radar return as I walk toward sensor. Bright sinusoid at freq. F is from motion of small retroreflector swinging in my right hand. Lower level ripple at freq. 2*F from overall walking motion. Yes PLL could work, but could post-processing improve a fit by considering full dataset at once (as you do visually looking at the plot) instead of a causal process like PLL that uses only prior time data? Was hoping for some method to detect different signals too, not necessarily periodic ones, similar to how you can see them in a spectrogram by eye. Commented Jan 4, 2023 at 16:14
• @user3732008 This is not a comment :) It is an answer to the question you are stating above. The PLL will detect any kind of frequency modulation, not necessarily periodic. In other words $u(t)$ would have had any sort of waveform and the PLL would still track it. There are ways to fit models to "image" data but they require you to have an idea of the model beforehand (e.g. a formula). So in this case, we would fit the model of a modulated sinusoid to the image data.
– A_A
Commented Jan 5, 2023 at 11:10
• If my idea is only a bandwidth-limited, single-valued function of time, is the PLL optimal in some sense? Say I make a 2D array of magnets, 1 magnet per pixel in the spectrogram, with magnetic strength in proportion to each pixel's brightness. Then drop a springy, flexible metal strip on top, shake it around and find the average position. I imagine the lowest-energy shape of that strip traces out the curve I'm looking for. In general, it may be different from a PLL output that cannot look ahead in time. Might that form a better fit to data, without using a more specific model? Commented Jan 5, 2023 at 15:55
• @user3732008 Please see edit. It would be useful to talk a little bit more about what the objective is. For example, are you building a classifier to figure out what sort of activity takes place in front of the radar?
– A_A
Commented Jan 5, 2023 at 23:34
• Yes classifier would be interesting, but I was trying for a best-effort quantitative speed vs. time curve, per detected target. I added another image to the bottom of my original question with another (more typical) signal, that may benefit from the adaptive contour idea. Commented Jan 6, 2023 at 2:38