CRTM_K_Matrix_Module.f90

Go to the documentation of this file.

!
! CRTM_K_Matrix_Module
!
! Module containing the CRTM K-matrix model function.
!
!
! CREATION HISTORY:
!       Written by:     Paul van Delst, 28-Jan-2005
!                       paul.vandelst@noaa.gov
!

MODULE crtm_k_matrix_module


  ! ------------
  ! Module usage
  ! ------------
  USE type_kinds,                 ONLY: fp
  USE message_handler,            ONLY: success, failure, warning, display_message
  USE crtm_parameters,            ONLY: set, not_set, zero, one, &
                                        max_n_layers           , &
                                        max_n_phase_elements   , &
                                        max_n_legendre_terms   , &
                                        max_n_stokes           , &
                                        max_n_angles           , &
                                        max_n_azimuth_fourier  , &
                                        max_source_zenith_angle, &
                                        max_n_streams          , &
                                        scattering_albedo_threshold
  USE crtm_spccoeff,              ONLY: sc, &
                                        spccoeff_isinfraredsensor , &
                                        spccoeff_ismicrowavesensor, &
                                        spccoeff_isvisiblesensor
  USE crtm_atmosphere_define,     ONLY: crtm_atmosphere_type, &
                                        crtm_atmosphere_destroy, &
                                        crtm_atmosphere_isvalid, &
                                        crtm_atmosphere_addlayercopy, &
                                        crtm_get_pressurelevelidx
  USE crtm_surface_define,        ONLY: crtm_surface_type, &
                                        crtm_surface_isvalid
  USE crtm_geometry_define,       ONLY: crtm_geometry_type, &
                                        crtm_geometry_isvalid
  USE crtm_channelinfo_define,    ONLY: crtm_channelinfo_type, &
                                        crtm_channelinfo_n_channels
  USE crtm_options_define,        ONLY: crtm_options_type, &
                                        crtm_options_isvalid
  USE crtm_atmosphere,            ONLY: crtm_atmosphere_addlayers, &
                                        crtm_atmosphere_addlayers_ad
  USE crtm_geometryinfo_define,   ONLY: crtm_geometryinfo_type, &
                                        crtm_geometryinfo_setvalue, &
                                        crtm_geometryinfo_getvalue
  USE crtm_geometryinfo,          ONLY: crtm_geometryinfo_compute
  USE crtm_predictor_define,      ONLY: crtm_predictor_type      , &
                                        crtm_predictor_associated, &
                                        crtm_predictor_destroy   , &
                                        crtm_predictor_create
  USE crtm_predictor,             ONLY: crtm_pvar_type => ivar_type, &
                                        crtm_compute_predictors    , &
                                        crtm_compute_predictors_ad
  USE crtm_atmabsorption,         ONLY: crtm_aavar_type => ivar_type , &
                                        crtm_compute_atmabsorption   , &
                                        crtm_compute_atmabsorption_ad
  USE crtm_atmoptics_define,      ONLY: crtm_atmoptics_type      , &
                                        crtm_atmoptics_associated, &
                                        crtm_atmoptics_create    , &
                                        crtm_atmoptics_destroy   , &
                                        crtm_atmoptics_zero
  USE crtm_aerosolscatter,        ONLY: crtm_compute_aerosolscatter   , &
                                        crtm_compute_aerosolscatter_ad
  USE crtm_cloudscatter,          ONLY: crtm_compute_cloudscatter   , &
                                        crtm_compute_cloudscatter_ad
  USE crtm_atmoptics,             ONLY: aovar_type  , &
                                        aovar_create, &
                                        crtm_no_scattering           , &
                                        crtm_include_scattering      , &
                                        crtm_compute_transmittance   , &
                                        crtm_compute_transmittance_ad, &
                                        crtm_combine_atmoptics       , &
                                        crtm_combine_atmoptics_ad
  USE crtm_sfcoptics_define,      ONLY: crtm_sfcoptics_type      , &
                                        crtm_sfcoptics_associated, &
                                        crtm_sfcoptics_create    , &
                                        crtm_sfcoptics_destroy
  USE crtm_sfcoptics,             ONLY: crtm_compute_surfacet   , &
                                        crtm_compute_surfacet_ad
  USE crtm_rtsolution,            ONLY: crtm_rtsolution_type      , &
                                        crtm_compute_nstreams     , &
                                        crtm_compute_rtsolution   , &
                                        crtm_compute_rtsolution_ad
  USE crtm_antennacorrection,     ONLY: crtm_compute_antcorr, &
                                        crtm_compute_antcorr_ad
  USE crtm_moleculescatter,       ONLY: crtm_compute_moleculescatter, &
                                        crtm_compute_moleculescatter_ad
  USE crtm_ancillaryinput_define, ONLY: crtm_ancillaryinput_type
  USE crtm_cloudcoeff,            ONLY: crtm_cloudcoeff_isloaded
  USE crtm_aerosolcoeff,          ONLY: crtm_aerosolcoeff_isloaded
  USE crtm_nltecorrection,        ONLY: nlte_predictor_type       , &
                                        nlte_predictor_isactive   , &
                                        compute_nlte_predictor    , &
                                        compute_nlte_predictor_ad , &
                                        compute_nlte_correction   , &
                                        compute_nlte_correction_ad
  USE accoeff_define,             ONLY: accoeff_associated
  USE nltecoeff_define,           ONLY: nltecoeff_associated
  USE crtm_planck_functions,      ONLY: crtm_planck_temperature   , &
                                        crtm_planck_temperature_ad

  ! Internal variable definition modules
  ! ...CloudScatter
  USE csvar_define, ONLY: csvar_type, &
                          csvar_associated, &
                          csvar_destroy   , &
                          csvar_create
  ! ...AerosolScatter
  USE asvar_define, ONLY: asvar_type, &
                          asvar_associated, &
                          asvar_destroy   , &
                          asvar_create
  ! ...Radiative transfer
  USE rtv_define,   ONLY: rtv_type      , &
                          rtv_associated, &
                          rtv_destroy   , &
                          rtv_create


  ! -----------------------
  ! Disable implicit typing
  ! -----------------------
  IMPLICIT NONE


  ! ------------
  ! Visibilities
  ! ------------
  ! Everything private by default
  PRIVATE
  ! Public procedures
  PUBLIC :: crtm_k_matrix
  PUBLIC :: crtm_k_matrix_version


  ! -----------------
  ! Module parameters
  ! -----------------
  ! Version Id for the module
  CHARACTER(*), PARAMETER :: module_version_id = &
  '$Id: CRTM_K_Matrix_Module.f90 60152 2015-08-13 19:19:13Z paul.vandelst@noaa.gov $'


