Skip to content

Commit

Permalink
Updating docs to handle_nan
Browse files Browse the repository at this point in the history
  • Loading branch information
davidorme committed Jul 24, 2023
1 parent 18451a3 commit aa6147b
Show file tree
Hide file tree
Showing 3 changed files with 30 additions and 27 deletions.
2 changes: 1 addition & 1 deletion CHANGES.txt
Original file line number Diff line number Diff line change
Expand Up @@ -142,4 +142,4 @@
approaches are now intended to be applied as penalties to GPP after P Model
fitting, allowing the two to be compared from the same P Model outputs.
0.10.1 - Addition of draft extention of subdaily memory_effect function that allows for
missing data arising from non-estimable Jmax and Vcmax.
missing data arising from non-estimable Jmax and Vcmax.
40 changes: 21 additions & 19 deletions docs/source/users/pmodel/subdaily_details/subdaily_calculations.md
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ steps used in the estimation process in order to show intermediates results but
practice, as shown in the [worked example](worked_example), most of these calculations
are handled internally by the model fitting in `pyrealm`.

```{code-cell} ipython3
```{code-cell}
:tags: [hide-input]
from importlib import resources
Expand All @@ -46,16 +46,15 @@ The code below uses half hourly data from 2014 for the [BE-Vie FluxNET
site](https://fluxnet.org/doi/FLUXNET2015/BE-Vie), which was also used as a
demonstration in {cite:t}`mengoli:2022a`.

```{code-cell} ipython3
```{code-cell}
data_path = resources.files("pyrealm_build_data") / "subdaily_BE_Vie_2014.csv"
data = np.genfromtxt(
data_path,
names=True,
delimiter=",",
dtype=None,
encoding="UTF8",
missing_values = "NA",
missing_values="NA",
)
# Extract the key half hourly timestep variables
Expand All @@ -74,7 +73,7 @@ This dataset can then be used to calculate the photosynthetic environment at the
subdaily timescale. The code below also estimates GPP under the standard P Model with no
slow responses for comparison.

```{code-cell} ipython3
```{code-cell}
# Calculate the photosynthetic environment
subdaily_env = PModelEnvironment(
tc=temp_subdaily,
Expand All @@ -100,7 +99,7 @@ best to sample those conditions. Typically those might be the observed environme
conditions at the observation closest to noon, or the mean environmental conditions in a
window around noon.

```{code-cell} ipython3
```{code-cell}
# Create the fast slow scaler
fsscaler = FastSlowScaler(datetime_subdaily)
Expand All @@ -114,19 +113,18 @@ fsscaler.set_window(
pmodel_fastslow = FastSlowPModel(
env=subdaily_env,
fs_scaler=fsscaler,
handle_nan=True,
ppfd=ppfd_subdaily,
fapar=fapar_subdaily,
)
```

```{code-cell} ipython3
```{code-cell}
:tags: [hide-input]
idx = np.arange(48 * 120, 48 * 130)
plt.figure(figsize=(10, 4))
plt.plot(
datetime_subdaily[idx], pmodel_subdaily.gpp[idx], label="Instantaneous model"
)
plt.plot(datetime_subdaily[idx], pmodel_subdaily.gpp[idx], label="Instantaneous model")
plt.plot(datetime_subdaily[idx], pmodel_fastslow.gpp[idx], "r-", label="Slow responses")
plt.ylabel = "GPP"
plt.legend(frameon=False)
Expand All @@ -145,7 +143,7 @@ The daily average conditions during the acclimation window can be sampled and us
inputs to the standard P Model to calculate the optimal behaviour of plants under those
conditions.

```{code-cell} ipython3
```{code-cell}
# Get the daily acclimation conditions for the forcing variables
temp_acclim = fsscaler.get_daily_means(temp_subdaily)
co2_acclim = fsscaler.get_daily_means(co2_subdaily)
Expand Down Expand Up @@ -176,7 +174,7 @@ temperatures so $J_{max}$ and $V_{cmax}$ must first be standardised to expected
at 25°C. This is acheived by multiplying by the reciprocal of the exponential part of
the Arrhenius equation ($h^{-1}$ in {cite}`mengoli:2022a`).

```{code-cell} ipython3
```{code-cell}
# Are these any of the existing values in the constants?
ha_vcmax25 = 65330
ha_jmax25 = 43900
Expand All @@ -191,18 +189,18 @@ jmax25_acclim = pmodel_acclim.jmax * (1 / calc_ftemp_arrh(tk_acclim, ha_jmax25))
The memory effect can now be applied to the three parameters with slow
responses to calculate realised values, here using the default 15 day window.

```{code-cell} ipython3
```{code-cell}
# Calculation of memory effect in xi, vcmax25 and jmax25
xi_real = memory_effect(pmodel_acclim.optchi.xi, alpha=1 / 15)
vcmax25_real = memory_effect(vcmax25_acclim, alpha=1 / 15)
jmax25_real = memory_effect(jmax25_acclim, alpha=1 / 15)
vcmax25_real = memory_effect(vcmax25_acclim, alpha=1 / 15, handle_nan=True)
jmax25_real = memory_effect(jmax25_acclim, alpha=1 / 15, handle_nan=True)
```

The plots below show the instantaneously acclimated values for $J_{max25}$,
$V_{cmax25}$ and $\xi$ in grey along with the realised slow reponses.
applied.

```{code-cell} ipython3
```{code-cell}
:tags: [hide-input]
fig, axes = plt.subplots(1, 3, figsize=(16, 5))
Expand Down Expand Up @@ -236,7 +234,7 @@ temperature at fast scales:
* These values are adjusted to the actual half hourly temperatures to give the fast
responses of $J_{max}$ and $V_{cmax}$.

```{code-cell} ipython3
```{code-cell}
tk_subdaily = subdaily_env.tc + pmodel_subdaily.const.k_CtoK
# Fill the realised jmax and vcmax from subdaily to daily
Expand All @@ -256,7 +254,7 @@ $\Gamma^\ast$ with the realised slow responses of $\xi$. The original implementa
$\Gamma^{\ast}$ and $c_a$, interpolated to the subdaily timescale and the actual
subdaily variation in VPD.

```{code-cell} ipython3
```{code-cell}
# Interpolate xi to subdaily scale
xi_subdaily = fsscaler.fill_daily_to_subdaily(xi_real)
Expand All @@ -273,7 +271,7 @@ Model, where $c_i$ includes the slow responses of $\xi$ and $V_{cmax}$ and $J_{m
include the slow responses of $V_{cmax25}$ and $J_{max25}$ and fast responses to
temperature.

```{code-cell} ipython3
```{code-cell}
# Calculate Ac
Ac_subdaily = (
vcmax_subdaily
Expand All @@ -300,3 +298,7 @@ GPP_subdaily = np.minimum(Ac_subdaily, Aj_subdaily) * PModelConst.k_c_molmass
diff = GPP_subdaily - pmodel_fastslow.gpp
print(np.nanmin(diff), np.nanmax(diff))
```

```{code-cell}
```
15 changes: 8 additions & 7 deletions docs/source/users/pmodel/subdaily_details/worked_example.md
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@ focal_datetime = np.where(datetimes == np.datetime64("2018-06-12 12:00:00"))[0]
# Plot the temperature data for an example timepoint and show the sites
focal_temp = ds["Tair"][focal_datetime] - 273.15
focal_temp.plot()
plt.plot(sites["lon"], sites["lat"], "xr");
plt.plot(sites["lon"], sites["lat"], "xr")
```

The WFDE data need some conversion for use in the PModel, along with the definition of
Expand Down Expand Up @@ -119,7 +119,12 @@ fsscaler.set_window(
half_width=np.timedelta64(1, "h"),
)
fs_pmod = FastSlowPModel(
env=pm_env, fs_scaler=fsscaler, fapar=fapar, ppfd=ppfd, alpha=1 / 15
env=pm_env,
fs_scaler=fsscaler,
handle_nan=True,
fapar=fapar,
ppfd=ppfd,
alpha=1 / 15,
)
```

Expand Down Expand Up @@ -161,7 +166,7 @@ ax2.set_title("Fast Slow")
# Add a colour bar
subfig2.colorbar(
im, ax=[ax1, ax2], shrink=0.55, label=r"GPP ($\mu g C\,m^{-2}\,s^{-1}$)"
);
)
```

## Time series predictions
Expand Down Expand Up @@ -200,7 +205,3 @@ for ax, st in zip(axes, sites["stid"].values):
axes[0].legend(loc="lower center", bbox_to_anchor=[0.5, 1], ncols=2, frameon=False)
plt.tight_layout()
```

```{code-cell}
```

0 comments on commit aa6147b

Please sign in to comment.