code: def to_phalf_from_pfull(arr, val_toa=0, val_sfc=0):
      
      phalf = np.zeros((arr.shape[0] + 1, arr.shape[1], arr.shape[2]))
      phalf[0] = val_toa
      phalf[-1] = val_sfc
      phalf[1:-1] = 0.5*(arr[:-1] + arr[1:])
      return phalf
docs: Compute data at half pressure levels from values at full levels.
  
      Could be the pressure array itself, but it could also be any other data
      defined at pressure levels.  Requires specification of values at surface
      and top of atmosphere.