dayalannair / FMCW_FFT_Radar

Implementation of the Xilinx FFT IP using a UART packetised data stream
MIT License
7 stars 1 forks source link

Fequency estimation algorithms in industry #37

Open dayalannair opened 2 years ago

dayalannair commented 2 years ago

Hi @jpt13653903 I would like to know what algorithms are used in industry for FMCW beat frequency extraction. So far I believe CFAR is a common method, and Im wondering if the root MUSIC algorithm provides better results (seems like it can provide super resolution through some interpolation). Not too sure on the maths but I think it determines which frequency is closest to the unit circle after decomposing a signal into eigenvalues.

Are there any alternatives to the two of these?

jpt13653903 commented 2 years ago

Right... Let's start at one side and work our way through...

Beat frequency extraction

The process is simple. You take the range-FFT and look for anything above a threshold. The index of the thing above the threshold is the beat frequency of the target.

You can perform various forms of interpolation to get a more accurate number (i.e. increasing range accuracy). This requires superior signal to noise ratio to work properly. The simplest option is to zero-pad your range-FFT and look for the local maximum.

CFAR detection

CFAR = Constant False Alarm Rate.

This has nothing to do with the beat frequency extraction. It is purely for target detection. Yes, you get hold of the beat-frequency by knowing at what index the detection occurred, but, fundamentally, you're doing detection, not frequency extraction. The same algorithm works equally well for pulsed systems, which does not have a "beat frequency".

There are many CFAR algorithms. One of the simpler ones is a range CFAR, also known as a cell-averaging CFAR. The general idea is to perform an estimate of the local noise floor, so that the detector can have an adjustable threshold so that the false alarm rate is kept constant even when the noise is range-dependent.

MUSIC

Read up all about it here. Especially pay attention to the phrase when the number of components is known in advance. In other words, you already know how many targets you have, before you've detected them. Not particularly practical for a general-purpose search radar...

This said, MUSIC does provide you with a super-resolution method that is more efficient than zero-padding the FFT, after you've done the detection.

CLEAN

This is a method by which you can suppress side-lobes by iterative subtraction. More info here.

It is useful when you have a small target in the presence of a much larger target. The side-lobes of the larger target shadows the smaller one, so if you can suppress the large target side-lobes, you can detect the small one.

Others

Have a look at Prony.

dayalannair commented 2 years ago

Great! Thanks for that. I seem to get good results from applying a Gaussian window and increasing the FFT point size. I then apply CFAR though need to do further research into the optimal algorithm/variant. The Gaussian window improves the Doppler estimation. I feel my method is kind of a hack because it doesnt use the algorithm associated with Gaussian windowing.

The MUSIC looks promising, though it seems to only work for simulated data which doesnt have a varying noise floor/clutter etc. From research it seems to be used mainly for DOA, though I have seen it used for beat extraction.

Im interested in finding an algorithm that produces the best results irrespective of computational cost. I will look into the information you provided.

jpt13653903 commented 2 years ago

A Gaussian window has a very wide main lobe. You'll do better with a Hann, Blackman-Harris or Kaizer.


I feel my method is kind of a hack because it doesn't use the algorithm associated with Gaussian windowing

I have absolutely no idea what you mean by this. A window is a window -- it is a function that you apply to the samples. There is no "algorithm" associated.


MUSIC is used for DOA because it is a super-resolution technique. It is a way to better the accuracy. It is not a frequency-extraction algorithm.


I then apply CFAR though need to do further research...

If you have a flat noise floor you can use a simple threshold above the mean (typically 13 dB). No need for complicated CFAR algorithms.


I should be asking -- what is your purpose in all this? What are you trying to optimise? The algorithm you choose depends on the purpose -- do you want a very low false alarm rate? Or do you want superior range accuracy? Or some other purpose?

dayalannair commented 2 years ago

I want better velocity resolution. Based on this pdf there is a method using an equation for Gaussian interpolation. My implementation multiplies a 200 point echo signal with a 200 point gauss win, and then I take a 1024 point fft. Not sure if this actually adds freq info, but I get good results: image Comparing the range and velocity estimates, I now get two 'clutter' frequencies which I filter out with an 'if' statement (if fdop> fdclutter). Also filter out negative doppler (targets moving away) using that condition.

