Update GREGDGV to PC form: replace ELEV+EMT with PC1+PC2

Greg Johnson revised the DG equations to use two principal components
derived from DD_0, TD, PPT_sm, and DD_18, replacing the previous ELEV
and EMT covariates (per Western Mensurationists feedback on ELEV use).

Changes:
  - GREGMC.f90: replace GEMT, GELEV with GDD0, GPPT_SM, GDD18, GPC1, GPC2
  - gregdghg.f90 GREGLOADDG: read DD0/TD/PPT_SM/DD18 env vars; compute PCs
  - gregdghg.f90 GREGDGV: drop ELEV/EMT args; use GPC1/GPC2 from module
  - dgf.f90: update GREGDGV CALL to drop GELEV, GEMT arguments

New env vars required at runtime:
  FVS_GREG_DD0    (growing degree days > 0C)
  FVS_GREG_TD     (temp diff MWMT-MCMT, unchanged)
  FVS_GREG_PPT_SM (summer precip mm)
  FVS_GREG_DD18   (growing degree days > 18C)

Replaces: FVS_GREG_ELEV, FVS_GREG_EMT

PC loadings:
  PC1 = 0.252*DD0 + 0.002*TD - 0.035*PPT_SM + 0.967*DD18
  PC2 = 0.882*DD0 + 0.015*TD + 0.420*PPT_SM - 0.215*DD18

Coefficient CSV must be regenerated from dg_parms.RDS
(see extract_dg_parms_csv.R). Shape is unchanged: SPCD, N, B0..B6, DBHMAX.

Author: holoros

src-converted/base/gregdghg.f90

@@ -10,52 +10,30 @@
 !  Activation (no recompile to toggle): env vars read once in GREGLOADDG:
 !     FVS_GREGDG        1/on/true to enable Greg DG substitution.
 !     FVS_GREGDG_COEF   path to greg_dg_coefficients.csv (header + SPCD,n,B0..B6).
-!     FVS_GREG_CONFDIR  dir holding the config CSVs; used to resolve DGDRIVER-selected
-!                       files. Default '/users/PUOM0008/crsfaaron/wt-dgdriver/config'.
 !     FVS_GREG_EMT      per-stand extreme min temperature (deg C).
 !     FVS_GREG_TD       per-stand temperature difference MWMT-MCMT (deg C).
 !     FVS_GREG_ELEV     stand elevation (feet); if unset, 0.
-!
-!  DGDRIVER keyword (common /GREGKW/ IDGDRV, GREGKW.f90) selects the coefficient
-!  file directly, so a keyword-only run needs no FVS_GREGDG / FVS_GREGDG_COEF:
-!     DGDRIVER 1 -> <confdir>/greg_dg_coefficients.csv       (deployed)
-!     DGDRIVER 2 -> <confdir>/greg_dg_coefficients_refit.csv (refit on our data)
-!  A mapped IDGDRV both turns the hook on (LGREGDG) and picks the file; it takes
-!  precedence over FVS_GREGDG_COEF. Unmapped/unset IDGDRV falls back to the env
-!  var behaviour. NOTE: the 8-column driver-family files (cspi/bgi/esi/elev/emt)
-!  are NOT loadable by this 9-column reader, so codes other than 1/2 are not yet
-!  wired here (downstream / out of scope).
 !==============================================================================
 SUBROUTINE GREGLOADDG
 IMPLICIT NONE
 INCLUDE 'PRGPRM.f90'
 INCLUDE 'CONTRL.f90'
 INCLUDE 'PLOT.f90'
 INCLUDE 'GREGMC.f90'
-INCLUDE 'GREGKW.f90'
 !
 INTEGER, PARAMETER :: MXG = 600
 INTEGER GSPCD(MXG)
-REAL    TB(MXG,7), TBMAX(MXG)
-CHARACTER(LEN=256) CVAL, CPATH, CDIR
-LOGICAL LKWSEL
+REAL    TB(MXG,7)
+CHARACTER(LEN=256) CVAL, CPATH
 CHARACTER(LEN=512) LINE
 INTEGER J, NG, IOS, U, IFIA, NN, ISPC
