Bug #2954

Correct handling of flux for cube sky models in csspec

Added by Martin Pierrick almost 5 years ago. Updated almost 5 years ago.

Status:ClosedStart date:07/11/2019
Priority:NormalDue date:
Assigned To:Martin Pierrick% Done:

100%

Category:-Estimated time:1.00 hour
Target version:1.6.2
Duration:

Description

When used with a cube sky model (GModelSpatialDiffuseCube), csspec currently returns normalization best-fit values and errors (multiplied by bin energy in MeV and bin energy in erg).

We would like to have physical fluxes instead, obtained from mutliplying the normalization by the cube flux integrated over that bin. Note: the upper limit computed by csspec are correct so everything is orrectly implemented in ctulimit. I attach an example csspec log file.

The code should include these bits:
  1. Make sure total cube spectrum is computed and cached
    if model.spatial().type() 'DiffuseMapCube’:
    model.spatial().set_mc_cone(GSkyDir(),180.0)
  2. Multiply normalization by flux in cube at this energy
    if model.spatial().type() 'DiffuseMapCube’:
    spectrum[i] *= model.spatial().spectrum().eval(energy[i],GTime())

csspec_binned_1-10TeV_LMC-Pion.log (6.08 KB) Martin Pierrick, 07/11/2019 03:45 PM

models.xml Magnifier (1.05 KB) Knödlseder Jürgen, 07/16/2019 10:10 AM

cube.fits (64.7 KB) Knödlseder Jürgen, 07/16/2019 10:10 AM

Recurrence

No recurrence.

History

#1 Updated by Knödlseder Jürgen almost 5 years ago

  • Target version set to 1.6.2

#2 Updated by Knödlseder Jürgen almost 5 years ago

I setup a test run using models.xml and cube.fits files. A 30 min simulation for a pointing at the Crab position gave the following csspec result:

2019-07-16T08:08:06: +===================+
2019-07-16T08:08:06: | Generate spectrum |
2019-07-16T08:08:06: +===================+
2019-07-16T08:08:06: === GEbounds ===
2019-07-16T08:08:06:  Number of intervals .......: 10
2019-07-16T08:08:06:  Energy range ..............: 100 GeV - 100 TeV
2019-07-16T08:08:46:  Bin 1 .....................: 1.602177e-10 +/- 4.421024e-12 erg/cm2/s
2019-07-16T08:08:46:  Bin 2 .....................: 2.742546e-08 +/- 7.623746e-09 [< 1.259704e-11] erg/cm2/s (TS = 16.325)
2019-07-16T08:08:46:  Bin 3 .....................: 2.062850e-07 +/- 1.437324e-08 [< 1.441020e-10] erg/cm2/s (TS = 609.284)
2019-07-16T08:08:46:  Bin 4 .....................: 1.841667e-07 +/- 1.140676e-08 [< 1.822600e-10] erg/cm2/s (TS = 1310.556)
2019-07-16T08:08:46:  Bin 5 .....................: 1.902318e-07 +/- 1.443690e-08 [< 1.410683e-10] erg/cm2/s (TS = 554.006)
2019-07-16T08:08:46:  Bin 6 .....................: 1.715272e-07 +/- 1.714087e-08 [< 9.521660e-11] erg/cm2/s (TS = 206.144)
2019-07-16T08:08:46:  Bin 7 .....................: 1.807496e-07 +/- 2.406497e-08 [< 8.060233e-11] erg/cm2/s (TS = 103.523)
2019-07-16T08:08:46:  Bin 8 .....................: 8.352029e-08 +/- 2.119759e-08 [< 3.329548e-11] erg/cm2/s (TS = 31.199)
2019-07-16T08:08:46:  Bin 9 .....................: 5.385808e-08 +/- 1.608778e-08 [< 1.990545e-11] erg/cm2/s (TS = 25.012)
2019-07-16T08:08:46:  Bin 10 ....................: 5.647644e-09 +/- 6.757822e-03 [< 9.737229e-12] erg/cm2/s (TS = 57.810)

#3 Updated by Knödlseder Jürgen almost 5 years ago

Adding the following code

# If the source model is a cube then multiply-in the cube
# spectrum
if source.spatial().classname() == 'GModelSpatialDiffuseCube':
    dir          = gammalib.GSkyDir()
    source.spatial().set_mc_cone(dir, 180.0)
    norm         = source.spatial().spectrum().eval(elogmean)
    fitted_flux *= norm
    e_flux      *= norm

gives the following result
2019-07-16T08:17:59: +===================+
2019-07-16T08:17:59: | Generate spectrum |
2019-07-16T08:17:59: +===================+
2019-07-16T08:17:59: === GEbounds ===
2019-07-16T08:17:59:  Number of intervals .......: 10
2019-07-16T08:17:59:  Energy range ..............: 100 GeV - 100 TeV
2019-07-16T08:18:38:  Bin 1 .....................: 3.078758e-14 +/- 8.495484e-16 erg/cm2/s
2019-07-16T08:18:38:  Bin 2 .....................: 5.328896e-12 +/- 1.481330e-12 [< 1.259704e-11] erg/cm2/s (TS = 16.325)
2019-07-16T08:18:38:  Bin 3 .....................: 4.045340e-11 +/- 2.818657e-12 [< 1.441020e-10] erg/cm2/s (TS = 609.284)
2019-07-16T08:18:38:  Bin 4 .....................: 3.611729e-11 +/- 2.237001e-12 [< 1.822600e-10] erg/cm2/s (TS = 1310.556)
2019-07-16T08:18:38:  Bin 5 .....................: 3.726838e-11 +/- 2.828337e-12 [< 1.410683e-10] erg/cm2/s (TS = 554.006)
2019-07-16T08:18:38:  Bin 6 .....................: 3.206621e-11 +/- 3.204407e-12 [< 9.521660e-11] erg/cm2/s (TS = 206.144)
2019-07-16T08:18:38:  Bin 7 .....................: 2.845845e-11 +/- 3.788953e-12 [< 8.060233e-11] erg/cm2/s (TS = 103.523)
2019-07-16T08:18:38:  Bin 8 .....................: 1.087201e-11 +/- 2.759333e-12 [< 3.329548e-11] erg/cm2/s (TS = 31.199)
2019-07-16T08:18:38:  Bin 9 .....................: 5.796319e-12 +/- 1.731400e-12 [< 1.990545e-11] erg/cm2/s (TS = 25.012)
2019-07-16T08:18:38:  Bin 10 ....................: 5.025189e-13 +/- 6.013009e-07 [< 9.737229e-12] erg/cm2/s (TS = 57.810)

#4 Updated by Knödlseder Jürgen almost 5 years ago

  • Status changed from New to Closed
  • % Done changed from 0 to 100

Merged into the 1.6.2 release.

Also available in: Atom PDF