clustertools.analysis.functions module

Functions to calculate StarCluster properties

Designed to accept StarCluster instance as in put to calculate key parameters

clustertools.analysis.functions.ckin(cluster, mmin=None, mmax=None, nmass=10, rmin=None, rmax=None, vmin=None, vmax=None, emin=None, emax=None, kwmin=0, kwmax=1, npop=None, indx=None, projected=False, **kwargs)

NAME: Find the kinematic concentration parameter ck

• see Bianchini et al. 2018, MNRAS, 475, 96

Parameters:

clusterclass

StarCluster instance

mmin/mmaxfloat

specific mass range

nmass

number of mass bins used to calculate alpha

rmin/rmax

specific radial range

vmin/vmaxfloat

specific velocity range

emin/emaxfloat

specific energy range

kwmin/kwmaxint

specific stellar evolution type range

npopint

population number

indxbool

specific subset of stars

projectedbool

use projected values (default: False)

Returns:

ckfloat

kinematic concentration

Other Parameters:

kwargsstr

key words for plotting

clustertools.analysis.functions.closest_star(cluster, projected=False, argument=False)

Find distance to closest star for each star - uses numba

Parameters:

clusterclass

positions of stars within the StarCluster

projectedbool

use projected values (default: False)

argumentbool

return argument of closest star as well (default: False)

Returns:

minimum_distancefloat

distance to closest star for each star

if argument:

argint

argument of closest star for each star

clustertools.analysis.functions.core_relaxation_time(cluster, coulomb=0.4, projected=False)

Calculate the core relaxation time (Stone & Ostriker 2015) of the cluster

• Stone, N.C. & Ostriker, J.P. 2015, ApJ, 806, 28

Parameters

clusterclass

StarCluster

coulombfloat

Coulomb parameter (default: 0.4)

projectedbool

use projected values (default: False)

methodstr

choose between Stone & Ostriker 2015 and other methods (in development)

Returns

trc (in Myr)

clustertools.analysis.functions.energies(cluster, specific=True, ids=None, projected=False, softening=0.0, full=True, parallel=False, **kwargs)

Calculate kinetic and potential energy of every star Parameters ———- cluster : class

StarCluster instance

specificbool

find specific energies (default: True)

ids: bool or int

if given, find the energues of a subset of stars defined either by an array of star ids, or a boolean array that can be used to slice the cluster. (default: None)

projectedbool

use projected values (default: False)

softeningfloat

Plummer softening length in cluster.units (default: 0.0)

fullbool

calculate distance of full array of stars at once with numbra (default: True)

parallelbool

calculate distances in parallel (default: False)

Returns

kin,potfloat

kinetic and potential energy of every star if the i_d argument is not used. If i_d argument is used, return an arrays with potential and kinetic energy in the same shape of i_d

History

2019 - Written - Webb (UofT) 2022 - Updated with support for multiple ids or an idexing array - Erik Gillis (UofT)

clustertools.analysis.functions.eta_function(cluster, mmin=None, mmax=None, nmass=10, rmin=None, rmax=None, vmin=None, vmax=None, emin=None, emax=None, kwmin=0, kwmax=1, npop=None, indx=None, projected=False, plot=False, meq=False, **kwargs)

NAME: Find power slope of velocity dispersion versus mass

• mass bins are set up so that there are an equal number of stars in each bin

Parameters:

clusterclass

StarCluster instance

mmin/mmaxfloat

specific mass range

nmass

number of mass bins used to calculate alpha

rmin/rmax

specific radial range

vmin/vmaxfloat

specific velocity range

emin/emaxfloat

specific energy range

kwmin/kwmaxint

specific stellar evolution type range

npopint

population number

indxbool

specific subset of stars

projectedbool

use projected values (default: False)

plotbool

plot the mass function

Returns:

m_meanfloat

mean mass in each bin

sigvmfloat

velocity dispersion of stars in each bin

etafloat

power-law slope of (sigvm ~ m^eta)

eetafloat

error in eta

yetafloat

y-intercept of fit to log(sigvm) vs log(m)

