Bug #1722

Prefactor offset for diffuse cube simulation

Added by Knödlseder Jürgen almost 9 years ago. Updated almost 9 years ago.

Status:ClosedStart date:03/01/2016
Priority:NormalDue date:
Assigned To:Knödlseder Jürgen% Done:

90%

Category:-
Target version:1.1.0
Duration:

Description

With 10 hours of simulation of a pi0 model, there is an offset in the fitted pi0 map cube normalisation:

2016-03-01T17:21:20: === GModelSky ===
2016-03-01T17:21:20:  Name ......................: pi0
2016-03-01T17:21:20:  Instruments ...............: all
2016-03-01T17:21:20:  Instrument scale factors ..: unity
2016-03-01T17:21:20:  Observation identifiers ...: all
2016-03-01T17:21:20:  Model type ................: DiffuseSource
2016-03-01T17:21:20:  Model components ..........: "MapCubeFunction" * "ConstantValue" * "Constant" 
2016-03-01T17:21:20:  Number of parameters ......: 3
2016-03-01T17:21:20:  Number of spatial par's ...: 1
2016-03-01T17:21:20:   Normalization ............: 1 [0.1,10]  (fixed,scale=1,gradient)
2016-03-01T17:21:20:  Number of spectral par's ..: 1
2016-03-01T17:21:20:   Value ....................: 0.914357 +/- 0.0300576 [0,10]  (free,scale=1,gradient)
2016-03-01T17:21:20:  Number of temporal par's ..: 1
2016-03-01T17:21:20:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)

Recurrence

No recurrence.

History

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

  • Status changed from New to In Progress
  • % Done changed from 0 to 10

Changing the log-log interpolation to log-lin (which appeared to me closer to the MC simulation) makes the fit even worse:

2016-03-01T17:41:35: === GModelSky ===
2016-03-01T17:41:35:  Name ......................: pi0
2016-03-01T17:41:35:  Instruments ...............: all
2016-03-01T17:41:35:  Instrument scale factors ..: unity
2016-03-01T17:41:35:  Observation identifiers ...: all
2016-03-01T17:41:35:  Model type ................: DiffuseSource
2016-03-01T17:41:35:  Model components ..........: "MapCubeFunction" * "ConstantValue" * "Constant" 
2016-03-01T17:41:35:  Number of parameters ......: 3
2016-03-01T17:41:35:  Number of spatial par's ...: 1
2016-03-01T17:41:35:   Normalization ............: 1 [0.1,10]  (fixed,scale=1,gradient)
2016-03-01T17:41:35:  Number of spectral par's ..: 1
2016-03-01T17:41:35:   Value ....................: 0.634118 +/- 0.0374409 [0,10]  (free,scale=1,gradient)
2016-03-01T17:41:35:  Number of temporal par's ..: 1
2016-03-01T17:41:35:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
2016-03-01T17:41:35: === GCTAModelIrfBackground ===
2016-03-01T17:41:35:  Name ......................: CTABackgroundModel
2016-03-01T17:41:35:  Instruments ...............: CTA
2016-03-01T17:41:35:  Instrument scale factors ..: unity
2016-03-01T17:41:35:  Observation identifiers ...: all
2016-03-01T17:41:35:  Model type ................: "PowerLaw" * "Constant" 
2016-03-01T17:41:35:  Number of parameters ......: 4
2016-03-01T17:41:35:  Number of spectral par's ..: 3
2016-03-01T17:41:35:   Prefactor ................: 1.0444 +/- 0.00424553 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
2016-03-01T17:41:35:   Index ....................: 0.0133982 +/- 0.00131418 [-5,5]  (free,scale=1,gradient)
2016-03-01T17:41:35:   PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
2016-03-01T17:41:35:  Number of temporal par's ..: 1
2016-03-01T17:41:35:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)

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

  • Assigned To set to Knödlseder Jürgen
  • % Done changed from 10 to 50

I investigated the GModelSpatialDiffuseCube class further and specifically replaced the former MC sampling using a cache array by a rejection method that is even faster than the previous method (that was always simulating the entire map cube, not only the simulation cone).

Now, only sky directions within the simulation cone are simulated, and the rejection method uses the same formulae to derive the sky map values as the eval() method. Hence the results should be 100% compliant.

