Reputation: 349
I am attempting to set a threshold for pixels in a FITS file using the code below. However, I get an error that says that:
IndexError: Index 2620 is out of bounds for axis 0 with size 2620
Any ideas on how to fix this?
This is the code:
from astropy.io import fits
import numpy as np
hdulist = fits.open("12ratio.fits")
origImData = hdulist[0].data
newImData = origImData*0
for x, y in np.nditer(origImData.shape):
curPixel = origImData[x, y]
if curPixel > 0.28 or curPixel < 3.11:
newImData[x, y] = curPixel
else:
newImData[x, y] = 0
newhdu = fits.PrimaryHDU(newImData)
newhdulist = fits.HDUList([newhdu])
newhdulist.writeto('modifiedratio12.fits')
Upvotes: 0
Views: 545
Reputation: 19205
The problem is your iteration scheme. You're iterating over the shape of the array, not the pixel indices.
In general there is a better approach for masking:
mask = (origImData > 0.28) & (origImData < 3.11)
newImData[mask] = origImData[mask]
will accomplish what you seem to be aiming for.
Or just as effective:
from astropy.io import fits
import numpy as np
hdulist = fits.open("12ratio.fits")
origImData = hdulist[0].data
mask = (origImData > 0.28) & (origImData < 3.11)
origImData[~mask] = 0
newhdu = fits.PrimaryHDU(data=origImData, header=hdulist[0].header)
newhdu.writeto('modifiedratio12.fits')
Here I've used the numpy array boolean not
operation (~
) to invert the mask. I've also skipped making an HDUList
since that is not needed for a single-extension FITS file.
There are two errors in the code you've shown:
for x, y in np.nditer(origImData.shape):
will assign x
and y
the values of origImData.shape
, i.e. x=2640
and y=(something else)
, it will not assign them the values 0,1,...,2640
as you intend.
The comparison operation:
if curPixel > 0.28 or curPixel < 3.11:
will always be true for real numbers: all numbers are either greater than 0.28 or less than 3.11.
Upvotes: 1