Time normalize kinematic data for gait analyses

Hello, I am having a bit of a problem with implementing my code. I have a set of data for determining the angles of the hip and the knee. its a lot of data and i need to time normalise it in order to graph it and analize. I have tried doing it 2 ways but my angle graph doesnt look right.

I would really apreciate if someody could check my code and point at any mistakes im making and suggest a better way of coding, any tips, any help.

Here is code n1

clear all
close all
[filename,folder]= uigetfile('*.xlsx');
cd(folder);
%f=1;
t=0:0.01:36.13;
%y=sin(2*pi*f*t);

Marker_trajectories=readmatrix(filename,'Sheet','BE302 Walk trajectories','Range','B6:AW3619');
A=readmatrix(filename,'sheet','anthropometrics','range','A1:A8');

p1 = Marker_trajectories(:,46:48);
p2 = Marker_trajectories(:,43:45);
p3 = Marker_trajectories(:,40:42);
p4= Marker_trajectories(:,37:39);
p5= Marker_trajectories(:,34:36);
p6= Marker_trajectories(:,31:33);
p8= Marker_trajectories(:,28:30);
p9= Marker_trajectories(:,25:27);
p10= Marker_trajectories(:,22:24);
p11= Marker_trajectories(:,19:21);
p12 = Marker_trajectories(:,16:18);
p13 = Marker_trajectories(:,13:15);
p7a = Marker_trajectories(:,10:12);
p14a = Marker_trajectories(:,7:9);
p7b = Marker_trajectories(:,4:6);
p14b = Marker_trajectories(:,1:3);
p7=(p7a+p7b)/2;
p14=(p14a+p14b)/2;
p15=(p14a+p7b)/2;
y=p2(:,2)-p15(:,1);

u_foot = (p1-p2)./vecnorm(p1-p2);
w_foot = cross(p1-p3,p2-p3)./vecnorm(cross(p1-p3,p2-p3));
v_foot= w_foot.*u_foot;

v_calf = (p3-p5)./vecnorm(p3-p5);
u_calf = cross(p4-p5,p3-p5)./vecnorm(cross(p4-p5,p3-p5));
w_calf= u_calf.*v_calf;
v_pelvis = (p14-p7)./vecnorm(p14-p7);
w_pelvis = cross(p7-p15,p14-p15)./vecnorm(cross(p7-p15,p14-p15));
u_pelvis= v_pelvis.*w_pelvis;

p_ankle= p3+ 0.0164*A(5)*u_foot+0.392*A(5)*v_foot+0.478*A(7)*w_pelvis;
p_knee=p5+0.5*A(4)*w_calf;
p_hip=p15+0.598*A(1)*u_pelvis-0.29*A(1)*w_pelvis;
i_pelvis=-v_pelvis;
j_pelvis=u_pelvis;
k_pelvis=w_pelvis;

k_thigh = (p_hip-p_knee)./vecnorm(p_hip-p_knee);
j_thigh = cross(p6-p_hip,p_knee-p_hip)./vecnorm(cross(p6-p_hip,p_knee-p_hip));
i_thigh=j_thigh.*k_thigh;
k_calf = (p_knee-p_ankle)./vecnorm(p_knee-p_ankle);
j_calf = cross(p5-p_knee,p_ankle-p_knee)./vecnorm(cross(p5-p_knee,p_ankle-p_knee));
i_calf=j_calf.*k_calf;

a_hip=-atan(norm(dot(k_thigh,j_pelvis,2))./(dot(k_thigh,k_pelvis,2)))*180/pi;
b_hip=asin(dot(k_thigh,i_pelvis,2));
y_hip=-atan(dot(j_thigh,i_pelvis,2)./dot(i_thigh,i_pelvis,2));
a_knee=-atan(dot(k_calf,j_thigh,2)./dot(k_calf,k_thigh,2));
b_knee=asin(dot(k_calf,i_thigh,2));
y_knee=-atan(dot(j_calf,i_thigh,2)./dot(i_calf,i_thigh,2));
[peaks,i_peaks]=findpeaks((-y),'MinPeakDistance',100);

%1. % points of the gait cycle
gc_perc=0:1:100;
% 3.i_peaks- heel strikes. number of full gait cycles i_peaks-1
%heel_strikes=size(i_peaks);
%gait_cycle=heel_strikes-1;
for i=1:size(i_peaks,2)-1
    i_gc=i_peaks(i):i_peaks(i+1);
    %  4b) number of rows
 
   num_rows=size(i_gc,2);
    gc=0:1/num_rows:100;
 
       gait_cycle=y(1,i_gc);
end       

  for j=1:size(i_peaks,2) -1 
 new_angle(:,ci) =interp1(a_hip,gc_perc,gc);
  end    
 
   figure(1);

plot(new_angle,gc_perc);

I used findpeaks to seperate the gait cycles and there are 33 cycles in my data. to separate the gait cycles we use a formula provided

