Image Combine with wcs offset

For simple and primitive image analysis for non-sidereal tracking data, it is quite common to use script like
```
imcombine @fits.list output=combined.fits combine="median" offset="wcs"
```

in IRAF (`fits.list` is the file containing the list of files), as explained [here](http://stsdas.stsci.edu/cgi-bin/gethelp.cgi?imcombine). *I want to request some simple features based on the difference between IRAF and python*.

The only module similar to the above IRAF task is `wcs_project` of `ccdproc`, AFAIK. However, the result from `ccdproc` differs from IRAF, which was unexpected for me.

For test, I combined 4 images using the first image's wcs as the `target_wcs`:

```python
import numpy as np
from ccdproc import CCDData, wcs_project, combine
from astropy.wcs import WCS

prefix = 'reproj_Tutorial/'
filelist = np.loadtxt(prefix+'fits.list', dtype=bytes).astype(str)

reproj = []
wcs_std = WCS((CCDData.read(prefix+filelist[0], unit='adu')).header)

for fname in filelist:
ccd = CCDData.read(prefix+fname, unit='adu')
reproj.append(wcs_project(ccd, wcs_std))

med = combine(reproj, output_file='test.fits', method='average')
```

Since I used the `target_wcs` as the wcs of the first input image, the `NAXIS` of `ccdproc`'s result is restricted to the field of view of the first image. I think most people who want to do the image combine with wcs offset want the whole image, not a "cropped" image at a certain FOV.

The following is the result from IRAF (left) and `ccdproc` (right):


The lower part and right part are "cropped" in `ccdproc`'s result, since they are out of the first wcs's FOV (`NAXIS=1024`). The `NAXIS=(1054, 1064)` for IRAF's result.

-----

So I want

1. Some options for `wcs_project` to reproduce similar result as IRAF (e.g., `restrict_NAXIS` or `fix_NAXIS` that has default value `True` for backward compatibility?)

2. And/or some functionality to use the wcs offset when combining images. I guess it will be very useful if some functions or modules, such as `median_combine`, have an option to take `offset={'wcs', 'world', 'physical', 'grid', }`, and do the reprojection and combination automatically based on the input. For example,

```python
combiner = Combiner([ccd1, ccd2, ...])
avg = combiner.median_combine(offset='wcs')
```

Start with the ccdproc.wcs_project and Combiner entry points described in the issue, and compare their current handling of target_wcs, NAXIS, and image combination. Define and test whether WCS-offset combination should expand beyond the target field of view, and whether reprojection should be performed automatically by the combiner.

Written by the indexing model from the issue text.