Laser Fields

Objects for creating laser beam objects, along with some common components.

Overview

laserBeam([kvec, s, pol, delta, phase, ...])

The base class for a single laser beam

infinitePlaneWaveBeam(kvec, pol, s, delta, ...)

Infinte plane wave beam

gaussianBeam(kvec, pol, s, delta, wb, **kwargs)

Collimated Gaussian beam

clippedGaussianBeam(kvec, pol, s, delta, wb, ...)

Clipped, collimated Gaussian beam

laserBeams([laserbeamparams, beam_type])

The base class for a collection of laser beams

conventional3DMOTBeams([k, pol, ...])

A collection of laser beams for 6-beam MOT

Details

class pylcp.laserBeam(kvec=None, s=None, pol=None, delta=None, phase=0.0, pol_coord='spherical', eps=1e-05)[source]

The base class for a single laser beam

Attempts to represent a laser beam as

\[\frac{1}{2}\hat{\boldsymbol{\epsilon}}(r, t) E_0(r, t) e^{i\mathbf{k}(r,t)\cdot\mathbf{r}-i \int dt\Delta(t) + i\phi(r, t)}\]

where \(\hat{\boldsymbol{\epsilon}}\) is the polarization, \(E_0\) is the electric field magnitude, \(\mathbf{k}(r,t)\) is the k-vector, \(\mathbf{r}\) is the position, \(\Delta(t)\) is the deutning, \(t\) is the time, and \(\phi\) is the phase.

Parameters
  • kvec (array_like with shape (3,) or callable) – The k-vector of the laser beam, specified as either a three-element list or numpy array or as callable function. If a callable, it must have a signature like (R, t), (R), or (t) where R is an array_like with shape (3,) and t is a float and it must return an array_like with three elements.

  • pol (int, float, array_like with shape (3,), or callable) – The polarization of the laser beam, specified as either an integer, float array_like with shape(3,), or as callable function. If an integer or float, if pol<0 the polarization will be left circular polarized relative to the k-vector of the light. If pol>0, the polarization will be right circular polarized. If array_like, polarization will be specified by the vector, whose basis is specified by pol_coord. If a callable, it must have a signature like (R, t), (R), or (t) where R is an array_like with shape (3,) and t is a float and it must return an array_like with three elements.

  • s (float or callable) – The intensity of the laser beam, normalized to the saturation intensity, specified as either a float or as callable function. If a callable, it must have a signature like (R, t), (R), or (t) where R is an array_like with shape (3,) and t is a float and it must return a float.

  • delta (float or callable) – Detuning of the laser beam. If a callable, it must have a signature like (t) where t is a float and it must return a float.

  • phase (float, optional) – Phase of laser beam. By default, zero.

  • pol_coord (string, optional) – Polarization basis of the input polarization vector: ‘cartesian’ or ‘spherical’ (default).

  • eps (float, optional) – Small distance to use in calculation of numerical derivatives. By default eps=1e-5.

eps

Small epsilon used for computing derivatives

Type

float

phase

Overall phase of the laser beam.

Type

float

cartesian_pol(R=array([0.0, 0.0, 0.0]), t=0)[source]

Returns the polarization in Cartesian coordinates.

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

pol – polarization of the laser beam at R and t in Cartesian basis.

Return type

array_like, size (3,)

delta(t=0.0)[source]

Returns the detuning of the laser beam at time t

Parameters

t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

delta – detuning of the laser beam at time t

Return type

float or array like

electric_field(R, t)[source]

The electric field at position R and t

Parameters
  • R (array_like, size (3,)) – vector of the position at which to return the kvector. By default, the origin.

  • t (float) – time at which to return the k-vector. By default, t=0.

Returns

Eq – electric field in the spherical basis.

Return type

array_like, shape (3,)

electric_field_gradient(R, t)[source]

The full derivative of electric field at position R and t

Parameters
  • R (array_like, size (3,)) – vector of the position at which to return the kvector. By default, the origin.

  • t (float) – time at which to return the k-vector. By default, t=0.

Returns

dEq – The full gradient of the electric field, in spherical coordinates.

