How to generate 2D colored noise

I'm unsuccessfully trying to generate a 2D array of pink noise. Can anyone explain me how to do this? Otherwise, is there a python library that allows to generate 2D (or higher dimensionality) colored noise (1/f noise)?

The answers given both here and in a related question (How to generate a pink noise image?) give you most of the picture - but I thought I'd lay out a full step-by-step to show how everything comes together.

Especially as it seems people run into issues understanding the "frequency matrix".

To convert white noise to pink noise

1. Generate some white noise: whitenoise = np.random.uniform(0, 1, (256, 256))
2. Fourier transform and shift low-frequencies to centre (this makes generating the "frequency matrix" more intuitive): ft_arr = np.fft.fftshift(np.fft.fft2(whitenoise))
3. Generate the "frequency matrix". In the 1-D case this would be an array of the actual frequencies which correspond to the amplitudes given by the transform. In the 2-D case, this is the distance from the center of our fftshifted fourier space, as the further out to the edges we go, the greater the frequency captured at that point. This can be generated via np.meshgrid and np.hypot :
   
   _x, _y = np.mgrid[0:ft_arr.shape[0], 0:ft_arr.shape[1]]
f = np.hypot(_x - ft_arr.shape[0] / 2, _y - ft_arr.shape[1] / 2)
   
   We then divide the fourier space by the frequencies.
   pink_ft_arr = ft_arr / f
4. Remove singularities that show up due to having a zero-frequency component. (Probably there's a smarter way to do this, but this gave me the expected output regardless):
   
   pink_ft_arr = np.nan_to_num(pink_ft_arr, nan=0, posinf=0, neginf=0)
   
5. Convert the fourier space pink noise back into image space:
   pinknoise = np.fft.ifft2(np.fft.ifftshift(pink_ft_arr)).real

You can then plot this using matplotlibs plt.imshow, or scale it into some sensible range and write it out as an image.

White noise, from np.random.uniform
Pink noise generated by the above.