The derivation of Boozer coordinates requires several steps of analysis. First, a general transformation is introduced that converts the familiar laboratory coordinate system into a set of arbitrary flux coordinates. Second, by using the relationships B ⋅ ∇ψ = J ⋅ ∇ψ = 0 and ∇ ⋅ B = ∇ ⋅ J = 0, both B and J can be cast into a cross-product form in flux coordinates, close to the desired form of Boozer coordinates. Third, by means of the relation ∇ × B = μ 0J, it is shown that B can also be written in a gradient form in flux coordinates, close to the desired form of Boozer coordinates. Fourth, it is shown how certain free functions appearing in the representation of B can be eliminated by means of an additional transformation of the angular flux coordinates χ, ζ. The new coordinates correspond to the actual Boozer coordinates. Fifth, the various free functions remaining in the expressions for B are rewritten in terms of physically recognizable quantities. Finally, the magnetic field expressed in Boozer coordinates is used to calculate the guiding center drifts of the particles.