Creating cutout fits files from one image with a list of RA/DEC positions

Astropy
1

I’m trying to open a large image, and then make cutout thumbnail fits files from that image for a list of RA/DEC positions, but I’m having an issue where things work for all filters for the first RA/DEC position, but the second RA/DEC position kicks an error:
astropy.nddata.utils.NoOverlapError: Arrays do not overlap.

This is because for some reason the original image is replaced with the cutout region when I try to make the fits file.

I’m following things at the bottom of the page here: Image Utilities — Astropy v5.1, where the suggestion is to use:
hdu.data = cutout.data, but I don’t want to overwrite the original image, since it’s very large survey data, and I’d rather just load it into memory and continue making thumbnails from that. But if you make a copy: hdu_for_cutout.data = hdu.data, then when you run hdu_for_cutout.data = cutout.data, it will actually change hdu.data as well. It’s weird that the original image data is inheriting the qualities of the copied variables, but maybe I’m missing something fundamental. Thank you!

2

When you do hdu_for_cutout.data = hdu.data you apparently don’t make a copy of the old hdu into a new hdu_for_cutout, you just make a second name that points to the same object in memory. That is often useful because it allows you to make shorter names for things that you frequently need (e.g. after my_source = hdu[123].data[11200:12000,3400:5300] you can say my_source in your code and don’t have to repeat the long hdu[123].data[11200:12000,3400:5300] every time), but it can be confusing if it happens by accident.

I can’t say for certain what happens exactly in your code, because you don’t show the entire code. How is hdu_for_cutout for cutout defined?

However, you don’t have to overwrite the data in the original hdu - it’s just what is done in that example. You can make a make new fits files from scratch or just work with the array you cut out directly in your code without saving it to a fits file first.

3

If you don’t mind using a different package then you can achieve that easily with gammapy:

from gammapy.maps import Map
from astropy.coordinates import SkyCoord
from astropy import units as u

postion = SkyCoord()
data = Map.read("your-data.fits")
cutout = data.cutout(position=SkyCoord("0d", "0d", frame="icrs"), width=0.1 * u.deg)

cutout.write("thumbnail.fits")

Maybe that helps…