Code Documentation
This page shows the documentation generated by sphinx automatically scanning our source code
Data Loaders
- class ParticlePhaseSpace.DataLoaders.Load_IAEA(data_schema: ~numpy.dtype, constants: dict, input_data: (<class 'str'>, <class 'pathlib.Path'>), n_records: int = -1, offset=0, **kwargs)
this loads a binary varian IAEA sent through the topas forums. Because this format is so arbitrary, uses are required to pass a data_schema variable indicating the order and types of the data in the phase space, because there is no general way for us to figure this out. Please see here for examples of how to use this data loader.
- Parameters:
data_schema (np.dtype) – the types and order of data, specified as an np.dtype
constants – any constants in the phase space
input_data – path to the a .phsp or .IAEAphsp file
n_records – specify how many rows of data to read in. By default, will read all rows.
offset – which row to start at. defaults to 0 (first row). Can be used in conjunction with n_records to read a large file in a series of small chunks
Example:
from ParticlePhaseSpace import PhaseSpace, DataLoaders from pathlib import Path import numpy as np file_name = Path(r'/home/brendan/Downloads/Varian_TrueBeam6MV_01.phsp') data_schema = np.dtype([ ('particle type', 'i1'), ('Ek', 'f4'), ('x', 'f4'), ('y', 'f4'), ('z', 'f4'), ('Cosine X', 'f4'), ('Cosine Y', 'f4') ]) constants = {'weight': np.int8(1)} ps_data = DataLoaders.Load_IAEA(data_schema=data_schema, constants=constants, input_data=file_name, n_records=int(1e5)) PS = PhaseSpace(ps_data)
- class ParticlePhaseSpace.DataLoaders.Load_PandasData(input_data: (<class 'str'>, <class 'pathlib.Path'>), particle_type: (<class 'str'>, None) = None, units: ~ParticlePhaseSpace.__unit_config__.UnitSet = <ParticlePhaseSpace.__unit_config__.UnitSet object>)
loads in pandas data of the format. This is used internally by ParticlePhaseSpace, and can also be used externally in cases where it is not desired to write a dedicated new data loader:
from ParticlePhaseSpace import DataLoaders import pandas as pd demo_data = pd.DataFrame( {'x [mm]': [0, 1, 2], 'y [mm]': [0, 1, 2], 'z [mm]': [0, 1, 2], 'px [MeV/c]': [0, 1, 2], 'py [MeV/c]': [0, 1, 2], 'pz [MeV/c]': [0, 1, 2], 'particle type [pdg_code]': [11, 11, 11], 'weight': [0, 1, 2], 'particle id': [0, 1, 2], 'time [ps]': [0, 1, 2]}) data = DataLoaders.Load_PandasData(demo_data)
- class ParticlePhaseSpace.DataLoaders.Load_TibarayData(input_data: (<class 'str'>, <class 'pathlib.Path'>), particle_type: (<class 'str'>, None) = None, units: ~ParticlePhaseSpace.__unit_config__.UnitSet = <ParticlePhaseSpace.__unit_config__.UnitSet object>)
Load ASCII data from tibaray of format x y z rxy Bx By Bz G t m q nmacro rmacro ID:
data_loc = Path(r'../tests/test_data/tibaray_test.dat') data = DataLoaders.Load_TibarayData(data_loc, particle_type=11) PS = PhaseSpace(data)
- class ParticlePhaseSpace.DataLoaders.Load_TopasData(input_data: (<class 'str'>, <class 'pathlib.Path'>), particle_type: (<class 'str'>, None) = None, units: ~ParticlePhaseSpace.__unit_config__.UnitSet = <ParticlePhaseSpace.__unit_config__.UnitSet object>)
DataLoader for Topas data. This data loader will read in both ascii and binary topas phase space (phsp) files. At present, we do not handle time or particle-id fields which may or may not be present in topas data. Behind the scenes, it relies on topas2numpy:
from ParticlePhaseSpace import DataLoaders from ParticlePhaseSpace import PhaseSpace from pathlib import Path data_loc = Path(r'../tests/test_data/coll_PhaseSpace_xAng_0.00_yAng_0.00_angular_error_0.0.phsp') data = DataLoaders.Load_TopasData(data_loc) PS = PhaseSpace(data)
- class ParticlePhaseSpace.DataLoaders.Load_p2sat_txt(input_data: (<class 'str'>, <class 'pathlib.Path'>), particle_type: (<class 'str'>, None) = None, units: ~ParticlePhaseSpace.__unit_config__.UnitSet = <ParticlePhaseSpace.__unit_config__.UnitSet object>)
Adapted from the p2sat ‘txt’ loader; loads csv data of format # weight x (um) y (um) z (um) px (MeV/c) py (MeV/c) pz (MeV/c) t (fs) Note that we use a hard coded seperator value “,”.
available_units = ParticlePhaseSpaceUnits() data_url = 'https://raw.githubusercontent.com/lesnat/p2sat/master/examples/ExamplePhaseSpace.csv' file_name = 'p2sat_txt_test.csv' request.urlretrieve(data_url, file_name) # read in ps_data = DataLoaders.Load_p2sat_txt(file_name, particle_type='electrons', units=available_units('p2_sat_UHI')) PS = PhaseSpace(ps_data)
PhaseSpace
- class ParticlePhaseSpace._ParticlePhaseSpace.PhaseSpace(data_loader)
This class holds phase space data in a pandas dataframe, and allowed users to utilise common libraries for plotting and analysis. It accepts data from any DataLoader Basic use is documented here.
- Parameters:
data_loader (_DataLoadersBase) – an instance of ParticlePhaseSpace.DataLoaders._DataLoadersBase
- assess_density_versus_r(Rvals=None, verbose: bool = True, beam_direction: str = 'z')
Assess how many particles are in a given radius
- Parameters:
Rvals – list of rvals to assess in mm, e.g. np.linspace(0, 2, 21)
verbose – prints data to screen if True
beam_direction (str, optional) – main direction in which beam is travelling. ‘x’, ‘y’, or ‘z’ (default)
- Return density_data:
a dataframe containing the roi vals and the proportion of particles inside
- calculate_twiss_parameters(beam_direction='z')
Calculate the RMS twiss parameters
- Parameters:
beam_direction (str, optional) – main direction in which beam is travelling. ‘x’, ‘y’, or ‘z’ (default)
- Returns:
None
- filter_by_boolean_index(boolean_index, in_place: bool = False, split: bool = False, verbose: bool = True)
filter data by input boolean index, keeping ‘True’ and discarding ‘False’
- Parameters:
boolean_index – an 1D array like structure of True and False values
in_place – if True, existing object is modified; if False a new object is returned
split – if True, will return two phase space objects: one where boolan_index=True and one where it equals False
- filter_by_time(t_start, t_finish, in_place: bool = False)
Generates a new PhaseSpace which only contains particles inside t_start and t_finish (inclusive). t_start and t_finish should be specfied in ps.
- Parameters:
t_start (float) – particles with t<t_start are removed.
t_finish (float) – particleswith t>t_finish are removed
- Returns:
new_instance: a new phase space object with data filtered by time
- get_downsampled_phase_space(downsample_factor: int = 10)
return a new phase space object which randomlt samples from the larger phase space. the new phase space has size ‘original data/downsample_factor’. the data is shuffled before being randomly sampled.
- Parameters:
downsample_factor (int) – the factor to downsample the phase space by
- merge(in_place=False)
merges identical data points by combining their weights. Typically, before performing a merge operation you will want to perform a ‘regrid’ operation. The underlying algorithm was developed by Leo Esnault for the p2sat code. Merged particles retain the particle ID from the first particle in the merged group.
- Parameters:
in_place – if True, self is operated on; if False, a new PhaseSpace is returned
- Returns:
new_PS if in_place is False.
- print_energy_stats(file_name=None)
Prints a summary of energy stats to the screen, which can optionally be saved to json
- Parameters:
file_name (Path, str, optional) – if specified, the data is saved as json in file_name
- print_twiss_parameters(file_name=None, beam_direction: str = 'z')
prints the twiss parameters if they exist they are always printed to the screen. if filename is specified, they are also saved to file as json
- Parameters:
file_name (str or Path, optional) – optional filename to write twiss data to. should be absolute path to an existing directory
beam_direction (str, optional) – the direction the beam is travelling in. “x”, “y”, or “z” (default)
- resample_via_gaussian_kde(n_new_particles_factor: int = 1, interpolate_weights: (<class 'bool'>, None) = None)
Generate a new phase space based on the existing data, by fitting a gaussian kernel density estimate to the 6-D space: x y z px py pz If interpolate_weights is set to True, we instead attempt to interpolate within a 7D space: x y z px py pz weight This method is fairly experimental and should be used with extreme caution!
- Parameters:
n_new_particles_factor – the returned Phase space will have size len(self)*n_new_particles_factor. In other words, when n_new_particles_factor=1, the new PhaseSpace will be the same size as the original.
- Returns:
new_PS: a new PhaseSpace object
- reset_phase_space()
reduce self._ps_data to only the required columns delete any other derived quantities such as twiss parameters this can be called whenever you want to reduce the memory footprint
- set_units(new_units: UnitSet)
converts ps_data to a new unit set. This will also reset the phase space to just the required columns
- Parameters:
new_units (UnitSet) – the new units to convert to
- sort(quantities_to_sort: (None, <class 'list'>) = None)
sort the data. Data will be sorted according quantities_to_sort, in order of quantity. Operates in place. example:
PS.sort(quantities_to_sort='x') PS.sort(quantities_to_sort=['x', 'y', 'z', 'px'])
- Parameters:
quantities_to_sort
PhaseSpace.plot
- class ParticlePhaseSpace._ParticlePhaseSpace._Plots(PS)
- energy_hist_1D(n_bins: int = 100, grid: bool = False)
generate a histogram plot of paritcle energies. Each particle spcies present in the phase space is overlaid on the same plot.
- Parameters:
n_bins (int, optional) – number of bins in histogram
grid (bool, optional) – turns grid on/off
- Returns:
None
- momentum_hist_1D(n_bins: int = 100, alpha: float = 0.5, grid: bool = False)
plot a histogram of particle momentum in x, y, z. a new histogram is generated for each particle species. histograms are overlaid.
- Parameters:
n_bins – number of bins in histogram
alpha – controls transparency of each histogram.
grid (bool, optional) – turns grid on/off
- n_particles_v_time(n_bins: int = 100, grid: bool = False)
basic plot of number of particles versus time; useful for quickly seperating out different bunches of electrons such that you can apply the ‘filter_by_time’ method
- Parameters:
n_bins (int) – number of bins for histogram
grid (bool, optional) – turns grid on/off
- particle_positions_hist_2D(beam_direction: str = 'z', quantity: str = 'intensity', grid: bool = True, log_scale: bool = False, bins: int = 100, normalize: bool = True, xlim=None, ylim=None, vmin=None, vmax=None)
plot a 2D histogram of data, either of accumulated number of particules or accumulated energy
- Parameters:
beam_direction (str, optional) – the direction the beam is travelling in. “x”, “y”, or “z” (default)
xlim (list, optional) – set the xlim for all plots, e.g. [-2,2]
ylim (list, optional) – set the ylim for all plots, e.g. [-2,2]
quantity (str) – quantity to accumulate; either ‘intensity’ or ‘energy
grid (bool, optional) – turns grid on/off
bins (int, optional) – number of bins in X/Y direction. n_pixels = bins ** 2
vmin (float, optional) – minimum color range
vmax (float, optional) – maximum color range
log_scale (bool, optional) – if True, log scale is used
normalize (bool, optional) – if True, data is normalized to 0-100 - otherwise raw values are used
- Returns:
None
- particle_positions_scatter_2D(beam_direction: str = 'z', weight_position_plot: bool = False, grid: bool = True, xlim=None, ylim=None)
produce a scatter plot of particle positions. one plot is produced for each unique species.
- Parameters:
beam_direction (str, optional) – the direction the beam is travelling in. “x”, “y”, or “z” (default)
weight_position_plot (bool) – if True, a gaussian kde is used to weight the particle positions. This can produce very informative and useful plots, but also be very slow. If it is slow, you could try downsampling the phase space first using get_downsampled_phase_space
grid (bool, optional) – turns grid on/off
xlim (list or None, optional) – set the xlim for all plots, e.g. [-2,2]
ylim (list or None, optional) – set the ylim for all plots, e.g. [-2,2]
- Returns:
None
- position_hist_1D(n_bins: int = 100, alpha: float = 0.5, grid: bool = False)
plot a histogram of particle positions in x, y, z. a new histogram is generated for each particle species. histograms are overlaid.
- Parameters:
n_bins – number of bins in histogram
alpha – controls transparency of each histogram.
grid (bool, optional) – turns grid on/off
- transverse_trace_space_hist_2D(beam_direction: str = 'z', plot_twiss_ellipse: bool = True, grid: bool = True, bins: int = 100, log_scale: bool = True, normalize: bool = True, xlim=None, ylim=None, vmin=None, vmax=None)
plot the intensity of the beam in trace space
- Parameters:
beam_direction (str, optional) – the direction the beam is travelling in. “x”, “y”, or “z” (default)
xlim (list, optional) – set the xlim for all plots, e.g. [-2,2]
ylim (list, optional) – set the ylim for all plots, e.g. [-2,2]
plot_twiss_ellipse (bool, optional) – if True, RMS ellipse from twiss parameters is overlaid.
grid (bool, optional) – turns grid on/off
log_scale (bool, optional) – if True, log scale is used
bins (int, optional) – number of bins in X/Y direction. n_pixels = bins ** 2
vmin (float, optional) – minimum color range
vmax (float, optional) – maximum color range
- transverse_trace_space_scatter_2D(beam_direction: str = 'z', plot_twiss_ellipse: bool = True, grid: bool = True, xlim=None, ylim=None)
Generate a scatter plot of x versus x’=px/pz and y versus y’=py/pz (these definitions are for beam_direction=’z’)
- Parameters:
beam_direction (str, optional) – main direction in which beam is travelling. ‘x’, ‘y’, or ‘z’ (default)
plot_twiss_ellipse (bool, optional) – if True, will overlay the RMS twiss ellipse onto the trace space
xlim (list, optional) – set xlim, e.g. [-2,2]
ylim (list, optional) – set ylim, e.g. [-2,2]
grid (bool, optional) – turns grid on/off
PhaseSpace.fill
- class ParticlePhaseSpace._ParticlePhaseSpace._Fill_Methods(PS)
Methods for calculating secondary quantities and adding them to ps_data
- beta_and_gamma()
add the relatavistic beta and gamma factors into self._PS._ps_data
- direction_cosines()
Calculate direction cosines, which are required for topas import: U (direction cosine of momentum with respect to X) V (direction cosine of momentum with respect to Y)
- kinetic_E()
- Uses energy-momementum relation to add
kinetic energy into self._ps_data
- relativistic_mass()
add relativistic mass to ps_data
- rest_mass()
add rest mass to self._PS._ps_data :return:
- velocity()
add velocities in m/s into self._PS._ps_data
PhaseSpace.transform
Data Exporters
- class ParticlePhaseSpace.DataExporters.CSV_Exporter(PhaseSpaceInstance: ~ParticlePhaseSpace._ParticlePhaseSpace.PhaseSpace, output_location: (<class 'str'>, <class 'pathlib.Path'>), output_name: str)
Export data to a csv format, in particular one which can be read by p2sat read text # weight x (um) y (um) z (um) px (MeV/c) py (MeV/c) pz (MeV/c) t (fs)
- class ParticlePhaseSpace.DataExporters.Topas_Exporter(PhaseSpaceInstance: ~ParticlePhaseSpace._ParticlePhaseSpace.PhaseSpace, output_location: (<class 'str'>, <class 'pathlib.Path'>), output_name: str, binary: bool = False)
output the phase space to topas ascii or binary format. the default output is ascii, the user can output binary by passing the flag binary as a boolean e.g. binary=True
- Note:
we do not handle any time features
every particle in the phase space is flagged as being a new history.