CONTAINS


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_K_Matrix
!
! PURPOSE:
!       Function that calculates the K-matrix of top-of-atmosphere (TOA)
!       radiances and brightness temperatures for an input atmospheric
!       profile or profile set and user specified satellites/channels.
!
! CALLING SEQUENCE:
!       Error_Status = CRTM_K_Matrix( Atmosphere       , &
!                                     Surface          , &
!                                     RTSolution_K     , &
!                                     Geometry         , &
!                                     ChannelInfo      , &
!                                     Atmosphere_K     , &
!                                     Surface_K        , &
!                                     RTSolution       , &
!                                     Options = Options  )
!
! INPUTS:
!       Atmosphere:     Structure containing the Atmosphere data.
!                       UNITS:      N/A
!                       TYPE:       CRTM_Atmosphere_type
!                       DIMENSION:  Rank-1 (n_Profiles)
!                       ATTRIBUTES: INTENT(IN)
!
!       Surface:        Structure containing the Surface data.
!                       UNITS:      N/A
!                       TYPE:       CRTM_Surface_type
!                       DIMENSION:  Same as input Atmosphere argument.
!                       ATTRIBUTES: INTENT(IN)
!
!       RTSolution_K:   Structure containing the RT solution K-matrix inputs.
!                       **NOTE: On EXIT from this function, the contents of
!                               this structure may be modified (e.g. set to
!                               zero.)
!                       UNITS:      N/A
!                       TYPE:       CRTM_RTSolution_type
!                       DIMENSION:  Rank-2 (n_Channels x n_Profiles)
!                       ATTRIBUTES: INTENT(IN OUT)
!
!       Geometry:       Structure containing the view geometry
!                       information.
!                       UNITS:      N/A
!                       TYPE:       CRTM_Geometry_type
!                       DIMENSION:  Same as input Atmosphere argument
!                       ATTRIBUTES: INTENT(IN)
!
!       ChannelInfo:    Structure returned from the CRTM_Init() function
!                       that contains the satellite/sesnor channel index
!                       information.
!                       UNITS:      N/A
!                       TYPE:       CRTM_ChannelInfo_type
!                       DIMENSION:  Rank-1 (n_Sensors)
!                       ATTRIBUTES: INTENT(IN)
!
! OPTIONAL INPUTS:
!       Options:        Options structure containing the optional forward model
!                       arguments for the CRTM.
!                       UNITS:      N/A
!                       TYPE:       CRTM_Options_type
!                       DIMENSION:  Same as input Atmosphere structure
!                       ATTRIBUTES: INTENT(IN), OPTIONAL
!
! OUTPUTS:
!       Atmosphere_K:   Structure containing the K-matrix Atmosphere data.
!                       **NOTE: On ENTRY to this function, the contents of
!                               this structure should be defined (e.g.
!                               initialized to some value based on the
!                               position of this function in the call chain.)
!                       UNITS:      N/A
!                       TYPE:       CRTM_Atmosphere_type
!                       DIMENSION:  Same as input RTSolution_K argument
!                       ATTRIBUTES: INTENT(IN OUT)
!
!       Surface_K:      Structure containing the tangent-linear Surface data.
!                       **NOTE: On ENTRY to this function, the contents of
!                               this structure should be defined (e.g.
!                               initialized to some value based on the
!                               position of this function in the call chain.)
!                       UNITS:      N/A
!                       TYPE:       CRTM_Surface_type
!                       DIMENSION:  Same as input RTSolution_K argument
!                       ATTRIBUTES: INTENT(IN OUT)
!
!       RTSolution:     Structure containing the solution to the RT equation
!                       for the given inputs.
!                       UNITS:      N/A
!                       TYPE:       CRTM_RTSolution_type
!                       DIMENSION:  Same as input RTSolution_K argument
!                       ATTRIBUTES: INTENT(IN OUT)
!
! FUNCTION RESULT:
!       Error_Status:   The return value is an integer defining the error status.
!                       The error codes are defined in the Message_Handler module.
!                       If == SUCCESS the computation was sucessful
!                          == FAILURE an unrecoverable error occurred
!                       UNITS:      N/A
!                       TYPE:       INTEGER
!                       DIMENSION:  Scalar
!
! SIDE EFFECTS:
!      Note that the input K-matrix arguments are modified upon exit, and
!      the output K-matrix arguments must be defined upon entry. This is
!      a consequence of the K-matrix formulation where, effectively, the
!      chain rule is being used and this funtion could reside anywhere
!      in the chain of derivative terms.
!
! COMMENTS:
!       - The Options optional structure arguments contain
!         spectral information (e.g. emissivity) that must have the same
!         spectral dimensionality (the "L" dimension) as the RTSolution
!         structures.
!
!:sdoc-:
!--------------------------------------------------------------------------------

  FUNCTION crtm_k_matrix( &
    Atmosphere  , &  ! FWD Input, M
    Surface     , &  ! FWD Input, M
    RTSolution_K, &  ! K   Input, L x M
    Geometry    , &  ! Input, M
    ChannelInfo , &  ! Input, n_Sensors
    Atmosphere_K, &  ! K   Output, L x M
    Surface_K   , &  ! K   Output, L x M
    RTSolution  , &  ! FWD Output, L x M
    Options     ) &  ! Optional FWD input,  M
  result( error_status )
    ! Arguments
    TYPE(crtm_atmosphere_type)       , INTENT(IN)     :: atmosphere(:)     ! M
    TYPE(crtm_surface_type)          , INTENT(IN)     :: surface(:)        ! M
    TYPE(crtm_rtsolution_type)       , INTENT(IN OUT) :: rtsolution_k(:,:) ! L x M
    TYPE(crtm_geometry_type)         , INTENT(IN)     :: geometry(:)       ! M
    TYPE(crtm_channelinfo_type)      , INTENT(IN)     :: channelinfo(:)    ! n_Sensors
    TYPE(crtm_atmosphere_type)       , INTENT(IN OUT) :: atmosphere_k(:,:) ! L x M
    TYPE(crtm_surface_type)          , INTENT(IN OUT) :: surface_k(:,:)    ! L x M
    TYPE(crtm_rtsolution_type)       , INTENT(IN OUT) :: rtsolution(:,:)   ! L x M
    TYPE(crtm_options_type), OPTIONAL, INTENT(IN)     :: options(:)        ! M
    ! Function result
    INTEGER :: error_status
    ! Local parameters
    CHARACTER(*), PARAMETER :: routine_name = 'CRTM_K_Matrix'
    ! Local variables
    CHARACTER(256) :: message
    LOGICAL :: options_present
    LOGICAL :: check_input
    LOGICAL :: user_emissivity, user_direct_reflectivity, user_n_streams
    LOGICAL :: user_antcorr, compute_antcorr
    LOGICAL :: apply_nlte_correction
    LOGICAL :: atmosphere_invalid, surface_invalid, geometry_invalid, options_invalid
    INTEGER :: rt_algorithm_id
    INTEGER :: ifov
    INTEGER :: nc, na
    INTEGER :: n, n_sensors,  sensorindex
    INTEGER :: l, n_channels, channelindex
    INTEGER :: m, n_profiles
    INTEGER :: j, ln
    INTEGER :: n_full_streams, mth_azi
    REAL(fp) :: source_za
    REAL(fp) :: wavenumber
    REAL(fp) :: transmittance, transmittance_k
    ! Local ancillary input structure
    TYPE(crtm_ancillaryinput_type) :: ancillaryinput
    ! Local options structure for default values
    TYPE(crtm_options_type) :: default_options
    ! Local atmosphere structure for extra layering
    TYPE(crtm_atmosphere_type) :: atm, atm_k
    ! Component variables
    TYPE(crtm_geometryinfo_type) :: geometryinfo
    TYPE(crtm_predictor_type)    :: predictor, predictor_k
    TYPE(crtm_atmoptics_type)    :: atmoptics, atmoptics_k
    TYPE(crtm_sfcoptics_type)    :: sfcoptics, sfcoptics_k
    ! Component variable internals
    TYPE(crtm_pvar_type)  :: pvar   ! Predictor
    TYPE(crtm_aavar_type) :: aavar  ! AtmAbsorption
    TYPE(csvar_type)      :: csvar  ! CloudScatter
    TYPE(asvar_type)      :: asvar  ! AerosolScatter
    TYPE(aovar_type)      :: aovar  ! AtmOptics
    TYPE(rtv_type)        :: rtv    ! RTSolution
    ! NLTE correction term predictors
    TYPE(nlte_predictor_type)   :: nlte_predictor, nlte_predictor_k

    ! ------
    ! SET UP
    ! ------
    error_status = success


    ! If no sensors or channels, simply return
    n_sensors  = SIZE(channelinfo)
    n_channels = sum(crtm_channelinfo_n_channels(channelinfo))
    IF ( n_sensors == 0 .OR. n_channels == 0 ) RETURN


    ! Check spectral arrays
    IF ( SIZE(rtsolution  ,dim=1) < n_channels .OR. &
         SIZE(atmosphere_k,dim=1) < n_channels .OR. &
         SIZE(surface_k   ,dim=1) < n_channels .OR. &
         SIZE(rtsolution_k,dim=1) < n_channels      ) THEN
      error_status = failure
      WRITE( message,'("RTSolution and K-matrix (Atm,Sfc,RT) structure arrays too small (",&
             &3(i0,","),i0,") for the number of requested channels (",i0,")")') &
             SIZE(rtsolution  ,dim=1), &
             SIZE(atmosphere_k,dim=1), &
             SIZE(surface_k   ,dim=1), &
             SIZE(rtsolution_k,dim=1), &
             n_channels
      CALL display_message( routine_name, message, error_status )
      RETURN
    END IF


    ! Check the number of profiles
    ! ...Number of atmospheric profiles.
    n_profiles = SIZE(atmosphere)
    ! ...Check the profile dimensionality of the other mandatory arguments
    IF ( SIZE(surface)            /= n_profiles .OR. &
         SIZE(rtsolution_k,dim=2) /= n_profiles .OR. &
         SIZE(geometry)           /= n_profiles .OR. &
         SIZE(atmosphere_k,dim=2) /= n_profiles .OR. &
         SIZE(surface_k   ,dim=2) /= n_profiles .OR. &
         SIZE(rtsolution  ,dim=2) /= n_profiles      ) THEN
      error_status = failure
      message = 'Inconsistent profile dimensionality for input arguments.'
      CALL display_message( routine_name, message, error_status )
      RETURN
    END IF
    ! ...Check the profile dimensionality of the other optional arguments
    options_present = .false.
    IF ( PRESENT( options ) ) THEN
      options_present = .true.
      IF ( SIZE( options ) /= n_profiles ) THEN
        error_status = failure
        message = 'Inconsistent profile dimensionality for Options optional input argument.'
        CALL display_message( routine_name, message, error_status )
        RETURN
      END IF
    END IF


    ! Allocate the profile independent surface optics local structure
    CALL crtm_sfcoptics_create( sfcoptics  , max_n_angles, max_n_stokes )
    CALL crtm_sfcoptics_create( sfcoptics_k, max_n_angles, max_n_stokes )
    IF ( (.NOT. crtm_sfcoptics_associated(sfcoptics  )) .OR. &
         (.NOT. crtm_sfcoptics_associated(sfcoptics_k)) ) THEN
      error_status = failure
      message = 'Error allocating SfcOptics data structures'
      CALL display_message( routine_name, message, error_status )
      RETURN
    END IF


    ! ------------
    ! PROFILE LOOP
    ! ------------
    profile_loop: DO m = 1, n_profiles


      ! Check the cloud and aerosol coeff. data for cases with clouds and aerosol
      IF( atmosphere(m)%n_Clouds > 0 .AND. .NOT. crtm_cloudcoeff_isloaded() )THEN
         error_status = failure
         WRITE( message,'("The CloudCoeff data must be loaded (with CRTM_Init routine) ", &
                &"for the cloudy case profile #",i0)' ) m
         CALL display_message( routine_name, message, error_status )
         RETURN
      END IF
      IF( atmosphere(m)%n_Aerosols > 0 .AND. .NOT. crtm_aerosolcoeff_isloaded() )THEN
         error_status = failure
         WRITE( message,'("The AerosolCoeff data must be loaded (with CRTM_Init routine) ", &
                &"for the aerosol case profile #",i0)' ) m
         CALL display_message( routine_name, message, error_status )
         RETURN
      END IF


      ! Copy over forward "non-variable" inputs to K-matrix outputs
      DO l = 1, n_channels
        ! ...Atmosphere
        atmosphere_k(l,m)%Climatology = atmosphere(m)%Climatology
        ! Loop over absorbers
        DO j = 1, atmosphere(m)%n_Absorbers
          atmosphere_k(l,m)%Absorber_ID(j)    = atmosphere(m)%Absorber_ID(j)
          atmosphere_k(l,m)%Absorber_Units(j) = atmosphere(m)%Absorber_Units(j)
        END DO
        ! Loop over and assign cloud types
        DO nc = 1, atmosphere(m)%n_Clouds
          atmosphere_k(l,m)%Cloud(nc)%Type = atmosphere(m)%Cloud(nc)%Type
        END DO
        ! Loop over and assign aerosol types
        DO na = 1, atmosphere(m)%n_Aerosols
          atmosphere_k(l,m)%Aerosol(na)%Type = atmosphere(m)%Aerosol(na)%Type
        END DO
        ! ...Surface
        surface_k(l,m)%Land_Coverage  = surface(m)%Land_Coverage
        surface_k(l,m)%Water_Coverage = surface(m)%Water_Coverage
        surface_k(l,m)%Snow_Coverage  = surface(m)%Snow_Coverage
        surface_k(l,m)%Ice_Coverage   = surface(m)%Ice_Coverage
        surface_k(l,m)%Land_Type  = surface(m)%Land_Type
        surface_k(l,m)%Water_Type = surface(m)%Water_Type
        surface_k(l,m)%Snow_Type  = surface(m)%Snow_Type
        surface_k(l,m)%Ice_Type   = surface(m)%Ice_Type
      END DO


      ! Check the optional Options structure argument
      ! ...Specify default actions
      check_input           = default_options%Check_Input
      user_emissivity       = default_options%Use_Emissivity
      user_antcorr          = default_options%Use_Antenna_Correction
      apply_nlte_correction = default_options%Apply_NLTE_Correction
      rt_algorithm_id       = default_options%RT_Algorithm_Id
      user_n_streams        = default_options%Use_N_Streams
      ! ...Check the Options argument
      IF (options_present) THEN
        ! Override input checker with option
        check_input = options(m)%Check_Input
        ! Check if the supplied emissivity should be used
        user_emissivity = options(m)%Use_Emissivity
        IF ( options(m)%Use_Emissivity ) THEN
          ! Are the channel dimensions consistent
          IF ( options(m)%n_Channels < n_channels ) THEN
            error_status = failure
            WRITE( message,'( "Input Options channel dimension (", i0, ") is less ", &
                   &"than the number of requested channels (",i0, ")" )' ) &
                   options(m)%n_Channels, n_channels
            CALL display_message( routine_name, message, error_status )
            RETURN
          END IF
          ! Check if the supplied direct reflectivity should be used
          user_direct_reflectivity = options(m)%Use_Direct_Reflectivity
        END IF
        ! Check if antenna correction should be attempted
        user_antcorr = options(m)%Use_Antenna_Correction
        ! Set NLTE correction option
        apply_nlte_correction = options(m)%Apply_NLTE_Correction

        ! Copy over ancillary input
        ancillaryinput%SSU    = options(m)%SSU
        ancillaryinput%Zeeman = options(m)%Zeeman
        ! Copy over surface optics input
        sfcoptics%Use_New_MWSSEM = .NOT. options(m)%Use_Old_MWSSEM
        ! Specify the RT algorithm
        rt_algorithm_id = options(m)%RT_Algorithm_ID
        ! Check if n_Streams should be used
        user_n_streams = options(m)%Use_N_Streams
        ! Check value for nstreams
        IF ( user_n_streams ) THEN
          IF ( options(m)%n_Streams <= 0 .OR. mod(options(m)%n_Streams,2) /= 0 .OR. &
               options(m)%n_Streams > max_n_streams ) THEN
              error_status = failure
              WRITE( message,'( "Input Options n_Streams (", i0, ") is invalid" )' ) &
                     options(m)%n_Streams
              CALL display_message( routine_name, message, error_status )
              RETURN
          END IF
        END IF
      END IF


      ! Check the input data if required
      IF ( check_input ) THEN
        ! ...Mandatory inputs
        atmosphere_invalid = .NOT. crtm_atmosphere_isvalid( atmosphere(m) )
        surface_invalid    = .NOT. crtm_surface_isvalid( surface(m) )
        geometry_invalid   = .NOT. crtm_geometry_isvalid( geometry(m) )
        IF ( atmosphere_invalid .OR. surface_invalid .OR. geometry_invalid ) THEN
          error_status = failure
          WRITE( message,'("Input data check failed for profile #",i0)' ) m
          CALL display_message( routine_name, message, error_status )
          RETURN
        END IF
        ! ...Optional input
        IF ( options_present ) THEN
          options_invalid = .NOT. crtm_options_isvalid( options(m) )
          IF ( options_invalid ) THEN
            error_status = failure
            WRITE( message,'("Options data check failed for profile #",i0)' ) m
            CALL display_message( routine_name, message, error_status )
            RETURN
          END IF
        END IF
      END IF


      ! Process geometry
      ! ...Compute derived geometry
      CALL crtm_geometryinfo_setvalue( geometryinfo, geometry=geometry(m) )
      CALL crtm_geometryinfo_compute( geometryinfo )
      ! ...Retrieve components into local variable
      CALL crtm_geometryinfo_getvalue( &
             geometryinfo, &
             ifov = ifov, &
             source_zenith_angle = source_za )


      ! Average surface skin temperature for multi-surface types
      CALL crtm_compute_surfacet( surface(m), sfcoptics )


      ! Add extra layers to current atmosphere profile
      ! if necessary to handle upper atmosphere
      error_status = crtm_atmosphere_addlayers( atmosphere(m), atm )
      IF ( error_status /= success ) THEN
        error_status = failure
        WRITE( message,'("Error adding FWD extra layers to profile #",i0)' ) m
        CALL display_message( routine_name, message, error_status )
        RETURN
      END IF
      ! ...Check the total number of Atm layers
      IF ( atm%n_Layers > max_n_layers ) THEN
        error_status = failure
        WRITE( message,'("Added layers [",i0,"] cause total [",i0,"] to exceed the ",&
               &"maximum allowed [",i0,"] for profile #",i0)' ) &
               atm%n_Added_Layers, atm%n_Layers, max_n_layers, m
        CALL display_message( routine_name, message, error_status )
        RETURN
      END IF
      ! ...Allocate the atmospheric optics structures based on Atm extension
      CALL crtm_atmoptics_create( atmoptics, &
                                  atm%n_Layers        , &
                                  max_n_legendre_terms, &
                                  max_n_phase_elements  )
      CALL crtm_atmoptics_create( atmoptics_k, &
                                  atm%n_Layers        , &
                                  max_n_legendre_terms, &
                                  max_n_phase_elements  )
      IF ( .NOT. crtm_atmoptics_associated( atmoptics ) .OR. &
           .NOT. crtm_atmoptics_associated( atmoptics_k ) ) THEN
        error_status = failure
        WRITE( message,'("Error allocating AtmOptics data structures for profile #",i0)' ) m
        CALL display_message( routine_name, message, error_status )
        RETURN
      END IF
      IF (options_present) THEN
        ! Set Scattering Switch
        atmoptics%Include_Scattering   = options(m)%Include_Scattering
        atmoptics_k%Include_Scattering = options(m)%Include_Scattering
      END IF
      ! ...Allocate the atmospheric optics internal structure
      CALL aovar_create( aovar, atm%n_Layers )


      ! Allocate the scattering internal variables if necessary
      ! ...Cloud
      IF ( atm%n_Clouds > 0 ) THEN
        CALL csvar_create( csvar, &
                           max_n_legendre_terms, &
                           max_n_phase_elements, &
                           atm%n_Layers        , &
                           atm%n_Clouds          )
      END IF
      ! ...Aerosol
      IF ( atm%n_Aerosols > 0 ) THEN
        CALL asvar_create( asvar, &
                           max_n_legendre_terms, &
                           max_n_phase_elements, &
                           atm%n_Layers        , &
                           atm%n_Aerosols        )
      END IF


      ! -----------
      ! SENSOR LOOP
      ! -----------
      ! Initialise channel counter for sensor(n)/channel(l) count
      ln = 0

      sensor_loop: DO n = 1, n_sensors


        ! Shorter name
        sensorindex = channelinfo(n)%Sensor_Index


        ! Check if antenna correction to be applied for current sensor
        IF ( user_antcorr                             .AND. &
             accoeff_associated( sc(sensorindex)%AC ) .AND. &
             ifov /= 0 ) THEN
          compute_antcorr = .true.
        ELSE
          compute_antcorr = .false.
        END IF


        ! Allocate the AtmAbsorption predictor structures
        CALL crtm_predictor_create( &
               predictor   , &
               atm%n_Layers, &
               sensorindex , &
               savefwv = 1   )
        CALL crtm_predictor_create( &
               predictor_k , &
               atm%n_Layers, &
               sensorindex   )
        IF ( (.NOT. crtm_predictor_associated(predictor)) .OR. &
             (.NOT. crtm_predictor_associated(predictor_k)) ) THEN
          error_status=failure
          WRITE( message,'("Error allocating predictor structures for profile #",i0, &
                 &" and ",a," sensor.")' ) m, sc(sensorindex)%Sensor_Id
          CALL display_message( routine_name, message, error_status )
          RETURN
        END IF


        ! ...Compute forward predictors
        CALL crtm_compute_predictors( sensorindex   , &  ! Input
                                      atm           , &  ! Input
                                      geometryinfo  , &  ! Input
                                      ancillaryinput, &  ! Input
                                      predictor     , &  ! Output
                                      pvar            )  ! Internal variable output


        ! Allocate the RTV structure if necessary
        IF( (atm%n_Clouds   > 0 .OR. &
            atm%n_Aerosols > 0 .OR. &
            spccoeff_isvisiblesensor( sc(sensorindex) ) ) .and. atmoptics%Include_Scattering ) THEN
          CALL rtv_create( rtv, max_n_angles, max_n_legendre_terms, atm%n_Layers )
          IF ( .NOT. rtv_associated(rtv) ) THEN
            error_status=failure
            WRITE( message,'("Error allocating RTV structure for profile #",i0, &
                   &" and ",a," sensor.")' ) m, trim(sc(sensorindex)%Sensor_Id)
            CALL display_message( routine_name, message, error_status )
            RETURN
          END IF
          ! Assign algorithm selector
          rtv%RT_Algorithm_Id = rt_algorithm_id
        END IF


        ! Compute NLTE predictors
        IF ( apply_nlte_correction ) THEN
          CALL compute_nlte_predictor( &
                 sc(sensorindex)%NC, &  ! Input
                 atm               , &  ! Input
                 geometryinfo      , &  ! Input
                 nlte_predictor      )  ! Output
        END IF


        ! ------------
        ! CHANNEL LOOP
        ! ------------
        channel_loop: DO l = 1, channelinfo(n)%n_Channels


          ! Channel setup
          ! ...Skip channel if requested
          IF ( .NOT. channelinfo(n)%Process_Channel(l) ) cycle channel_loop
          ! ...Shorter name
          channelindex = channelinfo(n)%Channel_Index(l)
          ! ...Increment the processed channel counter
          ln = ln + 1
          ! ...Assign sensor+channel information to output
          rtsolution(ln,m)%Sensor_Id        = channelinfo(n)%Sensor_Id
          rtsolution(ln,m)%WMO_Satellite_Id = channelinfo(n)%WMO_Satellite_Id
          rtsolution(ln,m)%WMO_Sensor_Id    = channelinfo(n)%WMO_Sensor_Id
          rtsolution(ln,m)%Sensor_Channel   = channelinfo(n)%Sensor_Channel(l)
          rtsolution_k(ln,m)%Sensor_Id        = rtsolution(ln,m)%Sensor_Id
          rtsolution_k(ln,m)%WMO_Satellite_Id = rtsolution(ln,m)%WMO_Satellite_Id
          rtsolution_k(ln,m)%WMO_Sensor_Id    = rtsolution(ln,m)%WMO_Sensor_Id
          rtsolution_k(ln,m)%Sensor_Channel   = rtsolution(ln,m)%Sensor_Channel


          ! Initialisations
          CALL crtm_atmoptics_zero( atmoptics )
          transmittance_k = zero


          ! Copy the input K-matrix atmosphere with extra layers if necessary
          atm_k = crtm_atmosphere_addlayercopy( atmosphere_k(ln,m), atm%n_Added_Layers )


          ! Determine the number of streams (n_Full_Streams) in up+downward directions
          IF ( user_n_streams ) THEN
            n_full_streams = options(m)%n_Streams
            rtsolution(ln,m)%n_Full_Streams = n_full_streams + 2
            rtsolution(ln,m)%Scattering_Flag = .true.
          ELSE
            n_full_streams = crtm_compute_nstreams( atm             , &  ! Input
                                                    sensorindex     , &  ! Input
                                                    channelindex    , &  ! Input
                                                    rtsolution(ln,m)  )  ! Output
          END IF
          ! ...Transfer stream count to scattering structures
          atmoptics%n_Legendre_Terms   = n_full_streams
          atmoptics_k%n_Legendre_Terms = n_full_streams


          ! Compute the gas absorption
          CALL crtm_compute_atmabsorption( sensorindex   , &  ! Input
                                           channelindex  , &  ! Input
                                           ancillaryinput, &  ! Input
                                           predictor     , &  ! Input
                                           atmoptics     , &  ! Output
                                           aavar           )  ! Internal variable output


          ! Compute the clear-sky atmospheric transmittance
          ! for use in FASTEM-X reflection correction
          CALL crtm_compute_transmittance(atmoptics,transmittance)


          ! Compute the molecular scattering properties
          ! ...Solar radiation
          IF( sc(sensorindex)%Solar_Irradiance(channelindex) > zero .AND. &
              source_za < max_source_zenith_angle) THEN
             rtv%Solar_Flag_true = .true.
          END IF
          ! ...Visible channel with solar radiation
          IF( spccoeff_isvisiblesensor( sc(sensorindex) ) .AND. rtv%Solar_Flag_true ) THEN
            rtv%Visible_Flag_true = .true.
            ! Rayleigh phase function has 0, 1, 2 components.
            IF( atmoptics%n_Legendre_Terms < 4 ) THEN
              atmoptics%n_Legendre_Terms = 4
              atmoptics_k%n_Legendre_Terms = atmoptics%n_Legendre_Terms
              rtsolution(ln,m)%Scattering_FLAG = .true.
              rtsolution(ln,m)%n_Full_Streams = atmoptics%n_Legendre_Terms + 2
            END IF
            rtv%n_Azi = min( atmoptics%n_Legendre_Terms - 1, max_n_azimuth_fourier )
            ! Get molecular scattering and extinction
            wavenumber = sc(sensorindex)%Wavenumber(channelindex)
            error_status = crtm_compute_moleculescatter( &
                             wavenumber, &
                             atm       , &
                             atmoptics   )
            IF ( error_status /= success ) THEN
              WRITE( message,'("Error computing MoleculeScatter for ",a,&
                     &", channel ",i0,", profile #",i0)') &
                     trim(channelinfo(n)%Sensor_ID), &
                     channelinfo(n)%Sensor_Channel(l), &
                     m
              CALL display_message( routine_name, message, error_status )
              RETURN
            END IF
          ELSE
            rtv%Visible_Flag_true = .false.
            rtv%n_Azi = 0
          END IF


          ! Compute the cloud particle absorption/scattering properties
          IF( atm%n_Clouds > 0 ) THEN
            error_status = crtm_compute_cloudscatter( atm         , &  ! Input
                                                      sensorindex , &  ! Input
                                                      channelindex, &  ! Input
                                                      atmoptics   , &  ! Output
                                                      csvar         )  ! Internal variable output
            IF (error_status /= success) THEN
              WRITE( message,'("Error computing CloudScatter for ",a,&
                     &", channel ",i0,", profile #",i0)' ) &
                     trim(channelinfo(n)%Sensor_ID), channelinfo(n)%Sensor_Channel(l), m
              CALL display_message( routine_name, message, error_status )
              RETURN
            END IF
          END IF


          ! Compute the aerosol absorption/scattering properties
          IF ( atm%n_Aerosols > 0 ) THEN
            error_status = crtm_compute_aerosolscatter( atm         , &  ! Input
                                                        sensorindex , &  ! Input
                                                        channelindex, &  ! Input
                                                        atmoptics   , &  ! In/Output
                                                        asvar         )  ! Internal variable output
            IF ( error_status /= success ) THEN
              WRITE( message,'("Error computing AerosolScatter for ",a,&
                     &", channel ",i0,", profile #",i0)' ) &
                     trim(channelinfo(n)%Sensor_ID), channelinfo(n)%Sensor_Channel(l), m
              CALL display_message( routine_name, message, error_status )
              RETURN
            END IF
          END IF


          ! Compute the combined atmospheric optical properties
          IF( atmoptics%Include_Scattering ) THEN
            CALL crtm_combine_atmoptics( atmoptics, aovar )
          END IF
          ! ...Save vertically integrated scattering optical depth for output
          rtsolution(ln,m)%SOD = atmoptics%Scattering_Optical_Depth


          ! Turn off FASTEM-X reflection correction for scattering conditions
          IF ( crtm_include_scattering(atmoptics) .AND. spccoeff_ismicrowavesensor( sc(sensorindex) ) ) THEN
            sfcoptics%Transmittance = -one
          ELSE
            sfcoptics%Transmittance = transmittance
          END IF


          ! Fill the SfcOptics structure for the optional emissivity input case.
          ! ...Indicate SfcOptics ARE to be computed
          sfcoptics%Compute = .true.
          ! Change SfcOptics emissivity/reflectivity contents/computation status
          IF ( user_emissivity ) THEN
            sfcoptics%Compute = .false.
            sfcoptics%Emissivity(1,1)       = options(m)%Emissivity(ln)
            sfcoptics%Reflectivity(1,1,1,1) = one - options(m)%Emissivity(ln)
            IF ( user_direct_reflectivity ) THEN
              sfcoptics%Direct_Reflectivity(1,1) = options(m)%Direct_Reflectivity(ln)
            ELSE
              sfcoptics%Direct_Reflectivity(1,1) = sfcoptics%Reflectivity(1,1,1,1)
            END IF

          END IF


          ! Fourier component loop for azimuth angles (VIS).
          ! mth_Azi = 0 is for an azimuth-averaged value (IR, MW)
          ! ...Initialise radiance
          rtsolution(ln,m)%Radiance = zero
          ! ...Initialise K-matrix atmospheric optics
          CALL crtm_atmoptics_zero( atmoptics_k )



          ! ###################################################
          ! TEMPORARY FIX : SENSOR-DEPENDENT AZIMUTH ANGLE LOOP
          ! ###################################################
          sensor_dependent_rtsolution: &
          IF ( spccoeff_isinfraredsensor( sc(sensorindex) ) .OR. &
               spccoeff_ismicrowavesensor( sc(sensorindex) ) ) THEN
            ! ------------------------------
            ! INFRARED and MICROWAVE sensors
            ! ------------------------------

            ! Set dependent component counters
            rtv%mth_Azi = 0
            sfcoptics%mth_Azi = 0

            ! Solve the forward radiative transfer problem
            error_status = crtm_compute_rtsolution( &
                             atm             , &  ! Input
                             surface(m)      , &  ! Input
                             atmoptics       , &  ! Input
                             sfcoptics       , &  ! Input
                             geometryinfo    , &  ! Input
                             sensorindex     , &  ! Input
                             channelindex    , &  ! Input
                             rtsolution(ln,m), &  ! Output
                             rtv               )  ! Internal variable output
            IF ( error_status /= success ) THEN
              WRITE( message,'( "Error computing RTSolution for ", a, &
                     &", channel ", i0,", profile #",i0)' ) &
                     trim(channelinfo(n)%Sensor_ID), channelinfo(n)%Sensor_Channel(l), m
              CALL display_message( routine_name, message, error_status )
              RETURN
            END IF

            ! Compute non-LTE correction to radiance if required
            IF ( apply_nlte_correction .AND. nlte_predictor_isactive(nlte_predictor) ) &
              CALL compute_nlte_correction( &
                     sc(sensorindex)%NC       , &  ! Input
                     channelindex             , &  ! Input
                     nlte_predictor           , &  ! Input
                     rtsolution(ln,m)%Radiance  )  ! In/Output

            ! Convert the radiance to brightness temperature
            CALL crtm_planck_temperature( &
                   sensorindex                            , & ! Input
                   channelindex                           , & ! Input
                   rtsolution(ln,m)%Radiance              , & ! Input
                   rtsolution(ln,m)%Brightness_Temperature  ) ! Output

            ! Compute Antenna correction to brightness temperature if required
            IF ( compute_antcorr ) THEN
              CALL crtm_compute_antcorr( &
                     geometryinfo    , &  ! Input
                     sensorindex     , &  ! Input
                     channelindex    , &  ! Input
                     rtsolution(ln,m)  )  ! Output
              CALL crtm_compute_antcorr_ad( &
                     geometryinfo      , &  ! Input
                     sensorindex       , &  ! Input
                     channelindex      , &  ! Input
                     rtsolution_k(ln,m)  )  ! Output
            END IF

            ! Compute the Planck temperature adjoijnt
            CALL crtm_planck_temperature_ad( &
                   sensorindex                              , & ! Input
                   channelindex                             , & ! Input
                   rtsolution(ln,m)%Radiance                , & ! Input
                   rtsolution_k(ln,m)%Brightness_Temperature, & ! Input
                   rtsolution_k(ln,m)%Radiance                ) ! Output
            rtsolution_k(ln,m)%Brightness_Temperature = zero

            ! Compute non-LTE correction adjoint if required
            IF ( apply_nlte_correction .AND. nlte_predictor_isactive(nlte_predictor) ) &
              CALL compute_nlte_correction_ad( &
                     sc(sensorindex)%NC         , &  ! Input
                     channelindex               , &  ! Input
                     rtsolution_k(ln,m)%Radiance, &  ! Input
                     nlte_predictor_k             )  ! Output

            ! The adjoint of the radiative transfer
            error_status = crtm_compute_rtsolution_ad( &
                             atm               , &  ! FWD Input
                             surface(m)        , &  ! FWD Input
                             atmoptics         , &  ! FWD Input
                             sfcoptics         , &  ! FWD Input
                             rtsolution(ln,m)  , &  ! FWD Input
                             rtsolution_k(ln,m), &  ! K   Input
                             geometryinfo      , &  ! Input
                             sensorindex       , &  ! Input
                             channelindex      , &  ! Input
                             atm_k             , &  ! K  Output
                             surface_k(ln,m)   , &  ! K  Output
                             atmoptics_k       , &  ! K  Output
                             sfcoptics_k       , &  ! K  Output
                             rtv                 )  ! Internal variable input
            IF ( error_status /= success ) THEN
              WRITE( message,'( "Error computing RTSolution_K for ", a, &
                     &", channel ", i0,", profile #",i0)' ) &
                     trim(channelinfo(n)%Sensor_ID), channelinfo(n)%Sensor_Channel(l), m
              CALL display_message( routine_name, message, error_status )
              RETURN
            END IF


          ELSE sensor_dependent_rtsolution
            ! --------------
            ! VISIBLE sensor
            ! --------------
            ! ...Fourier expansion over azimuth angle
            azimuth_fourier_loop: DO mth_azi = 0, rtv%n_Azi

              ! Set dependent component counters
              rtv%mth_Azi = mth_azi
              sfcoptics%mth_Azi = mth_azi

              ! Solve the forward radiative transfer problem
              error_status = crtm_compute_rtsolution( &
                               atm             , &  ! Input
                               surface(m)      , &  ! Input
                               atmoptics       , &  ! Input
                               sfcoptics       , &  ! Input
                               geometryinfo    , &  ! Input
                               sensorindex     , &  ! Input
                               channelindex    , &  ! Input
                               rtsolution(ln,m), &  ! Output
                               rtv               )  ! Internal variable output
              IF ( error_status /= success ) THEN
                WRITE( message,'( "Error computing RTSolution for ", a, &
                       &", channel ", i0,", profile #",i0)' ) &
                       trim(channelinfo(n)%Sensor_ID), channelinfo(n)%Sensor_Channel(l), m
                CALL display_message( routine_name, message, error_status )
                RETURN
              END IF

              ! The adjoint of the radiative transfer
              error_status = crtm_compute_rtsolution_ad( &
                               atm               , &  ! FWD Input
                               surface(m)        , &  ! FWD Input
                               atmoptics         , &  ! FWD Input
                               sfcoptics         , &  ! FWD Input
                               rtsolution(ln,m)  , &  ! FWD Input
                               rtsolution_k(ln,m), &  ! K   Input
                               geometryinfo      , &  ! Input
                               sensorindex       , &  ! Input
                               channelindex      , &  ! Input
                               atm_k             , &  ! K  Output
                               surface_k(ln,m)   , &  ! K  Output
                               atmoptics_k       , &  ! K  Output
                               sfcoptics_k       , &  ! K  Output
                               rtv                 )  ! Internal variable input
              IF ( error_status /= success ) THEN
                WRITE( message,'( "Error computing RTSolution_K for ", a, &
                       &", channel ", i0,", profile #",i0)' ) &
                       trim(channelinfo(n)%Sensor_ID), channelinfo(n)%Sensor_Channel(l), m
                CALL display_message( routine_name, message, error_status )
                RETURN
              END IF
            END DO azimuth_fourier_loop

            ! Still want to convert the final FORWARD radiance to brightness temperature
            CALL crtm_planck_temperature( &
                   sensorindex                            , & ! Input
                   channelindex                           , & ! Input
                   rtsolution(ln,m)%Radiance              , & ! Input
                   rtsolution(ln,m)%Brightness_Temperature  ) ! Output

          END IF sensor_dependent_rtsolution
          ! ###################################################
          ! TEMPORARY FIX : SENSOR-DEPENDENT AZIMUTH ANGLE LOOP
          ! ###################################################


          ! Compute the adjoint of the combined atmospheric optical properties
          IF( atmoptics%Include_Scattering ) THEN
            CALL crtm_combine_atmoptics_ad( atmoptics, atmoptics_k, aovar )
          END IF


          ! Compute the adjoint aerosol absorption/scattering properties
          IF ( atm%n_Aerosols > 0 ) THEN
            error_status = crtm_compute_aerosolscatter_ad( atm         , &  ! FWD Input
                                                           atmoptics   , &  ! FWD Input
                                                           atmoptics_k , &  ! K   Input
                                                           sensorindex , &  ! Input
                                                           channelindex, &  ! Input
                                                           atm_k       , &  ! K   Output
                                                           asvar         )  ! Internal variable input
            IF ( error_status /= success ) THEN
              WRITE( message,'("Error computing AerosolScatter_K for ",a,&
                     &", channel ",i0,", profile #",i0)' ) &
                     trim(channelinfo(n)%Sensor_ID), channelinfo(n)%Sensor_Channel(l), m
              CALL display_message( routine_name, message, error_status )
              RETURN
            END IF
          END IF


          ! Compute the adjoint cloud absorption/scattering properties
          IF ( atm%n_Clouds > 0 ) THEN
            error_status = crtm_compute_cloudscatter_ad( atm         , &  ! FWD Input
                                                         atmoptics   , &  ! FWD Input
                                                         atmoptics_k , &  ! K   Input
                                                         sensorindex , &  ! Input
                                                         channelindex, &  ! Input
                                                         atm_k       , &  ! K   Output
                                                         csvar         )  ! Internal variable input
            IF ( error_status /= success ) THEN
              WRITE( message,'("Error computing CloudScatter_K for ",a,&
                     &", channel ",i0,", profile #",i0)' ) &
                     trim(channelinfo(n)%Sensor_ID), channelinfo(n)%Sensor_Channel(l), m
              CALL display_message( routine_name, message, error_status )
              RETURN
            END IF
          END IF


          ! Compute the adjoint molecular scattering properties
          IF( rtv%Visible_Flag_true ) THEN
            wavenumber = sc(sensorindex)%Wavenumber(channelindex)
            error_status = crtm_compute_moleculescatter_ad( &
                             wavenumber , &
                             atmoptics_k, &
                             atm_k        )
            IF ( error_status /= success ) THEN
              WRITE( message,'("Error computing MoleculeScatter_K for ",a,&
                     &", channel ",i0,", profile #",i0)' ) &
                     trim(channelinfo(n)%Sensor_ID), &
                     channelinfo(n)%Sensor_Channel(l), &
                     m
              CALL display_message( routine_name, message, error_status )
              RETURN
            END IF
          END IF


          ! Compute the adjoint of the total atmospheric transmittance
          IF ( crtm_no_scattering(atmoptics) .AND. spccoeff_ismicrowavesensor(sc(sensorindex)) ) THEN
            transmittance_k = sfcoptics_k%transmittance
            sfcoptics_k%transmittance = zero
            CALL crtm_compute_transmittance_ad(atmoptics,transmittance_k,atmoptics_k)
          END IF


          ! Compute the adjoint gaseous absorption
          CALL crtm_compute_atmabsorption_ad( sensorindex , &  ! Input
                                              channelindex, &  ! Input
                                              predictor   , &  ! FWD Input
                                              atmoptics_k , &  ! K   Input
                                              predictor_k , &  ! K   Output
                                              aavar         )  ! Internal variable input


          ! K-matrix of the NLTE correction predictor calculations
          IF ( apply_nlte_correction ) THEN
            CALL compute_nlte_predictor_ad( &
                   nlte_predictor   , &  ! Input
                   nlte_predictor_k , &  ! Input
                   atm_k              )  ! Output
          END IF


          ! K-matrix of the predictor calculations
          CALL crtm_compute_predictors_ad( sensorindex   , &  ! Input
                                           atm           , &  ! FWD Input
                                           predictor     , &  ! FWD Input
                                           predictor_k   , &  ! K   Input
                                           ancillaryinput, &  ! Input
                                           atm_k         , &  ! K   Output
                                           pvar            )  ! Internal variable input


          ! Postprocess some input data
          ! ...K-matrix of average surface skin temperature for multi-surface types
          CALL crtm_compute_surfacet_ad( surface(m), sfcoptics_k, surface_k(ln,m) )
          ! ...K-matrix of the atmosphere layer addition
          error_status = crtm_atmosphere_addlayers_ad( atmosphere(m), atm_k, atmosphere_k(ln,m) )
          IF ( error_status /= success ) THEN
            error_status = failure
            WRITE( message,'("Error adding AD extra layers to profile #",i0)' ) m
            CALL display_message( routine_name, message, error_status )
            RETURN
          END IF

        END DO channel_loop


        ! Deallocate local sensor dependent data structures
        ! ...RTV structure
        IF ( rtv_associated(rtv) ) CALL rtv_destroy(rtv)
        ! ...Predictor structures
        CALL crtm_predictor_destroy( predictor )
        CALL crtm_predictor_destroy( predictor_k )

      END DO sensor_loop


      ! Deallocate local sensor independent data structures
      ! ...Atmospheric optics
      CALL crtm_atmoptics_destroy( atmoptics )
      CALL crtm_atmoptics_destroy( atmoptics_k )

    END DO profile_loop


    ! Destroy any remaining structures
    CALL crtm_sfcoptics_destroy( sfcoptics )
    CALL crtm_sfcoptics_destroy( sfcoptics_k )
    CALL crtm_atmosphere_destroy( atm_k )
    CALL crtm_atmosphere_destroy( atm )

  END FUNCTION crtm_k_matrix


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_K_Matrix_Version
!
! PURPOSE:
!       Subroutine to return the module version information.
!
! CALLING SEQUENCE:
!       CALL CRTM_K_Matrix_Version( Id )
!
! OUTPUTS:
!       Id:            Character string containing the version Id information
!                      for the module.
!                      UNITS:      N/A
!                      TYPE:       CHARACTER(*)
!                      DIMENSION:  Scalar
!                      ATTRIBUTES: INTENT(OUT)
!
!:sdoc-:
!--------------------------------------------------------------------------------

  SUBROUTINE crtm_k_matrix_version( Id )
    CHARACTER(*), INTENT(OUT) :: id
    id = module_version_id
  END SUBROUTINE crtm_k_matrix_version


END MODULE crtm_k_matrix_module