def METHOD_NAME(cellWidth, x, y, mesh_filename):
    """
    Add a ``meshDensity`` field into a planar MPAS mesh.  The mesh density is
    defined as:
      meshDensity = (minCellWidth / cellWidth)**4
    Parameters
    ----------
    cellWidth : ndarray
        m x n array of cell width in km
    x, y : ndarray
        Planar coordinates in meters
    mesh_filename : str
        The mesh file to add ``meshDensity`` to
    """
    minCellWidth = cellWidth.min()
    meshDensityVsXY = (minCellWidth / cellWidth)**4
    print('  minimum cell width in grid definition: {0:.0f} km'.format(minCellWidth))
    print('  maximum cell width in grid definition: {0:.0f} km'.format(cellWidth.max()))
    print('Open unstructured MPAS mesh file...')
    ds = nc4.Dataset(mesh_filename, 'r+')
    meshDensity = ds.variables['meshDensity']
    xCell = ds.variables['xCell'][:]
    yCell = ds.variables['xCell'][:]
    print('Interpolating and writing meshDensity...')
    mpasMeshDensity = interp_bilin(x, y, meshDensityVsXY, xCell, yCell)
    meshDensity[:] = mpasMeshDensity
    ds.close()