Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

MIRI Subarrays 390 Hz EMI #7729

Closed
stscijgbot-jp opened this issue Jul 14, 2023 · 168 comments · Fixed by #7857
Closed

MIRI Subarrays 390 Hz EMI #7729

stscijgbot-jp opened this issue Jul 14, 2023 · 168 comments · Fixed by #7857

Comments

@stscijgbot-jp
Copy link
Collaborator

Issue JP-3248 was created on JIRA by Michael Regan:

The majority of the MIRI subarrays have an 390 Hz EMI noise pattern in the raw data. This 390 Hz maps exactly to 256 pixels and has a constant amplitude of +- 4 DN. The effect of this is to imprint this into the rate images. For short integrations in LRSSLITLESS the correlated noise from this is quite apparent in the rate images. For longer integrations the net effect is to increase the read noise by ~20%.

The long term plan is a change to the sizes and locations of the subarrays to get the frame times to be in phase with the 390 Hz like the full frame images. For the previous and near term observations this can be fixed by adding a new step to the pipeline. A prototype has been implemented and shown to remove the 390 Hz to a fraction of a DN.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

For some context, I've posted a short powerpoint presentation given to the ramps-to-slopes WG last fall with some info on MIRI 10Hz and 390Hz noise. It has been determined that the period of the 390Hz waveform is precisely 256.0 pixel times, and a number of different subarray size/location/frametime combinations close to the existing subarray dimensions have been derived that will eliminate the 390Hz noise operationally.   

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Added a few more figures showing 390Hz noise in several of the existing MIRI subarrays, with some extracted waveforms and in/out of phase examples.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Misty Cracraft on JIRA:

Anton Koekemoer  Should this be presented at the next CALWG meeting? It should only affect MIRI, but it is a new addition to the pipeline being suggested.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Anton Koekemoer on JIRA:

Thanks very much Misty - indeed, I've been in discussion with with Mike and Eddie since they suggested this to me as a topic. (I also now added myself and also Greg Sloan as watchers on this ticket).

It might perhaps be efficiently scheduled as an 'out-of-cycle' CalWG mtg, ie not necessarily involving the monthly full CalWG since it's pretty MIRI-specific. though of course we'd advertise it to the full CalWG as always. Here are the dates I'd proposed to Mike and Eddie at the time when we first discussed this:

.. in terms of 'off-cycle' slots, the next one would be Jun 27 (since Jun 20 is actually taken already).

Then in July we have a full Cal WG meeting on Jul 11, and off-cycle slots available on Jul 18 and 25.

And finally in Aug we have a full Cal WG meeting on Aug 1, then with off-cycle slots available for any of the remaining Tuesdays until the next full Cal WG meeting on Sep 5, etc. 

So as soon as we hear from  Michael Regan and Eddie Bergeron about which date works best for them then we can schedule it for discussion.

 

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Anton Koekemoer on JIRA:

Mike emailed me a few weeks after the above and this has now been scheduled for discussion at JWST Cal WG Meeting 2023-07-18

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Anton Koekemoer on JIRA:

This was discussed at JWST Cal WG Meeting 2023-07-18 and approved as a candidate for implementation (following prioritization by the MIRI Instrument team as indicated in the criticality field for the ticket), please continue further discussion about implementation details/ questions in the comments here in this ticket.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Howard Bushouse on JIRA:

Eddie Bergeron can you provide SCSB with a copy of your IDL code that implements this correction?

Also, it was mentioned during the Cal WG discussion that Eddie's code could either derive a reference waveform from the exposure being processed or have it provided in a reference file. It seemed like the consensus for pipeline implementation was to use the reference file approach. Should we, however, keep the ability to generate a waveform from the data in the code that we implement, just in case it's useful in the future (or as an option for "power" users)?

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Anton Koekemoer on JIRA:

Howard Bushouse  wrote:

Also, it was mentioned during the Cal WG discussion that Eddie's code could either derive a reference waveform from the exposure being processed or have it provided in a reference file. It seemed like the consensus for pipeline implementation was to use the reference file approach. 

