Models#

The GCfit Models can be accessed through the core module at

>>> import gcfit

All models are based off of the single base class gcfit.core.Model.

To begin, we can start by exploring a model with some arbitrary default parameters. The gcfit.core.Model class gives default values for many arguments, which you may want to adjust yourself. See the documentation of said class for more explanation of the meaning of all available parameters.

>>> model = gcfit.Model(W0=6.3, M=5e5, rh=6.7, age=12, FeH=-0.7)
>>> model
<gcfit.core.data.Model object at 0x7f43609cb9a0

The model will automatically generate a number of mass bins, containing either stars or remnants of a certain type, which are used to solve the multimass version of the LIMEPY DF.

>>> # Mean masses per bin
>>> model.mj
<Quantity [ 0.11766481,  0.16234563,  0.22399308, ...,  7.10114089,
            8.94490091, 11.4071349] solMass>

>>> # Total mass per bin
>>> model.Mj
<Quantity [2.97895583e+04, 3.73135159e+04, 4.67358943e+04, ...,
           2.98231834e+01, 5.58230429e+01, 5.05028550e+01] solMass>

>>> # Stellar object types (MS, NS, WD, BH)
>>> model.star_types
array(['MS', 'MS', 'MS', 'MS', 'MS', 'MS', 'MS', 'MS', 'MS', 'WD', 'WD',
       'WD', 'WD', 'WD', 'WD', 'NS', 'BH', 'BH', 'BH', 'BH'], dtype='<U2')

>>> # Total mass and number of black holes, in their repective bins
>>> model.BH.Mj
<Quantity [12.26536265, 29.82318344, 55.82304292, 50.50285502] solMass>
>>> model.BH.Nj
<Quantity [2.07936342, 4.19977352, 6.24076705, 4.42730409]>

Notice that the majority of interesting quantities in gcfit.core.Model are stored as astropy.Quantity objects, with their respective units.

The radial profiles of a number of system properties, such as velocity dispersion, density and energy, are available for each mass bin, as well as a number of useful radii.

>>> # Radial profile domain
>>> model.r
<Quantity [0.00000000e+00, 2.08887713e-07, 2.29776485e-06, ...,
           6.00907685e+01, 6.00908280e+01, 6.00908306e+01] pc>

>>> # Density profile of the most massive main-sequence stars
>>> model.rhoj[model.nms - 1]
<Quantity [1.22615737e+05, 1.22615733e+05, 1.22615270e+05, ...,
           4.10890596e-21, 1.21353818e-25, 0.00000000e+0] solMass / pc3>

>>> # Half-mass radius of each mass bin
>>> model.rhj
<Quantity [9.16723465, 8.96575389, 8.66929104, ..., 0.0351056, 0.0252477,
           0.01440394] pc>

See gcfit.core.Model for further description of all available properties.

Models matching a number of historical DF formulations can also be created easily using the relevant generator functions. These functions mostly consist of setting a specific default value for the truncation parameter g.

>>> # Generate a King (1966) model
>>> king = gcfit.Model.king(6.3, 5e5, 6.7, age=12, FeH=-0.7)

>>> model.g, king.g
(1.5, 1)

Sampled Models#

These (multimass) models can also be sampled, in order to return a random distribution of stars matching the phase-space distribution of the models.

>>> sampled = model.sample()
>>> sampled
<gcfit.core.data.SampledModel object at 0x7f4360a3fa30>

>>> # Total number of stars in the system
>>> sampled.Nstars
<Quantity 1244880>

>>> # Cartesian coordinates of all stars, centred on the cluster centre
>>> sampled.pos.x
<Quantity [ 1.70500994e+01,  2.86122226e+00, -1.16756636e+00, ...,
            5.24831562e-03, -1.45917202e-03, -9.13248353e-03] pc>
>>> sampled.pos.z
<Quantity [-9.56794358e+00, -2.53661084e+00, -6.16515774e+00, ...,
            6.71300176e-03,  2.06363703e-03,  1.76499660e-03] pc>
>>> sampled.pos._fields
('x', 'y', 'z', 'r', 'theta', 'phi')

>>> # Radial and tangential velocities of each star
>>> sampled.vel.r
<Quantity [-0.71188986, -6.86505168,  3.51273965, ..., -3.84003514,
            3.12497992,  2.09300044] km / s>
>>> sampled.vel.t
<Quantity [6.83229216, 4.14966349, 4.9182554 , ..., 3.80989784, 4.72867768,
           1.41100603] km / s>
>>> sampled.vel._fields
('x', 'y', 'z', 'r', 't', 'theta', 'phi')

If a centre coordinate on the sky is given (as an astropy.SkyCoord with both position and velocity), the projected positions and velocities on the sky can also be computed.

>>> import astropy.units as u
>>> from astropy.coordinates import SkyCoord

>>> deg, masyr, kms = u.deg, u.Unit('mas/yr'), u.Unit('km/s')
>>> cen = SkyCoord(l=45. * deg, b=55. * deg,
>>>                pm_l_cosb=5 * masyr, pm_b=3 * masyr, radial_velocity=2 * kms,
>>>                frame='galactic')

>>> p_sampled = model.sample(centre=cen)

>>> p_sampled.galactic.lon
<Quantity [44.99073051, 45.27478485, 44.95426386, ..., 45.0009292 ,
           44.99996491, 45.00047616] deg>
>>> p_sampled.galactic.pm_b
<Quantity [2.59835011, 3.18921132, 3.36723663, ..., 3.06262506, 3.0400891 ,
           2.96432599] mas / yr>
>>> p_sampled.galactic._fields
('lat', 'lon', 'distance', 'pm_l_cosb', 'pm_b', 'v_los')

Observations#

Another useful class within GCfit is the gcfit.core.Observations class, which acts as a container for a number of observational datasets. These observations are key for all fitting (see below), but are also useful when working with individual models, as they contain a number of useful metadata fields about the cluster:

>>> obs = gcfit.Observations('NGC104')
Observations(cluster="NGC0104")

>>> model = gcfit.Model(W0=6.3, M=5e5, rh=6.7, observations=obs)
>>> model.age, model.FeH
(<Quantity 11.75 Gyr>, -0.72)

More information on the datafiles underlying this class, and how to create your own datafiles can be found at (TODO).

Models#

Sampled Models#

Observations#

This Page