This even works nicely when comparing the source model only, the resulting model scaling factor is basically one. However, once the background model is included, the same 10% underestimation is found. So the problem is possibly linked to the simulation of the background component and not the source component.

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

Here the model fitting result of a map cube for 100 hours. The fitting factor is very close to 1.0.

2016-03-02T14:27:29: === GModelSky ===
2016-03-02T14:27:29:  Name ......................: pi0
2016-03-02T14:27:29:  Instruments ...............: all
2016-03-02T14:27:29:  Instrument scale factors ..: unity
2016-03-02T14:27:29:  Observation identifiers ...: all
2016-03-02T14:27:29:  Model type ................: DiffuseSource
2016-03-02T14:27:29:  Model components ..........: "MapCubeFunction" * "ConstantValue" * "Constant" 
2016-03-02T14:27:29:  Number of parameters ......: 3
2016-03-02T14:27:29:  Number of spatial par's ...: 1
2016-03-02T14:27:29:   Normalization ............: 1 [0.1,10]  (fixed,scale=1,gradient)
2016-03-02T14:27:29:  Number of spectral par's ..: 1
2016-03-02T14:27:29:   Value ....................: 0.99671 +/- 0.00214087 [0,10]  (free,scale=1,gradient)
2016-03-02T14:27:29:  Number of temporal par's ..: 1
2016-03-02T14:27:29:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)

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

  • Target version set to 1.1.0
  • % Done changed from 50 to 90

I now implemented a rejection method also in GCTABackground3D::mc() so that the Monte Carlo simulations are identical to the model evaluation. The fits results are not prefect, but the remaining offset from 1 could be simply due to statistical fluctuations.

2016-03-03T09:14:56: === GModelSky ===
2016-03-03T09:14:56:  Name ......................: pi0
2016-03-03T09:14:56:  Instruments ...............: all
2016-03-03T09:14:56:  Instrument scale factors ..: unity
2016-03-03T09:14:56:  Observation identifiers ...: all
2016-03-03T09:14:56:  Model type ................: DiffuseSource
2016-03-03T09:14:56:  Model components ..........: "MapCubeFunction" * "ConstantValue" * "Constant" 
2016-03-03T09:14:56:  Number of parameters ......: 3
2016-03-03T09:14:56:  Number of spatial par's ...: 1
2016-03-03T09:14:56:   Normalization ............: 1 [0.1,10]  (fixed,scale=1,gradient)
2016-03-03T09:14:56:  Number of spectral par's ..: 1
2016-03-03T09:14:56:   Value ....................: 0.931227 +/- 0.0300119 [0,10]  (free,scale=1,gradient)
2016-03-03T09:14:56:  Number of temporal par's ..: 1
2016-03-03T09:14:56:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
2016-03-03T09:14:56: === GCTAModelIrfBackground ===
2016-03-03T09:14:56:  Name ......................: CTABackgroundModel
2016-03-03T09:14:56:  Instruments ...............: CTA
2016-03-03T09:14:56:  Instrument scale factors ..: unity
2016-03-03T09:14:56:  Observation identifiers ...: all
2016-03-03T09:14:56:  Model type ................: "PowerLaw" * "Constant" 
2016-03-03T09:14:56:  Number of parameters ......: 4
2016-03-03T09:14:56:  Number of spectral par's ..: 3
2016-03-03T09:14:56:   Prefactor ................: 1.00692 +/- 0.00412487 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
2016-03-03T09:14:56:   Index ....................: 0.00147384 +/- 0.00134355 [-5,5]  (free,scale=1,gradient)
2016-03-03T09:14:56:   PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
2016-03-03T09:14:56:  Number of temporal par's ..: 1
2016-03-03T09:14:56:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)

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

And here the results for 25 hours. There is still an offset of ~9%, but again, this corresponds to a tiny uncertainty in the background model:

2016-03-03T10:32:39: === GModelSky ===
2016-03-03T10:32:39:  Name ......................: pi0
2016-03-03T10:32:39:  Instruments ...............: all
2016-03-03T10:32:39:  Instrument scale factors ..: unity
2016-03-03T10:32:39:  Observation identifiers ...: all
2016-03-03T10:32:39:  Model type ................: DiffuseSource
2016-03-03T10:32:39:  Model components ..........: "MapCubeFunction" * "ConstantValue" * "Constant" 
2016-03-03T10:32:39:  Number of parameters ......: 3
2016-03-03T10:32:39:  Number of spatial par's ...: 1
2016-03-03T10:32:39:   Normalization ............: 1 [0.1,10]  (fixed,scale=1,gradient)
2016-03-03T10:32:39:  Number of spectral par's ..: 1
2016-03-03T10:32:39:   Value ....................: 0.91313 +/- 0.0190205 [0,10]  (free,scale=1,gradient)
2016-03-03T10:32:39:  Number of temporal par's ..: 1
2016-03-03T10:32:39:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
2016-03-03T10:32:39: === GCTAModelIrfBackground ===
2016-03-03T10:32:39:  Name ......................: CTABackgroundModel
2016-03-03T10:32:39:  Instruments ...............: CTA
2016-03-03T10:32:39:  Instrument scale factors ..: unity
2016-03-03T10:32:39:  Observation identifiers ...: all
2016-03-03T10:32:39:  Model type ................: "PowerLaw" * "Constant" 
2016-03-03T10:32:39:  Number of parameters ......: 4
2016-03-03T10:32:39:  Number of spectral par's ..: 3
2016-03-03T10:32:39:   Prefactor ................: 1.0103 +/- 0.00261397 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
2016-03-03T10:32:39:   Index ....................: 0.00279669 +/- 0.000848382 [-5,5]  (free,scale=1,gradient)
2016-03-03T10:32:39:   PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
2016-03-03T10:32:39:  Number of temporal par's ..: 1
2016-03-03T10:32:39:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)

#6 Updated by Knödlseder Jürgen almost 9 years ago

When fixing the index of the background the map cube model comes somewhat closer to the value of 1 (still for 25 hours):

2016-03-03T10:57:43: === GModelSky ===
2016-03-03T10:57:43:  Name ......................: pi0
2016-03-03T10:57:43:  Instruments ...............: all
2016-03-03T10:57:43:  Instrument scale factors ..: unity
2016-03-03T10:57:43:  Observation identifiers ...: all
2016-03-03T10:57:43:  Model type ................: DiffuseSource
2016-03-03T10:57:43:  Model components ..........: "MapCubeFunction" * "ConstantValue" * "Constant" 
2016-03-03T10:57:43:  Number of parameters ......: 3
2016-03-03T10:57:43:  Number of spatial par's ...: 1
2016-03-03T10:57:43:   Normalization ............: 1 [0.1,10]  (fixed,scale=1,gradient)
2016-03-03T10:57:43:  Number of spectral par's ..: 1
2016-03-03T10:57:43:   Value ....................: 0.943632 +/- 0.0167051 [0,10]  (free,scale=1,gradient)
2016-03-03T10:57:43:  Number of temporal par's ..: 1
2016-03-03T10:57:43:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
2016-03-03T10:57:43: === GCTAModelIrfBackground ===
2016-03-03T10:57:43:  Name ......................: CTABackgroundModel
2016-03-03T10:57:43:  Instruments ...............: CTA
2016-03-03T10:57:43:  Instrument scale factors ..: unity
2016-03-03T10:57:43:  Observation identifiers ...: all
2016-03-03T10:57:43:  Model type ................: "PowerLaw" * "Constant" 
2016-03-03T10:57:43:  Number of parameters ......: 4
2016-03-03T10:57:43:  Number of spectral par's ..: 3
2016-03-03T10:57:43:   Prefactor ................: 1.002 +/- 0.000669831 [0.001,1000] ph/cm2/s/MeV (free,scale=1,gradient)
2016-03-03T10:57:43:   Index ....................: 0 [-5,5]  (fixed,scale=1,gradient)
2016-03-03T10:57:43:   PivotEnergy ..............: 1e+06 [10000,1e+09] MeV (fixed,scale=1e+06,gradient)
2016-03-03T10:57:43:  Number of temporal par's ..: 1
2016-03-03T10:57:43:   Normalization ............: 1 (relative value) (fixed,scale=1,gradient)
That’s probably as good as it can be.

I merged the new code in devel and close this issue now.

#7 Updated by Knödlseder Jürgen almost 9 years ago

  • Status changed from In Progress to Closed

Also available in: Atom PDF