y=p2(:,2)-p15(:,1);

I am not sure how to apply the same principle to the rest of my data, i need to find the angles of hip and knee and i need to separate it into the gait cycles. But im getting an error Error in BE302_kinematics (line 86)

plot(new_angle,gc_perc); the program doesnt recognise the new_angle variable. I am not able to find what im doing wrong. So a wrote another code below:

F=a_hip;
x = linspace(0, 99, length(y));
z = 0:99;
y_z_normed = spline(x, y, z);
d=linspace(0,99,length(F));
a_hip_n=spline(d,F,z);

figure
subplot(2,1,1)
plot(y, 'b', 'LineWidth', 2)
xlim([0 length(y)])
title('Original Signal')
xlabel('Frames')
ylabel('a Hip')
subplot(2,1,2)
plot(y_z_normed, 'b', 'LineWidth', 2)
xlim([0 100])
title('Normalized Signal')
xlabel('Normalized Time (0-100%)')
ylabel('a Hip')
figure(3);
subplot(2,1,1)
plot(a_knee, 'b', 'LineWidth', 2)
xlim([0 length(a_knee)])
title('Original Signal')
xlabel('Frames')
ylabel('a Knee')
subplot(2,1,2)
plot(a_hip_n, 'b', 'LineWidth', 2)
xlim([0 100])
title('Normalized Signal')
xlabel('Normalized Time (0-100%)')
ylabel('a Knee')

with this i am getting a graph but it doesnt look right to me, I also need all of the gait cycles to be graphed together from 0 to 100% on top of each other. I tried doing it this way:

[pks, locs] = findpeaks(y, 'MinPeakDist',50);     % Determine Peaks & Indices
figure(1)
plot(a_hip, t)
hold on
plot(a_hip(locs), pks, '+r')
hold off
grid
for k1 = 1:numel(locs)-1
    apc{k1} = a_hip(locs(k1)-10:locs(k1+1)-10);                            % Define MUAP Frames
    tvc{k1} = t(locs(k1)-10:locs(k1+1)-10);
end
figure(2)
hold all
for k1 = 1:numel(apc)
    plot(tvc{k1}-tvc{k1}(1), apc{k1})                                   % Plot MUAP Frames
end
hold off
grid

I got this graph which is not right I also plotted it againced time not the gait cycle. My main issue i dont know how to tie the gait cycles i found from findpeaks to the rest of the data for finding the angles of the hip and the knee. I need to time normalise it for each gait cycle and graph the hip and knee angles for the average gait cycle.

Answer 1

[Star Strider]

Note that ‘i_peaks’ is a (33x1) column vector, That means that:

for j=1:size(i_peaks,2) -1

new_angle(:,ci) =interp1(a_hip,gc_perc,gc);

end

iterates from 1 to 0, with an implied increment of +1, so essentially does not iterate at all. That results in ‘new_angle’ being empty, and the error being thrown.

The next problem (after replacing the 2 with a 1, or using numel instead) is that ‘gc’ is not defined. A potential problem is that the interp1 function interpolates the first argument corresponding to the third argument, returning the interpolated value of the second argument. So if ‘gc’ is supposed to interpolate ‘gc_perc’, and return a corresponding value for ‘a_hip’, the first two arguments need to be reversed from their current placement.

.

[Neil Mcdonald]

Thank you for reply. I have fixed that but its still not working. Im getting this error Unrecognized function or variable 'gc'.

Error in BE302_kinematics (line 82)

new_angle(:,j) =interp1(gc,a_hip,gc_perc);

 for j=1:size(i_peaks,1) -1 
 new_angle(:,j) =interp1(gc,a_hip,gc_perc);
  end  
   figure(1);
   

[Neil Mcdonald]

I have changed my code and now I get gc, i_gc of the same size 1x114 but the size of a_hip is 3614x1. Can you please advise if interpr1 function will work with different lentgh arrays? Because right now im getting an error Error using interp1>reshapeAndSortXandV (line 445)

X and V must be of the same length.

Error in interp1 (line 128)

[X,V,orig_size_v] = reshapeAndSortXandV(X,V);

Error in BE302_kinematics (line 82)

new_angle(:,j) =interp1(gc,gc_perc,a_hip,"spline");

for i=1:size(i_peaks,1)-1
    i_gc=i_peaks(i):i_peaks(i+1);
    %  4b) number of rows
 
   num_rows=size(i_gc,2)-1;
    gc=0:100/num_rows:100;
 
       %gait_cycle=y(1,i_gc);
end       

  for j=1:size(i_peaks,1) -1 
 new_angle(:,j) =interp1(gc,a_hip,gc_perc,"spline");
  end  

[Star Strider]

I noted that in my Answer. You have to decide what ‘gc’ is supposed to be, however I believe it still needs to be the third argument. Specifically, since the other arguments are defined, if I understand correctly what you want to do, I believe it should be:

