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

Improve EnergyDispersion2D get_response and tests #911

Merged
merged 5 commits into from Feb 24, 2017

Conversation

@registerrier
Copy link
Contributor

@registerrier registerrier commented Feb 23, 2017

This introduces a different get_response() where integration is done with a given step in migration rather than with a fixed oversampling of e_reco bins.

  • The default value for the migra_step allows a significant decrease of computation time and allow for absolute precision <0.05.
  • Tests have been modified to test for a gaussian edisp where we can compute the desired value with satisfactory precision (<5% absolute)
  • Note that all tests on EnergyDispersion2D could be done with gaussian edisp rather than stored one.
@registerrier registerrier requested review from cdeil and joleroi Feb 23, 2017
@joleroi joleroi self-assigned this Feb 23, 2017
Copy link
Contributor

@joleroi joleroi left a comment

Thanks for this implementation and especially the test, it's always good to have tests from people that know what they're doing 😃 I left 2 comment mostly about code organization. I trust you that the gaussian edisp is coded correctly!


# Determine positions (bin indices) of e_reco bounds in migration array
# We ensure that no negative values are found
pos_mig = np.maximum(np.digitize(migra_e_reco, mig_array)-1 , 0)

This comment has been minimized.

@joleroi

joleroi Feb 23, 2017
Contributor

This line is quite involved I think. Could you split it up in 2? Like you already did with the comments

# Determine positions (bin indices) of e_reco bounds in migration array
pos_mig_temp = np.digitize(migra_e_reco, mig_array)-1

# We ensure that no negative values are found
pos_mig = np.maximum(pos_mig_temp , 0)
assert_equal(actual, desired)
# Here we test get_response with an expected gaussian shape for edisp

# First we define a fake EnergyDispersion2D

This comment has been minimized.

@joleroi

joleroi Feb 23, 2017
Contributor

The generation of the fake edisp is very nice! However, it should not happen in the test but be available as method on EnergyDispersion2D, similar to EnergyDispersion.from_gauss . The test should then just call this function.

This comment has been minimized.

@joleroi

joleroi Feb 23, 2017
Contributor

Note that all tests on EnergyDispersion2D could be done with gaussian edisp rather than stored one.

This would be also easier if the method you coded would live on the class

@registerrier
Copy link
Contributor Author

@registerrier registerrier commented Feb 23, 2017

I have added the _from_gauss function in the EnergyDispersion2D class.
I corrected a small issue that makes agreement between analytical solution and get_response better.
Use step of 1e-2 in the test gives result <0.05 in absolute error.

Copy link
Member

@cdeil cdeil left a comment

I left some more inline comments requesting some superficial code changes. I don't have time to review the formulas or method.

@@ -626,6 +626,48 @@ def offset(self):
return self.data.axis('offset')

@classmethod
def _from_gaus(cls, e_true, migra, bias, sigma, offset=None):

This comment has been minimized.

@cdeil

cdeil Feb 23, 2017
Member

The other edisp class in this file has from_gauss, you put _from_gaus.
Suggest to put gauss with 2 s in both.
I don't care if one or both are private or public.
The main point is that it's callable from tests or other places, like @joleroi said.

This comment has been minimized.

@joleroi

joleroi Feb 24, 2017
Contributor

👍 for the public method (like you put it). I can be useful e.g. for simulations

if offset is None:
offset = np.linspace(0.,3.,10)*u.deg

from scipy.special import erf

This comment has been minimized.

@cdeil

cdeil Feb 23, 2017
Member

We always put delayed imports at the very top of methods. Please do it here also for consistency with the rest of Gammapy.

This comment has been minimized.

@registerrier

registerrier Feb 24, 2017
Author Contributor

It is already imported in from_gauss from EnergyDispersion.
I think if we put it on the top, we have to add @requires_dependency('spicy') for EnergyDispersion test as well.
Do we want that?

Bin edges of offset
"""
if offset is None:
offset = np.linspace(0.,3.,10)*u.deg

This comment has been minimized.

@cdeil

cdeil Feb 23, 2017
Member

I know we have defaults in a lot of places that really only work well for HESS, like 3 deg FOV size here.
Probably we should remove those, no?
I.e. the caller must be explicit, but then we have no HESS-specific hardcoded stuff.
You know, be CTA friendly.

This comment has been minimized.

@registerrier

registerrier Feb 24, 2017
Author Contributor

Well, it really does not matter here since it is uniform in offset anyway. No?

@@ -702,7 +744,7 @@ def to_energy_dispersion(self, offset, e_true=None, e_reco=None):
e_reco_lo=ereco_lo, e_reco_hi=ereco_hi,
data=rm)

def get_response(self, offset, e_true, e_reco=None):
def get_response(self, offset, e_true, e_reco=None, mig_step=5e-3):

This comment has been minimized.

@cdeil

cdeil Feb 23, 2017
Member

I think we call it migra pretty consistenly everywhere (also in the spec, and in Gammapy).
So suggest to change mig_step to migra_step here.

This comment has been minimized.

@registerrier

registerrier Feb 24, 2017
Author Contributor

Done

assert_equal(actual, desired)
# Here we test get_response with an expected gaussian shape for edisp
from scipy.special import erf
# First we define a fake EnergyDispersion2D

This comment has been minimized.

@cdeil

cdeil Feb 23, 2017
Member

Suggest to split this whole part creating the test edisp object out into a separate method:

@staticmethod
def make_test_edisp2d():
    ...
    return edisp

Then if one ever needs to test it or if we want to call it from somewhere else, it's possible.
E.g.

from gammapy.irf.tests.test_energy_dispersion import TestEnergyDispersion2D
edisp = TestEnergyDispersion2D.make_test_edisp2d()
# play with it interactively or use it from other tests

This comment has been minimized.

@registerrier

registerrier Feb 24, 2017
Author Contributor

Since we also use mu, sigma, trues and migras array later it is simpler to keep it inside the call and not use an external make function

@cdeil
cdeil approved these changes Feb 24, 2017
@registerrier
Copy link
Contributor Author

@registerrier registerrier commented Feb 24, 2017

Modified test_data_store to add absolute tolerance of 2% to result.

@joleroi
Copy link
Contributor

@joleroi joleroi commented Feb 24, 2017

As far as I can see, we merge this and close #894 ?

@joleroi
Copy link
Contributor

@joleroi joleroi commented Feb 24, 2017

This looks like it was introduced in this PR
https://travis-ci.org/gammapy/gammapy/builds/204946729
I will have a look

@joleroi joleroi merged commit f1f9f55 into gammapy:master Feb 24, 2017
1 of 2 checks passed
1 of 2 checks passed
continuous-integration/travis-ci/pr The Travis CI build could not complete due to an error
Details
continuous-integration/appveyor/pr AppVeyor build succeeded
Details
@cdeil cdeil added the bug label Apr 28, 2017
@cdeil cdeil added this to the 0.6 milestone Apr 28, 2017
@cdeil cdeil changed the title Improved EnergyDispersion2D.get_response and associated tests Improve EnergyDispersion2D get_response and tests Apr 28, 2017
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Linked issues

Successfully merging this pull request may close these issues.

None yet

3 participants