iter Package

Package containing algorithms and cost functions for iterative optimisation. This section is split into two parts, a description of the optimisation methods and a description of the objective function classes.

Iterative optimisation methods

The package contains a series of iterative optimisation algorithms. There are currently two types of iterative optimisation algorithms, methods that attempt to approximate some target far-field and methods which attempt to combine a series of SLM patterns. The former can also be used to combine a series of patterns by first generating a target far-field, for example using otslm.tools.combine():

lin = otslm.simple.linear(sz, 10);
lg = otslm.simple.linear(sz, -10) + otslm.simple.lgmode(sz, 3, 0);

% Convert to complex amplitudes (could use finalize)
lin = exp(1i*2*pi*lin);
lg = exp(1i*2*pi*lg);

target = otslm.tools.combine({lin, lg}, 'method', 'farfield');

The methods which inherit from IterBase have an objective function property. For some methods the objective is required for the method to work, for other methods the objective is optional and can be used to track progress of the method. The objective can be set on construction or by setting the objective property. See the Objective functions section for available objectives. For an example of how to use these iterative methods, see examples.iterative, examples.iter_combine and Gerchberg-Saxton examples. A minimal example for methods which attempt to generate a particular far-field is shown below:

sz = [512, 512];
incident = ones(sz);

prop = otslm.tools.prop.FftForward.simpleProp(zeros(sz));
vismethod = @(U) prop.propagate(U .* incident);

target = otslm.simple.aperture(sz, sz(1)/20);
gs = otslm.iter.GerchbergSaxton(target, 'adaptive', 1.0, ...
    'vismethod', vismethod);

Table 1 compares the run-time and required number of iterations for some of the iterative optimisation methods. This table is based on Di Leonardo et al. 2007, a more detailed discussion can be found in the reference. This is only a guide, some methods may work better than other methods under certain circumstances. For instance, the direct search method can be used for fine tuning the output of other methods but takes too long for practical use when given a bad initial guess. The combination algorithm and 2-D optimisation algorithms have been combined, actual performance will be different but similar. There are a range of different extensions to the described methods which may improve performance for particular problems, such as using a super-pixel style approach with the Direct Search algorithm.

Table 1 Comparison of iterative methods
Iterative methods Num. Iterations Typical Run-time
Gerchberg-Saxton 30 5 s
Weighted GS 30 5 s
Adaptive-adaptive 30 5 s
Bowman 2017 < 200 2 m
Simulated Annealing \(10^4\) 10 m
Direct Search \(10^9\) days

GerchbergSaxton

The Gerchberg-Saxton algorithm is an iterative algorithm that involves iterating between the near-field and far-field and applying constraints to the fields after each iteration. The constraints could include a particular incident illumination or a desired far-field intensity or phase pattern. Components that are not constrained are free to change. The algorithm was originally described in

R. W. Gerchberg, O. A Saxton W., A practical algorithm for the determination of phase from image and diffraction plane pictures, Optik 35(1971) 237-250 (Nov 1971).

Details about the algorithm can be found on the Wikipedia page. A sketch of the algorithm for generating a target amplitude pattern using a phase-only device is shown below:

  1. Generate initial guess for the SLM phase pattern \(P\).
  2. Calculate output for phase pattern: \(\text{Proj}(P) \rightarrow O\).
  3. Multiply output phase by target amplitude: \(|T|\frac{O}{|O|} \rightarrow Q\).
  4. Calculate the complex amplitude required to generate \(Q\): \(\text{Inv}(Q) \rightarrow I\).
  5. Calculate new guess from the phase of \(I\): \(\text{Angle}(I) \rightarrow P\).
  6. Goto step 2 until converged.

\(\text{Proj}\) and \(\text{Inv}\) are the forward and inverse propagation methods, these could be, for example, the forward and inverse Fourier transforms. A constraint for the incident illumination can be added to the forward propagator or the constraint can be added at another step. There are other variants for generating a target phase field or applying other constraints to the far-field. If this guess is symmetric, these symmetries will influence the final output, this can be useful for generating symmetric target fields.

GerchbergSaxton also implements the adaptive-adaptive algorithm, which we can enable by setting the adaptive parameter to a non-unity value. The adaptive-adaptive algorithm is similar to the above except step 3 mixes the propagator amplitude output with the target amplitude instead of replacing it

\[ \begin{align}\begin{aligned}t = \alpha |T| + (1 - \alpha) |O|\\Q = t \frac{O}{|O|}\end{aligned}\end{align} \]

