parameters: (images,deformation,transform,referenceMesh)
return_statement: return imagesOut
docstring: images: array of n volumes with shape (x,y,z,n) (in subject space)
      deformation: the estimated deformation in reference space
      transform: the transform matrix from reference to subject space
      referenceMesh: original reference mesh
docstring_summary: images: array of n volumes with shape (x,y,z,n) (in subject space)
      deformation: the estimated deformation in reference space
      transform: the transform matrix from reference to subject space
      referenceMesh: original reference mesh
docstring_tokens: images : array of n volumes with shape ( x y z n ) ( in subject space ) deformation : the estimated deformation in reference space transform : the transform matrix from reference to subject space referenceMesh : original reference mesh
function: def imagesInAtlasSpace(images,deformation,transform,referenceMesh):
      """
      images: array of n volumes with shape (x,y,z,n) (in subject space)
      deformation: the estimated deformation in reference space
      transform: the transform matrix from reference to subject space
      referenceMesh: original reference mesh
      """
      imageSize = list((referenceMesh.points.max(axis=0)+1).astype(np.int))
      denseDeformation = np.zeros(imageSize + [3], dtype=np.double)
      mesh = referenceMesh #just to be clear that it's no longer the "reference" mesh
      mesh.points = mesh.points+deformation
      if np.max(np.absolute(deformation)) > 0:
          maxDeformation = np.max(deformation)
          minDeformation = np.min(deformation)
          deltaDeformation = maxDeformation - minDeformation
          mesh.alphas = (deformation - minDeformation) / deltaDeformation
          tmp = mesh.rasterize_warp(imageSize, -1)
          # Unvisited voxels are marked by zeroes in all three coordinates - except possibly for the origin
          # which has all three coordinates zero as its natural state
          validMask = np.absolute(tmp)
          validMask = np.sum(validMask, axis=3)  ## only zero where all three coordinates are zero
          denseDeformation = (tmp) / (2 ** 16 - 1) * (maxDeformation - minDeformation) + minDeformation
  
          # Due to tetrahedral inside/outside checking performed in the rasterizer, some voxels on the boundary of the
          # image grid are never visited. There also seems to be non-visited voxels inside the image grid (bug in rasterizer).
          # We can resolve both issues by giving to each non-visited voxel the deformation field value of the closest voxel
          # that has been visited, augmented by knowledge that certain deformation components of voxels on the image grid
          # boundaries are zero by theory (sliding boundary conditions)
          closestIndices = scipy.ndimage.morphology.distance_transform_edt(
              validMask == 0, return_indices=True, return_distances=False)
          denseDeformation = denseDeformation[closestIndices[0], closestIndices[1], closestIndices[2]]
  
          # Boundaries are known a priori
          denseDeformation[0,:,:,0] = 0
          denseDeformation[-1,:,:,0] = 0
          denseDeformation[:,0,:,1] = 0
          denseDeformation[:,-1,:,1] = 0
          denseDeformation[:,:,0,2] = 0
          denseDeformation[:,:,-1,2] = 0
  
      linspaces = [np.linspace(0,i-1,i) for i in imageSize]
      X,Y,Z = np.meshgrid(linspaces[0],linspaces[1],linspaces[2],indexing='ij')
      densePositions = np.array([X,Y,Z]).transpose(1,2,3,0)
      densePositions = densePositions + denseDeformation
      densePositions = densePositions.reshape(-1,densePositions.shape[-1])
      atlas = samseg.ProbabilisticAtlas()
      densePositionsInSubjectSpace = atlas.mapPositionsFromTemplateToSubjectSpace(densePositions,transform)
      imagesOut = np.zeros(imageSize + [images.shape[-1]])
      for i in range(images.shape[-1]):
          denseIntensitiesInTemplateSpace = map_coordinates(images[:,:,:,i],densePositionsInSubjectSpace.T,order=1,mode='nearest')
          imagesOut[:,:,:,i] = denseIntensitiesInTemplateSpace.reshape(imageSize)
      return imagesOut
function_tokens: def imagesInAtlasSpace ( images , deformation , transform , referenceMesh ) : imageSize = list ( ( referenceMesh . points . max ( axis = 0 ) + 1 ) . astype ( np . int ) ) denseDeformation = np . zeros ( imageSize + [ 3 ] , dtype = np . double ) mesh = referenceMesh #just to be clear that it's no longer the "reference" mesh mesh . points = mesh . points + deformation if np . max ( np . absolute ( deformation ) ) > 0 : maxDeformation = np . max ( deformation ) minDeformation = np . min ( deformation ) deltaDeformation = maxDeformation - minDeformation mesh . alphas = ( deformation - minDeformation ) / deltaDeformation tmp = mesh . rasterize_warp ( imageSize , - 1 ) # Unvisited voxels are marked by zeroes in all three coordinates - except possibly for the origin # which has all three coordinates zero as its natural state validMask = np . absolute ( tmp ) validMask = np . sum ( validMask , axis = 3 ) ## only zero where all three coordinates are zero denseDeformation = ( tmp ) / ( 2 ** 16 - 1 ) * ( maxDeformation - minDeformation ) + minDeformation # Due to tetrahedral inside/outside checking performed in the rasterizer, some voxels on the boundary of the # image grid are never visited. There also seems to be non-visited voxels inside the image grid (bug in rasterizer). # We can resolve both issues by giving to each non-visited voxel the deformation field value of the closest voxel # that has been visited, augmented by knowledge that certain deformation components of voxels on the image grid # boundaries are zero by theory (sliding boundary conditions) closestIndices = scipy . ndimage . morphology . distance_transform_edt ( validMask == 0 , return_indices = True , return_distances = False ) denseDeformation = denseDeformation [ closestIndices [ 0 ] , closestIndices [ 1 ] , closestIndices [ 2 ] ] # Boundaries are known a priori denseDeformation [ 0 , : , : , 0 ] = 0 denseDeformation [ - 1 , : , : , 0 ] = 0 denseDeformation [ : , 0 , : , 1 ] = 0 denseDeformation [ : , - 1 , : , 1 ] = 0 denseDeformation [ : , : , 0 , 2 ] = 0 denseDeformation [ : , : , - 1 , 2 ] = 0 linspaces = [ np . linspace ( 0 , i - 1 , i ) for i in imageSize ] X , Y , Z = np . meshgrid ( linspaces [ 0 ] , linspaces [ 1 ] , linspaces [ 2 ] , indexing = 'ij' ) densePositions = np . array ( [ X , Y , Z ] ) . transpose ( 1 , 2 , 3 , 0 ) densePositions = densePositions + denseDeformation densePositions = densePositions . reshape ( - 1 , densePositions . shape [ - 1 ] ) atlas = samseg . ProbabilisticAtlas ( ) densePositionsInSubjectSpace = atlas . mapPositionsFromTemplateToSubjectSpace ( densePositions , transform ) imagesOut = np . zeros ( imageSize + [ images . shape [ - 1 ] ] ) for i in range ( images . shape [ - 1 ] ) : denseIntensitiesInTemplateSpace = map_coordinates ( images [ : , : , : , i ] , densePositionsInSubjectSpace . T , order = 1 , mode = 'nearest' ) imagesOut [ : , : , : , i ] = denseIntensitiesInTemplateSpace . reshape ( imageSize ) return imagesOut