numpy.ogrid

Here are the examples of the python api numpy.ogrid taken from open source projects. By voting up you can indicate which examples are most useful and appropriate.

43 Examples 7

Example 1

Project: scikit-image Source File: draw.py
def _ellipse_in_shape(shape, center, radiuses):
    """Generate coordinates of points within ellipse bounded by shape."""
    r_lim, c_lim = np.ogrid[0:float(shape[0]), 0:float(shape[1])]
    r_o, c_o = center
    r_r, c_r = radiuses
    distances = ((r_lim - r_o) / r_r) ** 2 + ((c_lim - c_o) / c_r) ** 2
    return np.nonzero(distances < 1)

Example 2

Project: scikit-image Source File: test_denoise.py
def test_denoise_tv_chambolle_3d():
    """Apply the TV denoising algorithm on a 3D image representing a sphere."""
    x, y, z = np.ogrid[0:40, 0:40, 0:40]
    mask = (x - 22)**2 + (y - 20)**2 + (z - 17)**2 < 8**2
    mask = 100 * mask.astype(np.float)
    mask += 60
    mask += 20 * np.random.rand(*mask.shape)
    mask[mask < 0] = 0
    mask[mask > 255] = 255
    res = restoration.denoise_tv_chambolle(mask.astype(np.uint8), weight=0.1)
    assert_(res.dtype == np.float)
    assert_(res.std() * 255 < mask.std())

Example 3

Project: scikit-image Source File: test_unwrap.py
def test_unwrap_2d():
    x, y = np.ogrid[:8, :16]
    image = 2 * np.pi * (x * 0.2 + y * 0.1)
    yield check_unwrap, image
    mask = np.zeros(image.shape, dtype=np.bool)
    mask[4:6, 4:8] = True
    yield check_unwrap, image, mask

Example 4

Project: scikit-image Source File: test_unwrap.py
def test_unwrap_3d():
    x, y, z = np.ogrid[:8, :12, :16]
    image = 2 * np.pi * (x * 0.2 + y * 0.1 + z * 0.05)
    yield check_unwrap, image
    mask = np.zeros(image.shape, dtype=np.bool)
    mask[4:6, 4:6, 1:3] = True
    yield check_unwrap, image, mask

Example 5

Project: dask Source File: reductions.py
def _arg_combine(data, axis, argfunc, keepdims=False):
    """Merge intermediate results from ``arg_*`` functions"""
    axis = None if len(axis) == data.ndim or data.ndim == 1 else axis[0]
    vals = data['vals']
    arg = data['arg']
    if axis is None:
        local_args = argfunc(vals, axis=axis, keepdims=keepdims)
        vals = vals.ravel()[local_args]
        arg = arg.ravel()[local_args]
    else:
        local_args = argfunc(vals, axis=axis)
        inds = np.ogrid[tuple(map(slice, local_args.shape))]
        inds.insert(axis, local_args)
        vals = vals[inds]
        arg = arg[inds]
        if keepdims:
            vals = np.expand_dims(vals, axis)
            arg = np.expand_dims(arg, axis)
    return arg, vals

Example 6

Project: mayavi Source File: test_mlab_source.py
    def setUp(self):
        x, y, z = N.ogrid[-10:10:11j,
                          -10:10:12j,
                          -10:10:20j]
        self.x, self.y, self.z = x, y, z
        dims = (x.shape[0], y.shape[1], z.shape[2])
        self.v = v = N.ones(dims + (3,), float)
        v[...,2] = 2
        v[...,2] = 3
        self.s = s = N.ones(dims, float)
        src = sources.MArraySource()
        src.reset(x=x, y=y, z=z, u=v[...,0], v=v[...,1], w=v[...,2], scalars=s)
        self.src = src

Example 7

Project: mayavi Source File: helper_functions.py
def test_contour3d():
    x, y, z = numpy.ogrid[-5:5:64j, -5:5:64j, -5:5:64j]

    scalars = x * x * 0.5 + y * y + z * z * 2.0

    obj = contour3d(scalars, contours=4, transparent=True)
    return obj

Example 8

Project: mayavi Source File: helper_functions.py
@animate
def test_contour3d_anim(obj=None):
    obj = obj if obj is not None else test_contour3d()
    x, y, z = numpy.ogrid[-5:5:64j, -5:5:64j, -5:5:64j]
    # Now animate the contours.
    ms = obj.mlab_source
    for i in range(1, 10):
        ms.scalars = x * x * 0.5 + y * x * 0.1 * (i + 1) + z * z * 0.25
        yield

Example 9

Project: mayavi Source File: helper_functions.py
Function: test_volume_slice
def test_volume_slice():
    x, y, z = numpy.ogrid[-5:5:64j, -5:5:64j, -5:5:64j]

    scalars = x * x * 0.5 + y * y + z * z * 2.0

    obj = volume_slice(scalars, plane_orientation='x_axes')
    return obj

Example 10

Project: mayavi Source File: helper_functions.py
@animate
def test_volume_slice_anim(obj=None):
    obj = obj if obj is not None else test_volume_slice()
    x, y, z = numpy.ogrid[-5:5:64j, -5:5:64j, -5:5:64j]
    # Now animate the contours.
    ms = obj.mlab_source
    for i in range(1, 10):
        ms.scalars = x * x * 0.5 + y * x * 0.1 * (i + 1) + z * z * 0.25
        yield

Example 11

Project: neural-network-animation Source File: test_backend_pgf.py
Function: test_mixed_mode
@switch_backend('pgf')
def test_mixedmode():
    if not check_for('xelatex'):
        raise SkipTest('xelatex + pgf is required')

    rc_xelatex = {'font.family': 'serif',
                  'pgf.rcfonts': False}
    mpl.rcParams.update(rc_xelatex)

    Y, X = np.ogrid[-1:1:40j, -1:1:40j]
    plt.figure()
    plt.pcolor(X**2 + Y**2).set_rasterized(True)
    compare_figure('pgf_mixedmode.pdf')

Example 12