eetafloat

error in yeta

Other Parameters:

kwargsstr

key words for plotting

clustertools.analysis.functions.find_centre(cluster, xstart=0.0, ystart=0.0, zstart=0.0, vxstart=0.0, vystart=0.0, vzstart=0.0, indx=None, nsigma=1.0, nsphere=100, density=True, rmin=0.1, rmax=None, nmax=100, method='harfst', nneighbour=6)

Find the cluster’s centre

• The default assumes the cluster’s centre is the centre of density, calculated via the find_centre_of_density function.

• For density=False, the routine first works to identify a sphere of nsphere stars around the centre in which to perform a centre of mass calculation (similar to NBODY6). Stars beyond nsigma standard deviations are removed from the calculation until only nsphere stars remain. This step prevents long tidal tails from affecting the calculation

Parameters:

clusterclass

StarCluster

xstart,ystart,zstartfloat

starting position for centre

vxstart,vystart,vzstart

starting velocity for centre

indxbool

subset of stars to use when finding center

nsigmaint

number of standard deviations to within which to keep stars

nsphereint

number of stars in centre sphere (default:100)

densitybool

use Yohai Meiron’s centre of density calculator instead (default: True)

rminfloat

minimum radius to start looking for stars

rmaxfloat

maxmimum radius of sphere around which to estimate density centre (default: None cluster.units, uses maximum r)

nmaxint

maximum number of iterations to find centre

methodstr

method for finding the centre of density (‘harfst’ (default), ‘casertano’)

if method==’casertano’

nneighbourint

number of neighbours for calculation local densities

Returns:

xc,yc,zc,vxc,vyc,vzc - coordinates of centre of mass

clustertools.analysis.functions.find_centre_of_density(cluster, xstart=0.0, ystart=0.0, zstart=0.0, vxstart=0.0, vystart=0.0, vzstart=0.0, indx=None, nsphere=100, rmin=0.1, rmax=None, nmax=100, method='harfst', nneighbour=6)

Find cluster’s centre of density

• Find cluster’s centre of density:

  if method==’harfst’ use (Harfst, S., Gualandris, A., Merritt, D., et al. 2007, NewA, 12, 357) courtesy of Yohai Meiron if method==’casertano’ use (Casertano, S., Hut, P. 1985, ApJ, 298, 80)

Parameters:

clusterclass

StarCluster

xstart,ystart,zstartfloat

starting position for centre (default: 0,0,0)

vxstart,vystart,vzstartfloat

starting velocity for centre (default: 0,0,0)

indx: bool

subset of stars to perform centre of density calculation on (default: None)

nsphereint

number of stars in centre sphere (default:100)

rminfloat

minimum radius of sphere around which to estimate density centre (default: 0.1 cluster.units)

rmaxfloat

maxmimum radius of sphere around which to estimate density centre (default: None cluster.units, uses maximum r)

nmaxfloat

maximum number of iterations (default:100)

methodstr

method for finding the centre of density (‘harfst’ (default), ‘casertano’)

if method==’casertano’

nneighbourint

number of neighbours for calculation local densities

Returns:

xc,yc,zc,vxc,vyc,vzcfloat

coordinates of centre of mass

clustertools.analysis.functions.find_centre_of_mass(cluster)

Find the centre of mass of the cluster

Parameters:

clusterclass

StarCluster

Returns:

xc,yc,zc,vxc,vyc,vzcfloat

coordinates of centre of mass

clustertools.analysis.functions.half_mass_relaxation_time(cluster, coulomb=0.4, projected=False)

Calculate the half-mass relaxation time (Spitzer 1987) of the cluster - Spitzer, L. 1987, Dynamical evolution of globular clusters - Need to adjust amplitude for different input units

Parameters:

clusterclass

StarCluster

coulombfloat

Coulomb parameter (default: 0.4)

projectedbool

use projected values (default: False)

Returns:

trhfloat

half-mass relaxation time within radius rm (in Myr)

