Refactored Ppg for frequency based algorithm.

[Ceimour]

Changed driver settings in HRS3300.cpp for better signal to noise at 10Hz. Changed measurementStarted delay to 100ms in HeartRateTask.cpp. Removed unused ppg methods from HeartRateTask.cpp. Added ambient light sensor acquisition and handling in HeartRateTask.cpp. Added FFT and interpolation Ppg support code from Arduino FOSS implementations. Removed currently unused time domain based code.
Changed all doubles in PPG/Arduino code to float type so that FPU is used.

[Ceimour]

Complete refactoring the heart rate code. This has been tested against a predicate device and looks to be typically within +/-1bpm when heart rate is stable. This is my first pull request so I would appreciate any guidance on moving forward.

[minacode]

There are files under the Apache License 2.0. Are they compatible with this project?

Edit: seems like adding Apache 2.0 code to GPL code is compatible, but I am no lawyer.
Link

[minacode]

Since this PR contains a lot, any sort of documentation helps you to potentially get this through.

For the start it think you should describe a little why you did what you did to get an overview of the different changes. The maintainers don't have much time so presenting your ideas and reasons well is very helpful.

Then we could discuss if we like the changes and if they could maybe separated into smaller pull requests.

Also please tell us what this PR improves.

[Ceimour]

Thanks for taking a look at this! According to Apache.org (first 2 paragraphs) we are good to use Apache V2.0 within GPLV3, just not the other way around: https://www.apache.org/licenses/GPL-compatibility.html

[minacode]

Ok, cool.

Are you aware that your replies contain some email stuff? No offense, just wanted to tell you in case you didn't notice.

[Avamander]

If replying over email, it's better not to quote the previous message. I've also manually truncated it from people's comments so the thread would look nicer.

[minacode]

Your description sounds good as far as I read it. More reliable heart rate data would be a very good thing! Can you restore the formatting that it seemed to have when you wrote it? ;) You can just edit the comment if that is possible for you.

[Ceimour]

Motivation:
After receiving my new watch I was very impressed with the look, feel and quality of InfiniTime. I expected a hobby project but this is something different, something I actually want to use! I experimented with the various features and I noticed the heart rate monitor accuracy issue (like issue 532 below). I set up a dev environment and examined the incoming data. Much of the time it looked really good after Daniel Thompson’s filtering was applied, but then there would sporadic instances of high amplitude noise within the heart rate frequency range. So it’s the old garbage in, garbage out problem. Since I’ve had some experience working in the frequency domain I thought it might be easier to discriminate against the noise using frequency analysis.

Changes:

• No fundamental architecture changes.
• Previous code in the Ppg methods has been replaced to implement a sample buffer with an overlap window (currently at 8 samples), frequency analysis, a bandwidth filter that covers a wider heart rate range and various small methods for peak search, averaging, etc.
• All floating point types are 32bit float to take advantage of the FPU.
• No dynamic allocations.
• Adds 1-2% to stack and flash. Not sure of the exact amount.
• Changed driver settings in HRS3300.cpp for better signal to noise at 10Hz.

Improvements:

• First HR result in as little as 6.4 seconds
• Updates every 800ms
• Good stable HR accuracy (when HR is stable compare with a predicate device and should be within +/-1bpm)
• Tracks HR fairly well when the rate is changing slowly (experiment by waiting for a stable HR, take a breath, hold and HR should rapidly decrease)
• Uses spectral averaging to increase signal to noise ratio
• Uses sample averaging to hang on to previous HR (about 9 second history and drops to zero when the history buffer is all zeros)
• Spectral and sample averaging allow the algorithm to ignore quite a few noise events without going to zero at each noise event

Cons:

• Changing HR will lag a bit behind some predicate devices
• It is still possible to get errant readings due to noise although these generally go away quickly as the HR signal takes over

Pull Request:
I think this needs to be a single pull request because all the code in the algorithm is integrated to produce a single output. Although it looks like it addresses 532 and coincidentally 1405 and 1092 below.

Applied Algorithm (Ppg::HeartRate):
The algorithm is basically similar to the same processing chains found in many papers for extracting information from PPG:

1. Acquire samples from the HRS3300
   * If acquisition is reset or initial, acquire all 64 samples before processing (before going to step 2)
   * If acquisition is not reset, acquire 8 (overlapWindow) samples (before going to step 2)
