DM-54948: Switch source_injection to use empirical variance#71
Conversation
09281b5 to
33f39b7
Compare
taranu
left a comment
There was a problem hiding this comment.
Looks fine, just a couple of questions to consider.
| mask_plane_dict = mask.getMaskPlaneDict() | ||
| if bad_mask_names is None: | ||
| bad_mask_names = [] | ||
| bad_mask_names = [name for name in bad_mask_names if name in mask_plane_dict] |
There was a problem hiding this comment.
This could be in an else: block. More importantly, though, should the user be warned if they try to use a mask that doesn't exist in the visits they're injecting in?
There was a problem hiding this comment.
There is a standard list of "bad" mask names that are used by default in the standard SSI workflow. My hope that non-specialist users should never have to tweak this or even think much about it, much as a typical user does not need to know the nuts and bolts of ISR. Logging a missing mask warning would have users thinking that they should go back and fiddle with this list.
As far as the actual impact goes: The standard bad SSI mask list is generically present in every major data release so far, and presumably into the future, though we have seen that the Roman-DESC sims are an exception (which motivated this bug fix).
There was a problem hiding this comment.
Logging a missing mask warning would have users thinking that they should go back and fiddle with this list.
Shouldn't they, though, if it's a mask plane that doesn't exist in any of the inputs? I would hope for a warning if I overrode it and made a typo, for example.
I grant it might be slightly annoying if we, for example, add the relatively new SPIKE mask to the list and someone runs injection on an older set of inputs without that mask, but unfortunately there's no real way to tell whether a missing mask is a user error or intentional.
There was a problem hiding this comment.
I've added a logger to the new infer_gain_from_image function which reports when the user is configuring the task with a missing mask plane, and updated the unit test to check that this function is behaving as expected.
| try: | ||
| amplifiers = exposure.getDetector().getAmplifiers() | ||
| bboxes = [amplifier.getBBox() for amplifier in amplifiers] | ||
| except AttributeError: |
There was a problem hiding this comment.
Hmm... is this the best way to test if this is a single CCD? It should work for the foreseeable future, I suppose.
There was a problem hiding this comment.
I think the most important thing for this code is that the exposure has a .getAmplifiers() method to get a list of amplifiers, each of which has a defined BBox. I'm not aware of any plans to change the handling of CCD exposures to not have these, but if this is a possibility then it's certainly worth thinking about.
If a CCD exposure were to be supplied without a .getAmplifiers() method, that currently puts us in a situation of considerable ignorance. I'm guessing the best we can do is to treat the entire exposure as a single BBox, with a single overall effective gain. The remaining question is whether we should model the variance as Poisson or Gaussian. I could go either way, but I'm inclined to use the Gaussian distribution to reflect the fact that a CCD without defined amplifier regions is outside our known scope.
With these prescriptions, the code works as intended, although is_single_CCD would be False. I'm inclined to leave that name as-is, unless you think this potential case is a serious concern.
There was a problem hiding this comment.
Ok, as I can't really come up with a likely scenario where this would happen, I guess it can stay as is.
There was a problem hiding this comment.
Not quite the same thing, but I have had to guard against exposure.detector.getAmplifiers() returning None (it happened in a unittest that used mocks): https://github.com/lsst/pipe_tasks/blob/main/python/lsst/pipe/tasks/calibrateImage.py#L2013-L2014
There was a problem hiding this comment.
Thanks @laurenam, that's useful input.
To guard against possible None amplifier info, I've modified the bboxes constructor to use the full bbox if getAmplifiers() returns None.
Additionally, it no-longer made sense to calculate is_single_CCD within the gain map construction function, so instead I determine whether or not this is a coadd-level dataset by checking the contents of info.getCoaddInputs(). If there's a more optimal way to do this with just the exposure object loaded, let me know.
This commit guards against getAmplifiers() returning None, as may be the case for mocks. Here we also remove the is_single_CCD check from the get_gain_map function and instead compute is_coadd in the primary injection function by checking the content of getCoaddInputs().
No description provided.