where \(\alpha\) is the adaptive-adaptive factor.

class otslm.iter.GerchbergSaxton(target, varargin)

Implementation of Gerchberg-Saxton and Adaptive-Adaptive algorithms Inherits from IterBase.

Methods
  • run() – Run the iterative method
Properties
  • adaptive – Adaptive-adaptive factor (1 for Gerchberg-Saxton)
Inherited properties
  • guess – Best guess at hologram pattern
  • target – Target pattern the method tries to approximate
  • vismethod – Method used to do the visualisation
  • invmethod – Method used to calculate initial guess/inverse-visualisation
  • objective – Objective function used to evaluate fitness
  • fitness – Fitness evaluated after every iteration

See also GerchbergSaxton.

GerchbergSaxton(target, varargin)

Construct a new instance of the GerchbergSaxton iterative method

Usage
mtd = GerchbergSaxton(target, …)
Parameters
  • target – target pattern
Optional named arguments
  • adaptive num Adaptive-Adaptive factor. Default: 1.0, i.e. the method is Gerchberg-Saxton.
  • guess im Initial guess at complex amplitude pattern. If not image is supplied, a guess is created using invmethod.
  • vismethod fcn Function to calculate far-field. Takes one argument: the complex amplitude near-field. Default: @otslm.tools.prop.FftForward.simpleProp.propagate
  • invmethod fcn Function to calculate near-field. Takes one argument: the complex amplitude far-field. Default: @otslm.tools.prop.FftInverse.simpleProp.propagate
  • objective fcn Objective function to measure fitness. Default: @otslm.iter.objectives.FlatIntensity

DirectSearch

The direct search algorithm involves choosing a pixel, trying a range of possible values for that pixel, and keeping the choice which maximises some objective function. This is a expensive procedure, on a device with 512x512 pixels and 256 values per pixel, cycling over each pixel requires 67 million Fourier transforms. The process is further complicated since the optimal value for any pixel is not independent of every other pixel. However, this method can be useful for further improving a good guess, such as the output of one of the other methods.

A rough outline for the procedure is

  1. Choose an initial guess, \(P\)
  2. Randomly select a pixel to modify
  3. Generate a set of patterns \(P_i\) with a set \(\{i\}\) of different pixel values.
  4. Propagate these patterns and calculate the fitness \(F_i\)
  5. Choose the pattern which maximises the fitness \(P_j \rightarrow P\) where \(j = \text{argmax}_i F_i\).
  6. Go to 2 until converged
class otslm.iter.DirectSearch(target, varargin)

Optimiser to search through each pixel value to optimise hologram Inherits from IterBase.

This method randomly selects a pixel in the pattern and then tries every available level. The pixel value kept is the pixel value whic gives the best fitness.

The algorithm is described in Di Leonardo, et al., Opt. Express 15, 1913-1922 (2007)

Methods
  • run() – Run the iterative method
Properties
  • levels – Discrete levels that will be search in optimisation
Inherited properties
  • guess – Best guess at hologram pattern
  • target – Target pattern the method tries to approximate
  • vismethod – Method used to do the visualisation
  • invmethod – Method used to calculate initial guess/inverse-visualisation
  • objective – Objective function used to evaluate fitness
  • fitness – Fitness evaluated after every iteration

See also DirectSearch and SimulatedAnnealing.

DirectSearch(target, varargin)

Construct a new instance of the DirectSearch iterative method

Usage
mtd = DirectSearch(target, …) attempts to produce the target using the Direct Search algorithm.
Optional named arguments:
  • levels num – Number of discrete phase levels or array of levels between -pi and pi. Default: 256.
  • guess im – Initial guess at complex amplitude pattern. If not image is supplied, a guess is created using invmethod.
  • vismethod fcn – Function to calculate far-field. Takes one argument: the complex amplitude near-field. Default: @otslm.tools.prop.FftForward.simpleProp.propagate
  • invmethod fcn – Function to calculate near-field. Takes one argument: the complex amplitude far-field. Default: @otslm.tools.prop.FftInverse.simpleProp.propagate
  • objective fcn – Objective function to measure fitness. Default: @otslm.iter.objectives.FlatIntensity

SimulatedAnnealing

