Last active
January 28, 2016 14:17
-
-
Save edigiacomo/edcc4e135b0442fe42d6 to your computer and use it in GitHub Desktop.
Allerta comuni
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| # encoding: utf-8 | |
| # | |
| # Crea un GeoJSON dei comuni a partire dal raster di precipitazione in cui | |
| # aggiunge i seguenti attributi alle feature: | |
| # - preci_count: numero di celle nel comune | |
| # - preci_sum: somma della precipitazione nel comune | |
| # - preci_max: massima precipitazione nel comune | |
| # - preci_threshold_50mm_count: numero di celle nel comune che superano i 50mm | |
| # - preci_threshold_50mm_alert: "green" se nessuna cella supera i 50mm, | |
| # "yellow" se la superano meno di dieci, "red" se la superano 10 o più | |
| # celle. | |
| import json | |
| import numpy as np | |
| from rasterstats import zonal_stats | |
| def threshold_count(x): | |
| c = (x > 50.0).sum() | |
| # If all the values are masked return 0 | |
| if c is np.ma.masked: | |
| return 0 | |
| else: | |
| return int(c) | |
| def threshold_level(x): | |
| c = threshold_count(x) | |
| if c == 0: | |
| return "green" | |
| elif c < 10: | |
| return "yellow" | |
| else: | |
| return "red" | |
| vector = "comuni.json" | |
| raster = "preci.nc" | |
| outfile = "result.json" | |
| stats = zonal_stats( | |
| vector, raster, | |
| copy_properties=True, geojson_out=True, prefix="preci_", | |
| stats=("sum", "count", "max"), add_stats={ | |
| "threshold_50mm_count": threshold_count, | |
| "threshold_50mm_alert": threshold_level, | |
| }, | |
| ) | |
| with open(outfile, "w") as fp: | |
| json.dump({ | |
| "type": "FeatureCollection", | |
| "features": stats | |
| }, fp) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment