Projectile Motion of a particle with drag force

ProjectilewithDrag.m

%Projectile motion
clear all
close all
clc

global g
global h
global m
g=9.81;
h=0.3;
m=2;
v0=10;
x0=0;
y0=0;
alpha=pi/3;

X0=[0 v0*cos(alpha) 0 v0*sin(alpha)];
T=0;
X=X0;
t=2;
tsamp=0.01;

runtime=t/tsamp;
tspan=[0 tsamp];
    
for i=1:runtime
    t=i*tsamp;
    [t,y]=ode45(@ProjectileMotion,tspan,X0);
    X0=y(length(y),:);
    T=[T;i*tsamp];
    X=[X;X0]; 
end
 for i=2:size(X,1)
hold off
plot(X(i,1),X(i,3),'bo','markerfacecolor','b') 
axis([0,12,-4,8])
hold on
pause(0.05)
% drawnow
end

function dx=ProjectileMotion(t,x)
global g
global h
global m
dx(1,1)=x(2);
dx(2,1)=-(h/m)*x(2)*sqrt(x(2)^2+x(4)^2);
dx(3,1)=x(4);
dx(4,1)=-g-(h/m)*x(4)*sqrt(x(2)^2+x(4)^2);
end