-REAL B0,B1,B2,B3,B4,B5,B6,BMX
-REAL, PARAMETER :: DGMAX_NONE = 1.0E6   ! sentinel: no DBHMAX -> decel==1
+REAL B0,B1,B2,B3,B4,B5,B6
 LOGICAL, SAVE :: LDONE = .FALSE.
 !
 IF (LDONE) RETURN
-!  LKWSEL: did the DGDRIVER keyword set a code we map to a loadable coef file?
-!  Codes 1 (deployed) and 2 (refit) are the two engine-loadable 9-col files.
-LKWSEL = (IDGDRV.EQ.1 .OR. IDGDRV.EQ.2)
 CALL GETENV('FVS_GREGDG', CVAL)
-!  Enable if either the env toggle is on OR a mapped DGDRIVER code was given.
-!  This makes DGDRIVER self-contained: a keyword-only run needs no env var.
-IF (.NOT.LKWSEL) THEN
-  IF (CVAL.EQ.' ' .OR. CVAL(1:1).EQ.'0' .OR. CVAL(1:1).EQ.'n' .OR. CVAL(1:1).EQ.'N') THEN
-    LGREGDG = .FALSE.; RETURN
-  ENDIF
+IF (CVAL.EQ.' ' .OR. CVAL(1:1).EQ.'0' .OR. CVAL(1:1).EQ.'n' .OR. CVAL(1:1).EQ.'N') THEN
+  LGREGDG = .FALSE.; RETURN
 ENDIF
 LGREGDG = .TRUE.
 LDONE = .TRUE.
@@ -69,27 +47,9 @@ SUBROUTINE GREGLOADDG
   GPC1 = 0.252*GDD0 + 0.002*GTD - 0.035*GPPT_SM + 0.967*GDD18
   GPC2 = 0.882*GDD0 + 0.015*GTD + 0.420*GPPT_SM - 0.215*GDD18
 !
-!  ---- Coefficient-file resolution ------------------------------------------
-!  Precedence: a mapped DGDRIVER code selects the file (keyword self-contained);
-!  otherwise fall back to FVS_GREGDG_COEF (the original env-var behaviour).
-CPATH = ' '
-IF (LKWSEL) THEN
-  !  Resolve config dir robustly: FVS_GREG_CONFDIR if set, else a known abs dir.
-  CALL GETENV('FVS_GREG_CONFDIR', CDIR)
-  IF (CDIR.EQ.' ') CDIR = '/users/PUOM0008/crsfaaron/wt-dgdriver/config'
-  !  DGDRIVER code -> coefficient filename mapping.
-  IF (IDGDRV.EQ.1) THEN
-    CPATH = TRIM(CDIR)//'/greg_dg_coefficients.csv'        ! deployed
-  ELSE IF (IDGDRV.EQ.2) THEN
-    CPATH = TRIM(CDIR)//'/greg_dg_coefficients_refit.csv'  ! refit on our data
-  ENDIF
-  WRITE(JOSTND,*) 'GREGDG: DGDRIVER code ', IDGDRV, ' selected ', TRIM(CPATH)
-ELSE
-  CALL GETENV('FVS_GREGDG_COEF', CPATH)
-ENDIF
+CALL GETENV('FVS_GREGDG_COEF', CPATH)
 IF (CPATH.EQ.' ') THEN
-  WRITE(JOSTND,*) 'GREGDG: no coefficient file (DGDRIVER unset/unmapped and ', &
-       'FVS_GREGDG_COEF empty); NOT enabled.'
+  WRITE(JOSTND,*) 'GREGDG: FVS_GREGDG set but FVS_GREGDG_COEF empty; NOT enabled.'
   LGREGDG = .FALSE.; RETURN
 ENDIF
 U = 68
@@ -101,33 +61,20 @@ SUBROUTINE GREGLOADDG
 READ(U,'(A)',IOSTAT=IOS) LINE                ! header
 NG = 0
 10 CONTINUE