Project: neural-network-animation Source File: test_backend_pgf.py
@switch_backend('pgf')
def test_bbox_inches():
    if not check_for('xelatex'):
        raise SkipTest('xelatex + pgf is required')

    rc_xelatex = {'font.family': 'serif',
                  'pgf.rcfonts': False}
    mpl.rcParams.update(rc_xelatex)

    Y, X = np.ogrid[-1:1:40j, -1:1:40j]
    fig = plt.figure()
    ax1 = fig.add_subplot(121)
    ax1.plot(range(5))
    ax2 = fig.add_subplot(122)
    ax2.plot(range(5))
    plt.tight_layout()

    bbox = ax1.get_window_extent().transformed(fig.dpi_scale_trans.inverted())
    compare_figure('pgf_bbox_inches.pdf', savefig_kwargs={'bbox_inches': bbox})

Example 13

Project: openmc Source File: lattice.py
Function: indices
    @property
    def indices(self):
        if self.ndim == 2:
            return list(np.broadcast(*np.ogrid[
                :self.shape[1], :self.shape[0]]))
        else:
            return list(np.broadcast(*np.ogrid[
                :self.shape[2], :self.shape[1], :self.shape[0]]))

Example 14

Project: lyman Source File: mixedfx.py
    def _peaks_to_disks(self, peaks, r=4):

        zstat_img = nib.load(self.inputs.zstat_file)

        x, y = np.ogrid[-r:r + 1, -r:r + 1]
        disk = x ** 2 + y ** 2 <= r ** 2
        dilator = np.dstack([disk, disk, np.zeros_like(disk)])
        disk_data = self._dilate_peaks(zstat_img.shape, peaks, dilator)

        disk_img = nib.Nifti1Image(disk_data, zstat_img.get_affine())
        return disk_img

Example 15

Project: lyman Source File: mixedfx.py
    def _peaks_to_spheres(self, peaks, r=4):

        zstat_img = nib.load(self.inputs.zstat_file)

        x, y, z = np.ogrid[-r:r + 1, -r:r + 1, -r:r + 1]
        dilator = x ** 2 + y ** 2 + z ** 2 <= r ** 2
        sphere_data = self._dilate_peaks(zstat_img.shape, peaks, dilator)

        sphere_img = nib.Nifti1Image(sphere_data, zstat_img.get_affine())
        return sphere_img

Example 16

Project: pylearn-parsimony Source File: utils.py
Function: init
    def __init__(self, center, size, shape, **kwargs):
        super(Square, self).__init__(**kwargs)
        self.size = size
        self.center = center
        self.x_grid, self.y_grid, self.z_grid = np.ogrid[0:shape[0],
                                                         0:shape[1],
                                                         0:shape[2]]

Example 17

Project: pylearn-parsimony Source File: utils.py
Function: init
    def __init__(self, center, size, shape, **kwargs):
        super(Dot, self).__init__(**kwargs)
        self.size = size
        self.center = center
        self.x_grid, self.y_grid, self.z_grid = np.ogrid[0:shape[0],
                                                         0:shape[1],
                                                         0:shape[2]]

Example 18

Project: nideep Source File: fcn_solve.py
def upsample_filt(size):
    factor = (size + 1) // 2
    if size % 2 == 1:
        center = factor - 1
    else:
        center = factor - 0.5
    og = np.ogrid[:size, :size]
    return (1 - abs(og[0] - center) / factor) * \
           (1 - abs(og[1] - center) / factor)

Example 19

Project: dipy Source File: test_interpolate.py
def test_NearestNeighborInterpolator():
    # Place integers values at the center of every voxel
    l, m, n, o = np.ogrid[0:6.01, 0:6.01, 0:6.01, 0:4]
    data = l + m + n + o

    nni = NearestNeighborInterpolator(data, (1, 1, 1))
    a, b, c = np.mgrid[.5:6.5:1.6, .5:6.5:2.7, .5:6.5:3.8]
    for ii in xrange(a.size):
        x = a.flat[ii]
        y = b.flat[ii]
        z = c.flat[ii]
        expected_result = int(x) + int(y) + int(z) + o.ravel()
        assert_array_equal(nni[x, y, z], expected_result)
        ind = np.array([x, y, z])
        assert_array_equal(nni[ind], expected_result)
    assert_raises(OutsideImage, nni.__getitem__, (-.1, 0, 0))
    assert_raises(OutsideImage, nni.__getitem__, (0, 8.2, 0))

Example 20

Project: matplotlib2tikz Source File: heat.py
Function: plot
def plot():
    import matplotlib
    from matplotlib import pyplot as plt
    import numpy as np

    fig = plt.figure()
    x, y = np.ogrid[-10:10:100j, -10:10:100j]
    extent = (x.min(), x.max(), y.min(), y.max())
    cmap = matplotlib.cm.get_cmap('gray')
    plt.imshow(x*y, extent=extent, cmap=cmap)
    plt.colorbar()

    return fig

Example 21

Project: numdifftools Source File: extrapolation.py
    def _r_matrix(self, num_terms):
        step = self.step
        i, j = np.ogrid[0:num_terms + 1, 0:num_terms]
        r_mat = np.ones((num_terms + 1, num_terms + 1))
        r_mat[:, 1:] = (1.0 / self.step_ratio) ** (i * (step * j + self.order))
        return r_mat

Example 22

