syms d L R x y phi theta
assume (L,'positive')
eq(1) = -L*cos(theta) + x == 0;
eq(2) = -L*sin(theta) - L/2 + y == 0;
eq(3) = -0.5*L*cos(theta) + d*cos(phi) == 0;
eq(4) = -0.5*L*sin(theta) - d*sin(phi) + y == 0;
S = solve([eq(1),eq(2),eq(3),eq(4)],[x y phi theta])
x = simplify(subs(S.x,d,L/2))
y = simplify(subs(S.y,d,L/2))
phi = simplify(subs(S.phi,d,L/2))
theta = simplify(subs(S.theta,d,L/2))
gamma = asin(L/(4*R))
beta = pi/2 - phi - gamma