Um, what do you expect?
syms theta0
x=90; g=9.81; v0=30; y0=1.8; y=1;
th = -y+tan(theta0)*x-g*x^2/(2*v0^2*cos(theta0)^2)+y0;
pretty(th)
Do you expect an analytical solution, something nice and simple and very pretty?
thsol = solve(th)
You should be able to effectively reduce this to effectively a 4th degree polynomial in terms of sin(theta0), so 4 roots.
vpa(thsol)
There are probably infinitely many solutions since trig functions are periodic. I don't see how much better than that you can rationally expect? Honestly, I'd say you are pretty lucky, in that solutions are easily generated\. Just pick the one that makes sense to you.