htof

This repo contains htof, the package for parsing intermediate data from the Gaia and Hipparcos satellites, and reproducing five, seven, and nine (or higher) parameter fits to their astrometry.

.. image:: https://coveralls.io/repos/github/gmbrandt/HTOF/badge.svg?branch=main :target: https://coveralls.io/github/gmbrandt/HTOF?branch=main

.. image:: https://app.travis-ci.com/gmbrandt/HTOF.svg?branch=main :target: https://app.travis-ci.com/gmbrandt/HTOF

Parallax is handled by the :code:sky_path module which was written by Anthony Brown as a part of his astrometric-sky-path package: https://github.com/agabrown/astrometric-sky-path/

Installation

htof can be installed with PyPi and pip, by running

.. code-block:: bash

pip install htof

or by running

.. code-block:: bash

pip install .

while in the root directory of this repo. It can also be installed directly from github using

.. code-block:: bash

pip install git+https://github.com/gmbrandt/htof

Usage

HTOF has a rich variety of usages. We encourage the reader to consult the examples.ipynb jupyter notebook for a set of usage examples (e.g., fitting the standard astrometric model to data, combining astrometric missions). However, we also go into a few basic and specific use cases in this readme. Also see examples_recalibrating_hip2 and GenerateSyntheticGaiaAstrometry for more uses of htof.

