Python scipy.asarray函数代码示例

    def test3_raster_data_from_array(self):
        # A test based on this info;
        # http://en.wikipedia.org/wiki/Esri_grid
        # Let's hope no one edits the data....
        raster = [[-9999, -9999, 5, 2], [-9999, 20, 100, 36],
                  [3, 8, 35, 10], [32, 42, 50, 6],
                  [88, 75, 27, 9], [13, 5, 1, -9999]]
        upper_left_x = 0.
        upper_left_y = 300.
        cell_size = 50.0
        no_data_value = -9999

        # Just outside the midpoint of all sides
        lon = asarray([125, 125, 125, 125, 125, 125])
        lat = asarray([275, 225, 175, 125, 75, 25])

        raster = Raster.from_array(raster, upper_left_x, upper_left_y,
                                   cell_size, no_data_value)
        self.assertEqual(raster.ul_x, 0)
        self.assertEqual(raster.ul_y, 300)
        self.assertEqual(raster.x_pixel, 50)
        self.assertEqual(raster.y_pixel, -50)
        self.assertEqual(raster.x_size, 4)
        self.assertEqual(raster.y_size, 6)

        data = raster.raster_data_at_points(lon, lat)
        self.assertTrue(allclose(data, asarray([5.0, 100.0, 35.0,
                                                50.0, 27.0, 1.0])))

        # testing extent
        min_long, min_lat, max_long, max_lat = raster.extent()
        self.assertEqual(min_long, 0)
        self.assertEqual(min_lat, 0)
        self.assertEqual(max_long, 200)
        self.assertEqual(max_lat, 300)

def unique_rows(arr):
    """Returns a copy of arr with duplicate rows removed.
    
    From Stackoverflow "Find unique rows in numpy.array."
    
    Parameters
    ----------
    arr : :py:class:`Array`, (`m`, `n`). The array to find the unique rows of.
    
    Returns
    -------
    unique : :py:class:`Array`, (`p`, `n`) where `p` <= `m`
        The array `arr` with duplicate rows removed.
    """
    b = scipy.ascontiguousarray(arr).view(
        scipy.dtype((scipy.void, arr.dtype.itemsize * arr.shape[1]))
    )
    try:
        dum, idx = scipy.unique(b, return_index=True)
    except TypeError:
        # Handle bug in numpy 1.6.2:
        rows = [_Row(row) for row in b]
        srt_idx = sorted(range(len(rows)), key=rows.__getitem__)
        rows = scipy.asarray(rows)[srt_idx]
        row_cmp = [-1]
        for k in xrange(1, len(srt_idx)):
            row_cmp.append(rows[k-1].__cmp__(rows[k]))
        row_cmp = scipy.asarray(row_cmp)
        transition_idxs = scipy.where(row_cmp != 0)[0]
        idx = scipy.asarray(srt_idx)[transition_idxs]
    return arr[idx]

 def frombounds(
     cls, func, lbound, ubound, npop, crossover_rate=0.5, scale=None, strategy=("rand", 2, "bin"), eps=1e-6
 ):
     lbound = sp.asarray(lbound)
     ubound = sp.asarray(ubound)
     pop0 = rand(npop, len(lbound)) * (ubound - lbound) + lbound
     return cls(func, pop0, crossover_rate=crossover_rate, scale=scale, strategy=strategy, eps=eps)

def triplot(vertices, indices, labels=False):
    """
    Plot a 2D triangle mesh
    """
    
    vertices,indices = asarray(vertices),asarray(indices)

    #3d tensor [triangle index][vertex index][x/y value]
    triangles = vertices[numpy.ravel(indices),:].reshape((indices.shape[0],3,2))
    
    col = matplotlib.collections.PolyCollection(triangles)
    col.set_facecolor('grey')
    col.set_alpha(0.5)
    col.set_linewidth(1)

    #sub =  subplot(111)
    sub = matplotlib.pylab.gca()
    sub.add_collection(col,autolim=True)
    matplotlib.pylab.axis('off')
    sub.autoscale_view()

    if labels:
        barycenters = numpy.average(triangles,axis=1)
        for n,bc in enumerate(barycenters):
            matplotlib.pylab.text(bc[0], bc[1], str(n), {'color' : 'k', 'fontsize' : 8,
                                                         'horizontalalignment' : 'center',
                                                         'verticalalignment' : 'center'
                                                         })