2. Detrend the dataset
   * Remove the slope
   * Convert to the difference between each sample
   * Turns the problematic large, DC shifts into single sample spikes
3. Bandpass filter the data 30 to 240Hz
   * Removes low frequency data that can interfere with low HR detection
   * Removes the high frequency spikes from the detrending step
4. Apply a window function
   * Make the signal appear continuous to the FFT
5. Calculate the FFT (Fast Fourier Transform)
6. Calculate the power spectrum (ComplexToMagnitude())
7. Apply result to spectral average
   * Helps to enhance the more frequently occurring heart rate frequency while diminishing random noise
8. Check signal to noise ratio (SNR)
   * Nothing complex, just compare max to mean
9. Search for peaks in the spectrum only if SNR is good,
   * Look for peaks above the peak detection threshold (currently 60% of max)
   * If only one peak above threshold, return the peak position in peakLocation
10. Check peak width
    * Wide peaks can occur from noise or rapidly changing HR
    * If peak is too wide, set peak location to 0.0
11. Check if HR is within limits (currently 40-230bpm)
    * If out of range set HR to 0.0
12. If HR == 0.0, reset spectral averaging for next pass
13. Set Ambient Light Sensor threshold
    * Used in next pass to determine if ambient light was detected
    * HeartRateTask will reset the acquisition if ALS threshold is exceeded
14. Add current peak location (zero or otherwise) to Ppg::HeartRateAverage and assign the result back to peakLocation
15. Return result based on value of peakLocation
    * > 0.0 is valid HR
    * -1 tells HeartRateTask to reset acquisition
    * 0.0 is no valid peak found this pass
16. Goto 1 via HeartRateTask

Execution time ~20ms

Relevant feature requests:
Heart rate measurement accuracy
#532

Increase the heart rate data logging interval
#1405

Ring Buffer for Heart Rate Data
#1092

Potential feature requests that can be addressed in future pull requests:
Zero out HRM when stopped
#682

Bug: Heart rate monitoring continues while charging
#594

Sleep tracking
#307

Heart rate measurement is stops, when the screen is turned off
#183

[Ceimour]

Reformatted pull request description. Not sure what the SOP is here, should I delete the old one and this comment as well?

[minacode]

I don't think it is necessary to delete the comments. Thank you for posting it again! 😊 This will help a lot.

[yehoshuapw]

I didn't quite have as much time to test this as I wanted too, so this will be brief (making the hrs better was something I never quite got to, so thanks!)

it seems to have some issue - every about 10 seconds or so, it resets back to Not enough data.

[Ceimour]

Thanks for testing it! I know the signal gets polluted with quite a lot of noise from movement, heaving breathing, fast heart rate, etc., basically anything that adds a lot of energy. I have also noticed that I have hard time getting a reading when it is cold out and there is a large temperature difference between the top of the watch and my wrist. Anyway, if it continues to be an issue the averaging buffer size can be increased to mitigate. I hope that you and others will have time to try it for a longer period. I've only tried it on a couple of other people, so really looking forward to the feedback. We might also want to think of different instructions to the user like, "Not enough data, please remain still..." or something like that.

[Riksu9000]

From simple testing this seems to work very well. The fact that it resets to zero when I'm moving around and presumably can't get a clean signal is giving me confidence that the reported values are more trustworthy, and the faster reporting of the heartrate is a much better experience.

Can you explain how you can get better results with less frequent data? Is more data just not necessary? Will the code work as well if the frequency was increased, or would it require tweaking?

How is the average power consumption compared to the previous code?

I made some comments on the code. The suggestions should be applied in all applicable places. Some of these suggestions are from clang-tidy.

Libs should be places in the libs directory. Can you add the libs in one commit, and make the changes in another, so we can keep track of the changes? Libs also should not be reformatted.

[JF002]

I tested this PR, and it looks very promising. See the following video (the sealed pinetime on the left is running InfiniTime 1.11, the wired devkit on the right is running this PR):

new-hr-algo-2.mp4

The value is updated more often, it detects very quickly that I removed my finger, and it doesn't sense any heart beat in the breadboard, which is nice! 👍
I do not have the equipment needed to check the accuracy of the measurements, though.