new_angle(:,j) = interp1(gc_perc,a_hip,,gc);

That will return a value for ‘a_hip’ corresponding to the value of ‘gc_perc’ interpolated by ‘gc’ when you define ‘gc’ and since I do not know what it refers to, I leave that to you.

I just now happened to notice some atan calls. Consider using atan2 instead. It’s more robust.

.

[Star Strider]

‘Can you please advise if interpr1 function will work with different lentgh arrays?’

The first argument must be a vector. The second argument can be either a vector or a matrix, however it must have a corresponding dimension that matches the first argument, and the tthird argument must be a scalar or a vector.

X = linspace(0,10,11).'

X = 11×1

0 1 2 3 4 5 6 7 8 9

Y = randn(11,3)

Y = 11×3

0.3324 -1.2041 0.0289 -1.4620 0.5125 0.7410 1.2912 0.6179 -0.3363 0.5872 0.9363 0.7856 -0.3178 1.0000 1.5058 -0.8859 -1.0873 1.0522 1.2948 -0.2098 0.1146 -1.8857 -0.1262 1.4377 0.4107 0.8968 -0.4403 0.7492 0.6203 -1.1425

q = pi*[1 2]

q = 1×2

3.1416 6.2832

Out = interp1(X, Y, q)

Out = 2×3

0.4591 0.9453 0.8876 0.3941 -0.1862 0.4893

Also, please see my earlier Comment.

.

[Neil Mcdonald]

gc is the variable to create a % cycle based on number of rows.

it supposed to be the same length as a_hip but its not. I need to identify the rows in the hip data which correspond to the gait cycle im interested in, put it into a variable and that should be the same length as gc,

Can you advise what can be done? im struggling with implementing this

[Star Strider]

I am now lost. I have no idea what you are interpolating with respect to.

[Neil Mcdonald]

im interpolating onto gc_perc

[Star Strider]

O.K. To clarify, what are the vectors you have, and what are the results you want?

I’m still not clear on that.

[Neil Mcdonald]

I need to time normalise data for each gate cycle from 0 to 100% for hip and ankle angles.

[Neil Mcdonald]

gc is the vector based on the number of rows gc=0:100/num_rows:100, gc_perc -defines the percentage pointes of the gait cycle 0:1:100. a_hip is the angle for the hip. I need to time normalise the hip angles and plot it

[Star Strider]

If I understand correctly what you want to do, the rescale function could be the most efficient approach. I am not certain what the percentage is to be applied to, so it may be necessary to change the argument order and definition. This assumes that the angles are to be normalised between 0% and 100%.

Perhaps something like this, then, normalising on time —

t = linspace(0, 15, 50);                                                        % Time
a_hip = 25*exp(-(t-7).^2/10);                                                   % Angle
tr = rescale(t, 0, 100, 'InputMin',min(t), 'InputMax',max(t));     % Normalised Angle

figure
tiledlayout(2,1)
nexttile
plot(t, a_hip)
xlabel('Time')
ylabel('Angle')
nexttile
plot(tr, a_hip)
xlabel('Time')
ylabel('Normalised Angle')

EDIT — (2 Apr 2023 at 21:07)

Normalised on ‘t’.

.

[Neil Mcdonald]

I have 33 gait cycles, so they all have a lot of spikes and they need to be normalised in order to be able to analyze it. I need to separate each cycle and identify it by percentage points. Then use interp1 function to interpolate onto percentage, so when 0 cycle starts, when it reaches 100 cycle ends and so on.

[Neil Mcdonald]

yes from 0 to 100 % is related to time, 0 is initial contact with the ground and 100% foot off. 1 gait cycle contains 101 rows of data and I have 33 columns of such data. if i plot it raw it will all look messy, so we time normalising it so you can analize the mean angle and make sense of it. So the angles get normalised in resect of time.

[Star Strider]

So you are normalising on time, so that all the times are between 0 and 100. I was confused with respect to the interp1 call, since it was not obvious to me what the variables were, especially since one of them did not exist when I first ran that part of the code.

Normalising on time is straightforward to do. See my edited Comment.

EDIT — (2 Apr 2023 at 21:53)

Our latest comments collided in time. Please see my edited previous Comment. You can even create an anonymous function out of it, so all that is necessary is to provide an appropriate name for the result —

rescaleTime = @(t) rescale(t, 0, 100, 'InputMin',min(t), 'InputMax',max(t));    % Normalised Time
t1 = 0:10

t1 = 1×11

0 1 2 3 4 5 6 7 8 9 10

t2 = 0:12

t2 = 1×13

0 1 2 3 4 5 6 7 8 9 10 11 12

t1r = rescaleTime(t1)

t1r = 1×11

0 10 20 30 40 50 60 70 80 90 100

t2r = rescaleTime(t2)

t2r = 1×13

0 8.3333 16.6667 25.0000 33.3333 41.6667 50.0000 58.3333 66.6667 75.0000 83.3333 91.6667 100.0000

