yt.analysis_modules.cosmological_observation.light_ray.light_ray module

LightRay class and member functions.

class yt.analysis_modules.cosmological_observation.light_ray.light_ray.LightRay(parameter_filename, simulation_type=None, near_redshift=None, far_redshift=None, use_minimum_datasets=True, max_box_fraction=1.0, deltaz_min=0.0, minimum_coherent_box_fraction=0.0, time_data=True, redshift_data=True, find_outputs=False, load_kwargs=None)[source]

Bases: yt.analysis_modules.cosmological_observation.cosmology_splice.CosmologySplice

A 1D object representing the path of a light ray passing through a simulation. LightRays can be either simple, where they pass through a single dataset, or compound, where they pass through consecutive datasets from the same cosmological simulation. One can sample any of the fields intersected by the LightRay object as it passed through the dataset(s).

For compound rays, the LightRay stacks together multiple datasets in a time series in order to approximate a LightRay’s path through a volume and redshift interval larger than a single simulation data output. The outcome is something akin to a synthetic QSO line of sight.

Once the LightRay object is set up, use LightRay.make_light_ray to begin making rays. Different randomizations can be created with a single object by providing different random seeds to make_light_ray.

Parameters:
  • parameter_filename (string or Dataset) – For simple rays, one may pass either a loaded dataset object or the filename of a dataset. For compound rays, one must pass the filename of the simulation parameter file.
  • simulation_type (optional, string) – This refers to the simulation frontend type. Do not use for simple rays. Default: None
  • near_redshift (optional, float) – The near (lowest) redshift for a light ray containing multiple datasets. Do not use for simple rays. Default: None
  • far_redshift (optional, float) – The far (highest) redshift for a light ray containing multiple datasets. Do not use for simple rays. Default: None
  • use_minimum_datasets (optional, bool) – If True, the minimum number of datasets is used to connect the initial and final redshift. If false, the light ray solution will contain as many entries as possible within the redshift interval. Do not use for simple rays. Default: True.
  • max_box_fraction (optional, float) – In terms of the size of the domain, the maximum length a light ray segment can be in order to span the redshift interval from one dataset to another. If using a zoom-in simulation, this parameter can be set to the length of the high resolution region so as to limit ray segments to that size. If the high resolution region is not cubical, the smallest side should be used. Default: 1.0 (the size of the box)
  • deltaz_min (optional, float) – Specifies the minimum \(\Delta z\) between consecutive datasets in the returned list. Do not use for simple rays. Default: 0.0.
  • minimum_coherent_box_fraction (optional, float) – Use to specify the minimum length of a ray, in terms of the size of the domain, before the trajectory is re-randomized. Set to 0 to have ray trajectory randomized for every dataset. Set to np.inf (infinity) to use a single trajectory for the entire ray. Default: 0.
  • time_data (optional, bool) – Whether or not to include time outputs when gathering datasets for time series. Do not use for simple rays. Default: True.
  • redshift_data (optional, bool) – Whether or not to include redshift outputs when gathering datasets for time series. Do not use for simple rays. Default: True.
  • find_outputs (optional, bool) – Whether or not to search for datasets in the current directory. Do not use for simple rays. Default: False.
  • load_kwargs (optional, dict) – If you are passing a filename of a dataset to LightRay rather than an already loaded dataset, then you can optionally provide this dictionary as keywords when the dataset is loaded by yt with the “load” function. Necessary for use with certain frontends. E.g. Tipsy using “bounding_box” Gadget using “unit_base”, etc. Default : None
create_cosmology_splice(near_redshift, far_redshift, minimal=True, max_box_fraction=1.0, deltaz_min=0.0, time_data=True, redshift_data=True)

Create list of datasets capable of spanning a redshift interval.

For cosmological simulations, the physical width of the simulation box corresponds to some Delta z, which varies with redshift. Using this logic, one can stitch together a series of datasets to create a continuous volume or length element from one redshift to another. This method will return such a list

