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.

Sign up for GitHub

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

Jump to bottom

DM-5532 #4

Merged
merged 8 commits into from
Mar 30, 2016
Merged

DM-5532 #4

Show file tree
Hide file tree
Changes from 6 commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
9bc95d9
Update unit test for new star selector API
r-owen Mar 18, 2016
f0f302c
Clean up pyflakes warnings in unit test
r-owen Mar 18, 2016
bc39e70
Rename catalog to sourceCat in selectStars
r-owen Mar 18, 2016
f7e74c1
Modernize psfEx star selector
r-owen Mar 18, 2016
0d7cb93
Add documentation to psfEx star selector
r-owen Mar 18, 2016
101c2ec
Use StarSelector.run to make PSF candidates
r-owen Mar 22, 2016
ae95302
PsfExStarSelector: fix flag names
r-owen Mar 29, 2016
e6ea9ab
Fix flake8 picky errors
r-owen Mar 29, 2016
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
196 changes: 108 additions & 88 deletions python/lsst/meas/extensions/psfex/psfexStarSelector.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -29,27 +29,16 @@
except ImportError:
plt = None

from lsst.afw.table import SourceCatalog
from lsst.pipe.base import Struct
import lsst.pex.config as pexConfig
import lsst.pex.logging as pexLogging
import lsst.afw.display.ds9 as ds9
import lsst.afw.math as afwMath
import lsst.meas.algorithms as measAlg
from lsst.meas.algorithms.starSelectorRegistry import starSelectorRegistry
from lsst.meas.algorithms import StarSelectorTask
import lsst.meas.extensions.psfex as psfex

class PsfexStarSelectorConfig(pexConfig.Config):
badFlags = pexConfig.ListField(
doc="List of flags which cause a source to be rejected as bad",
dtype=str,
default=["base_PixelFlags_flag_edge",
"base_PixelFlags_flag_saturated.center",
"base_PixelFlags_flag_cr.center",
"base_PixelFlags_flag_bad",
"base_PixelFlags_flag_suspect.center",
"base_PsfFlux_flag",
#"parent", # actually this is a test on deblend_nChild
],
)
__all__ = ["PsfexStarSelectorConfig", "PsfexStarSelectorTask"]

class PsfexStarSelectorConfig(StarSelectorTask.ConfigClass):
fluxName = pexConfig.Field(
dtype=str,
doc="Name of photometric flux key ",
Expand Down Expand Up @@ -98,16 +87,6 @@ class PsfexStarSelectorConfig(pexConfig.Config):
default=100,
check = lambda x: x >= 0.0,
)
kernelSize = pexConfig.Field(
dtype=int,
doc = "size of the Psf kernel to create",
default = 21,
)
borderWidth = pexConfig.Field(
doc = "number of pixels to ignore around the edge of PSF candidate postage stamps",
dtype = int,
default = 0,
)

def validate(self):
pexConfig.Config.validate(self)
Expand All @@ -121,6 +100,17 @@ def validate(self):
if self.minFwhm > self.maxFwhm:
raise pexConfig.FieldValidationError("minFwhm (%f) > maxFwhm (%f)" % (self.minFwhm, self.maxFwhm))

def setDefaults(self):
self.badFlags = [
"base_PixelFlags_flag_edge",
"base_PixelFlags_flag_saturated.center",
"base_PixelFlags_flag_cr.center",
"base_PixelFlags_flag_bad",
"base_PixelFlags_flag_suspect.center",
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm pretty sure these .center strings aren't correct anymore.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Great catch.

"base_PsfFlux_flag",
#"parent", # actually this is a test on deblend_nChild
]