That should do what you want, and with reasonable efficiency.

.

[Neil Mcdonald]

if i want to identify the rows in the a_hip data which correspond to the gait cycle and create a new variable can i do it this way. i_gc is my gait cycle, gc_data is my new variable of a_hip data

for i=1:size(i_peaks,1)-1
    i_gc=i_peaks(i):i_peaks(i+1);
    %  number of rows
  num_rows=size(i_gc,2)-1;
    gc=0:100/num_rows:100;
    gc_data=a_hip(i_gc);
end       

[Neil Mcdonald]

This is the graph i was able to make, the top is raw data, below a time normalised knee angle.

Im still checking my code because some of the graphs dont look quite right, but thank you for your help i really appreciate.

[Star Strider]

My pleasure!

If my Answer helped you solve your problem, please Accept it!

.

[Neil Mcdonald]

Can you please have a look at the previous message:

if i want to identify the rows in the a_hip data which correspond to the gait cycle and create a new variable can i do it this way. i_gc is my gait cycle, gc_data is my new variable of a_hip data

ThemeCopy

for i=1:size(i_peaks,1)-1

i_gc=i_peaks(i):i_peaks(i+1);

% number of rows

num_rows=size(i_gc,2)-1;

gc=0:100/num_rows:100;

gc_data=a_hip(i_gc);

end

[Star Strider]

Yes, however you need to subscript them in a cell array to save them —

for i=1:size(i_peaks,1)-1

i_gc=i_peaks(i):i_peaks(i+1);

% number of rows

num_rows=size(i_gc,2)-1;

gc{i,:}=0:100/num_rows:100;

gc_data{i,:}=a_hip(i_gc);

end

See the result below.

.

% clear all
% close all
% [filename,folder]= uigetfile('*.xlsx');
% cd(folder);
%f=1;
filename = 'https://www.mathworks.com/matlabcentral/answers/uploaded_files/1343009/BE302%20walk%20data%202023.xlsx';

t=0:0.01:36.13;
%y=sin(2*pi*f*t);

Marker_trajectories=readmatrix(filename,'Sheet','BE302 Walk trajectories','Range','B6:AW3619');
A=readmatrix(filename,'sheet','anthropometrics','range','A1:A8');

p1 = Marker_trajectories(:,46:48);
p2 = Marker_trajectories(:,43:45);
p3 = Marker_trajectories(:,40:42);
p4= Marker_trajectories(:,37:39);
p5= Marker_trajectories(:,34:36);
p6= Marker_trajectories(:,31:33);
p8= Marker_trajectories(:,28:30);
p9= Marker_trajectories(:,25:27);
p10= Marker_trajectories(:,22:24);
p11= Marker_trajectories(:,19:21);
p12 = Marker_trajectories(:,16:18);
p13 = Marker_trajectories(:,13:15);
p7a = Marker_trajectories(:,10:12);
p14a = Marker_trajectories(:,7:9);
p7b = Marker_trajectories(:,4:6);
p14b = Marker_trajectories(:,1:3);
p7=(p7a+p7b)/2;
p14=(p14a+p14b)/2;
p15=(p14a+p7b)/2;
y=p2(:,2)-p15(:,1);

u_foot = (p1-p2)./vecnorm(p1-p2);
w_foot = cross(p1-p3,p2-p3)./vecnorm(cross(p1-p3,p2-p3));
v_foot= w_foot.*u_foot;

v_calf = (p3-p5)./vecnorm(p3-p5);
u_calf = cross(p4-p5,p3-p5)./vecnorm(cross(p4-p5,p3-p5));
w_calf= u_calf.*v_calf;
v_pelvis = (p14-p7)./vecnorm(p14-p7);
w_pelvis = cross(p7-p15,p14-p15)./vecnorm(cross(p7-p15,p14-p15));
u_pelvis= v_pelvis.*w_pelvis;

p_ankle= p3+ 0.0164*A(5)*u_foot+0.392*A(5)*v_foot+0.478*A(7)*w_pelvis;
p_knee=p5+0.5*A(4)*w_calf;
p_hip=p15+0.598*A(1)*u_pelvis-0.29*A(1)*w_pelvis;
i_pelvis=-v_pelvis;
j_pelvis=u_pelvis;
k_pelvis=w_pelvis;

k_thigh = (p_hip-p_knee)./vecnorm(p_hip-p_knee);
j_thigh = cross(p6-p_hip,p_knee-p_hip)./vecnorm(cross(p6-p_hip,p_knee-p_hip));
i_thigh=j_thigh.*k_thigh;
k_calf = (p_knee-p_ankle)./vecnorm(p_knee-p_ankle);
j_calf = cross(p5-p_knee,p_ankle-p_knee)./vecnorm(cross(p5-p_knee,p_ankle-p_knee));
i_calf=j_calf.*k_calf;

