Merge pull request #1 from JuliaImageRecon/tk/extensions

Redesign using extensions

Author: tknopp

.github/workflows/CI.yml

@@ -18,8 +18,7 @@ jobs:
       fail-fast: false
       matrix:
         version:
-          - '1.6'
-          - '1.9'
+          - '1' # Leave this line unchanged. '1' will automatically expand to the latest stable 1.x release of Julia.
           - 'nightly'
         os:
           - ubuntu-latest

Project.toml

@@ -4,23 +4,33 @@ authors = ["Tobias Knopp <tobias.knopp@tuhh.de> and contributors"]
 version = "1.0.0"
 
 [deps]
-FFTW = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 LinearOperators = "5c8ed15e-5a4c-59e4-a42b-c7e8811fb125"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+InteractiveUtils = "b77e0a4c-d291-57a0-90e8-8db25a27a240"
+
+[weakdeps]
 Wavelets = "29a6e085-ba6d-5f35-a997-948ac2efa89a"
+FFTW = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341"
+NFFT = "efe261a4-0d2b-5849-be55-fc731d526b0d"
+
+[extensions]
+LinearOperatorWaveletExt = "Wavelets"
+LinearOperatorFFTWExt = "FFTW"
+LinearOperatorNFFTExt = ["NFFT", "FFTW"]
 
 [compat]
-julia = "1.6"
-FFTW = "0.2, 1.0"
+julia = "1.9"
+FFTW = "1.0"
 LinearOperators = "2.3.3"
-Reexport = "0.2, 1.0"
-Wavelets = "0.8, 0.9"
+Reexport = "1.0"
+Wavelets = "0.9"
+NFFT = "0.13"
 
 [extras]
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [targets]
-test = ["Test"]
+test = ["Test", "FFTW", "Wavelets", "NFFT"]

README.md

