Fermi-LAT data analysis is performed in a coordinate system fixed with respect to distant stars. This poses problems for the inclusion of emission from the Sun and the Moon that are moving with respect to this coordinate system. This tutorial shows how to use the Solar System Tools to create a template of the Solar and Lunar emission for likelihood analysis.
Before going through this tutorial, you should go through the binned likelihood tutorial as some of the initial steps are very similar. A description of the Solar System Tools can be found in Johannesson & Orlando 2013. These tools are designed solely to account for constant emission from the moon and the quiescent sun. The templates produced will not account for emission due to solar flares.
You will need a cleaned event data file (as explained in the first step of the binned likelihood tutorial), a spacecraft data file, (also referred to as the "pointing and livetime history" file), and models for both the Solar and Lunar intensity (detailed model descriptions). You may use the files provided within this tutorial or select your own data files that can be retrieved from the LAT Data Server.
The exposure depends on the exact event selection and the Solar and Lunar templates will therefore depend on the event selection. Any changes made to the event selection requires regenerating the templates from scratch.
Precomputing the livetime for the dataset speeds up the exposure calculation. This is needed for calculating the exposure map in the next step.
This accounts for exposure as a function of energy, based on the cuts made. It is used to create an average intensity map for the Solar and Lunar emission.
To calculate the expected counts from the Sun and the Moon the livetime binned in distance from the Sun and the Moon is needed.
This accounts for exposure as a function of energy and distance from the Sun and the Moon, based on the cuts made.
This will create templates as a function of energy and sky-position (map cube) of the emission from the Sun and the Moon. These can be used in the xml file of the likelihood analysis in the same way as the diffuse background models.
Review the tutorial on binned likelihood analysis. For the rest of this tutorial the subselected event file and the corresponding spacecraft file from the binned tutorial will be used.
We use the livetime cube generated by the tool gtltcube in the binned likelihood tutorial.
Like in step 7 of the binned likelihood tutorial, we use the tool gtexpcube2 to calculate the exposure from the livetime cube. Since we are creating an all-sky template for the Solar and Lunar emission, we will create an all-sky exposure cube. Here it is important to select the proper instrument response function that matches the event file selection. In our case we used Pass 8 source class events, so we use the P8R2_SOURCE_V6 IRF. We generate the exposure map as follows:
prompt> gtexpcube2
Livetime cube file[] 3C279_binned_ltcube.fits
Counts map file[] 3C279_binned_srcmaps.fits
Output file name[] 3C279_binned_template_expcube.fits
Response functions to use[] P8R2_SOURCE_V6
Computing binned exposure map....................!
Here we have used a 0.5 degree pixel size, which should give sufficient resolution for most cases, because the Sun moves about 1 degree per day, making its average emission smeared over the sky. Since the Moon moves even faster, increasing the resolution is unnecessary. The output file generated by gtexpcube2 is available here.
NOTE: It is important to have the energy binning identical when running all the solar tools (gtexpcube2, gtexphpsun, and gtsuntemp). The tool gtsuntemp will fail with an error if that is not the case.
NOTE: If you get a "File not found" error, the tool just needs you to specify the IRF name and counts map file. In this case we need to specify P8R2_SOURCE_V6 IRF and 3C279_binned_srcmaps.fits respectively.
We use the tool gtltcubesun to calculate the livetime binned not only in sky-position and instrument angle, but also in distance from either the Sun or the Moon. This is needed to speed up the calculation of the emission from the Sun and the Moon that is assumed to be spherically symmetric. The usage of this tool is similar to that of gtltcube. Two additional parameters need to be specified: the solar system body, either SUN or MOON and the maximum angle for the distance from the solar system body. This angle should be 180 degrees for the Sun to account for its extended emission, while a value of 1 degree is more than enough for the Moon (indeed, a value of 0.5 degree can be safely used).
The tool gtltcubesun takes a long time to execute, and due to the nature of the memory storage, it is more efficient to split the calculation up into smaller time steps. A suitable time step for the calculation for the Sun is between a week and a month, while for the Moon it is appropriate to have month or year long periods. To facilitate this the hidden parameters tmin and tmax can be used to split the time range into smaller periods. The expected time system is MET and the start and stop values for the event selection can be found in the header of the event data fits file. For our dataset it is TSTART=239557417 and TSTOP=302572802. For the Sun we will split this up into 24 bins, while for the Moon we will split it up into 6 bins. In the end we will sum the livetime cubes with gtltsumsun.
While it is possible to change the pixel size, we do not recommend changing it from the default 0.5 degrees. Larger values will result in inaccurate calculation of the livetime, while smaller values significantly increase the memory needed to run the tools.
Here is an example run of gtltcubesun for the Sun:
prompt> gtltcubesun tmin=239557417 tmax=242183058
Event data file[] 3C279_binned_gti.fits
Spacecraft data file[] L130829131108BD489E7F38_SC00.fits
Output file[expCube.fits] 3C279_binned_ltcube_sun_01.fits
Name of solar system body [SUN,MOON][SUN] SUN
Step size in cos(theta) (0.:1.) [0.025]
Maximum angle in theta_sun (0.:180.) [180]
Pixel size (degrees)[0.5]
Working on file L130829131108BD489E7F38_SC00.fits
.....................!
And similarly for the Moon:
prompt> gtltcubesun tmin=239557417 tmax=250059981
Event data file[] 3C279_binned_gti.fits
Spacecraft data file[] L130829131108BD489E7F38_SC00.fits
Output file[expCube.fits] 3C279_binned_ltcube_moon_01.fits
Name of solar system body [SUN,MOON][SUN] MOON
Step size in cos(theta) (0.:1.) [0.025]
Maximum angle in theta_sun (0.:180.) [180] 1
Pixel size (degrees)[0.5]
Working on file L130829131108BD489E7F38_SC00.fits
.....................!
After we have created the files for each time bin, we have to sum them up using gtltsumsun. It requires two arguments: a file with a list of livetime cube filenames and the output filename. To create the list of files for gtltsumsun we issue the commands
prompt> ls 3C279_binned_ltcube_sun_*.fits > ltcube_files_sun.txt
prompt> ls 3C279_binned_ltcube_moon_*.fits > ltcube_files_moon.txt
Then it is a simple matter of running gtltsumsun to get the final livetime cubes:
prompt> gtltsumsun
Livetime cube 1 or list of files[] ltcube_files_sun.txt
Output file[] 3C279_binned_ltcube_sun.fits
Adding livetime cubes....................... Done
Merging cuts Done
Writing file Done
prompt> gtltsumsun
Livetime cube 1 or list of files[] ltcube_files_moon.txt
Output file[] 3C279_binned_ltcube_moon.fits
Adding livetime cubes..... Done
Merging cuts Done
Writing file Done
NOTE: Summing up the Solar livetime cubes for this example required 13 GB of memory. Make sure you run it on a machine with plenty of RAM and SWAP space if needed. The memory requirements will be smaller for observing periods shorter than a year. The Lunar livetime cubes are smaller because the distance from the Moon is cut of at 1 degree.
The livetime cubes generated in this step are available for the Sun and for the Moon.We use the tool gtexphpsun to calculate the exposure from the livetime cubes calculated with gtltcubesun and gtltsumsun in step 4. Unlike gtexpcube2, this tool calculates the exposure on a HEALPix grid and is therefore always full-sky. It uses the same pixel size as used in gtltcubesun but the energy binning has to be specified. It is important to use the same energy binning here as was used in step 3 when calculating the exposure map with gtexpcube2.
To generate the exposure map for the Sun:
prompt> gtexphpsun
Solar system tools livetime cube file[] 3C279_binned_ltcube_sun.fits
Output file name[] 3C279_binned_expcube_sun.fits
Response functions to use[] P8R2_SOURCE_V6
Start energy (MeV) of first bin[] 100
Stop energy (MeV) of last bin[] 100000
Number of logarithmically-spaced energy bins[] 30
Computing binned exposure map.....................!
And for the Moon:
prompt> gtexphpsun
Solar system tools livetime cube file[] 3C279_binned_ltcube_moon.fits
Output file name[] 3C279_binned_expcube_moon.fits
Response functions to use[] P8R2_SOURCE_V6
Start energy (MeV) of first bin[] 100
Stop energy (MeV) of last bin[] 100000
Number of logarithmically-spaced energy bins[] 30
Computing binned exposure map.....................!
The generated exposure cubes generated in this step are available for the Sun and for the Moon.
We use gtsuntemp to create an average template of the Solar and Lunar emission suitable for inclusion in likelihood analysis. It uses the exposure maps created in steps 3 and 5 along with a model for the intensity of the Sun and Moon (detailed model descriptions). Here is how we create the template for the Sun:
prompt> gtsuntemp
Exposure binned in healpix and solar angles[] 3C279_binned_expcube_sun.fits
Binned exposure[] 3C279_binned_template_expcube.fits
Fits file containing solar intensity profile[] solar_profile_v2r0.fits
Counts map file[] none
Output file name[] 3C279_binned_template_sun.fits
Size of the X axis in pixels[] 720
Size of the Y axis in pixels[] 360
Image scale (in degrees/pixel)[] 0.5
Coordinate system (CEL - celestial, GAL -galactic) (CEL|GAL) [] CEL
First coordinate of image center in degrees (RA or galactic l)[] 0
Second coordinate of image center in degrees (DEC or galactic b)[] 0
Rotation angle of image axis, in degrees[] 0
Projection method e.g. AIT|ARC|CAR|GLS|MER|NCP|SIN|STG|TAN[] CAR
Start energy (MeV) of first bin[] 100
Stop energy (MeV) of last bin[] 100000
Number of logarithmically-spaced energy bins[] 30
Computing binned solarTemplate map....................!
And similarly for the Moon:
prompt> gtsuntemp
Exposure binned in healpix and solar angles[] 3C279_binned_expcube_moon.fits
Binned exposure[] 3C279_binned_template_expcube.fits
Fits file containing solar intensity profile[] lunar_profile_v2r0.fits
Counts map file[] none
Output file name[] 3C279_binned_template_moon.fits
Size of the X axis in pixels[] 720
Size of the Y axis in pixels[] 360
Image scale (in degrees/pixel)[] 0.5
Coordinate system (CEL - celestial, GAL -galactic) (CEL|GAL) [] CEL
First coordinate of image center in degrees (RA or galactic l)[] 0
Second coordinate of image center in degrees (DEC or galactic b)[] 0
Rotation angle of image axis, in degrees[] 0
Projection method e.g. AIT|ARC|CAR|GLS|MER|NCP|SIN|STG|TAN[] CAR
Start energy (MeV) of first bin[] 100
Stop energy (MeV) of last bin[] 100000
Number of logarithmically-spaced energy bins[] 30
Computing binned solarTemplate map....................!
The templates generated in this step are available for the Sun and for the Moon. The templates can then be used in likelihood analysis by including them in the xml model, by adding the following instructions for the Sun:
<source name="sun" type="DiffuseSource">
<spectrum type="PowerLaw">
<parameter free="1" max="10" min="0" name="Prefactor" scale="1" value="1"/>
<parameter free="0" max="1" min="-1" name="Index" scale="1.0" value="0"/>
<parameter free="0" max="2e2" min="5e1" name="Scale" scale="1.0" value="1e2"/>
</spectrum>
<spatialModel file="3C279_binned_template_sun.fits" type="MapCubeFunction">
<parameter free="0" max="1e3" min="1e-3" name="Normalization" scale="1.0" value="1.0"/>
</spatialModel>
</source>
and for the Moon:
<source name="sun" type="DiffuseSource">
<spectrum type="PowerLaw">
<parameter free="1" max="10" min="0" name="Prefactor" scale="1" value="1"/>
<parameter free="0" max="1" min="-1" name="Index" scale="1.0" value="0"/>
<parameter free="0" max="2e2" min="5e1" name="Scale" scale="1.0" value="1e2"/>
</spectrum>
<spatialModel file="3C279_binned_template_moon.fits" type="MapCubeFunction">
<parameter free="0" max="1e3" min="1e-3" name="Normalization" scale="1.0" value="1.0"/>
</spatialModel>
</source>
Note that the templates provided with this tutorial cannot be used in analysis with a different data selection, because the templates depend on the data.
Last updated by: Elizabeth Ferrara 4/30/2015