-  READ(U,'(A)',IOSTAT=IOS) LINE
+  READ(U,*,IOSTAT=IOS) IFIA, NN, B0, B1, B2, B3, B4, B5, B6
   IF (IOS.NE.0) GO TO 20
-  IF (LINE.EQ.' ') GO TO 10
   IF (NG.GE.MXG) GO TO 20
-  !  Try 10-field record (SPCD,n,B0..B6,DBHMAX); if that fails, fall back to the
-  !  legacy 9-field record and leave DBHMAX at the sentinel (decel==1).
-  BMX = DGMAX_NONE
-  READ(LINE,*,IOSTAT=IOS) IFIA, NN, B0, B1, B2, B3, B4, B5, B6, BMX
-  IF (IOS.NE.0) THEN
-    BMX = DGMAX_NONE
-    READ(LINE,*,IOSTAT=IOS) IFIA, NN, B0, B1, B2, B3, B4, B5, B6
-    IF (IOS.NE.0) GO TO 10
-  ENDIF
-  IF (BMX.LE.0.0) BMX = DGMAX_NONE
   NG = NG + 1
   GSPCD(NG) = IFIA
   TB(NG,1)=B0; TB(NG,2)=B1; TB(NG,3)=B2; TB(NG,4)=B3
   TB(NG,5)=B4; TB(NG,6)=B5; TB(NG,7)=B6
-  TBMAX(NG)=BMX
   GO TO 10
 20 CONTINUE
 CLOSE(U)
 NGREGDG = NG
 !
 DO ISPC=1,MAXSP
   GHAVE_DG(ISPC) = .FALSE.
-  GDGMAX(ISPC) = DGMAX_NONE
   IFIA = -1
   IF (FIAJSP(ISPC).NE.' ') THEN
     READ(FIAJSP(ISPC),*,IOSTAT=IOS) IFIA
@@ -138,7 +85,6 @@ SUBROUTINE GREGLOADDG
       IF (GSPCD(J).EQ.IFIA) THEN
         GDG(ISPC,1)=TB(J,1); GDG(ISPC,2)=TB(J,2); GDG(ISPC,3)=TB(J,3)
         GDG(ISPC,4)=TB(J,4); GDG(ISPC,5)=TB(J,5); GDG(ISPC,6)=TB(J,6); GDG(ISPC,7)=TB(J,7)
-        GDGMAX(ISPC) = TBMAX(J)
         GHAVE_DG(ISPC) = .TRUE.
         EXIT
       ENDIF
@@ -147,7 +93,6 @@ SUBROUTINE GREGLOADDG
 ENDDO
 WRITE(JOSTND,*) 'GREGDG enabled: ', NG, ' species; DD0=', GDD0, ' TD=', GTD, &
        '  PPT_SM=', GPPT_SM, ' DD18=', GDD18, '  PC1=', GPC1, ' PC2=', GPC2
-IF (IDGDRV.GE.0) WRITE(JOSTND,*) 'GREGDG keyword DGDRIVER code=', IDGDRV
 RETURN
 END
 
@@ -159,7 +104,6 @@ SUBROUTINE GREGDGV(ISPC, DBH, CR, HT, BAL, G)
 INCLUDE 'GREGMC.f90'
 INTEGER ISPC
 REAL DBH, CR, HT, BAL, G, Z, CRC, HTC, BALC, ARGNUM, ARGDEN
-REAL DECEL, DBHMX, XR
 CRC = CR;  IF (CRC.LT.1.0E-4) CRC = 1.0E-4
 HTC = HT;  IF (HTC.LT.0.0) HTC = 0.0
 BALC = BAL; IF (BALC.LT.0.0) BALC = 0.0
@@ -171,34 +115,5 @@ SUBROUTINE GREGDGV(ISPC, DBH, CR, HT, BAL, G)
 IF (Z.GT.5.0)   Z = 5.0
 IF (Z.LT.-30.0) Z = -30.0
 G = EXP(Z); IF (G.LT.0.0) G = 0.0