Before gauss win with nfft = 200 the speed resolution is coarse:

image

After G win and 1024 nfft and doppler filtering:

image Real data collected from cars passing by on UCT road. My noise floor is not constant, so CFAR is needed. Currently have a script to compare SOCA, GOCA, SO and CA CFAR in MATLAB for various guard and training cell and false alarm rates. Still need to study these things, so far just testing them out.

Do you recommend a specific CFAR and window for 24 GHz (triangle) FMCW automotive radar? Using a 4x4 patch antenna with each row of patches connected (not 4x4 independent)

jpt13653903 commented 2 years ago

Two main problems before we get to the actual velocity estimate...

  1. You have more than 1 target
  2. You're using a window with sidelobes only 25 dB down from the main lobe.

You're not doing Doppler processing, so, to solve the first, you'll need to include a tracker, so that you can perform the velocity estimate using some form of predictive filter. You can probably get away with a simple α-β-filter, but I think you'll get better results with a predictive one, such as:

The point is, you need to filter one target at a time, so you need a tracker to separate the targets.

For the second problem, play around with different window functions.

jpt13653903 commented 2 years ago

Now... For the velocity estimate itself.

I can do better with my naked eye than whatever you're using, so do something else.

I think you might be trying to perform velocity estimation based only on one detection at a time (i.e. one triangle of up and down sweep). You'll struggle to get anything worth while out of that. Better to use the information from multiple sweeps (i.e. track and filter)


Not sure if this actually adds freq info, but I get good results

You're performing interpolation, so you can get a better estimate of the frequency by searching for the peak. The fact that you do not have a point-target makes this not work too great -- the car is too big, so you actually have a distributed target.


Comparing the range and velocity estimates...

I don't know what you mean by that, so I can't comment.


Be careful with what you call "resolution". Resolution is the concept of "how close together can two things get to each other before I perceive the two as a single detection?"

You are performing interpolation, not increasing the resolution.

jpt13653903 commented 2 years ago

Now for the CFAR.

