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 Background3D energy axis integration #1894

Merged
merged 4 commits into from Oct 30, 2018

Conversation

Projects
None yet
4 participants
@adonath
Member

adonath commented Oct 29, 2018

This is a change that I started at the coding sprint in Madrid, when debugging #1842. It changes
Background3D.integrate_on_energy_range() to use the _trapz_loglog method for the energy integration. This simplifies the existing code a lot, while improving accuracy. Note how the test values in gammapy/cube/tests/test_make.py now finally make sense. This change should also resolve #1830.

@adonath adonath self-assigned this Oct 29, 2018

@adonath adonath requested review from registerrier, JouvinLea and AtreyeeS Oct 29, 2018

@adonath adonath added the feature label Oct 29, 2018

@adonath adonath added this to the 0.9 milestone Oct 29, 2018

@cdeil

This seems good to me.

@registerrier and @AtreyeeS - I think you want to have oversampling (possibly even adaptive) for the energy binning, so maybe for you this PR partly goes in the wrong direction (the part that reduces the n_integration_bins option to a integrate flag option)?

@AtreyeeS

This comment has been minimized.

Contributor

AtreyeeS commented Oct 29, 2018

Thanks @adonath !
This does indeed resolve #1830.

I checked for varying number of energy bins, and now, the background counts are correctly estimated within fluctuations. Infact, making finer energy bins does not improve the computation, so this PR even obliterates the need for MapMaker to decide on optimal energy bins!

@cdeil

This comment has been minimized.

Member

cdeil commented Oct 29, 2018

@AtreyeeS - The background estimate will still be highly incorrect if the power-law approximation within the bin doesn't hold. So you can't just ask for the background e.g. in the energy bin 0.1 to 1 TeV, where the energy threshold might be at 0.1 TeV and the shape isn't a power-law in the energy bin.
It's a fundamental question at which level the energy grids get set up, and with this PR the assumption is that it's done by the user or in high-level code. I think that is OK, it's also what we currently assume in the sky model evaluation.

fov_lon = map_coord.skycoord.separation(pointing)
fov_lat = Angle(np.zeros_like(fov_lon), fov_lon.unit)
if n_integration_bins == 0:
energy_reco = map_coord[energy_axis.name] * energy_axis.unit
if not integrate:

This comment has been minimized.

@AtreyeeS

AtreyeeS Oct 29, 2018

Contributor

Is there any use case for integrate = False?

This comment has been minimized.

@registerrier

registerrier Oct 30, 2018

Contributor

Good remark.
A priori, background maps are always integrated quantities defined between energy edges and not at bin centers. In that sense, the log-log trapezoidal integral is the most relevant case.
The integrate = False is therefore probably irrelevant here.

This comment has been minimized.

@adonath

adonath Oct 30, 2018

Member

Done.

@registerrier

Maybe remove integrated`option as it is most likely useless here.

fov_lon = map_coord.skycoord.separation(pointing)
fov_lat = Angle(np.zeros_like(fov_lon), fov_lon.unit)
if n_integration_bins == 0:
energy_reco = map_coord[energy_axis.name] * energy_axis.unit
if not integrate:

This comment has been minimized.

@registerrier

registerrier Oct 30, 2018

Contributor

Good remark.
A priori, background maps are always integrated quantities defined between energy edges and not at bin centers. In that sense, the log-log trapezoidal integral is the most relevant case.
The integrate = False is therefore probably irrelevant here.

n_integration_bins=n_integration_bins,
)
energies = energy_axis.edges * energy_axis.unit
fov_lon, fov_lat, energy_reco = np.broadcast_arrays(

This comment has been minimized.

@registerrier

registerrier Oct 30, 2018

Contributor

+1 to move the broadcasting outside evaluate since this change has already been done with evaluate_at_coord.

n_integration_bins : int
Number of bins per energy bin in integration
integrate : bool
Whether to integrate the background model over the energy axes.

This comment has been minimized.

@registerrier

registerrier Oct 30, 2018

Contributor

I am not sure we should remove completely the possibility to have better integration but this is probably not the right place to pass the argument. We have two options:

  • we leave it to the user to have sufficient e_reco bins in the geometry
  • we perform the proper integral using the bkg nodes
    I tend to favor option 2. But in the meantime, let's go with option 1.

@adonath adonath dismissed stale reviews from registerrier and cdeil via d1f045a Oct 30, 2018

@adonath

This comment has been minimized.

Member

adonath commented Oct 30, 2018

@cdeil @AtreyeeS I think this change does not affect the need for the re-sampling API for Maps in general, as the latter is useful in many other contexts. It is just a better way to integrate the background model, which works for large energy bins as well. If users want highest precision, they can still set up a fine binned grid and re-sample the map later, as we have discussed in Madrid.

@@ -70,7 +70,7 @@ def geom(ebounds):
"counts": 34366,
"exposure": 5.971096e11,

This comment has been minimized.

@AtreyeeS

AtreyeeS Oct 30, 2018

Contributor

This is definitely not a part of this PR, but I see that the exposure values also vary a lot depending on the number of (true) energy bins. Perhaps changing from the evaluate() to _trapz_loglog integral might help in this case as well?

@adonath

This comment has been minimized.

Member

adonath commented Oct 30, 2018

Travis fails are unrelated. Merging now.

@adonath adonath merged commit a8ea355 into gammapy:master Oct 30, 2018

3 of 4 checks passed

continuous-integration/travis-ci/pr The Travis CI build could not complete due to an error
Details
Codacy/PR Quality Review Up to standards. A positive pull request.
Details
Scrutinizer Analysis: 2 updated code elements – Tests: passed
Details
continuous-integration/appveyor/pr AppVeyor build succeeded
Details

@adonath adonath deleted the adonath:change_bkg_3d_integration branch Nov 20, 2018

@cdeil cdeil changed the title from Change Background3D integration to _trapz_loglog to Improve Background3D energy axis integration Nov 23, 2018

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment