Add SFD map code to the DESI utilities #116

geordie666 · 2018-08-27T21:21:23Z

This PR transfers the dust map code from desitarget.mock.sfdmap to desiutil.

The context for porting this across is that the dust maps are now used in several places in desitarget (not just the mocks) and will likely also be needed for future updates to desispec. Some notes:

Changes the default environment variable SFD_DIR to DUST_DIR
Adds unit tests, which didn't exist in desitarget

…ts (desitarget in two places, and desispec, too, ultimately)

…es consistent with imaging

…ently produce E(B-V) values

moustakas · 2018-08-27T21:36:25Z

This looks outstanding. My only comment is that you consider renaming the sfdmap module something more generic like dust, so we can (in a separate PR) add one or more extinction curves, which will be needed in both desitarget.mock and desispec (to deredden spectra for Galactic extinction).

weaverba137 · 2018-08-27T22:02:36Z

There are some very significant problems with this PR. But there is clearly still some active development going on, so please let me know when things are stable so I can review.

geordie666 · 2018-08-27T22:32:58Z

@weaverba137: Looks like the critical unit tests are now passing, so feel free to elucidate the very significant problems.

weaverba137 · 2018-08-27T23:24:41Z

The last change came in just as I was about to start a review.

geordie666 · 2018-08-27T23:26:09Z

Go for it with the review. That's a completely cosmetic change (hence the [ci skip]). I was just improving coding style, but it's nothing significant.

weaverba137 · 2018-08-27T23:28:01Z

A GitHub review needs to be attached to a particular commit. If a new commit comes in while the review is going on, you have to refresh the files changed view and sometimes all your review comments can be lost.

geordie666 · 2018-08-27T23:31:00Z

Noted. In future I won't make commits while reviews are pending, or I'll wait until I'm done with all commits before requesting a review. Hopefully that didn't happen in this case.

weaverba137

In summary, the code should be rewritten to match DESI standards before it can be accepted.

weaverba137 · 2018-08-27T23:30:05Z

py/desiutil/dust.py

@@ -0,0 +1,384 @@
+# Copied on Nov/20/2016 from https://github.com/kbarbary/sfdmap/ commit: bacdbbd


Desiutil already has a license, and we need to ensure that all code is compatible with that license. Maybe we were sloppy about that when the code was in desitarget, but desiutil is typically held to a much higher standard.

I remain uncertain as to what to do about reusing code licensed in this manner. Now that you've brought this to my attention I certainly wouldn't merge this until we know the implications. I'm hoping @sbailey has some bright ideas.

weaverba137 · 2018-08-27T23:30:39Z

py/desiutil/dust.py

+desiutil.dust
+=============
+
+Get E(B-V) values from the Schlegel, Finkbeiner & Davis (1998) dust map.


It would be helpful to add a link to the paper itself, or the abstract on ADS.

weaverba137 · 2018-08-27T23:31:16Z

py/desiutil/dust.py

+    except ImportError:
+        raise ImportError("could not import fitsio or astropy.io.fits")
+
+__all__ = ['SFDMap', 'ebv']


No one in DESI should be doing from desiutil.dust import *, so this line is not necessary.

I have removed this line.

weaverba137 · 2018-08-27T23:31:48Z

py/desiutil/dust.py

+        raise ImportError("could not import fitsio or astropy.io.fits")
+
+__all__ = ['SFDMap', 'ebv']
+__version__ = "0.1.0"


Desiutil already has a version, desiutil.__version__, so this really shouldn't be here.

I've removed the version.

weaverba137 · 2018-08-27T23:32:10Z

py/desiutil/dust.py

+__version__ = "0.1.0"
+
+
+def _isiterable(obj):


There has to be a function in the Python standard library that does this.

I used np.atleast_1d() instead

weaverba137 · 2018-08-27T23:33:13Z

py/desiutil/dust.py

+
+# Create rotation matrix about a given axis (x, y, z)
+
+def zrotmat(angle):


All functions need documentation strings, so they can be processed by Sphinx.

I have gone through the code and added sphinx-compliant doc strings. So, this should be solved, If unit tests pass.

weaverba137 · 2018-08-27T23:34:08Z

py/desiutil/dust.py

+
+
+# constant ICRS --> FK5J2000 (See USNO Circular 179, section 3.5)
+try:


Why is this computation being performed at the module-level? There is even a comment below that notes that!

I have removed all of this code and use astropy.coordinates.SkyCoord objects instead.

weaverba137 · 2018-08-27T23:34:54Z

py/desiutil/dust.py

+
+
+# -----------------------------------------------------------------------------
+# bilinear interpolation (because scipy.ndimage.map_coordinates is really slow


This assertion should be justified, and if it is no longer true, we should use numpy/scipy/astropy.

I haven't actually checked that (as I have a stable and working new version of the code and don't want to mess around with it much further), but I have removed that assertion. If you feel really strongly that we should use scipy then I'll assess this in detail.

weaverba137 · 2018-08-27T23:37:06Z

py/desiutil/test/test_dust.py

+from .. import dust
+
+
+class TestDust(unittest.TestCase):