If you use HTOF, please cite the zenodo reference (https://doi.org/10.5281/zenodo.4104383) and the source paper (https://ui.adsabs.harvard.edu/abs/2021AJ....162..230B/abstract)

Usage: Fits without Parallax

The following examples show how one would both load in and fit a line to the intermediate astrometric data (IAD) from either Hipparcos data reduction or Gaia.

HTOF will download any missing IAD files for you (for Hip1, Hip2 (Java tool data), and Gaia). You should provide a valid directory though for htof to save the file into for future use. Currently, for the automatic download to work, you must provide a hipparcos name for the source (e.g., 27321). Include ONLY the numeric part of the name. For an example of downloading Gaia data, see Example 6 in examples.ipynb

For Gaia:

If the automatic download of the GOST scanning law does not work, or the source does not have a hipparcos ID. then you will have to download the GOST file manually with the user interface. Download a .csv of the predicted scans and scan epochs from GOST (https://gaia.esac.esa.int/gost/). In particular, using the 'submit for events forecast' feature on the website. One should select the widest range of dates possible because htof automatically restricts the predicted epochs of observations to the desired data release range (e.g., EDR3) and removes any astrometric gaps.

Let ra_vs_epoch, dec_vs_epoch be 1d arrays of ra and dec positions. Assume we want to fit to data from GaiaDR2 on the star with hip id 027321. The following lines parse the intermediate data and fit a line.

.. code-block:: python

from htof.main import Astrometry
import numpy as np
astro = Astrometry('GaiaDR2', '027321', 'htof/test/data_for_tests/GaiaDR2/IntermediateData', format='jyear')  # parse
# note that if you do not have a GOST csv file with 027321 in the name, inside of
# 'htof/test/data_for_tests/GaiaDR2/IntermediateData' , then htof will download it for you automatically!
ra_vs_epoch = dec_vs_epoch = np.zeros(len(astro.data), dtype=float) # dummy set of ra and dec to fit.
ra0, dec0, mu_ra, mu_dec = astro.fit(ra_vs_epoch, dec_vs_epoch)

:code:GaiaeDR3 will select all data corresponding to the eDR3 data interval and exclude eDR3 deadtimes. :code:GaiaDR2 will select all data corresponding to the DR2 data interval (excluding dead times). Finally, :code:Gaia will select all the data present in the GOST predicted observation file.

Note, selecting :code:Gaia and allowing HTOF to download the GOST scanning law may result in a wierd coverage of the scanning law. We recommend sticking with defined data release intervals, so :code:GaiaeDR3, etc.

In the :code:astro.fit part, ra_vs_epoch and dec_vs_epoch are the positions in right ascension and declination of the object. These arrays must have the same shape as astro.data.julian_day_epoch(), which are the epochs in the intermediate data. :code:format='jyear' specifies the time units of the output best fit parameters. The possible choices of format are the same as the choices for format in astropy.time.Time(val, format=format). If :code:format='jyear', then the output :code:mu_ra would have units of mas/year. If :code:jd then the output is mas/day. Both Hipparcos and Gaia catalogs list parallaxes in milli-arcseconds (mas), and so positional units are always in mas for HTOF.

For Hipparcos :

:code:Hip2 refers to the DVD IAD which is now obsolete. :code:Hip21 refers to the Java Tool Intermediate Astrometric Data (IAD) and best fit parameters. This is the preferred set of data to use with the 2007 re-reduction (preferred over the DVD IAD). The Hipparcos Java Tool data parser is meant for the 2014 Java tool data (Java tool first released at https://www.cosmos.esa.int/web/hipparcos/java-tools/intermediate-data, in 2014). As of 2021, there has not been an update to the Java tool.

Remember, Hip21 will be automatically downloaded per source. But if you want to avoid downloading the data on the fly, then the full Java Tool Intermediate Astrometric Data can be downloaded from https://www.cosmos.esa.int/web/hipparcos/hipparcos-2 and extracted (ignore the _MACOSX folder if there is one). One would then point any HTOF parser to the ResRec_JavaTool folder that contains the H00 etc. subfolders of the individual IAD files. So:

.. code-block:: python

from htof.main import Astrometry
astro = Astrometry('Hip21', star_id='027321', '/home/user/Downloads/ResRec_JavaTool_2014/ResRec_JavaTool_2014', format='jyear')  # parse
ra0, dec0, mu_ra, mu_dec = astro.fit(ra_vs_epoch, dec_vs_epoch)

We discuss enabling fits with parallax later. By default, the fit is a four-parameter fit: it returns the parameters to the line of best fit to the sky path ra_vs_epoch, dec_vs_epoch. If you want a 6 parameter or 8 parameter fit, specify fit_degree = 2 or fit_degree = 3 respectively. E.g.

.. code-block:: python

from htof.main import Astrometry
astro = Astrometry('GaiaDR2', '027321', 'htof/test/data_for_tests/GaiaDR2/IntermediateData', format='jyear',
                   fit_degree=2)
ra0, dec0, mu_ra, mu_dec, acc_ra, acc_dec = astro.fit(ra_vs_epoch, dec_vs_epoch)

If fit_degree = 3, then the additional last two parameters would be the jerk in right ascension and declination, respectively. The sky path in RA (for instance) should be reconstructed by ra0 + mu_ra*t + 1/2*acc_ra*t**2 where t are the epochs from astro.fitter.epoch_times minus the central epoch for RA (if provided).

HTOF allows fits of arbitrarily high degree. E.g. setting fit_degree=3 would give a 9 parameter fit (if using parallax as well). One should specify a central epoch for the fit, typically choosing the central epoch from the catalog (e.g. 2015.5 for GaiaDR2, 2016 for GaiaEDR3, 1991.25 for Hipparcos). You can specify the central epoch by:

.. code-block:: python

from htof.main import Astrometry

astro = Astrometry('GaiaDR2', '027321', 'htof/test/data_for_tests/GaiaDR2/IntermediateData',
                   central_epoch_ra=2015.5, central_epoch_dec=2015.5, format='jyear')
ra0, dec0, mu_ra, mu_dec = astro.fit(ra_vs_epoch, dec_vs_epoch)

The format of the central epochs must be specified along with the central epochs. The best fit sky path in right ascension would then be :code:ra0 + mu_ra * (epochs - centra_epoch_ra). The central epoch matters for numerical stability and covariances. E.g., dont choose a central epoch like the year 1200 for GaiaDR2. One should almost always choose the central epoch from the catalog.

One can access the BJD epochs with

.. code-block:: python

astro.data.julian_day_epoch()

If you want the standard (1-sigma) errors on the parameters, set :code:return_all=True when fitting:

.. code-block:: python

from htof.main import Astrometry

astro = Astrometry('GaiaDR2', '027321', 'htof/test/data_for_tests/GaiaDR2/IntermediateData',
                    central_epoch_ra=2015.5, central_epoch_dec=2015.5, format='jyear')
solution_vector, errors, chisq, residuals = astro.fit(ra_vs_epoch, dec_vs_epoch, return_all=True)

errors is an array the same shape as solution_vector, where each entry is the 1-sigma error for the parameter at the same location in the solution_vector array. chisq is the formal chisquared of the fit to the data, and residuals are the data - model residuals (given as a Nx2 shaped array, where N is the number of transits. The first column are the ra residuals and the second are the declination residuals). For a simple refit to the catalog IAD, residuals (converted to the AL basis) will equal (up to round off) the residuals given in the IAD. One could convert the residuals to the along scan basis by doing:

.. code-block:: python

from htof.special_parse import to_along_scan_basis
residuals = to_along_scan_basis(ra_decresiduals[:, 0], ra_decresiduals[:, 1], astro.data.scan_angle.values)
# now residuals will be a one dimensional array of length N (number of transits), giving the residuals along the
# scan.

For Hip1 and Hip21, HTOF loads in the real catalog errors and so the parameter error estimates (errors) should match those given in the catalog. For Hip2, the along scan errors are automatically inflated or deflated in accordance with D. Michalik et al. 2014. For Gaia we do not have the error estimates from the GOST tool. The AL errors are set to 1 mas by default and so the best-fit parameter errors to Gaia will not match those reported by the catalog.

`chi