From 47f8f1ab8136e2ec9d3d3bd7465674ddf4beff5f Mon Sep 17 00:00:00 2001 From: nHackel Date: Fri, 31 May 2024 18:15:29 +0200 Subject: [PATCH 1/9] Init working GPU reconstruction --- .../SinglePatchAlgorithms/SinglePatchAlgorithm.jl | 1 + src/LeastSquares.jl | 4 ++-- src/SystemMatrix/SystemMatrix.jl | 6 +++++- 3 files changed, 8 insertions(+), 3 deletions(-) diff --git a/src/Algorithms/SinglePatchAlgorithms/SinglePatchAlgorithm.jl b/src/Algorithms/SinglePatchAlgorithms/SinglePatchAlgorithm.jl index 2799548..9928070 100644 --- a/src/Algorithms/SinglePatchAlgorithms/SinglePatchAlgorithm.jl +++ b/src/Algorithms/SinglePatchAlgorithms/SinglePatchAlgorithm.jl @@ -44,6 +44,7 @@ function process(algo::SinglePatchReconstructionAlgorithm, params::Union{A, Proc @warn "System matrix and measurement have different element data type. Mapping measurment data to system matrix element type." result = map(eltype(algo.S),result) end + result = copyto!(similar(algo.S, size(result)...), result) return result end diff --git a/src/LeastSquares.jl b/src/LeastSquares.jl index a383360..d252928 100644 --- a/src/LeastSquares.jl +++ b/src/LeastSquares.jl @@ -51,7 +51,7 @@ Base.propertynames(params::RecoPlan{ElaborateSolverParameters}) = union([:solver getSolverKwargs(::Type{SL}) where SL <: AbstractLinearSolver = intersect(union(Base.kwarg_decl.(methods(SL))...), fieldnames(ElaborateSolverParameters)) -function process(t::Type{<:AbstractMPIRecoAlgorithm}, params::LeastSquaresParameters{SL}, u::Array) where SL +function process(t::Type{<:AbstractMPIRecoAlgorithm}, params::LeastSquaresParameters{SL}, u::AbstractArray) where SL N = size(params.S, 2) M = div(length(params.S), N) @@ -79,7 +79,7 @@ function process(t::Type{<:AbstractMPIRecoAlgorithm}, params::LeastSquaresParame if !isnothing(params.op) d[:] = params.op*d end - c[:, l] = real(d) + c[:, l] = Array(real(d)) end return c diff --git a/src/SystemMatrix/SystemMatrix.jl b/src/SystemMatrix/SystemMatrix.jl index 8a2febc..7e9f061 100644 --- a/src/SystemMatrix/SystemMatrix.jl +++ b/src/SystemMatrix/SystemMatrix.jl @@ -52,11 +52,12 @@ export AbstractSystemMatrixLoadingParameter abstract type AbstractSystemMatrixLoadingParameter <: AbstractSystemMatrixParameter end export DenseSystemMatixLoadingParameter -Base.@kwdef struct DenseSystemMatixLoadingParameter{F<:AbstractFrequencyFilterParameter, G<:AbstractSystemMatrixGriddingParameter} <: AbstractSystemMatrixLoadingParameter +Base.@kwdef struct DenseSystemMatixLoadingParameter{F<:AbstractFrequencyFilterParameter, matT <: AbstractArray, G<:AbstractSystemMatrixGriddingParameter} <: AbstractSystemMatrixLoadingParameter freqFilter::F gridding::G bgCorrection::Bool = false loadasreal::Bool = false + arrayType::Type{matT} = Array end function process(t::Type{<:AbstractMPIRecoAlgorithm}, params::DenseSystemMatixLoadingParameter, sf::MPIFile) # Construct freqFilter @@ -66,6 +67,9 @@ end function process(t::Type{<:AbstractMPIRecoAlgorithm}, params::DenseSystemMatixLoadingParameter, sf::MPIFile, frequencies::Vector{CartesianIndex{2}}) S, grid = getSF(sf, frequencies, nothing; toKwargs(params)...) @info "Loading SM" + if !isa(S, params.arrayType) + S = params.arrayType(S) + end return S, grid end From 8887c67c11caaa53e7c30cb54bdd8059497968d0 Mon Sep 17 00:00:00 2001 From: nHackel Date: Thu, 11 Jul 2024 16:40:57 +0200 Subject: [PATCH 2/9] Add default plans to MPIReco --- config/MultiPatch.toml | 42 +++++++++++++++ config/MultiPatchMotion.toml | 26 ++++++++++ config/SinglePatch.toml | 98 +++++++++++++++++++++++++++++++++++ config/SinglePatchSparse.toml | 70 +++++++++++++++++++++++++ src/AlgorithmInterface.jl | 19 ++++++- src/MPIReco.jl | 4 +- 6 files changed, 256 insertions(+), 3 deletions(-) create mode 100644 config/MultiPatch.toml create mode 100644 config/MultiPatchMotion.toml create mode 100644 config/SinglePatch.toml create mode 100644 config/SinglePatchSparse.toml diff --git a/config/MultiPatch.toml b/config/MultiPatch.toml new file mode 100644 index 0000000..2d12da3 --- /dev/null +++ b/config/MultiPatch.toml @@ -0,0 +1,42 @@ +_type = "RecoPlan{MultiPatchReconstructionAlgorithm}" +_module = "MPIReco" + +[parameter] +_type = "RecoPlan{MultiPatchParameters}" +_module = "MPIReco" + + [parameter.reco] + _type = "RecoPlan{MultiPatchReconstructionParameter}" + _module = "MPIReco" + + [parameter.reco.ffPos] + _type = "RecoPlan{DefaultFocusFieldPositions}" + _module = "MPIReco" + + [parameter.reco.freqFilter] + _type = "RecoPlan{SNRThresholdFrequencyFilterParameter}" + _module = "MPIReco" + + [parameter.reco.opParams] + _type = "RecoPlan{RegularMultiPatchOperatorParameter}" + _module = "MPIReco" + + [parameter.reco.ffPosSF] + _type = "RecoPlan{DefaultFocusFieldPositions}" + _module = "MPIReco" + + [parameter.reco.solverParams] + _type = "RecoPlan{SimpleSolverParameters}" + _module = "MPIReco" + + [parameter.pre] + _type = "RecoPlan{CommonPreProcessingParameters}" + _module = "MPIReco" + + [parameter.pre.bgParams] + _type = "RecoPlan{NoBackgroundCorrectionParameters}" + _module = "MPIReco" + + [parameter.post] + _type = "RecoPlan{NoPostProcessing}" + _module = "MPIReco" diff --git a/config/MultiPatchMotion.toml b/config/MultiPatchMotion.toml new file mode 100644 index 0000000..9d6aadd --- /dev/null +++ b/config/MultiPatchMotion.toml @@ -0,0 +1,26 @@ +_type = "RecoPlan{MultiPatchReconstructionAlgorithm}" +_module = "MPIReco" + +[parameter] +_type = "RecoPlan{MultiPatchParameters}" +_module = "MPIReco" + + [parameter.reco] + _type = "RecoPlan{PeriodicMotionReconstructionParameter}" + _module = "MPIReco" + + [parameter.reco.freqFilter] + _type = "RecoPlan{SNRThresholdFrequencyFilterParameter}" + _module = "MPIReco" + + [parameter.reco.solverParams] + _type = "RecoPlan{SimpleSolverParameters}" + _module = "MPIReco" + + [parameter.pre] + _type = "RecoPlan{PeriodicMotionPreProcessing}" + _module = "MPIReco" + + [parameter.post] + _type = "RecoPlan{NoPostProcessing}" + _module = "MPIReco" diff --git a/config/SinglePatch.toml b/config/SinglePatch.toml new file mode 100644 index 0000000..61cd8b1 --- /dev/null +++ b/config/SinglePatch.toml @@ -0,0 +1,98 @@ +_module = "MPIReco" +_type = "RecoPlan{SinglePatchReconstructionAlgorithm}" + +[parameter] +_module = "MPIReco" +_type = "RecoPlan{SinglePatchParameters}" + + [parameter.reco] + _module = "MPIReco" + _type = "RecoPlan{SinglePatchReconstructionParameter}" + + [parameter.reco.sfLoad] + _module = "AbstractImageReconstruction" + _type = "RecoPlan{ProcessResultCache}" + + [parameter.reco.sfLoad.param] + _module = "MPIReco" + _type = "RecoPlan{DenseSystemMatixLoadingParameter}" + + [parameter.reco.sfLoad.param.freqFilter] + _module = "MPIReco" + _type = "RecoPlan{SNRThresholdFrequencyFilterParameter}" + + [parameter.reco.sfLoad.param.gridding] + _module = "MPIReco" + _type = "RecoPlan{SystemMatrixGriddingParameter}" + + + [[parameter.reco._listener.sf]] + field = "fov" + _module = "AbstractImageReconstruction" + _type = "LinkedFieldListener" + plan = ["parameter", "reco", "sfLoad", "param", "gridding"] + + [parameter.reco._listener.sf.fn] + _module = "MPIReco" + _type = "defaultParameterCalibFov" + [[parameter.reco._listener.sf]] + field = "center" + _module = "AbstractImageReconstruction" + _type = "LinkedFieldListener" + plan = ["parameter", "reco", "sfLoad", "param", "gridding"] + + [parameter.reco._listener.sf.fn] + _module = "MPIReco" + _type = "defaultParameterCalibCenter" + [[parameter.reco._listener.sf]] + field = "gridsize" + _module = "AbstractImageReconstruction" + _type = "LinkedFieldListener" + plan = ["parameter", "reco", "sfLoad", "param", "gridding"] + + [parameter.reco._listener.sf.fn] + _module = "MPIReco" + _type = "defaultParameterGridSize" + [[parameter.reco._listener.sf]] + field = "minFreq" + _module = "AbstractImageReconstruction" + _type = "LinkedFieldListener" + plan = ["parameter", "reco", "sfLoad", "param", "freqFilter"] + + [parameter.reco._listener.sf.fn] + _module = "MPIReco" + _type = "defaultParameterMinFreq" + [[parameter.reco._listener.sf]] + field = "maxFreq" + _module = "AbstractImageReconstruction" + _type = "LinkedFieldListener" + plan = ["parameter", "reco", "sfLoad", "param", "freqFilter"] + + [parameter.reco._listener.sf.fn] + _module = "MPIReco" + _type = "defaultParameterMaxFreq" + [[parameter.reco._listener.sf]] + field = "recChannels" + _module = "AbstractImageReconstruction" + _type = "LinkedFieldListener" + plan = ["parameter", "reco", "sfLoad", "param", "freqFilter"] + + [parameter.reco._listener.sf.fn] + _module = "MPIReco" + _type = "defaultParameterRecChannels" + + [parameter.reco.solverParams] + _module = "MPIReco" + _type = "RecoPlan{ElaborateSolverParameters}" + + [parameter.pre] + _module = "AbstractImageReconstruction" + _type = "RecoPlan{ProcessResultCache}" + + [parameter.pre.param] + _module = "MPIReco" + _type = "RecoPlan{CommonPreProcessingParameters}" + + [parameter.post] + _module = "MPIReco" + _type = "RecoPlan{NoPostProcessing}" diff --git a/config/SinglePatchSparse.toml b/config/SinglePatchSparse.toml new file mode 100644 index 0000000..33d8ad3 --- /dev/null +++ b/config/SinglePatchSparse.toml @@ -0,0 +1,70 @@ +_module = "MPIReco" +_type = "RecoPlan{SinglePatchReconstructionAlgorithm}" + +[parameter] +_module = "MPIReco" +_type = "RecoPlan{SinglePatchParameters}" + + [parameter.reco] + _module = "MPIReco" + _type = "RecoPlan{SinglePatchReconstructionParameter}" + + [parameter.reco.sfLoad] + _module = "AbstractImageReconstruction" + _type = "RecoPlan{ProcessResultCache}" + + [parameter.reco.sfLoad.param] + _module = "MPIReco" + _type = "RecoPlan{SparseSystemMatrixLoadingParameter}" + + [parameter.reco.sfLoad.param.freqFilter] + _module = "MPIReco" + _type = "RecoPlan{SNRThresholdFrequencyFilterParameter}" + + [parameter.reco.sfLoad.param.gridding] + _module = "MPIReco" + _type = "RecoPlan{SystemMatrixGriddingParameter}" + + [[parameter.reco._listener.sf]] + field = "minFreq" + _module = "AbstractImageReconstruction" + _type = "LinkedFieldListener" + plan = ["parameter", "reco", "sfLoad", "param", "freqFilter"] + + [parameter.reco._listener.sf.fn] + _module = "MPIReco" + _type = "defaultParameterMinFreq" + [[parameter.reco._listener.sf]] + field = "maxFreq" + _module = "AbstractImageReconstruction" + _type = "LinkedFieldListener" + plan = ["parameter", "reco", "sfLoad", "param", "freqFilter"] + + [parameter.reco._listener.sf.fn] + _module = "MPIReco" + _type = "defaultParameterMaxFreq" + [[parameter.reco._listener.sf]] + field = "recChannels" + _module = "AbstractImageReconstruction" + _type = "LinkedFieldListener" + plan = ["parameter", "reco", "sfLoad", "param", "freqFilter"] + + [parameter.reco._listener.sf.fn] + _module = "MPIReco" + _type = "defaultParameterRecChannels" + + [parameter.reco.solverParams] + _module = "MPIReco" + _type = "RecoPlan{ElaborateSolverParameters}" + + [parameter.pre] + _module = "AbstractImageReconstruction" + _type = "RecoPlan{ProcessResultCache}" + + [parameter.pre.param] + _module = "MPIReco" + _type = "RecoPlan{CommonPreProcessingParameters}" + + [parameter.post] + _module = "MPIReco" + _type = "RecoPlan{NoPostProcessing}" diff --git a/src/AlgorithmInterface.jl b/src/AlgorithmInterface.jl index a7c15d0..963f5c6 100644 --- a/src/AlgorithmInterface.jl +++ b/src/AlgorithmInterface.jl @@ -42,9 +42,26 @@ struct MixedAlgorithm <: ReconstructionAlgorithmType end # TODO recoAlgorithmType # TODO undefined for certain "Algorithm" components #recoAlgorithmTypes(::Type{ConcreteRecoAlgorithm}) = SystemMatrixBasedAlgorithm() -export plandir +export planpath, plandir plandir() = abspath(homedir(), ".mpi", "RecoPlans") +function planpath(name::AbstractString) + for dir in [joinpath(@__DIR__, "..", "config"),plandir()] + filename = joinpath(dir, string(name, ".toml")) + if isfile(filename) + return filename + end + end + throw(ArgumentError("Could not find a suitable MPI reconstruction plan with name $name.\nCustom plans can be stored in $(plandir()).")) +end + +export reconstruct +function reconstruct(name::AbstractString, data::MPIFile; kwargs...) + plan = loadPlan(MPIReco, name, [AbstractImageReconstruction, MPIFiles, MPIReco, RegularizedLeastSquares]) + setAll!(plan; kwargs...) + return reconstruct(build(plan), data) +end + # Check if contains isSystemMatrixBased(::T) where T <: AbstractImageReconstructionAlgorithm = recoAlgorithmTypes(T) isa SystemMatrixBasedAlgorithm isXSpaceBased(::T) where T <: AbstractImageReconstructionAlgorithm = recoAlgorithmTypes(T) isa XSpaceBasedAlgorithm diff --git a/src/MPIReco.jl b/src/MPIReco.jl index 209d6b6..9d80cb8 100644 --- a/src/MPIReco.jl +++ b/src/MPIReco.jl @@ -16,10 +16,10 @@ module MPIReco using DistributedArrays # using TensorDecompositions using IniFile - import LinearAlgebra: ldiv!, \ + import LinearAlgebra: ldiv!, \, mul! # TODO sort out import for Base and AbstractImageReconstruction to avoid boiler plate import Base: put!, take! - import AbstractImageReconstruction: process, parameter + import AbstractImageReconstruction: process, parameter, reconstruct using FFTW using LinearOperatorCollection From 55c8af76ec30bcce5d66ab339fa9e81ab6aaefc3 Mon Sep 17 00:00:00 2001 From: nHackel Date: Thu, 11 Jul 2024 16:41:19 +0200 Subject: [PATCH 3/9] Adapt to kwarg changes in RegLS --- src/LeastSquares.jl | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/src/LeastSquares.jl b/src/LeastSquares.jl index d252928..17d3c96 100644 --- a/src/LeastSquares.jl +++ b/src/LeastSquares.jl @@ -28,15 +28,17 @@ end export ElaborateSolverParameters Base.@kwdef mutable struct ElaborateSolverParameters{SL} <: AbstractSolverParameters{SL} solver::Type{SL} - iterations::Int64=10 - enforceReal::Bool=true - enforcePositive::Bool=true + iterations::Int64 = 10 + enforceReal::Bool = true + enforcePositive::Bool = true + kwargWarning::Bool = true # Union of all kwargs normalizeReg::AbstractRegularizationNormalization = SystemMatrixBasedNormalization() - randomized::Union{Nothing, Bool} = false + randomized::Union{Nothing, Bool, Symbol} = nothing + seed::Union{Nothing, Int64} = nothing shuffleRows::Union{Nothing, Bool} = false rho::Union{Nothing, Float64} = nothing - normalize_rho::Union{Nothing, Bool} = nothing + vary_rho::Union{Nothing, Symbol} = nothing theta::Union{Nothing, Float64} = nothing restart::Union{Nothing, Symbol} = nothing regTrafo::Union{Nothing, Vector{Union{AbstractArray, AbstractLinearOperator}}} = nothing From a715ab3842f056f61c462facd81232352f4bb0f4 Mon Sep 17 00:00:00 2001 From: nHackel Date: Thu, 11 Jul 2024 16:42:45 +0200 Subject: [PATCH 4/9] Adapt tests to default plans --- test/Cartesian.jl | 6 +- test/MotionCompensation.jl | 33 +++++------ test/MultiGradient.jl | 30 +++++----- test/MultiPatch.jl | 40 ++++++------- test/Reconstruction.jl | 103 ++++++++++++++++----------------- test/SMExtrapolation.jl | 70 +++++++++++----------- test/Solvers.jl | 29 +++++----- test/Sparse.jl | 87 ++++++++++++---------------- test/recoPlans/Motion.toml | 26 --------- test/recoPlans/MultiPatch.toml | 42 -------------- test/recoPlans/Single.toml | 34 ----------- test/recoPlans/Sparse.toml | 34 ----------- test/runtests.jl | 2 - 13 files changed, 189 insertions(+), 347 deletions(-) delete mode 100644 test/recoPlans/Motion.toml delete mode 100644 test/recoPlans/MultiPatch.toml delete mode 100644 test/recoPlans/Single.toml delete mode 100644 test/recoPlans/Sparse.toml diff --git a/test/Cartesian.jl b/test/Cartesian.jl index 8378363..ea30800 100644 --- a/test/Cartesian.jl +++ b/test/Cartesian.jl @@ -14,7 +14,7 @@ using MPIReco numPeriodGrouping = 1 maxMixingOrder = -1 - plan = getPlan("Single") + plan = loadPlan(MPIReco, "SinglePatch", [MPIReco, RegularizedLeastSquares, MPIFiles, AbstractImageReconstruction]) setAll!(plan, :sf, bSF) setAll!(plan, :frames, 1:10) setAll!(plan, :numAverages, 10) @@ -49,8 +49,8 @@ using MPIReco setAll!(plan, :sf, bSFProc1) setAll!(plan, :numPeriodGrouping, 1) setAll!(plan, :numPeriodAverages, 1) - plan.parameter.pre.numPeriodGrouping = 100 # SM is preprocessed, only need to do it for meas data - plan.parameter.pre.numPeriodAverages = 65 + setAll!(plan.parameter.pre, :numPeriodGrouping, 100) # SM is preprocessed, only need to do it for meas data + setAll!(plan.parameter.pre, :numPeriodAverages, 65) @time c3 = reconstruct(build(plan), b) exportImage(joinpath(imgdir, "Cartesian3.png"), Array(c3[1,:,:,1,1])) diff --git a/test/MotionCompensation.jl b/test/MotionCompensation.jl index f339fe1..4e9444f 100644 --- a/test/MotionCompensation.jl +++ b/test/MotionCompensation.jl @@ -36,26 +36,25 @@ using MPIReco lambda = 0.01 iterations = 1 - plan = getPlan("Motion") - setAll!(plan, :sf, bSF) - setAll!(plan, :frames, 1:2) - setAll!(plan, :alpha, alpha) - setAll!(plan, :choosePeak, choosePeak) - setAll!(plan, :windowType, windowType) - setAll!(plan, :higherHarmonic, choosePeak) - setAll!(plan, :lambda, lambda) - setAll!(plan, :iterations, iterations) # reconstruction parameter - setAll!(plan, :SNRThresh, 10) - setAll!(plan, :minFreq, 80e3) - setAll!(plan, :recChannels, 1:3) - setAll!(plan, :λ, 0.0f0) - - c = reconstruct(build(plan), bMeas) + params = Dict{Symbol, Any}() + params[:sf] = bSF + params[:frames] = 1:2 + params[:alpha] = alpha + params[:choosePeak] = choosePeak + params[:windowType] = windowType + params[:higherHarmonic] = choosePeak + params[:lambda] = lambda + params[:iterations] = iterations # reconstruction parameter + params[:SNRThresh] = 10 + params[:minFreq] = 80e3 + params[:recChannels] = 1:3 + params[:λ] = 0.0f0 + + c = reconstruct("MultiPatchMotion", bMeas; params...) exportImage(joinpath(imgdir, "Motion1.png"), maximum(Array(c[1,:,:,:,1]),dims=3)[:,:,1]) @test compareImg("Motion1.png") - setAll!(plan, :bgParams, SimpleExternalBackgroundCorrectionParameters(;emptyMeas = bBG, bgFrames = bgFrames)) - c = reconstruct(build(plan), bMeas) + c = reconstruct("MultiPatchMotion", bMeas; params..., bgParams = SimpleExternalBackgroundCorrectionParameters(;emptyMeas = bBG, bgFrames = bgFrames)) exportImage(joinpath(imgdir, "Motion2.png"), maximum(Array(c[1,:,:,:,1]),dims=3)[:,:,1]) @test compareImg("Motion2.png") diff --git a/test/MultiGradient.jl b/test/MultiGradient.jl index 4cc325b..9a31280 100644 --- a/test/MultiGradient.jl +++ b/test/MultiGradient.jl @@ -8,18 +8,18 @@ using MPIReco bSFs = MultiMPIFile([joinpath(datadir, "calibrations", "13.mdf"), joinpath(datadir, "calibrations", "14.mdf")]) names = names = (:color, :x, :y, :z, :time) - plan = getPlan("MultiPatch") - setAll!(plan, :SNRThresh, 2) - setAll!(plan, :frames, [1]) - setAll!(plan, :minFreq, 80e3) - setAll!(plan, :recChannels, 1:2) - setAll!(plan, :iterations, 3) - setAll!(plan, :sf, bSFs) - setAll!(plan, :λ, 0.003) - setAll!(plan, :tfCorrection, false) - setAll!(plan, :roundPatches, false) - - c1 = reconstruct(build(plan), b) + params = Dict{Symbol, Any}() + params[:SNRThresh] = 2 + params[:frames] = [1] + params[:minFreq] = 80e3 + params[:recChannels] = 1:2 + params[:iterations] = 3 + params[:sf] = bSFs + params[:λ] = 0.003 + params[:tfCorrection] = false + params[:roundPatches] = false + + c1 = reconstruct("MultiPatch", b; params...) @test axisnames(c1) == names @test axisvalues(c1) == (1:1, -26.0u"mm":1.0u"mm":26.0u"mm", -26.75u"mm":1.0u"mm":26.25u"mm", 0.0u"mm":1.0u"mm":0.0u"mm", 0.0u"ms":0.6528u"ms":0.0u"ms") im1 = reverse(c1[1,:,:,1,1]',dims=1) @@ -31,7 +31,7 @@ using MPIReco dirs = ["1.mdf", "2.mdf", "3.mdf", "4.mdf"] b = MultiMPIFile(joinpath.(datadir, "measurements", "20211226_203921_MultiGradient", dirs)) - c2 = reconstruct(build(plan), b) + c2 = reconstruct("MultiPatch", b; params...) @test axisnames(c2) == names # TODO wo kommt dieser Floatingpoint Fehler rein? Lässt der sich vermeiden? @test axisvalues(c2) == (1:1, -26.0u"mm":1.0u"mm":26.0u"mm", -26.500000000000004u"mm":1.0u"mm":30.499999999999996u"mm", 0.0u"mm":1.0u"mm":0.0u"mm", 0.0u"ms":0.6528u"ms":0.0u"ms") @@ -45,8 +45,8 @@ using MPIReco b = MultiMPIFile([joinpath(datadir, "measurements", "20211226_203921_MultiGradient", "1.mdf")]) bSFs = MultiMPIFile([joinpath(datadir, "calibrations", "13.mdf")]) - setAll!(plan, :sf, bSFs) - c3 = reconstruct(build(plan), b) + params[:sf] = bSFs + c3 = reconstruct("MultiPatch", b; params...) im3 = reverse(c3[1,:,:,1,1]',dims=1) exportImage(joinpath(imgdir, "MultiGradient3.png"), im3) @test compareImg("MultiGradient3.png") diff --git a/test/MultiPatch.jl b/test/MultiPatch.jl index 164db8c..03b6be3 100644 --- a/test/MultiPatch.jl +++ b/test/MultiPatch.jl @@ -14,26 +14,26 @@ using MPIReco -19.375u"mm":1.25u"mm":19.375u"mm", values1[4:5]...) # basic multi-patch reconstruction - plan = getPlan("MultiPatch") - setAll!(plan, :SNRThresh, 5) - setAll!(plan, :frames, [1]) - setAll!(plan, :minFreq, 80e3) - setAll!(plan, :recChannels, 1:2) - setAll!(plan, :iterations, 1) - setAll!(plan, :spectralLeakageCorrection, false) - setAll!(plan, :sf, bSF) - setAll!(plan, :λ, 0) - setAll!(plan, :tfCorrection, false) - c1 = reconstruct(build(plan), b) + params = Dict{Symbol, Any}() + params[:SNRThresh] = 5 + params[:frames] = [1] + params[:minFreq] = 80e3 + params[:recChannels] = 1:2 + params[:iterations] = 1 + params[:spectralLeakageCorrection] = false + params[:sf] = bSF + params[:λ] = 0 + params[:tfCorrection] = false + c1 = reconstruct("MultiPatch", b; params...) @test axisnames(c1) == names @test axisvalues(c1) == values1 exportImage(joinpath(imgdir, "MultiPatch1.png"), Array(c1[1,:,:,1,1])) @test compareImg("MultiPatch1.png") # TODO test description - setAll!(plan, :roundPatches, true) - c2 = reconstruct(build(plan), b) - setAll!(plan, :roundPatches, false) + params[:roundPatches] = true + c2 = reconstruct("MultiPatch", b; params...) + params[:roundPatches] = false @test axisnames(c2) == names @test axisvalues(c2) == values1 exportImage(joinpath(imgdir, "MultiPatch2.png"), Array(c2[1,:,:,1,1])) @@ -42,8 +42,8 @@ using MPIReco # multi-patch reconstruction using multiple system matrices dirs = ["8.mdf", "9.mdf", "10.mdf", "11.mdf"] bSFs = MultiMPIFile(joinpath.(datadir, "calibrations", dirs)) - setAll!(plan, :sf, bSFs) - c3 = reconstruct(build(plan), b) + params[:sf] = bSFs + c3 = reconstruct("MultiPatch", b; params...) @test axisnames(c3) == names @test axisvalues(c3) == values2 exportImage(joinpath(imgdir, "MultiPatch3.png"), Array(c3[1,:,:,1,1])) @@ -57,16 +57,16 @@ using MPIReco S = [getSF(SF,freq,nothing,Kaczmarz, bgcorrection=false)[1] for SF in bSFs] SFGridCenter = zeros(3,4) opParams = ExplicitMultiPatchParameter(;tfCorrection = false, systemMatrices = S, SFGridCenter = SFGridCenter, mapping = mapping) - setAll!(plan, :opParams, opParams) + params[:opParams] = opParams FFPos = zeros(3,4) FFPos[:,1] = [-0.008, 0.008, 0.0] FFPos[:,2] = [-0.008, -0.008, 0.0] FFPos[:,3] = [0.008, 0.008, 0.0] FFPos[:,4] = [0.008, -0.008, 0.0] ffPos = CustomFocusFieldPositions(FFPos) - setAll!(plan, :ffPos, ffPos) - setAll!(plan, :ffPosSF, ffPos) - c4 = reconstruct(build(plan), b) + params[:ffPos] = ffPos + params[:ffPosSF] = ffPos + c4 = reconstruct("MultiPatch", b; params...) @test axisnames(c4) == names @test axisvalues(c4) == values2 exportImage(joinpath(imgdir, "MultiPatch4.png"), Array(c4[1,:,:,1,1])) diff --git a/test/Reconstruction.jl b/test/Reconstruction.jl index a0738ca..b699d2a 100644 --- a/test/Reconstruction.jl +++ b/test/Reconstruction.jl @@ -11,63 +11,62 @@ using MPIReco 0.0u"ms":0.6528u"ms":0.0u"ms") # standard reconstruction - plan = getPlan("Single") - setAll!(plan, :SNRThresh, 5) - setAll!(plan, :frames, 1:1) - setAll!(plan, :minFreq, 80e3), - setAll!(plan, :recChannels, 1:2) - setAll!(plan, :iterations, 1) - setAll!(plan, :spectralLeakageCorrection, true) - setAll!(plan, :sf, bSF) - setAll!(plan, :reg, [L2Regularization(0.0f0)]) - setAll!(plan, :solver, Kaczmarz) - setAll!(plan, :gridding, SystemMatrixGriddingParameter(;gridsize=calibSize(bSF), fov = calibFov(bSF))) + params = Dict{Symbol, Any}() + params[:SNRThresh] = 5 + params[:frames] = 1:1 + params[:minFreq] = 80e3 + params[:recChannels] = 1:2 + params[:iterations] = 1 + params[:spectralLeakageCorrection] = true + params[:sf] = bSF + params[:reg] = [L2Regularization(0.0f0)] + params[:solver] = Kaczmarz - c1 = reconstruct(build(plan), b) + c1 = reconstruct("SinglePatch", b; params...) @test axisnames(c1) == names @test axisvalues(c1) == values exportImage(joinpath(imgdir, "Reconstruction1.png"), Array(c1[1,:,:,1,1])) @test compareImg("Reconstruction1.png") # fused lasso - #setAll!(plan, :iterations, 100) - #setAll!(plan, :loadasreal, true) - #setAll!(plan, :solver, FusedLasso) - #setAll!(plan, :reg, [TVRegularization(0.01f0), L1Regularization(0.01f0)]) - #setAll!(plan, :normalizeReg, NoNormalization()) + #params[:iterations, 100) + #params[:loadasreal, true) + #params[:solver, FusedLasso) + #params[:reg, [TVRegularization(0.01f0), L1Regularization(0.01f0)]) + #params[:normalizeReg, NoNormalization()) #c2 = reconstruct(build(plan), b) - #setAll!(plan, :iterations, 1) - #setAll!(plan, :loadasreal, false) - #setAll!(plan, :solver, Kaczmarz) - #setAll!(plan, :reg, [L2Regularization(0.0f0)]) - #setAll!(plan, :normalizeReg, SystemMatrixBasedNormalization()) + #params[:iterations, 1) + #params[:loadasreal, false) + #params[:solver, Kaczmarz) + #params[:reg, [L2Regularization(0.0f0)]) + #params[:normalizeReg, SystemMatrixBasedNormalization()) #@test axisnames(c2) == names #@test axisvalues(c2) == values #exportImage(joinpath(imgdir, "Reconstruction2.png"), Array(c2[1,:,:,1,1])) #@test compareImg("Reconstruction2.png") # with interpolation - setAll!(plan, :gridding, SystemMatrixGriddingParameter(;gridsize=[100,100,1], fov = calibFov(bSF))) - c3 = reconstruct(build(plan), b) + params[:gridding] = SystemMatrixGriddingParameter(;gridsize=[100,100,1], fov = calibFov(bSF)) + c3 = reconstruct("SinglePatch", b; params...) @test axisnames(c3) == names @test axisvalues(c3) == (values[1], -27.8u"mm":0.4u"mm":11.8u"mm", -11.8u"mm":0.4u"mm":27.8u"mm", values[4:5]...) exportImage(joinpath(imgdir, "Reconstruction3.png"), Array(c3[1,:,:,1,1])) @test compareImg("Reconstruction3.png") # with fov adpation and center shift - setAll!(plan, :gridding, SystemMatrixGriddingParameter(;gridsize=[50,50,1], fov = [0.04,0.04,1], center = [0.0,-0.01,0])) - c4 = reconstruct(build(plan), b) + params[:gridding] = SystemMatrixGriddingParameter(;gridsize=[50,50,1], fov = [0.04,0.04,1], center = [0.0,-0.01,0]) + c4 = reconstruct("SinglePatch", b; params...) @test axisnames(c4) == names # TODO Tobi: does this make sense? @test axisvalues(c4) == (values[1], -27.6u"mm":0.8u"mm":11.6u"mm", -11.6u"mm":0.8u"mm":27.6u"mm", 499.5u"mm":1000.0u"mm":499.5u"mm", values[5]) exportImage(joinpath(imgdir, "Reconstruction4.png"), Array(c4[1,:,:,1,1])) @test compareImg("Reconstruction4.png") - setAll!(plan, :gridding, SystemMatrixGriddingParameter(;gridsize=calibSize(bSF), fov = calibFov(bSF))) + params[:gridding] = SystemMatrixGriddingParameter(;gridsize=calibSize(bSF), fov = calibFov(bSF)) # multi colored reconstruction (new interface) - setAll!(plan, :sf, MultiContrastFile([bSF,bSF])) - c5 = reconstruct(build(plan), b) + params[:sf] = MultiContrastFile([bSF,bSF]) + c5 = reconstruct("SinglePatch", b; params...) @test axisnames(c5) == names @test axisvalues(c5) == (1:2,values[2:end]...) exportImage(joinpath(imgdir, "Reconstruction5a.png"), Array(c5[1,:,:,1,1])) @@ -91,44 +90,44 @@ using MPIReco @test compareImg("Reconstruction6.png") =# - setAll!(plan, :sf, bSF) - setAll!(plan, :reg, [L2Regularization(0.1f0)]) - setAll!(plan, :gridding, SystemMatrixGriddingParameter(;gridsize=calibSize(bSF), fov = calibFov(bSF))) + params[:sf] = bSF + params[:reg] = [L2Regularization(0.1f0)] + params[:gridding] = SystemMatrixGriddingParameter(;gridsize=calibSize(bSF), fov = calibFov(bSF)) # channel weighting - setAll!(plan, :weightingParams, ChannelWeightingParameters(channelWeights = [1.0, 1.0, 1.0])) - c7a = reconstruct(build(plan), b) + params[:weightingParams] = ChannelWeightingParameters(channelWeights = [1.0, 1.0, 1.0]) + c7a = reconstruct("SinglePatch", b; params...) exportImage(joinpath(imgdir, "Reconstruction7a.png"), Array(c7a[1,:,:,1,1])) @test compareImg("Reconstruction7a.png") - setAll!(plan, :weightingParams, ChannelWeightingParameters(channelWeights = [1.0,0.001,1.0])) - c7b = reconstruct(build(plan), b) + params[:weightingParams] = ChannelWeightingParameters(channelWeights = [1.0,0.001,1.0]) + c7b = reconstruct("SinglePatch", b; params...) exportImage(joinpath(imgdir, "Reconstruction7b.png"), Array(c7b[1,:,:,1,1])) @test compareImg("Reconstruction7b.png") - setAll!(plan, :weightingParams, ChannelWeightingParameters(channelWeights = [0.001,1.0,1.0])) - c7c = reconstruct(build(plan), b) + params[:weightingParams] = ChannelWeightingParameters(channelWeights = [0.001,1.0,1.0]) + c7c = reconstruct("SinglePatch", b; params...) exportImage(joinpath(imgdir, "Reconstruction7c.png"), Array(c7c[1,:,:,1,1])) @test compareImg("Reconstruction7c.png") - setAll!(plan, :weightingParams, ChannelWeightingParameters(channelWeights = [1.0, 0.0, 1.0])) - c7d = reconstruct(build(plan), b) - setAll!(plan, :weightingParams, NoWeightingParameters()) - setAll!(plan, :recChannels, 1:1) - c7e = reconstruct(build(plan), b) - setAll!(plan, :recChannels, 1:2) + params[:weightingParams] = ChannelWeightingParameters(channelWeights = [1.0, 0.0, 1.0]) + c7d = reconstruct("SinglePatch", b; params...) + params[:weightingParams] = NoWeightingParameters() + params[:recChannel] = 1:1 + c7e = reconstruct("SinglePatch", b; params...) + params[:recChannel] = 1:2 @test isapprox(arraydata(c7d), arraydata(c7e)) - setAll!(plan, :weightingParams, ChannelWeightingParameters(channelWeights = [0.0, 1.0, 1.0])) - c7f = reconstruct(build(plan), b) - setAll!(plan, :weightingParams, NoWeightingParameters()) - setAll!(plan, :recChannels, 2:2) - c7g = reconstruct(build(plan), b) - setAll!(plan, :recChannels, 1:2) + params[:weightingParams] = ChannelWeightingParameters(channelWeights = [0.0, 1.0, 1.0]) + c7f = reconstruct("SinglePatch", b; params...) + params[:weightingParams] = NoWeightingParameters() + params[:recChannel] = 2:2 + c7g = reconstruct("SinglePatch", b; params...) + params[:recChannel] = 1:2 @test isapprox(arraydata(c7f), arraydata(c7g)) - setAll!(plan, :weightingParams, WhiteningWeightingParameters(whiteningMeas = bSF)) - c8 = reconstruct(build(plan), b) + params[:weightingParams] = WhiteningWeightingParameters(whiteningMeas = bSF) + c8 = reconstruct("SinglePatch", b; params...) exportImage(joinpath(imgdir, "Reconstruction8.png"), Array(c8[1,:,:,1,1])) @test compareImg("Reconstruction8.png") end diff --git a/test/SMExtrapolation.jl b/test/SMExtrapolation.jl index ae35717..987d73f 100644 --- a/test/SMExtrapolation.jl +++ b/test/SMExtrapolation.jl @@ -17,30 +17,28 @@ using MPIReco @test prod(SMsizeExtr) == size(SMextr2,1) b = MPIFile(joinpath(datadir, "measurements", "20211226_203916_MultiPatch", "1.mdf")) - plan = getPlan("Single") - setAll!(plan, :SNRThresh, 5) - setAll!(plan, :frames, 1:1) - setAll!(plan, :minFreq, 80e3), - setAll!(plan, :recChannels, 1:2) - setAll!(plan, :iterations, 1) - setAll!(plan, :spectralLeakageCorrection, true) - setAll!(plan, :sf, bSF) - setAll!(plan, :reg, [L2Regularization(0.0f0)]) - setAll!(plan, :solver, Kaczmarz) - setAll!(plan, :gridsize, calibSize(bSF)) - setAll!(plan, :fov, calibFov(bSF)) - c1 = reconstruct(build(plan), b) + params = Dict{Symbol, Any}() + params[:SNRThresh] = 5 + params[:frames] = 1:1 + params[:minFreq] = 80e3 + params[:recChannels] = 1:2 + params[:iterations] = 1 + params[:spectralLeakageCorrection] = true + params[:sf] = bSF + params[:reg] = [L2Regularization(0.0f0)] + params[:solver] = Kaczmarz + c1 = reconstruct("SinglePatch", b; params...) - setAll!(plan, :gridsize, [36,36,1]) - setAll!(plan, :fov, calibFov(bSF).+[0.006,0.006,0]) - c_extr = reconstruct(build(plan), b) + params[:gridsize] = [36,36,1] + params[:fov] = calibFov(bSF).+[0.006,0.006,0] + c_extr = reconstruct("SinglePatch", b; params...) @test size(c1[1,:,:,:,1]) .+ (4,4,0) == size(c_extr[1,:,:,:,1]) exportImage(joinpath(imgdir, "Extrapolated1.png"), Array(c_extr[1,:,:,1,1])) @test compareImg("Extrapolated1.png") - setAll!(plan, :gridsize, calibSize(bSF)) - setAll!(plan, :fov, calibFov(bSF).+[0.006,0.006,0]) - c_extr2 = reconstruct(build(plan), b) + params[:gridsize] = calibSize(bSF) + params[:fov] = calibFov(bSF).+[0.006,0.006,0] + c_extr2 = reconstruct("SinglePatch", b; params...) @test size(c1[1,:,:,:,1]) == size(c_extr2[1,:,:,:,1]) exportImage(joinpath(imgdir, "Extrapolated2.png"), Array(c_extr2[1,:,:,1,1])) @test compareImg("Extrapolated2.png") @@ -50,27 +48,27 @@ using MPIReco bdirs = ["1.mdf", "2.mdf", "3.mdf", "4.mdf"] b = MultiMPIFile(joinpath.(datadir, "measurements", "20211226_203916_MultiPatch", bdirs)) - plan = getPlan("MultiPatch") - setAll!(plan, :SNRThresh, 5) - setAll!(plan, :frames, 1:1) - setAll!(plan, :minFreq, 80e3) - setAll!(plan, :recChannels, 1:2) - setAll!(plan, :iterations, 1) - setAll!(plan, :spectralLeakageCorrection, false) - setAll!(plan, :sf, bSFs) - setAll!(plan, :λ, 0) - setAll!(plan, :tfCorrection, false) - c2 = reconstruct(build(plan), b) + params = Dict{Symbol, Any}() + params[:SNRThresh] = 5 + params[:frames] = 1:1 + params[:minFreq] = 80e3 + params[:recChannels] = 1:2 + params[:iterations] = 1 + params[:spectralLeakageCorrection] = false + params[:sf] = bSFs + params[:λ] = 0 + params[:tfCorrection] = false + c2 = reconstruct("MultiPatch", b; params...) calibsize=hcat((calibSize.(bSFs)...)).+[6, 6, 0]*ones(Int,4)' fov=hcat((calibFov.(bSFs)...)).+[0.006,0.006,0]*ones(Int,4)' - setAll!(plan, :opParams, RecoPlan(ExplicitMultiPatchParameter)) - setAll!(plan, :tfCorrection, false) - setAll!(plan, :gridsize, calibsize) - setAll!(plan, :fov, fov) - setAll!(plan, :mapping, 1:4) - c_extr3 = reconstruct(build(plan), b) + params[:opParams] = RecoPlan(ExplicitMultiPatchParameter) + params[:tfCorrection] = false + params[:gridsize] = calibsize + params[:fov] = fov + params[:mapping] = 1:4 + c_extr3 = reconstruct("MultiPatch", b; params...) @test size(c2[1,:,:,:,1]) .+ (6,6,0) == size(c_extr3[1,:,:,:,1]) exportImage(joinpath(imgdir, "ExtrapolatedMultiPatch1.png"), Array(c_extr3[1,:,:,1,1])) @test compareImg("ExtrapolatedMultiPatch1.png") diff --git a/test/Solvers.jl b/test/Solvers.jl index d0dd980..6828683 100644 --- a/test/Solvers.jl +++ b/test/Solvers.jl @@ -11,25 +11,24 @@ using MPIReco 0.0u"ms":0.6528u"ms":0.0u"ms") # standard reconstruction - plan = getPlan("Single") - plan.parameter.reco.solverParams = RecoPlan(ElaborateSolverParameters) - setAll!(plan, :SNRThresh, 5) - setAll!(plan, :frames, 1:1) - setAll!(plan, :minFreq, 80e3), - setAll!(plan, :recChannels, 1:2) - setAll!(plan, :spectralLeakageCorrection, true) - setAll!(plan, :sf, bSF) - setAll!(plan, :gridding, SystemMatrixGriddingParameter(;gridsize=calibSize(bSF), fov = calibFov(bSF))) - setAll!(plan, :weightingParams, WhiteningWeightingParameters(whiteningMeas = bSF)) + params = Dict{Symbol, Any}() + params[:SNRThresh] = 5 + params[:frames] = 1:1 + params[:minFreq] = 80e3 + params[:recChannels] = 1:2 + params[:spectralLeakageCorrection] = true + params[:sf] = bSF + params[:gridding] = SystemMatrixGriddingParameter(;gridsize=calibSize(bSF), fov = calibFov(bSF)) + params[:weightingParams] = WhiteningWeightingParameters(whiteningMeas = bSF) - setAll!(plan, :reg, [L2Regularization(0.1f0)]) - setAll!(plan, :iterations, 100) + params[:reg] = [L2Regularization(0.1f0)] + params[:iterations] = 100 for solver in [Kaczmarz, CGNR, FISTA, OptISTA, POGM] # , ADMM, SplitBregman] - setAll!(plan, :solver, solver) - setAll!(plan, :rho, 0.3f0) - c = reconstruct(build(plan), b) + params[:solver] = solver + params[:rho] = 0.3f0 + c = reconstruct("SinglePatch", b; params...) @test axisnames(c) == names @test axisvalues(c) == values exportImage(joinpath(imgdir, "Solver_$solver.png"), Array(c[1,:,:,1,1])) diff --git a/test/Sparse.jl b/test/Sparse.jl index 1e985e9..d22e58e 100644 --- a/test/Sparse.jl +++ b/test/Sparse.jl @@ -16,87 +16,72 @@ using MPIReco -0.5u"mm":1.0u"mm":-0.5u"mm", 0.0u"ms":65.28u"ms":0.0u"ms") - plan = getPlan("Single") - setAll!(plan, :SNRThresh, 2) - setAll!(plan, :frames, 1:100) - setAll!(plan, :numAverages, 100) - setAll!(plan, :minFreq, 80e3), - setAll!(plan, :recChannels, 1:2) - setAll!(plan, :iterations, 3) - setAll!(plan, :spectralLeakageCorrection, false) - setAll!(plan, :sf, bSF) - setAll!(plan, :reg, [L2Regularization(0.1f0)]) - setAll!(plan, :solver, Kaczmarz) - setAll!(plan, :gridding, SystemMatrixGriddingParameter(;gridsize=calibSize(bSF), fov = calibFov(bSF))) - c1 = reconstruct(build(plan), b) + params = Dict{Symbol, Any}() + params[:SNRThresh] = 2 + params[:frames] = 1:100 + params[:numAverages] = 100 + params[:minFreq] = 80e3 + params[:recChannels] = 1:2 + params[:iterations] = 3 + params[:spectralLeakageCorrection] = false + params[:sf] = bSF + params[:reg] = [L2Regularization(0.1f0)] + params[:solver] = Kaczmarz + params[:gridding] = SystemMatrixGriddingParameter(;gridsize=calibSize(bSF), fov = calibFov(bSF)) + c1 = reconstruct("SinglePatch", b; params...) @test axisnames(c1) == names @test axisvalues(c1) == values exportImage(joinpath(imgdir, "Sparse1.png"), Array(c1[1,:,:,1,1])) @test compareImg("Sparse1.png") - plan = getPlan("Sparse") - setAll!(plan, :SNRThresh, 2) - setAll!(plan, :frames, 1:100) - setAll!(plan, :minFreq, 80e3) - setAll!(plan, :recChannels, 1:2) - setAll!(plan, :iterations, 3) - setAll!(plan, :spectralLeakageCorrection, false) - setAll!(plan, :sf, bSF) - setAll!(plan, :reg, [L2Regularization(0.1f0)]) - setAll!(plan, :numAverages, 100) - setAll!(plan, :solver, Kaczmarz) - #setAll!(plan, :tfCorrection, false) - setAll!(plan, :redFactor, redFactor) - - setAll!(plan, :sparseTrafo, "FFT") - setAll!(plan, :useDFFoV, false) - c2 = reconstruct(build(plan), b) + params = Dict{Symbol, Any}() + params[:SNRThresh] = 2 + params[:frames] = 1:100 + params[:minFreq] = 80e3 + params[:recChannels] = 1:2 + params[:iterations] = 3 + params[:spectralLeakageCorrection] = false + params[:sf] = bSF + params[:reg] = [L2Regularization(0.1f0)] + params[:numAverages] = 100 + params[:solver] = Kaczmarz + #params[ :tfCorrection, false + params[:redFactor] = redFactor + + c2 = reconstruct("SinglePatchSparse", b; params..., sparseTrafo = "FFT", useDFFoV = false) @test axisnames(c2) == names @test axisvalues(c2) == values exportImage(joinpath(imgdir, "Sparse2.png"), Array(c2[1,:,:,1,1])) @test compareImg("Sparse2.png") - setAll!(plan, :sparseTrafo, "FFT") - setAll!(plan, :useDFFoV, true) - c3 = reconstruct(build(plan), b) + c3 = reconstruct("SinglePatchSparse", b; params..., sparseTrafo = "FFT", useDFFoV = true) @test axisnames(c3) == names @test axisvalues(c3) == valuesDF exportImage(joinpath(imgdir, "Sparse3.png"), Array(c3[1,:,:,1,1])) @test compareImg("Sparse3.png") - setAll!(plan, :SNRThresh, 3) - setAll!(plan, :iterations, 1) - setAll!(plan, :reg, [L2Regularization(0.01f0)]) - setAll!(plan, :sparseTrafo, "DCT-IV") - setAll!(plan, :useDFFoV, false) - c4 = reconstruct(build(plan), b) + params[:SNRThresh] = 3 + params[:iterations] = 1 + params[:reg] = [L2Regularization(0.01f0)] + c4 = reconstruct("SinglePatchSparse", b; params..., sparseTrafo = "DCT-IV", useDFFoV = false) @test axisnames(c4) == names @test axisvalues(c4) == values exportImage(joinpath(imgdir, "Sparse4.png"), Array(c4[1,:,:,1,1])) @test compareImg("Sparse4.png") - setAll!(plan, :spectralLeakageCorrection, true) - setAll!(plan, :sparseTrafo, "DCT-IV") - setAll!(plan, :useDFFoV, true) - c5 = reconstruct(build(plan), b) - setAll!(plan, :spectralLeakageCorrection, false) + c5 = reconstruct("SinglePatchSparse", b; params..., sparseTrafo = "DCT-IV", useDFFoV = true, spectralLeakageCorrection = true) @test axisnames(c5) == names @test axisvalues(c5) == valuesDF exportImage(joinpath(imgdir, "Sparse5.png"), Array(c5[1,:,:,1,1])) @test compareImg("Sparse5.png") - setAll!(plan, :sparseTrafo, "DST") - setAll!(plan, :useDFFoV, false) - c6 = reconstruct(build(plan), b) + c6 = reconstruct("SinglePatchSparse", b; params..., sparseTrafo = "DST", useDFFoV = false) @test axisnames(c6) == names @test axisvalues(c6) == values exportImage(joinpath(imgdir, "Sparse6.png"), Array(c6[1,:,:,1,1])) @test compareImg("Sparse6.png") - setAll!(plan, :spectralLeakageCorrection, true) - setAll!(plan, :sparseTrafo, "DST") - setAll!(plan, :useDFFoV, true) - c7 = reconstruct(build(plan), b) + c7 = reconstruct("SinglePatchSparse", b; params..., sparseTrafo = "DST", useDFFoV = true, spectralLeakageCorrection = true) @test axisnames(c7) == names @test axisvalues(c7) == valuesDF exportImage(joinpath(imgdir, "Sparse7.png"), Array(c7[1,:,:,1,1])) diff --git a/test/recoPlans/Motion.toml b/test/recoPlans/Motion.toml deleted file mode 100644 index 9d6aadd..0000000 --- a/test/recoPlans/Motion.toml +++ /dev/null @@ -1,26 +0,0 @@ -_type = "RecoPlan{MultiPatchReconstructionAlgorithm}" -_module = "MPIReco" - -[parameter] -_type = "RecoPlan{MultiPatchParameters}" -_module = "MPIReco" - - [parameter.reco] - _type = "RecoPlan{PeriodicMotionReconstructionParameter}" - _module = "MPIReco" - - [parameter.reco.freqFilter] - _type = "RecoPlan{SNRThresholdFrequencyFilterParameter}" - _module = "MPIReco" - - [parameter.reco.solverParams] - _type = "RecoPlan{SimpleSolverParameters}" - _module = "MPIReco" - - [parameter.pre] - _type = "RecoPlan{PeriodicMotionPreProcessing}" - _module = "MPIReco" - - [parameter.post] - _type = "RecoPlan{NoPostProcessing}" - _module = "MPIReco" diff --git a/test/recoPlans/MultiPatch.toml b/test/recoPlans/MultiPatch.toml deleted file mode 100644 index 2d12da3..0000000 --- a/test/recoPlans/MultiPatch.toml +++ /dev/null @@ -1,42 +0,0 @@ -_type = "RecoPlan{MultiPatchReconstructionAlgorithm}" -_module = "MPIReco" - -[parameter] -_type = "RecoPlan{MultiPatchParameters}" -_module = "MPIReco" - - [parameter.reco] - _type = "RecoPlan{MultiPatchReconstructionParameter}" - _module = "MPIReco" - - [parameter.reco.ffPos] - _type = "RecoPlan{DefaultFocusFieldPositions}" - _module = "MPIReco" - - [parameter.reco.freqFilter] - _type = "RecoPlan{SNRThresholdFrequencyFilterParameter}" - _module = "MPIReco" - - [parameter.reco.opParams] - _type = "RecoPlan{RegularMultiPatchOperatorParameter}" - _module = "MPIReco" - - [parameter.reco.ffPosSF] - _type = "RecoPlan{DefaultFocusFieldPositions}" - _module = "MPIReco" - - [parameter.reco.solverParams] - _type = "RecoPlan{SimpleSolverParameters}" - _module = "MPIReco" - - [parameter.pre] - _type = "RecoPlan{CommonPreProcessingParameters}" - _module = "MPIReco" - - [parameter.pre.bgParams] - _type = "RecoPlan{NoBackgroundCorrectionParameters}" - _module = "MPIReco" - - [parameter.post] - _type = "RecoPlan{NoPostProcessing}" - _module = "MPIReco" diff --git a/test/recoPlans/Single.toml b/test/recoPlans/Single.toml deleted file mode 100644 index 6a16ee4..0000000 --- a/test/recoPlans/Single.toml +++ /dev/null @@ -1,34 +0,0 @@ -_type = "RecoPlan{SinglePatchReconstructionAlgorithm}" -_module = "MPIReco" - -[parameter] -_type = "RecoPlan{SinglePatchParameters}" -_module = "MPIReco" - - [parameter.reco] - _type = "RecoPlan{SinglePatchReconstructionParameter}" - _module = "MPIReco" - - [parameter.reco.sfLoad] - _type = "RecoPlan{DenseSystemMatixLoadingParameter}" - _module = "MPIReco" - - [parameter.reco.sfLoad.freqFilter] - _type = "RecoPlan{SNRThresholdFrequencyFilterParameter}" - _module = "MPIReco" - - [parameter.reco.sfLoad.gridding] - _type = "RecoPlan{SystemMatrixGriddingParameter}" - _module = "MPIReco" - - [parameter.reco.solverParams] - _type = "RecoPlan{SimpleSolverParameters}" - _module = "MPIReco" - - [parameter.pre] - _type = "RecoPlan{CommonPreProcessingParameters}" - _module = "MPIReco" - - [parameter.post] - _type = "RecoPlan{NoPostProcessing}" - _module = "MPIReco" diff --git a/test/recoPlans/Sparse.toml b/test/recoPlans/Sparse.toml deleted file mode 100644 index 0db0cad..0000000 --- a/test/recoPlans/Sparse.toml +++ /dev/null @@ -1,34 +0,0 @@ -_type = "RecoPlan{SinglePatchReconstructionAlgorithm}" -_module = "MPIReco" - -[parameter] -_type = "RecoPlan{SinglePatchParameters}" -_module = "MPIReco" - - [parameter.reco] - _type = "RecoPlan{SinglePatchReconstructionParameter}" - _module = "MPIReco" - - [parameter.reco.sfLoad] - _type = "RecoPlan{SparseSystemMatrixLoadingParameter}" - _module = "MPIReco" - - [parameter.reco.sfLoad.freqFilter] - _type = "RecoPlan{SNRThresholdFrequencyFilterParameter}" - _module = "MPIReco" - - [parameter.reco.solverParams] - _type = "RecoPlan{SimpleSolverParameters}" - _module = "MPIReco" - - [parameter.pre] - _type = "RecoPlan{CommonPreProcessingParameters}" - _module = "MPIReco" - - [parameter.pre.bgParams] - _type = "RecoPlan{NoBackgroundCorrectionParameters}" - _module = "MPIReco" - - [parameter.post] - _type = "RecoPlan{NoPostProcessing}" - _module = "MPIReco" diff --git a/test/runtests.jl b/test/runtests.jl index 0d49bf9..92df5f8 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -26,8 +26,6 @@ function compareSSIM(expected, given; limit::Float64=SSIM_LIMIT, kwargs...) return ssim >= limit end -getPlan(plan::String) = loadPlan(joinpath("recoPlans", "$(plan).toml"), [MPIReco, RegularizedLeastSquares, MPIFiles]) - function exportImage(filename, I::AbstractMatrix) Iabs = abs.(I) Icolored = colorview(Gray, Iabs./maximum(Iabs)) From 265afe83015ccef73f7e96c828b353eddeff5f13 Mon Sep 17 00:00:00 2001 From: nHackel Date: Thu, 11 Jul 2024 16:43:23 +0200 Subject: [PATCH 5/9] Fix cached sparse matrix loading --- src/Algorithms/SinglePatchAlgorithms/SinglePatchAlgorithm.jl | 2 +- src/SystemMatrix/SystemMatrix.jl | 3 +++ 2 files changed, 4 insertions(+), 1 deletion(-) diff --git a/src/Algorithms/SinglePatchAlgorithms/SinglePatchAlgorithm.jl b/src/Algorithms/SinglePatchAlgorithms/SinglePatchAlgorithm.jl index 9928070..8f7c192 100644 --- a/src/Algorithms/SinglePatchAlgorithms/SinglePatchAlgorithm.jl +++ b/src/Algorithms/SinglePatchAlgorithms/SinglePatchAlgorithm.jl @@ -68,5 +68,5 @@ function getLinearOperator(algo::SinglePatchReconstructionAlgorithm, params::Sin end function getLinearOperator(algo::SinglePatchReconstructionAlgorithm, params::SinglePatchReconstructionParameter{<:SparseSystemMatrixLoadingParameter, S}) where {S} - return createLinearOperator(params.sfLoad.sparseTrafo, eltype(algo.S); shape=tuple(shape(algo.grid)...)) + return process(algo, params.sfLoad, eltype(algo.S), tuple(shape(algo.grid)...)) end \ No newline at end of file diff --git a/src/SystemMatrix/SystemMatrix.jl b/src/SystemMatrix/SystemMatrix.jl index 7e9f061..1b9375a 100644 --- a/src/SystemMatrix/SystemMatrix.jl +++ b/src/SystemMatrix/SystemMatrix.jl @@ -92,6 +92,9 @@ function process(t::Type{<:AbstractMPIRecoAlgorithm}, params::SparseSystemMatrix S, grid = getSF(sf, frequencies, params.sparseTrafo; toKwargs(params)...) return S, grid end +function process(t::Type{<:AbstractMPIRecoAlgorithm}, params::SparseSystemMatrixLoadingParameter, elType::Type{<:Number}, shape::NTuple{N, Int64}) where N + return createLinearOperator(params.sparseTrafo, elType; shape) +end function converttoreal(S::AbstractArray{Complex{T}},f) where T N = prod(calibSize(f)) From 89ae525c507be6381a760922a2bb602964bea3cd Mon Sep 17 00:00:00 2001 From: nHackel Date: Thu, 11 Jul 2024 16:46:07 +0200 Subject: [PATCH 6/9] Incr. dependencies --- Project.toml | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/Project.toml b/Project.toml index cc8b4d4..725503c 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "MPIReco" uuid = "e4246700-6248-511e-8146-a1d1f47669d2" authors = ["Tobias Knopp "] -version = "0.5.4" +version = "0.6.0" [deps] DSP = "717857b8-e6f2-59f4-9121-6e50c889abd2" @@ -32,12 +32,12 @@ FFTW = "1.3" ImageUtils = "0.2" IniFile = "0.5" LinearAlgebra = "1" -LinearOperators = "2.3.3" -LinearOperatorCollection = "1.2" +LinearOperators = "2.3" +LinearOperatorCollection = "2" MPIFiles = "0.13, 0.14, 0.15, 0.16" ProgressMeter = "1.2" Reexport = "1.0" -RegularizedLeastSquares = "0.14" +RegularizedLeastSquares = "0.16" SparseArrays = "1" Statistics = "1" ThreadPools = "2.1.1" From faeaf3e2dbf28e4ddb8eb8e42dec380fd1b4fcd9 Mon Sep 17 00:00:00 2001 From: nHackel Date: Thu, 11 Jul 2024 16:56:50 +0200 Subject: [PATCH 7/9] Add GPU tests --- Project.toml | 4 +++- test/Reconstruction.jl | 2 -- test/ReconstructionGPU.jl | 34 ++++++++++++++++++++++++++++++++++ test/gpu/cuda.jl | 5 +++++ test/gpu/rocm.jl | 5 +++++ test/runtests.jl | 33 +++++++++++++++++++++------------ 6 files changed, 68 insertions(+), 15 deletions(-) create mode 100644 test/ReconstructionGPU.jl create mode 100644 test/gpu/cuda.jl create mode 100644 test/gpu/rocm.jl diff --git a/Project.toml b/Project.toml index 725503c..b4ded33 100644 --- a/Project.toml +++ b/Project.toml @@ -31,6 +31,7 @@ DistributedArrays = "0.6" FFTW = "1.3" ImageUtils = "0.2" IniFile = "0.5" +JLArrays = "0.1" LinearAlgebra = "1" LinearOperators = "2.3" LinearOperatorCollection = "2" @@ -47,6 +48,7 @@ julia = "1.9" [extras] FileIO = "5789e2e9-d7fb-5bc7-8068-2c6fae9b9549" HTTP = "cd3eb016-35fb-5094-929b-558a96fad6f3" +JLArrays = "27aeb0d3-9eb9-45fb-866b-73c2ecf80fcb" LazyArtifacts = "4af54fe1-eca0-43a8-85a7-787d91b784e3" Scratch = "6c6a2e73-6563-6170-7368-637461726353" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" @@ -55,4 +57,4 @@ ImageQualityIndexes = "2996bd0c-7a13-11e9-2da2-2f5ce47296a9" Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" [targets] -test = ["Test", "HTTP", "FileIO", "LazyArtifacts", "Scratch", "ImageMagick", "ImageQualityIndexes", "Unitful"] +test = ["Test", "HTTP", "FileIO", "LazyArtifacts", "Scratch", "ImageMagick", "ImageQualityIndexes", "Unitful", "JLArrays"] diff --git a/test/Reconstruction.jl b/test/Reconstruction.jl index b699d2a..b0b1a84 100644 --- a/test/Reconstruction.jl +++ b/test/Reconstruction.jl @@ -1,5 +1,3 @@ -using MPIReco - @testset "single- and multi-channel in-memory reconstruction" begin bSF = MPIFile(joinpath(datadir, "calibrations", "12.mdf")) b = MPIFile(joinpath(datadir, "measurements", "20211226_203916_MultiPatch", "1.mdf")) diff --git a/test/ReconstructionGPU.jl b/test/ReconstructionGPU.jl new file mode 100644 index 0000000..489a0b4 --- /dev/null +++ b/test/ReconstructionGPU.jl @@ -0,0 +1,34 @@ +for arrayType in arrayTypes + @testset "Reconstruction $arrayType" begin + bSF = MPIFile(joinpath(datadir, "calibrations", "12.mdf")) + b = MPIFile(joinpath(datadir, "measurements", "20211226_203916_MultiPatch", "1.mdf")) + names = (:color, :x, :y, :z, :time) + values = (1:1, + -27.375u"mm":1.25u"mm":11.375u"mm", + -11.375u"mm":1.25u"mm":27.375u"mm", + 0.0u"mm":1.0u"mm":0.0u"mm", + 0.0u"ms":0.6528u"ms":0.0u"ms") + + # standard reconstruction + params = Dict{Symbol, Any}() + params[:SNRThresh] = 5 + params[:frames] = 1:1 + params[:minFreq] = 80e3 + params[:recChannels] = 1:2 + params[:iterations] = 1 + params[:spectralLeakageCorrection] = true + params[:sf] = bSF + params[:reg] = [L2Regularization(0.0f0)] + params[:solver] = CGNR + + c1 = reconstruct("SinglePatch", b; params...) + @test axisnames(c1) == names + @test axisvalues(c1) == values + + c2 = reconstruct("SinglePatch", b; params..., arrayType = arrayType) + @test axisnames(c2) == names + @test axisvalues(c2) == values + + @test isapprox(arraydata(c1), arraydata(c2)) + end +end \ No newline at end of file diff --git a/test/gpu/cuda.jl b/test/gpu/cuda.jl new file mode 100644 index 0000000..e7d0eb2 --- /dev/null +++ b/test/gpu/cuda.jl @@ -0,0 +1,5 @@ +using CUDA + +arrayTypes = [CuArray] + +include(joinpath(@__DIR__(), "..", "runtests.jl")) \ No newline at end of file diff --git a/test/gpu/rocm.jl b/test/gpu/rocm.jl new file mode 100644 index 0000000..ebf32fb --- /dev/null +++ b/test/gpu/rocm.jl @@ -0,0 +1,5 @@ +using AMDGPU + +arrayTypes = [ROCArray] + +include(joinpath(@__DIR__(), "..", "runtests.jl")) \ No newline at end of file diff --git a/test/runtests.jl b/test/runtests.jl index 92df5f8..5a8c4de 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,3 +1,4 @@ +using MPIReco using AbstractImageReconstruction using Test using Scratch @@ -32,18 +33,26 @@ function exportImage(filename, I::AbstractMatrix) save(filename, Icolored ) end +areTypesDefined = @isdefined arrayTypes +arrayTypes = areTypesDefined ? arrayTypes : [JLArray] + @testset "MPIReco" begin - #include("LoadSaveMDF.jl") - include("Reconstruction.jl") # FussedLasso causes segfault atm - include("Solvers.jl") - include("Cartesian.jl") - if !Sys.iswindows() - include("MotionCompensation.jl") + if !areTypesDefined + #include("LoadSaveMDF.jl") + include("Reconstruction.jl") # FussedLasso causes segfault atm + include("Solvers.jl") + include("Cartesian.jl") + if !Sys.iswindows() + include("MotionCompensation.jl") + end + include("MultiPatch.jl") + include("MultiGradient.jl") + include("Sparse.jl") + include("SMExtrapolation.jl") + include("ReconstructionGPU.jl") + #include("SMCenter.jl") + #include("SMRecovery.jl") # will be moved to SMTools + else + include("ReconstructionGPU.jl") end - include("MultiPatch.jl") - include("MultiGradient.jl") - include("Sparse.jl") - include("SMExtrapolation.jl") - #include("SMCenter.jl") - #include("SMRecovery.jl") # will be moved to SMTools end From 7bb5c2ce56b42fa8ff49b88b0211ed0fbc419409 Mon Sep 17 00:00:00 2001 From: nHackel Date: Thu, 11 Jul 2024 16:58:14 +0200 Subject: [PATCH 8/9] Add buildkite pipeline for CUDA and ROCm --- .buildkite/pipeline.yml | 37 +++++++++++++++++++++++++++++++++++++ 1 file changed, 37 insertions(+) create mode 100644 .buildkite/pipeline.yml diff --git a/.buildkite/pipeline.yml b/.buildkite/pipeline.yml new file mode 100644 index 0000000..44c9a3f --- /dev/null +++ b/.buildkite/pipeline.yml @@ -0,0 +1,37 @@ +steps: + - label: "Nvidia GPUs -- MPIReco.jl" + plugins: + - JuliaCI/julia#v1: + version: "1.10" + agents: + queue: "juliagpu" + cuda: "*" + command: | + julia --color=yes --project -e ' + using Pkg + Pkg.add("TestEnv") + using TestEnv + TestEnv.activate(); + Pkg.add("CUDA") + Pkg.instantiate() + include("test/gpu/cuda.jl")' + timeout_in_minutes: 30 + + - label: "AMD GPUs -- MPIReco.jl" + plugins: + - JuliaCI/julia#v1: + version: "1.10" + agents: + queue: "juliagpu" + rocm: "*" + rocmgpu: "*" + command: | + julia --color=yes --project -e ' + using Pkg + Pkg.add("TestEnv") + using TestEnv + TestEnv.activate(); + Pkg.add("AMDGPU") + Pkg.instantiate() + include("test/gpu/rocm.jl")' + timeout_in_minutes: 30 \ No newline at end of file From 13b977a524a12d2bfc0fce543637cefc90daf3ad Mon Sep 17 00:00:00 2001 From: Niklas Hackelberg Date: Fri, 12 Jul 2024 08:47:02 +0200 Subject: [PATCH 9/9] Update runtests.jl --- test/runtests.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/test/runtests.jl b/test/runtests.jl index 5a8c4de..e587fe0 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,5 +1,6 @@ using MPIReco using AbstractImageReconstruction +using JLArrays using Test using Scratch using ImageMagick