clustertools.analysis.functions.mass_function(cluster, mmin=None, mmax=None, nmass=10, rmin=None, rmax=None, vmin=None, vmax=None, emin=None, emax=None, kwmin=0, kwmax=1, npop=None, indx=None, projected=False, mcorr=None, plot=False, **kwargs)

Find mass function over a given mass range

• mass bins are set up so that there are an equal number of stars in each bin

Parameters:

clusterclass

StarCluster instance

mmin/mmaxfloat

specific mass range

nmass

number of mass bins used to calculate alpha

rmin/rmax

specific radial range

vmin/vmaxfloat

specific velocity range

emin/emaxfloat

specific energy range

kwmin/kwmaxint

specific stellar evolution type range

npopint

population number

indxbool

specific subset of stars

projectedbool

use projected values (default: False)

mcorrfloat

completeness correction for masses

plotbool

plot the mass function

Returns:

m_meanfloat

mean mass in each bin

m_histfloat

number of stars in each bin

dmfloat

dN/dm of each bin

alphafloat

power-law slope of the mass function (dN/dm ~ m^alpha)

ealphafloat

error in alpha

yalphafloat

y-intercept of fit to log(dN/dm) vs log(m)

eyalphafloat

error in yalpha

Other Parameters:

kwargsstr

key words for plotting

clustertools.analysis.functions.meq_function(cluster, mmin=None, mmax=None, nmass=10, rmin=None, rmax=None, vmin=None, vmax=None, emin=None, emax=None, kwmin=0, kwmax=1, npop=1, indx=None, projected=False, plot=False, **kwargs)

NAME: Find meq from velocity dispersion versus mass

• mass bins are set up so that there are an equal number of stars in each bin

• As per Bianchini, P. et al. 2016, MNRAS, 458, 3644, velocity dispersion versus mass is fit with the following: sigma(m)= sigma e^(-1/2 m/meq) if m<= meq

  = sigma0 e^(-1/2) (m/meq)^-1/2 if m > meq

Parameters:

clusterclass

StarCluster instance

mmin/mmaxfloat

specific mass range

nmass

number of mass bins used to calculate alpha

rmin/rmax

specific radial range

vmin/vmaxfloat

specific velocity range

emin/emaxfloat

specific energy range

kwmin/kwmaxint

specific stellar evolution type range

npopint

population number

indxbool

specific subset of stars

projectedbool

use projected values (default: False)

plotbool

plot the mass function

Returns:

m_meanfloat

mean mass in each bin

sigvmfloat

velocity dispersion of stars in each bin

meqfloat

Bianchini fit to sigvm vs m

emeqfloat

error in Bianchini fit to sigvm vs m

sigma0float

Bianchini fit to sigvm vs m

esigma0float

error in Bianchini fit to sigvm vs m

Other Parameters:

kwargsstr

key words for plotting

clustertools.analysis.functions.rcore(cluster, method='casertano', nneighbour=6, mfrac=0.1, projected=False, plot=False, **kwargs)

Calculate core radius of the cluster – The default method (method=’casertano’) follows Casertano, S., Hut, P. 1985, ApJ, 298, 80 to find the core – An alternative metrhod (method==’isothermal’) assumes the cluster is an isothermal sphere the core radius is where density drops to 1/3 central value — For projected core radius, the core radius is where the surface density profile drops to 1/2 the central value — Note that the inner mass fraction of stars used to calculate central density is set by mfrac (default 0.1 = 10%)

Parameters:

clusterclass

StarCluster instance

methodstring

method of calculating the core radius of a star cluster (default ‘casertano’)

if method ==’casertano’:

nneighbourint

number of neighbours for calculation local densities

if method==’isothermal’:

mfracfloat

inner mass fraction to be used to establish the central density

projectedbool

use projected values (default: False)

plotbool

plot the density profile and mark the core radius of the cluster (default: False)

Returns

——-

rcfloat

core radius

Other Parameters:

None

clustertools.analysis.functions.relaxation_time(cluster, rad=None, coulomb=0.4, projected=False, method='spitzer')

Calculate the relaxation time (Spitzer & Hart 1971) within a given radius of the cluster

