Iftach Sadeh / ctools

Commit dc4df07ca1bc85ad3e69bec35d233f4ee03cbed1

Authored by Jürgen Knödlseder
1 parent 6da50f30

Integrate code and add exponential spectral model test to science verification (#2454)

Inline Side-by-side
Showing 6 changed files with 186 additions and 22 deletions
ChangeLog
  View file @dc4df07
1 -2018-06-20 1 +2018-07-03
2 2
3 * Version 1.6.0 released 3 * Version 1.6.0 released
4 ======================== 4 ========================
@@ -6,6 +6,7 @@ @@ -6,6 +6,7 @@
6 Transparently update parameter files (#2513) 6 Transparently update parameter files (#2513)
7 Change formula for number of true energy bins in csphagen (#2490) 7 Change formula for number of true energy bins in csphagen (#2490)
8 Remove ctools::is_valid_filename() method 8 Remove ctools::is_valid_filename() method
  9 + Add Exponential spectral model to science verification (#2454)
9 Fix bug when using precomputed model in csresspec (#2448) 10 Fix bug when using precomputed model in csresspec (#2448)
10 Keep existing library path in installation script on Mac OS X (#2420) 11 Keep existing library path in installation script on Mac OS X (#2420)
11 csphagen produces Xspec compliant data (#2404) 12 csphagen produces Xspec compliant data (#2404)
1 New Features and Important Changes in ctools 1.6.0 1 New Features and Important Changes in ctools 1.6.0
2 2
3 -20 June 2018 3 +3 July 2018
4 4
5 5
6 Introduction 6 Introduction
doc/source/users/user_manual/models_spectral.rst
  View file @dc4df07
@@ -503,34 +503,49 @@ Multiplicative model @@ -503,34 +503,49 @@ Multiplicative model
503 model components that are multiplied. 503 model components that are multiplied.
504 504
505 505
506 -Exponential model: EBL absorption  
507 -^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^  
508 -  
509 - This spectral model component implements the exponential of an  
510 - arbitrary spectral model. It can be used to model a spectrum with  
511 - EBL absorption based on a tabulated model of opacity as a function  
512 - of photon energy. 506 +Exponential model
  507 +^^^^^^^^^^^^^^^^^
513 508
514 .. code-block:: xml 509 .. code-block:: xml
515 510
516 - <spectrum type="Multiplicative">  
517 - <spectrum type="PowerLaw">  
518 - <parameter name="Prefactor" scale="1e-16" value="5.7" min="1e-07" max="1000.0" free="1"/>  
519 - <parameter name="Index" scale="-1" value="2.48" min="0.0" max="+5.0" free="1"/>  
520 - <parameter name="PivotEnergy" scale="1e6" value="0.3" min="0.01" max="1000.0" free="0"/>  
521 - </spectrum>  
522 <spectrum type="Exponential"> 511 <spectrum type="Exponential">
523 <spectrum type="FileFunction" file="opacity.txt"> 512 <spectrum type="FileFunction" file="opacity.txt">
524 - <parameter scale="-1.0" name="Normalization" min="0." max="100.0" value="1.0" free="1"/> 513 + <parameter scale="-1.0" name="Normalization" min="0.0" max="100.0" value="1.0" free="1"/>
525 </spectrum> 514 </spectrum>
526 </spectrum> 515 </spectrum>
527 - </spectrum>  
528 516
529 - This example corresponds to a spectrum 517 + This spectral model component implements the exponential of an arbitrary
  518 + spectral model and computes
  519 +
  520 + .. math::
  521 + M_{\rm spectral}(E) = \exp \left( M_{\rm spectral}(E) \right)
  522 +
  523 + where :math:`M_{\rm spectral}(E)` is any spectral model component.
  524 +
  525 + The model can be used to describe a spectrum with EBL absorption based on a
  526 + tabulated model of opacity as a function of photon energy. The corresponding
  527 + XML file structure for such a model is shown below:
  528 +
  529 + .. code-block:: xml
530 530
531 - .. math:: 531 + <spectrum type="Multiplicative">
  532 + <spectrum type="PowerLaw">
  533 + <parameter name="Prefactor" scale="1e-16" value="5.7" min="1e-07" max="1000.0" free="1"/>
  534 + <parameter name="Index" scale="-1" value="2.48" min="0.0" max="+5.0" free="1"/>
  535 + <parameter name="PivotEnergy" scale="1e6" value="0.3" min="0.01" max="1000.0" free="0"/>
  536 + </spectrum>
  537 + <spectrum type="Exponential">
  538 + <spectrum type="FileFunction" file="opacity.txt">
  539 + <parameter scale="-1.0" name="Normalization" min="0.0" max="100.0" value="1.0" free="1"/>
  540 + </spectrum>
  541 + </spectrum>
  542 + </spectrum>
  543 +
  544 + This example corresponds to the function
  545 +
  546 + .. math::
532 M_{\rm spectral}(E) = k_0 \left( \frac{E}{E_0} \right)^{\gamma} 547 M_{\rm spectral}(E) = k_0 \left( \frac{E}{E_0} \right)^{\gamma}
533 - \times \exp\left( -\alpha \tau\left( E \right) \right) 548 + \times \exp\left( -\alpha \, \tau(E) \right)
534 549
535 where 550 where
536 551
@@ -540,3 +555,9 @@ Exponential model: EBL absorption @@ -540,3 +555,9 @@ Exponential model: EBL absorption
540 as first column and opacity :math:`\tau` as second column, respectively; 555 as first column and opacity :math:`\tau` as second column, respectively;
541 * the parameter :math:`\alpha` = ``Normalization`` represents an 556 * the parameter :math:`\alpha` = ``Normalization`` represents an
542 opacity scaling factor. 557 opacity scaling factor.
  558 +
  559 + .. note::
  560 + The ``Exponential`` model implements the function :math:`y=\exp(x)`,
  561 + hence in the example the ``scale`` attribute of the ``Normalization``
  562 + parameter was set to ``-1`` to implement the form
  563 + :math:`y=\exp(-x)`.
test/data/sciver/crab_exponential.xml 0 → 100755
  View file @dc4df07
  1 +<?xml version="1.0" standalone="no"?>
  2 +<source_library title="source library">
  3 + <source name="Crab" type="PointSource">
  4 + <spectrum type="Multiplicative">
  5 + <spectrum type="PowerLaw">
  6 + <parameter name="Prefactor" scale="1e-16" value="5.7" min="1e-07" max="1000.0" free="1"/>
  7 + <parameter name="Index" scale="-1" value="2.48" min="0.0" max="+5.0" free="1"/>
  8 + <parameter name="PivotEnergy" scale="1e6" value="0.3" min="0.01" max="1000.0" free="0"/>
  9 + </spectrum>
  10 + <spectrum type="Exponential">
  11 + <spectrum type="FileFunction" file="opacity.txt">
  12 + <parameter scale="-1.0" name="Normalization" min="0.0" max="100.0" value="1.0" free="1"/>
  13 + </spectrum>
  14 + </spectrum>
  15 + </spectrum>
  16 + <spatialModel type="PointSource">
  17 + <parameter name="RA" scale="1.0" value="83.6331" min="-360" max="360" free="0"/>
  18 + <parameter name="DEC" scale="1.0" value="22.0145" min="-90" max="90" free="0"/>
  19 + </spatialModel>
  20 + </source>
  21 + <source name="Background" type="CTAIrfBackground" instrument="CTA">
  22 + <spectrum type="PowerLaw">
  23 + <parameter name="Prefactor" scale="1.0" value="1.0" min="1e-3" max="1e+3" free="1"/>
  24 + <parameter name="Index" scale="1.0" value="0.0" min="-5.0" max="+5.0" free="1"/>
  25 + <parameter name="PivotEnergy" scale="1e6" value="1.0" min="0.01" max="1000.0" free="0"/>
  26 + </spectrum>
  27 + </source>
  28 +</source_library>
test/data/sciver/opacity.txt 0 → 100644
  View file @dc4df07
  1 +1.000000000000000000e+03 2.130514478920406932e-04
  2 +1.149756995397735636e+03 2.518904850653095901e-04
  3 +1.321941148466028608e+03 2.978098346395106704e-04
  4 +1.519911082952933157e+03 3.521002295303738675e-04
  5 +1.747528400007682876e+03 4.162877017994022770e-04
  6 +2.009233002565047855e+03 4.921764774211230195e-04
  7 +2.310129700083160515e+03 5.818996907177256166e-04
  8 +2.656087782946686730e+03 6.879793439775886213e-04
  9 +3.053855508833415570e+03 8.133971976442146796e-04
  10 +3.511191734215130964e+03 9.616785837061616292e-04
  11 +4.037017258596553802e+03 1.136991498172844079e-03
  12 +4.641588833612781855e+03 1.344263757996223352e-03
  13 +5.336699231206313016e+03 1.589321515566356410e-03
  14 +6.135907273413175972e+03 1.879053024242312688e-03
  15 +7.054802310718645458e+03 2.221602258153513362e-03
  16 +8.111308307896872066e+03 2.626597828671136097e-03
  17 +9.326033468832200015e+03 3.105423632092475100e-03
  18 +1.072267222010323167e+04 3.671538836091017399e-03
  19 +1.232846739442065882e+04 4.340856199333244930e-03
  20 +1.417474162926804820e+04 5.132189358332243760e-03
  21 +1.629750834620643400e+04 6.067781654187125867e-03
  22 +1.873817422860383158e+04 7.173931363836940321e-03
  23 +2.154434690031886566e+04 8.481730910262627252e-03
  24 +2.477076355991714081e+04 1.002794082986985247e-02
  25 +2.848035868435804878e+04 1.185602306313406941e-02
  26 +3.274549162877731578e+04 1.401736261295742636e-02
  27 +3.764935806792470976e+04 1.657271190995783169e-02
  28 +4.328761281083061476e+04 1.959389848391105901e-02
  29 +4.977023564332113892e+04 2.316584394175889502e-02
  30 +5.722367659350220492e+04 2.738895100301691130e-02
  31 +6.579332246575682075e+04 3.238192568902821533e-02
  32 +7.564633275546290679e+04 3.828511399411545918e-02
  33 +8.697490026177834079e+04 4.526444682809729891e-02
  34 +1.000000000000000000e+05 5.351610411734870898e-02
  35 +1.149756995397735591e+05 6.327202916619174833e-02
  36 +1.321941148466028681e+05 7.480644828010984959e-02
  37 +1.519911082952933211e+05 8.844357890887546980e-02
  38 +1.747528400007683085e+05 1.045667429754224553e-01
  39 +2.009233002565045899e+05 1.236291415542298683e-01
  40 +2.310129700083158095e+05 1.461665937613474131e-01
  41 +2.656087782946684165e+05 1.728125979296164794e-01
  42 +3.053855508833412314e+05 2.043161384190416152e-01
  43 +3.511191734215127653e+05 2.415627385884853706e-01
  44 +4.037017258596549509e+05 2.855993517002112791e-01
  45 +4.641588833612782182e+05 3.376637894080791180e-01
  46 +5.336699231206312543e+05 3.992195150257371594e-01
  47 +6.135907273413175717e+05 4.719967795681307421e-01
  48 +7.054802310718644876e+05 5.580412568466853740e-01
  49 +8.111308307896872284e+05 6.597715446871550427e-01
  50 +9.326033468832198996e+05 7.800471485542282979e-01
  51 +1.072267222010323079e+06 9.222488585137500472e-01
  52 +1.232846739442065824e+06 1.090373778823942263e+00
  53 +1.417474162926804740e+06 1.289147681313261629e+00
  54 +1.629750834620646900e+06 1.524157840651543472e+00
  55 +1.873817422860386781e+06 1.802010085340306489e+00
  56 +2.154434690031886566e+06 2.130514478920409349e+00
  57 +2.477076355991714168e+06 2.518904850653099992e+00
  58 +2.848035868435804732e+06 2.978098346395110951e+00
  59 +3.274549162877731491e+06 3.521002295303742535e+00
  60 +3.764935806792471092e+06 4.162877017994027540e+00
  61 +4.328761281083061360e+06 4.921764774211230709e+00
  62 +4.977023564332113601e+06 5.818996907177255373e+00
  63 +5.722367659350220114e+06 6.879793439775887087e+00
  64 +6.579332246575682424e+06 8.133971976442145646e+00
  65 +7.564633275546290912e+06 9.616785837061614828e+00
  66 +8.697490026177834719e+06 1.136991498172844395e+01
  67 +1.000000000000000000e+07 1.344263757996222175e+01
  68 +1.149756995397735573e+07 1.589321515566354392e+01
  69 +1.321941148466028646e+07 1.879053024242310599e+01
  70 +1.519911082952933200e+07 2.221602258153511400e+01
  71 +1.747528400007683039e+07 2.626597828671133072e+01
  72 +2.009233002565046027e+07 3.105423632092470854e+01
  73 +2.310129700083158165e+07 3.671538836091013280e+01
  74 +2.656087782946689427e+07 4.340856199333251197e+01
  75 +3.053855508833418787e+07 5.132189358332250606e+01
  76 +3.511191734215134382e+07 6.067781654187132290e+01
  77 +4.037017258596558124e+07 7.173931363836948094e+01
  78 +4.641588833612781763e+07 8.481730910262618295e+01
  79 +5.336699231206312776e+07 1.002794082986983994e+02
  80 +6.135907273413176090e+07 1.185602306313405592e+02
  81 +7.054802310718645155e+07 1.401736261295741031e+02
  82 +8.111308307896871865e+07 1.657271190995780898e+02
  83 +9.326033468832199275e+07 1.959389848391103612e+02
  84 +1.072267222010325342e+08 2.316584394175892214e+02
  85 +1.232846739442063421e+08 2.738895100301680827e+02
  86 +1.417474162926807702e+08 3.238192568902825315e+02
  87 +1.629750834620640278e+08 3.828511399411532921e+02
  88 +1.873817422860386968e+08 4.526444682809735696e+02
  89 +2.154434690031886697e+08 5.351610411734878880e+02
  90 +2.477076355991714001e+08 6.327202916619185089e+02
  91 +2.848035868435804844e+08 7.480644828010995298e+02
  92 +3.274549162877731919e+08 8.844357890887558824e+02
  93 +3.764935806792471409e+08 1.045667429754226077e+03
  94 +4.328761281083061695e+08 1.236291415542300228e+03
  95 +4.977023564332113862e+08 1.461665937613475762e+03
  96 +5.722367659350219965e+08 1.728125979296166861e+03
  97 +6.579332246575682163e+08 2.043161384190418630e+03
  98 +7.564633275546290874e+08 2.415627385884856722e+03
  99 +8.697490026177834272e+08 2.855993517002116278e+03
  100 +1.000000000000000000e+09 3.376637894080786737e+03
test/science_verification.py
  View file @dc4df07
@@ -11,7 +11,7 @@ @@ -11,7 +11,7 @@
11 # 11 #
12 # -------------------------------------------------------------------------- 12 # --------------------------------------------------------------------------
13 # 13 #
14 -# Copyright (C) 2015-2017 Juergen Knoedlseder 14 +# Copyright (C) 2015-2018 Juergen Knoedlseder
15 # 15 #
16 # This program is free software: you can redistribute it and/or modify 16 # This program is free software: you can redistribute it and/or modify
17 # it under the terms of the GNU General Public License as published by 17 # it under the terms of the GNU General Public License as published by
@@ -251,6 +251,7 @@ class sciver(gammalib.GPythonTestSuite): @@ -251,6 +251,7 @@ class sciver(gammalib.GPythonTestSuite):
251 self.append(self.spec_gauss, 'Test Gaussian model') 251 self.append(self.spec_gauss, 'Test Gaussian model')
252 self.append(self.spec_filefct, 'Test file function model') 252 self.append(self.spec_filefct, 'Test file function model')
253 self.append(self.spec_nodes, 'Test nodes model') 253 self.append(self.spec_nodes, 'Test nodes model')
  254 + self.append(self.spec_exponential, 'Test exponential model')
254 255
255 # Append spatial tests 256 # Append spatial tests
256 self.append(self.spat_ptsrc, 'Test point source model') 257 self.append(self.spat_ptsrc, 'Test point source model')
@@ -339,7 +340,7 @@ class sciver(gammalib.GPythonTestSuite): @@ -339,7 +340,7 @@ class sciver(gammalib.GPythonTestSuite):
339 mean_max = +lim_mean 340 mean_max = +lim_mean
340 std_min = 1.0-lim_std 341 std_min = 1.0-lim_std
341 std_max = 1.0+lim_std 342 std_max = 1.0+lim_std
342 - 343 +
343 # Test mean 344 # Test mean
344 mean = self.results[name]['mean'] 345 mean = self.results[name]['mean']
345 valid = (mean >= mean_min) and (mean <= mean_max) 346 valid = (mean >= mean_min) and (mean <= mean_max)
@@ -529,6 +530,19 @@ class sciver(gammalib.GPythonTestSuite): @@ -529,6 +530,19 @@ class sciver(gammalib.GPythonTestSuite):
529 self.test('Pull_Background_Index') 530 self.test('Pull_Background_Index')
530 return 531 return
531 532
  533 + # Test exponential model
  534 + def spec_exponential(self):
  535 + """
  536 + Test exponential model
  537 + """
  538 + self.pull('data/sciver/crab_exponential')
  539 + self.test('Pull_Crab_1:Prefactor')
  540 + self.test('Pull_Crab_1:Index')
  541 + self.test('Pull_Crab_2:Normalization')
  542 + self.test('Pull_Background_Prefactor')
  543 + self.test('Pull_Background_Index')
  544 + return
  545 +
532 # Test point source model 546 # Test point source model
533 def spat_ptsrc(self): 547 def spat_ptsrc(self):
534 """ 548 """