Abstract Operators

Linear operators

Basic Operators

AbstractOperators.Eye — Type.

Eye([domainType=Float64::Type,] dim_in::Tuple)

Eye([domainType=Float64::Type,] dims...)

Create the identity operator.

julia> op = Eye(Float64,(4,))
I  ℝ^4 -> ℝ^4

julia> op = Eye(2,3,4)
I  ℝ^(2, 3, 4) -> ℝ^(2, 3, 4)

julia> op*ones(2,3,4) == ones(2,3,4)
true

source

AbstractOperators.Zeros — Type.

Zeros(domainType::Type, dim_in::Tuple, [codomainType::Type,] dim_out::Tuple)

Create a LinearOperator which, when multiplied with an array x of size dim_in, returns an array y of size dim_out filled with zeros.

For convenience Zeros can be constructed from any AbstractOperator.

julia> Zeros(Eye(10,20))
0  ℝ^(10, 20) -> ℝ^(10, 20)

julia> Zeros([Eye(10,20) Eye(10,20)])
[0,0]  ℝ^(10, 20)  ℝ^(10, 20) -> ℝ^(10, 20)

source

AbstractOperators.DiagOp — Type.

DiagOp(domainType::Type, dim_in::Tuple, d::AbstractArray)

DiagOp(d::AbstractArray)

Creates a LinearOperator which, when multiplied with an array x, returns the elementwise product d.*x.

julia> D = DiagOp(Float64, (2, 2,), [1. 2.; 3. 4.])
╲  ℝ^(2, 2) -> ℝ^(2, 2)

julia> D*ones(2,2)
2×2 Array{Float64,2}:
 1.0  2.0
 3.0  4.0

source

AbstractOperators.GetIndex — Type.

GetIndex([domainType=Float64::Type,] dim_in::Tuple, idx...)

GetIndex(x::AbstractArray, idx::Tuple)

Creates a LinearOperator which, when multiplied with x, returns x[idx].

julia> x = collect(linspace(1,10,10));

julia> G = GetIndex(Float64,(10,), 1:3)
↓  ℝ^10 -> ℝ^3 

julia> G*x
3-element Array{Float64,1}:
 1.0
 2.0
 3.0

julia> GetIndex(randn(10,20,30),(1:2,1:4))
↓  ℝ^(10, 20, 30) -> ℝ^(2, 4)

source

AbstractOperators.MatrixOp — Type.

MatrixOp(domainType=Float64::Type, dim_in::Tuple, A::AbstractMatrix)

MatrixOp(A::AbstractMatrix)

MatrixOp(A::AbstractMatrix, n_colons)

Creates a LinearOperator which, when multiplied with a vector x::AbstractVector, returns the product A*x.

The input x can be also a matrix: the number of columns must be given either in the second entry of dim_in::Tuple or using the constructor MatrixOp(A::AbstractMatrix, n_colons).

julia> MatrixOp(Float64,(10,),randn(20,10))
▒  ℝ^10 -> ℝ^20 

julia> MatrixOp(randn(20,10))
▒  ℝ^10 -> ℝ^20

julia> MatrixOp(Float64,(10,20),randn(20,10))
▒  ℝ^(10, 20) -> ℝ^(20, 20)

julia> MatrixOp(randn(20,10),4)
▒  ℝ^(10, 4) -> ℝ^(20, 4)

source

AbstractOperators.LMatrixOp — Type.

LMatrixOp(domainType=Float64::Type, dim_in::Tuple, b::Union{AbstractVector,AbstractMatrix})

LMatrixOp(b::AbstractVector, number_of_rows::Int)

Creates a LinearOperator which, when multiplied with a matrix X::AbstractMatrix, returns the product X*b.

julia> op = LMatrixOp(Float64,(3,4),ones(4))
(⋅)b  ℝ^(3, 4) -> ℝ^3 

julia> op = LMatrixOp(ones(4),3)
(⋅)b  ℝ^(3, 4) -> ℝ^3

julia> op*ones(3,4)
3-element Array{Float64,1}:
 4.0
 4.0
 4.0

