Python numpy.broadcast函数代码示例

 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]]))

    def _validate_mask_and_error(self):
        """
        Raises ValueError if they don't match (using ~numpy.broadcast)
        """

        try:
            if self.mask is not None:
                np.broadcast(self.data, self.mask)
            maskmatch = True
        except ValueError:
            maskmatch = False

        try:
            if self.error is not None:
                np.broadcast(self.data, self.error)
            errmatch = True
        except ValueError:
            errmatch = False

        if not errmatch and not maskmatch:
            raise ValueError('NDData error and mask do not match data')
        elif not errmatch:
            raise ValueError('NDData error does not match data')
        elif not maskmatch:
            raise ValueError('NDData mask does not match data')

    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]))

def check_mean_sigma_keepdims(a, axis):
    mu1, sigma1 = mean_sigma(a, axis, keepdims=False)
    mu2, sigma2 = mean_sigma(a, axis, keepdims=True)

    assert_array_equal(mu1.ravel(), mu2.ravel())
    assert_array_equal(sigma1.ravel(), sigma2.ravel())

    assert_array_equal(np.broadcast(a, mu2).shape, a.shape)
    assert_array_equal(np.broadcast(a, sigma2).shape, a.shape)

 def _set_error(self, value):
     if value is not None:
         try:
             np.broadcast(self.data, value)
         except ValueError:
             raise ValueError("dimensions of `error` do not match data")
         else:
             self._error = value
     else:
         self._error = value

def logsubtrexp(minuend, subtrahend, sign_minuend = None, sign_subtrahend = None):

    if sign_minuend is None:
        sign_minuend = np.ones(minuend.shape)
    if sign_subtrahend is None:
        sign_subtrahend = np.ones(subtrahend.shape)
    if not (minuend.shape == sign_minuend.shape and subtrahend.shape == sign_subtrahend.shape):
        raise ValueError("sign arguments expected be of same shape as corresponding log-matrices")
    if not (np.abs(sign_minuend).all() and np.abs(sign_subtrahend).all()):
        raise ValueError("sign arguments expected to contain only +1 or -1 elements")
        
    b = np.broadcast(minuend, subtrahend)
    s_b = np.broadcast(sign_minuend, sign_subtrahend)
    abs_res = np.empty(b.shape)
    sign_res = np.empty(b.shape)
    
    for i in range(b.size):
        (m, s) = b.next()
        (sign_m, sign_s) = s_b.next()
        if sign_m > sign_s: # sign_m == 1 and sign_s == -1
            # this is equivalent to logsumexp(m, s)
            #print("sign_m > sign_s")
            sign_res.flat[i] = 1
            abs_res.flat[i] = logsumexp((m,s))
        elif sign_m < sign_s: # sign_m == -1 and sign_s == 1
            #print("sign_m < sign_s")
            sign_res.flat[i] = -1
            abs_res.flat[i] = logsumexp((m,s))
        else:
            #signs are eqal
            if m == s:                
                sign_res.flat[i] = 1
                abs_res.flat[i] = log(0)
            else:
                if sign_m == -1:
                    if m > s:
                        #print("m >= s")
                        sign_res.flat[i] = -1
                        abs_res.flat[i] = log(1 - exp(s - m)) + m
                    elif m < s:
                        #print("m < s")
                        sign_res.flat[i] = 1
                        abs_res.flat[i] = log(1 - exp(m - s)) + s
                else:# sign_m == 1
                    if m > s:
                        #print("m >= s")
                        sign_res.flat[i] = 1
                        abs_res.flat[i] = log(1 - exp(s - m)) + m
                    elif m < s:
                        #print("m < s")
                        sign_res.flat[i] = -1
                        abs_res.flat[i] = log(1 - exp(m - s)) + s
        #print(sign_m*exp(m),  sign_s*exp(s),  sign_m*exp(m) - sign_s*exp(s), sign_res.flat[i] * exp(abs_res.flat[i]))
    
    return (abs_res, sign_res)

    def getPsfAtPoints(self, bandnum, x, y):
        '''
        Reconstruct the SDSS model PSF from KL basis functions.

        x,y can be scalars or 1-d numpy arrays.

        Return value:
        if x,y are scalars: a PSF image
        if x,y are arrays:  a list of PSF images
        '''
        rtnscalar = np.isscalar(x) and np.isscalar(y)
        x = np.atleast_1d(x)
        y = np.atleast_1d(y)
        psf = fits_table(self.hdus[bandnum+1].data)
        psfimgs = None
        (outh, outw) = (None,None)

        # From the IDL docs:
        # http://photo.astro.princeton.edu/photoop_doc.html#SDSS_PSF_RECON
        #   acoeff_k = SUM_i{ SUM_j{ (0.001*ROWC)^i * (0.001*COLC)^j * C_k_ij } }
        #   psfimage = SUM_k{ acoeff_k * RROWS_k }
        for k in range(len(psf)):
            nrb = psf.nrow_b[k]
            ncb = psf.ncol_b[k]
            c = psf.c[k].reshape(5, 5)
            c = c[:nrb,:ncb]
            (gridi,gridj) = np.meshgrid(range(nrb), range(ncb))

            if psfimgs is None:
                psfimgs = [np.zeros_like(psf.rrows[k]) for xy
                        in np.broadcast(x,y)]
                (outh,outw) = (psf.rnrow[k], psf.rncol[k])
            else:
                assert(psf.rnrow[k] == outh)
                assert(psf.rncol[k] == outw)

            for i,(xi,yi) in enumerate(np.broadcast(x,y)):
                #print 'xi,yi', xi,yi
                acoeff_k = np.sum(((0.001 * xi)**gridi * (0.001 * yi)**gridj * c))
                if False: # DEBUG
                    print 'coeffs:', (0.001 * xi)**gridi * (0.001 * yi)**gridj
                    print 'c:', c
                    for (coi,ci) in zip(((0.001 * xi)**gridi * (0.001 * yi)**gridj).ravel(), c.ravel()):
                        print 'co %g, c %g' % (coi,ci)
                    print 'acoeff_k', acoeff_k

                #print 'acoeff_k', acoeff_k.shape, acoeff_k
                #print 'rrows[k]', psf.rrows[k].shape, psf.rrows[k]
                psfimgs[i] += acoeff_k * psf.rrows[k]

        psfimgs = [img.reshape((outh,outw)) for img in psfimgs]
        if rtnscalar:
            return psfimgs[0]
        return psfimgs