a_hip=-atan(norm(dot(k_thigh,j_pelvis,2))./(dot(k_thigh,k_pelvis,2)))*180/pi;
b_hip=asin(dot(k_thigh,i_pelvis,2));
y_hip=-atan(dot(j_thigh,i_pelvis,2)./dot(i_thigh,i_pelvis,2));
a_knee=-atan(dot(k_calf,j_thigh,2)./dot(k_calf,k_thigh,2));
b_knee=asin(dot(k_calf,i_thigh,2));
y_knee=-atan(dot(j_calf,i_thigh,2)./dot(i_calf,i_thigh,2));
[peaks,i_peaks]=findpeaks((-y),'MinPeakDistance',100);

%1. % points of the gait cycle
gc_perc=0:1:100;
% 3.i_peaks- heel strikes. number of full gait cycles i_peaks-1
%heel_strikes=size(i_peaks);
%gait_cycle=heel_strikes-1;
% for i=1:size(i_peaks,2)-1
%     i_gc=i_peaks(i):i_peaks(i+1);
%     %  4b) number of rows
% 
%    num_rows=size(i_gc,2);
%     gc=0:1/num_rows:100;
% 
%        gait_cycle=y(1,i_gc);
% end       
% 
%   for j=1:size(i_peaks,2) -1 
%  new_angle(:,ci) =interp1(a_hip,gc_perc,gc);
%   end    
% 
%    figure(1);
% 
% plot(new_angle,gc_perc);

for i=1:size(i_peaks,1)-1
    i_gc=i_peaks(i):i_peaks(i+1);
    %  number of rows
    num_rows=size(i_gc,2)-1;
    gc{i,:}=0:100/num_rows:100;
    gc_data{i,:}=a_hip(i_gc);
end       

gc

gc = 32×1 cell array
    {[ 0 0.9009 1.8018 2.7027 3.6036 4.5045 5.4054 6.3063 7.2072 8.1081 9.0090 9.9099 10.8108 11.7117 12.6126 13.5135 14.4144 15.3153 16.2162 17.1171 18.0180 18.9189 19.8198 20.7207 21.6216 … ]}
    {[0 0.9174 1.8349 2.7523 3.6697 4.5872 5.5046 6.4220 7.3394 8.2569 9.1743 10.0917 11.0092 11.9266 12.8440 13.7615 14.6789 15.5963 16.5138 17.4312 18.3486 19.2661 20.1835 21.1009 22.0183 … ]}
    {[0 0.9346 1.8692 2.8037 3.7383 4.6729 5.6075 6.5421 7.4766 8.4112 9.3458 10.2804 11.2150 12.1495 13.0841 14.0187 14.9533 15.8879 16.8224 17.7570 18.6916 19.6262 20.5607 21.4953 22.4299 … ]}
    {[0 0.9346 1.8692 2.8037 3.7383 4.6729 5.6075 6.5421 7.4766 8.4112 9.3458 10.2804 11.2150 12.1495 13.0841 14.0187 14.9533 15.8879 16.8224 17.7570 18.6916 19.6262 20.5607 21.4953 22.4299 … ]}
    {[  0 0.9091 1.8182 2.7273 3.6364 4.5455 5.4545 6.3636 7.2727 8.1818 9.0909 10 10.9091 11.8182 12.7273 13.6364 14.5455 15.4545 16.3636 17.2727 18.1818 19.0909 20 20.9091 21.8182 22.7273 … ]}
    {[ 0 0.8929 1.7857 2.6786 3.5714 4.4643 5.3571 6.2500 7.1429 8.0357 8.9286 9.8214 10.7143 11.6071 12.5000 13.3929 14.2857 15.1786 16.0714 16.9643 17.8571 18.7500 19.6429 20.5357 21.4286 … ]}
    {[0 0.9174 1.8349 2.7523 3.6697 4.5872 5.5046 6.4220 7.3394 8.2569 9.1743 10.0917 11.0092 11.9266 12.8440 13.7615 14.6789 15.5963 16.5138 17.4312 18.3486 19.2661 20.1835 21.1009 22.0183 … ]}
    {[0 0.9434 1.8868 2.8302 3.7736 4.7170 5.6604 6.6038 7.5472 8.4906 9.4340 10.3774 11.3208 12.2642 13.2075 14.1509 15.0943 16.0377 16.9811 17.9245 18.8679 19.8113 20.7547 21.6981 22.6415 … ]}
    {[0 0.9615 1.9231 2.8846 3.8462 4.8077 5.7692 6.7308 7.6923 8.6538 9.6154 10.5769 11.5385 12.5000 13.4615 14.4231 15.3846 16.3462 17.3077 18.2692 19.2308 20.1923 21.1538 22.1154 23.0769 … ]}
    {[0 0.9346 1.8692 2.8037 3.7383 4.6729 5.6075 6.5421 7.4766 8.4112 9.3458 10.2804 11.2150 12.1495 13.0841 14.0187 14.9533 15.8879 16.8224 17.7570 18.6916 19.6262 20.5607 21.4953 22.4299 … ]}
    {[0 0.9259 1.8519 2.7778 3.7037 4.6296 5.5556 6.4815 7.4074 8.3333 9.2593 10.1852 11.1111 12.0370 12.9630 13.8889 14.8148 15.7407 16.6667 17.5926 18.5185 19.4444 20.3704 21.2963 22.2222 … ]}
    {[ 0 0.8929 1.7857 2.6786 3.5714 4.4643 5.3571 6.2500 7.1429 8.0357 8.9286 9.8214 10.7143 11.6071 12.5000 13.3929 14.2857 15.1786 16.0714 16.9643 17.8571 18.7500 19.6429 20.5357 21.4286 … ]}
    {[ 0 0.9009 1.8018 2.7027 3.6036 4.5045 5.4054 6.3063 7.2072 8.1081 9.0090 9.9099 10.8108 11.7117 12.6126 13.5135 14.4144 15.3153 16.2162 17.1171 18.0180 18.9189 19.8198 20.7207 21.6216 … ]}
    {[ 0 0.8929 1.7857 2.6786 3.5714 4.4643 5.3571 6.2500 7.1429 8.0357 8.9286 9.8214 10.7143 11.6071 12.5000 13.3929 14.2857 15.1786 16.0714 16.9643 17.8571 18.7500 19.6429 20.5357 21.4286 … ]}
    {[  0 0.9091 1.8182 2.7273 3.6364 4.5455 5.4545 6.3636 7.2727 8.1818 9.0909 10 10.9091 11.8182 12.7273 13.6364 14.5455 15.4545 16.3636 17.2727 18.1818 19.0909 20 20.9091 21.8182 22.7273 … ]}
    {[ 0 0.8850 1.7699 2.6549 3.5398 4.4248 5.3097 6.1947 7.0796 7.9646 8.8496 9.7345 10.6195 11.5044 12.3894 13.2743 14.1593 15.0442 15.9292 16.8142 17.6991 18.5841 19.4690 20.3540 21.2389 … ]}