def compute_mean_vector(category_name, labellist, layer = 'fc8'):
    print category_name
    featurefile_list = glob.glob('%s/%s/*.mat' %(featurefilepath, category_name))
    
    # gather all the training samples for which predicted category
    # was the category under consideration
    correct_features = []
    for featurefile in featurefile_list:
        try:
            img_arr = loadmat(featurefile)
            predicted_category = labellist[img_arr['scores'].argmax()]
            if predicted_category == category_name:
                correct_features += [img_arr[layer]]
        except TypeError:
            continue
    
    # Now compute channel wise mean vector
    channel_mean_vec = []
    for channelid in range(correct_features[0].shape[0]):
        channel = []
        for feature in correct_features:
            channel += [feature[channelid, :]]
        channel = sp.asarray(channel)
        assert len(correct_features) == channel.shape[0]
        # Gather mean over each channel, to get mean channel vector
        channel_mean_vec += [sp.mean(channel, axis=0)]

    # this vector contains mean computed over correct classifications
    # for each channel separately
    channel_mean_vec = sp.asarray(channel_mean_vec)
    savemat('%s.mat' %category_name, {'%s'%category_name: channel_mean_vec})

def ZYFF(Te, EIJ):
    """Computes `ZY` and `FF`, used in other functions.
    
    If `EIJ` is a scalar, the output has the same shape as `Te`. If `EIJ` is an
    array, the output has shape `EIJ.shape` + `Te.shape`. This should keep the
    output broadcastable with `Te`.
    
    Parameters
    ----------
    Te : array of float
        Electron temperature. Shape is arbitrary.
    EIJ : scalar float or array of float
        Energy difference.
    """
    # Expand the dimensions of EIJ to produce the desired output shape:
    Te = scipy.asarray(Te, dtype=float)
    EIJ = scipy.asarray(EIJ, dtype=float)
    for n in xrange(Te.ndim):
        EIJ = scipy.expand_dims(EIJ, axis=-1)
    
    ZY = EIJ / (1e3 * Te)
    
    FF = scipy.zeros_like(ZY)
    mask = (ZY >= 1.5)
    FF[mask] = scipy.log((ZY[mask] + 1) / ZY[mask]) - (0.36 + 0.03 * scipy.sqrt(ZY[mask] + 0.01)) / (ZY[mask] + 1)**2
    mask = ~mask
    FF[mask] = scipy.log((ZY[mask] + 1) / ZY[mask]) - (0.36 + 0.03 / scipy.sqrt(ZY[mask] + 0.01)) / (ZY[mask] + 1)**2
    
    return ZY, FF

def computeOpenMaxProbability(openmax_fc8, openmax_score_u):
    """ Convert the scores in probability value using openmax
    
    Input:
    ---------------
    openmax_fc8 : modified FC8 layer from Weibull based computation
    openmax_score_u : degree

    Output:
    ---------------
    modified_scores : probability values modified using OpenMax framework,
    by incorporating degree of uncertainity/openness for a given class
    
    """
    prob_scores, prob_unknowns = [], []
    for channel in range(NCHANNELS):
        channel_scores, channel_unknowns = [], []
        for category in range(NCLASSES):
            channel_scores += [sp.exp(openmax_fc8[channel, category])]
                    
        total_denominator = sp.sum(sp.exp(openmax_fc8[channel, :])) + sp.exp(sp.sum(openmax_score_u[channel, :]))
        prob_scores += [channel_scores/total_denominator ]
        prob_unknowns += [sp.exp(sp.sum(openmax_score_u[channel, :]))/total_denominator]
        
    prob_scores = sp.asarray(prob_scores)
    prob_unknowns = sp.asarray(prob_unknowns)

    scores = sp.mean(prob_scores, axis = 0)
    unknowns = sp.mean(prob_unknowns, axis=0)
    modified_scores =  scores.tolist() + [unknowns]
    assert len(modified_scores) == 1001
    return modified_scores

    def test_Epicentral(self):
        dist = Distances(None,None,None,None,None,None,None,None,None,None,None)
        
        distance_type='Epicentral'

        rupture_centroid_lat=asarray((-31.0))
        rupture_centroid_lon=asarray((116.0))
        
        site_lat=asarray((-31,-32,-33,-34))        
        site_lon=asarray((116.0,116.0,116.0,116.0))

        distance=dist.raw_distances(site_lat,
                                    site_lon,
                                    rupture_centroid_lat,
                                    rupture_centroid_lon,
                                    lengths,
                                    azimuths,
                                    widths,
                                    dips,
                                    depths,
                                    depths_to_top,
                                    distance_type,
                                    projection)


        d=asarray((0,1,2,3))*(1.852*60)
        d=d[:,newaxis]
        d2=As_The_Cockey_Flies(rupture_centroid_lat,
                               rupture_centroid_lon,
                               site_lat,
                               site_lon)

        assert allclose(d,distance,rtol=0.001)
        assert allclose(d,d2)

    def test_Epicentral2(self):
        dist = Distances(None,None,None,None,None,None,None,None,None,None,None)
        
        distance_type='Epicentral'
        rupture_centroid_lat=asarray((-31.0,-33,-36))
        rupture_centroid_lon=asarray((116.0,118,222))
        

        site_lat=asarray((-31,-32,-33,-34))        
        site_lon=asarray((116.0,116.0,116.0,116.0))

        distance=dist.raw_distances(site_lat,
                                    site_lon,
                                    rupture_centroid_lat,
                                    rupture_centroid_lon,
                                    lengths,
                                    azimuths,
                                    widths,
                                    dips,
                                    depths,
                                    depths_to_top,
                                    distance_type,
                                    projection)
        

        
        d = As_The_Cockey_Flies(rupture_centroid_lat,
                                rupture_centroid_lon,
                                site_lat,
                                site_lon)

        assert allclose(d,distance,rtol=0.1)

