1function write(obj,filename,create_signature)
2%WRITE Write sequence to file.
3% WRITE(seqObj, filename) Write the sequence data to the given
4% filename using the open file format for MR sequences.
5%
6% Examples:
7% Write the sequence file to the my_sequences directory
8%
9% write(seqObj,'my_sequences/gre.seq')
10%
11% See also read
13if (nargin<3)
14 create_signature=true;
15end
17fid=fopen(filename, 'w');
18assert(fid ~= -1, 'Cannot open file: %s', filename);
19fprintf(fid, '# Pulseq sequence file\n');
20fprintf(fid, '# Created by MATLAB mr toolbox\n\n');
22% we always write files in the default current version, which may be
23% differen to one, loaded (and stored in the seq object)
24[version_major, version_minor, version_revision]=mr.aux.version('output');
25fprintf(fid, '[VERSION]\n');
26fprintf(fid, 'major %s\n', num2str(version_major));
27fprintf(fid, 'minor %s\n', num2str(version_minor));
28fprintf(fid, 'revision %s\n', num2str(version_revision));
29fprintf(fid, '\n');
31% handle RequiredExtensions definition
32if ~isempty(obj.rotationLibrary.keys)
33 RD=obj.getDefinition('RequiredExtensions');
34 if isempty(RD) || isempty(strfind(RD,'ROTATIONS'))
35 RD=mr.aux.strstrip([mr.aux.strstrip(RD) ' ROTATIONS']);
36 obj.setDefinition('RequiredExtensions', RD);
37 end
38end
40if ~isempty(obj.definitions)
41 fprintf(fid, '[DEFINITIONS]\n');
42 keys = obj.definitions.keys;
43 values = obj.definitions.values;
44 for i=1:length(keys)
45 fprintf(fid, '%s ', keys{i});
46 if (ischar(values{i}))
47 fprintf(fid, '%s ', values{i});
48 else
49 fprintf(fid, '%.9g ', values{i});
50 end
51 fprintf(fid, '\n');
52 end
53 fprintf(fid, '\n');
54end
56fprintf(fid, '# Format of blocks:\n');
57fprintf(fid, '# NUM DUR RF GX GY GZ ADC EXT\n');
58fprintf(fid, '[BLOCKS]\n');
59idFormatWidth = length(num2str(length(obj.blockEvents)));
60idFormatStr = ['%' num2str(idFormatWidth) 'd'];
61for i = 1:length(obj.blockEvents)
62 %fprintf(fid,[idFormatStr ' %2d %2d %3d %3d %3d %2d 0\n'],[i obj.blockEvents(i,:)]);
63 %fprintf(fid,[idFormatStr ' %2d %2d %3d %3d %3d %2d 0\n'],[i obj.blockEvents{i}]);
64 bd=obj.blockDurations(i)/obj.blockDurationRaster;
65 bdr=round(bd);
66 assert(abs(bdr-bd)<1e-6); % this may still trigger false alarms for very long delays due to the limited accuracy of the double
67 fprintf(fid,[idFormatStr ' %3d %3d %3d %3d %3d %2d %2d\n'], ...
68 [i bdr obj.blockEvents{i}(2:end)]);
69end
70fprintf(fid, '\n');
72if ~isempty(obj.rfLibrary.keys)
73 fprintf(fid, '# Format of RF events:\n');
74 fprintf(fid, '# id ampl. mag_id phase_id time_shape_id center delay freqPPM phasePPM freq phase use\n');
75 fprintf(fid, '# .. Hz .. .. .. us us ppm rad/MHz Hz rad ..\n');
76 fprintf(fid,['# Field ''use'' is the initial of: \n# ' ...
77 strtrim(cell2mat(cellfun(@(x) [x ' '], mr.getSupportedRfUse(), 'UniformOutput', false))) ...
78 '\n']);
79 fprintf(fid, '[RF]\n');
80 keys = obj.rfLibrary.keys;
81 for k = keys
82 libData1 = obj.rfLibrary.data(k).array(1:4);
83 libData2 = obj.rfLibrary.data(k).array(7:10);
84 center = obj.rfLibrary.data(k).array(5)*1e6; % us
85 delay = round(obj.rfLibrary.data(k).array(6)/obj.rfRasterTime)*obj.rfRasterTime*1e6; % a bit of a hack: round the delay
86 fprintf(fid, '%d %12g %d %d %d %g %g %g %g %g %g %c\n', [k libData1 center delay], libData2, obj.rfLibrary.type(k));
87 end
88 fprintf(fid, '\n');
89end
91arbGradMask = obj.gradLibrary.type == 'g';
92trapGradMask = obj.gradLibrary.type == 't';
94if any(arbGradMask)
95 fprintf(fid, '# Format of arbitrary gradients:\n');
96 fprintf(fid, '# time_shape_id of 0 means default timing (stepping with grad_raster starting at 1/2 of grad_raster)\n');
97 fprintf(fid, '# id amplitude first last amp_shape_id time_shape_id delay\n');
98 fprintf(fid, '# .. Hz/m Hz/m Hz/m .. .. us\n');
99 fprintf(fid, '[GRADIENTS]\n');
100 keys = obj.gradLibrary.keys;
101 for k = keys(arbGradMask)
102 fprintf(fid, '%d %12g %12g %12g %d %d %d\n', ...
103 [k obj.gradLibrary.data(k).array(1:5) ...
104 round(obj.gradLibrary.data(k).array(6)*1e6)]);
105 end
106 fprintf(fid, '\n');
107end
109if any(trapGradMask)
110 fprintf(fid, '# Format of trapezoid gradients:\n');
111 fprintf(fid, '# id amplitude rise flat fall delay\n');
112 fprintf(fid, '# .. Hz/m us us us us\n');
113 fprintf(fid, '[TRAP]\n');
114 keys = obj.gradLibrary.keys;
115 for k = keys(trapGradMask)
116 data = obj.gradLibrary.data(k).array;
117 data(2:end) = round(1e6*data(2:end));
118 fprintf(fid, '%2d %12g %3d %4d %3d %3d\n', [k data]);
119 end
120 fprintf(fid, '\n');
121end
123if ~isempty(obj.adcLibrary.keys)
124 fprintf(fid, '# Format of ADC events:\n');
125 fprintf(fid, '# id num dwell delay freqPPM phasePPM freq phase phase_id\n');
126 fprintf(fid, '# .. .. ns us ppm rad/MHz Hz rad ..\n');
127 fprintf(fid, '[ADC]\n');
128 keys = obj.adcLibrary.keys;
129 for k = keys
130 data = obj.adcLibrary.data(k).array.*[1 1e9 1e6 1 1 1 1 1];
131 fprintf(fid, '%d %d %.0f %.0f %g %g %g %g %d\n', [k data]);
132 end
133 fprintf(fid, '\n');
134end
136%if ~isempty(obj.delayLibrary.keys)
137% fprintf(fid, '# Format of delays:\n');
138% fprintf(fid, '# id delay (us)\n');
139% fprintf(fid, '[DELAYS]\n');
140% keys = obj.delayLibrary.keys;
141% for k = keys
142% fprintf(fid, '%d %d\n', ...
143% [k round(1e6*obj.delayLibrary.data(k).array)]);
144% end
145% fprintf(fid, '\n');
146%end
148if ~isempty(obj.extensionLibrary.keys)
149 fprintf(fid, '# Format of extension lists:\n');
150 fprintf(fid, '# id type ref next_id\n');
151 fprintf(fid, '# next_id of 0 terminates the list\n');
152 fprintf(fid, '# Extension list is followed by extension specifications\n');
153 fprintf(fid, '[EXTENSIONS]\n');
154 keys = obj.extensionLibrary.keys;
155 for k = keys
156 fprintf(fid, '%d %d %d %d\n', ...
157 [k round(obj.extensionLibrary.data(k).array)]);
158 end
159 fprintf(fid, '\n');
160end
162if ~isempty(obj.trigLibrary.keys)
163 fprintf(fid, '# Extension specification for digital output and input triggers:\n');
164 fprintf(fid, '# id type channel delay (us) duration (us)\n');
165% fprintf(fid, 'extension TRIGGERS 1\n'); % fixme: extension ID 1 is hardcoded here for triggers
166 fprintf(fid, ['extension TRIGGERS ',num2str(obj.getExtensionTypeID('TRIGGERS')),'\n']);
168 keys = obj.trigLibrary.keys;
169 for k = keys
170 fprintf(fid, '%d %d %d %d %d\n', ...
171 [k round(obj.trigLibrary.data(k).array.*[1 1 1e6 1e6])]);
172 end
173 fprintf(fid, '\n');
174end
176if ~isempty(obj.labelsetLibrary.keys) || ~isempty(obj.labelincLibrary.keys)
177 lbls=mr.getSupportedLabels();
179 if ~isempty(obj.labelsetLibrary.keys)
180 fprintf(fid, '# Extension specification for setting labels:\n');
181 fprintf(fid, '# id set labelstring\n');
182 tid=obj.getExtensionTypeID('LABELSET');
183 fprintf(fid, ['extension LABELSET ',num2str(tid),'\n']);
184 keys = obj.labelsetLibrary.keys;
185 for k = keys
186 fprintf(fid, '%d %d %s\n', ...
187 k, obj.labelsetLibrary.data(k).array(1),lbls{obj.labelsetLibrary.data(k).array(2)});
188 end
189 fprintf(fid, '\n');
190 end
191 if ~isempty(obj.labelincLibrary.keys)
192 fprintf(fid, '# Extension specification for increasing labels:\n');
193 fprintf(fid, '# id inc labelstring\n');
194 tid=obj.getExtensionTypeID('LABELINC');
195 fprintf(fid, ['extension LABELINC ',num2str(tid),'\n']);
196 lbls=mr.getSupportedLabels();
197 keys = obj.labelincLibrary.keys;
198 for k = keys
199 fprintf(fid, '%d %d %s\n', ...
200 k, obj.labelincLibrary.data(k).array(1),lbls{obj.labelincLibrary.data(k).array(2)});
201 end
202 fprintf(fid, '\n');
203 end
204end
206if ~isempty(obj.softDelayLibrary.keys)
207 fprintf(fid, '# Extension specification for soft delays:\n');
208 fprintf(fid, '# id num offset factor hint\n');
209 fprintf(fid, '# .. .. us .. ..\n');
210 fprintf(fid, ['extension DELAYS ',num2str(obj.getExtensionTypeID('DELAYS')),'\n']);
212 keys = obj.softDelayLibrary.keys;
213 for k = keys
214 fprintf(fid, '%d %d %g %g %s\n', ...
215 k, obj.softDelayLibrary.data(k).array(1), obj.softDelayLibrary.data(k).array(2)*1e6, obj.softDelayLibrary.data(k).array(3), obj.softDelayHints2{obj.softDelayLibrary.data(k).array(4)});
216 end
217 fprintf(fid, '\n');
218end
220if ~isempty(obj.rfShimLibrary.keys)
221 fprintf(fid, '# Extension specification for RF shimming:\n');
222 fprintf(fid, '# id num_chan magn_c1 phase_c1 magn_c2 phase_c2 ...\n');
223 fprintf(fid, ['extension RF_SHIMS ',num2str(obj.getExtensionTypeID('RF_SHIMS')),'\n']);
225 keys = obj.rfShimLibrary.keys;
226 for k = keys
227 fprintf(fid, '%d %d', [k length(obj.rfShimLibrary.data(k).array)/2]);
228 fprintf(fid, ' %g', obj.rfShimLibrary.data(k).array);
229 fprintf(fid, '\n');
230 end
231 fprintf(fid, '\n');
232end
234if ~isempty(obj.rotationLibrary.keys)
235 fprintf(fid, '# Extension specification for rotation events:\n');
236 fprintf(fid, '# id RotQuat0 RotQuatX RotQuatY RotQuatZ\n');
237 fprintf(fid, ['extension ROTATIONS ',num2str(obj.getExtensionTypeID('ROTATIONS')),'\n']);
239 keys = obj.rotationLibrary.keys;
240 for k = keys
241 fprintf(fid, '%d ', k );
242 fprintf(fid, ' %g', obj.rotationLibrary.data(k).array);
243 fprintf(fid, '\n');
244 end
245 fprintf(fid, '\n');
246end
248if ~isempty(obj.shapeLibrary.keys)
249 fprintf(fid, '# Sequence Shapes\n');
250 fprintf(fid, '[SHAPES]\n\n');
251 keys = obj.shapeLibrary.keys;
252 for k = keys
253 shape_dat = obj.shapeLibrary.data(k).array;
254 fprintf(fid, 'shape_id %d\n', k);
255 fprintf(fid, 'num_samples %d\n', shape_dat(1));
256 fprintf(fid, '%.9g\n', shape_dat(2:end));
257 fprintf(fid, '\n');
258 end
259end
261fclose(fid);
263if create_signature
264 % 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)
266 % re-open and read in the file
267 fid=fopen(filename, 'r');
268 buf=fread(fid);
269 fclose(fid);
271 % calculate the digest
272 md5hash=mr.aux.md5(buf);
273 %fprintf('%s\n',md5hash);
275 % store the signature in the object
276 obj.signatureType='md5';
277 obj.signatureFile='text';
278 obj.signatureValue=md5hash;
280 % re-open the file for appending
281 fid=fopen(filename, 'a');
282 fprintf(fid, '\n[SIGNATURE]\n'); % the preceding new line BELONGS to the signature (and needs to be sripped away to recalculate the signature)
283 fprintf(fid, '# This is the hash of the Pulseq file, calculated right before the [SIGNATURE]\n');
284 fprintf(fid, '# section was added. It can be reproduced/verified with md5sum if the file\n');
285 fprintf(fid, '# trimmed to the position right above [SIGNATURE]. The new line character\n');
286 fprintf(fid, '# preceding [SIGNATURE] BELONGS to the signature (and needs to be sripped away\n');
287 fprintf(fid, '# for recalculating/verification)\n');
288 fprintf(fid, 'Type md5\n');
289 fprintf(fid, 'Hash %s\n', md5hash);
290 fclose(fid);
291end
293end