Rev 554 | Rev 562 | Go to most recent revision | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed
| Rev 554 | Rev 561 | ||
|---|---|---|---|
| Line 352... | Line 352... | ||
| 352 | {
|
352 | {
|
| 353 | double vpd_mbar = climate->vpd * 10.; // convert from kPa to mbar |
353 | double vpd_mbar = climate->vpd * 10.; // convert from kPa to mbar |
| 354 | double temperature = climate->temperature; // average temperature of the day (°C) |
354 | double temperature = climate->temperature; // average temperature of the day (°C) |
| 355 | double daylength = daylength_h * 3600.; // daylength in seconds (convert from length in hours) |
355 | double daylength = daylength_h * 3600.; // daylength in seconds (convert from length in hours) |
| 356 | double rad = climate->radiation / daylength * 1000000; //convert from MJ/m2 (day sum) to average radiation flow W/m2 [MJ=MWs -> /s * 1,000,000 |
356 | double rad = climate->radiation / daylength * 1000000; //convert from MJ/m2 (day sum) to average radiation flow W/m2 [MJ=MWs -> /s * 1,000,000 |
| - | 357 | ||
| - | 358 | // the radiation: based on linear empirical function
|
|
| - | 359 | const double qa = -90.; |
|
| - | 360 | const double qb = 0.8; |
|
| - | 361 | double net_rad = qa + qb*rad; |
|
| 357 | 362 | ||
| 358 | //: Landsberg original: const double e20 = 2.2; //rate of change of saturated VP with T at 20C
|
363 | //: Landsberg original: const double e20 = 2.2; //rate of change of saturated VP with T at 20C
|
| 359 | const double VPDconv = 0.000622; //convert VPD to saturation deficit = 18/29/1000 |
364 | const double VPDconv = 0.000622; //convert VPD to saturation deficit = 18/29/1000 |
| 360 | const double latent_heat = 2460000.; // Latent heat of vaporization. Energy required per unit mass of water vaporized [J kg-1] |
365 | const double latent_heat = 2460000.; // Latent heat of vaporization. Energy required per unit mass of water vaporized [J kg-1] |
| 361 | 366 | ||
| Line 367... | Line 372... | ||
| 367 | // current response: see calculateSoilAtmosphereResponse(). This is basically a weighted average of min(water_response, vpd_response) for each species
|
372 | // current response: see calculateSoilAtmosphereResponse(). This is basically a weighted average of min(water_response, vpd_response) for each species
|
| 368 | double gC = mAvgMaxCanopyConductance * combined_response; |
373 | double gC = mAvgMaxCanopyConductance * combined_response; |
| 369 | 374 | ||
| 370 | 375 | ||
| 371 | double defTerm = mAirDensity * latent_heat * (vpd_mbar * VPDconv) * gBL; |
376 | double defTerm = mAirDensity * latent_heat * (vpd_mbar * VPDconv) * gBL; |
| 372 | // saturation vapor pressure (Running 1988, Eq. 1) in mbar
|
- | |
| 373 | double svp = 6.1078 * exp((17.269 * temperature) / (237.3 + temperature) ); |
- | |
| 374 | // the slope of svp is, thanks to http://www.wolframalpha.com/input/?i=derive+y%3D6.1078+exp+((17.269x)/(237.3%2Bx))
|
- | |
| 375 | double svp_slope = svp * ( 17.269/(237.3+temperature) - 17.269*temperature/((237.3+temperature)*(237.3+temperature)) ); |
- | |
| - | 377 | ||
| - | 378 | // with temperature-dependent slope of vapor pressure saturation curve
|
|
| - | 379 | // (following Allen et al. (1998), http://www.fao.org/docrep/x0490e/x0490e07.htm#atmospheric%20parameters)
|
|
| - | 380 | // svp_slope in mbar.
|
|
| - | 381 | double svp_slope = 4098. * (6.1078 * exp(17.269 * temperature / (temperature + 237.3))) / ((237.3+temperature)*(237.3+temperature)); |
|
| 376 | 382 | ||
| 377 | double div = (1. + svp_slope + gBL / gC); |
383 | double div = (1. + svp_slope + gBL / gC); |
| 378 | double Etransp = (svp_slope * rad + defTerm) / div; |
- | |
| - | 384 | double Etransp = (svp_slope * net_rad + defTerm) / div; |
|
| 379 | double canopy_transpiration = Etransp / latent_heat * daylength; |
385 | double canopy_transpiration = Etransp / latent_heat * daylength; |
| 380 | 386 | ||
| 381 | // calculate PET
|
387 | // calculate PET
|
| 382 | double div_evap = 2. + svp_slope; |
- | |
| 383 | double pet_day = (svp_slope*rad + defTerm) / div_evap / latent_heat * daylength; |
- | |
| - | 388 | double div_evap = 1. + svp_slope; |
|
| - | 389 | double pet_day = (svp_slope*net_rad + defTerm) / div_evap / latent_heat * daylength; |
|
| 384 | mPET[climate->month-1] += pet_day; |
390 | mPET[climate->month-1] += pet_day; |
| 385 | 391 | ||
| 386 | if (mInterception>0.) { |
392 | if (mInterception>0.) { |
| 387 | // we assume that for evaporation from leaf surface gBL/gC -> 0
|
393 | // we assume that for evaporation from leaf surface gBL/gC -> 0
|
| 388 | pet_day = qMin(pet_day, mInterception); |
394 | pet_day = qMin(pet_day, mInterception); |