Suppose I have a two dimensional Hilbert space $\{ |0 \rangle,|1\rangle \}$ with these states being orthonormal. Now suppose I have the Hamiltonian $H=|1\rangle \langle 0|+|0\rangle \langle 1| .$
It is clear that the (normalised) eigenstates of this are:
$|\phi_0\rangle=\frac{1}{\sqrt{2}}\Big(|0\rangle+|1\rangle\Big)$ and $|\phi_1\rangle=\frac{1}{\sqrt{2}}\Big(|0\rangle-|1\rangle\Big)$ with eigenvalues $+1$ and $-1$ respectively.
Now what I am confused about is calculating probabilities: I know that given the state $| \psi(t)\rangle$ of the system, if I want to calculate the probability of this state being in an arbitrary state $|\phi\rangle$, then I just need to calculate their overlap $\big($ie. $\mathcal{P}(\phi)=|\langle \phi | \psi(t) \rangle |^2\big).$ Suppose I chose the state to be $|0\rangle$, then I would need to find its overlap.
What is the actual grounding for this calculation $\dots$ From the principles of QM the only outcomes from measurement are eigenvalues of the operator in question, and the probability is given by the overlap of the eigenstate with the time evolving state of the system. So I should only be able to find the probabilities of being in states $|\phi_0 \rangle$ or $|\phi_1 \rangle$ by this axiom (obviously not true though!)
So how can I find the probability of being in state $|0\rangle$ if it isn't an eigenstate of the Hamiltonian operator. Is there an operator corresponding to the original basis states which we can measure (though this would be the identity)? What is the intuition behind this definition so that I can some appreciation as to why it makes sense to compute the overlap of states?