TNO Intern

Commit e7b744c4 authored by Florian Knappers's avatar Florian Knappers
Browse files

remove coarsen dataset

parent b93723a1
Loading
Loading
Loading
Loading
Loading
+3 −26
Original line number Diff line number Diff line
@@ -8,15 +8,14 @@ from pythermogis.workflow.utc.configuration import UTCConfiguration
from pythermogis.workflow.utc.doublet import DoubletInput, calculate_doublet_performance
from pythermogis.transmissivity.calculate_thick_perm_trans import _tpkt_kernel

def run_utc_workflow(config: UTCConfiguration, coarsen_dataset_by: int | None = None ) -> xr.DataTree:
def run_utc_workflow(config: UTCConfiguration) -> xr.DataTree:
    root = xr.DataTree()
    base = Path(config.input_data_dir)

    temperature_model = xr.open_dataset(base / "NetherlandsTemperatureModel_extended.nc")
    temperature_model = temperature_model.transpose("y", "x", "z")

    if coarsen_dataset_by is not None:
        temperature_model = temperature_model.coarsen(x=coarsen_dataset_by, y=coarsen_dataset_by, boundary="trim").mean()
    chunk_size = {'x': 200, 'y': 200}


    for aquifer in config.aquifers:
@@ -45,9 +44,6 @@ def run_utc_workflow(config: UTCConfiguration, coarsen_dataset_by: int | None =

            aquifer_ds[gridinfo.name] = ds[varname]

        if coarsen_dataset_by is not None:
            aquifer_ds = aquifer_ds.coarsen(x=coarsen_dataset_by, y=coarsen_dataset_by, boundary="trim").mean()

        root[aquifer] = xr.DataTree(dataset=aquifer_ds, name=aquifer)


@@ -100,8 +96,6 @@ def compute_results_for_aquifer(
    ]

    def cell_calculation(
            x,
            y,
            unknown_input_value,
            thickness,
            transmissivity,
@@ -141,22 +135,9 @@ def compute_results_for_aquifer(
        try:
            result = calculate_doublet_performance(config, doublet_input)
        except ValueError:
            print(
                f"Failure at x={x:.1f}, y={y:.1f} | "
                f"thick={thickness:.3f}, temp={temperature:.2f}, "
                f"trans={transmissivity:.3f}, trans_ntg={transmissivity_with_ntg:.3e}"
                f"ntg={ntg}, depth={depth:.3f}, porosity={porosity:.3f}"
            )
            return tuple(np.nan for _ in output_names)

        if result.utc > 1000:
            print(
                f"Hmmmmm at x={x:.1f}, y={y:.1f} | "
                f"thick={thickness:.3f}, temp={temperature:.2f}, "
                f"trans={transmissivity:.3f}, trans_ntg={transmissivity_with_ntg:.3e}, "
                f"ntg={ntg}, depth={depth:.3f}, porosity={porosity:.3f}, "
                f"utc:{result.utc}"
            )
            result = result._replace(utc=1000)


@@ -186,8 +167,6 @@ def compute_results_for_aquifer(
    )

    aquifer_ds["temperature"] = temperature
    y_coord, x_coord = xr.broadcast(aquifer_ds["y"], aquifer_ds["x"])


    valid_mask = (
            aquifer_ds["thickness"].notnull()
@@ -218,8 +197,6 @@ def compute_results_for_aquifer(

    results = xr.apply_ufunc(
        cell_calculation,
        x_coord,
        y_coord,
        1.0e30,
        aquifer_ds["thickness"],
        aquifer_ds["transmissivity"],
@@ -228,7 +205,7 @@ def compute_results_for_aquifer(
        aquifer_ds["depth"],
        aquifer_ds["porosity"],
        aquifer_ds["temperature"],
        input_core_dims=[[]] * 10,
        input_core_dims=[[]] * 8,
        output_core_dims=[[] for _ in output_names],
        vectorize=True,
        dask="parallelized",