code: def diagonal(a, offset=0, axis1=0, axis2=1):
      name = "diagonal-" + tokenize(a, offset, axis1, axis2)
  
      if a.ndim < 2:
          # NumPy uses `diag` as we do here.
          raise ValueError(
              "diag requires an array of at least two dimensions"
          )
  
      def _axis_fmt(axis, name, ndim):
          if axis < 0:
              t = ndim + axis
              if t < 0:
                  msg = "{}: axis {} is out of bounds for array of dimension {}"
                  raise AxisError(msg.format(name, axis, ndim))
              axis = t
          return axis
  
      def pop_axes(chunks, axis1, axis2):
          chunks = list(chunks)
          chunks.pop(axis2)
          chunks.pop(axis1)
          return tuple(chunks)
  
      axis1 = _axis_fmt(axis1, "axis1", a.ndim)
      axis2 = _axis_fmt(axis2, "axis2", a.ndim)
  
      if axis1 == axis2:
          raise ValueError("axis1 and axis2 cannot be the same")
  
      a = asarray(a)
      k = offset
      if axis1 > axis2:
          axis1, axis2 = axis2, axis1
          k = -offset
  
      free_axes = set(range(a.ndim)) - {axis1, axis2}
      free_indices = list(
          product(*(range(a.numblocks[i]) for i in free_axes))
      )
      ndims_free = len(free_axes)
  
      # equation of diagonal: i = j - k
      kdiag_row_start = max(0, -k)
      kdiag_col_start = max(0, k)
      kdiag_row_stop = min(a.shape[axis1], a.shape[axis2] - k)
      len_kdiag = kdiag_row_stop - kdiag_row_start
  
      if len_kdiag <= 0:
          xp = np
  
          if is_cupy_type(a._meta):
              import cupy
  
              xp = cupy
  
          out_chunks = pop_axes(a.chunks, axis1, axis2) + ((0,),)
          dsk = dict()
          for free_idx in free_indices:
              shape = tuple(
                  out_chunks[axis][free_idx[axis]]
                  for axis in range(ndims_free)
              )
              dsk[(name,) + free_idx + (0,)] = (
                  partial(xp.empty, dtype=a.dtype),
                  shape + (0,),
              )
  
          meta = meta_from_array(a, ndims_free + 1)
          return Array(dsk, name, out_chunks, meta=meta)
  
      # compute row index ranges for chunks along axis1:
      row_stops_ = np.cumsum(a.chunks[axis1])
      row_starts = np.roll(row_stops_, 1)
      row_starts[0] = 0
  
      # compute column index ranges for chunks along axis2:
      col_stops_ = np.cumsum(a.chunks[axis2])
      col_starts = np.roll(col_stops_, 1)
      col_starts[0] = 0
  
      # locate first chunk containing diagonal:
      row_blockid = np.arange(a.numblocks[axis1])
      col_blockid = np.arange(a.numblocks[axis2])
  
      row_filter = (row_starts <= kdiag_row_start) & (
          kdiag_row_start < row_stops_
      )
      col_filter = (col_starts <= kdiag_col_start) & (
          kdiag_col_start < col_stops_
      )
      (I,) = row_blockid[row_filter]
      (J,) = col_blockid[col_filter]
  
      # follow k-diagonal through chunks while constructing dask graph:
      dsk = dict()
      i = 0
      kdiag_chunks = ()
      while (
          kdiag_row_start < a.shape[axis1]
          and kdiag_col_start < a.shape[axis2]
      ):
          # localize block info:
          nrows, ncols = a.chunks[axis1][I], a.chunks[axis2][J]
          kdiag_row_start -= row_starts[I]
          kdiag_col_start -= col_starts[J]
          k = (
              -kdiag_row_start
              if kdiag_row_start > 0
              else kdiag_col_start
          )
          kdiag_row_end = min(nrows, ncols - k)
          kdiag_len = kdiag_row_end - kdiag_row_start
  
          # increment dask graph:
          for free_idx in free_indices:
              input_idx = (
                  free_idx[:axis1]
                  + (I,)
                  + free_idx[axis1 : axis2 - 1]
                  + (J,)
                  + free_idx[axis2 - 1 :]
              )
              output_idx = free_idx + (i,)
              dsk[(name,) + output_idx] = (
                  np.diagonal,
                  (a.name,) + input_idx,
                  k,
                  axis1,
                  axis2,
              )
  
          kdiag_chunks += (kdiag_len,)
          # prepare for next iteration:
          i += 1
          kdiag_row_start = kdiag_row_end + row_starts[I]
          kdiag_col_start = min(ncols, nrows + k) + col_starts[J]
          I = I + 1 if kdiag_row_start == row_stops_[I] else I
          J = J + 1 if kdiag_col_start == col_stops_[J] else J
  
      out_chunks = pop_axes(a.chunks, axis1, axis2) + (kdiag_chunks,)
      graph = HighLevelGraph.from_collections(
          name, dsk, dependencies=[a]
      )
      meta = meta_from_array(a, ndims_free + 1)
      return Array(graph, name, out_chunks, meta=meta)
  
imports: from __future__ import annotations
  import itertools
  from collections.abc import Sequence
  from functools import partial
  from itertools import product
  from numbers import Integral, Number
  import numpy
  from tlz import sliding_window
  from dask.array import chunk
  from dask.array.backends import array_creation_dispatch
  from dask.array.core import Array, asarray, block, blockwise, broadcast_arrays, broadcast_to, concatenate, normalize_chunks, stack
  from dask.array.numpy_compat import AxisError
  from dask.array.ufunc import greater_equal, rint
  from dask.array.utils import meta_from_array
  from dask.array.wrap import empty, full, ones, zeros
  from dask.base import tokenize
  from dask.highlevelgraph import HighLevelGraph
  from dask.utils import cached_cumsum, derived_from, is_cupy_type
label_desc: Use image node_id 13 for calling a global function with example usage: diagonal(a, offset, axis1, axis2) and returns: Array, axis, tuple, Array