ocelot.cpbd.beam
¶
definition of particles, beams and trajectories
Module Contents¶
Classes¶
class - container for twiss parameters |
|
particle |
|
array of particles of fixed size; for optimized performance |
|
Functions¶
|
|
|
Function to calculate twiss parameters form the ParticleArray |
|
Function calculates beam current from particleArray. |
|
|
|
|
|
|
|
|
|
|
|
Beam matching function, the beam is centered in the phase space |
|
|
|
|
|
|
|
Function to calculate beam current |
|
returns: |
|
|
|
|
|
|
|
Function to calculate slice parameters |
|
Function to calculate slice parameters |
|
reads ParticleArray() |
|
Method to generate ParticleArray with gaussian distribution. |
|
Generates BeamArray object |
-
ocelot.cpbd.beam.
_logger
¶
-
ocelot.cpbd.beam.
ne_flag
= True¶
-
ocelot.cpbd.beam.
nb_flag
= True¶
-
class
ocelot.cpbd.beam.
Twiss
(beam=None)¶ class - container for twiss parameters
-
__str__
(self)¶ Return str(self).
-
-
class
ocelot.cpbd.beam.
Particle
(x=0.0, y=0.0, px=0.0, py=0.0, s=0.0, p=0.0, tau=0.0, E=0.0, q=0.0)¶ particle to be used for tracking
-
__str__
(self)¶ Return str(self).
-
-
class
ocelot.cpbd.beam.
Beam
(x=0, xp=0, y=0, yp=0)¶ -
properties
= ['g', 'dg', 'emit_xn', 'emit_yn', 'p', 'pz', 'px', 'py']¶
-
property
g
(self)¶
-
property
dg
(self)¶
-
property
emit_xn
(self)¶
-
property
emit_yn
(self)¶
-
property
p
(self)¶
-
property
pz
(self)¶
-
property
px
(self)¶
-
property
py
(self)¶
-
len
(self)¶
-
to_array
(self, nslice=100, window_len=None)¶
-
sizes
(self)¶
-
-
class
ocelot.cpbd.beam.
BeamArray
(nslice=0)¶ Bases:
ocelot.cpbd.beam.Beam
-
idx_max
(self)¶
-
len
(self)¶
-
charge
(self)¶
-
params
(self)¶
-
sort
(self)¶
-
equidist
(self)¶
-
smear
(self, sw)¶
-
get_s
(self, s)¶
-
get_E
(self)¶
-
set_E
(self, E_GeV)¶
-
center
(self, s)¶
-
cut_empty_I
(self, thresh=0.01)¶
-
start_at_0
(self)¶
-
cleanup
(self)¶
-
__getitem__
(self, index)¶
-
__delitem__
(self, index)¶
-
pk
(self)¶
-
add_chirp
(self, chirp=0, order=1, s0=None)¶ adds energy chirp of given “order” to the beam around “center” position [m] chirp = dE/E0/ds[um]**order = dg/g0/s[um]**order so that g_new = g0 + dg = g0 + g0 (s-s0)**order * chirp
-
add_chirp_poly
(self, coeff, s0=None)¶ The method adds a polynomial energy chirp to the beam object.
coeff — coefficients for the chirp s0 — the point with respect to which the chirp will be introduced
The expression for the chirp:
E = E0((g0 + coeff[0])/g0 +
coeff[1]*(s - s0))**1 / 1! / ((speed_of_light * 1e-15)**1 * g0) +
coeff[2]*(s - s0))**2 / 2! / ((speed_of_light * 1e-15)**2 * g0) +
coeff[3]*(s - s0))**3 / 3! / ((speed_of_light * 1e-15)**3 * g0) + …
… + coeff[n]*(s - s0))**n / n! / ((speed_of_light * 1e-15)**n * g0))
where coeff[n] is represented in [1/fs**n] The convention for the coeff is introduced for convenient treatment this with respect to a radiation chirp in order to easily satisfy the resonant condition along the whole bunch in the case of linear electron bunch chirp. Here is the expresion:
2*dw/dt = (w0/g0) * dg/dt
@author: Andrei Trebushinin
-
add_wake
(self, tube_radius=0.005, tube_len=1, conductivity=36600000.0, tau=7.1e-15, roughness=6e-07, d_oxid=5e-09)¶
-
abstract
to_array
(self, *args, **kwargs)¶
-
-
class
ocelot.cpbd.beam.
ParticleArray
(n=0)¶ array of particles of fixed size; for optimized performance (x, x’ = px/p0),(y, y’ = py/p0),(ds = c*tau, p = dE/(p0*c)) p0 - momentum
-
class
LostParticleRecorder
(n)¶ Stores information about particles that are getting lost.
- Attributes:
lost_particles: List of indices of deleted particle. Notes: The indices of the initial ParticleArray are used. lp_to_pos_hist: Histogram of number of lost particles to position s. [(pos, num of lost particles), …]
-
add
(self, inds, position)¶
-
rm_tails
(self, xlim, ylim, px_lim, py_lim)¶ Method removes particles outside range [-xlim, +xlim], [-px_lim, +px_lim] …
-
__getitem__
(self, idx)¶
-
__setitem__
(self, idx, p)¶
-
list2array
(self, p_list)¶
-
array2list
(self)¶
-
array2ex_list
(self, p_list)¶
-
size
(self)¶
-
x
(self)¶
-
px
(self)¶
-
y
(self)¶
-
py
(self)¶
-
tau
(self)¶
-
p
(self)¶
-
property
t
(self)¶
-
property
n
(self)¶
-
thin_out
(self, nth=10, n0=0)¶ Method to thin out the particle array in n-th times. Means every n-th particle will be saved in new Particle array
- Parameters
nth – 10, every n-th particle will be taken to new Particle array
n0 – start from n0 particle
- Returns
New ParticleArray
-
rm_particle
(self, index)¶ Method removes a “bad” particle with particular “index”. :param index: :return:
-
rescale2energy
(self, energy)¶ Method to rescale beam coordinates with new energy
- Parameters
energy – new energy
- Returns
-
__str__
(self)¶ Return str(self).
-
delete_particles
(self, inds, record=True)¶ Deletes particles from the particle array via index. :param inds: Indices that will be removed from the particle array. :param record: If record is true the deleted particles will be saved in self.lost_particle_recorder. :return:
-
class
-
ocelot.cpbd.beam.
recalculate_ref_particle
(p_array)¶
-
ocelot.cpbd.beam.
get_envelope
(p_array, tws_i=Twiss(), bounds=None)¶ Function to calculate twiss parameters form the ParticleArray
- Parameters
p_array – ParticleArray
tws_i – optional, design Twiss,
bounds – optional, [left_bound, right_bound] - bounds in units of std(p_array.tau())
- Returns
Twiss()
-
ocelot.cpbd.beam.
get_current
(p_array, num_bins=200, **kwargs)¶ Function calculates beam current from particleArray.
- Parameters
p_array – particleArray
charge –
None, OBSOLETE, charge of the one macro-particle.
If None, charge of the first macro-particle is used
num_bins – number of bins
:return s, I - (np.array, np.array) - beam positions [m] and currents in [A]
-
ocelot.cpbd.beam.
gauss_from_twiss
(emit, beta, alpha)¶
-
ocelot.cpbd.beam.
waterbag_from_twiss
(emit, beta, alpha)¶
-
ocelot.cpbd.beam.
ellipse_from_twiss
(emit, beta, alpha)¶
-
ocelot.cpbd.beam.
moments
(x, y, cut=0)¶
-
ocelot.cpbd.beam.
m_from_twiss
(Tw1, Tw2)¶
-
ocelot.cpbd.beam.
beam_matching
(particles, bounds, x_opt, y_opt, remove_offsets=True)¶ Beam matching function, the beam is centered in the phase space
- Parameters
particles – ParticleArray
bounds – [start, stop] in rms of sigmas in longitudinal direction
x_opt – [alpha, beta, mu (phase advance)]
y_opt – [alpha, beta, mu (phase advance)]
remove_offsets – True, remove offsets in transverse planes
- Returns
transform ParticleArray (the same object)
-
ocelot.cpbd.beam.
sortrows
(x, col)¶
-
ocelot.cpbd.beam.
convmode_py
(A, B, mode)¶
-
ocelot.cpbd.beam.
convmode
¶
-
ocelot.cpbd.beam.
s2cur_auxil_py
(A, xiA, C, N, I)¶
-
ocelot.cpbd.beam.
s2cur_auxil
¶
-
ocelot.cpbd.beam.
s_to_cur
(A, sigma, q0, v)¶ Function to calculate beam current
- Parameters
A – s-coordinates of particles
sigma – smoothing parameter
q0 – bunch charge
v – mean velocity
- Returns
[s, I]
-
ocelot.cpbd.beam.
slice_analysis_py
(z, x, xs, M, to_sort)¶ returns: <x>, <xs>, <x^2>, <x*xs>, <xs^2>, np.sqrt(<x^2> * <xs^2> - <x*xs>^2) based on M particles in moving window Sergey, check please
-
ocelot.cpbd.beam.
slice_analysis
¶
-
ocelot.cpbd.beam.
simple_filter
(x, p, iter)¶
-
ocelot.cpbd.beam.
interp1
(x, y, xnew, k=1)¶
-
ocelot.cpbd.beam.
slice_analysis_transverse
(parray, Mslice, Mcur, p, iter)¶
-
class
ocelot.cpbd.beam.
SliceParameters
¶
-
ocelot.cpbd.beam.
global_slice_analysis_extended
(parray, Mslice, Mcur, p, iter)¶ Function to calculate slice parameters
- Parameters
parray – ParticleArray
Mslice – 5000, nparticles in the slice
Mcur – 0.01, smoothing parameters to calculate the beam current: smooth_param = m_std * np.std(p_array.tau())
p – 2, filter parameter in the func: simple_filter
iter – 2, filter parameter in the func: simple_filter
- Returns
s, I, ex, ey, me, se, gamma0, emitxn, emityn
-
ocelot.cpbd.beam.
global_slice_analysis
(parray, nparts_in_slice=5000, smooth_param=0.01, filter_base=2, filter_iter=2)¶ Function to calculate slice parameters
- Parameters
parray – ParticleArray
nparts_in_slice – 5000, nparticles in the slice
smooth_param – 0.01, smoothing parameters to calculate the beam current: smooth_param = m_std * np.std(p_array.tau())
filter_base – 2, filter parameter in the func: simple_filter
filter_iter – 2, filter parameter in the func: simple_filter
- Returns
SliceParameters,
-
ocelot.cpbd.beam.
parray2beam
(parray, step=1e-07)¶ reads ParticleArray() returns BeamArray() step [m] - long. size ob bin to calculate distribution parameters
-
ocelot.cpbd.beam.
generate_parray
(sigma_x=0.0001, sigma_px=2e-05, sigma_y=None, sigma_py=None, sigma_tau=0.001, sigma_p=0.0001, chirp=0.01, charge=5e-09, nparticles=200000, energy=0.13, tau_trunc=None, tws=None)¶ Method to generate ParticleArray with gaussian distribution.
- Note: in ParticleArray, {x, px} and {y, py} are canonical coordinates. {tau, p} is not, to make it canonical
the sign of “tau” should be flipped.
- Parameters
sigma_x – std(x), x is horizontal cartesian coordinate.
sigma_px – std(px), ‘px’ is conjugate momentum canonical momentum px/p0.
sigma_y – std(y), y is vertical cartesian coordinate.
sigma_py – std(py), ‘py’ is canonical momentum py/p0.
sigma_tau – std(tau), “tau” = c*t
sigma_p – std(p), ‘p’ is canonical momentum E/(c*p0)
chirp – energy chirp [unitless], linear correlation - p_i += chirp * tau_i/sigma_tau
charge – beam charge in [C], 5e-9 by default
nparticles – namber of particles, 200k by default
energy – beam energy in [GeV], 0.13 [GeV]
tau_trunc – None, if not [float] - truncated gauss distribution in “tau” direction.
tws – None, if Twiss obj - the beam is matched to twiss params.
- Returns
ParticleArray
-
ocelot.cpbd.beam.
generate_beam
(E, I=5000, l_beam=3e-06, **kwargs)¶ Generates BeamArray object accepts arguments with the same names as BeamArray().parameters() I - current in Amps E - beam ebergy in GeV
dE - rms energy spread in GeV emit_x, emit_n(both normalized), emit_xn, etc. shape - beam shape (‘gaussian’ of ‘flattop’) l_beam [m] - beam length in meters l_window [m] - window length in um
- by default: l_beam * 2 if flattop,
l_beam * 6 if gaussian,
nslice - number of slices in the beam