%% Exercise 3 Throwing a ball
% file: Exercise3_Balltrajectory.m HA 4.12.2018
%% (a)

clear; close all;format compact
g=9.81;
v=10; 
al=35; % Use sind,cosd, or al=(pi/180)*al and sin, cos

%
%% Find t-value, where ball hits the ground

Tend=10;   % Wild guess
t=linspace(0,Tend,1000);
h=v*t*sind(al)-1/2*g*t.^2; % sind(x), x degrees
h(end)<0  % If yes, then Tend is sufficiently large.

% This is the numerical way, solving the zeros by hand is even easier, see below. (One student did this way :-))
%

negind=find(h<0,1) % Find first meeting condition (h<0)
% Could use logical indexing also:
% lastposind=sum(h>=0)  % sum([1 1 1 1,...,1,0 0 0 ...]
%
Tend=t(negind)             % Update Tend
t=linspace(0,Tend,1000);   % Update t
h=v*t*sind(al)-1/2*g*t.^2; % Update h
negind=find(h<0,1)         % Update negind
% Check
hsignchange=h([negind-1 negind])
% Yes, indeed 
maxerr=t(2)-t(1)  % 
sameas=Tend/1000         % Same as this.

%

plot(t,h,t,h,'.');grid on
xlabel('t');ylabel('h');shg
%
% Using zoom suffciently many times, one gets: h=0, for t=1.16937
% Note: plot uses linear interpolation, which gives very good 
% approximation on small intervals, so the accuracy obtained by 
% zooming is much better than the distance of t-points suggests.
%

%
%% Solve the 2^{nd} degree equation numerically
% Remember Matlab's representation of polynomials.

coeff=[-g/2 v*sind(al) 0]
format long
hzero=roots(coeff)
format short
%{
coeff =
   -4.9050    5.7358         0
hzero =
                   0
   1.169370920185619
%}
% Hence all the 6 digits of the graphical method were correct.

%
%% The easiest way of course is to write the solution by hand:
% h=0 \iff v\sin(\alpha) = g t/2,
% hence:
% t_0=2v\sin(\alpha)/g

format long
Tend=2*v*sind(al)/g
format short   % back t default
% Showing that all digits given by "roots" are correct.
% Everything works fine!!
% One could also practice the symbolic toolbox, but let this be enough.

%
%% Max height:

t=linspace(0,Tend,1000);
h=v*t*sind(al)-1/2*g*t.^2;
[maxh,maxind]=max(h)  % Could use find also, but this is easiest,I guess.
max_t=t(maxind)
x=v*cosd(al)*t;
max_x=x(maxind)    % x-value for max(h) (perhaps was not asked)

%
%% (b)
% Same values as in a) Plot the ball's trajectory, i.e. the curve (x(t),h(t)).

figure(2);clf
t=linspace(0,Tend); % Enough to use 100 points.
x=v*cosd(al)*t;
h=v*sind(al)*t-1/2*g*t.^2;
plot(x,h);grid on
title('b) Trajectory of ball')
xlabel('x');ylabel('h');shg

%
%% (c) Trajectories for various angles(al)
% Same v and various \alpha

figure(3);clf
hold on
colors=['ygrbk'];  % In case you want to control the colors
                   % "yellow,green,red,blue,black"
c=0;
for al=[20 30 45 60 70]
    coeff=[-g/2 v*sind(al) 0];
    Tend=max(roots(coeff));
    t=linspace(0,Tend);
    x=v*cosd(al)*t;
    h=v*sind(al)*t-1/2*g*t.^2;
    c=c+1;
    plot(x,h,colors(c))
end   
title('Trajectories for various angles')
legend('\alpha=20','\alpha=30','\alpha=45','\alpha=60','\alpha=70')
grid on; xlabel('x');ylabel('h');shg

%
%% (d) Trajectories for \alpha=45, with various initial velocities
% Same v and various al

figure(4);clf
g=9.81;
al=45;
hold on
colors=['mgrbk'];  % "magenta,green,red,blue,black"
c=0;
for v=10:2:18
    coeff=[-g/2 v*sind(al) 0];
    Tmax=max(roots(coeff));
    t=linspace(0,Tmax);
    x=v*cosd(al)*t;
    h=v*sind(al)*t-1/2*g*t.^2;
    c=c+1;
    plot(x,h,colors(c))
end   
title('Trajectories for \alpha=45, various init. velocities')
legend('v=10','v=12','v=14','v=16','v=18')
grid on;xlabel('x');ylabel('h');shg

%