new version

VIIRSdownloader.py

@@ -7,84 +7,100 @@ Description of data displayed:
 https://www.earthdata.nasa.gov/learn/find-data/near-real-time/firms/vnp14imgtdlnrt#ed-viirs-375m-attributes
 @author: faraway
 """
-import sys
 import argparse
+import getpass
 import re
 import os.path
 import pickle
 import datetime as dt
 import webbrowser
 import requests
-import pandas as pd
+import geopandas as gpd
 import folium
+from shapely.geometry import Point
+from tqdm import tqdm
+from sentinelsat.sentinel import SentinelAPI, geojson_to_wkt
 from geopy.geocoders import Nominatim
 from geopy.extra.rate_limiter import RateLimiter
+from branca.element import Element, JavascriptLink, CssLink
 
 
 def check_file_format(args_inputfile):
-    pattern = re.compile(r'local_VIIRS_data-\d{4}-\d{2}-\d{2}.csv')
+    pattern = re.compile(r'local_VIIRS_data-\d{4}-\d{2}-\d{2}.zip')
     if not re.match(pattern, args_inputfile):
         raise argparse.ArgumentTypeError
     return args_inputfile
 
 
-def get_address(data_frame, date):
+def get_address(geo_data_frame, date):
     address_file = 'adressdump-' + date
     if not os.path.exists(address_file):
         address_list = []
         geolocator = Nominatim(user_agent='blub browser')
         reverse = RateLimiter(geolocator.reverse, min_delay_seconds=5)
-        for r in data_frame.iterrows():
+        p_bar = tqdm(total=len(geo_data_frame))
-            location = (r[1]['latitude'], r[1]['longitude'])
+        for r in geo_data_frame.iterrows():
-            address_list.append(reverse(location, language='en'))
+            address_list.append(
+                reverse(r[1]['geometry'].coords[0], language='en'))
+            p_bar.update(1)
+        p_bar.close()
+        try:
             with open(address_file, 'wb') as file:
                 pickle.dump(address_list, file)
+        except Exception as e:
+            SystemExit(e)
     else:
+        try:
             with open(address_file, 'rb') as file:
                 address_list = pickle.load(file)
+        except Exception as e:
+            SystemExit(e)
     return address_list
 
 
-def render_map(in_ukraine, address_list, date):
+def render_map(args, in_ukraine, address_list, date):
     src1 = (r'https://gibs.earthdata.nasa.gov/wmts-webmerc/VIIRS_NOAA20_CorrectedReflectance_BandsM11-I2-I1/' +
             'default/{time}/GoogleMapsCompatible_Level9/{z}/{y}/{x}.jpg')
     src2 = (r'https://gibs.earthdata.nasa.gov/wmts-webmerc/VIIRS_NOAA20_CorrectedReflectance_TrueColor/' +
             'default/{time}/GoogleMapsCompatible_Level9/{z}/{y}/{x}.jpg')
     src3 = (r'https://gibs.earthdata.nasa.gov/wmts-webmerc/VIIRS_NOAA20_CorrectedReflectance_BandsM3-I3-M11/' +
             'default/{time}/GoogleMapsCompatible_Level9/{z}/{y}/{x}.jpg')
+    fire_map = folium.Map(tiles='OpenStreetMap', min_zoom=5,
+                          no_wrap=True, control_scale=(True), crs='EPSG3857')
+    active_fire = folium.FeatureGroup('Thermal anomalies detected')
     if not in_ukraine.empty:
-        fire_map = folium.Map(location=in_ukraine[['latitude', 'longitude']].mean(
-        ).to_list(), zoom_start=4, max_bounds=True, crs='EPSG3857')
         for r in in_ukraine.iterrows():
-            location = (r[1]['latitude'], r[1]['longitude'])
+            location = (r[1]['geometry'].coords[0])
             folium.Marker(location=location, tooltip=(
-                r[1].to_string())).add_to(fire_map)
+                r[1].to_string())).add_to(active_fire)
         if address_list:
             for address in address_list:
                 if address is not None:
                     location = (address.raw['lat'], address.raw['lon'])
-                    folium.CircleMarker(location=location,
+                    folium.CircleMarker(location=location, tooltip=(
-                                        tooltip=(
                         str(address.raw['display_name'])),
-                                        radius=10).add_to(fire_map)
+                        radius=10).add_to(active_fire)
-    else:
+    fire_map.fit_bounds([(44.184598, 22.137059), (52.3797464, 40.2275801)])
-        fire_map = folium.Map(location=[(44.184598 + 52.3797464)/2, (22.137059 + 40.2275801)/2],
+    ukraine_borders = gpd.read_file('data.zip')
-                              zoom_start=4, max_bounds=True, crs='EPSG3857')
+    folium.GeoJson(
+        data=ukraine_borders['geometry'], name='Borders of Ukraine').add_to(fire_map)
     folium.raster_layers.TileLayer(
         tiles=src1,
+        min_zoom=5,
         subdomains='abc',
         name='VIIRS CorrectedReflectance_BandsM11-I2-I1',
         attr='NASA VIIRS',
         overlay=True,
-        layer='VIIRS_SNPP_CorrectedReflectance_BandsM11-I2-I1',
+        layer='VIIRS_NOAA20_CorrectedReflectance_BandsM11-I2-I1',
         tileMatrixSet='GoogleMapsCompatible_Level9',
         time=str(dt.date.fromisoformat(date) - dt.timedelta(days=1)),
         tileSize=256,
     ).add_to(fire_map)
     folium.raster_layers.TileLayer(
         tiles=src2,
+        min_zoom=5,
         subdomains='abc',
-        name='VIIRS_NOAA20_CorrectedReflectance_TrueColor',
+        name='VIIRS CorrectedReflectance TrueColor',
         attr='NASA VIIRS',
         overlay=True,
         layer='VIIRS_NOAA20_CorrectedReflectance_TrueColor',
@@ -94,8 +110,9 @@ def render_map(in_ukraine, address_list, date):
     ).add_to(fire_map)
     folium.raster_layers.TileLayer(
         tiles=src3,
+        min_zoom=5,
         subdomains='abc',
-        name='VIIRS_NOAA20_CorrectedReflectance_BandsM3-I3-M11',
+        name='VIIRS CorrectedReflectance BandsM3-I3-M11',
         attr='NASA VIIRS',
         overlay=True,
         layer='VIIRS_NOAA20_CorrectedReflectance_BandsM3-I3-M11',
@@ -103,21 +120,35 @@ def render_map(in_ukraine, address_list, date):
         time=str(dt.date.fromisoformat(date) - dt.timedelta(days=1)),
         tileSize=256,
     ).add_to(fire_map)
-    folium.LayerControl().add_to(fire_map)
+    active_fire.add_to(fire_map)
+    if args.username:
+        sar_area = download_sars_data(args, date)
+        sar_area.add_to(fire_map)
+    folium.map.LayerControl().add_to(fire_map)
+    #fire_map = playground(fire_map)
     return fire_map
 
 
-def main(args):
+def filter_coords(data_frame, bounding):
-    bounding = [44.184598, 52.3797464, 22.137059, 40.2275801]
+    data_frame = data_frame.astype({'CONFIDENCE': 'string'})
+    in_ukraine = data_frame[data_frame['LATITUDE'].between(bounding[0], bounding[1])
+                            & data_frame['LONGITUDE'].between(bounding[2],
+                                                              bounding[3])
+                            & (data_frame['CONFIDENCE'].str.contains('nominal')
+                               | data_frame['CONFIDENCE'].str.contains('high'))]
+    return in_ukraine
+
+
+def load_inputfile(args):
     if args.inputfile:
         local_file = args.inputfile
         date = re.search(r'\d{4}-\d{2}-\d{2}', args.inputfile).group(0)
     else:
         date = str(dt.date.today())
-        local_file = 'local_VIIRS_data-' + date + '.csv'
+        local_file = 'local_VIIRS_data-' + date + '.zip'
     if not os.path.exists(local_file):
         try:
-            remote_url = r'https://firms.modaps.eosdis.nasa.gov/data/active_fire/noaa-20-viirs-c2/csv/J1_VIIRS_C2_Europe_24h.csv'
+            remote_url = r'https://firms.modaps.eosdis.nasa.gov/data/active_fire/noaa-20-viirs-c2/shapes/zips/J1_VIIRS_C2_Europe_24h.zip'
             data = requests.get(remote_url, allow_redirects=True)
             try:
                 with open(local_file, 'wb') as file:
@@ -127,22 +158,93 @@ def main(args):
         except requests.exceptions.RequestException as e:
             raise SystemExit(e)
     html_file = local_file.split('.')[0] + '.html'
-    data_frame = pd.read_csv(local_file)
+    try:
-    data_frame = data_frame.astype({'confidence': 'string'})
+        geo_data_frame = gpd.read_file(local_file)
-    in_ukraine = data_frame[data_frame['latitude'].between(bounding[0], bounding[1])
+    except Exception as e:
-                            & data_frame['longitude'].between(bounding[2],
+        SystemExit(e)
-                                                              bounding[3])
+    geo_data_frame['geometry'] = geo_data_frame['geometry'].apply(
-                            & (data_frame['confidence'].str.contains('nominal')
+        lambda row: Point(row.y, row.x))
-                               | data_frame['confidence'].str.contains('high'))]
+    return local_file, html_file, geo_data_frame, date
 
+
+def download_sars_data(args, date):
+    if not args.password:
+        args.password = getpass.getpass('Password: ')
+    try:
+        api = SentinelAPI(args.username, args.password)
+    except Exception as e:
+        SystemExit(e)
+    bounding_geojson = {
+        "type": "FeatureCollection",
+        "features": [
+            {
+                "type": "Feature",
+                "properties": {},
+                "geometry": {
+                    "type": "Polygon",
+                    "coordinates": [
+                        [
+                            [
+                                22.137059,
+                                44.184598
+                            ],
+                            [
+                                40.2275801,
+                                44.184598
+                            ],
+                            [
+                                40.2275801,
+                                52.3797464
+                            ],
+                            [
+                                22.137059,
+                                52.3797464
+                            ],
+                            [
+                                22.137059,
+                                44.184598
+                            ]
+                        ]
+                    ]
+                }
+            }
+        ]
+    }
+    footprint = geojson_to_wkt(bounding_geojson)
+    products = api.query(
+        area=footprint,
+        date=(dt.date.fromisoformat(date) - dt.timedelta(days=1),
+              dt.date.fromisoformat(date)),
+        platformname='Sentinel-1',
+        producttype='GRD')
+    geo_data_frame = api.to_geodataframe(products)
+    geo_data_frame['beginposition'] = geo_data_frame['beginposition'].astype(
+        str)
+    geo_data_frame['endposition'] = geo_data_frame['endposition'].astype(str)
+    geo_data_frame['ingestiondate'] = geo_data_frame['ingestiondate'].astype(
+        str)
+    return folium.GeoJson(geo_data_frame.to_json(), name=('Areas covered by SAR'))
+
+
+def playground(map):
+    figure = map.get_root()
+    element = Element('alert("hi there");')
+    figure.script.add_child(element)
+    return map
+
+
+def main(args):
+    bounding = [44.184598, 52.3797464, 22.137059, 40.2275801]
+    local_file, html_file, geo_data_frame, date = load_inputfile(args)
+    in_ukraine = filter_coords(geo_data_frame, bounding)
     if not args.noreversegeolocation:
         address_list = get_address(in_ukraine, date)
     else:
         address_list = []
-    fire_map = render_map(in_ukraine, address_list, date)
+    fire_map = render_map(args, in_ukraine, address_list, date)
     fire_map.save(html_file)
     webbrowser.open_new(html_file)
-    sys.exit(0)
+    SystemExit(0)
 
 
 if __name__ == '__main__':
@@ -156,5 +258,9 @@ if __name__ == '__main__':
         help='Specify the cvs file containing the satellite data to be displayed. INPUTFILE must match local_VIIRS_data-YYYY-MM-DD.csv')
     parser.add_argument('-nr', '--noreversegeolocation', action='store_true',
                         help='Disable reverse geolocation')
+    parser.add_argument('-u', '--username', type=str,
+                        help='Username for https://scihub.copernicus.eu/')
+    parser.add_argument('-p', '--password', type=str,
+                        help='Password for https://scihub.copernicus.eu/')
     args = parser.parse_args()
     main(args)