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

Migrate SpectrumAnalysisIACT to the high-level interface #2330

Merged
merged 15 commits into from Sep 2, 2019

Conversation

@Bultako
Copy link
Member

commented Aug 28, 2019

This PR aims to add all features needed in the high-level interface to perform a 1D spectrum analysis.

In [1]: import yaml                                                                                                                                                                   

In [2]: from gammapy.scripts import Analysis                                                                                                                                          

In [3]: config = """
model:
    components:
        - name: Crab
          spectral:
            type: PowerLaw
            parameters:
            - name: index
              value: 2.6
              factor: 2.2
              scale: 1.0
              unit: ''
              min: .nan
              max: .nan
              frozen: false
            - name: amplitude
              value: 5.0e-11
              factor: 2.3
              scale: 1.0e-12
              unit: cm-2 s-1 TeV-1
              min: .nan
              max: .nan
              frozen: false
            - name: reference
              value: 1.0
              factor: 1.0
              scale: 1.0
              unit: TeV
              min: .nan
              max: .nan
              frozen: true
flux:
    fp_binning:
        lo_bnd: 1
        hi_bnd: 10
        nbin: 11
        unit: TeV
        interp: log
"""                                                                                                                                                                       

In [4]: config = yaml.load(config)                                                                                                                                                        
/Users/jer/anaconda/envs/gammapy-dev/bin/ipython:1: YAMLLoadWarning: calling yaml.load() without Loader=... is deprecated, as the default Loader is unsafe. Please read https://msg.pyyaml.org/load for full details.
  #!/Users/jer/anaconda/envs/gammapy-dev/bin/python

In [5]: a = Analysis(config)   
INFO:gammapy.scripts.analysis:Setting logging parameters (INFO).

In [6]: print(a.config)                                                                                                                                                               
fit: {}
flux:
  fp_binning:
    hi_bnd: 10
    interp: log
    lo_bnd: 1
    nbin: 11
    unit: TeV
general:
  logging:
    level: INFO
  outdir: .
geometry:
  axes: {}
model:
  components:
  - name: Crab
    spectral:
      parameters:
      - factor: 2.2
        frozen: false
        max: .nan
        min: .nan
        name: index
        scale: 1.0
        unit: ''
        value: 2.6
      - factor: 2.3
        frozen: false
        max: .nan
        min: .nan
        name: amplitude
        scale: 1.0e-12
        unit: cm-2 s-1 TeV-1
        value: 5.0e-11
      - factor: 1.0
        frozen: true
        max: .nan
        min: .nan
        name: reference
        scale: 1.0
        unit: TeV
        value: 1.0
      type: PowerLaw
observations:
  datastore: $GAMMAPY_DATA/hess-dl3-dr1
  filters:
  - filter_type: par_value
    value_param: Crab
    variable: TARGET_NAME
reduction:
  background:
    background_estimator: reflected
    on_region:
      center:
      - 83.633 deg
      - 22.014 deg
      frame: icrs
      radius: 0.1 deg
  containment_correction: true
  data_reducer: 1d


In [7]: a.get_observations()                                                                                                                                                          
INFO:gammapy.scripts.analysis:Fetching observations.
INFO:gammapy.scripts.analysis:Info for OBS_ID = 23592
- Start time: 53347.91
- Pointing pos: RA 82.01 deg / Dec 22.01 deg
- Observation duration: 1686.0 s
- Dead-time fraction: 6.212 %

INFO:gammapy.scripts.analysis:Info for OBS_ID = 23523
- Start time: 53343.92
- Pointing pos: RA 83.63 deg / Dec 21.51 deg
- Observation duration: 1687.0 s
- Dead-time fraction: 6.240 %

INFO:gammapy.scripts.analysis:Info for OBS_ID = 23526
- Start time: 53343.95
- Pointing pos: RA 83.63 deg / Dec 22.51 deg
- Observation duration: 1683.0 s
- Dead-time fraction: 6.555 %

