How to calculate Berry connection of a Hamiltonian numerically?

I hope there is a method to calculate berry connection of a Hamiltonian numerically, in order to calculate $theta$-term like:
$$theta=frac{1}{4 pi} int d^{3} k epsilon^{i j k} operatorname{Tr}left[mathcal{A}_{i} partial_{j} mathcal{A}_{k}+i frac{2}{3} mathcal{A}_{i} mathcal{A}_{j} mathcal{A}_{k}right]$$
my obstacle is how to avoid the random phase in numerical calculation.