Creating tiled TIFF files with IDL code

The code example shown at the bottom, writes a TIFF image using tiled image format instead of the normal TIFF image format. Standard TIFF images use strips that span all the way across the X dimension, and is only split in the Y dimension. According to the TIFF specification document, section 15:

For low-resolution to medium-resolution images, the standard TIFF method of breaking the image into strips is adequate. However high-resolution images can be accessed more efficiently—and compression tends to work better—if the image is broken into roughly square tiles instead of horizontally-wide but vertically narrow strips.

 

This is a quick test to check the data type support:

IDL> data = double(dist(800,600))

IDL> tiff6_write, 'test.tif', data & help, read_tiff('test.tif')

<Expression>    DOUBLE   = Array[800, 600]

IDL> tiff6_write, 'test.tif', float(data) & help, read_tiff('test.tif')

<Expression>    FLOAT    = Array[800, 600]

IDL> tiff6_write, 'test.tif', ulong(data) & help, read_tiff('test.tif')

<Expression>    ULONG    = Array[800, 600]

IDL> tiff6_write, 'test.tif', long(data) & help, read_tiff('test.tif')

<Expression>    LONG     = Array[800, 600]

IDL> tiff6_write, 'test.tif', uint(data) & help, read_tiff('test.tif')

<Expression>    UINT     = Array[800, 600]

IDL> tiff6_write, 'test.tif', fix(data) & help, read_tiff('test.tif')

<Expression>    INT      = Array[800, 600]

IDL> tiff6_write, 'test.tif', byte(data) & help, read_tiff('test.tif')

<Expression>   BYTE      = Array[800, 600]

 

This is the code used to produce the tiled TIFF files:

 

function tiff6_Ifd, tag, type, count,valueOffset

 compile_opt idl2,logical_predicate

 typelen = [1,1,1,2,4,8,1,1,2,4,8,4,8]

 ifd = { tiff6ifd, tag: 0us, type: 0us, count: 0ul, valueOffset: 0ul }

 ifd.tag = tag

 ifd.type = type

 ifd.count = count

 len = typelen[type]*count

 if len le 4 then begin

   case type of

     3: ifd.valueOffset= ulong([uint(valueOffset), 0us],0)

     4: ifd.valueOffset= ulong(valueOffset)

     8: ifd.valueOffset= ulong([fix(valueOffset), 0s],0)

     9: ifd.valueOffset= ulong(long(valueOffset),0)

     11: ifd.valueOffset= ulong(float(valueOffset),0)

     else: if max(type eq [1,2,6,7]) then begin

          ifd.valueOffset = ulong([byte(valueOffset), 0b,0b,0b],0)

       endif else message, 'Unsupported type: '+strtrim(type,2)

   endcase

 endif else begin

   ifd.valueOffset = ulong(valueOffset)

   if arg_present(valueOffset) then begin

     valueOffset += len

   endif

 endelse

 return, ifd

end

 

pro tiff6_write, filename, data

 compile_opt idl2,logical_predicate

 ;header

 byteord = byte(1us,0) ? 'II' : 'MM'

 header = { byteord: byteord, magic: 42us, ifdOffset: 8ul }

 numIfd = 18us

 

 ;determine length and type

 sz = size(data, /struct)

 sz.dimensions >= 1

 h = hash('BYTE',[1,1,1],'INT',[8,2,2],'UINT',[3,2,1],'LONG',[9,4,2], $

   'ULONG',[4,4,1],'FLOAT',[11,4,3],'DOUBLE',[12,8,3])

 if h.HasKey(sz.type_name)then begin

   info = (orderedhash(list('type','len','frmt'), $

     list(h[sz.type_name],/extract))).ToStruct()

 endif else message, 'Unsupported type: '+sz.type_name

 

 ;strings need to be 0 terminated and even length

 software = [byte('IDL 8.3'),0b]

 datetime = [byte(string(systime(/julian), format= $

   '(c(cyi4.4,":",cmoi2.2,":",cdi2.2,"",cHi2.2,":",cmi2.2,":",csi2.2))')),0b]

 

 ;compute tile info

 tileSize = 128

 tileByteSize = info.len*tileSize^2

 n = sz.dimensions/tileSize + (sz.dimensions mod tileSize ne 0)

 ntiles = n[0]*n[1]

 

 ptr = 8 + 2 + numIfd*12 + 4

 

 ;fill in fields

 ifds = list()

 ifds.Add, tiff6_ifd(256, 4, 1, sz.dimensions[0]) ; width

 ifds.Add, tiff6_ifd(257, 4, 1, sz.dimensions[1]) ; length

 ifds.Add, tiff6_ifd(258, 3, 1, info.len*8) ; bits

 ifds.Add, tiff6_ifd(259, 3, 1, 1) ; compression

 ifds.Add, tiff6_ifd(262, 3, 1, 1) ; photometric interpretation

 ifds.Add, tiff6_ifd(274, 3, 1, 1) ; orientation

 ifds.Add, tiff6_ifd(277, 3, 1, 1) ; samples per pixel

 ifds.Add, tiff6_ifd(282, 5, 1, ptr) ; xresolution

 ifds.Add, tiff6_ifd(283, 5, 1, ptr) ; yresolution

 ifds.Add, tiff6_ifd(284, 3, 1, 1) ; planar configuration

 ifds.Add, tiff6_ifd(296, 3, 1, 2) ; resolution unit

 ifds.Add, tiff6_ifd(305, 2, n_elements(software), ptr) ; software

 ifds.Add, tiff6_ifd(306, 2, n_elements(datetime), ptr) ; datetime

 ifds.Add, tiff6_ifd(322, 4, 1, tileSize) ; tile width

 ifds.Add, tiff6_ifd(323, 4, 1, tileSize) ; tile length

 ifds.Add, tiff6_ifd(324, 4, ntiles, ptr) ; tile offsets

 ifds.Add, tiff6_ifd(325, 4, ntiles, ptr) ; tile byte counts

 ifds.Add, tiff6_ifd(339, 3, 1, info.frmt) ; sample format

 

 ;fields that require more data

 aux = list()

 aux.Add, ulong([100,1])

 aux.Add, ulong([100,1])

 aux.Add, software

 aux.Add, datetime

 aux.Add, ulindgen(ntiles, start=ptr, increment=tileByteSize)

 tmp = replicate(tileByteSize, ntiles)

 aux.Add, tmp

 

 if numifd ne n_elements(ifds) then message, 'IFD size mismatch'

 ;write to file

 openw, lun, filename,/get_lun

 writeu, lun, header

 writeu, lun, numIfd

 writeu, lun, ifds.ToArray()

 writeu, lun, 0ul ; next IFD, formulti image

 foreach x, aux do writeu, lun, x

 tmp = make_array([tileSize,tileSize], type=size(data,/type))

 for j=0, n[1]-1 do begin

   ys = j*tileSize

   ye = (j*tileSize + tileSize-1) < (sz.dimensions[1]-1)

   for i=0, n[0]-1 do begin

     xs = i*tileSize

     xe = (i*tileSize +tileSize-1) < (sz.dimensions[0]-1)

     tmp[0,0] = data[xs:xe,ys:ye]

     writeu, lun, tmp

   endfor

 endfor

 free_lun, lun

end