#include <math.h>
/*--------------------------------------------------------------------------*/
double opp_befa( double chl,
double irr,
double sst,
double dayL ) {
/*
!C--------------------------------------------------------------------------*\
!Description: opp_befa - computes daily primary productivity using
the Behrenfeld-Falkowski (BeFa) algorithm. The BeFa
algorithm estimates productivity using surface chl
(mg m-3), surface irradiance (Einsteins m-2 d-1),
sea surface temperature (C), and day length (hours).
Pb_opt is modelled as a polynomial function of SST.
!Input Parameters:
chl Chlorophyll_a surface concentration in milligrams
chlorophyl per cubic meter
irr Photosynthetically available radiation in Einsteins per
day per square meter
sst Sea surface temperature in degrees Centigrade
dayL Length day in decimal hours.
!Output Parameters:
<return> Primary productivity in milligrams Carbon per square meter
per day
!Revision History:
First programmed up by Monica Chen at Rutgers
(1996)
Revised by K. Turpie at NASA
(August 1997)
Maintained by Don Shea at NASA
Now maintained by Robert O'Malley at Oregon State University
(April, 2005 - present)
changing z_eu calc to match updated Morel equation
O'Malley 02.13.2023
!References
Behrenfeld,M.J; Falkowski,P.G.; 1997. Photosynthetic Rates Derived
from Satellite-Based Chlorophyll Concentration. Limnology and
Oceanography, Volume 42, Number 1
!END------------------------------------------------------------------------*\
*/
double chl_tot,
z_eu,
pb_opt,
irrFunc,
npp;
double X,
log10z_eu;
/* Calculate euphotic depth (z_eu) with Morel's 2007 model. */
/* Previously had used Morel's 1988 model */
X = log10(chl);
log10z_eu = 1.524 - 0.436*X - 0.0145*X*X + 0.0186*X*X*X;
z_eu = pow(10,log10z_eu);
/* Calculate the Pb_opt from satellite sea surface temperature (sst). */
if (sst < -10.0L)
pb_opt = 0.00L;
else if (sst < -1.0L)
pb_opt = 1.13L;
else if (sst > 28.5L)
pb_opt = 4.00L;
else {
pb_opt = 1.2956 + 2.749e-1*sst + 6.17e-2*pow(sst,2) - 2.05e-2*pow(sst, 3)
+ 2.462e-3*pow(sst,4) - 1.348e-4*pow(sst,5) + 3.4132e-6*pow(sst,6)
- 3.27e-8*pow(sst,7);
}
/* calculate the irradiance function */
irrFunc = 0.66125L * irr / ( irr + 4.1L );
/* Return the primary production calculation. */
npp = pb_opt * chl * dayL * irrFunc * z_eu;
return npp;
}
|