TEAP (Toolbox for Emotion Analysis using Physiological Signals) doc
iread.m
Go to the documentation of this file.
1 function [data,HDR] = iread(HDR,CHAN,StartPos)
2 % IREAD reads image data
3 %
4 % [S,HDR] = iread(HDR)
5 %
6 % See also: fread, SREAD, SWRITE, SCLOSE, SSEEK, SREWIND, STELL, SEOF
7 
8 % This program is free software; you can redistribute it and/or
9 % modify it under the terms of the GNU General Public License
10 % as published by the Free Software Foundation; either version 3
11 % of the License, or (at your option) any later version.
12 
13 % $Id: iread.m 2205 2009-10-27 12:18:15Z schloegl $
14 % (C) 2005,2006,2008 by Alois Schloegl <a.schloegl@ieee.org>
15 % This is part of the BIOSIG-toolbox http://biosig.sf.net/
16 
17 
18 if 0,
19 
20 
21 elseif strcmp(HDR.TYPE,'IMAGE:BMP'),
22  fseek(HDR.FILE.FID,HDR.HeadLen,'bof');
23  nc = ceil((HDR.bits*HDR.IMAGE.Size(1))/32)*4;
24 
25  if (HDR.bits==1)
26  [data,c] = fread(HDR.FILE.FID,[nc,HDR.IMAGE.Size(2)*8],'ubit1');
27  %[data,c] = fread(HDR.FILE.FID,HDR.IMAGE.Size([2,1]),'ubit1');
28 size(data),
29 % data = reshape(data(1:nc,:)',HDR.IMAGE.Size);
30 
31  elseif (HDR.bits==4)
32  palr = [ 0,128, 0,128, 0,128, 0,192,128,255, 0,255, 0,255, 0,255];
33  palg = [ 0, 0,128,128, 0, 0,128,192,128, 0,255,255, 0, 0,255,255];
34  palb = [ 0, 0, 0, 0,128,128,128,192,128, 0, 0, 0,255,255,255,255];
35  tmp = uint8(fread(HDR.FILE.FID,[nc,HDR.IMAGE.Size(2)*2],'ubit4'));
36  data = palr(tmp(1:HDR.IMAGE.Size(1),:)'+1);
37  data(:,:,2) = palg(tmp(1:HDR.IMAGE.Size(1),:)'+1);
38  data(:,:,3) = palb(tmp(1:HDR.IMAGE.Size(1),:)'+1);
39  data = data(HDR.IMAGE.Size(2):-1:1,:,:);
40 
41  elseif (HDR.bits==8)
42  pal = uint8(colormap*256);
43  tmp = fread(HDR.FILE.FID,[nc,HDR.IMAGE.Size(2)],'uint8');
44  data = pal(tmp(1:HDR.IMAGE.Size(1),:)'+1,1);
45  data(:,:,2) = pal(tmp(1:HDR.IMAGE.Size(1),:)'+1,2);
46  data(:,:,3) = pal(tmp(1:HDR.IMAGE.Size(1),:)'+1,3);
47  data = reshape(data,[HDR.IMAGE.Size([2,1]),3]);
48  data = data(HDR.IMAGE.Size(2):-1:1,:,:);
49 
50  elseif (HDR.bits==24)
51  data = uint8(fread(HDR.FILE.FID,[nc,HDR.IMAGE.Size(2)],'uint8'));
52  HDR.BMP.Red = data((1:HDR.IMAGE.Size(1))*3,:)';
53  HDR.BMP.Green = data((1:HDR.IMAGE.Size(1))*3-1,:)';
54  HDR.BMP.Blue = data((1:HDR.IMAGE.Size(1))*3-2,:)';
55  data = HDR.BMP.Red;
56  data(:,:,2) = HDR.BMP.Green;
57  data(:,:,3) = HDR.BMP.Blue;
58  data = data(HDR.IMAGE.Size(2):-1:1,:,:);
59  else
60 
61  end;
62 
63 
64 elseif strcmp(HDR.TYPE,'IMAGE:FITS'),
65  if nargin>1,
66  KK = CHAN;
67  else
68  [tmp, KK] = max(HDR.IMAGE_Size); % select block
69  end;
70 
71  for KK = 1:length(HDR.FITS),
72 
73  status = fseek(HDR.FILE.FID,HDR.HeadLen(KK),'bof');
74 
75  HDR.AS.bps = abs(HDR.FITS{KK}.BITPIX)/8;
76  if HDR.FITS{KK}.BITPIX == 8,
77  HDR.GDFTYP = 'uint8';
78  elseif HDR.FITS{KK}.BITPIX == 16,
79  HDR.GDFTYP = 'int16';
80  elseif HDR.FITS{KK}.BITPIX == 32,
81  HDR.GDFTYP = 'int32';
82  elseif HDR.FITS{KK}.BITPIX == 64,
83  HDR.GDFTYP = 'int64';
84  elseif HDR.FITS{KK}.BITPIX == -32,
85  HDR.GDFTYP = 'float32';
86  elseif HDR.FITS{KK}.BITPIX == -64,
87  HDR.GDFTYP = 'float64';
88  else
89  warning('IREAD (FITS{KK})');
90  end;
91 
92  if isfield(HDR.FITS{KK},'BZERO') HDR.Off = HDR.FITS{KK}.BZERO;
93  else HDR.Off = 0; end;
94  if isfield(HDR.FITS{KK},'BSCALE') HDR.Cal = HDR.FITS{KK}.BSCALE;
95  else HDR.Cal = 1; end;
96  if isfield(HDR.FITS{KK},'BUNIT'), HDR.PhysDim = HDR.FITS{KK}.BUNIT;
97  else HDR.PhysDim = '[1]'; end;
98  if isfield(HDR.FITS{KK},'DATAMAX'), HDR.PhysMax = HDR.FITS{KK}.DATAMAX;
99  else HDR.PhysMax = NaN; end;
100  if isfield(HDR.FITS{KK},'DATAMIN'), HDR.PhysMin = HDR.FITS{KK}.DATAMIN;
101  else HDR.PhysMin = NaN; end;
102 
103  if ~isfield(HDR.FITS{KK},'XTENSION')
104  [data,c] = fread(HDR.FILE.FID, prod(HDR.IMAGE(KK).Size), HDR.GDFTYP);
105  data = reshape(data,HDR.IMAGE(KK).Size) * HDR.Cal + HDR.Off;
106 
107  elseif strncmp(HDR.FITS{KK}.XTENSION,'TABLE',5)
108  for k = 1:HDR.FITS{KK}.TFIELDS,
109  f = ['TTYPE',int2str(k)];
110  if isfield(HDR.FITS{KK},f)
111  HDR.TABLE{KK}.Label(k,:) = getfield(HDR.FITS{KK},f);
112  end;
113  tmp = getfield(HDR.FITS{KK},['TFORM',int2str(k)]);
114  typ(k) = tmp(1);
115  ix(k) = getfield(HDR.FITS{KK},['TBCOL',int2str(k)]);
116  end;
117  ix(k+1) = HDR.FITS{KK}.NAXIS1+1;
118 
119  sa = {};
120  [data,c] = fread(HDR.FILE.FID, HDR.IMAGE(KK).Size, HDR.GDFTYP);
121  data = data';
122  s = repmat(NaN, HDR.FITS{KK}.NAXIS2, HDR.FITS{KK}.TFIELDS);
123  for k = 1:HDR.FITS{KK}.TFIELDS,
124  [s(:,k), status,sa(:,k)] = str2double(char(data(:,ix(k):ix(k+1)-1)));
125  end;
126  HDR.TABLE{KK} = s;
127 
128  elseif strncmp(HDR.FITS{KK}.XTENSION,'BINTABLE',8)
129  sz = ones(1,HDR.FITS{KK}.TFIELDS);
130  for k = 1:HDR.FITS{KK}.TFIELDS,
131  f = ['TTYPE',int2str(k)];
132  if isfield(HDR.FITS{KK},f)
133  HDR.TABLE{KK}.Label(k,:) = getfield(HDR.FITS{KK},f);
134  end;
135  tmp = getfield(HDR.FITS{KK},['TFORM',int2str(k)]);
136  ix = min(find(tmp>'9'));
137  HDR.FITS{KK}.TYP(k) = tmp(ix);
138  tmp1 = str2double(tmp(1:ix-1));
139  if ~isempty(tmp1), sz(k)=tmp1; end;
140 
141  if 0,
142  elseif tmp(ix)=='L', GDFTYP{k} = 'uint8'; % char T, F
143  elseif tmp(ix)=='X', GDFTYP{k} = 'ubit1'; %ubit1
144  elseif tmp(ix)=='B', GDFTYP{k} = 'uint8'; % uint8
145  elseif tmp(ix)=='I', GDFTYP{k} = 'int16'; % int16
146  elseif tmp(ix)=='J', GDFTYP{k} = 'int32'; % int32
147  elseif tmp(ix)=='A', GDFTYP{k} = 'uchar'; % char
148  elseif tmp(ix)=='E', GDFTYP{k} = 'float32'; % float32
149  elseif tmp(ix)=='D', GDFTYP{k} = 'float64'; % double
150  elseif tmp(ix)=='C', sz(k) = sz(k)*2; GDFTYP{k} = 'float32'; % [1+i]*float32
151  elseif tmp(ix)=='M', sz(k) = sz(k)*2; GDFTYP{k} = 'float64'; % [1+i]*double
152  elseif tmp(ix)=='P', GDFTYP{k} = 'uint64'; % uint64, array desc ?
153  else
154  end;
155  end;
156 
157  for k2 = 1:HDR.FITS{KK}.NAXIS2,
158  for k1 = 1:HDR.FITS{KK}.TFIELDS,
159  [s,c] = fread(HDR.FILE.FID, [1,sz(k1)], GDFTYP{k1});
160 
161  if 0,
162  elseif HDR.FITS{KK}.TYP(k1)=='L',
163  sig{k1}(k2,:) = (s=='T');
164  elseif HDR.FITS{KK}.TYP(k1)=='A',
165  sig{k1}(k2,:) = char(s);
166  elseif any(HDR.FITS{KK}.TYP(k1)=='CM'),
167  sig{k1}(k2,:) = [1,i]*reshape(s,2,sz(k1)/2);
168  else
169  sig{k1}(k2,:) = s;
170  end;
171  end;
172  end;
173  HDR.TABLE{KK} = sig;
174 
175  elseif strncmp(HDR.FITS{KK}.XTENSION,'IMAGE',5)
176  [data,c] = fread(HDR.FILE.FID, prod(HDR.IMAGE(KK).Size), HDR.GDFTYP);
177 
178  [c,size(data),HDR.IMAGE(KK).Size]
179  data = reshape(data,HDR.IMAGE(KK).Size) * HDR.Cal + HDR.Off;
180 
181  else
182  [data,c] = fread(HDR.FILE.FID, prod(HDR.IMAGE(KK).Size), HDR.GDFTYP);
183  data = reshape(data,HDR.IMAGE(KK).Size) * HDR.Cal + HDR.Off;
184 
185  end;
186  end;
187 
188 
189 elseif strcmp(HDR.TYPE,'IMAGE:HGT'),
190  data = fread(HDR.FILE.FID,HDR.IMAGE.Size,'int16')';
191 
192 
193 elseif strcmp(HDR.TYPE,'IMAGE:TIFF'),
194 
195  if ~isfield(HDR.TIFF,'StripOffset')
196  if 1,
197  data = []; count = 0;
198  for k = 1:length(HDR.TIFF.TileOffset)
199  status = fseek(HDR.FILE.FID,HDR.TIFF.TileOffset(k),'bof');
200  [d,c] = fread(HDR.FILE.FID,HDR.TIFF.TileByteCount(k),HDR.GDFTYP);
201  data = [data; d]; count = count + c;
202  end;
203  else
204 
205  fprintf(HDR.FILE.stderr,'Warning IREAD (TIFF): TileOffset not supported yet.\n');
206  data = repmat(NaN,HDR.IMAGE.Size);
207  count = 0;
208  nr = HDR.IMAGE.Size(1)./HDR.TIFF.TileLength;
209  nc = HDR.IMAGE.Size(2)./HDR.TIFF.TileWidth;
210  for k1 = 1:nr
211  for k2 = 1:nc,
212  status = fseek(HDR.FILE.FID,HDR.TIFF.TileOffset(k1*nc-nc+k2),'bof');
213  [d,c] = fread(HDR.FILE.FID,HDR.TIFF.TileByteCount(k1*nc-nc+k2),HDR.GDFTYP);
214  d = permute(reshape(d,[HDR.TIFF.SamplesPerPixel,HDR.TIFF.TileWidth,HDR.TIFF.TileLength]),[3,2,1]);
215  [size(d),size(data)],
216  data ((k1-1)*HDR.TIFF.TileLength+1:k1*HDR.TIFF.TileLength,(k2-1)*HDR.TIFF.TileWidth+1:k2*HDR.TIFF.TileWidth,:) = d;
217  count = count + c;
218  end;
219  end;
220  end;
221 
222 
223  elseif ~any(HDR.TIFF.StripOffset(1:end-1)-HDR.TIFF.StripOffset(2:end)+HDR.TIFF.StripByteCounts(1:end-1))
224  status = fseek(HDR.FILE.FID,HDR.TIFF.StripOffset(1),'bof');
225  [data,count] = fread(HDR.FILE.FID,sum(HDR.TIFF.StripByteCounts)*length(HDR.Bits),HDR.GDFTYP);
226  else
227  data = []; count = 0;
228  for k = 1:length(HDR.TIFF.StripOffset)
229  status = fseek(HDR.FILE.FID,HDR.TIFF.StripOffset(k),'bof');
230  [d,c] = fread(HDR.FILE.FID,HDR.TIFF.StripByteCounts(k),HDR.GDFTYP);
231  data = [data; d]; count = count + c;
232  end;
233  end;
234 
235 
236  if HDR.TIFF.Compression==1,
237 
238  elseif HDR.TIFF.Compression==32773, % Unpack PackBits
239  D = data;
240  N = prod(HDR.IMAGE.Size);
241  data = zeros(N,1);
242  k1 = 1; k2 = 0;
243  while k2 < N; %length(D),
244  if (D(k1)<128),
245  data(k2+(1:D(k1)+1)) = D(k1+(1:D(k1)+1));
246  k2 = k2 + D(k1) + 1;
247  k1 = k1 + D(k1) + 2;
248  elseif (D(k1)==128),
249  k1 = k1 + 1;
250  elseif (D(k1)>128),
251  n = 1 - (D(k1) - 256);
252  data(k2+(1:n)) = D(k1+1);
253  k1 = k1 + 2;
254  k2 = k2 + n;
255  else
256  error('TIFF decompression error');
257  end;
258  end;
259 
260  if (k2~=length(data)) | (k1-1-length(D))
261  [k2,length(data),k1,length(D)],
262  warning('TIFF decompression not completed');
263  data = data(1:N);
264  end;
265  else
266  fprintf(HDR.FILE.stderr,'Error IREAD (TIFF): decompression Mode=%i not supported.\n',HDR.TIFF.Compression);
267  end;
268 
269  if prod(HDR.IMAGE.Size)<=length(data);
270  tmp = length(data)-prod(HDR.IMAGE.Size);
271  if tmp,
272  fprintf(HDR.FILE.stderr,'Warning IREAD (TIFF): %i bytes to much.\n',tmp);
273  end;
274 
275  data = permute(reshape(data(1:prod(HDR.IMAGE.Size)),HDR.IMAGE.Size([3,2,1])),[3,2,1]);
276  if HDR.IMAGE.Size(3)==1,
277  data= squeeze(data);
278  end;
279  else
280  fprintf(HDR.FILE.stderr,'Warning IREAD (TIFF): size of data does not fit.\n');
281  end;
282 
283  if isfield(HDR.TIFF,'ColorMap');
284  HDR.IMAGE.ColorMap = HDR.TIFF.ColorMap*2^-16;
285 
286  elseif isfield(HDR.TIFF,'SamplesPerPixel');
287  if HDR.TIFF.SamplesPerPixel==1,
288  data = data*2^(-HDR.Bits);
289  end;
290  end;
291 
292 
293 end;
294 
295 
iread
function iread(in HDR, in CHAN, in StartPos)