Profiles¶
- class plonk.Profile(snap, ndim=2, cmin=None, cmax=None, n_bins=100, aggregation='average', spacing='linear', coordinate='x', ignore_accreted=True)¶
Profiles.
A profile is a binning of particles in Cartesian slices, or cylindrical or spherical shells around the origin, i.e. (0, 0, 0). For cylindrical profiles the cylindrical cells are perpendicular to the xy-plane, i.e. the bins are averaged azimuthally and in the z-direction. Cartesian profiles can be in the x-, y-, and z-directions.
- Parameters
snap (SnapLike) – The Snap object.
ndim (optional) – The dimension of the profile. For ndim == 2, the radial binning is cylindrical in the xy-plane. For ndim == 3, the radial binning is spherical. For ndim == 1, the radial binning is Cartesian along the x-axis. Default is 2.
cmin (optional) – The minimum coordinate for binning. Can be a string, e.g. ‘10 au’, or a quantity with units, e.g. plonk.units(‘10 au’). Defaults to minimum on the particles.
cmax (optional) – The maximum coordinate for binning. Can be a string, e.g. ‘10 au’, or a quantity with units, e.g. plonk.units(‘10 au’). Defaults to the 99 percentile distance.
n_bins (optional) – The number of radial bins. Default is 100.
aggregation (optional) – The method to aggregate particle quantities in bins by. Options are ‘average’, ‘mean’, or ‘median’. Here ‘average’ is a mass-weighted average. Default is ‘average’.
spacing (optional) – The spacing of radial bins. Can be ‘linear’ or ‘log’. Default is ‘linear’.
coordinate (optional) – The coordinate (‘x’, ‘y’, or ‘z’) for Cartesian profiles only, i.e. when ndim==1. Default is ‘x’. For cylindrical and spherical profiles the coordinate is ‘radius’.
ignore_accreted (optional) – Ignore particles accreted onto sinks. Default is True.
Examples
Generate profile from snapshot.
>>> prof = plonk.load_profile(snap=snap) >>> prof = plonk.load_profile(snap=snap, n_bins=300) >>> prof = plonk.load_profile(snap=snap, cmin='10 au', cmax='300 au') >>> prof = plonk.load_profile(snap=snap, spacing='log')
To access a profile.
>>> prof['surface_density'] >>> prof['scale_height']
To set a new profile.
>>> prof['aspect_ratio'] = prof['scale_height'] / prof['radius']
Alternatively use the add_profile decorator.
>>> @prof.add_profile ... def mass(prof): ... M = prof.snap['mass'] ... return prof.particles_to_binned_quantity('sum', M)
Plot one or many quantities on the profile.
>>> prof.plot('radius', 'density') >>> prof.plot('radius', ['angular_momentum_x', 'angular_momentum_y'])
Plot a quantity on the profile with units.
>>> units = {'position': 'au', 'surface_density'='g/cm^2'} >>> prof.plot('radius', 'surface_density', units=units)
- add_alias(name, alias)¶
Add alias to array.
- add_profile(fn)¶
Decorate function to add profile to Profile.
- Parameters
fn (Callable) – A function that returns the profile as an array. The name of the function is the string with which to reference the array.
- Returns
The function which returns the array.
- Return type
Callable
- base_profile_name(name)¶
Get the base profile name from a string.
For example, ‘velocity_x’ returns ‘velocity’, ‘density’ returns ‘density’, ‘dust_fraction_001’ returns ‘dust_fraction’, ‘x’ returns ‘position’.
- particles_to_binned_quantity(aggregation, array)¶
Calculate binned quantities from particles.
This takes care of the bin indices and ignoring accreted particles (if requested in instantiating the profile).
- Parameters
aggregation (str) – The aggregation function that acts on particles in the radial bin.
array (pint.quantity.build_quantity_class.<locals>.Quantity) – The particle array.
- Return type
pint.quantity.build_quantity_class.<locals>.Quantity
- plot(x, y, units=None, std=None, label=None, ax=None, ax_kwargs={}, **kwargs)¶
Plot profile.
- Parameters
x (str) – The x axis to plot as a string.
y (Union[str, List[str]]) – The y axis to plot. Can be string or a list of strings.
units (Optional[Dict[str, Union[str, List[str]]]]) – The units of the plot as a dictionary. The keys correspond to quantities such as ‘position’, ‘density’, ‘velocity’, and so on. The values are strings representing units, e.g. ‘g/cm^3’ for density.
std (optional) – Add standard deviation on profile. Can be ‘shading’ or ‘errorbar’.
label (optional) – A label for the plot. Can be a string or a list of strings, one per y.
ax (optional) – A matplotlib Axes object to plot to.
ax_kwargs – Keyword arguments to pass to matplotlib Axes.
**kwargs – Keyword arguments to pass to Axes plot method.
- Returns
The matplotlib Axes object.
- Return type
ax
- set_units(**kwargs)¶
Set default unit for profiles.
- Parameters
kwargs – Keyword arguments with keys as the profile name, e.g. ‘pressure’, and with values as the unit as a string, e.g. ‘pascal’.
- Return type
Examples
Set multiple default units.
>>> profile.set_units(pressure='pascal', density='g/cm^3')
- to_dataframe(columns=None, units=None)¶
Convert Profile to DataFrame.
- Parameters
columns (optional) – A list of columns to add to the data frame. If None, add all loaded columns. Default is None.
units (optional) – A list of units corresponding to columns add to the data frame. Units must be strings, and must be base units. I.e. ‘cm’ not ‘10 cm’. If None, use default, i.e. cgs. Default is None.
- Returns
- Return type
DataFrame
- to_function(profile, **kwargs)¶
Create function via interpolation.
The function is of the coordinate of the profile, e.g. ‘radius’, and returns values of the selected profile, e.g. ‘scale_height’. The function is generated from the scipy.interpolate function interp1d.
- Parameters
profile (str) – The profile function to create as a string, e.g. ‘scale_height’.
- Returns
The function.
- Return type
Callable
Examples
Select all particles within a scale height in a disc.
>>> scale_height = prof.to_function('scale_height') >>> subsnap = snap[np.abs(snap['z']) < scale_height(snap['R'])]
- plonk.load_profile(snap, ndim=2, cmin=None, cmax=None, n_bins=100, aggregation='average', spacing='linear', coordinate='x', ignore_accreted=True)¶
Profiles.
A profile is a binning of particles in Cartesian slices, or cylindrical or spherical shells around the origin, i.e. (0, 0, 0). For cylindrical profiles the cylindrical cells are perpendicular to the xy-plane, i.e. the bins are averaged azimuthally and in the z-direction. Cartesian profiles can be in the x-, y-, and z-directions.
- Parameters
snap (SnapLike) – The Snap object.
ndim (optional) – The dimension of the profile. For ndim == 2, the radial binning is cylindrical in the xy-plane. For ndim == 3, the radial binning is spherical. For ndim == 1, the radial binning is Cartesian along the x-axis. Default is 2.
cmin (optional) – The minimum coordinate for binning. Can be a string, e.g. ‘10 au’, or a quantity with units, e.g. plonk.units(‘10 au’). Defaults to minimum on the particles.
cmax (optional) – The maximum coordinate for binning. Can be a string, e.g. ‘10 au’, or a quantity with units, e.g. plonk.units(‘10 au’). Defaults to the 99 percentile distance.
n_bins (optional) – The number of radial bins. Default is 100.
aggregation (optional) – The method to aggregate particle quantities in bins by. Options are ‘average’, ‘mean’, or ‘median’. Here ‘average’ is a mass-weighted average. Default is ‘average’.
spacing (optional) – The spacing of radial bins. Can be ‘linear’ or ‘log’. Default is ‘linear’.
coordinate (optional) – The coordinate (‘x’, ‘y’, or ‘z’) for Cartesian profiles only, i.e. when ndim==1. Default is ‘x’. For cylindrical and spherical profiles the coordinate is ‘radius’.
ignore_accreted (optional) – Ignore particles accreted onto sinks. Default is True.
- Return type
Examples
Generate profile from snapshot.
>>> prof = plonk.load_profile(snap=snap) >>> prof = plonk.load_profile(snap=snap, n_bins=300) >>> prof = plonk.load_profile(snap=snap, cmin='10 au', cmax='300 au') >>> prof = plonk.load_profile(snap=snap, spacing='log')
To access a profile.
>>> prof['surface_density'] >>> prof['scale_height']
To set a new profile.
>>> prof['aspect_ratio'] = prof['scale_height'] / prof['radius']
Alternatively use the add_profile decorator.
>>> @prof.add_profile ... def mass(prof): ... M = prof.snap['mass'] ... return prof.particles_to_binned_quantity('sum', M)
Plot one or many quantities on the profile.
>>> prof.plot('radius', 'density') >>> prof.plot('radius', ['angular_momentum_x', 'angular_momentum_y'])
Plot a quantity on the profile with units.
>>> units = {'position': 'au', 'surface_density'='g/cm^2'} >>> prof.plot('radius', 'surface_density', units=units)