.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery_tutorials/read_and_plot_maps.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. .. rst-class:: sphx-glr-example-title .. _sphx_glr_gallery_tutorials_read_and_plot_maps.py: .. _tutorial_read_plot_maps: Read and plot maps with MadcubaMap ################################## Introduction to using `madcubapy` to read and plot FITS maps. Usually, we would open and FITS files in python using `astropy` and plot them using `matplotlib`. This tutorial shows how this can be done with `madcubapy`, which uses the previous packages to do all the work, but in a more concise manner, simplifying this process for the user. We start by importing the necessary libraries for this tutorial. .. GENERATED FROM PYTHON SOURCE LINES 16-24 .. code-block:: Python from madcubapy.io import MadcubaMap from madcubapy.visualization import add_wcs_axes from madcubapy.visualization import add_manual_wcs_axes from madcubapy.visualization import insert_colorbar import matplotlib.pyplot as plt .. GENERATED FROM PYTHON SOURCE LINES 37-49 Reading the FITS map ==================== When MADCUBA exports a FITS map it also creates a :ref:`history file ` with the same name of the FITS file ended with *_hist.csv*. With the `~madcubapy.io.MadcubaMap` class we can open these FITS files alongside their history tables. We can read the FITS file with the `MadcubaMap.read() ` method, and the corresponding :ref:`history file ` will be loaded as well if present. .. GENERATED FROM PYTHON SOURCE LINES 49-53 .. code-block:: Python madcuba_map = MadcubaMap.read( "../examples/data/IRAS16293_SO_2-1_moment0_madcuba.fits") .. GENERATED FROM PYTHON SOURCE LINES 54-63 .. note:: Due to how MADCUBA saves some FITS header cards, several `astropy` warnings can pop up when reading FITS files. Usually these warnings are unexpected card names using non-standard conventions. The `~madcubapy.io.MadcubaMap` class resembles the `~astropy.nddata.CCDData` class from Astropy with ``data``, ``header``, ``wcs``, and ``unit`` attributes, with the addition of a ``hist`` attribute, the history table containing all the information from the :ref:`history file `: .. GENERATED FROM PYTHON SOURCE LINES 63-66 .. code-block:: Python madcuba_map.hist .. raw:: html
Table length=9
IndexFITSMacroTypeFROM_ROWTO_ROWRoi_RangeUserDate
int64str1str500str5int64int64int64str50str23
1C//ORIGIN=/home/madcubapy/data/MAD_CUB_member.uid___A001_X133d_X320d.IRAS_16293-2422_sci.spw29.cube.I.pbcor.fitsI------00
2C//INFO: Import CASA FITS FILE, ESTIMATED_SIGMA=0.005020903444186465P------madcubapy2022-10-14T14:39:05.876
3Crun("Import CASA FITS FILE","select='/home/madcubapy/data/member.uid___A001_X133d_X320d.IRAS_16293-2422_sci.spw29.cube.I.pbcor.fits' changehead='Vrad$2.5$km/s'");P------madcubapy2022-10-14T14:39:05.876
4Crun("Open Virtual Cube","select='/home/madcubapy/data/MAD_CUB_member.uid___A001_X133d_X320d.IRAS_16293-2422_sci.spw29.cube.I.pbcor.fits' ");P------madcubapy2023-05-02T13:12:32.467
5C//ORIGIN=CUBE MAD_CUB_member.uid___A001_X133d_X320d.IRAS_16293-2422_sci.spw29.cube.I.pbcor.fitsI------00
6C//selectWindow("CUBE MAD_CUB_member.uid___A001_X133d_X320d.IRAS_16293-2422_sci.spw29.cube.I.pbcor.fits");I------00
7Crun("INTEGRATED INTENSITY PLUGIN","ranges=500.0$5500.0$# axisunit=m s-1 interpolate=0");P------madcubapy2023-05-02T13:12:34.442
8Crun("Save Cube","select='/home/madcubapy/data/IRAS16293_SO_2-1_moment0_madcuba.fits'");P------madcubapy2023-05-02T13:12:47.239
9C//PYTHON: Open cube: '../examples/data/IRAS16293_SO_2-1_moment0_madcuba.fits'Py------docs2026-03-29T12:55:38.133


