`ADI.ADIAlgorithm`

— Type`ADI.ADIAlgorithm`

This abstract type is used for defining ADI algorithms. Algorithms are stateful objects that define the options for a given algorithm (e.g. the number of components used in PCA decomposition). The most direct usage of an algorithm is to use it to fit HCI data; that is, given a sample of pixels, apply the algorithm using the given options and return an object containing all the necessary information to reconstruct the input data.

See the extended help (`??ADIAlgorithm`

) for interface details.

**Extended help**

**Interface**

To extend `ADIAlgorithm`

you may implement the following

`ADI.fit(::Alg, data::AbstractMatrix; kwargs...)`

Fit the data (flattened into a matrix). To support RDI, ensure the `ref`

keyword argument is usable (`ref`

is also a flattened matrix). This is the only method you *need* to implement for a new `ADIAlgorithm`

, along with a suitable `ADIDesign`

.

ADI.jl automatically coaxes the `cube`

input into a matrix for use with `fit`

, appropriately handling the various geometries. When available, this input is a view, so if the algorithm requires dense arrays, make sure to call `collect`

when appropriate. If a given algorithm doesn't support the default operations, all that needs to be done is override the default behavior (for an example, see the `GreeDS`

implementation).

`reconstruct(::Alg, cube; kwargs...)`

Fit the data using the algorithm and return a cube with the estimate of the PSF. By default uses the reconstruction from the `ADIDesign`

fit to the data.

`subtract(::Alg, cube; kwargs...)`

Fit the data using the algorithm and return a cube that has had the PSF estimate subtracted. By default, calls `reconstruct`

and subtracts it from `cube`

.

```
process(::ADIAlgorithm, cube; kwargs...)
(::ADIAlgorithm)(cube; kwargs...)
```

Fully process the data (estimate, subtract, collapse). By default, derotates and collapses output from `subtract`

. You only need to define `process`

, since the functor version is supplied automatically.

`ADI.ADIAlgorithm`

— Method`ADI.ADIDesign`

— Type`ADI.ADIDesign`

An abstract type used for holding the output of `ADI.fit`

. The purpose of these types is to hold the *minimum* information required to reconstruct the PSF estimate (e.g. weights and principal components from PCA). ADI.jl includes two designs- `ADI.ClassicDesign`

and `ADI.LinearDesign`

.

**Interface**

When implementing a new algorithm, if your outputs do not fit into either of those designs, you will have to create your own `ADIDesign`

with the following methods-

`ADI.design(::Design)`

accessor for the pertinent data to approximate the PSF

`reconstruct(::Design)`

return the approximate PSF estimate as a matrix

`ADI.Classic`

— Type```
Classic(;method=median)
Classic(method)
```

Classic PSF subtraction using the median of entire data cube. If another statistic is desired, like the `mean`

, it can be passed as an argument as long as it supports slicing multi-dimensional arrays.

**References**

- Marois et al. 2006 Angular Differential Imaging: A Powerful High-Contrast Imaging Technique

`ADI.ClassicDesign`

— Type`ADI.ClassicDesign(n, frame)`

Output for the `Classic`

algorithm which contains the static frame unrolled into a vector (with size `Npx`

). `reconstruct`

will tile this vector in a non-allocating way `n`

times to appear like a flattened cube. `ADI.design`

will return the static frame .

`ADI.DoubleSDI`

— Type```
DoubleSDI(alg)
DoubleSDI(alg_spec, alg_temp)
```

A wrapper algorithm for spectral differential imaging (SDI) data reduced in two passes. The first pass uses `alg_spec`

to reduce each spectral cube slice in the SDI tensor. Then, the spectral residual frames will be reduced using `alg_temp`

, which will include the derotation and final combination.

The difference between `DoubleSDI`

and `SliceSDI`

is that `DoubleSDI`

does its first pass in the spectral slice, effectively collapsing the slice before performing ADI on the residual cube. `SliceSDI`

does its first pass in the temporal slice, collapsing it first before performing ADI on the residual cube.

**Examples**

```
julia> data, angles, scales = # load data...
# Median subtraction for each spectral slice,
# GreeDS{PCA} subtraction on spectral residual cube
julia> res = DoubleSDI(Classic(), GreeDS(15))(data, angles, scales)
```

`ADI.Framewise`

— Method```
Framewise(alg; limit=Inf, delta_rot=nothing)
Framewise(algs::AbstractVector; limit=Inf, delta_rot=nothing)
```

Wrap an algorithm such that the underlying data will be processed frame by frame. For each frame a reference library is created from the data. This reference can be filtered by rejecting frames which have not rotated a sufficient parallactic angle. `delta_rot`

sets the required arc length for rotation in units of the FWHM. If `delta_rot`

is `nothing`

there will be now temporal filtering. The number of frames retained can be specified with `limit`

, e.g. the 4 closest frames in time with the target frame.

The following keyword arguments must be provided to `reconstruct`

or `subtract`

`angles`

- the measured parallactic angles for each frame

If `delta_rot`

is provided, the following additional keyword arguments must be provided to `reconstruct`

or `subtract`

.

`fwhm`

