RTV_Define.f90

Go to the documentation of this file.

!
! RTV_Define
!
! Module containing the intermediate variables for RTSolution module.
!
!
! NOTE: Modified as generic "bucket" for all RT-related algorithm,
!       ADA, AMOM, and SOI.
!       This is initial step in truly separating out the algorithms
!       into their own modules. Currently each algorithm ties into 
!       the same RTV components (not good.)
!       
!
!
! CREATION HISTORY:
!       Written by:     Quanhua Liu,    QSS at JCSDA;    Quanhua.Liu@noaa.gov 
!                       Yong Han,       NOAA/NESDIS;     Yong.Han@noaa.gov
!                       Paul van Delst, paul.vandelst@noaa.gov
!                       08-Jun-2004

MODULE rtv_define

  ! ------------------
  ! Environment set up
  ! ------------------
  ! Module use statements
  USE type_kinds,            ONLY: fp
  USE message_handler,       ONLY: success, failure, display_message
  USE crtm_parameters,       ONLY: set, zero, one, two, pi, &
                                   max_n_layers, max_n_angles, max_n_legendre_terms, &
                                   degrees_to_radians, &
                                   secant_diffusivity, &
                                   scattering_albedo_threshold, &
                                   optical_depth_threshold, &
                                   rt_ada
  USE sensorinfo_parameters, ONLY: invalid_sensor
  USE crtm_sfcoptics,        ONLY: sovar_type => ivar_type
  ! Disable all implicit typing
  IMPLICIT NONE


  ! --------------------
  ! Default visibilities
  ! --------------------
  ! Everything private by default
  PRIVATE
  ! Parameters
  PUBLIC :: angle_threshold
  PUBLIC :: phase_threshold
  PUBLIC :: delta_optical_depth
  PUBLIC :: max_albedo
  PUBLIC :: small_od_for_sc
  PUBLIC :: max_n_doubling
  PUBLIC :: max_n_soi_iterations
  ! Datatypes
  PUBLIC :: aircraft_rt_type
  PUBLIC :: rtv_type
  ! Procedures
  PUBLIC :: rtv_associated
  PUBLIC :: rtv_destroy
  PUBLIC :: rtv_create
  
  ! -----------------
  ! Module parameters
  ! -----------------
  ! Version Id for the module
  CHARACTER(*),  PARAMETER :: module_rcs_id = &
  '$Id: RTV_Define.f90 60152 2015-08-13 19:19:13Z paul.vandelst@noaa.gov $'

  ! Threshold for determing if an additional stream
  ! angle is required for the satellite zenith angle
  REAL(fp), PARAMETER :: angle_threshold = 1.0e-7_fp

  ! Small positive value used to replace negative
  ! values in the computed phase function
  REAL(fp), PARAMETER :: phase_threshold = 1.0e-7_fp
  
  REAL(fp), PARAMETER :: delta_optical_depth = 1.0e-8_fp
  REAL(fp), PARAMETER :: max_albedo = 0.999999_fp
  
  ! Threshold layer optical depth for single scattering
  REAL(fp), PARAMETER :: small_od_for_sc = 1.e-5_fp

  ! The maximum number of doubling processes in the 
  ! the doubling-adding scheme. 
  INTEGER,  PARAMETER :: max_n_doubling = 55
  
  ! The maximum number of iterations for the SOI solution method.
  INTEGER,  PARAMETER :: max_n_soi_iterations = 75
  
  ! ---------------------
  ! Structure definitions
  ! ---------------------
  ! ...Aircraft model structure
  TYPE :: aircraft_rt_type
    ! The switch
    LOGICAL :: rt = .false.
    ! The output level index
    INTEGER :: idx
  END TYPE aircraft_rt_type  

  ! --------------------------------------
  ! Structure definition to hold forward
  ! variables across FWD, TL, and AD calls
  ! --------------------------------------
  TYPE :: rtv_type
  
    ! Type of sensor
    INTEGER :: sensor_type = invalid_sensor
    
    ! Type of RT algorithm
    INTEGER :: rt_algorithm_id = rt_ada
    
    ! Dimension information
    INTEGER :: n_layers         = 0       ! Total number of atmospheric layers
    INTEGER :: n_added_layers   = 0       ! Number of layers appended to TOA
    INTEGER :: n_angles         = 0       ! Number of angles to be considered
    INTEGER :: n_soi_iterations = 0       ! Number of SOI iterations
    

    REAL(fp):: cos_sun = zero           ! Cosine of sun zenith angle
    REAL(fp):: solar_irradiance = zero  ! channel solar iiradiance at TOA
    REAL(fp):: cosmic_background_radiance = zero ! For background temp=2.7253 
            
    ! Variable to hold the various portions of the
    ! radiance for emissivity retrieval algorithms
    ! Passed to FWD RTSolution structure argument output
    REAL(fp) :: up_radiance         = zero
    REAL(fp) :: down_solar_radiance = zero

    REAL(fp) :: secant_down_angle = 0

    ! Overcast radiance for cloud detection
    REAL(fp), DIMENSION( 0:MAX_N_LAYERS ) :: e_cloud_radiance_up= zero
    REAL(fp), DIMENSION( 0:MAX_N_LAYERS ) :: e_source_up        = zero
    REAL(fp), DIMENSION( 0:MAX_N_LAYERS ) :: e_level_trans_up   = one

    ! Emission model variables
    REAL(fp) :: total_od  = zero
    REAL(fp), DIMENSION(   MAX_N_LAYERS ) :: e_layer_trans_up   = zero
    REAL(fp), DIMENSION(   MAX_N_LAYERS ) :: e_layer_trans_down = zero
    REAL(fp), DIMENSION( 0:MAX_N_LAYERS ) :: e_level_rad_up     = zero
    REAL(fp), DIMENSION( 0:MAX_N_LAYERS ) :: e_level_rad_down   = zero

    ! Planck radiances
    REAL(fp)                               :: planck_surface    = zero
    REAL(fp), DIMENSION(  0:MAX_N_LAYERS ) :: planck_atmosphere = zero

    ! Quadrature information
    REAL(fp), DIMENSION( MAX_N_ANGLES ) :: cos_angle  = zero  ! Gaussian quadrature abscissa
    REAL(fp), DIMENSION( MAX_N_ANGLES ) :: cos_weight = zero  ! Gaussian quadrature weights

    ! Logical switches
    LOGICAL :: diffuse_surface = .true.
    LOGICAL :: is_solar_channel = .false.
    
    ! Aircraft model RT information
    TYPE(aircraft_rt_type) :: aircraft

    ! Scattering, visible model variables    
    INTEGER :: n_streams         = 0       ! Number of *hemispheric* stream angles used in RT    
    INTEGER :: mth_azi                     ! mth fourier component
    INTEGER :: n_azi                       ! number of fourier components
    LOGICAL :: solar_flag_true   = .false.
    LOGICAL :: visible_flag_true = .false. 
    LOGICAL :: scattering_rt     = .false.

    !-----------------------------------
    ! Variables used in the ADA routines
    !-----------------------------------
    ! Flag to indicate the following arrays have all been allocated
    LOGICAL :: is_allocated = .false.
     
    ! Phase function variables
    ! Forward and backward scattering phase matrices
    REAL(fp), ALLOCATABLE :: pff(:,:,:)  ! MAX_N_ANGLES, MAX_N_ANGLES+1, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: pbb(:,:,:)  ! MAX_N_ANGLES, MAX_N_ANGLES+1, MAX_N_LAYERS

    ! Positive and negative cosine angle Legendre phase functions
    REAL(fp), ALLOCATABLE :: pplus(:,:)  ! 0:MAX_N_LEGENDRE_TERMS, MAX_N_ANGLES
    REAL(fp), ALLOCATABLE :: pminus(:,:) ! 0:MAX_N_LEGENDRE_TERMS, MAX_N_ANGLES
    REAL(fp), ALLOCATABLE :: pleg(:,:)   ! 0:MAX_N_LEGENDRE_TERMS, MAX_N_ANGLES+1
    
    ! Original forward and backward scattering phase matrices.
    ! These may be slightly negative and, if so, need to be made
    ! positive and thus adjusted to ensure energy conservation
    REAL(fp), ALLOCATABLE :: off(:,:,:)  ! MAX_N_ANGLES, MAX_N_ANGLES+1, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: obb(:,:,:)  ! MAX_N_ANGLES, MAX_N_ANGLES+1, MAX_N_LAYERS
    
    ! Normalisation factor and intermediate sum used for original
    ! phase matrix energy conservation.
    REAL(fp), ALLOCATABLE :: n_factor(:,:)  ! MAX_N_ANGLES, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: sum_fac(:,:)   ! 0:MAX_N_ANGLES, MAX_N_LAYERS

    ! Adding-Doubling model variables
    REAL(fp), ALLOCATABLE :: inv_gamma(:,:,:)         ! MAX_N_ANGLES, MAX_N_ANGLES, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: inv_gammat(:,:,:)        ! MAX_N_ANGLES, MAX_N_ANGLES, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: refl_trans(:,:,:)        ! MAX_N_ANGLES, MAX_N_ANGLES, MAX_N_LAYERS

    REAL(fp), ALLOCATABLE :: s_layer_trans(:,:,:)     ! MAX_N_ANGLES, MAX_N_ANGLES, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: s_layer_refl(:,:,:)      ! MAX_N_ANGLES, MAX_N_ANGLES, MAX_N_LAYERS

    REAL(fp), ALLOCATABLE :: s_level_refl_up(:,:,:)   ! MAX_N_ANGLES, MAX_N_ANGLES, 0:MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: s_level_rad_up(:,:)      ! MAX_N_ANGLES, 0:MAX_N_LAYERS
     
    REAL(fp), ALLOCATABLE :: s_layer_source_up(:,:)   ! MAX_N_ANGLES, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: s_layer_source_down(:,:) ! MAX_N_ANGLES, MAX_N_LAYERS


    !------------------------------------
    ! Variables used in the AMOM routines
    !------------------------------------
    ! dimensions, MAX_N_ANGLES, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: thermal_c(:,:)
    REAL(fp), ALLOCATABLE :: eigva(:,:)
    REAL(fp), ALLOCATABLE :: exp_x(:,:)
    REAL(fp), ALLOCATABLE :: eigvalue(:,:)
    
    ! dimensions, MAX_N_ANGLES, MAX_N_ANGLES, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: hh(:,:,:)
    REAL(fp), ALLOCATABLE :: pm(:,:,:)
    REAL(fp), ALLOCATABLE :: pp(:,:,:)
    REAL(fp), ALLOCATABLE :: ppm(:,:,:)
    REAL(fp), ALLOCATABLE :: ppp(:,:,:)
    REAL(fp), ALLOCATABLE :: i_ppm(:,:,:)
    REAL(fp), ALLOCATABLE :: i_ppp(:,:,:)
    REAL(fp), ALLOCATABLE :: eigve(:,:,:)
    REAL(fp), ALLOCATABLE :: gm(:,:,:)
    REAL(fp), ALLOCATABLE :: i_gm(:,:,:)
    REAL(fp), ALLOCATABLE :: gp(:,:,:)
    REAL(fp), ALLOCATABLE :: eigvef(:,:,:)
    REAL(fp), ALLOCATABLE :: eigveva(:,:,:)
    REAL(fp), ALLOCATABLE :: a1(:,:,:)
    REAL(fp), ALLOCATABLE :: a2(:,:,:)
    REAL(fp), ALLOCATABLE :: a3(:,:,:)
    REAL(fp), ALLOCATABLE :: a4(:,:,:)
    REAL(fp), ALLOCATABLE :: a5(:,:,:)
    REAL(fp), ALLOCATABLE :: a6(:,:,:)
    REAL(fp), ALLOCATABLE :: gm_a5(:,:,:)
    REAL(fp), ALLOCATABLE :: i_gm_a5(:,:,:)


    !-----------------------------------
    ! Variables used in the SOI routines
    !-----------------------------------
    INTEGER :: number_soi_iter = 0
    REAL(fp), ALLOCATABLE :: e_layer_trans(:,:)           ! MAX_N_ANGLES, MAX_N_LAYERS
    REAL(fp), ALLOCATABLE :: s_level_iterrad_down(:,:,:)  ! MAX_N_ANGLES, 0:MAX_N_LAYERS, MAX_N_SOI_ITERATIONS
    REAL(fp), ALLOCATABLE :: s_level_iterrad_up(:,:,:)    ! MAX_N_ANGLES, 0:MAX_N_LAYERS, MAX_N_SOI_ITERATIONS

    INTEGER , ALLOCATABLE :: number_doubling(:)  ! n_Layers
    REAL(fp), ALLOCATABLE :: delta_tau(:)        ! n_Layers
    REAL(fp), ALLOCATABLE :: refl(:,:,:,:)       ! n_Angles, n_Angles, 0:MAX_N_DOUBLING, n_Layers
    REAL(fp), ALLOCATABLE :: trans(:,:,:,:)      ! n_Angles, n_Angles, 0:MAX_N_DOUBLING, n_Layers
    REAL(fp), ALLOCATABLE :: inv_bet(:,:,:,:)    ! n_Angles, n_Angles, 0:MAX_N_DOUBLING, n_Layers
    REAL(fp), ALLOCATABLE :: c1(:,:)             ! n_Angles, n_Layers
    REAL(fp), ALLOCATABLE :: c2(:,:)             ! n_Angles, n_Layers


    ! The surface optics forward variables
    TYPE(sovar_type) :: sov

  END TYPE rtv_type


