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
tickets/DM-12404: Update stack with new deblender API #48
Conversation
While the new deblender was being tested a large number of changes were made to the user branch u/fred/deblender. The vast majority of those commits contain code that is no longer necessary or usable in master and it was thought that merging the entire branch, with dozens of outdated commits, is undesirable. Instead the cumulative changes made on the user branch are contained in this commit, which consists almost entirely of the work done in DM-11329, which was already review by pschella in the user branch.
@@ -32,6 +33,8 @@ | |||
import lsst.afw.detection as afwDet | |||
import lsst.afw.table as afwTable | |||
|
|||
logger = lsst.log.Log.getLogger("meas.deblender.deblend") |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Adding logging would have been nice to see as a separate commit.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
yes, but the module was added to add a log, which was included in the same commit.
minTranslation = pexConfig.Field(dtype=float, default=1e-8, | ||
|
||
# Blend Configuration options | ||
minTranslation = pexConfig.Field(dtype=float, default=1e-3, |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Why the enormous difference in value?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Test revealed that there was no need to waste CPU cycles on higher precision positions.
translationMethod = pexConfig.Field(dtype=str, default="default", | ||
doc=("Method to use for fitting translations." | ||
"Currently 'default' is the only available option," | ||
"which performs a linear fit, but it is possible that we" | ||
"will use galsim or some other method as a future option")) | ||
edgeFluxThresh = pexConfig.Field(dtype=float, default=1.0, | ||
doc=("Boxes are resized when the flux at an edge is " | ||
"> edgeFluxThresh * bg_rms")) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Sad mixing of variable naming conventions :(
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes, that is what happens when a code using a non PEP8 naming convention uses a 3rd party code that plays by the rules. I updated the document to explicitly state "background RMS."
doc=("Calculate exact Lipschitz constant in every step" | ||
"(exact_lipschitz is True) or only calculate the" | ||
"Lipschitz constant with significant changes in A,S" | ||
"(exact_lipschitz is False)")) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Just say True
or False
and leave out the exact_lipschitz
part (which is confusing because of the name change).
"(exact_lipschitz is False)")) | ||
stepSlack = pexConfig.Field(dtype=float, default=0.2, | ||
doc=("Slack to use when updating the step size in each iteration of" | ||
"the bSDMM minimization")) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
In what units? Pixels?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I updated the description to be more clear.
"This is used to initialize the model for each source" | ||
"and to set the size of the bounding box for each source" | ||
"every `refinementSkip` iteration.")) | ||
usePsfConvolution = pexConfig.Field(dtype=bool, default=False, doc=("Perform PSF kernel matching?")) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Convolve images to a common PSF?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I clarified this, and made it True
by default, since this will be necessary for both LSST and HSC and is now fast enough to be run in production.
# If the default constraints are not used, set the constraints for | ||
# all of the sources | ||
constraints = None | ||
if (sorted(self.config.constraints) != ['+', '1', 'M', 'S'] |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ah, so it should be a list above? Or at least specify that ordering doesn't matter.
maxIter=self.config.maxIter, | ||
bg_scale=self.config.bgScale, | ||
relErr=self.config.relativeError | ||
) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Why the strange mix between snake case and camel case? and why abbreviate relativeError
?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This was an artifact of using a 3rd party package that adheres to PEP8 naming conventions (rel_err
is the parameter name in proxmin
and scarlet
. I fixed both names in meas_deblender
.
1<<maskPlane["CR"] | | ||
1<<maskPlane["NO_DATA"] | | ||
1<<maskPlane["SAT"] | | ||
1<<maskPlane["SUSPECT"]) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The user may want to specify the mask planes to flag.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Oh, and are you hoping the mask is the same for all masked images here?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The mask doesn't have to be the same in all bands for this to work, but the mask plane does. I do like your idea to let the user specify the keys to use for bad pixels in the mask plane and I think it's safe to assume that each image will have the same mask plane, that will be implemented on the next commit. @TallJimbo , please shout if you think there is a chance that a set of coadds might have different mask planes when passed to the deblender and I'll update the code.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
You should write this as:
badPixels = maskedImages[0].mask.getPlaneBitMask(self.config.badMask)
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think it's safe to assume that coadds from different bands will have the same mask planes, but adding an assert
to at least check that they have the same number of mask planes would probably be wise.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ok, just to be safe I calculate the mask for each maskedImage independently:
if badMask is None:
badMask = ["BAD", "CR", "NO_DATA", "SAT", "SUSPECT"]
fpMask = afwImage.Mask(bbox)
debResult.footprint.spans.setMask(fpMask, 1)
fpMask = ~fpMask.getArray().astype(bool)
mask = np.zeros(weights.shape, dtype=bool)
for fidx, mimg in enumerate(maskedImages):
badPixels = mimg.mask.getPlaneBitMask(badMask)
mask[fidx] = (mimg.getMask().array & badPixels) | fpMask
weights[mask] = 0
cx = peak.x+xmin | ||
cy = peak.y+ymin | ||
cx = src.center[1]+xmin | ||
cy = src.center[0]+ymin |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
(y, x)? Just checking ;)
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
yes, unfortunately I begrudgingly agreed to (y,x) pixel coordinates
No description provided.