- the FWHM of the instrument in pixels. Will be set to the width of a`MultiAnnulusView`

`r`

- The radius of the arc to calculate the parallactic angle threshold. Will be set automatically if using`AnnulusView`

or`MultiAnnulusView`

.

In addition, `Framewise`

versions of algorithms do not implement `ADI.fit`

and do not currently support RDI.

**Annular reduction**

In the case of reducing a `MultiAnnulusView`

, a vector of algorithms can be used, each one corresponding to each annulus of the view. In this case, too, `delta_rot`

can be given as a vector or as a tuple. If it is given as a tuple, `delta_rot`

will increase linearly from the first value to the last value across each annulus.

**Examples**

```
julia> cube, angles = # load data
julia> alg = PCA(10) # the algorithm to use on each reference
julia> res = Framewise(alg)(cube, angles);
julia> mav = MultiAnnulusView(cube, 5; inner=5);
julia> res_ann = Framewise(alg, delta_rot=(0.1, 1))(mav, angles);
```

`ADI.GreeDS`

— Type```
GreeDS(alg=PCA(); threshold=0)
GreeDS(ncomps; threshold=0, options...)
```

Performs the greedy disk subtraction (GreeDS) algorithm.

This method is an iterative approach to standard ADI reduction which seeks to minimize over-subtraction by constraining the low-rank matrix approximation from `alg`

. By default, uses `PCA`

. If `ncomps`

or other PCA options are provided, they will be passed to the constructor.

The GreeDS algorithm requires fully reconstructing a cube at each iteration, which requires knowing the geometry of the input (full-frame, annulus, etc.) and the corresponding parallactic angles. These angles must be passed as a keyword argument `angles`

. In the case of reducing data, e.g. `GreeDS()(cube, angles)`

the angles will be passed automatically. It is important to clarify, *these angles should correspond to the reference data in the case of RDI*, e.g. `GreeDS()(cube, angles; ref=ref_cube, angles=ref_angles)`

**Algorithms**

The following algorithms work natively with GreeDS: `PCA`

and `NMF`

**References**

- Pairet et al. 2018 "Reference-less algorithm for circumstellar disks imaging"
- Pairet et al. 2020 "MAYONNAISE: a morphological components analysis pipeline for circumstellar disks and exoplanets imaging in the near infrared"

`ADI.LOCI`

— Type`LOCI(; dist_threshold=nothing, metric=Cityblock())`

Local optimal combination of images (LOCI).

If provided, the frames used for the reference are filtered to only include the frames whose pairwise distances are within the `dist_threshold`

quantile. The distances are measured using Distances.jl and the metric used for measuring the distacnes can be specified with `metric`

(by default uses the Manhattan/cityblock metric).

**References**

- Lafreniere et al. (2007) "A New Algorithm for Point-Spread Function Subtraction in High-Contrast Imaging: A Demonstration with Angular Differential Imaging"

`ADI.LinearDesign`

— Type`ADI.LinearDesign(coeffs, basis)`

A "linear" design implies the use of some linear basis for reconstructing data along with a set of coefficients or weights. The reconstruction will be the matrix product of the basis and the weights, `Z * w`

.

`ADI.design`

will return `(basis, ceoffs)`

, and you can also extract them via iteration, like `Z, w = design`

.

`ADI.NMF`

— Type```
NMF(;ncomps=nothing)
NMF(ncomps)
```

Use non-negative matrix factorization (NMF, NNMF) to form a non-negative, low-rank, and orthonormal basis of the input. The implementation of the underlying NMF is provided by NMF.jl. The implementation uses a non-negative SVD for initialization and a coordinate-descent solver to fit.

If `ncomps`

is `nothing`

, it will be set to the number of frames in the reference cube when processed.

NMF is not designed to fit negative values. This algorithm will warn you (but will not error) if a target or reference cube contains negative values. The output may seem reasonable, but it is not well-defined with respect to the NMF algorithm. To overcome this, rescaling the data by its minimum before processing is necessary

```
target = cube .- minimum(cube)
S = reconstruct(NMF(), target, angles)
```

When doing full-frame reduction (e.g. `NMF()(cube, angles)`

) *this is handled automatically*, so this constraint only applies to the lower-level API and methods which rely on those, such as `GreeDS`

. **In general, if you see warnings, heed them.**

**References**

- Ren et al. 2018 Non-negative Matrix Factorization: Robust Extraction of Extended Structures

`ADI.PCA`

— Type```
PCA(ncomps; options...)
PCA(;ncomps=nothing, options...)
```

Use principal components analysis (PCA) to form a low-rank orthonormal basis of the input. Uses deterministic singular-value decomposition (SVD) to decompose data.

If `ncomps`

is `nothing`

, the basis will not be truncated (i.e. `ncomps`

is equal to the number of frames). `ncomps`

can be set to `:noise`

or `:pratio`

to automatically choose the number of components using the residual frame noise or principal ratio, respectively. For more information, see the extended help.

**References**

- Soummer, Pueyo, and Larkin (2012) "Detection and Characterization of Exoplanets and Disks Using Projections on Karhunen-Loève Eigenimages"