Simulated annealing is a stochastic method that can be useful for optimising systems with many degrees of freedom (such as patterns with many non-independent pixels). A description of the method can be found on the wikipedia page. The algorithm is analogous to cooling (annealing) of solids and chooses new state probabilistically depending on a temperature parameter. An outline follows

  1. Starting with an initial pattern \(P\) and temperature \(T\)
  2. Pick a random pattern which is similar to the current pattern
  3. Compare fitness of two patterns \(F_1\) and \(F_2\)
  4. Accept the new pattern if \(P(F_1, F_2, T) > \text{rand}(0, 1)\)
  5. Goto 2 until converged, gradually reducing temperature

There are several parameters that can be chosen which strongly affect the performance and convergence of the algorithm. The implementation currently only supports the following function

\[P(F_1, F_2, T) = \exp{-(F_2-F_1)/T}\]

The change in temperature can be controlled via the temperatureFcn optional parameter.

This implementation could be improved and we welcome suggestions.

class otslm.iter.SimulatedAnnealing(target, varargin)

Optimise the pattern using simulated annealing. Inherits from IterBase.

Methods
  • run() – Run the iterative method
Properties
  • levels – Discrete levels that will be search in optimisation
  • temperature – Current temperature of the system
  • maxTemperature – Scaling factor for new pattern guesses
  • temperatureFcn – Function used to calculate temperature in iteration
  • lastFitness – The fitness associated with the current guess
  • guess – Best guess at hologram pattern
  • target – Target pattern the method tries to approximate
  • vismethod – Method used to do the visualisation
  • invmethod – Method used to calculate initial guess/inverse-visualisation
  • objective – Objective function used to evaluate fitness
  • fitness – Fitness evaluated after every iteration

See also SimulatedAnnealing

SimulatedAnnealing(target, varargin)

Construct a new instance of the SimulatedAnnealing iterative method

mtd = SimulatedAnnealing(target, …) attempts to produce the target using the Simulated Annealing algorithm.

Optional named arguments:
  • levels num Number of discrete levels or array of levels between -pi and pi. Default: 256.
  • temperature num Initial temperature of the solver.
  • guess im Initial guess at complex amplitude pattern. If not image is supplied, a guess is created using invmethod.
  • vismethod fcn Function to calculate far-field. Takes one argument: the complex amplitude near-field. Default: @otslm.tools.prop.FftForward.simpleProp.propagate
  • invmethod fcn Function to calculate near-field. Takes one argument: the complex amplitude far-field. Default: @otslm.tools.prop.FftInverse.simpleProp.propagate
  • objective fcn Objective function to measure fitness. Default: @otslm.iter.objectives.FlatIntensity
static simpleTemperatureFcn(scale, decay)

Returns a exponentially decaying temperature function

fcn = simpleTemperatureFcn(scale, decay) creates a exponentially decaying temperature function. Scale is the initial temperature and decay is the exponential decay rate.

GerchbergSaxton3d

This function implements the 3-D analog of the Gerchberg-Saxton method. The method is described in

Hao Chen et al 2013 J. Opt. 15 035401

and

Graeme Whyte and Johannes Courtial 2005 New J. Phys. 7 117

For an outline of the Gerchberg-Saxton algorithm, see GerchbergSaxton.

class otslm.iter.GerchbergSaxton3d(target, varargin)

Implementation of 3-D Gerchberg-Saxton and Adaptive-Adaptive algorithms Inherits from GerchbergSaxton and IterBaseEwald.

This algorithm attempts to recreate the target volume using the 3-D analog of the Gerchberg-Saxton algorithm.

See Hao Chen et al 2013 J. Opt. 15 035401 and Graeme Whyte and Johannes Courtial 2005 New J. Phys. 7 117

Methods
  • run() – Run the iterative method
Properties
  • adaptive – Adaptive-adaptive factor (1 for Gerchberg-Saxton)
Inherited properties
  • guess – Best guess at hologram pattern
  • target – Target pattern the method tries to approximate
  • vismethod – Method used to do the visualisation
  • invmethod – Method used to calculate initial guess/inverse-visualisation
  • objective – Objective function used to evaluate fitness
  • fitness – Fitness evaluated after every iteration

See also GerchbergSaxton3d and GerchbergSaxton.

GerchbergSaxton3d(target, varargin)

Construct a new instance of the GerchbergSaxton3d iterative method

USage
mtd = GerchbergSaxton3d(target, …)
Parameters
  • target – target pattern to try and generate
