clear all

global B k m
B = 6000;       %Ns/m
k = 86064.685*9.81;  %N/m
m = 3084.43;    %kg

v_initial = 1.4; % m/s
pstart = m*v_initial;

start = [pstart; 0];  % start = [p_initial; x_initial]

t0 = 0;
tf = 1;
N = 1000;

[t,output] = rk4fixed('carcrash',[t0 tf],start,N);

%Discovering Important Data Points
fprintf('RESULTS for %f m/s, and %f damping\n', v_initial, B);
Maximum_Crush = max(output(:,2))
Initial_Velocity = v_initial

i=50;
while 0 < output(i,2)
    i=i+1;
%     if 0.001 > output(i,2)
%         End_Impact = t(i)
%         Exit_Velocity = velocity_usingP(i)
%     end
    End_Impact = t(i);
end
%END

%Plotting
subplot(4,1,1);
plot(t(1:i),output(1:i,1));
title('Momentum vs. Time'); xlabel('Time (sec)'); ylabel('Momentum');

subplot(4,1,2);
plot(t(1:i),output(1:i,2));
title('Position vs. Time'); xlabel('Time (sec)'); ylabel('Position');

% i=2;
% velocity(1) = v_initial;
% while i<=length(output)
%     velocity(i) = (output(i,2)-output(i-1,2))/t(i);
%     i=i+1;
% end
% 
% velocity
% i
% 
% subplot(4,1,3);
% plot(t,velocity);
% title('Velocity vs. Time'); xlabel('Time (sec)'); ylabel('Velocity');

velocity_usingP = output(:,1)/m;
subplot(4,1,3);
plot(t(1:i),velocity_usingP(1:i));
title('Velocity vs. Time'); xlabel('Time (sec)'); ylabel('Velocity');

ii=2;
acceleration(1) = 0;
while ii<=length(output)
    acceleration(ii) = (velocity_usingP(ii)-velocity_usingP(ii-1))/(t(ii)-t(ii-1));
    ii=ii+1;
end

% acceleration
% i

subplot(4,1,4);
plot(t(1:i),acceleration(1:i));
title('Acceleration vs. Time'); xlabel('Time (sec)'); ylabel('Acceleration');

Exit_Velocity = velocity_usingP(i)
Energy_Dissipated = 0.5*m*(Initial_Velocity)^2 - 0.5*m*(Exit_Velocity)^2