Polar Formulation of the Fourier series for 2d DFT

Question

I know that there is a complex and a polar formulation of the Fourier series. The polar form of can be written as a reformulation of the real valued Fourier series for the 1D case: $$f(x)=a_{0}+\sum_{k}^{\infty} a_{k}\cos(k\omega_{0}x)+b_{k} \sin(k\omega_{0}x)= d_{0}+\sum_{k=1}^{\infty}d_{k}cos(k\omega_{0}x-\phi_{k})$$ with $\omega_{0} := 2\pi /T$ and $$d_{0}=a_{0}, d_k:=\sqrt{a^2_k +b^2_k}\quad\text{as magnitude and the phase}\quad\phi_{k}:=\arctan{\frac{b_{k}}{a_{k}}}$$ , which can be calculated as follows for a dft: $$d_k \cos\phi_k = \sqrt{2/N}\sum^{N}_{j=1}z_j \cos(2\pi k(j-1)/N)$$ $$d_k \sin\phi_k = \sqrt{2/N}\sum^{N}_{j=1}z_j \sin(2\pi k(j-1)/N)$$ It's essentially taking the real and the imaginary parts separately and using trigonometric relationships to calculate the coefficients. Now for the 2D case I have seen this approach in a paper which delivers consistent results: $$d_k \cos\phi_k = 1/N\sum^{N}_{j=1}\left(x_j\cos\left(\frac{2\pi (k+1)(j-1)}{N}\right)-y_j\sin\left(\frac{2\pi (k+1)(j-1)}{N}\right)\right)$$ $$d_k \sin\phi_k = 1/N\sum^{N}_{j=1}\left(x_j\sin\left(\frac{2\pi (k+1)(j-1)}{N}\right)+y_j\cos\left(\frac{2\pi (k+1)(j-1)}{N}\right)\right)$$ And I just do not understand, where this is coming from. I also didn't encounter any book on DFT, which explains this. There were rare ones which mentioned at least the polar formulation. I know how the complex case works. Does someone know how the above formulas for the 2d and 1d case are derived or any literature about the polar form? I would really appreciate any information.