Optional named arguments
  • adaptive num – Adaptive-Adaptive factor. Default: 1.0, i.e. the method is Gerchberg-Saxton.
  • guess im – Initial guess at complex amplitude pattern. If not image is supplied, a guess is created using invmethod.
  • vismethod fcn – Function to calculate far-field. Takes one argument: the complex amplitude near-field. Default: @otslm.tools.prop.FftEwaldForward.simpleProp.propagate
  • invmethod fcn – Function to calculate near-field. Takes one argument: the complex amplitude far-field. Default: @otslm.tools.prop.FftEwaldInverse.simpleProp.propagate
  • objective fcn – Optional objective function to measure fitness. Default: @otslm.iter.objectives.FlatIntensity

CombineGerchbergSaxton

This function implements the Gerchberg-Saxton algorithm and similar iterative optimisers for generating point traps. The method can be used to combine a set of SLM patterns \(\phi_m\) into a single pattern in a similar way to otslm.tools.combine(). Starting with an initial guess at the phase pattern \(\phi^0\) the method proceeds as

\[\phi^{j+1} = \sum_n e^{i \phi_n} \eta_n^j \frac{V_n^j}{|V_n^j|}\]

where

\[V_m^j = \sum_{x,y} e^{i (\phi^j(x, y) - \phi_m(x, y))}\]

and \(x, y\) are the SLM pixel coordinates and \(\eta_n^j\) is an optional parameter for Adaptive-Adaptive or weighted versions of the algorithm (for Gerchberg-Saxton \(\eta = 1\)). To calculate the pattern we simply need to iterative the above equation for a few steps.

There are two relatively simple extensions to this algorithm. First is the Adaptive-Adaptive algorithm which involves setting

\[\eta = \alpha + \frac{1-\alpha}{|V_n^j|}\]

where \(\alpha\) is a factor between 0 and 1. The second extension is the weighted Gerchberg-Saxton algorithm which involves setting

\[\eta^{j+1} = \eta^j \frac{\langle V_n^j \rangle}{V_n^j}\]

where \(\langle \cdot \rangle\) denotes the average and we re-calculate \(\eta\) at each iteration starting with an initial value of 1.

To use the method we need to pass in a set of patterns to combine. For instance, we could have a set of 2 traps:

lin1 = otslm.simple.linear(sz, 10);
lin2 = otslm.simple.linear(sz, -5);

components = zeros([sz, 2]);
components(:, :, 1) = lin1;
components(:, :, 2) = lin1;

Then to use the iterative method we would run

mtd = otslm.iter.CombineGerchbergSaxton(2*pi*components, ...
   'weighted', true, 'adaptive', 1.0);
mtd.run(10);
imagesc(mtd.phase);

For a more complete example see examples.iter_combine. A more detailed discussion of these algorithms can be found in

R. Di Leonardo, et al., Opt. Express 15 (4) (2007) 1913-1922. https://doi.org/10.1364/OE.15.001913
class otslm.iter.CombineGerchbergSaxton(components, varargin)

Implementation of Gerchberg-Saxton type combination algorithms. Inherits from IterCombine.

This includes Gerchberg-Saxton, Adaptive-Adaptive and weighted GerchbergSaxton algorithms.

For details about these algorithms, see R. Di Leonardo, et al., Opt. Express 15 (4) (2007) 1913-1922. https://doi.org/10.1364/OE.15.001913

Properties
  • adaptive (numeric) – adaptive-adaptive factor.
  • weighted (logical) – if the method is weighted Gerchberg-Saxton.
Methods (inherited)
  • run() – Run the iterative method
  • showFitness() – Show the fitness graph
Properties (inherited)
  • components (real: 0, 2*pi) – NxMxD matrix of D patterns to be combined.
  • guess – Best guess at hologram pattern (complex)
  • target – Target pattern for estimating fitness (complex, optional)
  • vismethod – Method used to do the visualisation
  • invmethod – Method used to calculate initial guess/inverse-visualisation
  • phase – Phase of the best guess (real: 0, 2*pi)
  • amplitude – Amplitude of the best guess (real)
  • objective – Objective function used to evaluate fitness or []
  • fitness – Fitness evaluated after every iteration or []

See also CombineGerchbergSaxton and GerchbergSaxton.

CombineGerchbergSaxton(components, varargin)

Construct a new Gerchberg-Saxton combination iterative method.

Usage
mtd = IterCombine(components, …)
Parameters
  • components (real: 0, 2*pi) – NxMxD array of D phase patterns to be combined. Phase patterns should have range [0, 2*pi] or equivalent.
