.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "tutorials/details/theta_square_plot.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code or to run this example in your browser via Binder. .. rst-class:: sphx-glr-example-title .. _sphx_glr_tutorials_details_theta_square_plot.py: Make a theta-square plot ======================== This tutorial explains how to make such a plot, that is the distribution of event counts as a function of the squared angular distance, to a test position. .. GENERATED FROM PYTHON SOURCE LINES 13-16 Setup ----- .. GENERATED FROM PYTHON SOURCE LINES 16-27 .. code-block:: Python # %matplotlib inline import matplotlib.pyplot as plt from astropy.coordinates import SkyCoord from astropy import units as u from gammapy.data import DataStore from gammapy.maps import MapAxis from gammapy.makers.utils import ThetaSquaredTable from gammapy.visualization import plot_theta_squared_table .. GENERATED FROM PYTHON SOURCE LINES 28-33 Get some data ------------- Some data taken on the Crab by H.E.S.S. are used. .. GENERATED FROM PYTHON SOURCE LINES 33-38 .. code-block:: Python data_store = DataStore.from_dir("$GAMMAPY_DATA/hess-dl3-dr1") observations = data_store.get_observations([23523, 23526]) .. GENERATED FROM PYTHON SOURCE LINES 39-44 Define a test position ----------------------- Here we define the position of Crab .. GENERATED FROM PYTHON SOURCE LINES 44-49 .. code-block:: Python position = SkyCoord(83.6324, 22.0174, unit="deg") print(position) .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 50-60 Creation of the theta2 plot --------------------------- By default, the distribution of the OFF counts in squared angular distance is calculated from the mirror reflected coordinates of the test position, assuming therefore a single OFF position. However, one can set manually both the coordinates of the ``off_position``. It is worth to note that overlapping regions are always forbidden to avoid correlated OFF counts, therefore the user should take care in the choice of the ``off_position``. .. GENERATED FROM PYTHON SOURCE LINES 60-68 .. code-block:: Python separation = position.separation(observations[0].pointing.fixed_icrs) position_angle = 180 * u.deg off_position = position.directional_offset_by( position_angle=position_angle, separation=separation * 2 ) print(off_position) .. rst-class:: sphx-glr-script-out .. code-block:: none .. GENERATED FROM PYTHON SOURCE LINES 69-82 .. code-block:: Python theta2_axis = MapAxis.from_bounds(0, 0.2, nbin=20, interp="lin", unit="deg2") theta2_table_maker = ThetaSquaredTable( observations=observations, position=position, theta_squared_axis=theta2_axis, position_off=off_position, ) theta2_table = theta2_table_maker.run() plt.figure(figsize=(10, 5)) plot_theta_squared_table(theta2_table) plt.show() .. image-sg:: /tutorials/details/images/sphx_glr_theta_square_plot_001.png :alt: theta square plot :srcset: /tutorials/details/images/sphx_glr_theta_square_plot_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 83-88 Alternatively, it can be requested a number of reflected OFF positions that will be calculated through the `~gammapy.makers.WobbleRegionsFinder`. As mentioned before, the user should be cautious that the regions (ON and OFF) do not overlap, otherwise only the mirror reflected region will be adopted as OFF. .. GENERATED FROM PYTHON SOURCE LINES 88-104 .. code-block:: Python off_regions_number = 3 theta2_axis = MapAxis.from_bounds(0, 0.1, nbin=20, interp="lin", unit="deg2") theta2_table_maker_offreg = ThetaSquaredTable( observations=observations, position=position, theta_squared_axis=theta2_axis, off_regions_number=off_regions_number, ) theta2_table_offreg = theta2_table_maker_offreg.run() plt.figure(figsize=(10, 5)) plot_theta_squared_table(theta2_table_offreg) plt.show() .. image-sg:: /tutorials/details/images/sphx_glr_theta_square_plot_002.png :alt: theta square plot :srcset: /tutorials/details/images/sphx_glr_theta_square_plot_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 105-111 Making a theta2 plot for a given energy range --------------------------------------------- With the function `~gammapy.makers.utils.ThetaSquaredTable`, one can also select a fixed energy range. .. GENERATED FROM PYTHON SOURCE LINES 111-125 .. code-block:: Python theta2_table_maker_en = ThetaSquaredTable( observations=observations, position=position, theta_squared_axis=theta2_axis, energy_edges=[1.2, 11] * u.TeV, ) theta2_table_en = theta2_table_maker_en.run() plt.figure(figsize=(10, 5)) plot_theta_squared_table(theta2_table_en) plt.show() .. image-sg:: /tutorials/details/images/sphx_glr_theta_square_plot_003.png :alt: theta square plot :srcset: /tutorials/details/images/sphx_glr_theta_square_plot_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 126-139 Statistical significance of a detection --------------------------------------- To get the significance of a signal, the usual method consists of using the reflected background method (see the maker tutorial: :doc:`/user-guide/makers/reflected`) to compute the WStat statistics (see :ref:`wstat`, :doc:`/user-guide/stats/fit_statistics`). This is the well-known method of [LiMa1983]_ using ON and OFF regions. The following tutorials show how to get an excess significance: - :doc:`/tutorials/analysis-1d/spectral_analysis` - :doc:`/tutorials/analysis-1d/extended_source_spectral_analysis` .. _sphx_glr_download_tutorials_details_theta_square_plot.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: binder-badge .. image:: images/binder_badge_logo.svg :target: https://mybinder.org/v2/gh/gammapy/gammapy-webpage/v2.1?urlpath=lab/tree/notebooks/2.1/tutorials/details/theta_square_plot.ipynb :alt: Launch binder :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: theta_square_plot.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: theta_square_plot.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: theta_square_plot.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_