hand-merge PR #1124

Author: sdrabenh

GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSmoist_GridComp/Process_Library.F90

@@ -33,6 +33,13 @@ module GEOSmoist_Process_Library
     module procedure ICE_FRACTION_1D
     module procedure ICE_FRACTION_SC
   end interface ICE_FRACTION
+
+  ! SRF_TYPE constants
+  integer, parameter :: SRF_TYPE_LAND  = 1
+  integer, parameter :: SRF_TYPE_SNOW  = 2
+  integer, parameter :: SRF_TYPE_ICE   = 3
+  integer, parameter :: SRF_TYPE_OCEAN = 0
+
   ! ICE_FRACTION constants
    ! In anvil/convective clouds
    real, parameter :: aT_ICE_ALL = 252.16
@@ -525,7 +532,8 @@ function ICE_FRACTION_SC (TEMP,CNV_FRACTION,SRF_TYPE) RESULT(ICEFRCT)
       ICEFRCT_C = MAX(ICEFRCT_C,0.00)
       ICEFRCT_C = ICEFRCT_C**aICEFRPWR
      ! Sigmoidal functions like figure 6b/6c of Hu et al 2010, doi:10.1029/2009JD012384
-      if (SRF_TYPE >= 2.0) then
+      select case (nint(SRF_TYPE))
+      case (SRF_TYPE_SNOW, SRF_TYPE_ICE)
         ! Over snow (SRF_TYPE == 2.0) and ice (SRF_TYPE == 3.0)
         ICEFRCT_M  = 0.00
         if ( TEMP <= iT_ICE_ALL ) then
@@ -536,8 +544,8 @@ function ICE_FRACTION_SC (TEMP,CNV_FRACTION,SRF_TYPE) RESULT(ICEFRCT)
         ICEFRCT_M = MIN(ICEFRCT_M,1.00)
         ICEFRCT_M = MAX(ICEFRCT_M,0.00)
         ICEFRCT_M = ICEFRCT_M**iICEFRPWR
-      else if (SRF_TYPE == 1.0) then
-        ! Over Land
+      case (SRF_TYPE_LAND)
+        ! Over Land (SRF_TYPE == 1)
         ICEFRCT_M  = 0.00
         if ( TEMP <= lT_ICE_ALL ) then
            ICEFRCT_M = 1.000
@@ -547,8 +555,8 @@ function ICE_FRACTION_SC (TEMP,CNV_FRACTION,SRF_TYPE) RESULT(ICEFRCT)
         ICEFRCT_M = MIN(ICEFRCT_M,1.00)
         ICEFRCT_M = MAX(ICEFRCT_M,0.00)
         ICEFRCT_M = ICEFRCT_M**lICEFRPWR
-      else
-        ! Over Oceans
+      case (SRF_TYPE_OCEAN)
+        ! Over Oceans (SRF_TYPE == 0)
         ICEFRCT_M  = 0.00
         if ( TEMP <= oT_ICE_ALL ) then
            ICEFRCT_M = 1.000
@@ -558,7 +566,11 @@ function ICE_FRACTION_SC (TEMP,CNV_FRACTION,SRF_TYPE) RESULT(ICEFRCT)
         ICEFRCT_M = MIN(ICEFRCT_M,1.00)
         ICEFRCT_M = MAX(ICEFRCT_M,0.00)
         ICEFRCT_M = ICEFRCT_M**oICEFRPWR
-      endif
+      case default
+        ! You should not be here
+        print *, 'ICE_FRACTION_SC: Unknown SRF_TYPE = ',SRF_TYPE
+        error stop
+      end select
       ! Combine the Convective and MODIS functions
         ICEFRCT  = ICEFRCT_M*(1.0-CNV_FRACTION) + ICEFRCT_C*(CNV_FRACTION)
 #endif