Added new flameh_opt to use both FRP and User set value.

README.md

@@ -93,7 +93,7 @@ Current Canopy-App components:
     | lamdars          | Value representing influence of roughness sublayer (Massman et al., 2017)          |
     | dx_opt           | 0=Calculation of dx resolution/distance from lon; 1=user set dx grid resolution    |
     | dx_set           | user set real value of grid resolution (m) only if dx_opt=1                        |
-    | flameh_opt       | 0=Calculation of flame height from FRP (Byram, 1959); 1=user set flameh            |
+    | flameh_opt       | 0=Calculation of flame height from FRP (Byram, 1959); 1=user set flameh; 2=FRP calculation where available (active fires), otherwise user set flameh    |
     | flameh_set       | user set real value of flame height (m) only if flame_opt=1                        |
     | pai_opt          | integer (0=PAI fixed from Katul et al. 2004 veg types-->default; 1=PAI Massman et al. 2017 Eq. 19 calc; 2=PAI from model LAI+WAI; 3=user set PAI value) |
     | pai_set          | user set real value of PAI (default=4.0; only used if pai_opt=3)                   |

canopy_alloc.F90

@@ -42,8 +42,9 @@ SUBROUTINE canopy_alloc
         if(.not.allocated(dx))            allocate(dx(nlat*nlon))
         if(.not.allocated(dx_2d))         allocate(dx_2d(nlon,nlat))
         if(.not.allocated(waf))           allocate(waf(nlat*nlon))
-        if(.not.allocated(waf_2d))        allocate(waf_2d(nlon,nlat))
-
+        if(.not.allocated(waf_2d))        allocate(waf_2d(nlon,nlon))
+        if(.not.allocated(flameh))        allocate(flameh(nlat*nlon))
+        if(.not.allocated(flameh_2d))     allocate(flameh_2d(nlon,nlon))
     end if
 
 !-------------------------------------------------------------------------------

canopy_calcs.F90

@@ -103,10 +103,10 @@ SUBROUTINE canopy_calcs
 
 ! ... determine midflamepoint and flame height from user or FRP calculation
                                 call canopy_flameh(flameh_opt, flameh_set, dx_2d(i,j), modres, &
-                                    frpref, midflamepoint, flameh)
+                                    frpref, midflamepoint, flameh_2d(i,j))
 
-                                if (flameh .gt. 0.0) then !only calculate WAF when flameh > 0
-                                    call canopy_waf(hcmref, lamdars, rsl_opt, hgtref, flameh, &
+                                if (flameh_2d(i,j) .gt. 0.0) then !only calculate WAF when flameh > 0
+                                    call canopy_waf(hcmref, lamdars, rsl_opt, hgtref, flameh_2d(i,j), &
                                         firetype, d_h, zo_h, canBOT(midflamepoint), &
                                         canTOP(midflamepoint), waf_2d(i,j))
                                 end if
@@ -208,10 +208,10 @@ SUBROUTINE canopy_calcs
                             end do
 ! ... determine midflamepoint and flame height from user or FRP calculation
                             call canopy_flameh(flameh_opt, flameh_set, dx(loc), modres, &
-                                frpref, midflamepoint, flameh)
+                                frpref, midflamepoint, flameh(loc))
 
-                            if (flameh .gt. 0.0) then !only calculate WAF when flameh > 0
-                                call canopy_waf(hcmref, lamdars, rsl_opt, hgtref, flameh, &
+                            if (flameh(loc) .gt. 0.0) then !only calculate WAF when flameh > 0
+                                call canopy_waf(hcmref, lamdars, rsl_opt, hgtref, flameh(loc), &
                                     firetype, d_h, zo_h, canBOT(midflamepoint), &
                                     canTOP(midflamepoint), waf(loc))
                             end if
@@ -314,10 +314,10 @@ SUBROUTINE canopy_calcs
                             end do
 ! ... determine midflamepoint and flame height from user or FRP calculation
                             call canopy_flameh(flameh_opt, flameh_set, dx(loc), modres, &
-                                frpref, midflamepoint, flameh)
+                                frpref, midflamepoint, flameh(loc))
 
