CRTM_AtmOptics_Define.f90

Go to the documentation of this file.

!
! CRTM_AtmOptics_Define
!
! Module defining the CRTM AtmOptics object.
!
!
! CREATION HISTORY:
!       Written by:     Paul van Delst, 11-Oct-2011
!                       paul.vandelst@noaa.gov
!

MODULE crtm_atmoptics_define

  ! -----------------
  ! Environment setup
  ! -----------------
  ! Module use
  USE type_kinds           , ONLY: fp, long, double
  USE message_handler      , ONLY: success, failure, information, display_message
  USE compare_float_numbers, ONLY: default_n_sigfig, &
                                   OPERATOR(.equalto.), &
                                   compares_within_tolerance
  USE file_utility         , ONLY: file_open, file_exists
  USE binary_file_utility  , ONLY: open_binary_file      , &
                                   writegatts_binary_file, &
                                   readgatts_binary_file
  ! Disable implicit typing
  IMPLICIT NONE


  ! ------------
  ! Visibilities
  ! ------------
  ! Everything private by default
  PRIVATE
  ! Datatypes
  PUBLIC :: crtm_atmoptics_type
  ! Operators
  PUBLIC :: OPERATOR(==)
  PUBLIC :: OPERATOR(-)
  ! Procedures
  PUBLIC :: crtm_atmoptics_associated
  PUBLIC :: crtm_atmoptics_destroy
  PUBLIC :: crtm_atmoptics_create
  PUBLIC :: crtm_atmoptics_zero
  PUBLIC :: crtm_atmoptics_inspect
  PUBLIC :: crtm_atmoptics_validrelease
  PUBLIC :: crtm_atmoptics_info
  PUBLIC :: crtm_atmoptics_defineversion
  PUBLIC :: crtm_atmoptics_compare
  PUBLIC :: crtm_atmoptics_inquirefile
  PUBLIC :: crtm_atmoptics_readfile
  PUBLIC :: crtm_atmoptics_writefile


  ! ---------------------
  ! Procedure overloading
  ! ---------------------
  INTERFACE OPERATOR(==)
    MODULE PROCEDURE crtm_atmoptics_equal
  END INTERFACE OPERATOR(==)

  INTERFACE OPERATOR(-)
    MODULE PROCEDURE crtm_atmoptics_subtract
  END INTERFACE OPERATOR(-)

  INTERFACE crtm_atmoptics_inspect
    MODULE PROCEDURE scalar_inspect
    MODULE PROCEDURE rank1_inspect
  END INTERFACE crtm_atmoptics_inspect

  ! -----------------
  ! Module parameters
  ! -----------------
  CHARACTER(*), PARAMETER :: module_version_id = &
    '$Id: CRTM_AtmOptics_Define.f90 60152 2015-08-13 19:19:13Z paul.vandelst@noaa.gov $'
  ! Release and version
  INTEGER, PARAMETER :: atmoptics_release = 4  ! This determines structure and file formats.
  ! Close status for write errors
  CHARACTER(*), PARAMETER :: write_error_status = 'DELETE'
  ! Literal constants
  REAL(fp), PARAMETER :: zero = 0.0_fp
  REAL(fp), PARAMETER :: one  = 1.0_fp
  ! String lengths
  INTEGER,  PARAMETER :: ml = 256 ! Message length
  INTEGER,  PARAMETER :: sl =  80 ! String length


  ! ------------------------------
  ! AtmOptics data type definition
  ! ------------------------------
  !:tdoc+:
  TYPE :: crtm_atmoptics_type
    ! Allocation indicator
    LOGICAL :: is_allocated = .false.
    ! Release information
    INTEGER :: release = atmoptics_release
    ! Dimensions
    INTEGER :: n_layers           = 0  ! K dimension
    INTEGER :: n_legendre_terms   = 0  ! Ic dimension
    INTEGER :: n_phase_elements   = 0  ! Ip dimension
    ! ...Max. values to avoid unnecessary reallocation
    INTEGER :: max_layers         = 0  ! K-Max dimension
    INTEGER :: max_legendre_terms = 0  ! Ic-Max dimension
    INTEGER :: max_phase_elements = 0  ! Ip-Max dimension
    ! Scalar components
    LOGICAL  :: include_scattering = .true.
    INTEGER  :: loffset = 0   ! Start position in array for Legendre coefficients
    REAL(fp) :: scattering_optical_depth = zero
    ! Array components
    REAL(fp), ALLOCATABLE :: optical_depth(:)         ! K-Max
    REAL(fp), ALLOCATABLE :: single_scatter_albedo(:) ! K-Max
    REAL(fp), ALLOCATABLE :: asymmetry_factor(:)      ! K-Max
    REAL(fp), ALLOCATABLE :: delta_truncation(:)      ! K-Max
    REAL(fp), ALLOCATABLE :: phase_coefficient(:,:,:) ! 0:Ic-Max x Ip-Max x K-Max
  END TYPE crtm_atmoptics_type
  !:tdoc-:


CONTAINS


