Thanks to the discussion in the comments (Anthony's argument), we already know that this holds for smooth functions $f(z)=\sum a_nz^n$ (equivalently, it holds for rapidly decaying $a_n$'s).

So it now suffices to show that an arbitrary $f(e^{ix})=e^{ig(x)}\in H^{1/2}$ can be approximated in $H^{1/2}$ by smooth functions $e^{ih_k(x)}$: we then have that $\sum n|a_n|^2 = \lim_{k\to\infty} \sum n|a_n^{(k)}|^2$, and these latter sums are integers.

To do this, we recall that $H^{1/2}$ functions (in dimension $1$) can be characterized as those $L^2$ functions $u$ that also satisfy $$ \int\!\!\!\int dxdy\, \frac{|u(y)-u(x)|^2}{|x-y|^2} < \infty . \quad\quad\quad\quad (1) $$ The $L^2$ condition is automatic here since our functions satisfy $|f|=1$, and (1) for $u=e^{ig}$ implies the same condition for $g$ after some precautions associated with the non-uniqueness of $g$. More specifically, let's assume that $|u^{-1}(-1)|=0$, and let's then take $|g(x)|<\pi$. Then $|u(y)-u(x)|$ is comparable to $|g(x)-g(y)|$ on every set $A_n=\{ x: |g(x)\mp\pi|\ge 1/n \}$, and these exhaust our domain.

Thus $g\in H^{1/2}$ also, as claimed, and if we now approximate $h_k\to g$ in $H^{1/2}$, then $e^{ih_k}\to e^{ig}$ also, by a similar argument: (1) for $u=e^{ih_k}-e^{ig}$ together with $\|u\|_2$ essentially computes the $H^{1/2}$ norm, but this is controlled by (1) for $u=h_k-g$.

As correctly pointed by Anthony in the comments below, the argument I attempted here originally did not hold water.

These issues apparently have been studied quite extensively recently. See for example J. Bourgain, One cannot hear the winding number. In particular, in the introduction to this paper, Bourgain mentions that the answer to your question is yes; the result is attributed to Boutet de Monvel and Gabber.