kammerje / spaceKLIP

Pipeline for reducing JWST high-contrast imaging data. Published in Kammerer et al. 2022 and Carter et al. 2022.
https://ui.adsabs.harvard.edu/abs/2022SPIE12180E..3NK/abstract
MIT License
16 stars 9 forks source link

Redundancy in recenter_frames, update_nircam_centers, and align_frames()? #153

Closed AarynnCarter closed 1 month ago

AarynnCarter commented 3 months ago

Not clear to me what the difference is between these functions. Typically I have only run align_frames() and that has been all I needed.

asahooexo commented 3 months ago

I get the folllowing error when using recenter_frames() function

[spaceKLIP.imagetools:INFO] --> Concatenation JWST_NIRCAM_NRCALONG_F335M_MASKRND_MASK335R_SUB320A335R
[spaceKLIP.imagetools:INFO]   --> Recenter frames: jw04454003001_03106_00001_nrcalong_calints.fits
[py.warnings:WARNING] /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/pysynphot/spectrum.py:1251: DeprecationWarning: `alltrue` is deprecated as of NumPy 1.25.0, and will be removed in NumPy 2.0. Please use `all` instead.
  self.validate_wavetable()

[spaceKLIP.imagetools:INFO]   --> Recenter frames: generating WebbPSF image for absolute centering (this might take a while)
[py.warnings:WARNING] /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/pysynphot/spectrum.py:2910: DeprecationWarning: `alltrue` is deprecated as of NumPy 1.25.0, and will be removed in NumPy 2.0. Please use `all` instead.
  self.validate_wavetable()