Assuming that DUST_DIR is defined, do the tests actually cover the entire module? In other words, if we accept the module into desiutil, the total test coverage of desiutil should at least not decrease.

I've added a unit test and coverage is now ~94%.

weaverba137 · 2018-08-27T23:39:58Z

py/desiutil/dust.py

+import numpy as np
+
+# require a FITS reader of some sort.
+try:


This is not a common construction that I have seen in DESI software. In fact, desiutil already requires astropy, so there is no reason to add an additional dependency on fitsio.

I have explicitly imported astropy.io.fits instead of this construction.

…out the license, though

…emove version number

weaverba137 · 2018-08-28T15:54:30Z

An additional note as part of the previous review: as currently configured, no tests of this module are performed on Travis, because DUST_DIR is not currently defined in the Travis environment. This also means that Travis will report to Coveralls that this module is not covered.

Data files are required to test this module, and would need to be brought into the Travis environment. There are a few ways this could happen:

Download files as needed from a public site.
Include a "trimmed" version of the files in the desiutil test suite.
Use unittest.mock to simulate reading the files and inject appropriate data.

…ent variable

…only ~1.5x.

geordie666 · 2018-08-28T22:36:55Z

@weaverba137: Regarding your additional note...I chopped the dust maps down to the first 10 columns to reduce them in size to ~150K and use these chopped dust maps as the data for explicit unit tests. Your option (2), in other words.

weaverba137 · 2018-08-28T22:47:19Z

The code changes all look good to me know. But let's make sure tests pass, & get feedback on license before merging.

geordie666 · 2018-08-28T22:48:18Z

I'm going to clean up my code style, so give me 10 minutes before anyone looks at this.

weaverba137 · 2018-08-28T22:52:06Z

OK, I'll double-check when you're done.

geordie666 · 2018-08-28T23:54:22Z

I'm done, now. There are a handful of other recommended code style alterations, but it's for a case where I think it's easier to read as-is.

weaverba137 · 2018-08-29T00:04:24Z

I'm trying to build the sphinx documentation on my own system and the dust module is not showing up. Watch for additional commits to fix that.

moustakas · 2018-08-30T16:49:13Z

@weaverba137 @geordie666 Is this ready to be merged?

If so, @geordie666 can you update desitarget.mock.sfdmap to use the dust utilities in desiutil and add a deprecation warning?

geordie666 · 2018-08-30T16:53:19Z

@moustakas: We're waiting for guidance on what to do about reconciling the license for this code with the DESI BSD license from @sbailey.

moustakas · 2018-09-02T12:23:35Z

While this PR is still open, should we take the time to add one or more extinction curves?

Kyle Barbary has a very nice package which wraps on Cython for speed --
https://extinction.readthedocs.io/en/latest/

This package is MIT licensed, as was his implementation of reading the SFD dust maps, so we'd have to resolve the same licensing question(s). In addition, these curves would naturally go into the same desiutil.dust module.

dkirkby · 2018-09-02T15:04:24Z

+1 for adding at least one extinction curve. Schlaflay 2016 is also an option (code here). I doubt this would be a bottleneck for our applications, so straight numpy (w/o cython) should be sufficient if we impement this ourselves.

djschlegel · 2018-09-02T15:21:45Z

Let's not pull in another package just to read the SFD maps. That should quite literally be no more than 5 lines of code. Convert to Galactic coordinates, linearly interpolate the 4 nearest pixels from
x = 2048 SQRT {1 - n sin(b)} cos(l) + 2047.5
y = - 2048 n SQRT{1 - n sin(b)} sin(l) + 2047.5
where n=+1 for the NGP, and n=-1 for the SGP. The May 1999 update of the files added the (then-new) FITS header keywords to do the above, so maybe it's only 1 line if your FITS reader supports?

Extinction curves are a different question, and should be a different issue. One needs to de-extinction the standard stars during calibration, but I'd argue to keep the spectra in as-observed fluxes (no de-extinction), analogous to how photometric fluxes are reported. Please include Schlafly on that issue.

sbailey · 2018-09-05T17:28:35Z

I'm joining this PR comment stream late, so apologies if I miss the point on something:

I think license situation is fine: the main desiutil package is BSD licensed, and the dust.py file includes a comment that some portions were copied from https://github.com/kbarbary/sfdmap/ with the following MIT license, which seems to match both the spirit and the legalities of the open source licenses. It would be slightly better if the comment specifically mentioned which functions / classes were copied (or is it the whole file?).
Travis tests: looks like that has been resolved now, but in the future it is ok if Travis only runs a subset of the full tests if that is for a good reason, e.g. to avoid the overhead of large test files. The full suite runs nightly at NERSC so we do get regular coverage.
re: "Let's not pull in another package just to read the SFD maps". I agree that we shouldn't pull in a new external dependency that has to be separately installed just to read the dust maps, but it is good to have a standardized convenience way to do that within our own code (e.g. here in desiutil rather than N>1 nearly identical implementations of the same thing in desitarget, desispec, desisurvey, desisim, ...)

weaverba137 · 2018-09-05T17:36:14Z

Thanks @sbailey.

we do get regular coverage

