TNO Intern

Commit 234306ba authored by Hen Brett's avatar Hen Brett 🐔
Browse files

updating the examples and the documentation

parent edd7b295
Loading
Loading
Loading
Loading
+20 −21
Original line number Diff line number Diff line
@@ -32,10 +32,9 @@ This example corresponds to test case `test_example5` in `test_doc_examples.py`
directory of the repository.

```python
from pythermogis import calculate_doublet_performance, calculate_pos, calculate_pos_pvalues_singlelocation
from pythermogis import calculate_doublet_performance, calculate_pos
from pygridsio import read_grid, plot_grid, resample_xarray_grid
from matplotlib import pyplot as plt
import numpy as np
import xarray as xr
from pathlib import Path
from os import path
@@ -48,7 +47,7 @@ output_data_path = Path(path.dirname(__file__), "resources") / "test_output" / "
output_data_path.mkdir(parents=True, exist_ok=True)

# if set to True then simulation is always run, otherwise pre-calculated results are read (if available)
run_simulation = True
run_simulation = False

# grids can be in .nc, .asc, .zmap or .tif format: see pygridsio package
new_cellsize = 5000  # in m, this sets the resolution of the model; to speed up calcualtions you can increase the cellsize
@@ -80,9 +79,11 @@ plt.tight_layout() # ensure there is enough spacing
plt.savefig(output_data_path / "maps.png")

# plot net-present value at a single location as a function of p-value and find the probability of success
results["pos"] = calculate_pos(results) # calculate probability of success across whole aquifer
x, y = 125e3, 525e3  # define location
results_loc = results.sel(x=x, y=y, method="nearest")
pos = calculate_pos_pvalues_singlelocation(results_loc.npv, results_loc.p_value)
results_loc = results.sel(x=x,y=y,method="nearest") # select only the location of interest
pos = results_loc.pos.data   # get probability of success at location of interest as a single value

# plot npv versus p-value and a map showing the location of interest
fig, axes = plt.subplots(ncols=2, figsize=(10,5))
results_loc.npv.plot(y="p_value", ax=axes[0])
@@ -93,13 +94,11 @@ axes[0].legend()
axes[0].set_xlabel("net-present-value [Million €]")
axes[0].set_ylabel("p-value [%]")

pos = calculate_pos(results)
plot_grid(pos, axes=axes[1], add_netherlands_shapefile=True)
plot_grid(results.pos,axes=axes[1], add_netherlands_shapefile=True)
axes[1].scatter(x,y,marker="x",color="tab:red")

plt.tight_layout()  # ensure there is enough spacing
plt.savefig(output_data_path / "pos.png")

```

---
+64 −70
Original line number Diff line number Diff line
@@ -35,12 +35,11 @@ directory of the repository.

```python
from pythermogis import calculate_doublet_performance, calculate_pos_pvalues_singlelocation
from pygridsio import read_grid, plot_grid, resample_xarray_grid
from pygridsio import read_grid
from matplotlib import pyplot as plt
import xarray as xr
from numpy import np
from pathlib import Path
from os import path

input_data_path = Path(path.dirname(__file__), "resources") / "test_input" / "example_data"
# create a directory to write the output files to
output_data_path = Path(path.dirname(__file__), "resources") / "test_output" / "example_data"
@@ -54,10 +53,10 @@ reservoir_properties["depth"] = read_grid(input_data_path / "ROSL_ROSLU__top.zma
reservoir_properties["ln_permeability_mean"] = np.log( read_grid(input_data_path / "ROSL_ROSLU__perm.zmap"))
reservoir_properties["ln_permeability_sd"] =  read_grid(input_data_path / "ROSL_ROSLU__ln_perm_sd.zmap")

x, y = 125e3, 525e3  # define location
portfolioloc = [(125e3, 525e3), (100e3, 525e3), (110e3, 525e3), (125e3, 515e3),
              (125e3, 520e3)]  # define location


fig, axes = plt.subplots(ncols=2, nrows=2, figsize=(10, 10))
colors = plt.cm.tab10.colors
ic = 0
@@ -74,12 +73,10 @@ for x, y in portfolioloc:
                                                 results_loc.p_value)  # The order of the npv data has to be reversed as np.interp requires values to increase to function properly
   
   # plot npv versus p-value and a map showing the location of interest

   axe = axes[0, 0]
   results_loc.npv.plot(y="p_value", ax=axe, color=colors[ic])
   axe.set_title(f"Aquifer: ROSL_ROSLU\n  POS ")
   axe.axhline(pos, label=f"POS: {pos:.1f}% [ {x:.0f}m, {y:.0f}m]", ls="--", c=colors[ic])
    # axes.axhline(pos, label=f"probability of success: {pos:.1f}%",ls="--", c="tab:orange")
   axe.axvline(0.0, ls="--", c=colors[ic])
   axe.legend()
   axe.set_xlabel("net-present-value [Million €]")
@@ -94,7 +91,6 @@ for x, y in portfolioloc:
   prd_temp = results_loc.prd_temp.sel(p_value=50, method="nearest")
   axe.axhline(50.0, label=f"POS: {pos:.1f}%  TEMP res: {temp:.1f} inj: {inj_temp:.1f} prd: {prd_temp:.1f} ",
              ls="--", c=colors[ic])
    # axes.axhline(pos, label=f"probability of success: {pos:.1f}%",ls="--", c="tab:orange")
   kh50 = kh.sel(p_value=50, method="nearest")
   axe.axvline(kh50, ls="--", c=colors[ic])
   axe.legend()
@@ -105,7 +101,6 @@ for x, y in portfolioloc:
   results_loc.power.plot(y="p_value", ax=axe, color=colors[ic])
   axe.set_title(f"Aquifer: ROSL_ROSLU\n Power")
   axe.axhline(pos, label=f"POS: {pos:.1f}% [ {x:.0f}m, {y:.0f}m]", ls="--", c=colors[ic])
    # axes.axhline(pos, label=f"probability of success: {pos:.1f}%",ls="--", c="tab:orange")
   valp50 = results_loc.power.sel(p_value=50, method="nearest")
   axe.axvline(valp50, ls="--", c=colors[ic])
   axe.legend()
@@ -116,7 +111,6 @@ for x, y in portfolioloc:
   results_loc.cop.plot(y="p_value", ax=axe, color=colors[ic])
   axe.set_title(f"Aquifer: ROSL_ROSLU\n COP")
   axe.axhline(pos, label=f"POS: {pos:.1f}% [ {x:.0f}m, {y:.0f}m]", ls="--", c=colors[ic])
    # axes.axhline(pos, label=f"probability of success: {pos:.1f}%",ls="--", c="tab:orange")
   valp50 = results_loc.cop.sel(p_value=50, method="nearest")
   axe.axvline(valp50, ls="--", c=colors[ic])
   axe.legend()
−524 KiB
Loading image diff...
−213 KiB
Loading image diff...
−3.29 MiB

File deleted.

Loading