diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 0e7ae494..b830376e 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -45,8 +45,8 @@ jobs: shell: julia --color=yes --project=. {0} run: | using Pkg - Pkg.add(PackageSpec(path=pwd(), subdir="MRIBase", rev="master")) - Pkg.add(PackageSpec(path=pwd(), subdir="MRIFiles", rev="master")) + Pkg.add(PackageSpec(path=pwd(), subdir="MRIBase")) + Pkg.add(PackageSpec(path=pwd(), subdir="MRIFiles")) - uses: julia-actions/julia-buildpkg@v1 - uses: julia-actions/julia-runtest@v1 - uses: julia-actions/julia-processcoverage@v1 @@ -65,8 +65,8 @@ jobs: shell: julia --color=yes --project=. {0} run: | using Pkg - Pkg.add(PackageSpec(path=pwd(), subdir="MRIBase", rev="master")) - Pkg.add(PackageSpec(path=pwd(), subdir="MRIFiles", rev="master")) + Pkg.add(PackageSpec(path=pwd(), subdir="MRIBase")) + Pkg.add(PackageSpec(path=pwd(), subdir="MRIFiles")) - run: | julia --project=docs -e ' using Pkg diff --git a/MRIFiles/src/ISMRMRD/HDF5Types.jl b/MRIFiles/src/ISMRMRD/HDF5Types.jl index 37a541d1..e9829d9b 100644 --- a/MRIFiles/src/ISMRMRD/HDF5Types.jl +++ b/MRIFiles/src/ISMRMRD/HDF5Types.jl @@ -219,6 +219,10 @@ struct HDF5_Acquisition data::HVL_T{Float32} #for some reason this is stored as float end +struct HDF5_AcquisitionHeaderOnly + head::AcquisitionHeaderImmutable +end + function get_hdf5type_acquisition() off = i -> Cint(fieldoffset(HDF5_Acquisition,i)) datatype = h5t_create(H5T_COMPOUND, sizeof(HDF5_Acquisition) ) @@ -240,6 +244,16 @@ function get_hdf5type_acquisition() return datatype end +function get_hdf5type_acquisition_header_only() + off = i -> Cint(fieldoffset(HDF5_AcquisitionHeaderOnly,i)) + datatype = h5t_create(H5T_COMPOUND, sizeof(HDF5_AcquisitionHeaderOnly) ) + d1 = get_hdf5type_acquisitionheader() + h5t_insert(datatype, "head", off(1), d1) + h5t_close(d1) + + return datatype +end + function writeProfiles(file, dataset, profiles::Vector{Profile}) diff --git a/MRIFiles/src/ISMRMRD/ISMRMRD.jl b/MRIFiles/src/ISMRMRD/ISMRMRD.jl index 8c916472..96f02fa1 100644 --- a/MRIFiles/src/ISMRMRD/ISMRMRD.jl +++ b/MRIFiles/src/ISMRMRD/ISMRMRD.jl @@ -12,30 +12,44 @@ end reads the `ISMRMRDFile` f and stores the result in a `RawAcquisitionDataObject` """ -function MRIBase.RawAcquisitionData(f::ISMRMRDFile, dataset="dataset") +function MRIBase.RawAcquisitionData(f::ISMRMRDFile, dataset="dataset"; + slice=nothing, repetition=nothing, contrast=nothing) h5open(f.filename) do h headerStr = read(h["/$(dataset)/xml"]) params = GeneralParameters(headerStr[1]) - d = read(h["/$(dataset)/data"]) - M = length(d) + M = size(h["/$(dataset)/data"]) # for some reason this is a matrix + + # In the next lines we only read the headers of the profiles such that + # we can later filter on them + dTypeHeader = get_hdf5type_acquisition_header_only() + headerData = Array{HDF5_AcquisitionHeaderOnly}(undef, M) + HDF5.API.h5d_read(h["/$(dataset)/data"], dTypeHeader, H5S_ALL, H5S_ALL, H5P_DEFAULT, headerData) + profiles = Profile[] - for m=1:M - head = read_header(d[m].head) - D = Int(head.trajectory_dimensions) - chan = Int(head.active_channels) + for m in CartesianIndices(M) - traj = isempty(d[m].traj) ? Matrix{Float32}(undef,0,0) : reshape(d[m].traj, D, :) + if (repetition == nothing || headerData[m].head.idx.repetition in repetition) && + (slice == nothing || headerData[m].head.idx.slice in slice) && + (contrast == nothing || headerData[m].head.idx.contrast in contrast) - if !isempty(d[m].data) - dat = reshape(reinterpret(ComplexF32,d[m].data), :, chan) - else - dat = Matrix{ComplexF32}(undef,0,0) - end + d = h["/$(dataset)/data"][m] # Here we read the header again. This could/should be avoided + head = read_header(d.head) + D = Int(head.trajectory_dimensions) + chan = Int(head.active_channels) - push!(profiles, Profile(head,traj,dat) ) + traj = isempty(d.traj) ? Matrix{Float32}(undef,0,0) : reshape(d.traj, D, :) + + if !isempty(d.data) + dat = reshape(reinterpret(ComplexF32,d.data), :, chan) + else + dat = Matrix{ComplexF32}(undef,0,0) + end + + push!(profiles, Profile(head,traj,dat) ) + end end return RawAcquisitionData(params, profiles) end diff --git a/test/testISMRMRD.jl b/test/testISMRMRD.jl index d8e36d7e..3637f83a 100644 --- a/test/testISMRMRD.jl +++ b/test/testISMRMRD.jl @@ -1,6 +1,5 @@ using MRIReco using Test -using HTTP @testset "ISMRMRD" begin @@ -37,6 +36,14 @@ IrecoCopy = reconstruction(AcquisitionData(acqCopy), params) @test IrecoCopy == Ireco +### test partial read + +acqPartial = RawAcquisitionData(f, slice=1) # slice 1 is not contained in the file + +@test length(acqPartial.profiles) == 0 + +### test spiral data + filename = joinpath(datadir, "simple_spiral.h5") f = ISMRMRDFile(filename)