.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "generated/gallery/loop-bundle-rtv.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_generated_gallery_loop-bundle-rtv.py: Simulating AIA 171 Emission from a Bundle of Strands ===================================================== This example shows how to model AIA emission from a bundle of semi-circular strands for two different viewpoints. .. GENERATED FROM PYTHON SOURCE LINES 9-25 .. code-block:: Python import astropy.time import astropy.units as u import matplotlib.pyplot as plt from astropy.coordinates import SkyCoord from astropy.visualization import quantity_support from sunpy.map import pixelate_coord_path, sample_at_coords import synthesizAR from synthesizAR.instruments.sdo import InstrumentSDOAIA from synthesizAR.interfaces import RTVInterface from synthesizAR.models import semi_circular_bundle # sphinx_gallery_thumbnail_number = -1 .. GENERATED FROM PYTHON SOURCE LINES 26-30 First, we calculate a list of coordinates comprising a semi-circular bundle of strands. The bundle has a length of 50 Mm and a cross-sectional radius of 1 Mm and a total of 500 strands. The strands are uniformly distributed over the bundle. .. GENERATED FROM PYTHON SOURCE LINES 30-35 .. code-block:: Python obstime = astropy.time.Time.now() pos = SkyCoord( lon=0*u.deg, lat=0*u.deg, radius=1*u.AU, obstime=obstime, frame='heliographic_stonyhurst') bundle_coords = semi_circular_bundle(50 * u.Mm, 1*u.Mm, 500, observer=pos) .. GENERATED FROM PYTHON SOURCE LINES 36-38 As in other examples, we then use the coordinates of our strands to construct the `~synthesizAR.Skeleton` object. .. GENERATED FROM PYTHON SOURCE LINES 38-42 .. code-block:: Python strands = [synthesizAR.Strand(f'strand{i}', c) for i, c in enumerate(bundle_coords)] bundle = synthesizAR.Skeleton(strands) bundle.peek(observer=pos) .. image-sg:: /generated/gallery/images/sphx_glr_loop-bundle-rtv_001.png :alt: loop bundle rtv :srcset: /generated/gallery/images/sphx_glr_loop-bundle-rtv_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 43-45 We can also look at our bundle of strands from the side to confirm it has the desired geometry. .. GENERATED FROM PYTHON SOURCE LINES 45-48 .. code-block:: Python side_on_view = SkyCoord(lon=0*u.deg, lat=-90*u.deg, radius=1*u.AU, frame=pos.frame) bundle.peek(observer=side_on_view, grid_kwargs={'grid_spacing': 2*u.deg}) .. image-sg:: /generated/gallery/images/sphx_glr_loop-bundle-rtv_002.png :alt: loop bundle rtv :srcset: /generated/gallery/images/sphx_glr_loop-bundle-rtv_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 49-52 We will again use a simple RTV loop model to compute the thermal structure of each strand. This assigns a single temperature density to the entire loop based on the loop length and the specified heating rate. .. GENERATED FROM PYTHON SOURCE LINES 52-55 .. code-block:: Python rtv = RTVInterface(heating_rate=1e-4*u.Unit('erg cm-3 s-1')) bundle.load_loop_simulations(rtv) .. GENERATED FROM PYTHON SOURCE LINES 56-58 We can then compute the emission as observed by the 171 channel of AIA as viewed from an observer at 1 AU directly above the loop. .. GENERATED FROM PYTHON SOURCE LINES 58-63 .. code-block:: Python aia = InstrumentSDOAIA([0, 1]*u.s, pos, pad_fov=(40, 40)*u.pixel) maps = aia.observe(bundle, channels=aia.channels[2:3]) m_171 = maps['171'][0] m_171.peek() .. image-sg:: /generated/gallery/images/sphx_glr_loop-bundle-rtv_003.png :alt: AIA $171 \; \mathrm{\mathring{A}}$ 2026-04-15 06:12:21 :srcset: /generated/gallery/images/sphx_glr_loop-bundle-rtv_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none INFO: Apparent body location accounts for 499.00 seconds of light travel time [sunpy.coordinates.ephemeris] .. GENERATED FROM PYTHON SOURCE LINES 64-66 Additionally, we can look at intensity profiles along and across the loop axis using `sunpy.map.pixelate_coord_path` and `sunpy.map.sample_at_coords`. .. GENERATED FROM PYTHON SOURCE LINES 66-75 .. code-block:: Python coord_axis = SkyCoord(Tx=[-30, 30]*u.arcsec, Ty=0*u.arcsec, frame=m_171.coordinate_frame) coord_axis = pixelate_coord_path(m_171, coord_axis) profile_axis = sample_at_coords(m_171, coord_axis) coord_xs = SkyCoord(Tx=0*u.arcsec, Ty=[-10, 10]*u.arcsec, frame=m_171.coordinate_frame) coord_xs = pixelate_coord_path(m_171, coord_xs) profile_xs = sample_at_coords(m_171, coord_xs) .. GENERATED FROM PYTHON SOURCE LINES 76-80 Note that the intensity is highest at the footpoints because we are integrating through the most amount of plasma. Additionally, note that the cross-sectional profile has a width consistent with the spatial radius we specified when constructing our bundle of strands. .. GENERATED FROM PYTHON SOURCE LINES 80-92 .. code-block:: Python fig = plt.figure(figsize=(11, 5)) ax = fig.add_subplot(111, projection=m_171) m_171.plot(axes=ax) ax.plot_coord(coord_axis) ax.plot_coord(coord_xs) with quantity_support(): plt.figure(figsize=(11, 5)) plt.subplot(121) plt.plot(coord_axis.separation(coord_axis[0]).to('arcsec'), profile_axis, color='C0') plt.subplot(122) plt.plot(coord_xs.separation(coord_xs[0]).to('arcsec'), profile_xs, color='C1') .. rst-class:: sphx-glr-horizontal * .. image-sg:: /generated/gallery/images/sphx_glr_loop-bundle-rtv_004.png :alt: AIA $171 \; \mathrm{\mathring{A}}$ 2026-04-15 06:12:21 :srcset: /generated/gallery/images/sphx_glr_loop-bundle-rtv_004.png :class: sphx-glr-multi-img * .. image-sg:: /generated/gallery/images/sphx_glr_loop-bundle-rtv_005.png :alt: loop bundle rtv :srcset: /generated/gallery/images/sphx_glr_loop-bundle-rtv_005.png :class: sphx-glr-multi-img .. GENERATED FROM PYTHON SOURCE LINES 93-95 Finally, we can also compute the AIA 171 intensity as viewed from the side in order to see the semi-circular geometry of the loop bundle. .. GENERATED FROM PYTHON SOURCE LINES 95-98 .. code-block:: Python aia = InstrumentSDOAIA([0, 1]*u.s, side_on_view, pad_fov=(40, 40)*u.pixel) maps = aia.observe(bundle, channels=aia.channels[2:3]) maps['171'][0].peek() .. image-sg:: /generated/gallery/images/sphx_glr_loop-bundle-rtv_006.png :alt: AIA $171 \; \mathrm{\mathring{A}}$ 2026-04-15 06:12:21 :srcset: /generated/gallery/images/sphx_glr_loop-bundle-rtv_006.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none INFO: Apparent body location accounts for 499.00 seconds of light travel time [sunpy.coordinates.ephemeris] .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 19.393 seconds) .. _sphx_glr_download_generated_gallery_loop-bundle-rtv.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: loop-bundle-rtv.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: loop-bundle-rtv.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: loop-bundle-rtv.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_