gc_data

gc_data = 32×1 cell array
    {112×1 double}
    {110×1 double}
    {108×1 double}
    {108×1 double}
    {111×1 double}
    {113×1 double}
    {110×1 double}
    {107×1 double}
    {105×1 double}
    {108×1 double}
    {109×1 double}
    {113×1 double}
    {112×1 double}
    {113×1 double}
    {111×1 double}
    {114×1 double}

.

[Neil Mcdonald]

Thats giving me an error :Error using interp1>reshapeAndSortXandV (line 435)

X must be a vector of type double or single.

Error in interp1 (line 128)

[X,V,orig_size_v] = reshapeAndSortXandV(X,V);

Error in improved (line 80)

new_angle(:,j) =interp1(gc,gc_data,gc_perc,"spline");

Whast do i do?

[Star Strider]

I wasn’t certain where it was to go in the code.

It’s necessary to index the cell arrays:

new_angle(:,j) =interp1(gc{j},gc_data{j},gc_perc,"spline");

Try this —

% clear all
% close all
% [filename,folder]= uigetfile('*.xlsx');
% cd(folder);
%f=1;
filename = 'https://www.mathworks.com/matlabcentral/answers/uploaded_files/1343009/BE302%20walk%20data%202023.xlsx';

t=0:0.01:36.13;
%y=sin(2*pi*f*t);

Marker_trajectories=readmatrix(filename,'Sheet','BE302 Walk trajectories','Range','B6:AW3619');
A=readmatrix(filename,'sheet','anthropometrics','range','A1:A8');

p1 = Marker_trajectories(:,46:48);
p2 = Marker_trajectories(:,43:45);
p3 = Marker_trajectories(:,40:42);
p4= Marker_trajectories(:,37:39);
p5= Marker_trajectories(:,34:36);
p6= Marker_trajectories(:,31:33);
p8= Marker_trajectories(:,28:30);
p9= Marker_trajectories(:,25:27);
p10= Marker_trajectories(:,22:24);
p11= Marker_trajectories(:,19:21);
p12 = Marker_trajectories(:,16:18);
p13 = Marker_trajectories(:,13:15);
p7a = Marker_trajectories(:,10:12);
p14a = Marker_trajectories(:,7:9);
p7b = Marker_trajectories(:,4:6);
p14b = Marker_trajectories(:,1:3);
p7=(p7a+p7b)/2;
p14=(p14a+p14b)/2;
p15=(p14a+p7b)/2;
y=p2(:,2)-p15(:,1);

u_foot = (p1-p2)./vecnorm(p1-p2);
w_foot = cross(p1-p3,p2-p3)./vecnorm(cross(p1-p3,p2-p3));
v_foot= w_foot.*u_foot;