Project: pystruct Source File: edge_feature_graph_crf.py
Function: joint_feature
    def joint_feature(self, x, y):
        """Feature vector associated with instance (x, y).

        Feature representation joint_feature, such that the energy of the configuration
        (x, y) and a weight vector w is given by np.dot(w, joint_feature(x, y)).

        Parameters
        ----------
        x : tuple
            Input representation.

        y : ndarray or tuple
            Either y is an integral ndarray, giving
            a complete labeling for x.
            Or it is the result of a linear programming relaxation. In this
            case, ``y=(unary_marginals, pariwise_marginals)``.

        Returns
        -------
        p : ndarray, shape (size_joint_feature,)
            Feature vector associated with state (x, y).

        """
        self._check_size_x(x)
        features, edges = self._get_features(x), self._get_edges(x)
        n_nodes = features.shape[0]
        edge_features = self._get_edge_features(x)

        if isinstance(y, tuple):
            # y is result of relaxation, tuple of unary and pairwise marginals
            unary_marginals, pw = y
            unary_marginals = unary_marginals.reshape(n_nodes, self.n_states)

        else:
            y = y.reshape(n_nodes)
            gx = np.ogrid[:n_nodes]

            #make one hot encoding
            unary_marginals = np.zeros((n_nodes, self.n_states), dtype=np.int)
            gx = np.ogrid[:n_nodes]
            unary_marginals[gx, y] = 1

            ## pairwise
            pw = np.zeros((edges.shape[0], self.n_states ** 2))
            class_pair_ind = (y[edges[:, 1]] + self.n_states *
                              y[edges[:, 0]])
            pw[np.arange(len(edges)), class_pair_ind] = 1

        pw = np.dot(edge_features.T, pw)
        for i in self.symmetric_edge_features:
            pw_ = pw[i].reshape(self.n_states, self.n_states)
            pw[i] = (pw_ + pw_.T).ravel() / 2.

        for i in self.antisymmetric_edge_features:
            pw_ = pw[i].reshape(self.n_states, self.n_states)
            pw[i] = (pw_ - pw_.T).ravel() / 2.

        unaries_acc = np.dot(unary_marginals.T, features)

        joint_feature_vector = np.hstack([unaries_acc.ravel(), pw.ravel()])
        return joint_feature_vector

Example 23

Project: pystruct Source File: graph_crf.py
Function: joint_feature
    def joint_feature(self, x, y):
        """Feature vector associated with instance (x, y).

        Feature representation joint_feature, such that the energy of the configuration
        (x, y) and a weight vector w is given by np.dot(w, joint_feature(x, y)).

        Parameters
        ----------
        x : tuple
            Unary evidence.

        y : ndarray or tuple
            Either y is an integral ndarray, giving
            a complete labeling for x.
            Or it is the result of a linear programming relaxation. In this
            case, ``y=(unary_marginals, pariwise_marginals)``.

        Returns
        -------
        p : ndarray, shape (size_joint_feature,)
            Feature vector associated with state (x, y).

        """
        self._check_size_x(x)
        features, edges = self._get_features(x), self._get_edges(x)
        n_nodes = features.shape[0]

        if isinstance(y, tuple):
            # y is result of relaxation, tuple of unary and pairwise marginals
            unary_marginals, pw = y
            unary_marginals = unary_marginals.reshape(n_nodes, self.n_states)
            # accuemulate pairwise
            pw = pw.reshape(-1, self.n_states, self.n_states).sum(axis=0)
        else:
            y = y.reshape(n_nodes)
            gx = np.ogrid[:n_nodes]

            #make one hot encoding
            unary_marginals = np.zeros((n_nodes, self.n_states), dtype=np.int)
            gx = np.ogrid[:n_nodes]
            unary_marginals[gx, y] = 1

            ##accuemulated pairwise
            pw = np.dot(unary_marginals[edges[:, 0]].T,
                        unary_marginals[edges[:, 1]])

        unaries_acc = np.dot(unary_marginals.T, features)
        if self.directed:
            pw = pw.ravel()
        else:
            pw = compress_sym(pw)

        joint_feature_vector = np.hstack([unaries_acc.ravel(), pw])
        return joint_feature_vector

Example 24

Project: pystruct Source File: latent_node_crf.py
Function: joint_feature
    def joint_feature(self, x, y):
        """Feature vector associated with instance (x, y).

        Feature representation joint_feature, such that the energy of the
        configuration (x, y) and a weight vector w is given by
        np.dot(w, joint_feature(x, y)).

        Parameters
        ----------
        x : tuple
            Unary evidence.

        y : ndarray or tuple
            Either y is an integral ndarray, giving
            a complete labeling for x.
            Or it is the result of a linear programming relaxation. In this
            case, ``y=(unary_marginals, pariwise_marginals)``.

        Returns
        -------
        p : ndarray, shape (size_joint_feature,)
            Feature vector associated with state (x, y).

        """
        self._check_size_x(x)
        features, edges = self._get_features(x), self._get_edges(x)

        if isinstance(y, tuple):
            # y is result of relaxation, tuple of unary and pairwise marginals
            unary_marginals, pw = y
            # accuemulate pairwise
            pw = pw.reshape(-1, self.n_states, self.n_states).sum(axis=0)
        else:
            n_nodes = y.size
            gx = np.ogrid[:n_nodes]

            #make one hot encoding
            unary_marginals = np.zeros((n_nodes, self.n_states), dtype=np.int)
            gx = np.ogrid[:n_nodes]
            unary_marginals[gx, y] = 1

            ##accuemulated pairwise
            pw = np.dot(unary_marginals[edges[:, 0]].T,
                        unary_marginals[edges[:, 1]])
        n_visible = features.shape[0]
        unaries_acc = np.dot(unary_marginals[:n_visible,
                                             :self.n_input_states].T, features)

        joint_feature_vector = np.hstack([unaries_acc.ravel(),
                                          compress_sym(pw)])
        return joint_feature_vector

Example 25

Project: pystruct Source File: latent_node_crf.py
Function: joint_feature
    def joint_feature(self, x, y):
        """Feature vector associated with instance (x, y).

        Feature representation joint_feature, such that the energy of the
        configuration (x, y) and a weight vector w is given by
        np.dot(w, joint_feature(x, y)).

        Parameters
        ----------
        x : tuple
            Unary evidence.

        y : ndarray or tuple
            Either y is an integral ndarray, giving
            a complete labeling for x.
            Or it is the result of a linear programming relaxation. In this
            case, ``y=(unary_marginals, pariwise_marginals)``.

        Returns
        -------
        p : ndarray, shape (size_joint_feature,)
            Feature vector associated with state (x, y).

        """
        self._check_size_x(x)
        features, edges = self._get_features(x), self._get_edges(x)
        n_nodes = features.shape[0]
        edge_features = x[2]

        if isinstance(y, tuple):
            # y is result of relaxation, tuple of unary and pairwise marginals
            unary_marginals, pw = y
        else:
            n_nodes = y.size
            gx = np.ogrid[:n_nodes]

            # make one hot encoding
            unary_marginals = np.zeros((n_nodes, self.n_states), dtype=np.int)
            gx = np.ogrid[:n_nodes]
            unary_marginals[gx, y] = 1

            # pairwise
            pw = [np.outer(unary_marginals[edge[0]].T,
                           unary_marginals[edge[1]]).ravel()
                  for edge in edges]
            pw = np.vstack(pw)

        pw = np.dot(edge_features.T, pw)
        for i in self.symmetric_edge_features:
            pw_ = pw[i].reshape(self.n_states, self.n_states)
            pw[i] = (pw_ + pw_.T).ravel() / 2.

        for i in self.antisymmetric_edge_features:
            pw_ = pw[i].reshape(self.n_states, self.n_states)
            pw[i] = (pw_ - pw_.T).ravel() / 2.

        n_visible = features.shape[0]
        unaries_acc = np.dot(unary_marginals[:n_visible,
                                             :self.n_input_states].T, features)

        joint_feature_vector = np.hstack([unaries_acc.ravel(), pw.ravel()])
        return joint_feature_vector

