Checks for NDData and data/subshape dimension length agreement in subtract_psf

[Onoddil]

This PR fixes a few minor check/edge cases in subtract_psf.

1. The docstrings say that data can be an astropy.nddata.NDData object, but the first line previously checked for data.ndim which would throw an AttributeError for an NDData object.

2. Currently the only check on whether the PSF can be subtracted directly from the image or if a cutout is required is if subshape is None, but equally true would be the case where data and subshape are equal in both dimensions, which bypasses unnecessary code in extract_array and add_array.

Closes #818 (at least the case where both dimension lengths are equal; for the final edge case of one equal one non-equal dimension length add_array would need amending)