Fix Missing Aerosol Scattering Jacobians and Restore Backscatter Logic

.gitignore

@@ -9,3 +9,4 @@ build.bash
 rebuild.bash
 .gitignore
 *~
+conductor/

src/AtmScatter/AerosolScatter/ASvar_Define.f90

@@ -241,6 +241,11 @@ SUBROUTINE ASvar_Inspect( self)
       WRITE(*,'(5x,"g Aerosol index #",i0)') i4
       WRITE(*,'(5(1x,es22.15,:))') self%g(:,i4)
     END DO
+    WRITE(*,'(3x,"Backscatter coefficient (kb) :")')
+    DO i4 = 1, self%n_Aerosols
+      WRITE(*,'(5x,"kb Aerosol index #",i0)') i4
+      WRITE(*,'(5(1x,es22.15,:))') self%kb(:,i4)
+    END DO
     WRITE(*,'(3x,"Phase coefficients (pcoeff) :")')
     DO i4 = 1, self%n_Aerosols
       WRITE(*,'(5x,"pcoeff Aerosol index #",i0)') i4
@@ -588,6 +593,13 @@ FUNCTION ASvar_ReadFile( &
       msg = 'Error reading asymmetry factor - '//TRIM(io_msg)
       CALL Read_Cleanup(); RETURN
     END IF
+    ! ...Backscatter coefficient
+    READ( fid, IOSTAT=io_stat, IOMSG=io_msg ) &
+      ASvar%kb
+    IF ( io_stat /= 0 ) THEN
+      msg = 'Error reading backscatter coefficient - '//TRIM(io_msg)
+      CALL Read_Cleanup(); RETURN
+    END IF
     ! ...Phase coefficients
     READ( fid, IOSTAT=io_stat, IOMSG=io_msg ) &
       ASvar%pcoeff
@@ -763,6 +775,13 @@ FUNCTION ASvar_WriteFile( &
       msg = 'Error writing asymmetry factor - '//TRIM(io_msg)
       CALL Write_Cleanup(); RETURN
     END IF
+    ! ...Backscatter coefficient
+    WRITE( fid, IOSTAT=io_stat, IOMSG=io_msg ) &
+      ASvar%kb
+    IF ( io_stat /= 0 ) THEN
+      msg = 'Error writing backscatter coefficient - '//TRIM(io_msg)
+      CALL Write_Cleanup(); RETURN
+    END IF
     ! ...Phase coefficients
     WRITE( fid, IOSTAT=io_stat, IOMSG=io_msg ) &
       ASvar%pcoeff
@@ -842,6 +861,7 @@ ELEMENTAL FUNCTION ASvar_Equal( x, y ) RESULT( is_equal )
     IF ( ALL(x%ke       .EqualTo. y%ke       ) .AND. &
          ALL(x%w        .EqualTo. y%w        ) .AND. &
          ALL(x%g        .EqualTo. y%g        ) .AND. &
+         ALL(x%kb       .EqualTo. y%kb       ) .AND. &
          ALL(x%pcoeff   .EqualTo. y%pcoeff   ) .AND. &
          ALL(x%total_bs .EqualTo. y%total_bs ) ) &
       is_equal = .TRUE.

src/Atmosphere/CRTM_Atmosphere.f90

@@ -389,6 +389,7 @@ FUNCTION CRTM_Atmosphere_AddLayers( &
     IF ( Atm_In%n_Aerosols > 0 ) THEN
       DO i = 1, Atm_In%n_Aerosols
         Atm_Out%Aerosol(i)%Effective_Radius(1:n) = ZERO
+        Atm_Out%Aerosol(i)%Effective_Variance(1:n) = ZERO
         Atm_Out%Aerosol(i)%Concentration(1:n)    = ZERO
       END DO
     END IF
@@ -622,6 +623,8 @@ FUNCTION CRTM_Atmosphere_AddLayers_AD( &
                                                    Atm_Out_AD%Aerosol(i)%Concentration(n+1:nt)
         Atm_In_AD%Aerosol(i)%Effective_Radius(1:no) = Atm_In_AD%Aerosol(i)%Effective_Radius(1:no) + &
                                                       Atm_Out_AD%Aerosol(i)%Effective_Radius(n+1:nt)
+        Atm_In_AD%Aerosol(i)%Effective_Variance(1:no) = Atm_In_AD%Aerosol(i)%Effective_Variance(1:no) + &
+                                                        Atm_Out_AD%Aerosol(i)%Effective_Variance(n+1:nt)
         Atm_In_AD%Aerosol(i)%Type = Atm_Out_AD%Aerosol(i)%Type
       END DO
     END IF
@@ -645,6 +648,9 @@ FUNCTION CRTM_Atmosphere_AddLayers_AD( &
     END DO
     ! ...Temperature data
     Atm_In_AD%Temperature(1:no) = Atm_In_AD%Temperature(1:no) + Atm_Out_AD%Temperature(n+1:nt)
+    ! ...Relative humidity data
+    Atm_In_AD%Relative_Humidity(1:no) = Atm_In_AD%Relative_Humidity(1:no) + &
+                                        Atm_Out_AD%Relative_Humidity(n+1:nt)
     ! ...Pressure data
     Atm_In_AD%Pressure(1:no)       = Atm_In_AD%Pressure(1:no) + Atm_Out_AD%Pressure(n+1:nt)
     Atm_In_AD%Level_Pressure(0:no) = Atm_In_AD%Level_Pressure(0:no) + Atm_Out_AD%Level_Pressure(n:nt)
@@ -748,6 +754,7 @@ FUNCTION CRTM_Atmosphere_ClearSkyCopy( atm, atm_clear ) RESULT( err_stat )
     atm_clear%Level_Pressure = atm%Level_Pressure(0:k)
     atm_clear%Pressure       = atm%Pressure(1:k)
     atm_clear%Temperature    = atm%Temperature(1:k)
+    atm_clear%Relative_Humidity = atm%Relative_Humidity(1:k)
     atm_clear%Absorber       = atm%Absorber(1:k,:)
     atm_clear%Cloud_Fraction = atm%Cloud_Fraction(1:k)
     ! ...Aerosol components
@@ -884,6 +891,7 @@ FUNCTION CRTM_Atmosphere_ClearSkyCopy_TL( &
     atm_clear_TL%Level_Pressure = atm_TL%Level_Pressure(0:k)
     atm_clear_TL%Pressure       = atm_TL%Pressure(1:k)
     atm_clear_TL%Temperature    = atm_TL%Temperature(1:k)
+    atm_clear_TL%Relative_Humidity = atm_TL%Relative_Humidity(1:k)
     atm_clear_TL%Absorber       = atm_TL%Absorber(1:k,:)
     atm_clear_TL%Cloud_Fraction = atm_TL%Cloud_Fraction(1:k)
     ! ...Aerosol components
@@ -1019,6 +1027,7 @@ FUNCTION CRTM_Atmosphere_ClearSkyCopy_AD( &
     atm_AD%Level_Pressure(0:k) = atm_AD%Level_Pressure(0:k) + atm_clear_AD%Level_Pressure(0:k)
     atm_AD%Pressure(1:k)       = atm_AD%Pressure(1:k)       + atm_clear_AD%Pressure(1:k)      
     atm_AD%Temperature(1:k)    = atm_AD%Temperature(1:k)    + atm_clear_AD%Temperature(1:k)   
+    atm_AD%Relative_Humidity(1:k) = atm_AD%Relative_Humidity(1:k) + atm_clear_AD%Relative_Humidity(1:k)
     atm_AD%Absorber(1:k,:)     = atm_AD%Absorber(1:k,:)     + atm_clear_AD%Absorber(1:k,:)      
     atm_AD%Cloud_Fraction(1:k) = atm_AD%Cloud_Fraction(1:k) + atm_clear_AD%Cloud_Fraction(1:k)
 

src/Coefficients/AerosolCoeff/AerosolCoeff_Inspect/AerosolCoeff_Inspect.f90

@@ -27,20 +27,26 @@ PROGRAM AerosolCoeff_Inspect
   ! Parameters
   ! ----------
   CHARACTER(*), PARAMETER :: PROGRAM_NAME = 'AerosolCoeff_Inspect'
-  CHARACTER(*), PARAMETER :: PROGRAM_RCS_ID = &
+  CHARACTER(*), PARAMETER :: PROGRAM_RCS_ID = '$Id$'
 
   ! ---------
   ! Variables
   ! ---------
   INTEGER :: Error_Status
   CHARACTER(256) :: Filename
+  CHARACTER(256) :: Aerosol_Model
   TYPE(AerosolCoeff_type) :: A
 
   ! Output program header
   CALL Program_Message( PROGRAM_NAME, &
                         'Program to display the contents of a CRTM Binary format '//&
                         'AerosolCoeff file to stdout.', &
-                        '$Revision$' )
+                        PROGRAM_RCS_ID )
+
+  ! Get the aerosol model
+  WRITE( *,FMT='(/5x,"Enter the Aerosol Model (CRTM, CMAQ, GOCART-GEOS5, or NAAPS): ")',ADVANCE='NO' )
+  READ( *,'(a)' ) Aerosol_Model
+  Aerosol_Model = ADJUSTL( Aerosol_Model )
 
   ! Get the filename
   WRITE( *,FMT='(/5x,"Enter the Binary AerosolCoeff filename: ")',ADVANCE='NO' )
@@ -54,7 +60,7 @@ PROGRAM AerosolCoeff_Inspect
   END IF
 
   ! Read the binary data file
-  Error_Status = AerosolCoeff_Binary_ReadFile( Filename, A )
+  Error_Status = AerosolCoeff_Binary_ReadFile( Aerosol_Model, Filename, A )
   IF ( Error_Status /= SUCCESS ) THEN
     CALL Display_Message( PROGRAM_NAME, &
                           'Error reading Binary AerosolCoeff file '//&

src/Coefficients/AerosolCoeff/AerosolCoeff_netCDF_IO.f90

@@ -1152,16 +1152,21 @@ FUNCTION AerosolCoeff_netCDF_ReadFile( &
     ! ...kb variable, only for GOCART-GEOS5 scheme
     IF ( Aerosol_Model == 'GOCART-GEOS5' ) THEN
       NF90_Status = NF90_INQ_VARID( FileId,KB_VARNAME,VarId )
-      IF ( NF90_Status /= NF90_NOERR ) THEN
-        msg = 'Error inquiring '//TRIM(Filename)//' for '//KB_VARNAME//&
-              ' variable ID - '//TRIM(NF90_STRERROR( NF90_Status ))
-        CALL Read_Cleanup(); RETURN
-      END IF
-      NF90_Status = NF90_GET_VAR( FileId,VarID,AerosolCoeff%kb )
-      IF ( NF90_Status /= NF90_NOERR ) THEN
-        msg = 'Error reading '//KB_VARNAME//' from '//TRIM(Filename)//&
-              ' - '//TRIM(NF90_STRERROR( NF90_Status ))
-        CALL Read_Cleanup(); RETURN
+      ! Older GOCART-GEOS5 files may not include kb, so treat it as optional.
+      IF ( NF90_Status .EQ. NF90_ENOTVAR ) THEN
+        AerosolCoeff%kb = 0.0_fp
+      ELSE
+        IF ( NF90_Status /= NF90_NOERR ) THEN
+          msg = 'Error inquiring '//TRIM(Filename)//' for '//KB_VARNAME//&
+                ' variable ID - '//TRIM(NF90_STRERROR( NF90_Status ))
+          CALL Read_Cleanup(); RETURN
+        END IF
+        NF90_Status = NF90_GET_VAR( FileId,VarID,AerosolCoeff%kb )
+        IF ( NF90_Status /= NF90_NOERR ) THEN
+          msg = 'Error reading '//KB_VARNAME//' from '//TRIM(Filename)//&
+                ' - '//TRIM(NF90_STRERROR( NF90_Status ))
+          CALL Read_Cleanup(); RETURN
+        END IF
       END IF
     END IF
     ! ...ke variable

test/CMakeLists.txt

@@ -308,6 +308,24 @@ add_test(NAME test_Unit_Aerosol_Bypass
          COMMAND $<TARGET_FILE:Unit_Aerosol_Bypass>)
 set_tests_properties(test_Unit_Aerosol_Bypass PROPERTIES ENVIRONMENT "OMP_NUM_THREADS=$ENV{OMP_NUM_THREADS}")
 
+add_executable(Unit_AerosolScatter_TL mains/unit/Unit_Test/test_AerosolScatter_TL.f90)
+target_link_libraries(Unit_AerosolScatter_TL PRIVATE crtm)
+add_test(NAME test_Unit_AerosolScatter_TL
+         COMMAND $<TARGET_FILE:Unit_AerosolScatter_TL>)
+set_tests_properties(test_Unit_AerosolScatter_TL PROPERTIES ENVIRONMENT "OMP_NUM_THREADS=$ENV{OMP_NUM_THREADS}")
+
+add_executable(Unit_AerosolScatter_AD mains/unit/Unit_Test/test_AerosolScatter_AD.f90)
+target_link_libraries(Unit_AerosolScatter_AD PRIVATE crtm)
+add_test(NAME test_Unit_AerosolScatter_AD
+         COMMAND $<TARGET_FILE:Unit_AerosolScatter_AD>)
+set_tests_properties(test_Unit_AerosolScatter_AD PROPERTIES ENVIRONMENT "OMP_NUM_THREADS=$ENV{OMP_NUM_THREADS}")
+
+add_executable(Unit_AerosolScatter_K mains/unit/Unit_Test/test_AerosolScatter_K.f90)
+target_link_libraries(Unit_AerosolScatter_K PRIVATE crtm)
+add_test(NAME test_Unit_AerosolScatter_K
+         COMMAND $<TARGET_FILE:Unit_AerosolScatter_K>)
+set_tests_properties(test_Unit_AerosolScatter_K PROPERTIES ENVIRONMENT "OMP_NUM_THREADS=$ENV{OMP_NUM_THREADS}")
+
 add_executable(Unit_Aerosol_Bypass_TL mains/unit/Unit_Test/test_Aerosol_Bypass_TL.f90)
 target_link_libraries(Unit_Aerosol_Bypass_TL PRIVATE crtm)
 add_test(NAME test_Unit_Aerosol_Bypass_TL
@@ -378,6 +396,15 @@ foreach(testtype IN LISTS CCoef_Utils)
 	target_link_libraries(${testtype} PRIVATE crtm)
 endforeach()
 
+#AerosolCoeff utilities
+list (APPEND ACoeff_Utils
+  AerosolCoeff_Inspect
+)
+foreach(testtype IN LISTS ACoeff_Utils)
+	add_executable(${testtype} ${CRTM_SOURCE_DIR}/src/Coefficients/AerosolCoeff/${testtype}/${testtype}.f90)
+	target_link_libraries(${testtype} PRIVATE crtm)
+endforeach()
+
 list (APPEND SEcategoryCoef_Utils
   SEcategory_NC2BIN
   SEcategory_BIN2NC
@@ -546,6 +573,10 @@ set_tests_properties(test_TL_convergence_active_sensor PROPERTIES ENVIRONMENT OM
 
 list( APPEND crtm_test_input
 AerosolCoeff/Little_Endian/AerosolCoeff.bin
+AerosolCoeff/Little_Endian/AerosolCoeff.GOCART-GEOS5.bin
+AerosolCoeff/Little_Endian/AerosolCoeff.CMAQ.bin
+AerosolCoeff/netCDF/AerosolCoeff.GOCART-GEOS5.nc4
+AerosolCoeff/netCDF/AerosolCoeff.GOCART-GEOS5.BRC.kb.v2.nc
 AerosolCoeff/netCDF/AerosolCoeff.nc4
 CloudCoeff/Little_Endian/CloudCoeff.bin
 CloudCoeff/netCDF/CloudCoeff.nc4

test/mains/regression/forward/test_AOD/test_AOD.f90

@@ -12,6 +12,7 @@ PROGRAM test_AOD
   !
   ! Module usage
   USE CRTM_Module
+  USE File_Utility, ONLY: File_Exists
   ! Disable all implicit typing
   IMPLICIT NONE
   ! ============================================================================
@@ -45,10 +46,13 @@ PROGRAM test_AOD
   CHARACTER(256) :: Message
   CHARACTER(256) :: Version
   CHARACTER(256) :: Sensor_Id
+  CHARACTER(256) :: AerosolCoeff_File
   INTEGER :: Error_Status
   INTEGER :: Allocate_Status
   INTEGER :: n_Channels
   INTEGER :: l, m
+  LOGICAL :: has_new_coeff
+  LOGICAL :: has_old_coeff
   ! Declarations for RTSolution comparison
   INTEGER :: n_l, n_m
   CHARACTER(256) :: rts_File
@@ -92,12 +96,28 @@ PROGRAM test_AOD
   ! 2a. Initialise for the requested sensor
   ! ---------------------------------------
   WRITE( *,'(/5x,"Initializing the CRTM...")' )
+  has_new_coeff = File_Exists(COEFFICIENTS_PATH//'AerosolCoeff.GOCART-GEOS5.BRC.kb.v2.nc')
+  has_old_coeff = File_Exists(COEFFICIENTS_PATH//'AerosolCoeff.GOCART-GEOS5.nc4')
+  IF ( has_new_coeff ) THEN
+    AerosolCoeff_File = 'AerosolCoeff.GOCART-GEOS5.BRC.kb.v2.nc'
+  ELSE
+    AerosolCoeff_File = 'AerosolCoeff.GOCART-GEOS5.nc4'
+  END IF
   Error_Status = CRTM_Init( (/Sensor_Id/), &
                             ChannelInfo, &
                             Aerosol_Model = 'GOCART-GEOS5', &
                             AerosolCoeff_Format = 'netCDF', &
-                            AerosolCoeff_File = 'AerosolCoeff.GOCART-GEOS5.nc4', &
+                            AerosolCoeff_File = TRIM(AerosolCoeff_File), &
                             File_Path=COEFFICIENTS_PATH)
+  IF ( Error_Status /= SUCCESS .AND. has_old_coeff .AND. TRIM(AerosolCoeff_File) /= 'AerosolCoeff.GOCART-GEOS5.nc4' ) THEN
+    AerosolCoeff_File = 'AerosolCoeff.GOCART-GEOS5.nc4'
+    Error_Status = CRTM_Init( (/Sensor_Id/), &
+                              ChannelInfo, &
+                              Aerosol_Model = 'GOCART-GEOS5', &
+                              AerosolCoeff_Format = 'netCDF', &
+                              AerosolCoeff_File = TRIM(AerosolCoeff_File), &
+                              File_Path=COEFFICIENTS_PATH)
+  END IF
   IF ( Error_Status /= SUCCESS ) THEN
     Message = 'Error initializing CRTM'
     CALL Display_Message( PROGRAM_NAME, Message, FAILURE )

test/mains/regression/k_matrix/test_AOD/test_AOD.f90

@@ -12,6 +12,7 @@ PROGRAM test_AOD
   !
   ! Module usage
   USE CRTM_Module
+  USE File_Utility, ONLY: File_Exists
   ! Disable all implicit typing
   IMPLICIT NONE
   ! ============================================================================
@@ -43,10 +44,13 @@ PROGRAM test_AOD
   CHARACTER(256) :: Message
   CHARACTER(256) :: Version
   CHARACTER(256) :: Sensor_Id
+  CHARACTER(256) :: AerosolCoeff_File
   INTEGER :: Error_Status
   INTEGER :: Allocate_Status
   INTEGER :: n_Channels
   INTEGER :: l, m
+  LOGICAL :: has_new_coeff
+  LOGICAL :: has_old_coeff
   ! Declarations for Jacobian comparisons
   INTEGER :: n_l, n_m
   CHARACTER(256) :: atmk_File
@@ -95,12 +99,28 @@ PROGRAM test_AOD
   ! 2a. Initialise for the requested sensor
   ! ---------------------------------------
   WRITE( *,'(/5x,"Initializing the CRTM...")' )
+  has_new_coeff = File_Exists(COEFFICIENTS_PATH//'AerosolCoeff.GOCART-GEOS5.BRC.kb.v2.nc')
+  has_old_coeff = File_Exists(COEFFICIENTS_PATH//'AerosolCoeff.GOCART-GEOS5.nc4')
+  IF ( has_new_coeff ) THEN
+    AerosolCoeff_File = 'AerosolCoeff.GOCART-GEOS5.BRC.kb.v2.nc'
+  ELSE
+    AerosolCoeff_File = 'AerosolCoeff.GOCART-GEOS5.nc4'
+  END IF
   Error_Status = CRTM_Init( (/Sensor_Id/), &
                             ChannelInfo, &
                             Aerosol_Model = 'GOCART-GEOS5', &
                             AerosolCoeff_Format = 'netCDF', &
-                            AerosolCoeff_File = 'AerosolCoeff.GOCART-GEOS5.nc4', &
+                            AerosolCoeff_File = TRIM(AerosolCoeff_File), &
                             File_Path=COEFFICIENTS_PATH)
+  IF ( Error_Status /= SUCCESS .AND. has_old_coeff .AND. TRIM(AerosolCoeff_File) /= 'AerosolCoeff.GOCART-GEOS5.nc4' ) THEN
+    AerosolCoeff_File = 'AerosolCoeff.GOCART-GEOS5.nc4'
+    Error_Status = CRTM_Init( (/Sensor_Id/), &
+                              ChannelInfo, &
+                              Aerosol_Model = 'GOCART-GEOS5', &
+                              AerosolCoeff_Format = 'netCDF', &
+                              AerosolCoeff_File = TRIM(AerosolCoeff_File), &
+                              File_Path=COEFFICIENTS_PATH)
+  END IF
   IF ( Error_Status /= SUCCESS ) THEN
     Message = 'Error initializing CRTM'
     CALL Display_Message( PROGRAM_NAME, Message, FAILURE )