**Extended help**

**Optimizing ncomps**

There are a few ways to optimize `ncomps`

using the input data. Additional options for the optimization are listed below

`ncomps=:noise`

- residual noise optimization`ncomps=:pratio`

- principal ratio optimization

**Residual noise optimization**

This technique progressively increases `ncomps`

at each step measuring the pixel-to-pixel noise (standard deviation) in the residual data. Iteration will stop when the noise is not improving beyond a threshold. This is suited for data with similar statistical characteristics, such as an annulus more so than a full-frame cube.

`collapse=false`

- if true, the temporal median of the residual data will be used for measuring the noise.`noise_error=1e-3`

- the threshold for the minimal noise improvement looking back 2 iterations

**Principal ratio optimization**

This technique chooses the number of components required to explain some ratio of the total variance in the data. This is known as the *principal ratio* or the *explained variance ratio*. The explained variance is measured by transforming the singular values of the SVD decomposition (`Λ = @. S^2 / (n - 1)`

).

`pratio=0.9`

- the target principal ratio (between 0 and 1)

`ADI.SDIAlgorithm`

— Type`ADI.SDIAlgorithm <: ADI.ADIAlgorithm`

Spectral differential imaging (SDI) algorithms. These work on 4-D SDI tensors. To use these algorithms, simply treat them like functions (or call `process`

)

`(::SDIAlgorithm)(data::AbstractArray{T,4}, angles, scales; [ref] kwargs...)`

The data is expected to be laid out in `(nx, ny, nλ, nf)`

format, so you may need to `permutedims`

before processing the data. The `scales`

correspond to the relative wavelength scales for each spectrum, and can be retrieved with `HCIToolbox.scale_list`

.

**Algorithms**

The current SDI implementations are

`ADI.SingleSDI`

— Type`SingleSDI(alg)`

A wrapper algorithm for spectral differential imaging (SDI) data reduced in a single pass. This means that each channel will be scaled and then concatenated together, so an SDI tensor `(nx, ny, nλ, nf)`

becomes a stack `(nx, ny, nλ * nf)`

which is processed by the underlying `alg`

like ADI data.

`SingleSDI`

is the default SDI mode. This means instead of writing

`SingleSDI(PCA(15))(data, angles, scales)`

you can just write

`PCA(15)(data, angles, scales)`

`ADI.SliceSDI`

— Type```
SliceSDI(alg)
SliceSDI(alg_spec, alg_temp)
```

A wrapper algorithm for spectral differential imaging (SDI) data reduced in two passes. The first pass uses `alg_temp`

to reduce each temporal cube slice in the SDI tensor. These residuals will be rescaled and stacked into a new cube. Then, the temporal residual frames will be reduced using `alg_spec`

, which will include the derotation and final combination.

The difference between `SliceSDI`

and `DoubleSDI`

is that `DoubleSDI`

does its first pass in the spectral slice, effectively collapsing the slice before performing ADI on the residual cube. `SliceSDI`

does its first pass in the temporal slice, collapsing it first before performing ADI on the residual cube.

**Examples**

```
julia> data, angles, scales = # load data...
# Median subtraction for each spectral slice,
# GreeDS{PCA} subtraction on spectral residual cube
julia> res = SliceSDI(Classic(), GreeDS(15))(data, angles, scales)
```

`ADI.design`

— Function`ADI.design(::ADIDesign)`

Return the pertinent data required to form the PSF approximation. For example, weights and components for PCA/NMF or the median frame for classic ADI.

`ADI.expand_rotate`

— Method`expand_rotate(frame, angles, threshold; kwargs...)`

Takes a frame, expands it into a cube, rotates it clockwise by `angles`

, and min-clips at `threshold`

. Keyword arguments will be passed to `HCIToolbox.derotate`

.

`ADI.fit`

— Method`ADI.fit(::ADIAlgorithm, cube; [ref], kwargs...)`

Given the description of an algorithm and the appropriate options, take the pixels from `cube`

and fit them, returning an (`ADIDesign`

) containing the necessary information from the fit (e.g. the principal components from PCA decomposition).

If the algorithm supports reference differential imaging (RDI), the reference cube can be passed by the keyword argument `ref`

.

`ADI.process`

— Method`ADI.reconstruct`

— Method`reconstruct(alg, cube; [ref], kwargs...)`

Reconstruct the PSF approximation for the given algorithm, using `ref`

as the reference cube if given and supported by the algorithm.

**Examples**

```
julia> cube, angles = # load data
julia> S = reconstruct(PCA(10), cube);
julia> size(S) == size(cube)
true
julia> flat_res = collapse(cube .- S, angles); # form resid, derotate, and combine
```

`ADI.subtract`

— Method`subtract(alg, cube; [ref], kwargs...)`

Reconstruct the PSF approximation for the given algorithm and subtract it from `cube`

, using `ref`

as the reference cube if given and supported by the algorithm.

**Examples**

```
julia> cube, angles = # load data
julia> R = subtract(PCA(10), cube);
julia> size(R) == size(cube)
true
julia> flat_res = collapse(R, angles); # derotate, and combine
```