Example 26

Project: scikit-image Source File: texture.py
def greycoprops(P, prop='contrast'):
    """Calculate texture properties of a GLCM.

    Compute a feature of a grey level co-occurrence matrix to serve as
    a compact summary of the matrix. The properties are computed as
    follows:

    - 'contrast': :math:`\\sum_{i,j=0}^{levels-1} P_{i,j}(i-j)^2`
    - 'dissimilarity': :math:`\\sum_{i,j=0}^{levels-1}P_{i,j}|i-j|`
    - 'humogeneity': :math:`\\sum_{i,j=0}^{levels-1}\\frac{P_{i,j}}{1+(i-j)^2}`
    - 'ASM': :math:`\\sum_{i,j=0}^{levels-1} P_{i,j}^2`
    - 'energy': :math:`\\sqrt{ASM}`
    - 'correlation':
        .. math:: \\sum_{i,j=0}^{levels-1} P_{i,j}\\left[\\frac{(i-\\mu_i) \\
                  (j-\\mu_j)}{\\sqrt{(\\sigma_i^2)(\\sigma_j^2)}}\\right]


    Parameters
    ----------
    P : ndarray
        Input array. `P` is the grey-level co-occurrence histogram
        for which to compute the specified property. The value
        `P[i,j,d,theta]` is the number of times that grey-level j
        occurs at a distance d and at an angle theta from
        grey-level i.
    prop : {'contrast', 'dissimilarity', 'humogeneity', 'energy', \
            'correlation', 'ASM'}, optional
        The property of the GLCM to compute. The default is 'contrast'.

    Returns
    -------
    results : 2-D ndarray
        2-dimensional array. `results[d, a]` is the property 'prop' for
        the d'th distance and the a'th angle.

    References
    ----------
    .. [1] The GLCM Tutorial Home Page,
           http://www.fp.ucalgary.ca/mhallbey/tutorial.htm

    Examples
    --------
    Compute the contrast for GLCMs with distances [1, 2] and angles
    [0 degrees, 90 degrees]

    >>> image = np.array([[0, 0, 1, 1],
    ...                   [0, 0, 1, 1],
    ...                   [0, 2, 2, 2],
    ...                   [2, 2, 3, 3]], dtype=np.uint8)
    >>> g = greycomatrix(image, [1, 2], [0, np.pi/2], levels=4,
    ...                  normed=True, symmetric=True)
    >>> contrast = greycoprops(g, 'contrast')
    >>> contrast
    array([[ 0.58333333,  1.        ],
           [ 1.25      ,  2.75      ]])

    """
    assert_nD(P, 4, 'P')

    (num_level, num_level2, num_dist, num_angle) = P.shape
    assert num_level == num_level2
    assert num_dist > 0
    assert num_angle > 0

    # create weights for specified property
    I, J = np.ogrid[0:num_level, 0:num_level]
    if prop == 'contrast':
        weights = (I - J) ** 2
    elif prop == 'dissimilarity':
        weights = np.abs(I - J)
    elif prop == 'humogeneity':
        weights = 1. / (1. + (I - J) ** 2)
    elif prop in ['ASM', 'energy', 'correlation']:
        pass
    else:
        raise ValueError('%s is an invalid property' % (prop))

    # compute property for each GLCM
    if prop == 'energy':
        asm = np.apply_over_axes(np.sum, (P ** 2), axes=(0, 1))[0, 0]
        results = np.sqrt(asm)
    elif prop == 'ASM':
        results = np.apply_over_axes(np.sum, (P ** 2), axes=(0, 1))[0, 0]
    elif prop == 'correlation':
        results = np.zeros((num_dist, num_angle), dtype=np.float64)
        I = np.array(range(num_level)).reshape((num_level, 1, 1, 1))
        J = np.array(range(num_level)).reshape((1, num_level, 1, 1))
        diff_i = I - np.apply_over_axes(np.sum, (I * P), axes=(0, 1))[0, 0]
        diff_j = J - np.apply_over_axes(np.sum, (J * P), axes=(0, 1))[0, 0]

        std_i = np.sqrt(np.apply_over_axes(np.sum, (P * (diff_i) ** 2),
                                           axes=(0, 1))[0, 0])
        std_j = np.sqrt(np.apply_over_axes(np.sum, (P * (diff_j) ** 2),
                                           axes=(0, 1))[0, 0])
        cov = np.apply_over_axes(np.sum, (P * (diff_i * diff_j)),
                                 axes=(0, 1))[0, 0]

        # handle the special case of standard deviations near zero
        mask_0 = std_i < 1e-15
        mask_0[std_j < 1e-15] = True
        results[mask_0] = 1

        # handle the standard case
        mask_1 = mask_0 == False
        results[mask_1] = cov[mask_1] / (std_i[mask_1] * std_j[mask_1])
    elif prop in ['contrast', 'dissimilarity', 'humogeneity']:
        weights = weights.reshape((num_level, num_level, 1, 1))
        results = np.apply_over_axes(np.sum, (P * weights), axes=(0, 1))[0, 0]

    return results

Example 27