Yes, the pipeline implementation would use the reference file approach.

Should we, however, keep the ability to generate a waveform from the data in the code that we implement, just in case it's useful in the future (or as an option for "power" users)?

From Eddie's presentation it sounded as though this may be implemented as two separate pieces, ie one to generate the waveform (either from the current data or some other set of data) which could remain as a stand-alone script if not in the pipeline, while the other would be the pipeline modification to remove the waveform from the data (regardless of whether it's being supplied as a reference file or generated on-the-fly from the same dataset).

So this question should be discussed with more input from the MIRI team (tagging Sarah Kendrew and the other MIRI people on this ticket), with a few considerations being:

  • is this issue urgent enough that the team would like to first add this correction to the pipeline to remove the waveform from the data, assuming that it's provided in the form of a reference file? (in which case the waveform reference files could continue to be generated using the current IDL script that Eddie already has in hand).
  • what would be the cost difference in effort (ie time) from the SCSB side for the two scenarios, ie just implementing the pipeline correction to remove the waveform, vs the additional option of coding this in such a way that it's also able to generate it from the data?

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Anton Koekemoer on JIRA:

hi Howard Bushouse just checking whether you've received the IDL code from Eddie Bergeron  and Michael Regan  yet - if not, can you touch base with them please and let us know here when they will be planning to send it to you?

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Howard Bushouse on JIRA:

Anton Koekemoer No, have not received the code yet from either Eddie Bergeron or Michael Regan . Pinging here and via non-Jira channels.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Sorry for the delay, I'm working on it. Hopefully my cleanup and additional comments will make it easier for you to implement. Will aim to post it here by the end of the day tomorrow...ready or not.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Almost there. Just need a few more tests after the cleanup to make sure I didn't break anything. Will have it done today.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Ok, here's the main idl routine. There are some function calls. I'll post the other files, some test datasets and a little more description tomorrow. There's a brief algorithm overview at the top of the routine.

[^fix_miri_emi.pro]

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Maria Pena-Guerrero on JIRA:

Started to work on this. The step name I am going with for now is "emicorr" is everyone OK with this? 

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Anton Koekemoer on JIRA:

thanks, fine with me

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Maria Pena-Guerrero on JIRA:

what are all the EXP_TYPE values this step will be applied to?

@stscijgbot-jp
Copy link
Collaborator Author

stscijgbot-jp commented Sep 11, 2023

Comment by Howard Bushouse on JIRA:

The description mentions LRS slitless as having this problem, and I'm assuming Imaging mode would too (it has the most number of subarrays allowed), and perhaps subarrays used for coronagraphy as well?  So that would include MIR_IMAGE, MIR_LRS-SLITLESS, MIR_LYOT, and MIR_4QPM.  Eddie Bergeron Michael Regan can you confirm?

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Michael Regan on JIRA:

In general, this could run for any MIRI detector in any mode. 390 Hz is only in subarrays but there is a significant effect at 10 Hz that is in all read out modes.

 

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Ok, I've made a few mods to the fix_miri_emi.pro routine, fixing some bugs and adding a few features. The algorithm itself is unchanged. I've bundled up the main routine with all of the called subroutines and some example reference wave fits files into the attached tar file:

[^fix_miri_emi_lib.tar]

"tar xvf fix_miri_emi_lib.tar" will create a subdir called fix_miri_emi_lib with all the necessary routines in it. There are three main tools, and the rest are just dependencies. The example reference wave fits files are also in this subdir. The main tools to be concerned with are:

[^fix_miri_emi.pro]   ; this is it! (I've deleted the older version that was posted in this ticket).

[^miri_fft.pro]    ;  A utility tool for making a noise power spectrum of any miri exposure

[^miri_frametime.pro]    ; A tool for calculating the miri frametime of any fast or slow fullframe or subarray, based on Mike Ressler's frame timing doc

 

No need to code the latter two for the pipeline, but they are useful for other things as shown in the example calls in fix_miri_emi.pro. In particular, fix_miri_emi.pro needs some basic timing numbers for each subarray like rowclocks and frameclocks, however I've hardcoded those values into the case statements at the top of the main routine. Those should not change for the existing miri full/subroutines, but if new subarrays are defined you can use miri_frametime.pro to calculate the new values (or use Mike's own python routine for doing that - see doc reference in miri_frametime.pro).

 

