library(SimSpin)The purpose of this package is to measure the observable spin parameter \(\lambda_R\) for simulated galaxy models with an effort to provide a method of communication between observers and computational astrophysicists. The aim is to address a simple question: is the true spin parameter, \(\lambda\), and the one that we measure observationally, \(\lambda_R\), comparable?
This notebook will walk through a few of the basic steps that we can use to analyse a galaxy model. Here, we will use the SimSpin_example.hdf5 file included within the SimSpin package to explore the three basic functions available in this package. This is the Gadget snapshot for a \(10^{10} M_{\odot}\), S0 galaxy containing 50,000 disk and 75,000 bulge particles of equal mass. We read in this simulation file using the sim_data() function shown below. This function takes a SimSpin HDF5 file and writes it out as an R list of data frames that are easy to manipulate within the rest of the package. At this stage, you can specify the type of particle we want to analyse ptype = NA (all), 0 (gas), 1 (dark matter), 2 (disc), 3 (bulge), 4 (stars), 5 (boundaries) or any present combination, c(2,3), and the mass-to-light ratio of each of those particle types.
galaxy_file = sim_data(system.file("extdata", 'SimSpin_example.hdf5', package="SimSpin"),
ptype = NA,
m2l_disc = 2,
m2l_bulge = 1) # reading in the snapshot dataTo begin, we can use the sim_analysis() function to calculate a few useful quantities for our model. Using this function, we can assess several spatial and kinematic properties such as:
log(density) distributionsim_data = sim_analysis(simdata = galaxy_file,
bin_type = "r",
rmax = 200, # maximum radius considered, kpc
rbin = 200, # number of radial bins considered
DM_profile = list(profile="Hernquist", DM_mass=184.97, DM_a=34.51))When using sim_analysis(), we specify the snapshot filename, the maximum radius rmax in kpc we want to consider and the number of radial bins that we wish to divide the galaxy into, rbin. In the above example, all particles are grouped into bins of width rmax / rbin = 1 kpc. Because we have specified binddir = "r" our bins are spherical shells. If we wished to study the distributions of certain components in cylindrical coordinates, we would specify either bindir = "cr" (cylindrical radial coordinated) or ="z" (for 1D coordinates off the surface of the disc).
Below, a few of the possible distributions are plotted using the data produced by sim_analysis().
Alternatively, you can consider single components of the model. For example, if we wished to examine the properties of the galaxy disk, we specify ptype = 2. The same could be done to examine the bulge by specifying ptype = 3. If you wish to examine a combination of parameters, this is simply specified by ptype = c(2,3) to give, for example, the dark matter and disc components - though this vignette galaxy model contains no dark matter so this will result in an error as shown.
disk_file = sim_data(system.file("extdata", 'SimSpin_example.hdf5', package="SimSpin"),
ptype = 2) # considering just the disk particles
disk_data = sim_analysis(simdata = disk_file,
rmax = 200,
rbin = 200,
DM_profile = list(profile="Hernquist", DM_mass=184.97, DM_a=34.51))
bulge_file = sim_data(system.file("extdata", 'SimSpin_example.hdf5', package="SimSpin"),
ptype = 3) # considering just the disk particles
bulge_data = sim_analysis(simdata = bulge_file,
rmax = 200,
rbin = 200,
DM_profile = list(profile="Hernquist", DM_mass=184.97, DM_a=34.51))
dm_file = sim_data(system.file("extdata", 'SimSpin_example.hdf5', package="SimSpin"),
ptype = 1) # considering just the disk particles## Particles of PartType1 are missing in this model.
## Error in sim_data(system.file("extdata", "SimSpin_example.hdf5", package = "SimSpin"), : Npart Error
An important part of this package is to take mock integral field unit (IFU) observations of galaxy models. The output of such an observation is a kinematic data cube, in which the galaxy is binned in 3-dimensions; binned spatially in the projected x-y plane and binned via velocity in the z-direction. The function build_datacube() can be used to construct such a data product for simulated galaxy models.
SAMI_cube_phys = build_datacube(simdata = galaxy_file,
r200 = 200,
z = 0.06,
fov = 15,
ap_shape = "circular",
central_wvl = 4800,
lsf_fwhm = 2.56,
pixel_sscale = 0.5,
pixel_vscale = 1.04,
inc_deg = 70,
threshold = 25) Demontrating the cube produced by the ifu_cube() function, moving through each velocity plane.
The animation above demonstrates the sequential planes of the IFU data cube that is produced by the function. The $datacube that is produced by the function is a 3D array of counts; each frame represents the particles that have velocities that fall into the specified vbin and displays the spatial arrangement of particles that fall into each sbin grid cell.
The spatial and velocity midpoints of each bin are stored in $xbin_labels, $zbin_labels and $vbin_labels. The axis ratio of the galaxy observed is also calculated and output as $axis_ratio$a and $axis_ratio$b as the major and minor axes respectively.
The method emplotyed by this function is to assume that each particle in the simulation has some inherent uncertainty in its velocity. This mimics the idea that, when astronomers observe emission lines, those lines have a width representing an uncertianty in the true speed in the host environment at which the line originated. This is the line spread function (LSF), as specified by lsf_fwhm and stems from a spectral response of the observing telescope to a point like source. We use the LSF supplied by the user to dictate the width of the particle’s possible velocities.
As the LSF of IFU instruments can be well approximated by a Gaussian, we model each particles velocity as a Gaussian centred on the true value but with a full-width half-maximum specified by the characteristics of the mock observation telescope. Each particle’s associated Gaussian can then be summed and binned in both spatial projection and velocity space in order to construct the kinematic data cube.
In order to mimic the effects of beam smearing and observational seeing in our mock kinematic data cube, we can optionally specify the blur argument. In this, each projected sptial plane in the kinematic data cube will be convolved with a PSF kernel of either “Moffat” or “Gaussian” shape with a FWHM specified in arcseconds.
SAMI_cube_blur = build_datacube(simdata = galaxy_file,
r200 = 200,
z = 0.06,
fov = 15,
ap_shape = "circular",
central_wvl = 4800,
lsf_fwhm = 2.56,
pixel_sscale = 0.5,
pixel_vscale = 1.04,
inc_deg = 70,
threshold = 25,
blur = list("psf" = "Moffat", "fwhm" = 2)) The output of this function is in the same format as before, but with the data cube convolved such that the images produced are blurred by the form (“Moffat”) and degree (2’’) specified by list("psf" = "Moffat", "fwhm" = 2).
The main purpose of the SimSpin package is to measure the observable spin parameter, \(\lambda_R\). The final basic function, find_lambda(), starts by building an IFU data cube, as done above, and then analyses that cube to calculate \(\lambda_R\). Each cube is collapsed along each spaxel to create three images of the observed galaxy in projection: a flux image, a line-of-sight (LOS) velocity image and a LOS velocity dispersion image. The function then uses these images to calculate the the observed \(\lambda_R\).
First we calculate the effective radius of the observed galaxy. The observed axis ratio is calculated by examining the covariance of the spatial arrangement of the particles in projection. An ellipse with this axis ratio is then grown from the centre of the simulation until half the total number of particles are contained. This is one effective radius, R\(_{eff}\).
Within the red ellipse, we sum the contributions from each pixel using the equation:
\[\lambda_R = \frac{\sum_{i=1}^{n_p} F_i R_i |V_i|}{\sum_{i=1}^{n_p} F_i R_i \sqrt{V_i^2 + \sigma_i^2}},\] where \(F_i\) is the observed flux taken from the flux image, \(R_i\) is the radial position, \(V_i\) is the LOS velocity taken from the LOS velocity image, \(\sigma_i\) is the LOS velocity dispersion taken from the LOS dispersion image per pixel, \(i\), and summed across the total number of pixels, \(n_p\).
SAMI_lambda = find_lambda(simdata = galaxy_file,
r200 = 200,
z = 0.06,
fov = 15,
ap_shape = "circular",
central_wvl = 4800,
lsf_fwhm = 2.56,
pixel_sscale = 0.5,
pixel_vscale = 1.04,
inc_deg = 70,
measure_type = list(type = "fit", fac = 1),
threshold = 25,
blur = list("psf" = "Moffat", "fwhm" = 2))SAMI_lambda$lambda_r # the observed lambda_R spin parameter## [1] 0.1568817
This function will return the observed spin parameter ($obs_lambdar), three IFU images (the $counts_img, $velocity_img and $dispersion_img) and the coordinates of an ellipse marking the radius at which the \(\lambda_R\) parameter has been measured ($reff_ellipse).
If the measurement radius specified is larger than the aperture, the function will return a warning.
SAMI_lambda_warning = find_lambda(simdata = galaxy_file,
r200 = 200,
z = 0.06,
fov = 15,
ap_shape = "circular",
central_wvl = 4800,
lsf_fwhm = 2.56,
pixel_sscale = 0.5,
pixel_vscale = 1.04,
inc_deg = 70,
measure_type = list(type = "fit", fac = 5),
threshold = 25,
blur = list("psf" = "Moffat", "fwhm" = 2))## WARNING: reff > aperture, the value of $obs_lambdar produced will not be the true value evaluated at reff.
The level of blurring applied and the measurement radius within which \(\lambda_R\) is measured can be optionally altered to study the effects that these choices have on the measurement. If the measurement_radius or blur arguments are missing, no atmospheric seeing conditions will be applied and \(\lambda_R\) will, by default, be measured within 1 effective radius.
We can then compare the value of \(\lambda_R\) that is output from this calculation to the one measured using sim_analysis() from the simulated model.
Once we have calculated all of these parameters, we can see how the measured \(\lambda_R\) value compares to the true spin parameter we measure from the simulation. Obviously this is very dependent on the inclination angle at which we measure the galaxy and the accuracy to which we can measure the inclination observationally.
This package is designed to help astronomers investigate the observational limitations of measuring the spin parameter of a galaxy. This vignette has shown a series of possible steps to do so.