GammaLib

I guess you would also have to think about the XML interface. I am not familiar with these models, but it should be possible to access each model parameter from the XML file, too. I therefore think it might be worth to have individual classes for each model.
Is it common to fit the extension or position of these models? If not wouldn’t it be easier to just use FITS templates?

I was hoping to keep the profile and interaction values as model parameters (and not private variables), though always fixed, to make sure they are naturally saved to/read from the xml files.

Its not common to fit the position/halodist, since these are usually known beforehand through measuring star or galactic rotation velocities.

The extension (scale radius/density,exponents) is usually what’s being searched for, to improve DM upper limits, so it’d be nice to give the likelihood engine the freedom to play with those parameters, rather than using a static fits image that has to be changed every time you want to test a different scale parameter.

I guess a better example would be:

model = gammalib.GModelSpatialRadialDarkMatterHalo('zhao','selfannihilation') # sets and fixes these model parameters
model['profile'      ].string('einasto') # though you can still change them to something else
model['interaction'  ].string('decay')
model['halodist'     ].real(2.5e20) # distance from Earth to halo center, in meters (8kpc to Galactic Center = 2.5e20meters)
model['ra'           ].real(44.44 ) # deg
model['dec'          ].real(44.44 ) # deg
model['scale_radius' ].real(2e14  ) # r_s, meters
model['scale_density'].real(1e-5  ) # \rho_s, density at the scale radius, g/cm^3 ?
model['alpha'        ].real(1     ) # float
model['beta'         ].real(3     ) # float, only zhao
model['gamma'        ].real(2     ) # float, only zhao

I guess I’d be worried about having 6 separate model cpp/hpp files with huge names (that are only used in DM studies), when their only significant difference is in the eval() function.

I got today this request from Gonzalo Rodriguez Fdez. I guess it is related to this thread hence I post it here. I asked Gonzalo to sign-up and participate in the discussion.

Dear Jurgen,

In order to study the expected signal from dark matter within CTOOLS
we have to created several ascii files with the expected flux for
differences wimp mass and for difference annihilation channels. In
this way the analysis and the code con be a little messy.
I would like to suggest  the solution adopted by the Fermi
collaboration in their science tools. They use the DMFit tool to
calculated the gamma spectrum, then from a xml file they just model
the spectrum using 4 parameters, normalization, channel, branching
ratio and mass.

Do you think that this solutions can be implemented in CTOOLS?

Regards,
Gonzalo.

Hi all,

few comments:
1/ I prefer one unique function with parameters for DM density
2/ The integration of NFW profiles (divergent at the center) can be tricky, and needs some careful attention
3/ There is a public code (Clumpy) with many such features, including intrinsic boost from substructures, parameterization of DM ann/decay spectrum, etc, (c++), I wonder whether it can be linked to CTOOLS or a standalone code is better
4/ In the end of 2015/begin of 2016, several DM papers from CTA are expected. It would be good to have a complete cross-check with CTOOLS.

Thanks for the support!
Michele

Here, Fermi source model, you can find how the dark matter energy distribution is defined in the sciences tools definition files.

Thanks,
Gonzalo.

Thanks for the link.

So you are talking in fact about a spectral model (the thread was more for a spatial model). I guess we need both, spectral and spatial model components.

The spectral model from DMFit looks similar to a file function, just a bit more sophisticated. Is the format of the data file defined somewhere? Is it a fixed format (number of rows and columns fixed) or can this evolve?

Hi,

the Cirelli et al DM spectral have become quite standard for now. They have fixed columns and rows. And they should not change in the future often. Possibly new channels/new models may appear or some corrections.

Michele

Doro Michele wrote:

Hi,

the Cirelli et al DM spectral have become quite standard for now. They have fixed columns and rows. And they should not change in the future often. Possibly new channels/new models may appear or some corrections.

Michele

Wouldn’t it make sense to transform this ASCII file into a FITS file?

In that way metadata can be added to the header, and one would also easily cope with a possible evolution of the model.

Hi Juergen,

I am not accustomed to .fits format, so I cannot judge the complexity of doing it. But, being just a table of rows and columns with some headers probably, I guess it could be easily done.

M

Kelley-Hoskins Nathan wrote:

Well, it turns out the integral you need to calculate the J factor does not have an analytical solution. I think this means every time eval() is called, its going to have to numerically integrate that function, which is going to be slow. I see a GFunction and GIntegral, I imagine I can use these to do the integration in this model’s eval()?

Yes, you can use GFunction and GIntegral for that.

Is it possible to pre-compute things and put them into a cache (that’s how it works for some of the other models)?

I think there is no fundamental problem, you just need to get the units right.

The trick is to get a spatial template that normalizes to unity when integrated over the solid angle, so that the spectral component provides the flux density.

You can achieve this for any spatial template (in whatever units it is) by computing

My flux equation is this:

The I don’t think there is a component of that equation that can be split off that just has the units 1/sr.

get a spatial template that normalizes to unity when integrated over the solid angle

I’m not sure what is meant by this, or what the advantage is? Doesn’t this mess with the physical meaning of the above flux equation? And what would the limits of integration be for the solid angle integral?

Ah, I think I’m starting to understand. I’ll try this.

Is there a 2D integrator hidden somewhere in gammalib? It looks like GFunction/GIntegral only take one input parameter for integration. The integrand has no analytical integral with respect to or ...

where

Kelley-Hoskins Nathan wrote:

Is there a 2D integrator hidden somewhere in gammalib? It looks like GFunction/GIntegral only take one input parameter for integration. The integrand has no analytical integral with respect to or ...

where

For 2D you just nest a 1D integral within another (using the GFunction/GIntegral classes).