Optional named arguments
  • adaptive num Adaptive-Adaptive factor. Default: 1.0, i.e. the method is Gerchberg-Saxton.
  • weighted (logical) – If the method should use weighted Gerchberg-Saxton. Default: false.
  • target (complex) – approximate pattern for the target. This is only used for estimating the current fitness. Default: otslm.tools.combine(components, 'method', 'farfield').
  • guess (complex) – Initial guess at combination of patterns. Default: exp(2*pi*1i*random_super) where random_super = tools.combine(components, 'method', 'rsuper')
  • vismethod fcn Function to calculate far-field. Takes one argument: the complex amplitude near-field. Optional, only used for fitness evaluation. Default: @otslm.tools.prop.FftForward.simpleProp.propagate
  • invmethod fcn Function to calculate near-field. Takes one argument: the complex amplitude far-field. Optional, not used. Default: @otslm.tools.prop.FftInverse.simpleProp.propagate
  • objective fcn Objective function to measure fitness. Default: @otslm.iter.objectives.FlatIntensity

IterBase

This is the base class for iterative methods. It is an abstract class and cannot be directly instantiated. To implement your own iterative method class, inherit from this class and implement the abstract methods/properties.

class otslm.iter.IterBase(varargin)

Base class for iterative algorithm classes. Inherits from handle.

Methods
  • run() – Run the iterative method
  • showFitness() – Show the fitness graph
Properties
  • guess – Best guess at hologram pattern (complex)
  • target – Target pattern the method tries to approximate (complex)
  • vismethod – Method used to do the visualisation
  • invmethod – Method used to calculate initial guess/inverse-visualisation
  • phase – Phase of the best guess (real: 0, 2*pi)
  • amplitude – Amplitude of the best guess (real)
  • objective – Objective function used to evaluate fitness or []
  • fitness – Fitness evaluated after every iteration or []
Abstract methods
  • iteration() – run a single iteration of the method
IterBase(varargin)

Constructor for iterative algorithm (abstract) base class

Usage
mtd = IterBase(target, …)
Parameters
  • target – target pattern to generate
Optional named arguments
  • guess im Initial guess at complex amplitude pattern. If not image is supplied, a guess is created using invmethod.
  • vismethod fcn Function to calculate far-field. Takes one argument: the complex amplitude near-field. Default: @otslm.tools.prop.FftForward.simpleProp.propagate
  • invmethod fcn Function to calculate near-field. Takes one argument: the complex amplitude far-field. Default: @otslm.tools.prop.FftInverse.simpleProp.propagate
  • objective fcn Objective function to measure fitness. Default: @otslm.iter.objectives.FlatIntensity
evaluateFitness(mtd, varargin)

Evaluate the fitness of the current guess

Usage

score = mtd.evaluateFitness() visualises the current guess and evaluate the fitness.

score = mtd.evaluateFitness(guess) evaluate the fitness of the given guess. If guess is a stack of matrices, the returned score is a vector with size(trial, 3) elements. Guess should be a complex amplitude.

run(mtd, num_iterations, varargin)

Run the method for a specified number of iterations

Usage
result = mtd.run(num_iterations, …) run for the specified number of iterations.
Parameters
  • num_iterations (numeric) – Number of iterations
Optional named arguments
  • show_progress bool display a figure with optimisation progress
stopIterations(mtd, src, event)

Callback for the stop button in showFitness

Usage
mtd.stopIterations(…) arguments are ignored.

IterCombine

This is the base class for iterative methods which combine multiple input pattern. It is an abstract class inheriting from IterBase however not all properties are needed/used by classes inheriting from this method. For instance, the CombineGerchbergSaxton class only uses the vismethod to calculate the fitness when an objective function is supplied. If the objective is omitted the method doesn’t calculate the fitness and doesn’t need vismethod or invmethod.

class otslm.iter.IterCombine(components, varargin)

Base class for iterative combination algorithms. Inherits from IterBase.

Iterative methods that inherit from this class attempt to combine a set of SLM phase patterns \(\phi_m\) into a single phase pattern which generates a far-field phase pattern approximating the combination of each input phase pattern.

The target field is a optional and is only used for estimating fitness of the generated pattern.

Methods (inherited)
  • run() – Run the iterative method
  • showFitness() – Show the fitness graph
Properties
  • components (real: 0, 2*pi) – NxMxD matrix of D patterns to be combined.
