diff --git a/.build/build.spec b/.build/build.spec index a702bec1..301d6ce5 100644 --- a/.build/build.spec +++ b/.build/build.spec @@ -31,7 +31,7 @@ a = Analysis( pathex=[Path(build_dir), Path(project_root, "src")], binaries=[], datas=[], - hiddenimports=["fiat_build_time"], + hiddenimports=["fiat_build_time", "fiat.methods"], hookspath=[build_dir.as_posix()], hooksconfig={}, runtime_hooks=[Path(build_dir, 'runtime_hooks.py')], diff --git a/.github/workflows/docs.yml b/.github/workflows/docs.yml index 0236627f..f5439805 100644 --- a/.github/workflows/docs.yml +++ b/.github/workflows/docs.yml @@ -77,7 +77,7 @@ jobs: sudo apt install -y curl sudo apt install -y gdebi-core mkdir tmp - curl -L https://github.com/quarto-dev/quarto-cli/releases/download/v1.3.450/quarto-1.3.450-linux-amd64.deb --output tmp/quarto.deb + curl -L https://github.com/quarto-dev/quarto-cli/releases/download/v1.5.57/quarto-1.5.57-linux-amd64.deb --output tmp/quarto.deb chmod +x tmp/quarto.deb sudo gdebi -n tmp/quarto.deb diff --git a/.testdata/create_test_data.py b/.testdata/create_test_data.py index 4791e959..e3194d1d 100644 --- a/.testdata/create_test_data.py +++ b/.testdata/create_test_data.py @@ -31,8 +31,8 @@ def create_dbase_stucture(): def create_exposure_dbase(): """_summary_.""" with open(Path(p, "exposure", "spatial.csv"), "w") as f: - f.write("Object ID,Extraction Method,Ground Floor Height,Ground Elevation,") - f.write("Damage Function: Structure,Max Potential Damage: Structure\n") + f.write("object_id,extract_method,ground_flht,ground_elevtn,") + f.write("fn_damage_structure,max_damage_structure\n") for n in range(5): if (n + 1) % 2 != 0: dmc = "struct_1" @@ -63,13 +63,13 @@ def create_exposure_geoms(): ) field = ogr.FieldDefn( - "Object ID", + "object_id", ogr.OFTInteger, ) layer.CreateField(field) field = ogr.FieldDefn( - "ObjectName", + "object_name", ogr.OFTString, ) field.SetWidth(50) @@ -78,8 +78,8 @@ def create_exposure_geoms(): for idx, geom in enumerate(geoms): geom = ogr.CreateGeometryFromWkt(geom) ft = ogr.Feature(layer.GetLayerDefn()) - ft.SetField("Object ID", idx + 1) - ft.SetField("ObjectName", f"fp_{idx+1}") + ft.SetField("object_id", idx + 1) + ft.SetField("object_name", f"fp_{idx+1}") ft.SetGeometry(geom) layer.CreateFeature(ft) @@ -110,13 +110,13 @@ def create_exposure_geoms_2(): ) field = ogr.FieldDefn( - "Object ID", + "object_id", ogr.OFTInteger, ) layer.CreateField(field) field = ogr.FieldDefn( - "ObjectName", + "object_name", ogr.OFTString, ) field.SetWidth(50) @@ -124,8 +124,8 @@ def create_exposure_geoms_2(): geom = ogr.CreateGeometryFromWkt(geoms[0]) ft = ogr.Feature(layer.GetLayerDefn()) - ft.SetField("Object ID", 5) - ft.SetField("ObjectName", f"fp_{5}") + ft.SetField("object_id", 5) + ft.SetField("object_name", f"fp_{5}") ft.SetGeometry(geom) layer.CreateFeature(ft) @@ -158,13 +158,13 @@ def create_exposure_geoms_3(): ) field = ogr.FieldDefn( - "Object ID", + "object_id", ogr.OFTInteger, ) layer.CreateField(field) field = ogr.FieldDefn( - "ObjectName", + "object_name", ogr.OFTString, ) field.SetWidth(50) @@ -172,8 +172,8 @@ def create_exposure_geoms_3(): geom = ogr.CreateGeometryFromWkt(geoms[0]) ft = ogr.Feature(layer.GetLayerDefn()) - ft.SetField("Object ID", 5) - ft.SetField("ObjectName", f"fp_{5}") + ft.SetField("object_id", 5) + ft.SetField("object_name", f"fp_{5}") ft.SetGeometry(geom) layer.CreateFeature(ft) @@ -182,8 +182,8 @@ def create_exposure_geoms_3(): geom = ogr.CreateGeometryFromWkt(geoms[1]) ft = ogr.Feature(layer.GetLayerDefn()) - ft.SetField("Object ID", 6) - ft.SetField("ObjectName", f"fp_{6}") + ft.SetField("object_id", 6) + ft.SetField("object_name", f"fp_{6}") ft.SetGeometry(geom) layer.CreateFeature(ft) @@ -226,7 +226,7 @@ def create_exposure_grid(): for x, y in product(oneD, oneD): data[x, y] = 2000 + ((x + y) * 100) band.WriteArray(data) - band.SetMetadataItem("damage_function", "struct_1") + band.SetMetadataItem("fn_damage", "struct_1") band.FlushCache() src.FlushCache() @@ -325,7 +325,6 @@ def create_settings_geom(): doc = { "global": { "crs": "EPSG:4326", - "keep_temp_files": True, }, "output": { "path": "output/geom_event", @@ -403,7 +402,6 @@ def create_settings_grid(): doc = { "global": { "crs": "EPSG:4326", - "keep_temp_files": True, }, "output": { "path": "output/grid_event", diff --git a/docs/PDF_Documentation.qmd b/docs/PDF_Documentation.qmd index 52d5a955..6221541d 100644 --- a/docs/PDF_Documentation.qmd +++ b/docs/PDF_Documentation.qmd @@ -218,7 +218,7 @@ damages per asset [**exposure.geom**] `csv`: This will create an exposure CSV file within the exposure folder ([folder strucute](index.qmd)) that contains the [required information](exposure.qmd) per asset. `crs`: The projection of the exposure vector file. -`file1`: This will create an exposure vector file within the exposure folder with the assets'geometry and Object ID. +`file1`: This will create an exposure vector file within the exposure folder with the assets'geometry and object_id. [**vulnerability**] `file`: This will create an vulnerability curves CSV file within the vulnerability folder ([folder strucute](index.qmd)) that contains the damage curves. @@ -254,9 +254,9 @@ If the user prefers to create a more advanced model or to utilize grid data inst | var_ as_band | True/False Read netCDF subdatasets as raster bands | No | False | | **[exposure.geom]** | | | | | csv | File path to exposure.csv file, that contains information about e.g. asset type, max. potential damage, aggregation and so forth) | Yes | 'exposure/exposure.csv' | -| file1 | File path of exposure vector file with Object ID column to enable linking exposure.csv output to vector file. | Yes | 'exposure/buildings.gpkg' | +| file1 | File path of exposure vector file with object_id column to enable linking exposure.csv output to vector file. | Yes | 'exposure/buildings.gpkg' | | crs | Projection of exposure data | Yes, if crs is unknown in dataset | 'EPSG:32617' | -| index | Define the name of the index column of the data to link the exposure CSV file with the exposure vector file. | No | 'Object ID' | +| index | Define the name of the index column of the data to link the exposure CSV file with the exposure vector file. | No | 'object_id' | | **[exposure.grid]** | | | | | file | File path to exposure.nc grid file, that contains the spatial information and information about the maximum potential damage per cell. | Yes, if netCDF is provided as input | 'exposure/raster.nc' | | crs | Output projection | Yes, if crs is unknown in dataset | 'EPSG:32617' | @@ -324,7 +324,7 @@ For users who would want to create their own exposure data, or modify existing e | Field | Description | Required | Example | |----------------------------------------|--------------------------------------------------------------------------------------------|------------------------------------------|----------------------------| -| Object ID | Unique numerical indentifier of the object | Yes | 1 | +| object_id | Unique numerical indentifier of the object | Yes | 1 | | Object Name | Unique name of the object | No | fp_1 | | Primary Object Type | Object type | No | RES1_1SNB | | Secondary Object Type | More specification about object type | No | Res 1,1 Story no basement | @@ -342,8 +342,8 @@ For users who would want to create their own exposure data, or modify existing e A more detailed description of the data fields in the *exposure.csv* can be found below; - **Object ID/Object name** - Object ID and Object name are administrative information, which the user is free to choose. Input must be unique for each object, if they are not unique, FIAT gives a warning and stops the model built-up. + **object_id/Object name** + object_id and Object name are administrative information, which the user is free to choose. Input must be unique for each object, if they are not unique, FIAT gives a warning and stops the model built-up. **Primary/Secondary object type** The primary object type describes the category of the asset (e.g. residential or commercial). The secondary object type allows for a more detailed profile of the object (e.g. single-story home, or grocery store). The developer of the exposure dataset is free to set their own categories of object types. (*Exception: FIAT requires **roads** to be assigned as **primary object type = ‘road**’, to summarize road damages separately from buildings and utilities*.) diff --git a/docs/_quarto.yml b/docs/_quarto.yml index 7c593dc4..7125a677 100644 --- a/docs/_quarto.yml +++ b/docs/_quarto.yml @@ -36,13 +36,20 @@ website: page-footer: left: | Made possible by: - - Deltares - + + Deltares + Deltares + navbar: logo: _static/fiat.svg search: true @@ -86,9 +93,9 @@ website: - "setup_guide/general/docker.qmd" - section: "FIAT package" contents: - - setup_guide/kernel/installation.qmd - - setup_guide/kernel/application.qmd - - setup_guide/kernel/linux.qmd + - setup_guide/kernel/install.qmd + - setup_guide/kernel/dev.qmd + - setup_guide/kernel/build.qmd - title: "User guide" collapse-level: 1 contents: @@ -108,8 +115,9 @@ website: - text: Exposure data file: user_guide/data/exposure.qmd contents: - - user_guide/data/exposure/csv.qmd + - user_guide/data/exposure/data.qmd - user_guide/data/exposure/geometries.qmd + - user_guide/data/exposure/csv.qmd - user_guide/data/vulnerability.qmd - user_guide/data/supported.qmd - title: Examples @@ -176,6 +184,14 @@ quartodoc: children: separate - name: GridModel children: separate + - subtitle: Methods + desc: The hazard functions + package: fiat.methods + contents: + - ead.risk_density + - ead.calc_ead + - flood.calculate_hazard + - flood.calculate_damage # Logging - title: Logging diff --git a/docs/_static/theme-dark.scss b/docs/_static/theme-dark.scss index 34089287..b67b2127 100644 --- a/docs/_static/theme-dark.scss +++ b/docs/_static/theme-dark.scss @@ -16,6 +16,8 @@ $navbar-hl: $links ; // Code blocks $code-block-bg-alpha: -.8; +$code-bg: $highlight; +$code-color: $text; .navbar-nav .dropdown-menu { background-color: $highlight; // Dark mode dropdown background color @@ -68,3 +70,18 @@ $code-block-bg-alpha: -.8; // border-color: rgb(255, 72, 0); // // border-width: 1px; // } + +// :root { +// --footer-image: url('/_static/images/deltares-white.svg'); +// } +// +// footer img { +// content: var(--footer-image); +// } + +.footer-image-light { + display: none; +} +.footer-image-dark { + display: inline; +} diff --git a/docs/_static/theme-light.scss b/docs/_static/theme-light.scss index c0995e67..6d1692e4 100644 --- a/docs/_static/theme-light.scss +++ b/docs/_static/theme-light.scss @@ -23,3 +23,10 @@ $navbar-hl: $links ; border-color: $links; border: 1px solid $links; } + +.footer-image-dark { + display: none; +} +.footer-image-light { + display: inline; +} diff --git a/docs/changelog.qmd b/docs/changelog.qmd index 30607ba8..29427ecc 100644 --- a/docs/changelog.qmd +++ b/docs/changelog.qmd @@ -6,38 +6,59 @@ title: "What's new?" This contains the unreleased changes to Delft-FIAT. ### Added +- Attribute `size` of `GeomSource` object for in situ size return; `count` becomes private +- Attribute `size` of `GridSource` object for in situ size return; `count` becomes private +- Attributes `dtypes` and `geom_type` of `GeomModel` - Build time to FIAT cli (when build with pyinstaller), viewed with `fiat --version` +- Different types of exposure (e.g. 'damage', 'affected', 'outage' etc.) - Docker file for docker image creation +- Extra arguments to `grid.reproject` - Function (`generate_jobs`) to generate jobs for parallelization - Function (`execute_pool`) to execute code in parallel using `multiprocessing` +- Flood hazard/ damage functions (`methods` submodule) - General method of creating output files in `GeomModel` (`_setup_output_files`) - Method `_create_model_dirs` of `ConfigReader` object for creating result directories in one go - Method `add_handler` of the `Log` object; user setting of custom stream +- Method `add_feature_with_map` of `GeomModel` to set features with extra info +- Method `create` of `GeomModel` to create an ogr.DataSource (new dataset) +- Method `create_equal_grids` of `GridModel` for making hazard and exposure grid spatially equal - Method `set` of `ConfigReader` object -- Method `size` of `GeomSource` object for in situ size return; `count` becomes private -- Method `size` of `GridSource` object for in situ size return; `count` becomes private +- Not stopping when exposure and hazard grids (`GridModel`), but instead make them spatially equal - Numpy >= 2.0.0 support - Python 3.12 support -- Settings toml file: global setting 'loglevel' +- Settings toml file: global setting 'global.loglevel'; default 'INFO' +- Settings toml file: exposure setting 'exposure.types'; default 'flood' - Setting return period as a variable in hazard map bands (risk) - Support for using pixi for binary creation (properly) ### Changed - Better version of `BufferHandler` +- Exposure data headers are now lower-/ snakecase ('object_id' -> 'objectId'), see [docs](./user_guide/data/exposure/data.qmd) - Fixed binary creation in general, but also specifically for `GDAL >= v3.9.1` - Made read methods of `BaseModel`, `GeomModel` and `GridModel` public (removed underscore) +- Made csv files (exposure data) optional +- Moved hazard/ damage calculation function to `methods` submodule +- Proper checking for duplicate columns in csv files +- Settings toml file: exposure setting 'exposure.csv.file' (becomes optional) - Testing of workers (not properly caught due to using `multiprocessing`) - Testing only based on integers ### Deprecated - Base object `_BaseHandler`; incompatible with Python 3.12 +- Function `open_exp` from `fiat.io`, superseded by general use of `open_csv` - Method `add_c_handler` in favour of the more generalized `add_handler` - Methods `_create_output_dir`, `_create_tmp_dir` and `_create_risk_dir` of the `ConfigReader` object +- Object `ExposureTable`, now done via `TableLazy` +- Resolve stage of `GeomModel`; now properly handled in `GeomModel.run` - Setting return period via the name of the hazard band (risk) +- Settings toml file: global setting 'global.keep_temp_files' +- Settings toml file: output setting 'output.geom.chunk'; superseded by 'global.geom.chunk' - Support of `Python` versions under `3.9.0` +- Temporary files (`GeomModel`) - `TextHandler` object; unused ### Documentation +- Added methods to the api - Cleaner home page - Getting started remade into `Information` diff --git a/docs/examples/single_event.ipynb b/docs/examples/single_event.ipynb index e54e6038..b36cd028 100644 --- a/docs/examples/single_event.ipynb +++ b/docs/examples/single_event.ipynb @@ -79,7 +79,7 @@ "from pathlib import Path\n", "\n", "# check the output\n", - "out = open_csv(Path(\"../../.testdata/output/geom_event\", \"output.csv\"), index=\"Object ID\")\n", + "out = open_csv(Path(\"../../.testdata/output/geom_event\", \"output.csv\"), index=\"object_id\")\n", "print(out.columns)" ] }, @@ -98,8 +98,8 @@ "metadata": {}, "outputs": [], "source": [ - "assert float(out[2, \"Total Damage\"]) == 740\n", - "assert float(out[3, \"Total Damage\"]) == 1038" + "assert float(out[2, \"total_damage\"]) == 740\n", + "assert float(out[3, \"total_damage\"]) == 1038" ] } ], @@ -119,7 +119,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.11.7" + "version": "3.12.4" } }, "nbformat": 4, diff --git a/docs/setup_guide/general/conda.qmd b/docs/setup_guide/general/conda.qmd index 84d8dcfe..eed9a04f 100644 --- a/docs/setup_guide/general/conda.qmd +++ b/docs/setup_guide/general/conda.qmd @@ -6,7 +6,8 @@ In order to develop on **FIAT** locally, the Python package manager **Miniforge3 Download and install [Miniforge3](https://github.com/conda-forge/miniforge#mambaforge) -Initialize conda by running the following in the Miniforge prompt: +Make sure the conda binary (and mamba) is added to PATH. In windows this is simply done via the 'set environment variables' screen, on linux one can append the 'PATH' variable via the `.bashrc` configurations file (or another rc file corresponding with the shell in use). +Initialize conda by running the following command in your shell. ```bash conda init diff --git a/docs/setup_guide/index.qmd b/docs/setup_guide/index.qmd index 0a5cfd4b..401a8ecd 100644 --- a/docs/setup_guide/index.qmd +++ b/docs/setup_guide/index.qmd @@ -8,18 +8,18 @@ title: "Setup guide" @@ -28,6 +28,4 @@ title: "Setup guide" ::: ## Overview -FIAT requires Python 3.9 or greater, a [package manager](general/conda.qmd), and Git to be installed. After that, FIAT can be [installed](kernel/installation.qmd) for use or for development, or it can be [freezed](kernel/application.qmd) as an application/executable. - -It's possible to [run FIAT on Linux](kernel/linux.qmd) but this is currently in the bèta phase. +FIAT requires Python 3.10 or greater, a [package manager](general/conda.qmd), and Git to be installed. After that, FIAT can be [installed](kernel/install.qmd) for use or for [development](kernel/dev.qmd), or it can be [built](kernel/build.qmd) as an application/executable. diff --git a/docs/setup_guide/kernel/application.qmd b/docs/setup_guide/kernel/application.qmd deleted file mode 100644 index 30ba0343..00000000 --- a/docs/setup_guide/kernel/application.qmd +++ /dev/null @@ -1,28 +0,0 @@ ---- -title: "Freeze FIAT as an application" ---- - -To freeze FIAT as an application/executable, it is required to have the FIAT repository cloned on your local machine. A different environment is needed to build the application: - -- Create a yml for a seperate **build** environment: - -```bash -python make_env.py build -``` - -- Create the environment with mamba. This time, FIAT will be automatically installed with the environment: - -```bash -mamba env create -f environment.yml -``` - -- Go to the .build/core directory and execute the pybuild.bat script: - -```bash -cd ./.build - -win64.bat -``` - -That's it. -A FIAT application will be located in the `{root}/bin/Release` folder. diff --git a/docs/setup_guide/kernel/build.qmd b/docs/setup_guide/kernel/build.qmd new file mode 100644 index 00000000..94226373 --- /dev/null +++ b/docs/setup_guide/kernel/build.qmd @@ -0,0 +1,67 @@ +--- +title: "Building FIAT" +--- + +To build FIAT as an executable/ binary, it is required to clone the FIAT repository to your local drive. +This required [git](https://git-scm.com/) to be installed on your local machine. + +```bash +git clone git@github.com:Deltares/Delft-FIAT.git +``` + +Besides git, a python installation/ environment with the necessary packages is needed. +It is recommended to use [miniforge3](../general/conda.qmd) for this purpose. + +::: {.panel-tabset} + +## Windows + +- Create a yml for a seperate **build** environment: + +```bash +python make_env.py build +``` + +- Create the environment with mamba. This time, FIAT will be automatically installed with the environment: + +```bash +mamba env create -f environment.yml +``` + +- Go to the .build/ directory and execute the win64.bat script: + +```bash +cd ./.build + +win64.bat +``` + +That's it. +A FIAT application will be located in the `{root}/bin/Release` folder. + +## Linux + +- Create a yml for a seperate **build** environment: + +```bash +python make_env.py build +``` + +- Create the environment with mamba. This time, FIAT will be automatically installed with the environment: + +```bash +mamba env create -f environment.yml +``` + +- Set the rights of the shell script and exexute it: + +```bash +chmod u+x .build/linux64.sh + +.build/linux64.sh +``` + +That's it. +A FIAT application will be located in the `{root}/bin/Release` folder + +::: diff --git a/docs/setup_guide/kernel/dev.qmd b/docs/setup_guide/kernel/dev.qmd new file mode 100644 index 00000000..c78c2c49 --- /dev/null +++ b/docs/setup_guide/kernel/dev.qmd @@ -0,0 +1,37 @@ +--- +title: "Development install" +--- + +This is for those who wish to contribute to the development of FIAT. + +- First, clone the FIAT repository on Github into a local directory of choice: + +```bash +cd ~/{your path} + +git clone https://github.com/Deltares/Delft-FIAT.git fiat +``` + +- Create a new development environment. Make sure you either have tomli or tomllib (build-in with Python 3.11) in your base enviroment. Go into your cloned FIAT repository folder and create the environment file by running the *make_env.py* script: + +```bash +cd ~/{your path}/fiat + +python make_env.py dev +``` + +- Then, create and activate the new environment in conda: + +```bash +conda env create -f environment.yml + +conda activate fiat_dev +``` + +- To install all the required dependencies, run: + +```bash +pip install -e . +``` + +There you go. FIAT is now installed on your local machine for development purposes. diff --git a/docs/setup_guide/kernel/installation.qmd b/docs/setup_guide/kernel/install.qmd similarity index 50% rename from docs/setup_guide/kernel/installation.qmd rename to docs/setup_guide/kernel/install.qmd index a888ff0a..6f73a175 100644 --- a/docs/setup_guide/kernel/installation.qmd +++ b/docs/setup_guide/kernel/install.qmd @@ -11,7 +11,7 @@ To create a new environment follow the steps below. - Create a new environment: ```bash -conda create -n fiat python=3.11.* +conda create -n fiat python=3.12.* ``` - Activate the environment: @@ -19,50 +19,31 @@ conda create -n fiat python=3.11.* conda activate fiat` ``` -- Install FIAT from Github. After creating the new environment, you need to install all dependencies from the Deltares Github repository. You can use **pip install** to do so: - -```bash -pip install git+https://github.com/Deltares/Delft-FIAT.git -``` - -### Existing environment -If you want to install FIAT into an existing environment, simply activate the desired environment and run: - -```bash -pip install git+https://github.com/Deltares/Delft-FIAT.git -``` - -## For Development -This is for those who wish to contribute to the development of FIAT. - -- First, clone the FIAT repository on Github into a local directory of choice: +- Install FIAT from [pypi](): ```bash -cd ~/{your path} - -git clone https://github.com/Deltares/Delft-FIAT.git fiat +pip install delft_fiat ``` -- Create a new development environment. Make sure you either have tomli or tomllib (build-in with Python 3.11) in your base enviroment. Go into your cloned FIAT repository folder and create the environment file by running the *make_env.py* script: +- Install FIAT from [conda-forge]() ```bash -cd ~/{your path}/fiat +# When conda-forge is added as a channel +conda install delft_fiat -python make_env.py dev +# When not +conda install delft_fiat -c conda-forge ``` -- Then, create and activate the new environment in conda: +- Install FIAT from Github. After creating the new environment, you need to install all dependencies from the Deltares Github repository. You can use **pip install** to do so: ```bash -conda env create -f environment.yml - -conda activate fiat_dev +pip install git+https://github.com/Deltares/Delft-FIAT.git ``` -- To install all the required dependencies, run: +### Existing environment +If you want to install FIAT into an existing environment, simply activate the desired environment and run: ```bash -pip install -e . +pip install git+https://github.com/Deltares/Delft-FIAT.git ``` - -There you go. FIAT is now installed on your local machine for development purposes. diff --git a/docs/setup_guide/kernel/linux.qmd b/docs/setup_guide/kernel/linux.qmd deleted file mode 100644 index d67bf023..00000000 --- a/docs/setup_guide/kernel/linux.qmd +++ /dev/null @@ -1,5 +0,0 @@ ---- -title: "FIAT on Linux" ---- - -It's possible to run FIAT on Linux but priorities are set on releasing and finalizing other parts of the documentation so this will be filled at a later stage. diff --git a/docs/user_guide/data/exposure.qmd b/docs/user_guide/data/exposure.qmd index 1f0d5e51..2bcc9a6b 100644 --- a/docs/user_guide/data/exposure.qmd +++ b/docs/user_guide/data/exposure.qmd @@ -15,7 +15,7 @@ Because FIAT consists of two models, the exposure data can be supplied in either FIAT can be used to assess flood damages to buildings, roads, utilities, and other types of structures of which data is available. These type of assets are often spatially represented with vector, or geometry data. Building a FIAT model with geometry exposure data requires two data types: - [Geometry file(s)](exposure/geometries.qmd) -- [CSV file](exposure/csv.qmd) +- [CSV file (optional)](exposure/csv.qmd) ### Raster Data diff --git a/docs/user_guide/data/exposure/csv.qmd b/docs/user_guide/data/exposure/csv.qmd index c70f5a6b..5965ce57 100644 --- a/docs/user_guide/data/exposure/csv.qmd +++ b/docs/user_guide/data/exposure/csv.qmd @@ -2,61 +2,9 @@ title: CSV --- -The exposure data CSV file contains information about each asset in the area of interest that is needed for the damage calculation. Each row represents one asset, such as a building, road segment, or utility, and each column represents an attribute of the asset, such as its location, elevation or maximum potential damage value. -For users that want to create their own exposure data, or modify existing exposure data, a description of the default fields (columns) in the exposure data CSV can be found in @tbl-exposure. +The exposure data CSV file contains information about each asset in the area of interest that is optional for the damage calculation. Each row represents one asset, such as a building, road segment, or utility, and each column represents an attribute of the asset, such as its location, elevation or maximum potential damage value. +For users that want to create their own exposure data, or modify existing exposure data, a description of the default fields (columns) in the exposure data CSV can be found [here](./data.qmd). ::: {.callout-tip title="Tip"} Users can add as many columns to the exposure CSV as they want, for example, to later identify the most damaged administrative areas. All data field will also be saved in the output of FIAT. ::: - - -| Field | Required | Multiple | Example | -|--------------------------------------------------------|----------|----------|---------------------------| -| [Object ID](#object-id) | Yes | - | 1 | -| [Object Name](#object-name) | No | - | fp_1 | -| [Primary Object Type](#primarysecondary-object-type) | No | - | RES1_1SNB | -| [Secondary Object Type](#primarysecondary-object-type) | No | - | Res 1,1 Story no basement | -| [Extraction Method](#extraction-method) | Yes | - | 'centroid' | -| [Damage Function: < ANY >](#damage-functions) | Yes | Yes | struct_2 | -| [Ground Floor Height](#ground-floor-height) | Yes | - | 4 | -| [Ground Elevation](#ground-elevation) | No | - | 10.11 | -| [Max Potential Damage: < ANY >](#max-potential-damage) | Yes | Yes | 193457.00 | -| [< Additional attributes >](#additional-attributes) | No | | 1205 | -: exposure.csv input {#tbl-exposure .hover} - -A more detailed description of the data fields in the *exposure.csv* can be found below. - -#### Object ID -The object ID is used to link the geometries to the information in the exposure CSV. The object ID must be unique for each object. At the moment of writing, FIAT does not check whether the object IDs are indeed unique and it is the responsibility of the user to ensure this. At a later stage, FIAT will have a built-in check for this. - -#### Object name -The object name field can be chosen freely and can serve as a field for identifying the damages assets. - -#### Primary/Secondary object type -The primary object type describes the occupancy or category of the asset (e.g. residential or commercial). The secondary object type allows for a more detailed profile of the object (e.g. single-story home or grocery store). The developer of the exposure dataset is free to set their own categories of asset types. - -::: {.callout-tip title="Tip"} -Defining primary/secondary object types can facilitate the assignment of damage functions to the assets for the user by creating automatic look-up tables. -::: - -#### Extraction Method -The extraction method refers to how the water level or water depth is sampled per asset. The options are (1) *centroid*, which samples the water level or depth at the estimated centroid inside of the asset, or (2) area, which considers the water level or depth over the entire polygon or line and takes either an average or maximum. The user can set the choice for the latter per damage function, in the [vulnerability curves file](../vulnerability.qmd). - -::: {.callout-important} - In case the user selects 'area' as extraction method for certain assets, the geometries related to those assets should be a line or polygon. -::: - -#### Damage Functions -The user can input as many damage functions and related max potential damages as required. However, at least one damage function and related max potential damage must be provides per asset. Damage function column are named as **"Damage Function: `damage function name`"**. This name can be chosen freely but it is recommended to give it a descriptive name, such as *structure* or *content*, and it is required to use the same name in the max potential damage column. The value that must be entered is the name of the damage function, relating to the names in the corresponding damage functions in the [vulnerability curves file](../vulnerability.qmd). Globally, continental damage functions can be obtained from [European Commission's Joint Research Centre](https://publications.jrc.ec.europa.eu/repository/handle/JRC105688), but it is recommended to use more location-specific functions when available. - -#### Ground Floor Height -The ground floor height column indicates the height of the ground floor of an asset above the ground elevation. For example, when a building is built on poles, the ground floor (finished floor in the US) is lifted off the ground with a certain height. Usually, buildings are at least a certain amount above ground elevation because of the road and sidewalk. This value is used to calculate the water depth inside of the asset, for buildings. - -#### Ground Elevation -The ground elevation is the value of a digital elevation model (DEM) at the location of the asset. - -#### Max Potential Damage -The maximum potential damage corresponds to the damage functions for each asset. For each damage function type that was assigned, a maximum potential damage must also be assigned. These values represent the maximum damage to, for example, the structure, content, or other (e.g. inventory). There are methods to derive these values, based on building type and area of the building. Globally, maximum damage values per country can be obtained from [European Commission's Joint Research Centre](https://publications.jrc.ec.europa.eu/repository/handle/JRC105688), but it is recommended to use more location-specific damage values. In the US, [FEMA Hazus](https://www.fema.gov/flood-maps/products-tools/hazusis) is an industry standard in how to derive these values. - -#### Additional attributes -Attitional attributes can be added to the exposure CSV file. The added columns can be named freely. The resulting FIAT output contains all data that was put in, therefore they can be used by the [FIAT Toolbox](https://github.com/Deltares/fiat_toolbox) to calculate user-specified output metrics. For example, a user could configure the output metrics to always show the number of people affected with a high social vulnerability, per administrative region. The user can add columns to the CSV file manually or with the help of the [HydroMT-FIAT model builder](https://deltares.github.io/hydromt_fiat/latest/#), in which the user is free to select a descriptive *'Label name'* for the additional attributes. If the user wishes to connect the damages per administrative or other aggregation area to a GIS file, the **FIAT Toolbox offers a post-processing tool to do so. diff --git a/docs/user_guide/data/exposure/data.qmd b/docs/user_guide/data/exposure/data.qmd new file mode 100644 index 00000000..d595ca68 --- /dev/null +++ b/docs/user_guide/data/exposure/data.qmd @@ -0,0 +1,55 @@ +--- +title: Exposure data +--- + +In the table down below, an overview is provided od the necessary inputs in either the [geometry file](./geometries.qmd) or the [csv file](./csv.qmd)". + +| Field | Required | Multiple | Example | +|--------------------------------------------------------|----------|----------|---------------------------| +| [object_id](#object_id) | Yes | - | 1 | +| [object_name](#object-name) | No | - | fp_1 | +| [extract_method](#extraction-method) | Yes | - | centroid | +| [fn_< any >_< any >](#damage-functions) | Yes | Yes | struct_2 | +| [ground_flht](#ground-floor-height) | No | - | 4 | +| [ground_elevtn](#ground-elevation) | No | - | 10.11 | +| [max_< any >_< any >](#max-potential-damage) | Yes | Yes | 193457.00 | +| [< Additional attributes >](#additional-attributes) | No | | 1205 | +: exposure.csv input {#tbl-exposure .hover} + +A more detailed description of the data fields can be found below. + +#### Object ID +The object ID is used to link the geometries to the information in the exposure CSV. The object ID must be unique for each object. At the moment of writing, FIAT does not check whether the object IDs are indeed unique and it is the responsibility of the user to ensure this. At a later stage, FIAT will have a built-in check for this. + +#### Object name +The object name field can be chosen freely and can serve as a field for identifying the damages assets. + +#### Extraction Method +The extraction method refers to how the water level or water depth is sampled per asset. The options are (1) *centroid*, which samples the water level or depth at the estimated centroid inside of the asset, or (2) area, which considers the water level or depth over the entire polygon or line and takes either an average or maximum. The user can set the choice for the latter per damage function, in the [vulnerability curves file](../vulnerability.qmd). + +::: {.callout-important} + In case the user selects 'area' as extraction method for certain assets, the geometries related to those assets should be a line or polygon. +::: + +#### Damage Functions +The user can input as many damage functions and related max potential damages as required. However, at least one damage function and related max potential damage must be provides per asset. Damage function column are named as `fn_< exposure_type >_< component >` (for exposure types, see [here](../../settings/optional.qmd#exposure)). This name can be chosen freely but it is recommended to give it a descriptive name, such as **damage** (exposure type) in combination with **structure** or **content** (components), and it is required to use the same name in the max potential damage column. (i.e. `max_< exposure-type >_< component >`) The value that must be entered is the name of the damage function, relating to the names in the corresponding damage functions in the [vulnerability curves file](../vulnerability.qmd). Globally, continental damage functions can be obtained from [European Commission's Joint Research Centre](https://publications.jrc.ec.europa.eu/repository/handle/JRC105688), but it is recommended to use more location-specific functions when available. + +#### Ground Floor Height +The ground floor height column indicates the height of the ground floor of an asset above the ground elevation. For example, when a building is built on poles, the ground floor (finished floor in the US) is lifted off the ground with a certain height. Usually, buildings are at least a certain amount above ground elevation because of the road and sidewalk. This value is used to calculate the water depth inside of the asset, for buildings. + +::: {.callout-important} + Required for flood damage calculation, see [global model types](../../settings/optional.qmd#global) +::: + +#### Ground Elevation +The ground elevation is the value of a digital elevation model (DEM) at the location of the asset. + +::: {.callout-important} + Required for flood damage calculation, see [global model types](../../settings/optional.qmd#global) +::: + +#### Max Potential Damage +The maximum potential damage corresponds to the damage functions for each asset. For each damage function type that was assigned, a maximum potential damage must also be assigned. These values represent the maximum damage to, for example, the structure, content, or other (e.g. inventory). There are methods to derive these values, based on building type and area of the building. Globally, maximum damage values per country can be obtained from [European Commission's Joint Research Centre](https://publications.jrc.ec.europa.eu/repository/handle/JRC105688), but it is recommended to use more location-specific damage values. In the US, [FEMA Hazus](https://www.fema.gov/flood-maps/products-tools/hazusis) is an industry standard in how to derive these values. + +#### Additional attributes +Attitional attributes can be added to the exposure CSV file. The added columns can be named freely. The resulting FIAT output contains all data that was put in, therefore they can be used by the [FIAT Toolbox](https://github.com/Deltares/fiat_toolbox) to calculate user-specified output metrics. For example, a user could configure the output metrics to always show the number of people affected with a high social vulnerability, per administrative region. The user can add columns to the CSV file manually or with the help of the [HydroMT-FIAT model builder](https://deltares.github.io/hydromt_fiat/latest/#), in which the user is free to select a descriptive *'Label name'* for the additional attributes. If the user wishes to connect the damages per administrative or other aggregation area to a GIS file, the **FIAT Toolbox offers a post-processing tool to do so. diff --git a/docs/user_guide/data/exposure/geometries.qmd b/docs/user_guide/data/exposure/geometries.qmd index 4668fb63..6d2d0aa4 100644 --- a/docs/user_guide/data/exposure/geometries.qmd +++ b/docs/user_guide/data/exposure/geometries.qmd @@ -1,10 +1,13 @@ --- title: Geometry --- -A user can specify one or multiple **geometry files** in the [settings.toml](../../settings/index.qmd). It is advisable to give descriptive names to the geometry files, e.g., *buildings.gpkg* or *roads.gpkg*. Most-used file types are GeoPackages or Shapefiles but more file types are accepted as can be seen at the bottom of this page. All geometry files must have an attribute `Object ID` that holds unique IDs over all geometry files, linking to the `Object ID` column in the exposure CSV file (see below). +A user can specify one or multiple **geometry files** in the [settings.toml](../../settings/index.qmd). It is advisable to give descriptive names to the geometry files, e.g., *buildings.gpkg* or *roads.gpkg*. Most-used file types are GeoPackages or Shapefiles but more file types are accepted as can be seen at the bottom of this page. All geometry files must have an attribute `object_id`. These need to be unique over all geometry files if a csv file with exposure data is provided and corresponding to the `object_is`'s in the provided csv file. + +If no csv file provided with the exposure information, that information needs to be present in the geometry files(s) themselves. +An overview of all the required input is shown [here](./data.qmd). ::: {.callout-tip title="Example"} -When using a GeoPackage file with roads containing `Object IDs` 1-100, a second GeoPackage file with building footprints should contain `Object ID`s from 101-n. Then, those IDs should link to the corresponding rows and IDs in the exposure CSV file. +When using a GeoPackage file with roads containing `object_ids` 1-100, a second GeoPackage file with building footprints should contain `object_id`s from 101-n. Then, if applicable, those IDs should link to the corresponding rows and IDs in the exposure CSV file. ::: When providing a polygon or line geometry file, a user can choose to use the *area extraction method* ('extraction method' = 'area', see [here](csv.qmd#tbl-exposure)) for considering water levels or depths over the whole area of the asset. This can be more appropriate than point estimates for example for large buildings. diff --git a/docs/user_guide/data/supported.qmd b/docs/user_guide/data/supported.qmd index 4089e522..7c2cbad7 100644 --- a/docs/user_guide/data/supported.qmd +++ b/docs/user_guide/data/supported.qmd @@ -21,19 +21,54 @@ The amount of columns are determined by length of the header or the first line o ## Geometry files Below a list is presented of supported drivers in regards to the vector files: +### Reading :::: {.columns} ::: {.center width="70%" style="height: 500px; overflow: auto;"} ```{python} #| echo: false -#| label: tbl-geom-drivers -#| tbl-cap: "Available drivers for vector data" +#| label: tbl-geom-read-drivers +#| tbl-cap: "Available drivers for reading vector data" import pandas as pd -from fiat.util import GEOM_DRIVER_MAP +from fiat.util import GEOM_READ_DRIVER_MAP from IPython.display import HTML from osgeo import gdal -dr_map = dict(sorted(GEOM_DRIVER_MAP.items())) +dr_map = dict(sorted(GEOM_READ_DRIVER_MAP.items())) + +long_name = [] +for _dr in dr_map.values(): + dr = gdal.GetDriverByName(_dr) + long_name.append(dr.LongName) + dr = None + +df = pd.DataFrame( + data={ + "File extension": dr_map.keys(), + "Driver": dr_map.values(), + "Long name": long_name, + } +) + +HTML(df.to_html(index=False)) +``` +::: +:::: + +### Writing +:::: {.columns} +::: {.center width="70%" style="height: 500px; overflow: auto;"} +```{python} +#| echo: false +#| label: tbl-geom-write-drivers +#| tbl-cap: "Available drivers for writing vector data" + +import pandas as pd +from fiat.util import GEOM_WRITE_DRIVER_MAP +from IPython.display import HTML +from osgeo import gdal + +dr_map = dict(sorted(GEOM_WRITE_DRIVER_MAP.items())) long_name = [] for _dr in dr_map.values(): @@ -55,7 +90,7 @@ HTML(df.to_html(index=False)) :::: ## Gridded data files -Below a list is presented of supported drivers in regards to the raster files: +Below a list is presented of supported drivers (both reading and writing) in regards to the raster files: :::: {.columns} ::: {.center width="70%" style="height: 500px; overflow: auto;"} diff --git a/docs/user_guide/settings/computation.qmd b/docs/user_guide/settings/computation.qmd index f1690851..ccaa1b62 100644 --- a/docs/user_guide/settings/computation.qmd +++ b/docs/user_guide/settings/computation.qmd @@ -7,7 +7,7 @@ The inputs give the user a bit more control over the computational aspects of FI | Entry | Type | Default | |:---------------------------------|---------|-------------| | **[global]** | | | -| [threads](#global) | integer | | +| [threads](#global) | integer | 1 | | **[global.geom]** | | | | [chunk](#global.geom) | integer | - | | **[global.grid]** | | | @@ -31,7 +31,7 @@ This input benefits from multiple threads. ::: ::: {.callout-note} -This input is only applicable to the [GeomModel]() +This input is only applicable to the [GeomModel](../../info/models.qmd#geommodel) ::: #### [global.grid] @@ -39,5 +39,5 @@ This input is only applicable to the [GeomModel]() - `chunk`: Set the chunk size for the gridded calculations. This will chunk the data in rectangles with the goal of reducing the memory foodprint. An example would be `[1024, 1024,]`. ::: {.callout-note} -This input is only applicable to the [GridModel]() +This input is only applicable to the [GridModel](../../info/models.qmd#gridmodel) ::: diff --git a/docs/user_guide/settings/index.qmd b/docs/user_guide/settings/index.qmd index 38108041..0f2c61a5 100644 --- a/docs/user_guide/settings/index.qmd +++ b/docs/user_guide/settings/index.qmd @@ -33,10 +33,12 @@ File paths in the settings can be relative to the settings.toml file or absolute | [file](#hazard) | string | Yes | | | [elevation_reference](#hazard) | string | Yes | | | [risk](#hazard) | boolean | No | false | +| **[exposure.csv]** | | | | +| [file](#exposure.csv) | string | No | - | | **[exposure.geom]** | | | | | [file[n]](#exposure.geom) | string | Yes | | -| **[exposure.csv]** | | | | -| [file](#exposure.csv) | string | Yes | | +| **[exposure.grid]** | | | | +| [file](#exposure.grid) | string | Yes | | | **[vulnerability]** | | | | | [file](#vulnerability) | string | Yes | | : Most basic settings file input {#tbl-toml .hover} @@ -69,17 +71,28 @@ If provided, the suffix is mandatory. The suffix should match the suffix of the - `risk`: In case of a risk calculation this must be set to "true", for a single events this must be set to "false". +#### [exposure.csv] + +- `file`: The path to the exposure CSV file (recommended to be within the [exposure folder](../data/index.qmd)) that contains the [required information](../data/exposure.qmd) per asset. There can only be one exposure CSV file. + #### [exposure.geom] -- `file[n]`: The path to the exposure vector file (recommended to be within the [exposure folder](../data/index.qmd)) with the assets' geometry and Object ID. The user can provide multiple vector files. Therefore the '[n]' suffix, as the user can create mulitple entries for vector files (e.g. `file1`, `file2` etc.). +- `file[n]`: The path to the exposure vector file (recommended to be within the [exposure folder](../data/index.qmd)) with the assets' geometry and object_id. The user can provide multiple vector files. Therefore the '[n]' suffix, as the user can create mulitple entries for vector files (e.g. `file1`, `file2` etc.). ::: {.callout-warning} The suffix is mandatory. So if only one file is provided, name it `file1`. ::: +::: {.callout-note} +Only required when running the geometry based model. +::: -#### [exposure.csv] +#### [exposure.grid] -- `file`: The path to the exposure CSV file (recommended to be within the [exposure folder](../data/index.qmd)) that contains the [required information](../data/exposure.qmd) per asset. There can only be one exposure CSV file. +- `file`: The path to the exposure raster file (recommended to be within the [exposure folder](../data/index.qmd)). + +::: {.callout-note} +Only required when running the raster based model. +::: #### [vulnerability] diff --git a/docs/user_guide/settings/optional.qmd b/docs/user_guide/settings/optional.qmd index a07c5376..00fa34c6 100644 --- a/docs/user_guide/settings/optional.qmd +++ b/docs/user_guide/settings/optional.qmd @@ -9,18 +9,23 @@ Here we provide an overview of all optional/ non essential settings that the use | **[global]** | | | | [crs](#global) | string | - | | [loglevel](#global) | string | INFO | -| [keep_temp_files](#global) | boolean | false | -| [gdal_cache](#global) | int | - | +| [type](#global) | string | flood | +| **[global.grid]** | | | +| [prefer](#global.grid) | string | exposure | | **[hazard]** | | | | [crs](#hazard) | string | - | | [return_periods](#hazard) | list | - | | **[hazard.settings]** | | | | [subset](#hazard.settings) | string | - | | [var_as_band](#hazard.settings) | boolean | false | -| **[exposure.csv]** | | | -| [index](#exposure.csv) | string | 'Object ID' | +| **[exposure]** | | | +| [types](#exposure) | list | ['damage'] | +| **[exposure.csv.settings]** | | | +| [index](#exposure.csv.settings) | string | object_id | | **[exposure.geom]** | | | | [crs](#exposure.geom) | string | - | +| **[exposure.geom.settings]** | | | +| [index](#exposure.geom.settings) | string | object_id | | **[exposure.grid]** | | | | [crs](#exposure.grid) | string | - | | **[exposure.grid.settings]** | | | @@ -37,13 +42,11 @@ Here we provide an overview of all optional/ non essential settings that the use - `loglevel`: Set the loglevel of the fiat logger from the settings file. Choose from 'DEBUG', 'INFO', 'WARNING', 'ERROR' or 'DEAD'. -- `keep_temp_files`: Permanently save the temporary file outputs. These are produced during the calculations. +- `type`: Type of hazard. Should be the same name as the file containing the functions (e.g. 'flood' -> `fiat.methods.flood`) -- `gdal_cache`: Set the cache used by GDAL. +#### [global.grid] -::: {.callout-warning} -The entry `gdal_cache` performs poorly and will probably be deprecated soon. -::: +- `prefer`: Whether to spatially prefer exposure data or hazard data. The other will be warped when they are not equal. Chose 'exposure' or 'hazard'. #### [hazard] @@ -61,14 +64,22 @@ The entry `gdal_cache` performs poorly and will probably be deprecated soon. The entry `var_as_band` is only applicable to netCDF files. ::: -#### [exposure.csv] +#### [exposure] + +- `types`: Types of exposure. This could be monetary damages ('damage'), affected people ('affected') etc. + +#### [exposure.csv.settings] -- `index`: Set the index column of the csv file. In case of the exposure csv, if no entry is provided then FIAT will default to 'Object ID'. +- `index`: Set the index column of the csv file. In case of the exposure csv, if no entry is provided then FIAT will default to 'object_id'. #### [exposure.geom] - `crs`: Projection of the exposure geometry data if it cannot be inferred from the dataset itself. +#### [exposure.geom.settings] + +- `index`: Set the index column of the geom file(s). In case nothing is provided, the default value 'object_id' is used. + #### [exposure.grid] - `crs`: Projection of the gridded exposure data if it cannot be inferred from the dataset itself. diff --git a/pyproject.toml b/pyproject.toml index eaf82c6a..44b229b7 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -17,7 +17,9 @@ dynamic = ['version'] authors = [ {name = "B.W. Dalmijn", email = "brencodeert@outlook.com"}, - {name = "F.C. de Groen", email = "frederique.degroen@deltares.nl"} + {name = "F.C. de Groen", email = "frederique.degroen@deltares.nl"}, + {name = "S.A. Rautenbach", email = "sarah.rautenbach@deltares.nl"}, + {name = "P. Athanasiou", email = "panos.athanasiou@deltares.nl"}, ] maintainers = [ {name = "B.W. Dalmijn", email = "brencodeert@outlook.com"}, @@ -41,7 +43,7 @@ classifiers = [ # 3 - Alpha # 4 - Beta # 5 - Production/Stable -'Development Status :: 3 - Alpha', +'Development Status :: 4 - Beta', # Indicate who your project is intended for 'Intended Audience :: End Users/Desktop', @@ -171,9 +173,8 @@ exclude = ["docs"] "test/conftest.py" = ["E402"] "src/fiat/__init__.py" = ["E402", "F401", "F403"] "src/fiat/cli/__init__.py" = ["F403"] -"src/fiat/cy/__init__.py" = ["F403"] "src/fiat/gis/__init__.py" = ["F403"] -"src/fiat/gui/__init__.py" = ["F403"] +"src/fiat/methods/__init__.py" = ["F403"] "src/fiat/models/__init__.py" = ["F403"] [tool.ruff.pydocstyle] diff --git a/src/fiat/cfg.py b/src/fiat/cfg.py index a83fa37c..d5015d9c 100644 --- a/src/fiat/cfg.py +++ b/src/fiat/cfg.py @@ -130,14 +130,6 @@ def _create_model_dirs( ) self.set("output.path", _p) - # Temporary files directory - _p = self._create_dir( - self.get("output.path"), - ".tmp", - hidden=False, - ) - self.set("output.tmp.path", _p) - # Damage directory for grid risk calculations if self.get("hazard.risk") and check_config_grid(self): _p = self._create_dir( diff --git a/src/fiat/check.py b/src/fiat/check.py index a1696f47..c3703aae 100644 --- a/src/fiat/check.py +++ b/src/fiat/check.py @@ -1,6 +1,5 @@ """Checks for the data of FIAT.""" -import fnmatch from pathlib import Path from osgeo import osr @@ -41,7 +40,6 @@ def check_config_geom( ): """_summary_.""" _req_fields = [ - "exposure.csv.file", "exposure.geom.crs", "exposure.geom.file1", ] @@ -121,7 +119,8 @@ def check_grid_exact( ): msg = f"CRS of hazard data ({get_srs_repr(haz.get_srs())}) does not match the \ CRS of the exposure data ({get_srs_repr(exp.get_srs())})" - raise FIATDataError(msg) + logger.warning(msg) + return False gtf1 = [round(_n, 2) for _n in haz.get_geotransform()] gtf2 = [round(_n, 2) for _n in exp.get_geotransform()] @@ -129,12 +128,16 @@ def check_grid_exact( if gtf1 != gtf2: msg = f"Geotransform of hazard data ({gtf1}) does not match geotransform of \ exposure data ({gtf2})" - raise FIATDataError(msg) + logger.warning(msg) + return False if haz.shape != exp.shape: msg = f"Shape of hazard ({haz.shape}) does not match shape of \ exposure data ({exp.shape})" - raise FIATDataError(msg) + logger.warning(msg) + return False + + return True def check_internal_srs( @@ -196,7 +199,7 @@ def check_hazard_band_names( ): """_summary_.""" if risk: - return [f"{n}Y" for n in rp] + return [f"{n}y" for n in rp] if count == 1: return [""] @@ -243,13 +246,12 @@ def check_hazard_subsets( ## Exposure def check_exp_columns( columns: tuple | list, + specific_columns: tuple | list = [], ): """_summary_.""" _man_columns = [ - "Object ID", - "Ground Elevation", - "Ground Floor Height", - ] + "object_id", + ] + specific_columns _check = [item in columns for item in _man_columns] if not all(_check): @@ -257,29 +259,24 @@ def check_exp_columns( msg = f"Missing mandatory exposure columns: {_missing}" raise FIATDataError(msg) - dmg = fnmatch.filter(columns, "Damage Function: *") - dmg_suffix = [item.split(":")[1].strip() for item in dmg] - mpd = fnmatch.filter(columns, "Max Potential Damage: *") - mpd_suffix = [item.split(":")[1].strip() for item in mpd] - - if not dmg: - msg = "No damage function were given in exposure data." - raise FIATDataError(msg) - - if not mpd: - msg = "No maximum potential damages were given in exposure data" - raise FIATDataError(msg) - _check = [item in mpd_suffix for item in dmg_suffix] - if not any(_check): - msg = "Damage function and maximum potential damage do not have a single match" +def check_exp_derived_types( + type: str, + found: tuple | list, + missing: tuple | list, +): + """_summary_.""" + # Error when no columns are found for vulnerability type + if not found: + msg = f"For type: '{type}' no matching columns were found for \ +fn_{type}_* and max_{type}_* columns." raise FIATDataError(msg) - if not all(_check): - _missing = [item for item, b in zip(dmg_suffix, _check) if not b] + # Log when combination of fn and max is missing + if missing: logger.warning( f"No every damage function has a corresponding \ -maximum potential damage: {_missing}" + maximum potential damage: {missing}" ) @@ -288,7 +285,7 @@ def check_exp_grid_dmfs( dmfs: tuple | list, ): """_summary_.""" - _ef = [_i.get_metadata_item("damage_function") for _i in exp] + _ef = [_i.get_metadata_item("fn_damage") for _i in exp] _i = None _check = [item in dmfs for item in _ef] @@ -298,4 +295,13 @@ def check_exp_grid_dmfs( raise FIATDataError(msg) +def check_exp_index_col( + obj: object, + index_col: type, +): + """_summary_.""" + if index_col not in obj.columns: + raise FIATDataError(f"Index column ('{index_col}') not found in {obj.path}") + + ## Vulnerability diff --git a/src/fiat/cli/util.py b/src/fiat/cli/util.py index 6c7c555c..5cad211b 100644 --- a/src/fiat/cli/util.py +++ b/src/fiat/cli/util.py @@ -27,6 +27,6 @@ def run_log( func() except BaseException: t, v, tb = sys.exc_info() - logger.error(",".join(v.args)) + logger.error(",".join([str(item) for item in v.args])) # Exit with code 1 sys.exit(1) diff --git a/src/fiat/cy/__init__.py b/src/fiat/cy/__init__.py deleted file mode 100644 index 30efdc45..00000000 --- a/src/fiat/cy/__init__.py +++ /dev/null @@ -1 +0,0 @@ -"""Useless cy stuff.""" diff --git a/src/fiat/cy/cy.pyx b/src/fiat/cy/cy.pyx deleted file mode 100644 index e8752e24..00000000 --- a/src/fiat/cy/cy.pyx +++ /dev/null @@ -1,36 +0,0 @@ -import re -from libc.stdio cimport printf -from libc.string cimport strcpy, strlen -from cpython.mem cimport PyMem_Malloc, PyMem_Realloc, PyMem_Free - -cdef extern from "Python.h": - char* PyUnicode_AsUTF8(object unicode) - -cpdef char * return_char(char* st): - print(st) - return st.decode("utf-8") - -cdef char ** create_string_array(size): - cdef char **ret = PyMem_Malloc(size * sizeof(char *)) - return ret - -cpdef tuple csv(int size, int skip, r, f): - oid = [None]*size - cdef int *linesize = PyMem_Malloc( - size * sizeof(int)) - cdef int c = 0 - cdef char* l - try: - while True: - line = f.readline().strip() - if not line: - break - z = r.split(line)[1] - oid[c] = z.decode("utf-8") - linesize[c] = len(line) - c += 1 - for x in linesize[:5]: - print(x) - return ([int(x) for x in linesize[:size]],oid) - finally: - PyMem_Free(linesize) diff --git a/src/fiat/error.py b/src/fiat/error.py index ed2c1d62..322eb4b6 100644 --- a/src/fiat/error.py +++ b/src/fiat/error.py @@ -6,8 +6,8 @@ class DriverNotFoundError(Exception): def __init__(self, gog, path): self.base = f"{gog} data" - self.msg = f"Extension of file: {path.name} is not recoqnized" - super(DriverNotFoundError, self).__init__(self.base) + self.msg = f"Extension of file: {path.name} not recoqnized" + super(DriverNotFoundError, self).__init__(f"{self.base} -> {self.msg}") def __str__(self): return f"{self.base} -> {self.msg}" diff --git a/src/fiat/gis/grid.py b/src/fiat/gis/grid.py index be9aa6d3..70a4d628 100644 --- a/src/fiat/gis/grid.py +++ b/src/fiat/gis/grid.py @@ -30,7 +30,10 @@ def clip( def reproject( gs: GridSource, - crs: str, + dst_crs: str, + dst_gtf: list | tuple = None, + dst_width: int = None, + dst_height: int = None, out_dir: Path | str = None, resample: int = 0, ) -> object: @@ -40,8 +43,16 @@ def reproject( ---------- gs : GridSource Input object. - crs : str + dst_crs : str Coodinates reference system (projection). An accepted format is: `EPSG:3857`. + dst_gtf : list | tuple, optional + The geotransform of the warped dataset. Must be defined in the same + coordinate reference system as dst_crs. When defined, its only used when + both 'dst_width' and 'dst_height' are defined. + dst_width : int, optional + The width of the warped dataset in pixels. + dst_height : int, optional + The height of the warped dataset in pixels. out_dir : Path | str, optional Output directory. If not defined, if will be inferred from the input object. resample : int, optional @@ -63,14 +74,33 @@ def reproject( fname = Path(out_dir, f"{gs.path.stem}_repr_fiat{gs.path.suffix}") out_srs = osr.SpatialReference() - out_srs.SetFromUserInput(crs) + out_srs.SetFromUserInput(dst_crs) out_srs.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER) + warp_kw = {} + if all([item is not None for item in [dst_gtf, dst_width, dst_height]]): + warp_kw.update( + { + "xRes": dst_gtf[1], + "yRes": dst_gtf[5], + "outputBounds": ( + dst_gtf[0], + dst_gtf[3] + dst_gtf[5] * dst_height, + dst_gtf[0] + dst_gtf[1] * dst_width, + dst_gtf[3], + ), + "width": dst_width, + "height": dst_height, + } + ) + dst_src = gdal.Warp( str(fname_int), gs.src, + srcSRS=gs.get_srs(), dstSRS=out_srs, resampleAlg=resample, + **warp_kw, ) out_srs = None @@ -85,6 +115,6 @@ def reproject( dst_src = None gc.collect() - os.remove(fname_int) + os.unlink(fname_int) return open_grid(fname, **_gs_kwargs) diff --git a/src/fiat/gui/__init__.py b/src/fiat/gui/__init__.py deleted file mode 100644 index 0cc76e71..00000000 --- a/src/fiat/gui/__init__.py +++ /dev/null @@ -1,5 +0,0 @@ -"""GUI module of FIAT.""" - -__all__ = ["gui"] - -from .main import gui diff --git a/src/fiat/gui/main.py b/src/fiat/gui/main.py deleted file mode 100644 index 2d9b1146..00000000 --- a/src/fiat/gui/main.py +++ /dev/null @@ -1,19 +0,0 @@ -"""GUI main window.""" - -import sys - -from PySide2.QtWidgets import QApplication, QMainWindow - - -class gui(QMainWindow): - """_summary_.""" - - pass - - -if __name__ == "__main__": - app = QApplication(sys.argv) - w = gui() - w.show() - sys.exit(app.exec_()) - pass diff --git a/src/fiat/gui/widgets.py b/src/fiat/gui/widgets.py deleted file mode 100644 index 5cab9b73..00000000 --- a/src/fiat/gui/widgets.py +++ /dev/null @@ -1 +0,0 @@ -"""Widgets of the GUI.""" diff --git a/src/fiat/io.py b/src/fiat/io.py index aed849cf..03175f62 100644 --- a/src/fiat/io.py +++ b/src/fiat/io.py @@ -19,9 +19,11 @@ from fiat.util import ( DD_NEED_IMPLEMENTED, DD_NOT_IMPLEMENTED, - GEOM_DRIVER_MAP, + GEOM_READ_DRIVER_MAP, + GEOM_WRITE_DRIVER_MAP, GRID_DRIVER_MAP, NEED_IMPLEMENTED, + NEWLINE_CHAR, NOT_IMPLEMENTED, DummyLock, _dtypes_from_string, @@ -29,6 +31,7 @@ _read_gridsrouce_layers, _text_chunk_gen, deter_type, + find_duplicates, regex_pattern, replace_empty, ) @@ -220,8 +223,8 @@ def __init__( self, file: str | Path, srs: osr.SpatialReference, - layer_meta: ogr.FeatureDefn, - buffer_size: int = 20000, # geometries + layer_defn: ogr.FeatureDefn = None, + buffer_size: int = 100000, # geometries lock: Lock = None, ): """_summary_. @@ -253,18 +256,23 @@ def __init__( # Set for later use self.srs = srs - self.in_layer_meta = layer_meta self.flds = {} self.n = 1 + if layer_defn is None: + with open_geom(self.file, mode="r") as _r: + layer_defn = _r.layer.GetLayerDefn() + _r = None + self.layer_defn = layer_defn + # Create the buffer self.buffer = open_geom(f"/vsimem/{file.stem}.gpkg", "w") self.buffer.create_layer( srs, - layer_meta.GetGeomType(), + layer_defn.GetGeomType(), ) self.buffer.set_layer_from_defn( - layer_meta, + layer_defn, ) # Set some check vars # TODO: do this based om memory foodprint @@ -274,6 +282,7 @@ def __init__( def __del__(self): self.buffer = None + self.layer_defn = None def __reduce__(self) -> str | tuple[Any, ...]: pass @@ -289,10 +298,10 @@ def _reset_buffer(self): # Re-create self.buffer.create_layer( self.srs, - self.in_layer_meta.GetGeomType(), + self.layer_defn.GetGeomType(), ) self.buffer.set_layer_from_defn( - self.in_layer_meta, + self.layer_defn, ) self.create_fields(self.flds) @@ -306,33 +315,21 @@ def close(self): self._clear_cache() self.buffer.close() - def add_feature( + def add_feature_with_map( self, ft: ogr.Feature, fmap: dict, ): """_summary_.""" - _local_ft = ogr.Feature(self.buffer.layer.GetLayerDefn()) - _local_ft.SetFID(ft.GetFID()) - _local_ft.SetGeometry(ft.GetGeometryRef()) - for num in range(ft.GetFieldCount()): - _local_ft.SetField( - num, - ft.GetField(num), - ) - - for key, item in fmap.items(): - _local_ft.SetField( - key, - item, - ) + self.buffer.add_feature_with_map( + ft, + fmap=fmap, + ) if self.size + 1 > self.max_size: self.to_drive() - self.buffer.layer.CreateFeature(_local_ft) self.size += 1 - _local_ft = None def create_fields( self, @@ -349,6 +346,7 @@ def create_fields( def to_drive(self): """_summary_.""" # Block while writing to the drive + # self.buffer.close() self.lock.acquire() merge_geom_layers( self.file.as_posix(), @@ -357,6 +355,7 @@ def to_drive(self): ) self.lock.release() + # self.buffer = self.buffer.reopen(mode="w") self._reset_buffer() @@ -422,6 +421,18 @@ def write(self, b): self.to_drive() BytesIO.write(self, b) + def write_iterable(self, *args): + """_summary_. + + _extended_summary_ + """ + by = b"" + for arg in args: + by += ("," + "{}," * len(arg)).format(*arg).rstrip(",").encode() + by = by.lstrip(b",") + by += NEWLINE_CHAR.encode() + self.write(by) + ## Parsing class CSVParser: @@ -441,6 +452,7 @@ def __init__( self.meta["index_col"] = -1 self.meta["index_name"] = None self.meta["delimiter"] = delimiter + self.meta["dup_cols"] = None self.index = None self.columns = None self._nrow = self.data.size @@ -541,6 +553,7 @@ def _parse_meta( break self.columns = [item.strip() for item in line.split(self.delim)] + self.meta["dup_cols"] = find_duplicates(self.columns) self._resolve_column_headers() self._ncol = len(self.columns) break @@ -552,6 +565,17 @@ def _parse_meta( def _resolve_column_headers(self): """_summary_.""" _cols = self.columns + dup = self.meta["dup_cols"] + if dup is None: + dup = [] + # Solve duplicate values first + count = dict(zip(dup, [0] * len(dup))) + for idx, item in enumerate(_cols): + if item in dup: + _cols[idx] += f"_{count[item]}" + count[item] += 1 + + # Solve unnamed column headers _cols = [_col if _col else f"Unnamed_{_i+1}" for _i, _col in enumerate(_cols)] self.columns = _cols @@ -575,15 +599,6 @@ def read( **self.meta, ) - def read_exp(self): - """_summary_.""" - return ExposureTable( - data=self.data, - index=self.index, - columns=self.columns, - **self.meta, - ) - ## Structs class Grid( @@ -805,12 +820,6 @@ class GeomSource(_BaseIO, _BaseStruct): ``` """ - _type_map = { - "int": ogr.OFTInteger64, - "float": ogr.OFTReal, - "str": ogr.OFTString, - } - def __new__( cls, file: str, @@ -831,20 +840,27 @@ def __init__( _BaseStruct.__init__(self) _BaseIO.__init__(self, file, mode) - if self.path.suffix not in GEOM_DRIVER_MAP: - raise DriverNotFoundError("") + if mode == "r": + _map = GEOM_READ_DRIVER_MAP + else: + _map = GEOM_WRITE_DRIVER_MAP + + if self.path.suffix not in _map: + raise DriverNotFoundError(gog="Geometry", path=self.path) - driver = GEOM_DRIVER_MAP[self.path.suffix] + driver = _map[self.path.suffix] self._driver = ogr.GetDriverByName(driver) + info = gdal.VSIStatL(self.path.as_posix()) - if self.path.exists() and not overwrite: + if (self.path.exists() or info is not None) and not overwrite: self.src = self._driver.Open(self.path.as_posix(), self._mode) else: if not self._mode: raise OSError(f"Cannot create {self.path} in 'read' mode.") - self.src = self._driver.CreateDataSource(self.path.as_posix()) + self.create(self.path) + info = None self._count = 0 self._cur_index = 0 @@ -933,12 +949,15 @@ def reduced_iter( yield ft _c += 1 - def reopen(self): + def reopen( + self, + mode: str = "r", + ): """Reopen a closed GeomSource.""" if not self._closed: return self - obj = GeomSource.__new__(GeomSource, self.path) - obj.__init__(self.path) + obj = GeomSource.__new__(GeomSource, self.path, mode=mode) + obj.__init__(self.path, mode=mode) return obj @property @@ -968,6 +987,16 @@ def columns(self): """ return tuple(self._columns.keys()) + @property + @_BaseIO._check_state + def dtypes(self): + """_summary_.""" + if self.layer is not None: + _flds = self.layer.GetLayerDefn() + dt = [_flds.GetFieldDefn(_i).type for _i in range(_flds.GetFieldCount())] + _flds = None + return dt + @property @_BaseIO._check_state def fields(self): @@ -980,12 +1009,18 @@ def fields(self): _flds = None return fh + @property + @_BaseIO._check_state + def geom_type(self): + """Return the geometry type.""" + return self.layer.GetGeomType() + @property @_BaseIO._check_state def size( self, ): - """_summary_.""" + """Return the size (geometry count).""" if self.layer is not None: count = self.layer.GetFeatureCount() self._count = count @@ -1001,6 +1036,9 @@ def add_feature( Only in write (`'w'`) mode. + Note! Everything needs to already be compliant with the created/ edited + dataset. + Parameters ---------- ft : ogr.Feature @@ -1008,6 +1046,32 @@ def add_feature( """ self.layer.CreateFeature(ft) + @_BaseIO._check_mode + @_BaseIO._check_state + def add_feature_with_map( + self, + in_ft: ogr.Feature, + fmap: zip, + ): + """Add a feature with extra field data. + + Parameters + ---------- + in_ft : ogr.Feature + The feature to be added. + fmap : zip + Extra fields data, i.e. a zip object of fields id's + and the correspondingv alues + """ + ft = ogr.Feature(self.layer.GetLayerDefn()) + ft.SetFrom(in_ft) + + for key, item in fmap: + ft.SetField(key, item) + + self.layer.CreateFeature(ft) + ft = None + @_BaseIO._check_mode @_BaseIO._check_state def add_feature_from_defn( @@ -1037,12 +1101,28 @@ def add_feature_from_defn( self.layer.CreateFeature(out_ft) + @_BaseIO._check_mode + @_BaseIO._check_state + def create( + self, + path: Path | str, + ): + """Create a data source. + + Parameters + ---------- + path : Path | str + Path to the data source. + """ + self.src = None + self.src = self._driver.CreateDataSource(path.as_posix()) + @_BaseIO._check_mode @_BaseIO._check_state def create_field( self, name: str, - type: object, + type: int, ): """Add a new field. @@ -1052,13 +1132,13 @@ def create_field( ---------- name : str Name of the new field. - type : object + type : int Type of the new field. """ self.layer.CreateField( ogr.FieldDefn( name, - GeomSource._type_map[type], + type, ) ) self._retrieve_columns() @@ -1080,10 +1160,7 @@ def create_fields( are the data types of the new field. """ self.layer.CreateFields( - [ - ogr.FieldDefn(key, GeomSource._type_map[item]) - for key, item in fmap.items() - ] + [ogr.FieldDefn(key, item) for key, item in fmap.items()] ) self._retrieve_columns() @@ -1435,10 +1512,7 @@ def shape_xy(self): @property @_BaseIO._check_state def size(self): - """_summary_. - - _extended_summary_ - """ + """Return the number of bands.""" count = self.src.RasterCount self._count = count return self._count @@ -1507,7 +1581,7 @@ def deter_band_names( for n in range(self.size): name = self.get_band_name(n + 1) if not name: - _names.append(f"Band{n+1}") + _names.append(f"band{n+1}") continue _names.append(name) @@ -1623,7 +1697,6 @@ def __init__( ) -> object: """_summary_.""" # Declarations - self._dup_cols = None self.dtypes = () self.meta = kwargs self.index_col = -1 @@ -1647,11 +1720,7 @@ def __init__( columns = [f"col_{num}" for num in range(kwargs["ncol"])] # Some checking in regards to duplicates in column headers - if len(set(columns)) != len(columns): - _set = list(set(columns)) - _counts = [columns.count(elem) for elem in _set] - _dup = [elem for _i, elem in enumerate(_set) if _counts[_i] > 1] - self._dup_cols = _dup + self.columns_raw = columns # Create the column indexing self._columns = dict(zip(columns, range(kwargs["ncol"]))) @@ -1906,7 +1975,7 @@ class TableLazy(_Table): ---------- data : BufferHandler A stream. - index : strortuple, optional + index : str | tuple, optional The index column used as row indices. columns : list, optional The column headers of the table. @@ -2020,124 +2089,6 @@ def set_index( self.data = dict(zip(new_index, self.data.values())) -class ExposureTable(TableLazy): - """Create a table just for exposure data. - - Parameters - ---------- - data : BufferHandler - A stream. - index : strortuple, optional - The index column used as row indices. - columns : list, optional - The column headers of the table. - - Returns - ------- - object - An object containing a connection via a stream to a file. - """ - - def __init__( - self, - data: BufferHandler, - index: str | tuple = None, - columns: list = None, - **kwargs, - ): - TableLazy.__init__( - self, - data, - index, - columns, - **kwargs, - ) - - _dc_cols = ["Damage Function:", "Max Potential Damage"] - - for req in _dc_cols: - req_s = req.strip(":").lower().replace(" ", "_") - self.__setattr__( - req_s, - dict( - [ - (item.split(":")[1].strip(), self._columns[item]) - for item in self.columns - if item.startswith(req) - ] - ), - ) - - self._blueprint_columns = ( - ["Inundation Depth", "Reduction Factor"] - + [f"Damage: {item}" for item in self.damage_function.keys()] - + ["Total Damage"] - ) - - self._dat_len = len(self._blueprint_columns) - - def create_specific_columns(self, name: str): - """Generate new columns for output data. - - Parameters - ---------- - name : str - Exposure identifier. - - Returns - ------- - list - A list containing new columns. - """ - _out = [] - if name: - name = f"({name})" - for bp in self._blueprint_columns: - # _parts = bp.split(":") - - # if len(_parts) == 1: - # bp += f" {name}" - # _out.append(bp) - # continue - - bp += f" {name}" - _out.append(bp.strip()) - - return _out - - def create_all_columns( - self, - names: list, - extra: list = None, - ): - """Append existing columns of exposure database. - - Parameters - ---------- - names : list - Exposure identifiers. - extra : list, optional - Extra specified columns. - - Returns - ------- - list - List containing all columns. - """ - cols = [] - for n in names: - cols += self.create_specific_columns(n) - - if extra is not None: - cols += extra - - return cols - - def gen_dat_dtypes(self): - """Generate dtypes for the new columns.""" - return ",".join(["float"] * self._dat_len).encode() - - ## I/O mutating methods def merge_geom_layers( out_fn: Path | str, @@ -2239,42 +2190,6 @@ def open_csv( ) -def open_exp( - file: Path | str, - delimiter: str = ",", - header: bool = True, - index: str = None, -): - """_summary_. - - _extended_summary_ - - Parameters - ---------- - file : str - _description_ - header : bool, optional - _description_, by default True - index : str, optional - _description_, by default None - - Returns - ------- - _type_ - _description_ - """ - _handler = BufferHandler(file) - - parser = CSVParser( - _handler, - delimiter, - header, - index, - ) - - return parser.read_exp() - - def open_geom( file: Path | str, mode: str = "r", diff --git a/src/fiat/methods/__init__.py b/src/fiat/methods/__init__.py new file mode 100644 index 00000000..b6a34f59 --- /dev/null +++ b/src/fiat/methods/__init__.py @@ -0,0 +1,5 @@ +"""Logic submule of FIAT.""" + +__all__ = ["ead", "flood"] + +from . import ead, flood diff --git a/src/fiat/methods/ead.py b/src/fiat/methods/ead.py new file mode 100644 index 00000000..9a1b5871 --- /dev/null +++ b/src/fiat/methods/ead.py @@ -0,0 +1,88 @@ +"""EAD (Expected Annual Damages) related functionality.""" +import math + + +def calc_ead( + rp_coef: list, + dms: list, +) -> float: + """Calculate the EAD (risk). + + From a list of return periods and list of corresponding damages. + + Parameters + ---------- + rp_coef : list + List of return period coefficients. + dms : list + List of corresponding damages + (in the same order of the return periods coefficients). + + Returns + ------- + float + The Expected Annual Damage (EAD), or risk, as a log-linear integration over the + return periods. + """ + # Calculate the EAD + ead = sum([x * y for x, y in zip(rp_coef, dms)]) + return ead + + +def risk_density( + rp: list | tuple, +) -> list: + """Calculate the risk density factors from return periods values. + + Parameters + ---------- + rp : list | tuple + A list of return periods. + + Returns + ------- + list + List of risk density factors. + """ + # Step 1: Compute frequencies associated with T-values. + _rp = sorted(rp) + idxs = [_rp.index(n) for n in rp] + rp_u = sorted(rp) + rp_l = len(rp_u) + + f = [1 / n for n in rp_u] + lf = [math.log(1 / n) for n in rp_u] + + if rp_l == 1: + return f + + # Step 2: + c = [(1 / (lf[idx] - lf[idx + 1])) for idx in range(rp_l - 1)] + + # Step 3: + G = [(f[idx] * lf[idx] - f[idx]) for idx in range(rp_l)] + + # Step 4: + a = [ + ( + (1 + c[idx] * lf[idx + 1]) * (f[idx] - f[idx + 1]) + + c[idx] * (G[idx + 1] - G[idx]) + ) + for idx in range(rp_l - 1) + ] + b = [ + (c[idx] * (G[idx] - G[idx + 1] + lf[idx + 1] * (f[idx + 1] - f[idx]))) + for idx in range(rp_l - 1) + ] + + # Step 5: + alpha = [ + b[0] + if idx == 0 + else f[idx] + a[idx - 1] + if idx == rp_l - 1 + else a[idx - 1] + b[idx] + for idx in range(rp_l) + ] + + return [alpha[idx] for idx in idxs] diff --git a/src/fiat/methods/flood.py b/src/fiat/methods/flood.py new file mode 100644 index 00000000..5812cac2 --- /dev/null +++ b/src/fiat/methods/flood.py @@ -0,0 +1,135 @@ +"""Functions specifically for flood risk calculation.""" +import math + +from numpy import isnan +from osgeo import ogr + +from fiat.io import Table +from fiat.methods.util import AREA_METHODS + +MANDATORY_COLUMNS = ["ground_flht", "ground_elevtn"] +MANDATORY_ENTRIES = {"ref": "hazard.elevation_reference"} +NEW_COLUMNS = ["inun_depth"] + + +def calculate_hazard( + hazard: list, + reference: str, + ground_flht: float, + ground_elevtn: float = 0, + method: str = "mean", +) -> float: + """Calculate the hazard value for flood hazard. + + Parameters + ---------- + hazard : list + Raw hazard values. + reference : str + Reference, either 'dem' or 'datum'. + ground_flht : float + The height of the floor of an object (.e.g the door elevation). + ground_elevtn : float, optional + Ground height in reference to e.g. the ocean. + (Needed when 'reference' is 'datum') + method : str, optional + Chose 'max' or 'mean' for either the maximum value or the average, + by default 'mean'. + + Returns + ------- + float + A representative hazard value. + """ + _ge = 0 + if reference.lower() == "datum" and not math.isnan(ground_elevtn): + # The hazard data is referenced to a Datum + # (e.g., for flooding this is the water elevation). + _ge = ground_elevtn + + # Remove the negative hazard values to 0. + raw_l = len(hazard) + hazard = [n - _ge for n in hazard if (n - _ge) > 0.0001] + + if not hazard: + return math.nan, math.nan + + redf = 1 + + if method.lower() == "mean": + redf = len(hazard) / raw_l + + if len(hazard) > 1: + hazard = AREA_METHODS[method.lower()](hazard) + else: + hazard = hazard[0] + + # Subtract the Ground Floor Height from the hazard value + hazard -= ground_flht + + return hazard, redf + + +def calculate_damage( + hazard_value: float | int, + red_fact: float | int, + ft: ogr.Feature | list, + type_dict: dict, + vuln: Table, + vul_min: float | int, + vul_max: float | int, + vul_round: int, +) -> tuple: + """Calculate the damage corresponding with the hazard value. + + Parameters + ---------- + hazard_value : float | int + The representative hazard value. + red_fact : float | int + The reduction factor. How much to compensate for the lack of touching the grid + by an object (geometry). + ft : ogr.Feature | list + A feature or feature info (whichever has to contain the exposure data). + See docs on running FIAT with an without csv. + type_dict : dict + The exposure types and corresponding column id's. + vuln : Table + Vulnerability data. + vul_min : float | int + Minimum value of the index of the vulnerability data. + vul_max : float | int + Maximum value of the index of the vulnerability data. + vul_round : int + Significant decimals to be used. + + Returns + ------- + tuple + Damage values. + """ + # unpack type_dict + fn = type_dict["fn"] + maxv = type_dict["max"] + + # Define outgoing list of values + out = [0] * (len(fn) + 1) + + # Calculate the damage per catagory, and in total (_td) + total = 0 + idx = 0 + for key, col in fn.items(): + if isnan(hazard_value) or ft[col] is None or ft[col] == "nan": + val = "nan" + else: + hazard_value = max(min(vul_max, hazard_value), vul_min) + f = vuln[round(hazard_value, vul_round), ft[col]] + val = f * ft[maxv[key]] * red_fact + val = round(val, 2) + total += val + out[idx] = val + idx += 1 + + out[-1] = round(total, 2) + + return out diff --git a/src/fiat/methods/util.py b/src/fiat/methods/util.py new file mode 100644 index 00000000..8c42e096 --- /dev/null +++ b/src/fiat/methods/util.py @@ -0,0 +1,8 @@ +"""Calculation utility.""" + +from fiat.util import mean + +AREA_METHODS = { + "max": max, + "mean": mean, +} diff --git a/src/fiat/models/__init__.py b/src/fiat/models/__init__.py index 6a57c7fe..0df18f80 100644 --- a/src/fiat/models/__init__.py +++ b/src/fiat/models/__init__.py @@ -1,6 +1,7 @@ """Entry point for models.""" -__all__ = ["GeomModel", "GridModel"] +__all__ = ["GeomModel", "GridModel", "worker_geom", "worker_grid"] +from . import worker_geom, worker_grid from .geom import GeomModel from .grid import GridModel diff --git a/src/fiat/models/base.py b/src/fiat/models/base.py index a3ce62b3..a61c0f8a 100644 --- a/src/fiat/models/base.py +++ b/src/fiat/models/base.py @@ -1,5 +1,6 @@ """Base model of FIAT.""" +import importlib from abc import ABCMeta, abstractmethod from multiprocessing import Manager, get_context from os import cpu_count @@ -19,7 +20,6 @@ from fiat.gis.crs import get_srs_repr from fiat.io import open_csv, open_grid from fiat.log import spawn_logger -from fiat.models.calc import calc_rp_coef from fiat.util import NEED_IMPLEMENTED, deter_dec logger = spawn_logger("fiat.model") @@ -44,29 +44,25 @@ def __init__( self.exposure_grid = None self.hazard_grid = None self.vulnerability_data = None + # Type of calculations + self.type = self.cfg.get("global.type", "flood") + self.module = importlib.import_module(f"fiat.methods.{self.type}") + self.cfg.set("global.type", self.type) # Vulnerability data self._vul_step_size = 0.01 self._rounding = 2 self.cfg.set("vulnerability.round", self._rounding) - # Temporay files - self._keep_temp = False # Threading stuff self._mp_ctx = get_context("spawn") self._mp_manager = Manager() - self.max_threads = 1 - self.nthreads = 1 - self.chunk = None + self.threads = 1 self.chunks = [] - self.nchunk = 0 - self._set_max_threads() self._set_model_srs() + self._set_num_threads() self.read_hazard_grid() self.read_vulnerability_data() - if "global.keep_temp_files" in self.cfg: - self._keep_temp = self.cfg.get("global.keep_temp_files") - @abstractmethod def __del__(self): self.srs = None @@ -75,10 +71,6 @@ def __del__(self): def __repr__(self): return f"<{self.__class__.__name__} object at {id(self):#018x}>" - @abstractmethod - def _clean_up(self): - raise NotImplementedError(NEED_IMPLEMENTED) - def _set_model_srs(self): """_summary_.""" _srs = self.cfg.get("global.crs") @@ -115,12 +107,19 @@ def _set_model_srs(self): # Clean up gm = None - @abstractmethod - def _set_num_threads( - self, - ): + def _set_num_threads(self): """_summary_.""" - raise NotImplementedError(NEED_IMPLEMENTED) + max_threads = cpu_count() + user_threads = self.cfg.get("global.threads") + if user_threads is not None: + if user_threads > max_threads: + logger.warning( + f"Given number of threads ('{user_threads}') \ +exceeds machine thread count ('{max_threads}')" + ) + self.threads = min(max_threads, user_threads) + + logger.info(f"Using number of threads: {self.threads}") @abstractmethod def _setup_output_files( @@ -139,7 +138,7 @@ def read_hazard_grid(self): self.cfg.generate_kwargs("hazard.settings"), ) kw.update( - self.cfg.generate_kwargs("global.grid"), + self.cfg.generate_kwargs("global.grid.chunk"), ) data = open_grid(path, **kw) ## checks @@ -188,9 +187,6 @@ def read_hazard_grid(self): path, ) self.cfg.set("hazard.return_periods", rp) - # Directly calculate the coefficients - rp_coef = calc_rp_coef(rp) - self.cfg.set("hazard.rp_coefficients", rp_coef) # Information for output ns = check_hazard_band_names( @@ -219,7 +215,7 @@ def read_vulnerability_data(self): logger.info("Executing vulnerability checks...") # Column check - check_duplicate_columns(data._dup_cols) + check_duplicate_columns(data.meta["dup_cols"]) # upscale the data (can be done after the checks) if "vulnerability.step_size" in self.cfg: @@ -235,20 +231,6 @@ def read_vulnerability_data(self): # When all is done, add it self.vulnerability_data = data - def _set_max_threads(self): - """_summary_.""" - self.max_threads = cpu_count() - _max_threads = self.cfg.get("global.threads") - if _max_threads is not None: - if _max_threads > self.max_threads: - logger.warning( - f"Given number of threads ('{_max_threads}') \ -exceeds machine thread count ('{self.max_threads}')" - ) - self.max_threads = min(self.max_threads, _max_threads) - - logger.info(f"Available number of threads: {self.max_threads}") - @abstractmethod def run( self, diff --git a/src/fiat/models/calc.py b/src/fiat/models/calc.py deleted file mode 100644 index 18d67c70..00000000 --- a/src/fiat/models/calc.py +++ /dev/null @@ -1,203 +0,0 @@ -"""Logic of FIAT.""" - -import math - -from fiat.util import mean - -_inun_calc = { - "mean": mean, - "max": max, -} - - -## Calculates coefficients used to compute the EAD as a linear function of -## the known damages -# Args: -# T (list of ints): return periods T1 … Tn for which damages are known -# Returns: -# alpha [list of floats]: coefficients a1, …, an (used to compute the AED as -# a linear function of the known damages) -# In which D(f) is the damage, D, as a function of the frequency of exceedance, f. -# In order to compute this EAD, function D(f) needs to be known for -# the entire range of frequencies. Instead, D(f) is only given for the n -# frequencies as mentioned in the table above. So, in order to compute the integral -# above, some assumptions need to be made for function D(h): -# (i) For f > f1 the damage is assumed to be equal to 0 -# (ii) For f float: - """_summary_. - - Parameters - ---------- - haz : float - _description_ - idx : tuple - _description_ - values : tuple - _description_ - sig : int - significant figures - - Returns - ------- - float - Damage factor - """ - if math.isnan(haz): - return 0.0 - - # Clip based on min and max vulnerability values - haz = max(min(haz, idx[-1]), idx[0]) - - return values[idx.index(round(haz, sig))] - - -def calc_haz( - haz: list, - ref: str, - gfh: float, - ge: float = 0, - method: str = "mean", -) -> float: - """_summary_. - - Parameters - ---------- - haz : list - _description_ - ref : str - _description_ - gfh : float - _description_ - ge : float, optional - Ground Elevation, by default 0 - method : str, optional - _description_, by default "mean" - - Returns - ------- - float - _description_ - """ - _ge = 0 - if ref.lower() == "datum" and not math.isnan(ge): - # The hazard data is referenced to a Datum - # (e.g., for flooding this is the water elevation). - _ge = ge - - # Remove the negative hazard values to 0. - raw_l = len(haz) - haz = [n - _ge for n in haz if (n - _ge) > 0.0001] - - if not haz: - return math.nan, math.nan - - redf = 1 - - if method.lower() == "mean": - redf = len(haz) / raw_l - - if len(haz) > 1: - haz = _inun_calc[method.lower()](haz) - else: - haz = haz[0] - - # Subtract the Ground Floor Height from the hazard value - haz = haz - gfh - - return haz, redf - - -def calc_risk( - rp_coef: list, - dms: list, -) -> float: - """Calculate the EAD (risk). - - From a list of return periods and list of corresponding damages. - - Parameters - ---------- - rp_coef : list - List of return period coefficients. - dms : list - List of corresponding damages - (in the same order of the return periods coefficients). - - Returns - ------- - float - The Expected Annual Damage (EAD), or risk, as a log-linear integration over the - return periods. - """ - # Calculate the EAD - ead = sum([x * y for x, y in zip(rp_coef, dms)]) - return ead diff --git a/src/fiat/models/geom.py b/src/fiat/models/geom.py index b8ad9c11..28f18e1c 100644 --- a/src/fiat/models/geom.py +++ b/src/fiat/models/geom.py @@ -1,13 +1,19 @@ """Geom model of FIAT.""" +import copy import os +import re import time from pathlib import Path +from osgeo import ogr + from fiat.cfg import ConfigReader from fiat.check import ( check_duplicate_columns, check_exp_columns, + check_exp_derived_types, + check_exp_index_col, check_geom_extent, check_internal_srs, check_vs_srs, @@ -15,22 +21,20 @@ from fiat.gis import geom, overlay from fiat.gis.crs import get_srs_repr from fiat.io import ( - open_exp, + open_csv, open_geom, ) from fiat.log import setup_mp_log, spawn_logger +from fiat.models import worker_geom from fiat.models.base import BaseModel from fiat.models.util import ( - GEOM_MIN_CHUNK, - GEOM_MIN_WRITE_CHUNK, + EXPOSURE_FIELDS, + GEOM_DEFAULT_CHUNK, csv_def_file, - csv_temp_file, execute_pool, generate_jobs, - geom_threads, ) -from fiat.models.worker import geom_resolve, geom_worker -from fiat.util import create_1d_chunk +from fiat.util import create_1d_chunk, discover_exp_columns, generate_output_columns logger = spawn_logger("fiat.model.geom") @@ -60,22 +64,67 @@ def __init__( ): super().__init__(cfg) + # Set/ declare some variables + self.exposure_types = self.cfg.get("exposure.types", ["damage"]) + # Setup the geometry model - self.read_exposure_data() - self.read_exposure_geoms() + self.read_exposure() + self.get_exposure_meta() self._set_chunking() - self._set_num_threads() self._queue = self._mp_manager.Queue(maxsize=10000) def __del__(self): BaseModel.__del__(self) - def _clean_up(self): - """_summary_.""" - _p = self.cfg.get("output.tmp.path") - for _f in _p.glob("*"): - os.unlink(_f) - os.rmdir(_p) + def _discover_exposure_meta( + self, + columns: dict, + meta: dict, + index: int, + ): + """Simple method for sorting out the exposure meta.""" # noqa: D401 + # check if set from the csv file + if -1 not in meta: + meta[index] = {} + # Check the exposure column headers + check_exp_columns( + list(columns.keys()), + specific_columns=getattr(self.module, "MANDATORY_COLUMNS"), + ) + + # Check the found columns + types = {} + for t in self.exposure_types: + types[t] = {} + found, found_idx, missing = discover_exp_columns(columns, type=t) + check_exp_derived_types(t, found, missing) + types[t] = found_idx + meta[index].update({"types": types}) + + ## Information for output + extra = [] + if self.cfg.get("hazard.risk"): + extra = ["ead"] + new_fields, len1, total_idx = generate_output_columns( + getattr(self.module, "NEW_COLUMNS"), + types, + extra=extra, + suffix=self.cfg.get("hazard.band_names"), + ) + meta[index].update( + { + "new_fields": new_fields, + "slen": len1, + "total_idx": total_idx, + } + ) + else: + meta[index] = copy.deepcopy(meta[-1]) + new_fields = meta[index]["new_fields"] + + # Set the indices for the outgoing columns + idxs = list(range(len(columns), len(columns) + len(new_fields))) + meta[index].update({"idxs": idxs}) def _set_chunking(self): """_summary_.""" @@ -83,120 +132,104 @@ def _set_chunking(self): max_geom_size = max( [item.size for item in self.exposure_geoms.values()], ) - # Set calculations chunk size - self.chunk = max_geom_size - _chunk = self.cfg.get("global.geom.chunk") - if _chunk is not None: - self.chunk = max(GEOM_MIN_CHUNK, _chunk) - - # Set cache size for outgoing data - _out_chunk = self.cfg.get("output.geom.settings.chunk") - if _out_chunk is None: - _out_chunk = GEOM_MIN_WRITE_CHUNK - self.cfg.set("output.geom.settings.chunk", _out_chunk) - - # Determine amount of threads - self.nchunk = max_geom_size // self.chunk - if self.nchunk == 0: - self.nchunk = 1 - # Constrain by max threads - if self.max_threads < self.nchunk: - logger.warning( - f"Less threads ({self.max_threads}) available than \ -calculated chunks ({self.nchunk})" - ) - self.nchunk = self.max_threads - # Set the 1D chunks self.chunks = create_1d_chunk( max_geom_size, - self.nchunk, - ) - - def _set_num_threads(self): - """_summary_.""" - self.nthreads = geom_threads( - self.max_threads, - self.hazard_grid.size, - self.nchunk, + self.threads, ) + # Set the write size chunking + chunk_int = self.cfg.get("global.geom.chunk", GEOM_DEFAULT_CHUNK) + self.cfg.set("global.geom.chunk", chunk_int) def _setup_output_files(self): """_summary_.""" - # Create the temp file plus header - _nms = self.cfg.get("hazard.band_names") - for idx, _ in enumerate(_nms): - csv_temp_file( - self.cfg.get("output.tmp.path"), - idx + 1, - self.exposure_data.meta["index_name"], - self.exposure_data.create_specific_columns(_nms[idx]), - ) - - # Define the outgoing file - out_csv = "output.csv" - if "output.csv.name" in self.cfg: - out_csv = self.cfg.get("output.csv.name") - self.cfg.set("output.csv.name", out_csv) - - # Create an empty csv file for the separate thread to till - csv_def_file( - Path(self.cfg.get("output.path"), out_csv), - self.exposure_data.columns + tuple(self.cfg.get("output.new_columns")), - ) - - # Do the same for the geometry files - for key in self.exposure_geoms.keys(): - _add = key[-1] + # Setup the geometry output files + for key, gm in self.exposure_geoms.items(): # Define outgoing dataset - out_geom = f"spatial{_add}.gpkg" - if f"output.geom.name{_add}" in self.cfg: - out_geom = self.cfg.get(f"output.geom.name{_add}") - self.cfg.set(f"output.geom.name{_add}", out_geom) + out_geom = f"spatial{key}.fgb" + if f"output.geom.name{key}" in self.cfg: + out_geom = self.cfg.get(f"output.geom.name{key}") + self.cfg.set(f"output.geom.name{key}", out_geom) + # Open and write a layer with the necessary fields with open_geom( Path(self.cfg.get("output.path"), out_geom), mode="w", overwrite=True ) as _w: - pass + _w.create_layer(self.srs, gm.geom_type) + _w.create_fields(dict(zip(gm.fields, gm.dtypes))) + new = self.cfg.get("_exposure_meta")[key]["new_fields"] + _w.create_fields(dict(zip(new, [ogr.OFTReal] * len(new)))) + _w = None + + # Check whether to do the same for the csv + if self.exposure_data is not None: + out_csv = self.cfg.get("output.csv.name", "output.csv") + self.cfg.set("output.csv.name", out_csv) + + # Create an empty csv file for the separate thread to till + csv_def_file( + Path(self.cfg.get("output.path"), out_csv), + self.exposure_data.columns + + tuple(self.cfg.get("_exposure_meta")[-1]["new_fields"]), + ) + + def get_exposure_meta(self): + """Get the exposure meta regarding the data itself (fields etc.).""" + # Get the relevant column headers + meta = {} + if self.exposure_data is not None: + self._discover_exposure_meta( + self.exposure_data._columns, + meta, + -1, + ) + for key, gm in self.exposure_geoms.items(): + columns = gm._columns + self._discover_exposure_meta(columns, meta, key) + self.cfg.set("_exposure_meta", meta) + + def read_exposure(self): + """Read all the exposure files.""" + self.read_exposure_geoms() + csv = self.cfg.get("exposure.csv.file") + if csv is not None: + self.read_exposure_data() def read_exposure_data(self): - """_summary_.""" + """Read the exposure data file (csv).""" path = self.cfg.get("exposure.csv.file") logger.info(f"Reading exposure data ('{path.name}')") # Setting the keyword arguments from settings file - kw = {"index": "Object ID"} + kw = {"index": "object_id"} kw.update( self.cfg.generate_kwargs("exposure.csv.settings"), ) - data = open_exp(path, **kw) + data = open_csv(path, large=True, **kw) ##checks logger.info("Executing exposure data checks...") - # Check for mandatory columns - check_exp_columns(data.columns) - # Check for duplicate columns - check_duplicate_columns(data._dup_cols) - - ## Information for output - _ex = None - if self.cfg.get("hazard.risk"): - _ex = ["Risk (EAD)"] - cols = data.create_all_columns( - self.cfg.get("hazard.band_names"), - _ex, - ) - self.cfg.set("output.new_columns", cols) + check_duplicate_columns(data.meta["dup_cols"]) ## When all is done, add it self.exposure_data = data def read_exposure_geoms(self): - """_summary_.""" + """Read the exposure geometries.""" + # Discover the files _d = {} + # TODO find maybe a better solution of defining this in the settings file _found = [item for item in list(self.cfg) if "exposure.geom.file" in item] + _found = [item for item in _found if re.match(r"^(.*)file(\d+)", item)] + + # First check for the index_col + index_col = self.cfg.get("exposure.geom.settings.index", "object_id") + self.cfg.set("exposure.geom.settings.index", index_col) + + # For all that is found, try to read the data for file in _found: path = self.cfg.get(file) + suffix = int(re.findall(r"\d+", file.rsplit(".", 1)[1])[0]) logger.info( f"Reading exposure geometry '{file.split('.')[-1]}' ('{path.name}')" ) @@ -204,6 +237,9 @@ def read_exposure_geoms(self): ## checks logger.info("Executing exposure geometry checks...") + # check for the index column + check_exp_index_col(data, index_col=index_col) + # check the internal srs of the file _int_srs = check_internal_srs( data.get_srs(), @@ -227,45 +263,10 @@ def read_exposure_geoms(self): ) # Add to the dict - _d[file.rsplit(".", 1)[1]] = data + _d[suffix] = data # When all is done, add it self.exposure_geoms = _d - def resolve( - self, - ): - """Create permanent output. - - This is done but reading, loading and sorting the temporary output within - the `.tmp` folder within the output folder. \n - - - This method might become private. - """ - # Setup the locks for the worker - csv_lock = self._mp_manager.Lock() - geom_lock = self._mp_manager.Lock() - - # Setting up the jobs for resolving - jobs = generate_jobs( - { - "cfg": self.cfg, - "exp": self.exposure_data, - "exp_geom": self.exposure_geoms, - "chunk": self.chunks, - "csv_lock": csv_lock, - "geom_lock": geom_lock, - }, - ) - - # Execute the jobs - logger.info("Busy...") - execute_pool( - ctx=self._mp_ctx, - func=geom_resolve, - jobs=jobs, - threads=self.nthreads, - ) - def run( self, ): @@ -288,34 +289,37 @@ def run( # Start the receiver (which is in a seperate thread) _receiver.start() + # Exposure fields get function + field_func = EXPOSURE_FIELDS[self.exposure_data is None] + # Setup the jobs # First setup the locks - locks = [self._mp_manager.Lock() for _ in range(self.hazard_grid.size)] + lock1 = self._mp_manager.Lock() + lock2 = self._mp_manager.Lock() jobs = generate_jobs( { "cfg": self.cfg, "queue": self._queue, "haz": self.hazard_grid, - "idx": range(1, self.hazard_grid.size + 1), "vul": self.vulnerability_data, - "exp": self.exposure_data, + "exp_func": field_func, + "exp_data": self.exposure_data, "exp_geom": self.exposure_geoms, "chunk": self.chunks, - "lock": locks, + "lock1": lock1, + "lock2": lock2, }, - tied=["idx", "lock"], + # tied=["idx", "lock"], ) - logger.info(f"Using number of threads: {self.nthreads}") - # Execute the jobs in a multiprocessing pool _s = time.time() logger.info("Busy...") execute_pool( ctx=self._mp_ctx, - func=geom_worker, + func=worker_geom.worker, jobs=jobs, - threads=self.nthreads, + threads=self.threads, ) _e = time.time() - _s @@ -333,13 +337,5 @@ def run( Path(self.cfg.get("output.path"), "missing.log"), ) - logger.info("Producing model output from temporary files") - # Patch output from the seperate processes back together - self.resolve() logger.info(f"Output generated in: '{self.cfg.get('output.path')}'") - - if not self._keep_temp: - logger.info("Deleting temporary files...") - self._clean_up() - logger.info("Geom calculation are done!") diff --git a/src/fiat/models/grid.py b/src/fiat/models/grid.py index 24c2c3fc..bceba258 100644 --- a/src/fiat/models/grid.py +++ b/src/fiat/models/grid.py @@ -6,11 +6,13 @@ check_exp_grid_dmfs, check_grid_exact, ) +from fiat.gis import grid +from fiat.gis.crs import get_srs_repr from fiat.io import open_grid from fiat.log import spawn_logger +from fiat.models import worker_grid from fiat.models.base import BaseModel -from fiat.models.util import execute_pool, generate_jobs -from fiat.models.worker import grid_worker_exact, grid_worker_risk +from fiat.models.util import GRID_PREFER, execute_pool, generate_jobs logger = spawn_logger("fiat.model.grid") @@ -34,16 +36,62 @@ def __init__( ): super().__init__(cfg) + # Declare + self.equal = True + + # Setup the model self.read_exposure_grid() + self.create_equal_grids() def __del__(self): BaseModel.__del__(self) - def _clean_up(self): + def _setup_output_files(self): pass + def create_equal_grids(self): + """Make the hazard and exposure grid equal spatially if necessary.""" + if self.equal: + return + + # Get which way is preferred to reproject + prefer = self.cfg.get("global.grid.prefer", "exposure") + if prefer not in ["hazard", "exposure"]: + raise ValueError( + f"Preference value {prefer} not known. Chose from \ +'hazard' or 'exposure'." + ) + prefer_bool = prefer == "exposure" + + # Setup the data sets + data = self.exposure_grid + data_warp = self.hazard_grid + if not prefer_bool: + data = (self.hazard_grid,) + data_warp = self.exposure_grid + + # Reproject the data + logger.info( + f"Reprojecting {GRID_PREFER[not prefer_bool]} \ +data to {prefer} data" + ) + data_warped = grid.reproject( + data_warp, + get_srs_repr(data.get_srs()), + data.get_geotransform(), + *data.shape_xy, + ) + + # Set the output + if prefer_bool: + self.hazard_grid = data_warped + self.cfg.set("hazard.file", data_warped.path) + else: + self.exposure_grid = data_warped + self.cfg.set("exposure.grid.file", data_warped.path) + def read_exposure_grid(self): - """_summary_.""" + """Read the exposure grid.""" file = self.cfg.get("exposure.grid.file") logger.info(f"Reading exposure grid ('{file.name}')") # Set the extra arguments from the settings file @@ -52,13 +100,13 @@ def read_exposure_grid(self): self.cfg.generate_kwargs("exposure.grid.settings"), ) kw.update( - self.cfg.generate_kwargs("global.grid"), + self.cfg.generate_kwargs("global.grid.chunk"), ) data = open_grid(file, **kw) ## checks logger.info("Executing exposure data checks...") # Check exact overlay of exposure and hazard - check_grid_exact(self.hazard_grid, data) + self.equal = check_grid_exact(self.hazard_grid, data) # Check if all damage functions are correct check_exp_grid_dmfs( data, @@ -67,12 +115,6 @@ def read_exposure_grid(self): self.exposure_grid = data - def _set_num_threads(self): - pass - - def _setup_output_files(self): - pass - def resolve(self): """Create EAD output from the outputs of different return periods. @@ -86,7 +128,7 @@ def resolve(self): # Time the function _s = time.time() - grid_worker_risk( + worker_grid.worker_ead( self.cfg, self.exposure_grid.chunk, ) @@ -110,15 +152,13 @@ def run(self): } ) - logger.info(f"Using number of threads: {self.nthreads}") - # Execute the jobs _s = time.time() logger.info("Busy...") - pcount = min(self.max_threads, self.hazard_grid.size) + pcount = min(self.threads, self.hazard_grid.size) execute_pool( ctx=self._mp_ctx, - func=grid_worker_exact, + func=worker_grid.worker, jobs=jobs, threads=pcount, ) diff --git a/src/fiat/models/util.py b/src/fiat/models/util.py index 40197167..0f86581d 100644 --- a/src/fiat/models/util.py +++ b/src/fiat/models/util.py @@ -6,10 +6,54 @@ from pathlib import Path from typing import Callable, Generator -from fiat.util import NEWLINE_CHAR +from osgeo import ogr -GEOM_MIN_CHUNK = 50000 -GEOM_MIN_WRITE_CHUNK = 20000 +from fiat.io import TableLazy +from fiat.util import NEWLINE_CHAR, replace_empty + +GEOM_DEFAULT_CHUNK = 50000 +GRID_PREFER = { + False: "hazard", + True: "exposure", +} + + +def exposure_from_geom( + ft: ogr.Feature, + exp: TableLazy, + oid: int, + idxs_haz: list | tuple, + pattern: object, +): + """_summary_.""" + method = ft.GetField("extract_method") + haz = [ft.GetField(idx) for idx in idxs_haz] + return ft, method, haz + + +def exposure_from_csv( + ft: ogr.Feature, + exp: TableLazy, + oid: int, + idxs_haz: list | tuple, + pattern: object, +): + """_summary_.""" + ft_info_raw = exp[ft.GetField(oid)] + if ft_info_raw is None: + return None, None, None + + ft_info = replace_empty(pattern.split(ft_info_raw)) + ft_info = [x(y) for x, y in zip(exp.dtypes, ft_info)] + method = ft_info[exp._columns["extract_method"]].lower() + haz = [ft_info[idx] for idx in idxs_haz] + return ft_info, method, haz + + +EXPOSURE_FIELDS = { + True: exposure_from_geom, + False: exposure_from_csv, +} def csv_temp_file( @@ -47,7 +91,6 @@ def csv_def_file( def geom_threads( cpu_count: int, - haz_layers: int, chunks: int, ): """_summary_. @@ -57,7 +100,7 @@ def geom_threads( n = 1 if chunks == 0: chunks = 1 - n = chunks * haz_layers + n = chunks n = min(cpu_count, n) return n diff --git a/src/fiat/models/worker.py b/src/fiat/models/worker.py deleted file mode 100644 index 9a917774..00000000 --- a/src/fiat/models/worker.py +++ /dev/null @@ -1,527 +0,0 @@ -"""Workers of Delft-FIAT.""" - -import os -from math import floor, isnan, nan -from multiprocessing.synchronize import Lock -from pathlib import Path - -from numpy import full, ravel, unravel_index, where -from osgeo import osr - -from fiat.gis import geom, overlay -from fiat.io import ( - BufferedGeomWriter, - BufferedTextWriter, - GridSource, - open_csv, - open_grid, -) -from fiat.log import LogItem, Sender -from fiat.models.calc import calc_haz, calc_risk -from fiat.util import NEWLINE_CHAR, create_windows, regex_pattern, replace_empty - - -def geom_resolve( - cfg: object, - exp: object, - exp_geom: dict, - chunk: tuple | list, - csv_lock: Lock = None, - geom_lock: Lock = None, -): - """_summary_.""" - # pid - os.getpid() - - # Numerical stuff - risk = cfg.get("hazard.risk") - rp_coef = cfg.get("hazard.rp_coefficients") - sig_decimals = cfg.get("vulnerability.round") - - # Set srs as osr object - srs = osr.SpatialReference() - srs.SetFromUserInput(cfg.get("global.crs")) - - # Reverse the _rp_coef to let them coincide with the acquired - # values from the temporary files - if rp_coef: - rp_coef.reverse() - new_cols = cfg.get("output.new_columns") - - # For the temp files - _files = {} - _paths = Path(cfg.get("output.tmp.path")).glob("*.dat") - - # Open the temporary files lazy - for p in sorted(_paths): - _d = open_csv(p, index=exp.meta["index_name"], large=True) - _files[p.stem] = _d - _d = None - - # Open stream to output csv file - writer = BufferedTextWriter( - Path(cfg.get("output.path"), cfg.get("output.csv.name")), - mode="ab", - buffer_size=100000, - lock=csv_lock, - ) - - # Loop over all the geometry source files - for key, gm in exp_geom.items(): - # Get output filename - _add = key[-1] - out_geom = Path(cfg.get(f"output.geom.name{_add}")) - - # Setup the geometry writer - geom_writer = BufferedGeomWriter( - Path(cfg.get("output.path"), out_geom), - srs, - gm.layer.GetLayerDefn(), - buffer_size=cfg.get("output.geom.settings.chunk"), - lock=geom_lock, - ) - geom_writer.create_fields(zip(new_cols, ["float"] * len(new_cols))) - - # Loop again over all the geometries - for ft in gm.reduced_iter(*chunk): - row = b"" - - oid = ft.GetField(0) - ft_info = exp[oid] - - # If no data is found in the temporary files, write None values - if ft_info is None: - geom_writer.add_feature( - ft, - dict(zip(new_cols, [None] * len(new_cols))), - ) - row += f"{oid}".encode() - row += b"," * (len(exp.columns) - 1) - row += NEWLINE_CHAR.encode() - writer.write(row) - continue - - row += ft_info.strip() - vals = [] - - # Loop over all the temporary files (loaded) to - # get the damage per object - for item in _files.values(): - row += b"," - _data = item[oid].strip().split(b",", 1)[1] - row += _data - _val = [float(num.decode()) for num in _data.split(b",")] - vals += _val - - if risk: - ead = round( - calc_risk(rp_coef, vals[-1 :: -exp._dat_len]), - sig_decimals, - ) - row += f",{ead}".encode() - vals.append(ead) - row += NEWLINE_CHAR.encode() - writer.write(row) - geom_writer.add_feature( - ft, - dict(zip(new_cols, vals)), - ) - - geom_writer.to_drive() - geom_writer = None - - writer.to_drive() - writer = None - - # Clean up gdal objects - srs = None - - # Clean up the opened temporary files - for _d in _files.keys(): - _files[_d] = None - _files = None - - -def geom_worker( - cfg: object, - queue: object, - haz: GridSource, - idx: int, - vul: object, - exp: object, - exp_geom: dict, - chunk: tuple | list, - lock: Lock = None, -): - """_summary_.""" - # Extract the hazard band as an object - band = haz[idx] - # Setup some metadata - _pat = regex_pattern(exp.delimiter) - _ref = cfg.get("hazard.elevation_reference") - _rnd = cfg.get("vulnerability.round") - vul_min = min(vul.index) - vul_max = max(vul.index) - - # Setup the write and write the header - writer = BufferedTextWriter( - Path(cfg.get("output.tmp.path"), f"{idx:03d}.dat"), - mode="ab", - buffer_size=100000, - lock=lock, - ) - - # Setup connection with the main process for missing values: - _sender = Sender(queue=queue) - - # Loop over all the datasets - for _, gm in exp_geom.items(): - # Check if there actually is data for this chunk - if chunk[0] > gm._count: - continue - - # Loop over all the geometries in a reduced manner - for ft in gm.reduced_iter(*chunk): - row = b"" - - # Acquire data from exposure database - ft_info_raw = exp[ft.GetField(0)] - if ft_info_raw is None: - _sender.emit( - LogItem( - 2, - f"Object with ID: {ft.GetField(0)} -> \ -No data found in exposure database", - ) - ) - continue - ft_info = replace_empty(_pat.split(ft_info_raw)) - ft_info = [x(y) for x, y in zip(exp.dtypes, ft_info)] - row += f"{ft_info[exp.index_col]}".encode() - - # Get the hazard data from the exposure geometrie - if ft_info[exp._columns["Extraction Method"]].lower() == "area": - res = overlay.clip(band, haz.get_srs(), haz.get_geotransform(), ft) - else: - res = overlay.pin(band, haz.get_geotransform(), geom.point_in_geom(ft)) - - res[res == band.nodata] = nan - - # Calculate the inundation - inun, redf = calc_haz( - res.tolist(), - _ref, - ft_info[exp._columns["Ground Floor Height"]], - ft_info[exp._columns["Ground Elevation"]], - ) - row += f",{round(inun, 2)},{round(redf, 2)}".encode() - - # Calculate the damage per catagory, and in total (_td) - _td = 0 - for key, col in exp.damage_function.items(): - if isnan(inun) or str(ft_info[col]) == "nan": - _d = "nan" - else: - inun = max(min(vul_max, inun), vul_min) - _df = vul[round(inun, _rnd), ft_info[col]] - _d = _df * ft_info[exp.max_potential_damage[key]] * redf - _d = round(_d, 2) - _td += _d - - row += f",{_d}".encode() - - row += f",{round(_td, 2)}".encode() - - # Write this to the buffer - row += NEWLINE_CHAR.encode() - writer.write(row) - - # Flush the buffer to the drive and close the writer - writer.to_drive() - writer = None - - -def geom_worker_no_csv( - cfg: object, - queue: object, - haz: GridSource, - idx: int, - vul: object, - exp: dict, - chunk: tuple | list, - lock: Lock, -): - """_summary_.""" - for _, gm in exp.items(): - for ft in gm.reduced_iter(*chunk): - pass - - -def grid_worker_exact( - cfg: object, - haz: GridSource, - idx: int, - vul: object, - exp: GridSource, -): - """_summary_.""" - # Set some variables for the calculations - exp_bands = [] - write_bands = [] - exp_nds = [] - dmfs = [] - band_n = "" - - # Check the band names - if haz.size > 1: - band_n = "_" + cfg.get("hazard.band_names")[idx - 1] - - # Extract the hazard band as an object - haz_band = haz[idx] - # Set the output directory - _out = cfg.get("output.path") - if cfg.get("hazard.risk"): - _out = cfg.get("output.damages.path") - - # Create the outgoing netcdf containing every exposure damages - out_src = open_grid( - Path(_out, f"output{band_n}.nc"), - mode="w", - ) - out_src.create( - exp.shape_xy, - exp.size, - exp.dtype, - options=["FORMAT=NC4", "COMPRESS=DEFLATE"], - ) - out_src.set_srs(exp.get_srs()) - out_src.set_geotransform(exp.get_geotransform()) - # Create the outgoing total damage grid - td_out = open_grid( - Path( - _out, - f"total_damages{band_n}.nc", - ), - mode="w", - ) - td_out.create( - exp.shape_xy, - 1, - exp.dtype, - options=["FORMAT=NC4", "COMPRESS=DEFLATE"], - ) - # Set the neccesary attributes - td_out.set_geotransform(exp.get_geotransform()) - td_out.set_srs(exp.get_srs()) - td_band = td_out[1] - td_noval = -0.5 * 2**128 - td_band.src.SetNoDataValue(td_noval) - - # Prepare some stuff for looping - for idx in range(exp.size): - exp_bands.append(exp[idx + 1]) - write_bands.append(out_src[idx + 1]) - exp_nds.append(exp_bands[idx].nodata) - write_bands[idx].src.SetNoDataValue(exp_nds[idx]) - dmfs.append(exp_bands[idx].get_metadata_item("damage_function")) - - # Going trough the chunks - for _w, h_ch in haz_band: - td_ch = td_band[_w] - - # Per exposure band - for idx, exp_band in enumerate(exp_bands): - e_ch = exp_band[_w] - - # See if there is any exposure data - out_ch = full(e_ch.shape, exp_nds[idx]) - e_ch = ravel(e_ch) - _coords = where(e_ch != exp_nds[idx])[0] - if len(_coords) == 0: - write_bands[idx].src.WriteArray(out_ch, *_w[:2]) - continue - - # See if there is overlap with the hazard data - e_ch = e_ch[_coords] - h_1d = ravel(h_ch) - h_1d = h_1d[_coords] - _hcoords = where(h_1d != haz_band.nodata)[0] - - if len(_hcoords) == 0: - write_bands[idx].src.WriteArray(out_ch, *_w[:2]) - continue - - # Do the calculations - _coords = _coords[_hcoords] - e_ch = e_ch[_hcoords] - h_1d = h_1d[_hcoords] - h_1d = h_1d.clip(min(vul.index), max(vul.index)) - - dmm = [vul[round(float(n), 2), dmfs[idx]] for n in h_1d] - e_ch = e_ch * dmm - - idx2d = unravel_index(_coords, *[exp._chunk]) - out_ch[idx2d] = e_ch - - # Write it to the band in the outgoing file - write_bands[idx].write_chunk(out_ch, _w[:2]) - - # Doing the total damages part - # Checking whether it has values or not - td_1d = td_ch[idx2d] - td_1d[where(td_1d == td_noval)] = 0 - td_1d += e_ch - td_ch[idx2d] = td_1d - - # Write the total damages chunk - td_band.write_chunk(td_ch, _w[:2]) - - # Flush the cache and dereference - for _w in write_bands[:]: - write_bands.remove(_w) - _w.close() - _w = None - - # Flush and close all - exp_bands = None - td_band.close() - td_band = None - td_out = None - - out_src.close() - out_src = None - - haz_band = None - - -def grid_worker_loose(): - """_summary_.""" - pass - - -def grid_worker_risk( - cfg: object, - chunk: tuple, -): - """_summary_.""" - _rp_coef = cfg.get("hazard.rp_coefficients") - _out = cfg.get("output.path") - _chunk = [floor(_n / len(_rp_coef)) for _n in chunk] - td = [] - rp = [] - - # TODO this is really fucking bad; fix in the future - # Read the data from the calculations - for _name in cfg.get("hazard.band_names"): - td.append( - open_grid( - Path(cfg.get("output.damages.path"), f"total_damages_{_name}.nc"), - chunk=_chunk, - mode="r", - ) - ) - rp.append( - open_grid( - Path(cfg.get("output.damages.path"), f"total_damages_{_name}.nc"), - chunk=_chunk, - mode="r", - ) - ) - - # Create the estimatied annual damages output file - exp_bands = {} - write_bands = [] - exp_nds = [] - ead_src = open_grid( - Path(_out, "ead.nc"), - mode="w", - ) - ead_src.create( - rp[0].shape_xy, - rp[0].size, - rp[0].dtype, - options=["FORMAT=NC4", "COMPRESS=DEFLATE"], - ) - ead_src.set_srs(rp[0].get_srs()) - ead_src.set_geotransform(rp[0].get_geotransform()) - - # Gather and set information before looping through windows. - for idx in range(rp[0].size): - exp_bands[idx] = [obj[idx + 1] for obj in rp] - write_bands.append(ead_src[idx + 1]) - exp_nds.append(rp[0][idx + 1].nodata) - write_bands[idx].src.SetNoDataValue(exp_nds[idx]) - - # Do the calculation for the EAD - for idx, rpx in exp_bands.items(): - for _w in create_windows(rp[0].shape, _chunk): - ead_ch = write_bands[idx][_w] - # check for one - d_ch = rpx[0][_w] - d_1d = ravel(d_ch) - _coords = where(d_1d != exp_nds[0])[0] - - # Check if something is there - if len(_coords) == 0: - continue - - data = [_data[_w] for _data in rpx] - data = [ravel(_data)[_coords] for _data in data] - data = calc_risk(_rp_coef, data) - idx2d = unravel_index(_coords, *[_chunk]) - ead_ch[idx2d] = data - write_bands[idx].write_chunk(ead_ch, _w[:2]) - - rpx = None - - # Do some cleaning - exp_bands = None - for _w in write_bands[:]: - write_bands.remove(_w) - _w.close() - _w = None - ead_src.close() - ead_src = None - - # Create ead total outgoing dataset - td_src = open_grid( - Path(_out, "ead_total.nc"), - mode="w", - ) - td_src.create( - td[0].shape_xy, - 1, - td[0].dtype, - options=["FORMAT=NC4", "COMPRESS=DEFLATE"], - ) - td_src.set_srs(td[0].get_srs()) - td_src.set_geotransform(td[0].get_geotransform()) - td_band = td_src[1] - td_noval = -0.5 * 2**128 - td_band.src.SetNoDataValue(td_noval) - - # Do the calculations for total damages - for _w in create_windows(td[0].shape, _chunk): - # Get the data - td_ch = td_band[_w] - data = [_data[1][_w] for _data in td] - d_1d = ravel(data[0]) - _coords = where(d_1d != td[0][1].nodata)[0] - - # Check whether there is data to begin with - if len(_coords) == 0: - continue - - # Get data, calc risk and write it. - data = [ravel(_i)[_coords] for _i in data] - data = calc_risk(_rp_coef, data) - idx2d = unravel_index(_coords, *[_chunk]) - td_ch[idx2d] = data - td_band.write_chunk(td_ch, _w[:2]) - - # Cleaning up afterwards - td = None - td_band.close() - td_band = None - td_src.close() - td_src = None diff --git a/src/fiat/models/worker_geom.py b/src/fiat/models/worker_geom.py new file mode 100644 index 00000000..8f427206 --- /dev/null +++ b/src/fiat/models/worker_geom.py @@ -0,0 +1,168 @@ +"""Worker function for the geometry model (no csv).""" + +import importlib +from math import nan +from multiprocessing.queues import Queue +from multiprocessing.synchronize import Lock +from pathlib import Path +from typing import Callable + +from fiat.gis import geom, overlay +from fiat.io import ( + BufferedGeomWriter, + BufferedTextWriter, + GridSource, + Table, + TableLazy, +) +from fiat.log import LogItem, Sender +from fiat.methods.ead import calc_ead, risk_density +from fiat.util import DummyWriter, regex_pattern + + +def worker( + cfg, + queue: Queue, + haz: GridSource, + vul: Table, + exp_func: Callable, + exp_data: TableLazy, + exp_geom: dict, + chunk: tuple | list, + lock1: Lock, + lock2: Lock, +): + """_summary_.""" + # Setup the hazard type module + sender = Sender(queue=queue) + module = importlib.import_module(f"fiat.methods.{cfg.get('global.type')}") + func_hazard = getattr(module, "calculate_hazard") + func_damage = getattr(module, "calculate_damage") + man_columns = getattr(module, "MANDATORY_COLUMNS") + man_entries = tuple(getattr(module, "MANDATORY_ENTRIES").values()) + + # Get the bands to prevent object creation while looping + bands = [(haz[idx + 1], idx + 1) for idx in range(haz.size)] + + # More meta data + cfg_entries = [cfg.get(item) for item in man_entries] + index_col = cfg.get("exposure.geom.settings.index") + risk = cfg.get("hazard.risk", False) + rounding = cfg.get("vulnerability.round") + vul_min = min(vul.index) + vul_max = max(vul.index) + + if risk: + rp_coef = risk_density(cfg.get("hazard.return_periods")) + rp_coef.reverse() + + pattern = None + out_text_writer = DummyWriter() + if exp_data is not None: + man_columns = [exp_data.columns.index(item) for item in man_columns] + pattern = regex_pattern(exp_data.delimiter) + out_text_writer = BufferedTextWriter( + Path(cfg.get("output.path"), cfg.get("output.csv.name")), + mode="ab", + buffer_size=100000, + lock=lock1, + ) + + # Loop through the different files + for idx, gm in exp_geom.items(): + # Check if there actually is data for this chunk + if chunk[0] > gm._count: + continue + + # Get the object id column index + oid = gm.fields.index(index_col) + + # Some meta for the specific geometry file + field_meta = cfg.get("_exposure_meta")[idx] + slen = field_meta["slen"] + total_idx = field_meta["total_idx"] + types = field_meta["types"] + idxs = field_meta["idxs"] + + # Setup the dataset buffer writer + out_geom = Path(cfg.get(f"output.geom.name{idx}")) + out_writer = BufferedGeomWriter( + Path(cfg.get("output.path"), out_geom), + gm.get_srs(), + buffer_size=cfg.get("global.geom.chunk"), + lock=lock2, + ) + + # Loop over all the geometries in a reduced manner + for ft in gm.reduced_iter(*chunk): + out = [] + info, method, haz_kwargs = exp_func( + ft, + exp_data, + oid, + man_columns, + pattern, + ) + if info is None: + sender.emit( + LogItem( + 2, + f"Object with ID: {ft.GetField(oid)} -> \ +No data found in exposure database", + ) + ) + continue + for band, bn in bands: + # How to get the hazard data + if method == "area": + res = overlay.clip(band, haz.get_srs(), haz.get_geotransform(), ft) + else: + res = overlay.pin( + band, haz.get_geotransform(), geom.point_in_geom(ft) + ) + + res[res == band.nodata] = nan + + haz_value, red_fact = func_hazard( + res.tolist(), + *cfg_entries, + *haz_kwargs, + ) + out += [haz_value, red_fact] + for key, item in types.items(): + out += func_damage( + haz_value, + red_fact, + info, + item, + vul, + vul_min, + vul_max, + rounding, + ) + + # At last do (if set) risk calculation + if risk: + i = 0 + for ti in total_idx: + ead = round( + calc_ead(rp_coef, out[ti - i :: -slen]), + rounding, + ) + out.append(ead) + i += 1 + + # Write the feature to the in memory dataset + out_writer.add_feature_with_map( + ft, + zip( + idxs, + out, + ), + ) + out_text_writer.write_iterable(info, out) + + pass + out_writer.close() + out_writer = None + pass diff --git a/src/fiat/models/worker_grid.py b/src/fiat/models/worker_grid.py new file mode 100644 index 00000000..de817102 --- /dev/null +++ b/src/fiat/models/worker_grid.py @@ -0,0 +1,282 @@ +"""Worker functions for grid model.""" +from math import floor +from pathlib import Path + +from numpy import full, ravel, unravel_index, where + +from fiat.io import ( + GridSource, + open_grid, +) +from fiat.methods.ead import calc_ead, risk_density +from fiat.util import create_windows + + +def worker( + cfg: object, + haz: GridSource, + idx: int, + vul: object, + exp: GridSource, +): + """_summary_.""" + # Set some variables for the calculations + exp_bands = [] + write_bands = [] + exp_nds = [] + dmfs = [] + band_n = "" + + # Check the band names + if haz.size > 1: + band_n = "_" + cfg.get("hazard.band_names")[idx - 1] + + # Extract the hazard band as an object + haz_band = haz[idx] + # Set the output directory + _out = cfg.get("output.path") + if cfg.get("hazard.risk"): + _out = cfg.get("output.damages.path") + + # Create the outgoing netcdf containing every exposure damages + out_src = open_grid( + Path(_out, f"output{band_n}.nc"), + mode="w", + ) + out_src.create( + exp.shape_xy, + exp.size, + exp.dtype, + options=["FORMAT=NC4", "COMPRESS=DEFLATE"], + ) + out_src.set_srs(exp.get_srs()) + out_src.set_geotransform(exp.get_geotransform()) + # Create the outgoing total damage grid + td_out = open_grid( + Path( + _out, + f"total_damages{band_n}.nc", + ), + mode="w", + ) + td_out.create( + exp.shape_xy, + 1, + exp.dtype, + options=["FORMAT=NC4", "COMPRESS=DEFLATE"], + ) + # Set the neccesary attributes + td_out.set_geotransform(exp.get_geotransform()) + td_out.set_srs(exp.get_srs()) + td_band = td_out[1] + td_noval = -0.5 * 2**128 + td_band.src.SetNoDataValue(td_noval) + + # Prepare some stuff for looping + for idx in range(exp.size): + exp_bands.append(exp[idx + 1]) + write_bands.append(out_src[idx + 1]) + exp_nds.append(exp_bands[idx].nodata) + write_bands[idx].src.SetNoDataValue(exp_nds[idx]) + dmfs.append(exp_bands[idx].get_metadata_item("fn_damage")) + + # Going trough the chunks + for _w, h_ch in haz_band: + td_ch = td_band[_w] + + # Per exposure band + for idx, exp_band in enumerate(exp_bands): + e_ch = exp_band[_w] + + # See if there is any exposure data + out_ch = full(e_ch.shape, exp_nds[idx]) + e_ch = ravel(e_ch) + _coords = where(e_ch != exp_nds[idx])[0] + if len(_coords) == 0: + write_bands[idx].src.WriteArray(out_ch, *_w[:2]) + continue + + # See if there is overlap with the hazard data + e_ch = e_ch[_coords] + h_1d = ravel(h_ch) + h_1d = h_1d[_coords] + _hcoords = where(h_1d != haz_band.nodata)[0] + + if len(_hcoords) == 0: + write_bands[idx].src.WriteArray(out_ch, *_w[:2]) + continue + + # Do the calculations + _coords = _coords[_hcoords] + e_ch = e_ch[_hcoords] + h_1d = h_1d[_hcoords] + h_1d = h_1d.clip(min(vul.index), max(vul.index)) + + dmm = [vul[round(float(n), 2), dmfs[idx]] for n in h_1d] + e_ch = e_ch * dmm + + idx2d = unravel_index(_coords, *[exp._chunk]) + out_ch[idx2d] = e_ch + + # Write it to the band in the outgoing file + write_bands[idx].write_chunk(out_ch, _w[:2]) + + # Doing the total damages part + # Checking whether it has values or not + td_1d = td_ch[idx2d] + td_1d[where(td_1d == td_noval)] = 0 + td_1d += e_ch + td_ch[idx2d] = td_1d + + # Write the total damages chunk + td_band.write_chunk(td_ch, _w[:2]) + + # Flush the cache and dereference + for _w in write_bands[:]: + write_bands.remove(_w) + _w.close() + _w = None + + # Flush and close all + exp_bands = None + td_band.close() + td_band = None + td_out = None + + out_src.close() + out_src = None + + haz_band = None + + +def worker2(): + """_summary_.""" + pass + + +def worker_ead( + cfg: object, + chunk: tuple, +): + """_summary_.""" + _rp_coef = risk_density(cfg.get("hazard.return_periods")) + _out = cfg.get("output.path") + _chunk = [floor(_n / len(_rp_coef)) for _n in chunk] + td = [] + rp = [] + + # TODO this is really fucking bad; fix in the future + # Read the data from the calculations + for _name in cfg.get("hazard.band_names"): + td.append( + open_grid( + Path(cfg.get("output.damages.path"), f"total_damages_{_name}.nc"), + chunk=_chunk, + mode="r", + ) + ) + rp.append( + open_grid( + Path(cfg.get("output.damages.path"), f"total_damages_{_name}.nc"), + chunk=_chunk, + mode="r", + ) + ) + + # Create the estimatied annual damages output file + exp_bands = {} + write_bands = [] + exp_nds = [] + ead_src = open_grid( + Path(_out, "ead.nc"), + mode="w", + ) + ead_src.create( + rp[0].shape_xy, + rp[0].size, + rp[0].dtype, + options=["FORMAT=NC4", "COMPRESS=DEFLATE"], + ) + ead_src.set_srs(rp[0].get_srs()) + ead_src.set_geotransform(rp[0].get_geotransform()) + + # Gather and set information before looping through windows. + for idx in range(rp[0].size): + exp_bands[idx] = [obj[idx + 1] for obj in rp] + write_bands.append(ead_src[idx + 1]) + exp_nds.append(rp[0][idx + 1].nodata) + write_bands[idx].src.SetNoDataValue(exp_nds[idx]) + + # Do the calculation for the EAD + for idx, rpx in exp_bands.items(): + for _w in create_windows(rp[0].shape, _chunk): + ead_ch = write_bands[idx][_w] + # check for one + d_ch = rpx[0][_w] + d_1d = ravel(d_ch) + _coords = where(d_1d != exp_nds[0])[0] + + # Check if something is there + if len(_coords) == 0: + continue + + data = [_data[_w] for _data in rpx] + data = [ravel(_data)[_coords] for _data in data] + data = calc_ead(_rp_coef, data) + idx2d = unravel_index(_coords, *[_chunk]) + ead_ch[idx2d] = data + write_bands[idx].write_chunk(ead_ch, _w[:2]) + + rpx = None + + # Do some cleaning + exp_bands = None + for _w in write_bands[:]: + write_bands.remove(_w) + _w.close() + _w = None + ead_src.close() + ead_src = None + + # Create ead total outgoing dataset + td_src = open_grid( + Path(_out, "ead_total.nc"), + mode="w", + ) + td_src.create( + td[0].shape_xy, + 1, + td[0].dtype, + options=["FORMAT=NC4", "COMPRESS=DEFLATE"], + ) + td_src.set_srs(td[0].get_srs()) + td_src.set_geotransform(td[0].get_geotransform()) + td_band = td_src[1] + td_noval = -0.5 * 2**128 + td_band.src.SetNoDataValue(td_noval) + + # Do the calculations for total damages + for _w in create_windows(td[0].shape, _chunk): + # Get the data + td_ch = td_band[_w] + data = [_data[1][_w] for _data in td] + d_1d = ravel(data[0]) + _coords = where(d_1d != td[0][1].nodata)[0] + + # Check whether there is data to begin with + if len(_coords) == 0: + continue + + # Get data, calc risk and write it. + data = [ravel(_i)[_coords] for _i in data] + data = calc_ead(_rp_coef, data) + idx2d = unravel_index(_coords, *[_chunk]) + td_ch[idx2d] = data + td_band.write_chunk(td_ch, _w[:2]) + + # Cleaning up afterwards + td = None + td_band.close() + td_band = None + td_src.close() + td_src = None diff --git a/src/fiat/util.py b/src/fiat/util.py index 43b03eb9..236ca56e 100644 --- a/src/fiat/util.py +++ b/src/fiat/util.py @@ -1,6 +1,7 @@ """Base FIAT utility.""" import ctypes +import fnmatch import math import os import platform @@ -13,7 +14,7 @@ from types import FunctionType, ModuleType import regex -from osgeo import gdal +from osgeo import gdal, ogr # Define the variables for FIAT BLACKLIST = type, ModuleType, FunctionType @@ -45,6 +46,12 @@ "str": str, } +_fields_type_map = { + "int": ogr.OFTInteger64, + "float": ogr.OFTReal, + "str": ogr.OFTString, +} + def regex_pattern( delimiter: str, @@ -152,9 +159,69 @@ def flatten_dict(d: MutableMapping, parent_key: str = "", sep: str = "."): # Exposure specific utility -def gen_new_columns(): +def discover_exp_columns( + columns: dict, + type: str, +): """_summary_.""" - pass + dmg_idx = {} + + # Get column values + column_vals = list(columns.keys()) + + # Filter the current columns + dmg = fnmatch.filter(column_vals, f"fn_{type}_*") + dmg_suffix = [item.split("_")[-1].strip() for item in dmg] + mpd = fnmatch.filter(column_vals, f"max_{type}_*") + mpd_suffix = [item.split("_")[-1].strip() for item in mpd] + + # Check the overlap + _check = [item in mpd_suffix for item in dmg_suffix] + + # Determine the missing values + missing = [item for item, b in zip(dmg_suffix, _check) if not b] + for item in missing: + dmg_suffix.remove(item) + + fn = {} + maxv = {} + for val in dmg_suffix: + fn.update({val: columns[f"fn_{type}_{val}"]}) + maxv.update({val: columns[f"max_{type}_{val}"]}) + dmg_idx.update({"fn": fn, "max": maxv}) + + return dmg_suffix, dmg_idx, missing + + +def generate_output_columns( + specific_columns: tuple | list, + exposure_types: dict, + extra: tuple | list = [], + suffix: tuple | list = [""], +): + """_summary_.""" + default = specific_columns + ["red_fact"] + total_idx = [] + + # Loop over the exposure types + for key, value in exposure_types.items(): + default += [f"{key}_{item}" for item in value["fn"].keys()] + total_idx.append(len(default)) + default += [f"total_{key}"] + + total_idx = [item - len(default) for item in total_idx] + + out = [] + if len(suffix) == 1 and not suffix[0]: + out = default + else: + for name in suffix: + add = [f"{item}_{name}" for item in default] + out += add + + out += [f"{x}_{y}" for x, y in product(extra, exposure_types.keys())] + + return out, len(default), total_idx # GIS related utility @@ -201,7 +268,9 @@ def _read_gridsource_layers_from_info( pass -def _create_geom_driver_map(): +def _create_geom_driver_map( + write: bool = False, +): """_summary_.""" geom_drivers = {} _c = gdal.GetDriverCount() @@ -209,6 +278,9 @@ def _create_geom_driver_map(): for idx in range(_c): dr = gdal.GetDriver(idx) if dr.GetMetadataItem(gdal.DCAP_VECTOR): + edit = dr.GetMetadataItem(gdal.DCAP_DELETE_FIELD) + if write and edit is None: + continue if dr.GetMetadataItem(gdal.DCAP_CREATE) or dr.GetMetadataItem( gdal.DCAP_CREATE_LAYER ): @@ -230,8 +302,9 @@ def _create_geom_driver_map(): return geom_drivers -GEOM_DRIVER_MAP = _create_geom_driver_map() -GEOM_DRIVER_MAP[""] = "Memory" +GEOM_READ_DRIVER_MAP = _create_geom_driver_map() +GEOM_WRITE_DRIVER_MAP = _create_geom_driver_map(write=True) +GEOM_WRITE_DRIVER_MAP[""] = "Memory" def _create_grid_driver_map(): @@ -339,7 +412,48 @@ def generic_path_check( return path +# Logging utility +def progressbar( + iteration: int, + total: int, + prefix: str = "", + suffix: str = "", + decimals: int = 1, + bar_length: int = 50, +): + """Call in a loop to create terminal progress bar. + + @params: + iteration - Required : current iteration (Int) + total - Required : total iterations (Int) + prefix - Optional : prefix string (Str) + suffix - Optional : suffix string (Str) + decimals - Optional : positive number of decimals in percent complete (Int) + bar_length - Optional : character length of bar (Int) + """ + str_format = "{0:." + str(decimals) + "f}" + percents = str_format.format(100 * (iteration / float(total))) + filled_length = int(round(bar_length * iteration / float(total))) + bar = "█" * filled_length + "-" * (bar_length - filled_length) + + (sys.stdout.write("\r%s |%s| %s%s %s" % (prefix, bar, percents, "%", suffix)),) + + if iteration == total: + sys.stdout.write("\n") + sys.stdout.flush() + + # Misc. +def find_duplicates(elements: tuple | list): + """Find duplicate elements in an iterable object.""" + uni = list(set(elements)) + counts = [elements.count(elem) for elem in uni] + dup = [elem for _i, elem in enumerate(uni) if counts[_i] > 1] + if not dup: + return None + return dup + + def object_size(obj): """Calculate the actual size of an object (bit overestimated). @@ -392,6 +506,25 @@ def release(self): pass +class DummyWriter: + """Mimic the behaviour of an object that is capable of writing.""" + + def __init__(self, *args, **kwargs): + pass + + def close(self): + """Call dummy close.""" + pass + + def write(self, *args): + """Call dummy write.""" + pass + + def write_iterable(self, *args): + """Call dummy write iterable.""" + pass + + # Typing related stuff def deter_type( e: bytes, diff --git a/test/test_config.py b/test/test_config.py index 8c66ed74..76e58086 100644 --- a/test/test_config.py +++ b/test/test_config.py @@ -7,9 +7,6 @@ def test_settings(settings_files): # Assert path to itself assert cfg.filepath.name == "geom_risk.toml" - # Assert a global setting - assert cfg.get("global.keep_temp_files") - # Assert generated kwargs functionality haz_kw = cfg.generate_kwargs("hazard.settings") assert "var_as_band" in haz_kw diff --git a/test/test_gis.py b/test/test_gis.py index fc10a009..6b44ce67 100644 --- a/test/test_gis.py +++ b/test/test_gis.py @@ -34,13 +34,13 @@ def test_reproject(tmp_path, geom_data, grid_event_data): new_gm = geom.reproject( geom_data, dst_crs, - str(tmp_path), + out_dir=str(tmp_path), ) new_gr = grid.reproject( grid_event_data, dst_crs, - str(tmp_path), + out_dir=str(tmp_path), ) assert new_gm.get_srs().GetAuthorityCode(None) == "3857" diff --git a/test/test_io.py b/test/test_io.py index b437c470..a71901f8 100644 --- a/test/test_io.py +++ b/test/test_io.py @@ -13,19 +13,19 @@ def test_bufferedgeom(tmp_path, geom_data): ) assert writer.size == 0 - writer.add_feature( + writer.add_feature_with_map( geom_data.layer.GetFeature(1), {}, ) assert writer.size == 1 - writer.add_feature( + writer.add_feature_with_map( geom_data.layer.GetFeature(2), {}, ) assert writer.size == 2 - writer.add_feature( + writer.add_feature_with_map( geom_data.layer.GetFeature(3), {}, ) diff --git a/test/test_logic.py b/test/test_logic.py index b1a34cf2..ffca226e 100644 --- a/test/test_logic.py +++ b/test/test_logic.py @@ -1,20 +1,45 @@ -from fiat.models.calc import calc_haz, calc_risk, calc_rp_coef +from fiat.methods.ead import calc_ead, risk_density +from fiat.methods.flood import calculate_hazard def test_calc_haz(): - dmg, red_f = calc_haz([2.5, 5, 10], ref="dem", gfh=1.0, ge=0, method="mean") + dmg, red_f = calculate_hazard( + [2.5, 5, 10], + reference="dem", + ground_flht=1.0, + ground_elevtn=0, + method="mean", + ) assert int(dmg * 100) == 483 assert int(red_f) == 1 - dmg, red_f = calc_haz([0, 2.5, 5, 10], ref="dem", gfh=1.0, ge=0, method="mean") + dmg, red_f = calculate_hazard( + [0, 2.5, 5, 10], + reference="dem", + ground_flht=1.0, + ground_elevtn=0, + method="mean", + ) assert int(dmg * 100) == 483 assert int(red_f * 100) == 75 - dmg, red_f = calc_haz([0, 2.5, 5, 10], ref="datum", gfh=1.0, ge=1.0, method="mean") + dmg, red_f = calculate_hazard( + [0, 2.5, 5, 10], + reference="datum", + ground_flht=1.0, + ground_elevtn=1.0, + method="mean", + ) assert int(dmg * 100) == 383 assert int(red_f * 100) == 75 - dmg, red_f = calc_haz([0, 1.5, 5, 10], ref="datum", gfh=1.0, ge=1.0, method="mean") + dmg, red_f = calculate_hazard( + [0, 1.5, 5, 10], + reference="datum", + ground_flht=1.0, + ground_elevtn=1.0, + method="mean", + ) assert int(dmg * 100) == 350 assert int(red_f * 100) == 75 @@ -23,8 +48,8 @@ def test_calc_risk(): rps = [1, 2, 5, 25, 50, 100] dms = [5, 10, 50, 300, 1200, 3000] - coef = calc_rp_coef(rps) - ead = calc_risk(coef, dms) + coef = risk_density(rps) + ead = calc_ead(coef, dms) assert int(round(ead, 1) * 100) == 9850 @@ -33,8 +58,8 @@ def test_calc_risk_order(): rps = [50, 2, 100, 25, 1, 5] dms = [1200, 10, 3000, 300, 5, 50] - coef = calc_rp_coef(rps) - ead = calc_risk(coef, dms) + coef = risk_density(rps) + ead = calc_ead(coef, dms) assert int(round(ead, 1) * 100) == 9850 @@ -43,7 +68,7 @@ def test_calc_risk_one(): rps = [10] dms = [5] - coef = calc_rp_coef(rps) - ead = calc_risk(coef, dms) + coef = risk_density(rps) + ead = calc_ead(coef, dms) assert int(ead * 100) == 50 diff --git a/test/test_model.py b/test/test_model.py index 14aebf8d..2a9b1eae 100644 --- a/test/test_model.py +++ b/test/test_model.py @@ -17,9 +17,9 @@ def test_geom_event(tmp_path, configs): run_model(configs["geom_event"], tmp_path) # Check the output for this specific case - out = open_csv(Path(str(tmp_path), "output.csv"), index="Object ID") - assert int(float(out[2, "Total Damage"])) == 740 - assert int(float(out[3, "Total Damage"])) == 1038 + out = open_csv(Path(str(tmp_path), "output.csv"), index="object_id") + assert int(float(out[2, "total_damage"])) == 740 + assert int(float(out[3, "total_damage"])) == 1038 def test_geom_missing(tmp_path, configs): @@ -37,10 +37,10 @@ def test_geom_risk(tmp_path, configs): run_model(configs["geom_risk"], tmp_path) # Check the output for this specific case - out = open_csv(Path(str(tmp_path), "output.csv"), index="Object ID") - assert int(float(out[2, "Damage: Structure (5.0Y)"])) == 1804 - assert int(float(out[4, "Total Damage (10.0Y)"])) == 3840 - assert int(float(out[3, "Risk (EAD)"]) * 100) == 102247 + out = open_csv(Path(str(tmp_path), "output.csv"), index="object_id") + assert int(float(out[2, "damage_structure_5.0y"])) == 1804 + assert int(float(out[4, "total_damage_10.0y"])) == 3840 + assert int(float(out[3, "ead_damage"]) * 100) == 102247 def test_grid_event(tmp_path, configs): diff --git a/test/test_struct.py b/test/test_struct.py index da7de725..990d79d4 100644 --- a/test/test_struct.py +++ b/test/test_struct.py @@ -11,7 +11,7 @@ def test_geomsource(geom_data): bounds = [round(item * 10000) for item in bounds] assert bounds == [43550, 44400, 519600, 520450] - assert geom_data.fields == ["Object ID", "ObjectName"] + assert geom_data.fields == ["object_id", "object_name"] srs = geom_data.get_srs() assert srs.GetAuthorityCode(None) == "4326" diff --git a/test/test_worker.py b/test/test_worker.py deleted file mode 100644 index e098072d..00000000 --- a/test/test_worker.py +++ /dev/null @@ -1,98 +0,0 @@ -from pathlib import Path - -import pytest -from fiat.models.worker import ( - geom_resolve, - geom_worker, - grid_worker_exact, - grid_worker_risk, -) - - -@pytest.mark.dependency() -def test_geom_worker(tmp_path, geom_risk): - # Set model output directory - model = geom_risk - model.cfg.set_output_dir(Path(str(tmp_path), "..", "worker_geom")) - - # Create the files: - model._setup_output_files() - - # Invoke the worker directly - for idx in range(model.hazard_grid.size): - geom_worker( - model.cfg, - None, - model.hazard_grid, - idx + 1, - model.vulnerability_data, - model.exposure_data, - model.exposure_geoms, - model.chunks[0], - None, - ) - - # Assert output - files = list(model.cfg.get("output.tmp.path").iterdir()) - assert len(files) == 4 - # List of files check - expected_files = [f"{idx:03d}.dat" for idx in range(1, 5)] - output_files = [_f.name for _f in files] - assert sorted(output_files) == expected_files - - -@pytest.mark.dependency(depends=["test_geom_worker"]) -def test_geom_resolve(tmp_path, geom_risk): - # Set model output directory - model = geom_risk - model.cfg.set_output_dir(Path(str(tmp_path), "..", "worker_geom")) - - # Invoke the worker directly - geom_resolve( - model.cfg, - model.exposure_data, - model.exposure_geoms, - model.chunks[0], - None, - None, - ) - - # Assert the output - assert Path(model.cfg.get("output.path"), "output.csv").exists() - assert Path(model.cfg.get("output.path"), "spatial.gpkg").exists() - - -@pytest.mark.dependency() -def test_grid_worker(tmp_path, grid_risk): - # Set model output directory - model = grid_risk - model.cfg.set_output_dir(Path(str(tmp_path), "..", "worker_grid")) - - for idx in range(model.hazard_grid.size): - grid_worker_exact( - model.cfg, - model.hazard_grid, - idx + 1, - model.vulnerability_data, - model.exposure_grid, - ) - - # Assert the output - files = list(model.cfg.get("output.damages.path").iterdir()) - assert len(files) == 8 - - -@pytest.mark.dependency(depends=["test_grid_worker"]) -def test_grid_risk_worker(tmp_path, grid_risk): - # Set model output directory - model = grid_risk - model.cfg.set_output_dir(Path(str(tmp_path), "..", "worker_grid")) - - grid_worker_risk( - model.cfg, - model.exposure_grid.chunk, - ) - - # Assert the output - assert Path(model.cfg.get("output.path"), "ead.nc").exists() - assert Path(model.cfg.get("output.path"), "ead_total.nc").exists()