-
Notifications
You must be signed in to change notification settings - Fork 21
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
MRD Data opened in Python appears corrupted #214
Comments
Sorry, I believe I have confused some versions. The UI reported 0.7.6, which I believe is just for the core. I installed the latest release of the package, which overall is v0.7.5, I think (but contains core v0.7.6). |
There seems to be more beyond XML: Unexpectedly, the Also some other fields have unexpected values (e.g. flags, number_of_samples,...)
|
Hi @gabuzi, thanks for noticing this. I was aware that there could be a problem with the generated .mrd and Python's MRD reader. As when I tried it it also did not work. This helps us fix it. The problem with the Could you send us the code to open an .mrd in Python? I assume you opened the UI, pressed simulate, and tried opening the resulting .mrd located in the temp folder. Just to try to reproduce this. |
Hi @cncastillo import ismrmrd
with ismrmrd.File("Koma_signal.mrd") as f:
# hdr = f['dataset'].header # <- this fails because of the xml error
acquisitions = f['dataset'].acquisitions[:]
acq_headers = [a.getHead() for a in acquisitions]
# now you can look at the individual acq headers
# the encodingCounters of the n-th readout can e.g. be inspected via print(acq_headers[n].idx) |
@johannesmayer Regarding the binary acquisition headers, I can confirm that they do look ok with hdfview for data that is problematic to read via Python! |
Hi could any of you attach a valid MRD to compare? Edit: I generated one in MATLAB with https://github.com/ismrmrd/ismrmrd/blob/master/examples/matlab/create_dataset.m |
I found where the problem is. The ismrmrd spec is quite strict regarding the binary layout of the data. See https://ismrmrd.readthedocs.io/en/latest/mrd_raw_data.html#acquisitionheader. It defines the size of the fields and, crucially, their offset. However, hdf5 doesn't necessarily require strict binary layouts if the data is described well, and it seems like the files generated by KomaMRI do insert some padding in the binary output, and a header struct as it appears in hdf as seen from python ends up being 352 bytes instead of the required 340. The following code illustrates this (the second file is one from mridata.org that loads fine in python): import ismrmrd
def analyze_dtypes(ismrmrd_ds):
data_dtype = ds._file['dataset']['data'].dtype
print("data dtype (includes head)")
print_dtype(data_dtype)
head_dtype = ds._file['dataset']['data'][0]['head'].dtype
print("head dtype")
print_dtype(head_dtype)
def print_dtype(dtype):
print(dtype.descr)
print(dtype.str)
ds = ismrmrd.Dataset("Koma_0.7.5_default_signal.mrd", create_if_needed=False)
print("KomaMRI data")
analyze_dtypes(ds)
ds.close()
ds = ismrmrd.Dataset("52c2fd53-d233-4444-8bfd-7c454240d314.h5", create_if_needed=False)
print("mridata.org data")
analyze_dtypes(ds)
ds.close() which for me gives
We can find once |
It could be a bug in MRIFiles. I took a look to the function Code: using MRIFiles
structinfo(T) = [((i!=1 ? sum(sizeof(fieldtype(T,j)) for j=1:i-1) : 0), sizeof(fieldtype(T,i)), fieldname(T,i), fieldtype(T,i)) for i = 1:fieldcount(T)]
myfieldoffset(T, i) = i!=1 ? sum(sizeof(fieldtype(T,j)) for j=1:i-1) : 0
a = fieldoffset.(Ref(MRIFiles.AcquisitionHeaderImmutable),1:fieldcount(MRIFiles.AcquisitionHeaderImmutable)) .|> Int
b = myfieldoffset.(Ref(MRIFiles.AcquisitionHeaderImmutable),1:fieldcount(MRIFiles.AcquisitionHeaderImmutable)) .|> Int
[b a [0; diff(a.-b; dims=1)]] Output: julia> [b a [0; diff(a.-b; dims=1)]]
24×3 Matrix{Int64}:
0 0 0
2 8 6 # <- |V6
10 16 0
14 20 0
18 24 0
22 28 0
34 40 0
36 42 0
38 44 0
40 48 2 # <- |V2
168 176 0
170 178 0
172 180 0
174 182 0
176 184 0
178 188 2 # <- |V2
182 192 0
194 204 0
206 216 0
218 228 0
230 240 0
242 252 0
276 288 2 # <- |V2
308 320 0 The byte offsets i calculate correspond to the ones described in the official MRD docs (AcquisitionHeader) # byte offset, size in bytes, field, fieldtype
julia> structinfo(MRIFiles.AcquisitionHeaderImmutable)
24-element Vector{Tuple{Int64, Int64, Symbol, DataType}}:
(0, 2, :version, UInt16)
(2, 8, :flags, UInt64)
(10, 4, :measurement_uid, UInt32)
(14, 4, :scan_counter, UInt32)
(18, 4, :acquisition_time_stamp, UInt32)
(22, 12, :physiology_time_stamp, Tuple{UInt32, UInt32, UInt32})
(34, 2, :number_of_samples, UInt16)
(36, 2, :available_channels, UInt16)
(38, 2, :active_channels, UInt16)
(40, 128, :channel_mask, NTuple{16, UInt64})
(168, 2, :discard_pre, UInt16)
(170, 2, :discard_post, UInt16)
(172, 2, :center_sample, UInt16)
(174, 2, :encoding_space_ref, UInt16)
(176, 2, :trajectory_dimensions, UInt16)
(178, 4, :sample_time_us, Float32)
(182, 12, :position, Tuple{Float32, Float32, Float32})
(194, 12, :read_dir, Tuple{Float32, Float32, Float32})
(206, 12, :phase_dir, Tuple{Float32, Float32, Float32})
(218, 12, :slice_dir, Tuple{Float32, Float32, Float32})
(230, 12, :patient_table_position, Tuple{Float32, Float32, Float32})
(242, 34, :idx, EncodingCountersImmutable)
(276, 32, :user_int, NTuple{8, Int32})
(308, 32, :user_float, NTuple{8, Float32}) I will keep looking at it tomorrow :) |
Oh, this is really bad. I have not seen https://ismrmrd.readthedocs.io/en/latest/mrd_raw_data.html#acquisitionheader and assumed that the offsets correspond to the in-memory layout of the C struct (which in Julia is the immutable header type). And when I worked on this, I actually have seen the padded data in some public MRM files. So am I reading https://ismrmrd.readthedocs.io/en/latest/mrd_raw_data.html#acquisitionheader correctly, that the layout is always tight, i.e. that do padding is involved? In that case One issue though is that I don't understand how we can present the "unpadded" struct in Julia. Does somebody know how this is handled in the C bindings? The C struct will also have such padding. |
Looking at https://github.com/ismrmrd/ismrmrd/blob/master/libsrc/serialization.cpp#L11 it seems that the C struct is also serialized directly. I don't get that. Does it mean that the C struct is not padded? To solve this issue it would be great to ping somebody, how has more experience with the C/C++ implementation and can clarify these questions. |
@tknopp Thanks for tuning in! Looking at https://github.com/ismrmrd/ismrmrd/blob/d364e03d3faa3ca516da7807713b5acc72218a37/include/ismrmrd/ismrmrd.h#L259, especially including up to L315, I'm quite confident that the C/C++ datastructures are not padded. At least not within the AcquisitionHeader struct. They explicitly check at compile time ( |
C structs are padded by default unless you tell the compiler not to do so. But have found the line where this is actually done: I already did a little bit of research how we can solve this. Unfortunately we can't disable structural padding. But I think we will simply convert the struct to an |
And for reference: In our unit test we use the ISMRMRD data: And that seems to be padded. Would be great if someone could have a look and inspect that data. I can't load it with the python ismrmrd library but don't know where the issue is. So what I would definitely also need to fix the issue is correct data. Ideally again, some cartesian and spiral data. What is also important is that the data is small, so that we can use it in the unit tests. Edit: I was wrong. The data is not padded but I misinterpreted the way Julia reads it. Everything is fine. |
Ok, here is the fix: MagneticResonanceImaging/MRIReco.jl#167 Would be great if someone could test it. |
Hi @tknopp, Thanks for looking into this. I have tested the fix with cartesian and spiral data generated with Koma, and I think it does in fact solve the problem, as I do not see the data being corrupted in Python anymore. I will do further testing, but I believe I will close this issue soon :). |
Great, thanks so much for reporting. I did, by the way, already do the releases (including the move of the Limit... types) so you can already rely on that in Koma. |
Thanks! we will start using the Limit type in MRIBase now. |
Hi all!
I've just taken KomaMRI.jl for a quick test-drive and I noticed that I can't open the resulting ismrmrd files with the ismrmrd-python package. Reason seems to be a second nested level of the
<UserParameters>
element (see below). I also spotted a potential problem with the gpu user param value (annotated in the snippet below).I used v0.7.6 with julia 1.9.1
EDIT 1-Dec-2023 (@cncastillo):
Hi I think this issue can be separated in two:
<userParameters>
tag within XML header of exported MRD #236AcquisitionHeader
in exported MRDThe text was updated successfully, but these errors were encountered: