Real-time heart rate calculation
Show older comments
I have a real-time(phonocardiograph) plot which displays data from an electronic stethoscope via the COM port. I am trying to calculate the heart rate but having little success. I have created a uicontrol on the figure to display the BPM and have included the formulas to count the peak signals(beat_count). At the momement 'BPM:' displays on the plot figure but without a value so I am unable to get any BPM_avg reading.
a = arduino('COM4','Due');
y=0;
hPlot = plot(NaN);
intervalSize = 200;
currentInterval = 200;
t = 1; % number of samples
atInterval = 1;
beat_count = 0;
quitbutton = uicontrol('style','pushbutton',...
'string','Quit', ...
'fontsize',12, ...
'position',[10,2,50,20], ...
'callback','quitit=1;close');
quitit = 0;
bpmtext = uicontrol('style', 'text',...
'string', ['BPM: '],...
'fontsize', 12,...
'position', [80, 2, 100, 20]);
while(1)
k = 1;
while(t<currentInterval)
b=readVoltage(a, 'A0');
y=[y,b];
if ishandle(hPlot)
set(hPlot, 'YData', y);
else
break; % break out of the loop
end
xlabel('Samples')
ylabel('Voltage')
title('Phonocardiogram')
axis([currentInterval - intervalSize,currentInterval,0,3]);
%grid
t=t+k;
pause(0.002)
end
for m = 2 : length(b)-1
if(b(m) > b(m-1) & b(m) > b(m+1) & b(m) > 2.4)
%disp('Prominant peak found');
beat_count = beat_count + 1;
set(bpmtext, 'string', ['BPM: ',...
num2str(BPM_avg,4)]);
end
end
currentInterval = currentInterval + intervalSize;
atInterval = atInterval + 1;
if ~ishandle(hPlot)
break;
end
fs = 500;
N = length(b);
duration_in_seconds = N/fs;
duration_in_minutes = duration_in_seconds/60;
BPM_avg = beat_count/duration_in_minutes;
end
3 Comments
Star Strider
on 22 Mar 2017
I cannot follow your code. You likely need to be doing real-time signal processing of the sounds to isolate the ones you need from the ones you do not.
Note that the Phonocardiogram consists of as many as 4 usually distinct sounds, and varies over the respiratory cycle. It is normal in healthy people for the second heart sound, ‘S2’ (closure of the aortic and pulmonic valves), to be ‘physiologically split’, so that they do not occur at exactly the same time. I refer you to any textbook on cardiac physical diagnosis for the details and a discussion of the relevant physiology.
Iqra Chaudhary
on 14 Dec 2019
How to find heart beat rate of an audio wav file without using arduino?
Manju
on 25 May 2024
can u give me a code the shows the variation in the heart rate with respect to stress level
Accepted Answer
More Answers (1)
I'm going to put this as an Answer to emphasize the point again...your timing code of just a single read in a tight loop shows that whatever is going on in the Matlab interface to the Arduino readVoltage routine at least on your system it's just killing performance. I again STRONGLY suggest you read the blog of another experimenter with the Due for some very informative background on using the Due and on data collection in general.
It's clearly possible when running code directly on the board to run at much higher rates; as said above I don't understand how there can be as much latency using the Matlab interface, but until you can find some resolution there what you're trying to do just isn't going to work.
PS. I wonder if since there are 16 channels on the Due if they're all active the sample rate gets divided down by polling them all and just returning the one asked for, maybe? I didn't see where it says whether the A/D sampling is simultaneous or sequential, only that is 1MHz capable. As have emphasized off the board in our sidebar conversation, you need to dig into the guts of the board and the Matlab interface much more deeply to figure out what's going on here.
I don't suppose it could be some sort of communication issue slowing it down, is it? Can you test that somehow?
Mayhaps ask on one of the demo blogs from TMW what you should expect in a tight loop or, if you or your university has support, submit a support request and see if TMW will respond.
ADDENDUM
Btw, the last timing you sent is most curious...I pasted into a variable tacq and one gets the following interesting figure--