The easiest way to set things up in idl is to just add the fix_miri_emi_lib dir to your idl path after starting idl. If you want to run this code and have not used idl before, just follow the steps here. You can install idl on your machine using the self service portal, and you will get a network license as long as you are on the staff network or vpn. Start idl and set the path:

 

IDL> !path='/path/to/where/you/put/fix_miri_emi_lib:'+!path

 

Be sure to put this directory FIRST in your path. Idl uses precedence, and I've modified a version of fits_write.pro to correct an issue with the default version in astronlib that was making the output fits files incompatible with pipeline input. If you get errors trying to use the output in the pipeline, check to make sure you have this order correct, or at least make sure you are using my version of fits_write.pro and not the astron or other library versions.

 

There are a set of calling examples at the top of fix_miri_emi.pro. Give those a try. I've also included an example of how to create a reference wave file from any given dataset.

 

Changes from my previous version:

  • Fixed an indexing bug when running on single-int exposures
  • Removed the "default" reference wave file. You now have to specify either a desired frequency to phase (for a self-correction) or reference_wave_filename pointing at a fits file with the frequency in the header keyword "FREQ". The routine will exit with a message if neither is specified (If both a reference_wave_file and a frequency are specified, the reference wave takes precedence). User-specified frequencies should have very high precision. See examples in fix_miri_emi.pro.
  • Switched to a "per-integration with phase offset" calculation of the pixel phases to reduce the risk of overflowing the datatype on exposures with larger numbers of frames/ints. Gives exactly the same result as the old way. I left the old code in but commented out. You should still use the largest unsigned integer datatype you can for this (ul64 in this code).
  • Added a correction for the phase of the last reference pixel in a row (if there is one). Didn't affect the current subarrays since none of those reach the last refpix, but for generalized use with fullframes (e.g. 10 hz) this is needed. (Although its a tiny correction - leaving it out wouldn't matter, but to be strictly correct...)
  • Added optional scaling of the reference wave to match the amplitude of the self-wave (pa="phased amplitude"). I've noticed variations in the amplitudes of the phased waveforms in some datasets. Some of this is expected, such as with the 10Hz noise, but some is not. This scaling option seems very robust. I recommend that this be used by default with reference waves. It is ignored for self-correction of course.
  • Changed the ref-to-self vector cross correlation routine from brute-force to the classic ffta*conj(fftb) method. They give identical results. I left the old code in there, but commented out. Originally added a parabolic fit for sub-pixel accuracy but that's an unnecessary complication so took it back out. Use whatever python equivalent you want. 
  • Added an optional param to limit the number of integrations used to phase the self-wave. This can be useful for large datasets, since you can get a self-wave with sufficient S/N with just a few frames/ints and save processing time. By default this nints_to_phase starts with the first int in the exposure. The correction is still applied to all integrations in the exposure regardless of how many are used to phase.
  • Added a ",/check_phase" option to overplot the phase-matched reference wave and self-wave (with optional scaling) for your viewing pleasure :)  

 

Give it a try and let me know if you run into problems. Anyone else, feel free to try the idl routine on your favorite _uncal.fits  datasets and push them through the pipeline. In a following post to this ticket I'll include a few example plots and discussion about how to identify frequencies for correction and how to get them at the necessary precision to stay in phase across large datasets.

 

