code: def fl2norm2(xf, axis=(0, 1)):
      r
  
      xfs = xf.shape
      return (np.linalg.norm(xf)**2) / np.prod(np.array([xfs[k] for k in axis]))
docs: r"""
      Compute the squared :math:`\ell_2` norm in the DFT domain, taking
      into account the unnormalised DFT scaling, i.e. given the DFT of a
      multi-dimensional array computed via :func:`fftn`, return the
      squared :math:`\ell_2` norm of the original array.
  
      Parameters
      ----------
      xf : array_like
        Input array
      axis : sequence of ints, optional (default (0,1))
        Axes on which the input is in the frequency domain
  
      Returns
      -------
      x : float
        :math:`\|\mathbf{x}\|_2^2` where the input array is the result of
        applying :func:`fftn` to the specified axes of multi-dimensional
        array :math:`\mathbf{x}`