• Spitzer, L. Jr, Hart, M.H. 1971, ApJ, 164, 399 (Equation 5)

• Need to adjust amplitude for different input units

Parameters

clusterclass

StarCluster

radfloat

radius within which to calculate the relaxation time (defult: cluster.rm)

coulombfloat

Coulomb parameter (default: 0.4)

projectedbool

use projected values (default: False)

methodstr

choose between Spitzer & Hart 1971 and other methods (in development)

Returns:

trelaxfloat

relaxation time within radius rad

clustertools.analysis.functions.rlagrange(cluster, nlagrange=10, mfrac=None, projected=False)

Calculate lagrange radii of the cluster by mass

Parameters:

clusterclass

StarCluster

nlagrangeint

number of lagrange radii bins (default: 10)

mfracfloat

Exact masss fraction to calculate radius. Will supersede nlagrange if not None (default : None)

projectedbool

calculate projected lagrange radii (default: False)

Returns:

rnfloat

lagrange radii

clustertools.analysis.functions.rlimiting(cluster, pot=None, rgc=None, zgc=None, nrad=20, projected=False, plot=False, from_centre=False, **kwargs)

Calculate limiting radius of the cluster

• The limiting radius is defined to be where the cluster’s density reaches the local background density of the host galaxy

• for cases where the cluster’s orbital parameters are not set, it is possible to manually set rgc which is assumed to be in kpc.

Parameters:

clusterclass

StarCluster

potclass

GALPY potential used to calculate actions

rgc

Manually set galactocentric distance in kpc at which the tidal radius is to be evaluated (default: None)

zgcfloat

For non-spherically symmetric potentials, manually set distance in kpc above disk at which the tidal radius is to be evaluated. When set, rgc becomes radius in cylindrical coordinates (default: None)

nradint

number of radial bins used to calculate density profile (Default: 20)

projectedbool

use projected values (default: False)

plotbool

plot the density profile and mark the limiting radius of the cluster (default: False)

from_centrebool

calculate tidal radius based on location of cluster’s exact centre instead of its assigned galactocentric coordinates (default: False)

verbosebool

Returns

——-

rlfloat

limiting radius

Other Parameters:

kwargsstr

key words for plotting

clustertools.analysis.functions.rtidal(cluster, pot=None, rtiterate=0, rtconverge=0.9, indx=None, rgc=None, zgc=None, from_centre=False, plot=False, verbose=False, **kwargs)

