|
| 1 | +import numpy as np |
| 2 | +from ..wrapstruct import LabeledWrapStruct |
| 3 | +from ..volumeutils import Recoder |
| 4 | + |
| 5 | +transform_codes = Recoder(( |
| 6 | + (0, 'LINEAR_VOX_TO_VOX'), |
| 7 | + (1, 'LINEAR_RAS_TO_RAS'), |
| 8 | + (2, 'LINEAR_PHYSVOX_TO_PHYSVOX'), |
| 9 | + (14, 'REGISTER_DAT'), |
| 10 | + (21, 'LINEAR_COR_TO_COR')), |
| 11 | + fields=('code', 'label')) |
| 12 | + |
| 13 | + |
| 14 | +class StringBasedStruct(LabeledWrapStruct): |
| 15 | + def __init__(self, |
| 16 | + binaryblock=None, |
| 17 | + endianness=None, |
| 18 | + check=True): |
| 19 | + if binaryblock is not None and getattr(binaryblock, 'dtype', |
| 20 | + None) == self.dtype: |
| 21 | + self._structarr = binaryblock.copy() |
| 22 | + return |
| 23 | + super(StringBasedStruct, self).__init__(binaryblock, endianness, check) |
| 24 | + |
| 25 | + def __array__(self): |
| 26 | + return self._structarr |
| 27 | + |
| 28 | + |
| 29 | +class VolumeGeometry(StringBasedStruct): |
| 30 | + template_dtype = np.dtype([ |
| 31 | + ('valid', 'i4'), # Valid values: 0, 1 |
| 32 | + ('volume', 'i4', (3, 1)), # width, height, depth |
| 33 | + ('voxelsize', 'f4', (3, 1)), # xsize, ysize, zsize |
| 34 | + ('xras', 'f4', (3, 1)), # x_r, x_a, x_s |
| 35 | + ('yras', 'f4', (3, 1)), # y_r, y_a, y_s |
| 36 | + ('zras', 'f4', (3, 1)), # z_r, z_a, z_s |
| 37 | + ('cras', 'f4', (3, 1)), # c_r, c_a, c_s |
| 38 | + ('filename', 'U1024')]) # Not conformant (may be >1024 bytes) |
| 39 | + dtype = template_dtype |
| 40 | + |
| 41 | + def as_affine(self): |
| 42 | + affine = np.eye(4) |
| 43 | + sa = self.structarr |
| 44 | + A = np.hstack((sa['xras'], sa['yras'], sa['zras'])) * sa['voxelsize'] |
| 45 | + b = sa['cras'] - A.dot(sa['volume']) / 2 |
| 46 | + affine[:3, :3] = A |
| 47 | + affine[:3, [3]] = b |
| 48 | + return affine |
| 49 | + |
| 50 | + def to_string(self): |
| 51 | + sa = self.structarr |
| 52 | + lines = [ |
| 53 | + 'valid = {} # volume info {:s}valid'.format( |
| 54 | + sa['valid'], '' if sa['valid'] else 'in'), |
| 55 | + 'filename = {}'.format(sa['filename']), |
| 56 | + 'volume = {:d} {:d} {:d}'.format(*sa['volume'].flatten()), |
| 57 | + 'voxelsize = {:.15e} {:.15e} {:.15e}'.format( |
| 58 | + *sa['voxelsize'].flatten()), |
| 59 | + 'xras = {:.15e} {:.15e} {:.15e}'.format(*sa['xras'].flatten()), |
| 60 | + 'yras = {:.15e} {:.15e} {:.15e}'.format(*sa['yras'].flatten()), |
| 61 | + 'zras = {:.15e} {:.15e} {:.15e}'.format(*sa['zras'].flatten()), |
| 62 | + 'cras = {:.15e} {:.15e} {:.15e}'.format(*sa['cras'].flatten()), |
| 63 | + ] |
| 64 | + return '\n'.join(lines) |
| 65 | + |
| 66 | + @classmethod |
| 67 | + def from_image(klass, img): |
| 68 | + volgeom = klass() |
| 69 | + sa = volgeom.structarr |
| 70 | + sa['valid'] = 1 |
| 71 | + sa['volume'][:, 0] = img.shape[:3] # Assumes xyzt-ordered image |
| 72 | + sa['voxelsize'][:, 0] = img.header.get_zooms()[:3] |
| 73 | + A = img.affine[:3, :3] |
| 74 | + b = img.affine[:3, [3]] |
| 75 | + cols = A * (1 / sa['voxelsize']) |
| 76 | + sa['xras'] = cols[:, [0]] |
| 77 | + sa['yras'] = cols[:, [1]] |
| 78 | + sa['zras'] = cols[:, [2]] |
| 79 | + sa['cras'] = b + A.dot(sa['volume']) / 2 |
| 80 | + try: |
| 81 | + sa['filename'] = img.file_map['image'].filename |
| 82 | + except: |
| 83 | + pass |
| 84 | + |
| 85 | + return volgeom |
| 86 | + |
| 87 | + @classmethod |
| 88 | + def from_string(klass, string): |
| 89 | + volgeom = klass() |
| 90 | + sa = volgeom.structarr |
| 91 | + lines = string.splitlines() |
| 92 | + for key in ('valid', 'filename', 'volume', 'voxelsize', |
| 93 | + 'xras', 'yras', 'zras', 'cras'): |
| 94 | + label, valstring = lines.pop(0).split(' = ') |
| 95 | + assert label.strip() == key |
| 96 | + |
| 97 | + val = np.genfromtxt([valstring.encode()], |
| 98 | + dtype=klass.dtype[key]) |
| 99 | + sa[key] = val.reshape(sa[key].shape) if val.size else '' |
| 100 | + |
| 101 | + return volgeom |
| 102 | + |
| 103 | + |
| 104 | +class LinearTransform(StringBasedStruct): |
| 105 | + template_dtype = np.dtype([ |
| 106 | + ('mean', 'f4', (3, 1)), # x0, y0, z0 |
| 107 | + ('sigma', 'f4'), |
| 108 | + ('m_L', 'f4', (4, 4)), |
| 109 | + ('m_dL', 'f4', (4, 4)), |
| 110 | + ('m_last_dL', 'f4', (4, 4)), |
| 111 | + ('src', VolumeGeometry), |
| 112 | + ('dst', VolumeGeometry), |
| 113 | + ('label', 'i4')]) |
| 114 | + dtype = template_dtype |
| 115 | + |
| 116 | + def __getitem__(self, idx): |
| 117 | + val = super(LinearTransform, self).__getitem__(idx) |
| 118 | + if idx in ('src', 'dst'): |
| 119 | + val = VolumeGeometry(val) |
| 120 | + return val |
| 121 | + |
| 122 | + def to_string(self): |
| 123 | + sa = self.structarr |
| 124 | + lines = [ |
| 125 | + 'mean = {:6.4f} {:6.4f} {:6.4f}'.format( |
| 126 | + *sa['mean'].flatten()), |
| 127 | + 'sigma = {:6.4f}'.format(float(sa['sigma'])), |
| 128 | + '1 4 4', |
| 129 | + ('{:18.15e} ' * 4).format(*sa['m_L'][0]), |
| 130 | + ('{:18.15e} ' * 4).format(*sa['m_L'][1]), |
| 131 | + ('{:18.15e} ' * 4).format(*sa['m_L'][2]), |
| 132 | + ('{:18.15e} ' * 4).format(*sa['m_L'][3]), |
| 133 | + 'src volume info', |
| 134 | + self['src'].to_string(), |
| 135 | + 'dst volume info', |
| 136 | + self['dst'].to_string(), |
| 137 | + ] |
| 138 | + return '\n'.join(lines) |
| 139 | + |
| 140 | + @classmethod |
| 141 | + def from_string(klass, string): |
| 142 | + lt = klass() |
| 143 | + sa = lt.structarr |
| 144 | + lines = string.splitlines() |
| 145 | + for key in ('mean', 'sigma'): |
| 146 | + label, valstring = lines.pop(0).split(' = ') |
| 147 | + assert label.strip() == key |
| 148 | + |
| 149 | + val = np.genfromtxt([valstring.encode()], |
| 150 | + dtype=klass.dtype[key]) |
| 151 | + sa[key] = val.reshape(sa[key].shape) |
| 152 | + assert lines.pop(0) == '1 4 4' |
| 153 | + val = np.genfromtxt([valstring.encode() for valstring in lines[:4]], |
| 154 | + dtype='f4') |
| 155 | + sa['m_L'] = val |
| 156 | + lines = lines[4:] |
| 157 | + assert lines.pop(0) == 'src volume info' |
| 158 | + sa['src'] = np.asanyarray(VolumeGeometry.from_string('\n'.join(lines[:8]))) |
| 159 | + lines = lines[8:] |
| 160 | + assert lines.pop(0) == 'dst volume info' |
| 161 | + sa['dst'] = np.asanyarray(VolumeGeometry.from_string('\n'.join(lines))) |
| 162 | + return lt |
| 163 | + |
| 164 | + |
| 165 | +class LinearTransformArray(StringBasedStruct): |
| 166 | + template_dtype = np.dtype([ |
| 167 | + ('type', 'i4'), |
| 168 | + ('nxforms', 'i4'), |
| 169 | + ('subject', 'U1024'), |
| 170 | + ('fscale', 'f4')]) |
| 171 | + dtype = template_dtype |
| 172 | + _xforms = None |
| 173 | + |
| 174 | + def __init__(self, |
| 175 | + binaryblock=None, |
| 176 | + endianness=None, |
| 177 | + check=True): |
| 178 | + super(LinearTransformArray, self).__init__(binaryblock, endianness, check) |
| 179 | + self._xforms = [LinearTransform() |
| 180 | + for _ in range(self.structarr['nxforms'])] |
| 181 | + |
| 182 | + def __getitem__(self, idx): |
| 183 | + if idx == 'xforms': |
| 184 | + return self._xforms |
| 185 | + if idx == 'nxforms': |
| 186 | + return len(self._xforms) |
| 187 | + return super(LinearTransformArray, self).__getitem__(idx) |
| 188 | + |
| 189 | + def to_string(self): |
| 190 | + code = int(self['type']) |
| 191 | + header = [ |
| 192 | + 'type = {} # {}'.format(code, transform_codes.label[code]), |
| 193 | + 'nxforms = {}'.format(self['nxforms'])] |
| 194 | + xforms = [xfm.to_string() for xfm in self._xforms] |
| 195 | + footer = [ |
| 196 | + 'subject {}'.format(self['subject']), |
| 197 | + 'fscale {:.6f}'.format(float(self['fscale']))] |
| 198 | + return '\n'.join(header + xforms + footer) |
| 199 | + |
| 200 | + @classmethod |
| 201 | + def from_string(klass, string): |
| 202 | + lta = klass() |
| 203 | + sa = lta.structarr |
| 204 | + lines = string.splitlines() |
| 205 | + for key in ('type', 'nxforms'): |
| 206 | + label, valstring = lines.pop(0).split(' = ') |
| 207 | + assert label.strip() == key |
| 208 | + |
| 209 | + val = np.genfromtxt([valstring.encode()], |
| 210 | + dtype=klass.dtype[key]) |
| 211 | + sa[key] = val.reshape(sa[key].shape) if val.size else '' |
| 212 | + for _ in range(sa['nxforms']): |
| 213 | + lta._xforms.append( |
| 214 | + LinearTransform.from_string('\n'.join(lines[:25]))) |
| 215 | + lines = lines[25:] |
| 216 | + for key in ('subject', 'fscale'): |
| 217 | + # Optional keys |
| 218 | + if not lines[0].startswith(key): |
| 219 | + continue |
| 220 | + label, valstring = lines.pop(0).split(' ') |
| 221 | + assert label.strip() == key |
| 222 | + |
| 223 | + val = np.genfromtxt([valstring.encode()], |
| 224 | + dtype=klass.dtype[key]) |
| 225 | + sa[key] = val.reshape(sa[key].shape) if val.size else '' |
| 226 | + |
| 227 | + assert len(lta._xforms) == sa['nxforms'] |
| 228 | + |
| 229 | + return lta |
| 230 | + |
| 231 | + @classmethod |
| 232 | + def from_fileobj(klass, fileobj, check=True): |
| 233 | + return klass.from_string(fileobj.read()) |
0 commit comments