INFO:gammapy.scripts.analysis:Info for OBS_ID = 23559
- Start time: 53345.96
- Pointing pos: RA 85.25 deg / Dec 22.01 deg
- Observation duration: 1686.0 s
- Dead-time fraction: 6.398 %


In [8]: a.reduce()                                                                                                                                                                    
INFO:gammapy.scripts.analysis:Reducing data sets.
INFO:gammapy.spectrum.reflected:Found 46 reflected regions for the Obs #23592
INFO:gammapy.spectrum.reflected:Found 14 reflected regions for the Obs #23523
INFO:gammapy.spectrum.reflected:Found 14 reflected regions for the Obs #23526
INFO:gammapy.spectrum.reflected:Found 46 reflected regions for the Obs #23559
INFO:gammapy.spectrum.extract:Running <gammapy.spectrum.extract.SpectrumExtraction object at 0x10332bcc0>
INFO:gammapy.spectrum.extract:Process observation
 Info for OBS_ID = 23592
- Start time: 53347.91
- Pointing pos: RA 82.01 deg / Dec 22.01 deg
- Observation duration: 1686.0 s
- Dead-time fraction: 6.212 %

INFO:gammapy.spectrum.extract:Update observation meta info
INFO:gammapy.spectrum.extract:Offset : 1.501568530069493 deg

INFO:gammapy.spectrum.extract:Fill events
INFO:gammapy.spectrum.extract:Extract IRFs
INFO:gammapy.spectrum.extract:Apply containment correction
/Users/jer/git/gammapy/gammapy/gammapy/spectrum/extract.py:232: RuntimeWarning: invalid value encountered in true_divide
  self.containment = new_aeff.data.data.value / self._aeff.data.data.value
INFO:gammapy.spectrum.extract:Process observation
 Info for OBS_ID = 23523
- Start time: 53343.92
- Pointing pos: RA 83.63 deg / Dec 21.51 deg
- Observation duration: 1687.0 s
- Dead-time fraction: 6.240 %

INFO:gammapy.spectrum.extract:Update observation meta info
INFO:gammapy.spectrum.extract:Offset : 0.4995557435558715 deg

INFO:gammapy.spectrum.extract:Fill events
INFO:gammapy.spectrum.extract:Extract IRFs
INFO:gammapy.spectrum.extract:Apply containment correction
INFO:gammapy.spectrum.extract:Process observation
 Info for OBS_ID = 23526
- Start time: 53343.95
- Pointing pos: RA 83.63 deg / Dec 22.51 deg
- Observation duration: 1683.0 s
- Dead-time fraction: 6.555 %

INFO:gammapy.spectrum.extract:Update observation meta info
INFO:gammapy.spectrum.extract:Offset : 0.5004444451150709 deg

INFO:gammapy.spectrum.extract:Fill events
INFO:gammapy.spectrum.extract:Extract IRFs
INFO:gammapy.spectrum.extract:Apply containment correction
INFO:gammapy.spectrum.extract:Process observation
 Info for OBS_ID = 23559
- Start time: 53345.96
- Pointing pos: RA 85.25 deg / Dec 22.01 deg
- Observation duration: 1686.0 s
- Dead-time fraction: 6.398 %

INFO:gammapy.spectrum.extract:Update observation meta info
INFO:gammapy.spectrum.extract:Offset : 1.5021898826765052 deg

INFO:gammapy.spectrum.extract:Fill events
INFO:gammapy.spectrum.extract:Extract IRFs
INFO:gammapy.spectrum.extract:Apply containment correction

In [9]: a.fit()                                                                                                                                                                       
INFO:gammapy.scripts.analysis:Reading model.
INFO:gammapy.scripts.analysis:PowerLaw

Parameters: 

	   name     value   error      unit      min max frozen
	--------- --------- ----- -------------- --- --- ------
	    index 2.600e+00   nan                nan nan  False
	amplitude 5.000e-11   nan cm-2 s-1 TeV-1 nan nan  False
	reference 1.000e+00   nan            TeV nan nan   True