Calculate tidal radius of the cluster - The calculation uses Galpy (Bovy 2015_, which takes the formalism of Bertin & Varri 2008 to calculate the tidal radius – Bertin, G. & Varri, A.L. 2008, ApJ, 689, 1005 – Bovy J., 2015, ApJS, 216, 29 - riterate = 0 corresponds to a single calculation of the tidal radius based on the cluster’s mass (np.sum(cluster.m)) – Additional iterations take the mass within the previous iteration’s calculation of the tidal radius and calculates the tidal

radius again using the new mass until the change is less than 90%

• for cases where the cluster’s orbital parameters are not set, it is possible to manually set rgc which is assumed to be in kpc.

Parameters:

clusterclass

StarCluster instance

potclass

GALPY potential used to calculate tidal radius (default: None)

rtiterateint

how many times to iterate on the calculation of r_t (default: 0)

rtconvergefloat

criteria for tidal radius convergence within iterations (default 0.9)

indxbool

subset of stars to use when calculate the tidal radius (default: None)

rgcfloat

Manually set galactocentric distance in kpc at which the tidal radius is to be evaluated (default: None)

zgcfloat

For non-spherically symmetric potentials, manually set distance in kpc above disk at which the tidal radius is to be evaluated. When set, rgc becomes radius in cylindrical coordinates (default: None)

rofloat

GALPY radius scaling parameter

vofloat

GALPY velocity scaling parameter

from_centrebool

calculate tidal radius based on location of cluster’s exact centre instead of its assigned galactocentric coordinates (default: False)

plotbool

plot the x and y coordinates of stars and mark the tidal radius of the cluster (default: False)

verbosebool

Print information about iterative calculation of rt

Returns:

rtfloat

tidal radius

Other Parameters:

kwargsstr

key words for plotting

clustertools.analysis.functions.tapered_mass_function(cluster, mmin=None, mmax=None, nmass=10, rmin=None, rmax=None, vmin=None, vmax=None, emin=None, emax=None, kwmin=0, kwmax=1, npop=None, indx=None, projected=False, mcorr=None, plot=False, **kwargs)

Find a tapered mass function over a given mass range

• mass bins are set up so that there are an equal number of stars in each bin

• functional form of the tapered mass function is taken from De Marchi, Paresce & Portegies Zwart 2010

Parameters

clusterclass

StarCluster instance

mmin/mmaxfloat

specific mass range

nmass

number of mass bins used to calculate alpha

rmin/rmax

specific radial range

vmin/vmaxfloat

specific velocity range

emin/emaxfloat

specific energy range

kwmin/kwmaxint

specific stellar evolution type range

npopint

population number

indxbool

specific subset of stars

projectedbool

use projected values (default: False)

mcorrfloat

completeness correction for masses

plotbool

plot the mass function

Returns:

m_meanfloat

mean mass in each bin

m_histfloat

number of stars in each bin

dmfloat

dN/dm of each bin

alphafloat

power-law slope of the mass function (dN/dm ~ m^alpha)

ealphafloat

error in alpha

yalphafloat

y-intercept of fit to log(dN/dm) vs log(m)

eyalphafloat

error in yalpha

Other Parameters:

kwargsstr

key words for plotting

clustertools.analysis.functions.virial_radius(cluster, method='inverse_distance', full=True, H=70.0, Om=0.3, overdens=200.0, projected=False, plot=False, **kwargs)

Calculate virial radius of the cluster - Virial radius is calculated using either: – the average inverse distance between particles, weighted by their masses (default) – the radius at which the density is equal to the critical density of the Universe at the redshift of the system, multiplied by an overdensity constant Parameters ———- cluster : class

StarCluster

methodstr

method for calculating virial radius (default: ‘inverse_distance’, alternative: ‘critical_density’)

Returns:

rvfloat

virial radius

Other Parameters:

fullbool

Use Numba to calculate average inverse distance between stars (default:True)

Hfloat

Hubble constant

Omfloat

density of matter

overdensfloat

overdensity constant

projectedbool

calculate projected virial radius (default: False)

plotbool

plot cluster density profile and illustrate virial radius calculation

kwargsstr

key word arguments for plotting function

clustertools.analysis.functions.virial_radius_critical_density(cluster, H=70.0, Om=0.3, overdens=200.0, projected=False, plot=False, **kwargs)

Calculate virial radius of the cluster

• Virial radius is defined as the radius at which the density is equal to the critical density of the Universe at the redshift of the system, multiplied by an overdensity constant

• Note that this a quick method that is a bit of an approximation as it interpolates the cluster’s density profile. A more accurate (but expensive)

approach would be to subtract the product of the critical density and the overdensity constant from the density profile and find the root (in development)

Parameters:

clusterclass

StarCluster

Hfloat

Hubble constant

Omfloat

density of matter

overdensfloat

overdensity constant

projectedbool

calculate projected virial radius (default: False)

plotbool

plot cluster density profile and illustrate virial radius calculation

Returns:

r_vfloat

virial radius

Other Parameters:

kwargsstr

key word arguments for plotting function

clustertools.analysis.functions.virial_radius_inverse_distance(cluster, projected=False, full=True)

Calculate virial radius of the cluster - Virial radius is defined as the inverse of the average inverse distance between particles, weighted by their masses – Definition taken from AMUSE (www.amusecode.org) – Portegies Zwart S., McMillan S., 2018, Astrophysical Recipes; The art ofAMUSE, doi:10.1088/978-0-7503-1320-9

Parameters:

clusterclass

StarCluster

projectedbool

calculate projected virial radius (default: False)

fullbool

Use Numba to calculate average inverse distance between stars (default:True)

Returns:

r_vfloat

virial radius