Project: scipy Source File: measurements.py
Function: center_of_mass
def center_of_mass(input, labels=None, index=None):
    """
    Calculate the center of mass of the values of an array at labels.

    Parameters
    ----------
    input : ndarray
        Data from which to calculate center-of-mass. The masses can either
        be positive or negative.
    labels : ndarray, optional
        Labels for objects in `input`, as generated by `ndimage.label`.
        Only used with `index`.  Dimensions must be the same as `input`.
    index : int or sequence of ints, optional
        Labels for which to calculate centers-of-mass. If not specified,
        all labels greater than zero are used.  Only used with `labels`.

    Returns
    -------
    center_of_mass : tuple, or list of tuples
        Coordinates of centers-of-mass.

    Examples
    --------
    >>> a = np.array(([0,0,0,0],
    ...               [0,1,1,0],
    ...               [0,1,1,0],
    ...               [0,1,1,0]))
    >>> from scipy import ndimage
    >>> ndimage.measurements.center_of_mass(a)
    (2.0, 1.5)

    Calculation of multiple objects in an image

    >>> b = np.array(([0,1,1,0],
    ...               [0,1,0,0],
    ...               [0,0,0,0],
    ...               [0,0,1,1],
    ...               [0,0,1,1]))
    >>> lbl = ndimage.label(b)[0]
    >>> ndimage.measurements.center_of_mass(b, lbl, [1,2])
    [(0.33333333333333331, 1.3333333333333333), (3.5, 2.5)]

    Negative masses are also accepted, which can occur for example when
    bias is removed from measured data due to random noise.

    >>> c = np.array(([-1,0,0,0],
    ...               [0,-1,-1,0],
    ...               [0,1,-1,0],
    ...               [0,1,1,0]))
    >>> ndimage.measurements.center_of_mass(c)
    (-4.0, 1.0)

    If there are division by zero issues, the function does not raise an
    error but rather issues a RuntimeWarning before returning inf and/or NaN.

    >>> d = np.array([-1, 1])
    >>> ndimage.measurements.center_of_mass(d)
    (inf,)
    """
    normalizer = sum(input, labels, index)
    grids = numpy.ogrid[[slice(0, i) for i in input.shape]]

    results = [sum(input * grids[dir].astype(float), labels, index) / normalizer
               for dir in range(input.ndim)]

    if numpy.isscalar(results[0]):
        return tuple(results)

    return [tuple(v) for v in numpy.array(results).T]

Example 28

Project: VIP Source File: shapes.py
Function: get_circle
def get_circle(array, radius, output_values=False, cy=None, cx=None):           
    """Returns a centered circular region from a 2d ndarray. All the rest 
    pixels are set to zeros. 
    
    Parameters
    ----------
    array : array_like
        Input 2d array or image. 
    radius : int
        The radius of the circular region.
    output_values : {False, True}
        Sets the type of output.
    cy, cx : int
        Coordinates of the circle center.
        
    Returns
    -------
    values : array_like
        1d array with the values of the pixels in the circular region.
    array_masked : array_like
        Input array with the circular mask applied.
    """
    if not array.ndim == 2:
        raise TypeError('Input array is not a frame or 2d array.')
    sy, sx = array.shape
    if not cy and not cx:
        cy, cx = frame_center(array, verbose=False)
         
    yy, xx = np.ogrid[:sy, :sx]                                                 # ogrid is a multidim mesh creator (faster than mgrid)
    circle = (yy - cy)**2 + (xx - cx)**2                                        # eq of circle. squared distance to the center                                        
    circle_mask = circle < radius**2                                            # mask of 1's and 0's                                       
    if output_values:
        values = array[circle_mask]
        return values
    else:
        array_masked = array*circle_mask
        return array_masked

Example 29

Project: VIP Source File: shapes.py
def mask_circle(array, radius):                                      
    """ Masks (sets pixels to zero) a centered circle from a frame or cube. 
    
    Parameters
    ----------
    array : array_like
        Input frame or cube.
    radius : int
        Radius of the circular aperture.
    
    Returns
    -------
    array_masked : array_like
        Masked frame or cube.
        
    """
    if len(array.shape) == 2:
        sy, sx = array.shape
        cy = sy/2
        cx = sx/2
        xx, yy = np.ogrid[:sy, :sx]
        circle = (xx - cx)**2 + (yy - cy)**2               # squared distance to the center                                        
        hole_mask = circle > radius**2                                             
        array_masked = array*hole_mask
        
    if len(array.shape) == 3:
        n, sy, sx = array.shape
        cy = sy/2
        cx = sx/2
        xx, yy = np.ogrid[:sy, :sx]
        circle = (xx - cx)**2 + (yy - cy)**2               # squared distance to the center                                        
        hole_mask = circle > radius**2      
        array_masked = np.empty_like(array)
        for i in xrange(n):
            array_masked[i] = array[i]*hole_mask
        
    return array_masked    

Example 30

Project: mayavi Source File: test_mlab_integration.py
Function: test_volume_slice
    def test_volume_slice(self):
        x, y, z = np.ogrid[-5:5:20j, -5:5:20j, -5:5:20j]
        scalars = x * x * 0.5 + y * y + z * z * 2.0
        ipw = mlab.volume_slice(scalars, plane_orientation='y_axes')
        self.assertEqual(ipw.ipw.plane_orientation, 'y_axes')

Example 31

Project: mayavi Source File: test_mlab_source.py
    def test_reset1(self):
        "Test the reset method."
        x, y, z, v, s, src = self.get_data()
        self.check_traits()
        self.check_dataset()
        # Call reset again with just a few things changed to see if it
        # works correctly.

        x, y, z = N.ogrid[-10:10:11j,
                          -10:10:12j,
                          -10:10:20j]
        self.x, self.y, self.z = x, y, z

        dims = (x.shape[0], y.shape[1], z.shape[2])
        self.v = v = N.ones(dims + (3,), float)
        v[...,2] = 2
        v[...,2] = 3
        self.s = s = N.ones(dims, float)
        src = sources.MArraySource()
        src.reset(x=x, y=y, z=z, u=v[...,0], v=v[...,1], w=v[...,2], scalars=s,vectors=v)
        self.check_traits()
        self.check_dataset()

Example 32

Project: hyperspy Source File: test_samfire.py
def generate_test_model():

    # import hyperspy.api as hs
    from hyperspy.signals import Signal1D
    from hyperspy.components1d import (Gaussian, Lorentzian)
    import numpy as np
    from scipy.ndimage import gaussian_filter
    total = None
