<jupyter_start><jupyter_text>Use Astropy to analyze FITS images  Based on a tutorial <NAME><jupyter_code>%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from astropy.io import fits<jupyter_output><empty_output><jupyter_text>Open the FITS file<jupyter_code>fname = "HorseHead.fits"
hdu_list = fits.open(fname)
hdu_list.info()<jupyter_output><empty_output><jupyter_text>Generally, the image information is located in the PRIMARY block.  The blocks are numbered and can be accessed by intexting hdu_list.<jupyter_code>image_data = hdu_list[0].data
header = hdu_list[0].header   # contains some info about the file 
                            # thats contained in the header, in the beggining of the file 
print(header)<jupyter_output><empty_output><jupyter_text>Our data is now stored as a 2-D numpy array.  But how do we know the dimenstions of the image?  We can simply look at the shape of the array.<jupyter_code>print(type(image_data))
print(image_data.shape)<jupyter_output><empty_output><jupyter_text>At this point, we can close the FITS file because we've stored everything we wated to a variable.<jupyter_code>hdu_list.close()<jupyter_output><empty_output><jupyter_text>Shortcut:  use "getdata()" to just read in the image data and close the file.<jupyter_code>image_data = fits.getdata(fname)
print(type(image_data))
print(image_data.shape)<jupyter_output><empty_output><jupyter_text>Let's show the data<jupyter_code>plt.imshow(image_data, cmap='gray')
plt.colorbar()
plt.imshow(image_data, cmap='inferno')
plt.colorbar()<jupyter_output><empty_output><jupyter_text>Let's get some basic statistics about our image:<jupyter_code>print('Min: ', np.min(image_data))
print('Max: ', np.max(image_data))
print('Mean: ', np.mean(image_data))
print('Stdev: ', np.std(image_data))<jupyter_output><empty_output><jupyter_text>Plotting a histogram  To make a histogram with matplotlib.pyplot.hist(), we'll need to cast the data from a 2-D array to something one dimensional. In this case, let's use the ndarray.flatten() to return a 1-D numpy array.<jupyter_code>histogram = plt.hist(image_data.flatten(), bins='auto')<jupyter_output><empty_output><jupyter_text>Display the image with a logarthmic scale What if we want to use a logrithmic color scale ?  To do so, we can load the LogNorm object from matplotlib.<jupyter_code>from matplotlib.colors import LogNorm
plt.imshow(image_data, cmap='bone', norm=LogNorm())

# Choose the tick marks based on the histogram above 
cbar = plt.colorbar(ticks=[5.e3, 1.e3, 2.e4])
cbar.ax.set_yticklabels(['5,000', '10,000', '20,000'])<jupyter_output><empty_output><jupyter_text>Stacking Images Since the noise in an image results from a random process, we use stacking of separate images to improve the signal to noise ratio of objects we observe.  Here we are going to stack 5 images of M13 taken with a 10 inch telescope.<jupyter_code># make a list of filenames
image_list = ['M13_blue_0001.fits', 'M13_blue_0002.fits', 'M13_blue_0003.fits', 'M13_blue_0004.fits', 'M13_blue_0005.fits']
print(image_list)
# make an array of images from the list of images
image_concat = [fits.getdata(image) for image in image_list]
# sum the images together
final_image = np.sum(image_concat, axis=0)
# plot a histogram of the image pixel values
image_hist = plt.hist(final_image.flatten(), bins='auto')
image_hist = plt.hist(image_concat[0].flatten(), bins='auto')
plt.xlim([0,3500])
plt.imshow(image_concat[0])<jupyter_output><empty_output><jupyter_text>We'll use the keywords vmin and vmax to set limits on the color scaling for imshow.<jupyter_code>plt.imshow(final_image, cmap='gray', vmin=2E3, vmax=3E3)
plt.colorbar()<jupyter_output><empty_output><jupyter_text>Writing a new FITS file We can easily do this with the writeto() method.Warning: you'll receive an error if the file you are trying to write already exists.  That's why we've set clobber=True.<jupyter_code>outfile = 'stacked_M13_blue.fits'
hdu = fits.PrimaryHDU(final_image)
hdu.writeto(outfile, overwrite=True)
hdu_list= fits.open(outfile)
header = hdu_list[0].header
data = hdu_list[0].data
plt.imshow(data, cmap='gray', vmin=2E3, vmax=3E3)
print(header)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.10589077033150489, 'fraction_numerical': 0.019079418077748628, 'mean_word_length': 6.653284671532846, 'pattern_counts': {'":': 0, '<': 54, '<?xml version=': 0, '>': 54, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '13', 'dir': 'jupyter-structured-clean-dedup', 'id': '183592', 'license': 'MIT', 'n_tokens_mistral': 1529, 'n_tokens_neox': 1443, 'n_words': 406, 'path': 'astropy_tutorial.ipynb', 'repo_name': 'Veronika131/astr-119-session18', 'size': '4199'}