Parameters:
  • near_redshift (float) – The nearest (lowest) redshift in the cosmology splice list.
  • far_redshift (float) – The furthest (highest) redshift in the cosmology splice list.
  • minimal (bool) – If True, the minimum number of datasets is used to connect the initial and final redshift. If false, the list will contain as many entries as possible within the redshift interval. Default: True.
  • max_box_fraction (float) – In terms of the size of the domain, the maximum length a light ray segment can be in order to span the redshift interval from one dataset to another. If using a zoom-in simulation, this parameter can be set to the length of the high resolution region so as to limit ray segments to that size. If the high resolution region is not cubical, the smallest side should be used. Default: 1.0 (the size of the box)
  • deltaz_min (float) – Specifies the minimum delta z between consecutive datasets in the returned list. Default: 0.0.
  • time_data (bool) – Whether or not to include time outputs when gathering datasets for time series. Default: True.
  • redshift_data (bool) – Whether or not to include redshift outputs when gathering datasets for time series. Default: True.

Examples

>>> co = CosmologySplice("enzo_tiny_cosmology/32Mpc_32.enzo", "Enzo")
>>> cosmo = co.create_cosmology_splice(1.0, 0.0)
make_light_ray(seed=None, periodic=True, left_edge=None, right_edge=None, min_level=None, start_position=None, end_position=None, trajectory=None, fields=None, setup_function=None, solution_filename=None, data_filename=None, get_los_velocity=None, use_peculiar_velocity=True, redshift=None, field_parameters=None, njobs=-1)[source]
make_light_ray(seed=None, periodic=True,
left_edge=None, right_edge=None, min_level=None, start_position=None, end_position=None, trajectory=None, fields=None, setup_function=None, solution_filename=None, data_filename=None, use_peculiar_velocity=True, redshift=None, njobs=-1)

Create a light ray and get field values for each lixel. A light ray consists of a list of field values for cells intersected by the ray and the path length of the ray through those cells. Light ray data must be written out to an hdf5 file.

Parameters:
  • seed (optional, int) – Seed for the random number generator. Default: None.
  • periodic (optional, bool) – If True, ray trajectories will make use of periodic boundaries. If False, ray trajectories will not be periodic. Default : True.
  • left_edge (optional, iterable of floats or YTArray) – The left corner of the region in which rays are to be generated. If None, the left edge will be that of the domain. If specified without units, it is assumed to be in code units. Default: None.
  • right_edge (optional, iterable of floats or YTArray) – The right corner of the region in which rays are to be generated. If None, the right edge will be that of the domain. If specified without units, it is assumed to be in code units. Default: None.
  • min_level (optional, int) – The minimum refinement level of the spatial region in which the ray passes. This can be used with zoom-in simulations where the high resolution region does not keep a constant geometry. Default: None.
  • start_position (optional, iterable of floats or YTArray.) – Used only if creating a light ray from a single dataset. The coordinates of the starting position of the ray. If specified without units, it is assumed to be in code units. Default: None.
  • end_position (optional, iterable of floats or YTArray.) – Used only if creating a light ray from a single dataset. The coordinates of the ending position of the ray. If specified without units, it is assumed to be in code units. Default: None.
  • trajectory (optional, list of floats) – Used only if creating a light ray from a single dataset. The (r, theta, phi) direction of the light ray. Use either end_position or trajectory, not both. Default: None.
  • fields (optional, list) – A list of fields for which to get data. Default: None.
  • setup_function (optional, callable, accepts a ds) – This function will be called on each dataset that is loaded to create the light ray. For, example, this can be used to add new derived fields. Default: None.
  • solution_filename (optional, string) – Path to a text file where the trajectories of each subray is written out. Default: None.
  • data_filename (optional, string) – Path to output file for ray data. Default: None.
  • use_peculiar_velocity (optional, bool) – If True, the peculiar velocity along the ray will be sampled for calculating the effective redshift combining the cosmological redshift and the doppler redshift. Default: True.
  • redshift (optional, float) – Used with light rays made from single datasets to specify a starting redshift for the ray. If not used, the starting redshift will be 0 for a non-cosmological dataset and the dataset redshift for a cosmological dataset. Default: None.
  • njobs (optional, int) – The number of parallel jobs over which the segments will be split. Choose -1 for one processor per segment. Default: -1.