.. GENERATED FROM PYTHON SOURCE LINES 67-71 As a failsafe, a `~astropy.nddata.CCDData` object is also present inside the `~madcubapy.io.MadcubaMap` object in the ``ccddata`` attribute. This way, we can do anything a `~astropy.nddata.CCDData` object can with `madcubapy` through this attribute. .. GENERATED FROM PYTHON SOURCE LINES 73-269 Plotting the FITS map ===================== To plot FITS maps we can use the buiilt-in visualization functions of `madcubapy`, or the classical option: plotting the data array with `matplotlib` (and `astropy` for coordinates). .. tab-set:: .. tab-item:: Method 1: Using ``madcubapy`` The :func:`~madcubapy.visualization.add_wcs_axes()` method of `madcubapy` lets the user plot a `~madcubapy.io.MadcubaMap` or `~astropy.nddata.CCDData` quickly using only two statements. A `~matplotlib.figure.Figure` must be set before calling :func:`~madcubapy.visualization.add_wcs_axes`, and it must be the first input parameter of the function. The second mandatory parameter is ``fitsmap``, the map object to plot, which can be a `~madcubapy.io.MadcubaMap` or `~astropy.nddata.CCDData`. Addionally, the user can set a series of ``kwargs`` parameters, which are passed to :func:`matplotlib.pyplot.imshow`, like ``vmin`` or ``vmax``, for example. .. code-block:: python # Create empty figure fig = plt.figure(figsize=(6, 5)) # Use add_wcs_axes() to plot the map. ax, img = add_wcs_axes(fig=fig, fitsmap=madcuba_map, vmin=0, vmax=300) plt.show() .. figure:: read_and_plot_fig_wcsaxes1.png :alt: Image with pcolor :figclass: align-center We can very quickly add a colorbar using :func:`~madcubapy.visualization.add_colorbar` or :func:`~madcubapy.visualization.insert_colorbar`. Check the :ref:`Colorbar page ` to know more about how these two functions work. .. code-block:: python # Create empty figure fig = plt.figure(figsize=(6, 5)) # Use add_wcs_axes() to plot the map. ax, img = add_wcs_axes(fig=fig, fitsmap=madcuba_map, vmin=0, vmax=300) # Add colorbar cbar = insert_colorbar(ax) plt.show() .. figure:: read_and_plot_fig_wcsaxes2.png :alt: Image with pcolor :figclass: align-center The visualization functions of `madcubapy` automatically parse information from the FITS header like the units for the data. The same figure is recreated in Method #2 using the functions of `matplotlib` and `astropy` that `madcubapy` uses under the hood. .. note:: For this tutorial we are calling this function using the parameter names explicitly (i.e ``add_wcs_axes(fig=fig_object, fitsmap=map_object)``), and not by using positional arguments. If one opts to use positional arguments, to quickly plot one map the call to the function should be ``add_wcs_axes(fig_object, 1, 1, 1, map_object)``. This is due to the functionality of :func:`~madcubapy.visualization.add_wcs_axes()` and its manual version :func:`~madcubapy.visualization.add_manual_wcs_axes()` to plot several maps in one figure. This is explained thoroughly in :ref:`Plotting image with coordinates `. .. tab-item:: Method 2: Using ``matplotlib`` and ``astropy`` The classical option: plotting the data array of the FITS files using using `matplotlib` and `astropy`. We can first plot the data array without bothering with any coordinate system with the `matplotlib` function :func:`~matplotlib.pyplot.imshow`, and add the title, axis labels and a colorbar using other `matplotlib` functions. Keep in mind that usually we need to do some slicing to get rid of the spectral and polarization axes if present. .. code-block:: python # Create a figure and plot the image fig, ax = plt.subplots(1, 1, figsize=(6, 6)) map_data = madcuba_map.data[0, 0, :, :] # apply slicing img = ax.imshow(map_data, cmap='viridis', origin='lower', vmin=0, vmax=300) # Add a colorbar fig.colorbar(img) # Set axis labels and title ax.set_title('Madcuba Map') ax.set_xlabel('X axis') ax.set_ylabel('Y axis') plt.show() .. figure:: read_and_plot_fig_mpl1.png :alt: Image with imshow :figclass: align-center This figure uses pixel coordinates for the X and Y axes, and the colorbar does not have the same height as the image, resulting in an ugly astronomical map. We can correct the placement of the colorbar by manually creating a new axes with the size of the image, and use it for the colorbar. .. code-block:: python # New imports needed to create the colorbar import matplotlib.axes as maxes from mpl_toolkits.axes_grid1 import make_axes_locatable # Create a figure and plot the image fig, ax = plt.subplots(1, 1, figsize=(6, 5)) map_data = madcuba_map.data[0, 0, :, :] # apply slicing img = ax.imshow(map_data, cmap='viridis', origin='lower', vmin=0, vmax=300) # Add a colorbar with the same height as the image divider = make_axes_locatable(ax) cax = divider.append_axes("right", size="5%", pad=0.08, axes_class=maxes.Axes) fig.colorbar(img, cax=cax) # Set axis labels and title ax.set_title('Madcuba Map') ax.set_xlabel('X axis') ax.set_ylabel('Y axis') plt.show() .. figure:: read_and_plot_fig_mpl2.png :alt: Image with imshow :figclass: align-center The colorbar is correctly placed, but we still have pixel coordinates. We can use `Astropy `_ to create the plot using the WCS coordinates from the FITS file. .. code-block:: python # Create a figure and add a subplot with WCS projection fig = plt.figure(figsize=(5, 5)) ax = fig.add_subplot(1, 1, 1, projection=madcuba_map.wcs.celestial) map_data = madcuba_map.data[0, 0, :, :] # apply slicing img = ax.imshow(map_data, cmap='viridis', origin='lower', vmin=0, vmax=300) # Add a colorbar with the same height as the image divider = make_axes_locatable(ax) cax = divider.append_axes("right", size="5%", pad=0.05, axes_class=maxes.Axes) cbar = fig.colorbar(img, cax=cax, orientation='vertical') # Set title ax.set_title('Madcuba Map') plt.show() .. figure:: read_and_plot_fig_mpl3.png :alt: Image with imshow :figclass: align-center As we can see, astropy created a `~astropy.visualization.wcsaxes.WCSAxes` object and set the X and Y labels to the correct values read from the FITS file, but their representation looks kind of ugly, and we still do not have a label for the colorbar. We must set the labels manually, and we have to look into the FITS header to get the units of the data, located in the ``BUNIT`` card. For this FITS file, the units of the data are :math:`{\rm Jy \ beam^{-1} \ m \ s^{-1}}`. We can set the labels with: .. code-block:: python # Set correct labels ax.coords[0].set_axislabel("RA (ICRS)") ax.coords[1].set_axislabel("DEC (ICRS)") cbar.set_label(r'$I \ {\rm (Jy \ beam^{-1} \ m \ s^{-1})}$') .. figure:: read_and_plot_fig_mpl4.png :alt: Image with imshow :figclass: align-center .. GENERATED FROM PYTHON SOURCE LINES 271-281 Fix units ========= Some programs export the BUNIT fits card incorrectly with more than one slash and `astropy` has problems parsing the units from those strings. The :meth:`MadcubaMap.fix_units() ` method tries to fix this problem and correctly parse the units. For example, CARTA exports units like this. When we read a CARTA map, `madcubapy` warns us that no :ref:`history file ` has been found, which we can ignore when working with FITS files not created by MADCUBA. .. GENERATED FROM PYTHON SOURCE LINES 281-284 .. code-block:: Python carta_map = MadcubaMap.read("../examples/data/IRAS16293_SO2c_moment0_carta.fits") .. rst-class:: sphx-glr-script-out .. code-block:: none WARNING: Default history file not found. .. GENERATED FROM PYTHON SOURCE LINES 285-289 If we check the units of the map we can see that they are incorrect, because the BUNIT card has a string with multiple slashes: 'Jy/beam.km/s'. When `astropy` tries to parse this, it incorrectly places the :math:`{\rm km}` on the denominator: .. GENERATED FROM PYTHON SOURCE LINES 289-293 .. code-block:: Python print(f"Astropy parsed unit: {carta_map.unit}") print(f"FITS header BUNIT: {carta_map.header["BUNIT"]}") .. rst-class:: sphx-glr-script-out .. code-block:: none Astropy parsed unit: Jy / (beam km s) FITS header BUNIT: Jy/beam.km/s .. GENERATED FROM PYTHON SOURCE LINES 294-299 We can run the :meth:`MadcubaMap.fix_units() ` method of the MadcubaMap object to fix this. It is important to note that the `~astropy.nddata.CCDData` object inside `~madcubapy.io.MadcubaMap` also gets the correct units when they are fixed. .. GENERATED FROM PYTHON SOURCE LINES 299-305 .. code-block:: Python carta_map.fix_units() print(f"MadcubaMap unit: {carta_map.unit}") print(f"CCDData unit: {carta_map.ccddata.unit}") .. rst-class:: sphx-glr-script-out .. code-block:: none MadcubaMap unit: Jy km / (beam s) CCDData unit: Jy km / (beam s) .. GENERATED FROM PYTHON SOURCE LINES 306-316 Even though a CARTA exported map does not have a *_hist.csv* file, it is encouraged to read it as a MadcubaMap, even when the user wants to work with a `~astropy.nddata.CCDData` object if some functionality isn't working with a `~madcubapy.io.MadcubaMap` object. This is the reason why we have the failsafe `~astropy.nddata.CCDData` object in the `MadcubaMap.ccddata` attribute. We can work with it just as if it was an object directly read with the `CCDData.read()` method from `astropy` (because it actually is read this way). .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.080 seconds) .. _sphx_glr_download_gallery_tutorials_read_and_plot_maps.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: read_and_plot_maps.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: read_and_plot_maps.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: read_and_plot_maps.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_