petitRADTRANS.radtrans
======================

.. py:module:: petitRADTRANS.radtrans

.. autoapi-nested-parse::

   Stores the Radtrans object.



Classes
-------

.. autoapisummary::

   petitRADTRANS.radtrans.Radtrans


Module Contents
---------------

.. py:class:: Radtrans(pressures: numpy.typing.NDArray[numpy.floating] = None, wavelength_boundaries: numpy.typing.NDArray[numpy.floating] = None, line_species: list[str] = None, gas_continuum_contributors: list[str] = None, rayleigh_species: list[str] = None, cloud_species: list[str] = None, line_opacity_mode: str = 'c-k', line_by_line_opacity_sampling: int = 1, scattering_in_emission: bool = False, emission_angle_grid: numpy.typing.NDArray[numpy.floating] = None, anisotropic_cloud_scattering: bool = 'auto', path_input_data: str = None)

   Calculate spectra using a given set of opacities.

   Args:
       pressures (Optional):
           (bar) Array defining the pressure grid to be used for the spectral calculations.
       wavelength_boundaries (Optional):
           List containing left and right border of wavelength region to be considered, in micron. If nothing else
           is specified, it will be equal to ``[0.05, 300]``, hence using the full petitRADTRANS wavelength range
           (0.11 to 250 microns for ``'c-k'`` mode, 0.3 to 30 microns for the ``'lbl'`` mode). The larger the
           range the longer the computation time.
       line_species (Optional):
           List of strings, denoting which line absorber species to include.
       gas_continuum_contributors (Optional):
           List of strings, denoting which continuum absorber species to include.
       rayleigh_species (Optional):
           List of strings, denoting which Rayleigh scattering species to include.
       cloud_species (Optional):
           List of strings, denoting which cloud opacity species to include.
       line_opacity_mode (Optional[string]):
           If equal to ``'c-k'``: use low-resolution mode, at :math:`\\lambda/\\Delta \\lambda = 1000`, with the
           correlated-k assumption. if equal to ``'lbl'``: use high-resolution mode, at
           :math:`\\lambda/\\Delta \\lambda = 10^6`, with a line-by-line treatment.
       line_by_line_opacity_sampling (Optional[int]):
           If ``mode = 'lbl'``, then this will only load every line_by_line_opacity_sampling-nth point of the
           high-resolution opacities. This can be used to save time and RAM.
           The user should verify whether this leads to solutions which are identical to the results of the
           non down-sampled :math:`10^6` resolution.
       scattering_in_emission (Optional[bool]):
           Will be ``False`` by default.
           If ``True`` scattering will be included in the emission spectral calculations. Note that this increases
           the runtime of pRT!
       emission_angle_grid (Optional):
           Array defining the cosines of the angle grid to be used for the emission spectrum calculations,
           and their weights. The array is of shape (2, n_angles), with emission_angle_grid[0] being the cosines of
           the angles, and emission_angle_grid[1] being the weights.
           A dictionary of array with keys 'cos_angles' and 'weights' can also be used.
           If None, a default set of values and weights are used.
       anisotropic_cloud_scattering (Optional[bool, str]):
           If True, anisotropic cloud scattering opacities are used for the spectral calculations.
           If False, isotropic cloud scattering opacities are used for the spectral calculations.
           If 'auto' (recommended), anisotropic_cloud_scattering is set to True for emission spectrum
           calculations, and to False for transmission spectrum calculations.
       path_input_data (Optional[str]):
           Path to the input_data folder, containing the files to be loaded by petitRADTRANS.


   .. py:attribute:: __dat_opacity_files_warning_message
      :value: Multiline-String

      .. raw:: html

         <details><summary>Show Value</summary>

      .. code-block:: python

         """loading opacities from .dat files is discouraged, the HDF5 format offer better performances at for a lower memory usage
         
         Converting petitRADTRANS .dat opacity files into HDF5 can be done by executing:
         >>> from petitRADTRANS.__file_conversion import convert_all
         >>> convert_all()
         
         Alternatively, the petitRADTRANS HDF5 files can be downloaded (see https://petitradtrans.readthedocs.io/en/latest/content/available_opacities.html)"""

      .. raw:: html

         </details>




   .. py:attribute:: __line_opacity_property_setting_warning_message
      :value: Multiline-String

      .. raw:: html

         <details><summary>Show Value</summary>

      .. code-block:: python

         """setting a Radtrans line opacity property should be avoided
         These properties are loaded from the opacity data in the input_data directory and are inter-dependent (they need to be updated for consistency)
         It is recommended to create a new Radtrans instance instead"""

      .. raw:: html

         </details>




   .. py:attribute:: __property_setting_warning_message
      :value: Multiline-String

      .. raw:: html

         <details><summary>Show Value</summary>

      .. code-block:: python

         """setting a Radtrans property directly is not recommended
         Create a new Radtrans instance (recommended) or re-do all the setup steps necessary for the modification to be taken into account"""

      .. raw:: html

         </details>




   .. py:attribute:: _frequency_grid_misalignment_tolerance
      :type:  float
      :value: 0.002



   .. py:attribute:: _pressures


   .. py:attribute:: _line_opacity_mode
      :value: 'c-k'



   .. py:attribute:: _line_by_line_opacity_sampling
      :value: 1



   .. py:attribute:: _scattering_in_emission
      :value: False



   .. py:attribute:: _anisotropic_cloud_scattering
      :value: 'auto'



   .. py:attribute:: _path_input_data
      :value: None



   .. py:attribute:: __scattering_in_transmission
      :value: False



   .. py:attribute:: __sum_opacities
      :value: False



   .. py:attribute:: __absorber_present


   .. py:attribute:: __use_smart_clear_opacities
      :value: False



   .. py:attribute:: _lines_loaded_opacities


   .. py:attribute:: _cias_loaded_opacities


   .. py:attribute:: _clouds_loaded_opacities


   .. py:attribute:: _emission_angle_grid


   .. py:method:: __getattr__(name)

      Override of the object base __getattr__ method, in order to hint towards pRT3 names when pRT2 names are used.
              



   .. py:property:: anisotropic_cloud_scattering


   .. py:property:: cias_loaded_opacities


   .. py:property:: clouds_loaded_opacities


   .. py:property:: cloud_species


   .. py:property:: emission_angle_grid


   .. py:property:: frequencies


   .. py:property:: frequency_bins_edges


   .. py:property:: gas_continuum_contributors


   .. py:property:: line_by_line_opacity_sampling


   .. py:property:: lines_loaded_opacities


   .. py:property:: line_opacity_mode


   .. py:property:: line_species


   .. py:property:: path_input_data


   .. py:property:: pressures
      :type: numpy.typing.NDArray[numpy.floating]



   .. py:property:: rayleigh_species


   .. py:property:: scattering_in_emission


   .. py:property:: wavelength_boundaries


   .. py:method:: __check_anisotropic_cloud_scattering(mode)
      :staticmethod:



   .. py:method:: __check_line_opacity_mode(mode)
      :staticmethod:



   .. py:method:: __check_input_data_file_existence(path)
      :staticmethod:



   .. py:method:: __check_path_input_data(path)
      :staticmethod:



   .. py:method:: __check_pressures(pressures)
      :staticmethod:



   .. py:method:: __check_return_clear_spectrum_relevance(return_clear_spectrum, clouds_have_effect, opaque_cloud_top_pressure, cloud_fraction)
      :staticmethod:



   .. py:method:: __check_wavelength_boundaries(boundaries)
      :staticmethod:



   .. py:method:: __clouds_have_effect(mass_fractions, cloud_fraction)

      Check if the clouds have any effect, i.e. if the cloud species MMR is greater than 0.

      Args:
          mass_fractions: atmospheric mass mixing ratios
          cloud_fraction: cloud coverage fraction



   .. py:method:: __compute_additional_continuum_opacities(frequencies, pressures, additional_absorption_opacities_function, additional_scattering_opacities_function, cloud_photosphere_median_optical_depth)
      :staticmethod:



   .. py:method:: __compute_cloud_continuum_opacities(frequencies, frequency_bins_edges, pressures, temperatures, mass_fractions, mean_molar_masses, reference_gravity, scattering_in_transmission, sum_opacities, opaque_cloud_top_pressure, power_law_opacity_350nm, power_law_opacity_coefficient, clouds_have_effect, cloud_species, clouds_particles_porosity_factor, clouds_loaded_opacities, anisotropic_cloud_scattering, cloud_particles_mean_radii, cloud_particle_radius_distribution_std, cloud_f_sed, eddy_diffusion_coefficients, cloud_particles_radius_distribution, cloud_hansen_a, cloud_hansen_b, cloud_particle_number_density_grid, cloud_photosphere_median_optical_depth, complete_coverage_clouds, return_cloud_contribution, target_particle_radii, target_pressure)
      :staticmethod:



   .. py:method:: __compute_continuum_opacities(frequencies, frequency_bins_edges, pressures, temperatures, mass_fractions, mean_molar_masses, cias_loaded_opacities, gas_continuum_contributors, gray_opacity)
      :staticmethod:



   .. py:method:: __compute_continuum_scattering_opacities(frequencies, pressures, temperatures, mass_fractions, mean_molar_masses, rayleigh_species, haze_factor)
      :staticmethod:



   .. py:method:: __compute_combined_opacities(opacities, continuum_opacities, line_species, line_species_mass_fractions, line_opacity_mode, sum_opacities, g_gauss, weights_gauss)
      :staticmethod:



   .. py:method:: __get_base_species_mass_fractions(mass_fractions)
      :staticmethod:



   .. py:method:: __get_cia_contributors(gas_continuum_contributors)
      :staticmethod:



   .. py:method:: __get_base_opacities(use_smart_clear_opacities, complete_coverage_clouds, opacities, continuum_opacities_scattering, cloud_anisotropic_scattering_opacities=None, cloud_absorption_opacities=None)
      :staticmethod:


      Return the opacities that should be returned when not using cloud coverage fraction.



   .. py:method:: __get_cia_aliases(name: str) -> str
      :staticmethod:



   .. py:method:: __get_frequency_grids_intersection_indices(frequencies: numpy.typing.NDArray[numpy.floating], file_frequencies: numpy.typing.NDArray[numpy.floating], tolerate_grid_misalignment: bool = False) -> tuple[numpy.typing.NDArray[numpy.bool_], numpy.typing.NDArray[numpy.bool_]]
      :staticmethod:


      Get the indices to fill frequencies from a file within the Radtrans frequency grid.

      For example, if the Radtrans frequency grid is [0.1, ..., 0.3, ..., 3] and some loaded opacity frequency grid
      is in the interval [0.3, ..., 3, ..., 28], then the output indices will be the indices corresponding to
      [0.3, ..., 3] on the Radtrans frequency grid, and on the file frequency grid.

      Args:
          frequencies: the Radtrans frequency grid
          file_frequencies: the file frequency grid
          tolerate_grid_misalignment: if True, do not raise an error if the intersection indices sizes differ by at
          most 1.

      Returns:
          The indices of the intersection between both grids, for the Radtrans grid and the file grid0



   .. py:method:: __get_line_opacity_file(path_input_data, species, category)
      :staticmethod:



   .. py:method:: __get_non_cia_gas_continuum_contributions()
      :staticmethod:



   .. py:method:: __get_opacity_sources_dict()
      :staticmethod:



   .. py:method:: __init_clouds_particles_porosity_factor(clouds_particles_porosity_factor, cloud_species)
      :staticmethod:



   .. py:method:: __init_spectral_function(is_emission, reference_gravity, complete_coverage_clouds, return_clear_spectrum, mass_fractions, cloud_fraction, opaque_cloud_top_pressure)

      Initializations and checks common to calculate_flux and calculate_transit_radii.



   .. py:method:: __set_use_smart_clear_opacities(cloud_fraction, complete_coverage_clouds, return_clear_spectrum)


   .. py:method:: __set_sum_opacities(emission)


   .. py:method:: _calculate_flux(temperatures, reference_gravity, opacities, continuum_opacities_scattering, emission_geometry, star_irradiation_cos_angle, stellar_intensity, reflectances, emissivities, cloud_f_sed, photospheric_cloud_optical_depths, cloud_photosphere_wavelength_boundaries, cloud_anisotropic_scattering_opacities, cloud_absorption_opacities, adaptive_feautrier_iterations, return_contribution=False, return_rosseland_opacities=False)

      Calculate the flux.
      TODO complete docstring

      Args:
          temperatures:
          reference_gravity:
          opacities:
          continuum_opacities_scattering:
          emission_geometry:
          star_irradiation_cos_angle:
          stellar_intensity:
          reflectances:
          emissivities:
          cloud_f_sed:
          photospheric_cloud_optical_depths:
          cloud_photosphere_wavelength_boundaries:
          cloud_anisotropic_scattering_opacities:
          cloud_absorption_opacities:
          adaptive_feautrier_iterations:
          return_contribution:
          return_rosseland_opacities:

      Returns:




   .. py:method:: _calculate_opacities(temperatures, mass_fractions, mean_molar_masses, reference_gravity, opaque_cloud_top_pressure=None, cloud_particles_mean_radii=None, cloud_particle_radius_distribution_std=None, cloud_particles_radius_distribution='lognormal', cloud_particle_number_density_grid=None, cloud_hansen_a=None, cloud_hansen_b=None, clouds_particles_porosity_factor=None, cloud_f_sed=None, eddy_diffusion_coefficients=None, haze_factor=1.0, power_law_opacity_350nm=None, power_law_opacity_coefficient=None, gray_opacity=None, cloud_photosphere_median_optical_depth=None, cloud_fraction=1.0, complete_coverage_clouds=None, return_cloud_contribution=False, additional_absorption_opacities_function=None, additional_scattering_opacities_function=None, target_particle_radii=None, target_pressure=None)

      Combine total line opacities, according to mass fractions (abundances), also add continuum opacities,
      i.e. clouds, CIA...
      TODO complete docstring

      Args:
          temperatures:
          mass_fractions:
          mean_molar_masses:
          reference_gravity:
          cloud_particle_radius_distribution_std:
          cloud_f_sed:
          eddy_diffusion_coefficients:
          cloud_particles_mean_radii:
          cloud_particles_radius_distribution:
          cloud_hansen_a:
          cloud_hansen_b:
          cloud_particle_number_density_grid:
          return_cloud_contribution:
          additional_absorption_opacities_function:
          additional_scattering_opacities_function:

      Returns:




   .. py:method:: _calculate_transit_radii(temperatures, mean_molar_masses, reference_gravity, reference_pressure, planet_radius, variable_gravity, opacities, continuum_opacities_scattering, copy_opacities, return_contributions)

      Calculate the transit radii.
      TODO complete docstring

      Args:
          temperatures:
          mean_molar_masses:
          reference_gravity:
          reference_pressure:
          planet_radius:
          variable_gravity:
          opacities:
          continuum_opacities_scattering:
          copy_opacities:

      Returns:




   .. py:method:: _combine_ck_opacities(opacities, g_gauss, weights_gauss)
      :staticmethod:



   .. py:method:: _combine_opacities(line_species_mass_fractions, opacities, continuum_opacities)
      :staticmethod:



   .. py:method:: _compute_cia_opacities(cia_dicts, mass_fractions, pressures, temperatures, frequencies, mean_molar_masses)
      :staticmethod:


      Wrapper to _interpolate_cia, calculating each collision's combined mass fraction.



   .. py:method:: _compute_ck_flux(frequencies, temperatures, weights_gauss, emission_cos_angle_grid, emission_cos_angle_grid_weights, optical_depths, return_contribution)
      :staticmethod:



   .. py:method:: _integrate_cloud_opacities_from_dn_dr(dn_dr, atmosphere_densities, clouds_particles_densities, cloud_particles_radii_bins, cloud_particles_radii, clouds_absorption_opacities, clouds_scattering_opacities, clouds_particles_asymmetry_parameters, pressures, target_particle_radii, target_pressure)
      :staticmethod:


      Shared integration backend for both cloud opacity methods.
      Operates purely on dn/dr and handles all integration logic.



   .. py:method:: _compute_cloud_log_normal_particles_distribution_opacities(atmosphere_densities, clouds_particles_densities, clouds_mass_fractions, cloud_particles_mean_radii, cloud_particles_distribution_std, cloud_particles_radii_bins, cloud_particles_radii, clouds_absorption_opacities, clouds_scattering_opacities, clouds_particles_asymmetry_parameters, pressures, target_particle_radii, target_pressure)
      :staticmethod:


      This function reimplements calc_cloud_opas from fortran_radtrans_core.f90.
      For some reason it runs faster in python than in fortran, so we'll use this from now on.
      This function integrates the cloud opacity through the different layers of the atmosphere to get the total
      optical depth, scattering and anisotropic fraction.
      # TODO optical depth or opacity?
      # TODO complete docstring and vectorize on/off code

      Author: Francois Rozet

      Args:
          atmosphere_densities:
              Density of the atmosphere at each of its layer
          clouds_particles_densities:
              Density of each cloud particles
          clouds_mass_fractions:
              Mass fractions of each cloud at each atmospheric layer
          cloud_particles_mean_radii:
              Mean radius of each cloud particles at each atmospheric layer
          cloud_particles_distribution_std:
              Standard deviation of the log-normal cloud particles distribution
          cloud_particles_radii_bins:
              Bins of the particles cloud radii grid
          cloud_particles_radii:
              Particles cloud radii grid
          clouds_absorption_opacities:
              Cloud absorption opacities (radius grid, wavelength grid, clouds)
          clouds_scattering_opacities:
              Cloud scattering opacities (radius grid, wavelength grid, clouds)
          clouds_particles_asymmetry_parameters:
              Cloud particles asymmetry parameters (radius grid, wavelength grid, clouds)

      Returns:




   .. py:method:: _compute_cloud_opacities_using_particle_size_grid(atmosphere_densities, clouds_mass_fractions, cloud_particle_number_density_grid, clouds_particles_densities, cloud_particles_radii_bins, cloud_particles_radii, clouds_absorption_opacities, clouds_scattering_opacities, clouds_particles_asymmetry_parameters, pressures, target_particle_radii, target_pressure)
      :staticmethod:


      # TODO vectorize on/off code
      Args:
          atmosphere_densities:
              Density of the atmosphere at each of its layer
          clouds_mass_fractions:
              Mass fractions of each cloud at each atmospheric layer
          cloud_particle_number_density_grid:
              Output from cloud lognormal EXPERIMENTAL
          clouds_particles_densities:
              Density of each cloud particles
          cloud_particles_radii_bins:
              Bins of the particles cloud radii grid
          cloud_particles_radii:
              Particles cloud radii grid
          clouds_absorption_opacities:
              Cloud absorption opacities (radius grid, wavelength grid, clouds)
          clouds_scattering_opacities:
              Cloud scattering opacities (radius grid, wavelength grid, clouds)
          clouds_particles_asymmetry_parameters:
              Cloud particles asymmetry parameters (radius grid, wavelength grid, clouds)
          pressures:
              Pressure array
          target_particle_radii:
              Target particle radii array
          target_pressure:
              Target pressure array

      Returns:




   .. py:method:: _compute_cloud_opacities(pressures, temperatures, frequency_bins_edges, cloud_species_mass_fractions, mean_molar_masses, reference_gravity, cloud_particle_radius_distribution_std, clouds_loaded_opacities, sum_opacities, anisotropic_cloud_scattering, cloud_f_sed=None, eddy_diffusion_coefficients=None, cloud_particles_mean_radii=None, cloud_particles_radius_distribution='lognormal', cloud_hansen_a=None, cloud_hansen_b=None, cloud_particle_number_density_grid=None, clouds_particles_porosity_factor=None, photospheric_cloud_optical_depths=None, continuum_opacities=None, return_cloud_contribution=False, target_particle_radii=None, target_pressure=None)
      :staticmethod:


      Calculate cloud opacities for a defined atmospheric structure.
      # TODO complete docstring

      Args:
          temperatures:
          cloud_species_mass_fractions:
          mean_molar_masses:
          reference_gravity:
          cloud_particle_radius_distribution_std:
          cloud_f_sed:
          eddy_diffusion_coefficients:
          cloud_particles_mean_radii:
          cloud_particles_radius_distribution:
          cloud_hansen_a:
          cloud_hansen_b:
          cloud_particle_number_density_grid:
          return_cloud_contribution:

      Returns:




   .. py:method:: _compute_species_optical_depths(reference_gravity, pressures, cloud_opacities)
      :staticmethod:


      Calculate the optical depth of the clouds as function of
      frequency and pressure. The array with the optical depths is set to the
      ``tau_cloud`` attribute. The optical depth is calculated from the top of
      the atmosphere (i.e. the smallest pressure). Therefore, below the cloud
      base, the optical depth is constant and equal to the value at the cloud
      base.
      # TODO complete docstring

          Args:
              reference_gravity (float):
                  Surface gravity in cgs. Vertically constant for emission
                  spectra.



   .. py:method:: _compute_feautrier_radiative_transfer(frequency_bins_edges, temperatures, weights_gauss, emission_cos_angle_grid, emission_cos_angle_grid_weights, optical_depths, photon_destruction_probabilities, emission_geometry, stellar_intensity, star_irradiation_cos_angle, reflectances, emissivities, adaptive_feautrier_iterations, return_contribution)
      :staticmethod:



   .. py:method:: _compute_h_minus_free_free_xsec(wavelengths, temperatures, electron_partial_pressure)
      :staticmethod:


      Calculate the H- free-free cross-section in units of cm^2 per H per e- pressure (in cgs).
      Source: "The Observation and Analysis of Stellar Photospheres" by David F. Gray, p. 156
      # TODO complete docstring
      Args:
          wavelengths: (angstrom)
          temperatures:
          electron_partial_pressure:

      Returns:




   .. py:method:: _compute_h_minus_bound_free_xsec(wavelengths_bin_edges)
      :staticmethod:


      Calculate the H- bound-free cross-section in units of cm^2 per H-, as defined on page 155 of
      "The Observation and Analysis of Stellar Photospheres" by David F. Gray
      # TODO complete docstring
      Args:
          wavelengths_bin_edges: (angstrom)

      Returns:




   .. py:method:: _compute_h_minus_opacities(mass_fractions, pressures, temperatures, frequencies, frequency_bins_edges, mean_molar_masses, **kwargs)
      :staticmethod:


      Calculate the H- opacity.



   .. py:method:: _compute_non_cia_gas_induced_continuum_opacities(gas_continuum_contributors, mass_fractions, pressures, temperatures, frequencies, frequency_bins_edges, mean_molar_masses, **kwargs)
      :staticmethod:



   .. py:method:: _compute_optical_depths(pressures, reference_gravity, opacities, continuum_opacities_scattering, scattering_in_emission)
      :staticmethod:



   .. py:method:: _compute_optical_depths_wrapper(pressures, reference_gravity, opacities, continuum_opacities_scattering, scattering_in_emission, sum_opacities, photospheric_cloud_optical_depths=None, absorber_present=True, **custom_cloud_parameters)
      :staticmethod:



   .. py:method:: _compute_optical_depths_with_photospheric_cloud(pressures, reference_gravity, opacities, continuum_opacities_scattering, frequencies, weights_gauss, cloud_wavelengths, cloud_f_sed, cloud_anisotropic_scattering_opacities, cloud_absorption_opacities, photospheric_cloud_optical_depths, scattering_in_emission)
      :staticmethod:



   .. py:method:: _compute_power_law_opacities(power_law_opacity_350nm, power_law_opacity_coefficient, frequencies, n_layers)
      :staticmethod:



   .. py:method:: _compute_radius_hydrostatic_equilibrium(pressures, temperatures, mean_molar_masses, reference_gravity, reference_pressure, planet_radius, variable_gravity=True)
      :staticmethod:



   .. py:method:: _compute_rayleigh_scattering_opacities(rayleigh_species, pressures, temperatures, mass_fractions, mean_molar_masses, frequencies, haze_factor=1.0)
      :staticmethod:


      Add Rayleigh scattering opacities to scattering continuum opacities.

      Args:
          temperatures: temperatures in each atmospheric layer
          mass_fractions: dictionary of the Rayleigh scattering species mass fractions



   .. py:method:: _compute_rosseland_opacities(frequency_bins_edges, temperatures, weights_gauss, opacities, continuum_opacities_scattering, scattering_in_emission)
      :staticmethod:



   .. py:method:: _compute_transit_radii(opacities, continuum_opacities_scattering, pressures, temperatures, weights_gauss, mean_molar_masses, reference_gravity, reference_pressure, planet_radius, variable_gravity, summed_species_opacities, copy_opacities)
      :staticmethod:


      Calculate the planetary transmission spectrum.
      # TODO complete docstring
      Args:
          opacities:
          continuum_opacities_scattering:
          pressures:
          temperatures:
          weights_gauss:
          mean_molar_masses:
              Mean molecular weight in units of amu.
              (1-d numpy array, same length as pressure array).
          reference_gravity (float):
              Atmospheric gravitational acceleration at reference pressure and radius in units of
              dyne/cm^2
          reference_pressure (float):
              Reference pressure in bar.
          planet_radius (float):
              Planet radius in cm.
          variable_gravity (bool):
              If True, gravity in the atmosphere will vary proportional to 1/r^2, where r is the planet
              radius.
          summed_species_opacities (bool):
              If True, it is assumed that the given opacities contains the summed opacities of all species.
              This is the case in line-by-line.
          copy_opacities (bool):
              If True, a new array will be created for the transparencies, distinct from the opacities. This
              increases memory usage.

      Returns:
          * transmission radius in cm (1-d numpy array, as many elements as wavelengths)
          * planet radius as function of atmospheric pressure (1-d numpy array, as many elements as atmospheric
          layers)



   .. py:method:: _compute_transmission_spectrum_contribution(transit_radii, pressures, temperatures, mean_molar_masses, reference_gravity, reference_pressure, planet_radius, weights_gauss, opacities, continuum_opacities_scattering, scattering_in_transmission, variable_gravity)
      :staticmethod:



   .. py:method:: _handle_grid_misalignment(indices: numpy.typing.NDArray[numpy.bool_], frequencies_reference: numpy.typing.NDArray[numpy.floating], frequencies_test: numpy.typing.NDArray[numpy.floating]) -> numpy.typing.NDArray[numpy.bool_]
      :staticmethod:



   .. py:method:: _init_cia_loaded_opacities(cia_contributors)
      :staticmethod:



   .. py:method:: _init_frequency_grid()

      Initialize the Radtrans frequency grid, used to calculate the spectra.
      The frequency grid comes from the requested opacity files, in the following priority order:
          1. lines,
          2. CIA,
          3. clouds.
      If not opacities are provided, the mean of the wavelength boundaries is used. The frequency grid in that case
      has only 1 element.

      Returns:
          frequencies:
              (Hz) frequencies (center of bin) of the line opacities, also use for spectral calculations, of size N
          frequency_bins_edges:
              (Hz) edges of the frequencies bins, of size N+1
              for correlated-k only, number of points used to sample the g-space (1 in the case lbl is used)



   .. py:method:: _init_frequency_grid_from_frequency_grid(frequency_grid, wavelength_boundaries, sampling=1)
      :staticmethod:



   .. py:method:: _init_frequency_grid_from_lines()


   .. py:method:: _interpolate_cia(collision_dict, combined_mass_fractions, pressures, temperatures, frequencies, mean_molar_masses)
      :staticmethod:


      Interpolate CIA cross-sections onto the Radtrans (wavelength, temperature) grid and convert it into
      opacities.

      Args:
          combined_mass_fractions: combined mass fractions of the colliding species
              e.g., for H2-He and an atmosphere with H2 and He MMR of respectively 0.74 and 0.24,
              combined_mas_fractions = 0.74 * 0.24
              combined_mas_fractions is divided by the combined weight (e.g. for H2 and He, 2 * 4 AMU^2), so there is
              no units issue.

      Returns:
          A (wavelength, temperature) array containing the CIA opacities.



   .. py:method:: _interpolate_species_opacities(pressures, temperatures, n_g, n_frequencies, line_opacities_grid, line_opacities_temperature_pressure_grid, line_opacities_temperature_grid_size, line_opacities_pressure_grid_size)
      :staticmethod:



   .. py:method:: compute_bins_edges(middle_bin_points: numpy.typing.NDArray[numpy.floating]) -> numpy.typing.NDArray[numpy.floating]
      :staticmethod:


      Calculate bin edges for middle bin points.

      Args:
          middle_bin_points: array of size N containing the middle bin points to calculate the bin edges of

      Returns:
          array of size N+1 containing the bins edges



   .. py:method:: calculate_contribution_spectra(mode, opaque_cloud_top_pressure=None, power_law_opacity_350nm=None, power_law_opacity_coefficient=None, gray_opacity=None, cloud_photosphere_median_optical_depth=None, additional_absorption_opacities_function=None, additional_scattering_opacities_function=None, stellar_intensities=None, **kwargs)


   .. py:method:: calculate_flux(temperatures: numpy.typing.NDArray[numpy.floating], mass_fractions: dict[str, numpy.typing.NDArray[numpy.floating]], mean_molar_masses: numpy.typing.NDArray[numpy.floating], reference_gravity: float, planet_radius: float = None, opaque_cloud_top_pressure: float = None, cloud_particles_mean_radii: dict[str, numpy.typing.NDArray[numpy.floating]] = None, cloud_particle_radius_distribution_std: float = None, cloud_particles_radius_distribution: str = 'lognormal', cloud_hansen_a: dict[str, numpy.typing.NDArray[numpy.floating]] = None, cloud_hansen_b: dict[str, numpy.typing.NDArray[numpy.floating]] = None, cloud_particle_number_density_grid: dict[str, numpy.typing.NDArray[numpy.floating]] = None, clouds_particles_porosity_factor: dict[str, float] = None, cloud_f_sed: float = None, eddy_diffusion_coefficients: numpy.typing.NDArray[numpy.floating] = None, haze_factor: float = 1.0, power_law_opacity_350nm: float = None, power_law_opacity_coefficient: float = None, gray_opacity: float = None, cloud_photosphere_wavelength_boundaries: numpy.typing.NDArray[numpy.floating] = None, cloud_photosphere_median_optical_depth: float = None, cloud_fraction: float = 1.0, complete_coverage_clouds: list[str] = None, irradiation_geometry: str = 'dayside_ave', emission_geometry: str = None, stellar_intensities: numpy.typing.NDArray[numpy.floating] = None, star_effective_temperature: float = None, star_radius: float = None, orbit_semi_major_axis: float = None, star_irradiation_angle: float = 0.0, reflectances: numpy.typing.NDArray[numpy.floating] = None, emissivities: numpy.typing.NDArray[numpy.floating] = None, additional_absorption_opacities_function: Callable = None, additional_scattering_opacities_function: Callable = None, adaptive_feautrier_iterations: bool = False, frequencies_to_wavelengths: bool = True, return_contribution: bool = False, return_clear_spectrum: bool = False, return_photosphere_radius: bool = False, return_rosseland_optical_depths: bool = False, return_cloud_contribution: bool = False, return_opacities: bool = False, return_abundances: bool = False, return_particles_radii_bins: bool = False, return_particles_radii: bool = False, target_particle_radii: numpy.typing.NDArray[numpy.floating] = None, target_pressure: numpy.typing.NDArray[numpy.floating] = None) -> tuple[numpy.typing.NDArray[numpy.floating], numpy.typing.NDArray[numpy.floating], dict[str, Any]]

      Method to calculate the atmosphere's emitted flux (emission spectrum).

      Args:
          temperatures:
              The atmospheric temperature in K, at each atmospheric layer
              (1-d numpy array, same length as pressure array).
          mass_fractions:
              Dictionary of mass fractions for all atmospheric absorbers.
              Dictionary keys are the species names. Every mass fraction array has same length as pressure array.
          mean_molar_masses:
              The atmospheric mean molecular weight in amu, at each atmospheric layer
              (1-d numpy array, same length as pressure array).
          reference_gravity (float):
              Surface gravity in cgs. Vertically constant for emission spectra.
          planet_radius: planet radius at maximum pressure in cm. Only used to calculate the planet's changing
              Photospheric radius as function of wavelength, if return_photosphere_radius is True.
          opaque_cloud_top_pressure (Optional[float]):
              Pressure, in bar, where opaque cloud deck is added to the absorption opacity.
          cloud_particles_mean_radii (Optional):
              Dictionary of mean particle radii for all cloud species.
              Dictionary keys are the cloud species names. Every radius array has same length as pressure array.
          cloud_particle_radius_distribution_std (Optional[float]):
              Width of the log-normal cloud particle size distribution
          cloud_particles_radius_distribution (Optional[string]):
              The cloud particle size distribution to use.
              Can be either 'lognormal' (default) or 'hansen'.
              If hansen, the cloud_hansen_b parameters must be used.
          cloud_hansen_a (Optional[dict]):
              A dictionary of the 'a' parameter values for each included cloud species and for each atmospheric
              layer, formatted as the kzz argument. Equivalent to cloud_particles_mean_radii.
              If cloud_hansen_a is not included and dist is "hansen", then it will be computed using Kzz and fsed
              (recommended).
          cloud_hansen_b (Optional[dict]):
              A dictionary of the 'b' parameter values for each included cloud species and for each atmospheric
              layer, formatted as the kzz argument. This is the width of the hansen distribution normalized by
              the particle area (1/cloud_hansen_a^2)
          cloud_particle_number_density_grid (Optional[dict]):
              EXPERIMENTAL
          clouds_particles_porosity_factor (Optional[dict]):
              A dictionary of porosity factors depending on the cloud species. This can be useful when opacities
              are calculated using the Distribution of Hollow Spheres (DHS) method.
          cloud_f_sed (Optional[float]):
              Cloud settling parameter.
          eddy_diffusion_coefficients (Optional[float]):
              The atmospheric eddy diffusion coefficient in cgs (i.e. :math:`\rm cm^2/s`), at each atmospheric
              layer (1-d numpy array, same length as pressure array).
          haze_factor (Optional[float]):
              Scalar factor, increasing the gas Rayleigh scattering cross-section.
          power_law_opacity_350nm (Optional[float]):
              Scattering opacity at 0.35 micron, in cgs units (cm^2/g).
          power_law_opacity_coefficient (Optional[float]):
              Has to be given if kappa_zero is defined, this is the
              wavelength powerlaw index of the parametrized scattering
              opacity.
          gray_opacity (Optional[float]):
              Gray opacity value, to be added to the opacity at all pressures and wavelengths
              (units :math:`\rm cm^2/g`)
          cloud_photosphere_median_optical_depth (Optional[float]):
              Median optical depth (across ``wavelength_boundaries``) of the clouds from the top of the
              atmosphere down to the gas-only photosphere. This parameter can be used for enforcing the presence
              of clouds in the photospheric region.
          cloud_photosphere_wavelength_boundaries (Optional[NDArray]):
              Min and max wavelength boundaries when calculating the median of the photospheric cloud optical
              depth. By default, the whole spectral range is used.
          cloud_fraction:
              Mix a column without cloud with a column with cloud to the requested fraction
              (0 = clear, 1 = full cloud cover). By default, assume a full cloud cover.
          complete_coverage_clouds:
              List of clouds that have full coverage, regardless of the value of cloud_fraction.
          irradiation_geometry (Optional[string]):
              if equal to ``'dayside_ave'``: use the dayside average geometry.
              If equal to ``'planetary_ave'``: use the planetary average geometry.
              If equal to ``'non-isotropic'``: use the non-isotropic geometry.
          emission_geometry (Optional[string]):
              Deprecated, same as irradiation_geometry.
          stellar_intensities (Optional[array]):
              The stellar intensity to use. If None, it will be calculated using a PHOENIX model.
          star_effective_temperature (Optional[float]):
              The temperature of the host star in K, used only if the
              scattering is considered. If not specified, the direct light contribution is not calculated.
          star_radius (Optional[float]):
              The radius of the star in cm. If specified, used to scale the to scale the stellar flux,
              otherwise it uses PHOENIX radius.
          orbit_semi_major_axis (Optional[float]):
              The distance of the planet from the star. Used to scale the stellar flux when the scattering of the
              direct light is considered.
          star_irradiation_angle (Optional[float]):
              Inclination angle of the direct light with respect to the normal to the atmosphere. Used only in
              the non-isotropic geometry scenario.
          reflectances (Optional):
              Reflectances of the surface (layer with the highest pressure).
          emissivities (Optional):
              Emissivities of the surface (layer with the highest pressure).
          additional_absorption_opacities_function (Optional[function]):
              A python function that takes wavelength arrays in microns and pressure arrays in bars
              as input, and returns an absorption opacity matrix in units of cm^2/g, in the shape of
              number of wavelength points x number of pressure points.
              This opacity will then be added to the atmospheric absorption opacity.
              This must not be used to add atomic / molecular line opacities in low-resolution mode (c-k),
              because line opacities require a proper correlated-k treatment.
              It may be used to add simple cloud absorption laws, for example, which
              have opacities that vary only slowly with wavelength, such that the current
              model resolution is sufficient to resolve any variations.
          additional_scattering_opacities_function (Optional[function]):
              A python function that takes wavelength arrays in microns and pressure arrays in bars
              as input, and returns an isotropic scattering opacity matrix in units of cm^2/g, in the shape of
              number of wavelength points x number of pressure points.
              This opacity will then be added to the atmospheric absorption opacity.
              It may be used to add simple cloud absorption laws, for example, which
              have opacities that vary only slowly with wavelength, such that the current
              model resolution is sufficient to resolve any variations.
          adaptive_feautrier_iterations (Optional[bool]):
              If True, adapt the number of iterations of the Feautrier scattering-method
              (scattering_in_emission=True) based on the photon destruction probability.
          frequencies_to_wavelengths (Optional[bool]):
              if True, convert the frequencies (Hz) output to wavelengths (cm),
              and the flux per frequency output (erg.s-1.cm-2/Hz) to flux per wavelength (erg.s-1.cm-2/cm)
          return_contribution (Optional[bool]):
              If True, calculate the emission contribution function.
          return_clear_spectrum (Optional[bool]):
              If True, return the clear spectrum in addition to a cloudy spectrum.
          return_photosphere_radius (Optional[bool]):
              If True, calculate and return the photosphere radius.
          return_rosseland_optical_depths (Optional[bool]):
              If True, calculate and return the Rosseland opacities and optical depths.
          return_cloud_contribution (Optional[bool]):
              If True, calculate the cloud contribution.
          return_opacities (Optional[bool]):
              If True, return the absorption opacities and scattering opacities for species and clouds, as well
              as optical depths
          return_abundances (Optional[bool]):
              If True, return the mass fractions.
          return_particles_radii_bins (Optional[bool]):
              If True, return the particle radii bins.
          return_particles_radii (Optional[bool]):
              If True, return the particle radii.
          target_particle_radii (Optional[npt.NDArray[np.floating]]):
              # TODO complete docstring
          target_pressure (Optional[npt.NDArray[np.floating]]):
              # TODO complete docstring



   .. py:method:: calculate_photosphere_radius(temperatures: numpy.typing.NDArray[numpy.floating], mean_molar_masses: numpy.typing.NDArray[numpy.floating], reference_gravity: float, planet_radius: float, opacities: numpy.typing.NDArray[numpy.floating], continuum_opacities_scattering: numpy.typing.NDArray[numpy.floating], cloud_f_sed: float, cloud_photosphere_median_optical_depth: float, cloud_anisotropic_scattering_opacities: numpy.typing.NDArray[numpy.floating], cloud_absorption_opacities: numpy.typing.NDArray[numpy.floating], optical_depths: numpy.typing.NDArray[numpy.floating] = None, cloud_photosphere_wavelength_boundaries: numpy.typing.NDArray[numpy.floating] = None) -> numpy.typing.NDArray[numpy.floating]

      Calculate the photosphere radius.
      TODO complete docstring
      Args:
          temperatures:
          mean_molar_masses:
          reference_gravity:
          planet_radius:
          opacities:
          continuum_opacities_scattering:
          cloud_f_sed:
          cloud_photosphere_wavelength_boundaries:
          cloud_photosphere_median_optical_depth:
          cloud_anisotropic_scattering_opacities:
          cloud_absorption_opacities:
          optical_depths:
      Returns:




   .. py:method:: calculate_rosseland_planck_opacities(temperatures: numpy.typing.NDArray[numpy.floating], mass_fractions: dict[str, numpy.typing.NDArray[numpy.floating]], mean_molar_masses: numpy.typing.NDArray[numpy.floating], reference_gravity: float, opaque_cloud_top_pressure: float = None, cloud_particles_mean_radii: dict[str, numpy.typing.NDArray[numpy.floating]] = None, cloud_particle_radius_distribution_std: float = None, cloud_particles_radius_distribution: str = 'lognormal', cloud_hansen_a: float = None, cloud_hansen_b: float = None, clouds_particles_porosity_factor: dict[str, float] = None, cloud_f_sed: float = None, eddy_diffusion_coefficients: float = None, haze_factor: float = 1.0, power_law_opacity_350nm: float = None, power_law_opacity_coefficient: float = None, gray_opacity: float = None) -> tuple[numpy.typing.NDArray[numpy.floating], numpy.typing.NDArray[numpy.floating], dict[str, Any]]

      Method to calculate the atmosphere's Rosseland and Planck mean opacities.

      Args:
          temperatures:
              the atmospheric temperature in K, at each atmospheric layer
              (1-d numpy array, same length as pressure array).
          mass_fractions:
              dictionary of mass fractions for all atmospheric absorbers.
              Dictionary keys are the species names.
              Every mass fraction array
              has same length as pressure array.
          mean_molar_masses:
              the atmospheric mean molecular weight in amu,
              at each atmospheric layer
              (1-d numpy array, same length as pressure array).
          reference_gravity (float):
              Surface gravity in cgs. Vertically constant for emission
              spectra.
          opaque_cloud_top_pressure (Optional[float]):
              Pressure, in bar, where opaque cloud deck is added to the
              absorption opacity.
          cloud_particles_mean_radii (Optional):
              dictionary of mean particle radii for all cloud species.
              Dictionary keys are the cloud species names.
              Every radius array has same length as pressure array.
          cloud_particle_radius_distribution_std (Optional[float]):
              width of the log-normal cloud particle size distribution
          cloud_particles_radius_distribution (Optional[string]):
              The cloud particle size distribution to use.
              Can be either 'lognormal' (default) or 'hansen'.
              If hansen, the cloud_hansen_b parameters must be used.
          cloud_hansen_a (Optional[dict]):
              A dictionary of the 'a' parameter values for each
              included cloud species and for each atmospheric layer,
              formatted as the kzz argument. Equivalent to radius arg.
              If cloud_hansen_a is not included and dist is "hansen", then it will
              be computed using Kzz and fsed (recommended).
          cloud_hansen_b (Optional[dict]):
              A dictionary of the 'b' parameter values for each
              included cloud species and for each atmospheric layer,
              formatted as the kzz argument. This is the width of the hansen
              distribution normalized by the particle area (1/cloud_hansen_a^2)
          clouds_particles_porosity_factor (Optional[dict]):
              A dictionary containing the particle porosity factor for each cloud species.
          cloud_f_sed (Optional[float]):
              cloud settling parameter
          eddy_diffusion_coefficients (Optional):
              the atmospheric eddy diffusion coefficient in cgs
              (i.e. :math:`\rm cm^2/s`),
              at each atmospheric layer
              (1-d numpy array, same length as pressure array).
          haze_factor (Optional[float]):
              Scalar factor, increasing the gas Rayleigh scattering
              cross-section.
          power_law_opacity_350nm (Optional[float]):
              Scattering opacity at 0.35 micron, in cgs units (cm^2/g).
          power_law_opacity_coefficient (Optional[float]):
              Has to be given if kappa_zero is defined, this is the
              wavelength powerlaw index of the parametrized scattering
              opacity.
          gray_opacity (Optional[float]):
              Gray opacity value, to be added to the opacity at all
              pressures and wavelengths (units :math:`\rm cm^2/g`)



   .. py:method:: calculate_transit_radii(temperatures: numpy.typing.NDArray[numpy.floating], mass_fractions: dict[str, numpy.typing.NDArray[numpy.floating]], mean_molar_masses: numpy.typing.NDArray[numpy.floating], reference_gravity: float, reference_pressure: float, planet_radius: float, variable_gravity: bool = True, opaque_cloud_top_pressure: float = None, cloud_particles_mean_radii: dict[str, numpy.typing.NDArray[numpy.floating]] = None, cloud_particle_radius_distribution_std: float = None, cloud_particles_radius_distribution: str = 'lognormal', cloud_hansen_a: float = None, cloud_hansen_b: float = None, cloud_particle_number_density_grid: dict[str, numpy.typing.NDArray[numpy.floating]] = None, clouds_particles_porosity_factor: dict[str, float] = None, cloud_f_sed: float = None, eddy_diffusion_coefficients: float = None, haze_factor: float = 1.0, power_law_opacity_350nm: float = None, power_law_opacity_coefficient: float = None, gray_opacity: float = None, cloud_fraction: float = 1.0, complete_coverage_clouds: list[str] = None, additional_absorption_opacities_function: Callable = None, additional_scattering_opacities_function: Callable = None, frequencies_to_wavelengths: bool = True, return_contribution: bool = False, return_clear_spectrum: bool = False, return_cloud_contribution: bool = False, return_radius_hydrostatic_equilibrium: bool = False, return_opacities: bool = False, return_abundances: bool = False, return_particles_radii_bins: bool = False, return_particles_radii: bool = False, return_particles_densities: bool = False, target_particle_radii: numpy.typing.NDArray[numpy.floating] = None, target_pressure: numpy.typing.NDArray[numpy.floating] = None) -> tuple[numpy.typing.NDArray[numpy.floating], numpy.typing.NDArray[numpy.floating], dict[str, Any]]

      Method to calculate the atmosphere's transmission radius
      (for the transmission spectrum).

          Args:
              temperatures:
                  the atmospheric temperature in K, at each atmospheric layer
                  (1-d numpy array, same length as pressure array).
              mass_fractions:
                  dictionary of mass fractions for all atmospheric absorbers.
                  Dictionary keys are the species names.
                  Every mass fraction array
                  has same length as pressure array.
              mean_molar_masses:
                  the atmospheric mean molecular weight in amu,
                  at each atmospheric layer
                  (1-d numpy array, same length as pressure array).
              reference_gravity (float):
                  Surface gravity in cgs at reference radius and pressure.
              reference_pressure (float):
                  Reference pressure P0 in bar where R(P=P0) = R_pl,
                  where R_pl is the reference radius (parameter of this
                  method), and g(P=P0) = gravity, where gravity is the
                  reference gravity (parameter of this method)
              planet_radius (float):
                  Reference radius R_pl, in cm.
              variable_gravity (Optional[bool]):
                  Standard is ``True``. If ``False`` the gravity will be
                  constant as a function of pressure, during the transmission
                  radius calculation.
              opaque_cloud_top_pressure (Optional[float]):
                  Pressure, in bar, where opaque cloud deck is added to the
                  absorption opacity.
              cloud_particles_mean_radii (Optional):
                  dictionary of mean particle radii for all cloud species.
                  Dictionary keys are the cloud species names.
                  Every radius array has same length as pressure array.
              cloud_particle_radius_distribution_std (Optional[float]):
                  width of the log-normal cloud particle size distribution
              cloud_particles_radius_distribution (Optional[string]):
                  The cloud particle size distribution to use.
                  Can be either 'lognormal' (default) or 'hansen'.
                  If hansen, the cloud_hansen_b parameters must be used.
              cloud_hansen_a (Optional[dict]):
                  A dictionary of the 'a' parameter values for each
                  included cloud species and for each atmospheric layer,
                  formatted as the kzz argument. Equivalent to radius arg.
                  If cloud_hansen_a is not included and dist is "hansen", then it will
                  be computed using Kzz and fsed (recommended).
              cloud_hansen_b (Optional[dict]):
                  A dictionary of the 'b' parameter values for each
                  included cloud species and for each atmospheric layer,
                  formatted as the kzz argument. This is the width of the hansen
                  distribution normalized by the particle area (1/cloud_hansen_a^2)
              cloud_particle_number_density_grid (Optional[dict]):
                  EXPERIMENTAL
              clouds_particles_porosity_factor (Optional[dict]):
                  A dictionary of porosity factors depending on the cloud species. This can be useful when opacities
                  are calculated using the Distribution of Hollow Spheres (DHS) method.
              cloud_f_sed (Optional[float]):
                  cloud settling parameter
              eddy_diffusion_coefficients (Optional):
                  the atmospheric eddy diffusion coefficient in cgs
                  (i.e. :math:`\rm cm^2/s`),
                  at each atmospheric layer
                  (1-d numpy array, same length as pressure array).
              haze_factor (Optional[float]):
                  Scalar factor, increasing the gas Rayleigh scattering
                  cross-section.
              power_law_opacity_350nm (Optional[float]):
                  Scattering opacity at 0.35 micron, in cgs units (cm^2/g).
              power_law_opacity_coefficient (Optional[float]):
                  Has to be given if kappa_zero is defined, this is the
                  wavelength powerlaw index of the parametrized scattering
                  opacity.
              gray_opacity (Optional[float]):
                  Gray opacity value, to be added to the opacity at all
                  pressures and wavelengths (units :math:`\rm cm^2/g`)
              cloud_fraction:
                  Mix a column without cloud with a column with cloud to the requested fraction
                  (0 = clear, 1 = full cloud cover). By default, assume a full cloud cover.
              complete_coverage_clouds:
                  List of clouds that have full coverage, regardless of the value of cloud_fraction.
              additional_absorption_opacities_function (Optional[function]):
                  A python function that takes wavelength arrays in microns and pressure arrays in bars
                  as input, and returns an absorption opacity matrix in units of cm^2/g, in the shape of
                  number of wavelength points x number of pressure points.
                  This opacity will then be added to the atmospheric absorption opacity.
                  This must not be used to add atomic / molecular line opacities in low-resolution mode (c-k),
                  because line opacities require a proper correlated-k treatment.
                  It may be used to add simple cloud absorption laws, for example, which
                  have opacities that vary only slowly with wavelength, such that the current
                  model resolution is sufficient to resolve any variations.
              additional_scattering_opacities_function (Optional[function]):
                  A python function that takes wavelength arrays in microns and pressure arrays in bars
                  as input, and returns an isotropic scattering opacity matrix in units of cm^2/g, in the shape of
                  number of wavelength points x number of pressure points.
                  This opacity will then be added to the atmospheric absorption opacity.
                  It may be used to add simple cloud absorption laws, for example, which
                  have opacities that vary only slowly with wavelength, such that the current
                  model resolution is sufficient to resolve any variations.
              frequencies_to_wavelengths (Optional[bool]):
                  If True, convert the frequencies (Hz) output to wavelengths (cm).
              return_contribution (Optional[bool]):
                  If True the transmission and emission contribution function is calculated.
              return_clear_spectrum (Optional[bool]):
                  If True, return the clear spectrum in addition to a cloudy spectrum.
              return_cloud_contribution (Optional[bool]):
                  If True, calculate and return the cloud contribution.
              return_radius_hydrostatic_equilibrium (Optional[bool]):
                  If True, return the radius at hydrostatic equilibrium of the planet.
              return_opacities (Optional[bool]):
                  If True, return the absorption opacities and scattering opacities.
              return_abundances (Optional[bool]):
                  If True, return the mass fractions.
              return_particles_radii_bins (Optional[bool]):
                  If True, return the particles radii bins.
              return_particles_radii (Optional[bool]):
                  If True, return the particles radii.
              return_particles_densities (Optional[bool]):
                  If True, return the particles densities.
              target_particle_radii (Optional[npt.NDArray[np.floating]]):
                  # TODO complete docstring
              target_pressure (Optional[npt.NDArray[np.floating]]):
                  # TODO complete docstring



   .. py:method:: compute_pressure_hydrostatic_equilibrium(mean_molar_masses, reference_gravity, planet_radius, p0, temperature, radii, rk4=True)
      :staticmethod:



   .. py:method:: compute_star_spectrum(star_effective_temperature, orbit_semi_major_axis, frequencies, star_radius=None)
      :staticmethod:


      Method to get the PHOENIX spectrum of the star and rebin it
      to the wavelength points. If t_star is not explicitly written, the
      spectrum will be 0. If the distance is not explicitly written,
      the code will raise an error and break to urge the user to
      specify the value.

          Args:
              star_effective_temperature (float):
                  the stellar temperature in K.
              orbit_semi_major_axis (float):
                  the semi-major axis of the planet in cm.
              frequencies (float):
                  the frequencies on which to interpolate the spectrum in Hz.
              star_radius (float):
                  if specified, uses this radius in cm
                  to scale the flux, otherwise it uses PHOENIX radius.



   .. py:method:: get_wavelengths()

      Convert the frequency grid into wavelengths in centimeter.



   .. py:method:: load_all_opacities()


   .. py:method:: load_cia_opacities(path_input_data)


   .. py:method:: load_cia_opacity(file_path_hdf5, return_radtrans_opacities=True)
      :staticmethod:


      Load a CIA opacity file.



   .. py:method:: load_cloud_opacities(path_input_data)

      Load the cloud opacities.

      Args:
          path_input_data: path to petitRADTRANS input_data directory



   .. py:method:: load_cloud_opacity(file_path_hdf5, return_radtrans_opacities=True)
      :staticmethod:


      Load a cloud opacity file.



   .. py:method:: load_hdf5_ktables(file_path_hdf5, frequencies, g_size, temperature_pressure_grid_size, return_radtrans_opacities=True)
      :classmethod:


      Load k-coefficient tables in HDF5 format, based on the ExoMol setup.



   .. py:method:: load_hdf5_line_opacity_table(file_path_hdf5, frequencies, line_by_line_opacity_sampling=1, return_radtrans_opacities=True)
      :staticmethod:


      Load opacities (cm2.g-1) tables in HDF5 format, based on petitRADTRANS pseudo-ExoMol setup.



   .. py:method:: load_line_opacities(path_input_data)

      Load the line opacities for spectral calculation.
      The default pressure-temperature grid is a log-uniform (10, 13) grid.

      Args:
          path_input_data: path to petitRADTRANS input_data directory



   .. py:method:: load_line_opacities_pressure_temperature_grid(hdf5_file)
      :staticmethod:


      Load line opacities temperature grids.



   .. py:method:: rebin_star_spectrum(star_spectrum, star_wavelengths, wavelengths)
      :staticmethod:



