This file is indexed.

/usr/lib/python2.7/dist-packages/dipy/tracking/tests/test_propagation.py is in python-dipy 0.10.1-1.

This file is owned by root:root, with mode 0o644.

The actual contents of the file can be viewed below.

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
import os
import numpy as np
import numpy.testing

from dipy.data import get_data, get_sphere
from dipy.core.gradients import gradient_table
from dipy.reconst.gqi import GeneralizedQSamplingModel
from dipy.reconst.dti import TensorModel, quantize_evecs
from dipy.tracking import utils
from dipy.tracking.eudx import EuDX
from dipy.tracking.propspeed import ndarray_offset, eudx_both_directions
from dipy.tracking.metrics import length
from dipy.tracking.propspeed import map_coordinates_trilinear_iso

import nibabel as ni

from nose.tools import assert_true, assert_false, \
     assert_equal, assert_raises, assert_almost_equal

from numpy.testing import (assert_array_equal,
                           assert_array_almost_equal,
                           run_module_suite)


def stepped_1d(arr_1d):
    # Make a version of `arr_1d` which is not contiguous
    return np.vstack((arr_1d, arr_1d)).ravel(order='F')[::2]


def test_offset():
    # Test ndarray_offset function
    for dt in (np.int32, np.float64):
        index = np.array([1, 1], dtype=np.intp)
        A = np.array([[1,0,0],[0,2,0],[0,0,3]], dtype=dt)
        strides = np.array(A.strides, np.intp)
        i_size = A.dtype.itemsize
        assert_equal(ndarray_offset(index, strides, 2, i_size), 4)
        assert_equal(A.ravel()[4], A[1,1])
        # Index and strides arrays must be C-continuous. Test this is enforced
        # by using non-contiguous versions of the input arrays.
        assert_raises(ValueError, ndarray_offset,
                      stepped_1d(index), strides, 2, i_size)
        assert_raises(ValueError, ndarray_offset,
                      index, stepped_1d(strides), 2, i_size)


def test_trilinear_interp_cubic_voxels():
    A=np.ones((17,17,17))
    B=np.zeros(3)
    strides=np.array(A.strides, np.intp)
    A[7,7,7]=2
    points=np.array([[0,0,0],[7.,7.5,7.],[3.5,3.5,3.5]])
    map_coordinates_trilinear_iso(A,points,strides,3,B)
    assert_array_almost_equal(B,np.array([ 1. ,  1.5,  1. ]))
    # All of the input array, points array, strides array and output array must
    # be C-contiguous.  Check by passing in versions that aren't C contiguous
    assert_raises(ValueError, map_coordinates_trilinear_iso,
                  A.copy(order='F'), points, strides, 3, B)
    assert_raises(ValueError, map_coordinates_trilinear_iso,
                  A, points.copy(order='F'), strides, 3, B)
    assert_raises(ValueError, map_coordinates_trilinear_iso,
                  A, points, stepped_1d(strides), 3, B)
    assert_raises(ValueError, map_coordinates_trilinear_iso,
                  A, points, strides, 3, stepped_1d(B))


def test_eudx_further():
    """ Cause we love testin.. ;-)
    """

    fimg,fbvals,fbvecs=get_data('small_101D')

    img=ni.load(fimg)
    affine=img.get_affine()
    data=img.get_data()
    gtab = gradient_table(fbvals, fbvecs)
    tensor_model = TensorModel(gtab)
    ten = tensor_model.fit(data)
    x,y,z=data.shape[:3]
    seeds=np.zeros((10**4,3))
    for i in range(10**4):
        rx=(x-1)*np.random.rand()
        ry=(y-1)*np.random.rand()
        rz=(z-1)*np.random.rand()
        seeds[i]=np.ascontiguousarray(np.array([rx,ry,rz]),dtype=np.float64)

    sphere = get_sphere('symmetric724')

    ind = quantize_evecs(ten.evecs)
    eu=EuDX(a=ten.fa, ind=ind, seeds=seeds,
            odf_vertices=sphere.vertices, a_low=.2)
    T=[e for e in eu]

    #check that there are no negative elements
    for t in T:
        assert_equal(np.sum(t.ravel()<0),0)

    # Test eudx with affine
    def random_affine(seeds):
        affine = np.eye(4)
        affine[:3, :] = np.random.random((3, 4))
        seeds = np.dot(seeds, affine[:3, :3].T)
        seeds += affine[:3, 3]
        return affine, seeds

    # Make two random affines and move seeds
    affine1, seeds1 = random_affine(seeds)
    affine2, seeds2 = random_affine(seeds)

    # Make tracks using different affines
    eu1 = EuDX(a=ten.fa, ind=ind, odf_vertices=sphere.vertices,
               seeds=seeds1, a_low=.2, affine=affine1)
    eu2 = EuDX(a=ten.fa, ind=ind, odf_vertices=sphere.vertices,
               seeds=seeds2, a_low=.2, affine=affine2)

    # Move from eu2 affine2 to affine1
    eu2_to_eu1 = utils.move_streamlines(eu2, output_space=affine1,
                                        input_space=affine2)
    # Check that the tracks are the same
    for sl1, sl2 in zip(eu1, eu2_to_eu1):
        assert_array_almost_equal(sl1, sl2)


