Note
Go to the end to download the full example code.
Magnetic field along a field line#
How to get the value of the magnetic field along a field line traced through the PFSS solution.
import matplotlib.pyplot as plt
import astropy.constants as const
import astropy.units as u
from astropy.coordinates import SkyCoord
import sunpy.map
from sunkit_magex import pfss
Load a GONG magnetic field map.
gong_fname = pfss.sample_data.get_gong_map()
gong_map = sunpy.map.Map(gong_fname)
Files Downloaded: 0%| | 0/1 [00:00<?, ?file/s]
sunkit_magex.mrzqs200901t1304c2234_022.fits.gz: 0%| | 0.00/242k [00:00<?, ?B/s]
sunkit_magex.mrzqs200901t1304c2234_022.fits.gz: 0%| | 1.02k/242k [00:00<00:23, 10.2kB/s]
Files Downloaded: 100%|██████████| 1/1 [00:00<00:00, 3.57file/s]
Files Downloaded: 100%|██████████| 1/1 [00:00<00:00, 3.56file/s]
The PFSS solution is calculated on a regular 3D grid in (phi, s, rho), where rho = ln(r), and r is the standard spherical radial coordinate. We need to define the number of rho grid points, and the source surface radius.
From the boundary condition, number of radial grid points, and source surface, we now construct an Input object that stores this information
Now take a seed point, and trace a magnetic field line through the PFSS solution from this point
INFO: Missing metadata for solar radius: assuming the standard radius of the photosphere. [sunpy.map.mapbase]
From this field line we can extract the coordinates, and then use
.b_along_fline
to get the components of the magnetic field along the
field line.
From the plot we can see that the non-radial component of the mangetic field goes to zero at the source surface, as expected.
field_line = field_lines[0]
B = field_line.b_along_fline
r = field_line.coords.radius
fig, ax = plt.subplots()
ax.plot(r.to(const.R_sun), B[:, 0], label=r'$B_{r}$')
ax.plot(r.to(const.R_sun), B[:, 1], label=r'$B_{\theta}$')
ax.plot(r.to(const.R_sun), B[:, 2], label=r'$B_{\phi}$')
ax.legend()
ax.set_xlabel(r'r / $R_{\odot}$')
ax.set_ylabel(f'B ({B.unit})')
plt.show()
Total running time of the script: (0 minutes 3.737 seconds)