diff --git a/.github/workflows/build.yml b/.github/workflows/build.yml
index f5591067eb..cb1591fdeb 100644
--- a/.github/workflows/build.yml
+++ b/.github/workflows/build.yml
@@ -21,7 +21,7 @@ jobs:
build:
needs: [ check ]
if: always() && (needs.check.result == 'success' || needs.check.result == 'skipped')
- uses: OpenAstronomy/github-actions-workflows/.github/workflows/publish.yml@d68193b68216da64eafaa618f53c59f5d271c56e # v1.14.0
+ uses: OpenAstronomy/github-actions-workflows/.github/workflows/publish.yml@9f1f43251dde69da8613ea8e11144f05cdea41d5 # v1.15.0
with:
upload_to_pypi: ${{ (github.event_name == 'release') && (github.event.action == 'released') }}
targets: |
diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml
index 083d21625a..9dceff2eff 100644
--- a/.github/workflows/ci.yml
+++ b/.github/workflows/ci.yml
@@ -20,8 +20,8 @@ on:
type: string
required: false
default: ''
- crds_server:
- description: CRDS server
+ crds_server_url:
+ description: CRDS server URL
type: string
required: false
default: https://jwst-crds.stsci.edu
@@ -40,34 +40,23 @@ jobs:
python-version: '3.12'
- uses: pre-commit/action@2c7b3805fd2a0fd8c1884dcaebf91fc102a13ecd # v3.0.1
check:
- uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@d68193b68216da64eafaa618f53c59f5d271c56e # v1.14.0
+ uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@9f1f43251dde69da8613ea8e11144f05cdea41d5 # v1.15.0
with:
default_python: "3.12"
envs: |
- linux: check-dependencies
- linux: check-types
- latest_crds_contexts:
- uses: ./.github/workflows/contexts.yml
- crds_context:
- needs: [ latest_crds_contexts ]
- runs-on: ubuntu-latest
- steps:
- - id: context
- run: echo context=${{ github.event_name == 'workflow_dispatch' && (inputs.crds_context != '' && inputs.crds_context || needs.latest_crds_contexts.outputs.jwst) || needs.latest_crds_contexts.outputs.jwst }} >> $GITHUB_OUTPUT
- outputs:
- context: ${{ steps.context.outputs.context }}
test:
- uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@d68193b68216da64eafaa618f53c59f5d271c56e # v1.14.0
- needs: [ crds_context ]
+ uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@9f1f43251dde69da8613ea8e11144f05cdea41d5 # v1.15.0
with:
setenv: |
CRDS_PATH: /tmp/data/crds_cache
- CRDS_SERVER_URL: ${{ github.event_name == 'workflow_dispatch' && inputs.crds_server || 'https://jwst-crds.stsci.edu' }}
- CRDS_CONTEXT: ${{ needs.crds_context.outputs.context }}
+ CRDS_SERVER_URL: ${{ inputs.crds_server_url || 'https://jwst-crds.stsci.edu' }}
+ CRDS_CONTEXT: ${{ inputs.crds_context || 'jwst-edit' }}
CRDS_CLIENT_RETRY_COUNT: 3
CRDS_CLIENT_RETRY_DELAY_SECONDS: 20
cache-path: /tmp/data/crds_cache
- cache-key: crds-${{ needs.crds_context.outputs.context }}
+ cache-key: crds-${{ inputs.crds_context || 'jwst-edit' }}
envs: |
- linux: py310-oldestdeps-xdist-cov
pytest-results-summary: true
diff --git a/.github/workflows/ci_cron.yml b/.github/workflows/ci_cron.yml
index 3e7ab20618..1822318bbf 100644
--- a/.github/workflows/ci_cron.yml
+++ b/.github/workflows/ci_cron.yml
@@ -11,8 +11,8 @@ on:
type: string
required: false
default: ''
- crds_server:
- description: CRDS server
+ crds_server_url:
+ description: CRDS server URL
type: string
required: false
default: https://jwst-crds.stsci.edu
@@ -22,25 +22,14 @@ concurrency:
cancel-in-progress: true
jobs:
- latest_crds_contexts:
- uses: ./.github/workflows/contexts.yml
- crds_context:
- needs: [ latest_crds_contexts ]
- runs-on: ubuntu-latest
- steps:
- - id: context
- run: echo context=${{ github.event_name == 'workflow_dispatch' && (inputs.crds_context != '' && inputs.crds_context || needs.latest_crds_contexts.outputs.jwst) || needs.latest_crds_contexts.outputs.jwst }} >> $GITHUB_OUTPUT
- outputs:
- context: ${{ steps.context.outputs.context }}
test:
if: (github.repository == 'spacetelescope/jwst' && (github.event_name == 'schedule' || github.event_name == 'push' || github.event_name == 'workflow_dispatch' || contains(github.event.pull_request.labels.*.name, 'run scheduled tests')))
- uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@d68193b68216da64eafaa618f53c59f5d271c56e # v1.14.0
- needs: [ crds_context ]
+ uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@9f1f43251dde69da8613ea8e11144f05cdea41d5 # v1.15.0
with:
setenv: |
CRDS_PATH: /tmp/crds_cache
- CRDS_SERVER_URL: ${{ github.event_name == 'workflow_dispatch' && inputs.crds_server || 'https://jwst-crds.stsci.edu' }}
- CRDS_CONTEXT: ${{ needs.crds_context.outputs.context }}
+ CRDS_SERVER_URL: ${{ inputs.crds_server_url || 'https://jwst-crds.stsci.edu' }}
+ CRDS_CONTEXT: ${{ inputs.crds_context || 'jwst-edit' }}
CRDS_CLIENT_RETRY_COUNT: 3
CRDS_CLIENT_RETRY_DELAY_SECONDS: 20
cache-path: /tmp/crds_cache
diff --git a/.github/workflows/contexts.yml b/.github/workflows/contexts.yml
deleted file mode 100644
index b3047bdb63..0000000000
--- a/.github/workflows/contexts.yml
+++ /dev/null
@@ -1,26 +0,0 @@
-name: contexts
-
-on:
- workflow_call:
- outputs:
- jwst:
- value: ${{ jobs.contexts.outputs.jwst }}
- workflow_dispatch:
-
-jobs:
- contexts:
- name: retrieve latest CRDS contexts
- runs-on: ubuntu-latest
- outputs:
- jwst: ${{ steps.jwst_crds_context.outputs.pmap }}
- steps:
- - id: jwst_crds_context
- env:
- OBSERVATORY: jwst
- CRDS_SERVER_URL: https://jwst-crds.stsci.edu
- run: >
- echo "pmap=$(
- curl -s -X POST -d '{"jsonrpc": "1.0", "method": "get_default_context", "params": ["${{ env.OBSERVATORY }}", null], "id": 1}' ${{ env.CRDS_SERVER_URL }}/json/ --retry 8 --connect-timeout 10 |
- python -c "import sys, json; print(json.load(sys.stdin)['result'])"
- )" >> $GITHUB_OUTPUT
- - run: if [[ ! -z "${{ steps.jwst_crds_context.outputs.pmap }}" ]]; then echo ${{ steps.jwst_crds_context.outputs.pmap }}; else exit 1; fi
diff --git a/.github/workflows/tests_devdeps.yml b/.github/workflows/tests_devdeps.yml
index 57715057a6..5f79a9982d 100644
--- a/.github/workflows/tests_devdeps.yml
+++ b/.github/workflows/tests_devdeps.yml
@@ -20,8 +20,8 @@ on:
type: string
required: false
default: ''
- crds_server:
- description: CRDS server
+ crds_server_url:
+ description: CRDS server URL
type: string
required: false
default: https://jwst-crds.stsci.edu
@@ -31,25 +31,14 @@ concurrency:
cancel-in-progress: true
jobs:
- latest_crds_contexts:
- uses: ./.github/workflows/contexts.yml
- crds_context:
- needs: [ latest_crds_contexts ]
- runs-on: ubuntu-latest
- steps:
- - id: context
- run: echo context=${{ github.event_name == 'workflow_dispatch' && (inputs.crds_context != '' && inputs.crds_context || needs.latest_crds_contexts.outputs.jwst) || needs.latest_crds_contexts.outputs.jwst }} >> $GITHUB_OUTPUT
- outputs:
- context: ${{ steps.context.outputs.context }}
test:
if: (github.repository == 'spacetelescope/jwst' && (github.event_name == 'schedule' || github.event_name == 'push' || github.event_name == 'workflow_dispatch' || contains(github.event.pull_request.labels.*.name, 'run devdeps tests')))
- uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@d68193b68216da64eafaa618f53c59f5d271c56e # v1.14.0
- needs: [ crds_context ]
+ uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@9f1f43251dde69da8613ea8e11144f05cdea41d5 # v1.15.0
with:
setenv: |
CRDS_PATH: /tmp/data/crds_cache
- CRDS_SERVER_URL: ${{ github.event_name == 'workflow_dispatch' && inputs.crds_server || 'https://jwst-crds.stsci.edu' }}
- CRDS_CONTEXT: ${{ needs.crds_context.outputs.context }}
+ CRDS_SERVER_URL: ${{ inputs.crds_server_url || 'https://jwst-crds.stsci.edu' }}
+ CRDS_CONTEXT: ${{ inputs.crds_context || 'jwst-edit' }}
CRDS_CLIENT_RETRY_COUNT: 3
CRDS_CLIENT_RETRY_DELAY_SECONDS: 20
cache-path: /tmp/data/crds_cache
diff --git a/changes/8846.ami.rst b/changes/8846.ami.rst
new file mode 100644
index 0000000000..a965c823ff
--- /dev/null
+++ b/changes/8846.ami.rst
@@ -0,0 +1 @@
+Change how AMI observables are averaged: average fringe quantities before calculating additional observables. Update their error calculation: use covariance of amplitudes/phases (and derived quantities) and standard error of the mean. Code now expects an ASDF filename string for user-supplied affine2d and bandpass input arguments. Example file creation in documentation.
\ No newline at end of file
diff --git a/changes/8957.general.rst b/changes/8957.general.rst
new file mode 100644
index 0000000000..91a475dc0f
--- /dev/null
+++ b/changes/8957.general.rst
@@ -0,0 +1 @@
+remove ``okify_regtests`` script (move to ``ci_watson``)
diff --git a/changes/8978.assign_wcs.rst b/changes/8978.assign_wcs.rst
new file mode 100644
index 0000000000..2771c1d3c6
--- /dev/null
+++ b/changes/8978.assign_wcs.rst
@@ -0,0 +1 @@
+Fix negative SHUTTRID values under numpy 2
diff --git a/docs/jwst/ami_analyze/description.rst b/docs/jwst/ami_analyze/description.rst
index c0e2d27016..26b2f597ff 100644
--- a/docs/jwst/ami_analyze/description.rst
+++ b/docs/jwst/ami_analyze/description.rst
@@ -36,7 +36,7 @@ other options:
rotation search values. The default setting of '-3 3 1'
results in search values of [-3, -2, -1, 0, 1, 2, 3].
-:--bandpass: Synphot spectrum or suitable array to override filter/source
+:--bandpass: ASDF file containing suitable array to override filter/source
(default=None)
:--usebp: If True, exclude pixels marked DO_NOT_USE from fringe fitting
@@ -47,10 +47,73 @@ other options:
:--chooseholes: If not None, fit only certain fringes e.g. ['B4','B5','B6','C2']
(default=None)
-:--affine2d: User-defined Affine2d object (default=None)
+:--affine2d: ASDF file containing user-defined affine parameters (default=None)
:--run_bpfix: Run Fourier bad pixel fix on cropped data (default=True)
-
+
+
+Creating ASDF files
+^^^^^^^^^^^^^^^^^^^
+The optional arguments `bandpass` and `affine2d` must be written to `ASDF `_
+files to be used by the step. The step expects the contents to be stored with particular keys but the format is not currently
+enforced by a schema; incorrect ASDF file contents will cause the step to revert back to the defaults for each argument.
+
+Examples of how to create ASDF files containing the properly formatted information for each of the arguments follows.
+
+.. code-block:: python
+
+ # Create a F480M filter + Vega bandpass ASDF file
+
+ import asdf
+ from jwst.ami import utils
+ from stdatamodels.jwst import datamodels
+ from synphot import SourceSpectrum
+
+ # F480M throughput reference file from JWST CRDS
+ throughput_file = 'jwst_niriss_throughput_0012.fits'
+ nspecbin=19
+ throughput_model = datamodels.open(throughput_file)
+
+ filt_spec = utils.get_filt_spec(throughput_model)
+ src_spec = SourceSpectrum.from_vega()
+ bandpass = utils.combine_src_filt(filt_spec,
+ src_spec,
+ trim=0.01,
+ nlambda=nspecbin)
+
+ # This bandpass has shape (19, 2); each row is [throughput, wavelength]
+ asdf_name = 'bandpass_f480m_vega.asdf'
+ tree = {"bandpass": bandpass}
+ with open(asdf_name, 'wb') as fh:
+ af = asdf.AsdfFile(tree)
+ af.write_to(fh)
+ af.close()
+ throughput_model.close()
+
+
+.. code-block:: python
+
+ # Create an affine transform ASDF file to use for the model
+
+ import asdf
+ tree = {
+ 'mx': 1., # dimensionless x-magnification
+ 'my': 1., # dimensionless y-magnification
+ 'sx': 0., # dimensionless x shear
+ 'sy': 0., # dimensionless y shear
+ 'xo': 0., # x-offset in pupil space
+ 'yo': 0., # y-offset in pupil space
+ 'rotradccw': None
+ }
+
+ affineasdf = 'affine.asdf'
+
+ with open(affineasdf, 'wb') as fh:
+ af = asdf.AsdfFile(tree)
+ af.write_to(fh)
+ af.close()
+
+
Inputs
------
@@ -70,7 +133,7 @@ The ``ami_analyze`` step itself does not accept an ASN as input.
Outputs
-------
-The ``ami_analyze`` step produces three output files. The first two (``_ami-oi.fits`` and ``_amimulti-oi.fits``) contain the interferometric observables, and the third (``_amilg.fits``) contains the data, LG model, and residuals. These are described in more detail below.
+The ``ami_analyze`` step produces three output files. The first two (``_ami-oi.fits`` and ``_amimulti-oi.fits``) contain the interferometric observables, and the third (``_amilg.fits``) contains the data, LG model, and residuals. All products contain a ``PRIMARY`` HDU containing header keywords but no science data, as well as an ``ASDF`` extension. The files are described in more detail below.
The output file name syntax is exposure-based, using the input file name as the root, with
the addition of the association candidate ID and the "_ami-oi", "_amimulti-oi", or "amilg" product type suffix, e.g.
@@ -81,20 +144,19 @@ Interferometric observables
:Data model: `~jwst.datamodels.AmiOIModel`
:File suffix: _ami-oi.fits, _amimulti-oi.fits
-. The inteferometric observables are saved as OIFITS files, a registered FITS format
+The inteferometric observables are saved as OIFITS files, a registered FITS format
for optical interferometry, containing the following list of extensions:
1) ``OI_ARRAY``: AMI subaperture information
2) ``OI_TARGET``: target properties
-3) ``OI_T3``: extracted closure amplitudes, phases
+3) ``OI_T3``: extracted closure amplitudes, triple-product phases
4) ``OI_VIS``: extracted visibility (fringe) amplitudes, phases
5) ``OI_VIS2``: squared visibility (fringe) amplitudes
6) ``OI_WAVELENGTH``: filter information
For more information on the format and contents of OIFITS files, see the `OIFITS2 standard `_.
-The _ami-oi.fits file contains tables of median observables over all integrations of the input file. Errors
-are computed as the sigma-clipped standard deviation over integrations.
+The _ami-oi.fits file contains tables of observables averaged over all integrations of the input file. The error is taken to be the standard error of the mean, where the variance is the covariance between amplitudes and phases (e.g. fringe amplitudes and fringe phases, closure phases and triple-product amplitudes).
The _amimulti-oi.fits file contains observables for each integration, and does not contain error estimates. The
structure is the same as the _ami-oi.fits file, but the following data columns are 2D, with the second dimension being
the number of integrations: "PISTONS", "PIST_ERR", "VISAMP", "VISAMPERR", "VISPHI", "VISPHIERR", "VIS2DATA", "VIS2ERR", "T3AMP", "T3AMPERR", "T3PHI", "T3PHIERR".
diff --git a/docs/jwst/pipeline/calwebb_ami3.rst b/docs/jwst/pipeline/calwebb_ami3.rst
index c3f494861f..d8aeb69cbc 100644
--- a/docs/jwst/pipeline/calwebb_ami3.rst
+++ b/docs/jwst/pipeline/calwebb_ami3.rst
@@ -7,8 +7,7 @@ calwebb_ami3: Stage 3 Aperture Masking Interferometry (AMI) Processing
:Alias: calwebb_ami3
The stage 3 AMI pipeline is applied to associations of calibrated NIRISS AMI exposures.
-It computes fringe parameters for individual exposures, averages the fringe results from
-multiple exposures, and, optionally, corrects science target fringe parameters using the
+It computes fringe parameters for individual exposures, and, optionally, corrects science target fringe parameters using the
fringe results from reference PSF targets.
The steps applied by the ``calwebb_ami3`` pipeline are shown below.
@@ -83,15 +82,14 @@ Interferometric observables
:Data model: `~jwst.datamodels.AmiOIModel`
:File suffix: _ami-oi.fits
-
-For every input exposure, the fringe parameters and closure phases calculated
+For every input exposure, the interferometric observables calculated
by the :ref:`ami_analyze ` step are saved to an "_ami-oi.fits" product file,
-which is a FITS table of median observables over all integrations of the input file.
+which is a FITS table of averaged observables over all integrations of the input file.
Product names use the input "_calints" exposure-based file name, with the association candidate ID
included and the product type changed to "_ami-oi.fits", e.g.
"jw93210001001_03101_00001_nis_a0003_ami-oi.fits."
-Normalized Interferometric Observables
+Normalized interferometric observables
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
:Data model: `~jwst.datamodels.AmiOIModel`
:File suffix: _aminorm-oi.fits
@@ -102,3 +100,8 @@ via the :ref:`ami_normalize ` step, and will be saved to an
product file. This file has the same FITS table format as the "_ami-oi.fits" products.
The file name root uses the source-based output product name given in the ASN file,
e.g. "jw93210-a0003_t001_niriss_f480m-nrm_aminorm-oi.fits."
+
+.. note::
+
+ Users may wish to run the :ref:`ami_analyze step ` separately for access to flexible input arguments and to save additional diagnostic output products. See the step documentation for more details.
+
diff --git a/jwst/ami/ami_analyze_step.py b/jwst/ami/ami_analyze_step.py
index ddf6d3fd5f..3330e6aa25 100755
--- a/jwst/ami/ami_analyze_step.py
+++ b/jwst/ami/ami_analyze_step.py
@@ -2,6 +2,10 @@
from ..stpipe import Step
from . import ami_analyze
+from . import utils
+
+import numpy as np
+import asdf
__all__ = ["AmiAnalyzeStep"]
@@ -16,11 +20,11 @@ class AmiAnalyzeStep(Step):
rotation = float(default=0.0) # Rotation initial guess [deg]
psf_offset = string(default='0.0 0.0') # PSF offset values to use to create the model array
rotation_search = string(default='-3 3 1') # Rotation search parameters: start, stop, step
- bandpass = any(default=None) # Synphot spectrum or array to override filter/source
+ bandpass = string(default=None) # ASDF file containing array to override filter/source
usebp = boolean(default=True) # If True, exclude pixels marked DO_NOT_USE from fringe fitting
firstfew = integer(default=None) # If not None, process only the first few integrations
chooseholes = string(default=None) # If not None, fit only certain fringes e.g. ['B4','B5','B6','C2']
- affine2d = any(default=None) # None or user-defined Affine2d object
+ affine2d = string(default=None) # ASDF file containing user-defined affine parameters
run_bpfix = boolean(default=True) # Run Fourier bad pixel fix on cropped data
"""
@@ -32,6 +36,96 @@ def save_model(self, model, *args, **kwargs):
kwargs["suffix"] = ["ami-oi", "amimulti-oi", "amilg"][kwargs.pop("idx")]
return Step.save_model(self, model, *args, **kwargs)
+ def override_bandpass(self):
+ """
+ Read bandpass from asdf file and use it to override the default.
+
+ Expects an array of [effstims, wave_m]
+ (i.e. np.array((effstims,wave_m)).T) stored as 'bandpass' in asdf file,
+ where effstims are normalized countrates (unitless) and wave_m are the
+ wavelengths across the filter at which to compute the model (meters).
+
+ Returns
+ -------
+ bandpass: array
+ Array of [countrates, wavelengths]
+ """
+
+ try:
+ with asdf.open(self.bandpass, lazy_load=False) as af:
+ bandpass = np.array(af['bandpass'])
+
+ # assume it is an array of the correct shape
+ wavemin = np.min(bandpass[:,1])
+ wavemax = np.max(bandpass[:,1])
+ self.log.info('User-defined bandpass provided:')
+ self.log.info('\tOVERWRITING ALL NIRISS-SPECIFIC FILTER/BANDPASS VARIABLES')
+ self.log.info(f'Using {bandpass.shape[0]} wavelengths for fit.')
+ self.log.info(f'Wavelength min: {wavemin:.3e} \t Wavelength max: {wavemax:.3e}')
+
+ # update attribute and return
+ self.bandpass = bandpass
+ return bandpass
+
+ except FileNotFoundError:
+ message = f'File {self.bandpass} could not be found at the specified location.'
+ raise Exception(message)
+
+ except KeyError:
+ message1 = 'ASDF file does not contain the required "bandpass" key. '
+ message2 = 'See step documentation for info on creating a custom bandpass ASDF file.'
+ raise Exception((message1 + message2))
+
+ except (IndexError, ValueError):
+ message1 = f'Could not use bandpass from {self.bandpass}. It may have the wrong shape. '
+ message2 = 'See documentation for info on creating a custom bandpass ASDF file.'
+ raise Exception((message1 + message2))
+
+ def override_affine2d(self):
+ """
+ Read user-input affine transform from ASDF file.
+
+ Makes an Affine2d object (see utils.Affine2D class).
+ Input should contain mx,my,sx,sy,xo,yo,rotradccw.
+ """
+ try:
+ with asdf.open(self.affine2d, lazy_load=False) as af:
+ affine2d = utils.Affine2d(
+ mx = af['mx'],
+ my = af['my'],
+ sx = af['sx'],
+ sy = af['sy'],
+ xo = af['xo'],
+ yo = af['yo'],
+ rotradccw = af['rotradccw']
+ )
+ self.log.info(f'Using affine transform from ASDF file {self.affine2d}')
+ # now self.affine2d updated from string to object
+ self.affine2d = affine2d
+ return affine2d
+
+ except FileNotFoundError:
+ self.log.info(f'File {self.affine2d} could not be found at the specified location.')
+ self.log.info('\t **** DEFAULTING TO USE IDENTITY TRANSFORM ****')
+ affine2d = None
+
+ except KeyError:
+ message1 = 'ASDF file does not contain all of the required keys: mx, my, sx, sy ,xo, yo, rotradccw. '
+ message2 = 'See step documentation for info on creating a custom affine2d ASDF file.'
+ self.log.info((message1 + message2))
+ self.log.info('\t **** DEFAULTING TO USE IDENTITY TRANSFORM ****')
+ affine2d = None
+
+ except (IndexError, TypeError, ValueError):
+ message1 = f'Could not use affine2d from {self.affine2d}. '
+ message2 = 'See documentation for info on creating a custom bandpass ASDF file.'
+ self.log.info((message1 + message2))
+ self.log.info('\t **** DEFAULTING TO USE IDENTITY TRANSFORM ****')
+ affine2d = None
+
+ self.affine2d = affine2d
+ return affine2d
+
def process(self, input):
"""
Performs analysis of an AMI mode exposure by applying the LG algorithm.
@@ -83,6 +177,15 @@ def process(self, input):
raise RuntimeError("No THROUGHPUT reference file found. "
"ami_analyze cannot continue.")
+ # If there's a user-defined bandpass or affine, handle it
+ if bandpass is not None:
+ bandpass = self.override_bandpass()
+
+ if affine2d is not None:
+ # if it is None, handled in apply_LG_plus
+ affine2d = self.override_affine2d()
+
+
# Get the name of the NRM reference file to use
nrm_reffile = self.get_reference_file(input_model, 'nrm')
self.log.info(f'Using NRM reference file {nrm_reffile}')
diff --git a/jwst/ami/instrument_data.py b/jwst/ami/instrument_data.py
index 7b5e11614d..26d4f73d08 100755
--- a/jwst/ami/instrument_data.py
+++ b/jwst/ami/instrument_data.py
@@ -179,10 +179,11 @@ def read_data_model(self, input_model):
pscale_deg = np.mean([pscaledegx, pscaledegy])
self.pscale_rad = np.deg2rad(pscale_deg)
self.pscale_mas = pscale_deg * (60 * 60 * 1000)
- self.pav3 = input_model.meta.pointing.pa_v3
+
+ self.roll_ref = input_model.meta.wcsinfo.roll_ref
self.vparity = input_model.meta.wcsinfo.vparity
self.v3iyang = input_model.meta.wcsinfo.v3yangle
- self.parangh = input_model.meta.wcsinfo.roll_ref
+
self.crpix1 = input_model.meta.wcsinfo.crpix1
self.crpix2 = input_model.meta.wcsinfo.crpix2
self.pupil = input_model.meta.instrument.pupil
@@ -223,7 +224,7 @@ def read_data_model(self, input_model):
log.warning("All integrations will be analyzed")
self.nwav = scidata.shape[0]
[self.wls.append(self.wls[0]) for f in range(self.nwav - 1)]
- # Rotate mask hole centers by pav3 + v3i_yang to be in equatorial coordinates
+ # Rotate mask hole centers by roll_ref + v3i_yang to be in equatorial coordinates
ctrs_sky = self.mast2sky()
oifctrs = np.zeros(self.mask.ctrs.shape)
oifctrs[:, 0] = ctrs_sky[:, 1].copy() * -1
@@ -331,11 +332,19 @@ def reset_nwav(self, nwav):
def mast2sky(self):
"""
- Rotate hole center coordinates:
- Clockwise by the V3 position angle - V3I_YANG from north in degrees if VPARITY = -1
- Counterclockwise by the V3 position angle - V3I_YANG from north in degrees if VPARITY = 1
+ Rotate hole center coordinates.
+
+ Rotation of coordinates is:
+ Clockwise by the ROLL_REF + V3I_YANG from north in degrees if VPARITY = -1
+ Counterclockwise by the ROLL_REF + V3I_YANG from north in degrees if VPARITY = 1
Hole center coords are in the V2, V3 plane in meters.
+ Notes
+ -----
+ Nov. 2024 email discussion with Tony Sohn, Paul Goudfrooij confirmed V2/V3 coordinate
+ rotation back to "North up" equatorial orientation should use ROLL_REF + V3I_YANG
+ (= PA_APER).
+
Returns
-------
ctrs_rot: array
@@ -347,13 +356,13 @@ def mast2sky(self):
# NOT used for the fringe fitting itself
mask_ctrs = utils.rotate2dccw(mask_ctrs, np.pi / 2.0)
vpar = self.vparity # Relative sense of rotation between Ideal xy and V2V3
- rot_ang = self.pav3 - self.v3iyang # subject to change!
+ rot_ang = self.roll_ref + self.v3iyang
- if self.pav3 == 0.0:
+ if self.roll_ref == 0.0:
return mask_ctrs
else:
# Using rotate2sccw, which rotates **vectors** CCW in a fixed coordinate system,
- # so to rotate coord system CW instead of the vector, reverse sign of rotation angle. Double-check comment
+ # so to rotate coord system CW instead of the vector, reverse sign of rotation angle.
if vpar == -1:
# rotate clockwise
ctrs_rot = utils.rotate2dccw(mask_ctrs, np.deg2rad(-rot_ang))
diff --git a/jwst/ami/leastsqnrm.py b/jwst/ami/leastsqnrm.py
index d1c1a89524..79983d34e7 100644
--- a/jwst/ami/leastsqnrm.py
+++ b/jwst/ami/leastsqnrm.py
@@ -746,8 +746,8 @@ def populate_antisymmphasearray(deltaps, n=7):
deltaps: 1D float array
pistons between each pair of holes
- n: integer
- number of holes
+ n: integer, optional
+ number of holes (default=7)
Returns
-------
@@ -778,8 +778,8 @@ def populate_symmamparray(amps, n=7):
amps: 1D float array
fringe visibility between each pair of holes
- n: integer
- number of holes
+ n: integer, optional
+ number of holes (default=7)
Returns
-------
@@ -800,6 +800,43 @@ def populate_symmamparray(amps, n=7):
return arr
+def t3_amplitudes(amps, n=7):
+ """
+ Populate the triple-product amplitude array
+ (NOT closure amplitudes)
+
+ Parameters
+ ----------
+ amps: 1D float array
+ fringe visibility between each pair of holes
+
+ n: integer, optional
+ number of holes (default=7)
+
+ Returns
+ -------
+ cpamps: 1D float array
+ triple product amplitude array
+ """
+
+ arr = populate_symmamparray(amps, n=n)
+
+ cpamps = np.zeros(int(comb(n, 3)))
+
+ nn = 0
+ for kk in range(n - 2):
+ for ii in range(n - kk - 2):
+ for jj in range(n - kk - ii - 2):
+ cpamps[nn + jj] = (
+ arr[kk, ii + kk + 1]
+ * arr[ii + kk + 1, jj + ii + kk + 2]
+ * arr[jj + ii + kk + 2, kk]
+ )
+
+ nn += jj + 1
+
+ return cpamps
+
def redundant_cps(deltaps, n=7):
"""
@@ -810,8 +847,8 @@ def redundant_cps(deltaps, n=7):
deltaps: 1D float array
pistons between each pair of holes
- n: integer
- number of holes
+ n: integer, optional
+ number of holes (default=7)
Returns
-------
@@ -906,8 +943,8 @@ def closure_amplitudes(amps, n=7):
amps: 1D float array
fringe amplitudes
- N: integer
- number of holes
+ n: integer, optional
+ number of holes (default=7)
Returns
-------
@@ -934,3 +971,38 @@ def closure_amplitudes(amps, n=7):
nn = nn + ll + 1
return cas
+
+def q4_phases(deltaps, n=7):
+ """
+ Calculate phases for each set of 4 holes
+
+ Parameters
+ ----------
+ deltaps: 1D float array
+ pistons between each pair of holes
+
+ n: integer, optional
+ number of holes (default=7)
+
+ Returns
+ -------
+ quad_phases: 1D float array
+ quad phases
+ """
+ arr = populate_antisymmphasearray(deltaps, n=n) # fringe phase array
+ nn = 0
+ quad_phases = np.zeros(int(comb(n, 4)))
+
+ for ii in range(n - 3):
+ for jj in range(n - ii - 3):
+ for kk in range(n - jj - ii - 3):
+ for ll in range(n - jj - ii - kk - 3):
+ quad_phases[nn + ll] = (
+ arr[ii, jj + ii + 1]
+ + arr[ll + ii + jj + kk + 3, kk + jj + ii + 2]
+ - arr[ii, kk + ii + jj + 2]
+ - arr[jj + ii + 1, ll + ii + jj + kk + 3]
+ )
+ nn = nn + ll + 1
+
+ return quad_phases
diff --git a/jwst/ami/lg_model.py b/jwst/ami/lg_model.py
index 20ee86f33a..2981baaa8a 100755
--- a/jwst/ami/lg_model.py
+++ b/jwst/ami/lg_model.py
@@ -396,7 +396,9 @@ def fit_image(
self.fringeamp, self.fringephase = leastsqnrm.tan2visibilities(self.soln)
self.fringepistons = utils.fringes2pistons(self.fringephase, len(self.ctrs))
self.redundant_cps = leastsqnrm.redundant_cps(self.fringephase, n=self.N)
+ self.t3_amplitudes = leastsqnrm.t3_amplitudes(self.fringeamp, n=self.N) # RC 8/24
self.redundant_cas = leastsqnrm.closure_amplitudes(self.fringeamp, n=self.N)
+ self.q4_phases = leastsqnrm.q4_phases(self.fringephase, n=self.N) # RC 8/24
def create_modelpsf(self):
"""
diff --git a/jwst/ami/oifits.py b/jwst/ami/oifits.py
index e0bad2e233..b373babd24 100644
--- a/jwst/ami/oifits.py
+++ b/jwst/ami/oifits.py
@@ -8,16 +8,20 @@
CalibOifits class: takes two AmiOIModel datamodels and produces a final calibrated datamodel.
"""
import numpy as np
-
from scipy.special import comb
from astropy.stats import sigma_clipped_stats
from astropy.time.core import Time
+import logging
from stdatamodels.jwst import datamodels
+from . import leastsqnrm
+
+log = logging.getLogger(__name__)
+log.setLevel(logging.DEBUG)
class RawOifits:
- def __init__(self, fringefitter, method="median"):
+ def __init__(self, fringefitter, method="mean"):
"""
Class to store AMI data in the format required to write out to OIFITS files
Angular quantities of input are in radians from fringe fitting; converted to degrees for saving.
@@ -28,7 +32,7 @@ def __init__(self, fringefitter, method="median"):
Object containing nrm_list attribute (list of nrm objects)
and other info needed for OIFITS files
method: string
- Method to average observables: mean or median. Default median.
+ Method to average observables: mean or median. Default mean.
Notes
-----
@@ -41,44 +45,262 @@ def __init__(self, fringefitter, method="median"):
self.nslices = len(self.fringe_fitter.nrm_list) # n ints
self.n_baselines = int(comb(self.n_holes, 2)) # 21
self.n_closure_phases = int(comb(self.n_holes, 3)) # 35
- self.n_closure_amplitudes = int(comb(self.n_holes, 4))
+ self.n_closure_amplitudes = int(comb(self.n_holes, 4)) # also 35
self.method = method
+ if self.method not in ["mean", "median", "multi"]:
+ log.warning('"method" for saving OIFITS file must be one of ["mean","median","multi"]. Defaulting to use mean!')
+ self.method = "mean"
self.ctrs_eqt = self.fringe_fitter.instrument_data.ctrs_eqt
self.ctrs_inst = self.fringe_fitter.instrument_data.ctrs_inst
- self.pa = (
- self.fringe_fitter.instrument_data.pav3
- ) # header pav3, not including v3i_yang??
self.bholes, self.bls = self._makebaselines()
self.tholes, self.tuv = self._maketriples_all()
+ self.qholes, self.quads = self._makequads_all()
def make_obsarrays(self):
"""
Make arrays of observables of the correct shape for saving to datamodels
"""
- # arrays of observables, (nslices,nobservables) shape.
+ # empty arrays of observables, (nslices,nobservables) shape.
self.fringe_phases = np.zeros((self.nslices, self.n_baselines))
self.fringe_amplitudes = np.zeros((self.nslices, self.n_baselines))
self.closure_phases = np.zeros((self.nslices, self.n_closure_phases))
+ self.t3_amplitudes = np.zeros((self.nslices, self.n_closure_phases))
+ self.q4_phases = np.zeros((self.nslices, self.n_closure_amplitudes))
self.closure_amplitudes = np.zeros((self.nslices, self.n_closure_amplitudes))
self.pistons = np.zeros((self.nslices, self.n_holes))
# model parameters
self.solns = np.zeros((self.nslices, 44))
+ # populate with each integration's observables
for i, nrmslc in enumerate(self.fringe_fitter.nrm_list):
- self.fringe_phases[i, :] = np.rad2deg(nrmslc.fringephase) # FPs in degrees
+ self.fringe_phases[i, :] = nrmslc.fringephase # FPs in radians
self.fringe_amplitudes[i, :] = nrmslc.fringeamp
- self.closure_phases[i, :] = np.rad2deg(nrmslc.redundant_cps) # CPs in degrees
+ self.closure_phases[i, :] = nrmslc.redundant_cps # CPs in radians
+ self.t3_amplitudes[i, :] = nrmslc.t3_amplitudes
+ self.q4_phases[i, :] = nrmslc.q4_phases # quad phases in radians
self.closure_amplitudes[i, :] = nrmslc.redundant_cas
- self.pistons[i, :] = np.rad2deg(
- nrmslc.fringepistons
- ) # segment pistons in degrees
+ self.pistons[i, :] = nrmslc.fringepistons # segment pistons in radians
self.solns[i, :] = nrmslc.soln
self.fringe_amplitudes_squared = self.fringe_amplitudes ** 2 # squared visibilities
+ def rotate_matrix(self, cov_mat, theta):
+ """
+ Rotate a covariance matrix by an angle.
+
+ Parameters
+ ----------
+ cov_mat: array
+ The matrix to be rotated
+ theta: float
+ Angle by which to rotate the matrix (radians)
+
+ Returns
+ -------
+ cv_rotated: array
+ The rotated matrix
+
+ """
+ c, s = np.cos(theta), np.sin(theta)
+ r_mat = [[c, -s],
+ [s, c]]
+ # coordinate rotation from real/imaginary to absolute value/phase (modulus/argument)
+ return np.linalg.multi_dot([np.transpose(r_mat), cov_mat, r_mat])
+
+
+ def average_observables(self, averfunc):
+ """
+ Average all the observables.
+
+ Calculate covariance matrices between fringe amplitudes/fringe phases,
+ and between triple product amps/closure phases, and closure amplitudes/quad phases.
+ Convert r, theta (modulus, phase) to x,y. Calculate cov(x,y). Rotate resulting
+ 2x2 matrix back to r, theta. Take sqrt of relevant covariance element to be error.
+ This must be done with phases in radians.
+
+ Parameters
+ ----------
+ averfunc: function
+ np.mean (default) or np.median
+
+ Returns
+ -------
+ avg_sqv: array
+ Averaged squared visibilites
+ err_sqv: array
+ Standard error of the mean of averaged squared visibilities
+ avg_fa: array
+ Averaged fringe (visibility) amplitudes
+ err_fa: array
+ Standard error of the mean of averaged fringe (visibility) amplitudes
+ avg_fp: array
+ Averaged fringe phases (rad)
+ err_fp: array
+ Standard error of the mean of averaged fringe phases (rad)
+ avg_cp: array
+ Averaged closure phases (rad)
+ err_cp: array
+ Standard error of the mean of averaged closure phases (rad)
+ avg_t3amp: array
+ Averaged triple amplitudes
+ err_t3amp: array
+ Standard error of the mean of averaged triple amplitudes
+ avg_ca: array
+ Averaged closure amplitudes
+ err_ca: array
+ Standard error of the mean of averaged closure amplitudes
+ avg_q4phi: array
+ Averaged quad phases
+ err_q4phi: array
+ Standard error of the mean of averaged quad phases
+ avg_pist: array
+ Averaged segment pistons
+ err_pist: array
+ Standard error of the mean of averaged segment pistons
+
+ """
+ covmats_fringes, covmats_triples, covmats_quads = self.observable_covariances(averfunc)
+
+ if self.method == 'mean':
+ avg_fa, _, std_fa = sigma_clipped_stats(self.fringe_amplitudes, axis=0)
+ avg_fp, _, std_fp = sigma_clipped_stats(self.fringe_phases, axis=0)
+ avg_sqv, _, std_sqv = sigma_clipped_stats(self.fringe_amplitudes**2, axis=0)
+ avg_pist, _, err_pist = sigma_clipped_stats(self.pistons, axis=0)
+ else: # median
+ _, avg_fa, std_fa = sigma_clipped_stats(self.fringe_amplitudes, axis=0) # 21. std_fa is just for comparing to covariance
+ _, avg_fp, std_fp = sigma_clipped_stats(self.fringe_phases, axis=0) # 21
+ _, avg_sqv, std_sqv = sigma_clipped_stats(self.fringe_amplitudes**2, axis=0)
+ _, avg_pist, err_pist = sigma_clipped_stats(self.pistons, axis=0)
+
+ err_pist = err_pist/np.sqrt(self.nslices) # standard error of the mean
+ err_fa, err_fp = self.err_from_covmat(covmats_fringes)
+
+ # calculate squared visibility (fringe) amplitude uncertainties correctly
+ err_sqv = (2 * avg_fa * err_fa) / np.sqrt(self.nslices)
+
+ # calculate triple and quad quantities from **averaged** fringe amps and phases
+ avg_t3amp = leastsqnrm.t3_amplitudes(avg_fa, n=self.n_holes)
+ avg_cp = leastsqnrm.redundant_cps(avg_fp, n=self.n_holes)
+ err_t3amp, err_cp = self.err_from_covmat(covmats_triples)
+
+ avg_ca = leastsqnrm.closure_amplitudes(avg_fa, n=self.n_holes)
+ avg_q4phi = leastsqnrm.q4_phases(avg_fp, n=self.n_holes)
+ err_ca, err_q4phi = self.err_from_covmat(covmats_quads)
+
+ return (avg_sqv, err_sqv, avg_fa, err_fa,
+ avg_fp, err_fp, avg_cp, err_cp,
+ avg_t3amp, err_t3amp, avg_ca, err_ca,
+ avg_q4phi, err_q4phi, avg_pist, err_pist)
+
+ def err_from_covmat(self, covmatlist):
+ """
+ Derive observable errors from their covariance matrices.
+
+ Return sqrt of [0,0] and [1,1] elements of each of a list of covariance matrices,
+ divided by sqrt(N_ints), for use as observable errors (standard error of the mean).
+ If using median, error calculation is questionable because this is NOT the standard
+ error of the median.
+
+ Parameters
+ ----------
+ covmatlist: array
+ array of covariance matrices for each baseline/triple/quad
+ shape e.g. (21,2,2) or (35,2,2)
+
+ Returns
+ -------
+ err_00: array
+ standard errors of the mean of the first observable. shape e.g. (21)
+ err_11: array
+ standard errors of the mean of the second observable. shape e.g. (21)
+ """
+ err_00 = np.sqrt(np.array([covmat[0,0] for covmat in covmatlist]))/np.sqrt(self.nslices)
+ err_11 = np.sqrt(np.array([covmat[1,1] for covmat in covmatlist]))/np.sqrt(self.nslices)
+
+ return err_00, err_11
+
+ def observable_covariances(self, averfunc):
+ """
+ Calculate covariance matrices from each pair of observables.
+
+ For each baseline/triple/quad, calculate covariance between each
+ fringe amplitude/phase quantity.
+
+ Parameters
+ ----------
+ averfunc: function
+ np.mean (default) or np.median
+
+ Returns
+ -------
+ cov_mat_fringes: array
+ Array of 21 covariance matrices for fringes (amplitudes, phases)
+ cov_mat_triples: array
+ Array of 35 covariance matrices for triples (t3 amplitudes, closure phases)
+ cov_mat_quads: array
+ Array of 35 covariance matrices for quads (closure amplitudes, quad phases)
+
+ """
+ # loop over 21 baselines
+ cov_mat_fringes = []
+ # These operate on all slices at once, not on already-averaged integrations
+ for bl in np.arange(self.n_baselines):
+ fringeamps = self.fringe_amplitudes[:,bl]
+ fringephases = self.fringe_phases[:,bl]
+ covmat = self.cov_r_theta(fringeamps, fringephases, averfunc)
+ cov_mat_fringes.append(covmat)
+ # loop over 35 triples
+ cov_mat_triples = []
+ for triple in np.arange(self.n_closure_phases):
+ tripamp = self.t3_amplitudes[:,triple]
+ triphase = self.closure_phases[:,triple]
+ covmat = self.cov_r_theta(tripamp, triphase, averfunc)
+ cov_mat_triples.append(covmat)
+ # loop over 35 quads
+ cov_mat_quads = []
+ for quad in np.arange(self.n_closure_amplitudes):
+ quadamp = self.closure_amplitudes[:,quad]
+ quadphase = self.q4_phases[:,quad]
+ covmat = self.cov_r_theta(quadamp, quadphase, averfunc)
+ cov_mat_quads.append(covmat)
+
+ # covmats to be written to oifits. store in rawoifits object? TBD
+ # lists of cov mats have shape e.g. (21, 2, 2) or (35, 2, 2)
+
+ return np.array(cov_mat_fringes), np.array(cov_mat_triples), np.array(cov_mat_quads)
+
+
+ def cov_r_theta(self, rr, theta, averfunc):
+ """
+ Calculate covariance in polar coordinates.
+
+ Calculate covariance in x, y coordinates, then rotate covariance matrix
+ by **average** phase (over integrations) to get matrix in (r,theta).
+
+ Parameters
+ ----------
+ rr: array
+ complex number modulus
+ theta: array
+ complex number phase
+ averfunc: function
+ np.mean (default) or np.median
+
+ Returns
+ -------
+ cov_mat_r_theta: array (2,2)
+ Covariance matrix in r, theta coordinates
+
+ """
+ xx = rr * np.cos(theta)
+ yy = rr * np.sin(theta)
+ cov_mat_xy = np.cov(xx, yy)
+ return self.rotate_matrix(cov_mat_xy, averfunc(theta))
+
+
def make_oifits(self):
"""
Perform final manipulations of observable arrays, calculate uncertainties, and
@@ -114,16 +336,14 @@ def make_oifits(self):
flagVis = [False] * self.n_baselines
flagT3 = [False] * self.n_closure_phases
- # do the things done by populate_nrm here
- # average or don't, and get uncertainties
+ # Average observables (or don't), and get uncertainties
# Unwrap phases
shift2pi = np.zeros(self.closure_phases.shape)
shift2pi[self.closure_phases >= 6] = 2 * np.pi
shift2pi[self.closure_phases <= -6] = -2 * np.pi
self.closure_phases -= shift2pi
- if self.method not in ["mean", "median", "multi"]:
- self.method = "median"
+ # Now we are setting up the observables to be written out to OIFITS
# set these as attributes (some may exist and be overwritten)
if self.method == "multi":
self.vis2 = self.fringe_amplitudes_squared.T
@@ -134,35 +354,43 @@ def make_oifits(self):
self.e_visphi = np.zeros(self.visphi.shape)
self.closure_phases = self.closure_phases.T
self.e_cp = np.zeros(self.closure_phases.shape)
+ self.t3amp = self.t3_amplitudes.T
+ self.e_t3amp = np.zeros(self.t3amp.shape)
+ self.q4phi = self.q4_phases.T
+ self.e_q4phi = np.zeros(self.q4_phases.shape)
self.camp = self.closure_amplitudes.T
self.e_camp = np.zeros(self.camp.shape)
self.pist = self.pistons.T
self.e_pist = np.zeros(self.pist.shape)
- # apply sigma-clipping to uncertainties
- # sigma_clipped_stats returns mean, median, stddev. nsigma=3, niters=5
- elif self.method == "median":
- _, self.vis2, self.e_vis2 = sigma_clipped_stats(self.fringe_amplitudes_squared, axis=0)
- _, self.visamp, self.e_visamp = sigma_clipped_stats(self.fringe_amplitudes, axis=0)
- _, self.visphi, self.e_visphi = sigma_clipped_stats(self.fringe_phases, axis=0)
- _, self.closure_phases, self.e_cp = sigma_clipped_stats(self.closure_phases, axis=0)
- _, self.camp, self.e_camp = sigma_clipped_stats(self.closure_amplitudes, axis=0)
- _, self.pist, self.e_pist = sigma_clipped_stats(self.pistons, axis=0)
-
- else: # take the mean
- self.vis2, _, self.e_vis2 = sigma_clipped_stats(self.fringe_amplitudes_squared, axis=0)
- self.visamp, _, self.e_visamp = sigma_clipped_stats(self.fringe_amplitudes, axis=0)
- self.visphi, _, self.e_visphi = sigma_clipped_stats(self.fringe_phases, axis=0)
- self.closure_phases, _, self.e_cp = sigma_clipped_stats(self.closure_phases, axis=0)
- self.camp, _, self.e_camp = sigma_clipped_stats(self.closure_amplitudes, axis=0)
- self.pist, _, self.e_pist = sigma_clipped_stats(self.pistons, axis=0)
+
+ elif self.method == "mean":
+ (self.vis2, self.e_vis2, self.visamp, self.e_visamp,
+ self.visphi, self.e_visphi, self.closure_phases, self.e_cp,
+ self.t3amp, self.e_t3amp, self.camp, self.e_camp,
+ self.q4phi, self.e_q4phi, self.pist, self.e_pist) = self.average_observables(np.mean)
+
+ else: # take the median
+ (self.vis2, self.e_vis2, self.visamp, self.e_visamp,
+ self.visphi, self.e_visphi, self.closure_phases, self.e_cp,
+ self.t3amp, self.e_t3amp, self.camp, self.e_camp,
+ self.q4phi, self.e_q4phi, self.pist, self.e_pist) = self.average_observables(np.median)
+
+ # Convert angular quantities from radians to degrees
+ self.visphi = np.rad2deg(self.visphi)
+ self.e_visphi = np.rad2deg(self.e_visphi)
+ self.closure_phases = np.rad2deg(self.closure_phases)
+ self.e_cp = np.rad2deg(self.e_cp)
+ self.q4phi = np.rad2deg(self.q4phi)
+ self.e_q4phi = np.rad2deg(self.e_q4phi)
+ self.pist = np.rad2deg(self.pist)
+ self.e_pist = np.rad2deg(self.e_pist)
# prepare arrays for OI_ARRAY ext
self.staxy = instrument_data.ctrs_inst
- N_ap = len(self.staxy)
- tel_name = ["A%i" % x for x in np.arange(N_ap) + 1]
+ tel_name = ["A%i" % x for x in np.arange(self.n_holes) + 1]
sta_name = tel_name
- diameter = [0] * N_ap
+ diameter = [0] * self.n_holes
staxyz = []
for x in self.staxy:
@@ -170,105 +398,108 @@ def make_oifits(self):
line = [a[0], a[1], 0]
staxyz.append(line)
- sta_index = np.arange(N_ap) + 1
+ sta_index = np.arange(self.n_holes) + 1
pscale = instrument_data.pscale_mas / 1000.0 # arcsec
- isz = self.fringe_fitter.scidata.shape[
- 1
- ] # Size of the image to extract NRM data
- fov = [pscale * isz] * N_ap
- fovtype = ["RADIUS"] * N_ap
+ # Size of the image to extract NRM data
+ isz = self.fringe_fitter.scidata.shape[1]
+ fov = [pscale * isz] * self.n_holes
+ fovtype = ["RADIUS"] * self.n_holes
- m = datamodels.AmiOIModel()
- self.init_oimodel_arrays(m)
+ oim = datamodels.AmiOIModel()
+ self.init_oimodel_arrays(oim)
# primary header keywords
- m.meta.telescope = instrument_data.telname
- m.meta.origin = "STScI"
- m.meta.instrument.name = instrument_data.instrument
- m.meta.program.pi_name = instrument_data.pi_name
- m.meta.target.proposer_name = instrument_data.proposer_name
- m.meta.observation.date = observation_date.fits
- m.meta.oifits.array_name = instrument_data.arrname
- m.meta.oifits.instrument_mode = instrument_data.pupil
+ oim.meta.telescope = instrument_data.telname
+ oim.meta.origin = "STScI"
+ oim.meta.instrument.name = instrument_data.instrument
+ oim.meta.program.pi_name = instrument_data.pi_name
+ oim.meta.target.proposer_name = instrument_data.proposer_name
+ oim.meta.observation.date = observation_date.fits
+ oim.meta.oifits.array_name = instrument_data.arrname
+ oim.meta.oifits.instrument_mode = instrument_data.pupil
+
+ oim.meta.ami.roll_ref = instrument_data.roll_ref
+ oim.meta.ami.v3yangle = instrument_data.v3iyang
+ oim.meta.ami.vparity = instrument_data.vparity
# oi_array extension data
- m.array["TEL_NAME"] = tel_name
- m.array["STA_NAME"] = sta_name
- m.array["STA_INDEX"] = sta_index
- m.array["DIAMETER"] = diameter
- m.array["STAXYZ"] = staxyz
- m.array["FOV"] = fov
- m.array["FOVTYPE"] = fovtype
- m.array["CTRS_EQT"] = instrument_data.ctrs_eqt
- m.array["PISTONS"] = self.pist
- m.array["PIST_ERR"] = self.e_pist
+ oim.array["TEL_NAME"] = tel_name
+ oim.array["STA_NAME"] = sta_name
+ oim.array["STA_INDEX"] = sta_index
+ oim.array["DIAMETER"] = diameter
+ oim.array["STAXYZ"] = staxyz
+ oim.array["FOV"] = fov
+ oim.array["FOVTYPE"] = fovtype
+ oim.array["CTRS_EQT"] = instrument_data.ctrs_eqt
+ oim.array["PISTONS"] = self.pist
+ oim.array["PIST_ERR"] = self.e_pist
# oi_target extension data
- m.target['TARGET_ID'] = [1]
- m.target['TARGET'] = instrument_data.objname
- m.target['RAEP0'] = instrument_data.ra
- m.target['DECEP0'] = instrument_data.dec
- m.target['EQUINOX'] = [2000]
- m.target['RA_ERR'] = instrument_data.ra_uncertainty
- m.target['DEC_ERR'] = instrument_data.dec_uncertainty
- m.target['SYSVEL'] = [0]
- m.target['VELTYP'] = ['UNKNOWN']
- m.target['VELDEF'] = ['OPTICAL']
- m.target['PMRA'] = instrument_data.pmra
- m.target['PMDEC'] = instrument_data.pmdec
- m.target['PMRA_ERR'] = [0]
- m.target['PMDEC_ERR'] = [0]
- m.target['PARALLAX'] = [0]
- m.target['PARA_ERR'] = [0]
- m.target['SPECTYP'] = ['UNKNOWN']
+ oim.target['TARGET_ID'] = [1]
+ oim.target['TARGET'] = instrument_data.objname
+ oim.target['RAEP0'] = instrument_data.ra
+ oim.target['DECEP0'] = instrument_data.dec
+ oim.target['EQUINOX'] = [2000]
+ oim.target['RA_ERR'] = instrument_data.ra_uncertainty
+ oim.target['DEC_ERR'] = instrument_data.dec_uncertainty
+ oim.target['SYSVEL'] = [0]
+ oim.target['VELTYP'] = ['UNKNOWN']
+ oim.target['VELDEF'] = ['OPTICAL']
+ oim.target['PMRA'] = instrument_data.pmra
+ oim.target['PMDEC'] = instrument_data.pmdec
+ oim.target['PMRA_ERR'] = [0]
+ oim.target['PMDEC_ERR'] = [0]
+ oim.target['PARALLAX'] = [0]
+ oim.target['PARA_ERR'] = [0]
+ oim.target['SPECTYP'] = ['UNKNOWN']
# oi_vis extension data
- m.vis['TARGET_ID'] = 1
- m.vis['TIME'] = 0
- m.vis['MJD'] = observation_date.mjd
- m.vis['INT_TIME'] = instrument_data.itime
- m.vis['VISAMP'] = self.visamp
- m.vis['VISAMPERR'] = self.e_visamp
- m.vis['VISPHI'] = self.visphi
- m.vis['VISPHIERR'] = self.e_visphi
- m.vis['UCOORD'] = ucoord
- m.vis['VCOORD'] = vcoord
- m.vis['STA_INDEX'] = self._format_STAINDEX_V2(self.bholes)
- m.vis['FLAG'] = flagVis
+ oim.vis['TARGET_ID'] = 1
+ oim.vis['TIME'] = 0
+ oim.vis['MJD'] = observation_date.mjd
+ oim.vis['INT_TIME'] = instrument_data.itime
+ oim.vis['VISAMP'] = self.visamp
+ oim.vis['VISAMPERR'] = self.e_visamp
+ oim.vis['VISPHI'] = self.visphi
+ oim.vis['VISPHIERR'] = self.e_visphi
+ oim.vis['UCOORD'] = ucoord
+ oim.vis['VCOORD'] = vcoord
+ oim.vis['STA_INDEX'] = self._format_STAINDEX_V2(self.bholes)
+ oim.vis['FLAG'] = flagVis
# oi_vis2 extension data
- m.vis2['TARGET_ID'] = 1
- m.vis2['TIME'] = 0
- m.vis2['MJD'] = observation_date.mjd
- m.vis2['INT_TIME'] = instrument_data.itime
- m.vis2['VIS2DATA'] = self.vis2
- m.vis2['VIS2ERR'] = self.e_vis2
- m.vis2['UCOORD'] = ucoord
- m.vis2['VCOORD'] = vcoord
- m.vis2['STA_INDEX'] = self._format_STAINDEX_V2(self.bholes)
- m.vis2['FLAG'] = flagVis
+ oim.vis2['TARGET_ID'] = 1
+ oim.vis2['TIME'] = 0
+ oim.vis2['MJD'] = observation_date.mjd
+ oim.vis2['INT_TIME'] = instrument_data.itime
+ oim.vis2['VIS2DATA'] = self.vis2
+ oim.vis2['VIS2ERR'] = self.e_vis2
+ oim.vis2['UCOORD'] = ucoord
+ oim.vis2['VCOORD'] = vcoord
+ oim.vis2['STA_INDEX'] = self._format_STAINDEX_V2(self.bholes)
+ oim.vis2['FLAG'] = flagVis
# oi_t3 extension data
- m.t3['TARGET_ID'] = 1
- m.t3['TIME'] = 0
- m.t3['MJD'] = observation_date.mjd
- m.t3['T3AMP'] = self.camp
- m.t3['T3AMPERR'] = self.e_camp
- m.t3['T3PHI'] = self.closure_phases
- m.t3['T3PHIERR'] = self.e_cp
- m.t3['U1COORD'] = u1coord
- m.t3['V1COORD'] = v1coord
- m.t3['U2COORD'] = u2coord
- m.t3['V2COORD'] = v2coord
- m.t3['STA_INDEX'] = self._format_STAINDEX_T3(self.tholes)
- m.t3['FLAG'] = flagT3
+ oim.t3['TARGET_ID'] = 1
+ oim.t3['TIME'] = 0
+ oim.t3['MJD'] = observation_date.mjd
+ oim.t3['T3AMP'] = self.t3amp
+ oim.t3['T3AMPERR'] = self.e_t3amp
+ oim.t3['T3PHI'] = self.closure_phases
+ oim.t3['T3PHIERR'] = self.e_cp
+ oim.t3['U1COORD'] = u1coord
+ oim.t3['V1COORD'] = v1coord
+ oim.t3['U2COORD'] = u2coord
+ oim.t3['V2COORD'] = v2coord
+ oim.t3['STA_INDEX'] = self._format_STAINDEX_T3(self.tholes)
+ oim.t3['FLAG'] = flagT3
# oi_wavelength extension data
- m.wavelength["EFF_WAVE"] = wl
- m.wavelength["EFF_BAND"] = e_wl
+ oim.wavelength["EFF_WAVE"] = wl
+ oim.wavelength["EFF_BAND"] = e_wl
- return m
+ return oim
def init_oimodel_arrays(self, oimodel):
"""
@@ -381,29 +612,22 @@ def _maketriples_all(self):
-------
tarray: integer array
Triple hole indices (0-indexed),
- float array of two uv vectors in all triangles
+ float array of two uv vectors in all triangles
"""
tlist = []
+ uvlist = []
for i in range(self.n_holes):
for j in range(self.n_holes):
for k in range(self.n_holes):
if i < j and j < k:
tlist.append((i, j, k))
+ uvlist.append(
+ (
+ self.ctrs_eqt[i] - self.ctrs_eqt[j],
+ self.ctrs_eqt[j] - self.ctrs_eqt[k],
+ )
+ )
tarray = np.array(tlist).astype(int)
-
- tname = []
- uvlist = []
- # foreach row of 3 elts...
- for triple in tarray:
- tname.append("{0:d}_{1:d}_{2:d}".format(triple[0], triple[1], triple[2]))
-
- uvlist.append(
- (
- self.ctrs_eqt[triple[0]] - self.ctrs_eqt[triple[1]],
- self.ctrs_eqt[triple[1]] - self.ctrs_eqt[triple[2]],
- )
- )
-
return tarray, np.array(uvlist)
def _makebaselines(self):
@@ -425,6 +649,29 @@ def _makebaselines(self):
bllist.append(self.ctrs_eqt[i] - self.ctrs_eqt[j])
return np.array(blist).astype(int), np.array(bllist)
+ def _makequads_all(self):
+ """
+ Calculate all four-hole combinations (quads)
+
+ Returns
+ -------
+ qarray: int array of four-hole quads (0-based)
+ uvwlist: numpy array of u, v, w vectors for each quad
+ """
+ qlist = []
+ uvwlist = []
+ for i in range(self.n_holes):
+ for j in range(self.n_holes):
+ for k in range(self.n_holes):
+ for q in range(self.n_holes):
+ if i < j and j < k and k < q:
+ qlist.append((i, j, k, q))
+ uvwlist.append((self.ctrs_eqt[i] - self.ctrs_eqt[j],
+ self.ctrs_eqt[j] - self.ctrs_eqt[k],
+ self.ctrs_eqt[k] - self.ctrs_eqt[q]))
+ qarray = np.array(qlist).astype(int)
+ return qarray, np.array(uvwlist)
+
def _format_STAINDEX_T3(self, tab):
"""
Converts sta_index to save oifits T3 in the appropriate format
@@ -581,6 +828,6 @@ def calibrate(self):
# add calibrated header keywords
calname = self.caloimodel.meta.target.proposer_name # name of calibrator star
- self.calib_oimodel.meta.oifits.calib = calname
+ self.calib_oimodel.meta.ami.calibrator_object_id = calname
return self.calib_oimodel
diff --git a/jwst/ami/tests/test_ami_interface.py b/jwst/ami/tests/test_ami_interface.py
index 2da7d16eb3..23e765a343 100644
--- a/jwst/ami/tests/test_ami_interface.py
+++ b/jwst/ami/tests/test_ami_interface.py
@@ -48,7 +48,7 @@ def example_model(mock_nrm_reference_file):
model.meta.program.pi_name = "someone"
model.meta.target.catalog_name = ""
model.meta.visit.start_time = "2022-06-05 12:15:41.5020000"
- model.meta.pointing.pa_v3 = 171.8779402866089
+ model.meta.wcsinfo.roll_ref = 171.8779402866089
model.meta.wcsinfo.v3yangle = 0.56126717
model.meta.filename = "test_calints.fits"
model.meta.instrument.pupil = "NRM"
@@ -63,7 +63,6 @@ def test_ami_analyze_even_oversample_fail(example_model, oversample):
AmiAnalyzeStep.call(example_model, oversample=oversample)
-@pytest.mark.skip("reference files are not currently used")
def test_ami_analyze_no_reffile_fail(monkeypatch, example_model):
"""Make sure that ami_analyze fails if no throughput reffile is available"""
@@ -71,7 +70,7 @@ def mockreturn(input_model, reftype, observatory=None, asn_exptypes=None):
return "N/A"
monkeypatch.setattr(stpipe.crds_client, 'get_reference_file', mockreturn)
- with pytest.raises(RuntimeError, match="No throughput reference file found."):
+ with pytest.raises(RuntimeError, match="No THROUGHPUT reference file found."):
AmiAnalyzeStep.call(example_model)
diff --git a/jwst/ami/utils.py b/jwst/ami/utils.py
index 700f87a2c9..a9fa5e1712 100644
--- a/jwst/ami/utils.py
+++ b/jwst/ami/utils.py
@@ -1322,16 +1322,12 @@ def handle_bandpass(bandpass, throughput_model):
Array of weights, wavelengths used to generate model
"""
if bandpass is not None:
- log.info(
- "User-defined bandpass provided: OVERWRITING ALL NIRISS-SPECIFIC FILTER/BANDPASS VARIABLES"
- )
# bandpass can be user-defined synphot object or appropriate array
if isinstance(bandpass, synphot.spectrum.SpectralElement):
log.info("User-defined synphot spectrum provided")
wl, wt = bandpass._get_arrays(bandpass.waveset)
bandpass = np.array((wt, wl)).T
else:
- log.info("User-defined bandpass array provided")
bandpass = np.array(bandpass)
else:
diff --git a/jwst/assign_wcs/nirspec.py b/jwst/assign_wcs/nirspec.py
index 7902d6c3bc..26149cc262 100644
--- a/jwst/assign_wcs/nirspec.py
+++ b/jwst/assign_wcs/nirspec.py
@@ -726,8 +726,7 @@ def get_open_msa_slits(prog_id, msa_file, msa_metadata_id, dither_position,
quadrant = slitlet_rows[0]['shutter_quadrant']
ycen = j
xcen = slitlet_rows[0]['shutter_row'] # grab the first as they are all the same
- shutter_id = xcen + (ycen - 1) * 365 # shutter numbers in MSA file are 1-indexed
-
+ shutter_id = np.int64(xcen) + (np.int64(ycen) - 1) * 365 # shutter numbers in MSA file are 1-indexed
# Background slits all have source_id=0 in the msa_file,
# so assign a unique id based on the slitlet_id
source_id = slitlet_id
@@ -750,7 +749,7 @@ def get_open_msa_slits(prog_id, msa_file, msa_metadata_id, dither_position,
np.nan_to_num(s['estimated_source_in_shutter_x'], nan=0.5),
np.nan_to_num(s['estimated_source_in_shutter_y'], nan=0.5))
for s in slitlet_rows if s['background'] == 'N'][0]
- shutter_id = xcen + (ycen - 1) * 365 # shutter numbers in MSA file are 1-indexed
+ shutter_id = np.int64(xcen) + (np.int64(ycen) - 1) * 365 # shutter numbers in MSA file are 1-indexed
# y-size
jmin = min([s['shutter_column'] for s in slitlet_rows])
diff --git a/jwst/regtest/test_miri_mrs_spec3.py b/jwst/regtest/test_miri_mrs_spec3.py
index 1b71a4c984..6e27960713 100644
--- a/jwst/regtest/test_miri_mrs_spec3.py
+++ b/jwst/regtest/test_miri_mrs_spec3.py
@@ -100,10 +100,14 @@ def run_spec3_ifushort_extract1d(rtdata_module):
@pytest.mark.parametrize(
'output',
[
- 'jw01024-c1000_t002_miri_ch3-mediumlong_x1d.fits',
- 'jw01024-c1000_t002_miri_ch4-mediumlong_x1d.fits',
- 'jw01024-c1000_t002_miri_ch3-mediumlong_s3d.fits',
- 'jw01024-c1000_t002_miri_ch4-mediumlong_s3d.fits',
+ 'jw01024-c1000_t002_miri_ch3-long_x1d.fits',
+ 'jw01024-c1000_t002_miri_ch3-medium_x1d.fits',
+ 'jw01024-c1000_t002_miri_ch4-long_x1d.fits',
+ 'jw01024-c1000_t002_miri_ch4-medium_x1d.fits',
+ 'jw01024-c1000_t002_miri_ch3-long_s3d.fits',
+ 'jw01024-c1000_t002_miri_ch3-medium_s3d.fits',
+ 'jw01024-c1000_t002_miri_ch4-long_s3d.fits',
+ 'jw01024-c1000_t002_miri_ch4-medium_s3d.fits',
]
)
def test_spec3_ifulong(run_spec3_ifulong, fitsdiff_default_kwargs, output):
@@ -122,10 +126,14 @@ def test_spec3_ifulong(run_spec3_ifulong, fitsdiff_default_kwargs, output):
@pytest.mark.parametrize(
'output',
[
- 'jw01024-c1000_t002_miri_ch1-mediumlong_x1d.fits',
- 'jw01024-c1000_t002_miri_ch2-mediumlong_x1d.fits',
- 'jw01024-c1000_t002_miri_ch1-mediumlong_s3d.fits',
- 'jw01024-c1000_t002_miri_ch2-mediumlong_s3d.fits'
+ 'jw01024-c1000_t002_miri_ch1-long_x1d.fits',
+ 'jw01024-c1000_t002_miri_ch1-medium_x1d.fits',
+ 'jw01024-c1000_t002_miri_ch2-long_x1d.fits',
+ 'jw01024-c1000_t002_miri_ch2-medium_x1d.fits',
+ 'jw01024-c1000_t002_miri_ch1-long_s3d.fits',
+ 'jw01024-c1000_t002_miri_ch1-medium_s3d.fits',
+ 'jw01024-c1000_t002_miri_ch2-long_s3d.fits',
+ 'jw01024-c1000_t002_miri_ch2-medium_s3d.fits'
],
)
def test_spec3_ifushort(run_spec3_ifushort, fitsdiff_default_kwargs, output):
@@ -145,10 +153,14 @@ def test_spec3_ifushort(run_spec3_ifushort, fitsdiff_default_kwargs, output):
@pytest.mark.parametrize(
'output',
[
- 'miri_mrs_emsm_ch1-mediumlong_x1d.fits',
- 'miri_mrs_emsm_ch2-mediumlong_x1d.fits',
- 'miri_mrs_emsm_ch1-mediumlong_s3d.fits',
- 'miri_mrs_emsm_ch2-mediumlong_s3d.fits'
+ 'miri_mrs_emsm_ch1-long_x1d.fits',
+ 'miri_mrs_emsm_ch1-medium_x1d.fits',
+ 'miri_mrs_emsm_ch2-long_x1d.fits',
+ 'miri_mrs_emsm_ch2-medium_x1d.fits',
+ 'miri_mrs_emsm_ch1-long_s3d.fits',
+ 'miri_mrs_emsm_ch1-medium_s3d.fits',
+ 'miri_mrs_emsm_ch2-long_s3d.fits',
+ 'miri_mrs_emsm_ch2-medium_s3d.fits'
],
)
def test_spec3_ifushort_emsm(run_spec3_ifushort_emsm, fitsdiff_default_kwargs, output):
@@ -168,8 +180,10 @@ def test_spec3_ifushort_emsm(run_spec3_ifushort_emsm, fitsdiff_default_kwargs, o
@pytest.mark.parametrize(
'output',
[
- 'jw01024-c1000_t002_extract1dtest_miri_ch1-mediumlong_x1d.fits',
- 'jw01024-c1000_t002_extract1dtest_miri_ch2-mediumlong_x1d.fits'
+ 'jw01024-c1000_t002_extract1dtest_miri_ch1-long_x1d.fits',
+ 'jw01024-c1000_t002_extract1dtest_miri_ch1-medium_x1d.fits',
+ 'jw01024-c1000_t002_extract1dtest_miri_ch2-long_x1d.fits',
+ 'jw01024-c1000_t002_extract1dtest_miri_ch2-medium_x1d.fits'
],
)
def test_spec3_ifushort_extract1d(run_spec3_ifushort_extract1d, fitsdiff_default_kwargs, output):
diff --git a/jwst/scripts/tests/test_scripts.py b/jwst/scripts/tests/test_scripts.py
index b9ba4e8d08..eb133b15dc 100644
--- a/jwst/scripts/tests/test_scripts.py
+++ b/jwst/scripts/tests/test_scripts.py
@@ -9,7 +9,6 @@
'asn_make_pool',
'collect_pipeline_cfgs',
'create_data',
- 'okify_regtests',
'pointing_summary',
'schemadoc',
'set_telescope_pointing',
diff --git a/pyproject.toml b/pyproject.toml
index 6d61c69b5c..ef138c0999 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -74,7 +74,6 @@ collect_pipeline_cfgs = "jwst.scripts.collect_pipeline_cfgs:main"
create_data = "jwst.scripts.create_data:main"
csvconvert = "jwst.csv_tools.csvconvert:CSVConvertScript"
exp_to_source = "jwst.exp_to_source.main:Main"
-okify_regtests = "jwst.scripts.okify_regtests:main"
pointing_summary = "jwst.scripts.pointing_summary:main"
schemadoc = "jwst.scripts.schemadoc:main"
set_telescope_pointing = "jwst.scripts.set_telescope_pointing:main"
@@ -106,8 +105,6 @@ sdp = [
]
test = [
"ci-watson>=0.5.0",
- "colorama>=0.4.1",
- "readchar>=3.0",
"ruff",
"pysiaf>=0.13.0",
"pytest>=6.0.0",