code: def get_grid_data(xall, yall, zall, nbins=100, method='nearest'):
      
      from scipy.interpolate import griddata
      x, y = _np.meshgrid(
          _np.linspace(xall.min(), xall.max(), nbins),
          _np.linspace(yall.min(), yall.max(), nbins),
          indexing='ij')
      z = griddata(
          _np.hstack([xall[:,None], yall[:,None]]),
          zall, (x, y), method=method)
      return x, y, z
docs: Interpolate unstructured two-dimensional data.
  
      Parameters
      ----------
      xall : ndarray(T)
          Sample x-coordinates.
      yall : ndarray(T)
          Sample y-coordinates.
      zall : ndarray(T)
          Sample z-coordinates.
      nbins : int, optional, default=100
          Number of histogram bins used in x/y-dimensions.
      method : str, optional, default='nearest'
          Assignment method; scipy.interpolate.griddata supports the
          methods 'nearest', 'linear', and 'cubic'.
  
      Returns
      -------
      x : ndarray(nbins, nbins)
          The bins' x-coordinates in meshgrid format.
      y : ndarray(nbins, nbins)
          The bins' y-coordinates in meshgrid format.
      z : ndarray(nbins, nbins)
          Interpolated z-data in meshgrid format.