Skip to content

Commit

Permalink
Port #100 to master and commit.
Browse files Browse the repository at this point in the history
  • Loading branch information
tknopp committed Oct 11, 2022
1 parent 40508af commit 2e29cb9
Show file tree
Hide file tree
Showing 3 changed files with 48 additions and 16 deletions.
13 changes: 13 additions & 0 deletions MRIFiles/src/ISMRMRD/HDF5Types.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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) )
Expand All @@ -240,6 +244,15 @@ 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})

Expand Down
43 changes: 27 additions & 16 deletions MRIFiles/src/ISMRMRD/ISMRMRD.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,30 +12,41 @@ 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)
profiles = Profile[]

for m=1:M
head = read_header(d[m].head)
D = Int(head.trajectory_dimensions)
chan = Int(head.active_channels)
M = size(h["/$(dataset)/data"]) # for some reason this is a matrix

traj = isempty(d[m].traj) ? Matrix{Float32}(undef,0,0) : reshape(d[m].traj, D, :)
# 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)

if !isempty(d[m].data)
dat = reshape(reinterpret(ComplexF32,d[m].data), :, chan)
else
dat = Matrix{ComplexF32}(undef,0,0)
end
profiles = Profile[]

push!(profiles, Profile(head,traj,dat) )
for m in CartesianIndices(M)
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)
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)
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
Expand Down
8 changes: 8 additions & 0 deletions MRIFiles/test/testISMRMRD.jl
Original file line number Diff line number Diff line change
Expand Up @@ -36,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)
Expand Down

0 comments on commit 2e29cb9

Please sign in to comment.