INFO:gammapy.scripts.analysis:Fitting data sets to model.
INFO:gammapy.scripts.analysis:OptimizeResult

	backend    : minuit
	method     : minuit
	success    : True
	message    : Optimization terminated successfully.
	nfev       : 36
	total stat : 130.56


In [10]: a.get_flux_points()                                                                                                                                                           
INFO:gammapy.scripts.analysis:Calculating flux points.
INFO:gammapy.scripts.analysis:
      e_ref               ref_flux                 dnde                 dnde_ul                dnde_err        is_ul
       TeV              1 / (cm2 s)          1 / (cm2 s TeV)        1 / (cm2 s TeV)        1 / (cm2 s TeV)          
------------------ ---------------------- ---------------------- ---------------------- ---------------------- -----
1.1364636663857248  9.449329992294733e-12  2.816977541295253e-11  3.573691483792928e-11  3.505225434375799e-12 False
1.3768571648527583 3.4113700786811144e-12  1.828679888324125e-11  2.603817913822835e-11 3.4450407931548733e-12 False
1.6681005372000581  4.982516507677505e-12  1.365332939913724e-11  1.755910853352502e-11  1.799039326424294e-12 False
2.1544346900318834 3.2519752617892428e-12 6.1484585038692114e-12  8.309695131203611e-12   9.78342227704282e-13 False
2.6101572156825363 1.1740188049020575e-12  3.861118980517589e-12  6.052405645173305e-12  9.421227414399863e-13 False
3.1622776601683777 1.7147269105466929e-12 2.5897875508619773e-12  3.649825162635732e-12  4.731485364864076e-13 False
 3.831186849557287   6.19045803302273e-13  1.156618808270176e-12 2.1238966040340406e-12  3.955211988710733e-13 False
4.6415888336127775  9.041545956088386e-13  7.045665488041625e-13 1.1577572591321411e-12 1.9297021903699566e-13 False
 5.994842503189405  5.901211512741317e-13 3.6988307657563386e-13  6.485911472662804e-13 1.1511026416367203e-13 False
  7.26291750173621 2.1304384965865337e-13  4.204034950709673e-13  8.073308222764743e-13 1.5501752336572198e-13 False
  8.79922543569107 3.1116368887006885e-13 1.1711100306620158e-13  2.544714864537633e-13  5.304630090176535e-14 False

In [11]: a.flux_points_dataset.peek()                                                                                                                                                 
Out[11]: 
(<matplotlib.axes._subplots.AxesSubplot at 0x1111aec18>,
 <matplotlib.axes._subplots.AxesSubplot at 0x111b709b0>)
@cdeil cdeil self-assigned this Aug 29, 2019
@cdeil cdeil added the feature label Aug 29, 2019
@cdeil cdeil added this to the 0.14 milestone Aug 29, 2019
Copy link
Member

left a comment

@Bultako - Thanks!

Some comments inline.

You're adding a lot of things: config schema extension, observation handling, spectral analysis, WCS for 3D analysis.

I'm OK to add in any order, as long as there's a short-term merge point in sight and it doesn't become a long-running PR with 1000+ lines that I have to re-review many times for weeks.

One (very ambitious) goal could be to migrate the existing https://docs.gammapy.org/0.13/api/gammapy.scripts.SpectrumAnalysisIACT.html test case and get it to run with the new scheme (or even to delete the old class and test, and to update the docs and tutorials. If you don't want to do that much, I would suggest to finish up as quickly as possible (see inline comments) and then to merge.

gammapy/scripts/analysis.py Outdated Show resolved Hide resolved
gammapy/scripts/analysis.py Show resolved Hide resolved
gammapy/scripts/analysis.py Outdated Show resolved Hide resolved
gammapy/scripts/analysis.py Outdated Show resolved Hide resolved
gammapy/scripts/analysis.py Outdated Show resolved Hide resolved
gammapy/scripts/tests/test_analysis.py Outdated Show resolved Hide resolved
gammapy/scripts/config/schema.yaml Outdated Show resolved Hide resolved
@Bultako Bultako force-pushed the Bultako:HLI-1D branch from 0107db7 to 910711e Aug 29, 2019
@Bultako

