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

Reconstruction for 3D sequences #323

Open
katyabrui opened this issue Mar 12, 2024 · 14 comments · May be fixed by #324
Open

Reconstruction for 3D sequences #323

katyabrui opened this issue Mar 12, 2024 · 14 comments · May be fixed by #324
Labels
bug Something isn't working feature

Comments

@katyabrui
Copy link

What happened?

When using a 3D pulse sequence (MPRAGE) with brain_phantom3D, the recostruction is not performed propperly. Looks like a 3D k-space is considered as a set of 2D k-spaces or something like this. I attach the .seq file as an example (actually, I couldn't attach it as it is, so you will need to remove "txt" after you save it) MPRAGE_brain_25cmfov.seq.txt
mprage

Environment

OS x86_64-w64-mingw32
Julia 1.9.3
KomaMRIPlots 0.8.0
KomaMRIFiles 0.8.0
KomaMRI 0.8.0
KomaMRICore 0.8.0
KomaMRIBase 0.8.1
@katyabrui katyabrui added the bug Something isn't working label Mar 12, 2024
@curtcorum
Copy link
Contributor

curtcorum commented Mar 12, 2024

@katyabrui

Great that you are using KomaMRI, and thanks for reporting the issue!

An immediate work around might be to use:

your_obj = brain_phantom3D(; ss=your_ss, start_end=[177, 183]);

if using KomaUI:

obj_ui[] = brain_phantom3D(; ss=your_ss, start_end=[177, 183]);

This is a 3D phantom that is very thin along Z, will take a bit longer to simulate, but less so than the default brain_phantom3D(). ss defaults to 4 for both 2D and 3D.

Hope this helps while the issue is investigated!

=====================

This unfortunately does not fix the issue.

@cncastillo cncastillo changed the title Reconstruction for 2D sequences Reconstruction for 3D sequences Mar 12, 2024
@curtcorum
Copy link
Contributor

@katyabrui

I am trying the seq now. Nice job getting it through github!

@cncastillo
Copy link
Member

cncastillo commented Mar 12, 2024

This is related to #308, but instead of the "echo dimension," the problem is in the "slice dimension." Right now, we are guessing how to reconstruct stuff (the dimensions Nx, Ny, Nz); a proper way would be to specify the label counters in the "[EXTENSION]" section of the pulseq file. I believe previously, there was a warning that told the user we were doing this. We should put it back (@beorostica). For now, for the echoes, we could add a seq.DEF["Necho"] to reconstruct interleaved echoes.

In the meantime, you could specify the correct sequence dimensions in the REPL before reconstructing.

seq_ui[].DEF["Nx"] = 32
seq_ui[].DEF["Ny"] = 32
seq_ui[].DEF["Nz"] = 32

I can not test it right now, but I think it should work. For convenience, this also can be specified in the sequence file

[DEFINITIONS]
AdcRasterTime 6.25e-05 
BlockDurationRaster 1e-05 
GradientRasterTime 1e-05 
RadiofrequencyRasterTime 1e-06 
Nx 32
Ny 32
Nz 32

@curtcorum
Copy link
Contributor

@cncastillo @katyabrui

This seems to be the issue, the current seq_ui[].DEF below has Nz=1:

julia> seq_ui[].DEF
Dict{String, Any} with 12 entries:
  "signature"                => "9adc63c8e1f8331d67b3fb1c9bdeac98"
  "AdcRasterTime"            => 6.25e-5
  "GradientRasterTime"       => 1.0e-5
  "Nz"                       => 1
  "Nx"                       => 32
  "Ny"                       => 1024
  "PulseqVersion"            => 1004000
  "BlockDurationRaster"      => 1.0e-5
  "extension"                => [0, 0, 0, 0, 0, 0, 0, 0, 0, 0  …  0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
  "FileName"                 => "MPRAGE_brain_25cmfov.seq"
  "RadiofrequencyRasterTime" => 1.0e-6
  "additional_text"          => "# Format of extension lists:\n# id type ref next_id\n# next_id of 0 terminates the list\n# Extension list is followed by extension specification…

The Nz, Nx, Ny parts of DEF are generated by KomaMRI(?) and are not in the current seq file:

[DEFINITIONS]
AdcRasterTime 6.25e-05 
BlockDurationRaster 1e-05 
GradientRasterTime 1e-05 
RadiofrequencyRasterTime 1e-06 

I saw a similar issue with:

https://github.com/pulseq/pulseq/blob/master/matlab/demoSeq/writeGradientEcho3D.m

but had not gotten around to resolving yet. Will try the above fix, thanks!

@curtcorum
Copy link
Contributor

curtcorum commented Mar 12, 2024

@cncastillo @katyabrui

Changing the Nx-z parameters in seq.DEF:

julia> seq.DEF
Dict{String, Any} with 12 entries:
  "signature"                => "9adc63c8e1f8331d67b3fb1c9bdeac98"
  "AdcRasterTime"            => 6.25e-5
  "GradientRasterTime"       => 1.0e-5
  "Nz"                       => 32
  "Nx"                       => 32
  "Ny"                       => 32
  "PulseqVersion"            => 1004000
  "BlockDurationRaster"      => 1.0e-5
  "extension"                => [0, 0, 0, 0, 0, 0, 0, 0, 0, 0  …  0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
  "FileName"                 => "MPRAGE_brain_25cmfov.seq"
  "RadiofrequencyRasterTime" => 1.0e-6
  "additional_text"          => "# Format of extension lists:\n# id type ref next_id\n# next_id of 0 terminates the list\n# Extension list is followed by extension specifications\n[EXTENSIONS]\n\n#…

And running the simulation, results in:

julia> raw = simulate(obj, seq, sys);
┌ Info: Running simulation in the GPU (Quadro RTX 8000)
│   koma_version = v"0.8.1"
│   sim_method = Bloch()
│   spins = 78126
│   time_points = 81897
└   adc_points = 32768
Progress: 100%|███████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| Time: 0:01:04
  simulated_blocks:  2264
  rf_blocks:         1124
  acq_samples:       32768
 65.210778 seconds (65.12 M allocations: 3.968 GiB, 2.15% gc time)

which seems to be simulating everything:

MPRAGE_brain_raw

however:

julia> raw.params
Dict{String, Any} with 23 entries:
  "protocolName"                   => "NoName"
  "institutionName"                => "Pontificia Universidad Catolica de Chile"
  "reconSize"                      => [32, 32, 1]
  "enc_lim_repetition"             => Limit(0, 0, 0)
  "enc_lim_set"                    => Limit(0, 0, 0)
  "enc_lim_segment"                => Limit(0, 0, 0)
  "userParameters"                 => Dict{String, Any}("Nblocks"=>2264, "gpu"=>1, "gpu_device"=>0, "type_sim_parts"=>Bool[0, 1, 0, 1, 0, 1, 0, 1, 0, 1  …  1, 0, 1, 0, 1, 0, 1, 0, 1, 0], "precision…
  "enc_lim_phase"                  => Limit(0, 0, 0)
  "enc_lim_average"                => Limit(0, 0, 0)
  "enc_lim_slice"                  => Limit(0, 0, 0)
  "reconFOV"                       => Float32[258.065, 258.065, 1.0]
  "systemModel"                    => "v0.8.1"
  "H1resonanceFrequency_Hz"        => 63866203
  "patientName"                    => "brain3D"
  "enc_lim_contrast"               => Limit(0, 0, 0)
  "trajectory"                     => "other"
  "systemFieldStrength_T"          => 1.5
  "enc_lim_kspace_encoding_step_0" => Limit(0, 31, 16)
  "enc_lim_kspace_encoding_step_1" => Limit(0, 31, 16)
  "encodedFOV"                     => Float32[258.065, 258.065, 1.0]
  "enc_lim_kspace_encoding_step_2" => Limit(0, 0, 0)
  "systemVendor"                   => "KomaMRI.jl"
  "encodedSize"                    => [32, 32, 1]

with "encodedSize" => [32, 32, 1] it is not taking the values from seq.DEF as far as I can tell. and the acq parameters get generated correspondingly:


julia> acq = AcquisitionData( raw);
julia> acq.encodingSize
(32, 32)
julia> acq.fov
(258.06451416015625, 258.0649108886719, 1.0)

Trying to manually set the recon size fails with:

julia> Nx, Ny, Nz = [32, 32, 32];

julia> reconParams = Dict{Symbol,Any}(:reco=>"direct", :reconSize=>(Nx, Ny, Nz));

julia> image = reconstruction(acq, reconParams);
ERROR: ArgumentError: Nodes x have dimension 2 != 3

Next step is to fill in the encoded size correctly in raw.

@cncastillo
Copy link
Member

cncastillo commented Mar 12, 2024

with "encodedSize" => [32, 32, 1] it is not taking the values from seq.DEF as far as I can tell. and the acq parameters get generated correspondingly:

Ok, this should work. So it is an unrelated bug (@beorostica can you create an issue?).

julia> image = reconstruction(acq, reconParams);
ERROR: ArgumentError: Nodes x have dimension 2 != 3

Last time I tried using 3D k-spaces (nodes) with MRIReco it always gave an error, so we are just storing 2D k-spaces in the raw data (or at least that is the default ndims=2 of signal_to_raw_data). This could create this dimension mismatch. Maybe we just need to check if :reconSize has 3 components to save a 3D k-space, otherwise a 2D k-space. This is the main thing we need to explore to fix this issue.

@curtcorum
Copy link
Contributor

Manually editing raw.params:

julia> raw.params
Dict{String, Any} with 23 entries:
  "protocolName"                   => "NoName"
  "institutionName"                => "Pontificia Universidad Catolica de Chile"
  "reconSize"                      => [32, 32, 32]
  "enc_lim_repetition"             => Limit(0, 0, 0)
  "enc_lim_set"                    => Limit(0, 0, 0)
  "enc_lim_segment"                => Limit(0, 0, 0)
  "userParameters"                 => Dict{String, Any}("Nblocks"=>2264, "gpu"=>1, "gpu_device"=>0, "type_sim_parts"=>Bool[0, 1, 0, 1, 0, 1, 0, 1, 0, 1  …  1, 0, 1, 0, 1, 0, 1, 0, 1, 0], "precision…
  "enc_lim_phase"                  => Limit(0, 0, 0)
  "enc_lim_average"                => Limit(0, 0, 0)
  "enc_lim_slice"                  => Limit(0, 0, 0)
  "reconFOV"                       => [258.065, 258.065, 258.065]
  "systemModel"                    => "v0.8.1"
  "H1resonanceFrequency_Hz"        => 63866203
  "patientName"                    => "brain3D"
  "enc_lim_contrast"               => Limit(0, 0, 0)
  "trajectory"                     => "other"
  "systemFieldStrength_T"          => 1.5
  "enc_lim_kspace_encoding_step_0" => Limit(0, 31, 16)
  "enc_lim_kspace_encoding_step_1" => Limit(0, 31, 16)
  "encodedFOV"                     => [258.065, 258.065, 258.065]
  "enc_lim_kspace_encoding_step_2" => Limit(0, 31, 16)
  "systemVendor"                   => "KomaMRI.jl"
  "encodedSize"                    => [32, 32, 32]

gives mixed results:

julia> acq = AcquisitionData( raw);

julia> acq.encodingSize
(32, 32)

julia> acq.fov
(258.065, 258.065, 258.065)

julia> acq.traj
1-element Vector{MRIBase.Trajectory{Float32}}:
 MRIBase.Trajectory{Float32}("Custom", Float32[-0.20536378 -0.19211452 … 0.19211386 0.20536311; 0.21198764 0.21198764 … 0.21198764 0.21198764], Float32[0.0, 6.2f-5, 0.000124, 0.000186, 0.000248, 0.00031, 0.000372, 0.000434, 0.000496, 0.000558  …  0.001364, 0.001426, 0.001488, 0.00155, 0.001612, 0.001674, 0.001736, 0.001798, 0.00186, 0.001922], 0.0f0, 0.001f0, 32, 32, 1, false, true)

and unfortunately still does not reconstruct:

julia> Nx, Ny, Nz = raw.params["reconSize"][1:3]
3-element Vector{Int64}:
 32
 32
 32

julia> reconParams = Dict{Symbol,Any}(:reco=>"direct", :reconSize=>(Nx, Ny, Nz));

julia> reconParams
Dict{Symbol, Any} with 2 entries:
  :reconSize => (32, 32, 32)
  :reco      => "direct"

julia> image = reconstruction(acq, reconParams)
ERROR: ArgumentError: Nodes x have dimension 2 != 3
Stacktrace:
  [1] initParams(k::Matrix{Float64}, N::Tuple{Int64, Int64, Int64}, dims::UnitRange{Int64}; kargs::Base.Pairs{Symbol, Int64, Tuple{Symbol, Symbol}, NamedTuple{(:m, :σ), Tuple{Int64, Int64}}})
    @ NFFT ~/.julia/packages/NFFT/RT2hs/src/precomputation.jl:20
  [2] initParams
    @ ~/.julia/packages/NFFT/RT2hs/src/precomputation.jl:3 [inlined]
  [3] NFFT.NFFTPlan(k::Matrix{Float64}, N::Tuple{Int64, Int64, Int64}; dims::UnitRange{Int64}, fftflags::Nothing, kwargs::Base.Pairs{Symbol, Int64, Tuple{Symbol, Symbol}, NamedTuple{(:m, :σ), Tuple{Int64, Int64}}})
    @ NFFT ~/.julia/packages/NFFT/RT2hs/src/implementation.jl:78
  [4] NFFTPlan
    @ ~/.julia/packages/NFFT/RT2hs/src/implementation.jl:73 [inlined]
  [5] macro expansion
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:49 [inlined]
  [6] macro expansion
    @ ./timing.jl:393 [inlined]
  [7] #plan_nfft#130
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:48 [inlined]
  [8] plan_nfft
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:46 [inlined]
  [9] #plan_nfft#2
    @ ~/.julia/packages/AbstractNFFTs/Xd3qS/src/derived.jl:13 [inlined]
 [10] plan_nfft
    @ ~/.julia/packages/AbstractNFFTs/Xd3qS/src/derived.jl:13 [inlined]
 [11] samplingDensity(acqData::AcquisitionData{Float32, 2}, shape::Tuple{Int64, Int64, Int64})
    @ MRIBase ~/.julia/packages/MRIBase/oNfYy/src/Datatypes/AcqData.jl:288
 [12] setupDirectReco(acqData::AcquisitionData{Float32, 2}, recoParams::Dict{Symbol, Any})
    @ MRIReco ~/src/MRIReco.jl/src/Reconstruction/RecoParameters.jl:19
 [13] reconstruction(acqData::AcquisitionData{Float32, 2}, recoParams::Dict{Symbol, Any})
    @ MRIReco ~/src/MRIReco.jl/src/Reconstruction/Reconstruction.jl:37
 [14] top-level scope
    @ REPL[80]:1

probably because of acq.traj and other incorrect parameters or dimensions of trajectory and other related data?

@curtcorum
Copy link
Contributor

The data is getting there but in the wrong dimensionality/indexing:

julia> acq.kdata
1×32×1 Array{Matrix{ComplexF32}, 3}:
[:, :, 1] =
 [-5.97573-6.46294im; -21.4331+0.217442im; … ; 21.9411-4.22948im; -11.1961-2.9045im;;]  …  [22.9858-1.1932im; 5.25461-2.63809im; … ; 19.4627-5.74995im; -4.43083-8.30492im;;]

julia> size( acq.kdata[1, 1, 1])
(1024, 1)

julia> size( acq.kdata[1, 2, 1])
(1024, 1)

julia> size( acq.kdata[1, 32, 1])
(1024, 1)

@curtcorum
Copy link
Contributor

curtcorum commented Mar 12, 2024

@cncastillo @beorostica

Line 448 suggests it may work if the DEF items are placed in the seq file. Line 456 may be the problem when not present in 3D sequence? I will try the Nx-z DEF items in the seq file.

https://github.com/curtcorum/KomaMRI.jl/blob/0a8a93de7c30fec886d27a052b1eff34ceca3a45/KomaMRIFiles/src/Sequence/Pulseq.jl#L448

if !haskey(seq.DEF,"Nx")

https://github.com/curtcorum/KomaMRI.jl/blob/0a8a93de7c30fec886d27a052b1eff34ceca3a45/KomaMRIFiles/src/Sequence/Pulseq.jl#L456

seq.DEF["Nz"] = Nz #Number of unique RF frequencies, in a 3D acquisition this should not work

@curtcorum
Copy link
Contributor

curtcorum commented Mar 12, 2024

This is a test where the Nx-z parameters are in the [DEFINITIONS] block of the seq file:

[DEFINITIONS]
AdcRasterTime 6.25e-05 
BlockDurationRaster 1e-05 
GradientRasterTime 1e-05 
RadiofrequencyRasterTime 1e-06
Nx 32
Ny 32
Nz 32 
(base) curt@green:~/src/KomaMRI$ julia
               _
   _       _ _(_)_     |  Documentation: https://docs.julialang.org
  (_)     | (_) (_)    |
   _ _   _| |_  __ _   |  Type "?" for help, "]?" for Pkg help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 1.9.4 (2023-11-14)
 _/ |\__'_|_|_|\__'_|  |  Official https://julialang.org/ release
|__/                   |

julia> using Revise

julia> using NativeFileDialog, MAT, MRIReco

julia> using KomaMRIBase, KomaMRICore, KomaMRIFiles, KomaMRIPlots

julia> sys = Scanner();

julia> seq = read_seq( seq_file)
[ Info: Loading sequence gre3d_cac.seq ...
Sequence[ τ = 42992.22 ms | blocks: 6444 | ADC: 1024 | GR: 7290 | RF: 1074 | DEF: 15 ]

julia> seq.DEF
Dict{String, Any} with 15 entries:
  "Nz"                       => 32.0
  "Nx"                       => 32.0
  "signature"                => ""
  "Ny"                       => 32.0
  "AdcRasterTime"            => 1.0e-7
  "BlockDurationRaster"      => 1.0e-5
  "extension"                => [0, 0, 0, 0, 0, 0, 0, 0, 0, 0  …  0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
  "PulseqVersion"            => 1004001
  "FileName"                 => "gre3d_cac.seq"
  "GradientRasterTime"       => 1.0e-5
  "TotalDuration"            => 42.9922
  "FOV"                      => [0.19, 0.19, 0.19]
  "Name"                     => "gre3d_cac"
  "RadiofrequencyRasterTime" => 1.0e-6
  "additional_text"          => "# Format of extension lists:\n# id type ref next_id\n# next_id of 0 terminates the list\n# Extension list is followed by extension spec…

julia> get_flip_angles( seq[1])
1-element Vector{Float64}:
 8.0

julia> raw = simulate(obj, seq, sys, sim_params=sim_params);
┌ Info: Running simulation in the GPU (Quadro RTX 8000)
│   koma_version = v"0.8.1"
│   sim_method = Bloch()
│   spins = 78126
│   time_points = 74945
└   adc_points = 32768
Progress: 100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| Time: 0:00:10
  simulated_blocks:  2168
  rf_blocks:         1075
  acq_samples:       32768
 10.773805 seconds (9.39 M allocations: 498.809 MiB, 3.71% gc time)

julia> raw.params
Dict{String, Any} with 23 entries:
  "protocolName"                   => "gre3d_cac"
  "institutionName"                => "Pontificia Universidad Catolica de Chile"
  "reconSize"                      => [32, 32, 1]
  "enc_lim_repetition"             => Limit(0, 0, 0)
  "enc_lim_set"                    => Limit(0, 0, 0)
  "enc_lim_segment"                => Limit(0, 0, 0)
  "userParameters"                 => Dict{String, Any}("Nblocks"=>2168, "gpu"=>1, "gpu_device"=>0, "type_sim_parts"=>Bool[0, 1, 0, 1, 0, 1, 0, 1, 0, 1  …  1, 0, 1, 0, …
  "enc_lim_phase"                  => Limit(0, 0, 0)
  "enc_lim_average"                => Limit(0, 0, 0)
  "enc_lim_slice"                  => Limit(0, 0, 0)
  "reconFOV"                       => Float32[190.0, 190.0, 1.0]
  "systemModel"                    => "v0.8.1"
  "H1resonanceFrequency_Hz"        => 63866203
  "patientName"                    => "brain3D"
  "enc_lim_contrast"               => Limit(0, 0, 0)
  "trajectory"                     => "other"
  "systemFieldStrength_T"          => 1.5
  "enc_lim_kspace_encoding_step_0" => Limit(0, 30, 16)
  "enc_lim_kspace_encoding_step_1" => Limit(0, 30, 16)
  ⋮                                => ⋮

julia> acq = AcquisitionData( raw);

julia> reconParams = Dict{Symbol,Any}(:reco=>"direct", :reconSize=>(32, 32, 32));

julia> image = reconstruction(acq, reconParams);
ERROR: ArgumentError: Nodes x have dimension 2 != 3
Stacktrace:
  [1] initParams(k::Matrix{Float64}, N::Tuple{Int64, Int64, Int64}, dims::UnitRange{Int64}; kargs::Base.Pairs{Symbol, Int64, Tuple{Symbol, Symbol}, NamedTuple{(:m, :σ), Tuple{Int64, Int64}}})
    @ NFFT ~/.julia/packages/NFFT/RT2hs/src/precomputation.jl:20
  [2] initParams
    @ ~/.julia/packages/NFFT/RT2hs/src/precomputation.jl:3 [inlined]
  [3] NFFT.NFFTPlan(k::Matrix{Float64}, N::Tuple{Int64, Int64, Int64}; dims::UnitRange{Int64}, fftflags::Nothing, kwargs::Base.Pairs{Symbol, Int64, Tuple{Symbol, Symbol}, NamedTuple{(:m, :σ), Tuple{Int64, Int64}}})
    @ NFFT ~/.julia/packages/NFFT/RT2hs/src/implementation.jl:78
  [4] NFFTPlan
    @ ~/.julia/packages/NFFT/RT2hs/src/implementation.jl:73 [inlined]
  [5] macro expansion
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:49 [inlined]
  [6] macro expansion
    @ ./timing.jl:393 [inlined]
  [7] #plan_nfft#130
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:48 [inlined]
  [8] plan_nfft
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:46 [inlined]
  [9] #plan_nfft#2
    @ ~/.julia/packages/AbstractNFFTs/Xd3qS/src/derived.jl:13 [inlined]
 [10] plan_nfft
    @ ~/.julia/packages/AbstractNFFTs/Xd3qS/src/derived.jl:13 [inlined]
 [11] samplingDensity(acqData::AcquisitionData{Float32, 2}, shape::Tuple{Int64, Int64, Int64})
    @ MRIBase ~/.julia/packages/MRIBase/oNfYy/src/Datatypes/AcqData.jl:288
 [12] setupDirectReco(acqData::AcquisitionData{Float32, 2}, recoParams::Dict{Symbol, Any})
    @ MRIReco ~/src/MRIReco.jl/src/Reconstruction/RecoParameters.jl:19
 [13] reconstruction(acqData::AcquisitionData{Float32, 2}, recoParams::Dict{Symbol, Any})
    @ MRIReco ~/src/MRIReco.jl/src/Reconstruction/Reconstruction.jl:37
 [14] top-level scope
    @ REPL[24]:1

@curtcorum
Copy link
Contributor

curtcorum commented Mar 12, 2024

Here is some more testing:

julia> raw.params["encodedFOV"]
3-element Vector{Float32}:
 190.0
 190.0
   1.0

julia> raw.params["encodedFOV"]=[190.0, 190.0, 190.0]
3-element Vector{Float64}:
 190.0
 190.0
 190.0

julia> raw.params["encodedSize"]=[32, 32, 32]
3-element Vector{Int64}:
 32
 32
 32

julia> raw.params
Dict{String, Any} with 23 entries:
  "protocolName"                   => "gre3d_cac"
  "institutionName"                => "Pontificia Universidad Catolica de Chile"
  "reconSize"                      => [32, 32, 32]
  "enc_lim_repetition"             => Limit(0, 0, 0)
  "enc_lim_set"                    => Limit(0, 0, 0)
  "enc_lim_segment"                => Limit(0, 0, 0)
  "userParameters"                 => Dict{String, Any}("Nblocks"=>2168, "gpu"=>1, "gpu_device"=>0, "type_sim_parts"=>Bool[0, 1, 0, 1, 0, 1, 0, 1, 0, 1  …  1, 0, 1, 0, 1, 0, 1, 0, 1, 0…
  "enc_lim_phase"                  => Limit(0, 0, 0)
  "enc_lim_average"                => Limit(0, 0, 0)
  "enc_lim_slice"                  => Limit(0, 0, 0)
  "reconFOV"                       => [190.0, 190.0, 190.0]
  "systemModel"                    => "v0.8.1"
  "H1resonanceFrequency_Hz"        => 63866203
  "patientName"                    => "brain3D"
  "enc_lim_contrast"               => Limit(0, 0, 0)
  "trajectory"                     => "other"
  "systemFieldStrength_T"          => 1.5
  "enc_lim_kspace_encoding_step_0" => Limit(0, 30, 16)
  "enc_lim_kspace_encoding_step_1" => Limit(0, 30, 16)
  "encodedFOV"                     => [190.0, 190.0, 190.0]
  "enc_lim_kspace_encoding_step_2" => Limit(0, 30, 16)
  "systemVendor"                   => "KomaMRI.jl"
  "encodedSize"                    => [32, 32, 32]

julia> acq = AcquisitionData( raw);

julia> reconParams = Dict{Symbol,Any}(:reco=>"direct", :reconSize=>(32, 32, 32));

julia> image = reconstruction(acq, reconParams);
ERROR: ArgumentError: Nodes x have dimension 2 != 3
Stacktrace:
  [1] initParams(k::Matrix{Float64}, N::Tuple{Int64, Int64, Int64}, dims::UnitRange{Int64}; kargs::Base.Pairs{Symbol, Int64, Tuple{Symbol, Symbol}, NamedTuple{(:m, :σ), Tuple{Int64, Int64}}})
    @ NFFT ~/.julia/packages/NFFT/RT2hs/src/precomputation.jl:20
  [2] initParams
    @ ~/.julia/packages/NFFT/RT2hs/src/precomputation.jl:3 [inlined]
  [3] NFFT.NFFTPlan(k::Matrix{Float64}, N::Tuple{Int64, Int64, Int64}; dims::UnitRange{Int64}, fftflags::Nothing, kwargs::Base.Pairs{Symbol, Int64, Tuple{Symbol, Symbol}, NamedTuple{(:m, :σ), Tuple{Int64, Int64}}})
    @ NFFT ~/.julia/packages/NFFT/RT2hs/src/implementation.jl:78
  [4] NFFTPlan
    @ ~/.julia/packages/NFFT/RT2hs/src/implementation.jl:73 [inlined]
  [5] macro expansion
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:49 [inlined]
  [6] macro expansion
    @ ./timing.jl:393 [inlined]
  [7] #plan_nfft#130
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:48 [inlined]
  [8] plan_nfft
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:46 [inlined]
  [9] #plan_nfft#2
    @ ~/.julia/packages/AbstractNFFTs/Xd3qS/src/derived.jl:13 [inlined]
 [10] plan_nfft
    @ ~/.julia/packages/AbstractNFFTs/Xd3qS/src/derived.jl:13 [inlined]
 [11] samplingDensity(acqData::AcquisitionData{Float32, 2}, shape::Tuple{Int64, Int64, Int64})
    @ MRIBase ~/.julia/packages/MRIBase/oNfYy/src/Datatypes/AcqData.jl:288
 [12] setupDirectReco(acqData::AcquisitionData{Float32, 2}, recoParams::Dict{Symbol, Any})
    @ MRIReco ~/src/MRIReco.jl/src/Reconstruction/RecoParameters.jl:19
 [13] reconstruction(acqData::AcquisitionData{Float32, 2}, recoParams::Dict{Symbol, Any})
    @ MRIReco ~/src/MRIReco.jl/src/Reconstruction/Reconstruction.jl:37
 [14] top-level scope
    @ REPL[39]:1

julia> acq.
encodingSize      fov               kdata             sequenceInfo      subsampleIndices  traj
julia> size(acq.kdata)
(1, 32, 1)

julia> acq.kdata[1,1,1]
1024×1 Matrix{ComplexF32}:
 -20.856276f0 + 7.8534117f0im
 -23.674004f0 + 0.69720984f0im
 -31.270561f0 - 11.370958f0im
  -24.04995f0 + 42.965137f0im
   39.66914f0 + 48.67775f0im
  57.792076f0 + 25.581245f0im
  21.598516f0 - 4.7730045f0im
  -82.14332f0 - 29.625975f0im
   -40.8282f0 - 20.949701f0im
  13.959399f0 + 4.1346235f0im
 -11.851673f0 + 28.035168f0im
  24.585815f0 + 9.329832f0im
   56.90302f0 + 7.7491198f0im
  17.272808f0 + 13.794415f0im
   -1.95158f0 - 26.888247f0im
   8.668823f0 - 80.32206f0im
              ⋮
 -20.988113f0 - 13.841843f0im
   28.07649f0 - 22.968506f0im
  45.602043f0 + 26.476562f0im
   55.74817f0 + 22.142883f0im
  14.360283f0 + 28.204256f0im
   17.00085f0 + 1.0903955f0im
   5.857233f0 - 23.596966f0im
 -25.678976f0 + 40.5332f0im
  20.144512f0 + 9.556927f0im
  15.072937f0 - 12.6632595f0im
   9.912712f0 + 7.5607147f0im
 -23.415226f0 - 46.04351f0im
  6.5762672f0 + 27.526182f0im
 -19.309336f0 + 33.920307f0im
  -28.33278f0 - 93.11535f0im
  55.864494f0 + 14.333662f0im

julia> acq.kdata[1,32,1]
1024×1 Matrix{ComplexF32}:
 -23.373734f0 + 28.150126f0im
  34.316868f0 + 29.423656f0im
  1.5463567f0 - 3.6292593f0im
 -52.664368f0 - 23.4733f0im
 -17.954971f0 - 32.63766f0im
  7.2143764f0 - 14.636225f0im
 -15.502787f0 + 27.550005f0im
  -10.73413f0 + 2.437481f0im
 -19.018562f0 + 26.472218f0im
  -18.13644f0 + 36.484882f0im
   8.976894f0 + 4.1743774f0im
  41.983356f0 + 26.89516f0im
   46.48584f0 - 12.104265f0im
   41.22889f0 - 8.658623f0im
  40.798172f0 + 39.47073f0im
  5.5201626f0 + 19.247066f0im
              ⋮
  -2.829985f0 - 8.459925f0im
   32.78006f0 - 21.172215f0im
  38.438377f0 + 26.038582f0im
   58.80835f0 + 18.670595f0im
  20.287176f0 + 25.217157f0im
  14.318125f0 + 8.467522f0im
  16.699306f0 - 23.905882f0im
 -20.987934f0 + 35.578873f0im
   15.35615f0 + 21.750748f0im
   12.56637f0 - 13.1215515f0im
  14.982704f0 + 2.081315f0im
  -32.10189f0 - 36.650608f0im
 -4.2777576f0 + 27.28927f0im
 -5.1820927f0 + 39.759125f0im
 -43.646645f0 - 87.57974f0im
   39.00147f0 + 13.069844f0im

julia> size( acq.kdata[1,32,1])
(1024, 1)

It looks like fixing the parameters is not enough, the kdata is still packed up wrong?

@curtcorum
Copy link
Contributor

curtcorum commented Mar 13, 2024

Testing 3D non-Cartesian, apparently has the same issues as 3D Cartesian.

(base) curt@green:~/src/KomaMRI$ julia
               _
   _       _ _(_)_     |  Documentation: https://docs.julialang.org
  (_)     | (_) (_)    |
   _ _   _| |_  __ _   |  Type "?" for help, "]?" for Pkg help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 1.9.4 (2023-11-14)
 _/ |\__'_|_|_|\__'_|  |  Official https://julialang.org/ release
|__/                   |

julia> using Markdown

julia> using InteractiveUtils

julia> using Revise

julia> using NativeFileDialog, MAT, MRIReco

julia> using KomaMRIBase, KomaMRICore, KomaMRIFiles, KomaMRIPlots

julia> sys = Scanner()
Scanner
  B0: Float64 1.5
  B1: Float64 1.0e-5
  Gmax: Float64 0.06
  Smax: Int64 500
  ADC_Δt: Float64 2.0e-6
  seq_Δt: Float64 1.0e-5
  GR_Δt: Float64 1.0e-5
  RF_Δt: Float64 1.0e-6
  RF_ring_down_T: Float64 2.0e-5
  RF_dead_time_T: Float64 0.0001
  ADC_dead_time_T: Float64 1.0e-5


julia> sys.RF_dead_time_T=1e-5
1.0e-5

julia> sys.RF_ring_down_T=10e-6
1.0e-5

julia> sys
Scanner
  B0: Float64 1.5
  B1: Float64 1.0e-5
  Gmax: Float64 0.06
  Smax: Int64 500
  ADC_Δt: Float64 2.0e-6
  seq_Δt: Float64 1.0e-5
  GR_Δt: Float64 1.0e-5
  RF_Δt: Float64 1.0e-6
  RF_ring_down_T: Float64 1.0e-5
  RF_dead_time_T: Float64 1.0e-5
  ADC_dead_time_T: Float64 1.0e-5


julia> seq_file = pick_file( "/home/curt/src"; filterlist="seq")
Gtk-Message: 19:49:26.593: Failed to load module "canberra-gtk-module"
Gtk-Message: 19:49:26.593: Failed to load module "canberra-gtk-module"

(julia:37505): GLib-GIO-WARNING **: 19:49:26.671: Failed to create file monitor for /home/curt/.config/glib-2.0/settings/keyfile: Unable to find default local file monitor type
"/home/curt/src/pulseq_champaign_imaging_llc/matlab/pulseq_m/zte_petra_cac.seq"

julia> seq = read_seq( seq_file)
[ Info: Loading sequence zte_petra_cac.seq ...
Sequence[ τ = 4328.22 ms | blocks: 1052 | ADC: 515 | GR: 3145 | RF: 536 | DEF: 16 ]

julia> seq.DEF
Dict{String, Any} with 16 entries:
  "signature"                => "e4359a04de6485703a97f6d7c56107ae"
  "AdcRasterTime"            => 1.0e-7
  "GradientRasterTime"       => 1.0e-5
  "TotalDuration"            => 4.32822
  "FOV"                      => [0.256, 0.256, 0.256]
  "Name"                     => "petra"
  "Nz"                       => 1
  "Nx"                       => 200
  "Ny"                       => 515
  "SamplesPerShell"          => 515.0
  "PulseqVersion"            => 1004001
  "BlockDurationRaster"      => 1.0e-5
  "extension"                => [0, 0, 0, 0, 0, 0, 0, 0, 0, 0  …  0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
  "FileName"                 => "zte_petra_cac.seq"
  "additional_text"          => "# Format of extension lists:\n# id type ref next_id\n# next_id of 0 terminates the list\n# Extension list is followed by extension specifications\n[EXTENSIONS]\n\n# Exten…
  "RadiofrequencyRasterTime" => 1.0e-6

julia> seq_plot = plot_seq(seq; range=[0,100], darkmode=true, slider=true,)
[ Info: Listening on: 127.0.0.1:7953, thread id: 1

julia> get_flip_angles( seq[3])
1-element Vector{Float64}:
 3.996

julia> sim_params = KomaMRICore.default_sim_params()
Dict{String, Any} with 9 entries:
  "return_type" => "raw"
  "Nblocks"     => 20
  "gpu"         => true
  "Nthreads"    => 1
  "gpu_device"  => 0
  "sim_method"  => Bloch()
  "precision"   => "f32"
  "Δt"          => 0.001
  "Δt_rf"       => 5.0e-5

julia> sim_params.
age       count     idxfloor  keys      maxprobe  ndel      slots     vals
julia> sim_params["Δt_rf"]
5.0e-5

julia> sim_params["Δt_rf"]=1e-6
1.0e-6

julia> obj = brain_phantom3D(; ss=4, start_end=[178, 182]);

julia> sim_params
Dict{String, Any} with 9 entries:
  "return_type" => "raw"
  "Nblocks"     => 20
  "gpu"         => true
  "Nthreads"    => 1
  "gpu_device"  => 0
  "sim_method"  => Bloch()
  "precision"   => "f32"
  "Δt"          => 0.001
  "Δt_rf"       => 1.0e-6

julia> raw = simulate(obj, seq, sys, sim_params=sim_params);
┌ Info: Running simulation in the GPU (Quadro RTX 8000)
│   koma_version = v"0.8.1"
│   sim_method = Bloch()
│   spins = 13027
│   time_points = 117639
└   adc_points = 103000
Progress: 100%|█████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| Time: 0:00:13
  simulated_blocks:  1091
  rf_blocks:         536
  acq_samples:       103000
 13.879795 seconds (9.95 M allocations: 548.532 MiB, 0.90% gc time)

julia> raw.params
Dict{String, Any} with 23 entries:
  "protocolName"                   => "petra"
  "institutionName"                => "Pontificia Universidad Catolica de Chile"
  "reconSize"                      => [26, 26, 1]
  "enc_lim_repetition"             => Limit(0, 0, 0)
  "enc_lim_set"                    => Limit(0, 0, 0)
  "enc_lim_segment"                => Limit(0, 0, 0)
  "userParameters"                 => Dict{String, Any}("Nblocks"=>1091, "gpu"=>1, "gpu_device"=>0, "type_sim_parts"=>Bool[0, 1, 0, 1, 0, 1, 0, 1, 0, 1  …  1, 0, 1, 0, 1, 0, 1, 0, 1, 0], "precision"=>"f3…
  "enc_lim_phase"                  => Limit(0, 0, 0)
  "enc_lim_average"                => Limit(0, 0, 0)
  "enc_lim_slice"                  => Limit(0, 0, 0)
  "reconFOV"                       => Float32[256.0, 256.0, 1.0]
  "systemModel"                    => "v0.8.1"
  "H1resonanceFrequency_Hz"        => 63866203
  "patientName"                    => "brain3D"
  "enc_lim_contrast"               => Limit(0, 0, 0)
  "trajectory"                     => "other"
  "systemFieldStrength_T"          => 1.5
  "enc_lim_kspace_encoding_step_0" => Limit(0, 24, 13)
  "enc_lim_kspace_encoding_step_1" => Limit(0, 25, 13)
  "encodedFOV"                     => Float32[256.0, 256.0, 1.0]
  "enc_lim_kspace_encoding_step_2" => Limit(0, 0, 0)
  "systemVendor"                   => "KomaMRI.jl"
  "encodedSize"                    => [25, 26, 1]

julia> acq = AcquisitionData( raw);

julia> begin
               Nx, Ny = raw.params["reconSize"][1:2];
               reconParams = Dict{Symbol,Any}(:reco=>"direct", :reconSize=>(Nx, Ny));
               image = reconstruction(acq, reconParams);;
       end

6-dimensional AxisArray{ComplexF32,6,...} with axes:
    :x, (0.0:10.24:256.0) mm
    :y, (0.0:9.846153846153847:246.15384615384616) mm
    :z, (0.0:1.0:0.0) mm
    :echos, 1:1
    :coils, 1:1
    :repetitions, 1:1
And data, a 26×26×1×1×1×1 Array{ComplexF32, 6}:
[:, :, 1, 1, 1, 1] =
         0.0+0.0im                0.0+0.0im                0.0+0.0im                 0.0+0.0im        …          0.0+0.0im               0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im                0.0+0.0im                 0.0+0.0im                   0.0+0.0im               0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im                0.0+0.0im                 0.0+0.0im                   0.0+0.0im               0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im                0.0+0.0im          -0.0236505+0.0992867im             0.0+0.0im               0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im        -0.00777292-0.0254325im      0.012173+0.253888im         0.115967+0.127116im          0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im         -0.0793708+0.124642im       0.239155-0.130121im   …    0.0784993+0.198522im          0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im         -0.0569608+0.0734885im    -0.132711+0.00239755im  -0.00624282+0.0227152im       -0.401928-0.0373109im   -0.136928+0.164801im         0.0+0.0im              0.0+0.0im
  -0.0222453+0.124349im       0.14743+0.162326im     0.0811325+0.282061im      -0.248402+0.172537im          0.51513+0.188474im     0.265433+0.066368im    0.257239+0.142123im         0.0+0.0im
    -0.12643-0.103366im    -0.0969422-0.0783985im    0.0239488+0.0854764im     0.0508624-0.0679528im       0.0785211+0.757174im    -0.132184-0.129531im  -0.0516888-0.192597im         0.0+0.0im
   0.0145846+0.0472643im   -0.0715923-0.0530962im    0.0238834+0.221171im      -0.228209+0.151941im        -0.311803+1.07776im     0.0541115+0.721079im   0.0380341-0.126004im         0.0+0.0im
 -0.00197875+0.0977978im   -0.0290478+0.0666364im    0.0956486+0.0605195im     -0.020544+0.442739im   …   -0.0646548+1.2661im      -0.303796+1.2412im     -0.133189-0.060337im         0.0+0.0im
    0.037106+0.107172im      0.184824+0.0822097im    -0.133251-0.186752im       0.141473+0.988213im        -0.092847+1.7226im      -0.342548+1.21904im    0.0885918+0.442872im         0.0+0.0im
  -0.0738069-0.167435im     0.0514082+0.150021im     0.0903598-0.265952im      0.0595301+0.800169im         0.033736+2.22846im     0.0331009+1.57372im     0.261175+0.729681im   0.0829206-0.137986im
    0.129208+0.00561913im   -0.171928+0.324073im      0.142662-0.406718im       -0.24646+1.12177im         0.0174301+2.41918im    -0.0122723+2.06395im     0.196475+0.885627im         0.0+0.0im
   -0.194472-0.135396im     -0.157046+0.263632im      -0.23762-0.0996091im     -0.270387+1.47725im         -0.029394+2.2901im     -0.0107559+1.50242im    0.0865599+0.541747im         0.0+0.0im
    0.079245-0.0244986im     0.206436+0.0963664im    0.0358069-0.0455092im      0.218636+0.809432im   …     0.168469+2.19265im     -0.274536+1.46205im   -0.0206982-0.0357301im        0.0+0.0im
    0.111329-0.0116827im    0.0941668+0.186362im     0.0414217-0.0693443im      0.261088+0.25694im         0.0234841+1.6421im      -0.298419+1.49593im    0.0610674+0.0813188im        0.0+0.0im
   -0.099452+0.189673im     0.0295172+0.201996im    -0.0141953+0.182957im       0.151031+0.0302758im       -0.117989+1.24655im     0.0375125+0.718611im   0.0532909-0.0469774im        0.0+0.0im
         0.0+0.0im          -0.132238-0.202624im    -0.0268976+0.13931im       0.0730343-0.285966im      -0.00140481+0.758735im    0.0919741+0.324555im         0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im          0.0678213-0.0374102im      0.077156+0.133285im        -0.115163-0.118061im          0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im          0.0567113-0.04648im       0.0242621+0.280348im   …    -0.061529+0.135817im          0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im                0.0+0.0im           0.0612392+0.0382043im             0.0+0.0im               0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im                0.0+0.0im                 0.0+0.0im                   0.0+0.0im               0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im                0.0+0.0im                 0.0+0.0im                   0.0+0.0im               0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im                0.0+0.0im                 0.0+0.0im                   0.0+0.0im               0.0+0.0im              0.0+0.0im              0.0+0.0im
         0.0+0.0im                0.0+0.0im                0.0+0.0im                 0.0+0.0im        …          0.0+0.0im               0.0+0.0im              0.0+0.0im              0.0+0.0im

julia> image_plot = plot_image( abs.( image[:, :, 1]);)

(this does give a 3d to 2d projection image)

julia> reconParams = Dict{Symbol,Any}(:reco=>"direct", :reconSize=>(32, 32, 32));

julia> image = reconstruction(acq, reconParams);
ERROR: ArgumentError: Nodes x have dimension 2 != 3
Stacktrace:
  [1] initParams(k::Matrix{Float64}, N::Tuple{Int64, Int64, Int64}, dims::UnitRange{Int64}; kargs::Base.Pairs{Symbol, Int64, Tuple{Symbol, Symbol}, NamedTuple{(:m, :σ), Tuple{Int64, Int64}}})
    @ NFFT ~/.julia/packages/NFFT/RT2hs/src/precomputation.jl:20
  [2] initParams
    @ ~/.julia/packages/NFFT/RT2hs/src/precomputation.jl:3 [inlined]
  [3] NFFT.NFFTPlan(k::Matrix{Float64}, N::Tuple{Int64, Int64, Int64}; dims::UnitRange{Int64}, fftflags::Nothing, kwargs::Base.Pairs{Symbol, Int64, Tuple{Symbol, Symbol}, NamedTuple{(:m, :σ), Tuple{Int64, Int64}}})
    @ NFFT ~/.julia/packages/NFFT/RT2hs/src/implementation.jl:78
  [4] NFFTPlan
    @ ~/.julia/packages/NFFT/RT2hs/src/implementation.jl:73 [inlined]
  [5] macro expansion
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:49 [inlined]
  [6] macro expansion
    @ ./timing.jl:393 [inlined]
  [7] #plan_nfft#130
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:48 [inlined]
  [8] plan_nfft
    @ ~/.julia/packages/NFFT/RT2hs/src/NFFT.jl:46 [inlined]
  [9] #plan_nfft#2
    @ ~/.julia/packages/AbstractNFFTs/Xd3qS/src/derived.jl:13 [inlined]
 [10] plan_nfft
    @ ~/.julia/packages/AbstractNFFTs/Xd3qS/src/derived.jl:13 [inlined]
 [11] samplingDensity(acqData::AcquisitionData{Float32, 2}, shape::Tuple{Int64, Int64, Int64})
    @ MRIBase ~/.julia/packages/MRIBase/oNfYy/src/Datatypes/AcqData.jl:288
 [12] setupDirectReco(acqData::AcquisitionData{Float32, 2}, recoParams::Dict{Symbol, Any})
    @ MRIReco ~/src/MRIReco.jl/src/Reconstruction/RecoParameters.jl:19
 [13] reconstruction(acqData::AcquisitionData{Float32, 2}, recoParams::Dict{Symbol, Any})
    @ MRIReco ~/src/MRIReco.jl/src/Reconstruction/Reconstruction.jl:37
 [14] top-level scope
    @ REPL[32]:1

@curtcorum
Copy link
Contributor

@cncastillo @beorostica

Last time I tried using 3D k-spaces (nodes) with MRIReco it always gave an error, so we are just storing 2D k-spaces in the raw data (or at least that is the default ndims=2 of signal_to_raw_data). This could create this dimension mismatch. Maybe we just need to check if :reconSize has 3 components to save a 3D k-space, otherwise a 2D k-space. This is the main thing we need to explore to fix this issue.

KomaMRICore/src/rawdata/ISMRMRD.jl definitely is only doing 2D.

function signal_to_raw_data(
    signal, seq;
    phantom_name="Phantom", sys=Scanner(), sim_params=Dict{String,Any}(), ndims=2
)
    version = string(VersionNumber(Pkg.TOML.parsefile(joinpath(@__DIR__, "..", "..", "Project.toml"))["version"]))
    #Number of samples and FOV
    _, ktraj = get_kspace(seq) #kspace information
    mink = minimum(ktraj, dims=1)
    maxk = maximum(ktraj, dims=1)
    Wk = maxk .- mink
    Δx = 1 ./ Wk[1:2] #[m] Only x-y
    Nx = get(seq.DEF, "Nx", 1)
    Ny = get(seq.DEF, "Ny", 1)
    Nz = get(seq.DEF, "Nz", 1)
    if haskey(seq.DEF, "FOV")
        FOVx, FOVy, _ = seq.DEF["FOV"] #[m]
        if FOVx > 1 FOVx *= 1e-3 end #mm to m, older versions of Pulseq saved FOV in mm
        if FOVy > 1 FOVy *= 1e-3 end #mm to m, older versions of Pulseq saved FOV in mm
        Nx = round(Int64, FOVx / Δx[1])
        Ny = round(Int64, FOVy / Δx[2])
    else
        FOVx = Nx * Δx[1]
        FOVy = Ny * Δx[2]
    end

Is the idea to have it do 3D again, or fool it at least?

@curtcorum
Copy link
Contributor

curtcorum commented Mar 13, 2024

@curtcorum curtcorum linked a pull request Mar 13, 2024 that will close this issue
@cncastillo cncastillo linked a pull request May 17, 2024 that will close this issue
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
bug Something isn't working feature
Projects
None yet
Development

Successfully merging a pull request may close this issue.

3 participants