!################################################################################
!################################################################################
!##                                                                            ##
!##                         ## PUBLIC MODULE ROUTINES ##                       ##
!##                                                                            ##
!################################################################################
!################################################################################

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

  ELEMENTAL FUNCTION crtm_atmoptics_associated( self ) RESULT( status )
    TYPE(crtm_atmoptics_type), INTENT(IN) :: self
    LOGICAL :: status
    status = self%Is_Allocated
  END FUNCTION crtm_atmoptics_associated


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

  ELEMENTAL SUBROUTINE crtm_atmoptics_destroy( self )
    TYPE(crtm_atmoptics_type), INTENT(OUT) :: self
    self%Is_Allocated = .false.
    self%n_Layers           = 0
    self%n_Legendre_Terms   = 0
    self%n_Phase_Elements   = 0
    self%Max_Layers         = 0
    self%Max_Legendre_Terms = 0
    self%Max_Phase_Elements = 0
  END SUBROUTINE crtm_atmoptics_destroy


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_AtmOptics_Create
!
! PURPOSE:
!       Elemental subroutine to create an instance of an AtmOptics object.
!
! CALLING SEQUENCE:
!       CALL CRTM_AtmOptics_Create( AtmOptics       , &
!                                   n_Layers        , &
!                                   n_Legendre_Terms, &
!                                   n_Phase_Elements  )
!
! OBJECTS:
!       AtmOptics:          AtmOptics object structure.
!                           UNITS:      N/A
!                           TYPE:       CRTM_AtmOptics_type
!                           DIMENSION:  Scalar or any rank
!                           ATTRIBUTES: INTENT(IN OUT)
!
! INPUTS:
!       n_Layers:           Number of atmospheric layers.
!                           Must be > 0
!                           UNITS:      N/A
!                           TYPE:       INTEGER
!                           DIMENSION:  Conformable with the AtmOptics object.
!                           ATTRIBUTES: INTENT(IN)
!
!       n_Legendre_Terms:   The number of Legendre polynomial terms for the
!                           phase matrix.
!                           Must be > 0
!                           UNITS:      N/A
!                           TYPE:       INTEGER
!                           DIMENSION:  Same as n_Layers input.
!                           ATTRIBUTES: INTENT(IN)
!
!       n_Phase_Elements:   The number of phase elements for the phase matrix.
!                           Must be > 0
!                           UNITS:      N/A
!                           TYPE:       INTEGER
!                           DIMENSION:  Same as n_Layers input.
!                           ATTRIBUTES: INTENT(IN)
!
! COMMENTS:
!       Note the INTENT on the output AtmOptics argument is IN OUT rather than
!       just OUT. If the AtmOptics is already allocated to sufficient size for
!       the passed dimensions, no reallocation is performed, only the internal
!       dimension values are reset.
!
!:sdoc-:
!--------------------------------------------------------------------------------

  ELEMENTAL SUBROUTINE crtm_atmoptics_create( &
    self            , &  ! Output
    n_Layers        , &  ! Input
    n_Legendre_Terms, &  ! Input
    n_Phase_Elements  )  ! Input
    ! Arguments
    TYPE(crtm_atmoptics_type), INTENT(IN OUT) :: self
    INTEGER                  , INTENT(IN)     :: n_layers
    INTEGER                  , INTENT(IN)     :: n_legendre_terms
    INTEGER                  , INTENT(IN)     :: n_phase_elements
    ! Local variables
    INTEGER :: alloc_stat

    ! Check input
    IF ( n_layers         < 1 .OR. &
         n_legendre_terms < 1 .OR. &
         n_phase_elements < 1 ) THEN
      CALL crtm_atmoptics_destroy( self )
      RETURN
    END IF

    ! Check the allocation status
    IF ( crtm_atmoptics_associated( self ) ) THEN
      ! ...Does it need to be reallocated?
      IF ( self%Max_Layers         < n_layers         .OR. &
           self%Max_Legendre_Terms < n_legendre_terms .OR. &
           self%Max_Phase_Elements < n_phase_elements ) THEN
        CALL atmoptics_allocate( self, alloc_stat )
        IF ( alloc_stat /= 0 ) RETURN
      END IF
    ELSE
      ! ...Not allocated, so allocate it.
      CALL atmoptics_allocate( self, alloc_stat )
      IF ( alloc_stat /= 0 ) RETURN
    END IF

    ! Initialise dimensions (but not arrays!)
    self%n_Layers         = n_layers
    self%n_Legendre_Terms = n_legendre_terms
    self%n_Phase_Elements = n_phase_elements

    ! Set allocation indicator
    self%Is_Allocated = .true.

  CONTAINS

    PURE SUBROUTINE atmoptics_allocate(self,alloc_stat)
      TYPE(crtm_atmoptics_type), INTENT(OUT) :: self
      INTEGER                  , INTENT(OUT) :: alloc_stat
      ! Allocate object
      ALLOCATE( self%Optical_Depth( n_layers ), &
                self%Single_Scatter_Albedo( n_layers ), &
                self%Asymmetry_Factor( n_layers ), &
                self%Delta_Truncation( n_layers ), &
                self%Phase_Coefficient( 0:n_legendre_terms, n_phase_elements, n_layers ), &
                stat = alloc_stat )
      IF ( alloc_stat /= 0 ) RETURN
      ! Set maximum dimension values
      self%Max_Layers         = n_layers
      self%Max_Legendre_Terms = n_legendre_terms
      self%Max_Phase_Elements = n_phase_elements
    END SUBROUTINE atmoptics_allocate

  END SUBROUTINE crtm_atmoptics_create


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_AtmOptics_Zero
!
! PURPOSE:
!       Elemental subroutine to initialise the components of an AtmOptics
!       object to a value of zero.
!
! CALLING SEQUENCE:
!       CALL CRTM_AtmOptics_Zero( AtmOptics )
!
! OBJECTS:
!       AtmOptics:   AtmOptics object which is to have its components
!                    set to a zero value.
!                    UNITS:      N/A
!                    TYPE:       CRTM_AtmOptics_type
!                    DIMENSION:  Scalar or any rank
!                    ATTRIBUTES: INTENT(IN OUT)
!
!:sdoc-:
!--------------------------------------------------------------------------------

  ELEMENTAL SUBROUTINE crtm_atmoptics_zero( self )
    TYPE(crtm_atmoptics_type), INTENT(IN OUT) :: self
    self%Scattering_Optical_Depth = zero
    IF ( .NOT. crtm_atmoptics_associated( self ) ) RETURN
    self%Optical_Depth         = zero
    self%Single_Scatter_Albedo = zero
    self%Asymmetry_Factor      = zero
    self%Delta_Truncation      = zero
    self%Phase_Coefficient     = zero
  END SUBROUTINE crtm_atmoptics_zero


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_AtmOptics_Inspect
!
! PURPOSE:
!       Subroutine to print the contents of a AtmOptics object to stdout.
!
! CALLING SEQUENCE:
!       CALL CRTM_AtmOptics_Inspect( AtmOptics )
!
! OBJECTS:
!       AtmOptics:     AtmOptics object to display.
!                      UNITS:      N/A
!                      TYPE:       CRTM_AtmOptics_type
!                      DIMENSION:  Scalar
!                      ATTRIBUTES: INTENT(IN)
!
!:sdoc-:
!--------------------------------------------------------------------------------

  SUBROUTINE scalar_inspect(self)
    TYPE(CRTM_AtmOptics_type), INTENT(IN) :: self
    INTEGER :: ip, k
    WRITE(*,'(1x,"AtmOptics OBJECT")')
    ! Release/version info
    WRITE(*,'(3x,"Release          : ",i0)') self%Release
    ! Dimensions
    WRITE(*,'(3x,"n_Layers         : ",i0," (of max. ",i0,")")') self%n_Layers        , self%Max_Layers
    WRITE(*,'(3x,"n_Legendre_Terms : ",i0," (of max. ",i0,")")') self%n_Legendre_Terms, self%Max_Legendre_Terms
    WRITE(*,'(3x,"n_Phase_Elements : ",i0," (of max. ",i0,")")') self%n_Phase_Elements, self%Max_Phase_Elements
    WRITE(*,'(3x,"Scattering Optical Depth : ",es13.6)') self%Scattering_Optical_Depth
    IF ( .NOT. crtm_atmoptics_associated(self) ) RETURN
    ! Dimension arrays
    WRITE(*,'(3x,"Optical_Depth :")')
    WRITE(*,'(5(1x,es13.6,:))') self%Optical_Depth(1:self%n_Layers)
    WRITE(*,'(3x,"Single_Scatter_Albedo :")')
    WRITE(*,'(5(1x,es13.6,:))') self%Single_Scatter_Albedo(1:self%n_Layers)
    WRITE(*,'(3x,"Asymmetry_Factor :")')
    WRITE(*,'(5(1x,es13.6,:))') self%Asymmetry_Factor(1:self%n_Layers)
    WRITE(*,'(3x,"Delta_Truncation :")')
    WRITE(*,'(5(1x,es13.6,:))') self%Delta_Truncation(1:self%n_Layers)
    WRITE(*,'(3x,"Phase_Coefficient Legendre polynomial coefficients :")')
    DO k = 1, self%n_Layers
      DO ip = 1, self%n_Phase_Elements
        WRITE(*,'(5x,"Layer: ",i0,"; Phase element: ",i0)') k, ip
        WRITE(*,'(5(1x,es13.6,:))') self%Phase_Coefficient(0:self%n_Legendre_Terms,ip,k)
      END DO
      WRITE(*,*)
    END DO
  END SUBROUTINE scalar_inspect

  SUBROUTINE rank1_inspect( self )
    TYPE(CRTM_AtmOptics_type), INTENT(IN) :: self(:)
    INTEGER :: n, n_objects

    n_objects = SIZE(self)
    DO n = 1, n_objects
      WRITE(*, fmt='(1x,"OBJECT INDEX:",i0," - ")', advance='NO') n
      CALL scalar_inspect(self(n))
    END DO
  END SUBROUTINE rank1_inspect


