Zeeman
Zeeman geometry
import os
import matplotlib.pyplot as plt
import numpy as np
import pyarts3 as pyarts
# %% Get magnetic field at a sample position
ell = pyarts.arts.planets.Earth.ellipsoid
pos = [50e3, 0, 0]
mag = pyarts.arts.igrf(pos, ell, t="2000-03-11 14:39:37")
# %% Setup figure for multiple subplots
N = 3
M = 4
fig = plt.figure(figsize=(M * 8, N * 8))
# %% Store computed angles for later verification
angles = []
# %% loop over different LOS directions
for i in range(N):
for j in range(M):
los = [np.linspace(0, 180, N)[i], np.linspace(0, 360, M)[j]]
# Setup 3D subplot
ax = fig.add_subplot(N, M, i*M + j + 1, projection='3d')
# Plot and store angles
ang = pyarts.arts.zeeman.MagneticAngles(mag, los)
pyarts.plots.MagneticAngles.plot(ang, fig=fig, ax=ax, N=50)
ax.set_ylim(-1.2, 1.2)
ax.set_xlim(-1.2, 1.2)
ax.set_zlim(-1.2, 1.2)
angles.append([ang.eta, ang.theta])
if "ARTS_HEADLESS" not in os.environ:
plt.tight_layout()
plt.show()
assert np.allclose(angles,
np.array([[3.02347622, 1.07933684],
[0.92908112, 1.07933684],
[-1.16531399, 1.07933684],
[3.02347622, 1.07933684],
[0.21669967, 0.50432248],
[0.98174228, 2.12671818],
[-1.04334554, 1.92599378],
[0.21669967, 0.50432248],
[0.11811643, 2.06225581],
[2.21251154, 2.06225581],
[-1.97627867, 2.06225581],
[0.11811643, 2.06225581]])), "Angles do not match expected values."
(Source code, svg, pdf)
Zeeman
import os
import matplotlib.pyplot as plt
import numpy as np
import pyarts3 as pyarts
# Download catalogs
pyarts.data.download()
ws = pyarts.workspace.Workspace()
# %% Sampled frequency range
line_f0 = 118750348044.712
ws.freq_grid = np.linspace(-50e6, 50e6, 1001) + line_f0
# %% Species and line absorption
ws.abs_speciesSet(species=["O2-66"])
ws.ReadCatalogData()
ws.abs_bandsSelectFrequencyByLine(fmin=40e9, fmax=120e9)
ws.abs_bandsSetZeeman(species="O2-66", fmin=118e9, fmax=119e9)
ws.WignerInit()
# %% Use the automatic agenda setter for propagation matrix calculations
ws.spectral_propmat_agendaAuto()
# %% Grids and planet
ws.surf_fieldPlanet(option="Earth")
ws.surf_field[pyarts.arts.SurfaceKey("t")] = 295.0
ws.atm_fieldRead(
toa=100e3, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
)
ws.atm_fieldSchmidthFieldFromIGRF(time="2000-03-11 14:39:37")
# %% Checks and settings
ws.spectral_rad_transform_operatorSet(option="Tb")
# %% Core calculations
pos = [100e3, 0, 0]
los = [180.0, 0.0]
ws.ray_pathGeometric(pos=pos, los=los, max_stepsize=1000.0)
ws.rte_option = "constant"
ws.spectral_radClearskyEmission()
ws.spectral_radApplyUnitFromSpectralRadiance()
# %% Show results
fig, ax = pyarts.plot(ws.spectral_rad, freqs=(
ws.freq_grid - line_f0) / 1e6)
[a.set_xlabel("Frequency offset [MHz]") for a in ax.flatten()]
[a.set_ylabel("Spectral radiance [K]") for a in ax.flatten()]
fig.suptitle(f"Zeeman effect of {round(line_f0 / 1e6)} MHz O$_2$ line")
if "ARTS_HEADLESS" not in os.environ:
plt.show()
# %% Test
assert np.allclose(
ws.spectral_rad[::100],
np.array(
[[ 2.27651281e+02, 4.26102083e-04, 1.02747359e-04, 5.68792774e-02],
[ 2.30728865e+02, 6.60200934e-04, 1.59372202e-04, 7.04074781e-02],
[ 2.34671710e+02, 1.16870319e-03, 2.82624430e-04, 9.34017504e-02],
[ 2.40226519e+02, 2.61619541e-03, 6.34833393e-04, 1.40041579e-01],
[ 2.49649966e+02, 9.75456580e-03, 2.39140821e-03, 2.69997496e-01],
[ 2.09913635e+02, 2.41699652e+01, 1.73644537e+00, 5.52423749e-06],
[ 2.49649399e+02, 9.75810766e-03, 2.39230583e-03, -2.70074955e-01],
[ 2.40225428e+02, 2.61800091e-03, 6.35278864e-04, -1.40118138e-01],
[ 2.34670135e+02, 1.16990396e-03, 2.82918050e-04, -9.34777298e-02],
[ 2.30726817e+02, 6.61105654e-04, 1.59592427e-04, -7.04834283e-02],
[ 2.27648771e+02, 4.26832634e-04, 1.02924704e-04, -5.69551120e-02]]
),
), "Values have drifted from expected results in spectral radiance"
(Source code, svg, pdf)
Zeeman sensor
import os
import matplotlib.pyplot as plt
import numpy as np
import pyarts3 as pyarts
# Download catalogs
pyarts.data.download()
ws = pyarts.workspace.Workspace()
ws.rte_option = "constant"
# %% Sampled frequency range
line_f0 = 118750348044.712
ws.freq_grid = np.linspace(-50e6, 50e6, 1001) + line_f0
# %% Species and line absorption
ws.abs_speciesSet(species=["O2-66"])
ws.ReadCatalogData()
ws.abs_bandsSelectFrequencyByLine(fmin=40e9, fmax=120e9)
ws.abs_bandsSetZeeman(species="O2-66", fmin=118e9, fmax=119e9)
ws.WignerInit()
# %% Use the automatic agenda setter for propagation matrix calculations
ws.spectral_propmat_agendaAuto()
# %% Grids and planet
ws.surf_fieldPlanet(option="Earth")
ws.surf_field[pyarts.arts.SurfaceKey("t")] = 295.0
ws.atm_fieldRead(
toa=100e3, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
)
ws.atm_fieldIGRF(time="2000-03-11 14:39:37")
# %% Checks and settings
ws.spectral_rad_transform_operator = "Tb"
ws.ray_path_observer_agendaSetGeometric()
# %% Set up a sensor with Gaussian standard deviation channel widths on individual frequency ranges
pos = [100e3, 0, 0]
los = [180.0, 0.0]
ws.measurement_sensorSimpleGaussian(std=1e5, pos=pos, los=los, pol="RC")
# %% Core calculations
ws.measurement_vecFromSensor()
# %% Show results
fig, ax = pyarts.plot(ws.measurement_vec, xgrid=(
ws.freq_grid - line_f0) / 1e6)
ax.set_xlabel("Frequency offset [MHz]")
ax.set_ylabel("Spectral radiance [K]")
ax.set_title(
f"Zeeman effect of {round(line_f0 / 1e6)} MHz O$_2$ line with Gaussian channels on individual grids"
)
if "ARTS_HEADLESS" not in os.environ:
plt.show()
# %% Test
assert np.allclose(
ws.measurement_vec[::100],
np.array(
[227.72265043, 230.79931459, 234.76518855, 240.36673613,
249.92055137, 211.84612808, 249.3798804 , 240.08547841,
234.5767321 , 230.65637527, 227.60619589]
),
)
(Source code, svg, pdf)
Zeeman transmission
import os
import matplotlib.pyplot as plt
import numpy as np
import pyarts3 as pyarts
# Download catalogs
pyarts.data.download()
ws = pyarts.workspace.Workspace()
# %% Sampled frequency range
line_f0 = 118750348044.712
nf = 1001
ws.freq_grid = np.linspace(-50e6, 50e6, nf) + line_f0
# %% Species and line absorption
ws.abs_speciesSet(species=["O2-66"])
ws.ReadCatalogData()
ws.abs_bandsSelectFrequencyByLine(fmin=40e9, fmax=120e9)
ws.abs_bandsSetZeeman(species="O2-66", fmin=118e9, fmax=119e9)
ws.WignerInit()
# %% Use the automatic agenda setter for propagation matrix calculations
ws.spectral_propmat_agendaAuto()
# %% Grids and planet
ws.surf_fieldPlanet(option="Earth")
ws.surf_field[pyarts.arts.SurfaceKey("t")] = 295.0
ws.atm_fieldRead(
toa=100e3, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
)
ws.atm_fieldIGRF(time="2000-03-11 14:39:37")
# %% Checks and settings
ws.spectral_rad_space_agendaSet(option="Transmission")
ws.spectral_rad_surface_agendaSet(option="Transmission")
# %% Core calculations
pos = [100e3, 0, 0]
los = [180.0, 0.0]
ws.ray_pathGeometric(pos=pos, los=los, max_stepsize=1000.0)
ws.spectral_radClearskyTransmission()
# %% Show results
fig, ax = pyarts.plot(ws.spectral_rad, freqs=(
ws.freq_grid - line_f0) / 1e6, component='I')
ax.set_yscale('log')
ax.set_xlabel("Frequency offset [MHz]")
ax.set_ylabel("Spectral radiance [K]")
ax.set_title(f"Zeeman effect of {round(line_f0 / 1e6)} MHz O$_2$ line")
if "ARTS_HEADLESS" not in os.environ:
plt.show()
# %% Test
assert np.allclose(
ws.spectral_rad[::100],
np.array(
[[ 3.55416974e-06, -1.49626896e-10, -3.59945621e-11, -4.65675860e-08],
[ 1.75504614e-06, -1.14112664e-10, -2.74625108e-11, -2.83292389e-08],
[ 7.16577429e-07, -8.14996162e-11, -1.96228377e-11, -1.51461032e-08],
[ 2.07916647e-07, -5.19869800e-11, -1.25199070e-11, -6.43052931e-09],
[ 2.67056326e-08, -2.58535066e-11, -6.20359298e-12, -1.58692924e-09],
[ 7.87383592e-13, 7.80210615e-13, -6.46005354e-14, 8.84600141e-18],
[ 2.67295620e-08, -2.58964360e-11, -6.21391407e-12, 1.59197636e-09],
[ 2.08240338e-07, -5.21480127e-11, -1.25589127e-11, 6.46545336e-09],
[ 7.18110068e-07, -8.18625471e-11, -1.97105925e-11, 1.52574035e-08],
[ 1.75973088e-06, -1.14770061e-10, -2.76211969e-11, 2.85857722e-08],
[ 3.56540626e-06, -1.50680035e-10, -3.62484502e-11, 4.70615011e-08]]
),
), "Values have drifted from expected results in spectral radiance"
(Source code, svg, pdf)
Zeeman sun
import os
import matplotlib.pyplot as plt
import numpy as np
import pyarts3 as pyarts
# Download catalogs
pyarts.data.download()
ws = pyarts.workspace.Workspace()
ws.rte_option = "constant"
# %% Sampled frequency range
line_f0 = 118750348044.712
ws.freq_grid = [line_f0]
# %% Species and line absorption
ws.abs_speciesSet(species=["O2-66"])
ws.ReadCatalogData()
ws.abs_bandsSelectFrequencyByLine(fmin=40e9, fmax=120e9)
ws.abs_bandsSetZeeman(species="O2-66", fmin=118e9, fmax=119e9)
ws.WignerInit()
# %% Use the automatic agenda setter for propagation matrix calculations
ws.spectral_propmat_agendaAuto()
# %% Grids and planet
ws.surf_fieldPlanet(option="Earth")
ws.surf_field[pyarts.arts.SurfaceKey("t")] = 295.0
ws.atm_fieldRead(
toa=100e3, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
)
ws.atm_fieldIGRF(time="2000-03-11 14:39:37")
# %% Add a sun
ws.sunBlackbody()
ws.suns = [ws.sun]
# %% Checks and settings
ws.spectral_rad_transform_operatorSet(option="Tb")
ws.spectral_rad_space_agendaSet(option="SunOrCosmicBackground")
ws.spectral_rad_surface_agendaSet(option="Blackbody")
# %% Core calculations
pos = [90e3, 0, 0]
zens = np.linspace(0, 2, 21)
azis = np.linspace(-180, 180, 21)
res = np.empty((len(zens), len(azis)))
for izen in range(len(zens)):
for iazi in range(len(azis)):
los = [zens[izen], azis[iazi]]
ws.ray_pathGeometric(pos=pos, los=los, max_stepsize=1000.0)
ws.spectral_radClearskyEmission()
ws.spectral_radApplyUnitFromSpectralRadiance()
res[izen, iazi] = ws.spectral_rad[0][0]
# FIXME: Use some sort of Imager for measurement_vec for the above
r, theta = np.meshgrid(zens, np.rad2deg(azis))
fig, ax = plt.subplots(subplot_kw=dict(projection="polar"))
ax.contourf(theta, r, res.T)
if "ARTS_HEADLESS" not in os.environ:
plt.show()
assert np.allclose(
res[::3, ::7],
np.array(
[[5356.88051075, 5356.88051075, 5356.88051075],
[ 16.77380467, 16.66091577, 16.69389844],
[ 16.84063115, 16.61524188, 16.68151802],
[ 16.90725049, 16.56975243, 16.6696318 ],
[ 16.9736614 , 16.52444918, 16.65824055],
[ 17.03986264, 16.47933394, 16.64734506],
[ 17.10585295, 16.43440858, 16.63694619]])
)
(Source code, svg, pdf)
Zeeman sun scattering
import os
import matplotlib.pyplot as plt
import numpy as np
import pyarts3 as pyarts
# Download catalogs
pyarts.data.download()
ws = pyarts.workspace.Workspace()
# %% Sampled frequency range
ws.freq_grid = [pyarts.arts.convert.wavelen2freq(700e-9)]
# %% Species and line absorption
ws.abs_speciesSet(species=["O2-66"])
ws.ReadCatalogData()
ws.abs_bandsSelectFrequencyByLine(fmin=40e9, fmax=120e9)
ws.abs_bandsSetZeeman(species="O2-66", fmin=118e9, fmax=119e9)
ws.WignerInit()
# %% Use the automatic agenda setter for propagation matrix calculations
ws.spectral_propmat_agendaAuto()
# %% Grids and planet
ws.surf_fieldPlanet(option="Earth")
ws.surf_field[pyarts.arts.SurfaceKey("t")] = 295.0
ws.atm_fieldRead(
toa=100e3, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
)
ws.atm_fieldIGRF(time="2000-03-11 14:39:37")
# %% Add a sun
ws.sunBlackbody()
ws.suns = [ws.sun]
# %% Checks and settings
ws.spectral_rad_transform_operator = "Tb"
ws.spectral_rad_space_agendaSet(option="SunOrCosmicBackground")
ws.spectral_rad_surface_agendaSet(option="Blackbody")
ws.ray_path_observer_agendaSetGeometric()
ws.spectral_propmat_scat_agendaSet(option="AirSimple")
# %% Core calculations
pos = [90e3, 0, 0]
zens = np.linspace(0, 5, 21)
azis = np.linspace(-180, 180, 21)
res = np.empty((len(zens), len(azis)))
for izen in range(len(zens)):
for iazi in range(len(azis)):
los = [zens[izen], azis[iazi]]
ws.ray_pathGeometric(pos=pos, los=los, max_stepsize=1000.0)
ws.ray_path_suns_pathFromPathObserver(just_hit=1)
ws.spectral_radClearskyRayleighScattering()
ws.spectral_radApplyUnitFromSpectralRadiance()
res[izen, iazi] = ws.spectral_rad[0][0]
# FIXME: Use some sort of Imager for measurement_vec for the above
r, theta = np.meshgrid(zens, np.rad2deg(azis))
fig, ax = plt.subplots(subplot_kw=dict(projection="polar"))
ax.contourf(theta, r, res.T)
if "ARTS_HEADLESS" not in os.environ:
plt.show()
assert np.allclose(
res[1::3, 1::7],
np.array([[5771.9999155 , 5771.9999155 , 5771.9999155 ],
[ 642.44289971, 642.44289971, 642.44289971],
[ 642.44286517, 642.44286517, 642.44286524],
[ 642.4428158 , 642.4428158 , 642.44281587],
[ 642.44275627, 642.44275634, 642.44275648],
[ 642.44269338, 642.44269352, 642.44269372],
[ 642.44263584, 642.44263598, 642.44263625]]),
)
(Source code, svg, pdf)
Zeeman refractive
import os
import pyarts3 as pyarts
import numpy as np
import matplotlib.pyplot as plt
# Download catalogs
pyarts.data.download()
ws = pyarts.workspace.Workspace()
# %% Sampled frequency range
line_f0 = 53.0669e9
ws.freq_grid = np.linspace(-15e6, 15e6, 51) + line_f0
# %% Species and line absorption
ws.abs_speciesSet(species=["O2-66"])
ws.ReadCatalogData()
ws.abs_bandsSelectFrequencyByLine(fmin=40e9, fmax=120e9)
ws.abs_bandsSetZeeman(species="O2-66", fmin=52e9, fmax=54e9)
ws.WignerInit()
# %% Use the automatic agenda setter for propagation matrix calculations
ws.spectral_propmat_agendaAuto()
# %% Grids and planet
ws.surf_fieldPlanet(option="Earth")
ws.surf_field[pyarts.arts.SurfaceKey("t")] = 295.0
ws.atm_fieldRead(
toa=120e3, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
)
ws.atm_fieldIGRF()
# %% Checks and settings
ws.spectral_rad_transform_operatorSet(option="Tb")
@pyarts.arts_agenda(ws=ws, fix=True)
def single_propmat_agenda(ws):
ws.single_propmatInit()
ws.single_propmatAddVoigtLTE()
# %% Calculate and compare refractive and geometric ray paths
pos = [3571, 46, 7]
los = [20, 90]
res = []
ws.obs_pos = pos
ws.obs_los = los
ws.max_stepsize = 20000.0
# %% Show results
ws.ray_point_back_propagation_agendaSet(option="GeometricStepwise")
ws.spectral_radClearskyEmissionFrequencyDependentPropagation(max_tau=1e-2)
ws.spectral_radApplyUnitFromSpectralRadiance(ray_path=ws.spectral_ray_path[0])
geometric = ws.spectral_rad * 1.0
ws.ray_point_back_propagation_agendaSet(option="RefractiveStepwise")
ws.spectral_radClearskyEmissionFrequencyDependentPropagation(max_tau=1e-2)
ws.spectral_radApplyUnitFromSpectralRadiance(ray_path=ws.spectral_ray_path[0])
refractive = ws.spectral_rad * 1.0
# %% Show results
if "ARTS_HEADLESS" not in os.environ:
fig, ax = plt.subplots(2, 2, figsize=(10, 8))
freqs = ws.freq_grid / 1e9
pyarts.plots.StokvecVector.plot(
geometric - refractive, fig=fig, ax=ax, freqs=freqs)
for a in ax.flatten():
a.set_xlabel("Frequency offset [MHz]")
a.set_ylabel("Spectral radiance [K]")
ax[0, 0].set_title("Stokes I")
ax[0, 1].set_title("Stokes Q")
ax[1, 0].set_title("Stokes U")
ax[1, 1].set_title("Stokes V")
fig.suptitle(
f"Difference between geometric and refractive paths around the {round(line_f0 / 1e6)} MHz O$_2$ line")
fig.tight_layout()
plt.show()
# %% Tests
assert np.allclose(geometric.flatten()[::21], [1.25486103e+02, -5.56777144e-04, -2.63816009e-03, 6.12431001e-01,
1.38429184e+02, 1.04834522e-01, -2.21549573e-02, -5.18205673e-01,
1.29365186e+02, -4.48679095e-04])
assert np.allclose(refractive.flatten()[::21], [1.25552512e+02, -5.61580828e-04, -2.61219720e-03, 6.11932940e-01,
1.38488011e+02, 1.04448633e-01, -2.18946938e-02, -5.19985291e-01,
1.29439037e+02, -4.53301859e-04])
assert not np.allclose(geometric, refractive)
Zeeman adaptive path
import os
import matplotlib.pyplot as plt
import numpy as np
import pyarts3 as pyarts
from copy import deepcopy as copy
# Download catalogs
pyarts.data.download()
ws = pyarts.workspace.Workspace()
# %% Sampled frequency range
line_f0 = 118750348044.712
ws.freq_grid = np.linspace(-50e6, 50e6, 1001) + line_f0
# %% Species and line absorption
ws.abs_speciesSet(species=["O2-66"])
ws.ReadCatalogData()
ws.abs_bandsSelectFrequencyByLine(fmin=40e9, fmax=120e9)
ws.abs_bandsSetZeeman(species="O2-66", fmin=118e9, fmax=119e9)
ws.WignerInit()
# %% Use the automatic agenda setter for propagation matrix calculations
ws.spectral_propmat_agendaAuto()
# %% Grids and planet
ws.surf_fieldPlanet(option="Earth")
ws.surf_field[pyarts.arts.SurfaceKey("t")] = 295.0
ws.atm_fieldRead(
toa=100e3, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
)
ws.atm_fieldSchmidthFieldFromIGRF(time="2000-03-11 14:39:37")
# %% Checks and settings
ws.spectral_rad_transform_operatorSet(option="Tb")
# %% Core calculations
pos = [100e3, 0, 0]
los = [180.0, 0.0]
ws.ray_pathGeometric(pos=pos, los=los, max_stepsize=1000.0)
ws.ray_pointBackground()
ws.spectral_rad_bkgAgendasAtEndOfPath()
ws.atm_pathFromPath()
ws.freq_grid_pathFromPath()
ws.spectral_propmat_pathFromPath()
ws.spectral_radSetToBackground()
ws.spectral_radSinglePathEmissionFrequencyLoop()
ws.spectral_radApplyUnitFromSpectralRadiance()
srad0 = copy(ws.spectral_rad)
path0 = copy(ws.ray_path)
ws.spectral_propmat_pathAddAdaptiveHalfPath(
max_stepsize=100., max_tau=0.05, cutoff_tau=3.0)
ws.spectral_radSetToBackground()
ws.spectral_radSinglePathEmissionFrequencyLoop()
ws.spectral_radApplyUnitFromSpectralRadiance()
freqs = (ws.freq_grid - line_f0) / 1e6
f, a = pyarts.plot(srad0, freqs=freqs, label="Regular path")
[a.set_xlabel("Frequency offset [MHz]") for a in a.flatten()]
[a.set_ylabel("Spectral radiance [K]") for a in a.flatten()]
f.suptitle(f"Zeeman effect of {round(line_f0 / 1e6)} MHz O$_2$ line")
pyarts.plot(ws.spectral_rad, freqs=freqs, fig=f, ax=a, label="Adaptive path")
[a.legend() for a in a.flatten()]
fig, ax = plt.subplots(1, 1)
d0 = path0.distances(ws.surf_field.ellipsoid)
a0 = [p.pos[0] / 1e3 for p in path0]
ax.plot(d0, a0[:-1], label="Regular path")
d1 = ws.ray_path.distances(ws.surf_field.ellipsoid)
a1 = [p.pos[0] / 1e3 for p in ws.ray_path]
ax.plot(d1, a1[:-1], label="Adaptive path")
ax.legend()
ax.set_xlabel("Distance along path [m]")
ax.set_ylabel("Altitude [km]")
fig.suptitle("Ray paths")
if "ARTS_HEADLESS" not in os.environ:
plt.show()
assert len(ws.ray_path) > len(path0) # More points in adaptive path
assert np.allclose([np.sum(d0)], [np.sum(d1)])
assert np.allclose(sorted(a0, reverse=True), a0)
assert np.allclose(sorted(a1, reverse=True), a1)
Zeeman linear src
import os
import matplotlib.pyplot as plt
import numpy as np
import pyarts3 as pyarts
# Download catalogs
pyarts.data.download()
ws = pyarts.workspace.Workspace()
# %% Sampled frequency range
line_f0 = 118750348044.712
ws.freq_grid = np.linspace(-50e6, 50e6, 1001) + line_f0
# %% Species and line absorption
ws.abs_speciesSet(species=["O2-66"])
ws.ReadCatalogData()
ws.abs_bandsSelectFrequencyByLine(fmin=40e9, fmax=120e9)
ws.abs_bandsSetZeeman(species="O2-66", fmin=118e9, fmax=119e9)
ws.WignerInit()
# %% Use the automatic agenda setter for propagation matrix calculations
ws.spectral_propmat_agendaAuto()
# %% Grids and planet
ws.surf_fieldPlanet(option="Earth")
ws.surf_field[pyarts.arts.SurfaceKey("t")] = 295.0
ws.atm_fieldRead(
toa=100e3, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
)
ws.atm_fieldSchmidthFieldFromIGRF(time="2000-03-11 14:39:37")
# %% Checks and settings
ws.spectral_rad_transform_operatorSet(option="Tb")
# %% Core calculations
pos = [100e3, 0, 0]
los = [180.0, 0.0]
ws.ray_pathGeometric(pos=pos, los=los, max_stepsize=1000.0)
ws.rte_option = "lintau"
ws.spectral_radClearskyEmission()
ws.spectral_radApplyUnitFromSpectralRadiance()
# %% Show results
fig, ax = pyarts.plot(ws.spectral_rad, freqs=(
ws.freq_grid - line_f0) / 1e6)
[a.set_xlabel("Frequency offset [MHz]") for a in ax.flatten()]
[a.set_ylabel("Spectral radiance [K]") for a in ax.flatten()]
fig.suptitle(f"Zeeman effect of {round(line_f0 / 1e6)} MHz O$_2$ line")
if "ARTS_HEADLESS" not in os.environ:
plt.show()
# %% Test
assert np.allclose(
ws.spectral_rad[::100],
np.array(
[[ 2.27693775e+02, 4.25360933e-04, 1.02565822e-04, 5.68301097e-02],
[ 2.30768818e+02, 6.59206212e-04, 1.59126898e-04, 7.03542474e-02],
[ 2.34709071e+02, 1.16721699e-03, 2.82253340e-04, 9.33505189e-02],
[ 2.40261444e+02, 2.61334904e-03, 6.34105246e-04, 1.39994110e-01],
[ 2.49682208e+02, 9.74388124e-03, 2.38851950e-03, 2.69924367e-01],
[ 2.09884932e+02, 2.41682141e+01, 1.73723873e+00, 5.54512081e-06],
[ 2.49681641e+02, 9.74742248e-03, 2.38941674e-03, -2.70001907e-01],
[ 2.40260352e+02, 2.61515400e-03, 6.34550524e-04, -1.40070740e-01],
[ 2.34707495e+02, 1.16841714e-03, 2.82546778e-04, -9.34265633e-02],
[ 2.30766769e+02, 6.60110278e-04, 1.59346947e-04, -7.04302732e-02],
[ 2.27691263e+02, 4.26090913e-04, 1.02743016e-04, -5.69060525e-02]]
),
), "Values have drifted from expected results in spectral radiance"
(Source code, svg, pdf)
Zeeman linear prop
import os
import matplotlib.pyplot as plt
import numpy as np
import pyarts3 as pyarts
# Download catalogs
pyarts.data.download()
ws = pyarts.workspace.Workspace()
# %% Sampled frequency range
line_f0 = 118750348044.712
ws.freq_grid = np.linspace(-50e6, 50e6, 1001) + line_f0
# %% Species and line absorption
ws.abs_speciesSet(species=["O2-66"])
ws.ReadCatalogData()
ws.abs_bandsSelectFrequencyByLine(fmin=40e9, fmax=120e9)
ws.abs_bandsSetZeeman(species="O2-66", fmin=118e9, fmax=119e9)
ws.WignerInit()
# %% Use the automatic agenda setter for propagation matrix calculations
ws.spectral_propmat_agendaAuto()
# %% Grids and planet
ws.surf_fieldPlanet(option="Earth")
ws.surf_field[pyarts.arts.SurfaceKey("t")] = 295.0
ws.atm_fieldRead(
toa=100e3, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
)
ws.atm_fieldSchmidthFieldFromIGRF(time="2000-03-11 14:39:37")
# %% Checks and settings
ws.spectral_rad_transform_operatorSet(option="Tb")
# %% Core calculations
pos = [100e3, 0, 0]
los = [180.0, 0.0]
ws.ray_pathGeometric(pos=pos, los=los, max_stepsize=1000.0)
ws.rte_option = "linprop"
ws.spectral_radClearskyEmission()
ws.spectral_radApplyUnitFromSpectralRadiance()
# %% Show results
fig, ax = pyarts.plot(ws.spectral_rad, freqs=(
ws.freq_grid - line_f0) / 1e6)
[a.set_xlabel("Frequency offset [MHz]") for a in ax.flatten()]
[a.set_ylabel("Spectral radiance [K]") for a in ax.flatten()]
fig.suptitle(f"Zeeman effect of {round(line_f0 / 1e6)} MHz O$_2$ line")
if "ARTS_HEADLESS" not in os.environ:
plt.show()
# %% Test
assert np.allclose(
ws.spectral_rad[::100],
np.array(
[[ 2.27720431e+02, 1.92860863e-04, 4.80226318e-05, -2.39944534e-02],
[ 2.30795788e+02, 2.05916096e-04, 5.28115272e-05, -4.11877682e-02],
[ 2.34732775e+02, 6.07407437e-04, 1.53649762e-04, -3.35111787e-02],
[ 2.40285632e+02, 1.20481776e-03, 3.16497792e-04, -6.37947118e-02],
[ 2.49706649e+02, 5.90915798e-03, 1.64151310e-03, -9.03484330e-02],
[ 2.09884932e+02, 2.41682141e+01, 1.73723873e+00, 5.54512226e-06],
[ 2.49706087e+02, 5.90803897e-03, 1.64139792e-03, 9.04789486e-02],
[ 2.40284543e+02, 1.20399960e-03, 3.16336925e-04, 6.39391756e-02],
[ 2.34731202e+02, 6.06980643e-04, 1.53562668e-04, 3.36445494e-02],
[ 2.30793741e+02, 2.05389194e-04, 5.26935796e-05, 4.13367589e-02],
[ 2.27717921e+02, 1.92557235e-04, 4.79553616e-05, 2.41329191e-02]]
),
), "Values have drifted from expected results in spectral radiance"
(Source code, svg, pdf)
