Photometry with partially-saturated sources

[stscijgbot-jp]

Issue JP-2256 was created on JIRA by Mattia Libralato:

I have found an issue with stage-2 and stage-3 images with a partially-saturated source.

I simulated one source per exposure. The source properties were chosen to reach the same level of saturation in the central pixels when I "observe" the source with the 20 groups.

For each of the four filters analyzed (F560W, F770W, F1130W, F2100W), I simulated a 4-point dither pattern with 3 setups (always one integration):

• 10 groups. The brightest pixel of the star never saturates on the ramp. These data are called "NOT_SAT";
• 15 groups. Saturation occurs at group 11 and there are 4 groups in the "plateau". These data are called "MILD_PART_SAT";
• 20 groups. Saturation occurs at group 11 and there are 9 groups in the "plateau". These data are called "IMP_PART_SAT".

In "ramp.png" you can find an example of the setup.

In each filter, I performed aperture photometry on all stage-2 images, stage-3 images obtained with one stage-2 image at the time, and the stage-3 mosaic image. Aperture photometry was performed with aperture radii from 0 to 100 pixels with steps of 0.1 pixel. Sky background was evaluated in an annulus with inner and outer radii of 100 and 120 pixels. I then compared the magnitude of the source in all cases. The comparison between the magnitudes in the stage-3 mosaic image and each stage-2 image is shown in "zp.png".

Each black point corresponds to the magnitude difference with a given aperture radius. Red points are the average values of all points. The x-axis values correspond to the filters (1 = F560W, 2 = F770W, 3 = F1130W, 4 = F2100W). The large scatter is mainly due to photometry obtained with small aperture radii. With small radii (<10-15 pixels), the exact definition of the center matters in the estimate of the flux, especially at short wavelengths. In the "NOT SAT" case, the magnitudes of both sets of images are in agreement on average. However, when some groups along the ramp start to saturate, the magnitude difference between the two sets of images start to increase. Furthermore, the magnitude zero-point I see is filter dependent, and it is more important at short wavelengths (likely because of the sampling).

I attached a notebook with many more plots that can highlight the issue I just described.

It has to be something in the calwebb_image3 resample step. This is how I run this part:

fwhm = 2.5 # Filter-dependent value taken from Jdox
minobj = 1
snr = 100 # Filter-dependent value
sigma = 3
fit_geom = 'shift'
use2dhist = False

