Add bodemlagen_pwn_2024 v2.0.0 dataset#50
Conversation
Add mockup dataset bodemlagen_pwn_2024 (v2.0.0): many GeoJSON layers for dikte_aquitard (DS11..DS32) and top_aquitard (TS11..TS32) including interpolation points, masks, union_with_values files, and masks/clip (noordzee). Also add helper scripts (interpolate_layer_boundaries.py, interpolation_helper_functions.py, merge_masks.py) and a report PDF to support interpolation and mask-merge workflows.
|
The TODO's we colklected last week are listed in src/nhflodata/data/mockup/bodemlagen_pwn_2024/v2.0.0/interpolate_layer_boundaries.py of this branch |
- Updated `interpolate_layer_boundaries.py` to import helper functions directly and removed unused PyVista plotting code. - Modified `get_point_values` function in `interpolation_helper_functions.py` to accept a path parameter for improved flexibility. - Updated paths for reading shapefiles and GeoJSON in `get_point_values`. - Added new version 2.0.0 entry in `repository.yaml` for bodemlagen_pwn_2024, including metadata and changelog for Koster contour interpretations.
…, and ultimately compute the bottoms - Created geojson files for boundaries S13, S21, S22, S31, and S32. - Added Python script to convert boundary files and generate geojson outputs. - Included zip file containing necessary shapefiles for boundaries. - Implemented new geojson files for triwaco models (bergen and nhdz). - Added initialization and merging scripts for aquitard thickness layers.
|
The C values for both Bergen en NHDZ contain overlapping polygons. How to interpret this? This definitely the case for C1C.shp (Bergen resistance of 1C). Have a look at the by fill_boundary_with_polygons() returned The logic for dealing with overlapping polygons and filling up open areas within the boundaries is located in the I set the overlap priority arbitrairily to last, so in case of an overlap I am now using the data from the last row of the dataframe. |
|
It is stated correct in the report, but confusing and worth a check: Bergen uses the base of the aquitards and NHDZ the top of the aquitards. Also just to confirm, triwaco used both polygons and Kriging to compute the elevation and thickness of the resistance layers. My code to compute the base of the aquitards is located here: def _read_bergen_basis_aquitards() In je rapport noem je op p35 "Voor SDP 1G en 1H is de constante dikte uit DI1C.shp aangehouden." Maar DI1C zou ook met Kriging berekend moeten worden. Terwijl je voor diepere lagen wel naar de punten hebt gekeken. |
|
Contrary to the previous implementation logic that if either the layer had zero thickness or the hydraulic conductivity was missing, the clay layer is not present. I now want to use the thickness (configured via the bottoms) only to judge the presence of a layer. Hydraulic conductivities are defined everywhere. This follows more the modflow logic and makes it easier to interpret visuals where you expect a resistance if you see a clay layer with a thickness.
|
Introduce a new kh package for the bodemlagen_pwn_2024 mockup: adds an empty __init__.py, kh.py (loads Bergen/NHDZ K-value sources, assembles layer values and combines fills) and params_helper_functions.py (utilities to extract polygons, compute area-weighted mean/median, and fill a boundary with non-overlapping polygons). The fill_boundary_with_polygons function supports overlap resolution strategies, fill methods (fill_value/mean/median), optional override polygons, and logs/validates overlapping inputs. This enables creating complete, non-overlapping parameter coverage for the boundaries.
|
It must be somewhere in you report, but it is not completely clear to me. You describe 5 aquitards in your description of the Bergen model and you merge only 4 of them to NHDZ. In your Bergen description, you write that 1D contains almost the same data as 1C. Also looking at the resistance polygons, the 1D layer has 1000 days resistance for all polygons. Conclusion: We omit 1D from the Bergen model |
kh: Rework kh processing to compute C and KD layers from NHDZ and Bergen sources. Added mask-based clipping of KD inputs and assembled KD formulas into summed parts before filling boundaries; write overlap and fill GeoJSONs for C layers. Adjusted Bergen constants and name mapping and changed several fill/overlap priorities (use 'last' and 'sum' where appropriate). params_helper_functions: Add _compute_pairwise_overlaps and clip_polygons_to_mask helpers, strengthen geometry buffering/validation, and extend fill_boundary_with_polygons to return (result, overlaps_gdf, fill_gdf), detect/report overlaps more robustly, and support an overlap_priority 'sum' mode that sums values in overlapping regions. Overall improves handling of overlapping geometries and mask-based clipping for kh generation.
|
Let's discuss, right now we have stored in the shape files
MODFLOW requires Kh, Kv for aquifers and aquitards. Thus we compute the remaining as follows:
For Bergen, the data we got is:
And we compute:
Sadly, we cannot precompute those values and distribute those via NHFLO/data as we need the MODFLOW grid first for correct interpolation. Linear interpolation does not work well with the inverse relation of the resistance. As I am writing this, it would work with the KD of the nhdz aquitards and all terms that just require a multiplication with an anisotropy factor. It is strange to have a aquitards become thinner and therefor get a lower conductance. But the resistance is something that is directly derived from a pumping test.. To the point: I would like Kv values for the aquitards, that are either homogeneous in polygons or linear interpolated to the MODFLOW grid. Vincent, what is your opinion on this? |
Right now I am merging the two models, assuming Koster did a better job at estimating the parameters where both models overlap. This should be improved. |
Add multiple mock conductance GeoJSON datasets for bodemlagen_pwn_2024 v2.0.0 under src/nhflodata/data/mockup/.../conductances. Also move/rename the Python package from kh/ to conductances/ (moved __init__.py and renamed kh.py -> conductances.py and params_helper_functions.py) to reflect the new package layout for conductances data.
Add __init__.py modules (with brief docstrings) for botm, boundaries, conductances, dikte_aquitard and koster_drilling_interpretations to make these subpackages importable. Fix path casing in interpolation helper (Koster_drilling_interpretations -> koster_drilling_interpretations) to avoid file-not-found issues. Apply small formatting and readability improvements: wrap long path strings, reorder/import formatting in conversion scripts, add a trailing newline, and adjust dict formatting in botm merge code. These changes are maintenance-focused to improve package structure and correctness for conversion/interpolation workflows.
Expand mockup conductances data for bodemlagen_pwn_2024 v2.0.0: KDS12_NHDZ.geojson was populated with numerous Point features
Document open concerns from the Edinsi report (Vincent Post, Aug 2024) as TODO comments in the relevant source files: borehole location differences, North Sea layer extent, layer boundary gaps at the NHDZ/Bergen transition, conductance anomalies, and recommendations for extending deeper layers northward and investigating Eem clay splitting.
|
The check is failing due to an incompatibity between hatch and virtualenv. This should be fixed by hatch maintainers. We wait instead of implementing a hotfix |
Even though we made a lot of progress, there is still work to be done. The comments in this PR are ones that I would like to discuss with you. I you search TODO: in the code, you'll find many of your remarks from your report. Maybe the major one maybe to start with, and the most complex one is thinking about how do we make the move from c-values to kv-values for aquitards. I'd like to see an aquitard to become more permeable if it becomes thinner, and now those two are not really connected (but now the Kv is reduced iin order to keep the c the same). |
|
A comma was missing here in your script. Did it become part of the interpretation/conversion? |
|
There is additional data on conductances and depths available for the DWAT region. I don't think it is currently included.
|
2 participants
Add mockup dataset bodemlagen_pwn_2024 (v2.0.0): many GeoJSON layers for dikte_aquitard (DS11..DS32) and top_aquitard (TS11..TS32) including interpolation points, masks, union_with_values files, and masks/clip (noordzee). Also add helper scripts (interpolate_layer_boundaries.py, interpolation_helper_functions.py, merge_masks.py) and a report PDF to support interpolation and mask-merge workflows.