{{ message }}

# Consistent random number handling and improve sample_sphere #283

Merged
merged 10 commits into from Jul 30, 2015
Merged

# Consistent random number handling and improve sample_sphere#283

merged 10 commits into from Jul 30, 2015

## Conversation

### mapazarr commented Jul 21, 2015

 ~gammapy.utils.random.sample_sphere should use astropy Angle objects for better handling / less error prone. TODO: decide if we use Angle or Longitude and Latitude objects update/revise other functions using sample_sphere finish applying issue #238 (random_state parameters)
changed the title ~gammapy.utils.random.sample_sphere should use ~astropy.Qauntity objects ~gammapy.utils.random.sample_sphere should use ~astropy.Quantity objects Jun 9, 2015
mentioned this pull request Jun 9, 2015
added the cleanup label Jun 10, 2015
added this to the 0.3 milestone Jun 10, 2015
changed the title ~gammapy.utils.random.sample_sphere should use ~astropy.Quantity objects sample_sphere should use Quantity objects Jun 10, 2015

### cdeil commented Jun 10, 2015

 @mapazarr - I assigned this issue to you. And I shortened the title and removed the RST markup, which is not supported in Github issue titles. The goal should be that https://github.com/gammapy/gammapy/issues and https://gammapy.readthedocs.org/en/latest/changelog.html#id4 look nice if I just copy & paste the PR titles.

