Solving the same equations give different values

I am solving a few 3D physics motion equations. For the first function (Landing_Physics3D), I am solving for time, xf, and yf. In the second function (Landing_Angles), I take the xf and yf found from function 1 and enter it in function 2 to solve for phi2 and theta2 (previously given in function 1). Why are my phi2 and theta2 values in function 2 different than the theta and phi values entered in function 1 since all of the constants for both are exactly the same? I have the functions saved in separate files and then the code after them in a third file (had to compile it all here).

function [xf, yf, dt] = Landing_Physics3D(theta, phi, u, h1, h2, ax, ay, az)
syms dx dy dt
u_x = u*cosd(theta)*cosd(phi);
u_y = u*sind(theta)*cosd(phi);
u_z = u*sind(phi);
dz = h2 - h1;
f1 = xf == u_x*dt + 0.5*ax*dt^2;
f2 = yf == u_y*dt + 0.5*ay*dt^2;
f3 = dz == -u_z*dt + 0.5*az*dt^2;
[xf, yf, dt] = solve([f1, f2, f3], [xf, yf, dt]);
end 

function [d_t, phi2, theta2] = Landing_Angles(dx, dy, dz, u, ax, ay, az)
syms theta2 phi2 d_t
ux = u*cos(theta2*pi/180)*cos(phi2*pi/180);
uy = u*sin(theta2*pi/180)*cos(phi2*pi/180);
uz = u*sin(phi2*pi/180);
f1 = dx == ux*d_t + 0.5*ax*d_t^2;
f2 = dy == uy*d_t + 0.5*ay*d_t^2;
f3 = dz == -uz*d_t + 0.5*az*d_t^2;
[d_t, phi2, theta2] = solve([f1, f2, f3], [d_t, phi2, theta2]);
end

theta = 45; % deg
phi = 0; % deg
u = 16; % m/s
h1 = 60.96; % m
h2 = 0; 
ay = 0; % m/s^2
ax = 0;
az = -9.8;

[xf, yf, dt] = Landing_Physics3D(theta, phi, u, h1, h2, ax, ay, az);
xf = subs(xf);
yf = subs(yf);
a = double(xf);
xf = a(a>=0) % total distance traveled in x direction before impact, m
b = double(yf);
yf = b(b>=0) % total distance traveled in y direction before impact, m
c = double(dt);
tf = c(c>=0); % time before impact, s

dx = xf;
dy = yf;
dz = h2 - h1;
[d_t, phi2, theta2] = Landing_Angles(dx, dy, dz, u, ax, ay, az);
double(d_t)
double(phi2)
double(theta2)