\[\begin{split}\begin{pmatrix} \frac{dE_{-1}}{dx} & \frac{dE_0}{dx} & \frac{dE_{+1}}{dx} \\ \frac{dE_{-1}}{dy} & \frac{dE_0}{dy} & \frac{dE_{+1}}{dy} \\ \frac{dE_{-1}}{dz} & \frac{dE_0}{dz} & \frac{dE_{+1}}{dz} \\ \end{pmatrix}\end{split}\]

Return type

array_like, shape (3, 3)

intensity(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the intensity of the laser beam at position R and t

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

s – Saturation parameter of the laser beam at R and t.

Return type

float or array_like

jones_vector(xp, yp, R=array([0.0, 0.0, 0.0]), t=0)[source]

Returns the Jones vector at position

Parameters
  • xp (array_like, shape (3,)) – The x vector of the basis in which to calculate the Jones vector. Must be orthogonal to k.

  • yp (array_like, shape (3,)) – The y vector of the basis in which to calculate the Jones vector. Must be orthogonal to k and xp.

  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

pol – Jones vector of the laser beam at R and t in Cartesian basis.

Return type

array_like, size (2,)

kvec(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the k-vector of the laser beam

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

kvec – the k vector at position R and time t.

Return type

array_like, size(3,)

pol(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the polarization of the laser beam at position R and t

The polarization is returned in the spherical basis.

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

pol – polarization of the laser beam at R and t in spherical basis.

Return type

array_like, size (3,)

polarization_ellipse(xp, yp, R=array([0.0, 0.0, 0.0]), t=0)[source]

The polarization ellipse parameters of the laser beam at R and t

Parameters
  • xp (array_like, shape (3,)) – The x vector of the basis in which to calculate the polarization ellipse. Must be orthogonal to k.

  • yp (array_like, shape (3,)) – The y vector of the basis in which to calculate the polarization ellipse. Must be orthogonal to k and xp.

  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

  • psi (float) – \(\psi\) parameter of the polarization ellipse

  • chi (float) – \(\chi\) parameter of the polarization ellipse

project_pol(quant_axis, R=array([0.0, 0.0, 0.0]), t=0, treat_nans=False, calculate_norm=False, invert=False)[source]

Project the polarization onto a quantization axis.

Parameters
  • quant_axis (array_like, shape (3,)) – A normalized 3-vector of the quantization axis direction.

  • R (array_like, shape (3,), optional) – If polarization is a function of R is the 3-vectors at which the polarization shall be calculated.

  • calculate_norm (bool, optional) – If true, renormalizes the quant_axis. By default, False.

  • treat_nans (bool, optional) – If true, every place that nan is encoutnered, replace with the $hat{z}$ axis as the quantization axis. By default, False.

  • invert (bool, optional) – If true, invert the process to project the quantization axis onto the specified polarization.

Returns

projected_pol – The polarization projected onto the quantization axis.

Return type

array_like, shape (3,)

stokes_parameters(xp, yp, R=array([0.0, 0.0, 0.0]), t=0)[source]

The Stokes Parameters of the laser beam at R and t

Parameters
  • xp (array_like, shape (3,)) – The x vector of the basis in which to calculate the Stokes parameters. Must be orthogonal to k.

  • yp (array_like, shape (3,)) – The y vector of the basis in which to calculate the Stokes parameters. Must be orthogonal to k and xp.

  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

pol – Stokes parameters for the laser beam, [Q, U, V]

Return type

array_like, shape (3,)

class pylcp.infinitePlaneWaveBeam(kvec, pol, s, delta, **kwargs)[source]

Infinte plane wave beam

A beam which has spatially constant intensity, k-vector, and polarization.

\[\frac{1}{2}\hat{\boldsymbol{\epsilon}} E_0e^{i\mathbf{k}\cdot\mathbf{r}-i \int dt\Delta(t) + i\phi(r, t)}\]

where \(\hat{\boldsymbol{\epsilon}}\) is the polarization, \(E_0\) is the electric field magnitude, \(\mathbf{k}(r,t)\) is the k-vector, \(\mathbf{r}\) is the position, \(\Delta(t)\) is the deutning, \(t\) is the time, and \(\phi\) is the phase.

Parameters
  • kvec (array_like with shape (3,) or callable) – The k-vector of the laser beam, specified as either a three-element list or numpy array.

  • pol (int, float, array_like with shape (3,), or callable) – The polarization of the laser beam, specified as either an integer, float array_like with shape(3,). If an integer or float, if pol<0 the polarization will be left circular polarized relative to the k-vector of the light. If pol>0, the polarization will be right circular polarized. If array_like, polarization will be specified by the vector, whose basis is specified by pol_coord.

  • s (float or callable) – The intensity of the laser beam, specified as either a float or as callable function.

  • delta (float or callable) – Detuning of the laser beam. If a callable, it must have a signature like (t) where t is a float and it must return a float.

  • **kwargs – Additional keyword arguments to pass to laserBeam superclass.

Notes

This implementation is much faster, when it can be used, compared to the base laserBeam class.

electric_field_gradient(R, t)[source]

The full derivative of electric field at position R and t

Parameters
  • R (array_like, size (3,)) – vector of the position at which to return the kvector. By default, the origin.

  • t (float) – time at which to return the k-vector. By default, t=0.

Returns

dEq – The full gradient of the electric field, in spherical coordinates.

\[\begin{split}\begin{pmatrix} \frac{dE_{-1}}{dx} & \frac{dE_0}{dx} & \frac{dE_{+1}}{dx} \\ \frac{dE_{-1}}{dy} & \frac{dE_0}{dy} & \frac{dE_{+1}}{dy} \\ \frac{dE_{-1}}{dz} & \frac{dE_0}{dz} & \frac{dE_{+1}}{dz} \\ \end{pmatrix}\end{split}\]

Return type

array_like, shape (3, 3)

class pylcp.gaussianBeam(kvec, pol, s, delta, wb, **kwargs)[source]

Collimated Gaussian beam

A beam which has spatially constant k-vector and polarization, with a Gaussian intensity modulation. Specifically,

\[\frac{1}{2}\hat{\boldsymbol{\epsilon}} E_0 e^{-\mathbf{r}^2/w_b^2} e^{i\mathbf{k}\cdot\mathbf{r}-i \int dt\Delta(t) + i\phi(r, t)}\]

where \(\hat{\boldsymbol{\epsilon}}\) is the polarization, \(E_0\) is the electric field magnitude, \(\mathbf{k}(r,t)\) is the k-vector, \(\mathbf{r}\) is the position, \(\Delta(t)\) is the deutning, \(t\) is the time, and \(\phi\) is the phase. Note that because \(I\propto E^2\), \(w_b\) is the \(1/e^2\) radius.

Parameters
  • kvec (array_like with shape (3,) or callable) – The k-vector of the laser beam, specified as either a three-element list or numpy array.

  • pol (int, float, array_like with shape (3,), or callable) – The polarization of the laser beam, specified as either an integer, float array_like with shape(3,). If an integer or float, if pol<0 the polarization will be left circular polarized relative to the k-vector of the light. If pol>0, the polarization will be right circular polarized. If array_like, polarization will be specified by the vector, whose basis is specified by pol_coord.

  • s (float or callable) – The maximum intensity of the laser beam at the center, specified as either a float or as callable function.

  • delta (float or callable) – Detuning of the laser beam. If a callable, it must have a signature like (t) where t is a float and it must return a float.

  • wb (float) – The \(1/e^2\) radius of the beam.

  • **kwargs – Additional keyword arguments to pass to the laserBeam superclass.

intensity(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the intensity of the laser beam at position R and t

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

s – Saturation parameter of the laser beam at R and t.

Return type

float or array_like

class pylcp.clippedGaussianBeam(kvec, pol, s, delta, wb, rs, **kwargs)[source]

Clipped, collimated Gaussian beam

A beam which has spatially constant k-vector and polarization, with a Gaussian intensity modulation. Specifically,

\[\frac{1}{2}\hat{\boldsymbol{\epsilon}} E_0 e^{-\mathbf{r}^2/w_b^2} (|\mathbf{r}|<r_s) e^{i\mathbf{k}\cdot\mathbf{r}-i \int dt\Delta(t) + i\phi(r, t)}\]

where \(\hat{\boldsymbol{\epsilon}}\) is the polarization, \(E_0\) is the electric field magnitude, \(r_s\) is the radius of the stop, \(\mathbf{k}(r,t)\) is the k-vector, \(\mathbf{r}\) is the position, \(\Delta(t)\) is the deutning, \(t\) is the time, and \(\phi\) is the phase. Note that because \(I\propto E^2\), \(w_b\) is the \(1/e^2\) radius.

Parameters
  • kvec (array_like with shape (3,) or callable) – The k-vector of the laser beam, specified as either a three-element list or numpy array.

  • pol (int, float, array_like with shape (3,), or callable) – The polarization of the laser beam, specified as either an integer, float array_like with shape(3,). If an integer or float, if pol<0 the polarization will be left circular polarized relative to the k-vector of the light. If pol>0, the polarization will be right circular polarized. If array_like, polarization will be specified by the vector, whose basis is specified by pol_coord.

  • s (float or callable) – The maximum intensity of the laser beam at the center, specified as either a float or as callable function.

  • delta (float or callable) – Detuning of the laser beam. If a callable, it must have a signature like (t) where t is a float and it must return a float.

  • wb (float) – The \(1/e^2\) radius of the beam.

  • rs (float) – The radius of the stop.

  • **kwargs – Additional keyword arguments to pass to the laserBeam superclass.

intensity(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the intensity of the laser beam at position R and t

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

s – Saturation parameter of the laser beam at R and t.

Return type

float or array_like

class pylcp.laserBeams(laserbeamparams=None, beam_type=<class 'pylcp.fields.laserBeam'>)[source]

The base class for a collection of laser beams

Parameters
  • laserbeamparams (array_like of laserBeam or array_like of dictionaries) – If array_like contains laserBeams, the laserBeams in the array will be joined together to form a collection. If array_like is a list of dictionaries, the dictionaries will be passed as keyword arguments to beam_type

  • beam_type (laserBeam or laserBeam subclass, optional) – Type of beam to use in the collection of laserBeams. By default beam_type=laserBeam.

add_laser(new_laser)[source]

Add a laser to the collection

Parameters

new_laser (laserBeam or laserBeam subclass) –

cartesian_pol(R=array([0.0, 0.0, 0.0]), t=0)[source]

Returns the polarization of all laser beams in Cartesian coordinates.

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the polarization. By default, the origin.

  • t (float, optional) – time at which to return the polarization. By default, t=0.

Returns

pol – polarization of the laser beam at R and t in Cartesian basis.

Return type

array_like, shape (num_of_beams, 3)

delta(t=0)[source]

Returns the detuning of the laser beam at time t

Parameters

t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

delta – detuning of the laser beam at time t for all laser beams

Return type

float or array like

electric_field(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the electric field of the laser beams

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

E – the electric field vectors at position R and time t for each laser beam.

Return type

list of array_like, size(3,)

electric_field_gradient(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the gradient of the electric field of the laser beams

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

dE – the electric field gradient matrices at position R and time t for each laser beam.

Return type

list of array_like, size(3,)

intensity(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the intensity of the laser beam at position R and t

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

s – Saturation parameters of all laser beams at R and t.

Return type

list of float or array_like

jones_vector(xp, yp, R=array([0.0, 0.0, 0.0]), t=0)[source]

Jones vector at position R and time t

Parameters
  • xp (array_like, shape (3,)) – The x vector of the basis in which to calculate the Jones vector. Must be orthogonal to k.

  • yp (array_like, shape (3,)) – The y vector of the basis in which to calculate the Jones vector. Must be orthogonal to k and xp.

  • R (array_like, size (3,), optional) – vector of the position at which to evaluate the Jones vector. By default, the origin.

  • t (float, optional) – time at which to evaluate the Jones vector. By default, t=0.

Returns

pol – Jones vector of the laser beams at R and t in Cartesian basis.

Return type

array_like, size (num_of_beams, 2)

kvec(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the k-vector of the laser beam

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

kvec – the k vector at position R and time t for each laser beam.

Return type

list of array_like, size(3,)

pol(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the polarization of the laser beam at position R and t

The polarization is returned in the spherical basis.

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

pol – polarization of each laser beam at R and t in spherical basis.

Return type

list of array_like, size (3,)

polarization_ellipse(xp, yp, R=array([0.0, 0.0, 0.0]), t=0)[source]

The polarization ellipse parameters of the laser beam at R and t

Parameters
  • xp (array_like, shape (3,)) – The x vector of the basis in which to calculate the polarization ellipse. Must be orthogonal to k.

  • yp (array_like, shape (3,)) – The y vector of the basis in which to calculate the polarization ellipse. Must be orthogonal to k and xp.

  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

list of (psi, chi) – list of (\(\psi\), \(\chi\)) parameters of the polarization ellipses for each laser beam

Return type

list of tuples

project_pol(quant_axis, R=array([0.0, 0.0, 0.0]), t=0, **kwargs)[source]

Project the polarization onto a quantization axis.

Parameters
  • quant_axis (array_like, shape (3,)) – A normalized 3-vector of the quantization axis direction.

  • R (array_like, shape (3,), optional) – If polarization is a function of R is the 3-vectors at which the polarization shall be calculated.

  • calculate_norm (bool, optional) – If true, renormalizes the quant_axis. By default, False.

  • treat_nans (bool, optional) – If true, every place that nan is encoutnered, replace with the $hat{z}$ axis as the quantization axis. By default, False.

  • invert (bool, optional) – If true, invert the process to project the quantization axis onto the specified polarization.

Returns

projected_pol – The polarization projected onto the quantization axis for all laser beams

Return type

list of array_like, shape (3,)

stokes_parameters(xp, yp, R=array([0.0, 0.0, 0.0]), t=0)[source]

The Stokes Parameters of the laser beam at R and t

Parameters
  • xp (array_like, shape (3,)) – The x vector of the basis in which to calculate the Stokes parameters. Must be orthogonal to k.

  • yp (array_like, shape (3,)) – The y vector of the basis in which to calculate the Stokes parameters. Must be orthogonal to k and xp.

  • R (array_like, size (3,), optional) – vector of the position at which to calculate the Stokes parameters. By default, the origin.

  • t (float, optional) – time at which to calculate the Stokes parameters. By default, t=0.

Returns

pol – Stokes parameters for the laser beams, [Q, U, V]

Return type

array_like, shape (num_of_beams, 3)

total_electric_field(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the total electric field of the laser beams

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

E – the total electric field vector at position R and time t of all the laser beams

Return type

array_like, size(3,)

total_electric_field_gradient(R=array([0.0, 0.0, 0.0]), t=0.0)[source]

Returns the total gradient of the electric field of the laser beams

Parameters
  • R (array_like, size (3,), optional) – vector of the position at which to return the kvector. By default, the origin.

  • t (float, optional) – time at which to return the k-vector. By default, t=0.

Returns

dE – the total electric field gradient matrices at position R and time t of all laser beams.

Return type

array_like, size(3,)

class pylcp.conventional3DMOTBeams(k=1, pol=1, rotation_angles=[0.0, 0.0, 0.0], rotation_spec='ZYZ', beam_type=<class 'pylcp.fields.laserBeam'>, **kwargs)[source]

A collection of laser beams for 6-beam MOT

The standard geometry is to generate counter-progagating beams along all orthogonal axes \((\hat{x}, \hat{y}, \hat{z})\).

Parameters
  • k (float, optional) – Magnitude of the k-vector for the six laser beams. Default: 1

  • pol (int or float, optional) – Sign of the circular polarization for the beams moving along \(\hat{z}\). Default: +1. Orthogonal beams have opposite polarization by default.

  • rotation_angles (array_like) – List of angles to define a rotated MOT. Default: [0., 0., 0.]

  • rotation_spec (str) – String to define the convention of the Euler rotations. Default: ‘ZYZ’

  • beam_type (pylcp.laserBeam or subclass) – Type of beam to generate.

  • **kwargs – other keyword arguments to pass to beam_type