Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Setup SIWT data structures #51

Merged
merged 17 commits into from
Jun 7, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions .github/workflows/CI.yml
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,8 @@ on:
- test/**
- Project.toml
pull_request:
branches:
- master
paths:
- src/**
- test/**
Expand Down
2 changes: 1 addition & 1 deletion docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ makedocs(
"DWT" => "api/dwt.md",
"ACWT" => "api/acwt.md",
"SWT" => "api/swt.md",
"SIWPD" => "api/siwpd.md",
"SIWT" => "api/siwt.md",
"WaveMult" => "api/wavemult.md",
"Best Basis" => "api/bestbasis.md",
"Denoising" => "api/denoising.md",
Expand Down
4 changes: 0 additions & 4 deletions docs/src/api/bestbasis.md
Original file line number Diff line number Diff line change
Expand Up @@ -17,13 +17,11 @@ BestBasis.ShannonEntropyCost
BestBasis.LogEnergyEntropyCost
BestBasis.coefcost
BestBasis.tree_costs
BestBasis.tree_costs(::AbstractMatrix{T}, ::AbstractVector{BitVector}, ::SIBB) where T<:Number
```

# Best Basis Tree Selection
```@docs
BestBasis.bestbasis_treeselection
BestBasis.bestbasis_treeselection(::AbstractVector{Tc}, ::AbstractVector{Tt}) where {Tc<:AbstractVector{<:Union{Number,Nothing}}, Tt<:BitVector}
BestBasis.delete_subtree!
```

Expand All @@ -33,8 +31,6 @@ BestBasis.BestBasisType
BestBasis.LSDB
BestBasis.JBB
BestBasis.BB
BestBasis.SIBB
Wavelets.Threshold.bestbasistree
Wavelets.Threshold.bestbasistree(::AbstractMatrix{T}, ::Integer, ::SIBB) where T<:Number
BestBasis.bestbasistreeall
```
11 changes: 0 additions & 11 deletions docs/src/api/siwpd.md

This file was deleted.

49 changes: 49 additions & 0 deletions docs/src/api/siwt.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,49 @@
# [Shift Invariant Wavelet Packet Decomposition](@id siwt_api)

```@index
Modules = [SIWT]
```

## Data Structures
```@docs
SIWT.ShiftInvariantWaveletTransformNode
SIWT.ShiftInvariantWaveletTransformObject
```

## Signal Transform and Reconstruction
### Public API
```@docs
SIWT.siwpd
SIWT.isiwpd
```

### Private API
```@docs
SIWT.siwpd_subtree!
SIWT.isiwpd_subtree!
```

## Best Basis Search
### Public API
```@docs
SIWT.bestbasistree!
```

### Private API
```@docs
SIWT.bestbasis_treeselection!
```

## Single Step Transforms
### Private API
```@docs
SIWT.sidwt_step!
SIWT.isidwt_step!
```

## Other Utils
### Private API
```@docs
Wavelets.Util.isvalidtree(::ShiftInvariantWaveletTransformObject)
SIWT.delete_node!
```
2 changes: 1 addition & 1 deletion docs/src/api/utils.md
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,7 @@ Utils.finestdetailrange

## Tree traversing functions
```@docs
Wavelets.Util.isvalidtree
Wavelets.Util.isvalidtree(::AbstractMatrix,::BitVector)
Wavelets.Util.maketree
Utils.getchildindex
Utils.getparentindex
Expand Down
22 changes: 12 additions & 10 deletions docs/src/manual/bestbasis.md
Original file line number Diff line number Diff line change
Expand Up @@ -126,13 +126,15 @@ plot(p1, p2, p3, p4, p5, p6, layout=(3,2))

## [Best Basis of Shift-Invariant Wavelet Packet Decomposition](@id si_bestbasis)
One can think of searching for the best basis of the shift-invariant wavelet packet decomposition as a problem of finding ``\min_{b \in B} \sum_{x \in X} M_x(b)``, where ``X`` is all the possible shifted versions of an original signal ``y``. One can compute the best basis tree as follows:
```@example wt
xw = siwpd(x, wt)

# SIBB
tree = bestbasistree(xw, 7, SIBB());
nothing #hide
```

!!! warning
SIWPD is still undergoing large changes in terms of data structures and efficiency improvements. Syntax changes may occur in the next patch updates.
```@repl
x = [2,3,-4,5.0];
wt = wavelet(WT.haar);
xwObj = siwpd(x, wt, 1, 1);
xwObj.BestTree # Original tree (all decomposed nodes)

bestbasistree!(xwObj)
xwObj.BestTree # Best basis tree

x̂ = isiwpd(xwObj) # Reconstruction of signal
x̂ == x
```
13 changes: 7 additions & 6 deletions docs/src/manual/transforms.md
Original file line number Diff line number Diff line change
Expand Up @@ -142,14 +142,15 @@ plot(plot(img, title="Original"), p0, p1, p2, layout=(2,2))
The [Shift-Invariant Wavelet Decomposition (SIWPD)](https://israelcohen.com/wp-content/uploads/2018/05/ICASSP95.pdf) is developed by Cohen et. al.. While it is also a type of redundant transform, it does not follow the same methodology as the SWT and the ACWT. Cohen's main goal for developing this algorithm was to obtain a global minimum entropy from a signal and all its shifted versions. See [its best basis implementation](@ref si_bestbasis) for more information.

One can compute the SIWPD of a single signal as follows.
```@example wt
# decomposition
xw = siwpd(x, wt);
nothing # hide
```@repl
x = [2,3,-4,5.0];
wt = wavelet(WT.haar);
xwObj = siwpd(x, wt, 1, 1);
xwObj.Nodes[(0,0,0)] # Example of looking into a node
xwObj.BestTree # Looking into what's within the tree structure
```
For more information on the API, visit the [SIWT API page](@ref siwt_api).

!!! note
As of right now, there is not too many functions written based on the SIWPD, as it does not follow the conventional style of wavelet transforms. There is a lot of ongoing work to develop more functions catered for the SIWPD such as it's inverse transforms and group-implementations.



Expand Down
4 changes: 2 additions & 2 deletions src/WaveletsExt.jl
Original file line number Diff line number Diff line change
Expand Up @@ -4,9 +4,9 @@ module WaveletsExt

include("mod/Utils.jl")
include("mod/DWT.jl")
include("mod/SIWPD.jl")
include("mod/ACWT.jl")
include("mod/BestBasis.jl")
include("mod/SIWT.jl")
include("mod/SWT.jl")
include("mod/Denoising.jl")
include("mod/LDB.jl")
Expand All @@ -20,7 +20,7 @@ using Reexport
.Utils,
.LDB,
.SWT,
.SIWPD,
.SIWT,
.ACWT,
.Visualizations,
.WaveMult
Expand Down
125 changes: 2 additions & 123 deletions src/mod/BestBasis.jl
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,6 @@ export
LSDB,
JBB,
BB,
SIBB,
# best basis tree
bestbasistreeall

Expand All @@ -30,8 +29,7 @@ using

using
..Utils,
..DWT,
..SIWPD
..DWT

include("bestbasis/bestbasis_costs.jl")
include("bestbasis/bestbasis_tree.jl")
Expand Down Expand Up @@ -111,72 +109,6 @@ function bestbasis_treeselection(costs::AbstractVector{T},
return tree
end

# SIWPD tree selection
"""
bestbasis_treeselection(costs, tree)

Best basis tree selection on SIWPD.

# Arguments
- `costs::AbstractVector{Tc} where Tc<:AbstractVector{<:Union{Number, Nothing}}`: Cost of
each node.
- `tree::AbstractVector{Tt} where Tt<:BitVector`: SIWPD tree.

!!! warning
Current implementation works but is unstable, ie. we are still working on better
syntax/more optimized computations/better data structure.
"""
function bestbasis_treeselection(costs::AbstractVector{Tc},
tree::AbstractVector{Tt}) where
{Tc<:AbstractVector{<:Union{Number,Nothing}}, Tt<:BitVector}

@assert length(costs) == length(tree)
bt = deepcopy(tree)
bc = deepcopy(costs)
nn = length(tree)
for i in reverse(eachindex(bt))
if getchildindex(i,:left) > nn # current node is at bottom level
continue
end
level = floor(Int, log2(i))
for j in eachindex(bt[i]) # iterate through all available shifts
if !bt[i][j] # node of current shift does not exist
continue
end
@assert bt[getchildindex(i,:left)][j] ==
bt[getchildindex(i,:right)][j] ==
bt[getchildindex(i,:left)][j+1<<level] ==
bt[getchildindex(i,:right)][j+1<<level] == true || continue
# child cost of current shift of current node
cc = bc[getchildindex(i,:left)][j] + bc[getchildindex(i,:right)][j]
# child cost of circshift of current node
scc = bc[getchildindex(i,:left)][j+1<<level] +
bc[getchildindex(i,:right)][j+1<<level]
mincost = min(cc, scc, bc[i][j])
# mincost=parent, delete subtrees of both cc & scc
if mincost == bc[i][j]
delete_subtree!(bt, getchildindex(i,:left), j)
delete_subtree!(bt, getchildindex(i,:right), j)
delete_subtree!(bt, getchildindex(i,:left), j+1<<level)
delete_subtree!(bt, getchildindex(i,:right), j+1<<level)
# mincost=cc, delete subtrees of scc
elseif mincost == cc
bc[i][j] = mincost
delete_subtree!(bt, getchildindex(i,:left), j+1<<level)
delete_subtree!(bt, getchildindex(i,:right), j+1<<level)
# mincost=scc, delete subtrees of cc
else
bc[i][j] = mincost
delete_subtree!(bt, getchildindex(i,:left), j)
delete_subtree!(bt, getchildindex(i,:right), j)
end
end
end
# ensure only 1 node selected from each Ω(i,j)
@assert all(map(node -> sum(node), bt) .<= 1)
return bt
end

# Deletes subtree due to inferior cost
"""
delete_subtree!(bt, i, tree_type)
Expand Down Expand Up @@ -207,28 +139,6 @@ function delete_subtree!(bt::BitVector, i::Integer, tree_type::Symbol)
return bt
end

function delete_subtree!(bt::AbstractVector{BitVector}, i::Integer, j::Integer)
@assert 1 <= i <= length(bt)
level = floor(Int, log2(i))
bt[i][j] = false
if (getchildindex(i,:left)) < length(bt) # current node can contain subtrees
if bt[getchildindex(i,:left)][j] # left child of current shift
delete_subtree!(bt, getchildindex(i,:left), j)
end
if bt[getchildindex(i,:right)][j] # right child of current shift
delete_subtree!(bt, getchildindex(i,:right), j)
end
if bt[getchildindex(i,:left)][j+1<<level] # left child of added shift
delete_subtree!(bt, getchildindex(i,:left), j+1<<level)
end
if bt[getchildindex(i,:right)][j+1<<level] # right child of added shift
delete_subtree!(bt, getchildindex(i,:right), j+1<<level)
end
end
return nothing
end


## BEST BASIS TREES
# Best basis tree for LSDB
"""
Expand All @@ -237,8 +147,7 @@ end
Extension to the best basis tree function from Wavelets.jl. Given a set of decomposed
signals, returns different types of best basis trees based on the methods specified.
Available methods are the joint best basis ([`JBB`](@ref)), least statistically dependent
basis ([`LSDB`](@ref)), individual regular best basis ([`BB`](@ref)), and shift-invariant
best basis ([`SIBB`](@ref)).
basis ([`LSDB`](@ref)), and individual regular best basis ([`BB`](@ref)).

# Arguments
- `X::AbstractArray{T} where T<:AbstractFloat`: A set of decomposed signals, of sizes
Expand Down Expand Up @@ -299,37 +208,7 @@ function Wavelets.Threshold.bestbasistree(X::AbstractArray{T}, method::BB) where
tree = bestbasis_treeselection(costs, sz...)
return tree
end
# SIWPD Best basis tree
# TODO: find a way to compute bestbasis_tree without input d
"""
bestbasistree(y, d, method)

Computes the best basis tree for the shift invariant wavelet packet decomposition (SIWPD).

# Arguments
- `y::AbstractArray{T,2} where T<:Number`: A SIWPD decomposed signal.
- `d::Integer`: The number of depth computed for the decomposition.
- `method::SIBB`: The `SIBB()` method.

# Returns
- `Vector{BitVector}`: SIWPD best basis tree.

!!! warning
Current implementation works but is unstable, ie. we are still working on better
syntax/more optimized computations/better data structure.
"""
function Wavelets.Threshold.bestbasistree(y::AbstractArray{T,2},
d::Integer,
method::SIBB) where T<:Number

nn = size(y,2)
L = maxtransformlevels((nn+1)÷2)
ns = size(y,1)
tree = makesiwpdtree(ns, L, d)
costs = tree_costs(y, tree, method)
besttree = bestbasis_treeselection(costs, tree)
return besttree
end
# Default best basis tree search
function Wavelets.Threshold.bestbasistree(X::AbstractArray{T},
method::BestBasisType = JBB()) where
Expand Down
Loading