DM-33900: cp_pipe: calculate the gain using a pair of flats #125
Conversation
| @@ -237,6 +257,9 @@ def run(self, inputExp, inputDims): | |||
| inputDims : `list` | |||
| List of exposure IDs. | |||
|
|
|||
| taskMetadata : `list`[`lsst.pipe.base.TaskMetadata`] | |||
There was a problem hiding this comment.
I think a space is required between the back-tick and open-bracket.
| @@ -223,7 +243,7 @@ def runQuantum(self, butlerQC, inputRefs, outputRefs): | |||
| outputs = self.run(**inputs) | |||
| butlerQC.put(outputs, outputRefs) | |||
|
|
|||
| def run(self, inputExp, inputDims): | |||
| def run(self, inputExp, inputDims, taskMetadata): | |||
There was a problem hiding this comment.
Do the metadata items not need to be sorted the same as the exposures? I know the use case is for only a single flat pair to be supplied, but what happens if a full list is given? Is this not an issue because once sorted you can access everything via the exposure dimension?
|
|
||
| # Estimate the gain from the flat pair | ||
| # Overscan readnoise from post-ISR exposure metadata | ||
| readNoise = self.getReadNoiseFromMetadata(taskMetadata, ampName) |
There was a problem hiding this comment.
I see this just averages all the read noise values together (which is probably fine for PTC curves taken in sequence), which probably invalidates my previous question on sorting. That said, this is being done for each amplifier/exposure time pair. I think this needs to be pulled out of the exposure time loop, so it's done only once per amplifier.
| @@ -492,8 +526,10 @@ def measureMeanVarCov(self, exposure1, exposure2, region=None): | |||
| ---------- | |||
| exposure1 : `lsst.afw.image.exposure.ExposureF` | |||
| First exposure of flat field pair. | |||
|
|
|||
There was a problem hiding this comment.
Here and below: the extra lines between parameters isn't needed. See: https://developer.lsst.io/python/numpydoc.html#parameters
| 1/g = <(I1 - I2)^2/(I1 + I2)> - 1/mu(sigma^2 - 1/2g^2) | ||
|
|
||
| See below for the full solution. | ||
| https://www.wolframalpha.com/input/?i=solve+1%2Fy+%3D+c+-+(1%2Fm)*(s^2+-+1%2F(2y^2))+for+y |
There was a problem hiding this comment.
I think we should just insert the solution (it's just a simple quadratic).
| denom = (2*const*mu - 2*readNoise**2) | ||
|
|
||
| positiveSolution = (root + mu)/denom | ||
| negativeSolution = (mu - root)/denom # noqa: F841 unused, but the other solution |
There was a problem hiding this comment.
I think this can be dropped, left as a comment to avoid calculation, or moved to the docstring where we have the quadratic solution explained.
|
|
||
| readNoises = [] | ||
| for expMetadata in taskMetadata: | ||
| overscanNoise = expMetadata['isr'][expectedKey] |
There was a problem hiding this comment.
This probably needs to check that the key exists. My concern is appending a None to the list, and having a TypeError raised.
There was a problem hiding this comment.
I'll check for the key isr and if it doesn't exist, I'll continue the loop. If in the end the median readout noise can't be calculated, I'll return None (with a warining), and then in the function that estimates the gain, if the read noise is None but the correction type is not NONE, instead of raising (as it is right now), I'll issue a warning and default to correction type None, which does not require the noise, to at least return an estimation.
tests/test_ptc.py
Outdated
| for ampName in self.ampNames: | ||
| if exposurePair.gain[ampName] is np.nan: | ||
| continue | ||
| self.assertAlmostEqual(exposurePair.gain[ampName], inputGain, delta=0.1) |
There was a problem hiding this comment.
Is that large of a tolerance needed? I would think that the simulation would get better agreement.
There was a problem hiding this comment.
I see differences of ~0.065, especially with no correction.
Fix docstrings in other files
No description provided.