source

AbstractOperators.MyLinOp — Type.

MyLinOp(domainType::Type, dim_in::Tuple, [domainType::Type,] dim_out::Tuple, Fwd!::Function, Adj!::Function)

Construct a user defined LinearOperator by specifing its linear mapping Fwd! and its adjoint Adj!. The functions Fwd! and Adj must be in-place functions consistent with the given dimensions dim_in and dim_out and the domain and codomain types.

julia> n,m = 5,4;

julia> A = randn(n,m);

julia> op = MyLinOp(Float64, (m,),(n,), (y,x) -> mul!(y,A,x), (y,x) -> mul!(y,A',x))
A  ℝ^4 -> ℝ^5

julia> op = MyLinOp(Float64, (m,), Float64, (n,), (y,x) -> mul!(y,A,x), (y,x) -> mul!(y,A',x))
A  ℝ^4 -> ℝ^5

source

DSP

Transformations

AbstractOperators.DFT — Type.

DFT([domainType=Float64::Type,] dim_in::Tuple [,dims])

DFT(dim_in...)

DFT(x::AbstractArray [,dims])

Creates a LinearOperator which, when multiplied with an array x::AbstractArray{N}, returns the N-dimensional Discrete Fourier Transform over dimensions dims of x.

julia> DFT(Complex{Float64},(10,10))
ℱ  ℂ^(10, 10) -> ℂ^(10, 10)

julia> DFT(10,10)
ℱ  ℝ^(10, 10) -> ℂ^(10, 10)

julia> A = DFT(ones(3))
ℱ  ℝ^3 -> ℂ^3

julia> A*ones(3)
3-element Array{Complex{Float64},1}:
 3.0+0.0im
 0.0+0.0im
 0.0+0.0im

source

AbstractOperators.IDFT — Type.

IDFT([domainType=Float64::Type,] dim_in::Tuple [,dims])

IDFT(dim_in...)

IDFT(x::AbstractArray [,dims])

Creates a LinearOperator which, when multiplied with an array x::AbstractArray{N}, returns the N-dimensional Inverse Discrete Fourier Transform over dimensions dims of x.

julia> IDFT(Complex{Float64},(10,10))
ℱ⁻¹  ℂ^(10, 10) -> ℂ^(10, 10)

julia> IDFT(10,10)
ℱ⁻¹ ℝ^(10, 10) -> ℂ^(10, 10)

julia> A = IDFT(ones(3))
ℱ⁻¹  ℝ^3 -> ℂ^3

julia> A*ones(3)
3-element Array{Complex{Float64},1}:
 1.0+0.0im
 0.0+0.0im
 0.0+0.0im

source

AbstractOperators.RDFT — Type.

RDFT([domainType=Float64::Type,] dim_in::Tuple [,dims=1])

RDFT(dim_in...)

RDFT(x::AbstractArray [,dims=1])

Creates a LinearOperator which, when multiplied with a real array x, returns the DFT over the dimension dims, exploiting Hermitian symmetry.

julia> RDFT(Float64,(10,10))
ℱ  ℝ^(10, 10) -> ℂ^(6, 10)

julia> RDFT((10,10,10),2)
ℱ  ℝ^(10, 10, 10) -> ℂ^(10, 6, 10)

source

AbstractOperators.IRDFT — Type.

IRDFT([domainType=Float64::Type,] dim_in::Tuple, d::Int, [,dims=1])

IRDFT(x::AbstractArray, d::Int, [,dims=1])

Creates a LinearOperator which, when multiplied with a complex array x, returns the IDFT over the dimension dims, exploiting Hermitian symmetry. Like in the function BASE.irfft, d must satisfy div(d,2)+1 == size(x,dims).

julia> A = IRDFT(Complex{Float64},(10,),19)
ℱ⁻¹  ℂ^10 -> ℝ^19 

julia> A = IRDFT((5,10,8),19,2)
ℱ⁻¹  ℂ^(5, 10, 8) -> ℝ^(5, 19, 8)

source

AbstractOperators.DCT — Type.

DCT([domainType=Float64::Type,] dim_in::Tuple)

DCT(dim_in...)

DCT(x::AbstractArray)

Creates a LinearOperator which, when multiplied with an array x::AbstractArray{N}, returns the N-dimensional Inverse Discrete Cosine Transform of x.

julia> DCT(Complex{Float64},(10,10))
ℱc  ℂ^(10, 10) -> ℂ^(10, 10) 

julia> DCT(10,10)
ℱc  ℝ^(10, 10) -> ℂ^(10, 10) 

julia> A = DCT(ones(3))
ℱc  ℝ^3 -> ℝ^3

julia> A*ones(3)
3-element Array{Float64,1}:
 1.73205
 0.0
 0.0

source

AbstractOperators.IDCT — Type.

IDCT([domainType=Float64::Type,] dim_in::Tuple)

IDCT(dim_in...)

IDCT(x::AbstractArray)

Creates a LinearOperator which, when multiplied with an array x::AbstractArray{N}, returns the N-dimensional Discrete Cosine Transform of x.

julia> IDCT(Complex{Float64},(10,10))
ℱc⁻¹  ℂ^(10, 10) -> ℂ^(10, 10) 

julia> IDCT(10,10)
ℱc⁻¹  ℝ^(10, 10) -> ℂ^(10, 10) 

julia> A = IDCT(ones(3))
ℱc⁻¹  ℝ^3 -> ℝ^3

julia> A*[1.;0.;0.]
3-element Array{Float64,1}:
 0.57735
 0.57735
 0.57735

source

Convolution

AbstractOperators.Conv — Type.

Conv([domainType=Float64::Type,] dim_in::Tuple, h::AbstractVector)

Conv(x::AbstractVector, h::AbstractVector)

Creates a LinearOperator which, when multiplied with an array x::AbstractVector, returns the convolution between x and h. Uses conv and hence FFT algorithm.

source

AbstractOperators.Xcorr — Type.

Xcorr([domainType=Float64::Type,] dim_in::Tuple, h::AbstractVector)

Xcorr(x::AbstractVector, h::AbstractVector)

Creates a LinearOperator which, when multiplied with an array x::AbstractVector, returns the cross correlation between x and h. Uses xcross.

source

AbstractOperators.Filt — Type.

Filt([domainType=Float64::Type,] dim_in::Tuple, b::AbstractVector, [a::AbstractVector,])

Filt(x::AbstractVector, b::AbstractVector, [a::AbstractVector,])

Creates a LinearOperator which, when multiplied with an array x::AbstractVector, returns a vector y filtered by an IIR filter of coefficients b and a. If only b is provided a FIR is used to comute y instead.

source

AbstractOperators.MIMOFilt — Type.

MIMOFilt([domainType=Float64::Type,] dim_in::Tuple, B::Vector{AbstractVector}, [A::Vector{AbstractVector},])

MIMOFilt(x::AbstractMatrix, b::Vector{AbstractVector}, [a::Vector{AbstractVector},])

Creates a LinearOperator which, when multiplied with a matrix X, returns a matrix Y. Here a Multiple Input Multiple Output system is evaluated: the columns of X and Y represent the input signals and output signals respectively.

\[\mathbf{y}_i = \sum_{j = 1}^{M} \mathbf{h}_{i,j} * \mathbf{x}_j \]

where $\mathbf{y}_i$ and $\mathbf{x}_j$ are the $i$-th and $j$-th columns of the output Y and input X matrices respectively.

The filters $\mathbf{h}_{i,j}$ can be represented either by providing coefficients B and A (IIR) or B alone (FIR). These coefficients must be given in a Vector of Vectors.

For example for a 3 by 2 MIMO system (i.e. size(X,2) == 3 inputs and size(Y,2) == 2 outputs) B must be:

B = [b11, b12, b13, b21, b22, b23]

where bij are vector containing the filter coeffients of h_{i,j}.

julia> m,n = 10,3; #time samples, number of inputs

julia> B  = [[1.;0.;1.],[1.;0.;1.],[1.;0.;1.],[1.;0.;1.],[1.;0.;1.],[1.;0.;1.], ];
      #B = [   b11   ,     b12   ,    b13   ,   b21    ,   b22,       b23    , ]

julia> A  = [[1.;1.;1.],[2.;2.;2.],[      3.],[      4.],[      5.],[      6.], ];
      #A = [   a11   ,     a12   ,    a13   ,   a21    ,   a22,       a23    , ]

julia> op = MIMOFilt(Float64, (m,n), B, A)
※  ℝ^(10, 3) -> ℝ^(10, 2) 

julia> X = randn(m,n); #input signals

julia> Y = op*X;       #output signals

julia> Y[:,1] ≈ filt(B[1],A[1],X[:,1])+filt(B[2],A[2],X[:,2])+filt(B[3],A[3],X[:,3])
true

julia> Y[:,2] ≈ filt(B[4],A[4],X[:,1])+filt(B[5],A[5],X[:,2])+filt(B[6],A[6],X[:,3])
true

source

AbstractOperators.ZeroPad — Type.

ZeroPad([domainType::Type,] dim_in::Tuple, zp::Tuple)

ZeroPad(x::AbstractArray, zp::Tuple)

Create a LinearOperator which, when multiplied to an array x of size dim_in, returns an expanded array y of size dim_in .+ zp where y[1:dim_in[1], 1:dim_in[2] ... ] = x and zero elsewhere.

julia> Z = ZeroPad((2,2),(0,2))
[I;0]  ℝ^(2, 2) -> ℝ^(2, 4)

julia> Z*ones(2,2)
2×4 Array{Float64,2}:
 1.0  1.0  0.0  0.0
 1.0  1.0  0.0  0.0

source

Finite Differences

AbstractOperators.FiniteDiff — Type.

FiniteDiff([domainType=Float64::Type,] dim_in::Tuple, direction = 1)

FiniteDiff(x::AbstractArray, direction = 1)

Creates a LinearOperator which, when multiplied with an array x::AbstractArray{N}, returns the discretized gradient over the specified direction obtained using forward finite differences.

julia> FiniteDiff(Float64,(3,))
δx  ℝ^3 -> ℝ^2

julia> FiniteDiff((3,4),2)
δy  ℝ^(3, 4) -> ℝ^(3, 3)

julia> all(FiniteDiff(ones(2,2,2,3),1)*ones(2,2,2,3) .== 0)
true

source

AbstractOperators.Variation — Type.

Variation([domainType=Float64::Type,] dim_in::Tuple)

Variation(dims...)

Variation(x::AbstractArray)

Creates a LinearOperator which, when multiplied with an array x::AbstractArray{N}, returns a matrix with its ith column consisting of the vectorized discretized gradient over the ith `direction obtained using forward finite differences.

julia> Variation(Float64,(10,2))
Ʋ  ℝ^(10, 2) -> ℝ^(20, 2)

julia> Variation(2,2,2)
Ʋ  ℝ^(2, 2, 2) -> ℝ^(8, 3)

julia> Variation(ones(2,2))*[1. 2.; 1. 2.]
4×2 Array{Float64,2}:
 0.0  1.0
 0.0  1.0
 0.0  1.0
 0.0  1.0

source

Optimization

AbstractOperators.LBFGS — Type.

LBFGS(domainType::Type,dim_in::Tuple, M::Integer)

LBFGS(dim_in::Tuple, M::Integer)

LBFGS(x::AbstractArray, M::Integer)

Construct a Limited-Memory BFGS LinearOperator with memory M. The memory of LBFGS can be updated using the function update!, where the current iteration variable and gradient (x, grad) and the previous ones (x_prev and grad_prev) are needed:

julia> L = LBFGS(Float64,(4,),5)
LBFGS  ℝ^4 -> ℝ^4

julia> update!(L,x,x_prev,grad,grad_prev); # update memory

julia> d = L*grad; # compute new direction

Use reset!(L) to cancel the memory of the operator.

source