Yes, but that coverage is not reported to coveralls.io. Thus the coverage reported is lower than it should be, and there is not the immediate feedback of the Travis tests.

weaverba137 · 2018-09-05T17:37:31Z

PS, should we tag desiutil after this merge?

sbailey · 2018-09-05T20:13:54Z

Agreed that having tests included in Travis is better when possible without heroic efforts.

I'll merge this now.

Please go ahead and make a new tag and install at NERSC. I'll let @weaverba137 do that in order to pick the appropriate level of version number incrementing.

weaverba137 · 2018-09-05T20:29:22Z

1.9.12 is installed on cori and edison.

geordie666 added 7 commits August 27, 2018 13:03

moving dust map code to desiutil, as it's now used by multiple produc…

e98901e

…ts (desitarget in two places, and desispec, too, ultimately)

updating instances of SFD_DIR to DUST_DIR to make environment variabl…

2e213eb

…es consistent with imaging

updating default scaling to 1, as that's how the imaging surveys curr…

8760ed7

…ently produce E(B-V) values

adding some unit tests for the SFD map code

76eac90

adding SFDmap code to the read the docs api

83502fe

updating the changes docs [ci skip]

b2c0cfa

skipping test if the DUST directory doesn't exist

2e4f033

mocking up fitsio in sphinx config file

5b70eae

weaverba137 self-requested a review August 27, 2018 21:52

geordie666 added 2 commits August 27, 2018 14:59

changing 'sfdmap' name to 'dust' throughout

15b8947

changing 'sfdmap' name to 'dust' throughout

4f06846

geordie666 added 2 commits August 27, 2018 15:14

cosmetic changes to coding style

3c816c8

further cosmetic changes to coding style

4329fad

additional cosmetic changes [ci skip]

ab197a2

weaverba137 requested changes Aug 27, 2018

View reviewed changes

geordie666 added 3 commits August 28, 2018 08:46

adding link to SFD98 abstract, augmenting license...still not sure ab…

4fd5e10

…out the license, though

Use astropy as FITS reader, don't expect 'from x import *' imports; r…

defe0aa

…emove version number

use numpy.atleast_1d in order to remove the _isiterable function

c38b5ee

geordie666 added 4 commits August 28, 2018 09:21

should be importing fits reader from astropy.io not from astropy

ac7a4f2

explicitly restore scalar returns for backwards-compatability

f63871b

updating unit tests to use small fake map instead of rely on environm…

d307b4b

…ent variable

adding a new unit test

65dbc6f

geordie666 added 5 commits August 28, 2018 13:40

switch to using astropy for coordinate conversions. It's slower, but …

5ad3b10

…only ~1.5x.

handle tuples properly

8659180

additional unit test coverage

db3b7ec

improving doc strings for sphinx

2156bf7

finishing up doc strings

028c840

mock up astropy for sphinx

804cabf

weaverba137 approved these changes Aug 28, 2018

View reviewed changes

geordie666 added 2 commits August 28, 2018 16:01

cleaning up coding style

f81990a

a couple of additional coding style changes

ea609dc

fix some fiddly Sphinx formatting problems

9500148

sbailey merged commit 3fefb5e into master Sep 5, 2018

sbailey deleted the ADMsfdmap branch September 5, 2018 20:14

moustakas mentioned this pull request Sep 5, 2018

use new desiutil.dust.SFDMap module desihub/desitarget#366

Merged

		@@ -0,0 +1,384 @@
		# Copied on Nov/20/2016 from https://github.com/kbarbary/sfdmap/ commit: bacdbbd


		# Create rotation matrix about a given axis (x, y, z)

		def zrotmat(angle):



		# constant ICRS --> FK5J2000 (See USNO Circular 179, section 3.5)
		try:



		# -----------------------------------------------------------------------------
		# bilinear interpolation (because scipy.ndimage.map_coordinates is really slow

Add SFD map code to the DESI utilities #116

Add SFD map code to the DESI utilities #116

Conversation

geordie666 commented Aug 27, 2018

moustakas commented Aug 27, 2018

weaverba137 commented Aug 27, 2018

geordie666 commented Aug 27, 2018 • edited Loading

weaverba137 commented Aug 27, 2018

geordie666 commented Aug 27, 2018

weaverba137 commented Aug 27, 2018

geordie666 commented Aug 27, 2018

weaverba137 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

weaverba137 commented Aug 28, 2018

geordie666 commented Aug 28, 2018

weaverba137 commented Aug 28, 2018

geordie666 commented Aug 28, 2018

weaverba137 commented Aug 28, 2018

geordie666 commented Aug 28, 2018

weaverba137 commented Aug 29, 2018

moustakas commented Aug 30, 2018

geordie666 commented Aug 30, 2018

moustakas commented Sep 2, 2018

dkirkby commented Sep 2, 2018

djschlegel commented Sep 2, 2018

sbailey commented Sep 5, 2018

weaverba137 commented Sep 5, 2018

weaverba137 commented Sep 5, 2018

sbailey commented Sep 5, 2018

weaverba137 commented Sep 5, 2018

geordie666 commented Aug 27, 2018 •

edited

Loading