3. Disort

Clearsky radiance

 1import pyarts
 2import numpy as np
 3import matplotlib.pyplot as plt
 4
 5PLOT = False
 6toa = 100e3
 7lat = 0
 8lon = 0
 9NQuad = 40
10
11ws = pyarts.Workspace()
12
13# %% Sampled frequency range
14
15line_f0 = 118750348044.712
16ws.frequency_grid = [line_f0]
17ws.frequency_grid = np.linspace(-20e9, 2e6, 101) + line_f0
18
19# %% Species and line absorption
20
21ws.absorption_speciesSet(species=["O2-66"])
22ws.ReadCatalogData()
23
24# %% Use the automatic agenda setter for propagation matrix calculations
25ws.propagation_matrix_agendaAuto()
26
27# %% Grids and planet
28
29ws.surface_fieldSetPlanetEllipsoid(option="Earth")
30ws.surface_field[pyarts.arts.SurfaceKey("t")] = 295.0
31ws.atmospheric_fieldRead(
32    toa=toa, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
33)
34ws.atmospheric_fieldIGRF(time="2000-03-11 14:39:37")
35# ws.atmospheric_field[pyarts.arts.AtmKey.t] = 300.0
36
37# %% Checks and settings
38
39ws.spectral_radiance_unit = "Tb"
40ws.spectral_radiance_space_agendaSet(option="UniformCosmicBackground")
41ws.spectral_radiance_surface_agendaSet(option="Blackbody")
42
43# %% Core Disort calculations
44
45ws.disort_spectral_radiance_fieldClearsky(
46    longitude=lon,
47    latitude=lat,
48    disort_quadrature_dimension=NQuad,
49    disort_legendre_polynomial_dimension=1,
50    disort_fourier_mode_dimension=1,
51)
52
53# %% Equivalent ARTS calculations
54
55ws.ray_pathGeometricDownlooking(
56    latitude=lat,
57    longitude=lon,
58    max_step=1000.0,
59)
60ws.spectral_radianceClearskyEmission()
61
62# %% Plot results
63
64if PLOT:
65    f = ws.frequency_grid - line_f0
66
67    plt.semilogy(
68        f,
69        ws.disort_spectral_radiance_field[:, 0, 0, : (NQuad // 2)],
70        label="disort",
71    )
72    plt.semilogy(f, ws.spectral_radiance[:, 0], "k--", lw=3)
73    plt.semilogy(
74        f,
75        ws.disort_spectral_radiance_field[:, 0, 0, 0],
76        "g:",
77        lw=3,
78    )
79    plt.semilogy(
80        f,
81        ws.disort_spectral_radiance_field[:, 0, 0, (NQuad // 2) - 1],
82        "m:",
83        lw=3,
84    )
85    plt.ylabel("Spectral radiance [W sr$^{-1}$ m$^{-2}$ Hz$^{-1}$]")
86    plt.xlabel("Dirac frequency [index count]")
87    plt.title("Downlooking")
88
89# %% The last test should be that we are close to the correct values
90
91assert np.allclose(
92    ws.disort_spectral_radiance_field[:, 0, 0, NQuad // 2-1]
93    / ws.spectral_radiance[:, 0],
94    1,
95    rtol=1e-3,
96), "Bad results, clearsky calculations are not close between DISORT and ARTS"

Clearsky flux

 1import pyarts
 2import numpy as np
 3import matplotlib.pyplot as plt
 4
 5toa = 100e3
 6lat = 0
 7lon = 0
 8NQuad = 40
 9
10ws = pyarts.Workspace()
11
12# %% Sampled frequency range
13
14line_f0 = 118750348044.712
15ws.frequency_grid = [line_f0]
16ws.frequency_grid = np.linspace(-20e9, 2e6, 5) + line_f0
17
18# %% Species and line absorption
19
20ws.absorption_speciesSet(species=["O2-66"])
21ws.ReadCatalogData()
22
23# %% Use the automatic agenda setter for propagation matrix calculations
24ws.propagation_matrix_agendaAuto()
25
26# %% Grids and planet
27
28ws.surface_fieldSetPlanetEllipsoid(option="Earth")
29ws.surface_field[pyarts.arts.SurfaceKey("t")] = 295.0
30ws.atmospheric_fieldRead(
31    toa=toa, basename="planets/Earth/afgl/tropical/", missing_is_zero=1
32)
33ws.atmospheric_fieldIGRF(time="2000-03-11 14:39:37")
34# ws.atmospheric_field[pyarts.arts.AtmKey.t] = 300.0
35
36# %% Checks and settings
37
38ws.spectral_radiance_unit = "Tb"
39ws.spectral_radiance_space_agendaSet(option="UniformCosmicBackground")
40ws.spectral_radiance_surface_agendaSet(option="Blackbody")
41
42# %% Core Disort calculations
43
44ws.disort_settings_agendaSet(option="Clearsky")
45
46ws.ray_pathGeometricDownlooking(longitude=lon, latitude=lat, max_step=40_000)
47ws.ray_path_atmospheric_pointFromPath()
48ws.ray_path_frequency_gridFromPath()
49ws.ray_path_propagation_matrixFromPath()
50ws.ray_path_pointLowestFromPath()
51
52ws.disort_spectral_flux_fieldFromAgenda(
53    disort_quadrature_dimension=NQuad,
54    disort_legendre_polynomial_dimension=1,
55    disort_fourier_mode_dimension=1,
56)
57
58assert np.allclose(
59    ws.disort_spectral_flux_field[:, :2].flatten()
60    / np.array(
61        [
62            2.65924430e-15,
63            2.66162733e-15,
64            2.75440515e-15,
65            9.57958182e-18,
66            2.08076327e-17,
67            5.39749082e-16,
68            2.93074072e-15,
69            2.93357370e-15,
70            3.03918211e-15,
71            9.99616600e-18,
72            2.34511632e-17,
73            6.12773186e-16,
74            3.19264874e-15,
75            3.19665668e-15,
76            3.33784028e-15,
77            1.03768102e-17,
78            3.05943372e-17,
79            8.07853130e-16,
80            3.35251230e-15,
81            3.36075470e-15,
82            3.65036247e-15,
83            1.07209025e-17,
84            6.78926741e-17,
85            1.56163327e-15,
86            3.37235023e-15,
87            2.86992821e-15,
88            3.97673151e-15,
89            1.52446216e-15,
90            2.80896759e-15,
91            3.94566396e-15,
92        ]
93    ),
94    1,
95), "Mismatch from historical Disort spectral fluxes"