def draw_regions(original_image, mask_image, outline=["green", "blue", "red"],
                 fill=[None, None, None], mode_pass=[False, False, False]):
    heigth = original_image.size[1]
    regions = find_regions(mask_image)
    filtered_regions = filter_regions(regions,
                                      noise_heigth=NOISE_FACTOR * heigth)
    out = original_image.convert("RGB")
    for i, region in enumerate(filtered_regions):
        if mode_pass[i]:
            widths = scipy.asarray([f.width() for f in filtered_regions[i]])
            heigths = scipy.asarray([f.heigth() for f in filtered_regions[i]])
            mode_width = float(scipy.stats.mstats.mode(widths)[0])
            mode_heigth = float(scipy.stats.mstats.mode(heigths)[0])
        if outline[i] or fill[i]:
            draw = ImageDraw.Draw(out)
            for r in filtered_regions[i]:
                if (not mode_pass[i]
                    or (mode_pass[i]
                        and (mode_width - mode_pass[i] <= r.width() \
                             <= mode_width + mode_pass[i]
                             or mode_heigth - mode_pass[i] <= r.heigth() \
                             <= mode_heigth + mode_pass[i]))):
                    draw.rectangle(r.box(), outline=outline[i], fill=fill[i])
            del draw
    return out

 def test_As_The_Cockey_Flies(self):
     # Test data from GA website 
     # As_The_Cockey_Flies implements the Great Circle method
     # 
     # http://www.ga.gov.au/earth-monitoring/geodesy/geodetic-techniques/distance-calculation-algorithms.html
     rupture_centroid_lat = asarray((-30))
     rupture_centroid_lon = asarray((150))
     
     site_lat = asarray((-31,-31,-32,-33,-34,-35,-40,-50,-60,-70,-80))
     site_lon = asarray((150,151,151,151,151,151,151,151,151,151,151))
     
     expected = asarray([[111.120],
                         [146.677],
                         [241.787],
                         [346.556],
                         [454.351],
                         [563.438],
                         [1114.899],
                         [2223.978],
                         [3334.440],
                         [4445.247],
                         [5556.190]])
     
     d = As_The_Cockey_Flies(rupture_centroid_lat,
                             rupture_centroid_lon,
                             site_lat,
                             site_lon)
     
     assert allclose(d,expected)

    def ref_indicator(self, coord):
        """
        Return the value of the indicator for the reference coordinates if appropriate.

        NOTE: Currently we will simply implement the reference discretization as non-overlapping boxes.
        """
        # create a distance vector
        distancevec = sp.asarray(coord) - sp.asarray(self.ref_center)

        # if any collective variable is periodic, construct dr, the adjuct for minimum image convetion for the periodic cv's
        if self.wrapping is not None:

            # build dr
            dr = np.zeros(distancevec.shape)

            # add values to dr if the CV wraps
            for i in xrange(len(self.wrapping)):
                if self.wrapping[i] != 0.0:
                    # This is an old trick from MD codes to find the minimum distance between two points.
                    dr[i] = self.wrapping[i] * np.rint(distancevec[i]/self.wrapping[i])

            # add min image vector
            distancevec -= dr

        # We return 1.0 if all the distances are smaller than the width of the box from the center, 0.0 otherwise.
        return float(np.prod(self.ref_width > np.abs(distancevec)))

