/ concept-collection / seqlab
Sign in
concept-collection / seqlab
seqlab / src / engine / pulseq / +mr / @Sequence / writeBinary.m
288 lines · 10.9 KBBlameHistoryRaw
1function writeBinary(obj,filename,create_signature)
2%WRITEBINARY Write sequence to file in binary format.
3% WRITEBINARY(seqObj, filename) Write the sequence data to the given
4% filename using the binary version of the Pulseq open file format for MR
5% sequences. The file specification is available at
6% http://pulseq.github.io
7%
8% Examples:
9% Write the sequence file to the sequences directory
11% writeBinary(seqObj,'sequences/gre.bseq')
13% See also readBinary
15if (nargin<3)
16 create_signature=true;
17end
19% handle RequiredExtensions definition (same as write())
20if ~isempty(obj.rotationLibrary.keys)
21 RD=obj.getDefinition('RequiredExtensions');
22 if isempty(RD) || isempty(strfind(RD,'ROTATIONS'))
23 RD=mr.aux.strstrip([mr.aux.strstrip(RD) ' ROTATIONS']);
24 obj.setDefinition('RequiredExtensions', RD);
25 end
26end
28binaryCodes = obj.getBinaryCodes();
29fid=fopen(filename, 'w');
30fwrite(fid, binaryCodes.fileHeader, 'int64');
31fwrite(fid, int64(obj.version_major), 'int64');
32fwrite(fid, int64(obj.version_minor), 'int64');
33fwrite(fid, int64(obj.version_revision), 'int64');
35if ~isempty(obj.definitions)
36 fwrite(fid, binaryCodes.section.definitions, 'int64');
37 keys = obj.definitions.keys;
38 values = obj.definitions.values;
39 fwrite(fid, length(keys), 'int64');
40 for i = 1:length(keys)
41 fwrite(fid, length(keys{i}),'int32');
42 fwrite(fid, keys{i},'char');
43 val = values{i};
44 fwrite(fid, length(val), 'int32');
45 if ischar(val)
46 fwrite(fid, 'c', 'char');
47 fwrite(fid, val, 'char');
48 elseif isinteger(val)
49 fwrite(fid, 'i', 'char');
50 fwrite(fid, val, 'int32');
51 elseif isfloat(val)
52 fwrite(fid, 'f', 'char');
53 fwrite(fid, val, 'float64');
54 else
55 error(['unknown type of the value type for ' keys{i} ]);
56 end
57 end
58end
60% Blocks: write count, then per block: duration (int64, in blockDurationRaster units)
61% followed by 6 event IDs (int32): rf, gx, gy, gz, adc, ext
62fwrite(fid, binaryCodes.section.blocks, 'int64');
63fwrite(fid, length(obj.blockEvents), 'int64');
64for i = 1:length(obj.blockEvents)
65 bd = obj.blockDurations(i) / obj.blockDurationRaster;
66 bdr = round(bd);
67 assert(abs(bdr - bd) < 1e-6);
68 fwrite(fid, bdr, 'int64'); % block duration in raster units
69 fwrite(fid, obj.blockEvents{i}(2:end), 'int32'); % rf, gx, gy, gz, adc, ext
70end
72% RF: amp(f64) mag_id(i32) phase_id(i32) time_shape_id(i32) center(i64,us)
73% delay(i64,us) freqPPM(f64) phasePPM(f64) freq(f64) phase(f64) use(char)
74% array layout: [amp mag_id phase_id time_shape_id center delay freqPPM phasePPM freq phase]
75if ~isempty(obj.rfLibrary.keys)
76 keys = obj.rfLibrary.keys;
77 fwrite(fid, binaryCodes.section.rf, 'int64');
78 fwrite(fid, length(keys), 'int64');
79 for k = keys
80 data = obj.rfLibrary.data(k).array;
81 fwrite(fid, k, 'int32');
82 fwrite(fid, data(1), 'float64'); % amp
83 fwrite(fid, data(2:4), 'int32'); % mag_id, phase_id, time_shape_id
84 fwrite(fid, round(data(5)*1e12), 'int64'); % center (ps)
85 fwrite(fid, round(data(6)*1e12), 'int64'); % delay (ps)
86 fwrite(fid, data(7:10), 'float64'); % freqPPM, phasePPM, freq, phase
87 fwrite(fid, obj.rfLibrary.type(k), 'char'); % use
88 end
89end
91arbGradMask = obj.gradLibrary.type=='g';
92trapGradMask = obj.gradLibrary.type=='t';
94% Arbitrary gradients: amp(f64) first(f64) last(f64) amp_shape_id(i32)
95% time_shape_id(i32) delay(i32,us)
96% array layout: [amp first last amp_shape_id time_shape_id delay]
97if any(arbGradMask)
98 keys = obj.gradLibrary.keys;
99 fwrite(fid, binaryCodes.section.gradients, 'int64');
100 fwrite(fid, length(keys(arbGradMask)), 'int64');
101 for k = keys(arbGradMask)
102 data = obj.gradLibrary.data(k).array;
103 fwrite(fid, k, 'int32');
104 fwrite(fid, data(1:3), 'float64'); % amp, first, last
105 fwrite(fid, data(4:5), 'int32'); % amp_shape_id, time_shape_id
106 fwrite(fid, round(data(6)*1e12), 'int64'); % delay (ps)
107 end
108end
110% Trapezoid gradients: amp(f64) rise(i64,us) flat(i64,us) fall(i64,us) delay(i64,us)
111% array layout: [amp rise flat fall delay]
112if any(trapGradMask)
113 keys = obj.gradLibrary.keys;
114 fwrite(fid, binaryCodes.section.trapezoids, 'int64');
115 fwrite(fid, length(keys(trapGradMask)), 'int64');
116 for k = keys(trapGradMask)
117 data = obj.gradLibrary.data(k).array;
118 fwrite(fid, k, 'int32');
119 fwrite(fid, data(1), 'float64'); % amp
120 fwrite(fid, data(2:5)*1e12, 'int64'); % rise, flat, fall, delay (ps)
121 end
122end
124% ADC: num(i64) dwell(i64,ns) delay(i64,us) freqPPM(f64) phasePPM(f64)
125% freq(f64) phase(f64) phase_id(i32)
126% array layout: [num dwell delay freqPPM phasePPM freq phase phase_id]
127if ~isempty(obj.adcLibrary.keys)
128 keys = obj.adcLibrary.keys;
129 fwrite(fid, binaryCodes.section.adc, 'int64');
130 fwrite(fid, length(keys), 'int64');
131 for k = keys
132 data = obj.adcLibrary.data(k).array;
133 fwrite(fid, k, 'int32');
134 fwrite(fid, data(1), 'int64'); % num
135 fwrite(fid, round(data(2)*1e12), 'int64'); % dwell (ps)
136 fwrite(fid, round(data(3)*1e12), 'int64'); % delay (ps)
137 fwrite(fid, data(4:7), 'float64'); % freqPPM, phasePPM, freq, phase
138 fwrite(fid, data(8), 'int32'); % phase_id
139 end
140end
142if ~isempty(obj.shapeLibrary.keys)
143 keys = obj.shapeLibrary.keys;
144 fwrite(fid, binaryCodes.section.shapes, 'int64');
145 fwrite(fid, length(keys), 'int64');
146 for k = keys
147 shape = obj.shapeLibrary.data(k).array;
148 num_samples = shape(1);
149 data = shape(2:end);
150 fwrite(fid, k, 'int32');
151 fwrite(fid, num_samples, 'int64'); % num uncompressed
152 fwrite(fid, length(data), 'int64'); % num compressed
153 fwrite(fid, data, 'float32');
154 end
155end
157% Extensions: id(i32) type(i32) ref(i32) next_id(i32)
158if ~isempty(obj.extensionLibrary.keys)
159 keys = obj.extensionLibrary.keys;
160 fwrite(fid, binaryCodes.section.extensions, 'int64');
161 fwrite(fid, length(keys), 'int64');
162 for k = keys
163 fwrite(fid, k, 'int32');
164 fwrite(fid, round(obj.extensionLibrary.data(k).array), 'int32'); % type, ref, next_id
165 end
166end
168% Triggers: id(i32) type(i32) channel(i32) delay(i64,ps) duration(i64,ps)
169if ~isempty(obj.trigLibrary.keys)
170 keys = obj.trigLibrary.keys;
171 fwrite(fid, binaryCodes.section.triggers, 'int64');
172 fwrite(fid, obj.getExtensionTypeID('TRIGGERS'), 'int32'); % extension type ID
173 fwrite(fid, length(keys), 'int64');
174 for k = keys
175 data = obj.trigLibrary.data(k).array;
176 fwrite(fid, k, 'int32');
177 fwrite(fid, data(1:2), 'int32'); % type, channel
178 fwrite(fid, round(data(3:4)*1e12), 'int64'); % delay, duration (ps)
179 end
180end
182% Labels (LABELSET and LABELINC): id(i32) value(i32) label_index(i32)
183if ~isempty(obj.labelsetLibrary.keys)
184 keys = obj.labelsetLibrary.keys;
185 fwrite(fid, binaryCodes.section.labelset, 'int64');
186 fwrite(fid, obj.getExtensionTypeID('LABELSET'), 'int32'); % extension type ID
187 fwrite(fid, length(keys), 'int64');
188 for k = keys
189 data = obj.labelsetLibrary.data(k).array;
190 fwrite(fid, k, 'int32');
191 fwrite(fid, data(1:2), 'int32'); % value, label_index
192 end
193end
195if ~isempty(obj.labelincLibrary.keys)
196 keys = obj.labelincLibrary.keys;
197 fwrite(fid, binaryCodes.section.labelinc, 'int64');
198 fwrite(fid, obj.getExtensionTypeID('LABELINC'), 'int32'); % extension type ID
199 fwrite(fid, length(keys), 'int64');
200 for k = keys
201 data = obj.labelincLibrary.data(k).array;
202 fwrite(fid, k, 'int32');
203 fwrite(fid, data(1:2), 'int32'); % value, label_index
204 end
205end
207% Soft delays: id(i32) num(i32) offset(i64,ps) factor(f64) hint_len(i32) hint(chars)
208if ~isempty(obj.softDelayLibrary.keys)
209 keys = obj.softDelayLibrary.keys;
210 fwrite(fid, binaryCodes.section.softdelays, 'int64');
211 fwrite(fid, obj.getExtensionTypeID('DELAYS'), 'int32'); % extension type ID
212 fwrite(fid, length(keys), 'int64');
213 for k = keys
214 data = obj.softDelayLibrary.data(k).array;
215 hint_str = obj.softDelayHints2{data(4)};
216 fwrite(fid, k, 'int32');
217 fwrite(fid, data(1), 'int32'); % num
218 fwrite(fid, round(data(2)*1e12), 'int64'); % offset (ps)
219 fwrite(fid, data(3), 'float64'); % factor
220 fwrite(fid, length(hint_str), 'int32'); % hint string length
221 fwrite(fid, hint_str, 'char'); % hint string (no null terminator)
222 end
223end
225% RF shims: id(i32) num_chan(i32) mag_c1(f64) phase_c1(f64) ...
226if ~isempty(obj.rfShimLibrary.keys)
227 keys = obj.rfShimLibrary.keys;
228 fwrite(fid, binaryCodes.section.rfshims, 'int64');
229 fwrite(fid, obj.getExtensionTypeID('RF_SHIMS'), 'int32'); % extension type ID
230 fwrite(fid, length(keys), 'int64');
231 for k = keys
232 chan_data = obj.rfShimLibrary.data(k).array;
233 fwrite(fid, k, 'int32');
234 fwrite(fid, length(chan_data)/2, 'int32'); % num channels
235 fwrite(fid, chan_data, 'float64'); % mag/phase pairs
236 end
237end
239% Rotations: id(i32) q0(f64) qx(f64) qy(f64) qz(f64)
240if ~isempty(obj.rotationLibrary.keys)
241 keys = obj.rotationLibrary.keys;
242 fwrite(fid, binaryCodes.section.rotations, 'int64');
243 fwrite(fid, obj.getExtensionTypeID('ROTATIONS'), 'int32'); % extension type ID
244 fwrite(fid, length(keys), 'int64');
245 for k = keys
246 fwrite(fid, k, 'int32');
247 fwrite(fid, obj.rotationLibrary.data(k).array, 'float64'); % quaternion [q0 qx qy qz]
248 end
249end
251fclose(fid);
253if create_signature
254 % sign the file (this version with re/loading the file is a factor 2 faster than the sprintf() based one that kept a memory-copy of the data written)
256 % re-open and read in the file
257 fid=fopen(filename, 'r');
258 buf=fread(fid);
259 fclose(fid);
261 % calculate the digest
262 md5hash=mr.aux.md5(buf);
263 %fprintf('%s\n',md5hash);
265 % store the signature in the seq object
266 obj.signatureType='md5';
267 obj.signatureFile='bin';
268 obj.signatureValue=md5hash;
270 % re-open the file for appending
271 fid=fopen(filename, 'a');
272 fseek(fid, 0, 'eof'); % Octave seems to need this inspite of 'a'
273 fpos=ftell(fid);
274 fwrite(fid, binaryCodes.section.signature, 'int64');
275 % signature type: length,string
276 fwrite(fid, length(obj.signatureType), 'int32');
277 fwrite(fid, obj.signatureType, 'char');
278 % signature: length,data (as bytes, not characters)
279 fwrite(fid, length(obj.signatureValue)/2, 'int32');
280 for i=1:length(obj.signatureValue)/2
281 fwrite(fid, hex2dec(obj.signatureValue(i*2-1:i*2)), 'uint8');
282 end
283 % the original length of the file prior to adding the signature for easier signature validation : int64
284 fwrite(fid,fpos,'int64');
285 fclose(fid);
286end
288end
moveopenescclose