Mc Missile
Mc Missile

Reputation: 735

Making an Image Cube using multiple fits images

I have a bunch of fits files that can be read using the below script

from astropy.io import fits
hdu = fits.open('file.fits')
data = hdu[0].data

I am trying to make an image cube using the data read from multiple fits files. (An image cube is a 3D image that contains data from multiple fits file where the x and y axis is the 2D image dimension and the 3rd axis is time or frequency)

I believe it can be done usin spectral _cube module, however most of the documentation only talks about how to read an image cube and not how to make one using individual fits files.

So far I have tried the following script.

#In the below script data is a 3D numpy array
from spectral_cube import SpectralCube
cube = SpectralCube(data=data)
cube.write('new_cube.fits', format='fits')

However, the above script gives an error saying 3 arguments required while only 2 given.

Upvotes: 2

Views: 5182

Answers (2)

Iva Laginja
Iva Laginja

Reputation: 242

The simplest way to do this is to simply put the images you want to have in your cube into a numpy array and then save that array as a fits file. You can also save them into a numpy array directly, but appending lists is easier if you do it in a for loop, instead of doing it explicitly for each image, like I do it here.

import numpy as np
from astropy import fits

# Read the images you want to concatenate into a cube
img1 = fits.getdata('img1.fits')
img2 = fits.getdata('img2.fits')

# Make a list that will hold all your images
img_list = []
img_list.append(img1)
img_list.append(img2)

# Cast the list into a numpy array
img_array = np.array(img_list)

# Save the array as fits - it will save it as an image cube
fits.writeto('mycube.fits', img_array)

Upvotes: 3

Mc Missile
Mc Missile

Reputation: 735

Apparently, one does not need to use the spectral_cube module to make image cubes. It can be much easily done using the AstroPy python module. Below is the script.

from astropy.io import fits
import numpy as np

cube = np.zeros((50,1000,1000)) #Here 1000x1000 is the dimension of the individual fits images and 50 is the third perpendicular axis(time/freq)

for i in range(50):
    hdu = fits.open('image' + str(i) + '.fits') #The fits images that you want to combine have the name string 'image' + str(i) + '.fits'
    data = hud[0].data[:,:]
    cube[i,:,:] = data

hdu_new = fits.PrimaryHDU(cube)
hdu_new.writeto('cube.fits')

Upvotes: 2

Related Questions