[py.warnings:WARNING] /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/astropy/io/fits/file.py:446: AstropyUserWarning: File may have been truncated: actual file length (131072) is smaller than the expected size (285120)
  warnings.warn(

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
Cell In[14], line 1
----> 1 ImageTools.recenter_frames()

File ~/repos/spaceKLIP/spaceKLIP/imagetools.py:1690, in ImageTools.recenter_frames(self, method, subpix_first_sci_only, spectral_type, kwargs, subdir)
   1684 for k in range(data.shape[0]):
   1685     
   1686     # For the first SCI frame, get the star position
   1687     # and the shift between the star and coronagraphic
   1688     # mask position.
   1689     if j == ww_sci[0] and k == 0:
-> 1690         xc, yc, xshift, yshift = self.find_nircam_centers(data0=data[k].copy(),
   1691                                                           key=key,
   1692                                                           j=j,
   1693                                                           spectral_type=spectral_type,
   1694                                                           date=head_pri['DATE-BEG'],
   1695                                                           output_dir=output_dir)
   1697     # Apply the same shift to all SCI and REF frames.
   1698     shifts += [np.array([-(xc - data.shape[-1]//2), -(yc - data.shape[-2]//2)])]

File /opt/miniconda3/envs/spaceklip2/lib/python3.10/contextlib.py:79, in ContextDecorator.__call__.<locals>.inner(*args, **kwds)
     76 @wraps(func)
     77 def inner(*args, **kwds):
     78     with self._recreate_cm():
---> 79         return func(*args, **kwds)

File ~/repos/spaceKLIP/spaceKLIP/imagetools.py:1887, in ImageTools.find_nircam_centers(self, data0, key, j, spectral_type, date, output_dir, fov_pix, oversample, use_coeff)
   1879 APERNAME = self.database.obs[key]['APERNAME'][j]
   1880 kwargs = {
   1881     'fov_pix': fov_pix,
   1882     'oversample': oversample,
   (...)
   1885     'sp': spectrum
   1886 }
-> 1887 psf = JWST_PSF(APERNAME, FILTER, **kwargs)
   1889 # Get SIAF reference pixel position.
   1890 apsiaf = psf.inst_on.siaf_ap

File ~/repos/spaceKLIP/spaceKLIP/psf.py:198, in JWST_PSF.__init__(self, apername, filt, date, fov_pix, oversample, sp, use_coeff, **kwargs)
    196 # Load date-specific OPD files?
    197 if date is not None:
--> 198     inst_on.load_wss_opd_by_date(date=date, choice='closest', verbose=False, plot=False)
    199     inst_off.load_wss_opd_by_date(date=date, choice='closest', verbose=False, plot=False)
    201 # Renormalize spectrum to have 1 e-/sec within bandpass to obtain normalized PSFs

File /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/webbpsf/webbpsf_core.py:1613, in JWInstrument.load_wss_opd_by_date(self, date, choice, verbose, plot, **kwargs)
   1611     date = astropy.time.Time.now().isot
   1612 opd_fn = webbpsf.mast_wss.get_opd_at_time(date, verbose=verbose, choice=choice, **kwargs)
-> 1613 self.load_wss_opd(opd_fn, verbose=verbose, plot=plot, **kwargs)

File /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/webbpsf/webbpsf_core.py:1460, in JWInstrument.load_wss_opd(self, filename, output_path, backout_si_wfe, verbose, plot, save_ote_wfe)
   1458 if verbose:
   1459     print(f"Importing and format-converting OPD from {filename}")
-> 1460 opdhdu = webbpsf.mast_wss.import_wss_opd(filename)
   1462 # Mask out any pixels in the OPD array which are outside the OTE pupil.
   1463 # This is mostly cosmetic, and helps mask out some edge effects from the extrapolation + interpolation in
   1464 # resizing the OPDs
   1465 ote_pupil_mask = utils.get_pupil_mask() != 0

File /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/webbpsf/mast_wss.py:253, in import_wss_opd(filename, npix_out, verbose)
    226 """Import an OPD produced by the JWST WSS, and convert to the right format for use with WebbPSF.
    227 
    228 This includes:
   (...)
    248 
    249 """
    251 wasopd = fits.open(filename)
--> 253 inputOPD = wasopd["RESULT_PHASE"].data
    254 npix_in = inputOPD.shape[0]
    256 wasopd[0].header.add_history("OPD file retrieved from MAST for use by WebbPSF.")

File /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/astropy/utils/decorators.py:837, in lazyproperty.__get__(self, obj, owner)
    835         val = obj_dict.get(self._key, _NotFound)
    836         if val is _NotFound:
--> 837             val = self.fget(obj)
    838             obj_dict[self._key] = val
    839 return val

File /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/astropy/io/fits/hdu/image.py:250, in _ImageBaseHDU.data(self)
    247 if len(self._axes) < 1:
    248     return
--> 250 data = self._get_scaled_image_data(self._data_offset, self.shape)
    251 self._update_header_scale_info(data.dtype)
    253 return data

File /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/astropy/io/fits/hdu/image.py:809, in _ImageBaseHDU._get_scaled_image_data(self, offset, shape)
    802 """
    803 Internal function for reading image data from a file and apply scale
    804 factors to it.  Normally this is used for the entire image, but it
    805 supports alternate offset/shape for Section support.
    806 """
    807 code = BITPIX2DTYPE[self._orig_bitpix]
--> 809 raw_data = self._get_raw_data(shape, code, offset)
    810 raw_data.dtype = raw_data.dtype.newbyteorder(">")
    812 if self._do_not_scale_image_data or (
    813     self._orig_bzero == 0 and self._orig_bscale == 1 and self._blank is None
    814 ):
    815     # No further conversion of the data is necessary

File /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/astropy/io/fits/hdu/base.py:552, in _BaseHDU._get_raw_data(self, shape, code, offset)
    550     return np.ndarray(shape, dtype=code, buffer=self._buffer, offset=offset)
    551 elif self._file:
--> 552     return self._file.readarray(offset=offset, dtype=code, shape=shape)
    553 else:
    554     return None

File /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/astropy/io/fits/file.py:393, in _File.readarray(self, size, offset, dtype, shape)
    390             else:
    391                 raise
--> 393     return np.ndarray(
    394         shape=shape, dtype=dtype, offset=offset, buffer=self._mmap
    395     )
    396 else:
    397     count = reduce(operator.mul, shape)

TypeError: buffer is too small for requested array
kammerje commented 2 months ago

@AarynnCarter The alignment function should be enough to align all subsequent frames in a concatenation to the first frame in that concatenation. KLIP should work after that. However, you will be assuming that the star position (KLIP center) is the CRPIX position in the header which is what we have done for MIRI, but not for NIRCam. For NIRCam, the recentering function finds the position of the star using WebbPSF_ext models, and then centers the first frame of a concatenation on the position of the star, before running the alignment, so that the KLIP center is more accurate. The update centers function was only required while the CRPIX values in CRDS were still wrong in the beginning of the mission. Now it should be redundant, but let's not throw this function out because we might want to reprocess old data. All three functions have their right to exist.

kammerje commented 2 months ago

@JarronL Ananya's issue looks like a WebbPSF/WebbPSF_ext problem. Can you say anything about it?

kammerje commented 1 month ago

@AarynnCarter Can we close this after my explanations?

AarynnCarter commented 1 month ago

Okay gotcha! Yes, I'll close now.

JarronL commented 1 month ago

@JarronL Ananya's issue looks like a WebbPSF/WebbPSF_ext problem. Can you say anything about it?

@kammerje, This appears to be an issue when loading the date-specific OPD map through webbpsf. Following the error trace, I'm not sure why it would be failing. Maybe @mperrin has an idea?

File ~/repos/spaceKLIP/spaceKLIP/psf.py:198, in JWST_PSF.__init__(self, apername, filt, date, fov_pix, oversample, sp, use_coeff, **kwargs)
    196 # Load date-specific OPD files?
    197 if date is not None:
--> 198     inst_on.load_wss_opd_by_date(date=date, choice='closest', verbose=False, plot=False)
    199     inst_off.load_wss_opd_by_date(date=date, choice='closest', verbose=False, plot=False)
    201 # Renormalize spectrum to have 1 e-/sec within bandpass to obtain normalized PSFs
mperrin commented 1 month ago

In terms of the error message, right at the top of the error trace this seems pretty noteworthy to me:

[py.warnings:WARNING] /opt/miniconda3/envs/spaceklip2/lib/python3.10/site-packages/astropy/io/fits/file.py:446: AstropyUserWarning: File may have been truncated: actual file length (131072) is smaller than the expected size (285120)

It's saying the file was truncated, which presumably is likely due to some incomplete download or transient network issue. Likely not a broader problem in general. @asahooexo were you able to resolve your problem back in March? Perhaps by redownloading some file?

mperrin commented 1 month ago

Going back for a moment to the question of the redundancy or not-redundancy of these functions. it would be nice if @kammerje's nice explanation for these was added somewhere into the docs, not just here in this GitHub issue. I'll file a separate issue for this.

asahooexo commented 1 month ago

@mperrin I no longer see this error with the latest version of the code on develop.