Add fire danger tool

docs/Utilities/Fire_danger.md

@@ -0,0 +1,73 @@
+[[/Utilities/Fire_danger]] -- Fire danger forecast data tool
+
+# Synopsis
+
+Shell:
+~~~
+$ gridlabd fire_danger [OPTIONS ...]
+	-f|--forecast=DAY1[,DAY2[,...] 
+	-d|--date=DATE1[,DATE2[,...]] 
+	-t|--type=TYPE1[,TYPE2[,...]]  
+~~~
+
+GLM:
+~~~
+#python -m fire_danger [OPTIONS ...]
+	-f|--forecast=DAY1[,DAY2[,...] 
+	-d|--date=DATE1[,DATE2[,...]] 
+	-t|--type=TYPE1[,TYPE2[,...]]
+~~~
+
+Python:
+~~~
+>>> import sys
+>>> sys.path.append("/usr/local/share/gridlabd/")
+>>> import fire_danger
+>>> 
+~~~
+
+# Options
+
+- `-h|--help|help`:   display this help
+- `-r|--refresh`:     force refresh of files in the cache
+- `-s|--show`:        show image instead of print image name
+
+# Description
+
+Downloads the USGS fire danger `DAY`-day forecast map for the `DATE` given. The
+map is stored in the `CACHEDIR` (by default `$GLD_ETC/usgs/firedanger/`) using
+the file name `TYPE_DATE_DAY.tif`. The output is the full pathname where the
+data is stored. `DAY`, `DATE`, and `TYPE` may be specified as comma delimited values,
+in which case, all the combination of the specified values are downloaded.
+
+DAY may be specified as 1 through 7.  DATE must be specified using the format
+YYYY-MM-DD. Valid TYPE values are:
+
+- `fpi` - Fire Potential Index
+- `lfp` - Large Fire Probability
+- `fsp` - Fire Spread Probability
+
+# Example
+
+The following example displays the 5-day fire potential index forecast map for August 20, 2021.
+
+~~~
+sh$ gridlabd fire_danger -f=5 -d=2021-08-20 -t=fpi --show
+~~~
+
+The following example obtains the image name the same data
+
+~~~
+sh$ gridlabd fire_danger -f=5 -d=2021-08-20 -t=fpi
+/usr/local/opt/gridlabd/current/share/gridlabd/usgs/firedanger/fpi_20210820_5.tif
+~~~
+
+The following example obtains the image name in Python:
+
+~~~
+>>> import sys
+>>> sys.path.append("/usr/local/share/gridlabd")
+>>> import fire_danger
+>>> fire_danger.get_data(5,'20210820','fpi')
+'/usr/local/opt/gridlabd/current/share/gridlabd/usgs/firedanger/fpi_20210820_5.tif'
+~~~

geodata/geodata_firerisk.py

@@ -0,0 +1,117 @@
+
+"""GridLAB-D Geodata Fire Danger Package
+
+This package gathers data on the 7 day wildfire forecast risk from the USGS Fire Danger Map
+
+INPUT:
+
+    latitude    the latitude of the location
+
+    longitude   the longitude of the location
+
+    day of forecast 1-7
+    date of forecast (cooresponds to day 1)
+    type of forecast:  fpi - Fire Potential Index, lfp - Large Fire Probability, fsp - Fire Spread Probability
+
+OUTPUT:
+
+    fireDanger  Firedanger Value for all listed lat,longs
+
+
+OPTIONS:
+
+"""
+
+import sys
+import os
+from click import option
+import pandas
+import rasterio as rio
+from rasterio.plot import show
+from pyproj import Transformer
+from datetime import datetime
+
+today= datetime.today().strftime('%Y%m%d')
+
+share = os.getenv("GLD_ETC")
+if not share:
+    share = "/usr/local/share/gridlabd"
+if share not in sys.path:
+    sys.path.append("/usr/local/share/gridlabd")
+import fire_danger
+
+#
+# Defaults
+#
+default_options = dict(
+    day= "1",
+    date= datetime.today().strftime('%Y%m%d'),
+    type ="fpi",
+    )
+
+default_config = dict(
+    cachedir = f"{share}/usgs/firedanger"
+    )
+
+#
+# Implementation of fireRisk package
+#
+def apply(data, options=default_options, config=default_config, warning=print):
+    """Gather fire danger data for various locations
+    ARGUMENTS:
+
+        data (pandas.DataFrame)
+
+            The data frame must contain `latitude` and `longitude` fields for which the firedanger 
+                will be calculated
+
+        options (dict)
+
+            'day' = # of days ahead of forecast date
+            'date'= date of forecast given in format 'YYYYMMDD'
+            'type'= type of forecast (fpi, lfp, fsp)
+
+        config (dict)
+
+            There are no configuration options
+
+    RETURNS:
+
+        list object of fire_danger values
+
+            The specified fire danger value for the geocodes provided.
+
+    """
+    USGSURL = "https://edcintl.cr.usgs.gov/downloads/sciweb1/shared/firedanger/download-tool/source_rasters/w{TYPE}-forecast-{DAYAHEAD}/emodis-w{TYPE}-forecast-{DAYAHEAD}_data_{DATE}_{DATE}.zip"
+
+    filePath= fire_danger.get_data(options['day'],options['date'],options['type'],USGSURL,config['cachedir'])
+    with rio.open(filePath) as usgsMap:
+        transformer= Transformer.from_crs("epsg:4326",usgsMap.crs)
+        coords= [transformer.transform(x, y) for x,y in data.itertuples(index=False)]
+        coordsRC= []
+        band= usgsMap.read(1)
+        for x, y in coords:            
+            row,col= usgsMap.index(x,y)
+            coordsRC.append(band[row,col])
+    return coordsRC
+
+#
+# Perform validation tests
+#
+if __name__ == '__main__':
+
+    import unittest
+
+    class TestFireRisk(unittest.TestCase):
+
+        def test_firerisk(self):
+            testData = pandas.DataFrame({
+                "lat":[34.8020,36.636],
+                "long":[-104.1632, -121.9298],
+                })
+            testOptions=dict(day="4",date='20220415',type="fpi")
+            result = apply(data=testData,options=testOptions)
+            # print(result) 
+            self.assertEqual(result,[114,0])
+
+    unittest.main()

tools/Makefile.mk

@@ -2,6 +2,7 @@ dist_pkgdata_DATA += tools/create_filter.py
 dist_pkgdata_DATA += tools/create_player.py
 dist_pkgdata_DATA += tools/create_poles.py
 dist_pkgdata_DATA += tools/eia_recs.py
+dist_pkgdata_DATA += tools/fire_danger.py
 dist_pkgdata_DATA += tools/fit_filter.py
 dist_pkgdata_DATA += tools/gridlabd-editor.png
 dist_pkgdata_DATA += tools/gridlabd-editor.py

tools/autotest/test_fire_danger.glm

@@ -0,0 +1,6 @@
+#python -m fire_danger -f=5 -d=2021-08-20 -t=fpi
+#system cp /usr/local/opt/gridlabd/current/share/gridlabd/usgs/firedanger/fpi_20210820_5.tif test_fire_danger.tif
+
+#ifexist ../test_fire_danger.tif
+#on_exit 0 diff -q ../test_fire_danger.tif test_fire_danger.tif
+#endif

tools/fire_danger.py

@@ -0,0 +1,134 @@
+"""Fire danger forecast data tool
+
+SYNTAX
+
+	$ gridlabd fire_danger [OPTIONS ...]
+		-f|--forecast=DAY1[,DAY2[,...] 
+		-d|--date=DATE1[,DATE2[,...]] 
+		-t|--type=TYPE1[,TYPE2[,...]]  
+
+OPTIONS:
+
+	-h|--help|help   display this help
+	-r|--refresh     force refresh of files in the cache
+	-s|--show        show image instead of print image name
+
+DESCRIPTION
+
+Downloads the USGS fire danger DAY ahead forecast map for the DATE given. The
+map is stored in the CACHEDIR (by default `$GLD_ETC/usgs/firedanger/`) using
+the file name `TYPE_DATE_DAY.tif`. The output is the full pathname where the
+data is stored. DAY, DATE, and TYPE may be specified as comma delimited values,
+in which case, all the combination of the specified values are downloaded.
+
+DAY may be specified as 1 through 7.  DATE must be specified using the format
+YYYY-MM-DD. Valid TYPE values are:
+
+	fpi - Fire Potential Index
+	lfp - Large Fire Probability
+	fsp - Fire Spread Probability
+"""
+
+import sys, os
+import requests
+from io import BytesIO
+import zipfile
+import datetime
+import rasterio
+from rasterio.plot import show
+
+CACHEDIR = "/usr/local/opt/gridlabd/current/share/gridlabd/usgs/firedanger"
+USGSURL = "https://edcintl.cr.usgs.gov/downloads/sciweb1/shared/firedanger/download-tool/source_rasters/w{TYPE}-forecast-{DAYAHEAD}/emodis-w{TYPE}-forecast-{DAYAHEAD}_data_{DATE}_{DATE}.zip"
+REFRESH = False
+SHOW = False
+
+class FireDangerBadRequest(Exception):
+	pass
+
+class FireDangerInvalidOption(Exception):
+	pass
+
+class FireDangerMissingOption(Exception):
+	pass
+
+def get_data(dayahead,date,maptype,url=USGSURL,cachedir=CACHEDIR):
+
+	filename = f"{cachedir}/{maptype}_{date}_{dayahead}.tif"
+	if not os.path.exists(filename) or REFRESH:
+		usgsurl = url.format(DAYAHEAD=dayahead,DATE=date,TYPE=maptype)
+		reply = requests.get(usgsurl,stream=True)
+		if 200 <= reply.status_code < 300:
+			archive = zipfile.ZipFile(BytesIO(reply.content),mode="r")
+			os.makedirs(cachedir,exist_ok=True)
+			cachename = f"emodis-w{maptype}_data_{date}_{date}.tiff"
+			archive.extract(cachename,f"{cachedir}")
+			os.rename(f"{cachedir}/{cachename}",filename)
+		else:		
+			raise FireDangerBadRequest(f"{usgsurl} (code {reply.status_code})")
+	return filename
+
+def main(args):
+
+	if not args:
+		print("Syntax: gridlabd fire_danger -f|--forecast=DAYAHEAD -d|--date=YYYY-MM-DD -t|--type=TYPE [OPTIONS ...]",file=sys.stderr)
+		return
+
+	DAYAHEAD = None
+	DATE = None
+	TYPE = None
+	for arg in args:
+		spec = arg.split("=")
+		if len(spec) == 1:
+			tag = arg
+			value = True
+		elif len(spec) == 2:
+			tag = spec[0]
+			value = spec[1]
+		else:
+			tag = spec[0]
+			value = "=".join(spec[1:])
+		if tag in ["-d","--date"]:
+			DATE = value
+		elif tag in ["-f","--forecast"]:
+			DAYAHEAD = value
+		elif tag in ["-t","--type"]:
+			TYPE = value
+		elif tag in ["-h","--help","help"]:
+			print(__doc__,file=sys.stdout)
+			return
+		elif tag in ["-r","--refresh"]:
+			global REFRESH
+			REFRESH = True
+		elif tag in ["-s","--show"]:
+			global SHOW
+			SHOW = True
+		else:
+			raise FireDangerInvalidOption(arg)
+
+	if not DAYAHEAD:
+		raise FireDangerMissingOption("-f|--forecast=DAYAHEAD")
+
+	if not DATE:
+		raise FireDangerMissingOption("-d|--date=YYYY-MM-DD")
+
+	if not TYPE:
+		raise FireDangerMissingOption("-t|--type=TYPE")
+
+	for day in DAYAHEAD.split(","):
+		for date in DATE.split(","):
+			for mtype in TYPE.split(","):
+				filename = get_data(day,datetime.datetime.strptime(date,"%Y-%m-%d").strftime("%Y%m%d"),mtype)		
+				if SHOW:
+					show(rasterio.open(filename))
+				else:
+					print(filename,file=sys.stdout)
+
+
+if __name__ == "__main__":
+
+	BASENAME = os.path.splitext(os.path.basename(sys.argv[0]))[0]
+	try:
+		main(sys.argv[1:])
+	except Exception as err:
+		etype,evalue,etrace = sys.exc_info()
+		print(f"ERROR [{BASENAME}]: {etype.__name__} {evalue}",file=sys.stderr)

tools/requirements.txt

@@ -1,3 +1,4 @@
 metar>=1.8.0
 boto3>=1.19.11
 control>=0.9.0
+rasterio>=1.2.10