-!  ---- Size-based deceleration (COR-style plateau) --------------------------
-!  Native NE-TWIGS DG plateaus via a size calibration this hook otherwise omits,
-!  so an unbounded compounding loop runs QMD away over multi-century horizons.
-!  Multiply the annual increment by a logistic that is ~1 until DBH nears the
-!  per-species maximum diameter GDGMAX (from the DBHMAX coef column) and ramps
-!  to 0 as DBH -> GDGMAX. Onset ~85% of max, half-width 4% of max, so trees far
-!  below max (all short remeasurement intervals) are essentially untouched.
-!  Sentinel GDGMAX (no DBHMAX column) leaves DECEL=1 -> exact old behaviour.
-DBHMX = GDGMAX(ISPC)
-IF (DBHMX .LT. 1.0E5) THEN
-  XR = (DBH - 0.85*DBHMX) / (0.04*DBHMX)
-  IF (XR .GT. 30.0) THEN
-    DECEL = 0.0
-  ELSE IF (XR .LT. -30.0) THEN
-    DECEL = 1.0
-  ELSE
-    DECEL = 1.0 / (1.0 + EXP(XR))
-  ENDIF
-  IF (DECEL .LT. 0.0) DECEL = 0.0
-  IF (DECEL .GT. 1.0) DECEL = 1.0
-  G = G * DECEL
-ENDIF
 RETURN
 END
-
-BLOCK DATA GREGKWBD
-!  Initialise keyword-selected driver codes to -1 (unset) before any keyword runs.
-IMPLICIT NONE
-INCLUDE 'GREGKW.f90'
-DATA IDGDRV /-1/, IMORTDRV /-1/
-END

src-converted/common/GREGKW.f90

@@ -1,6 +1,12 @@
-!  GREGKW -- keyword-selected Greg codes (DGDRIVER / MORTDRVR keywords).
-!    IDGDRV   -- DG site driver code, -1 unset; 0 none,1 elev,2 cspi,3 bgi,4 esi,5 emt
-!    IMORTDRV -- mort coefficient TIER code (NOT a driver family; additive
-!                gompit tiers), -1 unset; 0 crown-only,1 size,2 size+BGI
-INTEGER IDGDRV, IMORTDRV
-COMMON /GREGKW/ IDGDRV, IMORTDRV
+!  GREGKW -- keyword-layer control codes for Greg Johnson hook family.
+!  INCLUDE in any routine that needs to read or set hook-driver codes.
+!  Initialised to -1 (unset) by BLOCK DATA GREGKWBD in greghghg.f90;
+!  set by keyword handlers in initre.f90 (DGDRIVER, HGDRIVER, CROWNDRIVER).
+!
+!  IDGDRV   -- DG site driver code, -1 unset; 0 none,1 elev,2 cspi,3 bgi,4 esi,5 emt
+!  IMORTDRV -- mortality coefficient TIER code, -1 unset; 0 crown,1 size,2 size+BGI
+!  IHGDRV   -- HG arm code, -1 unset; 1 Bayes-calibrated regional,2 Greg ORGANON CONUS,
+!              3 CONUS species-dependent,4 CONUS species-independent
+!  ICRWDRV  -- crown recession driver code, -1 unset; 1 Greg CONUS recession model
+INTEGER IDGDRV, IMORTDRV, IHGDRV, ICRWDRV
+COMMON /GREGKW/ IDGDRV, IMORTDRV, IHGDRV, ICRWDRV

src-converted/common/GREGMC.f90

@@ -1,24 +1,55 @@
 !----------
-!  GREGMC -- Greg Johnson native diameter-growth (and future HG) substitution state.
-!  INCLUDE after PRGPRM.f90 so MAXSP/MAXTRE are defined. Mirrors GOMPMC.f90.
+!  GREGMC -- Greg Johnson DG, HG, and crown-recession substitution state.
+!  INCLUDE after PRGPRM.f90 so MAXSP is defined. Mirrors GOMPMC.f90.
+!
+!  ==== DG variables (PC-form, 9-column CSV B0..B6 + DBHMAX) ====
 !
 !    LGREGDG  -- master switch; .TRUE. when Greg DG substitution is active (set by
-!                GREGLOADDG from the FVS_GREGDG env var).
+!                GREGLOADDG from the FVS_GREGDG env var or DGDRIVER keyword).
 !    NGREGDG  -- number of FVS species with a fitted Greg DG row.
 !    GDG      -- per-FVS-species DG coefficient table, GDG(ISPC,1:7) = B0..B6,
 !                resolved from the FIA-keyed coefficient file via FIAJSP.