im3 = Image3Pipeline()    
im3.tweakreg.kernel_fwhm = fwhm
im3.tweakreg.snr_threshold = snr
im3.tweakreg.minobj = minobj
im3.tweakreg.sigma = sigma
im3.tweakreg.fitgeometry = fit_geom
im3.tweakreg.use2dhist = use2dhist
im3.source_catalog.kernel_fwhm = fwhm
im3.source_catalog.snr_threshold = snr
im3.skymatch.save_results = True
im3.outlier_detection.save_results = True
im3.resample.save_results = True
im3.resample.kernel = 'square'
im3.source_catalog.save_results = True
im3.save_results = True```
I tested this issue with the MIRI imager, I do not know if this is instrument-dependent or not.

 

[jdavies-st]

I will add here that if you're trying to test saturation and are using simulated data, then turn tweakreg off. It has known issues and will complicate your analysis. Simulated data uses the exact same distortion models as the pipeline, so there will be nothing to be corrected by tweakreg.

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

For completeness, the issue is present even when I compare a stage-3 image obtained with only one stage-2 image and the corresponding stage-2 image.

[stscijgbot-jp]

Comment by Alicia Canipe on JIRA:

Setting this to "calwg_medium" priority to get it on the sprint list, we can bump it up higher if other instruments find it's an issue for them as well. May need more investigating. Misty Cracraft FYI

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

I re-processed all files with the latest pipeline version (1.3.1) and re-uploaded the data on Box (still uploading the new version; it will be v5). For the tweakreg step, I used the 'shift' geom (since there is only 1 source, only a shift can be computed). I slightly improved the notebook as well.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Mattia Libralato when you say "I performed aperture photometry on all stage-2 images, stage-3 images" does that mean you did photometry on both the unrectified "cal" and rectified "i2d" products, or were the measurements only done on i2d products? What I'm getting at is whether this is a discrepancy between before and after resampling, or just some kind of dependency on the number of images combined during resampling.

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Sorry for the late reply (I was on leave).

Howard Bushouse : I did aperture photometry on both unrectified and rectified products (both the single "i2d" images obtained with only 1 "cal" image and on the combined "i2d" mosaic). The discrepancies between unrectified and rectified products are different if I compare the "cal" and single "i2d" files or "cal" and mosaic "i2d" files. Although on average the differences are similar, the shorter the wavelength, the more different is the discrepancy between the two "cal" vs. "i2d" sets compared. The plots in the notebook can better describe the issue. So, I would say that sampling is a factor, but there might be something else.

James Davies [X] : I re-run the pipeline with "im3.tweakreg.skip = True". I'm uploading the new files on box. The result is still the same.

I just checked two random images and I made a correction in the "Description" above. The partial saturation does not occur just in one pixel, but in more than one. Furthermore, I noticed that some pixels have a DQ flag of 4 (jump detected).

**

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Investigation so far has revealed the following:

1) For the tests I've run so far, the tweakreg step appears to have no effect on the results at all, because the step correctly recognizes that it has too few matches between the images to compute a solution and hence the step gets skipped.

2) For the F770W exposures with no saturation, I have been able to verify - using just quick-n-dirty box-like photometry - that the measured flux in the individual cal products is nearly an exact match to the flux measured from each corresponding individual i2d product (i.e. the i2d products created by calwebb_image2). But the flux measured from the combined i2d image created by calwebb_image3 is in fact ~4% lower than any of the individual cal and i2d images. It appears this is due to a few pixels near the core of the single source getting flagged as outliers by the outlier_detection step, which then causes them to be excluded from the combined i2d image. I did an experiment of changing the outlier_detection parameter "scale" from its default value of "0.5 0.4" to "1.0 0.8" (i.e. just doubled the thresholds) and that resulted in no pixels near the source getting flagged anymore and the measured flux from the resulting combined i2d image matched the individual cal and i2d images to ~0.03%. A vast improvement compared to 4%!

3) For the F770W exposures with high levels of saturation, both the individual i2d images and the combined i2d image have measured fluxes ~4% lower than the individual cal images, even when I use the same higher values for the outlier_detection "scale" parameter. The individual and combined i2d images both indicate lower WHT values for the pixels near the core of the source that went into saturation in the original ramp. This is expected, because the smaller number of good ramp groups that contributed to those pixels in turn leads to higher VAR_RNOISE values than surrounding pixels (which used all groups) and the weighting in the resample step is now based on inverse readnoise. But it's still a mystery (at least to me) how this could affect the results for the individual i2d images, where each output pixel only has 1 input pixel contributing to it and hence the weights should be irrelevant. So now it's time to investigate how the weights are applied when there's only 1 image being resampled.

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Hi Howard Bushouse, thanks for the update! I'm relieved that you found these discrepancies as well. About your point 2, it is indeed a nice improvement!

About your point 3, what you have found is indeed puzzling. I hope the solution for the individual i2d images will be able to improve the combined i2d images as well.

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Howard Bushouse : I have a question about how the resampling is done. When you have one image, the value of a pixel in the i2d image is a combination of the flux of multiple close-by pixels from the cal.fits image (because of the distortion). Is it correct? If so, does the combination include a weight somewhere? Maybe this could the reason of the 4% discrepancy even with only one image.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Yes, there's definitely something "fishy" with the way the inverse variance weights are being applied. In addition to that, I've also found a bug in the way the specification for the "weight_type" parameter is being handled. As an experiment, I tried running resample with "weight_type=exptime" (instead of "ivm"), to see what difference it would make and found that my user input for the parameter was being ignored and the code was still applying ivm weighting. Once fixing that, I have now verified that when I resample a single image, using "weight_type=exptime" so that all pixels have the same weight, the fluxes/photometry in the resample image come out correctly, i.e. they match the values in the individual cal images. So now we need to figure out why the "ivm" weighting is not being applied correctly.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Adding documentation here for Mihai Cara who is going to start looking at the lower-level drizzle routines in an effort to determine whether the weights are being applied correctly. First, background information for him:

1) Simulated MIRI images with a single bright source at various exposure levels. The first exposure does NOT have the source go into saturation. The second one here has a longer exposure time that allows 3-4 pixels in the core of the source to go into saturation for some number of groups in the original ramps. Those saturated groups get flagged and rejected from ramp fitting, so that good slopes are still determined. But because those pixels use fewer groups than their neighbors, and because the variance due to read noise depends on the number of groups used, the pixels that went into saturation have higher values of VAR_RNOISE than their neighbors (which did not saturate).

2) The resample step uses the inverse of the VAR_RNOISE to do IVM weighting. Hence those pixels that saturated get lower weight than their neighbors, because their VAR_RNOISE values are higher. This appears to result in the resampled values for those pixels being systematically lower than the same pixels in the simulation that did not go into saturation. This happens even for the case of resampling a single image, such that each output pixel has only 1 image contributing to it, but of course can have contributions from fractions of multiple pixels. The guess is that the output flux is being biased lower due to the fact that the input pixels with higher flux have lower weights than the neighboring pixels that have lower flux but higher weights.

3) Example images follow. First, snapshots of the source in the exposure that did NOT saturate, showing the resampled SCI, WHT, and VAR_RNOISE arrays. The blue circle in the WHT snapshot shows the central pixels of the source, which in this exposure did NOT saturate and hence the weights are all relatively constant within the core of the source. Second, snapshots of the source in the exposure that DID saturate, showing the resampled SCI, WHT, and VAR_RNOISE arrays. The blue circle in the WHT snapshot indicates the location of the central few pixels that went into saturation and hence have higher VAR_RNOISE values than their neighbors and, equivalently, lower WHT values than their neighbors.

!SCI no saturation.png|thumbnail! !WHT no saturation.png|thumbnail! !RNOISE no saturation.png|thumbnail!

!SCI saturation.png|thumbnail! !WHT saturation.png|thumbnail! !RNOISE saturation.png|thumbnail!

[stscijgbot-jp]

Comment by Karl Gordon on JIRA:

Hum. An unintended consequence of the use of VAR_RNOISE for the weighting. Michael Regan Thoughts?

[stscijgbot-jp]

Comment by James Davies [X] on JIRA:

No, I think we currently think that it's actually a bug in the C drizzle code in how it applies weights, independent of what type of weights they are. The math here is really simple, but the C code in drizzle is not.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Mihai Cara example data that I used for my experimentation is located in █████████████████████████████████████████████████████████████████████████████ which is the directory containing the images with moderately heavy saturation (and hence many groups excluded for pixels at the core of a bright star).

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

So, I did my own investigation and I can confirm this issue and I think it is to be expected and I did not find anything wrong with the "drizzle" itself.

Howard Bushouse wrote: ??This is expected, because the smaller number of good ramp groups that contributed to those pixels in turn leads to higher VAR_RNOISE values than surrounding pixels (which used all groups) and the weighting in the resample step is now based on inverse readnoise. But it's still a mystery (at least to me) how this could affect the results for the individual i2d images, where each output pixel only has 1 input pixel contributing to it and hence the weights should be irrelevant. So now it's time to investigate how the weights are applied when there's only 1 image being resampled.??

Your first sentence is exactly correct. However, it still applies to single drizzled images because multiple pixels from the input image contribute to a single output pixel. Your second sentence ("each output pixel only has 1 input pixel") is not correct in general unless there is no distortion AND mapping between input and output coordinates is "identity" (x', y' = x, y) (i.e., no sub-pixel shifts) AND pixscale is 1. OR, if you use "point" kernel (pixfrac->0).

So, in general, each output pixel is a linear combination of several input pixels. The coefficients of the linear combination are determined by the overlap area of each input pixel with the output pixel. It is also combined with the weight assigned to each input pixel. So, if you weigh the brightest pixel very low relative to its lower intensity neighbors, I think it should not be surprising that you will get lower value in the output pixel and hence lower flux measurement.

For reference, the exact formula is provided in Casertano et al., AJ 2000, Appendix A.2 ("Drizzled Images") - see equation A1 and I also recommend Appendix A.3.

So, if you were to do flux computation on input image as sum(pixels*weights) and the same for the output aperture (image) then you will get a good agreement which shows that "drizzle" does exactly what it was designed to do. Maybe that is not what you want to do but then probably you should either switch to using "point" kernel or constant weighting (like EXPTIME) or modify the prescription for the computation of the IVM. However, as is, the code is working as intended.

I suggest that you repeat your experiments with a "point" kernel or EXPTIME weighting and see if that meets your needs.

Maybe there is an issue with how VAR_RNOISE is estimated? Maybe the noise in the saturated pixel is not 5x larger than that of its neighboring pixels? This is probably not a completely new issue and the same problem must have been encountered in processing NICMOS or WFC3 IR data? How was this dealt with for those instruments?

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Hi Mihai Cara, thanks for the update. 

If the "drizzle" step is performing as expected, I believe the issue is then more "philosophical" than practical. Many programs rely on using partially-saturated ramps.

I run the two tests, one using weight_type=EXPTIME and kernel="SQUARE", and one with weight_type="EXPTIME" and kernel="POINT". In the former case, I get the same issue. In the latter test, the result is slightly different. When I compare the single-i2d images with the cal images, there is a simple zero-point of 0.015-0.02 mag at all aperture radii in all filters but F2100W, for which the mag difference increases with the aperture radius. I think the issue with the F2100W images is related to the MIRISim simulations (a mismatch between reference-file versions that creates an extended gradient in the image). [Note that using the kernel="POINT" option with only one image is not a viable solution because the distortion leaves some pixels "empty".] Anyhow, the comparison between combined-i2d and cal images is always the same.

I am surprised that with no weight the result is still not correct. I have a question that might be related (apologies if otherwise). I run the image3 pipeline with weight_type="EXPTIME", but when I opened the resulting fits file with ds9, the header said:

NDRIZ = 4 / Number of groups/pointings included in resample
TEXPTIME= 222.0032 / [s] Total exposure time for product 
RESWHT = 'ivm ' / Type of drizzle weighting to use in resampling 
PIXFRAC = 1.0 / Drizzle parameter describing pixel shrinking 
PXSCLRT = 1.0 / Pixel scale ratio relative to native detector s```
There is "IVM" again there. I compared the combined-i2d image made with the weight_type="EXPTIME" and kernel="POINT" and the one made with weight_type="IVM" and kernel="POINT", and the 2 images are identical, which is weird. I tried also with kernel="SQUARE", same result. Could it be that somewhere in the pipeline the parameter passed as input is ignored. In the output log of the image3 pipeline, the parameters are those I gave in input.

