The Bloch sphere is beautifully minimalist.
Conventionally, a qubit has four real parameters; $a e^{i\chi} |0\rangle + b e^{i\phi} |1\rangle.$ However, some quick insight reveals that the $a$-vs-$b$ tradeoff only has one degree of freedom due to the normalization $a^2 + b^2 = 1$ and some more careful insight reveals that, in the way we construct expectation values in QM, you cannot observe $\chi$ or $\phi$ themselves but only the difference $\chi - \phi$, which is $2\pi$-periodic. (This is covered further in the comments below but briefly: QM only predicts averages $\langle \psi|\hat A|\psi\rangle$ and shifting the overall phase of a wave function by some $|\psi\rangle\mapsto e^{i\theta}|\psi\rangle$ therefore cancels itself out in every prediction.)
So if you think at the most abstract about what you need, you just draw a line from 0 to 1 representing the $a$-vs-$b$ tradeoff: how much is this in one of these two states? Then you draw circles around it: how much is the phase difference? What stops it from being a cylinder is that the phase difference ceases to matter when $a=1$ or $b=1$, hence the circles must shrink down to points. Et voila, you have something which is topologically equivalent to a sphere. The sphere contains all of the information you need for experiments, and nothing else.
It’s also physical, a real sphere in 3D space.
This is the more shocking fact. Given only the simple picture above, you could be forgiven for thinking that this was all harmless mathematics: no! In fact the quintessential qubit is a spin-$\frac 12$ system, with the Pauli matrices indicating the way that the system is spinning around the $x$, $y$, or $z$ axes. This is a system where we identify $|0\rangle$ with $|\uparrow\rangle$, $|1\rangle$ with $|\downarrow\rangle$, and the phase difference comes in by choosing the $+x$-axis via $|{+x}\rangle = \sqrt{\frac 12} |0\rangle + \sqrt{\frac 12} |1\rangle.$
The orthogonal directions of space are not Hilbert-orthogonal in the QM treatment, because that’s just not how the physics of this system works. Hilbert-orthogonal states are incommensurate: if you’re in this state, you’re definitely not in that one. But this system has a spin with a definite total magnitude of $\sqrt{\langle L^2 \rangle} = \sqrt{3/4} \hbar$, but only $\hbar/2$ of it points in the direction that it is “most pointed along,” meaning that it must be distributed on some sort of “ring” around that direction. Accordingly, when you measure that it’s in the $+z$-direction it turns out that it’s also sort-of half in the $+x$, half in the $-x$ direction. (Here “sort-of” means: it is, if you follow up with an $x$-measurement.)
So let’s ask “which direction is the spin-$\frac12$ most spinning in?” This requires constructing an observable. To give an example, if the $+z$-direction is most-spun-in by a state $|\uparrow\rangle$ then the observable for $z$-spin is the Pauli matrix $\sigma_z = |\uparrow\rangle\langle\uparrow| - |\downarrow\rangle\langle\downarrow|,$ $+1$ in that state, $-1$ in the Hilbert-perpendicular state $\langle \downarrow | \uparrow \rangle = 0.$ Similarly if you look at $\sigma_x = |\uparrow\rangle \langle \downarrow | + |\downarrow \rangle\langle \uparrow |$ you will see that the $|{+x}\rangle$ state defined above is an eigenvector with eigenvalue +1 and similarly there should be a $|{-x}\rangle \propto |\uparrow\rangle - |\downarrow\rangle$ satisfying $\langle {+x}|{-x}\rangle = 0,$ and you can recover $\sigma_x = |{+x}\rangle\langle{+x}| - |{-x}\rangle\langle{-x}|.$
Then the state orthogonal to $|\psi\rangle = \alpha |0\rangle + \beta |1\rangle$ is $|\bar \psi\rangle = \beta^*|0\rangle - \alpha^* |1\rangle,$ so the observable which is +1 in that state or -1 in the opposite state is:$$
\begin{align}
|\psi\rangle\langle\psi| - |\bar\psi\rangle\langle\bar\psi| &= \begin{bmatrix}\alpha\\\beta\end{bmatrix}\begin{bmatrix}\alpha^*&\beta^*\end{bmatrix} - \begin{bmatrix}\beta^*\\-\alpha^*\end{bmatrix} \begin{bmatrix}\beta & -\alpha\end{bmatrix}\\
&=\begin{bmatrix}|\alpha|^2 - |\beta|^2 & 2 \alpha\beta^*\\
2\alpha^*\beta & |\beta|^2 - |\alpha|^2\end{bmatrix}
\end{align}$$Writing this as $v_i \sigma_i$ where the $\sigma_i$ are the Pauli matrices we get:$$v_z = |\alpha|^2 - |\beta|^2,\\
v_x + i v_y = 2 \alpha^* \beta.$$
Now letting $\alpha = \cos(\theta/2)$ and $\beta = \sin(\theta/2) e^{i\phi}$ we find out that these are:$$\begin{align} v_z &= \cos^2(\theta/2) - \sin^2(\theta/2) &=&~ \cos \theta,\\
v_x &= 2 \cos(\theta/2)\sin(\theta/2) ~\cos(\phi) &=&~ \sin \theta~\cos\phi, \\
v_y &= 2 \cos(\theta/2)\sin(\theta/2) ~\sin(\phi) &=&~ \sin \theta~\sin\phi.
\end{align}$$So the Bloch prescription uses a $(\theta, \phi)$ which are simply the spherical coordinates of the point on the sphere which such a $|\psi\rangle$ is “most spinning in the direction of.”
So instead of being a purely theoretical visualization, we can say that the spin-$\frac 12$ system, the prototypical qubit, actually spins in the direction given by the Bloch sphere coordinates! (At least, insofar as a spin-up system spins up.) It is ruthlessly physical: you want to wave it away into a mathematical corner and it says, “no, for real systems I’m pointed in this direction in real 3D space and you have to pay attention to me.”
How these answer your questions.
Yes, N and S are spatially parallel but in the Hilbert space they are orthogonal. This Hilbert-orthogonality means that a system cannot be both spin-up and spin-down. Conversely the lack of Hilbert-orthogonality between, say, the $z$ and $x$ directions means that when you measure the $z$-spin you can still have nonzero measurements of the spin in the $x$-direction, which is a key feature of such systems. It is indeed a little confusing to have two different notions of “orthogonal,” one for physical space and one for the Hilbert space, but it comes from having two different spaces that you’re looking at.
One way to see why the angles are physically very useful is given above. But as mentioned in the first section, you can also view it as a purely mathematical exercise of trying to describe the configuration space with a sphere: then you naturally have the polar angle as the phase difference, which is $2\pi$-periodic, so that is a naturally ‘azimuthal’ coordinate; therefore the way that the coordinate lies along 0/1 should be a ‘polar’ coordinate with $0$ mapping to $|0\rangle$ and $\pi$ mapping to $|1\rangle$. The obvious way to do this is with $\cos(\theta/2)$ mapping from 1 to 0 along this range, as the amplitude for the $|0\rangle$ state; the fact that $\cos^2 + \sin^2 = 1$ means that the $|1\rangle$ state must pick up a $\sin(\theta/2)$ amplitude to match it.
Since the action of a unitary matrix $U$ on the Bloch sphere is defined by
$$\sigma_j\mapsto U\,\sigma_j\,U^\dagger = U\,\sigma_j\,U^{-1}\tag{1}$$
with $\sigma_j$ standing for the Pauli matrices, any scalar phase factor $e^{i\,\phi}$ is quotiented away by the map (which is the big A Adjoint representation of the unitary group). Unitary matrices always act on the Bloch sphere by rotation, never reflexion. This is because $U(2)$, as opposed to $O(3)$, is a topologically connected group and a reflexion matrix, say $\left(\begin{array}{cc}-1&0\\0&1\end{array}\right)$, is simply the $SU(2)$ member $i\,\sigma_z=\left(\begin{array}{cc}i&0\\0&-i\end{array}\right)$ multipied by the phase factor $i$. So the action of the two are the same rotation.
Given 1, the Bloch vector $S =\sum\limits_j s_j\,\sigma_j$ becomes $S^\prime = \sum\limits_j s_j^\prime\,\sigma_j$ where:
$$s_j^\prime = \frac{1}{2}\,\mathrm{tr}(\sigma_j\,\sum\limits_k\,(s_k\,U\,\sigma_k\,U^\dagger))\tag{2}$$
So the matrix of your Bloch sphere rotator is the $3\times 3$ matrix $\left[\frac{1}{2}\,\mathrm{tr}(\sigma_j\,U\,\sigma_k\,U^\dagger))\right]$.
You can find the axis of rotation by finding the eigenvector of the rotation matrix that corresponds to the eigenvalue 1, which is usually the only real eigenvector (unless the rotation is trivial).
Alternatively you should first divide out any phase factor to set your matrix's determinant to 1, so that you're dealing with a member of $SU(2)$. This is not essential for the method above (because the phase is cancelled by the Adjoint representation) but it is here. Now calculate the components of the Pauli matrices in your operator:
$$u_j = -\frac{i}{2}\,\mathrm{tr}(\sigma_j\,U)$$
The vector with the $u_j$ as components is the axis of rotation, where the $j^{th}$ co-ordinate axis is the axis of rotation of the Pauli matrix $\sigma_j$. The magnitude $\sqrt{\sum_j\,u_j^2}$ is the sine of half the rotation angle.
Question from OP
Anyway, why is it possible to find a rotation axis for any unitary matrix?It is possible to decompose a unitary matrix in the product of the three Operators x,y,z , so there shoul be 3 axes? –
Answer: (1) The image of the group $U(2)$ under the Adjoint representation is the group $SO(3)$; then (2) Euler's rotation theorem shows that all members of the latter group are rotations.
I tried to find the rotation axis of the matrix .... This gives me the $e_x$
as rotation axis which is not correct.
This is what I get for your matrix. Applying my second method, the matrix is, modulo phase, equivalent to the unimodular (member of $SU(2)$) matrix $\frac{1}{\sqrt{2}}\left(\begin{array}{cc}1&i\\i&1\end{array}\right)$ (multiply everything by $-i$). This, by inspection, is the matrix $\exp\left(i\,\frac{\pi}{4}\,\sigma_x\right) = \cos\left(\frac{\pi}{4}\right)\,\mathrm{id}+ i\,\sin\left(\frac{\pi}{4}\right)\,\sigma_x$. So, depending on how you define your rotation axes, it's a rotation of $2\times \frac{\pi}{4}$, i.e. a right angled turn about the axis corresponding to $\sigma_x$. So, unless I have done something dreadfully wrong, I agree with your calculations and there must be a typo somewhere if the answer is meant to be a rotation about $(1,\,1,\,1)$.
For a matrix equivalent (modulo a phase) to $\frac{1}{\sqrt{2}} \begin{pmatrix} 1&-i\\1 &i\end{pmatrix}$, we multiply by $\exp(-\pi\,i/4)$ to get an $SU(2)$ member. Then my second method entails the following calculations (I checked the matrix was unitary first):
So indeed the axis is along $(1,\,1,\,1)$. Alternatively, we can simply calculate the matrix of the the transformation and work out which eigenvector corresponds to the eigenvalue 1; the Frobenius (trace norm) of the matrix logarithm gives us the angle of rotation. So its a rotation of through an angle of $2\pi/3$ about the axis $(1,1,1)$. These are the relevant calculations:
Best Answer
This is a very good question that many people ask when they learn about the Bloch sphere. There are two important reasons
Physical reason
The most famous two state system in quantum mechanics is the spin 1/2 particle. When you draw the spin 1/2 states on the Bloch sphere something very nice happens, the points on the Bloch sphere correspond to the physical directions of the states. For example the $|\uparrow\,\rangle$ state corresponds to a point that points in the $+z$ direction. Similarly the state $|x+\rangle=\tfrac{ 1}{\sqrt{2}}(|\uparrow\,\rangle+|\downarrow\,\rangle)$ points in the $+x$ direction on the Bloch sphere. To be more precise for some state $|\psi\rangle$ the vector of expectation values given by $\langle\psi|\vec S|\psi\rangle=\langle \vec S\rangle =(\langle S_x\rangle,\langle S_y\rangle,\langle S_z\rangle)^T$ points in the same direction as the Bloch vector. See also this picture:
Mathematical reason
Another reason that the Bloch representation is good is that states which are physically the same correspond to the same point on the Bloch sphere. This is also mentioned in the comments. In your diagram $|0\rangle$ and $-|0\rangle$ are antipodal points but they represent the same state up to an arbitrary phase. Remember that $|\psi\rangle$ and $e^{i\alpha}|\psi\rangle$ represent the same state. The Bloch sphere intuitively leaves out any global phase. To represent the full state on a Bloch sphere you would have to specify an additional phase for any point. In Steane's "Introduction to Spinors" this is represented by drawing a state as a little flag. The direction of the flag gives the point on the Bloch sphere and the angle that flag makes gives the global phase. So in this representation the states $|0\rangle$ and $-|0\rangle$ both point in the z-direction but the first one has its flag rotated to $\alpha=0$ and the second one has it rotated to $\alpha=\pi$ radians. Interestingly when you rotate a spin 1/2 particle by $360^\circ$ its state will pick up a minus sign and you would have to rotate it by another $360^\circ$ before it becomes exactly the same state.