This section is dedicated to a few tutorials on how to use ixpeobssim.
Defining a point source¶
The heart of the code lies in how we define the source that we wish to simulate
and ultimately study. In ixpeobssim we do that by writing a configuration file
in python. Several example sources can be found in the
config directory, here we
will describe how to simulate a single point source with a 50 % constant
polarization degree, constant polarization angle of 30 degrees and a power-law
spectrum with index 2 and normalization of 10. Below is a snapshot of a
configuration file defining this single point source:
import numpy from ixpeobssim.config import file_path_to_model_name from ixpeobssim.srcmodel.roi import xPointSource, xROIModel from ixpeobssim.srcmodel.spectrum import power_law from ixpeobssim.srcmodel.polarization import constant # Define the source properties and create the actual source. ra, dec = 45., 45. PL_NORM = 10. PL_INDEX = 2. #Constant polarization degree and angle pol_deg = constant(0.5) pol_ang = constant(numpy.radians(30.)) #Create the spectrum with a power-law shape with index 2 and normalization 10 spec = power_law(PL_NORM, PL_INDEX) #Create the point source using the xPointSource class src = xPointSource('Point source', ra, dec, spec, pol_deg, pol_ang) # Create the complete region of interest. ROI_MODEL = xROIModel(ra, dec) ROI_MODEL.add_source(src)
We create our point source by passing our parameters (ra, dec, spectrum, ect) to
ixpeobssim.srcmodel.roi.xPointSource class. The first value passed
to the xPointSource class is the name of the source, in this case we simply
called it ‘Point source’. You can also specify the column density (value used to
calculate the Galactic absorption) and the redshift of your source by passing it
to the class.
When running a simulation in ixpeobssim we define our region of interest (ROI)
and populate it with the sources that we wish to simulate. This is what is seen
in the above configuration file. We define the point source and then define our
ROI_MODEL by calling the
The ra and dec values that we pass to this class define the center of the ROI.
We then proceed to add the source we just defined to the ROI. Note that we can
add any number of sources to our ROI and each source can have their own positions
(different ra and dec values than those used to define the center of the ROI).
This is all you need to create a single point source in ixpeobssim. Once you have created your configuration file you can run a simulation of the duration you wish as described in the Quick start section under the xpobssim section of the documentation.
Defining an extended source¶
Now that we know how to define a simple point source we can move on to a more complicated example of how to build an extended source. There are two ways to go about doing that in ixpeobssim both of which rely on observation fits files:
create the extended source using the
create the extended source using the
Here we will illustrate how to use both cases.
ixpeobssim.srcmodel.roi.xExtendedSource class is very similar
to what we did in the case of a single point source in the sense that we will
need to define the ROI by calling the xROIModel class and passing it the ra and
dec values for the center of your region of interest. You will also need to
provide the parameters describing the source itself (just as was done for the
point source). The major difference lies in the fact that you need to provide
an image fits file (for example a Chandra observation) that describes the
morphology of the extended source. So for example if we were to use the same
parameters used for the single point source but decide to assign them to the
morphology of CasA, then you would do the following:
import numpy import os from ixpeobssim.srcmodel.roi import xExtendedSource, xROIModel from ixpeobssim.srcmodel.polarization import constant from ixpeobssim.srcmodel.spectrum import power_law from ixpeobssim.core.spline import xInterpolatedUnivariateSplineLinear from ixpeobssim import IXPEOBSSIM_CONFIG # Define the source properties and create the actual source. #This is the center of the extended source (CasA) ra, dec = 350.8664167, 58.8117778 #Parameters to define the power-law spectrum of the source #with index 2 and normalization 10 PL_NORM = 10. PL_INDEX = 2. spec = power_law(PL_NORM, PL_INDEX) #Constant polarization degree and angle pol_deg = constant(0.5) pol_ang = constant(numpy.radians(30.)) #Here is the path to the Chandra observation in the energy range 4-6 keV #of CasA. img_file_path = os.path.join(IXPEOBSSIM_CONFIG, 'fits', 'casa_4p0_6p0_keV.fits') #Create the extended source calling it CasA using the xExtendedSource class src = xExtendedSource('CasA', img_file_path, spec, pol_deg, pol_ang) # Create the complete region of interest by calling the xROIModel class and #adding CasA to the ROI. ROI_MODEL = xROIModel(ra, dec) ROI_MODEL.add_source(src)
ixpeobssim.srcmodel.roi.xChandraObservation class allows you to use
a Chandra photon list to describe the source that you wish to simulate,
conserving the temporal, energetic and spatial information. How this is
implemented in ixpeobssim is described in the Implementation details section.
The only parameters that you need to provide are those specifying the polarization
degree and angle and of course the observation fits file. When building a
configuration file using a Chandra photon list you need to use the
ixpeobssim.srcmodel.roi.xChandraROIModel class instead of the
ixpeobssim.srcmodel.roi.xROIModel class that was used in the previous
examples. It is also possible to pass a map describing the polarization degree
and angle for the extended source that you want to simulate. To do this you need
to use the
Here is an example of how to build a configuration file to simulate CasA using
the morphology from a Chandra observation and describing the polarization
properties with fit maps:
import os from ixpeobssim.srcmodel.roi import xChandraObservation, xChandraROIModel from ixpeobssim.srcmodel.polarization import xStokesSkyMap from ixpeobssim import IXPEOBSSIM_CONFIG #Path to the Chandra observation file chandra_evt_file_path = 'path/to/evt/file/casa.fits' #Here we call the xChandraROIModel class instead of the XROIModel class #as was done for the previous examples and we pass the path to the #Chandra observation file. You need to specidy the acis detector of the #Chandra photon list being used for the simulation. ROI_MODEL = xChandraROIModel(chandra_evt_file_path, acis='S') #Here you define the polarization degree and angle you want to assign #to the source. In this case the polarization map is a simple #geometrical, radially-symmetric, model. pol_mapx_path = os.path.join(IXPEOBSSIM_CONFIG, 'fits', 'casa_pmax050_reg000_x.fits') pol_mapy_path = os.path.join(IXPEOBSSIM_CONFIG, 'fits', 'casa_pmax050_reg000_y.fits') #We call the xStokesSkyMap class that handles the polarization maps. #The "rotate" argument is intended to deal with input maps where the #reference system is rotated by 90 degrees (counter-clockwise) with #respect to the standard system that we use throughout, i.e., with the #x-axis horizontal and the angles measured from the x-axis itself #counter clockwise. In this case the two components of the polarization #vector are properly swapped, i.e., x -> y and -y -> x. polarization_map = xStokesSkyMap.load_from_xy(pol_mapx_path, pol_mapy_path, rotate=True) #We define the functions for the polarization degree and angle def polarization_angle(E, t, ra, dec): return polarization_map.polarization_angle(ra, dec) def polarization_degree(E, t, ra, dec): return polarization_map.polarization_degree(ra, dec) #Here we create the source called CasA using the XChandraObservation #class and pass the name of the source, the polarization degree and #angle you defined above. casa = xChandraObservation('CasA', polarization_degree, polarization_angle) #Add the source CasA to the ROI model ROI_MODEL.add_source(casa)
In this example we used a home-made map to define the polarization degree and
angle of the source. In ixpeobssim we have the capabilities of accepting maps
(in fits format) of the polarization parameters (Q/U, x/y polarization
components and polarization degree and angle) through the
In this case we have created a simple model to describe the polarization as a
geometrical, radially-symmetric toy-model. Once we have created the polarization
maps, the xChandraObservation uses the spatial and spectral information of the
source from the Chandra observation and combines it with the polarization model
to build the source to simulate.
You can pass an arbitrary number of layers in different energy bands to the
ixpeobssim.srcmodel.polarization.xStokesSkyMap class, each layer is
essentially a Stokes sky map and we interpolate the polarimetric information in
three dimensions (Ra, Dec, energy). For an example of how to implement a more
complex simulation of an extended source see the ixpeobssim/examples/tycho_3d.py
Defining a periodic source¶
We have the
ixpeobssim.srcmodel.roi.xPeriodicSource class that allows
us to easily build a point-like periodic source to simulate. The basic
ingredients to describe the point source are the same as in the simple point
source configuration but this time we need to include also the information on
the ephemeris of the pulsar. The ephemeris are handled by the
ixpeobssim.srcmodel.roi.xEphemeris class. Here is an example of how to
write a configuration file describing a pulsar with a power-law photon spectrum
with normalization that varies as a function of the pulse phase. To illustrate
the flexibility of the code, here we will also make the polarization degree vary
as a function of both energy and pulsar phase.
import numpy from ixpeobssim.config import file_path_to_model_name from ixpeobssim.srcmodel.roi import xPeriodicPointSource, xEphemeris, xROIModel from ixpeobssim.srcmodel.spectrum import power_law from ixpeobssim.srcmodel.polarization import constant ra, dec = 45., 45. #Method that returns a power-law normalization as a function of phase. def pl_norm(phase): """Photon spectrum: power-law normalization as a function of the pulse phase. """ return 1.25 + numpy.cos(4 * numpy.pi * phase) #Constant power-law index pl_index = 2. #Use the function power_law to build the photon spectrum of the pulsar. spec = power_law(pl_norm, pl_index) #Method that returns the polarization degree as a function of energy and phase def pol_deg(E, phase, ra=None, dec=None): """Polarization degree as a function of the dynamical variables. Since we're dealing with a point source the sky direction (ra, dec) is irrelevant and, as they are not used, defaulting the corresponding arguments to None allows to call the function passing the energy and phase only. """ norm = numpy.clip(E / 10., 0., 1.) return norm * (0.5 + 0.25 * numpy.cos(4 * numpy.pi * (phase - 0.25))) #Constant polarization angle pol_ang = constant(numpy.radians(30.)) #Use the xEphemeris class to pass to the xPeriodicPointSource class ephemeris = xEphemeris(0., 1.) src = xPeriodicPointSource('Periodic source', ra, dec, spec, pol_deg, pol_ang, ephemeris) #Add the source to the ROI ROI_MODEL = xROIModel(ra, dec, src)
The configuration file for this particular example is included as one of the so-called toy_models (toy_periodic_source.py) in the examples folder of the ixpeobssim package.
Running a simulation with xpobssim is very simple and all you really need to provide in input is the configuration file and the duration of the simulation. Here is the command line that you need to type (assuming that the current working directory is the ixpeobssim root folder) to run a simulation of the periodic source (described above) for a duration of 1000 seconds:
xpobssim.py --config ixpeobssim/config/toy_periodic_source.py --duration 1000
This will produce three event (FITS) files (one for each Detector Unit) in your $IXPEOBSSIM_DATA folder.
By default the energy range for the simulation is the same as the one where the response functions are defined (currently between 1 and 12 keV). The emin and emax command-line switches allow you to override the default, should that be needed for whatever reason (e.g., if your spectral model is only defined in a smaller energy range). However this feature should be used with caution, as some padding on both ends of the nominal 2–8 keV IXPE energy band is generally needed due to the energy dispersion of the detector.
Once the ixpe observation-simulation has been run, ixpeobssim has a few tools to help in the analysis/manipulation and visualization of the data. These are:
More details on these tools can be found in the Application reference section of the documentation. For the sake of this tutorial we will go through a few examples of how to use these tools.
All the possible binning algorithms that this tool provides are described in the Binned data products and Application reference sections of the documentation. Here we will provide a few examples on what can be done with this tool.
We will go through a few options using the output from the single point source example provided above. Once you have the simulated data of the single point source and you want to create a count map all you need to do is:
xpbin.py $IXPEOBSSIM_DATA/single_point_source_du*.fits --algorithm CMAP
This will create your count map, and if you want to visualize your count map you can use the xpbinview tool. xpbinview is capable of visualizing any of the binned outputs of xpbin without having to specify which binning algorithm you used.
And this is what you will see:
In both of the examples illustrated here, we made use of the wild card (*) to pass all three of the event files to the tool. It is also possible to pass a list of the event files, if you prefer.
To produce a polarization cube, all you need to do is run the xpbin script
PCUBE option to the algorithm argument.
xpbin.py $IXPEOBSSIM_DATA/single_point_source_du*.fits --algorithm PCUBE
If you run the
PCUBE algorithm the default is to bin the data in 4 energy bins.
If you would like to have an
PCUBE in one single energy bin, you can specify it
via this command line option:
xpbin.py $IXPEOBSSIM_DATA/single_point_source_du*.fits --algorithm PCUBE --ebins 1
The xpbin tool also has the
--ebinning option that allows you to pass a list
containing the bin edges for the energy binning you want for your analysis.
There are several other options and binning algorithms that the xpbin script
can do, for a complete list all you need to do make use of the help option or
read the Application reference of the documentation.
As the name implies, this tool allows you to select subsets of your simulated data. Say for example that you just simulated the MSH15-52 nebula + pulsar and that you want to analyze the pulsar. To do that all you need to do is to run the xpeselect script passing the radius of the region around the pulsar and the RA and DEC of the pulsar.
xpselect.py $IXPEOBSSIM_DATA/msh1552_simulated_data_du*.fits --ra RA --dec DEC --rad radius_region
xpeselect also has the options of selecting on energy, time, phase and source id. This last option is for cases where you have more than one source in your ROI.
This particular tool is needed when you are simulating periodic sources, in particular xpphase creates an extra column in the event file for the phase of the periodic source. All this tool needs is the file list (output from the simulation) and the configuration file used to run the simulation.
The following is obsolete and must be edited.
xpphase.py $IXPEOBSSIM_DATA/toy_periodic_source_du*.fits --config ixpeobssim/config/toy_periodic_source.py
This will create three new fits files in your
$IXPEOBSSIM_DATA directory that
are the same as the output of the simulation but now with an extra column for
the phase of the pulsar. The name of the output fits files will be the same as
the input but with an added suffix ‘phase’. You can also specify the suffix for
the output files by including the
--suffix option when running xpophase.
Below is a snapshot of the fits file columns after running the above command.
This section is searching for a writer :-)
Each ixpeobssim application is wrapped in the
module so that it can be called from within a generic Python script with the
exact same arguments that one would pass from command line. Several examples of
how to do this are described in the Analysis pipelines section of the ixpeobssim
documentation. There are also several Python scripts in the example folder of
the ixpeobssim package that illustrate how to perform a full analysis chain
using the pipeline facility.