I think there are a few things that must be improved in the code, but I really hope we'll complete the review, validate the results and approve those changes :)

I fully agree with @Riksu9000's review :

• Document how this signal processing work and why it's better than the previous one (do not go into much details, but give a bit of info for those like me who do not know much about signal processing)
• Move the 3rd party libraries in the folder libs.

I also added a few comments in my review. I know there are quite a few things to change, but they are mostly minor changes needed to conform with our coding standards :)

We must keep in mind memory usage : as I wrote in a comment : this implementation uses 4x time amount of RAM memory than before (from 200B to 1024+B). This seems to be fine as the compiler does not complain about overlapping regions, but we know that we are already quite tight in RAM memory.

[JF002]

Here is the current measurements while running this branch:

When the sensor is not running, with the medium brightness, the watch use ~10-20mA (15mA avg).
When the sensor is running, the watch draws ~12-30mA (20mA avg). So the sensor roughly draws 5mA in average.

How is the average power consumption compared to the previous code?

The original code draws 24mA avg (12 - 40mA) when the sensor is running, so this new implementation draws 4mA less than the original code, which is nice!

[Ceimour]

Basic Algorithm How and Why

How
Since we are working in frequency space we are primarily concerned with Nyquist and frequency resolution. The highest analyzed frequency is 4Hz (240bpm) so we need to acquire data at a minimum 8Hz. Frequency resolution is 1/TotalTime(seconds). I wanted to keep the total acquisition time to somewhere around 5 seconds for good UX. So I ended up at 10Hz with a 64 sample buffer resulting in a frequency resolution of 1/6.4s or 0.15625Hz, multiplied by 60 it is 9.375bpm. Yes that's poor resolution for pulse rate results! Thankfully the pulse rate frequency data has a nice Gaussian distribution and we can interpolate fairly accurate values. That coupled with pulse rate result averaging gives values that agree nicely with my off the shelf pulse oximeter. So we can certainly run at higher rates keeping in mind that we still need to acquire for at least 6.4s to keep the existing frequency resolution. It looks like it's taking ~20ms to run the existing algorithm, so the hard limit would be ~50Hz. Keep in mind that to keep the existing frequency resolution we would need to acquire for 6.4s resulting in a larger data buffer.

Why
We can get better results using this method because the metrics for determining good versus bad pulse rate data are simpler. We can just look for a single peak above a defined threshold in a power spectrum that has low noise characteristics (simple max to average in this case) and assume that peak is the pulse. Then an output averaging ring buffer is used to store the good pulse rate values and also store failed pulse rates as zero. When the averaging buffer becomes all zeroes, reset the acquisition as if we had just pressed the "Start" button. So basically we just look at the values we think are good results. There are times one may experience incorrect readings but they should be short lived after physically remaining still and allowing good values to fill the buffer.

[lman0]

any news ?
it's the most precise heart rate counter since infinitme creation.
it make the pinetime nearly on par with an spo2 meter (for the heart rate part).
when i say nearly on par , in fact it's about 99% of the time the same value that on spo2 meter .
it's my feeling and when i compare with a owned spo2meter.
i use it for 1 week and the value are way better.

sometime it's restart measuring ...
is there any way to avoid this @Ceimour ?
(is there any movement that triger this behaviour....?)
or maybe it will be resolved with #1092 ?

[Ceimour]

@lman0 , thanks for testing this!
You are correct that movement can trigger a reset in the data acquisition. For a future improvement I have implemented an adaptive filter that will mitigate this behavior somewhat. But since this pull request is already quite large I hesitate to add it at this point.
A reset is produced when the data is noisy (caused primarily from motion, but the sensor sometimes just outputs bad data) for a long period. By long period, I mean more than 10 seconds. The reset text should say "Acquiring, please remain still..." or something like that (another future pull request). The basic philosophy is to attempt to only display the truth. Although it is averaged so it could be the true pulse rate up to 10 seconds in the past. But I think to display values older than that when bad data has been present for over ~10 seconds is misleading.
One more thing that can cause a reset is ambient light changes. Try removing the watch when getting a good reading and you should see an immediate reset.
I think that in the future with wording and display changes it will provide a better user experience. It's all about user expectations!

[Riksu9000]

There are still unaddressed review comments and the PR needs to be cleaned up. Someone from the dev team may end up creating a new PR based on these changes, so we can get it merged soon.

[Ceimour]

@lman0
Thanks for sharing the perception that the pulse rate measurement going into reset while in active measurement mode seems confusing. In this recent commit I have removed that bit of code that invokes the reset while measuring.
Please give a try and let me know if it feels like the correct behavior.
Thanks again!

[lman0]

@Ceimour can you check why the action crashed?
I can't test your new commit since there us no artifact...

[Ceimour]

@lman0
I see that this repeats until exit with failure:
Retrieving check runs named Compare build size on 25fe200...
Retrieved 0 check runs named Compare build size
No completed checks named Compare build size, waiting for 10 seconds...

[Riksu9000]

The PR needs to be rebased and commits that fix mistakes squashed. The libs should be added in one commit, and modified in another, so we can see the changes made.

[Riksu9000]

Mark reviews as resolved when you've resolved them, please.

[JF002]

@Ceimour Amazing job 🥇 , thank you!

I've just checked the flash memory usage : these changes add a bit less than 2KB, which seems much more reasonable than before !

Regarding RAM usage, I couldn't find any dynamic allocation (which is nice), and the whole HeartRateTask class uses less than 1KB. It's probably more than the previous implementation, but these improvements are probably worth it.

This PR looks good to me. I'll ask @Riksu9000 and @FintasticMan to check whether all their conversations (review) are now solved :)

[Ceimour]

@JF002 Thanks!

[Riksu9000]

The last commit contains fixes for the first commit, which should be applied to the first commit to meet the commit requirements, however the entire commit should not necessarily be squashed. Ideally each commit would contain much fewer changes, because the descriptions are long.

[FintasticMan]

Other than the points Riksu has raised, I have no issues with this code.

[FintasticMan]

This looks good to me! I've opened a PR to the FFT library that would allow us to specify the square root function to use without modifying the library: kosme/arduinoFFT#83. I don't think we should wait for that before merging this. We can update the library and put it in a git submodule in a later PR.

[JF002]

@Ceimour Thank you very much for your work and also for you patience and courage while we reviewed this PR!
I'm really happy to finally merge those changes which greatly improve the measurements and usability of the HR sensor in InfiniTime!

[Ceimour]

@JF002 Thanks! And thanks to everyone who tested and helped guide me through the process. You all made the code and implementation much better!

[Anasdec]

I noticed the HR readings drops a lot to zero especially during motion. I was wondering if you have got the time to implement the adaptive filtering above?

Otherwise, is changing some of the parameters such as peakDetectionThreshold, maxPeakWidth, and signalToNoiseThreshold would help?

Thanks.

[Ceimour]

@Anasdec Thanks for your interest in in improving this. Adaptive filtering unfortunately would require the use of the math libraries mentioned above (adds ~8KB). Not sure it would help that much since the fundamental issue is the poor signal quality from the sensor. I came to the conclusion that the ONLY use case for this current implementation is stationary, spot measurements. I still hold out hope that someone out there will find a way to get better signal to noise from the sensor.

You may be able to get a better motion related experience by increasing the spectrum averaging and/or results averaging. The downside is more averaging, more lag in the reading.

If you are really keen on playing with peak and signal metrics, I recommend using a dev unit wired to a debugger such that the watch heart rate sensor can be in contact with the wrist. Then make a small change and observe the effect. Iterate to get desired end result. Warning, it will take up a lot of your time!

[Anasdec]

@Ceimour Thanks for coming back. I'm keen to look into this, but I don't have much experience.

To increase the spectrum averaging as you suggest. Does that mean I need to do the following :

Change: datalength from 64 to 128 and the Hanning Coefficients (via python - numpy.hanning(128) and only use the first 64 values)?

and increasing dataAverage array from 20 to, say, 30?

am I correct?

[Ceimour]

@Anasdec
In Ppg.h:

To increase spectrum averaging adjust:
static constexpr uint16_t spectralAvgMax = 2;

To increase heart rate result history, adjust the array size as you mentioned above:
std::array<float, 20> dataAverage;

Note that as the spectrum average is increased, the peak may spread out (peak width) as the heart rate varies. If you find that you are suddenly unable to get a result, the max peak width may need to be increased:
static constexpr float maxPeakWidth = 2.5f;

Cheers!