Короче говоря: я хочу представить данные астрономии Гайи в изображениях TESS в Python. Как это возможно?
Ниже приведена детальная версия.
У меня 64x64 пикселя TESS изображения звезды с Gaia ID 4687500098271761792 . На странице 8 Руководства обсерватории TESS говорится, что 1 пиксель составляет ~ 21 угловую секунду. Используя Gaia Archive , я ищу эту звезду (ниже главные функции , нажимаю Поиск .) И отправляю запрос, чтобы увидеть звезды в пределах 1000 угловых секунд, примерно радиус нам нужен. Имя, которое я использую для поиска: Gaia DR2 4687500098271761792, как показано ниже:
Отправить запрос, и я получаю список из 500 звезд с RA и DEC координатами. Выберите CSV и Download results, я получу список звезд вокруг 4687500098271761792 . Этот результирующий файл также можно найти здесь . Это входные данные Gaia , которые мы хотим использовать.
Из TESS у нас есть 4687500098271761792_med.fits , файл изображения. Мы строим это, используя:
from astropy.io import fits
from astropy.wcs import WCS
import matplotlib.pyplot as plt
hdul = fits.open("4687500098271761792_med.fits")[0]
wcs = WCS(hdul.header)
fig = plt.figure(figsize=(12,12))
fig.add_subplot(111, projection=wcs)
plt.imshow(hdul.data)
и получите хорошую картинку:
и кучу предупреждений, большая часть которых была любезно объяснена здесь (предупреждения в Q, объяснение в комментариях) .
Обратите внимание, что действительно хорошо, что мы используем проекцию WCS . Чтобы проверить, давайте просто нанесем данные в hdul.data, не заботясь о проекции:
plt.imshow(hdul.data)
Результат:
Почти так же, как и раньше, но теперь метки осей - это просто числа пикселей, а не RA и DEC , как было бы предпочтительно. Значения DEC и RA на первом графике составляют около -72 ° и 16 ° соответственно, что хорошо, учитывая, что каталог Gaia дал нам звезды в окрестности 4687500098271761792 примерно с этими координатами , Таким образом, прогноз, кажется, достаточно хорошо.
Теперь давайте попробуем нанести звезды Гайи на график imshow(). Мы читаем в файле CSV, который мы скачали ранее, и извлекаем из него значения RA и DEC объектов:
import pandas as pd
df=pd.read_csv("4687500098271761792_within_1000arcsec.csv")
ralist=df['ra'].tolist()
declist=df['dec'].tolist()
Участок для проверки:
plt.scatter(ralist,declist,marker='+')
Форма не является кругом, как ожидалось. Это может быть индикатором будущих неприятностей.
Давайте попробуем преобразовать эти значения RA и DEC в WCS и построить их таким образом:
for index, each in enumerate(ralist):
ra, dec = wcs.all_world2pix([each], [declist[index]], 1)
plt.scatter(ra, dec, marker='+', c='k')
Результат:
Функция all_world2pix от здесь . Параметр 1 просто устанавливает, где находится источник. Описание all_world2pix гласит:
Здесь начало координат - это координата в верхнем левом углу изображения.
В стандартах FITS и Fortran это 1. В стандартах Numpy и C
это 0.
Тем не менее, форма распределения точек, которую мы получаем, вовсе не обещает. Давайте соединим данные TESS и Gaia:
hdul = fits.open("4687500098271761792_med.fits")[0]
wcs = WCS(hdul.header)
fig = plt.figure(figsize=(12,12))
fig.add_subplot(111, projection=wcs)
plt.imshow(hdul.data)
for index, each in enumerate(ralist):
ra, dec = wcs.all_world2pix([each], [declist[index]], 1)
plt.scatter(ra, dec, marker='+', c='k')
Получаем:
, которого нет рядом с тем, что было бы идеально. Я ожидаю, что на нем будет картинка imshow() с множеством маркеров, а маркеры должны находиться там, где звезды находятся на изображении TESS. Блокнот Jupyter, в котором я работал, доступен здесь .
Какие шаги я пропускаю или что я делаю не так?
Дальнейшие разработки
В ответе на другой вопрос , keflavich любезно предложено использовать аргумент transform для построения в мировых координатах . Пробовал это с некоторыми примерами точек (изогнутый крест на графике ниже). Также нанесены данные Gaia на график без их обработки, в итоге они оказались сосредоточены в очень узком пространстве. Применив к ним метод transform, получил, казалось бы, очень похожий результат, чем раньше. Код (а также здесь ):
import pandas as pd
df=pd.read_csv("4687500098271761792_within_1000arcsec.csv")
ralist=df['ra'].tolist()
declist=df['dec'].tolist()
from astropy.io import fits
from astropy.wcs import WCS
import matplotlib.pyplot as plt
hdul = fits.open("4687500098271761792_med.fits")[0]
wcs = WCS(hdul.header)
fig = plt.figure(figsize=(12,12))
fig.add_subplot(111, projection=wcs)
plt.imshow(hdul.data)
ax = fig.gca()
ax.scatter([16], [-72], transform=ax.get_transform('world'))
ax.scatter([16], [-72.2], transform=ax.get_transform('world'))
ax.scatter([16], [-72.4], transform=ax.get_transform('world'))
ax.scatter([16], [-72.6], transform=ax.get_transform('world'))
ax.scatter([16], [-72.8], transform=ax.get_transform('world'))
ax.scatter([16], [-73], transform=ax.get_transform('world'))
ax.scatter([15], [-72.5], transform=ax.get_transform('world'))
ax.scatter([15.4], [-72.5], transform=ax.get_transform('world'))
ax.scatter([15.8], [-72.5], transform=ax.get_transform('world'))
ax.scatter([16.2], [-72.5], transform=ax.get_transform('world'))
ax.scatter([16.6], [-72.5], transform=ax.get_transform('world'))
ax.scatter([17], [-72.5], transform=ax.get_transform('world'))
for index, each in enumerate(ralist):
ax.scatter([each], [declist[index]], transform=ax.get_transform('world'),c='k',marker='+')
for index, each in enumerate(ralist):
ax.scatter([each], [declist[index]],c='b',marker='+')
и полученный сюжет:
Ожидается этот изгибный крест, так как TESS не выровнен с линиями постоянной широты и долготы (т.е. плечи креста не должны быть параллельны сторонам изображения TESS, нанесенного с помощью imshow()). Теперь давайте попробуем построить постоянные линии RA и DEC (или, скажем, постоянные линии широты и долготы), чтобы лучше понять, почему точки данных из Gaia смещены. Разверните код выше несколькими строками:
ax.coords.grid(True, color='green', ls='solid')
overlay = ax.get_coords_overlay('icrs')
overlay.grid(color='red', ls='dotted')
Результат обнадеживает:
(см. Тетрадь здесь .)