# blurs = [0., 0.5, 1., 2.,5.]
    blurs = [1.5]
    radius = 5
    domain = 15
# do circle/domain
    cent = (domain // 2, domain // 2)
    y, x = np.ogrid[-cent[0]:domain - cent[0], -cent[1]:domain - cent[1]]
    mask = x * x + y * y <= radius * radius
    lor_map = None
    for blur in blurs:

        s = Signal1D(np.ones((domain, domain, 1024)))
        cent = tuple([int(0.5 * i) for i in s.data.shape[:-1]])
        m0 = s.create_model()

        gs01 = Lorentzian()
        m0.append(gs01)
        gs01.gamma.map['values'][:] = 50
        gs01.gamma.map['is_set'][:] = True
        gs01.centre.map['values'][:] = 300
        gs01.centre.map['values'][mask] = 400
        gs01.centre.map['values'] = gaussian_filter(
            gs01.centre.map['values'],
            blur)
        gs01.centre.map['is_set'][:] = True
        gs01.A.map['values'][:] = 100 * \
            np.random.random((domain, domain)) + 300000
        gs01.A.map['values'][mask] *= 0.75
        gs01.A.map['values'] = gaussian_filter(gs01.A.map['values'], blur)
        gs01.A.map['is_set'][:] = True

        gs02 = Gaussian()
        m0.append(gs02)
        gs02.sigma.map['values'][:] = 15
        gs02.sigma.map['is_set'][:] = True
        gs02.centre.map['values'][:] = 400
        gs02.centre.map['values'][mask] = 300
        gs02.centre.map['values'] = gaussian_filter(
            gs02.centre.map['values'],
            blur)
        gs02.centre.map['is_set'][:] = True
        gs02.A.map['values'][:] = 50000
        gs02.A.map['is_set'][:] = True

        gs03 = Lorentzian()
        m0.append(gs03)
        gs03.gamma.map['values'][:] = 20
        gs03.gamma.map['is_set'][:] = True
        gs03.centre.map['values'][:] = 100
        gs03.centre.map['values'][mask] = 900
        gs03.centre.map['is_set'][:] = True
        gs03.A.map['values'][:] = 100 * \
            np.random.random((domain, domain)) + 50000
        gs03.A.map['values'][mask] *= 0.
        gs03.A.map['is_set'][:] = True

        s11 = m0.as_signal(show_progressbar=False)
        if total is None:
            total = s11.data.copy()
            lor_map = gs01.centre.map['values'].copy()
        else:
            total = np.concatenate((total, s11.data), axis=1)
            lor_map = np.concatenate(
                (lor_map, gs01.centre.map['values'].copy()), axis=1)

    s = Signal1D(total)
    s.add_poissonian_noise()
    s.data += 0.1
    s.estimate_poissonian_noise_variance()

    m = s.inav[:, :7].create_model()
    g = Gaussian()
    l1 = Lorentzian()
    l2 = Lorentzian()
    g.sigma.value = 50
    g.centre.value = 400
    g.A.value = 50000
    l1.gamma.value = 40
    l1.centre.value = 300
    l1.A.value = 300000
    l2.gamma.value = 15
    l2.centre.value = 100
    l2.A.value = 50000
    l2.centre.bmin = 0
    l2.centre.bmax = 200
    l2.A.bmin = 30000
    l2.A.bmax = 100000
    l2.gamma.bmin = 0
    l2.gamma.bmax = 60
    m.extend([g, l1, l2])
    m.assign_current_values_to_all()
    l2.active_is_multidimensional = True
    return m, gs01, gs02, gs03

Example 33

Project: invesalius3 Source File: cursor_actors.py
    def _build_actor(self):
        """
        Function to plot the circle
        """
        r = self.radius
        sx, sy, sz = self.spacing
        if self.orientation == 'AXIAL':
            xi = math.floor(-r/sx)
            xf = math.ceil(r/sx) + 1
            yi = math.floor(-r/sy)
            yf = math.ceil(r/sy) + 1
            zi = 0
            zf = 1
        elif self.orientation == 'CORONAL':
            xi = math.floor(-r/sx)
            xf = math.ceil(r/sx) + 1
            yi = 0
            yf = 1
            zi = math.floor(-r/sz)
            zf = math.ceil(r/sz) + 1
        elif self.orientation == 'SAGITAL':
            xi = 0
            xf = 1
            yi = math.floor(-r/sy)
            yf = math.ceil(r/sy) + 1
            zi = math.floor(-r/sz)
            zf = math.ceil(r/sz) + 1

        z,y,x = numpy.ogrid[zi:zf,yi:yf, xi:xf]

        circle_m = (z*sz)**2 + (y*sy)**2 + (x*sx)**2 <= r**2
        circle_i = to_vtk(circle_m.astype('uint8'),
                                          self.spacing, 0, self.orientation)
        circle_ci = self._set_colour(circle_i, self.colour)

        if self.mapper is None:
            self.actor.SetInputData(circle_ci)
            self.actor.InterpolateOff()
            self.actor.PickableOff()
            self.actor.SetDisplayExtent(circle_ci.GetExtent())
        else:
            self.mapper.SetInputData(circle_ci)
            self.mapper.BorderOn()

            self.mapper.SetOrientation(ORIENTATION[self.orientation])

Example 34

Project: invesalius3 Source File: cursor_actors.py
    def _calculate_area_pixels(self):
        """
        Return the cursor's pixels.
        """
        r = self.radius
        if self.orientation == 'AXIAL':
            sx = self.spacing[0]
            sy = self.spacing[1]
        elif self.orientation == 'CORONAL':
            sx = self.spacing[0]
            sy = self.spacing[2]
        elif self.orientation == 'SAGITAL':
            sx = self.spacing[1]
            sy = self.spacing[2]

        xi = math.floor(-r/sx)
        xf = math.ceil(r/sx) + 1
        yi = math.floor(-r/sy)
        yf = math.ceil(r/sy) + 1

        y,x = numpy.ogrid[yi:yf, xi:xf]

        index = (y*sy)**2 + (x*sx)**2 <= r**2
        self.points = index

Example 35

Project: neupy Source File: sofm.py
def neuron_neighbours(neurons, center, radius):
    """
    Function find all neighbours neurons by radius and coords.

    Parameters
    ----------
    neurons : arary-like
        Array element with neurons.
    center : tuple
        Index of the main neuron for which function must find neighbours.
    radius : int
        Radius indetify which neurons hear the main one are neighbours.

    Returns
    -------
    array-like
        Return matrix with the same dimension as ``neurons`` where center
        element and it neighbours positions filled with value ``1``
        and other as ``0``.

    Examples
    --------
    >>> import numpy as np
    >>> from neupy.algorithms.competitive.sofm import neuron_neighbours
    >>>
    >>> neuron_neighbours(np.zeros((3, 3)), (0, 0), 1)
    array([[ 1.,  1.,  0.],
           [ 1.,  0.,  0.],
           [ 0.,  0.,  0.]])
    >>> neuron_neighbours(np.zeros((5, 5)), (2, 2), 2)
    array([[ 0.,  0.,  1.,  0.,  0.],
           [ 0.,  1.,  1.,  1.,  0.],
           [ 1.,  1.,  1.,  1.,  1.],
           [ 0.,  1.,  1.,  1.,  0.],
           [ 0.,  0.,  1.,  0.,  0.]])
    """
    center_x, center_y = center
    nx, ny = neurons.shape

    y, x = np.ogrid[-center_x:nx - center_x, -center_y:ny - center_y]
    mask = (x * x + y * y) <= radius ** 2
    neurons[mask] = 1

    return neurons

Example 36

Project: pylinac Source File: mask.py
def sector_mask(shape, center, radius, angle_range=(0, 360)):
    """Return a circular arc-shaped boolean mask.

    Parameters
    ----------
    shape : tuple
        Shape of the image matrix. Usually easiest to pass something like array.shape
    center : Point, iterable
        The center point of the desired mask.
    radius : int, float
        Radius of the mask.
    angle_range : iterable
        The angle range of the mask. E.g. the default (0, 360) creates an entire circle.
        The start/stop angles should be given in clockwise order. 0 is right (0 on unit circle).

    References
    ----------
    https://stackoverflow.com/questions/18352973/mask-a-circular-sector-in-a-numpy-array/18354475#18354475
    """

    x, y = np.ogrid[:shape[0], :shape[1]]
    cy, cx = center.x, center.y
    # tmin, tmax = np.deg2rad(angle_range)
    tmin, tmax = angle_range

    # ensure stop angle > start angle
    if tmax < tmin:
        tmax += 2 * np.pi

    # convert cartesian --> polar coordinates
    r2 = (x - cx) * (x - cx) + (y - cy) * (y - cy)
    theta = np.arctan2(x - cx, y - cy) - tmin

    # wrap angles between 0 and 2*pi
    theta %= (2 * np.pi)

    # circular mask
    circmask = r2 <= radius * radius

    # angular mask
    anglemask = theta <= (tmax - tmin)

    return circmask * anglemask

Example 37

Project: artiq Source File: fit_image.py
Function: build
    def build(self, data, meta):
        super(Fit2DGaussParabola, self).build(data, meta)
        self.xy = np.ogrid[:data.shape[0], :data.shape[1]]

Example 38

Project: minivect Source File: lazy_numpy.py
Function: test_2
def test2():
    N = 200
    i, j, k = np.ogrid[:N, :N, :N]
    dtype = np.double
    i, j, k = i.astype(dtype), j.astype(dtype), k.astype(dtype)

    numpy_result = np.empty((N, N, N), dtype=dtype)

    t = time.time()
    numpy_result[...] = i * j * k
    print time.time() - t

    our_result = np.empty((N, N, N), dtype=dtype)

    # Lazy evaluation slice assignment
#    t = time.time()
    lazy_dst = lazy_array(our_result)
    lazy_i, lazy_j, lazy_k = lazy_array(i), lazy_array(j), lazy_array(k)
#    lazy_dst[...] = lazy_i * lazy_j * lazy_k
#    print time.time() - t
#
#    assert np.all(numpy_result == our_result)

    # Lazy evaluation with compilation separate from timing
    f, a = lazy_dst.slice_assign(lazy_i * lazy_j * lazy_k)
    t = time.time()
    for i in range(10):
        f(*a)
    print time.time() - t

Example 39

Project: openmc Source File: triso.py
Function: classify
    def classify(self, lattice):
        """Determine lattice element indices which might contain the TRISO particle.

        Parameters
        ----------
        lattice : openmc.RectLattice
            Lattice to check

        Returns
        -------
        list of tuple
            (z,y,x) lattice element indices which might contain the TRISO
            particle.

        """

        ll, ur = self.region.bounding_box
        if lattice.ndim == 2:
            (i_min, j_min), p = lattice.find_element(ll)
            (i_max, j_max), p = lattice.find_element(ur)
            return list(np.broadcast(*np.ogrid[
                j_min:j_max+1, i_min:i_max+1]))
        else:
            (i_min, j_min, k_min), p = lattice.find_element(ll)
            (i_max, j_max, k_max), p = lattice.find_element(ur)
            return list(np.broadcast(*np.ogrid[
                k_min:k_max+1, j_min:j_max+1, i_min:i_max+1]))

Example 40

Project: openmc Source File: triso.py
def create_triso_lattice(trisos, lower_left, pitch, shape, background):
    """Create a lattice containing TRISO particles for optimized tracking.

    Parameters
    ----------
    trisos : list of openmc.model.TRISO
        List of TRISO particles to put in lattice
    lower_left : Iterable of float
        Lower-left Cartesian coordinates of the lattice
    pitch : Iterable of float
        Pitch of the lattice elements in the x-, y-, and z-directions
    shape : Iterable of float
        Number of lattice elements in the x-, y-, and z-directions
    background : openmc.Material
        A background material that is used anywhere within the lattice but
        outside a TRISO particle

    Returns
    -------
    lattice : openmc.RectLattice
        A lattice containing the TRISO particles

    """

    lattice = openmc.RectLattice()
    lattice.lower_left = lower_left
    lattice.pitch = pitch

    indices = list(np.broadcast(*np.ogrid[:shape[2], :shape[1], :shape[0]]))
    triso_locations = {idx: [] for idx in indices}
    for t in trisos:
        for idx in t.classify(lattice):
            if idx in sorted(triso_locations):
                # Create copy of TRISO particle with materials preserved and
                # different cell/surface IDs
                t_copy = copy.deepcopy(t)
                t_copy.id = None
                t_copy.fill = t.fill
                t_copy._surface.id = None
                triso_locations[idx].append(t_copy)
            else:
                warnings.warn('TRISO particle is partially or completely '
                              'outside of the lattice.')

    # Create universes
    universes = np.empty(shape[::-1], dtype=openmc.Universe)
    for idx, triso_list in sorted(triso_locations.items()):
        if len(triso_list) > 0:
            outside_trisos = openmc.Intersection(*[~t.region for t in triso_list])
            background_cell = openmc.Cell(fill=background, region=outside_trisos)
        else:
            background_cell = openmc.Cell(fill=background)

        u = openmc.Universe()
        u.add_cell(background_cell)
        for t in triso_list:
            u.add_cell(t)
            iz, iy, ix = idx
            t.center = lattice.get_local_coordinates(t.center, (ix, iy, iz))

        if len(shape) == 2:
            universes[-1 - idx[0], idx[1]] = u
        else:
            universes[idx[0], -1 - idx[1], idx[2]] = u
    lattice.universes = universes

    # Set outer universe
    background_cell = openmc.Cell(fill=background)
    lattice.outer = openmc.Universe(cells=[background_cell])

    return lattice

Example 41

Project: AWS-Lambda-ML-Microservice-Skeleton Source File: measurements.py
Function: center_of_mass
def center_of_mass(input, labels=None, index=None):
    """
    Calculate the center of mass of the values of an array at labels.

    Parameters
    ----------
    input : ndarray
        Data from which to calculate center-of-mass.
    labels : ndarray, optional
        Labels for objects in `input`, as generated by `ndimage.label`.
        Only used with `index`.  Dimensions must be the same as `input`.
    index : int or sequence of ints, optional
        Labels for which to calculate centers-of-mass. If not specified,
        all labels greater than zero are used.  Only used with `labels`.

    Returns
    -------
    center_of_mass : tuple, or list of tuples
        Coordinates of centers-of-mass.

    Examples
    --------
    >>> a = np.array(([0,0,0,0],
    ...               [0,1,1,0],
    ...               [0,1,1,0],
    ...               [0,1,1,0]))
    >>> from scipy import ndimage
    >>> ndimage.measurements.center_of_mass(a)
    (2.0, 1.5)

    Calculation of multiple objects in an image

    >>> b = np.array(([0,1,1,0],
    ...               [0,1,0,0],
    ...               [0,0,0,0],
    ...               [0,0,1,1],
    ...               [0,0,1,1]))
    >>> lbl = ndimage.label(b)[0]
    >>> ndimage.measurements.center_of_mass(b, lbl, [1,2])
    [(0.33333333333333331, 1.3333333333333333), (3.5, 2.5)]

    """
    normalizer = sum(input, labels, index)
    grids = numpy.ogrid[[slice(0, i) for i in input.shape]]

    results = [sum(input * grids[dir].astype(float), labels, index) / normalizer
               for dir in range(input.ndim)]

    if numpy.isscalar(results[0]):
        return tuple(results)

    return [tuple(v) for v in numpy.array(results).T]

Example 42

Project: dipy Source File: test_interpolate.py
def test_TriLinearInterpolator():
    # Place (0, 0, 0) at the bottom left of the image
    l, m, n, o = np.ogrid[.5:6.51, .5:6.51, .5:6.51, 0:4]
    data = l + m + n + o
    data = data.astype("float32")

    tli = TriLinearInterpolator(data, (1, 1, 1))
    a, b, c = np.mgrid[.5:6.5:1.6, .5:6.5:2.7, .5:6.5:3.8]
    for ii in xrange(a.size):
        x = a.flat[ii]
        y = b.flat[ii]
        z = c.flat[ii]
        expected_result = x + y + z + o.ravel()
        assert_array_almost_equal(tli[x, y, z], expected_result, decimal=5)
        ind = np.array([x, y, z])
        assert_array_almost_equal(tli[ind], expected_result)

    # Index at 0
    expected_value = np.arange(4) + 1.5
    assert_array_almost_equal(tli[0, 0, 0], expected_value)
    # Index at shape
    expected_value = np.arange(4) + (6.5 * 3)
    assert_array_almost_equal(tli[7, 7, 7], expected_value)

    assert_raises(OutsideImage, tli.__getitem__, (-.1, 0, 0))
    assert_raises(OutsideImage, tli.__getitem__, (0, 7.01, 0))

Example 43

Project: nipy Source File: array_coords.py
Function: get_item
    def __getitem__(self, index):
        """
        Create an AffineTransform coordinate map with into self.coords with
        slices created as in np.mgrid/np.ogrid.
        """
        dtype = self.coordmap.function_domain.coord_dtype
        results = [a.ravel().astype(dtype) for a in np.ogrid[index]]
        if len(results) != len(self.coordmap.function_domain.coord_names):
            raise ValueError('the number of slice objects must match ' 
                             'the number of input dimensions')
        cmaps = []
        for i, result in enumerate(results):
            if result.shape[0] > 1:
                step = result[1] - result[0]
            else:
                step = 0
            start = result[0]
            cmaps.append(AffineTransform(
                    CoordinateSystem(['i%d' % i], coord_dtype=dtype),
                    CoordinateSystem([self.coordmap.function_domain.coord_names[i]],
                                     coord_dtype=dtype),
                    np.array([[step, start],[0,1]], dtype=dtype)))

        shape = [result.shape[0] for result in results]
        cmap = cmap_product(*cmaps)

        # Identify the origin in the range of cmap
        # with the origin in the domain of self.coordmap

        cmap = shifted_range_origin(cmap, 
                                    np.zeros(cmap.ndims[1]),
                                    self.coordmap.function_domain.name)

        return ArrayCoordMap(compose(self.coordmap, cmap), tuple(shape))