For the last comment, this is something I'm curious as well. It would be nice to see if this happens with WFC3/IR data as well.

 

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Maybe these last findings are a different issue. I opened another tickect (https://jira.stsci.edu/browse/JP-2367)

 

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

If you are resampling only one input image with kernel='point' then likely there will be no difference between output images obtained with different weighting ('ivm' or 'exptime'). That is because, normally, no more than one input pixel can contribute to a single output pixel (but you could change pixel scales - i.e., make output pixels very large - to get more than one input pixel in each output pixel, if you wish).

[stscijgbot-jp]

Comment by Nadia Dencheva on JIRA:

Mattia Libralato The last communication here was at the ned of Dec. Just checking if you feel the issue was resolved or if there's something still left to look into. Thanks!

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Hi Nadia Dencheva, I installed the latest version of the pipeline and done what Mihai Cara suggested above (resample step with "point" kernel and EXPTIME weight). The result is still the same. I'm following another ticket (JP-1570) and I understand that also DQ flags have a weight in the process. So, maybe this last result is also expected. That said, if this is what is expected by the drizzle algorithm, I don't think there is much we can do. Maybe Howard Bushouse and Karl Gordon have something to add.

 

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

??... and done what Mihai Cara suggested above (resample step with "point" kernel and EXPTIME weight). The result is still the same.??

When you say "the result is still the same" what do you exactly mean? what files are you comparing? single i2d? combined i2d?

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

I'm comparing single cal fits files with both single and combined i2d.

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

I guess my question is when you say "The result is still the same", do you mean that combined i2d file does not change (=you get identical images) when you modify weight or kernel types in resample?

Individual i2d images should not be affected by weight (maybe a little but expect differences to be negligible) but they should change if kernel changes (a stronger factor).

If you see no change at all in combined i2d images when you modify weighting or kernel - that is likely a bug in the software. There is one more possibility: your 'cal' images do not have VAR_RNOISE extension and then IVM is equivalent to EXPTIME ... I believe. (because all noise will be set to 1)

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

"The result s still the same" means the issue with the photometry on partially saturated sources is always present. 

The combined i2d files change if I change weight and kernel.

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

I assume that "the issue with the photometry on partially saturated sources" is that aperture photometry on input images ('cal') gives higher fluxes than in the resampled image, right? If so, have you repeated aperture photometry while taking into account the weights of input pixels as suggested in my [#comment-564361] above?

??... do flux computation on input image as sum(pixels*weights) and the same for the output aperture (image)??

https://jira.stsci.edu/browse/JP-1570 is similar to this issue in many ways but it is also fundamentally different. There, the aperture photometry on input images was computed fundamentally in a different way from what was inputed into the resample step - it guaranteed that input will be different from the output. In both tickets I suggest that you compare "same kind" of aperture photometries: since resample weighs each pixel - do the same on input apertures. So, my notebook posted in https://jira.stsci.edu/browse/JP-1570 could be useful for this ticket too.

Unlike JP-1570, you do get good agreement when sources are not saturated, so your test setup is better than in JP-1570. You are getting differences only for saturated sources and that is because resample uses weights. You should do the same on your input data.

Setting weighting to EXPTIME should make the agreement better with aperture photometry on input images when the later does not use weighting. But it will never perfectly agree because of partial pixels as described in the Fruchter & Hook paper (so weighting is still better for comparison) unless you also set kernel to "point".

When weighting is EXPTIME, there is no way for the drizzle algorithm employed by the resample step to discriminate between saturated or not saturated pixels.

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Mihai Cara Yes, the issue I'm referring to is that aperture photometry on partially-saturated sources in a cal fits file gives different results than if using i2d single or combined fits files.

I will try to compute the aperture photometry as sum(pixels*weights). However, if I run the pipeline using kernel=point and weight=exptime, which weight should I use? I expect that the only difference should be related to the distorted vs. undistorted frame where I'm measuring the flux. However, I modified the notebook to use the pixel area map (PAM) on the cal fits files prior to measure the flux of the source.

I will run again the pipeline, just to be sure that it is using the correct kernel and weight.

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

??... if I run the pipeline using kernel=point and weight=exptime, which weight should I use???

If you set weight='exptime' then +when performing aperture photometry on input images+ you could simply do sum(pixels), i.e., set weights=1 in sum(pixels*weights) - this is assuming that all input images have identical exposure time. If weight='ivm' then weights should be set to 1/VAR_RNOISE (again, all of this is related to input image's apertures).

For the output apertures (in 'i2d' images) you should use the output weight of the resampled image (from the 'WHT' extension of the resampled image). This is needed, in principle, even when weight='exptime' because input pixels can contribute to multiple output pixels and are weighted based on overlap area. This issue is less important when kernel is set to "point".

The kernel controls how the flux from input images is dropped onto output grid. A point kernel will minimize chances of multiple input pixels contributing to output pixels. So, in that case, while I still suggest using output weights computed by the resample step when performing aperture photometry, setting weights=1 should be "good enough".

To answer your question directly: if weight='exprime' and kernel='point' then it is OK to not use weights +assuming+ your images do not have any "bad pixels". If some input pixels in some input images are flagged - this complicates things and then again I recommend using weights.

By the way, while PAM can be used to account for area change due to distortions, I do not believe it would be a major factor in aperture photometry, as confirmed by your own simulations for non-saturated sources - you seem to get good agreement regardless of PAM.

You have complex simulations with multiple factors that affect your images (saturation/rdnoise, IVM weighting, DQ arrays, sky background, etc.) You have to account for all factors in both input and output images.

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

Just a note. I am not saying that users should be using any weighting on resampled image or that they shouldn't be using weighting. I am just arguing, that for the testing purposes, if one wants to confirm that measurements on input images give same results as on the output images, weighting, in general, should be used. If your images have some pixels flagged in DQ - then you should exclude them from input images because that is what resample would do. Corresponding output pixels will be filled with values from other input images but will have lower weights. and so on... All of this has to be taken into account when validating resample step. Whether or not to use weights in the end should be a decision by individual researchers. That is, I am not commenting on how to do science and how to interpret weights and results. My comments should be are strictly applied to validation procedure.

And again, maybe IVM is not what you want or maybe you do not want to do aperture photometry on resampled images ... (???)

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Thanks for the explanations! I will use the weight when using exptime and kernel. I will let you know the results.

Yes, the PAM has a negligible effect because the star is at the center of the MIRI imager where the distortion is small. I would expect to see larger differences when placing the star in a corner.

Finally, this test is not meant as a validation of the pipeline. The ticket was filled in light of the unexpected results of some tests I was running in preparation of commissioning. If these results are what it is expected by the pipeline, then I'm good. However, I think that a brief summary of what has been discussed in this ticket should be added to the documentation just to help users understanding a bit more what is happening.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

Mattia Libralato Similar to the investigation for JP-1570, I think the work done here so far has shown no obvious flaw in the implementation of the resample step in the pipeline. As outlined in my comment in that other ticket (https://jira.stsci.edu/browse/JP-1570?focusedCommentId=572609&page=com.atlassian.jira.plugin.system.issuetabpanels%3Acomment-tabpanel#comment-572609) I'm wondering whether the investigation for this particular ticket can be declared complete, and we create a new ticket focused on updating or creating appropriate documentation to better explain some of the issues that've been discovered here.

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Howard Bushouse I was able to only run part of a test on resampled images obtained with weight='exptime' and kernel='point'. With the suggestions from Mihai Cara, the differences between the photometry on cal and i2d fits files is smaller than before. Differences are on average the order of <0.01 mag. The centermost 10-40 pixels (depends on the filter) show larger mag offsets, but the larger the aperture radius, the smaller the offset. So, it could be related both to sampling (in the case of F560W) and pipeline resampling. Note that because of the kernel some points in the i2d single images have "filler" values, so this might also explain the mag offsets at small radii. Finally, there are some weird mag differences (0.05 mag) sometimes, but they might be related to the simulations themselves (we know about artificial image background gradients). Bottom line: the results are better than before, although still not perfect.

I would like to run other tests using weight='ivm' and kernel='square' as when I opened the ticket, but I can't assure I will be able to run them on a short timescale.

I leave the choice of whether close the ticket or not to you.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

I've created JP-2560 to address potential updates to documentation for any of the issues discussed here and am closing this ticket (investigation of the data is complete).

[stscijgbot-jp]

Comment by Anton Koekemoer on JIRA:

Re-opening this ticket since the ticket was closed apparently without having discussed the proposed solution to possibly implement exptime weighting instead of the default ivm weighting (thereby, closing the ticket meant it disappeared from our lists for all the ticket sthat we're monitoring).

Tagging Karl Gordon, Mattia Libralato,  Howard Bushouse,  Mihai Cara,  Swara Ravindranath  so that we can determine the best path forward for resolving this issue.

 

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Hi Anton Koekemoer . I read the email from Tyler Pauly . Can you post here the figures as well?

I believe (and this is what came out from the discussion in this ticket) that the pipeline does what it is told to do. It is the weight the problem. In general, I would prefer to just combine the images without any weight because I just want to obtain a resampled image. Of course, the DQ flag should be considered so not to include in the computation pixels that should not be considered. If I'm not wrong, "weight=exptime" basically does this if all images have the same exp time.

I don't remember how exp time is defined. If it is the default NgroupsNintsNexp, then what I said before is OK. If it is the real exp time (so the actual number of groups used in the ramp fit, not the groups in the header), then you might have a similar issue as with "weight=ivm".

Maybe the question about the wrong weight/noise for bright sources should be redirected to other people.

[stscijgbot-jp]

Comment by Anton Koekemoer on JIRA:

hi Mattia Libralato  thanks very much for your reply. Indeed, currently exptime would be simply the value for the entire image, not per pixel taking into account the actual number of ramps (which indeed would likely led to the same issue as with ivm). So as long as you were able to confirm that 'exptime' gave you photometry that was closer to what you expected, then I interpret that as the recommendation, ie to use this as the weighting for partially saturated sources.

We can then consider in the CalWG about whether to keep the current default 'ivm' weighting or change to something different, eg exptime, but at least your tests are helpful in showing that this behaves better for partially saturated sources.

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

I would discuss this with the CalWG.

I have to check the results with my old simulations, too much time has passed. However, it seems that Tyler did more or less the same test I did, so now we have two independent tests, with two different instruments, showing that in some cases the ivm weight can be problematic.

I would not mix weights within the same pipeline call, either we use ivm or exptime for all pixels, not only those that are partially saturated (maybe I misinterpreted you statement; in this case apologies). One could think a way to weight these pixels knowing that they are "special" (knowing they are so from the DQ flag), but again it would be an ad-hoc tweak.

[stscijgbot-jp]

Comment by Eddie Schlafly on JIRA:

A potential hack-ish solution here that preserves the benefits of inverse variance weighting (downweighting high-noise pixels & pixels affected by CRs) while not biasing saturated sources would be to lie about the read noise variance in saturated pixels.  i.e., compute the read noise variances in these pixels as if the ramps were unaffected by saturation.

This clearly would give a poor prediction for the amount of read noise, but the read noise on these heavily source-dominated pixels is not terribly relevant, and the pixels have a saturated flag anyway.

One will still not be optimal on saturated sources—you'd be weighting by the inverse read noise, and you'd probably just want constant weights for these signal dominated pixels.  But at least you wouldn't be biased.

[stscijgbot-jp]

Comment by Howard Bushouse on JIRA:

I'm temporarily reassigning this back to Anton Koekemoer just to confirm that the notes on this issue from the Cal WG discussion on 2022-08-02 ([https://outerspace.stsci.edu/display/JWSTCC/2022-08-02+Meeting+Notes)] are the final word on this, which is:

"... the main issue was that the readthedocs pages discouraged the use of EXPTIME weighting, which confused the user, so the consensus was that the readthedocs pages should be updated to clarify that EXPTIME weighting is also a viable option in some cases, but that case should be taken in all cases to ensure accurate results."

So all you want is an update to the RTD docs for resample and that'll be it?  Also could you clarify for me the phrase "but that case should be taken in all cases to ensure accurate results"? Does "that case" refer to the case of using EXPTIME weighting and hence you're suggesting that EXPTIME should be used in all cases (i.e. effectively changing the default parameter setting)?

[stscijgbot-jp]

Comment by Misty Cracraft on JIRA:

Pinging the responsibles of this ticket. (Anton Koekemoer  and Mattia Libralato )  Is this now down to a documentation issue, or was there more that needed to be done by SCSB? Is anyone recommending changing the weighting, or do we just update the documentation and close the ticket?

[stscijgbot-jp]

Comment by Anton Koekemoer on JIRA:

hi Misty, thanks for pinging.

First, in answer to the question from Howard Bushouse   I think that "case" was a typo, and it should have read "care" instead, ie it should say this (I've updated the meetings page as well):

 ... but that care should be taken in all cases to ensure accurate results.

 

I'm checking the rest of the discussion now to see what else may still be needed.

[stscijgbot-jp]

Comment by Anton Koekemoer on JIRA:

Ok, as for how to resolve this ticket – there have been a range of pipeline changes since this ticket was first opened, some of which might now yield improved results, and we also now have real data to investigate this with.

So I'm going to pass this question back to Mattia Libralato  to ask whether he can have a go at re-running his code, but now with real data if it's similar enough to his original simulated data, and report back here on whether there still seems to be a problem in the pipeline, or whether it's now resolved, once he has looked at the results.

[stscijgbot-jp]

Comment by Anton Koekemoer on JIRA:

tagging Mattia Libralato  to see if the current IVM weighting (using the inverse of VAR_RNOISE) is still an issue for photometry on partially saturated sources – if it is, then Nadia Dencheva  just mentioned in  a meeting that the folks on Roman have been looking into this and may be able to provide some alternative schemes for weighting, which we could discuss in a future JWST CalWG meeting if it seems useful.

So  Mattia Libralato  could you let us know please, whether this is still an issue, if you'd like us to schedule a JWST CalWG discussion of it with the Roman folks?

[stscijgbot-jp]

Comment by Mattia Libralato on JIRA:

Anton Koekemoer  I still have to test this, so I do not have an answer at the moment.

[stscijgbot-jp]

Comment by Melanie Clarke on JIRA:

Anton Koekemoer - The NIRSpec team has noticed significant problems with using IVM weighting in the resample_spec step.  The data quality is noticeably worse, compared to exptime weighting.  We would be interested in exploring other weighting schemes.

[stscijgbot-jp]

Comment by Anton Koekemoer on JIRA:

thanks Melanie Clarke  this is helpful to know. Would you be able to provide a link to a page or other collection of the test cases used, and the results, so that SCSB can start taking a look at these results as well?

[stscijgbot-jp]

Comment by Melanie Clarke on JIRA:

We're still working on the analysis, as part of a calibration assessment for MOS mode.  I'll pull something together and update here.

[stscijgbot-jp]

Comment by Mihai Cara on JIRA:

Melanie Clarke Could you, please, see if installing this proposed fix (will replace current JWST pipeline):

pip install git+https://github.com/mcara/jwst.git@median-flt-ivm ```
makes any difference in data quality?

[stscijgbot-jp]

Comment by Melanie Clarke on JIRA:

Mihai Cara - sure, I'll check.

[stscijgbot-jp]

Comment by Melanie Clarke on JIRA:

Mihai Cara - I see no difference between ivm and ivm-smed weights on my sample data.  Both are significantly noisier than exptime weights.

Anton Koekemoer - I have placed some sample data and a notebook demonstrating the problem, including testing Mihai's changes, on central store, at:

████████████████████████████████████████████████████

The jwst pipeline version installed is from Mihai's fork, on the median-flt-ivm branch.

[stscijgbot-jp]

Comment by Paul Goudfrooij on JIRA:

Dear Melanie Clarke - I took a quick look at your notebook since I was interested in your findings, given that we found similar results for NIRISS data of bright stars earlier this year, which led to the development of a new step within the detector1 pipeline called "charge_migration" which is now ready for testing (see JP-3299). To check whether your results are relevant within that framework, I checked your notebook and it seems like you are extracting the uppermost source visible on the rate.fits files, correct?Under that assumption, I looked at the count rates in that spectrum as well as the number of groups and the frame read time and it looks like the pixels near the peak of the PSF and near the peak of the spectrum along the dispersion (where you see the biggest differences between the results with ivm and exptime weights) are saturating at some point during the ramp, correct? (Or at least reaching well beyond about 25,000 ADU during the integration.) If so, my (somewhat educated) guess is that this may well be due to effects of charge migration (I will check this further tomorrow when I have more time). You may well want to test the "charge_migration" step for this dataset to see whether that mitigates the difference between ivm- and exptime-weighted results of the "resample_spec" step (it does work very well for NIRISS imaging data).