Your function loops over some values and for each one computes a row vector.
To store the row vector in a matrix, specify an integer row. I've also adjusted your function to return the value:
q = integrateQuaternions;
size(q)
function Qfinal = integrateQuaternions(BR,EA)
% Parameters
A = 0; % Initial Conditions (angles halved and in rad)
B = 0;
C = -45 * (pi/180);
p = 0 * (pi/180); % Body rates in radians (0)
q = 5.5 * (pi/180); % (5.N8 deg/s)
r = 3.5 *(pi/180); % ((1+0.5*N9) deg/s)
BR = [p q r];
EA = [0 0 (-45*(pi/180))];
QFinal = nan(length(1:.1:30),4);
row = 1;
for i = 1:0.1:30 % for 30 second period with 0.1s timestep
% Initial Attitude
q0 = i.*(cosd(C)*cosd(B)*cosd(A) + sind(C)*sind(B)*sind(A));
q1 = i.*(cosd(C)*cosd(B)*sind(A) - sind(C)*sind(B)*cosd(A));
q2 = i.*(cosd(C)*sind(B)*cosd(A) + sind(C)*cosd(B)*sind(A));
q3 = i.*(-cosd(C)*sind(B)*sind(A) + sind(C)*cosd(B)*cosd(A));
Q = [q0 q1 q2 q3];
%Euler Interation
q0dot = i.*(-0.5*((q1*p)+(q2*q)+(q3*r)));
q1dot = i.*(0.5*((q0*p)-(q3*q)+(q2*r)));
q2dot = i.*(0.5*((q3*p)+(q0*q)-(q1*r)));
q3dot = i.*(-0.5*((q2*p)-(q1*q)-(q0*r)));
% Normalizing
Qdot = [q0dot q1dot q2dot q3dot];
mu = sqrt(q0^2 + q1^2 + q2^2 + q3^2);
qnplus1 = Q + Qdot;
Qx = (qnplus1./mu);
Qfinal(row,:) = Qx;
row = row + 1;
end
end
Note that a more common MATLAB approach looks more like
% timesteps = 1:.1:30;
% for i = 1:length(timesteps)
% t = timesteps(i);
% q0 = t .* ...
% ...
% Qfinal(i,:) = ...
% end