-!    GHAVE_DG -- .TRUE. for FVS species that carry a fitted Greg DG row; Greg DG
-!                owns their diameter increment. Unfit species keep native NE-TWIGS.
-!    GEMT/GTD/GELEV -- per-stand climate (extreme min temp, temp diff) and elevation
-!                (feet), resolved once per stand by GREGLOADDG.
-!    GDGMAX   -- per-FVS-species maximum diameter (in) for the size-based DG
-!                deceleration in GREGDGV. Read from an optional 10th column
-!                (DBHMAX) of the coefficient file; if the column is absent the
-!                array is left at a large sentinel so decel==1 (old behaviour).
+!    GDGMAX   -- per-species maximum diameter (in) for size-based DG deceleration.
+!    GHAVE_DG -- .TRUE. for FVS species that carry a fitted Greg DG row.
+!    GDD0/GTD/GPPT_SM/GDD18 -- per-stand ClimateNA inputs for PC-form DG.
+!    GPC1/GPC2 -- principal component climate scalars derived from the above.
+!
+!  ==== HG variables (Chapman-Richards ORGANON form, 11-column CSV) ====
+!
+!    LGREGHG   -- master switch for Greg HG substitution.
+!    NGREGHG   -- number of FVS species with a fitted Greg HG row.
+!    GHG       -- GHG(ISPC,1:9): (1)=MX asymptote, (2:9)=B1..B8.
+!    GHAVE_HG  -- .TRUE. for species that carry a fitted Greg HG row.
+!    GMCW      -- GMCW(ISPC,1:4): (1)=FORM code, (2:4)=A,B,C for MCW quadratic.
+!    GHAVE_MCW -- .TRUE. for species that carry a fitted MCW row.
+!    GEMT/GTD/GELEV -- per-stand climate for HG (EMT deg C, TD deg C, elevation ft).
+!                      GTD is shared with DG (same FVS_GREG_TD env var).
+!
+!  ==== Crown recession variables (3-param per-species model) ====
+!
+!    LGREG_CRW  -- master switch for Greg crown recession substitution.
+!    NGREG_CRW  -- number of FVS species with a fitted crown recession row.
+!    GCRW       -- GCRW(ISPC,1:3) = b0, b1, b2.
+!    GHAVE_CRW  -- .TRUE. for species that carry a fitted crown recession row.
 !----------
+! -- DG --
 INTEGER NGREGDG
 REAL    GDG(MAXSP,7), GDGMAX(MAXSP), GDD0, GTD, GPPT_SM, GDD18, GPC1, GPC2
 LOGICAL LGREGDG, GHAVE_DG(MAXSP)
 COMMON /GREGMR/ GDG, GDGMAX, GDD0, GTD, GPPT_SM, GDD18, GPC1, GPC2
 COMMON /GREGMI/ NGREGDG
 COMMON /GREGML/ LGREGDG, GHAVE_DG
+! -- HG --
+INTEGER NGREGHG
+REAL    GHG(MAXSP,9), GMCW(MAXSP,4), GEMT, GELEV
+LOGICAL LGREGHG, GHAVE_HG(MAXSP), GHAVE_MCW(MAXSP)
+COMMON /GREGHR/ GHG, GMCW, GEMT, GELEV
+COMMON /GREGHC/ NGREGHG
+COMMON /GREGHL/ LGREGHG, GHAVE_HG, GHAVE_MCW
+! -- Crown recession --
+INTEGER NGREG_CRW
+REAL    GCRW(MAXSP,3)
+LOGICAL LGREG_CRW, GHAVE_CRW(MAXSP)
+COMMON /GREGCR/ GCRW
+COMMON /GREGCC/ NGREG_CRW
+COMMON /GREGCL/ LGREG_CRW, GHAVE_CRW