v_calf = (p3-p5)./vecnorm(p3-p5);
u_calf = cross(p4-p5,p3-p5)./vecnorm(cross(p4-p5,p3-p5));
w_calf= u_calf.*v_calf;
v_pelvis = (p14-p7)./vecnorm(p14-p7);
w_pelvis = cross(p7-p15,p14-p15)./vecnorm(cross(p7-p15,p14-p15));
u_pelvis= v_pelvis.*w_pelvis;

p_ankle= p3+ 0.0164*A(5)*u_foot+0.392*A(5)*v_foot+0.478*A(7)*w_pelvis;
p_knee=p5+0.5*A(4)*w_calf;
p_hip=p15+0.598*A(1)*u_pelvis-0.29*A(1)*w_pelvis;
i_pelvis=-v_pelvis;
j_pelvis=u_pelvis;
k_pelvis=w_pelvis;

k_thigh = (p_hip-p_knee)./vecnorm(p_hip-p_knee);
j_thigh = cross(p6-p_hip,p_knee-p_hip)./vecnorm(cross(p6-p_hip,p_knee-p_hip));
i_thigh=j_thigh.*k_thigh;
k_calf = (p_knee-p_ankle)./vecnorm(p_knee-p_ankle);
j_calf = cross(p5-p_knee,p_ankle-p_knee)./vecnorm(cross(p5-p_knee,p_ankle-p_knee));
i_calf=j_calf.*k_calf;

a_hip=-atan(norm(dot(k_thigh,j_pelvis,2))./(dot(k_thigh,k_pelvis,2)))*180/pi;
b_hip=asin(dot(k_thigh,i_pelvis,2));
y_hip=-atan(dot(j_thigh,i_pelvis,2)./dot(i_thigh,i_pelvis,2));
a_knee=-atan(dot(k_calf,j_thigh,2)./dot(k_calf,k_thigh,2));
b_knee=asin(dot(k_calf,i_thigh,2));
y_knee=-atan(dot(j_calf,i_thigh,2)./dot(i_calf,i_thigh,2));
[peaks,i_peaks]=findpeaks((-y),'MinPeakDistance',100);

%1. % points of the gait cycle
gc_perc=0:1:100;
% 3.i_peaks- heel strikes. number of full gait cycles i_peaks-1
%heel_strikes=size(i_peaks);
%gait_cycle=heel_strikes-1;
% for i=1:size(i_peaks,2)-1
%     i_gc=i_peaks(i):i_peaks(i+1);
%     %  4b) number of rows
% 
%    num_rows=size(i_gc,2);
%     gc=0:1/num_rows:100;
% 
%        gait_cycle=y(1,i_gc);
% end       
% 

for i=1:size(i_peaks,1)-1
    i_gc=i_peaks(i):i_peaks(i+1);
    %  number of rows
    num_rows=size(i_gc,2)-1;
    gc{i,:}=0:100/num_rows:100;
    gc_data{i,:}=a_hip(i_gc);
end       

gc