def xcorrv(a, b=None, lag=None, dtype=None):
    """vectorial cross correlation by taking the expectation over an outer product"""

    # checks
    a = sp.asarray(a)
    b = sp.asarray(b or a)
    if not (a.ndim == b.ndim):
        raise ValueError('a.ndim !== b.ndim')

    #if a.size != b.size:
    #    raise ValueError('a.size != b.size')
    #if a.size < 2:
    #    raise ValueError('a.size < 2')

    if lag is None:
        lag = int(a.shape[0] - 1)
    if lag > a.shape[0] - 1:
        raise ValueError('lag > vector len - 1')

    # init
    lag_range = xrange(int(-lag), int(lag) + 1)
    rval = sp.empty((a.shape[1], b.shape[1], len(lag_range)), dtype=dtype or a.dtype)

    # calc
    for tau in lag_range:
        prod = a.T[:, None, max(0, +tau):min(len(a), len(a) + tau)] * \
               b.T[None, :, max(0, -tau):min(len(b), len(b) - tau)].conj()
        rval[..., lag + tau] = prod.mean(axis=-1)

    # return
    return rval

    def test2_raster_data_at_points(self):
        # Write a file to test
        f = tempfile.NamedTemporaryFile(suffix='.aai',
                                        prefix='test_misc',
                                        delete=False)
        f.write('ncols 3   \r\n')
        f.write('nrows 2 \r\n')
        f.write('xllcorner +0.   \r\n')
        f.write('yllcorner +8. \r\n')
        f.write('cellsize 1   \r\n')
        f.write('NODATA_value -9999 \r\n')
        f.write('1 2 -9999   \r\n')
        f.write('4 5 6')
        f.close()
        # lon 0 - 3
        # lat 8 - 10

        # Just outside the midpoint of all sides
        lon = asarray([-0.0001, 1.5, 3.0001, 1.5])
        lat = asarray([9., 10.00001, 9.0, 7.99999])
        raster = Raster.from_file(f.name)
        data = raster.raster_data_at_points(lon, lat)
        self.assertTrue(numpy.all(numpy.isnan(data)))

        # Inside lower left corner of No data cell

        lon = asarray([2.0001])
        lat = asarray([9.000019])
        raster = Raster.from_file(f.name)
        data = raster.raster_data_at_points(lon, lat)
        self.assertTrue(numpy.all(numpy.isnan(data)))

        os.remove(f.name)

 def __call__(self, Xi, Xj, ni, nj, hyper_deriv=None, symmetric=False):
     """Evaluate the covariance between points `Xi` and `Xj` with derivative order `ni`, `nj`.
     
     Parameters
     ----------
     Xi : :py:class:`Matrix` or other Array-like, (`M`, `N`)
         `M` inputs with dimension `N`.
     Xj : :py:class:`Matrix` or other Array-like, (`M`, `N`)
         `M` inputs with dimension `N`.
     ni : :py:class:`Matrix` or other Array-like, (`M`, `N`)
         `M` derivative orders for set `i`.
     nj : :py:class:`Matrix` or other Array-like, (`M`, `N`)
         `M` derivative orders for set `j`.
     hyper_deriv : Non-negative int or None, optional
         The index of the hyperparameter to compute the first derivative
         with respect to. If None, no derivatives are taken. Hyperparameter
         derivatives are not supported at this point. Default is None.
     symmetric : bool, optional
         Whether or not the input `Xi`, `Xj` are from a symmetric matrix.
         Default is False.
     
     Returns
     -------
     Kij : :py:class:`Array`, (`M`,)
         Covariances for each of the `M` `Xi`, `Xj` pairs.
     
     Raises
     ------
     NotImplementedError
         If the `hyper_deriv` keyword is not None.
     """
     if hyper_deriv is not None:
         raise NotImplementedError("Hyperparameter derivatives have not been implemented!")
     n_cat = scipy.asarray(scipy.concatenate((ni, nj), axis=1), dtype=int)
     X_cat = scipy.asarray(scipy.concatenate((Xi, Xj), axis=1), dtype=float)
     n_cat_unique = unique_rows(n_cat)
     k = scipy.zeros(Xi.shape[0], dtype=float)
     # Loop over unique derivative patterns:
     if self.num_proc > 1:
         pool = multiprocessing.Pool(processes=self.num_proc)
     for n_cat_state in n_cat_unique:
         idxs = scipy.where(scipy.asarray((n_cat == n_cat_state).all(axis=1)).squeeze())[0]
         if (n_cat_state == 0).all():
             k[idxs] = self.cov_func(Xi[idxs, :], Xj[idxs, :], *self.params)
         else:
             if self.num_proc > 1 and len(idxs) > 1:
                 k[idxs] = scipy.asarray(
                     pool.map(_ArbitraryKernelEval(self, n_cat_state), X_cat[idxs, :]),
                     dtype=float
                 )
             else:
                 for idx in idxs:
                     k[idx] = mpmath.chop(mpmath.diff(self._mask_cov_func,
                                                      X_cat[idx, :],
                                                      n=n_cat_state,
                                                      singular=True))
     
     if self.num_proc > 0:
         pool.close()
     return k