!----------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_AtmOptics_ValidRelease
!
! PURPOSE:
!       Function to check the AtmOptics object Release value.
!
! CALLING SEQUENCE:
!       IsValid = CRTM_AtmOptics_ValidRelease( AtmOptics )
!
! OBJECTS:
!       AtmOptics:     AtmOptics object for which the Release component
!                      is to be checked.
!                      UNITS:      N/A
!                      TYPE:       CRTM_AtmOptics_type
!                      DIMENSION:  Scalar
!                      ATTRIBUTES: INTENT(IN)
!
! FUNCTION RESULT:
!       IsValid:       Logical value defining the release validity.
!                      UNITS:      N/A
!                      TYPE:       LOGICAL
!                      DIMENSION:  Scalar
!
!:sdoc-:
!----------------------------------------------------------------------------------

  FUNCTION crtm_atmoptics_validrelease( self ) RESULT( IsValid )
    ! Arguments
    TYPE(crtm_atmoptics_type), INTENT(IN) :: self
    ! Function result
    LOGICAL :: isvalid
    ! Local parameters
    CHARACTER(*), PARAMETER :: routine_name = 'CRTM_AtmOptics_ValidRelease'
    ! Local variables
    CHARACTER(ML) :: msg

    ! Set up
    isvalid = .true.


    ! Check release is not too old
    IF ( self%Release < atmoptics_release ) THEN
      isvalid = .false.
      WRITE( msg,'("An AtmOptics data update is needed. ", &
                  &"AtmOptics release is ",i0,". Valid release is ",i0,"." )' ) &
                  self%Release, atmoptics_release
      CALL display_message( routine_name, msg, information ); RETURN
    END IF


    ! Check release is not too new
    IF ( self%Release > atmoptics_release ) THEN
      isvalid = .false.
      WRITE( msg,'("An AtmOptics software update is needed. ", &
                  &"AtmOptics release is ",i0,". Valid release is ",i0,"." )' ) &
                  self%Release, atmoptics_release
      CALL display_message( routine_name, msg, information ); RETURN
    END IF

  END FUNCTION crtm_atmoptics_validrelease


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_AtmOptics_Info
!
! PURPOSE:
!       Subroutine to return a string containing version and dimension
!       information about an AtmOptics object.
!
! CALLING SEQUENCE:
!       CALL CRTM_AtmOptics_Info( AtmOptics, Info )
!
! OBJECTS:
!       AtmOptics:     AtmOptics object about which info is required.
!                      UNITS:      N/A
!                      TYPE:       CRTM_AtmOptics_type
!                      DIMENSION:  Scalar
!                      ATTRIBUTES: INTENT(IN)
!
! OUTPUTS:
!       Info:          String containing version and dimension information
!                      about the AtmOptics object.
!                      UNITS:      N/A
!                      TYPE:       CHARACTER(*)
!                      DIMENSION:  Scalar
!                      ATTRIBUTES: INTENT(OUT)
!
!:sdoc-:
!--------------------------------------------------------------------------------

  SUBROUTINE crtm_atmoptics_info( self, Info )
    ! Arguments
    TYPE(crtm_atmoptics_type), INTENT(IN)  :: self
    CHARACTER(*),              INTENT(OUT) :: info
    ! Parameters
    INTEGER, PARAMETER :: carriage_return = 13
    INTEGER, PARAMETER :: linefeed = 10
    ! Local variables
    CHARACTER(2000) :: long_string

    ! Write the required data to the local string
    WRITE( long_string, &
           '(a,1x,"AtmOptics RELEASE: ",i2,3x, &
           &"N_LAYERS=",i0,2x,&
           &"N_LEGENDRE_TERMS=",i0,2x,&
           &"N_PHASE_ELEMENTS=",i0 )' ) &
           achar(carriage_return)//achar(linefeed), &
           self%Release, &
           self%n_Layers, &
           self%n_Legendre_Terms, &
           self%n_Phase_Elements

    ! Trim the output based on the
    ! dummy argument string length
    info = long_string(1:min(len(info), len_trim(long_string)))

  END SUBROUTINE crtm_atmoptics_info


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_AtmOptics_DefineVersion
!
! PURPOSE:
!       Subroutine to return the module version information.
!
! CALLING SEQUENCE:
!       CALL CRTM_AtmOptics_DefineVersion( 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_atmoptics_defineversion( Id )
    CHARACTER(*), INTENT(OUT) :: id
    id = module_version_id
  END SUBROUTINE crtm_atmoptics_defineversion


!--------------------------------------------------------------------------------
!:sdoc+:
! NAME:
!       CRTM_AtmOptics_Compare
!
! PURPOSE:
!       Elemental function to compare two CRTM_AtmOptics objects to within
!       a user specified number of significant figures.
!
! CALLING SEQUENCE:
!       is_comparable = CRTM_AtmOptics_Compare( x, y, n_SigFig=n_SigFig )
!
! OBJECTS:
!       x, y:          Two CRTM AtmOptics objects to be compared.
!                      UNITS:      N/A
!                      TYPE:       CRTM_AtmOptics_type
!                      DIMENSION:  Scalar or any rank
!                      ATTRIBUTES: INTENT(IN)
!
! OPTIONAL INPUTS:
!       n_SigFig:      Number of significant figures to compare floating point
!                      components.
!                      UNITS:      N/A
!                      TYPE:       INTEGER
!                      DIMENSION:  Scalar or same as input
!                      ATTRIBUTES: INTENT(IN), OPTIONAL
!
! FUNCTION RESULT:
!       is_equal:      Logical value indicating whether the inputs are equal.
!                      UNITS:      N/A
!                      TYPE:       LOGICAL
!                      DIMENSION:  Same as inputs.
!:sdoc-:
!--------------------------------------------------------------------------------

  ELEMENTAL FUNCTION crtm_atmoptics_compare( &
    x, &
    y, &
    n_SigFig ) &
  result( is_comparable )
    TYPE(crtm_atmoptics_type), INTENT(IN) :: x, y
    INTEGER,         OPTIONAL, INTENT(IN) :: n_sigfig
    LOGICAL :: is_comparable
    ! Variables
    INTEGER :: ic, ip, k
    INTEGER :: n

    ! Set up
    is_comparable = .false.
    IF ( PRESENT(n_sigfig) ) THEN
      n = abs(n_sigfig)
    ELSE
      n = default_n_sigfig
    END IF

    ! Check the structure association status
    IF ( (.NOT. crtm_atmoptics_associated(x)) .OR. &
         (.NOT. crtm_atmoptics_associated(y))      ) RETURN

    ! Check dimensions
    IF ( (x%n_Layers         /= y%n_Layers        ) .OR. &
         (x%n_Legendre_Terms /= y%n_Legendre_Terms) .OR. &
         (x%n_Phase_Elements /= y%n_Phase_Elements)      ) RETURN

    ! Check scalar components
    IF ( (x%Include_Scattering .NEQV. y%Include_Scattering ) .OR. &
         (x%lOffset             /=    y%lOffset            ) .OR. &
         (.NOT. compares_within_tolerance(x%Scattering_Optical_Depth,&
                                          y%Scattering_Optical_Depth, n)) ) RETURN


    ! Check floating point arrays
    k  = x%n_Layers
    ip = x%n_Phase_Elements
    ic = x%n_Legendre_Terms
    IF ( (.NOT. all(compares_within_tolerance( &
                      x%Optical_Depth(1:k), &
                      y%Optical_Depth(1:k), &
                      n                      ))) .OR. &
         (.NOT. all(compares_within_tolerance( &
                      x%Single_Scatter_Albedo(1:k), &
                      y%Single_Scatter_Albedo(1:k), &
                      n                      ))) .OR. &
         (.NOT. all(compares_within_tolerance( &
                      x%Asymmetry_Factor(1:k), &
                      y%Asymmetry_Factor(1:k), &
                      n                      ))) .OR. &
         (.NOT. all(compares_within_tolerance( &
                      x%Delta_Truncation(1:k), &
                      y%Delta_Truncation(1:k), &
                      n                      ))) .OR. &
         (.NOT. all(compares_within_tolerance( &
                      x%Phase_Coefficient(0:ic,1:ip,1:k), &
                      y%Phase_Coefficient(0:ic,1:ip,1:k), &
                      n                      ))) ) RETURN


    ! If we get here, the structures are comparable
    is_comparable = .true.

  END FUNCTION crtm_atmoptics_compare


!------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_AtmOptics_InquireFile
!
! PURPOSE:
!       Function to inquire AtmOptics object files.
!
! CALLING SEQUENCE:
!       Error_Status = CRTM_AtmOptics_InquireFile( &
!                        Filename             , &
!                        n_Objects = n_Objects, &
!                        Release   = Release    )
!
! INPUTS:
!       Filename:           Character string specifying the name of the
!                           data file to inquire.
!                           UNITS:      N/A
!                           TYPE:       CHARACTER(*)
!                           DIMENSION:  Scalar
!                           ATTRIBUTES: INTENT(IN)
!
! OPTIONAL OUTPUTS:
!       n_Objects:          Number of AtmOptics objects contained in the file.
!                           UNITS:      N/A
!                           TYPE:       INTEGER
!                           DIMENSION:  Scalar
!                           ATTRIBUTES: INTENT(OUT), OPTIONAL
!
!       Release:            The data/file release number. Used to check
!                           for data/software mismatch.
!                           UNITS:      N/A
!                           TYPE:       INTEGER
!                           DIMENSION:  Scalar
!                           ATTRIBUTES: INTENT(OUT), OPTIONAL
!
! 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 file inquire was successful
!                              == FAILURE an unrecoverable error occurred.
!                           UNITS:      N/A
!                           TYPE:       INTEGER
!                           DIMENSION:  Scalar
!
!:sdoc-:
!------------------------------------------------------------------------------

  FUNCTION crtm_atmoptics_inquirefile( &
    Filename , &  ! Input
    n_Objects, &  ! Optional output
    Release  , &  ! Optional output
    Title    , &  ! Optional output
    History  , &  ! Optional output
    Comment  ) &  ! Optional output
  result( err_stat )
    ! Arguments
    CHARACTER(*),           INTENT(IN)  :: filename
    INTEGER     , OPTIONAL, INTENT(OUT) :: n_objects
    INTEGER     , OPTIONAL, INTENT(OUT) :: release
    CHARACTER(*), OPTIONAL, INTENT(OUT) :: title
    CHARACTER(*), OPTIONAL, INTENT(OUT) :: history
    CHARACTER(*), OPTIONAL, INTENT(OUT) :: comment
    ! Function result
    INTEGER :: err_stat
    ! Function parameters
    CHARACTER(*), PARAMETER :: routine_name = 'CRTM_AtmOptics_InquireFile'
    ! Function variables
    CHARACTER(ML) :: msg
    CHARACTER(ML) :: io_msg
    INTEGER :: io_stat
    INTEGER :: fid
    INTEGER :: n
    TYPE(crtm_atmoptics_type) :: atmoptics


    ! Setup
    err_stat = success
    ! ...Check that the file exists
    IF ( .NOT. file_exists( filename ) ) THEN
      msg = 'File '//trim(filename)//' not found.'
      CALL inquire_cleanup(); RETURN
    END IF


    ! Open the file
    err_stat = open_binary_file( filename, fid )
    IF ( err_stat /= success ) THEN
      msg = 'Error opening '//trim(filename)
      CALL inquire_cleanup(); RETURN
    END IF


    ! Read the release
    READ( fid, iostat=io_stat, iomsg=io_msg ) &
      atmoptics%Release
    IF ( io_stat /= 0 ) THEN
      msg = 'Error reading Release - '//trim(io_msg)
      CALL inquire_cleanup(); RETURN
    END IF
    IF ( .NOT. crtm_atmoptics_validrelease( atmoptics ) ) THEN
      msg = 'AtmOptics Release check failed.'
      CALL inquire_cleanup(); RETURN
    END IF


    ! Read the number of objects
    READ( fid, iostat=io_stat,iomsg=io_msg ) n
    IF ( io_stat /= 0 ) THEN
      msg = 'Error reading dimensions from '//trim(filename)//' - '//trim(io_msg)
      CALL inquire_cleanup(); RETURN
    END IF


    ! Read the global attributes
    err_stat = readgatts_binary_file( &
                 fid, &
                 title   = title  , &
                 history = history, &
                 comment = comment  )
    IF ( err_stat /= success ) THEN
      msg = 'Error reading global attributes'
      CALL inquire_cleanup(); RETURN
    END IF


    ! Close the file
    CLOSE( fid, iostat=io_stat, iomsg=io_msg )
    IF ( io_stat /= 0 ) THEN
      msg = 'Error closing '//trim(filename)//' - '//trim(io_msg)
      CALL inquire_cleanup(); RETURN
    END IF


    ! Assign the return arguments
    IF ( PRESENT(n_objects) ) n_objects = n
    IF ( PRESENT(release  ) ) release   = atmoptics%Release

  CONTAINS

    SUBROUTINE inquire_cleanup()
      IF ( file_open(fid) ) THEN
        CLOSE( fid, iostat=io_stat, iomsg=io_msg )
        IF ( io_stat /= 0 ) &
          msg = trim(msg)//'; Error closing input file during error cleanup - '//trim(io_msg)
      END IF
      err_stat = failure
      CALL display_message( routine_name, msg, err_stat )
    END SUBROUTINE inquire_cleanup

  END FUNCTION crtm_atmoptics_inquirefile


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_AtmOptics_ReadFile
!
! PURPOSE:
!       Function to read AtmOptics object files.
!
! CALLING SEQUENCE:
!       Error_Status = CRTM_AtmOptics_ReadFile( &
!                        AtmOptics          , &
!                        Filename           , &
!                        No_Close = No_Close, &
!                        Quiet    = Quiet     )
!
! OBJECTS:
!       AtmOptics:      AtmOptics object array containing the data read from file.
!                       UNITS:      N/A
!                       TYPE:       CRTM_AtmOptics_type
!                       DIMENSION:  Rank-1
!                       ATTRIBUTES: INTENT(OUT), ALLOCATABLE
!
! INPUTS:
!       Filename:       Character string specifying the name of a
!                       AtmOptics data file to read.
!                       UNITS:      N/A
!                       TYPE:       CHARACTER(*)
!                       DIMENSION:  Scalar
!                       ATTRIBUTES: INTENT(IN)
!
! OPTIONAL INPUTS:
!       No_Close:       Set this logical argument to *NOT* close the datafile
!                       upon exiting this routine. This option is required if
!                       the AtmOptics data is embedded within another file.
!                       If == .FALSE., File is closed upon function exit [DEFAULT].
!                          == .TRUE.,  File is NOT closed upon function exit
!                       If not specified, default is .FALSE.
!                       UNITS:      N/A
!                       TYPE:       LOGICAL
!                       DIMENSION:  Scalar
!                       ATTRIBUTES: INTENT(IN), OPTIONAL
!
!       Quiet:          Set this logical argument to suppress INFORMATION
!                       messages being printed to stdout
!                       If == .FALSE., INFORMATION messages are OUTPUT [DEFAULT].
!                          == .TRUE.,  INFORMATION messages are SUPPRESSED.
!                       If not specified, default is .FALSE.
!                       UNITS:      N/A
!                       TYPE:       LOGICAL
!                       DIMENSION:  Scalar
!                       ATTRIBUTES: INTENT(IN), OPTIONAL
!
! 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 file read was successful
!                          == FAILURE, an unrecoverable error occurred.
!                       UNITS:      N/A
!                       TYPE:       INTEGER
!                       DIMENSION:  Scalar
!
!:sdoc-:
!------------------------------------------------------------------------------

  FUNCTION crtm_atmoptics_readfile( &
    AtmOptics, &  ! Output
    Filename , &  ! Input
    No_Close , &  ! Optional input
    Quiet    , &  ! Optional input
    Title    , &  ! Optional output
    History  , &  ! Optional output
    Comment  , &  ! Optional output
    Debug    ) &  ! Optional input (Debug output control)
  result( err_stat )
    ! Arguments
    TYPE(crtm_atmoptics_type), ALLOCATABLE, INTENT(OUT) :: atmoptics(:)
    CHARACTER(*),                           INTENT(IN)  :: filename
    LOGICAL,                      OPTIONAL, INTENT(IN)  :: no_close
    LOGICAL,                      OPTIONAL, INTENT(IN)  :: quiet
    CHARACTER(*),                 OPTIONAL, INTENT(OUT) :: title
    CHARACTER(*),                 OPTIONAL, INTENT(OUT) :: history
    CHARACTER(*),                 OPTIONAL, INTENT(OUT) :: comment
    LOGICAL,                      OPTIONAL, INTENT(IN)  :: debug
    ! Function result
    INTEGER :: err_stat
    ! Function parameters
    CHARACTER(*), PARAMETER :: routine_name = 'CRTM_AtmOptics_ReadFile'
    ! Function variables
    CHARACTER(ML) :: msg
    CHARACTER(ML) :: count_msg
    CHARACTER(ML) :: io_msg
    LOGICAL :: close_file
    LOGICAL :: noisy
    INTEGER :: io_stat
    INTEGER :: alloc_stat
    INTEGER :: fid
    INTEGER :: n, n_objects
    TYPE(crtm_atmoptics_type) :: dummy

    ! Setup
    err_stat = success
    ! ...Check No_Close argument
    close_file = .true.
    IF ( PRESENT(no_close) ) close_file = .NOT. no_close
    ! ...Check Quiet argument
    noisy = .true.
    IF ( PRESENT(quiet) ) noisy = .NOT. quiet
    ! ...Override Quiet settings if debug set.
    IF ( PRESENT(debug) ) THEN
      IF ( debug ) noisy = .true.
    END IF


    ! Check if the file is open.
    IF ( file_open( filename ) ) THEN
      ! ...Inquire for the logical unit number
      INQUIRE( file=filename, number=fid )
      ! ...Ensure it's valid
      IF ( fid < 0 ) THEN
        msg = 'Error inquiring '//trim(filename)//' for its FileID'
        CALL read_cleanup(); RETURN
      END IF
    ELSE
      ! ...Open the file if it exists
      IF ( file_exists( filename ) ) THEN
        err_stat = open_binary_file( filename, fid )
        IF ( err_stat /= success ) THEN
          msg = 'Error opening '//trim(filename)
          CALL read_cleanup(); RETURN
        END IF
      ELSE
        msg = 'File '//trim(filename)//' not found.'
        CALL read_cleanup(); RETURN
      END IF
    END IF


    ! Read and check the release and version
    READ( fid, iostat=io_stat, iomsg=io_msg ) &
      dummy%Release
    IF ( io_stat /= 0 ) THEN
      msg = 'Error reading Release/Version - '//trim(io_msg)
      CALL read_cleanup(); RETURN
    END IF
    IF ( .NOT. crtm_atmoptics_validrelease( dummy ) ) THEN
      msg = 'Release check failed.'
      CALL read_cleanup(); RETURN
    END IF


    ! Read the number of objects
    READ( fid, iostat=io_stat,iomsg=io_msg ) n_objects
    IF ( io_stat /= 0 ) THEN
      msg = 'Error reading dimensions from '//trim(filename)//' - '//trim(io_msg)
      CALL read_cleanup(); RETURN
    END IF


    ! Allocate the output array
    ALLOCATE( atmoptics(n_objects), stat=alloc_stat )
    IF ( alloc_stat /= 0 ) THEN
      msg = 'Error allocating output object array'
      CALL read_cleanup(); RETURN
    END IF


    ! Read the global attributes
    err_stat = readgatts_binary_file( &
                 fid, &
                 title   = title  , &
                 history = history, &
                 comment = comment  )
    IF ( err_stat /= success ) THEN
      msg = 'Error reading global attributes'
      CALL read_cleanup(); RETURN
    END IF


    ! Loop over all the objects
    read_loop: DO n = 1, n_objects


      ! Generate count message for error output
      WRITE(count_msg,'("for object (",i0,")")') n


      ! Read the dimensions
      READ( fid, iostat=io_stat, iomsg=io_msg ) &
        dummy%n_Layers        , &
        dummy%n_Legendre_Terms, &
        dummy%n_Phase_Elements
      IF ( io_stat /= 0 ) THEN
        msg = 'Error reading dimensions '//trim(count_msg)//' from '//&
              trim(filename)//' - '//trim(io_msg)
        CALL read_cleanup(); RETURN
      END IF
      ! ...Allocate the object
      CALL crtm_atmoptics_create( &
             atmoptics(n)          , &
             dummy%n_Layers        , &
             dummy%n_Legendre_Terms, &
             dummy%n_Phase_Elements  )
      IF ( .NOT. crtm_atmoptics_associated( atmoptics(n) ) ) THEN
        msg = 'Allocation failed '//trim(count_msg)
        CALL read_cleanup(); RETURN
      END IF


      ! Read the data
      ! ...Read the scalar data
      READ( fid, iostat=io_stat, iomsg=io_msg ) &
        atmoptics(n)%Scattering_Optical_Depth
      IF ( io_stat /= 0 ) THEN
        msg = 'Error reading scalar data '//trim(count_msg)//' from '//&
              trim(filename)//' - '//trim(io_msg)
        CALL read_cleanup(); RETURN
      END IF
      ! ...Read the profile data
      READ( fid, iostat=io_stat, iomsg=io_msg ) &
        atmoptics(n)%Optical_Depth        , &
        atmoptics(n)%Single_Scatter_Albedo, &
        atmoptics(n)%Asymmetry_Factor     , &
        atmoptics(n)%Delta_Truncation
      IF ( io_stat /= 0 ) THEN
        msg = 'Error reading profile data '//trim(count_msg)//' from '//&
              trim(filename)//' - '//trim(io_msg)
        CALL read_cleanup(); RETURN
      END IF
      ! ...Read the scattering phase matrix coefficients
      READ( fid, iostat=io_stat, iomsg=io_msg ) &
        atmoptics(n)%Phase_Coefficient
      IF ( io_stat /= 0 ) THEN
        msg = 'Error reading phase matrix coefficients '//trim(count_msg)//' from '//&
              trim(filename)//' - '//trim(io_msg)
        CALL read_cleanup(); RETURN
      END IF

    END DO read_loop


    ! Close the file
    IF ( close_file ) THEN
      CLOSE( fid, iostat=io_stat, iomsg=io_msg )
      IF ( io_stat /= 0 ) THEN
        msg = 'Error closing '//trim(filename)//' - '//trim(io_msg)
        CALL read_cleanup(); RETURN
      END IF
    END IF


    ! Output an info message
    IF ( noisy ) THEN
      WRITE( msg,'("Number of objects read from ",a,": ",i0)' ) trim(filename), n_objects
      CALL display_message( routine_name, msg, information )
    END IF

  CONTAINS

    SUBROUTINE read_cleanup()
      IF ( file_open(filename) ) THEN
        CLOSE( fid, iostat=io_stat, iomsg=io_msg )
        IF ( io_stat /= 0 ) &
          msg = trim(msg)//'; Error closing input file during error cleanup - '//trim(io_msg)
      END IF
      IF ( ALLOCATED(atmoptics) ) THEN
        DEALLOCATE( atmoptics, stat=alloc_stat )
        IF ( alloc_stat /= 0 ) &
          msg = trim(msg)//'; Error deallocating object array during error cleanup'
      END IF
      err_stat = failure
      CALL display_message( routine_name, msg, err_stat )
    END SUBROUTINE read_cleanup

  END FUNCTION crtm_atmoptics_readfile


!--------------------------------------------------------------------------------
!:sdoc+:
!
! NAME:
!       CRTM_AtmOptics_WriteFile
!
! PURPOSE:
!       Function to write AtmOptics object files.
!
! CALLING SEQUENCE:
!       Error_Status = CRTM_AtmOptics_WriteFile( &
!                        AtmOptics          , &
!                        Filename           , &
!                        No_Close = No_Close, &
!                        Quiet    = Quiet     )
!
! OBJECTS:
!       AtmOptics:      AtmOptics object array containing the data to write to file.
!                       UNITS:      N/A
!                       TYPE:       CRTM_AtmOptics_type
!                       DIMENSION:  Rank-1
!                       ATTRIBUTES: INTENT(IN)
!
! INPUTS:
!       Filename:       Character string specifying the name of a
!                       AtmOptics format data file to write.
!                       UNITS:      N/A
!                       TYPE:       CHARACTER(*)
!                       DIMENSION:  Scalar
!                       ATTRIBUTES: INTENT(IN)
!
! OPTIONAL INPUTS:
!       No_Close:       Set this logical argument to *NOT* close the datafile
!                       upon exiting this routine. This option is required if
!                       the AtmOptics data is to be embedded within another file.
!                       If == .FALSE., File is closed upon function exit [DEFAULT].
!                          == .TRUE.,  File is NOT closed upon function exit
!                       If not specified, default is .FALSE.
!                       UNITS:      N/A
!                       TYPE:       LOGICAL
!                       DIMENSION:  Scalar
!                       ATTRIBUTES: INTENT(IN), OPTIONAL
!
!       Quiet:          Set this logical argument to suppress INFORMATION
!                       messages being printed to stdout
!                       If == .FALSE., INFORMATION messages are OUTPUT [DEFAULT].
!                          == .TRUE.,  INFORMATION messages are SUPPRESSED.
!                       If not specified, default is .FALSE.
!                       UNITS:      N/A
!                       TYPE:       LOGICAL
!                       DIMENSION:  Scalar
!                       ATTRIBUTES: INTENT(IN), OPTIONAL
!
! 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 file write was successful
!                          == FAILURE, an unrecoverable error occurred.
!                       UNITS:      N/A
!                       TYPE:       INTEGER
!                       DIMENSION:  Scalar
!
!:sdoc-:
!------------------------------------------------------------------------------

  FUNCTION crtm_atmoptics_writefile( &
    AtmOptics, &  ! Input
    Filename , &  ! Input
    No_Close , &  ! Optional input
    Quiet    , &  ! Optional input
    Title    , &  ! Optional input
    History  , &  ! Optional input
    Comment  , &  ! Optional input
    Debug    ) &  ! Optional input (Debug output control)
  result( err_stat )
    ! Arguments
    TYPE(crtm_atmoptics_type), INTENT(IN) :: atmoptics(:)
    CHARACTER(*),              INTENT(IN) :: filename
    LOGICAL,         OPTIONAL, INTENT(IN) :: no_close
    LOGICAL,         OPTIONAL, INTENT(IN) :: quiet
    CHARACTER(*),    OPTIONAL, INTENT(IN) :: title
    CHARACTER(*),    OPTIONAL, INTENT(IN) :: history
    CHARACTER(*),    OPTIONAL, INTENT(IN) :: comment
    LOGICAL,         OPTIONAL, INTENT(IN) :: debug
    ! Function result
    INTEGER :: err_stat
    ! Function parameters
    CHARACTER(*), PARAMETER :: routine_name = 'CRTM_AtmOptics_WriteFile'
    ! Function variables
    CHARACTER(ML) :: msg
    CHARACTER(ML) :: count_msg
    CHARACTER(ML) :: io_msg
    LOGICAL :: close_file
    LOGICAL :: noisy
    INTEGER :: io_stat
    INTEGER :: fid
    INTEGER :: n, n_objects


    ! Setup
    err_stat = success
    ! ...Check No_Close argument
    close_file = .true.
    IF ( PRESENT(no_close) ) close_file = .NOT. no_close
    ! ...Check Quiet argument
    noisy = .true.
    IF ( PRESENT(quiet) ) noisy = .NOT. quiet
    ! ...Override Quiet settings if debug set.
    IF ( PRESENT(debug) ) THEN
      IF ( debug ) noisy = .true.
    END IF
    ! ...Check there is data to write
    IF ( .NOT. all(crtm_atmoptics_associated( atmoptics )) ) THEN
      msg = 'Unassociated objects in input array.'
      CALL write_cleanup(); RETURN
    END IF
    n_objects= SIZE(atmoptics)


    ! Check if the file is open.
    IF ( file_open( filename ) ) THEN
      ! ...Inquire for the logical unit number
      INQUIRE( file=filename, number=fid )
      ! ...Ensure it's valid
      IF ( fid < 0 ) THEN
        msg = 'Error inquiring '//trim(filename)//' for its FileID'
        CALL write_cleanup(); RETURN
      END IF
    ELSE
      ! ...Open the file for output
      err_stat = open_binary_file( filename, fid, for_output=.true. )
      IF ( err_stat /= success ) THEN
        msg = 'Error opening '//trim(filename)
        CALL write_cleanup(); RETURN
      END IF
    END IF


    ! Write the release
    WRITE( fid, iostat=io_stat, iomsg=io_msg ) atmoptics_release
    IF ( io_stat /= 0 ) THEN
      msg = 'Error writing Release - '//trim(io_msg)
      CALL write_cleanup(); RETURN
    END IF


    ! Write the number of objects
    WRITE( fid, iostat=io_stat,iomsg=io_msg ) SIZE(atmoptics)
    IF ( io_stat /= 0 ) THEN
      msg = 'Error writing dimensions to '//trim(filename)//' - '//trim(io_msg)
      CALL write_cleanup(); RETURN
    END IF


    ! Write the global attributes
    err_stat = writegatts_binary_file( &
                 fid, &
                 write_module = module_version_id, &
                 title        = title  , &
                 history      = history, &
                 comment      = comment  )
    IF ( err_stat /= success ) THEN
      msg = 'Error writing global attributes'
      CALL write_cleanup(); RETURN
    END IF


    ! Loop over all the objects
    write_loop: DO n = 1, n_objects


      ! Generate count message for error output
      WRITE(count_msg,'("for object (",i0,")")') n


      ! Write the dimensions
      WRITE( fid, iostat=io_stat, iomsg=io_msg ) &
        atmoptics(n)%n_Layers        , &
        atmoptics(n)%n_Legendre_Terms, &
        atmoptics(n)%n_Phase_Elements
      IF ( io_stat /= 0 ) THEN
        msg = 'Error writing dimensions '//trim(count_msg)//' to '//&
              trim(filename)//' - '//trim(io_msg)
        CALL write_cleanup(); RETURN
      END IF


      ! Write the data
      ! ...Write the scalar data
      WRITE( fid, iostat=io_stat, iomsg=io_msg ) &
        atmoptics(n)%Scattering_Optical_Depth
      IF ( io_stat /= 0 ) THEN
        msg = 'Error writing scalar data '//trim(count_msg)//' to '//&
              trim(filename)//' - '//trim(io_msg)
        CALL write_cleanup(); RETURN
      END IF
      ! ...Write the profile data
      WRITE( fid, iostat=io_stat, iomsg=io_msg ) &
        atmoptics(n)%Optical_Depth(1:atmoptics(n)%n_Layers)        , &
        atmoptics(n)%Single_Scatter_Albedo(1:atmoptics(n)%n_Layers), &
        atmoptics(n)%Asymmetry_Factor(1:atmoptics(n)%n_Layers)     , &
        atmoptics(n)%Delta_Truncation(1:atmoptics(n)%n_Layers)
      IF ( io_stat /= 0 ) THEN
        msg = 'Error writing profile data '//trim(count_msg)//' to '//&
              trim(filename)//' - '//trim(io_msg)
        CALL write_cleanup(); RETURN
      END IF
      ! ...Write the scattering phase matrix coefficients
      WRITE( fid, iostat=io_stat, iomsg=io_msg ) &
        atmoptics(n)%Phase_Coefficient(0:atmoptics(n)%n_Legendre_Terms, &
                                       1:atmoptics(n)%n_Phase_Elements, &
                                       1:atmoptics(n)%n_Layers)
      IF ( io_stat /= 0 ) THEN
        msg = 'Error writing phase matrix coefficients '//trim(count_msg)//' to '//&
              trim(filename)//' - '//trim(io_msg)
        CALL write_cleanup(); RETURN
      END IF

    END DO write_loop


    ! Close the file
    IF ( close_file ) THEN
      CLOSE( fid, iostat=io_stat, iomsg=io_msg )
      IF ( io_stat /= 0 ) THEN
        msg = 'Error closing '//trim(filename)//' - '//trim(io_msg)
        CALL write_cleanup(); RETURN
      END IF
    END IF


    ! Output an info message
     IF ( noisy ) THEN
      WRITE( msg,'("Number of objects written to ",a,": ",i0)' ) trim(filename), n_objects
      CALL display_message( routine_name, msg, information )
     END IF

   CONTAINS

     SUBROUTINE write_cleanup()
       IF ( file_open(filename) ) THEN
         CLOSE( fid, iostat=io_stat, iomsg=io_msg )
         IF ( io_stat /= 0 ) &
           msg = trim(msg)//'; Error closing output file during error cleanup - '//trim(io_msg)
       END IF
       err_stat = failure
       CALL display_message( routine_name, msg, err_stat )
     END SUBROUTINE write_cleanup

  END FUNCTION crtm_atmoptics_writefile


!################################################################################
!################################################################################
!##                                                                            ##
!##                          ## PRIVATE PROCEDURES ##                          ##
!##                                                                            ##
!################################################################################
!################################################################################

!--------------------------------------------------------------------------------
!
! NAME:
!       CRTM_AtmOptics_Equal
!
! PURPOSE:
!       Elemental function to test the equality of two AtmOptics objects.
!       Used in OPERATOR(==) interface block.
!
! CALLING SEQUENCE:
!       is_equal = CRTM_AtmOptics_Equal( x, y )
!
!         or
!
!       IF ( x == y ) THEN
!         ...
!       END IF
!
! OBJECTS:
!       x, y:          Two AtmOptics objects to be compared.
!                      UNITS:      N/A
!                      TYPE:       CRTM_AtmOptics_type
!                      DIMENSION:  Scalar or any rank
!                      ATTRIBUTES: INTENT(IN)
!
! FUNCTION RESULT:
!       is_equal:      Logical value indicating whether the inputs are equal.
!                      UNITS:      N/A
!                      TYPE:       LOGICAL
!                      DIMENSION:  Same as inputs.
!
!--------------------------------------------------------------------------------

  ELEMENTAL FUNCTION crtm_atmoptics_equal( x, y ) RESULT( is_equal )
    TYPE(crtm_atmoptics_type), INTENT(IN) :: x, y
    LOGICAL :: is_equal

    ! Set up
    is_equal = .false.

    ! Check the object association status
    IF ( (.NOT. crtm_atmoptics_associated(x)) .OR. &
         (.NOT. crtm_atmoptics_associated(y))      ) RETURN

    ! Check contents
    ! ...Release/version info
    IF ( (x%Release /= y%Release) ) RETURN
    ! ...Dimensions
    IF ( (x%n_Layers         /= y%n_Layers        ) .OR. &
         (x%n_Legendre_Terms /= y%n_Legendre_Terms) .OR. &
         (x%n_Phase_Elements /= y%n_Phase_Elements) ) RETURN
    ! ...Scalar data
    IF ( x%Scattering_Optical_Depth .equalto. y%Scattering_Optical_Depth ) &
      is_equal = .true.
    ! ...Array data
    is_equal = is_equal .AND. &
               all(x%Optical_Depth(1:x%n_Layers)         .equalto. y%Optical_Depth(1:y%n_Layers)        ) .AND. &
               all(x%Single_Scatter_Albedo(1:x%n_Layers) .equalto. y%Single_Scatter_Albedo(1:y%n_Layers)) .AND. &
               all(x%Asymmetry_Factor(1:x%n_Layers)      .equalto. y%Asymmetry_Factor(1:y%n_Layers)     ) .AND. &
               all(x%Delta_Truncation(1:x%n_Layers)      .equalto. y%Delta_Truncation(1:y%n_Layers)     ) .AND. &
               all(x%Phase_Coefficient(0:x%n_Legendre_Terms, 1:x%n_Phase_Elements, 1:x%n_Layers) .equalto. &
                   y%Phase_Coefficient(0:y%n_Legendre_Terms, 1:y%n_Phase_Elements, 1:y%n_Layers) )

  END FUNCTION crtm_atmoptics_equal


!--------------------------------------------------------------------------------
!
! NAME:
!       CRTM_AtmOptics_Subtract
!
! PURPOSE:
!       Pure function to subtract two CRTM AtmOptics objects.
!       Used in OPERATOR(-) interface block.
!
! CALLING SEQUENCE:
!       aodiff = CRTM_AtmOptics_Subtract( ao1, ao2 )
!
!         or
!
!       aodiff = ao1 - ao2
!
!
! INPUTS:
!       ao1, ao2: The AtmOptics objects to difference.
!                   UNITS:      N/A
!                   TYPE:       CRTM_AtmOptics_type
!                   DIMENSION:  Scalar
!                   ATTRIBUTES: INTENT(IN OUT)
!
! RESULT:
!       aodiff:    AtmOptics object containing the differenced components.
!                   UNITS:      N/A
!                   TYPE:       CRTM_AtmOptics_type
!                   DIMENSION:  Scalar
!
!--------------------------------------------------------------------------------

  ELEMENTAL FUNCTION crtm_atmoptics_subtract( ao1, ao2 ) RESULT( aodiff )
    TYPE(crtm_atmoptics_type), INTENT(IN) :: ao1, ao2
    TYPE(crtm_atmoptics_type) :: aodiff
    INTEGER :: ic, ip, k

    ! Check input
    ! ...If input structures not allocated, do nothing
    IF ( (.NOT. crtm_atmoptics_associated(ao1)) .OR. &
         (.NOT. crtm_atmoptics_associated(ao2))      ) RETURN
    ! ...If input structure for different sizes, do nothing
    IF ( ao1%n_Layers         /= ao2%n_Layers         .OR. &
         ao1%n_Legendre_Terms /= ao2%n_Legendre_Terms .OR. &
         ao1%n_Phase_Elements /= ao2%n_Phase_Elements ) RETURN
    ! ...If input structure for different scattering setup, do nothing
    IF ( (ao1%Include_Scattering .NEQV. ao2%Include_Scattering ) .AND. &
         (ao1%lOffset             /=    ao2%lOffset            ) ) RETURN

    ! Copy the first structure
    aodiff = ao1

    ! And subtract the second one's components from it
    ! ...The scalar values
    aodiff%Scattering_Optical_Depth = aodiff%Scattering_Optical_Depth - ao2%Scattering_Optical_Depth
    ! ...The arrays
    k  = aodiff%n_Layers
    ip = aodiff%n_Phase_Elements
    ic = aodiff%n_Legendre_Terms
    aodiff%Optical_Depth(1:k)               = aodiff%Optical_Depth(1:k)               - ao2%Optical_Depth(1:k)
    aodiff%Single_Scatter_Albedo(1:k)       = aodiff%Single_Scatter_Albedo(1:k)       - ao2%Single_Scatter_Albedo(1:k)
    aodiff%Asymmetry_Factor(1:k)            = aodiff%Asymmetry_Factor(1:k)            - ao2%Asymmetry_Factor(1:k)
    aodiff%Delta_Truncation(1:k)            = aodiff%Delta_Truncation(1:k)            - ao2%Delta_Truncation(1:k)
    aodiff%Phase_Coefficient(0:ic,1:ip,1:k) = aodiff%Phase_Coefficient(0:ic,1:ip,1:k) - ao2%Phase_Coefficient(0:ic,1:ip,1:k)

  END FUNCTION crtm_atmoptics_subtract

END MODULE crtm_atmoptics_define