Examples

Make a light ray from multiple datasets:

>>> import yt
>>> from yt.analysis_modules.cosmological_observation.light_ray.api import         ...     LightRay
>>> my_ray = LightRay("enzo_tiny_cosmology/32Mpc_32.enzo", "Enzo",
...                   0., 0.1, time_data=False)
...
>>> my_ray.make_light_ray(seed=12345,
...                       solution_filename="solution.txt",
...                       data_filename="my_ray.h5",
...                       fields=["temperature", "density"],
...                       use_peculiar_velocity=True)

Make a light ray from a single dataset:

>>> import yt
>>> from yt.analysis_modules.cosmological_observation.light_ray.api import         ...     LightRay
>>> my_ray = LightRay("IsolatedGalaxy/galaxy0030/galaxy0030")
...
>>> my_ray.make_light_ray(start_position=[0., 0., 0.],
...                       end_position=[1., 1., 1.],
...                       solution_filename="solution.txt",
...                       data_filename="my_ray.h5",
...                       fields=["temperature", "density"],
...                       use_peculiar_velocity=True)
plan_cosmology_splice(near_redshift, far_redshift, decimals=3, filename=None, start_index=0)

Create imaginary list of redshift outputs to maximally span a redshift interval.

If you want to run a cosmological simulation that will have just enough data outputs to create a cosmology splice, this method will calculate a list of redshifts outputs that will minimally connect a redshift interval.

Parameters:
  • near_redshift (float) – The nearest (lowest) redshift in the cosmology splice list.
  • far_redshift (float) – The furthest (highest) redshift in the cosmology splice list.
  • decimals (int) – The decimal place to which the output redshift will be rounded. If the decimal place in question is nonzero, the redshift will be rounded up to ensure continuity of the splice. Default: 3.
  • filename (string) – If provided, a file will be written with the redshift outputs in the form in which they should be given in the enzo dataset. Default: None.
  • start_index (int) – The index of the first redshift output. Default: 0.

Examples

>>> from yt.analysis_modules.api import CosmologySplice
>>> my_splice = CosmologySplice('enzo_tiny_cosmology/32Mpc_32.enzo', 'Enzo')
>>> my_splice.plan_cosmology_splice(0.0, 0.1, filename='redshifts.out')
yt.analysis_modules.cosmological_observation.light_ray.light_ray.non_periodic_ray(ds, left_edge, right_edge, ray_length, max_iter=5000, min_level=None, my_random=None)[source]
yt.analysis_modules.cosmological_observation.light_ray.light_ray.periodic_adjust(p, left=None, right=None)[source]

Return the point p adjusted for periodic boundaries.

yt.analysis_modules.cosmological_observation.light_ray.light_ray.periodic_distance(coord1, coord2)[source]

Calculate length of shortest vector between to points in periodic domain.

yt.analysis_modules.cosmological_observation.light_ray.light_ray.periodic_ray(start, end, left=None, right=None)[source]

Break up periodic ray into non-periodic segments. Accepts start and end points of periodic ray as YTArrays. Accepts optional left and right edges of periodic volume as YTArrays. Returns a list of lists of coordinates, where each element of the top-most list is a 2-list of start coords and end coords of the non-periodic ray:

[[[x0start,y0start,z0start], [x0end, y0end, z0end]],
[[x1start,y1start,z1start], [x1end, y1end, z1end]], ...,]
yt.analysis_modules.cosmological_observation.light_ray.light_ray.vector_length(start, end)[source]

Calculate vector length.