Discrepancy between simulated counts and model for radial source
|Assigned To:||Knödlseder Jürgen||% Done:|
If I simulate a source with a
GModelSpatialRadialDisk spatial model with IRFs from a run of the H.E.S.S. public data release 1, I experience discrepancies in the number of events between simulated counts and modelled counts.
The simulated number of counts is 91389 and the number of model counts is 85285. I simulate only the source and no background.
The residual map (counts-model) does also not look empty (See attached image).
I also attached a script for reproducing the simulation.
#4 Updated by Knödlseder Jürgen over 3 years ago
Have a look at #2860, I just increased the integration precision for extended source models since it was apparently not sufficient for the broad H.E.S.S. PSF. Could your issue be related to this? New code is in the
devel branch since yesterday evening.
#8 Updated by Knödlseder Jürgen over 3 years ago
- Tracker changed from Support to Bug
- Status changed from New to In Progress
- Target version set to 1.7.0
- % Done changed from 0 to 10
I looked a bit in the problem and played around with the integration precisions in
GCTAResponseIrf::irf_radial. Below a summary of the results, obtained using a modified
simulation.py script where the cube energy range was limited to 0.5-90 TeV (to make sure to stay away from the threshold). I tried a variant where the maximum delta of the PSF was restricted to 99% of the PSF containment, but this led to some artefacts in the model cube at high energies. Setting the number of iterations in rho and phi to 6, as done in #2860, reduces the discrepancy between model and counts to about 3%. 1% is only reached when further increasing the integration precision in phi, but this brings the computation time to 10 min instead of 1 min when compared to the initial precision. In other words, the speed penalty would be a factor of 10, while with 6/6 it would be a factor of 4, which is already quite bad.
|rho||phi||Nobs||Nmodel||ctmodel time||delta_max||model cube aspect|
|5||5||84169||78538 (-6.7%)||1m12sec||99% of PSF||Artefacts in the model cube at high energies|
|6||6||84169||82332.2 (-2.2%)||4m09sec||99% of PSF||Artefacts in the model cube at high energies|
|6||6||84132||81336 (-3.3%)||5m03sec||Full range||Still faint ring-like structures at high energies|
|7||6||84132||81355.4 (-3.3%)||9m44sec||Full range||Model looks now very smooth, extremely little radial structure|
|6||7||84132||83359.9 (-0.9%)||9m37sec||Full range||Still faint ring-like structures at high energies|
|7||7||84132||83381.6 (-0.9%)||19m09sec||Full range||Almost perfectly smooth|
|5||7||84132||83625.4 (-0.6%)||4m53sec||Full range||Strong ring-like structures over all energies|
I prefer for the time being (and release 1.6) to fix the precision at 6/6, and defer it to the next release to think about a more precise but fast integration scheme. I keep this issue open for that purpose.
#10 Updated by Knödlseder Jürgen over 3 years ago
Tiziani Domenico wrote:
I agree, the computation time of ctmodel should be kept short.
Just to keep in track of it: We should also not forget about
GCTAResponseIrf::irf_elliptical, where I see similar problems when going to large semi-axes.
That’s the method where I discovered the issue. I changed the precision in both methods.