Skip to content
Snippets Groups Projects
Commit 48be21b2 authored by Paul McCarthy's avatar Paul McCarthy
Browse files

Voxels are no longer aligned at corners in id/pixdim mode.

parent 36ac1747
No related branches found
No related tags found
No related merge requests found
......@@ -359,6 +359,7 @@ class LabelAtlas(Atlas):
"""
voxelLoc = transform.transform([worldLoc], self.worldToVoxMat.T)[0]
voxelLoc = voxelLoc.round()
if voxelLoc[0] < 0 or \
voxelLoc[1] < 0 or \
......@@ -404,6 +405,7 @@ class ProbabilisticAtlas(Atlas):
location is out of bounds.
"""
voxelLoc = transform.transform([worldLoc], self.worldToVoxMat.T)[0]
voxelLoc = voxelLoc.round()
if voxelLoc[0] < 0 or \
voxelLoc[1] < 0 or \
......
......@@ -409,12 +409,18 @@ def loadClusterResults(featdir, settings, contrast):
# level, the coords are in mm, and need to
# be transformed to voxels.
for c in clusters:
c.zmaxx, c.zmaxy, c.zmaxz = transform.transform(
[[c.zmaxx, c.zmaxy, c.zmaxz]], coordXform)[0]
c.zcogx, c.zcogy, c.zcogz = transform.transform(
[[c.zcogx, c.zcogy, c.zcogz]], coordXform)[0]
c.copemaxx, c.copemaxy, c.copemaxz = transform.transform(
[[c.copemaxx, c.copemaxy, c.copemaxz]], coordXform)[0]
zmax = [c.zmaxx, c.zmaxy, c.zmaxz]
zcog = [c.zcogx, c.zcogy, c.zcogz]
copemax = [c.copemaxx, c.copemaxy, c.copemaxz]
zmax = transform.transform([zmax], coordXform)[0].round()
zcog = transform.transform([zcog], coordXform)[0].round()
copemax = transform.transform([copemax], coordXform)[0].round()
c.zmaxx, c.zmaxy, c.zmaxz = zmax
c.zcogx, c.zcogy, c.zcogz = zcog
c.copemax, c.copemaxy, c.copemaxz = copemax
return clusters
......
......@@ -544,7 +544,7 @@ class LocationPanel(fslpanel.FSLEyesPanel):
source, coords, target, xformed))
if target == 'display': self._displayCtx.location.xyz = xformed
elif target == 'voxel': self.voxelLocation.xyz = np.floor(xformed)
elif target == 'voxel': self.voxelLocation.xyz = np.round(xformed)
elif target == 'world': self.worldLocation.xyz = xformed
......@@ -599,7 +599,7 @@ class LocationPanel(fslpanel.FSLEyesPanel):
vloc = opts.transformCoords(
[self._displayCtx.location.xyz], 'display', 'voxel')[0]
vloc = tuple(map(int, np.floor(vloc)))
vloc = tuple(map(int, np.round(vloc)))
if overlay.is4DImage():
vloc = vloc + (opts.volume,)
......
......@@ -178,30 +178,28 @@ class ModelOpts(fsldisplay.DisplayOpts):
elif oldRefImage is None:
refOpts = self.displayCtx.getOpts(refImage)
newLoc = transform.transform(
[oldLoc],
refOpts.getTransform(self.coordSpace, 'display'))[0]
newLoc = refOpts.transformCoords([oldLoc],
self.coordSpace,
'display')[0]
elif refImage is None:
if oldRefImage is not None:
oldRefOpts = self.displayCtx.getOpts(oldRefImage)
newLoc = transform.transform(
[oldLoc],
oldRefOpts.getTransform('display', self.coordSpace))[0]
newLoc = oldRefOpts.transformCoords([oldLoc],
'display',
self.coordSpace)[0]
elif propName == 'coordSpace':
if self.refImage is not None:
refOpts = self.displayCtx.getOpts(self.refImage)
worldLoc = transform.transform(
worldLoc = refOpts.transformCoords(
[oldLoc],
refOpts.getTransform(
self.getLastValue('coordSpace'),
'world'))[0]
newLoc = transform.transform(
self.getLastValue('coordSpace'),
'world')[0]
newLoc = refOpts.transformCoords(
[worldLoc],
refOpts.getTransform(
'world',
self.coordSpace))[0]
'world',
self.coordSpace)[0]
elif propName == 'transform':
......
......@@ -53,7 +53,7 @@ in the display coordinate system. In other words, a voxel at location::
will be transformed such that, in the display coordinate system, it occupies
the space::
[x - x + 1, y - y + 1, z - z + 1]
[x-0.5 - x+0.5, y-0.5 - y+0.5, z-0.5 - z+0.5]
For example, the voxel::
......@@ -62,17 +62,11 @@ For example, the voxel::
is drawn such that it occupies the space::
[2 - 3, 3 - 4, 4 - 5]
A similar transformation is applied to image data which is displayed in
``pixdim`` space, scaled appropriately.
[1.5 - 2.5, 2.5 - 3.5, 3.5 - 4.5]
This convention was adopted so that multiple images would be aligned at the
bottom left corner when displayed in ``id`` or ``pixzim`` space. But this
convention is in contrast to the convention taken when images are displayed in
world, or ``affine`` space. The ``qform`` and ``sform`` transformation
This convention is the same as the convention taken when images are displayed
in world, or ``affine`` space. The ``qform`` and ``sform`` transformation
matrices in the ``NIFTI1`` specification assume that the voxel coordinates
``[x, y, z]`` correspond to the centre of a voxel. As an example, assuming
that our affine transformation is an identity matrix, the voxel::
......@@ -83,12 +77,6 @@ that our affine transformation is an identity matrix, the voxel::
for an image displayed in ``affine`` space would occupy the space::
[1.5 - 2.5, 2.5 - 3.5, 3.5 - 4.5]
.. note:: With the introduction of **GedMode** I am also going to change this
convention - integer voxel coordinates will map to the voxel centre
(as for the ``affine`` convention described above) regardless of the
coordinate system the image is displayed in.
"""
......@@ -185,13 +173,9 @@ class ImageOpts(fsldisplay.DisplayOpts):
:attr:`.DisplayOpts.bounds` property accordingly.
"""
if self.transform == 'affine': origin = 'centre'
else: origin = 'corner'
lo, hi = transform.axisBounds(
self.overlay.shape[:3],
self.getTransform('voxel', 'display'),
origin=origin)
self.getTransform('voxel', 'display'))
self.bounds[:] = [lo[0], hi[0], lo[1], hi[1], lo[2], hi[2]]
......@@ -276,96 +260,47 @@ class ImageOpts(fsldisplay.DisplayOpts):
def getTransformOffsets(self, from_, to_):
"""When an image is displayed in ``id``/``pixdim`` space, voxel
coordinates map to the voxel corner; i.e. a voxel at ``(0, 1, 2)``
occupies the space ``(0 - 1, 1 - 2, 2 - 3)``.
"""Returns a set of offsets which should be applied to coordinates
before/after applying a transfromation.
In contrast, when an image is displayed in affine space, voxel
coordinates map to the voxel centre, so our voxel from above will
When an image is displayed in ``id``, ``pixdim`` and ``affine`` space,
voxel coordinates map to the voxel centre, so our voxel from above will
occupy the space ``(-0.5 - 0.5, 0.5 - 1.5, 1.5 - 2.5)``. This is
dictated by the NIFTI specification.
dictated by the NIFTI specification. See the
:ref:`note on coordinate systems <volumeopts-coordinate-systems>`.
This function returns some offsets to ensure that the coordinate
transformation from the source space to the target space is valid,
given the above requirements.
A tuple containing two sets of offsets (each of which is a tuple of
three values). The first set is to be applied to the source coordinates
before transformation, and the second set to the target coordinates
after the transformation.
three values). The first set is to be applied to the source
coordinates (in the ``from_`` space) before transformation, and the
second set to the target coordinates (in the ``to_`` space) after the
transformation.
See also the :meth:`transformCoords` method, which will perform the
transformation correctly for you, without you having to worry about
these offsets.
.. note:: These offsets, and this method, will soon become obsolete -
see the note about **GedMode** in the
:ref:`note on coordinate systems
<volumeopts-coordinate-systems>`.
"""
displaySpace = self.transform
pixdim = np.array(self.overlay.pixdim[:3])
offsets = {
# world to voxel transformation
# (regardless of the display space):
#
# add 0.5 to the resulting voxel
# coords, so the _propagate method
# can just floor them to get the
# integer voxel coordinates
('world', 'voxel', displaySpace) : ((0, 0, 0), (0.5, 0.5, 0.5)),
# World to display transformation:
#
# if displaying in id/pixdim space,
# we add half a voxel so that the
# resulting coords are centered
# within a voxel, instead of being
# in the voxel corner
('world', 'display', 'id') : ((0, 0, 0), (0.5, 0.5, 0.5)),
('world', 'display', 'pixdim') : ((0, 0, 0), pixdim / 2.0),
# Display to voxel space:
# If we're displaying in affine space,
# we have the same situation as the
# world -> voxel transform above
('display', 'voxel', 'affine') : ((0, 0, 0), (0.5, 0.5, 0.5)),
# Display to world space:
#
# If we're displaying in id/pixdim
# space, voxel coordinates map to
# the voxel corner, so we need to
# subtract half the voxel width to
# the coordinates before transforming
# to world space.
('display', 'world', 'id') : ((-0.5, -0.5, -0.5), (0, 0, 0)),
('display', 'world', 'pixdim') : (-pixdim / 2.0, (0, 0, 0)),
# Voxel to display space:
#
# If the voxel location was changed,
# we want the display to be moved to
# the centre of the voxel If displaying
# in affine space, voxel coordinates
# map to the voxel centre, so we don't
# need to offset. But if in id/pixdim,
# we need to add 0.5 to the voxel coords,
# as otherwise the transformation will
# put us in the voxel corner.
('voxel', 'display', 'id') : ((0.5, 0.5, 0.5), (0, 0, 0)),
('voxel', 'display', 'pixdim') : ((0.5, 0.5, 0.5), (0, 0, 0)),
}
return offsets.get((from_, to_, displaySpace), ((0, 0, 0), (0, 0, 0)))
these silly offsets.
.. note:: This method was written during a crazy time when, in ``id``
or ``pixdim`` space, voxels ``(0, 0, 0)`` of images overlaid
on each other were aligned at the voxel corner, whereas in
``affine`` space, they were aligned at the voxel
centre. This is no longer the case, so this method is not
actually necessary. But it is still here, and still being
used, just in case we need to change these conventions again
in the future.
"""
return (0, 0, 0), (0, 0, 0)
def transformCoords(self, coords, from_, to_):
"""Transforms the given coordinates from ``from_`` to ``to_``, including
correcting for display space offsets (see :meth:`getTransformOffsets`).
"""Transforms the given coordinates from ``from_`` to ``to_``.
The ``from_`` and ``to_`` parameters must both be one of:
......
......@@ -518,8 +518,7 @@ class VoxelSelection(AnnotationObject):
yax,
zpos,
self.voxToDisplayMat,
self.displayToVoxMat,
origin='corner')
self.displayToVoxMat)
verts = np.array(verts, dtype=np.float32).ravel('C')
texs = np.array(texs, dtype=np.float32).ravel('C')
......
......@@ -149,9 +149,7 @@ def updateShaderState(self):
shape = np.array(list(self.image.shape[:3]) + [0], dtype=np.float32)
invShape = 1.0 / shape
modThreshold = [opts.modThreshold / 100.0, 0.0, 0.0, 0.0]
if opts.transform in ('id', 'pixdim'): offset = [0.0, 0.0, 0.0, 0.0]
else: offset = [0.5, 0.5, 0.5, 0.0]
offset = [0.5, 0.5, 0.5, 0.0]
# Vertex program inputs
shaders.setVertexProgramVector( 0, invShape)
......
......@@ -8,9 +8,7 @@
#
# program.local[0] - (first three components) inverse of image shape
# program.local[1] - (first three components) Offset to apply to transformed
# voxel coordinates before flooring to integers (this
# should be 0 in id/pixdim display mode, and 0.5 in
# affine display mode).
# voxel coordinates before flooring to integers.
#
# Outputs:
# result.position - the vertex position
......
......@@ -201,9 +201,12 @@ def updateShaderState(self):
d2vMat = opts.getTransform('display', 'voxel')
v2dMat = opts.getTransform('voxel', 'display')
if opts.transform in ('id', 'pixdim'): offset = [0, 0, 0]
else: offset = [0.5, 0.5, 0.5]
# The shader adds these offsets to
# transformed voxel coordinates, so
# it can floor them to get integer
# voxel coordinates
offset = [0.5, 0.5, 0.5]
offset = np.array(offset, dtype=np.float32)
imageDims = np.array(imageDims, dtype=np.float32)
d2vMat = np.array(d2vMat, dtype=np.float32).ravel('C')
......@@ -336,16 +339,12 @@ def hardwareDraw(self, zpos, xform=None):
elif opts.transform == 'pixdim':
resolution = map(lambda r, p: max(r, p), resolution, image.pixdim[:3])
if opts.transform == 'affine': origin = 'centre'
else: origin = 'corner'
vertices = glroutines.calculateSamplePoints(
image.shape,
resolution,
v2dMat,
self.xax,
self.yax,
origin)[0]
self.yax)[0]
vertices[:, self.zax] = zpos
......
......@@ -419,17 +419,13 @@ class GLImageObject(GLObject):
corresponding to each vertex
"""
if self.displayOpts.transform == 'affine': origin = 'centre'
else: origin = 'corner'
vertices, voxCoords, texCoords = glroutines.slice2D(
self.image.shape[:3],
self.xax,
self.yax,
zpos,
self.displayOpts.getTransform('voxel', 'display'),
self.displayOpts.getTransform('display', 'voxel'),
origin=origin)
self.displayOpts.getTransform('display', 'voxel'))
if xform is not None:
vertices = transform.transform(vertices, xform)
......
......@@ -99,6 +99,9 @@ def calculateSamplePoints(shape, resolution, xform, xax, yax, origin='centre'):
2).
:arg yax: The vertical display coordinate system axis (0, 1, or 2).
:arg origin: ``centre`` or ``corner``. See the
:func:`.transform.axisBounds` function.
"""
xres = resolution[xax]
......
......@@ -454,7 +454,7 @@ class OrthoEditProfile(orthoviewprofile.OrthoViewProfile):
opts = self._displayCtx.getOpts(self.__currentOverlay)
voxel = opts.transformCoords([canvasPos], 'display', 'voxel')[0]
return np.int32(np.floor(voxel))
return np.int32(np.round(voxel))
def __drawCursorAnnotation(self, canvas, voxel, blockSize=None):
......
......@@ -993,7 +993,7 @@ class TimeSeriesPanel(plotpanel.PlotPanel):
return None
vox = opts.transformCoords([[x, y, z]], 'display', 'voxel')[0]
vox = np.floor(vox)
vox = np.round(vox)
if vox[0] < 0 or \
vox[1] < 0 or \
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment