Seestar S50: Строим диаграмму Герцшпрунга-Рассела для шарового скопления M 13

[apximhd]

Как я уже писал ранее, смарт-телескоп Seestar S50 сохраняет изображения в формате FITS. Пайтоновская библиотека astropy умеет работать с этим форматом. А раз так, можно с игрушкой Seestar S50 делать вполне профессиональные вещи: например, построить диаграмму Герцшпрунга–Рассела для какого-нибудь шарового скопления.

Поскольку размер типичного шарового скопления мал по сравнению с расстоянием до него, можно считать все его звёзды находящимися примерно на одинаковом расстоянии от Земли, а значит, вместо абсолютной звёздной величины использовать видимую, учитывая, что для всех звёзд скопления первая преобразуется во вторую простым прибавлением одинаковой для всех звёзд константы.

Первоначально Герцшпрунг и Рассел в качестве отличительного признака звёзд использовали спектральный класс, что требовало получения и изучения спектров звёзд и делало обработку большого количества данных очень трудозатратной. Однако впоследствии было обнаружено, что вместо спектрального класса можно использовать показатель цвета, который вычисляется как (I1 - I2)/(I1 + I2), где I1 и I2 — интенсивности за соответствующими светофильтрами. Обычно берутся синий и жёлто-зелёный фильтры (обозначение B и V) или зелёный и красный. Для анализа фотографического изображения удобнее всего взять интенсивности в каналах G и R, поскольку при изменении температуры звезды разница в этих каналах меняется сильнее всего.

Вот изображение шарового скопления M 13 в созвездии «Геркулес» — здесь приведён скриншот из программы Siril, которая позволяет открывать, просматривать и обрабатывать FITS-файлы.

Загруженный в программу Siril файл M 13_crop_calibrated.fit

А вот полученная для этого скопления диаграмма Герцшпрунга–Рассела, построенная при помощи пайтоновского скрипта, текст которого приводится далее. На этой диаграмме видна очень хилая главная последовательность, что говорит о почтенном возрасте скопления, и хорошо выраженная отходящая от главной последовательности ветвь красных гигантов, переходящая вверху в ветвь красных сверхгигатнов то есть звёзд, находящихся на последней стадии своей эволюции.

Диаграмма Цветовой индекс — видимая звёздная величина для шарового скопления M 13

Точка, в которой ветвь красных гигантов отходит от главной последовательности, позволяет определить возраст скопления, поскольку чем дальше в красную область мы сдвигаемся по главной последовательности, тем больше время жизни звезды. В шаровых скоплениях не идёт процесс звездообразования, а значит, можно утверждать, что отсутствующая левая часть главной последовательности говорит нам о том, что все находившиеся там звёзды уже проэволюционировали, а все звёзды справа от точки ответвления гигантов — ещё не успели. Современные данные говорят о том, что возраст скопления М 13 порядка 12 миллиардов лет.

Вот код пайтоновского скрипта, при помощи которого построена приведённая выше диаграмма:

import numpy as np
import matplotlib.pyplot as plt
from astropy.io import fits
from photutils.detection import DAOStarFinder
from photutils.aperture import CircularAperture, aperture_photometry

# Load FITS file and check the header
def load_fits_image(filename):
    with fits.open(filename) as hdul:
        data = hdul[0].data # Main image
        header = hdul[0].header
        mag_zero_point = header.get("MAGZP", 30.0) # Get magnitude zero point if available
    return data, mag_zero_point
 
# Detect stars in the image
def find_stars(image, threshold=5):
    mean, std = np.mean(image), np.std(image)
    star_finder = DAOStarFinder(fwhm=3, threshold=threshold * std)
    stars = star_finder(image - mean)
    return stars
 
# Compute integrated intensity for each star
def compute_intensity(image, stars, radius=5):
    positions = np.transpose((stars['xcentroid'], stars['ycentroid']))
    apertures = CircularAperture(positions, r=radius)
    phot_table = aperture_photometry(image, apertures)
    return phot_table['aperture_sum']
 
# Main processing function
def process_fits(filename):
    image, mag_zero_point = load_fits_image(filename)
    if len(image.shape) == 3:
        stars = find_stars(image[1]) # Use the G channel to detect stars
        intensities = []
        for i in range(len(stars)):
            star_intensity = []
            for channel in range(3):
                intensity = compute_intensity(image[channel], stars)[i]
                star_intensity.append(intensity)
            intensities.append(tuple(star_intensity))
        return intensities, mag_zero_point
    else:
        print("Error: The image does not contain RGB color channels.")
        return None, None
 
# Run the program
fits_filename = '/Users/Alexey_Pasechnik/Pictures/Seestar/MyWorks/M 13/M 13_crop_calibrated.fit'
intensity_results, mag_zero_point = process_fits(fits_filename)
 
# Convert RGB to G-R and M_V
G_R = [2 * (g - r) / (g + r) for r, g, b in intensity_results]
M_V = [-2.5 * np.log10(r + g + b) + mag_zero_point for r, g, b in intensity_results]
 
# Define colormap based on G-R index
colors = plt.cm.coolwarm((np.array(G_R) - min(G_R)) / (max(G_R) - min(G_R)))
 
# Plot Hertzsprung-Russell diagram
plt.figure(figsize=(8, 6))
plt.scatter(G_R, M_V, c=colors, marker='*')
plt.gca().invert_yaxis() # Invert Y-axis as brighter stars have lower magnitudes
plt.xlabel("Color Index G-R")
plt.ylabel("Apparent magnitude $m_V$")
plt.title("Hertzsprung-Russell Diagram")
plt.grid()
plt.show()

Comments

[2born.livejournal.com]

Обалдеть! Хоть хватайся и делай из этого лабу для студентов по астрофизике!

[al-pas.livejournal.com]

Ага, уже обсудили с другом из университета Турку это. Десяток девайсов на группу вполне достаточно — студентов можно по 2-3 человека на девайс распределить и выставить на мороз снимать звёзды.

[2born.livejournal.com]

:))) Выставить на мороз — это ключевой момент:))))

[al-pas.livejournal.com]

А что делать? Осенью в Финляндии погоды нет, а в апреле уже белые ночи начинаются.