Regarding the reference files. As currently written, this code reads individual reference waves from separate fits files as an amplitude vector (only, evenly spaced in phase from 0-1.0) with the header keyword FREQ specifying the frequency. Since different full/subarrays will have different correction needs - only some subarrays need 390Hz correction, while all full/sub can benefit from a 10hz correction, and some subarrays have multiple noise components - I'd suggest a single reference file format, with either a table of variable length records for the wave vectors, or a multi-extension fits with an extension for each wave vector, and then keywords with each record/extension indicating which subarrays it should be applied to. Then the EMICORR step can recursively go though and apply each appropriate correction. New waves can be added by simply adding new records or extensions to the single reference file. If you feel it necessary, the phase at each amplitude could also be specified, but so far this seems to work with enough precision just assuming equal phase spacing in the reference wave vector. There is a reference_wave phase per bin calculation in the code ("pa_phase") but it is currently unused.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Maria Pena-Guerrero on JIRA:

thank you, I'll start working on this.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Maria Pena-Guerrero on JIRA:

I am having trouble locating Mike's python code to calculate the time frames, can you please add it in the ticket too, please?

@stscijgbot-jp
Copy link
Collaborator Author

stscijgbot-jp commented Sep 28, 2023

Comment by Eddie Bergeron on JIRA:

The reference is: DFM478_fps_timing_fpga2_initial.pdf (MIRI DFM 478 - 04.02 "MIRI FPS Exposure Time Calculations", (SCE FPGA2), 28 October 2014, M. Ressler, MIRI Project Scientist, JPL. Here is the doc, python code in the text. I believe the equivalent numbers you need for "rowclocks" and "frameclocks" in fix_miri_emi.pro from Mike's code are roi_row_clocks and frame_clocks, but check these against what miri_frametime.pro returns. Should be the same, but if different let me know. My idl version is not a translation of his python, it was developed directly from the prescription in the doc.

Note that I've hardcoded some of the current flight constants in miri_frametime.pro that are required params in calls to Mike's code, so look for those params in miri_frametime.pro. Technically some of these could change in the future, but that's unlikely.

edit: Maria, I see that this doc is marked itar. I'll e-mail you with the location of this doc in the itar area on central store.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Those other "flight constant" params, like the main subarray definitions themselves (colstart, colstop, rowstart, rowstop), should be in a database somewhere. 

@stscijgbot-jp
Copy link
Collaborator Author

stscijgbot-jp commented Oct 13, 2023

Comment by Maria Pena-Guerrero on JIRA:

I created the following dictionaries based on the docstring in the IDL script, to correct for the frequencies according to subarray. Eddie Bergeron and/or Michael Regan, can you please that the dictionary of frequencies to correct for according to subarray is correct or how should it be changed? Please look at the image frees2correct.png: 

!freqs2correct.png|thumbnail!

Thank you. 

 

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Maria Pena-Guerrero on JIRA:

I have completed the implementation of the step. We now need to request CRDS to add a rule to look for a reference file for the emicorr step. For now, the resulting value can be set to 'N/A'. Howard Bushouse should I file the Jira request or should MIRI? 

@stscijgbot-jp
Copy link
Collaborator Author

stscijgbot-jp commented Oct 26, 2023

Comment by Maria Pena-Guerrero on JIRA:

Another question Howard Bushouse, Michael Regan, Eddie Bergeron, is the plan to have 1 reference file with all frequencies to correct for OR 1 file per frequency?

(I am coding the file per frequency approach)

And, are we allowing the user to supply a reference file? (if so, as the code is right now, they would have to submit a file per frequency)

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Howard Bushouse on JIRA:

For the reference file(s) I would think that you would at least bundle things to the level of having one ref file per subarray, which would contain the list of freqs to correct for that type of subarray. Then the ref file would use SUBARRAY as a selector in CRDS. Or you could even bundle everything up into a single ref file that contains a table (or equivalent) listing freqs for each subarray type. The cal step would then have to find the matching subarray entry in the ref file and read the list of associated freqs. The only downside to that approach is that if you want to change the list for one subarray you have to redeliver the whole ref file. Oh, and I would encourage the use of a ref file in ASDF format, rather than the clunky FITS BINTABLEs.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Maria Pena-Guerrero on JIRA:

Hmmm... Howard Bushouse, in thinking about it a little bit more, from the instrument's team it would be painful to do any changes to the frequencies because then all the files would have to be modified. So I think it makes more sense to have all frequencies in one file and choose from them according to subarray. This would make it easier for users to create their own reference file as well. I am using the fits format, should I switch to ASDF and abandon fits entirely?

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Howard Bushouse on JIRA:

Yes, for this type of info, which is a list or table of entries, as opposed to something like an image, I would strongly advocate for a ref file in ASDF format. So much easier to work with than FITS BINTABLEs. There are examples of similar ASDF ref files already in CRDS.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Yes, the 10Hz "heater" noise is EMI pickup from the detector heater/temp-sense circuit and it affects all MIRI data, not just the subarrays. There had always been a plan to try to remove this with a pipeline correction. This general purpose emicorr step takes care of it, just like it does for the FPGA-induced 390 Hz in affected subarrays. It just has a different physical source and characteristic shape in the data. Two for one.

In full frame the 10hz noise looks a little like 1/f noise, but its periodic rather than a random walk. I looked at the emicorr.fits result that Misty just made and compared it to both the uncorrected uncal and the idl-corrected version. Looks good to me. Here's an animated gif of the group 45-44 CDS, with and without correction. Hard stretch here (+/- 20 DN), so the signal areas are all white.

!emicorr_full_cds_compare.gif!

@stscijgbot-jp
Copy link
Collaborator Author

stscijgbot-jp commented Dec 19, 2023

Comment by Misty Cracraft on JIRA:

Though my tests just showed that a SLOWR1 file was not corrected. This is part of my shower testing notebook, run in a local temp directory. Eddie Bergeron  should the 10Hz apply to slowr1 mode?

2023-12-19 11:36:56,016 - stpipe.Detector1Pipeline.emicorr - INFO - Using user-supplied reference file: /Users/cracraft/MIRI/temp/emicorr/jwst_miri_emicorr_0002.asdf

2023-12-19 11:36:56,050 - stpipe.Detector1Pipeline.emicorr - INFO - Using reference file to get subarray case.

2023-12-19 11:36:56,051 - stpipe.Detector1Pipeline.emicorr - INFO - With configuration: Subarray=FULL, Read_pattern=SLOWR1, Detector=MIRIMAGE

2023-12-19 11:36:56,051 - stpipe.Detector1Pipeline.emicorr - WARNING - No correction match for this configuration

2023-12-19 11:36:56,052 - stpipe.Detector1Pipeline.emicorr - WARNING - Step skipped

2023-12-19 11:36:56,056 - stpipe.Detector1Pipeline.emicorr - INFO - Step emicorr done

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Rosa Diaz on JIRA:

I was going to ask the same question. The reference file does not seem to cover SLOWR1 and FASTR1, only FAST and SLOW. If we need to cover those cases too, via the reference file, then we need to modify the file again.

@stscijgbot-jp
Copy link
Collaborator Author

stscijgbot-jp commented Dec 19, 2023

Comment by Maria Pena-Guerrero on JIRA:

The reference file is ok. There was another bug in the code that it could not find that configuration. This was fixed. Misty Cracraft please run your tests on this branch (and I hope this time we caught them all):

pip install git+https://github.com/penaguerrero/jwst@emicorrbugfix2

Rosa Diaz Mike said that the same correction applies for either FAST or FASTR1 and the same for SLOW and SLOWR1.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Yes, the 10Hz correction should be applied to ALL MIRI data: FAST*, SLOW*, all subarrays and all fullframes. It's in every MIRI exposure, including all ground test data that used the flight-like thermal control. The 10Hz noise amplitude (1-1.5 DN PtoP) is much lower than the 390Hz noise (9 DN PtoP) though.

@stscijgbot-jp
Copy link
Collaborator Author

stscijgbot-jp commented Dec 19, 2023

Comment by Rosa Diaz on JIRA:

Yes, I understood that, but not sure if the repeated info or mapping (FAST array = FASTR1 array, SLOW array = SLOWR1 array) needed to be in the reference file (duplicate entries) or was in the code. I assume this is in the code then.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Misty Cracraft on JIRA:

Latest test shows it running on SLOWR1 data.

 

2023-12-19 14:16:37,370 - stpipe.Detector1Pipeline.emicorr - INFO - Using reference file to get subarray case.

2023-12-19 14:16:37,372 - stpipe.Detector1Pipeline.emicorr - INFO - With configuration: Subarray=FULL, Read_pattern=SLOWR1, Detector=MIRIMAGE

2023-12-19 14:16:37,409 - stpipe.Detector1Pipeline.emicorr - INFO - Will correct data for the following 1 frequencies:

2023-12-19 14:16:37,409 - stpipe.Detector1Pipeline.emicorr - INFO - ['Hz10_slow_MIRIMAGE']

2023-12-19 14:16:37,410 - stpipe.Detector1Pipeline.emicorr - INFO - Correcting for frequency: 10.039216 Hz (1 out of 1)

2023-12-19 14:16:37,419 - stpipe.Detector1Pipeline.emicorr - INFO - Subtracting self-superbias from each group of each integration 2023-12-19 14:19:03,889 - stpipe.Detector1Pipeline.emicorr - INFO - Phase calculation per integration

2023-12-19 14:19:04,034 - stpipe.Detector1Pipeline.emicorr - INFO - Calculating the phase amplitude for 500 bins

2023-12-19 14:19:15,753 - stpipe.Detector1Pipeline.emicorr - INFO - Using reference file to measure phase shift

2023-12-19 14:19:15,807 - stpipe.Detector1Pipeline.emicorr - INFO - Creating phased-matched noise model to subtract from data

2023-12-19 14:19:15,936 - stpipe.Detector1Pipeline.emicorr - INFO - Subtracting EMI noise from data

2023-12-19 14:20:00,448 - stpipe.Detector1Pipeline.emicorr - INFO - Saved model in /ifs/jwst/wit/miri/cracraft/redcat_tests/jw01335002001_02101_00001_mirimage_emicorr.fits

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Yep, this looks good. Harder to see the difference visually in slow because it wraps to nearly even-odd striping and the scattered background is higher in this set. I checked it against the idl correction using the Hz10_slow_MIRIMAGE_reference wave and it matched nicely. Looks like the correct reference wave was selected. The 10Hz reference waves I provided need some work, but the code itself seems to be working as expected.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Maria Pena-Guerrero on JIRA:

ok, good to hear about the testing results! Misty Cracraft can you please let us know when you have run all your tests on this branch?

Regarding the limit in the number of integrations used to find the correction to limit the running time, I tested with a dataset with integrations=100, groups=100, y=224, x=288. This is what I found:

 
||integrations limit||390Hz (128 bins)||10 Hz (500 bins)||
|100|18 min|21 min|
|50|11 min|16 min|
|25|6 min|8 min|
|20|5 min|7 min|
|15|4 min|4 min|

Based on this tests, Michael Regan and Eddie Bergeron, what would be a good number of integrations to find the phase for integrations of about 100 and greater than 1000?

 

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Misty Cracraft on JIRA:

At this point, the code seems to be operating on all the modes as it should. Eddie has confirmed that the output data are looking as they should. My full notebook won't finish today (it tests a lot of data), but I'm not seeing any data skipping this step now (like I was earlier), so I think my part of the testing on this is probably about done.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

The important number for the top end limit is number of cycles that are being phased together to make the PA wave. You need enough S/N in that phased wave to get a good phase match to the reference wave. If it's too noisy that match could fail, but the phasing is a powerful technique - it's robust against all sorts of residual errors in in the processed input dataset such as the signal subtraction, CRs, bad pix etc.

So rather than looking only at the number of integrations, start by looking at the product of integrations and groups. That is, integrations*groups = total groups. You can't do partial integrations, but if we say we want a conservative minimum number of groups to phase we can then just round up to the next higher number of integrations. Any dataset with fewer total groups than this limit would just use all available integrations in the exp.

Using the MASK1550 example here, your 4 minute timings for the 15 integration case is for 15100= 1500 groups. The number of cycles of each frequency in this many groups is (frameclocksnints*ngroups_per_int)/period_in_clocks:

(23968L15100.)/256.0 = 140438 cycles of 390 Hz

(23968L15100.)/9960.94 = 3609.3 cycles of 10Hz

These are way more than enough. In this case, 1 integration of 100 groups would probably be enough - 9300 and 240 cycles each respectively.  The worst-case scenario would be for 10Hz in the smallest subarray, SUB64.

(8512L1100.)/256.0 = 3325 cycles of 390 Hz in 100 groups

(8512L1100.)/9960.94 = 85.5 cycles of 10 Hz in 100 groups

This is also probably enough, but on the low side for 10Hz. If we are conservative and set the limit at 4x this (min 400 total groups), that would probably be enough. This doesn't scale linearly...I worry about short integrations where measuring the signal rate might be problematic, but I just don't have a good feel for that. Since this is an upper limit we can be as conservative as we like, so it's just a matter of how much the extra execution time costs.

Maria, can you run a few more timing tests with this set, for 1, 5 and 10 integrations (100, 500 and 1000 total groups)? The trend seems to be flattening already. If there isn't much execution time difference between 400 groups and 1000 groups we could just take the higher number.

 

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

If you wanted to be more sophisticated and squeeze the execution efficiency you could choose an Ncycles limit and then just calculate how many integrations you need to get that many cycles for a given frequency (rounding up to the nearest Nintegrations). This would naturally scale with frame size, since a full frame will hold more cycles than a subarray. If we just set a single Ngroups limit, fullframes will take longer than subarray and be over-constrained. However if the difference in execution time is less than a factor of 2 this probably doesn't matter.

The right way to do all of this is either a S/N metric on the PA, or a "goodness of phase match" metric relative to the reference wave. The current code has neither of those (at least my idl code doesn't), but it does have a parameter for limiting the number of integrations to process, which is why we are talking about an Nintegrations upper limit.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Misty Cracraft on JIRA:

Latest testing shows the code working on all modes and readout patterns. It also is faster than yesterday (Howard asked me to re-test after an update last night). Eddie has confirmed that the results look good, so I'm saying this step passes MIRI testing. If Mike, Eddie or anyone else chooses to do any further testing, they can comment on those tests here, but for now, I'm checking this off of the MIRI pipeline testing tasks.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Greg Sloan on JIRA:

I've been watching the work on the EMI correction fly by the past few days on this ticket. Congratulations to all on the hard work and success!

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Maria Pena-Guerrero on JIRA:

further testing showed that subarray MASK1065 is not defined neither in the reference file nor in the default cases. I implemented code to avoid a crash and simply skip emicorr if the subarray is not defined, but wanted to let you know in case you want to add it to the defaults and/or the reference file.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Code to skip emicorr for unrecognized/undefined subarrays is a good thing in any case. However, MASK1065 should indeed be in the defaults and emicorr reference file. MASK1065 data does need 390Hz (and 10Hz) correction. Maria, what exactly is meant by "default cases" vs the reference file here? I'm unfamiliar.

Misty, is there a MASK1065 example in your test set?

@stscijgbot-jp
Copy link
Collaborator Author

stscijgbot-jp commented Dec 21, 2023

Comment by Maria Pena-Guerrero on JIRA:

the default cases is a dictionary in the pipeline code in case CRDS throws back N/A, then the code autocorrects following that dictionary.

note that this entry needs to be added to the reference file.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Rosa Diaz on JIRA:

The test cases that seem to be crashing are handled by DMS Mike Swam or Richard Spencer  should be able to provide with an example.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Maria Pena-Guerrero on JIRA:

I think this was my bad, I skipped the line while creating the subarray_cases. Here is the updated reference file with all the info Misty had requested to include: [^jwst_miri_emicorr_0003.asdf]

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Meaghan McDonald on JIRA:

The new file passes certification. I will deliver to CRDS TEST now.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Meaghan McDonald on JIRA:

The new emicorr file has been delivered to CRDS TEST. The context to now use for testing is jwst_1182.pmap.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Misty Cracraft on JIRA:

I just ran a test using the latest dev version of build 10.1 and crds-test and all looks good for the FASTR1 and SLOWR1 data in my test set. I do have one FAST mode coronagraph file (MASK1065), and when I run it through emicorr, it does appropriately apply the 390 and 10 Hz corrections.

2024-01-31 17:35:04,350 - stpipe.Detector1Pipeline.emicorr - INFO - With configuration: Subarray=MASK1065, Read_pattern=FAST, Detector=MIRIMAGE

2024-01-31 17:35:04,371 - stpipe.Detector1Pipeline.emicorr - INFO - Will correct data for the following 2 frequencies:

2024-01-31 17:35:04,372 - stpipe.Detector1Pipeline.emicorr - INFO - ['Hz390', 'Hz10']

But when I then run the shower code on it, I get large regions of NaNs, which is a problem that the EMICORR step corrected for all of the FASTR1 coron data. I suppose this does not appropriately treat FAST data.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Misty Cracraft on JIRA:

Eddie Bergeron Do you know what the emi patterns are for FAST mode? Is this something we should look into? The TA data is taken in FAST and FASTGRPAVG, but this may need to be a discussion in the MIRI team as to what we should do with that.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Mike had mentioned to me that only TA used FAST or FASTGRPAVG and so those were not considered for the emicorr. My understanding is that the difference between FAST and FASTR1 is that FASTR1 has one full reset frame between integrations and FAST does not. Is that correct? If so, it would be possible to adjust the pixel phase calculation to work for FAST exposures. It's a simple enough change.

FASTGRPAVG is another beast. The out-of-phase 390hz will be averaged into a different pattern in the averaged frames. This is easily calculable, but I'm not sure what effect that would have on the phase finding. It should work, but it's not as simple as the fix for FAST vs FASTR1.

So the desire here is to emi-correct the TA exposures in the archive for science use?

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Misty Cracraft on JIRA:

I've brought this up as an agenda item for next week's MIRI meeting. I'll relay your information and bring back the decision of the team as to what we want to do with TA data.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Eddie Bergeron on JIRA:

Ok, great. Again, the fix for FAST vs FASTR1 should be very simple - a one-line conditional statement in Maria's code. A fix for FASTGRPAVG would be more involved. Keep in mind that the group averaging will "beat down" the emi amplitude to some extent, so the correction is smaller to begin with. I suspect the 390 might remain somewhat coherent through the averaging, but this will depend on the subarray dimensions. Also note there would be no need to change the emicorr reference file for these fixes. The existing reference waves are applicable as-is. It would just be a python code change.

The decision depends on the goals to some extent. If the goal is to correct the data for science post-obs, then sure. But if the goal is to use the data to assess the TA itself maybe you'd want to leave the 390 in there since that's what the onboard TA was working with. Seems unlikely that the 390 would affect the TA, but if it was, the solution would be to change the onboard TA subarray dimensions to be in-phase for 390.

 

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Misty Cracraft on JIRA:

Ok, the MIRI team talked about this question today. We would like an emi correction for FAST mode, but we are not requesting one for any FASTGRPAVG data. Can we do that with this ticket, or should we file another ticket specifically asking for that?

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Howard Bushouse on JIRA:

Misty Cracraft Given that the baseline implementation is already in place for B10.1, I would suggest a new ticket for enhancements/additions to the original.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by Misty Cracraft on JIRA:

Ok. I opened a new ticket for the FAST mode updates. I have tested the code for FASTR1 and SLOWR1 and consider the step to work as initially requested. This ticket can be closed while updates happen in the new ticket.

@stscijgbot-jp
Copy link
Collaborator Author

Comment by David Law on JIRA:

Closing, with a link to https://jira.stsci.edu/browse/JP-3533 for ongoing work.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging a pull request may close this issue.

2 participants