.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "tutorials/analysis-3d/flux_profiles.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_analysis-3d_flux_profiles.py: Flux Profile Estimation ======================= Learn how to estimate flux profiles on a Fermi-LAT dataset. Prerequisites ------------- Knowledge of 3D data reduction and datasets used in Gammapy, see for instance the first analysis tutorial. Context ------- A useful tool to study and compare the spatial distribution of flux in images and data cubes is the measurement of flux profiles. Flux profiles can show spatial correlations of gamma-ray data with e.g. gas maps or other type of gamma-ray data. Most commonly flux profiles are measured along some preferred coordinate axis, either radially distance from a source of interest, along longitude and latitude coordinate axes or along the path defined by two spatial coordinates. Proposed Approach ----------------- Flux profile estimation essentially works by estimating flux points for a set of predefined spatially connected regions. For radial flux profiles the shape of the regions are annuli with a common center, for linear profiles it’s typically a rectangular shape. We will work on a pre-computed `~gammapy.datasets.MapDataset` of Fermi-LAT data, use `~regions.SkyRegion` to define the structure of the bins of the flux profile and run the flux profile extraction using the `~gammapy.estimators.FluxProfileEstimator` .. GENERATED FROM PYTHON SOURCE LINES 37-45 .. code-block:: Python import numpy as np from astropy import units as u from astropy.coordinates import SkyCoord # %matplotlib inline import matplotlib.pyplot as plt .. GENERATED FROM PYTHON SOURCE LINES 46-49 Setup ----- .. GENERATED FROM PYTHON SOURCE LINES 49-55 .. code-block:: Python from IPython.display import display from gammapy.datasets import MapDataset from gammapy.estimators import FluxPoints, FluxProfileEstimator from gammapy.maps import RegionGeom from gammapy.modeling.models import PowerLawSpectralModel .. GENERATED FROM PYTHON SOURCE LINES 56-58 Check setup ----------- .. GENERATED FROM PYTHON SOURCE LINES 58-67 .. code-block:: Python from gammapy.utils.check import check_tutorials_setup from gammapy.utils.regions import ( make_concentric_annulus_sky_regions, make_orthogonal_rectangle_sky_regions, ) check_tutorials_setup() .. rst-class:: sphx-glr-script-out .. code-block:: none System: python_executable : /Users/mregeard/Workspace/dev/code/gammapy/gammapy/.tox/build_docs/bin/python python_version : 3.11.10 machine : x86_64 system : Darwin Gammapy package: version : 1.3.dev1205+g00f44f94ac path : /Users/mregeard/Workspace/dev/code/gammapy/gammapy/.tox/build_docs/lib/python3.11/site-packages/gammapy Other packages: numpy : 1.26.4 scipy : 1.14.1 astropy : 5.2.2 regions : 0.10 click : 8.1.7 yaml : 6.0.2 IPython : 8.28.0 jupyterlab : not installed matplotlib : 3.9.2 pandas : not installed healpy : 1.17.3 iminuit : 2.30.1 sherpa : 4.16.1 naima : 0.10.0 emcee : 3.1.6 corner : 2.2.2 ray : 2.37.0 Gammapy environment variables: GAMMAPY_DATA : /Users/mregeard/Workspace/dev/code/gammapy/gammapy-data/ .. GENERATED FROM PYTHON SOURCE LINES 68-71 Read and Introduce Data ----------------------- .. GENERATED FROM PYTHON SOURCE LINES 71-77 .. code-block:: Python dataset = MapDataset.read( "$GAMMAPY_DATA/fermi-3fhl-gc/fermi-3fhl-gc.fits.gz", name="fermi-dataset" ) .. GENERATED FROM PYTHON SOURCE LINES 78-80 This is what the counts image we will work with looks like: .. GENERATED FROM PYTHON SOURCE LINES 80-85 .. code-block:: Python counts_image = dataset.counts.sum_over_axes() counts_image.smooth("0.1 deg").plot(stretch="sqrt") plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_001.png :alt: flux profiles :srcset: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 86-89 There are 400x200 pixels in the dataset and 11 energy bins between 10 GeV and 2 TeV: .. GENERATED FROM PYTHON SOURCE LINES 89-93 .. code-block:: Python print(dataset.counts) .. rst-class:: sphx-glr-script-out .. code-block:: none WcsNDMap geom : WcsGeom axes : ['lon', 'lat', 'energy'] shape : (400, 200, 11) ndim : 3 unit : dtype : >i4 .. GENERATED FROM PYTHON SOURCE LINES 94-108 Profile Estimation ------------------ Configuration ~~~~~~~~~~~~~ We start by defining a list of spatially connected regions along the galactic longitude axis. For this there is a helper function `~gammapy.utils.regions.make_orthogonal_rectangle_sky_regions`. The individual region bins for the profile have a height of 3 deg and in total there are 31 bins. Its starts from lon = 10 deg and goes to lon = 350 deg. In addition, we have to specify the `wcs` to take into account possible projections effects on the region definition: .. GENERATED FROM PYTHON SOURCE LINES 108-118 .. code-block:: Python regions = make_orthogonal_rectangle_sky_regions( start_pos=SkyCoord("10d", "0d", frame="galactic"), end_pos=SkyCoord("350d", "0d", frame="galactic"), wcs=counts_image.geom.wcs, height="3 deg", nbin=51, ) .. GENERATED FROM PYTHON SOURCE LINES 119-122 We can use the `~gammapy.maps.RegionGeom` object to illustrate the regions on top of the counts image: .. GENERATED FROM PYTHON SOURCE LINES 122-129 .. code-block:: Python geom = RegionGeom.create(region=regions) ax = counts_image.smooth("0.1 deg").plot(stretch="sqrt") geom.plot_region(ax=ax, color="w") plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_002.png :alt: flux profiles :srcset: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none /Users/mregeard/Workspace/dev/code/gammapy/gammapy/.tox/build_docs/lib/python3.11/site-packages/regions/shapes/rectangle.py:205: UserWarning: Setting the 'color' property will override the edgecolor or facecolor properties. return Rectangle(xy=xy, width=width, height=height, .. GENERATED FROM PYTHON SOURCE LINES 130-134 Next we create the `~gammapy.estimators.FluxProfileEstimator`. For the estimation of the flux profile we assume a spectral model with a power-law shape and an index of 2.3 .. GENERATED FROM PYTHON SOURCE LINES 134-143 .. code-block:: Python flux_profile_estimator = FluxProfileEstimator( regions=regions, spectrum=PowerLawSpectralModel(index=2.3), energy_edges=[10, 2000] * u.GeV, selection_optional=["ul"], ) .. GENERATED FROM PYTHON SOURCE LINES 144-146 We can see the full configuration by printing the estimator object: .. GENERATED FROM PYTHON SOURCE LINES 146-150 .. code-block:: Python print(flux_profile_estimator) .. rst-class:: sphx-glr-script-out .. code-block:: none FluxProfileEstimator -------------------- energy_edges : [ 10. 2000.] GeV fit : n_jobs : None n_sigma : 1 n_sigma_ul : 2 norm : Parameter(name='norm', value=1.0, factor=1.0, scale=1.0, unit=Unit(dimensionless), min=nan, max=nan, frozen=False, prior=None, id=0x143e28450) null_value : 0 parallel_backend : None reoptimize : False selection_optional : ['ul'] source : 0 spectrum : PowerLawSpectralModel sum_over_energy_groups : False .. GENERATED FROM PYTHON SOURCE LINES 151-156 Run Estimation ~~~~~~~~~~~~~~ Now we can run the profile estimation and explore the results: .. GENERATED FROM PYTHON SOURCE LINES 158-163 .. code-block:: Python profile = flux_profile_estimator.run(datasets=dataset) print(profile) .. rst-class:: sphx-glr-script-out .. code-block:: none FluxPoints ---------- geom : RegionGeom axes : ['lon', 'lat', 'energy', 'projected-distance'] shape : (1, 1, 1, 51) quantities : ['norm', 'norm_err', 'norm_ul', 'ts', 'npred', 'npred_excess', 'stat', 'stat_null', 'counts', 'success'] ref. model : pl n_sigma : 1 n_sigma_ul : 2 sqrt_ts_threshold_ul : 2 sed type init : likelihood .. GENERATED FROM PYTHON SOURCE LINES 164-173 We can see the flux profile is represented by a `~gammapy.estimators.FluxPoints` object with a `projected-distance` axis, which defines the main axis the flux profile is measured along. The `lon` and `lat` axes can be ignored. Plotting Results ~~~~~~~~~~~~~~~~ Let us directly plot the result using `~gammapy.estimators.FluxPoints.plot`: .. GENERATED FROM PYTHON SOURCE LINES 173-178 .. code-block:: Python ax = profile.plot(sed_type="dnde") ax.set_yscale("linear") plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_003.png :alt: flux profiles :srcset: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 179-183 Based on the spectral model we specified above we can also plot in any other sed type, e.g. energy flux and define a different threshold when to plot upper limits: .. GENERATED FROM PYTHON SOURCE LINES 183-192 .. code-block:: Python profile.sqrt_ts_threshold_ul = 2 plt.figure() ax = profile.plot(sed_type="eflux") ax.set_yscale("linear") plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_004.png :alt: flux profiles :srcset: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_004.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 193-197 We can also plot any other quantity of interest, that is defined on the `~gammapy.estimators.FluxPoints` result object. E.g. the predicted total counts, background counts and excess counts: .. GENERATED FROM PYTHON SOURCE LINES 197-209 .. code-block:: Python quantities = ["npred", "npred_excess", "npred_background"] fig, ax = plt.subplots() for quantity in quantities: profile[quantity].plot(ax=ax, label=quantity.title()) ax.set_ylabel("Counts") plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_005.png :alt: flux profiles :srcset: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_005.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 210-216 Serialisation and I/O ~~~~~~~~~~~~~~~~~~~~~ The profile can be serialised using `~gammapy.estimators.FluxPoints.write`, given a specific format: .. GENERATED FROM PYTHON SOURCE LINES 216-231 .. code-block:: Python profile.write( filename="flux_profile_fermi.fits", format="profile", overwrite=True, sed_type="dnde", ) profile_new = FluxPoints.read(filename="flux_profile_fermi.fits", format="profile") ax = profile_new.plot() ax.set_yscale("linear") plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_006.png :alt: flux profiles :srcset: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_006.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 232-235 The profile can be serialised to a `~astropy.table.Table` object using: .. GENERATED FROM PYTHON SOURCE LINES 235-240 .. code-block:: Python table = profile.to_table(format="profile", formatted=True) display(table) .. rst-class:: sphx-glr-script-out .. code-block:: none x_min x_max x_ref e_ref e_min e_max ... npred_excess stat stat_null is_ul counts success deg deg deg keV keV keV ... ------------------- ------------------ ------------------ ------------ ------------ ------------- ... ------------ --------- --------- ----- ------ ------- -0.1960784313725492 0.1960784313725492 0.0 70710678.119 10000000.000 500000000.000 ... 0.0 0.000 0.000 False 0.0 False 0.1960784313725492 0.5882352941176467 0.392156862745098 70710678.119 10000000.000 500000000.000 ... 0.0 0.000 0.000 False 0.0 False 0.5882352941176467 0.9803921568627443 0.7843137254901955 70710678.119 10000000.000 500000000.000 ... 17.448364 -818.848 -816.789 True 163.0 True 0.9803921568627443 1.3725490196078436 1.176470588235294 70710678.119 10000000.000 500000000.000 ... 192.60901 -3013.446 -2892.855 False 448.0 True 1.3725490196078436 1.764705882352942 1.5686274509803928 70710678.119 10000000.000 500000000.000 ... 320.16006 -4348.966 -4068.317 False 599.0 True 1.764705882352942 2.15686274509804 1.960784313725491 70710678.119 10000000.000 500000000.000 ... 74.95967 -2218.623 -2199.632 False 354.0 True 2.15686274509804 2.549019607843138 2.3529411764705888 70710678.119 10000000.000 500000000.000 ... 66.89662 -2452.717 -2438.503 False 367.0 True 2.549019607843138 2.941176470588236 2.745098039215687 70710678.119 10000000.000 500000000.000 ... 45.943016 -2174.005 -2166.994 False 339.0 True 2.941176470588236 3.3333333333333344 3.137254901960785 70710678.119 10000000.000 500000000.000 ... 247.27599 -3887.531 -3713.385 False 531.0 True 3.3333333333333344 3.725490196078432 3.529411764705883 70710678.119 10000000.000 500000000.000 ... 223.8438 -3183.973 -3014.276 False 458.0 True 3.725490196078432 4.117647058823529 3.9215686274509807 70710678.119 10000000.000 500000000.000 ... 293.2505 -4201.048 -3957.227 False 566.0 True 4.117647058823529 4.509803921568627 4.313725490196078 70710678.119 10000000.000 500000000.000 ... 163.63292 -3040.782 -2958.880 False 448.0 True 4.509803921568627 4.901960784313727 4.7058823529411775 70710678.119 10000000.000 500000000.000 ... 56.33687 -2313.852 -2303.609 False 356.0 True 4.901960784313727 5.294117647058824 5.098039215686276 70710678.119 10000000.000 500000000.000 ... 21.54317 -1979.035 -1977.449 True 315.0 True 5.294117647058824 5.686274509803923 5.4901960784313735 70710678.119 10000000.000 500000000.000 ... -15.065341 -1643.463 -1642.639 True 272.0 True 5.686274509803923 6.07843137254902 5.882352941176471 70710678.119 10000000.000 500000000.000 ... 31.023888 -1692.050 -1688.285 True 282.0 True 6.07843137254902 6.470588235294118 6.274509803921569 70710678.119 10000000.000 500000000.000 ... 9.594084 -2034.334 -2034.033 True 319.0 True 6.470588235294118 6.862745098039216 6.666666666666667 70710678.119 10000000.000 500000000.000 ... -51.3714 -1931.193 -1923.389 True 311.0 True 6.862745098039216 7.254901960784315 7.058823529411765 70710678.119 10000000.000 500000000.000 ... -29.780582 -1861.277 -1858.467 True 303.0 True 7.254901960784315 7.647058823529413 7.450980392156864 70710678.119 10000000.000 500000000.000 ... -60.758495 -1577.587 -1565.087 True 263.0 True 7.647058823529413 8.039215686274511 7.843137254901962 70710678.119 10000000.000 500000000.000 ... -62.7998 -1754.528 -1742.253 True 288.0 True 8.039215686274511 8.43137254901961 8.235294117647062 70710678.119 10000000.000 500000000.000 ... -91.55876 -2098.681 -2076.661 True 328.0 True ... ... ... ... ... ... ... ... ... ... ... ... ... 11.176470588235311 11.568627450980387 11.372549019607849 70710678.119 10000000.000 500000000.000 ... 67.51383 -3256.816 -3246.003 False 471.0 True 11.568627450980387 11.960784313725515 11.764705882352951 70710678.119 10000000.000 500000000.000 ... 53.64453 -3129.927 -3122.952 False 453.0 True 11.960784313725515 12.352941176470619 12.156862745098067 70710678.119 10000000.000 500000000.000 ... 17.3615 -2782.182 -2781.412 True 408.0 True 12.352941176470619 12.745098039215696 12.549019607843157 70710678.119 10000000.000 500000000.000 ... 45.757122 -2877.351 -2871.969 False 425.0 True 12.745098039215696 13.137254901960773 12.941176470588236 70710678.119 10000000.000 500000000.000 ... 51.08866 -2705.136 -2698.118 False 412.0 True 13.137254901960773 13.529411764705875 13.333333333333325 70710678.119 10000000.000 500000000.000 ... 99.729225 -3216.032 -3190.152 False 458.0 True 13.529411764705875 13.921568627451002 13.72549019607844 70710678.119 10000000.000 500000000.000 ... 63.567993 -2403.516 -2391.050 False 373.0 True 13.921568627451002 14.313725490196079 14.11764705882354 70710678.119 10000000.000 500000000.000 ... 39.272743 -2526.278 -2521.871 False 384.0 True 14.313725490196079 14.705882352941181 14.50980392156863 70710678.119 10000000.000 500000000.000 ... -25.449812 -1936.261 -1934.240 True 313.0 True 14.705882352941181 15.098039215686308 14.901960784313744 70710678.119 10000000.000 500000000.000 ... 32.01157 -2178.012 -2174.879 True 357.0 True 15.098039215686308 15.490196078431385 15.294117647058847 70710678.119 10000000.000 500000000.000 ... -8.61883 -1914.901 -1914.661 True 311.0 True 15.490196078431385 15.882352941176464 15.686274509803924 70710678.119 10000000.000 500000000.000 ... 23.800396 -2016.119 -2014.256 True 327.0 True 15.882352941176464 16.274509803921593 16.078431372549026 70710678.119 10000000.000 500000000.000 ... 44.554707 -1957.576 -1950.597 False 321.0 True 16.274509803921593 16.666666666666693 16.470588235294144 70710678.119 10000000.000 500000000.000 ... 79.45122 -2767.621 -2750.068 False 421.0 True 16.666666666666693 17.05882352941177 16.862745098039234 70710678.119 10000000.000 500000000.000 ... 70.58883 -2830.574 -2817.297 False 431.0 True 17.05882352941177 17.4509803921569 17.254901960784338 70710678.119 10000000.000 500000000.000 ... 26.727758 -2356.744 -2354.687 True 374.0 True 17.4509803921569 17.843137254902004 17.647058823529452 70710678.119 10000000.000 500000000.000 ... 18.242876 -2463.364 -2462.413 True 370.0 True 17.843137254902004 18.23529411764708 18.03921568627454 70710678.119 10000000.000 500000000.000 ... 62.1283 -2803.387 -2792.682 False 410.0 True 18.23529411764708 18.627450980392155 18.431372549019617 70710678.119 10000000.000 500000000.000 ... 32.074463 -2186.572 -2183.230 True 336.0 True 18.627450980392155 19.01960784313726 18.823529411764707 70710678.119 10000000.000 500000000.000 ... 0.995932 -877.032 -877.027 True 172.0 True 19.01960784313726 19.411764705882383 19.21568627450982 70710678.119 10000000.000 500000000.000 ... 0.0 0.000 0.000 False 0.0 False 19.411764705882383 19.803921568627494 19.60784313725494 70710678.119 10000000.000 500000000.000 ... 0.0 0.000 0.000 False 0.0 False Length = 51 rows .. GENERATED FROM PYTHON SOURCE LINES 241-244 No we can also estimate a radial profile starting from the Galactic center: .. GENERATED FROM PYTHON SOURCE LINES 244-252 .. code-block:: Python regions = make_concentric_annulus_sky_regions( center=SkyCoord("0d", "0d", frame="galactic"), radius_max="1.5 deg", nbin=11, ) .. GENERATED FROM PYTHON SOURCE LINES 253-255 Again we first illustrate the regions: .. GENERATED FROM PYTHON SOURCE LINES 255-264 .. code-block:: Python geom = RegionGeom.create(region=regions) gc_image = counts_image.cutout( position=SkyCoord("0d", "0d", frame="galactic"), width=3 * u.deg ) ax = gc_image.smooth("0.1 deg").plot(stretch="sqrt") geom.plot_region(ax=ax, color="w") plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_007.png :alt: flux profiles :srcset: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_007.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none /Users/mregeard/Workspace/dev/code/gammapy/gammapy/.tox/build_docs/lib/python3.11/site-packages/regions/core/compound.py:160: UserWarning: Setting the 'color' property will override the edgecolor or facecolor properties. patch = mpatches.PathPatch(path, **mpl_kwargs) .. GENERATED FROM PYTHON SOURCE LINES 265-267 This time we define two energy bins and include the fit statistic profile in the computation: .. GENERATED FROM PYTHON SOURCE LINES 267-274 .. code-block:: Python flux_profile_estimator = FluxProfileEstimator( regions=regions, spectrum=PowerLawSpectralModel(index=2.3), energy_edges=[10, 100, 2000] * u.GeV, selection_optional=["ul", "scan"], ) .. GENERATED FROM PYTHON SOURCE LINES 275-276 The configuration of the fit statistic profile is done throught the norm parameter: .. GENERATED FROM PYTHON SOURCE LINES 276-278 .. code-block:: Python flux_profile_estimator.norm.scan_values = np.linspace(-1, 5, 11) .. GENERATED FROM PYTHON SOURCE LINES 279-280 Now we can run the estimator, .. GENERATED FROM PYTHON SOURCE LINES 280-284 .. code-block:: Python profile = flux_profile_estimator.run(datasets=dataset) .. GENERATED FROM PYTHON SOURCE LINES 285-287 and plot the result: .. GENERATED FROM PYTHON SOURCE LINES 287-292 .. code-block:: Python profile.plot(axis_name="projected-distance", sed_type="flux") plt.show() .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_008.png :alt: flux profiles :srcset: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_008.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 293-296 However because of the powerlaw spectrum the flux at high energies is much lower. To extract the profile at high energies only we can use: .. GENERATED FROM PYTHON SOURCE LINES 296-301 .. code-block:: Python profile_high = profile.slice_by_idx({"energy": slice(1, 2)}) plt.show() .. GENERATED FROM PYTHON SOURCE LINES 302-304 And now plot the points together with the likelihood profiles: .. GENERATED FROM PYTHON SOURCE LINES 304-312 .. code-block:: Python fig, ax = plt.subplots() profile_high.plot(ax=ax, sed_type="eflux", color="tab:orange") profile_high.plot_ts_profiles(ax=ax, sed_type="eflux") ax.set_yscale("linear") plt.show() # sphinx_gallery_thumbnail_number = 2 .. image-sg:: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_009.png :alt: flux profiles :srcset: /tutorials/analysis-3d/images/sphx_glr_flux_profiles_009.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 10.712 seconds) .. _sphx_glr_download_tutorials_analysis-3d_flux_profiles.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/main?urlpath=lab/tree/notebooks/dev/tutorials/analysis-3d/flux_profiles.ipynb :alt: Launch binder :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: flux_profiles.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: flux_profiles.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: flux_profiles.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_