search

gammapy / gammapy-benchmarks

Performance benchmarks for Gammapy
BSD 3-Clause "New" or "Revised" License
4 stars 14 forks source link

Add check on small source extension fit #106

Closed registerrier closed 2 years ago

registerrier commented 3 years ago

This PR adds a new validation in the checks directory.

It performs a 3D fit of the PKS 2155 data during flare from the HESS DL3 DR1. The objective is to compare the fit behavior between pointlike and gaussian source component.

With current dev version, it shows that the two fit yield significantly different results in terms of flux.

adonath commented 3 years ago

Thanks a lot @registerrier! From a quick test to run the script locally I received the following error:

Traceback (most recent call last):
  File "make.py", line 165, in <module>
    cli()
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 1137, in __call__
    return self.main(*args, **kwargs)
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 1062, in main
    rv = self.invoke(ctx)
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 1668, in invoke
    return _process_result(sub_ctx.command.invoke(sub_ctx))
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 1404, in invoke
    return ctx.invoke(self.callback, **ctx.params)
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 763, in invoke
    return __callback(*args, **kwargs)
  File "make.py", line 68, in run_analyses
    conf_result = fit.confidence([dataset], "sigma", 3)
  File "/Users/adonath/github/adonath/gammapy/gammapy/modeling/fit.py", line 304, in confidence
    **kwargs,
  File "/Users/adonath/github/adonath/gammapy/gammapy/modeling/iminuit.py", line 132, in confidence_iminuit
    result = minuit.minos(var=var, sigma=sigma, maxcall=maxcall)
  File "src/iminuit/_libiminuit.pyx", line 1017, in iminuit._libiminuit.Minuit.minos
RuntimeError: Function minimum is not valid. Make sure MIGRAD converged first

Can you reproduce this locally? Or is this basically expected behaviour, because the Gaussian fit is expected to fail? Here is my environment for testing:

System:

    python_executable      : /Users/adonath/opt/anaconda3/envs/gammapy-dev/bin/python3.7 
    python_version         : 3.7.0      
    machine                : x86_64     
    system                 : Darwin     

Gammapy package:

    version                : 0.18.3.dev2034+ga6af2c26b 
    path                   : /Users/adonath/github/adonath/gammapy/gammapy 

Other packages:

    numpy                  : 1.21.1     
    scipy                  : 1.5.3      
    astropy                : 4.3        
    regions                : 0.5        
    click                  : 8.0.1      
    yaml                   : 5.4.1      
    IPython                : 7.24.1     
    jupyterlab             : 3.0.16     
    matplotlib             : 3.4.2      
    pandas                 : 1.2.4      
    healpy                 : 1.13.0     
    iminuit                : 1.5.1      
    sherpa                 : 4.13.1     
    naima                  : not installed 
    emcee                  : 2.2.1      
    corner                 : 2.1.0      

Gammapy environment variables:

    GAMMAPY_DATA           : /Users/adonath/github/gammapy/gammapy-data
registerrier commented 3 years ago

Thanks a lot @registerrier! From a quick test to run the script locally I received the following error:

Traceback (most recent call last):
  File "make.py", line 165, in <module>
    cli()
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 1137, in __call__
    return self.main(*args, **kwargs)
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 1062, in main
    rv = self.invoke(ctx)
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 1668, in invoke
    return _process_result(sub_ctx.command.invoke(sub_ctx))
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 1404, in invoke
    return ctx.invoke(self.callback, **ctx.params)
  File "/Users/adonath/opt/anaconda3/envs/gammapy-dev/lib/python3.7/site-packages/click/core.py", line 763, in invoke
    return __callback(*args, **kwargs)
  File "make.py", line 68, in run_analyses
    conf_result = fit.confidence([dataset], "sigma", 3)
  File "/Users/adonath/github/adonath/gammapy/gammapy/modeling/fit.py", line 304, in confidence
    **kwargs,
  File "/Users/adonath/github/adonath/gammapy/gammapy/modeling/iminuit.py", line 132, in confidence_iminuit
    result = minuit.minos(var=var, sigma=sigma, maxcall=maxcall)
  File "src/iminuit/_libiminuit.pyx", line 1017, in iminuit._libiminuit.Minuit.minos
RuntimeError: Function minimum is not valid. Make sure MIGRAD converged first

Can you reproduce this locally? Or is this basically expected behaviour, because the Gaussian fit is expected to fail? Here is my environment for testing:

It does not fail for me but with many warning messages, I suppose the re-optimization is causing problems. The obtained 3 sigma interval is meaningless (err ~3e-5 deg). What minimum do you have? I my case optimize converges to sigma ~2e-3 deg but the flux becomes ~3 times larger than for the point source.

I am not up to date on everything:

System:

    python_executable      : /Users/terrier/Code/anaconda3/envs/gammapy-dev/bin/python 
    python_version         : 3.7.0      
    machine                : x86_64     
    system                 : Darwin     

Gammapy package:

    version                : 0.18.3.dev1261+g3691edad3 
    path                   : /Users/terrier/Code/gammapy-dev/gammapy/gammapy 

Other packages:

    numpy                  : 1.19.1     
    scipy                  : 1.5.2      
    astropy                : 4.0.1.post1 
    regions                : 0.4        
    click                  : 7.1.2      
    yaml                   : 5.3.1      
    IPython                : 7.18.1     
    jupyterlab             : 2.2.8      
    matplotlib             : 3.3.4      
    pandas                 : 1.1.3      
    healpy                 : 1.13.0     
    iminuit                : 1.5.2      
    sherpa                 : 4.12.1     
    naima                  : 0.9.1      
    emcee                  : 2.2.1      
    corner                 : 2.1.0      

Gammapy environment variables:

    GAMMAPY_DATA           : /Users/terrier/Code/gammapy-dev/gammapy-datasets/
adonath commented 3 years ago

@registerrier Updating to iminuit==1.5.2 fixed the issue. However the optimiser tries the "wildest" positions:

WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (305.11332314, 20.52171332)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (305.11332314, 20.52171332)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (319.79523339, -62.56772905)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (319.79523339, -62.56772905)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (135.42739612, 87.0571366)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (135.42739612, 87.0571366)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (47.8200215, -62.29209085)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (47.8200215, -62.29209085)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (184.73428323, -87.63581342)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (184.73428323, -87.63581342)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (73.47988294, -66.0717863)> is outside valid IRF map range, using nearest IRF defined within
WARNING:gammapy.irf.core:Position <SkyCoord (ICRS): (ra, dec) in deg
    (73.47988294, -66.0717863)> is outside valid IRF map range, using nearest IRF defined within
...

For this study it's probably fine to just fix the position during the fit...

registerrier commented 2 years ago

Merging this now. Further changes will come in follow up PRs.