diff --git a/fsl/fsleyes/displaycontext/tensoropts.py b/fsl/fsleyes/displaycontext/tensoropts.py
index e3d55e6b9b7e714fdccfa7a3726e848b9629cb65..f8455697db7028cc88d3c1f3ca4181d19e78a3f4 100644
--- a/fsl/fsleyes/displaycontext/tensoropts.py
+++ b/fsl/fsleyes/displaycontext/tensoropts.py
@@ -25,7 +25,10 @@ class TensorOpts(vectoropts.VectorOpts):
# Enable/disable lighting effects
lighting = props.Boolean(default=False)
-
+
+
+ # Tensor ellipsoid resolution
+ tensorResolution = props.Int(default=10)
def __init__(self, *args, **kwargs):
diff --git a/fsl/fsleyes/gl/__init__.py b/fsl/fsleyes/gl/__init__.py
index fc892789d59f2216a91476427df7d624135cb3e0..253938ee5e44f88eb713d03d11acd66e56ac6b65 100644
--- a/fsl/fsleyes/gl/__init__.py
+++ b/fsl/fsleyes/gl/__init__.py
@@ -123,6 +123,7 @@ as:
~fsl.fsleyes.gl.gllinevector.GLLineVector
~fsl.fsleyes.gl.glrgbvector.GLRGBVector
~fsl.fsleyes.gl.glmodel.GLModel
+ ~fsl.fsleyes.gl.gltensor.GLTensor
These objects are created and destroyed automatically by :class:`.SliceCanvas`
instances, so application code does not need to worry about them too much.
diff --git a/fsl/fsleyes/gl/gl21/gllinevector_funcs.py b/fsl/fsleyes/gl/gl21/gllinevector_funcs.py
index b788a27dffc8991ad0593fa41e04624c2c14c16c..8617f60fe715dec76d99fa222eea16c113481e34 100644
--- a/fsl/fsleyes/gl/gl21/gllinevector_funcs.py
+++ b/fsl/fsleyes/gl/gl21/gllinevector_funcs.py
@@ -43,11 +43,10 @@ def init(self):
:class:`.Image` overlay.
"""
- self.shaders = None
- self.vertexBuffer = gl.glGenBuffers(1)
- self.texCoordBuffer = gl.glGenBuffers(1)
- self.vertexIDBuffer = gl.glGenBuffers(1)
- self.lineVertices = None
+ self.shaders = None
+ self.vertexBuffer = gl.glGenBuffers(1)
+ self.texCoordBuffer = gl.glGenBuffers(1)
+ self.vertexIDBuffer = gl.glGenBuffers(1)
self._vertexResourceName = '{}_{}_vertices'.format(
type(self).__name__, id(self.image))
diff --git a/fsl/fsleyes/gl/gl21/gltensor_funcs.py b/fsl/fsleyes/gl/gl21/gltensor_funcs.py
index 0e925ea1b3e9600f788dcd127c63e6afd1de3de2..300a95b8f0a1d4b45e220b64a126466f55d444c5 100644
--- a/fsl/fsleyes/gl/gl21/gltensor_funcs.py
+++ b/fsl/fsleyes/gl/gl21/gltensor_funcs.py
@@ -7,10 +7,15 @@
import logging
-import OpenGL.GL as gl
-import fsl.fsleyes.gl as fslgl
+import numpy as np
+import OpenGL.GL as gl
+import OpenGL.raw.GL._types as gltypes
+import OpenGL.GL.ARB.instanced_arrays as arbia
+
+import fsl.utils.transform as transform
import fsl.fsleyes.gl.resources as glresources
+import fsl.fsleyes.gl.routines as glroutines
import fsl.fsleyes.gl.textures as textures
import fsl.fsleyes.gl.shaders as shaders
@@ -20,7 +25,6 @@ log = logging.getLogger(__name__)
def init(self):
-
image = self.image
v1 = image.V1()
@@ -38,7 +42,7 @@ def init(self):
imgs = [ v1, v2, v3, l1, l2, l3]
for name, img in zip(names, imgs):
- texName = '{}_{}'.format(type(self).__name__, id(img))
+ texName = '{}_{}_{}'.format(type(self).__name__, name, id(img))
if name[0] == 'v':
prefilter = vPrefilter
@@ -58,17 +62,32 @@ def init(self):
setattr(self, '{}Texture'.format(name), tex)
+ self.vertexBuffer = gl.glGenBuffers(1)
+ self.voxelBuffer = gl.glGenBuffers(1)
self.shaders = None
compileShaders(self)
+ updateShaderState(self)
def destroy(self):
- pass
+ gl.glDeleteBuffers(1, gltypes.GLuint(self.vertexBuffer))
+ gl.glDeleteBuffers(1, gltypes.GLuint(self.voxelBuffer))
+ gl.glDeleteProgram(self.shaders)
+
+ names = ['v1', 'v2', 'v3', 'l1', 'l2', 'l3']
+
+ for name in names:
+ attrName = '{}Texture'.format(name)
+ tex = getattr(self, attrName)
+
+ glresources.delete(tex.getTextureName())
+ setattr(self, attrName, None)
def compileShaders(self):
+
if self.shaders is not None:
gl.glDeleteProgram(self.shaders)
@@ -77,18 +96,175 @@ def compileShaders(self):
self.shaders = shaders.compileShaders(vertShaderSrc, fragShaderSrc)
+ vertUniforms = ['v1Texture', 'v2Texture', 'v3Texture',
+ 'l1Texture', 'l2Texture', 'l3Texture',
+ 'v1ValXform', 'v2ValXform', 'v3ValXform',
+ 'l1ValXform', 'l2ValXform', 'l3ValXform',
+ 'voxToDisplayMat', 'imageShape', 'resolution']
+
+ vertAtts = ['voxel', 'vertex']
+
+ fragUniforms = ['imageTexture', 'modTexture', 'modThreshold',
+ 'xColourTexture', 'yColourTexture', 'zColourTexture',
+ 'voxValXform', 'cmapXform', 'imageShape',
+ 'useSpline']
+
+ for vu in vertUniforms:
+ loc = gl.glGetUniformLocation(self.shaders, vu)
+ setattr(self, '{}Pos'.format(vu), loc)
+
+ for fu in fragUniforms:
+ if hasattr(self, '{}Pos'.format(fu)):
+ continue
+ loc = gl.glGetUniformLocation(self.shaders, fu)
+ setattr(self, '{}Pos'.format(fu), loc)
+
+ for va in vertAtts:
+ loc = gl.glGetAttribLocation(self.shaders, va)
+ setattr(self, '{}Pos'.format(va), loc)
+
def updateShaderState(self):
- pass
+
+ gl.glUseProgram(self.shaders)
+
+ opts = self.displayOpts
+
+ # Textures used by the fragment shader
+ gl.glUniform1i(self.imageTexturePos, 0)
+ gl.glUniform1i(self.modTexturePos, 1)
+ gl.glUniform1i(self.xColourTexturePos, 2)
+ gl.glUniform1i(self.yColourTexturePos, 3)
+ gl.glUniform1i(self.zColourTexturePos, 4)
+
+ # Textures used by the vertex shader
+ gl.glUniform1i(self.v1TexturePos, 5)
+ gl.glUniform1i(self.v2TexturePos, 6)
+ gl.glUniform1i(self.v3TexturePos, 7)
+ gl.glUniform1i(self.l1TexturePos, 8)
+ gl.glUniform1i(self.l2TexturePos, 9)
+ gl.glUniform1i(self.l3TexturePos, 10)
+
+ # Voxel value transforms use by the fragment shader
+ cmapXform = self.xColourTexture.getCoordinateTransform()
+ voxValXform = self.imageTexture.voxValXform
+ v1ValXform = self.v1Texture.voxValXform
+ v2ValXform = self.v2Texture.voxValXform
+ v3ValXform = self.v3Texture.voxValXform
+ l1ValXform = self.l1Texture.voxValXform
+ l2ValXform = self.l2Texture.voxValXform
+ l3ValXform = self.l3Texture.voxValXform
+
+ voxValXform = np.array(voxValXform, dtype=np.float32).ravel('C')
+ cmapXform = np.array(cmapXform, dtype=np.float32).ravel('C')
+ v1ValXform = np.array(v1ValXform, dtype=np.float32).ravel('C')
+ v2ValXform = np.array(v2ValXform, dtype=np.float32).ravel('C')
+ v3ValXform = np.array(v3ValXform, dtype=np.float32).ravel('C')
+ l1ValXform = np.array(l1ValXform, dtype=np.float32).ravel('C')
+ l2ValXform = np.array(l2ValXform, dtype=np.float32).ravel('C')
+ l3ValXform = np.array(l3ValXform, dtype=np.float32).ravel('C')
+
+ gl.glUniformMatrix4fv(self.voxValXformPos, 1, False, voxValXform)
+ gl.glUniformMatrix4fv(self.cmapXformPos, 1, False, cmapXform)
+ gl.glUniformMatrix4fv(self.v1ValXformPos, 1, False, v1ValXform)
+ gl.glUniformMatrix4fv(self.v2ValXformPos, 1, False, v2ValXform)
+ gl.glUniformMatrix4fv(self.v3ValXformPos, 1, False, v3ValXform)
+ gl.glUniformMatrix4fv(self.l1ValXformPos, 1, False, l1ValXform)
+ gl.glUniformMatrix4fv(self.l2ValXformPos, 1, False, l2ValXform)
+ gl.glUniformMatrix4fv(self.l3ValXformPos, 1, False, l3ValXform)
+
+ # Other miscellaneous uniforms
+ imageShape = np.array(self.image.shape[:3], dtype=np.float32)
+ resolution = opts.tensorResolution
+ modThreshold = opts.modThreshold / 100.0
+ useSpline = 0
+
+ gl.glUniform3fv(self.imageShapePos, 1, imageShape)
+ gl.glUniform1f( self.resolutionPos, resolution)
+ gl.glUniform1f( self.modThresholdPos, modThreshold)
+ gl.glUniform1f( self.useSplinePos, useSpline)
+
+ # Vertices and indices
+ vertices = glroutines.unitSphere(resolution)
+
+ self.nVertices = len(vertices)
+
+ vertices = vertices.ravel('C')
+
+ gl.glBindBuffer(gl.GL_ARRAY_BUFFER, self.vertexBuffer)
+ gl.glBufferData(
+ gl.GL_ARRAY_BUFFER, vertices.nbytes, vertices, gl.GL_STATIC_DRAW)
+ gl.glVertexAttribPointer(
+ self.vertexPos, 3, gl.GL_FLOAT, gl.GL_FALSE, 0, None)
+
+ gl.glUseProgram(0)
def preDraw(self):
- pass
+ gl.glUseProgram(self.shaders)
+ self.v1Texture.bindTexture(gl.GL_TEXTURE5)
+ self.v2Texture.bindTexture(gl.GL_TEXTURE6)
+ self.v3Texture.bindTexture(gl.GL_TEXTURE7)
+ self.l1Texture.bindTexture(gl.GL_TEXTURE8)
+ self.l2Texture.bindTexture(gl.GL_TEXTURE9)
+ self.l3Texture.bindTexture(gl.GL_TEXTURE10)
+ gl.glEnableVertexAttribArray(self.voxelPos)
+ gl.glEnableVertexAttribArray(self.vertexPos)
+
def draw(self, zpos, xform=None):
- pass
+
+ image = self.image
+ opts = self.displayOpts
+ v2dMat = opts.getTransform('voxel', 'display')
+ d2vMat = opts.getTransform('display', 'voxel')
+
+ resolution = np.array([opts.resolution] * 3)
+
+ if opts.transform == 'id':
+ resolution = resolution / min(image.pixdim[:3])
+ elif opts.transform == 'pixdim':
+ resolution = map(lambda r, p: max(r, p), resolution, image.pixdim[:3])
+
+ voxels = glroutines.calculateSamplePoints(
+ image.shape,
+ resolution,
+ v2dMat,
+ self.xax,
+ self.yax)[0]
+
+ voxels[:, self.zax] = zpos
+
+ voxels = transform.transform(voxels, d2vMat)
+ nVoxels = len(voxels)
+
+ voxels = np.array(voxels, dtype=np.float32).ravel('C')
+
+ gl.glBindBuffer(gl.GL_ARRAY_BUFFER, self.voxelBuffer)
+ gl.glBufferData(
+ gl.GL_ARRAY_BUFFER, voxels.nbytes, voxels, gl.GL_STATIC_DRAW)
+ gl.glVertexAttribPointer(
+ self.voxelPos, 3, gl.GL_FLOAT, gl.GL_FALSE, 0, None)
+ arbia.glVertexAttribDivisorARB(self.voxelPos, 1)
+
+ if xform is None: xform = v2dMat
+ else: xform = transform.concat(v2dMat, xform)
+
+ xform = np.array(xform, dtype=np.float32).ravel('C')
+ gl.glUniformMatrix4fv(self.voxToDisplayMatPos, 1, False, xform)
+
+ gl.glDrawArraysInstanced(gl.GL_QUADS, 0, self.nVertices, nVoxels)
def postDraw(self):
- pass
+ gl.glUseProgram(0)
+ gl.glBindBuffer(gl.GL_ARRAY_BUFFER, 0)
+ self.v1Texture.unbindTexture()
+ self.v2Texture.unbindTexture()
+ self.v3Texture.unbindTexture()
+ self.l1Texture.unbindTexture()
+ self.l2Texture.unbindTexture()
+ self.l3Texture.unbindTexture()
+ gl.glDisableVertexAttribArray(self.voxelPos)
+ gl.glDisableVertexAttribArray(self.vertexPos)
diff --git a/fsl/fsleyes/gl/gl21/gltensor_vert.glsl b/fsl/fsleyes/gl/gl21/gltensor_vert.glsl
index 6874d1b7eae141c85a05d4ccfec3bfd83a86696d..438564dfed665d042594b99aa580e3a1d6ffcaa1 100644
--- a/fsl/fsleyes/gl/gl21/gltensor_vert.glsl
+++ b/fsl/fsleyes/gl/gl21/gltensor_vert.glsl
@@ -5,27 +5,6 @@
*/
#version 120
-#define PI 3.141592653589793
-#define PI_ON_2 1.5707963267948966
-
-
-/*
- Required inputs:
-
- - Textures containing V1, V2, V3 and L1, L2, L3
-
- - Voxel coordinates (these are the vertices)
- - Vertex index
- - Ellipsoid resolution
-
- **We use (index % resolution) to calculate:**
-
- - U and V angles
-
- -
-
- */
-
uniform sampler3D v1Texture;
uniform sampler3D v2Texture;
uniform sampler3D v3Texture;
@@ -46,53 +25,15 @@ uniform vec3 imageShape;
uniform float resolution;
-attribute float index;
attribute vec3 voxel;
+attribute vec3 vertex;
varying vec3 fragVoxCoord;
varying vec3 fragTexCoord;
-vec3 ellipsoidVertex(float l1, float l2, float l3, float u, float v) {
-
- float cosu = cos(u);
- float cosv = cos(v);
- float sinu = sin(u);
- float sinv = sin(v);
-
- float spcu = sign(cosu) * pow(abs(cosu), 2);
- float spcv = sign(cosv) * pow(abs(cosv), 2);
- float spsu = sign(sinu) * pow(abs(sinu), 2);
- float spsv = sign(sinv) * pow(abs(sinv), 2);
-
- vec3 x;
-
- x.x = spcu * spcv;
- x.y = spcu * spsv;
- x.z = spsu;
-
- return x;
-}
-
-
-
void main(void) {
- float umin = -PI_ON_2;
- float umax = PI_ON_2;
- float vmin = -PI;
- float vmax = PI;
- float ustep = (umax - umin) / resolution;
- float vstep = (vmax - vmin) / resolution;
-
- // Index of this vertex
- // within the ellipsoid
- float ellipsoidIndex = mod(index, resolution);
-
- // Ellipsoid angles for this vertex
- float u = umin + ustep * ellipsoidIndex;
- float v = vmin + vstep * ellipsoidIndex;
-
// Lookup the tensor parameters from the textures
vec3 texCoord = (voxel + 0.5) / imageShape;
@@ -112,25 +53,18 @@ void main(void) {
l2 = l2 * l2ValXform[0].x + l2ValXform[3].x;
l3 = l3 * l3ValXform[0].x + l3ValXform[3].x;
- // Calculate the position of
- // this vertex on the ellipsoid.
- // Vertices are grouped into quads -
- // figure out what corner we're on
- vec3 pos;
-
- float corner = mod(ellipsoidIndex, 4);
-
- if (corner == 0) pos = ellipsoidVertex(l1, l2, l3, u, v);
- else if (corner == 1) pos = ellipsoidVertex(l1, l2, l3, u + ustep, v);
- else if (corner == 2) pos = ellipsoidVertex(l1, l2, l3, u + ustep, v + vstep);
- else if (corner == 3) pos = ellipsoidVertex(l1, l2, l3, u, v + vstep);
-
- // Transform the vertex from
- // the fibre coordinate system
- // to the voxel coordinate system
- mat3 eigvecs = mat3(v1, v2, v3);
+ // Transform the sphere by the tensor
+ // eigenvalues and eigenvectors, and
+ // shift it by the voxel coordinates
+ // to transform it in to the image
+ // voxel coordinate system.
+ vec3 pos = vertex;
+ pos.x *= l1 * 150;
+ pos.y *= l2 * 150;
+ pos.z *= l3 * 150;
- pos = pos * eigvecs + voxel;
+ pos = mat3(v1, v2, v3) * pos;
+ pos += voxel;
// Transform from voxels into display
gl_Position = gl_ModelViewProjectionMatrix *
diff --git a/fsl/fsleyes/gl/routines.py b/fsl/fsleyes/gl/routines.py
index aac8bc38063b50607fd9c274fa0d2b1529ae37bf..61bc7738c08c3aa430fd9cfa76756b76864a92dd 100644
--- a/fsl/fsleyes/gl/routines.py
+++ b/fsl/fsleyes/gl/routines.py
@@ -83,6 +83,9 @@ def show2D(xax, yax, width, height, lo, hi):
elif zax == 1:
gl.glRotatef(270, 1, 0, 0)
+ gl.glDepthMask(gl.GL_FALSE)
+ gl.glDisable(gl.GL_DEPTH_TEST)
+
def calculateSamplePoints(shape, resolution, xform, xax, yax, origin='centre'):
"""Calculates a uniform grid of points, in the display coordinate system
@@ -635,3 +638,52 @@ def planeEquation(xyz1, xyz2, xyz3):
(x3 * ((y1 * z2) - (y2 * z1))))
return eq
+
+
+def unitSphere(res):
+ """Generates a unit sphere, as described in the *Sphere Generation*
+ article, on Paul Bourke's excellent website:
+
+ http://paulbourke.net/geometry/circlesphere/
+
+ :arg res: Resolution - the number of angles to sample.
+
+ :returns: A ``numpy.float32`` array of size ``(4 * res**2, 3)`` containing
+ the ``(x, y, z)`` vertices which can be used to draw a unit
+ sphere (using the ``GL_QUAD`` primitive type).
+ """
+
+ # TODO Eliminate duplicate vertices,
+ # and return an index array
+
+ ustep = np.pi / res
+ vstep = np.pi * 2.0 / res
+
+ u = np.linspace(-np.pi / 2, np.pi / 2 + ustep, res + 1, dtype=np.float32)
+ v = np.linspace(-np.pi, np.pi + vstep, res + 1, dtype=np.float32)
+
+ cosu = np.cos(u)
+ cosv = np.cos(v)
+ sinu = np.sin(u)
+ sinv = np.sin(v)
+
+ vertices = np.zeros((res * res * 4, 3), dtype=np.float32)
+
+ cucv = np.outer(cosu[:-1], cosv[:-1]).flatten()
+ cusv = np.outer(cosu[:-1], sinv[:-1]).flatten()
+ cu1cv = np.outer(cosu[1:], cosv[:-1]).flatten()
+ cu1sv = np.outer(cosu[1:], sinv[:-1]).flatten()
+ cu1cv1 = np.outer(cosu[1:], cosv[1:]) .flatten()
+ cu1sv1 = np.outer(cosu[1:], sinv[1:]) .flatten()
+ cucv1 = np.outer(cosu[:-1], cosv[1:]) .flatten()
+ cusv1 = np.outer(cosu[:-1], sinv[1:]) .flatten()
+
+ su = np.repeat(sinu[:-1], res)
+ s1u = np.repeat(sinu[1:], res)
+
+ vertices.T[:, ::4] = [cucv, cusv, su]
+ vertices.T[:, 1::4] = [cu1cv, cu1sv, s1u]
+ vertices.T[:, 2::4] = [cu1cv1, cu1sv1, s1u]
+ vertices.T[:, 3::4] = [cucv1, cusv1, su]
+
+ return vertices
diff --git a/fsl/fsleyes/gl/textures/texture.py b/fsl/fsleyes/gl/textures/texture.py
index 4ba6e0a9f07de22ba80d81a7480815fa58f7cb4d..38bf5af7f48553921e70a0b376f3dd3502ccad33 100644
--- a/fsl/fsleyes/gl/textures/texture.py
+++ b/fsl/fsleyes/gl/textures/texture.py
@@ -143,7 +143,6 @@ class Texture(object):
if textureUnit is not None:
gl.glActiveTexture(textureUnit)
- gl.glEnable(self.__ttype)
gl.glBindTexture(self.__ttype, self.__texture)
@@ -155,7 +154,6 @@ class Texture(object):
if self.__textureUnit is not None:
gl.glActiveTexture(self.__textureUnit)
- gl.glDisable(self.__ttype)
gl.glBindTexture(self.__ttype, 0)