class EventHandler(object):
"""A class to handle key strokes with matplotlib displays"""
def __init__(self, axes, xs, ys, x, y, frames=[0]):
Expand Down Expand Up @@ -193,27 +183,86 @@ def plot(mag, width, centers, clusterId, marker="o", markersize=2, markeredgewid

return fig

class PsfexStarSelector(object):
## \addtogroup LSST_task_documentation
## \{
## \page PsfexStarSelectorTask
## \ref PsfexStarSelectorTask_ "PsfexStarSelectorTask"
## \copybrief PsfexStarSelectorTask
## \}

class PsfexStarSelectorTask(StarSelectorTask):
"""!A star selector whose algorithm is not yet documented

@anchor PsfexStarSelectorTask_

@section meas_extensions_psfex_psfexStarSelectorStarSelector_Contents Contents

- @ref meas_extensions_psfex_psfexStarSelectorStarSelector_Purpose
- @ref meas_extensions_psfex_psfexStarSelectorStarSelector_Initialize
- @ref meas_extensions_psfex_psfexStarSelectorStarSelector_IO
- @ref meas_extensions_psfex_psfexStarSelectorStarSelector_Config
- @ref meas_extensions_psfex_psfexStarSelectorStarSelector_Debug

@section meas_extensions_psfex_psfexStarSelectorStarSelector_Purpose Description

A star selector whose algorithm is not yet documented

@section meas_extensions_psfex_psfexStarSelectorStarSelector_Initialize Task initialisation

@copydoc \_\_init\_\_

@section meas_extensions_psfex_psfexStarSelectorStarSelector_IO Invoking the Task

Like all star selectors, the main method is `run`.

@section meas_extensions_psfex_psfexStarSelectorStarSelector_Config Configuration parameters

See @ref PsfexStarSelectorConfig

@section meas_extensions_psfex_psfexStarSelectorStarSelector_Debug Debug variables

PsfexStarSelectorTask has a debug dictionary with the following keys:
<dl>
<dt>display
<dd>bool; if True display debug information
<dt>displayExposure
<dd>bool; if True display the exposure and spatial cells
<dt>plotFwhmHistogram
<dd>bool; if True plot histogram of FWHM
<dt>plotFlags
<dd>bool: if True plot the sources coloured by their flags
<dt>plotRejection
<dd>bool; if True plot why sources are rejected
</dl>

For example, put something like:
@code{.py}
import lsstDebug
def DebugInfo(name):
di = lsstDebug.getInfo(name) # N.b. lsstDebug.Info(name) would call us recursively
if name.endswith("objectSizeStarSelector"):
di.display = True
di.displayExposure = True
di.plotFwhmHistogram = True

return di

lsstDebug.Info = DebugInfo
@endcode
into your `debug.py` file and run your task with the `--debug` flag.
"""
ConfigClass = PsfexStarSelectorConfig
usesMatches = False # selectStars does not use its matches argument

def __init__(self, config):
"""Construct a star selector using psfex's algorithm
def selectStars(self, exposure, sourceCat, matches=None):
"""!Select stars from source catalog

@param[in] config: An instance of PsfexStarSelectorConfig
"""
self.config = config

def selectStars(self, exposure, catalog, matches=None):
"""Return a list of PSF candidates that represent likely stars
@param[in] exposure the exposure containing the sources
@param[in] sourceCat catalog of sources that may be stars (an lsst.afw.table.SourceCatalog)
@param[in] matches astrometric matches; ignored by this star selector

A list of PSF candidates may be used by a PSF fitter to construct a PSF.

@param[in] exposure: the exposure containing the sources
@param[in] catalog: a SourceCatalog containing sources that may be stars
@param[in] matches: astrometric matches; ignored by this star selector

@return psfCandidateList: a list of PSF candidates.
@return a Struct containing:
- starCat a subset of sourceCat containing the selected stars
"""
import lsstDebug
display = lsstDebug.Info(__name__).display
Expand All @@ -226,8 +275,6 @@ def selectStars(self, exposure, catalog, matches=None):
lsstDebug.Info(__name__).plotFlags # Plot the sources coloured by their flags
plotRejection = display and plt and \
lsstDebug.Info(__name__).plotRejection # Plot why sources are rejected
# create a log for my application
logger = pexLogging.Log(pexLogging.getDefaultLog(), "meas.extensions.psfex.psfexStarSelector")

#-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-
#
Expand All @@ -244,21 +291,21 @@ def selectStars(self, exposure, catalog, matches=None):
maxelong = (maxellip + 1.0)/(1.0 - maxellip) if maxellip < 1.0 else 100

# Unpack the catalogue
shape = catalog.getShapeDefinition()
ixx = catalog.get("%s.xx" % shape)
iyy = catalog.get("%s.yy" % shape)
shape = sourceCat.getShapeDefinition()
ixx = sourceCat.get("%s.xx" % shape)
iyy = sourceCat.get("%s.yy" % shape)

fwhm = 2*np.sqrt(2*np.log(2))*np.sqrt(0.5*(ixx + iyy))
elong = 0.5*(ixx - iyy)/(ixx + iyy)

flux = catalog.get(fluxName)
fluxErr = catalog.get(fluxErrName)
flux = sourceCat.get(fluxName)
fluxErr = sourceCat.get(fluxErrName)
sn = flux/np.where(fluxErr > 0, fluxErr, 1)
sn[fluxErr <= 0] = -psfex.psfex.cvar.BIG

flags = 0x0
for i, f in enumerate(self.config.badFlags):
flags = np.bitwise_or(flags, np.where(catalog.get(f), 1 << i, 0))
flags = np.bitwise_or(flags, np.where(sourceCat.get(f), 1 << i, 0))
#
# Estimate the acceptable range of source widths
#
Expand Down Expand Up @@ -327,7 +374,7 @@ def selectStars(self, exposure, catalog, matches=None):
ds9.mtv(mi, frame=frame, title="PSF candidates")

with ds9.Buffering():
for i, source in enumerate(catalog):
for i, source in enumerate(sourceCat):
if good[i]:
ctype = ds9.GREEN # star candidate
else:
Expand Down Expand Up @@ -370,7 +417,8 @@ def selectStars(self, exposure, catalog, matches=None):

if displayExposure:
global eventHandler
eventHandler = EventHandler(plt.axes(), imag, fwhm, catalog.getX(), catalog.getY(), frames=[frame])
eventHandler = EventHandler(plt.axes(), imag, fwhm, sourceCat.getX(), sourceCat.getY(),
frames=[frame])

if plotFlags or plotRejection:
while True:
Expand Down Expand Up @@ -399,39 +447,11 @@ def selectStars(self, exposure, catalog, matches=None):
else:
break

#
# Time to use that stellar classification to generate psfCandidateList
#
with ds9.Buffering():
psfCandidateList = []
if True:
catalog = [s for s,g in zip(catalog, good) if g]
else:
catalog = catalog[good]
starCat = SourceCatalog(sourceCat.schema)
for source, isGood in zip(sourceCat, good):
if isGood:
starCat.append(source)

for source in catalog:
try:
psfCandidate = measAlg.makePsfCandidate(source, exposure)
# The setXXX methods are class static, but it's convenient to call them on
# an instance as we don't know Exposure's pixel type
# (and hence psfCandidate's exact type)
if psfCandidate.getWidth() == 0:
psfCandidate.setBorderWidth(self.config.borderWidth)
psfCandidate.setWidth(self.config.kernelSize + 2*self.config.borderWidth)
psfCandidate.setHeight(self.config.kernelSize + 2*self.config.borderWidth)

im = psfCandidate.getMaskedImage().getImage()
vmax = afwMath.makeStatistics(im, afwMath.MAX).getValue()
if not np.isfinite(vmax):
continue
psfCandidateList.append(psfCandidate)

if display and displayExposure:
ds9.dot("o", source.getX() - mi.getX0(), source.getY() - mi.getY0(),
size=4, frame=frame, ctype=ds9.CYAN)
except Exception as err:
logger.logdebug("Failed to make a psfCandidate from source %d: %s" % (source.getId(), err))

return psfCandidateList

starSelectorRegistry.register("psfex", PsfexStarSelector)
return Struct(
starCat = starCat,
)
29 changes: 9 additions & 20 deletions tests/testPsfexPsf.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -31,16 +31,10 @@
python
>>> import testPsfexPsf; testPsfexPsf.run()
"""

import os, sys
from math import *
import math
import numpy as np
import unittest
import eups
import lsst.utils.tests as utilsTests
import lsst.pex.exceptions as pexExceptions
import lsst.pex.logging as logging
import lsst.pex.policy as pexPolicy
import lsst.afw.image as afwImage
import lsst.afw.coord as afwCoord
import lsst.afw.detection as afwDetection
Expand All @@ -49,14 +43,11 @@
import lsst.afw.table as afwTable
import lsst.afw.display.ds9 as ds9
import lsst.daf.base as dafBase
import lsst.afw.display.utils as displayUtils
import lsst.meas.algorithms as measAlg
import lsst.meas.algorithms.defects as defects
import lsst.meas.algorithms.utils as maUtils
import lsst.afw.cameraGeom as cameraGeom
import lsst.meas.extensions.psfex.psfexPsfDeterminer as psfexPsfDeterminer
# register the PSF determiner
import lsst.meas.extensions.psfex.psfexPsfDeterminer
assert lsst.meas.extensions.psfex.psfexPsfDeterminer # make pyflakes happy
from lsst.meas.base import SingleFrameMeasurementTask
from lsst.meas.algorithms.detection import SourceDetectionTask

try:
type(verbose)
Expand All @@ -77,8 +68,8 @@ def psfVal(ix, iy, x, y, sigma1, sigma2, b):
c, s = np.cos(theta), np.sin(theta)
u, v = c*dx - s*dy, s*dx + c*dy

return (exp(-0.5*(u**2 + (v*ab)**2)/sigma1**2) +
b*exp(-0.5*(u**2 + (v*ab)**2)/sigma2**2))/(1 + b)
return (math.exp(-0.5*(u**2 + (v*ab)**2)/sigma1**2) +
b*math.exp(-0.5*(u**2 + (v*ab)**2)/sigma2**2))/(1 + b)

class SpatialModelPsfTestCase(unittest.TestCase):
"""A test case for SpatialModelPsf"""
Expand All @@ -92,7 +83,6 @@ def measure(footprintSet, exposure):

measureSources = SingleFrameMeasurementTask(schema,config=config)

tab = afwTable.SourceTable.make(schema)
catalog = afwTable.SourceCatalog(schema)
if display:
ds9.mtv(exposure, title="Original", frame=0)
Expand Down Expand Up @@ -220,8 +210,7 @@ def tearDown(self):
@staticmethod
def setupDeterminer(exposure):
"""Setup the starSelector and psfDeterminer"""
starSelectorFactory = measAlg.starSelectorRegistry["objectSize"]
starSelectorConfig = starSelectorFactory.ConfigClass()
starSelectorConfig = measAlg.ObjectSizeStarSelectorTask.ConfigClass()

starSelectorConfig.sourceFluxField = "base_GaussianFlux_flux"
starSelectorConfig.badFlags = ["base_PixelFlags_flag_edge",
Expand All @@ -231,7 +220,7 @@ def setupDeterminer(exposure):
]
starSelectorConfig.widthStdAllowed = 0.5 # Set to match when the tolerance of the test was set

starSelector = starSelectorFactory(starSelectorConfig)
starSelector = measAlg.ObjectSizeStarSelectorTask(starSelectorConfig)

psfDeterminerFactory = measAlg.psfDeterminerRegistry["psfex"]
psfDeterminerConfig = psfDeterminerFactory.ConfigClass()
Expand Down Expand Up @@ -285,7 +274,7 @@ def testPsfexDeterminer(self):
starSelector, psfDeterminer = SpatialModelPsfTestCase.setupDeterminer(self.exposure)
metadata = dafBase.PropertyList()

psfCandidateList = starSelector.selectStars(self.exposure, self.catalog)
psfCandidateList = starSelector.run(self.exposure, self.catalog).psfCandidates
psf, cellSet = psfDeterminer.determinePsf(self.exposure, psfCandidateList, metadata)
self.exposure.setPsf(psf)

Expand Down