indicating when data reaches a certain threshold

Hi,

I have a data set with time in one column and force in the second column during a persons walking trial. I'm trying to establish heel contact by trying to find when rising force value crosses a specific threshold in. I would like to identify the time and force value for when the data crosses this threshold value. I'm still very knew to Matlab and havent't had any luck with using the find function or trying to threshold and could use some assistance for how I can go about establising heel strike and when the foor leaves the ground.

Thank you.

Answer 1

[Image Analyst]

Try this:

clc;    % Clear the command window.
fprintf('Beginning to run %s.m ...\n', mfilename);
close all;  % Close all figures (except those of imtool.)
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 18;
lineWidth = 3;
data = readmatrix('S1C1.txt');
t = data(:, 2);
force = data(:, 3);
plot(t, force, 'b-', 'LineWidth', lineWidth)
grid on;
title('Heel Strike Forces', 'FontSize', fontSize);
xlabel('Time', 'FontSize', fontSize);
ylabel('Force', 'FontSize', fontSize);
threshold = 5; % Whatever...
yline(threshold, 'Color', 'm', 'LineWidth', lineWidth)
heelDownIndexes = reshape(force > threshold, 1, []); % Threshold and make into a row vector.
startingIndexes = strfind(heelDownIndexes, [0, 1]) + 1;
endingIndexes = strfind(heelDownIndexes, [1, 0]) + 1;
% Put vertical green lines at the start of the force pulse.
xline(t(startingIndexes), 'Color', 'g', 'LineWidth', lineWidth)
% Put vertical red lines at the end of the force pulse.
xline(t(endingIndexes), 'Color', 'r', 'LineWidth', lineWidth)
g = gcf;
g.WindowState = 'maximized'

[William Rose]

@Perri Johnson, @Image Analyst, I never knew about the xline() (or yline()) command. Thanks!

[Perri Johnson]

@Image Analyst, Thank you for the help. I just got settled and about to work on my code some more from the recommendations you and @William Rose have offered.

Answer 2

[William Rose]

@Perri Johnson

FYI here's a paper on determining ground contact time.

https://pubmed.ncbi.nlm.nih.gov/20478742/

Will get back to you on the other part.

Answer 3

[William Rose]

@Perri Johnson

SInce you did not provide data I will provide some

data=load('grf_z_raw.txt');
t=data(:,1); y=data(:,2);
plot(t,y,'-r'); xlabel('Time (s)'); ylabel('GRF_z (N)');

It's upside down so lets flip it.

y=-y; plot(t,y,'-r'); xlabel('Time (s)'); ylabel('GRF_z (N)');

OK that's better.

tHS=min(t(y>20)); fprintf('Heel strike time=%.3f\n',tHS);

Heel strike time=0.434

OK

[Perri Johnson]

William,

Thank you, this really helps. I'm going to work on adapting this more to identify toe off as well and for more than one time when it occurs since I have multiple instances in one data set.

Answer 4

[Image Analyst]

Try this to find the first time the heel gets contact:

times = data(:, 1); % Times are in column 1
force = data(:, 2); % Forces are in column 2
threshold = 0.7; % Whatever it is....
% Find first time signal exceeds threshold:
indexOfFirstHeelContact = find(force > threshold, 1, 'first');
% Get the force there
force1 = force(indexOfFirstHeelContact)
% Get the time there.
time1 = times(indexOfFirstHeelContact)

If you have multiple times where if crosses the threshold, you can use strfind() instead of find():

heelDownIndexes = force > threshold;
startingIndexes = strfind(heelDownIndexes, [0, 1]) + 1;

Attach your actual data if you still need help.

[Perri Johnson]

Thank you for the advice. I was using the find function at first but couldn't quite get it how I needed. Lookes like I was just missing those last two lines of code.

[Image Analyst]

So do you want to find just the first one, or all of them? I gave both solutions. If this works, could you Accept the best answer, and Vote for any other working ones? Thanks in advance. 🙂

[Perri Johnson]

So I have ten steps worth of data and trying to find the instance for each step where the force value rises above/below the threshold I need. Right now both answers work or identifying the very first instance of the first step at heel strike.

the problem with just setting a variable for when the force is above a certain threshold is that it takes every instance throughout stance phse when the force is greater than the threshold. I'm just looking for the moment of heel contact. For the strfind command I'm getting an error that the "Input strings mush have one row" when trying to adjust to multiple steps.

[Image Analyst]

Try this:

heelDownIndexes = reshape(force > threshold, 1, []); % Threshold and make into a row vector.
startingIndexes = strfind(heelDownIndexes, [0, 1]) + 1;

If you need to get rid of small noise regions that are too short to be legitimate, you can use bwareaopen() to get rid of runs less than a certain number of indexes long. Again, if you need more help, attach your data.

[Image Analyst]

See my other answer for a full demo.

Answer 5

[William Rose]

@Perri Johnson,

Here is ground reaction force for 10 second sampled at 2.5 kHz. (Really one step, duplicated to make 10.)

clear;
data=load('grf_z_10step.txt');
t=data(:,1); y=data(:,2);
thresh=20;        %detection threshold, in Y-units
jh=0; jt=0;
for i=2:length(t)
    if y(i)>=thresh && y(i-1)<thresh
        jh=jh+1;
        ths(jh)=t(i);
    elseif y(i)<thresh && y(i-1)>=thresh
        jt=jt+1;
        tto(jt)=t(i);
    end
end
fprintf('Found %d heel strikes and %d toe offs.\n',jh,jt);

Found 10 heel strikes and 10 toe offs.

figure; plot(t,y,'-k',ths,thresh,'rx',tto,thresh,'bo');
xlabel('Time (s)'); ylabel('GRF_z (N)'); grid on

Looks good.

[Image Analyst]

You should really try the strfind() trick I mentioned rather than scanning the vector. It's a nice vectorized way of doing it rather than a for loop. Nice trick for your bag of tricks.

[William Rose]

@Perri Johnson, @Image Analyst is right that their approach is a lot more elegant than my simple approach. Which is why I gave it a thumbs up. When @Image Analyst says stuff, I pay attention, because they are really smart and really good at Matlab. The only virtue of my approach is that it is easy to understand. I like the plot better with red and green triangles for heelstrike and toe off: replace the figure; plot(); command above with

figure; plot(t,y,'-k'); hold on
plot(ths,thresh,'r^','MarkerFaceColor','r');
plot(tto,thresh,'gv','MarkerFaceColor','g');

For reasons I don't understand, the legend command keeps putting up a legend with the wrong toe-off symbol, which is why have not been including a legend on the plots above. For example if I include

legend('GRF','HS','TO');

I get the following, which is not right.

So I'm leaving the legend off for now. Good luck.

[Perri Johnson]

@William Rose, @Image Analyst, You've both given good points and ways to perform the following operation in Matlab. There's always more than one way to perform an operation. With being new to Matlab, it's still tough to wrap my mind around how to perform certain tasks and both of these were good learning points for me. Really grateful for you both taking the time out of your schedule to teach me.

[Image Analyst]

@Perri Johnson, you're quite welcome, and thanks for accepting and voting for the answers to award reputation points. If you have any more questions, feel free to ask. Finding and defining peaks and valleys can be very tricky as the signal gets noisier, and it's a whole order of magnitude more difficult if you're dealing with 2-D images rather than 1-D signals.