You do not have any problem with detection. Your signal to noise ratio is ENORMOUS (I'm estimating about 35 dB from your figure). Don't worry about it. You can probably use a fixed threshold with no worries about false alarm rates.

dayalannair commented 2 years ago

As suspected, I am lacking a fair amount of understanding 😢 .

You're not doing Doppler processing, so, to solve the first, you'll need to include a tracker, so that you can perform the velocity estimate using some form of predictive filter. You can probably get away with a simple α-β-filter, but I think you'll get better results with a predictive one, such as:

Yes, the multi target scenario is definitely causing issues. For the case of two closely related targets, is there a way of separating/assigning the up and down beat frequencies correctly? E.g. in the figure below, the outer two correspond to the same target, but the amplitudes change. My method of taking the target with the largest SNR will cause incorrect assignment of up and down beat frequencies. I will look into the tracking methods you provided. image

I think you might be trying to perform velocity estimation based only on one detection at a time (i.e. one triangle of up and down sweep). You'll struggle to get anything worth while out of that. Better to use the information from multiple sweeps (i.e. track and filter)

For this I assume it means filter out one target (same one each time) for multiple sweeps and integrate the return signal. I think my method of storing the received data (writing each detection to textfile as a string) is slowing the process down, though if I use a large array and store after recording it should improve. How will I know if the signals are coherent without knowing the radar schematic? Is it generally coherent? I think it must be.

Apologies. I meant from the range plot, I could see fixed targets with a very low velocity (corresponding to 100 and 200 Hz), so I filtered using an if statement: if f_doppler > 200 Hz, continue processing, else skip.

You do not have any problem with detection. Your signal to noise ratio is ENORMOUS (I'm estimating about 35 dB from your figure). Don't worry about it. You can probably use a fixed threshold with no worries about false alarm rates.

I think the case presented above is for only one sweep, which is not a consistent result. Will double check. I created a MATLAB script that can plot the threshold, CFAR detections, and FFT and see that detections can be missed if the training and guard cells are not well considered. Also depends on the CFAR algorithm. I set the transmitter feed through section to the last sample on the left (hold) in MATLAB, this should be fine?

EDIT: Im not sure where you are seeing a 35 dB SNR? Seems like you are referring to the transmitter feedthrough image

dayalannair commented 2 years ago

You have more than 1 target

The multi target case is shown below. From what Ive seen now on Kalman filtering, would the result not just eliminate the further away target as shown below? I guess its more for the Doppler estimate accuracy... image

dayalannair commented 2 years ago

The guassian window seems to offer the best ML width. SL ratio not that bad? Kaiser probably needs narrower shape factor image

dayalannair commented 2 years ago

DSP Block Diagram

I feel like something is missing... main issue is the multi target. I think some pairing algorithms exist for up and down beat frequencies though Im sure there will be exceptions where it fails. MSc_radar_processing .

jpt13653903 commented 2 years ago

Im not sure where you are seeing a 35 dB SNR? Seems like you are referring to the transmitter feedthrough

Ah -- I understand. You're "transmitter feed-through" side-lobes are killing your detection. You really need to get rid of it. The usual method is to match the coupling and LO delays, so that the transmitter feed-through mixes down to DC, so that you can remove it with a DC-block capacitor.

If you don't have the luxury of matching cable lengths, you need to add an SFC filter before the ADC. In other words, filter the transmitter coupling out before you take the FFT.

jpt13653903 commented 2 years ago

on Kalman filtering, would the result not just eliminate the further away target

No

A Kalman filter makes the assumption of a single target.

jpt13653903 commented 2 years ago

The guassian window seems to offer the best ML width. SL ratio not that bad? Kaiser probably needs narrower shape factor

Good to know -- no that you've done the investigation ;-)

I'm personally quite attached to the Hamming -- good compromise between side-lobe levels and main lobe width -- but specifically for the case of designing FIR filters. For detection we often favour a well-tuned Kaizer.

jpt13653903 commented 2 years ago

DSP Block Diagram

I had a quick chat with our radar guru at RRS. He feels the same way I do -- drop the up-down sweep idea and just do proper Doppler processing. This works especially well if you have control over the sweep parameters (i.e. sweep rate, sweep repetition rate, etc.). The up-down sweep thing works well if you have a single target, not when you have many targets in the presence of clutter.

dayalannair commented 2 years ago

If you don't have the luxury of matching cable lengths, you need to add an SFC filter before the ADC. In other words, filter the transmitter coupling out before you take the FFT.

Unfortunately I only have access to the raw IQ data. For interest, what does SFC stand for? Papers Ive found do not. Ive found an STC filter definition which is probably what you mean:

Sensitivity time control (STC), also known as swept-gain control, is a system used to attenuate the very strong signals returned from nearby ground clutter targets in the first few range gates of a radar receiver. Without this attenuation, the receiver would routinely saturate due to the strong signals.

I had a quick chat with our radar guru at RRS. He feels the same way I do -- drop the up-down sweep idea and just do proper Doppler processing. This works especially well if you have control over the sweep parameters (i.e. sweep rate, sweep repetition rate, etc.). The up-down sweep thing works well if you have a single target, not when you have many targets in the presence of clutter.

Fair point. My radar has onboard processing which is limited such that only range is obtainable from Sawtooth modulation. This means/could mean:

I will look into these. For interest, here is the datasheet and SDK manual

dayalannair commented 2 years ago

In terms of interpolation, is it worth the computational cost? Currently considering:

  1. Gaussian interpol
  2. parabolic interp
  3. Cubic spline - not sure if youve heard of this?

Resampling - is there merit in this application? Since sampling is at twice the max beat, would this cause oversampling/aliasing? seems like it only increases processing time due to larger sampling set

jpt13653903 commented 2 years ago

what does SFC stand for?

"Sensitivity frequency control". Same concept as STC, but all your range information is in the frequency of the signal, not in time. An SFC filter is a 40 dB / decade high-pass filter that does for FMCW radars the same thing as STC does for pulsed systems.

jpt13653903 commented 2 years ago

In terms of interpolation, is it worth the computational cost? Currently considering...

I'm not saying you "should". I'm saying that what you are doing at the moment is interpolation. Zero-padding an FFT is one of the easiest ways to perform interpolation.

And yes -- it's worth the cost, especially in your case.

jpt13653903 commented 2 years ago

My radar has onboard processing which is limited...

You don't need to store everything. Do the range FFT, then detection, and then store only the range samples you're interested in. After many sweeps (say 128 or 256 or so), do the Doppler FFT.

Also, the process only works if your target stays within one range-bin for the duration of a burst, so you cannot make the bursts too long. You can work out what your integration time should be based on the range bin size and expected target velocity. And yes, increased burst time gives you better Doppler resolution.

The device gives your the raw I/Q output, so you don't need to use their velocity estimation.

A higher sweep repetition rate gives you a higher blind velocity (i.e. the point at which the Doppler dimension Nyquest wraps around.) Think about it in terms of sampling theory -- every sweep performs one Doppler sample.

You're dealing with very slow moving targets, so you don't need a fast sweep repetition rate. Do work it out though.

Maximum range in an FMCW system is dependent on the sweep rate (i.e. the slope of the frequency sweep), in combination with your ADC sampling rate. Go through the exercise of drawing the pictures and working it out.

The radar is coherent if your triggering is stable. In other words, the spacing between consecutive sweeps must be constant. To test: use a test target of constant velocity and perform the Doppler FFT. The more "sharp" the spectrum is, the better the sweep-to-sweep coherence.

The above test / measurement is not the best one, but it is relatively easy to do and it will show you if you have sweep-to-sweep coherence or not.

dayalannair commented 2 years ago

The above test / measurement is not the best one, but it is relatively easy to do and it will show you if you have sweep-to-sweep coherence or not.

Ive just run the program with data stored in an array/list:

return_code, results, raw_results = uRAD_USB_SDK11.detection(ser)   
# Extract results from outputs
I.append(raw_results[0])
Q.append(raw_results[1])
t_i.append(time())
period = t_i[len(t_i)-1]-t_i[len(t_i)-2]
print(str(period))

The period is not consistent. This surely means the results are not coherent?

some results: 0.008167505264282227 0.015860795974731445 0.007972478866577148 0.008144617080688477 0.015987634658813477 0.00786447525024414 0.008001327514648438 0.016076087951660156 0.007973432540893555 0.015952348709106445 0.008043766021728516 0.008137226104736328 0.015816926956176758 0.008045434951782227 0.008123636245727539 0.01590704917907715 0.008077383041381836 0.008093833923339844 0.015837430953979492 0.008001327514648438

This is for triangle modulation, 200 + 200 samples (up, down), bw = 240 MHz. period is 1ms per sweep so 2ms. Theoretically, we should see 2 ms between results. Must be non coherent? #:cry: So no Doppler processing, though I will try it out anyway.

dayalannair commented 2 years ago

The radar is coherent if your triggering is stable. In other words, the spacing between consecutive sweeps must be constant. To test: use a test target of constant velocity and perform the Doppler FFT. The more "sharp" the spectrum is, the better the sweep-to-sweep coherence.

If I use a deterministic system to trigger, results may be coherent. This depends if the radar system itself i.e. controller sending the data is also coherent/deterministic.

jpt13653903 commented 2 years ago

return_code, results, raw_results = uRAD_USB_SDK11.detection(ser)

Yes -- If you software trigger the thing you can't do Doppler processing. Pity.

jpt13653903 commented 2 years ago

This said -- your timing results are very consistently 8 or 16 ms, so you might be able to get away with placing the samples on a 8 ms grid and filling the samples you don't have with zeros. It'll give you higher side-lobes (because of the gaps in the FFT), but it might not hinder you too much...

Worth investigating...

dayalannair commented 2 years ago

Yes! I subconsciously noticed that (why I sent that many outputs 😃 ) which is very interesting. I will look into that. Also will try to find the root cause, could it be OS or program related? or hardware...

Seems like its randomly 8 or 16. needs an if statement to see where padding is needed.

jpt13653903 commented 2 years ago

It's OS scheduling related. Could also be related to the resolution of your timer function. Or any one of a number of reasons.

dayalannair commented 2 years ago

FFT of Sawtooth data

Seems like some coherency/partial coherency as the phase matches up around the target frequency? or am i misinterpreting...did not apply the grid fill yet. image

dayalannair commented 2 years ago

Some time stamps for consecutive sweeps are the same, meaning I need to increase precision. However. here it says Windows can only achieve around 16 ms of precision, similar to what is said here. This could explain the results. Need a more precise timing function and also need to increase precision of the time to string conversion.

time_ns may work:

time_ns does not provide nanosecond resolution, it just returns the time in nanoseconds. This may or may not improve the resolution. Measured clock resolutions for Linux and Windows are documented in PEP 564. For Windows, the measured resolution of time_ns was 318 us. – wovano Jul 6,

Time resolution of 318 us should be fine.

dayalannair commented 2 years ago

timeit seems better though cant find a comparison. Will try both

EDIT: timeit seems to be for timing code and so may not be suitable for the real time process timing

jpt13653903 commented 2 years ago

Windows can do much better, you just need to find the Python library to link into that.

Side note -- be careful when you quote stuff -- if it includes an at-mention you're causing that person to receive a notification and pulling them into the conversation.

dayalannair commented 2 years ago

Ahh yes my bad. time_ns doesnt seem to offer much improvement in resolution:

0.00810393600000000 0.00797388800000000 0.0158850560000000 0.00819020800000000 0.00797798400000000 0.00800435200000000 0.00798003200000000 0.00799820800000000 0.00808422400000000 0.0158417920000000 0 0.0155599360000000 0 0.0156592640000000 0.0158120960000000 0 0.0158668800000000 0 0.0153858560000000

As you can see, increasing str conversion precision does nothing as predicted. My concern is the 0s which mean the time stamp was not updated. The values above are the differences between consecutive time stamps.

Windows can do much better, you just need to find the Python library to link into that.

I would have to run the C code in python/interleave it somehow or continue hunting for a better lib. Apparently Linux can improve timing precision.

jpt13653903 commented 2 years ago

Yes -- Linux implementations generally use the clock_ family of functions, which are much better resolution.

dayalannair commented 2 years ago

I think it may be best to just assume/not use the time data as it is just for info sake. From real data, see how many sweeps the target remains in a range bin and choose frame size accordingly.

I will try running it on Linux as well for interest.

dayalannair commented 2 years ago

Last thing for today

Does the phase below indicate (rough) coherency? To be it seems like it. I was walking away from the radar so there should be some phase shift for Doppler processing. Giving it a try now. image

jpt13653903 commented 2 years ago

Difficult to say...

Plot the resulting range-doppler map

dayalannair commented 2 years ago

Okay, so I set it up to work based on a specified number of sweeps per frame. I took 512 sweeps of me walking away from the radar. I then plot each frame per second and observe the target moving away. It seems that larger sweep sizes offer (other than incr dopp resol) clearer depictions of the scene. I think I should stick below the theoretical number of sweeps per frame, calculated to be aroud 60 based on the time a target stays in a range resolution cell (+-1.5m) while travelling at the max speed.

Range in meters on y axis, still figuring out conversion for x axis. 50 sweeps per frame, last frame: image

I think the transmitter feed through is contributing to/causing the smears. I will:

  1. Add mechanical separation
  2. filter

Setting it to zero just eliminates DC and not smearing

jpt13653903 commented 2 years ago

Please label your axes... I assume range on X and Doppler on Y?

The transmitter feed-through will show a large peak at 0 range, which I don't see. I only see the large peak at 26 Range and 0 Doppler, which is probably clutter...?

The smearing can be one of numerous things... The cars wheels do not move at the same speed as the car, to name but one example.

dayalannair commented 2 years ago

Range in meters on y axis, still figuring out conversion for x axis. 50 sweeps per frame, last frame:

Range is the y axis, so the bright point is the feedthrough. I would think there is less smearing due to multiple Dopplers as the radar was observing me walk away from it (single moving component). How would I convert the second fft dimension axis to velocity? Apparently phase is related to angular frequency which I will look into. Currently trying to get matlabs 2D CFAR working