How to find zero on plot?

21 次查看(过去 30 天)
Sam
Sam 2014-12-27
评论: Star Strider 2014-12-28
In attachment, the blue line represents the right heel marker walking up the stairs. The orange line is the derivative of the position of the heel, so it represents the speed of the heel marker. I need to know where the heel touches the ground, so when the speed of the heel is zero. How can I find these points?
  2 个评论
Azzi Abdelmalek
Azzi Abdelmalek 2014-12-27
On the plot, or on your data?
Sam
Sam 2014-12-27
data (thought it would be the same... sorry)

请先登录,再进行评论。

回答(3 个)

Image Analyst
Image Analyst 2014-12-27
编辑:Image Analyst 2014-12-27
Try this
zeroSpeedIndexes = orangeSignal = 0; % logical array of where speed is zero.
To be more robust, you might want to check if it's less than some small number (in case it never hits exactly zero), and then find the centroid, or first index, of those indexes.
zeroSpeedIndexes = abs(orangeSignal) <= someSmallNumber;
For example, let's say the speed is around 0 to 10 when the heel is moving, but sometimes the speed never quite hits zero but instead bounces around 0.05 or whatever. So someSmallNumber would be 0.05. This is more robust than just checking for an exact zero.
To get actual indexes ( in case you need them but probably don't) use find():
linearIndexes = find(zeroSpeedIndexes);
  3 个评论
Image Analyst
Image Analyst 2014-12-27
Sam, why do you make us work blind? See this and this.
I don't know what "define" means to you. I gave you code to identify elements where the speed falls below some low level or is exactly zero (if you're lucky enough to have that with actual real-world data).
Sam
Sam 2014-12-27
I have 5 subjects which are measured 5 times. I need to find where the Right Heel touches the ground, so I need to find where the speed of the heel is zero. I can find the points where it crosses the zero-line, but how can I select for all subject the correct point... Because sometimes I have 6 zero-points... (I'd rather don't want to use ginput command...)
for welke_pp=1:5 % for 5 subjects
for i_testen=1:5 %for 5 measurements per subject
RHEE = data_stair_rise(welke_pp,i_testen).VideoSignals(:, strcmp('RASI', data_stair_rise(welke_pp,i_testen).VideoSignals_headers)); %extract position data
RHEE = abs(RHEE);
T = 0:size(RHEE)-1;
dRHEE = gradient(RHEE,T); %calculate first derivative of data ==> speed
dRHEE = dRHEE*100;
zeroSpeedIndexes_down = dRHEE <= 0; %points below 0
zeroSpeedIndexes_up = dRHEE >=0; %point above 0
zeroSpeed_down = zeros;
zeroSpeed_up = zeros;
clear zeroSpeed_down;
clear zeroSpeed_up;
zeroSpeed_down(1,:) = zeroSpeedIndexes_down;
zeroSpeed_up(1,:) = zeroSpeedIndexes_up;
V1 = [0 diff(zeroSpeed_down)>0] ; %find where it cross the zero-line
welk1 = find(V1); %find where it crosses the zero-line
V2 = [0 diff(zeroSpeed_up)>0] ; %find where it crosses the zero-line
welk2 = find(V2); %find where it crosses the zero-line
end
end

请先登录,再进行评论。


Azzi Abdelmalek
Azzi Abdelmalek 2014-12-27
x=[1 2 0 3 4 5 3 2 0 4 5]
out=find(x==0)

Star Strider
Star Strider 2014-12-27
Going back to your Question earlier:
d = matfile('Sam_RHEE.mat');
RH = d.RHEE; % Signal Vector
T = 0:size(RH)-1; % Time Vector
dRH = gradient(RH,T); % Derivative (w.r.t. T, ‘T’ not required)
dRHcs = dRH .* circshift(dRH, [-1 0]);
dRHzx = find(dRHcs < 0);
for k1 = 1:size(dRHzx,1)-1
dRH0(k1) = interp1(dRH(dRHzx(k1):dRHzx(k1)+1), T(dRHzx(k1):dRHzx(k1)+1), 0);
end
figure(1)
plot(T, RH)
hold on
plot(T, dRH*100)
yl = ylim;
plot([dRH0; dRH0], [ones(1,size(dRH0,2))*yl(1); ones(1,size(dRH0,2))*yl(2)], '-g')
hold off
grid
This finds the points where the derivative ‘dRH’ crosses zero (which is what I believe you want) in the ‘dHR0’ vector, and then plots vertical green lines at those points, so you can match the derivative and the ‘RHEE’ signals at those times.
The plot:
  4 个评论
Sam
Sam 2014-12-28
Thanks! It helped a lot. But I still have a few short questions. Can you please explain the following parts of your code? (for example what circshift does, etc.)
dRHcs = dRH .* circshift(dRH, [-1 0]);
dRHzx = find(dRHcs < 0);
for k1 = 1:size(dRHzx,1)-1
dRH0(k1) = interp1(dRH(dRHzx(k1):dRHzx(k1)+1), T(dRHzx(k1):dRHzx(k1)+1), 0);
end
And this part:
plot([dRH0; dRH0], [ones(1,size(dRH0,2))*yl(1); ones(1,size(dRH0,2))*yl(2)], '-g')
Star Strider
Star Strider 2014-12-28
My pleasure!
Line-by-line:
1. dRHcs = dRH .* circshift(dRH, [-1 0]);
This takes ‘dRH’ and multiplies it by a one-position-circularly-shifted version of itself. (See the documentation for circshift for details.) All the products that have the same sign, whether positive or negative, will have positive products. If adjacent elements have opposite signs, that product will be negative.
2. dRHzx = find(dRHcs <= 0);
This provides the indices of the negative or zero elements of the multiplication in 1.
3. The for loop calculates the exact zero-crossings using adjacent elements of ‘dRHzx’ and the interp1 function. Note that the arguments to it are reversed from the usual, because the desired output is the value where ‘dRH’ is zero. It has to do this in a loop because it takes each estimated zero-crossing in turn and calculates the exact values for it. (It does not include the last estimated zero-crossing, because using circshift can create an artefactual zero-crossing at the end of the vector. It is best not to include it.)
4. The plot call looks more complicated than it is. The plot function requires at least two independent and two dependent variable values to plot a line. The usual way to plot a series of vertical lines (and the way I did it here), is to create a (2xN) matrix of different x-coordinates and a matching (2xN) matrix of essentially the same y-coordinates, so that here a vertical (green) line is plotted at every zero-crossing for the ‘dRH’ variable. If you want to see how it works in detail, create these:
xmat = [dRH0; dRH0];
ymat = [ones(1,size(dRH0,2))*yl(1); ones(1,size(dRH0,2))*yl(2)];
after the plot (because plot creates the ‘yl’ vector), and look at them. (The same idea — with the appropriate changes — works if you want to plot horizontal lines.)

请先登录,再进行评论。

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by