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 10 forks source link

Pixel scale in pyKLIP JWST.py #171

Open kammerje opened 4 months ago

kammerje commented 4 months ago

Unfortunately, another update of pyKLIP JWST.py is needed. Currently, the pixel scale is taken as the square root of the PIXAR_A2 header keyword, but this header keyword does not exist if the photometry step of the JWST pipeline has been skipped. Skipping this step is common practice for AMI and KPI data. Hence, doing KLIP RDI with KPI data is currently not possible. Tagging @semaphoreP and @AarynnCarter.

JarronL commented 4 months ago

Recommend grabbing XSciScale and YSciScale attributes from pysiaf aperture object if header keyword doesn't exist. The equivalent would then be np.sqrt(XSciScale * YSciScale).

kammerje commented 4 months ago

Agreed! I have implemented it already and can push it. It makes pysiaf a pyKLIP dependency though @semaphoreP.

semaphoreP commented 4 months ago

Isn't this just grabbing some header keyword? I just want to understand why we need a whole package for this, unless this is quite complicated.

kammerje commented 4 months ago

@semaphoreP Unfortunately not! Pixel scales are not in the headers by default. The pixel area is added when running the JWST pipeline photometry step, but if this step is skipped, then the current implementation will not work.

I have updated pyKLIP so that it only tries to import pysiaf when the photometry step has been skipped, so that the package is not required by default. See https://bitbucket.org/pyKLIP/pyklip/commits/0be44828dda8e7c22c465409d6a2decf18750eb6.

semaphoreP commented 4 months ago

Can this be a step in spaceKLIP to import pysiaf and insert the relevant keywords into the header? I'd like to both avoid forcing all pyklip uses to install jwst specific packages, but I also rather avoid not installing dependencies when possible.

JarronL commented 4 months ago

Jason's recommendation could work, but we should be careful about implementation. Jens' solution also looks feasible as it assumes the only reason you would be running this portion of the code is for JWST-specific tasks, so the user would already have pysiaf installed as it's a spaceKLIP dependency. As long as you put the import inside the JWST-specific function call, then there should not be any conflict with other portions of the code being used for other purposes. One potential hiccup might be for automated tests and coverage (but pysiaf could be a requirement in a separate requirement-tests.txt file).

semaphoreP commented 4 months ago

Both for testing and coverage, but also to keep JWST-specific things in spaceKLIP, it would be nice to do this in spaceKLIP. We could just have a step that checks for the platescale, does nothing if it's already there, or uses pysiaf to add it into the header. I'd like to keep the JWST interface in pyklip as simple as possible, so we don't have to remember if some pre-processing step is done in spaceKLIP or pyKLIP (it should always be done in spaceKLIP).

JarronL commented 4 months ago

Yeah, I agree. I had earlier advocated for moving much of the JWST stuff from pyklip to spaceklip, but it's been a while since I've looked at the pyklip code to see what makes sense.

kammerje commented 3 months ago

Ok, I have implemented a utility function in spaceKLIP which populates the PIXAR_A2 SCI header keyword if it is not already available before running pyKLIP (here). I have also reverted pyKLIP back to the previous state which did not need pysiaf as a dependency. On the way, I have also added a fallback option for when SVO is unavailable to the JWST class in pyKLIP (here). I hope this concludes our work here @semaphoreP and the pyKLIP jwst branch can now be merged into main.

kammerje commented 3 months ago

For completeness, I wanted to note that with NIRCam LW for instance, the pixel area in arcsec^2 is 0.003918172030450555 when populated by the JWST pipeline photometry step, but 0.003896778872914225 (i.e., slightly different) when using the average of the X and Y pixel scales from pysiaf. Any ideas why this could be the case @mperrin?

mperrin commented 3 months ago

My first guess for what might cause that difference in PIXAR_S2 values was probably sightly different choice of apertures, or different versions... But upon looking a bit closer at the code in the pipeline, I see that the JWST pipeline photometry step is not getting values directly from pysiaf, but instead is retrieving reference files from CRDS for pixel area. Looking at the history comments of those, it seems like the difference is due to including the geometric distortions, and in particular the slight non-orthogonality of the X and Y axes. Example history text from https://jwst-crds.stsci.edu/browse/jwst_nircam_area_0042.fits:

The nominal pixel area that is stored in PIXAR_A2 and PIXAR_SR was calculated by multiplying the nominal pixel scale in the x and y directions by the sine of the angle between the x and y axes. This angle was calculated using the arctangent of the Science to Ideal coefficients in the y and x directions. This file was created using create_pam.py

Follow up question for @JarronL and @juliengirard, maybe you can confirm if the above interpretation is correct? It looks like the pixel area reference files for NRC_CORON were all derived from commissioning data, so I presume may eventually get updated with more precise values as the geometric distortion knowledge improves for coronagraphy modes?