# Last edited on 2021-08-01 19:17:28 by jstolfi

def basis_2d_grid_choose(xmin,xmax,nx, ymin,ymax,ny, span):
  # Rturns a list {c} of /centroids/ and a list {r} of nominal /radii/ half-widths 
  # that define a set of {nx} by {ny} hump-like function basis elements on a square grid.
  # Both lists will have {nx*ny} elements, which will be 2-tuples of floats.
  # 
  # The centroids will comprise a regular grid of points (2-tuples of
  # floats) in the rectangle {[xmin _ xmax] × [ymin _ ymax]}. 
  # All elements will have the same {x}-radius and same {y}-radius.
  # 
  # If {span} is true, the {x}-range of the centroids will be that
  # interval. Otherwise, it will be displaced inwards by the nominal
  # {x}-radius at either end. In either case, the radius along each axis
  # will be equal to the centroid spacing along that axis.
  #
  # As a special case, if {nx} is 1 the centroid's {x} coordinate will
  # always be the midpoint of {[xmin _ xmax]}. If {span} s true, the nominal
  # {x}-radius will be infinite, meaning that the element will be
  # constant along the {x}-direction. If span is false, the nominal
  # {x}-radius will be half the width of that interval.
  # 
  # The analogous logic applies on the {y} axis.
 
  assert nx >= 1, "invalid {nx}"
  assert ny >= 1, "invalid {ny}"
 
  nb = nx*ny
  
  # Number of full steps in each coord interval:
  ndx = nx - 1 if span else nx + 1;
  ndy = ny - 1 if span else ny + 1;

  # Nominal radii along each axis:
  rbx = +inf if ndx == 0 else (xmax - xmin)/ndx
  rby = +inf if ndy == 0 else (ymax - ymin)/ndy

  c = []   # Centroids.
  r = [];  # Nominal radii.
  
  for kx in range(nx):
    tx = 0.5 if ndx == 0 else (kx if span else kx + 1)/ndx
    xck = (1-tx)*xmin + tx*xmax
    for ky in range(ny):
      ty = 0.5 if ndy == 0 else (ky if span else ky + 1)/ndy
      yck = (1-ty)*ymin + ty*ymax
      c.append((xck,yck))
      r.append((rbx,rby))
  return c, r
  # ----------------------------------------------------------------------