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.

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

Simple utility continuum baseline checker #538

Merged
merged 24 commits into from
Nov 1, 2019

Conversation

brechmos-stsci
Copy link

This PR will create a simple spectrum checker to see if the flux baseline is near zero as this is required for several of the other methods (e.g., find_lines which assumes the flux data is continuum subtracted).

Fixes #535

Copy link
Member

@eteq eteq left a comment

Choose a reason for hiding this comment

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

Several inline suggestions.

A few broader comments:

  1. I think the docstring needs to describe the flow of what happens given various values of eps and whether or not the uncertainty is present (but maybe you were waiting on a first pass?)
  2. This function doesn't seem to be using the mask. Presumably the user does not want masked values to be included, so I'd say the thing to do is just s.flux[~s.mask] if a mask is present.
  3. There needs to be some way to turn the decorator "off", which the warning message should explain. I suggest using the configuration item machinery (see http://docs.astropy.org/en/stable/config/#adding-new-configuration-items) with just a boolean config item like always_check_continuum that defaults to True but just causes the decorator to immediately return True if it's set to False. Then you tell the user to use the with version of the config item machinery to supress the check.
  4. Once the modification for the "no uncertainty and no unit" case is in, we should test this on some realistic spectrum. Maybe use the spectrum from the "getting started" example? The goal is to check whether or not the three cases actually yield sensible answers, and also do some timing tests relative to other analysis functions like a line flux measurement.

specutils/analysis/flux.py Outdated Show resolved Hide resolved
specutils/analysis/flux.py Outdated Show resolved Hide resolved
specutils/analysis/flux.py Outdated Show resolved Hide resolved
specutils/analysis/flux.py Outdated Show resolved Hide resolved
if hasattr(spectrum, 'uncertainty') and spectrum.uncertainty:
return np.median(spectrum.flux / spectrum.uncertainty.quantity) < eps
else:
raise Exception('Spectrum flux has units, eps does not, either include uncertainty or use units on eps')
Copy link
Member

Choose a reason for hiding this comment

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

I think this else case should instead use the median absolute deviation (in the "scaled to be like sigma" form - https://docs.astropy.org/en/stable/api/astropy.stats.mad_std.html ) as the "uncertainty" estimate, rather than raising an error. Then the decorator can default to using that case with an eps of .01 or something like that, and we can count on this function always running on a spectrum.

(This is what makes it a heuristic for continuum subtraction - not all spectra that have their continuum subtracted really will have this property, but it's a reasonable guess for many spectra.)

Copy link
Author

Choose a reason for hiding this comment

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

Right, but, this was my question in #535 (comment).

Copy link
Member

Choose a reason for hiding this comment

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

Oops, sorry - to finish out that conversation here, I don't really understand your point there. In the case of your example linked there, that looks to have done exactly what was expected: the fluxes are close together, and that closeness indicates the "uncertainty per pixel" is probably low, so therefore the fact that the average offset from zero is larger than the scatter is exactly the sign of "continuum not subtracted" that we are looking for.

@brechmos-stsci
Copy link
Author

  1. I think the docstring needs to describe the flow of what happens given various values of eps and whether or not the uncertainty is present (but maybe you were waiting on a first pass?)

I was waiting for the first pass to be done before finishing up the documentation.

  1. This function doesn't seem to be using the mask. Presumably the user does not want masked values to be included, so I'd say the thing to do is just s.flux[~s.mask] if a mask is present.

Right, I was waiting to make sure the rest was good, but I added it in now.

  1. There needs to be some way to turn the decorator "off", which the warning message should explain. I suggest using the configuration item machinery (see http://docs.astropy.org/en/stable/config/#adding-new-configuration-items) with just a boolean config item like always_check_continuum that defaults to True but just causes the decorator to immediately return True if it's set to False. Then you tell the user to use the with version of the config item machinery to supress the check.

I added this in, but you or someone will have to check it to make sure it is implemented correctly .I looked at astropy to see how it was used there, and mimiced it.

  1. Once the modification for the "no uncertainty and no unit" case is in, we should test this on some realistic spectrum. Maybe use the spectrum from the "getting started" example? The goal is to check whether or not the three cases actually yield sensible answers, and also do some timing tests relative to other analysis functions like a line flux measurement.

Sure.

@brechmos-stsci brechmos-stsci marked this pull request as ready for review October 28, 2019 20:28
Copy link
Member

@eteq eteq left a comment

Choose a reason for hiding this comment

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

Several inline suggestions. Also, as per discussion in #535, should change the name to is_continuum_below_threshold

Configuration parameters for specutils.
"""

always_check_continuum = _config.ConfigItem(
Copy link
Member

Choose a reason for hiding this comment

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

Instead of "always", perhaps do_continuum_function_check or something like that? Only because "always" could be taken as "for all functions" and it's really just "selected functions".

Copy link
Author

Choose a reason for hiding this comment

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

changed

spectrum : `~specutils.spectra.spectrum1d.Spectrum1D`
The spectrum object over which the width will be calculated.

eps: float
Copy link
Member

Choose a reason for hiding this comment

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

eps -> threshold. And also should say that it can be a flux-unit quantity or a float/dimensionless qunatity, and that determines if it's a flux threshold vs a "signal/noise" threshold

Copy link
Author

Choose a reason for hiding this comment

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

already done

if uncertainty:
return np.median(flux / uncertainty.quantity) < eps
else:
warnings.warn('Spectrum does not appear to be continuum subtracted.', AstropyUserWarning)
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 confused why this warning is here?

Copy link
Author

Choose a reason for hiding this comment

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

I think you saw an old PR.

if check:
spectrum = args[0]
if not is_continuum_near_zero(spectrum, eps):
warnings.warn('Spectrum does not appear to be continuum subtracted.',
Copy link
Member

Choose a reason for hiding this comment

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

This should also have some text that refers back to always_check_continuum and tells the user exactly what they need to do to turn it off (or links to a page in the docs that explains it).

Copy link
Author

Choose a reason for hiding this comment

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

I added a short message with a link back to the astropy configuration section.

@eteq
Copy link
Member

eteq commented Oct 28, 2019

Oops, I think my review there was against an old version, but I didn't realize there was a new commit until I had already done the review. So ignore anything that was fixed in 72446dc that's in my comments.

Copy link
Member

@eteq eteq left a comment

Choose a reason for hiding this comment

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

A couple other inline items, but overall looking pretty good now!

def actual_decorator(function):
@wraps(function)
def wrapper(*args, **kwargs):
if check:
Copy link
Member

Choose a reason for hiding this comment

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

Hmm. I might be confused but I think this isn't quite right - I think you want to do

Suggested change
if check:
if check():

because a configuration item can be changed at runtime.

Copy link
Member

Choose a reason for hiding this comment

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

Or alternatively,

Suggested change
if check:
if check() if callable(check) else check:

which would let you pass in True/False in addition to a configuration item

Copy link
Author

Choose a reason for hiding this comment

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

Based on our conversation, this changed to an alternative form (that is correct :) ).

@@ -97,6 +97,7 @@ def _consecutive(data, stepsize=1):
return np.split(data, np.where(np.diff(data) != stepsize)[0]+1)


@warn_continuum_below_threshold(threshold=0.01, check=True)
Copy link
Member

Choose a reason for hiding this comment

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

I think you don't want check=True here (or below), right? Because then the do_continuum_function_check option gets ignored?

Copy link
Author

@brechmos-stsci brechmos-stsci Oct 29, 2019

Choose a reason for hiding this comment

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

hah, good point After we talked about it, we realized the check parameter doesn't even need to be there.

Copy link
Member

@eteq eteq left a comment

Choose a reason for hiding this comment

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

A few more things from a more detailed review. Mostly just docstring wordings but I also realized a couple of subtleties in the implementation.

specutils/analysis/flux.py Outdated Show resolved Hide resolved
specutils/analysis/flux.py Outdated Show resolved Hide resolved
specutils/analysis/flux.py Show resolved Hide resolved
specutils/analysis/flux.py Outdated Show resolved Hide resolved
specutils/analysis/flux.py Outdated Show resolved Hide resolved
specutils/analysis/flux.py Outdated Show resolved Hide resolved
specutils/analysis/flux.py Show resolved Hide resolved
specutils/analysis/flux.py Outdated Show resolved Hide resolved
Craig Jones and others added 6 commits October 30, 2019 14:53
Co-Authored-By: Erik Tollerud <erik.tollerud@gmail.com>
Co-Authored-By: Erik Tollerud <erik.tollerud@gmail.com>
Co-Authored-By: Erik Tollerud <erik.tollerud@gmail.com>
Co-Authored-By: Erik Tollerud <erik.tollerud@gmail.com>
Copy link
Member

@eteq eteq left a comment

Choose a reason for hiding this comment

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

I had a one-line addition but it was easier just to push up so I did that. Thanks @brechmos-stsci , this now all looks good!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

Check for continuum subtracted data
4 participants