@@ -2,4 +2,22 @@
 
 [![Build Status](https://github.com/JuliaImageRecon/LinearOperatorCollection.jl/actions/workflows/CI.yml/badge.svg?branch=main)](https://github.com/JuliaImageRecon/LinearOperatorCollection.jl/actions/workflows/CI.yml?query=branch%3Amain)
 
-[![codecov.io](http://codecov.io/JuliaImageRecon/LinearOperatorCollection.jl/coverage.svg?branch=master)](http://codecov.io/github/JuliaImageRecon/LinearOperatorCollection.jl?branch=master)
+[![codecov.io](http://codecov.io/JuliaImageRecon/LinearOperatorCollection.jl/coverage.svg?branch=master)](http://codecov.io/github/JuliaImageRecon/LinearOperatorCollection.jl?branch=master)
+
+## Purpose
+
+This package contains a collection of linear operators that are in particular useful for multi-dimensional signal and image processing tasks. All operators are build using the LinearOperators.jl base type and derive from `AbstractLinearOperator`. For example this package
+provides operators for the FFT (Fast Fourier Transform) and its non-equidistant variant (NFFT), the DCT (Discrete Cosine Transform), and the Wavelet transform. This package, however, does not implement
+these transformation itself but uses established libraries for them. So in fact, LinearOperatorCollection's main purpose is to add a wrapper around low-level libraries like
+FFTW.jl and NFFT.jl, which allows to use the transformation as if they would be linear operators, i.e. implement `Op * x`, `adjoint(Op) * x` and the `mul!` based in-place variants of the former.
+
+## Installation
+
+Within Julia, use the package manager to install this package:
+```julia
+using Pkg
+Pkg.add("LinearOperatorCollection")
+```
+
+## Usage
+...

src/DCTOp.jl ext/LinearOperatorFFTWExt/DCTOp.jlsrc/DCTOp.jl renamed to ext/LinearOperatorFFTWExt/DCTOp.jl

@@ -1,6 +1,6 @@
-export DCTOp
+export DCTOpImpl
 
-mutable struct DCTOp{T} <: AbstractLinearOperator{T}
+mutable struct DCTOpImpl{T} <: DCTOp{T}
   nrow :: Int
   ncol :: Int
   symmetric :: Bool
@@ -20,19 +20,19 @@ mutable struct DCTOp{T} <: AbstractLinearOperator{T}
   dcttype::Int
 end
 
-LinearOperators.storage_type(op::DCTOp) = typeof(op.Mv5)
+LinearOperators.storage_type(op::DCTOpImpl) = typeof(op.Mv5)
 
 """
-  DCTOp(T::Type, shape::Tuple, dcttype=2)
+  DCTOpImpl(T::Type, shape::Tuple, dcttype=2)
 
-returns a `DCTOp <: AbstractLinearOperator` which performs a DCT on a given input array.
+returns a `DCTOpImpl <: AbstractLinearOperator` which performs a DCT on a given input array.
 
 # Arguments:
 * `T::Type`       - type of the array to transform
 * `shape::Tuple`  - size of the array to transform
 * `dcttype`       - type of DCT (currently `2` and `4` are supported)
 """
-function DCTOp(T::Type, shape::Tuple, dcttype=2)
+function LinearOperatorCollection.DCTOp(T::Type; shape::Tuple, dcttype=2)
 
   tmp=Array{Complex{real(T)}}(undef, shape) 
   if dcttype == 2
@@ -50,7 +50,7 @@ function DCTOp(T::Type, shape::Tuple, dcttype=2)
     error("DCT type $(dcttype) not supported")
   end
 
-  return DCTOp{T}(prod(shape), prod(shape), false, false,
+  return DCTOpImpl{T}(prod(shape), prod(shape), false, false,
                       prod!, nothing, tprod!,
                       0, 0, 0, true, false, true, T[], T[],
                       plan, dcttype)
@@ -68,6 +68,6 @@ function dct_multiply4(res::Vector{T}, plan::P, x::Vector{T}, tmp::Array{T,D}, f
   res .= factor.*vec(tmp)
 end
 
-function Base.copy(S::DCTOp)
-  return DCTOp(eltype(S), size(S.plan), S.dcttype)
+function Base.copy(S::DCTOpImpl)
+  return DCTOpImpl(eltype(S), size(S.plan), S.dcttype)
 end

src/DSTOp.jl ext/LinearOperatorFFTWExt/DSTOp.jlsrc/DSTOp.jl renamed to ext/LinearOperatorFFTWExt/DSTOp.jl

@@ -1,6 +1,6 @@
-export DSTOp
+export DSTOpImpl
 
-mutable struct DSTOp{T} <: AbstractLinearOperator{T}
+mutable struct DSTOpImpl{T} <: DSTOp{T}
   nrow :: Int
   ncol :: Int
   symmetric :: Bool
@@ -20,7 +20,7 @@ mutable struct DSTOp{T} <: AbstractLinearOperator{T}
   iplan
 end
 
-LinearOperators.storage_type(op::DSTOp) = typeof(op.Mv5)
+LinearOperators.storage_type(op::DSTOpImpl) = typeof(op.Mv5)
 
 """
   DSTOp(T::Type, shape::Tuple)
@@ -31,15 +31,15 @@ returns a `LinearOperator` which performs a DST on a given input array.
 * `T::Type`       - type of the array to transform
 * `shape::Tuple`  - size of the array to transform
 """
-function DSTOp(T::Type, shape::Tuple)
+function LinearOperatorCollection.DSTOp(T::Type; shape::Tuple)
   tmp=Array{Complex{real(T)}}(undef, shape)
 
   plan = FFTW.plan_r2r!(tmp,FFTW.RODFT10)
   iplan = FFTW.plan_r2r!(tmp,FFTW.RODFT01)
 
   w = weights(shape, T)
 
-  return DSTOp{T}(prod(shape), prod(shape), true, false
+  return DSTOpImpl{T}(prod(shape), prod(shape), true, false
             , (res,x) -> dst_multiply!(res,plan,x,tmp,w)
             , nothing
             , (res,x) -> dst_bmultiply!(res,iplan,x,tmp,w)
@@ -72,6 +72,6 @@ function dst_bmultiply!(res::Vector{T}, plan::P, x::Vector{T}, tmp::Array{T,D},
   res[:] .= vec(tmp)./(8*length(tmp))
 end
 
-function Base.copy(S::DSTOp)
-  return DSTOp(eltype(S), size(S.plan))
+function Base.copy(S::DSTOpImpl)
+  return DSTOpImpl(eltype(S), size(S.plan))
 end

src/FFTOp.jl ext/LinearOperatorFFTWExt/FFTOp.jlsrc/FFTOp.jl renamed to ext/LinearOperatorFFTWExt/FFTOp.jl

@@ -1,7 +1,7 @@
-export FFTOp
+export FFTOpImpl
 import Base.copy
 
-mutable struct FFTOp{T} <: AbstractLinearOperator{T}
+mutable struct FFTOpImpl{T} <: FFTOp{T}
   nrow :: Int
   ncol :: Int
   symmetric :: Bool
@@ -23,10 +23,10 @@ mutable struct FFTOp{T} <: AbstractLinearOperator{T}
   unitary::Bool
 end
 
-LinearOperators.storage_type(op::FFTOp) = typeof(op.Mv5)
+LinearOperators.storage_type(op::FFTOpImpl) = typeof(op.Mv5)
 
 """
-  FFTOp(T::Type, shape::Tuple, shift=true, unitary=true)
+  FFTOp(T::Type; shape::Tuple, shift=true, unitary=true)
 
 returns an operator which performs an FFT on Arrays of type T
 
@@ -36,7 +36,7 @@ returns an operator which performs an FFT on Arrays of type T
 * (`shift=true`)  - if true, fftshifts are performed
 * (`unitary=true`)  - if true, FFT is normalized such that it is unitary
 """
-function FFTOp(T::Type, shape::NTuple{D,Int64}, shift::Bool=true; unitary::Bool=true, cuda::Bool=false) where D
+function LinearOperatorCollection.FFTOp(T::Type; shape::NTuple{D,Int64}, shift::Bool=true, unitary::Bool=true, cuda::Bool=false) where D
 
   #tmpVec = cuda ? CuArray{T}(undef,shape) : Array{Complex{real(T)}}(undef, shape)
   tmpVec = Array{Complex{real(T)}}(undef, shape)
@@ -54,7 +54,7 @@ function FFTOp(T::Type, shape::NTuple{D,Int64}, shift::Bool=true; unitary::Bool=
   let shape_=shape, plan_=plan, iplan_=iplan, tmpVec_=tmpVec, facF_=facF, facB_=facB
 
   if shift
-    return FFTOp{T}(prod(shape), prod(shape), false, false
+    return FFTOpImpl{T}(prod(shape), prod(shape), false, false
               , (res, x) -> fft_multiply_shift!(res, plan_, x, shape_, facF_, tmpVec_) 
               , nothing
               , (res, x) -> fft_multiply_shift!(res, iplan_, x, shape_, facB_, tmpVec_) 
@@ -64,7 +64,7 @@ function FFTOp(T::Type, shape::NTuple{D,Int64}, shift::Bool=true; unitary::Bool=
               , shift
               , unitary)
   else
-    return FFTOp{T}(prod(shape), prod(shape), false, false
+    return FFTOpImpl{T}(prod(shape), prod(shape), false, false
             , (res, x) -> fft_multiply!(res, plan_, x, facF_, tmpVec_) 
             , nothing
             , (res, x) -> fft_multiply!(res, iplan_, x, facB_, tmpVec_)
@@ -91,6 +91,6 @@ function fft_multiply_shift!(res::AbstractVector{T}, plan::P, x::AbstractVector{
 end
 
 
-function Base.copy(S::FFTOp)
-  return FFTOp(eltype(S), size(S.plan), S.shift, unitary=S.unitary)
+function Base.copy(S::FFTOpImpl)
+  return FFTOpImpl(eltype(S), size(S.plan), S.shift, unitary=S.unitary)
 end

ext/LinearOperatorFFTWExt/LinearOperatorFFTWExt.jl

@@ -0,0 +1,9 @@
+module LinearOperatorFFTWExt
+
+using LinearOperatorCollection, FFTW
+
+include("FFTOp.jl")
+include("DCTOp.jl")
+include("DSTOp.jl")
+
+end

ext/LinearOperatorNFFTExt/LinearOperatorNFFTExt.jl

@@ -0,0 +1,7 @@
+module LinearOperatorNFFTExt
+
+using LinearOperatorCollection, NFFT, FFTW
+
+include("NFFTOp.jl")
+
+end