How is it possible to guess what state the qubit was in by measuring it?

Let's describe the one qubit pure stat in Bloch sphere notation (in order to avoid the global phase ambiguity):

$$|psi rangle = cosbig(frac{theta}{2}big) |0rangle + e^{ivarphi} sinbig(frac{theta}{2}big)$$

The problem can be solved with the quantum state tomography, but in this answer, I want to consider a slightly different approach for dealing with the pure states. This answer is a generalized version of this answer. Here we assume that the described state can be prepared as many times as we want.

We are going to execute three different experiments to estimate the $alpha = cosbig(frac{theta}{2}big)$ and $beta = e^{ivarphi} sinbig(frac{theta}{2}big)$, where $theta$ is in $[0, pi]$ range, $varphi$ is in $[-pi, pi)$ range (in the Bloch sphere formalism $varphi$ is in $[0, 2pi)$, but we can, without any problems, take $[-pi, pi)$ range for our convenience). The first experiment will give us $theta$ and the last two experiments will give us $varphi$ from which it will be straightforward to calculate $alpha$, $beta$, and thus the $|psi rangle$ pure state.

The first experiment: determining $theta$.

Execute $N$ times (bigger $N$ will give a better answer) $Z$ basis measurements. Note that the probabilities of measuring $|0rangle$ $left( P(0) right)$ and $|1rangle$ $left(P(1)right)$ states have these relations with the outcomes of the experiment:

$$ P(0) = lim_{N rightarrow infty} frac{N_{0}}{N} = |alpha|^2 = cos^2 big(frac{theta}{2}big) P(1) = lim_{N rightarrow infty} frac{N_{1}}{N} = |beta|^2 = sin^2big(frac{theta}{2}big) $$

So:

$$theta = 2 arccos big(sqrt{(P(0))}big) = 2 arcsin big(sqrt{(P(1))}big)$$

because values of $frac{theta}{2}$ is in the $[0, frac{pi}{2}]$ range.

The second experiment: the absolute value of the relative phase.

We have defined the relative phase $varphi$ in the range $[-pi, pi)$ and here we are going to find the $|varphi|$. For that, we apply a Hadamard gate:

$$H left( cosbig(frac{theta}{2}big) |0rangle + e^{ivarphi} sin big(frac{theta}{2} big) right) = = frac{1}{sqrt{2}} left( cosbig(frac{theta}{2}big) + e^{ivarphi} sin big(frac{theta}{2} big) right) |0rangle + frac{1}{sqrt{2}} left( cosbig(frac{theta}{2}big) - e^{ivarphi} sin big(frac{theta}{2} big) right) |1rangle$$

If we calculate the probabilities of measuring $|0rangle$ ($P'(0)$) and $|1rangle$ $P'(1)$ we will obtain that:

$$P'(0) = frac{1}{2} left| cosbig(frac{theta}{2}big) + e^{ivarphi} sin big(frac{theta}{2} big) right|^2 = frac{1 + cos(varphi) sin(theta)}{2} P'(1) = frac{1}{2} left| cosbig(frac{theta}{2}big) - e^{ivarphi} sin big(frac{theta}{2} big) right|^2 = frac{1 - cos(varphi) sin(theta)}{2} $$

That is why:

$$varphi = pm arccos left( frac{P'(0) - P'(1)}{sin(theta)} right)$$

because the range of usual principal value arccosine function is equal to $[0, pi]$. We know $theta$, we know how to calculate $P'(0)$ and $P'(1)$ from the new experiment with Hadamard gate, so we will be able to find $|varphi|$. If $sin(theta) = 0$, then we can just skip the second and the third experiments, because in that case we have either $|psirangle = |1rangle$ ($theta = pi$) or $|psirangle = |0rangle$ ($theta = 0$).

Also, note that:

$$langle X rangle = langle psi | X | psi rangle = langle psi |H Z H| psi rangle = P(0) - P(1)$$

So, the formula can be written in this way:

$$theta = pm arccos left( frac{langle X rangle}{sin(theta)} right)$$

The third experiment: determining the sign of the relative phase.

For this we will need to apply $S^{dagger}$ then $H$ gates to the initial state before $N$ measurements:

$$H S^{dagger} left( cosbig(frac{theta}{2}big) |0rangle + e^{ivarphi} sin big(frac{theta}{2} big) right) = = frac{1}{sqrt{2}} left( cosbig(frac{theta}{2}big) - i e^{ivarphi} sin big(frac{theta}{2} big) right) |0rangle + frac{1}{sqrt{2}} left( cosbig(frac{theta}{2}big) + i e^{ivarphi} sin big(frac{theta}{2} big) right) |1rangle $$

The probabilities:

$$P''(0) = frac{1}{2} left| cosbig(frac{theta}{2}big) - ie^{ivarphi} sin big(frac{theta}{2} big) right|^2 = frac{1 + sin(varphi) sin(theta)}{2} P''(1) = frac{1}{2} left| cosbig(frac{theta}{2}big) + ie^{ivarphi} sin big(frac{theta}{2} big) right|^2 = frac{1 - sin(varphi) sin(theta)}{2} $$

For the sign of the relative phase we have:

$$sign(varphi) = sign left( arcsin left(frac{P''(0) - P''(1)}{sin(theta)} right) right) = sign left( P''(0) - P''(1) right)$$

because the range of usual principal value of arcsine function is $[-frac{pi}{2}, frac{pi}{2}]$ and $sin(theta) > 0$ $left( theta in [0, pi] right)$.

Also, note that for the expectation value of the $Y$ operator (as can be seen from this answer) we have this relation:

$$langle Y rangle = langle psi| Y | psirangle = langle psi| S H Z H S^{dagger} | psirangle = P''(0) - P''(1)$$

By taking this into account and combining the last two experiments we can obtain the relative phase:

$$varphi = sign big( langle Y rangle big) arccos left( frac{langle X rangle}{sin(theta)} right)$$