### cdeil commented Jun 10, 2015

 It can also be useful to link to docs / source code in the issue description, in this case to gammapy.utils.random.sample_sphere, which I typed as [gammapy.utils.random.sample_sphere](https://gammapy.readthedocs.org/en/latest/api/gammapy.utils.random.sample_sphere.html)  To be more specific, this function should probably use astropy.coordinates.Longitude and astropy.coordinates.Latitude objects, which are Quantity sub-classes.

### mapazarr commented Jul 21, 2015

 Implemented. Merge?
assigned cdeil and unassigned mapazarr Jul 21, 2015

### cdeil commented Jul 21, 2015

 I'd prefer it if we use the more explicit Longitude and Latitude objects from astropy.coordinates. They have convenience methods for angle wrapping. I've never used them much, so if that doesn't work well for some reason, please at least use Angle objects.
 units : {'rad', 'deg'} Units of input range and returned angles lon_range : ~astropy.units.Quantity, optional Longitude range (min, max) in range (0, 360) deg

#### cdeil Jul 21, 2015 Member

Would it be easy to make the code just work for longitudes in the range (0, 360) and (-180, 180)?
I.e. not have this requirement of lon_range being in the range (0, 360)?

#### mapazarr Jul 22, 2015 Author Member

OK. Done.

 lon_range = np.radians(lon_range) if lon_range is not None: lon_unit = lon_range.unit lon_range.to('radian')

#### cdeil Jul 21, 2015 Member

@mapazarr – Here I'm not sure what the best API is.
Is it better to always return angles in the same units as the inputs?
What's more intuitive / convenient for callers?

#### mapazarr Jul 22, 2015 Author Member

It was intended for returning the same units as the input, but actually it's nonsense. After looking again to the function, I have simplified this part.

 lat_range = np.radians(lat_range) if lat_range is not None: lat_unit = lon_range.unit lat_range.to('radian')

#### cdeil Jul 21, 2015 Member

I think this statement has no effect:

In [3]: a = Angle(10, 'degree')

In [5]: a
Out[5]: <Angle 10.0 deg>


Also: the numpy trig functions work with Angle, so no need for unit conversion to radian on input?

#### mapazarr Jul 22, 2015 Author Member

As I said, this has been simplified. You're right no need to transform units: it was a relic of the old version. I removed it already.

### cdeil commented Jul 21, 2015

 Can you please look at the tests for this function and see if they cover the common and some edge cases well, and if not, please add a test or two?

### mapazarr commented Jul 22, 2015

 Sorry, I missed the Longitude Latitude part of the comment. I tried to use them, but it's not straightforward (the function call in the test returns (0, 0)), and would need deeper debugging, for which I have no time right now. Indeed, Angle would be much more appropriate as Quantity. I am using Angle right now.

### mapazarr commented Jul 22, 2015

 I added a few more tests.
added 4 commits Jul 21, 2015
 Using Quantity objects in sample_sphere. 
 369efbd 
 Added entry to CHANGES log. 
 2f07b12 
 sample_sphere using Angle, not Quantity 
 f779df7 
 Simplified sample_sphere, added support for lon in (-180, 180) deg ra… 
 37675c1 
…nge and added a few more tests.
force-pushed the mapazarr:sample_sphere_cleanup branch from 223711e to 37675c1 Jul 22, 2015
 assert_quantity_allclose(lat, Angle([0.20700192, 0.08988736], 'radian')) # test specify a limited range lon_min = Angle(40., 'degree')

#### cdeil Jul 22, 2015 Member

Just to save a few lines, please use

lon_range = Angle([40, 45], 'deg')

here and then pass this to sample_sphere directly and use lon_range[0] and lon_range[1] in the asserts below.

#### mapazarr Jul 23, 2015 Author Member

Done.

 lon, lat = sample_sphere(size=10, lon_range=Angle([lon_min, lon_max]), lat_range=Angle([lat_min, lat_max])) assert (lon_min <= lon).all() & (lon < lon_max).all()

#### cdeil Jul 22, 2015 Member

I find this a bit more readable (avoids the duplicated .all():

assert np.alltrue(lon_min <= lon) & (lon < lon_max))

I'm not sure if either of these spellings is the best way to write this assert (i.e. short and easy to read, but still give a good message when it fails).

There's also numpy.testing.assert_array_less, but I haven't used it.

So please take 10 minutes to thing / play with this some more, to make sure we have a good pattern for such asserts (we'll probably have 100s once we start adding lots of tests for Gammapy).
But in the end, yes, this is not very important and whatever you put, I'll merge.

#### mapazarr Jul 23, 2015 Author Member

This is shorter:

assert ((lon_range[0] <= lon) & (lon < lon_range[1])).all()

#### mapazarr Jul 23, 2015 Author Member

Drawback: it doesn't show the values in case of failure. For this reason I also tested numpy.testing.assert_array_less, but this has the same problem as #306, an there is no good method in Astropy. I also don't think we should start reimplementing all numpy assertions for quantities, So I'll use the above-mentioned syntax.

 return np.degrees(lon), np.degrees(lat) if lon_format == '-180_180_deg': # convert back to (-180, 180) deg lon -= Angle(180., 'degree')

#### cdeil Jul 22, 2015 Member

You're re-implementing longitude angle wrapping here, although it's available in Astropy:

I know you said you tried using Longitude objects and there was some issue.
But I'd prefer you work out that issue even if it takes an extra hour or two and we start using Longitude objects in Gammapy.
There's many places already where we need this (catalogs, observation lists) and using Longitude objects should make for easier / more readable code.

If you can't get this to work Ok with Longitude, can you post a short explanation why?
Could it be solved with a feature request for Angle or Longitude in Astropy?

In the end, if you conclude that it's just in this specific function that Longitude doesn't work for some reason,
OK, I'll still merge as-is.

#### mapazarr Jul 23, 2015 Author Member

Ok, I gave it some more thought, and revised the docs about Longitude and Latitude. I problem I see is that the Longitude class, because of the wrapping method, is not very convenient for setting longitude ranges. For instance, when setting the range (0, 360) deg, it will be wrapped to (0, 0) deg!!! One would have to define the 2 \pi interval as (0, 360-epsilon) deg, with epsilon = 1e-8 for instance. And to remember this every time one defines the longitude interval is kind of hard (it's easy to forget).

In my opinion, Longitude does not make things easier (at least for the case of this function). With Angle you need to wrap angles, but you don't need to think about the wraping of the intrerval length. With Longitude, you don't need to wrap angles, but you need to remember to check if the interval is wrapped correctly.

 @@ -54,6 +54,7 @@ Pull requests - Remove gammapy.shower package [#291] (Christoph Deil) - Add cube background model class [#299] (Manuel Paz Arribas) - Use assert_quantity_allclose from Astropy [#306] (Manuel Paz Arribas) - sample_sphere should use Quantity objects [#283] (Manuel Paz Arribas)

#### cdeil Jul 22, 2015 Member

Change PR title and this line in the changelog to:

Change sample_sphere to use Angle objects


#### mapazarr Jul 29, 2015 Author Member

Done (forgot to comment on this).

 sample_sphere using Longitude and Latitude 
 8c048ba 
changed the title sample_sphere should use Quantity objects sample_sphere should use Angle objects Jul 23, 2015
changed the title sample_sphere should use Angle objects Change sample_sphere to use Angle objects Jul 23, 2015
 Updated name of PR in CHANGES log. 
 557cb94 

### mapazarr commented Jul 23, 2015

 Ok, the new version uses Longitude and Latitude, but please look at my comments and the TODO list I just added at the beginning of the issue.
 # test lon range explicitly (0, 360) deg epsilon = 1.e-8 lon_range = Longitude([0., 360.-epsilon], 'degree')

#### mapazarr Jul 23, 2015 Author Member

@cdeil This way of defining the 360 deg boundary is necessary to avoid the Longitude wrapper to convert it to 0 deg. I'm not sure this is better than just using Angle objects, even if the angles might need to be wrapped.

BTW, I'm not even sure if we have to wrap them at all. Won't the method just work fine, provided the angle ranges are not broader than 360 deg? (I haven't tested it).

mentioned this pull request Jul 29, 2015

### mapazarr commented Jul 29, 2015

 I am applying the fixes for issue #238 here. I started adding random_state parameters to some of the functions, and added the entry in the high-level doc.
 catalog = make_catalog_random_positions_sphere(size=100, center='Milky Way') catalog['Source_Name'] = ['source_{:04d}'.format(_) for _ in range(len(catalog))] catalog['Association_Radius'] = Angle(10 * np.random.random(len(catalog)), unit='deg') #catalog['Association_Radius'] = Angle(rng.uniform(0, 10, len(catalog)), unit='deg') # TODO: not working with rng!!! (but why did it work before??!!!)

#### mapazarr Jul 29, 2015 Author Member

If I use check_random_state at this point, the result is not reproducible and the test at the end fails.

But the more strange thing is: why was it reproducible before? The test uses random functions: make_catalog_random_positions_sphere, which makes use of sample_sphere and sample_sphere_distance without passing the seed as a parameter.
Was the seed passed implicitly, when setting it at the beginning of the test?

Ideas? I guess I need to pass the random_state parameter to sample_sphere and sample_sphere_distance all the way through make_catalog_random_positions_sphere, right?

#### mapazarr Jul 29, 2015 Author Member

Alternativelly: I can adapt only the functions that use a random generator, and leave the tests unchanged. This would be faster.

#### cdeil Jul 29, 2015 Member

This was deterministic before because of the np.random.seed(0) call at the top of the test_catalog_xmatch_circle function. Because we were using the global numpy rng everywhere, it was deterministic.

As you say, now that we're using the pattern that all functions using rng calls take the random_state parameter, yes, you have to pass it around.

#### mapazarr Jul 30, 2015 Author Member

Done.

 @@ -126,21 +138,24 @@ def make_catalog_random_positions_sphere(size, center='Earth', distance=Quantity return table ''' def make_cat_gauss_random(n_sources=100, glon_sigma=30, glat_sigma=1,

#### mapazarr Jul 29, 2015 Author Member

@cdeil This function was placed between 2 sets of triple single quotes. Why? Was it an bug? Was the function commented out for some reason? Or are these quotes this some kind of decorator? I removed the quotes.

#### cdeil Jul 29, 2015 Member

It's commented out.
Don't update it here, we might remove it.

If you want, you can add a TODO or make an issue that these functions need to be removed or improved.
(no time to deal with it at the moment)

#### mapazarr Jul 30, 2015 Author Member

Done.

 Started to apply issue #238: random_state parameter. 
 bdc3520 
force-pushed the mapazarr:sample_sphere_cleanup branch from 554618c to bdc3520 Jul 29, 2015
 @@ -54,6 +54,7 @@ Pull requests - Remove gammapy.shower package [#291] (Christoph Deil) - Add cube background model class [#299] (Manuel Paz Arribas) - Use assert_quantity_allclose from Astropy [#306] (Manuel Paz Arribas) - Change sample_sphere to use Angle objects [#283] (Manuel Paz Arribas)

#### cdeil Jul 29, 2015 Member

Please change the title of this PR and the changelog entry to reflect the content.
Maybe "Consistent random number handling and improve sample_sphere"?

#### mapazarr Jul 30, 2015 Author Member

Done.

 Continuing with issue #238: random_state parameter. 
 381d67c 
changed the title Change sample_sphere to use Angle objects Consistent random number handling and improve sample_sphere Jul 30, 2015
added 2 commits Jul 30, 2015
 Finished issue #238: random_state parameter. 
 2385e88 
 Going back to use Angle in sample_sphere. 
 5105a5b 

### mapazarr commented Jul 30, 2015

 @cdeil I am done with this PR. Ready for final review and evtly. merge.

### cdeil commented Jul 30, 2015

 Thanks!
added a commit that referenced this pull request Jul 30, 2015
 Merge pull request #283 from mapazarr/sample_sphere_cleanup 
 4504ec0 
Consistent random number handling and improve sample_sphere
merged commit 4504ec0 into gammapy:master Jul 30, 2015
1 check passed
1 check passed
continuous-integration/travis-ci/pr The Travis CI build passed
Details
mentioned this pull request Jul 30, 2015

### cdeil commented Jul 30, 2015

 I'm continuing with this a bit in #310.