Python fitpack.bisplev函数代码示例

    def correct(self, pos):
        delta1 = fitpack.bisplev(pos[1], pos[0], [self.xSplineKnotsX,
                                                  self.xSplineKnotsY,
                                                  self.xSplineCoeff,
                                                  self.splineOrder,
                                                  self.splineOrder],
                                 dx=0, dy=0)

        delta0 = fitpack.bisplev(pos[1], pos[0], [self.ySplineKnotsX,
                                                  self.ySplineKnotsY,
                                                  self.ySplineCoeff,
                                                  self.splineOrder,
                                                  self.splineOrder],
                                 dx=0, dy=0)
        return delta0 + pos[0], delta1 + pos[1]

    def splineFuncY(self, x, y):
        """
        calculates the displacement matrix using fitpack for the Y
        direction

        @param x: points in the x direction
        @type x: ndarray
        @param y: points in the y direction
        @type y: ndarray

        @return: displacement matrix for the Y direction
        @rtype: ndarray
        """
        if x.ndim == 2:
            if abs(x[1:, :] - x[:-1, :] - numpy.zeros((x.shape[0] - 1, x.shape[1]))).max() < 1e-6:
                x = x[0]
                y = y[:, 0]
            elif abs(x[:, 1:] - x[:, :-1] - numpy.zeros((x.shape[0], x.shape[1] - 1))).max() < 1e-6:
                x = x[:, 0]
                y = y[0]

        yDispArray = fitpack.bisplev(
            x, y, [self.ySplineKnotsX,
                   self.ySplineKnotsY,
                   self.ySplineCoeff,
                   self.splineOrder,
                   self.splineOrder ],
            dx=0, dy=0).transpose()

        return yDispArray

    def get_derivative(self, x, dx=None):
        """Evaluate derivative at x"""
        if dx is None:
            dx = [0, 0]

        knots = self.rbs.get_knots()
        coeffs = self.rbs.get_coeffs()
        tck = [knots[0], knots[1], coeffs, self.kx0, self.kx1]
        return fitpack.bisplev(x[0], x[1], tck, dx[0], dx[1])

    def spline2array(self, timing=False):
        """
        Calculates the displacement matrix using fitpack
        bisplev(x, y, tck, dx = 0, dy = 0)

        :param timing: profile the calculation or not
        :type timing: bool

        :return: Nothing !
        :rtype: float or ndarray

        Evaluate a bivariate B-spline and its derivatives. Return a
        rank-2 array of spline function values (or spline derivative
        values) at points given by the cross-product of the rank-1
        arrays x and y. In special cases, return an array or just a
        float if either x or y or both are floats.
        """
        if self.xDispArray is None:
            x_1d_array = numpy.arange(self.xmin, self.xmax + 1)
            y_1d_array = numpy.arange(self.ymin, self.ymax + 1)
            startTime = time.time()
            self.xDispArray = fitpack.bisplev(
                x_1d_array,
                y_1d_array,
                [self.xSplineKnotsX, self.xSplineKnotsY, self.xSplineCoeff, self.splineOrder, self.splineOrder],
                dx=0,
                dy=0,
            ).transpose()
            intermediateTime = time.time()
            self.yDispArray = fitpack.bisplev(
                x_1d_array,
                y_1d_array,
                [self.ySplineKnotsX, self.ySplineKnotsY, self.ySplineCoeff, self.splineOrder, self.splineOrder],
                dx=0,
                dy=0,
            ).transpose()
            if timing:
                logger.info(
                    "Timing for: X-Displacement spline evaluation: %.3f sec,"
                    " Y-Displacement Spline evaluation:  %.3f sec."
                    % ((intermediateTime - startTime), (time.time() - intermediateTime))
                )

 def check_5(self,f=f2,kx=3,ky=3,xb=0,xe=2*pi,yb=0,ye=2*pi,Nx=20,Ny=20,s=0):
     x=xb+(xe-xb)*arange(Nx+1,dtype=float)/float(Nx)
     y=yb+(ye-yb)*arange(Ny+1,dtype=float)/float(Ny)
     xy=makepairs(x,y)
     tck=bisplrep(xy[0],xy[1],f(xy[0],xy[1]),s=s,kx=kx,ky=ky)
     tt=[tck[0][kx:-kx],tck[1][ky:-ky]]
     t2=makepairs(tt[0],tt[1])
     v1=bisplev(tt[0],tt[1],tck)
     v2=f2(t2[0],t2[1])
     v2.shape=len(tt[0]),len(tt[1])
     err = norm2(ravel(v1-v2))
     assert_(err < 1e-2, err)
     put(err)

 def test_rerun_lwrk2_too_small(self):
     # in this setting, lwrk2 is too small in the default run. Here we
     # check for equality with the bisplrep/bisplev output because there,
     # an automatic re-run of the spline representation is done if ier>10.
     x = np.linspace(-2, 2, 80)
     y = np.linspace(-2, 2, 80)
     z = x + y
     xi = np.linspace(-1, 1, 100)
     yi = np.linspace(-2, 2, 100)
     tck = bisplrep(x, y, z)
     res1 = bisplev(xi, yi, tck)
     interp_ = SmoothBivariateSpline(x, y, z)
     res2 = interp_(xi, yi)
     assert_almost_equal(res1, res2)

 def test_bispev(self):
     x_1d_array = numpy.arange(self.spline.xmin, self.spline.xmax + 1)
     y_1d_array = numpy.arange(self.spline.ymin, self.spline.ymax + 1)
     t0 = time.time()
     dx_ref = fitpack.bisplev(
         x_1d_array,
         y_1d_array,
         [
             self.spline.xSplineKnotsX,
             self.spline.xSplineKnotsY,
             self.spline.xSplineCoeff,
             self.spline.splineOrder,
             self.spline.splineOrder,
         ],
         dx=0,
         dy=0,
     )
     t1 = time.time()
     logger.debug(self.spline.xSplineKnotsX.dtype)
     logger.debug(self.spline.xSplineKnotsY.dtype)
     logger.debug(self.spline.xSplineCoeff.dtype)
     dx_loc = _bispev.bisplev(
         x_1d_array,
         y_1d_array,
         [
             self.spline.xSplineKnotsX,
             self.spline.xSplineKnotsY,
             self.spline.xSplineCoeff,
             self.spline.splineOrder,
             self.spline.splineOrder,
         ],
     )
     t2 = time.time()
     logger.debug("Scipy timings: %.3fs\t cython timings: %.3fs" % (t1 - t0, t2 - t1))
     logger.debug("%s, %s" % (dx_ref.shape, dx_loc.shape))
     logger.debug(dx_ref)
     logger.debug(dx_loc)
     logger.debug("delta = %s" % abs(dx_loc - dx_ref).max())
     if logger.getEffectiveLevel() == logging.DEBUG:
         fig = pylab.figure()
         ax1 = fig.add_subplot(121)
         ax2 = fig.add_subplot(122)
         ax1.imshow(dx_ref)
         ax2.imshow(dx_loc)
         fig.show()
         six.moves.input()
     self.assert_(abs(dx_loc - dx_ref).max() < 2e-5, "Result are similar")

    def splineFuncY(self, x, y, list_of_points=False):
        """
        calculates the displacement matrix using fitpack for the Y
        direction

        @param x: points in the x direction
        @type x: ndarray
        @param y: points in the y direction
        @type y: ndarray
        @param list_of_points: if true, consider the zip(x,y) instead of the of the square array
        @return: displacement matrix for the Y direction
        @rtype: ndarray
        """
        if x.ndim == 2:
            if abs(x[1:, :] - x[:-1, :] - numpy.zeros((x.shape[0] - 1, x.shape[1]))).max() < 1e-6:
                x = x[0]
                y = y[:, 0]
            elif abs(x[:, 1:] - x[:, :-1] - numpy.zeros((x.shape[0], x.shape[1] - 1))).max() < 1e-6:
                x = x[:, 0]
                y = y[0]

        if list_of_points and x.ndim == 1 and len(x) == len(y):
            lx = ly = len(x)
            x_order = x.argsort()
            y_order = y.argsort()
            x = x[x_order]
            y = y[y_order]
            x_unordered = numpy.zeros(lx, dtype=int)
            y_unordered = numpy.zeros(ly, dtype=int)
            x_unordered[x_order] = numpy.arange(lx)
            y_unordered[y_order] = numpy.arange(ly)

        yDispArray = fitpack.bisplev(
            x, y, [self.ySplineKnotsX,
                   self.ySplineKnotsY,
                   self.ySplineCoeff,
                   self.splineOrder,
                   self.splineOrder ],
            dx=0, dy=0)
        if list_of_points and x.ndim == 1:
            return yDispArray[x_unordered, y_unordered]
        else:
            return yDispArray.T

    def make_diff_func(self) :
        ''' everytime parameters are changed, the diffusion 
            function must be REMADE or it's data adapted !! 
        '''
        if  (self.difftype == 'MC_2D4pint_ext1D' 
                    or self.difftype == 'MC_2D4pint_func'):
                from GridUtils import GridFunc2D
                self.diff11   = GridFunc2D([self.param[0],self.param[1]],N.asarray(self.param[2]),False)
                self.diff12   = GridFunc2D([self.param[0],self.param[1]],N.asarray(self.param[3]),False)
                self.diff21   = GridFunc2D([self.param[0],self.param[1]],N.asarray(self.param[4]),False)
                self.diff22   = GridFunc2D([self.param[0],self.param[1]],N.asarray(self.param[5]),False)
                self.diff11_u = None;self.diff12_u = None;self.diff21_u = None;self.diff22_u = None;
                self.diff11_v = None;self.diff12_v = None;self.diff21_v = None;self.diff22_v = None;
                #grideval means a ugrid and vgrid is given, the grid made by this is evaluated
                #this is used in the plotting function
                self.diff11_grideval   = GridFunc2D([self.param[0],self.param[1]],N.asarray(self.param[2]),True)
                self.diff12_grideval   = GridFunc2D([self.param[0],self.param[1]],N.asarray(self.param[3]),True)
                self.diff21_grideval   = GridFunc2D([self.param[0],self.param[1]],N.asarray(self.param[4]),True)
                self.diff22_grideval   = GridFunc2D([self.param[0],self.param[1]],N.asarray(self.param[5]),True)


        elif self.difftype == 'MC_bspline_ext1D' or \
                  self.difftype == 'MC_bspline1storder_ext1D' : 
                self.arrayaware = False           
                weight = ones(self.nrparam)
                # standard dev on diffusion = 0.2, inverse is 5
                weight[:] =5 * weight[:] 
                #create an interp object
                from scipy.interpolate.fitpack import bisplrep,bisplev
                #param[0] should be list of u values
                #param[1] should be list of v values
                #param[2] should be list of D(u,v) values
                #create the data of the bicubic bspline:
                # no smoother so s = 0
                if self.difftype == 'MC_bspline1storder_ext1D' :
                    self.splinedataD11 = bisplrep(self.param[0], self.param[1],
                        self.param[2], w= weight, s=None,
                        kx = 1, ky=1, full_output = 1, quiet = 0)
                else :
                    #nknots = 11
                    #tx = range(nknots)
                    #tx = [self.umin + x/(1.*(nknots-1)) * (self.umax-self.umin) for x in tx]
                    #ty = range(nknots)
                    #ty = [self.vmin + x/(1.*(nknots-1)) * (self.vmax-self.vmin) for x in ty]
                    #tx = tx[0]*3 + tx + tx[-1]*3
                    #tx[:1] = [tx[0]]*4 ; tx[-1:] = [tx[-1]]*4
                    #ty[:1] = [ty[0]]*4 ; ty[-1:] = [ty[-1]]*4
                    #print 'tx,ty',tx, ty
                    self.splinedataD11 = bisplrep(self.param[0], self.param[1],
                        self.param[2], w= weight, s=2100, 
                        kx = 3, ky=3, full_output = 1, quiet = 0)
                #print 'D11 output :', self.splinedataD11[1:]
                if self.splinedataD11[2] > 0 :
                    print 'ERROR in spline interpolation'
                    sys.exit()
                self.splinedataD11 = self.splinedataD11[0]
                self.diff11   = lambda u,v: bisplev([u], [v], self.splinedataD11, \
                            dx = 0, dy=0)
                self.diff11_u = lambda u,v: bisplev([u], [v] , self.splinedataD11,\
                            dx = 1, dy=0)
                self.diff11_v = lambda u,v: bisplev([u], [v] , self.splinedataD11,\
                            dx = 0, dy=1)
                self.diff11_grideval   = lambda u,v: bisplev(u, v, self.splinedataD11, \
                            dx = 0, dy=0)
                self.diff11_u_grideval = lambda u,v: bisplev(u, v , self.splinedataD11,\
                            dx = 1, dy=0)
                self.diff11_v_grideval = lambda u,v: bisplev(u, v , self.splinedataD11,\
                            dx = 0, dy=1)
                if self.difftype == 'MC_bspline1storder_ext1D' :
                    self.splinedataD12 = bisplrep(self.param[0], self.param[1],
                        self.param[3], w= weight, s=None,
                        kx = 1, ky=1, full_output = 1, quiet = 1)
                else :
                    self.splinedataD12 = bisplrep(self.param[0], self.param[1],
                        self.param[3], w= weight, s=None,
                        kx = 3, ky=3, full_output = 1, quiet = 0)
                #print 'D12 output :', self.splinedataD12[1:]
                if self.splinedataD12[2] > 0 :
                    print 'ERROR in spline interpolation'
                    sys.exit()
                self.splinedataD12 = self.splinedataD12[0]
                self.diff12   = lambda u,v: bisplev([u], [v], self.splinedataD12, \
                            dx = 0, dy=0)
                self.diff12_u = lambda u,v: bisplev([u], [v] , self.splinedataD12,\
                            dx = 1, dy=0)
                self.diff12_v = lambda u,v: bisplev([u], [v] , self.splinedataD12,\
                            dx = 0, dy=1)
                self.diff12_grideval   = lambda u,v: bisplev(u, v, self.splinedataD12, \
                            dx = 0, dy=0)
                self.diff12_u_grideval = lambda u,v: bisplev(u, v , self.splinedataD12,\
                            dx = 1, dy=0)
                self.diff12_v_grideval = lambda u,v: bisplev(u, v , self.splinedataD12,\
                            dx = 0, dy=1)
                if self.difftype == 'MC_bspline1storder_ext1D' :
                    self.splinedataD21 = bisplrep(self.param[0], self.param[1],
                        self.param[4], w= weight, s=None,
                        kx = 1, ky=1, full_output = 1, quiet = 1)
                else :
                    self.splinedataD21 = bisplrep(self.param[0], self.param[1],
                        self.param[4], w= weight, s=None,
#.........这里部分代码省略.........