This comment has been minimized.

Copy link
Member Author

commented Aug 29, 2019

@cdeil

Ok to merge this PR soon and avoid having 1000+ lines.
I will add the features to do the fitting and then I think we would have cover SpectrumAnalysisIACT use case. I will leave the documentation and notebooks for another PR.

I've taken some of your suggestions and commented on the rest.

@Bultako Bultako force-pushed the Bultako:HLI-1D branch 2 times, most recently from 6902206 to 4590fd5 Aug 29, 2019
@Bultako Bultako changed the title Add high-level features to address 1D spectrum use case Migrate SpectrumAnalysisIACT to the high-level interface Aug 29, 2019
Copy link
Member

left a comment

Left some inline comments.
The most important one is to find a good solution to auto-create some types like Quantity, Angle, ... from the config file.
Looking at how they implemented it might be a useful example:
https://github.com/tardis-sn/tardis-setups/blob/4718593b8041aefd27e29e64919448e33353669f/2014/2014_kerzendorf_sim/appendix_A1/tardis_example.yml#L16

gammapy/scripts/analysis.py Show resolved Hide resolved
gammapy/scripts/analysis.py Outdated Show resolved Hide resolved
Copy link
Member

left a comment

Some more minor comments inline.

You can see how Tardis has quantity in the schema here and how they process it here.
Maybe we should do it the same way?

gammapy/scripts/config/schema.yaml Show resolved Hide resolved
gammapy/scripts/analysis.py Outdated Show resolved Hide resolved
@Bultako Bultako force-pushed the Bultako:HLI-1D branch 3 times, most recently from d26e574 to 7d10823 Aug 30, 2019
@Bultako Bultako force-pushed the Bultako:HLI-1D branch 2 times, most recently from 5d35a51 to 8ba510d Aug 31, 2019
@Bultako Bultako force-pushed the Bultako:HLI-1D branch from 8ba510d to c96f78d Aug 31, 2019
@Bultako Bultako force-pushed the Bultako:HLI-1D branch from 0028e10 to 13233b8 Aug 31, 2019
@Bultako Bultako force-pushed the Bultako:HLI-1D branch from 13233b8 to 556f84e Aug 31, 2019
@Bultako Bultako force-pushed the Bultako:HLI-1D branch from 5620b08 to 560695b Sep 1, 2019
@Bultako

This comment has been minimized.

Copy link
Member Author

commented Sep 2, 2019

@cdeil
This PR is ready for final review before merging.

@Bultako Bultako force-pushed the Bultako:HLI-1D branch from 560695b to b43268c Sep 2, 2019
@cdeil
cdeil approved these changes Sep 2, 2019
Copy link
Member

left a comment

@Bultako - Thanks!

I'm merging this in now.

@cdeil cdeil merged commit 126854a into gammapy:master Sep 2, 2019
6 of 9 checks passed
6 of 9 checks passed
Codacy/PR Quality Review Hang in there, Codacy is reviewing your Pull request.
Details
Scrutinizer Running
Details
continuous-integration/travis-ci/pr The Travis CI build is in progress
Details
gammapy.gammapy Build #20190902.2 succeeded
Details
gammapy.gammapy (DevDocs) DevDocs succeeded
Details
gammapy.gammapy (Lint) Lint succeeded
Details
gammapy.gammapy (Test Python36) Test Python36 succeeded
Details
gammapy.gammapy (Test Windows36) Test Windows36 succeeded
Details
gammapy.gammapy (Test Windows37) Test Windows37 succeeded
Details
@Bultako Bultako deleted the Bultako:HLI-1D branch Sep 2, 2019
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
2 participants
You can’t perform that action at this time.