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

Added kosmos trace algorithm #85

Merged
merged 18 commits into from
Mar 29, 2022
Merged

Conversation

ojustino
Copy link
Contributor

@ojustino ojustino commented Mar 4, 2022

I ported a trace algorithm from the KOSMOS repository's apextract.py. I refactored much of it and took on one of wishlist requests of using astropy.modeling's Gaussian1D instead of a custom function to create Gaussians. The fits end up being slightly different, but I have a local notebook where I've verified that the calculated inputs to the respective fitters are the same in this version and the original.

I wrote it as a descendant of Trace instead of as a class method of ArrayTrace as suggested in the tracing discussion. ArrayTrace requires an initial array, while the KOSMOS implementation can run with no guess or just a float. I didn't think requiring an array when it wasn't needed made sense, though we can change the format if I've missed something.

The original has arguments for specifying wavelength and spatial directions that I commented out since it seems that's still up for discussion. I also left out functionality for displaying the trace after the calculation.

Copy link
Member

@kecnry kecnry left a comment

Choose a reason for hiding this comment

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

I would still vote in favor of trace = ArrayTrace.kosmos_from_image(image, bins, guess, window) or trace = ArrayTrace.from_image(image, alg='kosmos', **kwargs) (which wouldn't require an array, it creates and sets the array) unless there will be additional methods we expect to add for this that won't be applicable to the general parent class (or if we expect that to be the case for other algorithms). My thinking is that the inheritance is not as obvious to the user that they really have the same functionality as an ArrayTrace object, just constructed differently. But we can defer that decision for a while, if necessary, as it really is just moving around existing code.

This also really needs a test case, which will also make it easier to make sure its behaving as expected and thoroughly step through the code.

specreduce/tracing.py Show resolved Hide resolved
specreduce/tracing.py Outdated Show resolved Hide resolved
specreduce/tracing.py Outdated Show resolved Hide resolved
Comment on lines +196 to +198
width_guess = (2 if width_guess < 2
else 25 if width_guess > 25
else width_guess)
Copy link
Member

Choose a reason for hiding this comment

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

this scares me a bit... should we at least add warnings?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

I copied these numbers as they were, so maybe we should learn where they come from and whether they're worth keeping before adding warnings. @jradavenport, do you remember how you came up with this range of valid guesses for the width of the trace's peak?

Copy link
Contributor

Choose a reason for hiding this comment

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

they were good values I arrived at experimentally to ensure the trace didn't go off the rails... 2 pixels being a minimum to (hopefully) ensure I wasn't locking onto a cosmic ray. 25 pixel maximum to keep me looking at a point source instead of just noise

Copy link
Contributor Author

Choose a reason for hiding this comment

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

What size images were you working with? Maybe we can turn the 25 into a fraction of the image's size in the cross-dispersion direction.

Copy link
Contributor

Choose a reason for hiding this comment

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

i think this works more or less well enough for now. DRAGONS has some more sophisticated algorithms for peak finding, centering, and width measurements that we'll likely want to port over. we can compare/contrast methods when we get to that point and adjust as necessary.

specreduce/tracing.py Outdated Show resolved Hide resolved
Copy link
Contributor

@duytnguyendtn duytnguyendtn left a comment

Choose a reason for hiding this comment

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

Substituting the FlatTrace with KosmosTrace in the boxcar extract notebook produced results that by-eye looked appropriate. I strongly suggest a scientific review from someone who can scientifically validate the results (see @PatrickOgle and of course @tepickering)

specreduce/tracing.py Outdated Show resolved Hide resolved
specreduce/tracing.py Outdated Show resolved Hide resolved
Comment on lines 193 to 198
width_guess = yy_above_half_max / 2.355

# enforce some (maybe sensible?) rules about trace peak width
width_guess = (2 if width_guess < 2
else 25 if width_guess > 25
else width_guess)
Copy link
Contributor

Choose a reason for hiding this comment

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

While I'm all for intelligent software, I share @kecnry's concern about these hardcoded values. I'm wondering if these values can be extracted from some attribute of the image that's coming in? But that's outside my expertise

Copy link
Contributor Author

@ojustino ojustino Mar 9, 2022

Choose a reason for hiding this comment

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

I agree on the if/else clause, but the 2.355 in the first line you've highlighted comes from this equation relating full width at half maximum to standard deviation.

specreduce/tracing.py Outdated Show resolved Hide resolved
ojustino and others added 2 commits March 9, 2022 16:27
Co-authored-by: Kyle Conroy <kyleconroy@gmail.com>
Co-authored-by: Duy Tuong Nguyen <dtn5ah@virginia.edu>
@PatrickOgle
Copy link

It seems to work on the example data in compare_extractions, replacing FlatTrace with KosmosTrace(img). However, I do not understand the output format of the trace. I would expect a 1D array of values with the trace location for each x value in the input image. Please add a notebook example of how KosmosTrace works and overplot the location of the computed trace on the input 2D image. Also, the fake dataset in that notebook would be more interesting if it had a curved trace.

@ojustino
Copy link
Contributor Author

@PatrickOgle the output shouldn't be different. For either FlatTrace or KosmosTrace, you should be able to get that 1D array of trace locations using the trace attribute of the finished product (e.g., result.trace).

@ojustino
Copy link
Contributor Author

ojustino commented Mar 15, 2022

I'm collecting current feedback in one comment. Some are noted as still being in progress, but please offer opinions and highlight those that are preventing a merge.

  1. Format. Is there any more feedback on whether we prefer this functionality to be a method of ArrayTrace or inherit from Trace as its own class?
  2. Tests. I'm writing some.
  3. Notebook example? Patrick has asked for a notebook example, but would one be compelling? Not as much goes into preparing for a trace as with the extractions.
  4. Magic numbers. Being discussed in this thread.
  5. Licensing. I'll copy the steps I mentioned here if there are no objections.
  6. disp_axis (EDIT) I will turn disp_axis into a private variable and add a check that its value is 1 before calculating the trace. Hopefully, we can make it an argument in the future.

@PatrickOgle
Copy link

PatrickOgle commented Mar 16, 2022

OK, knowing that I can access the trace via result.trace enabled me to overplot the trace on the image. Here is an example with bins=8. The result from the kosmos spline-fitting algorithm is somewhat more sensitive to variations in noise than I would like. We should consider (for a future PR) adding a trace class for a polynomial function of order n. In this particular example, I would likely fit a linear function.

kosmos_trace_8bins

Copy link

@PatrickOgle PatrickOgle left a comment

Choose a reason for hiding this comment

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

Works as advertised. I was a bit confused by the output format--that I had to retrieve the result by a method, and that the class carries along the input image. However, that is OK if it is properly documented and examples are eventually given. The result of the spline fit is a bit sensitive to the noise. We should consider adding another trace method that uses a polynomial fit rather than a spline.

@ojustino
Copy link
Contributor Author

Everything from my summary coment except the example notebook (marked as a future task) should either be resolved or commented on with # NOTE in the code from these newest commits.

Otherwise, the trace now works on NDData objects (previously assumed an array) and is more explicit about how masked values are handled. Results are still the same for non-masked images like the one Patrick tested earlier.

@tepickering
Copy link
Contributor

great! i'll give things a full review and do a notebook run-through myself later today so we can get this merged.

Copy link
Member

@kecnry kecnry left a comment

Choose a reason for hiding this comment

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

All my comments have either been addressed or postponed to future work. The code itself and tests look good. Whichever get merged first (this or #86), it would be useful to rebase the other and update the notebook to cover the case of an initial rough trace -> background subtraction -> kosmos trace -> extraction.

specreduce/tracing.py Show resolved Hide resolved
make window size error message more explicit

Co-authored-by: Duy Tuong Nguyen <dtn5ah@virginia.edu>
Copy link
Contributor

@tepickering tepickering left a comment

Choose a reason for hiding this comment

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

left a couple comments. it might be worth adding an upper limit check for bins, but it's also fine to leave it and make an issue for it. there's still more work to do on tracing and that fix can be done as part of that.

Comment on lines +196 to +198
width_guess = (2 if width_guess < 2
else 25 if width_guess > 25
else width_guess)
Copy link
Contributor

Choose a reason for hiding this comment

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

i think this works more or less well enough for now. DRAGONS has some more sophisticated algorithms for peak finding, centering, and width measurements that we'll likely want to port over. we can compare/contrast methods when we get to that point and adjust as necessary.

tepickering and others added 3 commits March 25, 2022 19:19
Fix error message for `self.window` to adhere to line length limits.
fix f-string in window error message
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.

6 participants