Haipeng Li Lab 2 python code (1)
pdf
keyboard_arrow_up
School
University of California, Irvine *
*We aren’t endorsed by this school
Course
115A
Subject
Industrial Engineering
Date
Jan 9, 2024
Type
Pages
8
Uploaded by DeanResolve98810
10/19/23, 11:05 PM
lab 2 - Jupyter Notebook
localhost:8888/notebooks/lab 2.ipynb
1/8
In [5]:
In [5]:
In [6]:
In [7]:
In [8]:
In [9]:
In [10]:
In [11]:
In [12]:
Filename: kb200915_00021.fits
No.
Name
Ver
Type
Cards
Dimensions
Format
0
PRIMARY
1 PrimaryHDU
417
(2300, 2200)
int16 (rescales t
o uint16)
1
Exposure Events
1 TableHDU
17
4R x 2C
[A23, A20]
Maximum count value:65535
Minimum count value:0
Dimensions of the image: 2200 x 2300
import
astropy.io.fits
as
fits
import
matplotlib.pyplot
as
plt
hdul
=
fits.open(
'kb200915_00021.fits'
)
hdul.info( )
fits_file
=
fits.open(
'kb200915_00021.fits'
)
image_data
=
fits_file[
0
].data
plt.imshow(image_data, cmap
=
'gray'
)
plt.show()
max_count
=
image_data.max()
min_count
=
image_data.min()
print
(
f"Maximum count value:
{max_count}
"
)
print
(
f"Minimum count value:
{min_count}
"
)
height, width
=
image_data.shape
print
(
f"Dimensions of the image:
{height}
x
{width}
"
)
10/19/23, 11:05 PM
lab 2 - Jupyter Notebook
localhost:8888/notebooks/lab 2.ipynb
2/8
In [13]:
In [14]:
X pixel positions defining the extent of the spectral slice: 1420 to 1470
Maximum count value in the spectral slice: 49539
Minimum count value in the spectral slice: 9946
start_column
=
1420
end_column
=
1470
spectral_slice
=
image_data[:, start_column:end_column]
x_start
=
start_column
x_end
=
end_column
max_count_slice
=
spectral_slice.max()
min_count_slice
=
spectral_slice.min()
print
(
f"X pixel positions defining the extent of the spectral slice:
{x_sta
print
(
f"Maximum count value in the spectral slice:
{max_count_slice}
"
)
print
(
f"Minimum count value in the spectral slice:
{min_count_slice}
"
)
import
numpy
as
np
spectrum
=
np.sum(spectral_slice, axis
=
1
)
y_positions
=
np.arange(
len
(spectrum))
plt.plot(y_positions, spectrum)
plt.xlabel(
'Y Pixel Position'
)
plt.ylabel(
'Counts'
)
plt.title(
'Spectral Slice Spectrum'
)
plt.grid(
True
)
plt.show()
10/19/23, 11:05 PM
lab 2 - Jupyter Notebook
localhost:8888/notebooks/lab 2.ipynb
3/8
In [17]:
In [ ]:
In [28]:
Standard deviation (scatter) of counts: 995.9909034665728
Median of the distribution: 43952.5
Mean of the distribution: 44007.4072
Standard deviation of the Poisson distribution: 209.40476072486985
x_start_roi
=
1420
x_end_roi
=
1470
y_start_roi
=
100
y_end_roi
=
150
roi
=
image_data[y_start_roi:y_end_roi, x_start_roi:x_end_roi]
counts_roi
=
roi.flatten()
plt.hist(counts_roi, bins
=
50
, color
=
'blue'
, alpha
=
0.7
)
plt.xlabel(
'Count Level'
)
plt.ylabel(
'Number of Pixels'
)
plt.title(
'Distribution of Counts in Region of Interest'
)
plt.grid(
True
)
plt.show()
std_count
=
np.std(counts_roi)
print
(
f"Standard deviation (scatter) of counts:
{std_count}
"
)
median_count
=
np.median(counts_roi)
mean_count
=
np.mean(counts_roi)
print
(
f"Median of the distribution:
{median_count}
"
)
print
(
f"Mean of the distribution:
{mean_count}
"
)
mean_of_dataset
=
np.mean(counts_roi)
num_values
=
len
(counts_roi)
poisson_distribution
=
np.random.poisson(mean_of_dataset, num_values)
std_of_poisson
=
np.std(poisson_distribution)
print
(
f"Standard deviation of the Poisson distribution:
{std_of_poisson}
"
)
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
10/19/23, 11:05 PM
lab 2 - Jupyter Notebook
localhost:8888/notebooks/lab 2.ipynb
4/8
In [29]:
In [30]:
In [19]:
Standard deviation of the Normal distribution: 212.242041706014
The distribution of counts per pixel in the subset of the image appears to
follow a...
Poisson distribution, as it matches the Poisson distribution more closely.
Gain for this image: 0.161
num_values
=
len
(counts_roi)
normal_distribution
=
np.random.normal(mean_of_dataset, np.sqrt(mean_of_dat
std_of_normal
=
np.std(normal_distribution)
print
(
f"Standard deviation of the Normal distribution:
{std_of_normal}
"
)
plt.hist(counts_roi, bins
=
50
, color
=
'blue'
, alpha
=
0.7
, label
=
'Data (Counts
plt.hist(poisson_distribution, bins
=
50
, color
=
'red'
, alpha
=
0.5
, label
=
'Pois
plt.hist(normal_distribution, bins
=
50
, color
=
'green'
, alpha
=
0.5
, label
=
'Nor
plt.xlabel(
'Count Level'
)
plt.ylabel(
'Number of Pixels'
)
plt.title(
'Distribution of Counts with Poisson and Normal Distributions'
)
plt.legend()
plt.grid(
True
)
plt.show()
print
(
"The distribution of counts per pixel in the subset of the image appe
if
std_of_poisson
<
std_of_normal:
print
(
"Poisson distribution, as it matches the Poisson distribution mor
else
:
print
(
"Normal distribution, as it matches the Normal distribution more
from
astropy.io
import
fits
image_data
=
fits.open(
'kb200915_00021.fits'
)[
0
].data
hdr
=
fits.getheader(
'kb200915_00021.fits'
,
0
)
gain
=
hdr[
'GAIN1'
]
image_electrons
=
image_data
*
gain
print
(
f"Gain for this image:
{gain}
"
)
10/19/23, 11:05 PM
lab 2 - Jupyter Notebook
localhost:8888/notebooks/lab 2.ipynb
5/8
In [37]:
In [ ]:
In [ ]:
In [ ]:
{'0.7"': 12, '0.8"': 71, '1.0"': 121, '1.2"': 173, '1.5"': 233}
from
scipy.signal
import
find_peaks
import
glob
file_pattern
=
'd0923_*.fits.gz'
file_names
=
glob.glob(file_pattern)
if
not
file_names:
print
(
"No files found matching the pattern!"
)
exit()
data
=
None
with
fits.open(file_names[
0
])
as
hdulist:
for
hdu
in
hdulist:
if
hdu.data
is
not
None
:
data
=
hdu.data
break
if
data
is
None
:
print
(
f"No data found in
{file_names[
0
]}
!"
)
exit()
profile
=
np.sum(data, axis
=
0
)
peaks, _
=
find_peaks(profile, distance
=
50
)
slit_widths
=
[
'0.7"'
,
'0.8"'
,
'1.0"'
,
'1.2"'
,
'1.5"'
]
slit_centers
=
{}
for
idx, width
in
enumerate
(slit_widths):
slit_centers[width]
=
peaks[idx]
print
(slit_centers)
10/19/23, 11:05 PM
lab 2 - Jupyter Notebook
localhost:8888/notebooks/lab 2.ipynb
6/8
In [44]:
file_pattern
=
'd0923_*.fits.gz'
file_names
=
glob.glob(file_pattern)
exposure_times
=
[]
counts_per_pixel_chip3
=
[]
counts_per_pixel_chip7
=
[]
for
file_name
in
file_names:
hdulist
=
fits.open(file_name)
if
hdulist[
2
].data
is
None
or
hdulist[
1
].data
is
None
:
print
(
"NONE"
)
header
=
hdulist[
0
].header
exposure_time
=
header.get(
'EXPTIME'
,
0.0
)
exposure_times.append(exposure_time)
image_data_chip3
=
hdulist[
1
].data
image_data_chip7
=
hdulist[
2
].data
mean_counts_chip3
=
np.mean(image_data_chip3)
mean_counts_chip7
=
np.mean(image_data_chip7)
counts_per_pixel_chip3.append(mean_counts_chip3)
counts_per_pixel_chip7.append(mean_counts_chip7)
hdulist.close()
median_counts_chip3
=
np.median(counts_per_pixel_chip3)
median_counts_chip7
=
np.median(counts_per_pixel_chip7)
errors_chip3
=
np.std(counts_per_pixel_chip3)
errors_chip7
=
np.std(counts_per_pixel_chip7)
counts_per_second_chip3
=
np.array(counts_per_pixel_chip3)
/
np.maximum(np.
counts_per_second_chip7
=
np.array(counts_per_pixel_chip7)
/
np.maximum(np.
plt.figure(figsize
=
(
10
,
6
))
plt.plot(exposure_times, [median_counts_chip3]
*
len
(exposure_times), label
plt.plot(exposure_times, [median_counts_chip7]
*
len
(exposure_times), label
plt.xlabel(
"Exposure Time (seconds)"
)
plt.ylabel(
"Median Counts per Pixel"
)
plt.legend()
plt.title(
"Median Counts per Pixel as a Function of Exposure Time"
)
plt.grid()
plt.show()
plt.figure(figsize
=
(
10
,
6
))
plt.plot(exposure_times, counts_per_second_chip3, label
=
"Chip 3 Counts per
plt.plot(exposure_times, counts_per_second_chip7, label
=
"Chip 7 Counts per
plt.xlabel(
"Exposure Time (seconds)"
)
plt.ylabel(
"Counts per Pixel per Second"
)
plt.legend()
plt.title(
"Counts per Pixel per Second as a Function of Exposure Time"
)
plt.grid()
plt.show()
plt.figure(figsize
=
(
10
,
6
))
plt.errorbar(exposure_times, counts_per_pixel_chip3, yerr
=
errors_chip3, lab
plt.errorbar(exposure_times, counts_per_pixel_chip7, yerr
=
errors_chip7, lab
plt.xlabel(
'Exposure Time (seconds)'
)
plt.ylabel(
'Counts per Pixel'
)
plt.legend()
plt.title(
'CCD Response with Error Bars'
)
plt.grid(
True
)
plt.show()
Your preview ends here
Eager to read complete document? Join bartleby learn and gain access to the full version
- Access to all documents
- Unlimited textbook solutions
- 24/7 expert homework help
10/19/23, 11:05 PM
lab 2 - Jupyter Notebook
localhost:8888/notebooks/lab 2.ipynb
7/8
10/19/23, 11:05 PM
lab 2 - Jupyter Notebook
localhost:8888/notebooks/lab 2.ipynb
8/8
In [ ]: