'PRODCURV.BAS--program to generate reasonable NPP, mortality,
'biomass, and CWD curves from long-term forest measurements
'including coarse woody debris
'S. Acker, Dept. of Forest Science, Oregon State Univ.
'Jan. 1999

CLS
INPUT "Enter name for output file (complete path) "; outfil$
OPEN outfil$ FOR OUTPUT AS #1
INPUT "Enter 1st NPP value, age"; p5, t5
INPUT "Enter biomass and age at first measurement"; biom0, t0
int1 = t5 - t0
INPUT "Enter lag (years) before appreciable bole production "; lag
INPUT "Enter lag (years) before appreciable bole mortality"; mlag
INPUT "Enter decay rate for coarse woody debris "; dkrate
INPUT "Enter target mortality rate and age ", targmort, ttarg
INPUT "Treat 1st NPP value as maximum? (Y/N)"; p5max$

endyr = t0
IF t0 MOD 5 <> 0 THEN
    endyr = t0 - (t0 MOD 5) + 5
END IF

IF UCASE$(p5max$) = "N" THEN
    INPUT "Enter maximum NPP value"; pmax
    INPUT "Enter year of max NPP"; tmax
ELSE
    pmax = p5
    tmax = t5
END IF
PRINT "Enter name of file with observed of NPPB, mortality, biomass"
INPUT "(complete path) ", obsfil$
PRINT "Indicate order of observed values in data file:"
PRINT "1 = age, production, mortality, biomass"
PRINT "2 = age, biomass, mortality, CWD, production"
INPUT infmt

'fit parabola for ascending limb of NPP time trend
a = pmax / ((2 * tmax * lag) - tmax ^ 2 - lag ^ 2)
b = (-2 * tmax * pmax) / ((2 * tmax * lag) - tmax ^ 2 - lag ^ 2)
c = -(a * lag ^ 2 + b * lag)

biom = 0
mort = 0
cwd = 0

'print parameters of model run
PRINT #1, "output of PRODCURV.BAS", DATE$
PRINT #1, ""
PRINT #1, "Observed data from " + obsfil$
PRINT #1, "Initial obs. NPP = " + STR$(p5) + " at age = " + STR$(t5)
PRINT #1, "Initial obs. biomass = " + STR$(biom0) + " at age = " + STR$(t0)
PRINT #1, "Lag before appreciable NPPB = " + STR$(lag)
PRINT #1, "Lag before appreciable mortality =" + STR$(mlag)
PRINT #1, "Decay rate of CWD = " + STR$(dkrate)
PRINT #1, "Target mortality rate = " + STR$(targmort) + " at age = " + STR$(ttarg)
PRINT #1, "Initial obs. NPPB used as maximum? " + p5max$
IF UCASE$(p5max$) = "N" THEN
    PRINT #1, "Maximum NPPB = " + STR$(pmax) + " at age = " + STR$(tmax)
END IF
PRINT #1,

PRINT "Year", "NPPB  "; "MORT  "; "BIOMASS "; "CWD"
PRINT #1, "Year", "NPPB  "; "MORT  "; "BIOMASS "; "CWD"

FOR i = 0 TO endyr STEP 5
REM if max is before 1st measurement, decrease nppb linearly
    IF i > tmax AND p5max$ = "N" THEN
        nppb = pmax - (pmax - p5) * (i - tmax) / (t5 - tmax)
    ELSE
        nppb = a * i ^ 2 + b * i + c
    END IF
    IF nppb < .00001 THEN nppb = 0
'no mortality until age=mlag, then linear increase to target mort rate
'at specified age, based on plot data
    mort = (targmort / (ttarg - mlag)) * (i - mlag)
    IF mort < 0 THEN mort = 0
    'accomodate initial yrs not divisible by 5
    IF i > t0 AND t0 MOD 5 <> 0 THEN
        intlen = t0 MOD 5
        biom = biom + (nppb * intlen) - (mort * intlen)
        cwd = (cwd * (1 - dkrate) ^ intlen) + (mort * intlen)
        year = t0
    ELSE
        biom = biom + nppb * 5 - mort * 5
        cwd = (cwd * (1 - dkrate) ^ 5) + (mort * 5)
        year = i
    END IF
    PRINT year, : PRINT USING "##.##"; nppb; : PRINT USING "##.##"; mort;
    PRINT USING "#####.#"; biom; : PRINT USING "#####.#"; cwd
    PRINT #1, year, : PRINT #1, USING "##.##"; nppb;
    PRINT #1, USING "##.##"; mort;
    PRINT #1, USING "#####.#"; biom; : PRINT #1, USING "#####.#"; cwd
NEXT i

'read in observed data, starting with first year with production
'and mortality rates (i.e., 2nd measurement year);
'format of file--descriptive text on lines 1 & 2,
'number of measurement years in file on line 3,
'age,NPPB,mortality,biomass one line per measurement year

OPEN obsfil$ FOR INPUT AS #2
LINE INPUT #2, description$
LINE INPUT #2, more.info$
INPUT #2, numobs
FOR i = 1 TO numobs
    IF infmt = 1 THEN
        INPUT #2, age(i), nppb(i), mort(i), biom(i)
    ELSEIF infmt = 2 THEN
        INPUT #2, age(i), biom(i), mort(i), dead, nppb(i)
    END IF
NEXT i
CLOSE #2
FOR i = 1 TO numobs
    IF i = 1 THEN
            begin = t0
    ELSE
            begin = age(i - 1)
    END IF
    intlen = age(i) - begin
    cwd = (cwd * (1 - dkrate) ^ intlen) + (mort(i) * intlen)
    PRINT age(i), : PRINT USING "##.##"; nppb(i);
    PRINT USING "##.##"; mort(i);
    PRINT USING "#####.#"; biom(i); : PRINT USING "#####.#"; cwd
    PRINT #1, age(i), : PRINT #1, USING "##.##"; nppb(i);
    PRINT #1, USING "##.##"; mort(i);
    PRINT #1, USING "#####.#"; biom(i); : PRINT #1, USING "#####.#"; cwd
NEXT i


CLOSE #1