def test_eudx_bad_seed():
    """Test passing a bad seed to eudx"""
    fimg, fbvals, fbvecs = get_data('small_101D')

    img = ni.load(fimg)
    affine = img.get_affine()
    data = img.get_data()
    gtab = gradient_table(fbvals, fbvecs)
    tensor_model = TensorModel(gtab)
    ten = tensor_model.fit(data)
    ind = quantize_evecs(ten.evecs)

    sphere = get_sphere('symmetric724')
    seed = [1000000., 1000000., 1000000.]
    eu = EuDX(a=ten.fa, ind=ind, seeds=[seed],
              odf_vertices=sphere.vertices, a_low=.2)
    assert_raises(ValueError, list, eu)

    print(data.shape)
    seed = [1., 5., 8.]
    eu = EuDX(a=ten.fa, ind=ind, seeds=[seed],
              odf_vertices=sphere.vertices, a_low=.2)
    track = list(eu)

    seed = [-1., 1000000., 1000000.]
    eu = EuDX(a=ten.fa, ind=ind, seeds=[seed],
              odf_vertices=sphere.vertices, a_low=.2)
    assert_raises(ValueError, list, eu)


def test_eudx_boundaries():
    """
    This test checks that the tracking will exclude seeds in both directions.
    Here we create a volume of shape (50, 60, 40) and we will add 2 seeds
    exactly at the volume's boundaries (49, 0, 0) and (0, 0, 0). Those should
    not generate any streamlines as EuDX does not interpolate on the boundary
    voxels. We also add 3 seeds not in the boundaries which should generate
    streamlines without a problem.
    """

    fa = np.ones((50, 60, 40))
    ind = np.zeros(fa.shape)
    sphere = get_sphere('repulsion724')

    seed = [49., 0, 0]
    seed2 = [0., 0, 0]
    seed3 = [48., 0, 0]
    seed4 = [1., 0, 0]
    seed5 = [5., 5, 5]

    eu = EuDX(a=fa, ind=ind, seeds=[seed, seed2, seed3, seed4, seed5],
              odf_vertices=sphere.vertices, a_low=.2,
              total_weight=0.)
    track = list(eu)

    assert_equal(len(track), 3)


def test_eudx_both_directions_errors():
    # Test error conditions for both directions function
    sphere = get_sphere('symmetric724')
    seed = np.zeros(3, np.float64)
    qa = np.zeros((4, 5, 6, 7), np.float64)
    ind = qa.copy()
    # All of seed, qa, ind, odf_vertices must be C-contiguous.  Check by
    # passing in versions that aren't C contiguous
    assert_raises(ValueError, eudx_both_directions,
                  stepped_1d(seed), 0, qa, ind, sphere.vertices, 0.5, 0.1,
                  1., 1., 2)
    assert_raises(ValueError, eudx_both_directions,
                  seed, 0, qa[..., ::2], ind, sphere.vertices, 0.5, 0.1,
                  1., 1., 2)
    assert_raises(ValueError, eudx_both_directions,
                  seed, 0, qa, ind[..., ::2], sphere.vertices, 0.5, 0.1,
                  1., 1., 2)
    assert_raises(ValueError, eudx_both_directions,
                  seed, 0, qa, ind, sphere.vertices[::2], 0.5, 0.1,
                  1., 1., 2)


if __name__ == '__main__':
    run_module_suite()