Converting an RGB image to fits (astropy) I am trying to save a NumPy RGB array (a 3d array with dimension 512x512x3), to a .fits file using astropy.io.
Saving this array to a .png file works fine, but when I'm trying to save it to a .fits file using HDUList.writeto, it produces just garbage colors. 
Also, I have to reshape the image like so:

img = img.reshape(img.shape[2], img.shape[0], img.shape[1])


or it will try to store just 3 column of the image instead of 512.

Each element in the 3d array contains a value between 0 and 255.
I haven't done much astronomical image processing before, but as this question is unanswered I'll give it a shot - hopefully to some avail. If the problem is more specific, a code sample/image sample would probably be useful for further diagnosis, but otherwise this example may help. It discusses the process of writing a 3-channel image to separate FITS images. I would try this first, and check if this output is also garbage - it might just be an issue with how you're using the HDUList.writeto() function. The relevant code is pasted below (with a couple of edits, as the example uses the Pillow Image class at first):

Split the three channels (RGB):

r, g, b = img[:, :, 0], img[:, :, 1], img[:, :, 2]


Write out the channels as separate FITS images

red = fits.PrimaryHDU(data=r)
red.writeto('red.fits')

green = fits.PrimaryHDU(data=g)
green.writeto('green.fits')

blue = fits.PrimaryHDU(data=b)
blue.writeto('blue.fits')


Also out of curiosity - why is your "numpy RGB array" storing values between 0 and 512? The standard for images is usually integers between 0-255 or floats between 0 and 1; but my guess is that astronomy requires a larger bit depth? In any case, you should also make sure the dtype of your numpy array is large enough for the desired bit depth - often images are stored with dtype="uint8" which ranges from 0 to 255.