Properties (inherited)
  • guess – Best guess at hologram pattern (complex)
  • target – Target pattern for estimating fitness (complex, optional)
  • vismethod – Method used to do the visualisation
  • invmethod – Method used to calculate initial guess/inverse-visualisation
  • phase – Phase of the best guess (real: 0, 2*pi)
  • amplitude – Amplitude of the best guess (real)
  • objective – Objective function used to evaluate fitness or []
  • fitness – Fitness evaluated after every iteration or []
Abstract methods
  • iteration() – run a single iteration of the method
IterCombine(components, varargin)

Constructor for iterative combination algorithms (abstract) base class

Usage
mtd = IterCombine(components, …)
Parameters
  • components (real: 0, 2*pi) – NxMxD array of D phase patterns to be combined. Phase patterns should have range [0, 2*pi] or equivalent.
Optional named arguments
  • target (complex) – approximate pattern for the target. This is only used for estimating the current fitness. Default: otslm.tools.combine(components, 'method', 'farfield').
  • guess (complex) – Initial guess at combination of patterns. Default: exp(2*pi*1i*random_super) where random_super = tools.combine(components, 'method', 'rsuper')
  • vismethod fcn Function to calculate far-field. Takes one argument: the complex amplitude near-field. Optional: only used for fitness evaluation. Default: @otslm.tools.prop.FftForward.simpleProp.propagate
  • invmethod fcn Function to calculate near-field. Takes one argument: the complex amplitude far-field. Optional: not used. Default: @otslm.tools.prop.FftInverse.simpleProp.propagate
  • objective fcn Objective function to measure fitness. Default: @otslm.iter.objectives.FlatIntensity

IterBaseEwald

This is the base class for iterative methods that 3-D Fourier transforms and an Ewald sphere far-field mapping. This is class can be combined with IterBase to provide the 3-D specialisation. Currently only used by GerchbergSaxton3d.

class otslm.iter.IterBaseEwald(target, varargin)

Abstract base class for 3-D Ewald iterative algorithm classes Inherits from IterBase.

Methods
  • run() – Run the iterative method
Properties
  • guess – Best guess at hologram pattern (complex, matrix)
  • target – Target pattern the method tries to approximate (volume)
  • vismethod – Method used to do the visualisation
  • invmethod – Method used to calculate initial guess/inverse-visualisation
  • phase – Phase of the best guess (real: 0, 2*pi)
  • amplitude – Amplitude of the best guess (real)
  • objective – Objective function used to evaluate fitness or []
  • fitness – Fitness evaluated after every iteration or []
Abstract methods
  • iteration() – run a single iteration of the method
IterBaseEwald(target, varargin)

Abstract constructor for 3-D iterative algorithm base class

Usage
mtd = IterBaseEwald(target, …) constructs a new instance. target should be a 3-D volume. Guess, if supplied, should be a 2-D matrix for the pattern on the SLM.
Optional named arguments:
  • guess im Initial guess at complex amplitude pattern. If no image is supplied, a guess is created using invmethod.
  • vismethod fcn Function to calculate far-field. Takes one argument: the complex amplitude near-field. Default: @otslm.tools.prop.FftEwaldForward.simpleProp.propagate
  • invmethod fcn Function to calculate near-field. Takes one argument: the complex amplitude far-field. Default: @otslm.tools.prop.FftEwaldInverse.simpleProp.propagate
  • objective fcn Objective function to measure fitness. Default: @otslm.iter.objectives.FlatIntensity

bsc

This function attempts to optimise the beam using vector spherical wave functions. The function may be unstable/change in future releases but demonstrates how OTT can be used with OTSLM.

otslm.iter.bsc(sz, target, varargin)

Optimisation in vector spherical wave function basis

Requires the optical tweezers toolbox (OTT).

Usage
[pattern, beam, coeffs] = bsc(target, …) attempt to produce target using a phase pattern. Returns the phase pattern matched to the beam (bsc) and optimised basis weighting coefficients.
Parameters
  • target – target pattern
Optional named parameters
  • ‘incident’ pattern – Incident illumination on SLM
  • ‘roi’ func – Region to optimise (default: roiAll)
  • ‘basis’ str – BSC basis to optimise in (default: vswf_lg)
  • ‘basis_size’ num – Number of basis functions to use
  • ‘polarisation’ [x y] – Polarisation of the basis functions
  • ‘wavelength’ num – Wavelength in medium [m]
  • ‘speed’ num – Speed in medium [m/s]
  • ‘NA’ num – Numberical aperture of objective
  • ‘pixel_size’ num – Size of pixels in target [m]
  • ‘method’ str – Optimisation method to use
  • ‘radius’ num – Radius for hologram unwrapping (default: 1.0)

