[ismrmrd] 54/281: Working on the matlab
Ghislain Vaillant
ghisvail-guest at moszumanska.debian.org
Wed Jan 14 20:00:55 UTC 2015
This is an automated email from the git hooks/post-receive script.
ghisvail-guest pushed a commit to annotated tag ismrmrd0.5
in repository ismrmrd.
commit 5d425b8a1871cc717cc4087f219b048d37a4df3b
Author: Souheil Inati <souheil.inati at nih.gov>
Date: Thu Sep 27 11:32:05 2012 -0400
Working on the matlab
---
ismrmrd_hdf5.cpp | 5 +-
matlab/+ismrmrd/Acquisition.m | 21 ++--
matlab/+ismrmrd/AcquisitionHeader.m | 23 +++-
matlab/+ismrmrd/file.m | 217 ++++++++++++++++++++++++++----------
matlab/testread.m | 66 -----------
matlab/testrecon.m | 10 --
schema/ismrmrd.xsd | 2 +-
7 files changed, 192 insertions(+), 152 deletions(-)
diff --git a/ismrmrd_hdf5.cpp b/ismrmrd_hdf5.cpp
index 4eec514..6f3c8af 100644
--- a/ismrmrd_hdf5.cpp
+++ b/ismrmrd_hdf5.cpp
@@ -499,9 +499,10 @@ boost::shared_ptr<std::string> IsmrmrdDataset::readHeader()
try {
boost::shared_ptr<DataSet> xml_dataset = boost::shared_ptr<DataSet>(new DataSet(file_->openDataSet(xml_header_path_.c_str())));
boost::shared_ptr<DataType> datatype(new StrType(0, H5T_VARIABLE));
- xml_dataset->read(*ret,*datatype);
+ DataType dt = xml_dataset->getDataType();
+ xml_dataset->read(*ret,dt);
} catch( Exception& e ) {
- std::cout << "Exception caught while readin XML Header to HDF5 file" << std::endl;
+ std::cout << "Exception caught while reading XML Header to HDF5 file" << std::endl;
std::cout << e.getDetailMsg() << std::endl;
return ret;
}
diff --git a/matlab/+ismrmrd/Acquisition.m b/matlab/+ismrmrd/Acquisition.m
index 5d15dc7..e9e4fe3 100644
--- a/matlab/+ismrmrd/Acquisition.m
+++ b/matlab/+ismrmrd/Acquisition.m
@@ -5,8 +5,8 @@ classdef Acquisition
properties
head = ismrmrd.AcquisitionHeader;
- traj = [];
- data = [];
+ traj = single([]);
+ data = single([]);
end % Properties
@@ -14,7 +14,16 @@ classdef Acquisition
methods
function obj = set.head(obj,v)
- obj.head = v;
+ if isa(v,'ismrmrd.AcquisitionHeader')
+ obj.head = v;
+ else
+ % not of the correct type, hope it's a struct
+ % and try to copy one element at a time
+ u = fieldnames(obj.head);
+ for p = 1:length(u)
+ obj.head = setfield(obj.head,u{p},getfield(v,u{p}));
+ end
+ end
end
function obj = set.traj(obj,v)
@@ -22,11 +31,7 @@ classdef Acquisition
end
function obj = set.data(obj,v)
- if isreal(v)
- obj.data = single(v(1:2:end) + 1j*v(2:2:end));
- else
- obj.data = single(v);
- end
+ obj.data = single(v);
end
function b = isFlagSet(obj,flag)
diff --git a/matlab/+ismrmrd/AcquisitionHeader.m b/matlab/+ismrmrd/AcquisitionHeader.m
index 1ad6853..f333c74 100644
--- a/matlab/+ismrmrd/AcquisitionHeader.m
+++ b/matlab/+ismrmrd/AcquisitionHeader.m
@@ -2,11 +2,11 @@ classdef AcquisitionHeader
properties
version = uint16(0); % First unsigned int indicates the version %
- flag = uint64(0); % bit field with flags %
+ flags = uint64(0); % bit field with flags %
measurement_uid = uint32(0); % Unique ID for the measurement %
scan_counter = uint32(0); % Current acquisition number in the measurement %
acquisition_time_stamp = uint32(0); % Acquisition clock %
- physiology_time_stamp = zeros(1,3,'uint32'); % Physiology time stamps, e.g. ecg, breating, etc. %
+ physiology_time_stamp = zeros(3,1,'uint32'); % Physiology time stamps, e.g. ecg, breating, etc. %
number_of_samples = uint16(0); % Number of samples acquired %
available_channels = uint16(0); % Available coils %
active_channels = uint16(0); % Active coils on current acquisiton %
@@ -31,8 +31,8 @@ classdef AcquisitionHeader
obj.version = uint16(v);
end
- function obj = set.flag(obj,v)
- obj.flag = uint64(v);
+ function obj = set.flags(obj,v)
+ obj.flags = uint64(v);
end
function obj = set.measurement_uid(obj,v)
@@ -70,7 +70,7 @@ classdef AcquisitionHeader
if (length(v)~=16)
error('AcquisitionHeader.channel_mask must have 16 elements')
end
- obj.channel_mask = uint16(v);
+ obj.channel_mask = uint64(v);
end
function obj = set.discard_pre(obj,v)
@@ -117,6 +117,19 @@ classdef AcquisitionHeader
end
obj.patient_table_position = single(v);
end
+
+ function obj = set.idx(obj,v)
+ if isa(v,'ismrmrd.EncodingCounters')
+ obj.idx = v;
+ else
+ % not of the correct type, hope it's a struct
+ % and try to copy one element at a time
+ u = fieldnames(obj.idx);
+ for p = 1:length(u)
+ obj.idx = setfield(obj.idx,u{p},getfield(v,u{p}));
+ end
+ end
+ end
function obj = set.user_int(obj,v)
if (length(v)~=8)
diff --git a/matlab/+ismrmrd/file.m b/matlab/+ismrmrd/file.m
index fb4aa5e..79d42ef 100644
--- a/matlab/+ismrmrd/file.m
+++ b/matlab/+ismrmrd/file.m
@@ -4,7 +4,6 @@ classdef file
properties
fid = -1;
filename = '';
- xml = ''; % remove this and replace with a get method
datapath = '';
xmlpath = '';
end
@@ -41,102 +40,200 @@ classdef file
H5G.close(group_id);
end
H5P.close(lapl_id);
-
- % Read the XML
- obj.readxml();
-
- % Check if the Data exists
- % if it does not exist, create it
- lapl_id=H5P.create('H5P_LINK_ACCESS');
- if (H5L.exists(obj.fid, obj.datapath, lapl_id) == 0)
- % Data does not exist
- % create with rank 2, unlimited, and set the chunk size
- dims = [1 1];
- maxdims = [H5ML.get_constant_value('H5S_UNLIMITED') 1];
- chunk = [1 1];
- space = H5S.create_simple(2, dims, maxdims);
- dcpl = H5P.create('H5P_DATASET_CREATE');
- H5P.set_chunk (dcpl, chunk);
- data_id = H5D.create(obj.fid, obj.datapath, ...
- ismrmrd.hdf5_datatypes.getType_Acquisition, ...
- space, dcpl);
- H5D.close(data_id);
- end
- H5P.close(lapl_id);
end
-
- % destructor
- % is this necessary?
- %function delete(obj)
- % obj.close();
- %end
function obj = close(obj)
H5F.close(obj.fid);
end
- function obj = readxml(obj)
+ function xmlstring = readxml(obj)
- % Set variable length string type
- xml_dtype = H5T.copy('H5T_C_S1');
- H5T.set_size(xml_dtype,'H5T_VARIABLE');
-
% Check if the XML header exists
- % if it does not exist, create it
% TODO: set it's value to the default
lapl_id=H5P.create('H5P_LINK_ACCESS');
if (H5L.exists(obj.fid,obj.xmlpath,lapl_id) == 0)
- xml_dtype = H5T.copy('H5T_C_S1');
- H5T.set_size(xml_dtype,'H5T_VARIABLE');
- xml_space_id = H5S.create_simple (1, 1, []);
- xml_id = H5D.create (obj.fid, obj.xmlpath, xml_dtype, xml_space_id, 'H5P_DEFAULT');
- H5D.close(xml_id);
+ error('No XML header found.');
end
H5P.close(lapl_id);
% Open
xml_id = H5D.open(obj.fid, obj.xmlpath);
-
+
+ % Get the type
+ xml_dtype = H5D.get_type(xml_id);
+
% Read the data
hdr = H5D.read(xml_id, xml_dtype, 'H5S_ALL', 'H5S_ALL', 'H5P_DEFAULT');
- obj.xml = hdr{1};
+
+ % Output depends on whether or not the stored string was variale length
+ if (H5T.is_variable_str(xml_dtype))
+ xmlstring = hdr{1};
+ else
+ xmlstring = hdr';
+ end
- % Close the XML
+ % Close the XML
+ H5T.close(xml_dtype);
H5D.close (xml_id);
end
- function obj = writexml(obj,xmlstring)
+ function writexml(obj,xmlstring)
% No validation is performed. You're on your own.
- % Set variable length string type
- xml_dtype = H5T.copy('H5T_C_S1');
- H5T.set_size(xml_dtype,'H5T_VARIABLE');
-
+ % TODO: add error checking on the write and return a status
+
% Check if the XML header exists
% if it does not exist, create it
+ % if it exists, modify the size appropriately
lapl_id=H5P.create('H5P_LINK_ACCESS');
- if (H5L.exists(obj.fid,obj.xmlpath,lapl_id) == 0)
- xml_space_id = H5S.create_simple (1, 1, []);
- xml_id = H5D.create (obj.fid, obj.xmlpath, xml_dtype, xml_space_id, 'H5P_DEFAULT');
- H5S.close(xml_space_id);
- H5D.close(xml_id);
+ if (H5L.exists(obj.fid,obj.xmlpath,lapl_id) == 1)
+ % Delete it
+ H5L.delete(obj.fid, obj.xmlpath,'H5P_DEFAULT');
end
H5P.close(lapl_id);
-
- % Open the data set and get the space
- xml_id = H5D.open(obj.fid, obj.xmlpath);
- xml_space_id = H5D.get_space(xml_id);
- %H5S.select_all(space);
+
+ % Set variable length string type
+ xml_dtype = H5T.copy('H5T_C_S1');
+ % Matlab is having trouble writing variable length strings
+ % that are longer that 512 characters. Switched to fixed
+ % length.
+ % H5T.set_size(xml_dtype,'H5T_VARIABLE');
+ H5T.set_size(xml_dtype, length(xmlstring));
+ xml_space_id = H5S.create_simple (1, 1, []);
+ xml_id = H5D.create (obj.fid, obj.xmlpath, xml_dtype, ....
+ xml_space_id, 'H5P_DEFAULT');
+ H5S.close(xml_space_id);
% Write the data
- H5D.write(xml_id, xml_dtype, 'H5S_ALL', xml_space_id, 'H5P_DEFAULT', xmlstring);
+ H5D.write(xml_id, xml_dtype, ...
+ 'H5S_ALL', 'H5S_ALL', 'H5P_DEFAULT', xmlstring);
% Close the XML
- H5S.close(xml_space_id);
H5D.close(xml_id);
end
+
+ function nacq = getNumberOfAcquisitions(obj)
+
+ dset = H5D.open(obj.fid, obj.datapath);
+ space = H5D.get_space(dset);
+ H5S.get_simple_extent_dims(space);
+ [~,dims,~] = H5S.get_simple_extent_dims(space);
+ nacq = dims(1);
+ H5S.close(space);
+ H5D.close(dset);
+
+ end
+
+ function acq = readAcquisition(obj, nacq)
+
+ % Open the data
+ dset_id = H5D.open(obj.fid, obj.datapath);
+
+ % Open the data space
+ file_space_id = H5D.get_space(dset_id);
+
+ % Initialize the return object
+ acq = ismrmrd.Acquisition;
+
+ % Create a mem_space for reading
+ dims = [1 1];
+ mem_space_id = H5S.create_simple(2,dims,[]);
+
+ % Read the desired acquisition
+ offset = [nacq-1 0];
+ H5S.select_hyperslab(file_space_id,'H5S_SELECT_SET',offset,[1 1],[1 1],[1 1]);
+ d = H5D.read(dset_id, ismrmrd.hdf5_datatypes.getType_Acquisition, ...
+ mem_space_id, file_space_id, 'H5P_DEFAULT');
+
+ % Set the structure bits
+ acq.head = d.head(1);
+ acq.traj = d.traj{1};
+ t = d.data{1};
+ acq.data = t(1:2:end) + 1j*t(2:2:end);
+
+ % Clean up
+ H5S.close(mem_space_id);
+ H5S.close(file_space_id);
+ H5D.close(dset_id);
+ end
+
+ function appendAcquisition(obj, acq)
+ % Append an acquisition
+
+ % TODO: Check the type of the input
+
+ % Check if the Data exists
+ % if it does not exist, create it
+ % if it does exist increase it's size by one
+ lapl_id=H5P.create('H5P_LINK_ACCESS');
+ if (H5L.exists(obj.fid, obj.datapath, lapl_id) == 0)
+ % Data does not exist
+ % create with rank 2, unlimited, and set the chunk size
+ dims = [1 1];
+ maxdims = [H5ML.get_constant_value('H5S_UNLIMITED') 1];
+ file_space_id = H5S.create_simple(2, dims, maxdims);
+
+ dcpl = H5P.create('H5P_DATASET_CREATE');
+ chunk = [1 1];
+ H5P.set_chunk (dcpl, chunk);
+ data_id = H5D.create(obj.fid, obj.datapath, ...
+ ismrmrd.hdf5_datatypes.getType_Acquisition, ...
+ file_space_id, dcpl);
+ H5P.close(dcpl);
+ H5S.close(file_space_id);
+
+ else
+ % Open the data
+ data_id = H5D.open(obj.fid, obj.datapath);
+
+ % Open the data space
+ file_space_id = H5D.get_space(data_id);
+
+ % Get the size, increment by one
+ H5S.get_simple_extent_dims(file_space_id);
+ [~,dims,~] = H5S.get_simple_extent_dims(file_space_id);
+ dims = [dims(1)+1, 1];
+ H5D.set_extent (data_id, dims);
+ H5S.close(file_space_id);
+
+ end
+ H5P.close(lapl_id);
+
+ % Get the file space (room for one more acq)
+ file_space_id = H5D.get_space(data_id);
+ [~,dims,~] = H5S.get_simple_extent_dims(file_space_id);
+
+ % Select the last block
+ offset = [dims(1)-1 0];
+ H5S.select_hyperslab(file_space_id,'H5S_SELECT_SET',offset,[],[],[]);
+
+ % Mem space
+ mem_space_id = H5S.create_simple(2,[1 1],[]);
+
+ % Pack the acquisition into the correct struct for writing
+ swarn = warning('query','MATLAB:structOnObject');
+ warning('off', 'MATLAB:structOnObject')
+ d.head(1) = struct(acq.head);
+ d.head(1).idx = struct(acq.head.idx);
+ warning(swarn.state, 'MATLAB:structOnObject')
+ d.traj{1} = acq.traj;
+ t = zeros(2*length(acq.data),1,'single');
+ t(1:2:end) = real(acq.data);
+ t(2:2:end) = imag(acq.data);
+ d.data{1} = t;
+
+ % Write
+ H5D.write(data_id, ismrmrd.hdf5_datatypes.getType_Acquisition(), ...
+ mem_space_id, file_space_id, 'H5P_DEFAULT', d);
+
+ % Clean up
+ H5S.close(mem_space_id);
+ H5S.close(file_space_id);
+ H5D.close(data_id);
+ end
+
end
end
diff --git a/matlab/testread.m b/matlab/testread.m
deleted file mode 100644
index 08e78bd..0000000
--- a/matlab/testread.m
+++ /dev/null
@@ -1,66 +0,0 @@
-function [xmlhdr data nacq] = testread(filename)
-
- % Open the HDF5 File
- file = H5F.open(filename, 'H5F_ACC_RDONLY', 'H5P_DEFAULT');
-
- % Open the xml
- dset_id = H5D.open(file, '/dataset/xml');
-
- % Get the dataspace
- file_space_id = H5D.get_space(dset_id);
-
- % Set variable length string type
- ismrm.datatypes.xmlhead = H5T.copy('H5T_C_S1');
- H5T.set_size(ismrm.datatypes.xmlhead,'H5T_VARIABLE');
-
- % Read the data
- hdr = H5D.read(dset_id, ismrm.datatypes.xmlhead, 'H5S_ALL', 'H5S_ALL', 'H5P_DEFAULT');
- xmlhdr = hdr{1};
-
- % Close the XML Header
- H5D.close (dset_id);
- H5S.close (file_space_id);
-
- % Open the data
- dset_id = H5D.open(file, '/dataset/data');
-
- % Get the number of acquisitions
- file_space_id = H5D.get_space(dset_id);
- H5S.get_simple_extent_dims(file_space_id);
- [~,dims,~] = H5S.get_simple_extent_dims(file_space_id);
- nacq = dims(1);
-
- % Initialize the return array of structs
- data(nacq) = ismrmrd.Acquisition;
-
- % Loop over the acquisitions and read them in
- dims = [1 1];
- mem_space_id = H5S.create_simple(2,dims,[]);
-
- for p = 1:nacq
-
- offset = [p-1 0];
- H5S.select_hyperslab(file_space_id,'H5S_SELECT_SET',offset,[1 1],[1 1],[1 1]);
- d = H5D.read(dset_id, ismrmrd.hdf5_datatypes.getType_Acquisition, ...
- mem_space_id, file_space_id, 'H5P_DEFAULT');
-
- data(p).head = d.head(1);
- data(p).traj = d.traj{1};
- data(p).data = d.data{1};
-
- end
-
- % Close the data
- H5S.close (file_space_id);
- H5D.close (dset_id);
-
- % Close the file
- H5F.close (file)
-
-end
-
- % Read all of the acquisitions at once
- %data = H5D.read(dset_id, ismrmrd.hdf5_datatypes.getType_Acquisition, 'H5S_ALL','H5S_ALL','H5P_DEFAULT');
-
-
-
\ No newline at end of file
diff --git a/matlab/testrecon.m b/matlab/testrecon.m
deleted file mode 100644
index f3bcfe8..0000000
--- a/matlab/testrecon.m
+++ /dev/null
@@ -1,10 +0,0 @@
-[xmlhdr data] = testread('~/Projects/ismrmrd_data/testdata.h5');
-
-nacq = size(data,2);
-k = zeros([data(1).head.number_of_samples nacq]);
-for p = 1:nacq;
- k(:,p) = data(p).data;
-end
-
-imagesc(abs(fftshift(ifft2(fftshift(k)))))
-axis image; axis off;
\ No newline at end of file
diff --git a/schema/ismrmrd.xsd b/schema/ismrmrd.xsd
index eac214a..386f513 100644
--- a/schema/ismrmrd.xsd
+++ b/schema/ismrmrd.xsd
@@ -7,7 +7,7 @@
<xs:element maxOccurs="1" minOccurs="0" name="subjectInformation" type="subjectInformationType"/>
<xs:element maxOccurs="1" minOccurs="0" name="acquisitionSystemInformation" type="acquisitionSystemInformationType"/>
<xs:element maxOccurs="1" minOccurs="1" name="experimentalConditions" type="experimentalConditionsType"/>
- <xs:element maxOccurs="65535" minOccurs="1" name="encoding">
+ <xs:element maxOccurs="unbounded" minOccurs="1" name="encoding">
<xs:complexType>
<xs:all>
<xs:element maxOccurs="1" minOccurs="1" name="encodedSpace" type="encodingSpaceType"/>
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/ismrmrd.git
More information about the debian-science-commits
mailing list