C C C C C A LINKED FOREST PRODUCTIVITY-SOIL WATER, CARBON, AND NITROGEN C MODEL WRITTEN BY JOHN PASTOR AND W.M. POST, ENVIRONMENTAL C SCIENCES DIVISION, OAK RIDGE NATIONAL LABORATORY, OAK RIDGE, C TENNESSEE 37831, USA. C WRITTEN IN FORTRAN TO RUN ON THE ORNL IBM 3033, C THE MODEL IS BASED ON THE PREVIOUS MODELS OF C BOTKIN ET AL. (1972. J. ECOL. 60:849-872) AND SHUGART AND WEST C (1977. J. ENVIRON. MGMT. 5: 161-179). FULL DOCUMENTATION IS C GIVEN IN PASTOR AND POST (1985. ORNL/TM - 9519). C C THE DECOMPOSITION SUBROUTINE HAS BEEN MODIFIED TO INCLUDE C THREE CLASSES OF ORGANIC MATTER SUCH THAT THE INFLUENCE C OF EACH SPECIES EXTENDS FARTHER INTO THE DECAY PROCESS. C THE OUTPUT SUBROUTINE HAS BEEN MODIFIED TO GIVE BIOMASS C OF UP TO 10 SPECIES SPECIFIED IN INPUT. C C JOHN PASTOR. APRIL 1986. UNIV. OF MINNESOTA. C C MAIN PROGRAM C C C THE MAIN PROGRAM ESTABLISHES ALL COMMON BLOCKS, C THE SEEDS FOR THE RANDOM NUMBER GENERATOR, C AND CALLS SUBROUTINES IN ORDER OF EXECUTION. C ALL CALCULATIONS OF BIOLOGICAL INTEREST ARE DONE IN THE C SUBROUTINES. C COMMON/WATER/T(12),VT(12),RT(12),R(12),VR(12),FC,DRY,BGS,EGS,PLAT, 1FJ,AET COMMON/FOREST/NTREES(100),DBH(1500),IAGE(1500),KSPRT(100), > NEWTR(100),SUMLA(700),NEW(100),SWTCH(5) COMMON/PARAM/AAA(100,6),DMAX(100),DMIN(100),B3(100),B2(100), > ITOL(100),AGEMX(100),G(100),SPRTND(100),SPRTMN(100), 1SPRTMX(100),SWITCH(100,5),MPLANT(100),D3(100),FROST(100), 1TL(100),CM1(100),CM2(100),CM3(100),CM4(100),CM5(100), 2FWT(100),SLTA(100),SLTB(100),RTST(100),FRT(100) COMMON/CONST/NSPEC,DEGD COMMON/DEAD/NOGRO(1500),NTEMP(1500) COMMON/PROD/AWP(1500) COMMON/COUNT/NTOT,NYEAR,KPRNT,NMAX,KLAST,NWRITE,KWRITE COMMON/TEMP/DTEMP(1500),ITEMP(1500) COMMON/SEED/USEED(15) COMMON/INTERP/IPOLAT,X(10) COMMON/LINEAR/TSAV(45,12),VTSAV(45,12),RSAV(45,12),VRSAV(45,12) COMMON/DCMP/AVAILN,TYL(20),C(100,15),FDAT(20,10),FF(20,3), 1 NCOHRT,HCN,BASESC,BASESN,SCO2,TYLN COMMON/GMLT/SMGF(100),SNGF(100),DEGDGF(100) COMMON/SPECIE/SPEC(10),BMSPEC INTEGER SPEC,BMSPEC INTEGER USEED LOGICAL SWITCH,SWTCH C..... C.....SEEDS FOR RANDOM NUMBER GENERATOR C..... C.....USEED(1) - KILL- AGE DEPENDENT MORTALITY C.....USEED(2) - KILL- SLOW GROWTH MORTALITY C.....USEED(3) - NOT CURRENTLY USED C.....USEED(4) - BIRTH- SELECT NUMBER OF TREES TO SPROUT C.....USEED(5) - NOT CURRENTLY USED C.....USEED(6) - NOT CURRENTLY USED C.....USEED(7) - BIRTH- USED IN SMALL MAMMAL SWITCH C.....USEED(8) - NOT CURRENTLY USED C.....USEED(9) - BIRTH- DETERMINES THE NUMBER OF SEEDLINGS TO PLANT C.....USEED(10) - NOT CURRENTLY USED C.....USEED(11) - BIRTH- USED TO CALCULATE DBH FOR SEEDLINGS C.....USEED(12) AND (13) - TEMPE- STOCHASTIC VARIATION OF TEMPERATURE C.....USEED(14) AND (15) - MOIST- STOCHASTIC VARIATION OF RAINFALL C..... C..... C.....INITIALIZE SEED FOR RANDOM NUMBER GENERATOR C..... USEED(1) = 75364 USEED(2) = 82625 USEED(3) = 79154 USEED(4) = 79324 USEED(5) = 31697 USEED(6) = 91917 USEED(7) = 89819 USEED(8) = 37517 USEED(9) = 17119 USEED(10) = 72641 USEED(11) = 53797 USEED(12) = 91712 USEED(13) = 73319 USEED(14) = 51291 USEED(15) = 92761 OPEN(UNIT=6,FILE='OUT.FIL',ACCESS='SEQUENTIAL',STATUS='UNKNOWN') C..... C.....IPLOT - CURRENT PLOT C.....KYR - CURRENT YEAR C.....NYEAR - NUMBER OF YEARS SIMULATED C.....KLAST - NUMBER OF PLOTS SIMULATED C..... IPLOT = 0 C..... C.....READ INPUT DATA C..... CALL INPUT C..... C.....BEGIN PLOT LOOP C..... DO 60 K=1,KLAST C..... C.....INITIALIZE ARRAYS TO START ON BARE PLOT C..... CALL PLOTIN(IPLOT) KYR=0 C..... C.....ANNUAL LOOP WITHIN EACH PLOT C..... DO 50 I=1,NYEAR KYR = I CALL TEMPE(DEGD,KYR) CALL MOIST(KYR) CALL DECOMP(KYR) CALL GMULT(KYR) CALL BIRTH(KYR) CALL GROW(KYR) CALL KILL(KYR) IF(KYR.EQ.1) CALL OUTPUT(KYR,IPLOT) IF(MOD(KYR,KPRNT) .EQ. 0) CALL OUTPUT(KYR,IPLOT) 50 CONTINUE 60 CONTINUE STOP END C C SUBROUTINE BIRTH C C C SUBROUTINE BIRTH CALCULATES SEEDLING AND SPROUT BIRTH BASED ON C SPECIES FECUNDITY, SEEDBED CONDITIONS, SUSCEPTIBILITY TO BROWSING, C AND THE DEGREE TO WHICH LIGHT, SOIL MOISTURE, AND DEGREE DAYS ARE C LESS THAN OPTIMUM FOR GROWTH. SOIL MOISTURE AND DEGREE DAY C MULTIPLIERS ARE SUPPLIED BY SUBROUTINE GMULT. C A SPECIES CAN HAVE SPROUTS IF AT LEAST ONE C TREE WITH DIAMETER BETWEEN SPRTMN AND SPRTMX DIED LAST C YEAR (KSPRT INCREMENTED BY 1 IN KILL). C RANDOM NUMBERS USED TO DETERMINE OCCURENCE OF BROWSING, NUMBERS C OF SEEDLINGS AND SPROUTS, AND DBH SUPPLIED BY URAND. C SUBROUTINE BIRTH(KYR) COMMON/WATER/T(12),VT(12),RT(12),R(12),VR(12),FC,DRY,BGS,EGS,PLAT, 1FJ,AET COMMON/FOREST/NTREES(100),DBH(1500),IAGE(1500),KSPRT(100), > NEWTR(100),SUMLA(700),NEW(100),SWTCH(5) COMMON/PARAM/AAA(100,6),DMAX(100),DMIN(100),B3(100),B2(100), > ITOL(100),AGEMX(100),G(100),SPRTND(100),SPRTMN(100), 1SPRTMX(100),SWITCH(100,5),MPLANT(100),D3(100),FROST(100), 2TL(100),CM1(100),CM2(100),CM3(100),CM4(100),CM5(100), 3FWT(100),SLTA(100),SLTB(100),RTST(100),FRT(100) COMMON/CONST/NSPEC,DEGD COMMON/DEAD/NOGRO(1500),NTEMP(1500) COMMON/COUNT/NTOT,NYEAR,KPRNT,NMAX,KLAST,NWRITE,KWRITE COMMON/TEMP/DTEMP(1500),ITEMP(1500) COMMON/SEED/USEED(15) COMMON/GMLT/SMGF(100),SNGF(100),DEGDGF(100) INTEGER USEED LOGICAL SWITCH,SWTCH C..... C.....INITIALIZE FOLIAGE BIOMASS (FOLW) AND FOLIAGE AREA (FOLA) C..... FOLW = 0. FOLA = 0. NL = 1 C..... C.....CALCULATE LEAF WEIGHT IN G/PLOT (FOLW; ABER ET AL. 1982. C.....FOR.SCI. 28:31-45) AND LEAF AREA INDEX (FOLA; SOLLINS ET AL. 1973. C.....EDFB-IBP-73-2). C..... DO 30 J=1,NSPEC IF (NTREES(J).EQ.0) GO TO 30 NU = NL+NTREES(J)-1 RET=FRT(J) DO 20 K=NL,NU AGE=IAGE(K) IF(AGE.LT.RET) RET=AGE FOLW=FOLW+(((SLTA(J)+SLTB(J)*DBH(K))/2.) & **2*3.14*FWT(J)*RET) FOLA=FOLA+1.9283295E-4*DBH(K)**2.129 20 CONTINUE NL = NL+NTREES(J) 30 CONTINUE C..... C.....CALCULATE AMOUNT OF LIGHT AT FOREST FLOOR (% FULL SUNLIGHT) C.....(ABER ET AL. 1982. FOREST SCI. 28:31-45; C.....PASTOR AND MCCLAUGHERTY, UNPUBL.) C..... AL=1.0*EXP(-FOLW/93750.) C..... C.....TOTAL NUMBER OF TREES IN STAND C..... NTOT = NL-1 C..... C.....DETERMINE WHICH SPECIES ARE ELIGIBLE FOR PLANTING THIS YEAR C..... C.....SWITCH 1 IS TRUE IF THE SPECIES REQUIRES LEAF LITTER FOR C.....SUCCESSFUL RECRUITMENT C.....SWITCH 2 IS TRUE IF THE SPECIES REQUIRES MINERAL SOIL C.....SWITCH 3 IS TRUE IF THE SPECIES RECRUITMENT IS REDUCED BY HOT YEAR C.....SWITCH 4 IS TRUE IF THE SPECIES IS A PREFERRED FOOD OF DEER C.....OR SMALL MAMMALS C.....SWITCH 5 REDUCES SEEDING RATE OF DESIRABLE MAST C..... DO 40 J=1,5 SWTCH(J) = .TRUE. 40 CONTINUE SWTCH(3) = .FALSE. C..... C.....SET SWITCHES BASED ON VALUE OF LF AREA, DEGD, AND RANDOM NUMBER C..... IF (FOLA.GE..1) SWTCH(1) = .FALSE. IF (FOLA.LE..2) SWTCH(2) = .FALSE. C..... C.....BROWSE - A RANDOM NUMBER SIMULATING THE OCCURENCE OF BROWSING C..... YFL=URAND(USEED(7)) BROWSE=YFL IF (BROWSE.GT..5) SWTCH(4) = .FALSE. IF(FOLA .LE. .05) SWTCH(5) = .FALSE. NW = 0 C..... C.....SELECT ELIGIBLE SPECIES C..... DO 60 J=1,NSPEC C..... C.....END RECRUITMENT OF ASPEN, PIN CHERRY AND MOST PINES IF AVAILABLE C.....LIGHT IS < 60% OF FULL SUNLIGHT AND RECRUITMENT OF PAPER C.....BIRCH AND WHITE PINE IF AVAILABLE LIGHT IS < 30% OF FULL SUNLIGHT C..... IF(AL.LT..60.AND.J.EQ.20) GO TO 60 IF(AL.LT..30.AND.J.EQ.25) GO TO 60 IF(AL.LT..60.AND.J.EQ.58) GO TO 60 IF(AL.LT..60.AND.J.EQ.59) GO TO 60 IF(AL.LT..60.AND.J.EQ.60) GO TO 60 IF(AL.LT..60.AND.J.EQ.30) GO TO 60 IF(AL.LT..60.AND.J.EQ.31) GO TO 60 IF(AL.LT..60.AND.J.EQ.32) GO TO 60 IF(AL.LT..30.AND.J.EQ.33) GO TO 60 IF(AL.LT..60.AND.J.EQ.34) GO TO 60 IF(AL.LT..60.AND.J.EQ.35) GO TO 60 IF(AL.LT..30.AND.J.EQ.55) GO TO 60 DO 50 K=3,5 IF (SWITCH(J,K).AND.SWTCH(K)) GO TO 60 50 CONTINUE C..... C.....ALLOW ONLY THOSE SPECIES WHOSE DEGREE DAY TOLERANCES SPAN C.....THE SIMULATED DEGREE DAYS THIS YEAR TO BE ELIGIBLE FOR SEEDING C..... IF(DEGD .LE. DMIN(J) .OR. DEGD .GE. DMAX(J)) GO TO 60 C..... C.....ALLOW ONLY THOSE SPECIES WHOSE FROST TOLERANCE IS LESS THAN THE C.....JANUARY MEAN TEMPERATURE TO BE ELIGIBLE FOR SEEDING C..... IF(FROST(J) .GT. RT(1)) GO TO 60 C..... C.....PLACE ELIGIBLE SPECIES NUMBERS IN ARRAY NEWTR C..... NW = NW+1 NEWTR(NW) = J 60 CONTINUE C..... C.....CHECK TO SEE IF THERE ARE ANY NEW TREES C..... IF (NW.EQ.0) GO TO 145 C..... C.....PLACE IAGE, DBH, AND NOGRO INTO TEMPORARY ARRAYS C..... DO 70 I=1,NTOT ITEMP(I) = IAGE(I) DTEMP(I) = DBH(I) NTEMP(I) = NOGRO(I) 70 CONTINUE C..... C.....BEGIN MAIN LOOP FOR PLANTING C..... DO 140 K=1,NW NSP=NEWTR(K) C..... C.....CALCULATE SEEDLING LIGHT MULTIPLIERS C.....(SHUGART AND WEST. 1977. J. ENV. MGMT. 5:161-179) C.....IF ITOL IS LESS THAN 2, SPECIES IS SHADE TOLERANT C..... SLITE=1.5*(1.-EXP(-1.136*(AL-.08))) IF(ITOL(NSP).LT.2) SLITE=1.0-EXP(-4.64*(AL-.05)) IF(SLITE.LE.0.) SLITE=0. C..... C.....REDUCE MAXIMUM NUMBER OF SEEDLINGS TO THE EXTENT THAT LIGHT, C.....SOIL MOISTURE, AND DEGREE DAYS ARE LESS THAN OPTIMUM FOR C.....GROWTH OF EACH SPECIES C..... YFL=URAND(USEED(9)) NPLANT=MPLANT(NSP)*SLITE*SMGF(NSP)*DEGDGF(NSP)*YFL C..... C.....SEE IF ANY STUMPS OF THIS SPECIES ARE AVAILABLE FOR SPROUTING C..... IF(KSPRT(NSP).LE.0) GO TO 75 C..... C.....SEE IF SPECIES CAN STUMP SPROUT C..... IF(SPRTND(NSP).LE.0) GO TO 75 C..... C.....IF AVAILABLE LIGHT IS GREATER THAN 50% OF FULL SUNLIGHT, C.....DETERMINE NUMBER OF STUMP SPROUTS AND ADD TO NPLANT C..... YFL=URAND(USEED(4)) IF(AL.GE..50) & NPLANT=NPLANT+(SPRTND(NSP)*SLITE*SMGF(NSP)*DEGDGF(NSP)* & KSPRT(NSP)*YFL) 75 CONTINUE NSUM = 0 DO 80 I=1,NSP 80 NSUM = NSUM+NTREES(I) C..... C.....PLANT SEEDLINGS AND SPROUTS C..... NL = NSUM+1 NUP = NTOT IF(NPLANT.EQ.0) GO TO 140 DO 90 J=1,NPLANT NTOT = NTOT+1 IF (NTOT.GT.1500) CALL ERR6 NSUM = NSUM+1 NTREES(NSP) = NTREES(NSP)+1 ITEMP(NSUM) = 0 C..... C.....CALCULATE DBH FOR NEW TREES. DBH = 1.42CM +/- RANDOM AMOUNT C..... SIZE=1.27 YFL=URAND(USEED(11)) DTEMP(NSUM) = SIZE+.30*(1.0-YFL)**3 NTEMP(NSUM) = 0 90 CONTINUE IF (NL.GT.NUP) GO TO 110 N1 = NSUM+1 DO 100 L=NL,NUP DTEMP(N1) = DBH(L) ITEMP(N1) = IAGE(L) NTEMP(N1) = NOGRO(L) N1 = N1+1 100 CONTINUE C..... C.....REINITIALIZE ORIGINAL DBH AND AGE ARRAYS - INCLUDING NEW TREES C..... C..... 110 DO 120 I=1,NTOT IAGE(I) = ITEMP(I) DBH(I) = DTEMP(I) NOGRO(I) = NTEMP(I) 120 CONTINUE 140 CONTINUE 145 CONTINUE C..... C.....INCREMENT AGES BY ONE YEAR C..... DO 150 I=1,NTOT IAGE(I) = IAGE(I)+1 150 CONTINUE C..... C.....REINITIALIZE ARRAY KSPRT C..... DO 160 I=1,NSPEC KSPRT(I)=0 160 CONTINUE RETURN END C C SUBROUTINE DECOMP C C SUBROUTINE DECOMP CALCULATES CARBON AND NITROGEN FLOWS THROUGH C SOIL. AVAILABLE N (AVAILN) IS USED IN GMULT TO CALCULATE C SOIL NITROGEN GROWTH MULTIPLIERS. AET IS FED IN FROM MOIST. C THIS YEAR'S LEAF, TWIG, ROOT, AND WOOD LITTER IS FED IN FROM KILL C (ARRAY TYL). THE SIMULATION STARTS ON BARE GROUND (ONLY HUMUS C PRESENT. BASESC AND BASESN ARE STARTING HUMUS WEIGHT AND N C CONTENTS READ IN INPUT). THREE TYPES OF SOIL ORGANIC MATTER ARE C RECOGNIZED: COHORTS EITHER IMMOBILIZING OR RAPIDLY MINERALIZING C NITROGEN AND A HOMOGENOUS HUMUS POOL SLOWLY MINERALIZING N. C C SUBROUTINE DECOMP(KYR) COMMON/WATER/T(12),VT(12),RT(12),R(12),VR(12),FC,DRY, 1BGS,EGS,PLAT,FJ,AET COMMON/DCMP/AVAILN,TYL(20),C(100,15),FDAT(20,10),FF(20,3), 1NCOHRT,HCN,BASESC,BASESN,SCO2,TYLN C..... C.....INITIALIZATION C.....ALL VARIABLES IN T/HA C.....FCO2 - CO2-C FROM LITTER IMMOBILIZING N C.....HCO2 - CO2-C FROM HUMUS C.....SCO2 - TOTAL SOIL CO2-C C.....FFW - NEW WELL DECAYED WOOD COHORT C.....TNIMOB - TOTAL N IMMOBILIZATION C.....TNMIN - TOTAL N MINERALIZATION C.....AVAILN - AVAILABLE N C.....ASHFREE - ASHFREE WEIGHT OF LEAF LITTER C.....TYLC - LEAF LITTER C CONTENT C.....TYLN - LEAF LITTER N CONTENT C.....TYLNC - LEAF LITTER N:C C.....FF - WEIGHTS AND N CONTENT OF FOREST FLOOR BY LITTER TYPE C.....AET - ACTUAL EVAPOTRANSPIRATION IN MM (FROM MOIST) C..... C..... FCO2=0. HCO2=0. SCO2=0. FFW=0. TNIMOB=0. AVAILN=0. FNMIN=0. HNMIN=0. TNMIN=0. TYLN=0. DO 9 I=1,20 DO 9 J=1,3 9 FF(I,J)=0. C..... C.....CALCULATE LITTER N C..... DO 8 I=1,12 8 TYLN=TYLN+(TYL(I)*FDAT(I,2)) C..... C.....CALCULATE AET MULTIPLIER C..... XAET=AET IF(XAET.GT.600.) XAET=600. AETM = (-1.*XAET)/(-1200. + XAET) C..... C......CREATE NEW COHORTS C.....THE FIRST ROW OF THE C ARRAY HOLDS DATA ON HUMUS. C.....ROWS 2 THROUGH NCOHRT HOLD DATA ON LITTER COHORTS C.....C(I,1) - WEIGHT (T/HA) C.....C(I,2) - N CONTENT (T/HA) C.....C(I,3) AND C(I,4) - N CHANGE PARAMETERS C.....C(I,5) - LITTER TYPE. 1 THROUGH 12 ARE LEAVES, 13 IS ROOT, C..... 14 AND 15 ARE WOOD, 16 IS TWIGS, 17 IS WELL DECAYED WOOD, C..... AND 18 IS HUMUS C.....C(I,6) - DESTINATION WHEN TRANSFERRED. 1=HUMUS, C..... 2=WELL DECAYED WOOD C.....C(I,7) - CURRENT % LIGNIN C.....C(I,8) AND C(I,9) - LIGNIN DECAY PARAMETERS C.....C(I,10) - ORIGINAL WEIGHT C.....C(I,11) - CURRENT % N C.....C(I,12) - FRACTION OF ORIGINAL WEIGHT WHICH WILL BECOME C..... HUMUS OR WELL DECAYED WOOD. WHEN THIS FRACTION C..... IS REACHED, THE COHORT IS TRANSFERRED TO THE C..... DESTINATION SPECIFIED BY C(I,6). THIS FRACTION C..... IS BASED ON THE LIGNIN CONTENT OF LEAVES, TAKEN C..... FROM DEHAAN (1977) SOIL ORGANIC MATTER STUDIES, C..... VOL.1, IEA. C..... DO 3 I=1,16 IF(TYL(I).EQ.0.) GO TO 3 NCOHRT=NCOHRT+1 IF(NCOHRT.GT.100) CALL ERR10 C(NCOHRT,1)=TYL(I)*FDAT(I,10) C(NCOHRT,2)=TYL(I)*FDAT(I,2) DO 2 J=3,9 C(NCOHRT,J)=FDAT(I,J) 2 CONTINUE C(NCOHRT,10)=TYL(I)*FDAT(I,10) C(NCOHRT,11)=FDAT(I,2) C(NCOHRT,12)=FDAT(I,7)*1.7039 + 0.0955 IF(C(NCOHRT,5).EQ.14) C(NCOHRT,12)=.30 IF(C(NCOHRT,5).EQ.15) C(NCOHRT,12)=.30 IF(C(NCOHRT,5).EQ.16) C(NCOHRT,12)=.30 3 CONTINUE C..... C......CALCULATE DECAY MULTIPLIER, SIMULATING EFFECT OF GAPS ON DECAY C..... TYLL=TYL(17) CCLL=1.54+.0457*(FC-DRY) C CCLL=1.0 IF(TYLL.GT.CCLL)TYLL=CCLL DECMLT=1.0+(-.50+.075*(FC-DRY))*(1.0-TYLL/CCLL) C..... C.....BYPASS FOREST FLOOR COHORT CALCULATIONS IF THERE IS NO FLOOR C..... IF(NCOHRT.EQ.1) GO TO 15 C..... C.....LOOP TO CALCULATE LITTER DECAY, N IMMOBILIZATION, LIGNIN DECAY, C.....AND LITTER CO2 EVOLUTION C..... DO 4 I=2,NCOHRT C..... C.....CALCULATE % WT LOSS BASED ON AET AND LIGNIN:N RATIO C..... PWTLOS=(.9804+.09352*AET)-((-.4956+.00193*AET)*(C(I,7)/C(I,11))) PWTLOS=(DECMLT*PWTLOS)/100. IF(PWTLOS.GT..99) PWTLOS=.99 C..... C.....WT LOSS OF LARGE WOOD (DBH>10CM) IS 3%; C.....WT LOSS OF SMALL WOOD IS 10%; WT LOSS C.....OF WELL DECAYED WOOD IS 5%; WT LOSS C.....OF TWIGS IS LESS THAN 20% C..... LT=C(I,5) IF(LT.EQ.14) PWTLOS=.10 IF(LT.EQ.15) PWTLOS=.03 IF(LT.EQ.17) PWTLOS=.05 IF(LT.EQ.16.AND.PWTLOS.GT..20) PWTLOS=.20 C..... C.....CALCULATE ACTUAL WT LOSS (T/HA) C..... WTLOSS=PWTLOS*C(I,1) C..... C.....CALCULATE FRACTION OF ORGANIC MATTER REMAINING C..... POMR = (C(I,1)-WTLOSS)/C(I,10) C..... C.....FIND NEW N CONCENTRATION IN COHORT C..... C(I,11)=C(I,3)-C(I,4)*POMR C..... C.....RETAIN COHORT FOR ANOTHER YEAR OF DECAY IF FRACTION C.....REMAINING IS GREATER THAN FRACTION WHICH WILL BECOME HUMUS C.....OR WELL DECAYED WOOD C..... IF(POMR.GT.C(I,12)) GO TO 7 C..... C.....IF COHRT IS TO BE TRANSFERRED TO HUMUS, RECALCULATE WTLOSS C.....AND N CONCENTRATION SO THAT THE TRANSFER OCCURS AT THE C.....FRACTION SPECIFIED BY THE INITIAL LIGNIN CONCENTRATION C..... WTLOSS = C(I,1)-C(I,12)*C(I,10) C(I,11) = C(I,3)-C(I,4)*C(I,12) C..... C.....CALCULATE ABSOLUTE CHANGE IN N CONTENT C..... DELTAN = C(I,2)-C(I,11)*(C(I,1)-WTLOSS) IF(DELTAN.LT.0.) TNIMOB=TNIMOB-DELTAN IF(DELTAN.GE.0.) FNMIN=FNMIN+DELTAN C..... C.....TRANSFER COHORTS C..... IF(C(I,6).NE.1.) GO TO 5 C(1,1)=C(1,1)+C(I,1)-WTLOSS C(1,2)=C(1,2)+C(I,11)*(C(I,1)-WTLOSS) C(I,1)=0. GO TO 7 C..... C.....FFW - TEMPORARY VARIABLE ASSIGNED TO WELL DECAYED WOOD COHORT C..... 5 FFW=FFW+C(I,1)-WTLOSS C(I,1)=0. C..... C.....UPDATE COHORTS C..... 7 IF(C(I,1).EQ.0.) GO TO 14 C(I,1)=C(I,1)-WTLOSS C(I,2)=C(I,1)*C(I,11) C(I,7)=C(I,8)-C(I,9)*(C(I,1)/C(I,10)) C..... C.....CALCULATE LITTER COHORT CO2 EVOLUTION C..... 14 FCO2=FCO2+(WTLOSS*0.48) 4 CONTINUE C..... C.....THROUGHFALL IS 16% OF LEAF LITTER N (COLE AND RAPP.1982. C.....PP.341-410 IN D.E.REICHLE, ED. DYNAMIC PROPERTIES OF C.....FOREST ECOSYSTEMS. IBP 23. CAMBRIDGE UNIV. PRESS) C.....AND SUPPLIES SOME OF IMMOBILIZATION DEMANDS C..... TNIMOB=TNIMOB-.16*TYLN C..... C.....CALCULATE HUMUS N MINERALIZATION C..... 15 HNMIN = C(1,2)*0.035*DECMLT*AETM C..... C.....SUBTRACT MINERALIZED N FROM HUMUS N POOL C.....AND CALCULATE HUMUS CO2 C..... HNNEW=C(1,2)-HNMIN HOMNEW=C(1,1)*(HNNEW/C(1,2)) HCO2=(C(1,1)-HOMNEW)*0.48 C(1,1)=HOMNEW C(1,2)=HNNEW C..... C.....HCN - HUMUS C:N RATIO C..... HCN=(0.48*C(1,1))/C(1,2) C..... C.....ADD ADD HUMUS N MINERALIZATION TO COHORT N C.....MINERALIZATION TO GET TOTAL N MINERALIZATION C..... TNMIN = FNMIN + HNMIN C..... C.....SUBTRACT IMMOBILIZATION FROM TOTAL MOINERALIZATION C.....TO GET AVAILABLE N TO TREES CC..... AVAILN = TNMIN - TNIMOB C..... C.....CALCULATE TOTAL SOIL RESPIRATION C..... SCO2=FCO2+HCO2 C..... C.....REMOVE TRANSFERRED COHORTS C..... IX=0 DO 20 I=1,NCOHRT IF(C(I,1).EQ.0.) GO TO 16 DO 12 J=1,12 12 C(I-IX,J)=C(I,J) GO TO 20 16 IX=IX+1 20 CONTINUE NCOHRT=NCOHRT-IX C..... C.....CREATE NEW WELL DECAYED WOOD COHORT C..... IF(FFW.EQ.0.) GO TO 30 NCOHRT=NCOHRT+1 IF(NCOHRT.GT.1500) CALL ERR10 C(NCOHRT,1)=FFW C(NCOHRT,2)=FFW*FDAT(17,2) DO 22 J=3,9 22 C(NCOHRT,J)=FDAT(17,J) C(NCOHRT,10)=FFW C(NCOHRT,11)=FDAT(17,2) C(NCOHRT,12)=.50 C..... C.....CALCULATE TOTAL WT AND N CONTENT BY FOREST FLOOR COMPARTMENT C..... 30 DO 35 I=1,NCOHRT LT=C(I,5) FF(LT,1)=C(I,5) FF(LT,2)=FF(LT,2)+C(I,1) FF(LT,3)=FF(LT,3)+C(I,2) 35 CONTINUE FF(19,1)=19. DO 36 LT=1,12 FF(19,2)=FF(19,2)+FF(LT,2) FF(19,3)=FF(19,3)+FF(LT,3) 36 CONTINUE FF(19,2)=FF(19,2)+FF(18,2)+FF(13,2) FF(19,3)=FF(19,3)+FF(18,3)+FF(13,3) 40 RETURN END C C SUBROUTINE ERR C C C THIS SUBROUTINE IS A LIST OF ERRORS THAT MAY OCCUR IN THE C PROGRAM. THE ERRORS FOCUS ON VARIABLES USED AS DO LOOP COUNTERS C AND WHICH MAY EXCEED DIMENSIONS OF CERTAIN ARRAYS. C SUBROUTINE ERR ENTRY ERR1 WRITE(6,10) 10 FORMAT(' STATISTICAL ARRAYS IN SUBROUTINE OUTPUT CAN ONLY &HANDLE 100 PLOTS.',/' REDUCE KLAST OR REDIMENSION ARRAYS.') GO TO 110 ENTRY ERR2 WRITE(6,20) 20 FORMAT(' STATISTICAL ARRAYS IN SUBROUTINE OUTPUT CAN ONLY &HANDLE 70 LINES',/' OF OUTPUT. REDUCE NYEAR, INCREASE KPRNT, &OR REDIMENSION ARRAYS.') GO TO 110 ENTRY ERR3 WRITE(6,30) 30 FORMAT(' NUMBER OF SPECIES, NSPEC, HAS EXCEEDED 100 IN INPUT.',/ &' REDUCE NSPEC OR REDIMENSION ALL ARRAYS IN COMMON BLOCKS',/ &' PARAM AND GMLT AND APPROPRIATE ARRAYS IN ',/ &' COMMON BLOCK FOREST.') GO TO 110 ENTRY ERR4 WRITE(6,40) 40 FORMAT(' NUMBER OF LITTER COHORTS, NCOHRT, HAS EXCEEDED 100',/ &' IN SUBROUTINE INPUT. REDUCE NCOHRT OR REDIMENSION C ARRAY.') GO TO 110 ENTRY ERR5 WRITE(6,50) 50 FORMAT(' NUMBER OF TREES, NTOT1, HAS EXCEEDED 1500 IN KILL.') GO TO 110 ENTRY ERR6 WRITE(6,60) 60 FORMAT(' NUMBER OF TREES, NTOT, HAS EXCEEDED 1500 IN BIRTH.') GO TO 110 ENTRY ERR7 WRITE(6,70) 70 FORMAT(' NUMBER OF TREES, NTOT, HAS EXCEEDED 1500 IN OUTPUT.') GO TO 110 ENTRY ERR8 WRITE(6,80) 80 FORMAT(' NUMBER OF TREES, NTOT, HAS EXCEEDED 1500 IN GROW.') GO TO 110 ENTRY ERR9 WRITE (6,90) 90 FORMAT(' IHT HAS EXCEEDED 700 IN GROW. AT LEAST ONE TREE',/ &' IS GREATER THAN 70M TALL. REDIMENSION ARRAY SUMLA.') GO TO 110 ENTRY ERR10 WRITE(6,100) 100 FORMAT(' THE NUMBER OF LITTER COHORTS, NCOHRT, HAS EXCEEDED',/ &' 100 IN SUBROUTINE DECOMP. ') ENTRY ERR11 WRITE(6,120) 120 FORMAT(' NUMBER OF SPECIES FOR BIOMASS ANALYSIS, BMSPEC',/ &' HAS EXCEEDED 10') 110 STOP END C C SUBROUTINE GGNORD C C C SUBROUTINE GGNORD CALCULATES NORMALLY DISTRIBUTED RANDOM NUMBERS C SUPPLIED BY URAND. IT IS CALLED FROM SUBROUTINES TEMPE AND MOIST. C SUBROUTINE GGNORD(NSEED1,NSEED2,Z) DIMENSION Z(1) DATA PI2/0.62831853E01/ K = 0 C..... C.....FIND RANDOM NUMBERS C..... A1 = URAND(NSEED1) A2 = URAND(NSEED2) K = K+1 C..... C.....CALCULATE NORMALLY DISTRIBUTED RANDOM NUMBERS. C.....(EMSHOFF AND SISSON.1970.DESIGN AND USE OF COMPUTER C.....SIMULATION MODELS. MACMILLAN) C..... Z(K) = SQRT(-.2E01*ALOG(A1))*SIN(PI2*A2) K = K+1 Z(K) = SQRT(-0.2E01*ALOG(A1))*COS(PI2*A2) RETURN END C C SUBROUTINE GMULT C C C SUBROUTINE GMULT CALCULATES DEGREE DAY, SOIL MOISTURE, AND SOIL C NITROGEN MULTIPLIERS USED IN SUBROUTINES BIRTH AND GROW BASED C ON ON DEGD (SUPPLIED BY TEMPE), FJ (SUPPLIED BY MOIST), AND C AVAILN (SUPPLIED BY DECOMP), RESPECTIVELY. C SUBROUTINE GMULT(KYR) COMMON/WATER/T(12),VT(12),RT(12),R(12),VR(12),FC,DRY,BGS,EGS,PLAT, 1FJ,AET COMMON/PARAM/AAA(100,6),DMAX(100),DMIN(100),B3(100),B2(100), 1ITOL(100),AGEMX(100),G(100),SPRTND(100),SPRTMN(100), 2SPRTMX(100),SWITCH(100,5),MPLANT(100),D3(100),FROST(100), 3TL(100),CM1(100),CM2(100),CM3(100),CM4(100),CM5(100), 4FWT(100),SLTA(100),SLTB(100),RTST(100),FRT(100) COMMON/CONST/NSPEC,DEGD COMMON/DCMP/AVAILN,TYL(20),C(100,15),FDAT(20,10),FF(20,3), 1NCOHRT,HCN,BASESC,BASESN,SCO2,TYLN COMMON/GMLT/SMGF(100),SNGF(100),DEGDGF(100) LOGICAL SWITCH, SWTCH C..... C.....CALCULATE TOTAL NUMBER OF GROWING SEASON DEGREE DAYS C..... TGS = EGS-BGS + 1 C..... C.....CALCULATE RELATIVE NITROGEN AVAILABILITY, C.....(ABER ET AL. 1979. CAN.J.FOR.RES. 9:10-14) C..... AVAILN=AVAILN+.005 C AVAILN=.024 AVLMC = -170.0 + 4.0*(AVAILN*1000.) C..... C.....CALCULATE GROWTH MULTIPLIERS C..... DO 40 I=1,NSPEC C..... C.....CALCULATE SPECIES DEGREE DAY MULTIPLIERS C.....(BOTKIN ET AL. 1972. J. ECOL. 60:849-872) C..... DEGDGF(I)=4.0*(DEGD-DMIN(I))*(DMAX(I)-DEGD)/ &(DMAX(I)-DMIN(I))**2 IF(DEGDGF(I).LT.0.) DEGDGF(I)=0. IF(DEGDGF(I).EQ.0.) GO TO 40 C..... C.....CALCULATE SPECIES SOIL MOISTURE MULTIPLIERS C..... DROUT = D3(I)*TGS IF(DROUT.LT.FJ) DROUT=FJ SMGF(I)=SQRT((DROUT-FJ)/DROUT) IF(SMGF(I).EQ.0.) GO TO 40 C..... C.....CALCULATE SPECIES SOIL NITROGEN GROWTH MULTIPLIER C.....(MITCHELL AND CHANDLER. 1939. BLACK ROCK FOREST BULLETIN 11; C.....ABER ET AL. 1979.CAN.J.FOR.RES. 9:10-14). C.....CONN - %N IN GREEN LEAVES C..... CONN=CM1(I)*(1.-10.**((-1.*CM3(I))*(AVLMC+CM2(I)))) SNGF(I)=CM4(I)+CM5(I)*CONN IF(SNGF(I).LT.0.) SNGF(I)=0. 40 CONTINUE RETURN END C C SUBROUTINE GROW C C SUBROUTINE GROW CALCULATES DIAMETER GROWTH FOR EACH TREE C BY DECREASING MAXIMAL GROWTH TO THE EXTENT THAT THE MOST C LIMITING RESOURCE IS LESS THAN OPTIMAL. C SUBROUTINE GROW(KYR) COMMON/WATER/T(12),VT(12),RT(12),R(12),VR(12),FC,DRY,BGS,EGS,PLAT, 1FJ,AET COMMON/FOREST/NTREES(100),DBH(1500),IAGE(1500),KSPRT(100), > NEWTR(100),SUMLA(700),NEW(100),SWTCH(5) COMMON/PARAM/AAA(100,6),DMAX(100),DMIN(100),B3(100),B2(100), > ITOL(100),AGEMX(100),G(100),SPRTND(100),SPRTMN(100), 1SPRTMX(100),SWITCH(100,5),MPLANT(100),D3(100),FROST(100), 2TL(100),CM1(100),CM2(100),CM3(100),CM4(100),CM5(100), 3FWT(100),SLTA(100),SLTB(100),RTST(100),FRT(100) COMMON/CONST/NSPEC,DEGD COMMON/DEAD/NOGRO(1500),NTEMP(1500) COMMON/PROD/AWP(1500) COMMON/COUNT/NTOT,NYEAR,KPRNT,NMAX,KLAST,NWRITE,KWRITE COMMON/TEMP/DTEMP(1500),ITEMP(1500) COMMON/GMLT/SMGF(100),SNGF(100),DEGDGF(100) LOGICAL SWITCH,SWTCH C..... C.....INITIALIZE WOOD PRODUCTION C..... DO 6 I=1,1500 AWP(I) = 0.0 6 CONTINUE C..... C.....CALCULATE TOTAL NUMBER OF TREES C..... NTOT = 0 DO 10 I=1,NSPEC 10 NTOT = NTOT+NTREES(I) IF (NTOT.EQ.0) RETURN IF(NTOT.GT.1500) CALL ERR8 C..... C.....INITIALIZE CANOPY LEAF BIOMASS PROFILE C..... DO 20 I=1,700 20 SUMLA(I) = 0. C..... C.....LOOP FOR CALCULATING CANOPY PROFILE C..... NL = 1 DO 40 J=1,NSPEC IF (NTREES(J).EQ.0) GO TO 40 NU = NL+NTREES(J) -1 RET=FRT(J) DO 30 K=NL,NU AGE=IAGE(K) IF(AGE.LT.RET) RET=AGE C..... C..... C.....HEIGHT PROFILE IS CALCULATED IN .1 METER UNITS C.....(BOTKIN ET AL. 1972. J. ECOL. 60:849-872) C..... IHT = (B2(J)*DBH(K)-B3(J)*DBH(K)**2)/10.+1. IF (IHT.GT.700) CALL ERR9 C..... C.....CALCULATE AND SUM LEAF BIOMASS FOR TREES OF APPROXIMATELY C.....THE SAME HEIGHT (ABER ET AL. 1982. FOREST SCI. 28:31-45) C..... SUMLA(IHT) = SUMLA(IHT)+(((SLTA(J)+SLTB(J)*DBH(K))/2.) & **2*3.14*FWT(J)*RET) 30 CONTINUE NL = NL+NTREES(J) 40 CONTINUE C..... C.....CALCULATE CUMULATIVE LEAF BIOMASS DOWN THROUGH THE CANOPY C..... DO 50 J=1,699 J1 = 700-J SUMLA(J1) = SUMLA(J1)+SUMLA(J1+1) 50 CONTINUE C..... C.....MAIN LOOP FOR CALCULATING DIAMETER INCREMENT C..... NL = 1 DO 80 I=1,NSPEC IF (NTREES(I).EQ.0) GO TO 80 NU = NL+NTREES(I)-1 DO 70 J=NL,NU C..... C.....CALCULATE LEAF BIOMASS OF ALL TALLER TREES (SLAR) C..... HT = B2(I)*DBH(J)-B3(I)*DBH(J)**2 IHT = HT/10.+2. SLAR = SUMLA(IHT) C..... C.....CALCULATE AVAILABLE LIGHT TO THIS TREE (% FULL SUNLIGHT) C.....(ABER ET AL. 1982. FOREST SCI. 28: 31-45; C.....PASTOR AND MCCLAUGHERTY UNPUBL.) C..... AL = 1.0*EXP(-SLAR/93750.) C..... C.....CALCULATE AVAILABLE LIGHT MULTIPLIER IF TREE IS SHADE INTOLERANT C.....(SHUGART AND WEST. 1977. J. ENV. MGMT. 5:161-179) C..... IF(ITOL(I).LT.2) GO TO 58 ALGF=2.24*(1.0-EXP(-1.136*(AL-.08))) GO TO 59 C..... C.....CALCULATE AVAILABLE LIGHT MULTIPLIER IF TREE IS SHADE TOLERANT C.....(SHUGART AND WEST. 1977) C..... 58 ALGF=1.0-EXP(-4.64*(AL-.05)) 59 IF(ALGF.LT.0.) ALGF=0. C..... C.....CALCULATE MAXIMUM TREE VOLUME C.....(BOTKIN ET AL. 1972. J. ECOL. 60:849-872) C..... GR = (137.+.25*B2(I)**2/B3(I))*(0.5*B2(I)/B3(I)) C..... C.....CALCULATE DIAMETER INCREMENT UNDER OPTIMAL CONDITIONS C..... DNCMAX = G(I)*DBH(J)*(1.0-(137.*DBH(J)+B2(I)*DBH(J)**2-B3(I) * > DBH(J)**3)/GR)/(274.+3.0*B2(I)*DBH(J)-4.0*B3(I)*DBH(J)**2) C..... C.....CHOOSE SMALLEST GROWTH MULTIPLIER FOR THIS TREE C.....(LIEBIG'S LAW OF THE MINIMUM) C..... GF=AMIN1(ALGF,SMGF(I),SNGF(I),DEGDGF(I)) C..... C.....REDUCE DIAMETER INCREMENT TO THE EXTENT THAT CONDITIONS ARE C.....LESS THAN OPTIMUM FOR GROWTH C..... DINC=DNCMAX*GF C..... C.....CHECK INCREMENT LESS THAN MINIMUM REQUIRED FOR GROWTH. C.....IF DINC LESS THAN 1.0 MM OR 10% OF DNCMAX, OR IF C.....JANUARY TEMPERATURE LESS THAN FROST TOLERANCE, FLAG TREE C.....IN NOGRO C..... IF(DINC.LT..1.OR.FROST(I).GT.RT(1)) DINC=0. IF(DINC .GE. .10*DNCMAX) NOGRO(J) = 0 IF(DINC .LT. .10*DNCMAX) NOGRO(J) =NOGRO(J) -1 C..... C.....CALCULATE WOODY BIOMASS (KG) BEFORE INCREMENTING DIAMETER C.....(SOLLINS ET AL. 1973. EDFB-IBP-73-2) C..... AB1 = .1193*DBH(J)**2.393 C..... C.....INCREMENT DIAMETER C..... 60 DBH(J) = DBH(J)+DINC C..... C.....CALCULATE WOODY BIOMASS AFTER INCREMENTING DIAMETER C..... AB2 = .1193*DBH(J)**2.393 C..... C.....CALCULATE NET INCREASE IN WOODY BIOMASS (ABOVEGROUND C.....WOODY PRODUCTION IN KG) C..... AWP(J) = AB2-AB1 70 CONTINUE NL = NL+NTREES(I) 80 CONTINUE RETURN END C C SUBROUTINE INPUT C C C INPUT READS RUN PARAMETERS, LATITUDE, LONGITUDE, DAYS OF THE YEAR C THE GROWING SEASON BEGINS AND ENDS, SOIL FIELD MOISTURE CAPACITY C AND WILTING POINTS, MONTHLY TEMPERATURE, PRECIPITATION, AND THEIR C STND DEV, SPECIES PARAMETERS, DECOMPOSITION PARAMETERS, AND C STARTING HUMUS WEIGHT AND N CONTENT C SUBROUTINE INPUT COMMON/INTERP/IPOLAT,X(10) COMMON/WATER/T(12),VT(12),RT(12),R(12),VR(12),FC,DRY,BGS,EGS,PLAT, 1FJ,AET COMMON/PARAM/AAA(100,6),DMAX(100),DMIN(100),B3(100),B2(100), > ITOL(100),AGEMX(100),G(100),SPRTND(100),SPRTMN(100), 1SPRTMX(100),SWITCH(100,5),MPLANT(100),D3(100),FROST(100), 2TL(100),CM1(100),CM2(100),CM3(100),CM4(100),CM5(100), 3FWT(100),SLTA(100),SLTB(100),RTST(100),FRT(100) COMMON/COUNT/NTOT,NYEAR,KPRNT,NMAX,KLAST,NWRITE,KWRITE COMMON/CONST/NSPEC,DEGD COMMON/LINEAR/TSAV(45,12),VTSAV(45,12),RSAV(45,12),VRSAV(45,12) COMMON/DCMP/AVAILN,TYL(20),C(100,15),FDAT(20,10),FF(20,3), 1 NCOHRT,HCN,BASESC,BASESN,SCO2,TYLN COMMON/SPECIE/SPEC(10),BMSPEC INTEGER BMSPEC,SPEC LOGICAL SWITCH, SWTCH C..... C.....READ RUN PARAMETERS C..... C.....KPRNT - PRINT INTERVAL C.....KLAST - NUMBER OF PLOTS SIMULATED C.....NYEAR - NUMBER OF YEARS SIMULATED C.....NMAX - NUMBER OF TIMES OUTPUT IS CALLED PER PLOT C.....NWRITE - NUMBER OF TIMES OUTPUT IS CALLED PER RUN C.....KWRITE - COUNTS NUMBER OF TIMES OUTPUT IS CALLED. C..... WHEN KWRITE=NWRITE, MEANS AND CONFIDENCE INTERVALS C..... OF SELECTED VARIABLES ARE CALCULATED IN OUTPUT. C..... OPEN(UNIT=5,FILE='LINKAGES.DAT',ACCESS='SEQUENTIAL',STATUS='OLD') READ(5,1005) KPRNT,KLAST,NYEAR C OPEN(UNIT=9,FILE='OUT.FIL',ACCESS='SEQUENTIAL',STATUS='UNKNOWN') 1005 FORMAT(6X,I3,7X,I3,7X,I5) IF(KLAST.GT.100) CALL ERR1 NMAX=NYEAR/KPRNT + 1 C IF(NMAX.GT.70) CALL ERR2 NWRITE=NMAX*KLAST KWRITE=0 C..... C.....IPOLAT - NUMBER OF BREAK POINTS IN LINEAR INTERPOLATIONS C..... (EQUALS NUMBER OF ENTRIES IN ARRAY X) C..... READ(5,1010) IPOLAT 1010 FORMAT(7X,I2) C..... C.....ARRAY X CONTAINS YEARS IN WHICH CLIMATE CHANGES AND BETWEEN C.....WHICH LINEAR INTERPOLATIONS MUST BE MADE. THE FIRST ENTRY C.....MUST BE EQUAL TO ZERO AND THE LAST EQUAL TO NYEAR. C..... READ(5,1015) (X(I),I=1,10) 1015 FORMAT(10F7.0) C..... C.....READ LATITUDE, LONGITUDE, DAYS GROWING SEASON BEGINS AND ENDS, C.....SOIL FIELD MOISTURE CAPACITY AND WILTING POINT C..... READ(5,1020) PLAT,PLONG,BGS,EGS,FC,DRY TGS = EGS-BGS+1 1020 FORMAT(5X,F5.2,7X,F5.1,5X,F4.0,5X,F4.0,4X,F4.1,5X,F4.1) C..... C.....READ MONTHLY TEMPERATURE, STND DEV, PRECIPITATION, STND DEV C..... READ(5,1035) ((TSAV(J,K),K=1,12),J=1,IPOLAT) READ(5,1035) ((VTSAV(J,K),K=1,12),J=1,IPOLAT) READ(5,1035) ((RSAV(J,K),K=1,12),J=1,IPOLAT) READ(5,1035) ((VRSAV(J,K),K=1,12),J=1,IPOLAT) 1035 FORMAT(8X,12F6.0) C..... C.....WRITE CLIMATE DATA TO PRINTER C..... WRITE(6,1025) PLAT,PLONG,BGS,EGS 1025 FORMAT(' LATITUDE=',F5.1,' LONGITUDE=',F5.1, 1' GROWING SEASON: BEGINS DAY ',F5.1, 2' ENDS DAY ',F5.1) WRITE(6,7036) 7036 FORMAT(' '/13X,'J',5X,'F',5X,'M',5X,'A',5X,'M',5X,'J',5X,'J',5X, 1'A',5X,'S',5X,'O',5X,'N',5X,'D') WRITE(6,7037) (TSAV(1,K),K=1,12) 7037 FORMAT(' TEMP (C)',12F6.1) WRITE(6,7038) (VTSAV(1,K),K=1,12) 7038 FORMAT(' STND DEV',12F6.1) WRITE(6,7039) (RSAV(1,K),K=1,12) 7039 FORMAT(' PPT (CM)',12F6.1) WRITE(6,7040) (VRSAV(1,K),K=1,12) 7040 FORMAT(' STND DEV',12F6.1/) C..... C.....READ NUMBER OF SPECIES (NSPEC) C..... READ(5,9000) NSPEC IF(NSPEC.GT.100) CALL ERR3 C..... C.....READ NUMBER OF SPECIES (BMSPEC) AND C.....SPECIES NUMBERS (SPEC) FOR BIOMASS ANALYSIS C..... READ(5,5000) BMSPEC IF(BMSPEC.GT.10) CALL ERR11 READ(5,5001) (SPEC(ISP), ISP=1,BMSPEC) 5000 FORMAT(7X, I3) 5001 FORMAT(10I3) C..... C.....INPUT INDIVIDUAL SPECIES INFORMATION C.....AAA - SPECIES NAME C.....DMAX - MAXIMUM GROWING DEGREE DAYS C.....DMIN - MINIMUM GROWING DEGREE DAYS C.....B3 - INDIVIDUAL SPECIES CONSTANT USED IN GROW C.....B2 - INDIVIDUAL SPECIES CONSTANT USED IN GROW C.....ITOL - LIGHT TOLERANCE CLASS C.....AGEMX - MAXIMUM AGE OF SPECIES C.....G - SCALES THE GROWTH RATE OF EACH SPECIES C.....SPRTND - TENDENCY TO STUMP SPROUT C.....SPRTMN - MINIMUM SIZE TREE THAT WILL SPROUT C.....SPRTMX - MAXIMUM SIZE TREE THAT WILL SPROUT C.....SWITCH - REPRODUCTION SWITCHES USED IN BIRTH C.....MPLANT - MAXIMUM NUMBER OF SEEDLINGS TO PLANT C.....NUM - SPECIES IDENTIFICATION NUMBER C.....D3 - PROPORTION OF GROWING SEASON SPECIES CAN WITHSTAND DROUGHT C.....FROST - MINIMUM JANUARY TEMPERATURE SPECIES CAN WITHSTAND C.....TL - LEAF LITTER TYPE C.....CM1 THROUGH CM5 - PARAMETERS TO CALCULATE NITROGEN GROWTH FACTORS C.....FWT - LEAF WEIGHT/UNIT CROWN AREA C.....SLTA,SLTB - PARAMETERS TO CALCULATE CROWN AREA C.....RTST - ROOT/SHOOT RATIO C.....FRT - FOLIAGE RETENTION TIME IN YEARS C..... WRITE(6,4000) 4000 FORMAT(' AAA',22X,'DMAX',2X,'DMIN',4X,'B3',4X,'B2',2X,'ITOL', & 1X,'AGEMX',3X,'G',2X,'SPRTND',1X,'SPRTMN',1X,'SPRTMX',1X, & 'SWITCH',1X,'MPLANT',1X,'NUM'/' D3',1X,'FROST',1X,'TL',2X, & 'CM1',4X,'CM2',4X,'CM3',4X,'CM4',1X,'CM5',2X,'FWT',2X,'SLTA', & 2X,'SLTB',1X,'RTST',1X,'FRT'/) DO 10 J=1,NSPEC READ(5,9001) (AAA(J,I),I=1,6),DMAX(J),DMIN(J),B3(J),B2(J), > ITOL(J),AGEMX(J),G(J),SPRTND(J),SPRTMN(J),SPRTMX(J), > (SWITCH(J,I),I=1,5),MPLANT(J),NUM,D3(J),FROST(J), > TL(J),CM1(J),CM2(J),CM3(J),CM4(J),CM5(J),FWT(J),SLTA(J), > SLTB(J),RTST(J),FRT(J) C..... C.....RESCALE NITROGEN GROWTH MULTIPLIER PARAMETERS SO THAT C.....GROWTH IS NOT N LIMITED WHEN N AVAILABILITIES ARE GREATER THAN C.....SATURATED LEVELS FOR EACH TYPE OF RESPONSE C..... IF(CM4(J).EQ.-5.0) CM4(J)=CM4(J)/3.0 IF(CM5(J).EQ.2.9) CM5(J)=CM5(J)/3.0 IF(CM4(J).EQ.-1.2) CM4(J)=CM4(J)/2.4 IF(CM5(J).EQ.1.3) CM5(J)=CM5(J)/2.4 IF(CM4(J).EQ.-0.6) CM4(J)=CM4(J)/1.75 IF(CM5(J).EQ.1.0) CM5(J)=CM5(J)/1.75 C..... C.....WRITE SPECIES PARAMETERS TO PRINTER C..... WRITE(6,9002) (AAA(J,I),I=1,6),DMAX(J),DMIN(J),B3(J),B2(J), > ITOL(J),AGEMX(J),G(J),SPRTND(J),SPRTMN(J),SPRTMX(J), > (SWITCH(J,I),I=1,5),MPLANT(J),NUM,D3(J),FROST(J), > TL(J),CM1(J),CM2(J),CM3(J),CM4(J),CM5(J),FWT(J),SLTA(J), > SLTB(J),RTST(J),FRT(J) 10 CONTINUE C..... C.....INPUT FOREST FLOOR DECOMPOSITION PARAMETERS C.....NLVAR - # OF VARIABLES USED TO CALCULATE DECOMPOSITION C.....NLT - # OF LITTER TYPES C..... READ(5,9003) NLVAR,NLT 9003 FORMAT(8X,I2,7X,I2) DO 15 I=1,NLT READ(5,9004) (FDAT(I,J),J=1,NLVAR) 9004 FORMAT(F3.0,3F6.4,F4.0,F3.0,F5.3,F7.4,F5.3,F4.2) 15 CONTINUE C..... C.....NCOHRT - NUMBER OF LITTER COHORTS INITIALLY PRESENT C.....IF ONLY HUMUS IS PRESENT, NCOHRT = 1 C..... READ(5,9005) NCOHRT 9005 FORMAT(9X,I1) IF(NCOHRT.GT.100) CALL ERR4 C..... C.....BASESC - STARTING HUMUS WEIGHT (T/HA) C.....BASESN - STARTING HUMUS N CONTENT (T/HA) C..... READ(5,9006) BASESC,BASESN 9006 FORMAT(F3.0,F6.4) C..... C.....WRITE INITIAL SOIL CONDITIONS TO PRINTER C..... WRITE(6,9007) FC,DRY 9007 FORMAT(' '/' FIELD CAPACITY (CM)=',F5.1,' WILTING POINT (CM)=', &F5.1) WRITE(6,9008) BASESC,BASESN 9008 FORMAT(' INITIAL SOIL O.M.=',F4.1,' INITIAL SOIL N=',F5.2) RETURN 9000 FORMAT(6X,I3) 9001 FORMAT(6A4,F6.0,F5.0,F4.4,F5.2,I1,F6.2,F5.1,F2.0,F4.1,F4.0, 1 5L1,I4,I5/F5.3,F5.0,F3.0,F5.2,F7.2,F7.5,F5.1,F4.1,F5.0, 2 2F5.3,F4.1,F3.0) 9002 FORMAT(' ',6A4,F6.0,F6.0,F7.4,F6.2,I5,F6.0,F7.2,F5.0,F6.1,F7.0,3X, 1 5L1,I4,I6/' ',F5.3,F5.0,F4.0,F5.2,F7.2,F8.5,F5.1,F4.1,F5.0, 2 2F6.3,F4.1,F5.0/) END C C SUBROUTINE KILL C C C SUBROUTINE KILL KILLS TREES BY AGE DEPENDENT MORTALITY (ONLY 1% C REACH MAXIMUM AGE) AND AGE INDEPENDENT MORTALITY (PROBABILITY OF C SURVIVING 10 CONSECUTIVE YEARS OF SLOW GROWTH (SEE GROW) = 1%). C DECISIONS ON WHETHER OR NOT TO KILL A TREE ARE PARTLY BASED ON C RANDOM NUMBERS SUPPLIED BY URAND. C KILL ALSO CALCULATES LITTER AMOUNTS, WHICH ARE DECAYED IN C SUBROUTINE DECOMP. C SUBROUTINE KILL(KYR) COMMON/FOREST/NTREES(100),DBH(1500),IAGE(1500),KSPRT(100), > NEWTR(100),SUMLA(700),NEW(100),SWTCH(5) COMMON/PARAM/AAA(100,6),DMAX(100),DMIN(100),B3(100),B2(100), > ITOL(100),AGEMX(100),G(100),SPRTND(100),SPRTMN(100), 1SPRTMX(100),SWITCH(100,5),MPLANT(100),D3(100),FROST(100), 2TL(100),CM1(100),CM2(100),CM3(100),CM4(100),CM5(100), 3FWT(100),SLTA(100),SLTB(100),RTST(100),FRT(100) COMMON/CONST/NSPEC,DEGD COMMON/DEAD/NOGRO(1500),NTEMP(1500) COMMON/PROD/AWP(1500) COMMON/COUNT/NTOT,NYEAR,KPRNT,NMAX,KLAST,NWRITE,KWRITE COMMON/TEMP/DTEMP(1500),ITEMP(1500) COMMON/SEED/USEED(15) COMMON/DCMP/AVAILN,TYL(20),C(100,15),FDAT(20,10),FF(20,3), 1NCOHRT,HCN,BASESC,BASESN,SCO2,TYLN INTEGER USEED LOGICAL SWITCH,SWTCH KNT = 0 C..... C......INITIALIZE LITTER C..... DO 5 L=1,20 TYL(L)=0.0 5 CONTINUE C..... C.....INITIALIZE PLOT BASAL AREA C..... BA=0.0 C..... C.....BEGIN MAIN KILLING LOOP C.... DO 30 I=1,NSPEC IF (NTREES(I).EQ.0) GO TO 30 NL = KNT+1 NU = NTREES(I)+KNT DO 20 K=NL,NU C..... C.....CALCULATE LEAF PRODUCTION (T/HA) C.....(ABER ET AL. 1982. FOREST SCI. 28:31-45) C..... FOLW=((SLTA(I)+SLTB(I)*DBH(K))/2.)**2*3.14*FWT(I)*.000012 C..... C.....CALCULATE BASAL AREA C..... BA=BA+.0314*(DBH(K)*0.5)**2 C..... C.....KILL TREES BASED ON PROBABILITY THAT ONLY 1% REACH MAXIMUM AGE C..... YFL=URAND(USEED(1)) IF (YFL.LE.(4.605/AGEMX(I))) GO TO 10 C..... C.....CHECK TO SEE IF THERE WAS ANY GROWTH FOR TREE C.....IF NOT, INCREASE PROBABILITY OF DIEING SO THAT THE TREE HAS C.....A 1% CHANCE OF SURVIVING 10 CONSECUTIVE YEARS OF SLOW GROWTH C..... IF (NOGRO(K).GT. -2) GO TO 19 YFL=URAND(USEED(2)) IF (YFL.GT.0.368) GO TO 19 10 CONTINUE NTREES(I) = NTREES(I)-1 C..... C.....CHECK TO SEE IF DEAD TREE CAN STUMP SPROUT. INCREMENT KSPRT C.....IF TREE CAN SPROUT C..... IF (DBH(K).GT.SPRTMN(I).AND.DBH(K).LT.SPRTMX(I)) KSPRT(I) = > KSPRT(I)+1 C..... C.....CALCULATE WOODY LITTER IN T/HA (SOLLINS ET AL. 1973. C.....EDFB-IBP-73-2) C..... BD=.60 IF(DBH(K).LE.10.) TYL(14)=TYL(14)+BD*(.00143*DBH(K)**2.393) IF(DBH(K).GT.10.) TYL(15)=TYL(15)+BD*(.00143*DBH(K)**2.393) C..... C.....FLAG DEAD TREES C..... DBH(K) = -1.0 C..... C.....CALCULATE LEAF LITTER BY QUALITY CLASS IN T/HA C.....IF THE TREE IS SLOW GROWING BUT DIDN'T DIE, LEAF LITTER C.....IS HALVED (BRAY AND GORHAM.1964.ADVANCES IN ECOLOGICAL C.....RESEARCH 2:101-157). C.....IF THE TREE DIED, TOTAL LEAF BIOMASS IS RETURNED TO THE SOIL. C..... 19 L=TL(I) IF(NOGRO(K).EQ.-2.AND.DBH(K).GT.-1.) FOLW=FOLW*0.5 IF(DBH(K).LT.0.) FOLW=FOLW*FRT(I) TYL(L)=TYL(L)+FOLW C..... C.....CALCULATE ROOT LITTER (T/HA) C..... TYL(13)=TYL(13)+1.3*FOLW*RTST(I) 20 CONTINUE KNT = NU 30 CONTINUE C..... C.....CALCULATE TOTAL LEAF LITTER (T/HA) C..... DO 35 L=1,12 TYL(17)=TYL(17) + TYL(L) 35 CONTINUE C..... C.....CALCULATE TWIG LITTER IN T/HA (CHRISTENSEN.1977.OIKOS 28:177-186) C..... TYL(16)=BA/333. C..... C.....CALCULATE TOTAL LITTER (T/HA) C..... TYL(18)=TYL(13)+TYL(14)+TYL(15)+TYL(16)+TYL(17) C..... C.....REWRITE DIAMETERS AND AGES TO ELIMINATE DEAD TREES C..... K = 0 DO 40 I=1,1500 IF (DBH(I).EQ.0.) GO TO 50 IF (DBH(I).LT.0.) GO TO 40 K = K+1 DBH(K) = DBH(I) IAGE(K) = IAGE(I) NOGRO(K) = NOGRO(I) 40 CONTINUE 50 NTOT = K IF (NTOT.EQ.0) RETURN NTOT1 = K+1 IF(NTOT1.GT.1500) CALL ERR5 C..... C.....ELIMINATE DEAD TREES C..... DO 60 I=NTOT1,NU DBH(I) = 0. IAGE(I) = 0 NOGRO(I) = 0 60 CONTINUE RETURN END C C SUBROUTINE LININT C C C SUBROUTINE LININT INTERPOLATES MONTHLY TEMPERATURES, PRECIPITATION C AND THEIR STND DEV FOR ALL YEARS BRACKETED BY TWO YEARS OF C DIFFERENT CLIMATES. THESE YEARS ARE SUPPLIED IN ARRAY X. C SUBROUTINE LININT(P1,P2,XX,NTYPE) DIMENSION P1(12),P2(12) COMMON/INTERP/IPOLAT,X(10) COMMON /LINEAR/TSAV(45,12),VTSAV(45,12),RSAV(45,12),VRSAV(45,12) NPTS = IPOLAT NPT1 = NPTS-1 C..... C.....FIND YEARS BETWEEN WHICH LINEAR INTERPOLATIONS SHOULD BE C.....MADE. XX - CURRENT YEAR. X(I) AND X(I+1) - BRACKETING YEARS C.....AS SPECIFIED IN ARRAY X C..... DO 250 I=1,NPT1 IF(XX .GT. X(I) .AND. XX .LE. X(I+1)) GO TO 300 250 CONTINUE 300 CONTINUE DO 500 K=1,12 C..... C.....IF TEMPE CALLS LININT, NTYPE = 1 C.....IF MOIST CALLS LININT, NTYPE = 2 C..... IF(NTYPE .EQ. 2) GO TO 400 C..... C.....INTERPOLATE MEAN MONTHLY TEMPERATURES (C) BETWEEN YEARS C.....OF DIFFERENT CLIMATES C..... P1(K) = TSAV(I,K)+((TSAV(I+1,K)-TSAV(I,K))/ 1(X(I+1)-X(I)))*(XX-X(I)) C..... C.....INTERPOLATE STND DEVS OF MONTHLY TEMPERATURES C..... P2(K) = VTSAV(I,K)+((VTSAV(I+1,K)-VTSAV(I,K))/ 1(X(I+1)-X(I)))*(XX-X(I)) GO TO 450 C..... C.....INTERPOLATE MEAN MONTHLY RAINFALL (CM) BETWEEN YEARS C.....OF DIFFERENT CLIMATES C..... 400 P1(K) = RSAV(I,K)+((RSAV(I+1,K)-RSAV(I,K))/ 1(X(I+1)-X(I)))*(XX-X(I)) C..... C.....INTERPOLATE STND DEVS OF MONTHLY RAINFALL C..... P2(K) = VRSAV(I,K)+((VRSAV(I+1,K)-VRSAV(I,K))/ 1(X(I+1)-X(I)))*(XX-X(I)) 450 CONTINUE 500 CONTINUE RETURN END C C SUBROUTINE MOIST C C C SUBROUTINE MOIST CALCULATES THE FRACTION OF THE C GROWING SEASON WITH UNFAVORABLE SOIL MOISTURE FOR GROWTH C (FJ) USED IN SUBROUTINE GMULT TO DETERMINE SOIL MOISTURE C GROWTH MULTIPLIERS, AND ACTUAL EVAPOTRANSPIRATION (AET) C USED IN SUBROUTINE DECOMP TO DETERMINE DECAY RATES. C THE SUBROUTINE SIMULATES THE METHOD OF THORNTHWAITE AND MATHER C (1957) AS MODIFIED BY PASTOR AND POST (1984). C TEMPERATURES ARE PROVIDED BY SUBROUTINE TEMPE. MONTHLY C PRECIPITATION IS CALCULATED THE SAME WAY AS TEMPERATURES WERE IN C TEMPE. C C DATA THAT IS REQUIRED AND READ BY SUBROUTINE INPUT: C T = AVERAGE MONTHLY TEMPERATURES (JAN-DEC) CENTIGRADE C R = AVERAGE MONTHLY RAINFALL (JAN-DEC) CENTIMETERS C VR = STANDARD DEVIATION OF AVERAGE MONTHLY RAINFALL C FC = CENTIMETERS OF WATER THE SOIL CAN HOLD AT FIELD CAPACITY C DRY = CENTIMETERS OF WATER BELOW WHICH TREE GROWTH STOPS C (-15 BARS) C BGS = YEAR DAY ON WHICH THE GROWING SEASON BEGINS C EGS = YEAR DAY ON WHICH THE GROWING SEASON ENDS C PLAT = LATITUDE OF PLOT (DEGREES NORTH) C SUBROUTINE MOIST(KYR) COMMON/WATER/T(12),VT(12),RT(12),R(12),VR(12),FC,DRY,BGS,EGS,PLAT, 1FJ,AET COMMON/COUNT/NTOT,NYEAR,KPRNT,NMAX,KLAST,NWRITE,KWRITE COMMON/CONST/NSPEC,DEGD COMMON/SEED/USEED(15) DIMENSION Z(2),CLAT(12,26) C DIMENSION Z(2),CLAT(12,26),DLAT(12,6) DIMENSION DAYS(12) C EQUIVALENCE (CLAT(1,21),DLAT(1,1)) INTEGER USEED C..... C.....MONTHLY CORRECTION FACTORS FOR 25-50 DEGREES LATITUDE NORTH C.....DAYS(K) = NUMBER OF DAYS BETWEEN MID-MONTH K-1 AND K C..... C DATA DLAT/.80,.81,1.02,1.13,1.28,1.29,1.31,1.21,1.04,.94,.79,.75, C 6 .79,.81,1.02,1.13,1.29,1.31,1.32,1.22,1.04,.94,.79,.74, C 7 .77,.80,1.02,1.14,1.30,1.32,1.32,1.22,1.04,.93,.78,.73, C 8 .76,.80,1.02,1.14,1.31,1.33,1.34,1.23,1.05,.93,.77,.72, C 9 .75,.79,1.02,1.14,1.32,1.34,1.35,1.24,1.05,.93,.76,.71, C * .74,.78,1.02,1.15,1.33,1.36,1.37,1.25,1.06,.92,.76,.70/ DATA CLAT/.93,.89,1.03,1.06,1.15,1.14,1.17,1.12,1.02,.99,.91,.91, 6 .92,.88,1.03,1.06,1.15,1.15,1.17,1.12,1.02,.99,.91,.91, 7 .92,.88,1.03,1.07,1.16,1.15,1.18,1.13,1.02,.99,.90,.90, 8 .91,.88,1.03,1.07,1.16,1.16,1.18,1.13,1.02,. 8,.90,. 0, 9 .91,.87,1.03,1.07,1.17,1.16,1.19,1.13,1.03,.98,.90,.89, * .90,.87,1.03,1.08,1.18,1.17,1.20,1.14,1.03,.98,.89,.88, 1 .90,.87,1.03,1.08,1.18,1.18,1.20,1.14,1.03,.98,.89,.88, 2 .89,.86,1.03,1.08,1.19,1.19,1.21,1.15,1.03,.98,.88,.87, 3 .88,.86,1.03,1.09,1.19,1.20,1.22,1.15,1.03,.97,.88,.86, 4 .88,.85,1.03,1.09,1.20,1.20,1.22,1.16,1.03,.97,.87,.86, 5 .87,.85,1.03,1.09,1.21,1.21,1.23,1.16,1.03,.97,.86,.85, 6 .87,.85,1.03,1.10,1.21,1.22,1.24,1.16,1.03,.97,.86,.84, 7 .86,.84,1.03,1.10,1.22,1.23,1.25,1.17,1.03,.97,.85,.83, 8 .85,.84,1.03,1.10,1.23,1.24,1.25,1.17,1.04,.96,.84,.83, 9 .85,.84,1.03,1.11,1.23,1.24,1.26,1.18,1.04,.96,.84,.82, * .84,.83,1.03,1.11,1.24,1.25,1.27,1.18,1.04,.96,.83,.81, 1 .83,.83,1.03,1.11,1.25,1.26,1.27,1.19,1.04,.96,.82,.80, 2 .82,.83,1.03,1.12,1.26,1.27,1.28,1.19,1.04,.95,.82,.79, 3 .81,.82,1.02,1.12,1.26,1.28,1.29,1.20,1.04,.95,.81,.77, 4 .81,.82,1.02,1.13,1.27,1.29,1.30,1.20,1.04,.95,.80,.76, 5 .80,.81,1.02,1.13,1.28,1.29,1.31,1.21,1.04,.94,.79,.75, 6 .79,.81,1.02,1.13,1.29,1.31,1.32,1.22,1.04,.94,.79,.74, 7 .77,.80,1.02,1.14,1.30,1.32,1.32,1.22,1.04,.93,.78,.73, 8 .76,.80,1.02,1.14,1.31,1.33,1.34,1.23,1.05,.93,.77,.72, 9 .75,.79,1.02,1.14,1.32,1.34,1.35,1.24,1.05,.93,.76,.71, * .74,.78,1.02,1.15,1.33,1.36,1.37,1.25,1.06,.92,.76,.70/ DATA DAYS/31.,28.,31.,30.,31.,30.,31.,31.,30.,31.,30.,30./ DATA NCT/0/ C..... C.....ADJUST LATITUDE POINTER C..... LAT=(PLAT+.5)-24 C..... C.....INITIALIZE WATER CONTENT OF SOIL IN JANUARY TO FC C..... XFC=10.*FC WATER=FC C..... C.....INITIALIZE THORNTHWAITE PARAMETERS C.....TE = TEMPERATURE EFFICIENCY C.....A = EXPONENT OF EVAPOTRANSPIRATION FUNCTION C.....U = POTENTIAL EVAPOTRANSPIRATION C.....AET = ACTUAL EVAPOTRANSPIRATION C.....ACCPWL = ACCUMULATED POTENTIAL WATER LOSS C..... AET=0.0 ACCPWL=0.0 TE=0. RSAVE = RT(1) DO 10 K=1,12 IF(RT(K) .LT. 0.0) RT(K) = 0.0 TE=TE+(.2*RT(K))**1.514 10 CONTINUE A=.675*TE**3-77.1*TE**2+17920.0*TE+492390.0 A=.000001*A C..... C.....INITIALIZE THE NUMBER OF DRY DAYS (DD), C.....AND CURRENT DAY OF YEAR (CDAY) C..... DD=0.0 CDAY=15.0 C..... C.....CALL LININT TO INTERPOLATE MONTHLY RAINFALL AND STND DEV C.....BETWEEN YEARS OF DIFFERENT CLIMATE C..... YR = FLOAT(KYR) CALL LININT(R,VR,YR,2) C..... C.....MAIN LOOP FOR YEARLY WATER BALANCE CALCULATION BY MONTH C..... DO 50 K=1,12 OWATER=WATER NCT=NCT+1 IF(NCT.EQ.2) GO TO 36 C..... C.....CALL GGNORD TO PROVIDE NORMALLY DISTRIBUTED RANDOM NUMBER C..... CALL GGNORD(USEED(14),USEED(15),Z) GO TO 38 36 Z(1)=Z(2) NCT=0 C..... C.....CALCULATE THIS MONTH'S RAINFALL (INTERPOLATED MEAN +/- NORMALLY C.....DISTRIBUTED RANDOM NUMBER TIMES THE INTERPOLATED STND DEV) C..... 38 RAIN=R(K)+VR(K)*Z(1) IF(RAIN.LT.0.0) RAIN=0.0 RAIN=AMAX1(0.0,RAIN) TTMP=RT(K) IF(TTMP.LT.0.0) TTMP=0.0 C..... C.....CALCULATE POTENTIAL EVAPOTRANSPIRATION (U) C.....(THORNTHWAITE AND MATHER.1957.PUBL. IN CLIMATOLOGY 10:83-311). C..... U=1.6*((10.0*TTMP/TE)**A)*CLAT(K,LAT) C..... C.....CALCULATE POTENTIAL WATER LOSS THIS MONTH C..... PWL=RAIN-U C..... C.....IF RAIN SATISFIES U THIS MONTH, DON'T DRAW ON SOIL WATER C..... IF(PWL.GE.0.0) GO TO 55 C..... C.....IF RAIN DOESN'T SATISFY U, ADD THIS MONTH'S POTENTIAL C.....WATER LOSS TO ACCUMULATED POTENTIAL WATER LOSS FROM SOIL C..... ACCPWL=ACCPWL + PWL XACPWL=ACCPWL*10. C..... C.....CALCULATE WATER RETAINED IN SOIL GIVEN SO MUCH ACCUMULATED C.....POTENTIAL WATER LOSS C.....(PASTOR AND POST.1984. CAN.J.FOR.RES. 14:466-467) C..... WATER=FC*(EXP((.000461-1.10559/XFC)*(-1.*XACPWL))) IF(WATER.LT.0.0) WATER=0.0 C..... C.....CSM - CHANGE IN SOIL MOISTURE DURING THIS MONTH C..... CSM = WATER-OWATER C..... C.....CALCULATE ACTUAL EVAPOTRANSPIRATION (AET) IF SOIL WATER C.....IS DRAWN DOWN THIS MONTH C..... AET=AET+(RAIN-CSM) GO TO 56 55 WATER = OWATER+PWL IF(WATER.GE.FC)WATER = FC CSM=WATER-OWATER C..... C.....IF SOIL IS PARTIALLY RECHARGED, REDUCE ACCUMULATED C.....POTENTIAL WATER LOSS ACCORDINGLY C..... ACCPWL=ACCPWL+CSM C..... C.....IF SOIL IS COMPLETELY RECHARGED, RESET ACCUMULATED C.....POTENTIAL WATER LOSS TO ZERO C..... IF(WATER.GE.FC)ACCPWL=0. C..... C.....IF SOIL WATER IS NOT DRAWN UPON, ADD U TO AET C..... AET=AET+U 56 CONTINUE OCDAY=CDAY CDAY=CDAY+DAYS(K) DDI=0.0 C..... C.....INCREMENT THE NUMBER OF DRY DAYS, INTERPOLATING C.....IF NECESSARY, TRUNCATE AT ENDS OF GROWING SEASON C..... IF(CDAY .LE. BGS) GO TO 50 IF(OCDAY .GE. EGS) GO TO 50 IF(OWATER .GE. DRY .AND. WATER .GE. DRY) GO TO 50 IF(OWATER .GT. DRY .AND. WATER .LT. DRY) GO TO 20 IF(OWATER .LT. DRY .AND. WATER .GT. DRY) GO TO 30 DDI=DAYS(K) IF(OCDAY .LT. BGS .AND. CDAY .GT. BGS) DDI=CDAY-BGS IF(OCDAY .LT. EGS .AND. CDAY .GT. EGS) DDI=EGS-OCDAY GO TO 40 20 DDI=DAYS(K)*(DRY-WATER)/(OWATER-WATER) IF(OCDAY .LT. BGS .AND. CDAY .GT. BGS) DDI=AMIN1(DDI,CDAY-BGS) IF(OCDAY .LT. EGS .AND. CDAY .GT. EGS) DDI=EGS-CDAY+DDI IF(DDI .LT. 0.0) DDI=0.0 GO TO 40 30 DDI=DAYS(K)*(DRY-OWATER)/(WATER-OWATER) IF(OCDAY .LT. BGS .AND. CDAY .GT. BGS) DDI=OCDAY+DDI-BGS IF(DDI .LT. 0.0) DDI=0.0 IF(OCDAY .LT. EGS .AND. CDAY .GT. EGS) DDI=AMIN1(DDI,EGS-OCDAY) 40 DD=DD+DDI 50 CONTINUE C..... C.....SAVE TOTAL NUMBER OF DRY DAYS FOR YEAR C..... FJ=DD RT(1) = RSAVE C..... C.....CONVERT AET FROM CM TO MM C..... AET=AET*10. RETURN END C C SUBROUTINE OUTPUT C C C SUBROUTINE OUTPUT CALCULATES AND STORES SPECIES BIOMASS, TOTAL C ABOVEGROUND BIOMASS AND NPP, TOTAL NUMBER OF STEMS, AND LEAF AREA. C IT ALSO STORES AVAILABLE N, SOIL C:N, SOIL O.M., LEAF LITTER, C TOTAL LITTER, AND AET. C MEANS AND 95% CONFIDENCE INTERVALS ARE CALCULATED ON ALL STAND C LEVEL VARIABLES AT THE END OF THE RUN. C ARRAY ST CONTAIN'S STUDENT'S T FOR N=1 TO GREATER THAN OR C EQUAL TO 30. C SUBROUTINE OUTPUT(KYR,IPLOT) COMMON/WATER/T(12),VT(12),RT(12),R(12),VR(12),FC,DRY,BGS,EGS,PLAT, 1FJ,AET COMMON/FOREST/NTREES(100),DBH(1500),IAGE(1500),KSPRT(100), > NEWTR(100),SUMLA(700),NEW(100),SWTCH(5) COMMON/PARAM/AAA(100,6),DMAX(100),DMIN(100),B3(100),B2(100), 1 ITOL(100),AGEMX(100),G(100),SPRTND(100),SPRTMN(100), 2 SPRTMX(100),SWITCH(100,5),MPLANT(100),D3(100),FROST(100), 3 TL(100),CM1(100),CM2(100),CM3(100),CM4(100),CM5(100), 4 FWT(100),SLTA(100),SLTB(100),RTST(100),FRT(100) COMMON/COUNT/NTOT,NYEAR,KPRNT,NMAX,KLAST,NWRITE,KWRITE COMMON/CONST/NSPEC,DEGD COMMON/PROD/AWP(1500) COMMON/DCMP/AVAILN,TYL(20),C(100,15),FDAT(20,10),FF(20,3), 1NCOHRT,HCN,BASESC,BASESN,SCO2,TYLN COMMON/SPECIE/SPEC(10),BMSPEC DIMENSION BAR(100) DIMENSION TSTEM(70,100),TAB(70,100),FL(70,100),TOTL(70,100) DIMENSION TNAP(70,100),AVLN(70,100),CN(70,100),SCO2C(70,100) DIMENSION SOM(70,100),ET(70,100),XSTEM(70),XTAB(70),XFL(70) DIMENSION XTL(70),XNAP(70),XAVLN(70),XCN(70),XSCO2C(70) DIMENSION XSOM(70),XET(70),NYR(70),AVG(70,11) DIMENSION SSTEM(70),STAB(70),SFL(70),STL(70),SNAP(70) DIMENSION SAVLN(70),SCN(70),SSCO2C(70),SSOM(70),SET(70) DIMENSION CONF(70,11),ST(30) DIMENSION XBMSPE(10,70),SSBMSP(10,70),AVGBM(70,10) DIMENSION CONFBM(70,10),XBAR(10,70,100) DIMENSION NFREQ(100), TSBAR(100) INTEGER BMSPEC, SPEC LOGICAL SWITCH, SWTCH DATA ST/12.706,4.303,3.182,2.776,2.571,2.447,2.365,2.306,2.262, &2.228,2.201,2.179,2.160,2.145,2.131,2.120,2.110,2.101,2.093, &2.086,2.080,2.074,2.069,2.064,2.060,2.056,2.052,2.048,2.045, &2.042/ C..... C.....INITIALIZATION C.....AREA - LEAF AREA C.....FOLW - LEAF BIOMASS C.....AVAILN - AVAILABLE NITROGEN (FROM GMULT) C.....TBAR - TOTAL ABOVEGROUND BIOMASS C.....TAWP - TOTAL ABOVEGROUND WOODY PRODUCTION C.....NTOT - NUMBER OF TREES C..... KWRITE = KWRITE + 1 KYR1 = KYR/KPRNT+1 AREA=0.0 FOLW=0.0 AVAILN=AVAILN*1000. TBAR=0.0 TAWP = 0.0 NYR(KYR1) = KYR YR=FLOAT(KYR) NTOT = 0 TYLN=TYLN*1000. C..... C.....CALCULATE SPECIES BIOMASS, TOTAL BIOMASS, TOTAL NUMBER OF C.....STEMS, LEAF AREA, AND TOTAL WOODY PRODUCTION C..... NL = 1 DO 20 I =1,NSPEC BAR(I)=0.0 IF (NTREES(I).EQ.0) GO TO 20 NU = NL+NTREES(I)-1 RET=FRT(I) DO 10 J=NL,NU AGE=IAGE(J) IF(AGE.LT.RET) RET=AGE C..... C.....CALCULATE LEAF BIOMASS (KG/TREE) C.....(ABER ET AL. 1982. FOREST SCI. 28:31-45) C..... FOLW=((SLTA(I)+SLTB(I)*DBH(J))/2.)**2*3.14*FWT(I)* & RET*.001 C..... C.....CALCULATE SPECIES BIOMASS (KG/PLOT) C.....(SOLLINS ET AL. 1973. EDFB-IBP-73-2) C..... BAR(I) = BAR(I)+.1193*DBH(J)**2.393 + FOLW C..... C.....CALCULATE LEAF AREA INDEX C.....(SOLLINS ET AL. 1973. EDFB-IBP-73-2) C..... AREA=AREA+1.9283295E-4*DBH(J)**2.129 C..... C.....CALCULATE WOODY PRODUCTION (KG/PLOT) C..... TAWP = TAWP + AWP(J) 10 CONTINUE C..... C.....CALCULATE TOTAL ABOVEGROUND BIOMASS (KG/PLOT) C..... TBAR=TBAR+BAR(I) NL = NU+1 C..... C.....CALCULATE TOTAL NUMBER OF TREES PER PLOT C..... NTOT = NTOT +NTREES(I) IF(NTOT.GT.1500) CALL ERR7 20 CONTINUE C..... C.....CONVERT NUMBER OF TREES PER PLOT TO NUMBER PER HA C..... ATOT=NTOT ATOT=ATOT*12. C..... C.....CONVERT TOTAL ABOVEGROUND BIOMASS AND WOODY PRODUCTION C.....TO T/HA. C..... TBAR=TBAR*.012 TAWP = TAWP*.012 C..... C.....CALCULATE TOTAL ABOVEGROUND PRODUCTION C..... TYNAP = TAWP+TYL(17) C..... C.....CONVERT SPECIES BIOMASS TO T/HA C..... DO 30 IV1 = 1,NSPEC BAR(IV1) = BAR(IV1)*.012 30 CONTINUE C..... C.....PLACE THIS YEAR'S RESULTS IN STORAGE ARRAYS FOR STATISTICAL C.....CALCULATIONS C.....TSTEM - NUMBER OF STEMS C.....TAB - TOTAL ABOVEGROUND BIOMASS C.....FL - LEAF LITTER C.....TOTL - LEAF LITTER N C.....TNAP - TOTAL NET ABOVEGROUND PRODUCTION C.....AVLN - AVAILABLE NITROGEN C.....CN - HUMUS C:N RATIO C.....SCO2C - SOIL CO2 EVOLUTION C.....SOM - SOIL ORGANIC MATTER C.....ET - ACTUAL EVAPOTRANSPIRATION C.....XBAR - SPECIES(I) BIOMASS C..... TSTEM(KYR1,IPLOT) = ATOT TAB(KYR1,IPLOT) = TBAR FL(KYR1,IPLOT) =TYL(17) TOTL(KYR1,IPLOT) = TYLN TNAP(KYR1,IPLOT) = TYNAP AVLN(KYR1,IPLOT) = AVAILN CN(KYR1,IPLOT) = HCN SCO2C(KYR1,IPLOT) = SCO2 SOM(KYR1,IPLOT) = FF(19,2) ET(KYR1,IPLOT) = AET DO 90 III = 1,BMSPEC XBAR(III,KYR1,IPLOT) = BAR(SPEC(III)) 90 CONTINUE C..... C.....WRITE STATEMENT TO STORE DATA ON TAPE C..... C WRITE(9) (BAR(IV1),IV1=1,NSPEC),(NTREES(IV1),IV1=1,NSPEC), C 1TBAR,ATOT,TYNAP,DEGD,FJ,AVAILN,((FF(K,L),K=1,18),L=1,3), C 2(TYL(J),J=1,20),HCN,SCO2,AET C..... IF(KYR .GT. 1) GO TO 9990 DO 9989 IV1=1,NSPEC NFREQ(IV1)=0 TSBAR(IV1)=0.0 9989 CONTINUE 9990 CONTINUE DO 9991 IV1 = 1,NSPEC IF(BAR(IV1) .LE. 0.0) GO TO 9991 NFREQ(IV1)=NFREQ(IV1)+1 TSBAR(IV1)=TSBAR(IV1)+BAR(IV1) 9991 CONTINUE IF(KYR .LT. NYEAR) GO TO 9992 DO 9993 IV1=1,NSPEC IF(NFREQ(IV1) .LE. 0) GO TO 9993 C WRITE(9,9994) IV1,(AAA(IV1,IV2),IV2=1,6),NFREQ(IV1),TSBAR(IV1) 9994 FORMAT(5X,I5,3X,6A4,I10,F15.3) 9993 CONTINUE 9992 CONTINUE C WRITE(9,999) KYR, C 1BAR(1),BAR(3),BAR(5),BAR(19),BAR(20),BAR(27),BAR(28),BAR(32), C 2BAR(33),BAR(44),BAR(55),BAR(58),BAR(60),BAR(65),BAR(67), C 3BAR(69),BAR(72) C 1(BAR(IV1),IV1=1,NSPEC), C 2TBAR,DEGD,FJ,AVAILN,AET C 999 FORMAT(I10,8F10.3) C 999 FORMAT(I7,/,10(8F10.3,/)) C..... C.....BYPASS STATISTICAL CALCULATIONS IF IT IS NOT THE LAST YEAR C.....OF THE LAST PLOT C..... IF(KWRITE.LT.NWRITE) GO TO 60 C..... C.....INITIALIZE SUMMATION ARRAYS C..... DO 35 I=1,70 XSTEM(I) = 0.0 XTAB(I) = 0.0 XFL(I) = 0.0 XTL(I) = 0.0 XNAP(I) = 0.0 XAVLN(I) = 0.0 XCN(I) = 0.0 XSCO2C(I) = 0.0 XSOM(I) = 0.0 XET(I) = 0.0 SSTEM(I) = 0.0 STAB(I) = 0.0 SFL(I) = 0.0 STL(I) = 0.0 SNAP(I) = 0.0 SAVLN(I) = 0.0 SCN(I) = 0.0 SSCO2C(I) = 0.0 SSOM(I) = 0.0 SET(I) = 0.0 DO 91 III = 1,BMSPEC XBMSPE(III,I) = 0.0 SSBMSP(III,I)= 0.0 91 CONTINUE 35 CONTINUE C..... C.....CHOOSE APPROPRIATE STUDENT'S T FOR CONFIDENCE INTERVALS C..... L = KLAST IF(L.GT.30) L=30 TS=ST(L) C..... C.....BEGIN MAIN STATISTICAL LOOP C..... DO 40 I = 1,NMAX C..... C.....ACCUMULATE SUMS C..... DO 50 J = 1,KLAST XSTEM(I) = XSTEM(I) + TSTEM(I,J) XTAB(I) = XTAB(I) + TAB(I,J) XFL(I) = XFL(I) + FL(I,J) XTL(I) = XTL(I) + TOTL(I,J) XNAP(I) = XNAP(I) + TNAP(I,J) XAVLN(I) = XAVLN(I) + AVLN(I,J) XCN(I) = XCN(I) + CN(I,J) XSCO2C(I) = XSCO2C(I) + SCO2C(I,J) XSOM(I) = XSOM(I) + SOM(I,J) XET(I) = XET(I) + ET(I,J) DO 93 III = 1,BMSPEC XBMSPE(III,I) = XBMSPE(III,I) + XBAR(III,I,J) 93 CONTINUE C..... C.....ACCUMULATE SUMS OF SQUARES C..... SSTEM(I) = SSTEM(I) + TSTEM(I,J)**2 STAB(I) = STAB(I) + TAB(I,J)**2 SFL(I) = SFL(I) + FL(I,J)**2 STL(I) = STL(I) + TOTL(I,J)**2 SNAP(I) = SNAP(I) + TNAP(I,J)**2 SAVLN(I) = SAVLN(I) + AVLN(I,J)**2 SCN(I) = SCN(I) + CN(I,J)**2 SSCO2C(I) = SSCO2C(I) + SCO2C(I,J)**2 SSOM(I) = SSOM(I) + SOM(I,J)**2 SET(I) = SET(I) + ET(I,J)**2 DO 94 III = 1,BMSPEC SSBMSP(III,I) = SSBMSP(III,I) + XBAR(III,I,J)**2 94 CONTINUE 50 CONTINUE PLOTS = KLAST C..... C.....FIND MEANS C..... AVG(I,1) = NYR(I) AVG(I,2) = XSTEM(I)/PLOTS AVG(I,3) = XTAB(I)/PLOTS AVG(I,4) = XFL(I)/PLOTS AVG(I,5) = XTL(I)/PLOTS AVG(I,6) = XNAP(I)/PLOTS AVG(I,7) = XAVLN(I)/PLOTS AVG(I,8) = XCN(I)/PLOTS AVG(I,9) = XSCO2C(I)/PLOTS AVG(I,10) = XSOM(I)/PLOTS AVG(I,11) = XET(I)/PLOTS DO 95 III = 1,BMSPEC AVGBM(I,III) = XBMSPE(III,I)/PLOTS 95 CONTINUE C..... C.....BYPASS CONFIDENCE INTERVAL CALCULATIONS IF ONLY ONE PLOT C..... IF(PLOTS.EQ.1.) GO TO 40 C..... C.....FIND 95% CONFIDENCE INTERVALS C..... 51 CONF(I,1) = NYR(I) 52 Y2=((SSTEM(I)-(XSTEM(I)**2)/PLOTS)/(PLOTS-1.))/ &PLOTS IF(Y2.LT..001) Y2=0.0 CONF(I,2)=TS*SQRT(Y2) 53 Y3=((STAB(I)-(XTAB(I)**2)/PLOTS)/(PLOTS-1.))/ &PLOTS IF(Y3.LT..001) Y3=0.0 CONF(I,3)=TS*SQRT(Y3) 54 Y4=((SFL(I)-(XFL(I)**2)/PLOTS)/(PLOTS-1.))/PLOTS IF(Y4.LT..001) Y4=0.0 CONF(I,4)=TS*SQRT(Y4) 55 Y5=((STL(I)-(XTL(I)**2)/PLOTS)/(PLOTS-1.))/PLOTS IF(Y5.LT..001) Y5=0.0 CONF(I,5)=TS*SQRT(Y5) 56 Y6=((SNAP(I)-(XNAP(I)**2)/PLOTS)/(PLOTS-1.))/ &PLOTS IF(Y6.LT..001) Y6=0.0 CONF(I,6)=TS*SQRT(Y6) 57 Y7=((SAVLN(I)-(XAVLN(I)**2)/PLOTS)/(PLOTS-1.))/ &PLOTS IF(Y7.LT..001) Y7=0.0 CONF(I,7)=TS*SQRT(Y7) 58 Y8=((SCN(I)-(XCN(I)**2)/PLOTS)/(PLOTS-1.))/PLOTS IF(Y8.LT..001) Y8=0.0 CONF(I,8)=TS*SQRT(Y8) 59 Y9=((SSCO2C(I)-(XSCO2C(I)**2)/PLOTS)/(PLOTS-1.))/ &PLOTS IF(Y9.LT..001) Y9=0.0 CONF(I,9)=TS*SQRT(Y9) 41 Y10=((SSOM(I)-(XSOM(I)**2)/PLOTS)/(PLOTS-1.))/ &PLOTS IF(Y10.LT..001) Y10=0.0 CONF(I,10)=TS*SQRT(Y10) 42 Y11=((SET(I)-(XET(I)**2)/PLOTS)/(PLOTS-1.))/PLOTS IF(Y11.LT.001) Y11=0.0 CONF(I,11)=TS*SQRT(Y11) DO 96 III = 1,BMSPEC YBM = ((SSBMSP(III,I)-(XBMSPE(III,I)**2)/PLOTS)/ & (PLOTS-1.))/PLOTS IF(YBM.LT.0.001) YBM=0.0 CONFBM(I,III)=TS*SQRT(YBM) 96 CONTINUE 40 CONTINUE C..... C.....WRITE STATISTICS TO PRINTER C..... WRITE(6,3003) WRITE(6,3001) 3001 FORMAT(' '/' YR',5X,'NUM',5X,'A.G.',4X,'LEAF',4X,'LEAF ',3X, &'A.G.',4X,'AVAIL',3X,'HUMUS',3X,'SOIL',4X,'SOIL',4X,'AET') WRITE(6,3002) 3002 FORMAT(' ',7X,'STEMS',2X,'BIOMASS',2X,'LITTER',2X,'LITR N',3X, &'NPP',7X,'N',6X,'C:N',4X,'CO2-C',3X,'O.M.') WRITE(6,3003) 3003 FORMAT(' '/' --------------------------------------------------', &'---------------------------------'/) WRITE(6,3004) KLAST 3004 FORMAT(' '/' ',26X,' AVERAGE ACROSS ',I3,' PLOTS') WRITE(6,3005) 3005 FORMAT(' ',8X,'------------------------------------------------', &'---------------------------'/) WRITE(6,3000) ((AVG(I,J),J=1,11),I=1,NMAX) IF(PLOTS.EQ.1.) GO TO 99 WRITE(6,3006) 3006 FORMAT(' '/' ',26X,' 95% CONFIDENCE INTERVALS') WRITE(6,3005) WRITE(6,3000) ((CONF(I,J),J=1,11),I=1,NMAX) 99 CONTINUE WRITE(6,3005) WRITE(6,3095) 3095 FORMAT(' '/' ',31X,' SPECIES BIOMASS') WRITE(6,3096) (SPEC(I), I=1,BMSPEC) 3096 FORMAT(' '/' YR',1X,I6,9(2X,I6)) WRITE(6,3003) WRITE(6,3004) KLAST WRITE(6,3005) DO 97 I = 1,NMAX WRITE(6,3090) AVG(I,1), (AVGBM(I,III), III=1,BMSPEC) 97 CONTINUE IF(PLOTS.EQ.1.) GO TO 60 3090 FORMAT(F5.0,10F8.0) WRITE(6,3006) WRITE(6,3005) DO 98 I = 1,NMAX WRITE(6,3090) AVG(I,1), (CONFBM(I,III), III=1,BMSPEC) 98 CONTINUE 60 CONTINUE 3000 FORMAT(F5.0,2F8.0,3F8.1,F8.0,F8.0,F8.0,F8.1,F8.0) RETURN END C C SUBROUTINE PLOTIN C C C SUBROUTINE PLOTIN INITIALIZES VARIABLES TO START SIMULATION ON C BARE PLOTS. C NTREES CONTAINS NUMBER OF TREES FOR EACH SPECIES C DBH CONTAINS DIAMETER AT BREAST HEIGHT FOR EACH TREE C KSPRT IS USED TO FLAG TREES ELIGIBLE TO SPROUT C NOGRO IS USED TO FLAG SLOWLY GROWING TREES C IAGE CONTAINS THE AGE OF EACH TREE C C CONTAINS DATA ON LITTER COHORTS C C(1,1) IS HUMUS WEIGHT (BASESC IS STARTING VALUE) C C(1,2) IS HUMUS N CONTENT (BASESN IS STARTING VALUE) C C(I,5) IS "LITTER TYPE" FOR HUMUS C NCOHRT IS NUMBER OF LITTER COHORTS PRESENT C (IF NOCOHRT = 1, ONLY HUMUS IS PRESENT) C TYL CONTAINS THIS YEAR'S LITTER C SUBROUTINE PLOTIN(IPLOT) COMMON/FOREST/NTREES(100),DBH(1500),IAGE(1500),KSPRT(100), > NEWTR(100),SUMLA(700),NEW(100),SWTCH(5) COMMON/CONST/NSPEC,DEGD COMMON/DEAD/NOGRO(1500),NTEMP(1500) COMMON/DCMP/AVAILN,TYL(20),C(100,15),FDAT(20,10),FF(20,3), 1NCOHRT,HCN,BASESC,BASESN,SCO2,TYLN IPLOT = IPLOT+1 DO 10 I=1,NSPEC NTREES(I) = 0 DBH(I) = 0. NOGRO(I) = 0 KSPRT(I) = 0 IAGE(I) = 0 10 CONTINUE NSPE1 = NSPEC+1 DO 20 I=NSPE1,1500 DBH(I) = 0. NOGRO(I) = 0 IAGE(I) = 0 20 CONTINUE DO 30 I=1,100 DO 30 J=1,15 C(I,J)=0.0 30 CONTINUE C(1,1)=BASESC C(1,2)=BASESN C(1,5)=18. NCOHRT=1 DO 40 I=1,20 TYL(I)=0.0 40 CONTINUE RETURN END C C SUBROUTINE TEMPE C C C TEMPE CALCULATES GROWING SEASON DEGREE DAYS (DEGD) BASED ON C MONTHLY TEMPERATURES NORMALLY DISTRIBUTED AROUND A SPECIFIED C MEAN WITH A SPECIFIED STND DEV. THE TEMPERATURES ARE SUPPLIED C BY SUBROUTINE GGNORD USING RANDOM NUMBER GENERATOR URAND AND C ARE LINERALY INTERPOLATED BETWEEN YEARS OF DIFFERENT CLIMATES C BY SUBROUTINE LININT C SUBROUTINE TEMPE(DEGD,KYR) COMMON/WATER/T(12),VT(12),RT(12),R(12),VR(12),FC,DRY,BGS,EGS,PLAT, 1FJ,AET COMMON/SEED/USEED(15) DIMENSION DAYS(12),Z(2) INTEGER USEED DATA DAYS/31.,28.,31.,30.,31.,30.,31.,31.,30.,31.,30.,31./ DATA NCT/0/ YR = FLOAT(KYR) C..... C.....DDBASE IS TEMPERATURE (C) ABOVE WHICH DEGREE DAYS (DEGD) C.....ARE COUNTED C..... DDBASE = 5.56 C..... C.....INITIALIZE DEGREE DAYS AND MONTHLY TEMPERATURES USED TO C.....CALCULATE THEM (IN ARRAY RT) C..... DEGD = 0. DO 3 I=1,12 3 RT(I) = 0. C..... C.....CALL SUBROUTINE LININT TO MAKE LINEAR INTERPOLATIONS OF MONTHLY C.....TEMPERATURES AND STND DEVS BETWEEN YEARS OF DIFFERENT CLIMATE. C..... CALL LININT(T,VT,YR,1) DO 1 I=1,12 NCT = NCT+1 IF(NCT .EQ. 2) GO TO 36 C..... C.....CALL GGNORD TO PROVIDE NORMALLY DISTRIBUTED RANDOM NUMBERS C..... CALL GGNORD(USEED(12),USEED(13),Z) GO TO 38 36 Z(1) = Z(2) NCT = 0 38 CONTINUE C..... C.....CALCULATE MONTHLY TEMPERATURES (INTERPOLATED MEAN +/- C.....NORMALLY DISTRIBUTED RANDOM NUMBER TIMES INTERPOLATED C.....STND DEV C..... RT(I) = T(I)+VT(I)*Z(1) IF(RT(I) .LE. DDBASE) GO TO 1 C..... C.....SUM DEGREE DAYS FOR CONSECUTIVE MONTHS C..... DEGD = DEGD+(RT(I)-DDBASE)*DAYS(I) 1 CONTINUE RETURN END C C FUNCTION URAND C C DUMMY CALL FOR TRANSITION TO CALLING RANGEN(R) FOR VAX C REAL FUNCTION URAND(IY) INTEGER IY REAL R DATA IFRST/1/ C R = 0.0 C IF(IFRST .EQ. 1) R = .1 IFRST = 0 C URAND = RANGEN(R) URAND = RAN(IY) RETURN END