code: def interpolate(self, lons, lats, zdata, order=1):
          
  
  
          shape = lons.shape
  
          lons, lats = self._check_integrity(lons, lats)
  
          if order not in [0,1,3]:
              raise ValueError("order must be 0, 1, or 3")
          if zdata.size != self.npoints:
              raise ValueError("data must be of size {}".format(self.npoints))
  
          zdata = self._shuffle_field(zdata)
  
          zi, zierr, ierr = _stripack.interp_n(order, lats, lons,\
                                        self._x, self._y, self._z, zdata,\
                                        self.lst, self.lptr, self.lend)
  
          if ierr != 0:
              print('Warning some points may have errors - check error array\n'.format(ierr))
              zi[zierr < 0] = np.nan
  
          return zi.reshape(shape), zierr.reshape(shape)
docs: Base class to handle nearest neighbour, linear, and cubic interpolation.
          Given a triangulation of a set of nodes on the unit sphere, along with data
          values at the nodes, this method interpolates (or extrapolates) the value
          at a given longitude and latitude.
  
          Parameters
          ----------
           lons : float / array of floats, shape (l,)
                  longitudinal coordinate(s) on the sphere
           lats : float / array of floats, shape (l,)
                  latitudinal coordinate(s) on the sphere
           zdata : array of floats, shape (n,)
                  value at each point in the triangulation
                  must be the same size of the mesh
           order : int (default=1)
                  order of the interpolatory function used
                   0 = nearest-neighbour
                   1 = linear
                   3 = cubic
  
          Returns
          -------
           zi : float / array of floats, shape (l,)
              interpolated value(s) at (lons, lats)
           err : int / array of ints, shape (l,)
              whether interpolation (0), extrapolation (1) or error (other)