CONTAINS



!################################################################################
!################################################################################
!##                                                                            ##
!##                        ## PRIVATE MODULE ROUTINES ##                       ##
!##                                                                            ##
!################################################################################
!################################################################################

!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       RTV_Associated
!
! PURPOSE:
!       Elemental function to test if the allocatable components of an
!       RTV object have been allocated.
!
! CALLING SEQUENCE:
!       Status = RTV_Associated( RTV )
!
! OBJECTS:
!       RTV:     RTV structure which is to have its member's
!                status tested.
!                UNITS:      N/A
!                TYPE:       RTV_type
!                DIMENSION:  Scalar or any rank
!                ATTRIBUTES: INTENT(IN)
!
! FUNCTION RESULT:
!       Status:  The return value is a logical value indicating the
!                status of the RTV members.
!                .TRUE.  - if ANY of the RTV allocatable or
!                          pointer members are in use.
!                .FALSE. - if ALL of the RTV allocatable or
!                          pointer members are not in use.
!                UNITS:      N/A
!                TYPE:       LOGICAL
!                DIMENSION:  Same as input RTV argument
!
!:sdoc-:
!--------------------------------------------------------------------------------

  ELEMENTAL FUNCTION rtv_associated( RTV ) RESULT( Status )
    ! Arguments
    TYPE(rtv_type), INTENT(IN) :: rtv
    ! Function result
    LOGICAL :: status

    ! Test the structure members
    status = rtv%Is_Allocated

  END FUNCTION rtv_associated


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       RTV_Destroy
! 
! PURPOSE:
!       Elemental subroutine to re-initialize RTV objects.
!
! CALLING SEQUENCE:
!       CALL RTV_Destroy( RTV )
!
! OBJECTS:
!       RTV:          Re-initialized RTV structure.
!                     UNITS:      N/A
!                     TYPE:       RTV_type
!                     DIMENSION:  Scalar OR any rank
!                     ATTRIBUTES: INTENT(OUT)
!
!:sdoc-:
!--------------------------------------------------------------------------------

  ELEMENTAL SUBROUTINE rtv_destroy( RTV )
    TYPE(rtv_type), INTENT(OUT) :: rtv

    ! Belts and braces
    rtv%Is_Allocated = .false.
    
  END SUBROUTINE rtv_destroy

  