gc = 32×1 cell array
    {[ 0 0.9009 1.8018 2.7027 3.6036 4.5045 5.4054 6.3063 7.2072 8.1081 9.0090 9.9099 10.8108 11.7117 12.6126 13.5135 14.4144 15.3153 16.2162 17.1171 18.0180 18.9189 19.8198 20.7207 21.6216 … ]}
    {[0 0.9174 1.8349 2.7523 3.6697 4.5872 5.5046 6.4220 7.3394 8.2569 9.1743 10.0917 11.0092 11.9266 12.8440 13.7615 14.6789 15.5963 16.5138 17.4312 18.3486 19.2661 20.1835 21.1009 22.0183 … ]}
    {[0 0.9346 1.8692 2.8037 3.7383 4.6729 5.6075 6.5421 7.4766 8.4112 9.3458 10.2804 11.2150 12.1495 13.0841 14.0187 14.9533 15.8879 16.8224 17.7570 18.6916 19.6262 20.5607 21.4953 22.4299 … ]}
    {[0 0.9346 1.8692 2.8037 3.7383 4.6729 5.6075 6.5421 7.4766 8.4112 9.3458 10.2804 11.2150 12.1495 13.0841 14.0187 14.9533 15.8879 16.8224 17.7570 18.6916 19.6262 20.5607 21.4953 22.4299 … ]}
    {[  0 0.9091 1.8182 2.7273 3.6364 4.5455 5.4545 6.3636 7.2727 8.1818 9.0909 10 10.9091 11.8182 12.7273 13.6364 14.5455 15.4545 16.3636 17.2727 18.1818 19.0909 20 20.9091 21.8182 22.7273 … ]}
    {[ 0 0.8929 1.7857 2.6786 3.5714 4.4643 5.3571 6.2500 7.1429 8.0357 8.9286 9.8214 10.7143 11.6071 12.5000 13.3929 14.2857 15.1786 16.0714 16.9643 17.8571 18.7500 19.6429 20.5357 21.4286 … ]}
    {[0 0.9174 1.8349 2.7523 3.6697 4.5872 5.5046 6.4220 7.3394 8.2569 9.1743 10.0917 11.0092 11.9266 12.8440 13.7615 14.6789 15.5963 16.5138 17.4312 18.3486 19.2661 20.1835 21.1009 22.0183 … ]}
    {[0 0.9434 1.8868 2.8302 3.7736 4.7170 5.6604 6.6038 7.5472 8.4906 9.4340 10.3774 11.3208 12.2642 13.2075 14.1509 15.0943 16.0377 16.9811 17.9245 18.8679 19.8113 20.7547 21.6981 22.6415 … ]}
    {[0 0.9615 1.9231 2.8846 3.8462 4.8077 5.7692 6.7308 7.6923 8.6538 9.6154 10.5769 11.5385 12.5000 13.4615 14.4231 15.3846 16.3462 17.3077 18.2692 19.2308 20.1923 21.1538 22.1154 23.0769 … ]}
    {[0 0.9346 1.8692 2.8037 3.7383 4.6729 5.6075 6.5421 7.4766 8.4112 9.3458 10.2804 11.2150 12.1495 13.0841 14.0187 14.9533 15.8879 16.8224 17.7570 18.6916 19.6262 20.5607 21.4953 22.4299 … ]}
    {[0 0.9259 1.8519 2.7778 3.7037 4.6296 5.5556 6.4815 7.4074 8.3333 9.2593 10.1852 11.1111 12.0370 12.9630 13.8889 14.8148 15.7407 16.6667 17.5926 18.5185 19.4444 20.3704 21.2963 22.2222 … ]}
    {[ 0 0.8929 1.7857 2.6786 3.5714 4.4643 5.3571 6.2500 7.1429 8.0357 8.9286 9.8214 10.7143 11.6071 12.5000 13.3929 14.2857 15.1786 16.0714 16.9643 17.8571 18.7500 19.6429 20.5357 21.4286 … ]}
    {[ 0 0.9009 1.8018 2.7027 3.6036 4.5045 5.4054 6.3063 7.2072 8.1081 9.0090 9.9099 10.8108 11.7117 12.6126 13.5135 14.4144 15.3153 16.2162 17.1171 18.0180 18.9189 19.8198 20.7207 21.6216 … ]}
    {[ 0 0.8929 1.7857 2.6786 3.5714 4.4643 5.3571 6.2500 7.1429 8.0357 8.9286 9.8214 10.7143 11.6071 12.5000 13.3929 14.2857 15.1786 16.0714 16.9643 17.8571 18.7500 19.6429 20.5357 21.4286 … ]}
    {[  0 0.9091 1.8182 2.7273 3.6364 4.5455 5.4545 6.3636 7.2727 8.1818 9.0909 10 10.9091 11.8182 12.7273 13.6364 14.5455 15.4545 16.3636 17.2727 18.1818 19.0909 20 20.9091 21.8182 22.7273 … ]}
    {[ 0 0.8850 1.7699 2.6549 3.5398 4.4248 5.3097 6.1947 7.0796 7.9646 8.8496 9.7345 10.6195 11.5044 12.3894 13.2743 14.1593 15.0442 15.9292 16.8142 17.6991 18.5841 19.4690 20.3540 21.2389 … ]}

gc_data

gc_data = 32×1 cell array
    {112×1 double}
    {110×1 double}
    {108×1 double}
    {108×1 double}
    {111×1 double}
    {113×1 double}
    {110×1 double}
    {107×1 double}
    {105×1 double}
    {108×1 double}
    {109×1 double}
    {113×1 double}
    {112×1 double}
    {113×1 double}
    {111×1 double}
    {114×1 double}

for j=1:size(i_peaks,1) -1 
    new_angle(:,j) =interp1(gc{j},gc_data{j},gc_perc,"spline");
end    

   figure(1);

plot(new_angle,gc_perc);

.

[Afet]

Hey,

I have a related question to this. I already have kinematic data calculated in OpenSim. I have 5 trials per subject and comparing simultanously recorded markerless and Vicon data. I have synchronized the systems however there is still a shift in the data and I need to time normalize it. This code is super helpful but I got a bit confused on which part exactly normalizes and finds the intances of the gait cycle. In my data there is two gait cycles and I will use only the first one. Is there a more simplified version of this code? Or any advice on how I can do it?

[Star Strider]

@Afet — I believe the part that finds the instances of the gait cycle is the is the last findpeaks call and the code following it. (That seems to be adapted from code I wrote for a different problem about calculating the ensembles and ensemble average of an EMG problem, since it refers to MUAPs.) The rest is not code I wrote and my involvement here (about five months ago) was limited to solving this particular problem (the details of which I do not now remember).

If you have a specific problem with your data, it would be best to post it as a new Question, including the necessary code and the relevant data files.