Queremos brindarte la mejor solución que descubrimos en internet. Nosotros deseamos que te resulte de mucha utilidad y si quieres aportar cualquier detalle que nos pueda ayudar a mejorar hazlo con libertad.
Solución:
Primero tengo que decir, ¡gran pregunta! Muy detallada y reproducible. Revisé su pregunta e intenté rehacer el ejercicio comenzando desde su repositorio de git y descargando el catálogo del archivo GAIA.
EDITAR
Programáticamente, su código está bien (ver PARTE VIEJA a continuación para un enfoque ligeramente diferente). El problema con los puntos que faltan es que solo se obtienen 500 puntos de datos al descargar el archivo csv del archivo GAIA. Por lo tanto, parece como si todos los puntos de la consulta estuvieran abarrotados en una forma extraña. Sin embargo, si restringe el radio de la búsqueda a un valor más pequeño, puede ver que hay puntos que se encuentran dentro de la imagen TESS:
compare con la versión que se muestra a continuación en la PARTE ANTIGUA. El código es el mismo que se muestra a continuación, solo que el archivo csv descargado es para un radio más pequeño. Por lo tanto, parece que acaba de descargar una parte de todos los datos disponibles del archivo GAIA al exportar a csv. La forma de eludir esto es hacer la búsqueda como lo hizo. Luego, en la página de resultados, haga clic en Show query in ADQL form en la parte inferior y en la consulta se muestra en cambio de formato SQL:
Select Top 500
para
Select
al principio de la consulta.
PARTE ANTIGUA (el código está bien y funciona, pero mi conclusión es incorrecta):
Para graficar usé aplpy – usa matplotlib en segundo plano – y terminó con el siguiente código:
from astropy.io import fits
from astropy.wcs import WCS
import aplpy
import matplotlib.pyplot as plt
import pandas as pd
from astropy.coordinates import SkyCoord
import astropy.units as u
from astropy.io import fits
fits_file = fits.open("4687500098271761792_med.fits")
central_coordinate = SkyCoord(fits_file[0].header["CRVAL1"],
fits_file[0].header["CRVAL2"], unit="deg")
figure = plt.figure(figsize=(10, 10))
fig = aplpy.FITSFigure("4687500098271761792_med.fits", figure=figure)
cmap = "gist_heat"
stretch = "log"
fig.show_colorscale(cmap=cmap, stretch=stretch)
fig.show_colorbar()
df = pd.read_csv("4687500098271761792_within_1000arcsec.csv")
# the epoch found in the dataset is J2015.5
df['coord'] = SkyCoord(df["ra"], df["dec"], unit="deg", frame="icrs",
equinox="J2015.5")
coords = df["coord"].tolist()
coords_degrees = [[coord.ra.degree, coord.dec.value] for coord in df["coord"]]
ra_values = [coord[0] for coord in coords_degrees]
dec_values = [coord[1] for coord in coords_degrees]
width = (40*u.arcmin).to(u.degree).value
height = (40*u.arcmin).to(u.degree).value
fig.recenter(x=central_coordinate.ra.degree, y=central_coordinate.dec.degree,
width=width, height=height)
fig.show_markers(central_coordinate.ra.degree,central_coordinate.dec.degree,
marker="o", c="white", s=15, lw=1)
fig.show_markers(ra_values, dec_values, marker="o", c="blue", s=15, lw=1)
fig.show_circles(central_coordinate.ra.degree,central_coordinate.dec.degree,
radius=(1000*u.arcsec).to(u.degree).value, edgecolor="black")
fig.save("GAIA_TESS_test.png")
Sin embargo, esto da como resultado una trama similar a la tuya:
Para comprobar mi sospecha de que las coordenadas del archivo GAIA se muestran correctamente, dibujo un círculo de 1000 segundos de arco desde el centro de la imagen TESS. Como puede ver, se alinea aproximadamente con la forma circular del lado exterior (visto desde el centro de la imagen) de la nube de puntos de datos de las posiciones de GAIA. Simplemente creo que estos son todos los puntos en el archivo GAIA DR2 que se encuentran dentro de la región que buscó. La nube de datos incluso parece tener un límite cuadrado en el interior, que podría provenir de algo así como un campo de visión cuadrado.
Muy buen ejemplo. Solo para mencionar que también puede integrar la consulta al archivo Gaia usando el módulo astroquery.gaia incluido en astropy
https://astroquery.readthedocs.io/en/latest/gaia/gaia.html
De esta manera, podrá ejecutar las mismas consultas que están dentro de la interfaz de usuario del archivo Gaia y cambiar a diferentes fuentes de una manera más fácil.
from astroquery.simbad import Simbad
import astropy.units as u
from astropy.coordinates import SkyCoord
from astroquery.gaia import Gaia
result_table = Simbad.query_object("Gaia DR2 4687500098271761792")
raValue = result_table['RA']
decValue = result_table['DEC']
coord = SkyCoord(ra=raValue, dec=decValue, unit=(u.hour, u.degree), frame='icrs')
query = """SELECT TOP 1000 * FROM gaiadr2.gaia_source
WHERE CONTAINS(POINT('ICRS',gaiadr2.gaia_source.ra,gaiadr2.gaia_source.dec),
CIRCLE('ICRS',ra,dec,0.2777777777777778))=1 ORDER BY random_index""".format(ra=str(coord.ra.deg[0]),dec=str(coord.dec.deg[0]))
job = Gaia.launch_job_async(query)
r = job.get_results()
ralist = r['ra'].tolist()
declist = r['dec'].tolist()
import matplotlib.pyplot as plt
plt.scatter(ralist,declist,marker='+')
plt.show()
Tenga en cuenta que he agregado el orden por random_index que eliminará este extraño comportamiento no circular. Este índice es bastante útil para no forzar la salida completa para las pruebas iniciales.
Además, he declarado las coordenadas de salida para la ascensión recta de Simbad como horas.
Por último, he utilizado la consulta asíncrona que tiene menos limitaciones en tiempo de ejecución y máximo filas en la respuesta.
También puede cambiar la consulta a
query = """SELECT * FROM gaiadr2.gaia_source
WHERE CONTAINS(POINT('ICRS',gaiadr2.gaia_source.ra,gaiadr2.gaia_source.dec),
CIRCLE('ICRS',ra,dec,0.2777777777777778))=1""".format(ra=str(coord.ra.deg[0]),dec=str(coord.dec.deg[0]))
(eliminando la limitación a 1000 filas) (en este caso, no es necesario el uso del índice aleatorio) para tener una respuesta completa del servidor.
Por supuesto, esta consulta tarda un tiempo en ejecutarse (alrededor de 1,5 minutos). La consulta completa devolverá 103574 filas.