bowman2017

This function provides an interface for Bowman, et al. Optics Express 25, 11692 (2017). This requires a suitable Python version and various libraries. The wrapper may be unstable and will hopefully be improved in future releases.

otslm.iter.bowman2017(target, varargin)

Wrapper for Bowman 2017 conjugate gradient implementation

See Bowman, et al. Optics Express 25, 11692 (2017) If you use this method, please consider citing Bowman 2017.

Warning

This wrapper may be unstable and may change in future releases.

Usage
pattern = bowman2017(target, …) attempt to generate the target using a phase pattern optimised using conjugate gradient method.
Parameters
  • target – target pattern to generate
Optional named parameters
  • ‘guess’ – Initial guess at the phase
  • ‘iterations’ – Number of iterations (default: 200)
  • ‘steepness’ – Steepness for Bowman cost function (default: 9.0)
  • ‘incident’ – Incident illumination (default: ones)
  • ‘roisize’ – Optimisation region size (default: min(size)/2)

Objective functions

Objective functions are contained in the otslm.iter.objectives sub-package. These functions are used with the above optimisation methods for both optimisation and diagnostics.

To evaluate the fitness (similarity) between a trial pattern and a target, we can construct a new objective instance and call the evaluate method. For example, using the Flatness objective:

% Setup the trial and target
sz = [256, 256];
target = ones(sz);
trial = randn(sz) + 1.0;

% Setup the objective
obj = otslm.iter.objectives.Flatness('target', target);

% Evaluate the fitness
fitness = obj.evaluate(trial);

It is possible to reuse the objective multiple times or test the trial pattern against a different target pattern when evaluate is called:

new_target = zeros(sz);
fitness = obj.evaluate(trial, new_target);

Objective classes support a region of interest mask. The region of interest can either be a logical mask or a function which selects a region of the image, for example:

% Select only half of the image with a function
obj.roi = @(pattern) pattern(1:end/2, :)
fitness = obj.evaluate(trial);

% Use a logical array
obj.roi = otslm.simple.aperture(sz, 128);
fitness = obj.evaluate(trial);

Objective base class

class otslm.iter.objectives.Objective(varargin)

Abstract base class for optimisation objective functions.

To use this class, you need to inherit from it and implement the evaluate_internal function.

Methods
  • evaluate() – Evaluate the fitness of the specified pattern
Properties
  • target (numeric) – Target pattern to compare with (or [] for no default). This is only used if no target is provided in evaluate().
  • type (enum) – Type of optimisation function (‘min’ or ‘max’) This property isn’t widely used (may change in future version). For now, most functions simply have this property set to ‘min’.
  • roi (logical|empty|function_handle) – Region of interest to apply to target and trial. Must be either a logical array the same size as target, an empty matrix for no roi, or a function handle. If roi is a function handle, the function should have the signature masked = roi(pattern). Calling the function should select elements from the pattern for comparison. The function is applied to both the target and the trial pattern.
Abstract methods
  • evaluate_internal() – Implementation called by evaluate(). Signature: fitness = obj.evaluate_internal(target, trial). The roi has already been applied to the trial and target.

See also Objective, Intensity and Flatness.

Objective(varargin)

Construct a new objective function instance

Usage
obj = Objective(…) construct a new objective function instance.
Optional named arguments
  • roi [] | logical | function_handle – specify the roi to use when evaluating the fitness function. Can be a logical array or a function handle. Default: []
  • target [] | matrix – specify the target pattern for this objective. If not supplied, the target must be supplied when the evaluate function is called. Default: []
evaluate(trial, target)

Evaluate the fitness of the specified trial pattern.

Usage
fitness = obj.evaluate(trial, [target]) evaluate the specified trial pattern. If target is not specified, uses the internal target pattern set during construction.
Parameters
  • trial (numeric) – pattern to compare to target
  • target (numeric) – pattern to compare to trial. Optional. Default target is obj.target.

Bowman2017

class otslm.iter.objectives.Bowman2017(varargin)

Cost function used in Bowman et al. 2017 paper. Inherits from Objective.

\[C = 10^d * (1.0 - \sum_{nm} \sqrt{I_nm T_nm} \cos(phi_nm - psi_nm))^2\]