!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       RTV_Create
! 
! PURPOSE:
!       Elemental subroutine to create an instance of the RTV object.
!
! CALLING SEQUENCE:
!       CALL RTV_Create( RTV, &
!                        n_Angles        , &
!                        n_Legendre_Terms, &
!                        n_Layers          )
!
! OBJECTS:
!       RTV:               RTV structure.
!                          UNITS:      N/A
!                          TYPE:       RTV_type
!                          DIMENSION:  Scalar or any rank
!                          ATTRIBUTES: INTENT(OUT)
!
! INPUTS:
!       n_Angles:          Number of 
!                          Must be > 0.
!                          UNITS:      N/A
!                          TYPE:       INTEGER
!                          DIMENSION:  Same as RTV object
!                          ATTRIBUTES: INTENT(IN)
!
!       n_Legendre_Terms:  Number of 
!                          Must be > 0.
!                          UNITS:      N/A
!                          TYPE:       INTEGER
!                          DIMENSION:  Same as RTV object
!                          ATTRIBUTES: INTENT(IN)
!
!       n_Layers:          Number of atmospheric layers.
!                          Must be > 0.
!                          UNITS:      N/A
!                          TYPE:       INTEGER
!                          DIMENSION:  Same as RTV object
!                          ATTRIBUTES: INTENT(IN)
!
!:sdoc-:
!--------------------------------------------------------------------------------

  ELEMENTAL SUBROUTINE rtv_create( &
    RTV, &
    n_Angles        , &
    n_Legendre_Terms, &
    n_Layers          )
    ! Arguments
    TYPE(rtv_type), INTENT(OUT) :: rtv
    INTEGER       , INTENT(IN)  :: n_angles        
    INTEGER       , INTENT(IN)  :: n_legendre_terms
    INTEGER       , INTENT(IN)  :: n_layers        
    ! Local variables
    INTEGER :: alloc_stat

    ! Check input
    IF ( n_angles < 1 .OR. n_legendre_terms < 1 .OR. n_layers < 1 ) RETURN
    
    ! Perform the allocation for phase function variables
    ALLOCATE( rtv%Pff(n_angles, n_angles+1, n_layers) , &
              rtv%Pbb(n_angles, n_angles+1, n_layers) , &
              rtv%Pplus( 0:n_legendre_terms, n_angles), &
              rtv%Pminus(0:n_legendre_terms, n_angles), &
              rtv%Pleg(0:n_legendre_terms, n_angles+1), &
              rtv%Off(n_angles, n_angles+1, n_layers) , &
              rtv%Obb(n_angles, n_angles+1, n_layers) , &
              rtv%n_Factor (n_angles, n_layers)       , &
              rtv%sum_fac(0:n_angles, n_layers)       , &
              stat = alloc_stat )
    IF ( alloc_stat /= 0 ) RETURN

    ! Perform the allocation for adding-doubling variables
    ALLOCATE( rtv%Inv_Gamma( n_angles, n_angles, n_layers)       , &
              rtv%Inv_GammaT(n_angles, n_angles, n_layers)       , &
              rtv%Refl_Trans(n_angles, n_angles, n_layers)       , &
              rtv%s_Layer_Trans(n_angles, n_angles, n_layers)    , &
              rtv%s_Layer_Refl( n_angles, n_angles, n_layers)    , &
              rtv%s_Level_Refl_UP(n_angles, n_angles, 0:n_layers), &
              rtv%s_Level_Rad_UP(n_angles, 0:n_layers)           , &
              rtv%s_Layer_Source_UP(  n_angles, n_layers)        , &
              rtv%s_Layer_Source_DOWN(n_angles, n_layers)        , &
              stat = alloc_stat )
    IF ( alloc_stat /= 0 ) RETURN

    ! Perform the allocation for AMOM variables
    ALLOCATE( rtv%Thermal_C(n_angles, n_layers)        , &
              rtv%EigVa(n_angles, n_layers)            , &
              rtv%Exp_x(n_angles, n_layers)            , &
              rtv%EigValue(n_angles, n_layers)         , &
              rtv%HH(n_angles, n_angles, n_layers)     , &
              rtv%PM(n_angles, n_angles, n_layers)     , &
              rtv%PP(n_angles, n_angles, n_layers)     , &
              rtv%PPM(n_angles, n_angles, n_layers)    , &
              rtv%PPP(n_angles, n_angles, n_layers)    , &
              rtv%i_PPM(n_angles, n_angles, n_layers)  , &
              rtv%i_PPP(n_angles,n_angles,n_layers)    , &
              rtv%EigVe(n_angles, n_angles, n_layers)  , &
              rtv%Gm(n_angles, n_angles, n_layers)     , &
              rtv%i_Gm(n_angles, n_angles, n_layers)   , &
              rtv%Gp(n_angles, n_angles, n_layers)     , &
              rtv%EigVeF(n_angles, n_angles, n_layers) , &
              rtv%EigVeVa(n_angles,n_angles,n_layers)  , &
              rtv%A1(n_angles,n_angles, n_layers)      , &
              rtv%A2(n_angles, n_angles, n_layers)     , &
              rtv%A3(n_angles, n_angles, n_layers)     , &
              rtv%A4(n_angles, n_angles, n_layers)     , &
              rtv%A5(n_angles, n_angles, n_layers)     , &
              rtv%A6(n_angles, n_angles, n_layers)     , &
              rtv%Gm_A5(n_angles, n_angles, n_layers)  , &
              rtv%i_Gm_A5(n_angles, n_angles, n_layers), & 
              stat = alloc_stat )
    IF ( alloc_stat /= 0 ) RETURN

    ! Perform the allocation for SOI variables
    ALLOCATE( rtv%e_Layer_Trans( n_angles, n_layers), &
              rtv%s_Level_IterRad_DOWN( n_angles, 0:n_layers, max_n_soi_iterations ), &
              rtv%s_Level_IterRad_UP( n_angles, 0:n_layers, max_n_soi_iterations ), &  
              rtv%Number_Doubling(n_layers), &
              rtv%Delta_Tau(n_layers), &      
              rtv%Refl(n_angles, n_angles, 0:max_n_doubling, n_layers), &
              rtv%Trans(n_angles, n_angles, 0:max_n_doubling, n_layers), &
              rtv%Inv_BeT(n_angles, n_angles, 0:max_n_doubling, n_layers), &
              rtv%C1(n_angles, n_layers), &
              rtv%C2(n_angles, n_layers), &
              stat = alloc_stat )
    IF ( alloc_stat /= 0 ) RETURN

    ! Set dimensions
    rtv%n_Layers         = n_layers
    rtv%n_Angles         = n_angles
    
    rtv%n_SOI_Iterations = 0

    ! Set the allocate flag
    rtv%Is_Allocated = .true.
    
  END SUBROUTINE rtv_create
  
END MODULE rtv_define