-                            if (flameh .gt. 0.0) then !only calculate WAF when flameh > 0
-                                call canopy_waf(hcmref, lamdars, rsl_opt, hgtref, flameh, &
+                            if (flameh(loc) .gt. 0.0) then !only calculate WAF when flameh > 0
+                                call canopy_waf(hcmref, lamdars, rsl_opt, hgtref, flameh(loc), &
                                     firetype, d_h, zo_h, canBOT(midflamepoint), &
                                     canTOP(midflamepoint), waf(loc))
                             end if

canopy_canvars_mod.F90

@@ -23,7 +23,7 @@ MODULE canopy_canvars_mod
     real(rk)       ::    sigma1        !Standard deviation of shape function below zcanmax (z/h)
     real(rk)       ::    d_h           !Zero plane displacement heights (z/h)
     real(rk)       ::    zo_h          !Surface (soil+veg) roughness lengths (z/h)
-    real(rk)       ::    flameh        !Flame Height (m)
+    !real(rk)       ::    flameh        !Flame Height (m)
 
 !-------------------------------------------------------------------------------
 ! Allocatable canopy variable arrays
@@ -46,6 +46,8 @@ MODULE canopy_canvars_mod
     real(rk), allocatable :: Kz_3d         ( : , : , : )  ! Eddy Diffusivities -- 3D (m2/s)
     real(rk), allocatable :: rjcf          ( :, : )       ! Photolysis Attenuation Correction Factors
     real(rk), allocatable :: rjcf_3d       ( : , : , : )  ! Photolysis Attenuation Correction Factors -- 3D
+    real(rk), allocatable :: flameh        ( : )          ! Flame Height (m)
+    real(rk), allocatable :: flameh_2d     ( : , : )      ! Flame Height -- 2D (m)
 
 !-------------------------------------------------------------------------------
 ! Canopy-App Program and version descriptors.
@@ -104,6 +106,7 @@ MODULE canopy_canvars_mod
 
     TYPE(fld2ddata), ALLOCATABLE, TARGET :: fld2dxyt ( : )
     TYPE(fld2ddata), POINTER     :: c_waf
+    TYPE(fld2ddata), POINTER     :: c_flameh
 
 !-------------------------------------------------------------------------------
 ! Time-varying 3d fields at cell centers for output NETCDF

canopy_dealloc.F90

@@ -44,6 +44,9 @@ SUBROUTINE canopy_dealloc
         if(allocated(dx_2d))         deallocate(dx_2d)
         if(allocated(waf))           deallocate(waf)
         if(allocated(waf_2d))        deallocate(waf_2d)
+        if(allocated(flameh))        deallocate(flameh)
+        if(allocated(flameh_2d))     deallocate(flameh_2d)
+
     end if
 
 !-------------------------------------------------------------------------------

canopy_ncf_io_mod.F90

@@ -340,6 +340,19 @@ SUBROUTINE canopy_outncf_init
         c_waf%iend(1) = nlon
         c_waf%iend(2) = nlat
 
+        c_flameh%fld = fillreal
+        c_flameh%fldname = 'FLAMEH'
+        c_flameh%long_name = 'flame height'
+        c_flameh%units = 'm'
+        c_flameh%fillvalue = fillreal
+        c_flameh%dimnames(1) = 'nlon'
+        c_flameh%dimnames(2) = 'nlat'
+        c_flameh%istart(1) = 1
+        c_flameh%istart(2) = 1
+        c_flameh%iend(1) = nlon
+        c_flameh%iend(2) = nlat
+
+
         !-------------------------------------------------------------------------------
         ! Time-varying 3d fields at cell centers.
         !-------------------------------------------------------------------------------
@@ -432,13 +445,16 @@ SUBROUTINE canopy_outncf_alloc
 
         nfld2dxyt = nfld2dxyt + 1  !WAF
 
+        nfld2dxyt = nfld2dxyt + 1  !FLAMEH
+
         ALLOCATE ( fld2dxyt ( nfld2dxyt ) )
 
         DO nn = 1, nfld2dxyt
             ALLOCATE ( fld2dxyt(nn)%fld(nlon,nlat) )
         ENDDO
 
-        c_waf    => fld2dxyt( 1 )
+        c_waf       => fld2dxyt( 1 )
+        c_flameh    => fld2dxyt( 2 )
 
         !-------------------------------------------------------------------------------
         ! Time-varying 3d fields at cell centers.
@@ -1309,6 +1325,7 @@ SUBROUTINE canopy_write_ncf (OUTPREFX)
             !-------------------------------------------------------------------------------
 
             c_waf%fld = waf_2d
+            c_flameh%fld = flameh_2d
 
             !-------------------------------------------------------------------------------
             ! Time-varying 3d fields at cell centers.

canopy_txt_io_mod.F90

@@ -79,9 +79,10 @@ SUBROUTINE write_txt(TXTPREFX)
                 write(*,*)  '-------------------------------'
                 open(11, file=TRIM(TXTPREFX)//'_output_waf.txt')
                 write(11, '(a30, f6.1)') 'Reference height, h: ', href_set, 'm'
-                write(11, '(a8, a9, a19, a11)') 'Lat', 'Lon', 'Canopy height (m)', 'WAF'
+                write(11, '(a8, a9, a19, a19, a11)') 'Lat', 'Lon', 'Canopy height (m)', 'Flame height (m)', 'WAF'
                 do loc=1, nlat*nlon
-                    write(11, '(f8.2, f9.2, f19.2, es15.7)')  variables(loc)%lat, variables(loc)%lon, hcmref, waf(loc)
+                    write(11, '(f8.2, f9.2, f19.2, f19.2, es15.7)')  variables(loc)%lat, variables(loc)%lon, &
+                        variables(loc)%fh, flameh(loc), waf(loc)
                 end do
             end if
 

canopy_waf_mod.F90

@@ -49,6 +49,12 @@ SUBROUTINE CANOPY_FLAMEH( FLAMEH_OPT, FLAMEH_SET, DX, MODRES, &
             end if
         else if (FLAMEH_OPT .eq. 1) then  !user set value
             FLAMEH = FLAMEH_SET
+        else if (FLAMEH_OPT .eq. 2) then  !both FRP calc and user set
+            if (FRP .gt. 0.0) then
+                FLAMEH = CalcFlameH(FRP,DX)
+            else
+                FLAMEH = FLAMEH_SET
+            end if
         else
             write(*,*)  'Wrong FLAMEH_OPT choice of ', FLAMEH_OPT, ' in namelist...exiting'
             call exit(2)

namelist.canopy

@@ -22,8 +22,8 @@
   lamdars    =  1.25
   dx_opt     =  0
   dx_set     =  12000.0
-  flameh_opt =  0
-  flameh_set =  2.0
+  flameh_opt =  2
+  flameh_set =  1.0
   pai_opt    =  0
   pai_set    =  4.0
   lu_opt     =  0