target and trial should be the complex field amplitudes.

Properties
  • scale – d scaling factor in cost function
  • field – ‘complex’, ‘phase’, or ‘amplitude’ for optimisation type
  • normalize – Normalize target/trial every evaluation

See also Bowman2017 and Intensity.

Bowman2017(varargin)

Construct a new objective function instance

Usage
obj = Bowman2017(…) construct a new objective function instance.
Optional named arguments
  • scale num – d scaling factor in cost function. Default: 0.5
  • field [char] – One of ‘complex’, ‘phase’, or ‘amplitude’ for optimisation type. Default: ‘complex’.
  • normalize bool – Normalize target/trial every evaluation. Default: true.
  • roi [] | logical | function_handle – specify the roi to use when evaluating the fitness function. Can be a logical array or a function handle. Default: []
  • target [] | matrix – specify the target pattern for this objective. If not supplied, the target must be supplied when the evaluate function is called. Default: []

FlatIntensity

class otslm.iter.objectives.FlatIntensity(varargin)

Objective function for pattern flatness and intensity. Inherits from Intensity and Flatness.

Evaluates the fitness using the Intensity and Flatness method and adds the result:

Properties
  • flatness – scaling factor for pattern flatness

See also FlatIntensity and Bowman2017.

FlatIntensity(varargin)

Construct a new objective function instance

Usage
obj = FlatIntensity(…)
Optional named arguments
  • roi [] | logical | function_handle – specify the roi to use when evaluating the fitness function. Can be a logical array or a function handle. Default: []
  • target [] | matrix – specify the target pattern for this objective. If not supplied, the target must be supplied when the evaluate function is called. Default: []
  • flatness num – Scaling factor for pattern flatness. Default: 0.5.

Flatness

class otslm.iter.objectives.Flatness(varargin)

Objective function for pattern flatness Inherits from Objective.

See also Flatness, Intensity and FlatIntensity.

Flatness(varargin)

Construct a new objective function instance

Usage
obj = Flatness(…)
Optional named arguments
  • roi [] | logical | function_handle specify the roi to use when evaluating the fitness function. Can be a logical array or a function handle. Default: []
  • target [] | matrix specify the target pattern for this objective. If not supplied, the target must be supplied when the evaluate function is called. Default: []

Goorden2014

class otslm.iter.objectives.Goorden2014(varargin)

Fidelity function from Goorden, et al. 2014 paper.

\[F = |\textrm{conj}(target) * trial|^2\]

Error is 1 - F.

Properties
  • normalize (logical) – True if the patterns should be normalized by the area (i.e., \(F = F/A^2\)).

See also Goorden2014, Flatness and :class`Bowman2017`.

Goorden2014(varargin)

Construct a new objective function instance

Usage
obj = Goorden2014(…)
Optional named arguments
  • normalize logical – If true, normalized the pattern by the number of pixels in the pattern. Default: true.
  • roi [] | logical | function_handle – specify the roi to use when evaluating the fitness function. Can be a logical array or a function handle. Default: []
  • target [] | matrix – specify the target pattern for this objective. If not supplied, the target must be supplied when the evaluate function is called. Default: []

Intensity

class otslm.iter.objectives.Intensity(varargin)

Objective function for pattern intensity.

See also Intensity, Flatness and FlatIntensity.

Intensity(varargin)

Construct a new objective function instance

Usage
obj = Intensity(…)
Optional named arguments
  • roi [] | logical | function_handle specify the roi to use when evaluating the fitness function. Can be a logical array or a function handle. Default: []
  • target [] | matrix specify the target pattern for this objective. If not supplied, the target must be supplied when the evaluate function is called. Default: []

RmsIntensity

class otslm.iter.objectives.RmsIntensity(varargin)

Objective function for pattern RMS intensity Inherits from Objective.

Evaluates the fitness according to

\[F = \sqrt{\textrm{mean}((|t|^2 - |T|^2)^2)}\]

Where \(t\) and \(T\) are the trial and target pattern complex amplitudes.

See also RmsIntensity, Intensity and Flatness.

RmsIntensity(varargin)

Construct a new objective function instance

Usage
obj = RmsIntensity(…)
Optional named arguments
  • roi [] | logical | function_handle specify the roi to use when evaluating the fitness function. Can be a logical array or a function handle. Default: []
  • target [] | matrix specify the target pattern for this objective. If not supplied, the target must be supplied when the evaluate function is called. Default: []