def test_mean_sigma_keepdims(axis):
    np.random.seed(0)
    a = np.random.random((4, 5, 6))
    mu1, sigma1 = mean_sigma(a, axis, keepdims=False)
    mu2, sigma2 = mean_sigma(a, axis, keepdims=True)

    assert_array_equal(mu1.ravel(), mu2.ravel())
    assert_array_equal(sigma1.ravel(), sigma2.ravel())

    assert_array_equal(np.broadcast(a, mu2).shape, a.shape)
    assert_array_equal(np.broadcast(a, sigma2).shape, a.shape)

 def _assert_incompatible_broadcast(self, shape1, shape2):
   if shape1.dims is not None and shape2.dims is not None:
     zeros1 = np.zeros(shape1.as_list())
     zeros2 = np.zeros(shape2.as_list())
     with self.assertRaises(ValueError):
       np.broadcast(zeros1, zeros2)
     with self.assertRaises(ValueError):
       np.broadcast(zeros2, zeros1)
   with self.assertRaises(ValueError):
     common_shapes.broadcast_shape(shape1, shape2)
   with self.assertRaises(ValueError):
     common_shapes.broadcast_shape(shape2, shape1)

def test_broadcastable():
    for ndim1 in range(1, 4):
        for ndim2 in range(1, 4):
            for shape1 in itertools.permutations(range(1, 4), ndim1):
                for shape2 in itertools.permutations(range(1, 4), ndim2):
                    try:
                        np.broadcast(np.zeros(shape1),
                                     np.zeros(shape2))
                        result = True
                    except ValueError:
                        result = False
                    assert result == broadcastable(shape1, shape2)

 def _set_flags(self, value):
     if value is not None:
         if isinstance(value, np.ndarray):
             try:
                 np.broadcast(self.data, value)
             except ValueError:
                 raise ValueError("dimensions of `flags` do not match data")
             else:
                 self._flags = value
         else:
             raise TypeError("`flags` should be a Numpy array")
     else:
         self._flags = value

    def __init__(self, w, mu, *args, **kwargs):
        _, sd = get_tau_sd(tau=kwargs.pop('tau', None),
                           sd=kwargs.pop('sd', None))

        distshape = np.broadcast(mu, sd).shape
        self.mu = mu = tt.as_tensor_variable(mu)
        self.sd = sd = tt.as_tensor_variable(sd)

        if not distshape:
            distshape = np.broadcast(mu.tag.test_value, sd.tag.test_value).shape

        super(NormalMixture, self).__init__(w, Normal.dist(mu, sd=sd, shape=distshape),
                                            *args, **kwargs)

 def near(x, l, h):
     #RESULT: array(3) & array([3, 4])
     a = broadcast(x,l,h)
     has_iterable = a.shape
     x,l,h = tuple(atleast_1d(i) for i in (x,l,h))
     b = broadcast(x,l,h)
     _x,_l,_h = (empty(b.shape), empty(b.shape), empty(b.shape))
     _x.flat = [i for i in x]
     _l.flat = [i for i in l]
     _h.flat = [i for i in h]
     result = asarray(map(_near, x, l, h))
     if not has_iterable:
         result = asarray(result[0])
     return result

    def _natural_indices(self):
        """Iterate over all possible (x,y) or (x,y,z) lattice element indices.

        This property is used when constructing distributed cell and material
        paths. Most importantly, the iteration order matches that used on the
        Fortran side.

        """
        if self.ndim == 2:
            nx, ny = self.shape
            return np.broadcast(*np.ogrid[:nx, :ny])
        else:
            nx, ny, nz = self.shape
            return np.broadcast(*np.ogrid[:nx, :ny, :nz])

 def _assert_broadcast(self, expected, shape1, shape2):
   if shape1.dims is not None and shape2.dims is not None:
     expected_np = expected.as_list()
     zeros1 = np.zeros(shape1.as_list())
     zeros2 = np.zeros(shape2.as_list())
     self.assertAllEqual(expected_np, np.broadcast(zeros1, zeros2).shape)
     self.assertAllEqual(expected_np, np.broadcast(zeros2, zeros1).shape)
     self.assertEqual(
         expected, common_shapes.broadcast_shape(shape1, shape2))
     self.assertEqual(
         expected, common_shapes.broadcast_shape(shape2, shape1))
   else:
     self.assertEqual(expected, common_shapes.broadcast_shape(shape1, shape2))
     self.assertEqual(expected, common_shapes.broadcast_shape(shape2, shape1))