Here is what ended up working for me:
def setSlicePoseFromSliceNormalAndPosition(sliceNode, sliceNormal, slicePosition, offsetOrCenter='offset'):
""" Set slice pose from the provided plane normal and position. The most similar canonical view (sagittal, coronal, or axial) is determined by the direction of the sliceNormal. For example, if the largest component of the sliceNormal is in the superior direction (or negative superior direction), then the most similar canonical view is the axial view. Care is taken to not allow misleading views, such as an axial-appearing view which has the patient's right on the right hand side (reversed from the standard axial view). If the supplied normal would naturally lead to such a misleading view, it is negated. :param offsetOrCenter should be either 'offset' or 'center', and is used to determine whether the view should be shifed to the given slicePosition by jumping by offsetting or by centering (using sliceNode.JumpSlicesByOffsetting or sliceNode.JumpSlicesByCentering, respectively). """
# Ensure sliceNormal is numpy array
sliceNormal = np.array(sliceNormal)
# Find largest and smallest components of sliceNormal
largestComponentAxis = np.argmax(np.abs(sliceNormal)) # used to determine most similar view
smallestComponentAxis = np.argmin(np.abs(sliceNormal)) # used to decide which of the slice X or Y axes to precisely align
if largestComponentAxis==smallestComponentAxis:
raise ValueError('Largest and smallest slice normal components are the same!!')
viewLookup = {0: 'sagittal',
1: 'coronal',
2: 'axial'}
mostSimilarView = viewLookup[largestComponentAxis]
#print mostSimilarView
# In order to make the adjusted view as similar as possible to a standard view, it may be necessary to use
# the negation of the given slice normal rather than the given slice normal.
# The cross product relationships among the 3 axes are
# Y = NxX, X = YxN, N = XxY (where x means vector cross product here, recall that order matters)
# In each of the mostSimilarView cases which follow, we
# 1) use the sign of the largest component to determine whether the sliceNormal needs to be negated
# to more closely conform to the corresponding canonical view,
# 2) use the smallest component to determine which of the remaining axes is closer to lying in the new slice view
# plane. The smaller the component in the direction of the sliceNormal, the larger the component must be in
# the orthogonal plane.
# 3) set the direction of the slice X axis so that the projection of the most in-plane of the R,A, or S axes
# points precisely in the correct direction. If it is the new Y axis which is to be aligned, the X axis direction
# is determined by finding the cross product of the Y axis with the sliceNormal.
#
# Note that we can handle the normalization to unit vectors later, so if we want the slice Y axis to align as closely
# as possible to the Superior axis, we can assign sliceAxisY to be [0,0,1] for now, even if the sliceNormal guarantees
# that the vector [0,0,1] is not in the slice plane. After the cross products and normalization, the sliceAxisY
if mostSimilarView=='sagittal':
# Canonical view is image up is superior, image right is posterior, and patient right points out of screen
if sliceNormal[largestComponentAxis]>0:
sliceNormal = -sliceNormal # invert the normal if it would lead the Right axis to point into the screen
if smallestComponentAxis==2:
# if superior is closest to in-plane, say that up should be superior
sliceAxisY = [0,0,1]
sliceAxisX = np.cross(sliceAxisY,sliceNormal)
#print 'up should be superior'
elif smallestComponentAxis==1:
# if right is closest to in-plane, say that right should be posterior
sliceAxisX = [0,-1,0]
#print 'right should be posterior'
elif mostSimilarView=='coronal':
# Canonical view is image up is superior, image right is patient left, and anterior points into the screen
if sliceNormal[largestComponentAxis]<0:
sliceNormal = -sliceNormal # invert the normal if it would lead the anterior axis to point out of the screen
if smallestComponentAxis==2:
# if superior is closest to in-plane, say that up should be superior
sliceAxisY = [0,0,1]
sliceAxisX = np.cross(sliceAxisY,sliceNormal)
#print('up should be superior')
elif smallestComponentAxis==0:
# if right is closest to in-plane, say that right should be left
sliceAxisX = [-1,0,0]
print 'right should be patient left'
elif mostSimilarView=='axial':
# Canonical view is image up is anterior and image right is left
if sliceNormal[largestComponentAxis]>0:
sliceNormal=-sliceNormal # invert the normal if it would lead the Superior axis to point into the screen
if smallestComponentAxis==1:
# if anterior is closest to in-plane, say that up should be anterior
sliceAxisY = [0,1,0]
sliceAxisX = np.cross(sliceAxisY,sliceNormal)
print 'up should be anterior'+' sliceAxisX='+str(sliceAxisX)
elif smallestComponentAxis==0:
# if right is closest to in-plane, say that right should be left
sliceAxisX = [-1,0,0]
print 'right should be patient left'+' sliceAxisX='+str(sliceAxisX)
else:
# Should throw error here
raise ValueError('mostSimilarView must be axial, coronal, or sagittal, but is "'+mostSimilarView+'"')
vtkSliceToRAS = sliceNode.GetSliceToRAS()
oldPosition = [vtkSliceToRAS.GetElement(0,3),
vtkSliceToRAS.GetElement(1,3),
vtkSliceToRAS.GetElement(2,3)]
# In the new SliceToRAS matrix, sliceAxisX should be the first column, sliceAxisY should be the second column,
# and the normal vector should be the third colum (sliceAxisZ). These should all be normalized and orthogonal.
# SetSliceToRASByNTP handles this by 1) normalizing the input N and T vectors (slice normal and sliceAxisX),
# 2) crossing them to generate sliceAxisY, and then 3) crossing the slice normal and sliceAxisY to generate
# an updated version of sliceAxisX which is guaranteed to be orthogonal to both sliceAxisY and the slice normal.
sliceNormal = np.array(sliceNormal)
magNormal = np.linalg.norm(sliceNormal)
sliceNormal = sliceNormal/magNormal # normalize
sliceAxisX = np.array(sliceAxisX)
magSliceAxisX = np.linalg.norm(sliceAxisX)
sliceAxisX = sliceAxisX/magSliceAxisX
sliceAxisY = np.cross(sliceNormal,sliceAxisX) # Y = NxX (recall that order matters for cross products)
magSliceAxisY = np.linalg.norm(sliceAxisY)
sliceAxisY = sliceAxisY/magSliceAxisY # normalize
sliceAxisX = np.cross(sliceAxisY,sliceNormal) # X = YxN (this order ensures X points in the same general direction as the original X, but is guaranteed to be orthogonal)
magSliceAxisX = np.linalg.norm(sliceAxisX)
# Set the elements of the SliceToRAS matrix
for rowInd in range(3):
vtkSliceToRAS.SetElement(rowInd,0,sliceAxisX[rowInd])
vtkSliceToRAS.SetElement(rowInd,1,sliceAxisY[rowInd])
vtkSliceToRAS.SetElement(rowInd,2,sliceNormal[rowInd])
if offsetOrCenter.lower()=='center':
# use the requested position as the slice center position
vtkSliceToRAS.SetElement(rowInd,3,slicePosition[rowInd])
sliceNode.UpdateMatrices() # This applies the changes to the SliceToRAS matrix
if offsetOrCenter.lower()=='offset':
# keep the original slice center position, but then adjust offset until it contains requested position
sliceNode.JumpSliceByOffsetting(slicePosition[0],slicePosition[1],slicePosition[2])
elif offsetOrCenter.lower()=='center':
sliceNode.JumpSliceByCentering(slicePosition[0],slicePosition[1],slicePosition[2]) [/code]