For some reason, the acquire time goes down with each loop; dramatically at the beginning and then looks like is approaching a steady-state rate, maybe. I don't know what would cause this behavior; seems very odd, indeed.
Why don't you try the above loop with nSamples at some sizable number and collect a waveform or two and see what you get? If you approach 60 hz as it appears this is, how long a time series do you want? Multiply N seconds by 60 and collect that many samples. Of course, do the preallocation and collect both the reading and the transpired time as shown and
plot(t,y)
to see what you get. Post those results back here. It would, of course also be informative to just collect the timing information for several trials of a large number to see if the above is consistent and what steady-state rate is. You might want to "pretrigger" by some number of samples before recording "real" data...
16 Comments
bilal malik
on 30 Mar 2017
TMW--The Mathworks, owners/publishers of Matlab.
How could taking a set of data with a loop similar to your above to see what the result is on timing be out of the scope of what you're trying to do?
If you're talking about the link, the importance there is the background insight into what the DUE board is doing and at least partially why speed isn't what you think it might/should be--altho there's still a big difference from what he gets directly from the device vis a vis what you're showing through the Matlab gateway (another level of distance from the board).
IF your objective is as was stated, the scope has to include figuring out a way to retrieve data fast enough and on a repeatable enough timing to be able to analyze it. So far it doesn't seem that either of those objectives may be achievable this way; does that change the scope or does it mean you've got to find another way to get to the end objective? Have you had the discussion with your advisor I suggested on the problems you're encountering yet?
Specifically, again, what I'm suggesting is put the following in a script file and run it...
a=arduino('COM4','Due');
n=input('How many seconds to run roughly)? ');
nSamples=60*n; % assume roughly 60 Hz from past experience for now
t=zeros(nSamples,1); y=zeros(size(t)); % preallocate for time, data
tic % start the timer
for i=1:nSamples % loop over the number set
y(i)=readVoltage(a, 'A0'); % read the A/D, save
t(i)=toc; % save elapsed time since tic accumulative
end
plot(t,y) % see what you got
save('Timing','t','y') % save to .mat file for further looking at
dpb
on 30 Mar 2017
a=arduino('COM4','Due');
Oh, what is COM4? Is it serial or USB port? If it's serial and set a a pretty low baud rate, that could be where all the time is being lost.
bilal malik
on 30 Mar 2017
dpb
on 30 Mar 2017
It's waiting for you...answer the question. :)
>> n=input('How many seconds to run roughly)? ')
How many seconds to run roughly)? 5
n =
5
>>
If it isn't showing the question, you've got a coding error somewhere...
Well, normally you create a serial port object via serial and then connect the device with fopen using that object handle. You can then set the properties of that port object. That somehow must all be inside the arduino code that creates that object that you address. I see there's then a read-only port property but don't see anywhere that is allows setting baud rate so I presume it must just be default. A gargle search finds references to 115200 baud, though, so if it's that high it's not up to USB speeds but wouldn't seem like it should be limiting to as slow a cycle rate as you're seeing. So, "I dunno"; seemed like a good question.
I did see here that you can load user-written libraries which apparently are way to package code to run on the Arduino; looks like that is the way you get to more in-depth control with Matlab-only interface.
I thought there was an optional USB terminal port; doesn't appear as it is set up that way on looking more; the USB port seems to be available to stream data over but guess it would be writing one of these libraries to be able to load that to the Due to execute.
But, try
a.port
after you execute the arduino command and see what it shows you...
I'm still interested in seeing what you get if you run the timing 'spearmint above a couple of times (note you'll have to change filename or you'll overwrite previous .mat file so be sure you're done with the first set of data or make that change before doing it again).
bilal malik
on 31 Mar 2017
<arduinoio setup-and-configuration>. See link to <propertylist> at that page. Looks like it's 'Port' not 'port'...need to be more aggressive in looking for minor errors such as typos, it seems. (Apologize if seems "preachy", but ... :) )
What happens if you just open a connection at the command line and type the object name used? Will it not echo all the properties for you? You're not taking much advantage of the interactive nature of Matlab if you don't use command line extensively for such exploration (or use the debugger inside a script/function the same way when the object/variable is in scope).
For the data, well, looks smoother than what was before; if as we presume it's only approximately 60 Hz sample rate and that's the beat rate, you'd expect to asynchronously hit/miss the peaks and that's pretty-much what it looks like it did.(*) Did you look at the timing information to see what the timing actually is/was? Look at diff(t) to see (since we captured the running time from trigger, that's the dt); a plot of those against observation number and 1/dt for instantaneous sampling frequency would be interesting to see how it behaved with time; does the same initial delay occur and it reaches a near steady state? If so, what's that final sample rate? How uniform is that sample rate? So many questions, so few answers (so far)...attach the .mat file.
Again, think about what that experiment can tell you and try to get the most out of it you can...it's probably going to confirm what we think we already know that unless you figure out where the overhead is that polling isn't going to work, but it would surely seem you'd want to know just what it is that it can and is doing...
It seems hard to have to coax such information/explorations out...
(*) Thought 'spearmint--What would happen if you had a pure sine at 60 Hz and sampled at 60Hz? Now, change the sample rate to 60Hz +/-few Hz. How does the sampled waveform change? Now repeat both of those with a changing phase shift of where the initial sample is taken and you should begin to understand what you're looking at.
"Not sure it makes much sense to me."
Well, besides the comments above, it's clear directly from the plot that the previous estimate of somewhere in the 60 Hz polling rate is maintained; the overall end time is 6+ sec in first; interestingly however it's 5+ sec in the second (presuming, of course, you asked for 5 sec in running both and entered that value without a typo); the internal assumption of needing 60X that number of samples to get the time span seemed to hold up moderately well. The question remains then why there still seems if both were 5x60-->300 samples, the 1 sec overall time difference. Also, there appears to be a sizable lag for the first acquisition; granted I put the initial tic before the for loop to avoid the need for a branch around it inside the loop so there should be a minor delay, but wouldn't expect what it appears is the length of that one. Ergo, one must presume there's something going on regarding the acquisition about timing for it. And, of course, what happens if you subsequently sample for longer time periods; how stable will if finally get and what is the final true steady-state rate as one would presume must be?
The second trace almost looks like lost connection to the device or somesuch; looks mostly like noise but still the timing is useful as that was the prime thing looking for here. ADDENDUM: Actually the above thought 'spearmint well may play into the answer here as well...
Oh! The other thing need to do...rearrange the test code so the Arduino connection is made first, before the timing. Then, pass the object to the test function. This way if you re run multiple times you won't have whatever initialization effects there are on the board impacting the beginning of the test; as it's written now you're going through that every time; only need to do that once a day or after a major crash or something, not routinely.
I didn't look to see; is there some cleanup code to close a connection; is there, by any chance, any problems arising with having issued a zillion arduino commands to create a connection without cleaning up after prior ones before making a new one even though you reused the Matlab variable? Don't know, just a thought...
dpb
on 31 Mar 2017
Another bunch of questions
- Is this device connected to the Arduino continuously?
- Is it externally powered or get its power from the Arduino?
- When you start a sample, are you turning it on the starting, is it already on and you just start acquiring, ...?
IOW, am wondering if there are startup transients owing to the device that explain the early spikes or whether that's a fignewton of the board somehow that the A/D has been running for ages with some spurious voltage on the input that has to be bled off before there's anything real or perhaps it's open and then a capacitive load applied or the like...what is the actual sensor in the device, what's it input impedance; could you be loading-down the A/D (not that think that would affect the timing altho I noted in the link I posted there's discussion there of there being settling times built in; if that were to be firmware based on input signal stability rather than just a fixed number of clocks such effects could be major. Not likely, but not out of realm of possibility, either.)
Again, more questions than answers...
bilal malik
on 31 Mar 2017
By "device" I mean the instrument you're reading with the Arduino.
You did read the comments of the previous link about stability of the reference power and all that if just use the USB port for the 3V reference, right? Have you tried a known input voltage from a reference and then from a signal generator to ensure you can actually measure a known input?
Reason for the last is apparent from attached which is from the last of the three files above--

The top of those is
plot(t,-y)
The middle is same except expanded by
xlim([15 30])
to observe the first peak/disturbance area more closely. The third is from
dt=diff(t); % the sampling interval
i1=find(t>15,1); % the index corresponding to 15 sec
i2=find(t<30,1,'last'); % the 30 sec mark
plot(t(i1:i2),dt(i1:i2) % and the sampling interval corresponding
What is most disturbing is that all the peaks in the sampled waveform are exactly correlated in time with an increase in sampling time as I had somewhat postulated might be going on earlier.
I think this shows the measurement is mostly noise induced from the PC and/or the Arduino board itself; NOT the monitor device. Is it, in fact, even hooked up or operating at all? Doesn't look like it to me from here...
Looks like time to break out your EE training... :)
ADDENDUM
Finishing up on some looking at the sampling rate...


First of those is just the instantaneous frequency (1/dt) from each sampling interval vs time which shows it is basically a random value about an approximately constant mean (the solid line is a the fitted linear regression; altho the slope isn't zero; the line is
fhat=0.0177*t+58.6
>> [fhat(2) fhat(end) fhat(end)-fhat(2)]
ans =
58.6928 60.8584 2.1656
>>
(The reason for f(2) for first is there's one less element in f=./diff(t); than in t so for convenience I augmented it with the predicted fhat(0)). Colors weren't easy but there's a barely discernible dotted white line above/below the predicted line that's the 1-sigma bound. Those values were
>> sd=std([fhat(1);f]-fhat); % sd from regression
>> [mean(f) sd sd/mean(f)*100]
ans =
59.7833 6.7612 11.3094
>>
The second is just histogram that does show essentially a normally distributed value around the mean with a small cluster of longer time intervals (ergo, lower sampling rate). These would be expected latencies inherent in the occasional task-switching, etc., under Windows when polling.
While if it were sufficient to see the information you're looking to capture you might be able to live with the jitter to simply display waveforms, it would be difficult to do any very precise timing. Since the actual beat rate is roughly the same rate as this sampling rate, it's pretty clear that this won't work as a data collection scheme as we had pretty much already concluded; just needed to really make sure that is so without other distractions.
As noted before, perhaps there are some ways to make the polling run faster and it does really seem incredible that there is such a delay going on, but unless you can figure out where all the clocks are going, this is "just the way it is"--you've pared the loop down to nothing but the request for a new sample and that's as fast as it ran and it didn't do anything different over a period of time. That indicates there isn't some initialization going on that the first couple of very short samples looked like might be and that the achievable rate might possibly be better. I can't for the life of me imagine why it should be so slow, but can't deny the results.
And, since the top two plots are -y, what are peaks there are actually dips; I first (as mentioned in passing above) thought you had likely inadvertently switched polarity on the device and so thought should be the inverse of the readings.
Then when comparing the sampling rate, observed the coincident nature of the peaks illustrated. But, what may actually be happening is that the board is loading the input during sampling and pulling down the input signal and the values are actually as recorded. That would indicate the USB power supply isn't able to supply the device with all the current it needs so the voltage drops.
Fixing that with external reference/power supply wouldn't solve your timing problems, but at least you might then get clean data. As mentioned before, looks like it's time to go back to ground zero and follow a very precise sequence of steps to uncover just what is actually happening on the board. The outline in that previous link of what that individual did to confirm his results is, as noted before, a pretty good tutorial on how to approach a new piece of hardware; particularly one such as this that is open source and while a neat piece of gear, NOT guaranteed by the manufacturer to do much of anything as, say, one can rely on a Tektronix or Agilent lab instrument to actually measure what it says it measures; you really don't know that here without confirming it independently.
PCs (and computers in general) are just so rife with noise issues that without careful attention to grounding, shielding, etc., etc., it's easy to measure something but not necessarily anything much related to what it is one is trying to measure. I'm sorry to say it looks a lot like that may be what's happening here.
What does that app use for its collection; you should model your input arrangement after it to at least get a clean signal; then the timing issue can be worked on.
Unless the investigations into what the arduino.m initialization code is doing bear fruit as to how to increase the response of Matlab via readVoltage or you get some feedback from TMW or another user of the Matlab interface that has solved rapid data collection with Matlab-only code, it would appear to write a user library and load it or use the Arduino-supplied IDE are your choices.
However, if the data being sampled are still just noise, how fast you can sample it won't change that it's still only noise, so need to solve the first problems first.
This is, however, a very good introduction into what "the real world" looks like it seems...things are rarely as simple as one would hope and learning to attack in an orderly fashion is probably the most important lesson you can learn from this project.
bilal malik
on 2 Apr 2017
dpb
on 2 Apr 2017
I was asking about how it is physically connected to read the device? How the data transfer occurs is immaterial to trying to solve your instrumentation issues.
What was in the arduino.m file out of curiosity? Is it m-code too or did it drop down into a mex file or does it show an Arduino "sketch" or what?
Again, I'm still curious about whether there was anything attached to the Arduino on those timing runs that was powered up???
bilal malik
on 2 Apr 2017
Categories
Find more on Real-Time Simulation and Testing in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!