It's common knowledge (and has been discussed in other questions on this site) that the standard BCS ground state $ \left|\Psi_{BCS}\right\rangle = \prod_k \left( u_k + v_k c_{k\uparrow}^{\dagger} c_{-k\downarrow}^{\dagger}\right) \left|0\right\rangle$ does not have a well-defined particle number and that this doesn't matter in bulk superconductors because the standard deviation is $\Delta N \propto \sqrt{N}$ and hence irrelevant for $N\rightarrow \infty$.
But I also read that you can arrive at a BCS state with well-defined particle number by first defining
$$\left|\Psi_{BCS}(\phi)\right\rangle = \left( |u_k| + e^{i\phi} |v_k| c_{k\uparrow}^{\dagger} c_{-k\downarrow}^{\dagger}\right) \left|0\right\rangle $$
and then "integrating out" the phase according to
$$\left|\Psi_{BCS}(N)\right\rangle = \int_{0}^{2\pi} \mathrm{d}\phi\, e^{-iN\phi/2} \left|\Psi_{BCS}(\phi)\right\rangle \,\, ,$$ which gives you a BCS state with precisely N particles at the cost of having a completely ill-defined phase.
If that is true, then surely this is the actual "physical" state of a superconductor and the original BCS state is merely used for convenience.
But that, in turn, would make the well-defined phase of the superconducting state a mere mathematical artifact, when every other textbook highlights it as something very fundamental (and if I remember correctly, it's very important for things like the Josephson effect as well).
So, can anyone point to the error in the logic above?