POSEIDON.emission ================= .. py:module:: POSEIDON.emission .. autoapi-nested-parse:: Radiative transfer calculations for generating emission spectra. Attributes ---------- .. autoapisummary:: POSEIDON.emission.cp POSEIDON.emission.block POSEIDON.emission.thread Functions --------- .. autoapisummary:: POSEIDON.emission.find_nearest POSEIDON.emission.planck_lambda_arr POSEIDON.emission.planck_lambda_arr_GPU POSEIDON.emission.emission_single_stream POSEIDON.emission.emission_single_stream_w_albedo POSEIDON.emission.emission_single_stream_GPU POSEIDON.emission.determine_photosphere_radii POSEIDON.emission.determine_photosphere_radii_GPU POSEIDON.emission.slice_gt POSEIDON.emission.setup_tri_diag POSEIDON.emission.tri_diag_solve POSEIDON.emission.emission_Toon POSEIDON.emission.numba_cumsum POSEIDON.emission.reflection_Toon POSEIDON.emission.emission_bare_surface POSEIDON.emission.reflection_bare_surface POSEIDON.emission.assign_assumptions_and_compute_single_stream_emission POSEIDON.emission.assign_assumptions_and_compute_thermal_scattering_emission POSEIDON.emission.assign_assumptions_and_compute_reflection Module Contents --------------- .. py:data:: cp .. py:data:: block .. py:data:: thread .. py:function:: find_nearest(array, value) .. py:function:: planck_lambda_arr(T, wl) Compute the Planck function spectral radiance for a range of model wavelengths and atmospheric temperatures. :param T: Array of temperatures in each atmospheric layer (K). :type T: np.array of float :param wl: Wavelength grid (μm). :type wl: np.array of float :returns: Planck function spectral radiance as a function of layer temperature and wavelength in SI units (W/m^2/sr/m). :rtype: B_lambda (2D np.array of float) .. py:function:: planck_lambda_arr_GPU(T, wl, B_lambda) GPU variant of the 'planck_lambda_arr' function. Compute the Planck function spectral radiance for a range of model wavelengths and atmospheric temperatures. :param T: Array of temperatures in each atmospheric layer (K). :type T: np.array of float :param wl: Wavelength grid (μm). :type wl: np.array of float :returns: Planck function spectral radiance as a function of layer temperature and wavelength in SI units (W/m^2/sr/m). :rtype: B_lambda (2D np.array of float) .. py:function:: emission_single_stream(T, dz, wl, kappa, Gauss_quad=2) Compute the emergent top-of-atmosphere flux from a planet or brown dwarf. This function considers only pure thermal emission (i.e. no scattering). :param T: Temperatures in each atmospheric layer (K). :type T: np.array of float :param dz: Vertical extent of each atmospheric layer (m). :type dz: np.array of float :param wl: Wavelength grid (μm). :type wl: np.array of float :param kappa: Extinction coefficient in each layer as a function of wavelength (m^-1). :type kappa: 2D np.array of float :param Gauss_quad: Gaussian quadrature order for integration over emitting surface (Options: 2 / 3). :type Gauss_quad: int :returns: Spectral surface flux in SI units (W/m^2/sr/m). :rtype: F (np.array of float) .. py:function:: emission_single_stream_w_albedo(T, dz, wl, kappa, Gauss_quad=2, surf_reflect=[], index_below_P_surf=0) Compute the emergent top-of-atmosphere flux from a planet or brown dwarf. This function considers only pure thermal emission (i.e. no scattering). :param T: Temperatures in each atmospheric layer (K). :type T: np.array of float :param dz: Vertical extent of each atmospheric layer (m). :type dz: np.array of float :param wl: Wavelength grid (μm). :type wl: np.array of float :param kappa: Extinction coefficient in each layer as a function of wavelength (m^-1). :type kappa: 2D np.array of float :param Gauss_quad: Gaussian quadrature order for integration over emitting surface (Options: 2 / 3). :type Gauss_quad: int :param surf_reflect: numpy.ndarray Surface reflectivity as a function of wavelength. :param index_below_P_surf: int Index below P_surf, so that the blackbody can be computed. (Note that P_surf can be between two pressure levels, we take the lower one) :returns: Spectral surface flux in SI units (W/m^2/sr/m). :rtype: F (np.array of float) .. py:function:: emission_single_stream_GPU(T, dz, wl, kappa, Gauss_quad=2) GPU variant of the 'emission_rad_transfer' function. Compute the emergent top-of-atmosphere flux from a planet or brown dwarf. This function considers only pure thermal emission (i.e. no scattering). :param T: Temperatures in each atmospheric layer (K). :type T: np.array of float :param dz: Vertical extent of each atmospheric layer (m). :type dz: np.array of float :param wl: Wavelength grid (μm). :type wl: np.array of float :param kappa: Extinction coefficient in each layer as a function of wavelength (m^-1). :type kappa: 2D np.array of float :param Gauss_quad: Gaussian quadrature order for integration over emitting surface (Options: 2 / 3). :type Gauss_quad: int :returns: Spectral surface flux in SI units (W/m^2/sr/m). :rtype: F (np.array of float) .. py:function:: determine_photosphere_radii(dtau, r_low, wl, photosphere_tau=2 / 3) Interpolate optical depth to find the radius corresponding to the photosphere (by default at tau = 2/3). :param dtau: Vertical optical depth across each layer (starting from the top of the atmosphere) as a function of layer and wavelength. :type dtau: 2D np.array of float :param r_low: Radius at the lower boundary of each layer (m). :type r_low: np.array of float :param wl: Wavelength grid (μm). :type wl: np.array of float :param photosphere_tau: Optical depth to determine photosphere radius. :type photosphere_tau: float :returns: Photosphere radius as a function of wavelength (m). :rtype: R_p_eff (np.array of float) .. py:function:: determine_photosphere_radii_GPU(tau_lambda, r_low, wl, R_p_eff, photosphere_tau=2 / 3) GPU variant of the 'determine_photosphere_radii' function. Interpolate optical depth to find the radius corresponding to the photosphere (by default at tau = 2/3). :param dtau: Vertical optical depth across each layer (starting from the top of the atmosphere) as a function of layer and wavelength. :type dtau: 2D np.array of float :param r_low: Radius at the lower boundary of each layer (m). :type r_low: np.array of float :param wl: Wavelength grid (μm). :type wl: np.array of float :param photosphere_tau: Optical depth to determine photosphere radius. :type photosphere_tau: float :returns: Photosphere radius as a function of wavelength (m). :rtype: R_p_eff (np.array of float) .. py:function:: slice_gt(array, lim) Function to replace values with upper or lower limit .. py:function:: setup_tri_diag(N_layer, N_wl, c_plus_up, c_minus_up, c_plus_down, c_minus_down, b_top, b_surface, surf_reflect, gamma, dtau, exptrm_positive, exptrm_minus) Before we can solve the tridiagonal matrix (See Toon+1989) section "SOLUTION OF THE TwO-STREAM EQUATIONS FOR MULTIPLE LAYERS", we need to set up the coefficients. :param N_layer: number of layers in the model :type N_layer: int :param N_wl: number of wavelength points :type N_wl: int :param c_plus_up: c-plus evaluated at the top of the atmosphere :type c_plus_up: array :param c_minus_up: c_minus evaluated at the top of the atmosphere :type c_minus_up: array :param c_plus_down: c_plus evaluated at the bottom of the atmosphere :type c_plus_down: array :param c_minus_down: c_minus evaluated at the bottom of the atmosphere :type c_minus_down: array :param b_top: The diffuse radiation into the model at the top of the atmosphere :type b_top: array :param b_surface: The diffuse radiation into the model at the bottom. Includes emission, reflection of the unattenuated portion of the direct beam :type b_surface: array :param surf_reflect: Surface reflectivity :type surf_reflect: array :param g1: table 1 toon et al 1989 :type g1: array :param g2: table 1 toon et al 1989 :type g2: array :param g3: table 1 toon et al 1989 :type g3: array :param lamba: Eqn 21 toon et al 1989 :type lamba: array :param gamma: Eqn 22 toon et al 1989 :type gamma: array :param dtau: Opacity per layer :type dtau: array :param exptrm_positive: Eqn 44, exponential terms needed for tridiagonal rotated layered, clipped at 35 :type exptrm_positive: array :param exptrm_minus: Eqn 44, exponential terms needed for tridiagonal rotated layered, clipped at 35 :type exptrm_minus: array :returns: coefficient of the positive exponential term :rtype: array .. py:function:: tri_diag_solve(l, a, b, c, d) ' Tridiagonal Matrix Algorithm solver, a b c d can be NumPy array type or Python list type. refer to this wiki_ and to this explanation_. .. _wiki: http://en.wikipedia.org/wiki/Tridiagonal_matrix_algorithm .. _explanation: http://www.cfd-online.com/Wiki/Tridiagonal_matrix_algorithm_-_TDMA_(Thomas_algorithm) A, B, C and D refer to: .. math:: A(I)*X(I-1) + B(I)*X(I) + C(I)*X(I+1) = D(I) This solver returns X. :param A: :type A: array or list :param B: :type B: array or list :param C: :type C: array or list :param C: :type C: array or list :returns: Solution, x :rtype: array .. py:function:: emission_Toon(P, T, wl, dtau_tot, kappa_Ray, kappa_cloud, kappa_tot, w_cloud, g_cloud, zone_idx, surf_reflect, kappa_cloud_seperate, hard_surface=0, tridiagonal=0, Gauss_quad=5, numt=1, T_surf=0) This function uses the source function method, which is outlined here : https://agupubs.onlinelibrary.wiley.com/doi/pdf/10.1029/JD094iD13p16287 The result of this routine is the top of the atmosphere thermal flux as a function of gauss and chebychev points across the disk. Everything here is in CGS units: Fluxes - erg/s/cm^3 Temperature - K Wave grid - cm-1 Pressure ; dyne/cm2 Reminder: Flux = pi * Intensity, so if you are trying to compare the result of this with a black body you will need to compare with pi * BB ! Adapted from get_thermal_1d :param nlevel: Number of levels which occur at the grid points (not to be confused with layers which are mid points) :type nlevel: int :param wno: Wavenumber grid in inverse cm :type wno: numpy.ndarray :param N_wl: Number of wavenumber points :type N_wl: int :param numg: Number of gauss points (think longitude points) :type numg: int :param numt: Number of chebychev points (think latitude points) :type numt: int :param tlevel: Temperature as a function of level (not layer) :type tlevel: numpy.ndarray :param dtau: This is a matrix of nlayer by nwave. This describes the per layer optical depth. :type dtau: numpy.ndarray :param w0: This is a matrix of nlayer by nwave. This describes the single scattering albedo of the atmosphere. Note this is free of any Raman scattering or any d-eddington correction that is sometimes included in reflected light calculations. :type w0: numpy.ndarray :param cosb: This is a matrix of nlayer by nwave. This describes the asymmetry of the atmosphere. Note this is free of any Raman scattering or any d-eddington correction that is sometimes included in reflected light calculations. :type cosb: numpy.ndarray :param plevel: Pressure for each level (not layer, which is midpoints). CGS units (dyne/cm2) :type plevel: numpy.ndarray :param ubar1: This is a matrix of ng by nt. This describes the outgoing incident angles and is generally computed in `picaso.disco` :type ubar1: numpy.ndarray :param surf_reflect: Surface reflectivity as a function of wavenumber. :type surf_reflect: numpy.ndarray :param hard_surface: 0 for no hard surface (e.g. Jupiter/Neptune), 1 for hard surface (terrestrial) :type hard_surface: int :param tridiagonal: 0 for tridiagonal, 1 for pentadiagonal :type tridiagonal: int :returns: Thermal flux in CGS units (erg/cm3/s) in a matrix that is numg x numt x nwno :rtype: numpy.ndarray .. py:function:: numba_cumsum(mat) Function to compute cumsum along axis=0 to bypass numba not allowing kwargs in cumsum .. py:function:: reflection_Toon(P, wl, dtau_tot, kappa_Ray, kappa_cloud, kappa_tot, w_cloud, g_cloud, zone_idx, surf_reflect, kappa_cloud_seperate, single_phase=3, multi_phase=0, frac_a=1, frac_b=-1, frac_c=2, constant_back=-0.5, constant_forward=1, Gauss_quad=5, numt=1, toon_coefficients=0, tridiagonal=0, b_top=0) Computes toon fluxes given tau and everything is 1 dimensional. This is the exact same function as `get_flux_geom_3d` but is kept separately so we don't have to do unecessary indexing for fast retrievals. :param nlevel: Number of levels in the model :type nlevel: int :param wno: Wave number grid in cm -1 :type wno: array of float :param nwno: Number of wave points :type nwno: int :param numg: Number of Gauss angles :type numg: int :param numt: Number of Chebyshev angles :type numt: int :param DTAU: This is the opacity contained within each individual layer (defined at midpoints of "levels") WITHOUT D-Eddington Correction Dimensions=# layer by # wave :type DTAU: ndarray of float :param TAU: This is the cumulative summed opacity WITHOUT D-Eddington Correction Dimensions=# level by # wave :type TAU: ndarray of float :param W0: This is the single scattering albedo, from scattering, clouds, raman, etc WITHOUT D-Eddington Correction Dimensions=# layer by # wave :type W0: ndarray of float :param COSB: This is the asymmetry factor WITHOUT D-Eddington Correction Dimensions=# layer by # wave :type COSB: ndarray of float :param GCOS2: Parameter that allows us to directly include Rayleigh scattering = 0.5*tau_rayleigh/(tau_rayleigh + tau_cloud) :type GCOS2: ndarray of float :param ftau_cld: Fraction of cloud extinction to total = tau_cloud/(tau_rayleigh + tau_cloud) :type ftau_cld: ndarray of float :param ftau_ray: Fraction of rayleigh extinction to total = tau_rayleigh/(tau_rayleigh + tau_cloud) :type ftau_ray: ndarray of float :param dtau_og: This is the opacity contained within each individual layer (defined at midpoints of "levels") WITHOUT the delta eddington correction, if it was specified by user Dimensions=# layer by # wave :type dtau_og: ndarray of float :param tau_og: This is the cumulative summed opacity WITHOUT the delta eddington correction, if it was specified by user Dimensions=# level by # wave :type tau_og: ndarray of float :param w0_og: Same as w0 but WITHOUT the delta eddington correction, if it was specified by user :type w0_og: ndarray of float :param cosb_og: Same as cosbar buth WITHOUT the delta eddington correction, if it was specified by user :type cosb_og: ndarray of float :param surf_reflect: Surface reflectivity :type surf_reflect: float :param ubar0: matrix of cosine of the incident angle from geometric.json :type ubar0: ndarray of float :param ubar1: matrix of cosine of the observer angles :type ubar1: ndarray of float :param cos_theta: Cosine of the phase angle of the planet :type cos_theta: float :param F0PI: Downward incident solar radiation :type F0PI: array :param single_phase: Single scattering phase function, default is the two-term henyey-greenstein phase function :type single_phase: str :param multi_phase: Multiple scattering phase function, defulat is N=2 Legendre polynomial :type multi_phase: str :param frac_a: (Optional), If using the TTHG phase function. Must specify the functional form for fraction of forward to back scattering (A + B * gcosb^C) :type frac_a: float :param frac_b: (Optional), If using the TTHG phase function. Must specify the functional form for fraction of forward to back scattering (A + B * gcosb^C) :type frac_b: float :param frac_c: (Optional), If using the TTHG phase function. Must specify the functional form for fraction of forward to back scattering (A + B * gcosb^C), Default is : 1 - gcosb^2 :type frac_c: float :param constant_back: (Optional), If using the TTHG phase function. Must specify the assymetry of back scatterer. Remember, the output of A & M code does not separate back and forward scattering. :type constant_back: float :param constant_forward: (Optional), If using the TTHG phase function. Must specify the assymetry of forward scatterer. Remember, the output of A & M code does not separate back and forward scattering. :type constant_forward: float :param tridiagonal: 0 for tridiagonal, 1 for pentadiagonal :type tridiagonal: int :param toon_coefficients: 0 for quadrature (default) 1 for eddington :type toon_coefficients: int :returns: * *intensity at the top of the atmosphere for all the different ubar1 and ubar2* * *To Do* * *-----* * *- F0PI Solar flux shouldn't always be 1.. Follow up to make sure that this isn't a bad* -- hardwiring to solar, despite "relative albedo" .. py:function:: emission_bare_surface(T_surf, wl, surf_reflect) Compute the emergent top-of-atmosphere flux from a bare rock This function considers only pure thermal emission (i.e. no scattering). Uses gauss weights of 5 to stay consistent with emission_Toon :param T: Temperatures in each atmospheric layer (K). :type T: np.array of float :param wl: Wavelength grid (μm). :type wl: np.array of float :param surf_reflect: numpy.ndarray Surface reflectivity as a function of wavelength. :returns: Spectral surface flux in SI units (W/m^2/sr/m). :rtype: F (np.array of float) .. py:function:: reflection_bare_surface(wl, surf_reflect, Gauss_quad=5) Compute the emergent top-of-atmosphere flux from a bare rock This function considers only pure thermal emission (i.e. no scattering). Uses gauss weights of 5 to stay consistent with emission_Toon :param wl: Wavelength grid (μm). :type wl: np.array of float :param surf_reflect: numpy.ndarray Surface reflectivity as a function of wavelength. :param Gauss_quad: Gaussian quadrature order for integration over emitting surface (Options: 2 / 3). :type Gauss_quad: int :returns: Spectral surface flux in SI units (W/m^2/sr/m). :rtype: F (np.array of float) .. py:function:: assign_assumptions_and_compute_single_stream_emission(P, T, dz, wl, kappa_tot, dtau_tot, kappa_gas, kappa_Ray, kappa_cloud, kappa_cloud_seperate, zone_idx, Gauss_quad, device, P_cloud, cloud_dim, aerosol_species, f_cloud, albedo_deck, disable_atmosphere, surface, surface_model, P_surf, T_surf, surf_reflect, surf_reflect_array, surface_component_percentages, surface_percentage_apply_to) Assigns assumptions (surfaces, clouds) and computes single-stream emission Maximum complexity allowed: patchy 2D clouds (1 aerosol max), surfaces (and decks) with albedos :param TO DO: :returns: Spectral surface flux in SI units (W/m^2/sr/m). :rtype: F (np.array of float) .. py:function:: assign_assumptions_and_compute_thermal_scattering_emission(P, T, dz, wl, kappa_tot, dtau_tot, kappa_gas, kappa_Ray, kappa_cloud, kappa_cloud_seperate, zone_idx, Gauss_quad, device, P_cloud, cloud_dim, aerosol_species, w_cloud, g_cloud, f_cloud, f_both, f_aerosol_1, f_aerosol_2, f_clear, albedo_deck, disable_atmosphere, surface, surface_model, P_surf, T_surf, surf_reflect, surf_reflect_array, surface_component_percentages, surface_percentage_apply_to) Assigns assumptions (surfaces, clouds) and computes thermal scattering emission Maximum complexity allowed: patchy 2D clouds (2 aerosols max), surfaces (and decks) with albedos :param TO DO: :returns: Spectral surface flux in SI units (W/m^2/sr/m). :rtype: F (np.array of float) .. py:function:: assign_assumptions_and_compute_reflection(P, T, dz, wl, kappa_tot, dtau_tot, kappa_gas, kappa_Ray, kappa_cloud, kappa_cloud_seperate, zone_idx, reflection_up_to_5um, P_cloud, cloud_dim, aerosol_species, w_cloud, g_cloud, f_cloud, f_both, f_aerosol_1, f_aerosol_2, f_clear, albedo_deck, disable_atmosphere, surface, surface_model, P_surf, T_surf, surf_reflect, surf_reflect_array, surface_component_percentages, surface_percentage_apply_to) Assigns assumptions (surfaces, clouds) and computes reflection Maximum complexity allowed: patchy 2D clouds (2 aerosols max), surfaces (and decks) with albedos 2 aerosols max for patchy 2d clouds available for gas giants (no surface or shiny deck) 1 aerosol max for patchy 2d clouds available for all models :param TO DO: :returns: albedo (np.array of float)