Temporal variability of climate time Series using MODIS satellite data

Introduction

The colab is avaiable here

Install planetary_computer and import other libraries

# Install libraries
!pip install planetary_computer
!pip install odc-stac
!pip install osmnx
!pip install rioxarray
!pip install zarr
!pip install geojson

# Import libraries
from urllib.parse import urlparse
from pystac.extensions.eo import EOExtension as eo

import matplotlib.pyplot as plt
import planetary_computer
import pystac
import pystac_client
import odc.stac
import rich.table
import osmnx as ox
import folium
import numpy as np
import rioxarray
import rasterio
import xarray as xr
import zarr
import os

Use pystac to import a stac catalog from microsoft planetary computer

catalog = pystac_client.Client.open(
    "https://planetarycomputer.microsoft.com/api/stac/v1",
    modifier=planetary_computer.sign_inplace,
)

Choose an area, time of interest and dataset from planetary computer

bbox = [3.771744,43.610725,4.335823,43.814463] # 738785.0528,6274273.2065,779106.8653,6298491.2641
year = "2021"
months = {
    "January": "01",
    "February": "02",
    "March": "03",
    "April": "04",
    "May": "05",
    "June": "06",
    "July": "07",
    "August": "08",
    "September": "09",
    "October": "10",
    "november": "11",
    "December": "12",
}
items = dict()

# Fetch the collection of interest and print available items
for name, number in months.items():
    datetime = f"{year}-{number}"
    search = catalog.search(
        collections=["modis-11A1-061"],
        bbox=bbox,
        datetime=datetime,
    )
    items[name] = search.item_collection()[0]

print(items)
print(len(items))
{'January': <Item id=MYD11A1.A2021031.h18v04.061.2021041154946>, 'February': <Item id=MYD11A1.A2021059.h18v04.061.2021062101022>, 'March': <Item id=MYD11A1.A2021090.h18v04.061.2021092202346>, 'April': <Item id=MYD11A1.A2021120.h18v04.061.2021122015755>, 'May': <Item id=MYD11A1.A2021151.h18v04.061.2021152213240>, 'June': <Item id=MYD11A1.A2021181.h18v04.061.2021183192128>, 'July': <Item id=MYD11A1.A2021212.h18v04.061.2021213191344>, 'August': <Item id=MYD11A1.A2021243.h18v04.061.2021244183922>, 'September': <Item id=MYD11A1.A2021273.h18v04.061.2021277221831>, 'October': <Item id=MYD11A1.A2021304.h18v04.061.2021306213324>, 'november': <Item id=MYD11A1.A2021334.h18v04.061.2021335171453>, 'December': <Item id=MYD11A1.A2021365.h18v04.061.2022004015033>}
12

Choose the right asset

Each item has several available assets, including the original HDF file and a Cloud-optimized GeoTIFF of each subdataset. In our case, our interest is the ` LST_Day_1km` asset, that corresponds to the temperature mesured during daytime

rich_table = rich.table.Table("Key", "Title")
for key, asset in items["March"].assets.items():
    rich_table.add_row(key, asset.title)
rich_table

┏━━━━━━━━━━━━━━━━━━┳━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━┓
┃ Key               Title                                                               ┃
┡━━━━━━━━━━━━━━━━━━╇━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━━┩
│ hdf              │ Source data containing all bands                                    │
│ QC_Day           │ Quality control for daytime LST and emissivity                      │
│ Emis_31          │ Band 31 emissivity                                                  │
│ Emis_32          │ Band 32 emissivity                                                  │
│ QC_Night         │ Quality control for nighttime LST and emissivity                    │
│ metadata         │ Federal Geographic Data Committee (FGDC) Metadata                   │
│ LST_Day_1km      │ Daily daytime 1km grid Land-surface Temperature                     │
│ Clear_day_cov    │ Day clear-sky coverage                                              │
│ Day_view_angl    │ View zenith angle of daytime Landsurface Temperature                │
│ Day_view_time    │ (local solar) Time of daytime Land-surface Temperature observation  │
│ LST_Night_1km    │ Daily nighttime 1km grid Land-surface Temperature                   │
│ Clear_night_cov  │ Night clear-sky coverage                                            │
│ Night_view_angl  │ View zenith angle of nighttime Landsurface Temperature              │
│ Night_view_time  │ (local solar) Time of nighttime Landsurface Temperature observation │
│ tilejson         │ TileJSON with default rendering                                     │
│ rendered_preview │ Rendered preview                                                    │
└──────────────────┴─────────────────────────────────────────────────────────────────────┘

The MODIS coordinate reference system is a sinusoidal grid, which means that views in a naïve XY raster look skewed. For visualization purposes, we reproject to a spherical Mercator projection for intuitive, north-up visualization.

data = odc.stac.load(
    items.values(),
    crs="EPSG:2154",
    bands="LST_Day_1km",
    resolution=500,
    bbox=bbox,
)
# Get data data of July as a sample
raster = items["July"].assets["LST_Day_1km"].extra_fields["raster:bands"]

# Use raster scale to redefine temperature in celsius
data = (data["LST_Day_1km"] * raster[0]["scale"]) - 273.15
data
<xarray.DataArray 'LST_Day_1km' (time: 12, y: 48, x: 92)>
array([[[-273.15, -273.15, -273.15, ...,   16.55,   16.71,   16.71],
        [-273.15, -273.15, -273.15, ...,   16.55,   16.71,   16.71],
        [-273.15, -273.15, -273.15, ...,   16.73,   16.69,   16.69],
        ...,
        [  15.51,   15.21,   15.21, ...,   11.71,   11.35,   11.35],
        [  15.47,   15.55,   15.55, ...,   11.71,   11.35,   11.35],
        [  15.47,   15.55,   15.55, ...,   13.47,   14.21,   14.21]],

       [[-273.15, -273.15, -273.15, ..., -273.15, -273.15, -273.15],
        [-273.15, -273.15, -273.15, ..., -273.15, -273.15, -273.15],
        [-273.15, -273.15, -273.15, ..., -273.15,   20.53,   20.53],
        ...,
        [  16.67,   17.83,   17.83, ...,   13.59,   14.23,   14.23],
        [-273.15,   17.83,   17.83, ...,   13.59,   14.23,   14.23],
        [-273.15,   17.83,   17.83, ...,   15.63,   15.97,   15.97]],

       [[  23.83,   23.27,   23.27, ...,   31.33,   31.65,   31.65],
        [  24.17,   23.93,   23.93, ...,   31.33,   31.65,   31.65],
        [  24.17,   23.93,   23.93, ...,   31.13,   30.55,   30.55],
        ...,
...
        ...,
        [-273.15, -273.15, -273.15, ..., -273.15, -273.15, -273.15],
        [-273.15, -273.15, -273.15, ..., -273.15, -273.15, -273.15],
        [-273.15, -273.15, -273.15, ..., -273.15, -273.15, -273.15]],

       [[   8.65,    8.65,    8.65, ...,    9.33,    9.35,    9.35],
        [   8.63,    8.63,    8.63, ...,    9.33,    9.35,    9.35],
        [   8.63,    8.63,    8.63, ...,    9.55,    9.47,    9.47],
        ...,
        [  10.17,   10.23,   10.23, ...,    6.01,    6.09,    6.09],
        [  10.21,   10.17,   10.17, ...,    6.01,    6.09,    6.09],
        [  10.21,   10.17,   10.17, ...,    6.91,    7.67,    7.67]],

       [[  16.45,   16.07,   16.07, ..., -273.15, -273.15, -273.15],
        [  16.59,   16.17,   16.17, ..., -273.15, -273.15, -273.15],
        [  16.59,   16.17,   16.17, ..., -273.15, -273.15, -273.15],
        ...,
        [  16.65,   16.99,   16.99, ..., -273.15, -273.15, -273.15],
        [  16.63,   16.79,   16.79, ..., -273.15, -273.15, -273.15],
        [  16.63,   16.79,   16.79, ..., -273.15, -273.15, -273.15]]])
Coordinates:
  * y            (y) float64 6.303e+06 6.302e+06 ... 6.28e+06 6.279e+06
  * x            (x) float64 7.622e+05 7.628e+05 ... 8.072e+05 8.078e+05
    spatial_ref  int32 2154
  * time         (time) datetime64[ns] 2021-01-31 2021-02-28 ... 2021-12-31

Interpolate data and fill missing values

There are many nodata values. Lets interpolate missing values using rioxarray. But me must first define our nodata value. In this case, we can define as -273.15, since Modis defines nodata value as 0

data_with_nodata = data.rio.write_nodata(-273.15, inplace=True)
print(f'No data value is {data_with_nodata.rio.nodata}')
No data value is -273.15

We plot the month of mars

data_with_nodata[2].plot()
<matplotlib.collections.QuadMesh at 0x798254418730>

png

With the nodata defined, we can interpolate na’s using the nearest neighbour method and fill missing data

# Interpolate missing values
data_interpolated=data_with_nodata.rio.interpolate_na(method='nearest')

# Plot march again
data_interpolated[2].plot()
<matplotlib.collections.QuadMesh at 0x798254a2e380>

png

Plot all the months together

g = data_interpolated.plot.imshow(cmap="Spectral_r", col="time", vmin=-10, vmax=60, size=4,col_wrap=4)
datetimes = data.time.to_pandas().dt.strftime("%B")

for ax, datetime in zip(g.axes.flat, datetimes):
    ax.set_title(datetime)
<ipython-input-39-d1854156f2f9>:4: DeprecationWarning: self.axes is deprecated since 2022.11 in order to align with matplotlibs plt.subplots, use self.axs instead.
  for ax, datetime in zip(g.axes.flat, datetimes):

png

Plot images with folium

Create bbox with Fiona

import fiona
import rasterio
import rasterio.mask
import json
import geojson
from fiona.crs import from_epsg


from shapely.geometry import mapping
import json

import requests
response = requests.get('http://app3.ouicodedata.com/collections/agropolis-zone/items')

with open('./response.json', 'w') as outfile:
  json.dump(response.json(), outfile)

# Load JSON data from a file
with open('./response.json','r') as infile:
    with open('./agropoilis_zone.geojson', 'w') as outfile:
      data = json.load(infile)
      geojson.dump(data,outfile)

with fiona.open("./agropoilis_zone.geojson", "r",crs=from_epsg(2154)) as geojson:
    shapes = [feature["geometry"] for feature in geojson]
    print(from_epsg(2154))


EPSG:2154
import os
from rasterio import windows
import random
from rasterio import features
from rasterio import warp
from rasterio.windows import get_data_window
import rasterio.mask

for i,d in enumerate(data_interpolated,start=1):
    if not os.path.exists('temp'): os.makedirs('temp')
    temp_file_name = '_'.join(['./temp/modis_tile',list(months.values())[i-1],'.tiff'])
    d.rio.to_raster(temp_file_name)
    with rasterio.open(temp_file_name) as src:
        out_image, out_transform = rasterio.mask.mask(src, shapes, crop=True)
        out_meta = src.meta
        out_meta.update({"driver": "GTiff",
                 "height": out_image.shape[1],
                 "width": out_image.shape[2],
                 "transform": out_transform})
        if not os.path.exists('cropped'): os.makedirs('cropped')
        croped_file_name = '_'.join(['./cropped/cropped_modis_tile',list(months.values())[i-1],'.tiff'])
        with rasterio.open(croped_file_name, "w", **out_meta) as dest:
            dest.write(out_image)

Export images and plot then together with Folium

from pyproj import Transformer
f = folium.Figure(width=700, height=500)
m = folium.Map(location=(43.7237837,4.0003508)).add_to(f)
from rasterio.windows import Window


for i,d in enumerate(data,start=1):
    temp_file_name = '_'.join(['./temp/modis_tile',list(months.values())[i-1],'.tiff'])
    dst_crs = 'EPSG:4326'
    with rasterio.open(temp_file_name) as src:

        img = src.read()

        src_crs = src.crs['init'].upper()
        min_lon, min_lat, max_lon, max_lat = src.bounds

    ## Conversion from UTM to WGS84 CRS
    bounds_orig = [[min_lat, min_lon], [max_lat, max_lon]]

    bounds_fin = []

    for item in bounds_orig:
        #converting to lat/lon
        lat = item[0]
        lon = item[1]

        proj = Transformer.from_crs(int(src_crs.split(":")[1]), int(dst_crs.split(":")[1]), always_xy=True)

        lon_n, lat_n = proj.transform(lon, lat)

        bounds_fin.append([lat_n, lon_n])

    # Finding the centre latitude & longitude
    centre_lon = bounds_fin[0][1] + (bounds_fin[1][1] - bounds_fin[0][1])/2
    centre_lat = bounds_fin[0][0] + (bounds_fin[1][0] - bounds_fin[0][0])/2


    layer = folium.raster_layers.ImageOverlay(
        image = img.transpose(1, 2, 0),
        name=list(months.keys())[i-1],    # add a name to the layer
        bounds = bounds_fin,
        opacity=0.7,
        interactive=False,
        cross_origin=False,
        zindex=1,
    )
    # Overlay the image
    m.add_child(layer)

folium.LayerControl().add_to(m)
m

<iframe srcdoc=”<!DOCTYPE html> <html> <head>

&lt;meta http-equiv=&quot;content-type&quot; content=&quot;text/html; charset=UTF-8&quot; /&gt;

    &lt;script&gt;
        L_NO_TOUCH = false;
        L_DISABLE_3D = false;
    &lt;/script&gt;

&lt;style&gt;html, body {width: 100%;height: 100%;margin: 0;padding: 0;}&lt;/style&gt;
&lt;style&gt;#map {position:absolute;top:0;bottom:0;right:0;left:0;}&lt;/style&gt;
&lt;script src=&quot;https://cdn.jsdelivr.net/npm/leaflet@1.9.3/dist/leaflet.js&quot;&gt;&lt;/script&gt;
&lt;script src=&quot;https://code.jquery.com/jquery-1.12.4.min.js&quot;&gt;&lt;/script&gt;
&lt;script src=&quot;https://cdn.jsdelivr.net/npm/bootstrap@5.2.2/dist/js/bootstrap.bundle.min.js&quot;&gt;&lt;/script&gt;
&lt;script src=&quot;https://cdnjs.cloudflare.com/ajax/libs/Leaflet.awesome-markers/2.0.2/leaflet.awesome-markers.js&quot;&gt;&lt;/script&gt;
&lt;link rel=&quot;stylesheet&quot; href=&quot;https://cdn.jsdelivr.net/npm/leaflet@1.9.3/dist/leaflet.css&quot;/&gt;
&lt;link rel=&quot;stylesheet&quot; href=&quot;https://cdn.jsdelivr.net/npm/bootstrap@5.2.2/dist/css/bootstrap.min.css&quot;/&gt;
&lt;link rel=&quot;stylesheet&quot; href=&quot;https://netdna.bootstrapcdn.com/bootstrap/3.0.0/css/bootstrap.min.css&quot;/&gt;
&lt;link rel=&quot;stylesheet&quot; href=&quot;https://cdn.jsdelivr.net/npm/@fortawesome/fontawesome-free@6.2.0/css/all.min.css&quot;/&gt;
&lt;link rel=&quot;stylesheet&quot; href=&quot;https://cdnjs.cloudflare.com/ajax/libs/Leaflet.awesome-markers/2.0.2/leaflet.awesome-markers.css&quot;/&gt;
&lt;link rel=&quot;stylesheet&quot; href=&quot;https://cdn.jsdelivr.net/gh/python-visualization/folium/folium/templates/leaflet.awesome.rotate.min.css&quot;/&gt;

        &lt;meta name=&quot;viewport&quot; content=&quot;width=device-width,
            initial-scale=1.0, maximum-scale=1.0, user-scalable=no&quot; /&gt;
        &lt;style&gt;
            #map_c411caddeadbebf70f7cebfd52330cd5 {
                position: relative;
                width: 100.0%;
                height: 100.0%;
                left: 0.0%;
                top: 0.0%;
            }
            .leaflet-container { font-size: 1rem; }
        &lt;/style&gt;


            &lt;style&gt;
                .leaflet-image-layer {
                    /* old android/safari*/
                    image-rendering: -webkit-optimize-contrast;
                    image-rendering: crisp-edges; /* safari */
                    image-rendering: pixelated; /* chrome */
                    image-rendering: -moz-crisp-edges; /* firefox */
                    image-rendering: -o-crisp-edges; /* opera */
                    -ms-interpolation-mode: nearest-neighbor; /* ie */
                }
            &lt;/style&gt;

</head> <body>

        &lt;div class=&quot;folium-map&quot; id=&quot;map_c411caddeadbebf70f7cebfd52330cd5&quot; &gt;&lt;/div&gt;

</body> <script>

        var map_c411caddeadbebf70f7cebfd52330cd5 = L.map(
            &quot;map_c411caddeadbebf70f7cebfd52330cd5&quot;,
            {
                center: [43.7237837, 4.0003508],
                crs: L.CRS.EPSG3857,
                zoom: 10,
                zoomControl: true,
                preferCanvas: false,
            }
        );





        var tile_layer_8ebb7bdf73674872bc467e44a609f99a = L.tileLayer(
            &quot;https://{s}.tile.openstreetmap.org/{z}/{x}/{y}.png&quot;,
            {&quot;attribution&quot;: &quot;Data by \u0026copy; \u003ca target=\&quot;_blank\&quot; href=\&quot;http://openstreetmap.org\&quot;\u003eOpenStreetMap\u003c/a\u003e, under \u003ca target=\&quot;_blank\&quot; href=\&quot;http://www.openstreetmap.org/copyright\&quot;\u003eODbL\u003c/a\u003e.&quot;, &quot;detectRetina&quot;: false, &quot;maxNativeZoom&quot;: 18, &quot;maxZoom&quot;: 18, &quot;minZoom&quot;: 0, &quot;noWrap&quot;: false, &quot;opacity&quot;: 1, &quot;subdomains&quot;: &quot;abc&quot;, &quot;tms&quot;: false}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var image_overlay_06ec3d8b1ee2472e0dee5d2551cc2776 = L.imageOverlay(
            &quot;data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFwAAAAwCAYAAACc2R15AAAJmUlEQVR42uXbBY8USxSG4f7DBIIEJ0hwd0Jwd3d3d3d3CAGChzA330neztnaap3ZAXJv8iULd5nuerq66lR1TzJ58uRGUYYNG9aYMWOGZd26dZYtW7bkZsmSJZkZM2ZMmgkTJjSdlStX5ubIkSOW2P/bu3dvh6xfv96ydu3aNCNGjCidBQsWWG7cuGH5/Pmz5du3b5YkD3r8+PGGXQW8LHQrwfPQwW4GvCx6z549LW0BrwrdDvSq2HXB+/btawH8x48faYDOBR84cGCawYMHWyZNmmTZvHmz5dy5c5bz589b6kB3FfiaNWssYG/durUT9vLlyyuB56FXAp87d24mdjPgZbFbCQ7awYMHLceOHbOE4MKuCs4dHgPv3r17h3z48CHN+/fvO8TASYgdguvEycmTJy379++3LF682FIF+k+Ag10X3KOH0JXBY+iAe+ws8DrYrQLnPHwuX75sOXHiRFPgHlsZN26cZeTIkYXgjx8/tjx9+tTSCTxEj2Eru3btsujkxo4da/lT4DFswIVdFTyvd4fg5PDhw5YLFy5Yrl27lmIXgvfp0yfNokWLugScO4e0Evzo0aMWzS2UZ/v27bNs2rQpzenTp9Mw/OiiTJs2rUPqgvtcvXrVUgpcCSsZykPh8ztloHVyIXZZ8BUrVkQTw24nONiVwbt169YB24OH1UwIzkSlxJBJDPvWrVuWssB5AXvHjh2W69evW4BljFZOnTqVhqFn+/btjSlTpnSISkwPjkEMnOoNYP5MqZpQwgm7CDxEP3DggMVjE4+ch+3BlTrIPkAXgW/YsKETuLA9OEgkCxwHhc/juJTNVHG54IMGDWosXLjQwkGXLl1qUeOYaHbv3p1m4sSJtbG7Ap3Fz5UrVyycs8AVkHXuHjzE9ujCloEA6ZClwXXQ+fPndwIXdlVwYcfAq2A3C65xXMDNgCvTp09P44c52h+Cg86ccenSJQvH37lzpyVhLKZsoqgfMmSIZdasWRYmHW7JZcuWpatKTj4GXhW7GXQmTjWQnw8dOmRhTKXH8ffbtm2z8BnNgCtdAi7sMuB1sZsBp4Ex8Nu3b6cQtGfPnj2dPmf06NGZARvgefPmWfjzxo0bLVw47jTdVVooJshT5Gt4UdiIETYzraJbj4vEvgqfsXr1asvQoUMtIfSZM2csrQanpBMgDaVaYUwVdhlwAP8YuCLkZsGrYJcFB5qoUYAI2/+/sN5mElNi4Fnw+uzZs2d3ypw5cyzhEOMvrFa+SbilyRjN3zMhsndCT1KAphRk9clnhdg3b94sDV6EHmKH4B63CFwJoZVRo0ZZevXqlSaG3WXg/tYtAw402K0A9w8VfFh4sIABWA1VuNWZo3iwImw9QCkLnoXOqhR4TFQa4qbJOhk+fHhDYXLkxBguwNVwwskDytMeej+1KOOnoM+ePWupAu7LsyLoZsH9U6s64OE2wF8Dfu/ePUu4x+DjocnUqVMLsf2kzvEAp+FULYzFanP4mBBoPVpU+vfvbwnBPXoWuN9bYRGkDpeCs1xlWUt9KuDjx49baCC9jovDJMoFoZfpQj18+NCSBx6D9uCxxLDbCc7YHwMPN7OErba/efOmOXB6eliX6yJQjt2/f99y9+5di9+/yINWVMYVwSt+D4eGAs85c1tr1czEyq4hxwCNYYdqSwmhfTw4x6cMVehsBg40twI9AHC/rPcrOaUZ8CJsLcM54SroXQVOG2PgoMewX758mebTp0+dwf1+MQdj1ckQQjXCrchKE3C/18xkyVUugs4Dz0NnHQEwS3g6Dm3yFRlVFu1gXO7du7eF9utuLgJXwueXpAM4m1OrVq2ycGLa5/5T4Gw0xcCz4FsFzsJP4L4c5bFcFrbG6SxwXgb6/v17IwEv3NsWOLfYgAEDLOwHZy1jNd7TYMDBo+Tjdg43+YnmihD8wYMHheDcobSDyomhRaDgcA7h2gNsUhU8C/2vBWdyjoEXoeeBq1LQ4g08j9oKcLCz0D9+/GjRi0EJdTa7akwuQmNMY2+8X79+FhYErEa5JQXOv7948aKF7UpOht8NeyjYVcA9urD9gxDWAjxjzAJXmgHX3lEIrnOtBE4vbRe4x64KTsLanjmDrQU+uy44G3NZCcGVt2/fWhhSfv782UiAplxTb2ARwcFYCLB7yKYUCwUdkMkJAODDZTsAulBMZF0NHn5+CAt6FngR9j8LzjlwfO4OVqtKWXCqFOaSOuBloGPod+7csbx7987Ci52/f/9uJOGSWOBs0jOh+jEL8HCjh9KQMHnxWTzioteyV6PwnLFV4GB7cI8e29ple5dUxa4FzgSTB64NLsbwLGwPzvBAWHFqxcrvMtGB9OzZMwu3ogfPQvdvx9IO3xa/pVCEXYTuj8U+ktrKu4zPnz+3fP361fLr1y+L/kvCyUUJeyWNYkcRcN5AUqqAs0VQF9y/H8hXPML3vz12EXgMOws8PA7YMXB6NuDWw6uCa2zz0GE8eIidBc5zP+YTwLklHz16lP57j50H7o/psauAh+isyn10LHZJmad4czYKzi9p7GN4YcxjvKWmZTIpA+4bnFdd+Jd2qoL36NEjBfdvHlQBz8MGfObMmZlR/Q84x6NKef36teXLly8WXQQDZ3Ly71Mr1NKcaNZLjTH0ImjCatWvWHnFl00fnTznwCt3wg7BPbpWmexUcqdSGJTFVtQBWwquHlwXnJPy2DSGiqBZcE6+CjjL+rrg/oGG2p8HrgAebo+EbdXdm2gSoowB3D+hV7gAeScajouA63M4IKtMj+xLOCoLJscnT56kWNTsrAUA540o1gHajw7B9TkKHSv89kYMGux/Epy0E1zQrPTY06BA4Dg8MMjDVlRI5IGzNcDGHNvUTKyU1LoDE30NgtqYA3CFmABZceo1tqzbsAw4e+7+bVQPztcyGLv1M/s6nIt/vzsEZ8igZ9cFD7FJHvY/A67qRz0YbMD52YOrh8fAw9KThx7Av3jxwsKdw9Ci3w1Xp1nYMXB/Lgovf/IVHB7gqLJjMk9UoDNJsTQFnAmJ/Wo1jvc5FLCF6PfCFUpJ/RxeQB/AFRYMwmGGD5/YhOjtBAc9hG4LOH9Hr2XiqALO7wBO/f3q1atccI8eW6azJcHnstYAmnWG368hDDcxbJIHzqsVCtva/k22RA2kB1CC8ayPW0ffkOCbyDzRzgLn31IB+K9jCAFsLSjKgIMegitZm1D/O3BVH/qZlSP4QLBUplyiIToXxm6tLn2qgLMSZY+GGp1JTQ0HPfzWGRWZf5knTBa4x1Z0DsxfDHf/AT3t4WDmuiQpAAAAAElFTkSuQmCC&quot;,
            [[43.60734015888762, 3.767730367838], [43.81775157192726, 4.3423021816889245]],
            {&quot;crossOrigin&quot;: false, &quot;interactive&quot;: false, &quot;opacity&quot;: 0.7, &quot;zindex&quot;: 1}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var image_overlay_a89a77ffe85cc350b6d5e694476e9cb5 = L.imageOverlay(
            &quot;data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFwAAAAwCAYAAACc2R15AAAE6ElEQVR42u3b10vsQBTH8f3/3/VFsT0oCnaxd8XesPeuiBUbKiJ7OQd+4Tg7k53JJpvZvQpf8Hr3pnySTGJukqmpqcmm1fr6uraFhQUvm5+f52ZnZ7nh4WFuaGiIm5iY4La3t7N7e3u/+v7+5jI+gvuKHhUc2KmCh2H7Ck/YLS0tQc3NzVx3dzeHDVDS4D6hS+wwcEqHXTLgtMf4AF5dXa2to6MjJ6CXLPjU1JSX2CbwgYGBHOzUwF2wAb6zs8PoacCHYevAOzs7GZyik2vJgxcTPh82RcBqEpxaWlriSho8afgo0N6Bn5+fxw6eBHwh2GUDboMeB3yh2FRfXx+HZSpr8ELg48D2BpywiwmuIhQL2ztwV3QJXii6aSPEie0FuMT2CVwWF/YfuGOFYqcOvru7mxp4FPSenh6ObrmWDbgLugqe9F4uwefm5rg4wGtra5MHJ+y0wV3RdeBbW1scLYvNNAhXV6YY2OUGLjs7O/uVCTpxcImdBHiS6IWA50NPFdwW/Q/cEft/Aw9DzxQD2xdwW/SSATdhJwWeFHoc4Cb0TDGww8Bt0P/AHbF9A8+HHhe4Dj0yeHt7e1Ca4Ens5V6AS2BX7Hzg+dCTAA9Djwv88PCQ04JjJv39/b8yQVOnp6f/Nfj4+Hh2cnIyFDsHvLe3N0tJ8DBkie0L+MHBQazoXV1dHFkAfHl5mcN/ChM2wKmTkxOOTCQ2dXFxERSA4wlQusOVD/vm5iZ28DB0G/Co6CZwHOEu4DApa3Bgx42ugs/MzHDT09OcCk7zBvDR0dGvtOCUzVDy8PDgNTieQ4zrikVixwpuix0FnJ4itQE3obuCu6KHgUtsEzhOqLqd6/r6mgvAV1dXs9Tg4GA2bG8Htiv48fFx8Kx0uYED2gQObCtwiS6xXcAJW4JHRTeBq9gSfG1trSB0dV1U8I2NjZwwbwmNsG4B+OjoaM41OD4UBRzYKrgNug04xscwcGQLbjM8qr/4kIXa5eUlh2tyuW5GcPkhV3CJTdGh5IKOZUISDwGcPp8PHGGZdeC2J38VXIduAudffOiaUj0kdAi24Co2wF3QbbEBLtFN2Pv7+8Eyy1ywTeBYhqurK059NALvAEUGv7u7C7YiVlwHrYKr6EBSf1GQex7GzsXFxSAdOABpui7gLpe2RQP//PzkooB/fHz8ApfoJnCJXgg41ouw4wDHEKGCj42NceptkoaGBq6uri4oI8+4mDC21MvLCwdw+h7o9/f3HD6L5FhGXwROn5fgQFVPeDr4fODqOLmyshJkC047hi22t+D4ArhEjwI+MjLC0ckdgBSOFroOltgSHFc2JmyKzhc22Ki1tTUI723izmtFRQVXWVnJ0TOK9fX1XIbGGCw8Zg40Ffz19TXniuXr64t7f3/n3t7eOB24ek2KO2zqMARsDFvABrhEjxNch65ecdBnab4lC44VonAEyJ+ZwAEEVBqzcZLC+zVh4LRB5UbG9OjeiTxiTeAUhji8fQzwqqoqDvfAGxsbs21tbVxGHp63t7ccnRSpx8dHjqCfn5+5p6cnDhuDQCm8h4gN8PPzw9H39O8pTA8bQIWNGxwvNeHn6nBlApfoJQFOezodEfg7FRx7OuZn2gA6cISVxdAhwXHSUt8gcwHHn3F00LkO09vc3ORwtGBZmpqaOFwG0o6BEyquuv4Bbij46rTk9coAAAAASUVORK5CYII=&quot;,
            [[43.60734015888762, 3.767730367838], [43.81775157192726, 4.3423021816889245]],
            {&quot;crossOrigin&quot;: false, &quot;interactive&quot;: false, &quot;opacity&quot;: 0.7, &quot;zindex&quot;: 1}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var image_overlay_9be1b876fc220c2b22bc75a44a970ca5 = L.imageOverlay(
            &quot;data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFwAAAAwCAYAAACc2R15AAAJCElEQVR42u2cZ29USRBF3///EUiIDyCMENE2OWeTc44mJwM2iKBZ3dKe2evr92BmbWu9mkG6kmde6O7T3VXV1T00hw4d6u3bt6+0c+fO0q5du0p79uxZor1795ZmZmZKt27dKt25c6d07969/neXLl1apCtXrpQuX75cOnfuXOns2bO9kydPlk6cOFG6ePHiIl29enWRVA569uxZ6c2bN4v08uXL0uPHj0vXr1/v1+XMmTMl2vHgwYPS8+fPS0+fPi09efKk//fdu3dLN2/eLGWdeeeFCxd6c3NzpZ8/f5b41xw9erR3+PDh0vT0dIkOkPgO0OogSZAkGn3//v2SKgRwFSydP3++BDwazTuOHTvWf+/BgwdLCV6wJDrWGw6I2dnZEtCAyDPqZMpE1JF73717V6LT9L4x8BUETt2oE6Dfvn1bevjw4RLQzLLTp0+Xjh8/XqKuN27c6Hfc9+/fSwsLC6WGXpJ48MiRI/2GY2a2b99eYgpSURrjgOlpTAb3olOnTpVUDqZq//79JerCvUDDLLi4huh0OlwNlwQJc0bdBgEu2MMCp8yRAA6YtPcSQHO28UyO7NevX/cePXpUSp8DaNoBN12jk5hlX79+LWFiGj2YwGm8BOjJyckSI4NKeKP4jl6+fft26dq1ayWuU55XlsoDho6U05JwgC9evOh/R+PGwFcQOKAB8urVqyWA8SN8xt5jRgSWejND0yfovRJ+QKJjc6ZS13SaYomD1cCQAE8kVU5TTkvS30xxAGNa6FUKZ4TIhCT8dFbYPqRKY3YQjaGzqKgDB8owwKnnagOHz5oErhGcjWA0JnBNbcnjagDzPpwkn5lR1E8CeFdcj/kQaOBTds5U3gVwPY9fATTP0nlcbxTq4Ri3bNnS27p1a4nvdu/eXWLq0znuGGkUoyfDP2BS4THwv2EvB7jH1Tl9gUrjPLTK1SffU0Eq7qYEAIBG+U6J7zAhHipKgP748WNJIdznz59LdDr3Ykqw2QBXJ6QJxIRglrjebNu2rbdjx47+cp4FDs5S1yWFcBKF+mjMhczvgGcsOwjwtN2rBZyYeRjgwM4ox6MqH/G/BT41NdVfaQKcKelRB0AzXMrlsyDRwHRwXWB8IcOURukbgNyWt2FQUH+eoTMFhVEOLOpCedkuDQg6ZWDgcoha1uMcCQcFe9SB52zkXZQv4AB+//79IrHS5DpqyGEAmigF0AcOHCgpgpG8sZgMACPAc50KC1quFtMsZAd4riSFaQEMn8m1SGnW0qHTsRKjcc0BpyLDAAfoMMAzOTUIcE9uJXCHhnQPoBEOj1VjW5IME8Gs+PDhQ6scfJkUgJIdzLw3DoIR4jYNUwHgzKEgNxVAIrTzXIwnnPQ3DcQMAQDlEt+X+SSyMvwbA18mcM9XS1zLWD2jFspR2WlKcIR8TyyNYxRIQGf+HJOSwOUnKixkRUm8DXCcJBDdWWH/clGRy3QP9XiWhuZyHBDc52YGiIAHLiDGwA14ghwEuDsyiXuw+0x5z1fkQiNnTkY4mmWYHUADCZhofn6+9O3bt/7fAMa0pHiX4Hs+pQF05k5wNoxSNyV4eu4BeOZb3HZmg9MBesw7Br4CwN02+6ZzRhp6P9M/6wDMXC5LTHvAp0nJvItvGAPnx48fJVacwCWn/evXr/49mVFMs+TZRzpW7W1YwmNKAA40XzTkHmWaEq5zf5spGVXgOP1lAc+Vpe8P+sLFG5rRCB3StjuDuAenSSPVmHR0mB2A8ywmS2VzDQdH7oQddqACWvBlTtykYDrSweK8VU9GtkdZS4ADMdOZEkAZeaMO3KMR3+WX3Dcw2tURS4ADk9HsSuCYEtKUdJLnMnLBkTkOpA7lfXQg5WJ+WNQwimQ6aEwejUiT4llCOp1nP336VAIanYZjdPg40gz/BgaeTnMM/J9MHzE0RxykvId38BkzorIzkyiT07CibItOcg8PcU/uYme87jG7b8n57pArox2+z8NGGn2Yl0GPTuh5nsGxAolR2wbcYQ8KnPqOLHDf+vKNXSKNzKUA1c8p5mYw8kiJM4m5XvD4vCGHwhaa79bwN9fy+BkwPW+SKUw/quByU0U5uc/pCSbJl+2Ek5lvobH5eU0DzxG3loFj392uSzTaDxJ1wfOjEal8H8KEZJTUtsvjm9INm8JtUzudWJoUvz8XPYRheSjHnSnXclcGR5fH3DztOgzwPEy0WsC7ttXWLPCMbNI8uInI5FdmFBnZDrQrG5imRvd0rVyRZw7JmXTJ8+NN7jv6kj13qRO8nzzKLSxCrzx56lt0ubwfOeCZI1kJ4AnGnR3X8vgD92I+MvJw05L35CaFH6fwvLmbHzoHEyCHmZ2RUVDmYSQObKLMOkoN09mX6gBmHxOzk8A5HidlhvB3wPNM4GoBz/MrgwAnQlkJ4CyolgU8j+vmee623/w4cE333EmhMb4sdoieN8+pndGIA89ceK5UM/vov57o2gXiXt6hzgHmly9fSuRbWoEDU/JjxFKmXdvOR48qcD478BzhpAa4Lpu/CHgCzvPeyHMgGbH4XqED98Z2pTZpXDpIgc73/Qm4/5oOM4eJwZSkiXEzkyFjdtK/Bu4mY7WBU+lhgDOyhwHOtUGBe3STTjeTYXz2qCfjbZ7NM5di3bCKlDAVvq0mZSU8vs6fY+QvGDz/QQO7FhGZ+VOjKDOB5ywAuMrKuuWuU5uDHQNfJnBPUeQKNlep/neX78kz6O5rsAgkAjnqjUh9F3A/dJmhHVB5sReSFckVZptNTDvLu7g3HZQvqSkHWLnEz5+WjIEPCTxzGDQ2Z5bvr/pslfKHABnp+JZf1++F6BzMkd5D5+av2DZv3lzasGFDadOmTaWNGzeWCrjnnVlaUxiA2pJHea4kR4inKP/0Q9Y8UeXiGiPr/wqcrcz/HHibo3IYeWQ4TwywONMZd/+5o5Q/8spy2kxI11lHB44TZg1ClLd+/frSunXrShMTEyX9hKcP3J1jeuL00K4cEblY8VHbBbwrMhgDHxI4I74tB9L23rYjb7m5QeM488iWIGfXHTgnEDxC8qyjm4qsU55Tp4NlMnCODt//KxEsBHXz//YE8/MX7G+51/YmauYAAAAASUVORK5CYII=&quot;,
            [[43.60734015888762, 3.767730367838], [43.81775157192726, 4.3423021816889245]],
            {&quot;crossOrigin&quot;: false, &quot;interactive&quot;: false, &quot;opacity&quot;: 0.7, &quot;zindex&quot;: 1}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var image_overlay_d1707f042f199bc480911a410925d0ed = L.imageOverlay(
            &quot;data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFwAAAAwCAYAAACc2R15AAABtUlEQVR42u3ZzYqDQAzA8bz/M/UderIexFIQxC+wFylU6C4OzKJbK3ZM4mTM4X+SufyETGDgdDr9SK4sS7TqunYqy7LVgWK7Y38DLR58T2wXaAV3wN4CLRqcGxsDWiw4JzYmtEhwF9iqqmbjhhYHjon9CZwS+hDgc6OjaRrT/X6fxIEtAhx7Vs+Bc2F7D459Ob5er7c4sb0Gx9pGLpfLJAUnxO77/g38er2a0jQ1HR58C+5cCo6M/QnaFkXRpDiOTUmSmG6321+HAaeAHhq2kTH2+XxeBOdAB2nYz+dzscfjYbLr3xjbtgROjQ4SsLuuW13btiZ71s7sMfjQEjglOviK/Q2yC7jtEzgVOvgE7YpMBU6BDntiYwBvBedGB25sCuSt4JzoLODUyBjgXOgQCjQHOAY6hAKNAc6BDiFhY4BTo0Mo0NzgrugQEjYWOCU6hAI9B57nucm++Nhv34aJDtKhLepSW8HX/BAycN9x9wBfyhlcCq5v4P8DH6AxcUWDS8QVCS4ZVxR4CLhiwEPB9R48NFxR4CHgeg8eKrKC71hRFCYFV3AFV3AFV3AFV3AFV3BPwH8B7rUxKEpomD0AAAAASUVORK5CYII=&quot;,
            [[43.60734015888762, 3.767730367838], [43.81775157192726, 4.3423021816889245]],
            {&quot;crossOrigin&quot;: false, &quot;interactive&quot;: false, &quot;opacity&quot;: 0.7, &quot;zindex&quot;: 1}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var image_overlay_047cf7a4389ddf495e7be954fc800a97 = L.imageOverlay(
            &quot;data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFwAAAAwCAYAAACc2R15AAAH10lEQVR42u2baU8UTRSF+///FxNRIwoouIviBuK+4K6gImBCQoK8uZU8k9OHut3AyBuQ+XASZuju6X5u1d2qurlz587OrVu3im7evNnSjRs3dq5du1Z09erVlubm5qq6e/fuQNevX68qruvif/x23Fdodna26N69e0X8zv3793cWFhaK5ufni/h87ty5oqdPn+7Sw4cPi+L8EJ9fvXpV9Pbt25Y+fvw40Js3b4qePHlS9Pjx4yJ+d3Fxsej58+c729vbRX/+/Cnic1N7eKBeuXJloJmZmZZ4cEAAWQ3n505PTxddvny5pYsXLw506dKlIgyNuCYGCGM8ePCgCGgIA7ghHj16NDiH+0ac6yxGwP8ScI6N70LA49mBquJ6XAPxPcfFzNjc3Cza2toq4nNTcyUBZ2pqqiWH5FM9AKiBQgBGgAdqXHdycrJoYmKiiN8DNIMAt6dGzoADUcX/Dhs4ruhEAAee+3s9h2Px6c+ePSsCKu5jaWlp8D9cBWDx/4gZFK4Et/Pt27eilZWVIj6nwAGBAA08RkIWVFXuHriWXh/Q7lK4J4BrUOaBfTSPgA8B/Pbt26mACQgHrsYh6+EzkICKNEiSuXiQdOAY5PXr1wPgnz9/Lvry5UsRnxtN/TTA+VQHkgdPFa6kD3gc22dAjH9cgAfsIwtc/Tvn+rHk5biSmoDnWYnGFYUewkjAe/HiRRFwP3z4MAD+8uXLItwP4DGsBk1Ao/fv31fVxEPy8DGlAXzhwoWWmPIO2QOjjtYR8A7gAD0o8K6R7xCBrnIjeeVZE67Fixjg6ncABhqg3717t0tkRvhmr2i5FoaIYBmGCuE6EKBxP03AjYLj/PnzVY2PjxcBXKF0Qc6A6yg/TOD+3WEBD9jaAnDgq6urLXUCj//hy4FH/h3TMqskMxej2Uw2Q/rckc4UdztahWalu6eQ7i5C/Db35v0iwNd6Lr3APUh6sDypwNUfe8+EcyL7+Pr1a0tejO0CDjSt/moVJsfplAaAG8uvwbVrVahXsFnmo/Lejw6CTG4AfCrwNHB6ADyWwLOyvwt4DfYwwDXQoXADwMINOPDw0yFy9TAO+bXLZxAVZlw/fmvgUnhgwGddu5oyI6lbyYD779V6L26ArNyvdQwZlZ5ZjIAfELiObOAqYKS9alX0TAAONAe+vLy8Sxz76dOnln79+tUS1x7k4TyoT2kebi85tBtrL8Azt6HXys4lHXRX8s8Bz3Lovw3cR7hnJYDW6tFBA5jMAxdCzyNSOSCsra0Vra+vt7SxsVHE/3/8+DEQXUBAYwiOBXgc06o0PZPISmyVp25ZWriXvot3G2vnj4AfAHjWus0yka6spNZJzBYjPEjiJhTS79+/WwK0K9K6nz9/tuQuBOAo3Em4rVAYufECp7ZynwHvK3y6euQnDXhkNhEzhgLubkLz6WyBIpsFajjPzTMXpZmSB1DvoajLoRwnDcSFOCyVFzAOWl1IiF6NqgD3lM7h1lZzMr980oDXwOOrdflOZ9gu4Exbna7Zph1XzWc78K68G8BZa6ArRvh9Y4Da6jwjjn63uwkfxQ67a6T3AvebP+nANRPRpTt6L74jzKtU+jFZTdD4NgRuWNfs+K5r9SWkYLOq1N1FrX/T18+p9dF5Dl+7JD3UUcZ3QKCIUdgAx1D7Aa6wTyTwCFYBgurSl+RU2o+JWADYTKR8KDIgn1Xa12l8Q4/vJtJ9GN6B4+GQVqPZxs1s65wWQ5nr0MDt18Xo3CvPAeQjDdxv4igCBzC/zT1RWRKw8NMBgyqzBhdlYD2ecR9cH0VsqGVGVL+NLxup3/LVacBjnK51Rg9egNH8PlurzPot6pr2A5zRd9jACcbHBjjHcG5tX3nWC+dcHlJ7KL5Vgf91Aef6/DZG94CIYbUD6UG5FTRrPjvbaep7q2s+NDuHz7rvz7c1jID/BeA+a7jRmHpMcQD4wyjsbLubzxLfYBkP6cABkwHv6vVkKZ8uQmfAS9DUqO7B0Xct+a7UAOzA/UZqwPF1ewXetb+wCzj3sB/gof0C933kvttWGewbuKY5Kh7OV0BUrO+56OD5FuIu0Blw1UGAZ30jql0Hrq+aZBuHWsBr++5837W+sjEM8GiF/p/Aa+rLq92lOHBPAkJZkuEeIc4rwLN3VTyp100ztQ00uqHdt+vqeqAD95eWkAfRWjDtextuGOBdxxwYuAa6owi8Nqr7gHeNcN2h1SV3p17kKXQ3Ap81BaboKyPcXyxiStOYAQwQwz1oG1LFLlJA6zW80YOr4RwHHqPDC46+1af4uwu29ou8VhgBHwK4jvQMtDbe+oDXagVt2oV49dG3evN97EJuNCh6K5Plou/fvxepOwAeW3wxAO4B4KR+9IlV/A8BPoTRvQYAYtdyn/vxvQB36McaeA12F3B950aBxz0ChBGv7w0hQHv25Rs5UUDzYiV+S9NFDYCA5h5osp09e7ZobGys6MyZM0XxGjp/N1GAsMbn++R8yclXREIeADESCiMB0oOwG4vj4m+Ae857GMA1h/bG0z8J3MGHqEa9UtVWcLYbLFte87eSURdwf/s4DIQh9TXL0OnTp4tOnTrVEgYINRGceGAgAbMGOltwdeGG4noYEJi11zxGwIcAHt+rv2eEuw9nZtVGPD1mb4YBvPZmXAa8TzXg3guqvcvvBQ7VKEGytlf+P1IB94zuRhRRAAAAAElFTkSuQmCC&quot;,
            [[43.60734015888762, 3.767730367838], [43.81775157192726, 4.3423021816889245]],
            {&quot;crossOrigin&quot;: false, &quot;interactive&quot;: false, &quot;opacity&quot;: 0.7, &quot;zindex&quot;: 1}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var image_overlay_35a70477cdb1aec56ce558f46abbb54b = L.imageOverlay(
            &quot;data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFwAAAAwCAYAAACc2R15AAAIBElEQVR42u2baU8UTRSF+///CmL4aNQQIRoRRXEXV8R93xF3UNzgzb3JMzlzrGp6XkAxMyYnGafpquqnqs69VdXTzM/Pb5w9ezZ15syZPp06dWpjbm6uT6dPn05duHAhFfer4jv+Ju5XUW7pOm24ePFin6iH6zdu3EhFW6anp1OHDx9O3b17t0+zs7OpmZmZnviOtpw7dy5FPZcuXUpRj+vmzZsb+/fvb9Xjx483vn//nvr582efmsuXL/8GAsVDOTQHXgLEQwAJ0Wk8tIpr3Iv4/uTJk6mFhYVUPDwQa23hXgWOeI54/tCVK1f6pPWMgG8DcGYHgK9fv55aXFxMBdjQkydPep9RCbDr48ePffrw4UMqgQPE4SpgV+0hQw7NO+Do0aOpsIFDhw6lsAXK984CyO3bt3u6c+dOChuotePYsWM9wEeOHEl1AR6wR8BbgBNH3CZi9jpoxD0OtU3Up1LYr1+/7unZs2cpOq+JBgFEAyYP7oHUg00JOJ1REw87MTFR1fHjx1NYCWVrUAQAYP8m8Pv37/fpnwV+7dq1FN6t/oqwgatXr6a8A6KdDhqdOHEi1QWwZ2QqB/7o0aOUW01aigcoDVIA9imO2hpR6gw65ODBg32amppKUS8PqcA9oHUBzrPtJHCHveuAxzUA491MeerDszUruXXrVqoGHNC0TWcoYB24lxWx4d69eymg0ZGlNUco/taBP3jwIAXonqVEIyKKayRXETxrgbDkZ12m4NADr025rQCvjXBGeYh6yFyYYYAGsoJmEcLf0DnYj87K8+fPp3zFWrKBUEB5+vRpn+hQyvKVeVyrAX/+/HnqxYsXqSYCUwD3hQGjHuBto3U3AgfKIMAZhYMAj+/Dlinj8+fPqbdv36ZevnyZImtpBR6jn8J96RtTDLC1bIQGlsQKk6nu+y10SIClTqZ4TVqvrxNIFek8bIMUMz4z/bnmVlXK1ii3M3B/cNewAvdnVQ9X4MvLy6k3b970CcDv37/vU8MUruXbAKJSFQ3wzMUzm5L3A8RB83/apdBqqRaivrgHe8GGWJUigGtwY5SWBtdQAPdR2gU495SAx0h2v9aARoDzncJS/h9t5b42sZcSnRDxoCHfrsEBeGkbE0up2cZmK84QluL77iHfI2dUPnz4sKgInFiGj2isY7uAh4YSuO6tlMBiGe73ZCAx7bEFdvSYMbV8n/+/evWqp1+/fqUImr5bGHaSK032LNT/SioB57tauqj3uB1xzYNmabeQDnT/9ROePw08AO8ocA8mXYB74OkC3GcIM8ltQsWIr/k8FsJUjykOYIdDSoc9ARqobVpbW0t5maGwvB5wT8dcDluB1xRBSxcuKkZiKb0bAd8G4D79fVdNBTy8G+A8eCm1q2UuBELSNexDLQQbWFlZKWoQ4JxdApn0kJn9m6WUwG4mtwvdUnXQ/yLwEP++fPmS+ivAKcR39jSlAgryYFZL8ULMFOrhaK20D1IDjpUAOazi3bt3KRYjQP306VOK60tLSyn9B3Cs49u3byk6ryRdHDW+hGd0qkqLHtWwAg/Qq6urKfZMVJTfClyh1XbjXLV3OEK6V1ESwErXKJ/6Sce807TjAE3ax5kilhIA3EoAo4e/mi4ymnVEo4GB+9J+2IGH7379+jWlL/DwYg/gGemURScpcGYXbUzgvrKsLWNDZBz6qkLI0yeVL60dqgZSv0ZDa4E35KAJkroKDMXDMhp9tAKPQDcCvkXgATG8GUiABiIi01hfX09p1sG9+Lzn2Ph9AKYtZF08Rx6x+amNvkPnoGtninoeOEj6t9nuH+BLCx3K99MZRheLF4DsauC1kb2bgDO7GNF+sgIYfFpTuR8/fhRVs48QdlQDTozQ9jMIaGMCJyChEvDaWaK+nANYT+18s4ggVhINQ4wMtzBVF+CeO3cBrrC7AvcBs+uBYwMIDy7l6Mw2D44AJ5cupXKMYEC7tcTfOHACaQ141EWne5s0iDa1l87bgHsnBXRA+CilMl14sPjgs18bAd8CcIKJ576q2pIY8OzwacrlD0U92Aerx/gMLGyhFjz5v8NWUS65uw4a2kAbsRYGR+6llBY6m8mBhyi8BmKrwH0V1wU4nwcB3ga7BpwZ6sBL7+MMDFwzBM0aGIV+Sg08Tb1IuzzlAogHKH0oFxbCw1Jv3F/b+aNcXwiV5G3TNvnIhgGsOgPXNM/BDjtwTwMjEfA0F36lk7NGAQOPjMJPoPle96wRf1OzBwIUsEvAa2AIgCV5ByM/T9RNJO8cLCcshXYC3O1vy8A1KO40cIXdFXgb7Brw0gHuZsD912aAp9MBrK9U+AKNmc+2BpDZIIyf1jS6N8LqkEJpGJW5V6m4RgNdGrxcTFcfTXEfUDYDrdc83SSI0UZsiOcaAd8G4NEmX1jRRn+uqMeBI9pImktZMdNZCQOY7WTOFngHPhTAx8bGNprSqtAXDyqHjxhVtfNBTcdq8rQtPnu9Xj+jtzQIdhPwgP3HgDvILsD5rPXqaQz5t7+mpnLgfpAc5eLzBEm3FL6nHt0NhBsvLwH3wIEDqX379qX27NmzMT4+nmp0f4OG+ql2F9HwNrg1S0GltMw3+EfAtwi89MKM/nKXGeR7zINIAStk9sk1J6/FE8qIDTpfj3B4A2h++sjvl+Lz3r17U41Gax8BXQC72oB7IK2VMQI+IHD1bZXvkbjUp90y/g9onqeU/fh+CvK9En3BnvSvdtTI75T01x38VmpycjL1H+2Zxlqr6oJTAAAAAElFTkSuQmCC&quot;,
            [[43.60734015888762, 3.767730367838], [43.81775157192726, 4.3423021816889245]],
            {&quot;crossOrigin&quot;: false, &quot;interactive&quot;: false, &quot;opacity&quot;: 0.7, &quot;zindex&quot;: 1}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var image_overlay_aa3f48d500ceccff387b5b78900cb83f = L.imageOverlay(
            &quot;data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFwAAAAwCAYAAACc2R15AAAJNUlEQVR42u2c127cShBE+f/fYMNPhh+c4JxzzjnnnHMO8BrVwFmU6pIrUpINXckPBYi75ITTM909M1w1V65cGZ05c6Z05MiR0sGDB0sHDhwYbd++vbRz587Srl27SlyfOHGidPny5bHOnj1bOn78eGnPnj0lyufzkydPjsU9lM916sKFC53y8rp0+PDh0c2bN0fPnz+fVrTTRX/R3r17xzp06FDp9OnTo1u3bpUo6+XLl6VGwLs6LW3ZsqW0devWEvfQADp76dKlkoAfPXq0tH///pIDE/A2EGnQPwUb4H2gtwGfBH3eAaeMLtiTgM817LmEDmj6q/qePHlSevfuXen169el5uLFi2OXQqP00ObNm0vr168vbdy4sbRv377SuXPnSlevXi3JcLgUd0nS7t27SwCn7Dbg3LNjx47Sn4Q9V9AdOHXOS+AoYTtw2jZT4PRpEvDpoD9+/LgTukRskluhXc+ePSvhSiirAZzkjaTDAMe10Hj5KYlnFSSxeAZFPu8CLjlsCbc0CfhMYQvAEOgCPgn6ggLOLKF+dO3atRmPbAHQrOwLfTrgSAxu375dwqVQxtOnT0uNg8ZKEukg4Pn81KlTpWPHjpXkuzKAyr1IGMNTp77AaVMbcMGeBNxhdwEfAh3gun8S7HkJfPny5aVJwB2eYG/atGksYHcB7wsb4H2gO3CJmMQsPH/+/OjevXulV69eld6+fVvKMhoFN1Kvbdu2jTZs2FBas2ZNiY4SLKmkLVpjDKzNvUBuUwJPgA5bUm7bF3gbSAfeF/ofAS7YswHu/laZj0T2wSxx0Hfu3JmiNuAsHlIAJ47Qjj4ZSRfwNuhApD4GElkY0iwG6Pv370sAJ3hiiEYBce3atWPQiJEHLFwMHfJsIQNcF3ApQbumA90XuO5pg602JvAXL150Qu8DXLDv3r07BvvmzZuSypWoZ5yHTwLuoJlGQLlx48Z/8m5WiYCX3OdOgi0xSpA64mI1i6E1lSVcmhtsCPCETnsJgCxsGFgYWn2nbb2BC4r7Yfx5juyZAFe5Q4ALZkJ38LMBDqyEDhhlEA5c9RL8SRCohxmsNtHGzLcRrmTsUnJ5ClhAA5HGAlmdoMO5onTgEh2RkdwALjUOoAn84cOHJQwDCEZaGwiHLSP0AY4Lod55A5yODAEu2F3A79+/P25gG/QhwJkNDjuBA5kpz34HUKgXA2Tb5WbIq3NFybXXISkn13MNO1xkIUAj4AGPzjkoLI8bImjRMb5Xh0kZXQKNfFQAHXj4UjrO98DkWkbhGYzSBlx1LDrgDjuBA70PcJXP94LDyAY49/JsrvwEOFO5dAMZaBX88h6eJWjms7isJpfUgGZlCDxSIj/ZwQiAZmpz7Rvz3AvEBK4R5cDdZfFsCoMuWOAOuy9wPwXpCzyzEsCyK+fuQfIAmdOfa/w/AORC+Bs4gP/06VMJVwPEz58/j0Waxz1fv34t/fz5s/T9+/eSysGtPHr0aNSwkmRjKE90GGVKw4jOKH1pnokq2OKyKI+OtwGX5Ov+AZ8D4HwH2AcPHpQITID3bAQ3Rt6t59Vo5O6B8tLf5xqA+jX9Afzx48cSkBBwgSl4gMYlffjwofTjx4/Sr1+/St++fSvpXu9zw8oS4ADKs0oJAIjDAUZgBs25Bi79X4Grj+rbrIB3bYnqM9+SlZjyLIIoEyBaHZJvE/AwHNu8169fL9FZz9XdYBLgiTfsw6jj7FXjShIshsBdfPnyZfwdn1GGvpMwEtcqlzYxqAq49lEkTnZIC3ElgpfZASN7sQPXLOFZQOPvNZuok9y8DiASOHkynRRgPgM01+mzyUjkXsjfAc6UBjwgErwEaAzIPRnc5EIISN4pCQPkBpiA45J4JoMlKR5lyy0AFODcMxj4qlWrRi4HzoheTMAFk2vg6XOAAhrhw2kbscRnos/IZt26dSNEoJGUheShMNfAzDeQ3Cj4eVwTMyb30GUc3FiucjEgriQDrpRwgAhUXIyeyy3TXIYzer1MRv0Q4On6WoFz00IErs4qm8nVYO6DABw3ob8BmvVwLwYhY/IsC2Mzu/RZw8kOnfJTlVxa535IvrwpkLkaxR0AMV8SEnSH7ate3BFwASToAM1NJK5zCT+vgeNL/xZwlNsMk4AzE30fReJeZgOG8T3rTO1wB8BMw0i5n4KLSVej73zVLAFebBtcwKRDgnQLfM6zbSf3fu7nwMn3ue4DHNACmUv2PsDzkKAPcIfdFzjfzRvg/jJPX+AedNT4tlWmlJtidJYZIWE4Vpa5Jez7I7mD6CdDEq6EsvRMmxvza7WnYcqTrulvQHou7gse9lL4XlAZjbkIktzF4N/8jdkEzohACxY4f/8N4G1K0Ck1HndAJ/gOF8MqFVfjq1GmfAbAtqDJPXlvGgkjYmAJ4+R+Ue2HA5kXMn2PBOUbVVxrVHJSxNlnvjAkySCzBU4HhwDPw+g+wMlQhgB32A6cWUcbp4zwvsD9fXLfIWQPJYUfZ1YocMwEOJ0AHoDzEITVr+rK18+6AiHfK/0jCOIyMADXPINrk/vKLGoicOD55hRBkkCa++CLGTjtIFZoJhH/MEJuN9AOzfCGLMPl755I5Oh+suKHAqgNOKtFz1pmA5yFGW3ygeDA8eESIBiJ1AcIdym5lPe98DkBrlE6F8AZ9fl6whDgk6A77C7ggqk2MTP7AvfRnruArDSJGZTpx4T5enbXYboyskZQ0q0A0/eq/fUxOivxbAIHgp/WzGaUO2zfS0k5cO7P89PcxliwwN39DAHeF7YD95OqNABtpD+ChhsgK8GVkGfni52KH8wQQAOYH6CtXLlyinTuUL/TzFMYXEounwU/07t8aTLfB8lfOaemA+6w24yR9frh9qIF3gW7Dfok4G2iHsUepm7+uJbYRBtxDwJMhkJ+z+Evh8J8zrMCSx9hQT2AXbJkyRStWLFi/AOHxk9FcrRhQdJESRX4ygljJHRBmA50G/ghsP8BnyVwCdgZG7rkv4Jgq5nOsepF9JMVoTIXYhH1Aphgma9FaA2QP7PE0ABfunRpadmyZSW1hU27RoByOpM/AlMjCeAOuwu45IFsCPC+oBc8cH/NYDrgbdlDX+C8VDRUDnz16tWlLuC+oQTw3CvJF5R8dZubXfkrC97tkfIfH/wGz3awyBSztDMAAAAASUVORK5CYII=&quot;,
            [[43.60734015888762, 3.767730367838], [43.81775157192726, 4.3423021816889245]],
            {&quot;crossOrigin&quot;: false, &quot;interactive&quot;: false, &quot;opacity&quot;: 0.7, &quot;zindex&quot;: 1}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var image_overlay_6610ffc3843fc286b5f922d0ff999b64 = L.imageOverlay(
            &quot;data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFwAAAAwCAYAAACc2R15AAAIcElEQVR42u2cSU+UTRSF3///G1wYTFgYiIkrF0ZBREUcGJzAWVDAEUVAROwv9yZP5/SxqoevNSHpXpyku3mHqqeq7lSlzfPnz1uPHj1K3blzJ7W0tNTW4uJi6tatW6nbt2+n+P3evXuptbW1th48eJBaXl7uEL+jlZWV1P3799vye7jm8ePHqSdPnrTFb7T/xYsXqdevX6devnyZWl9fT2nbeC7fuYZnbGxspN6+fdva3t5O8VswC3GP8grFe7j23bt3qQ8fPqSaN2/etP/48OHDFDcGyPn5+Q4BHiDcw8tDgOAankcnV1dXU/xdB3RhYSHFAHAvUBnUAO4DyPufPXuW4rsOEu2lDQiIAEJbW1tj4H8TOM/lWgDFzGZ2w+fp06cp+uGgdTXu7u6mvn37lvr69WuqiSXHiPHyAH337t3UjRs3UgDne2lZ+3KlIYClkzw7TBOgEaaKZ/EeQEeHMRUA0EHQa1G0jzZ5W7oB95ndD3B4jARwvvtK6wace1CseFY99tz9Sm3l6irHhHz58iX1+fPnVBOgfQnGkr5582YKwIglX2oAL2bZMoAuzEU4aOD7M9wRAiJsKmLm0e4x8L8AXJ1XKD4DmpkOGFYq/WGVBFgHznfuffXqVYoIJ7S5uZnyaIjnApzBi/sxQ+/fv08Bug1cZ0SIMA8Tcv369RS/o1IYR0N0ZqnUcfnz3PEBU4HzeRDgwPjXwLn/VAIPOQg6A7RSPIx6gUbxPFYbbeA7wpRoGIc5+PTpUwoz4z6KNkd/aAv34CT/MClxMfDCZFy9ejU1NzeXmp2dTTl4Be7hni7lEJ0bAxfg2Oj/C9wTCZWbi+ggAHxQuAaIKDoMaB8UVHKWunpVmADM1M7Ozh/Cb2BaPH9g8sRAYUL29/dTBwcHKcJCwDcBWcE66MuXL6cYCLW5jHA32A6cGTUIcGbYaQQObI9KXO04vBtwZnWIGgqQNf6lAdRiiDzc1GgW6LG/L32tf/g9Hoc76HgObXPQ3It5UgeJw0N+jSdU0SZmdt/AA1DMWMI9wkGFPYrAY2bzXk++dOV5VAJ4j06w7Q0mAuAedwMRmw0gDZN4hlcSAa4NphPAUkeqUtMCFBybF5GQhnaYIcyND5oPVgyCJlWnCrjHpf0A12XYL3A6Nghwri0BV3vPdQH248ePKU2uQv4M7omVwbXYewcOYBS/xSpqHHAto9RlWUsAuMcrfuostaCksbmbgICHs3JpiBjSvzEYwHJoY+BDAtdSql4DLJ/x/B2oCsWjE6/Z8J64hnswGQyamxJ9TrSj8YwSc+COr1Tv1rq5ymvXWhzyQlMtaYkZVQPtmwRj4ALcq3H9ACeS0LpGyGNnNQE+GDXg6uR82QMN8G4eNMzzevfe3l6K74iIh3Y2xN2A93o0cAMeIBkM7uEZ1LR1Q8BjV/f02Ecgq8bAhwCupkKrcy7dK/TYGfl+ocL37NDNAVAVtMPz+gd/V1NBIvP9+/eUp/Ionq+xfeOgCQM1EQH+zMxMhzwpwsSUqoOjDpz2DAUcmL6Zqw4TAQuwJbOBmfE43zNMri+B7uZIySaBjskADvBKAjrChPz8+bNDXB+DwszW9zfXrl1rhXB4zFI1KYDlGlL1MfC99szmOlZNrCqPiNKkOHDfJtPYuOYAvf7tpdiQp9buNOOzPw/5fiMd0TMjLhyf7xLFIGAWMB0Ap7TqCqB8BvTQwDElCtyTkVEADsCafOa73Sfj1GMWHU6TzJLIQ52Vn/fQiMXPfWACvNrnCY1nnPp8F9cyKBr60QnMAmaiJnV+LhxgL9j9AtdDRCMJHNgeQTj4kjk5PDzsEL/7AR/NSj0n0Oy3IaTDpChsOl6rWTsgNT81E+IhX1zjAwhonGQpafGkpARZdaqB+8w7jcBrVToHrjsu7vCA5dmimgxMiWeS3Mtg8X7NDTzbzVoK8HTJO0hAuzNTU+LhXy1NLzWiW+m0VAIdFDhw/jVwL6qdeuBeFQSq1i8ctFftfC9RnaIDr9U/FHQtHOSZuqOjzjFEP7Qu33jKXUo8/DxGSQwKD66VUPU7nz0jGxngtUxvGOB+9EDNgpsDr0fTAV36fHb57ribiRLg2rPimpqz9HAQExb98HbzHSaRYTdal/YzI+4sS4c3NZ0PMUtL5zyGAd4LUgl4rxldgz0ocPrh7fYTYjGxBwautRM9voYNxy47aKBqHcJ3tL0eTUe0Kse9HtJ5PUSdn8Nyh8d3DQN//PjRIUxMt90hTEjtsGkC73bo0v+hkRerRhG4bzgHWC/IwctPJgSvRmeyx8FefNETSLU4m9+9Tg08hebLGCAMjp7xqJmBmq0uhXncozvpeoYkgAL4+Pg49evXr9TR0VFqaOD671z+NXCfpf0AB8YgwGtxdQ14zPJeoHkW/SmdafRTDF7Ui92wRmsjXtlzB8gs1di5lkH6AXU6FnIH50uf63RXpleSoim3m4ya6dK0fwx8SOClA5Ve6dOB4XkAR/zuTlL7jsmgPMJJ40uXLqWmp6dT58+fbzVas/YddEyIxpXElG5K/HwfHtuLSLofWCv8A0T3FX2QfFAAEu2snWAdA+8C3Ddy43OsEtpS2kqj/Z50MSheAWTvUUU0ginhHt2mgw9JHyeOz549mzpz5kzq3LlzqampqTb0RvcZ/R8T+T6h7pzUCk5+7q+UyHgqXdpTrDm4MfAhgZd2XqLjtWNntI02szq1uMZq85o5gwlMdYY+GF5Lof/xTkwI/zrkypUrKaBOTEykJicnUxcuXGhdvHgx1cQDgOigNXOqVQFLO+l+KspnI863ttU1Bj4gcI3dFbTv9fkuTQyCVw6JmPw/JfAMV4ETlrmjBebJyUnr9+/fKZyjg/YaTQw+faVtXkfibE8ADun5eszOf0iadF2eUGWUAAAAAElFTkSuQmCC&quot;,
            [[43.60734015888762, 3.767730367838], [43.81775157192726, 4.3423021816889245]],
            {&quot;crossOrigin&quot;: false, &quot;interactive&quot;: false, &quot;opacity&quot;: 0.7, &quot;zindex&quot;: 1}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var image_overlay_ba0f0d14f5a326f6d0321ceb4226f921 = L.imageOverlay(
            &quot;data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAFwAAAAwCAYAAACc2R15AAAKXklEQVR42u2bZY8bSxBF5w/nSxQpicLMzMzMzMzMzMycSIkCfrolnVG5dsae8Xr3ORDpSll7oPt0961qcHLkyJHKrl27TLt37zYdOHDApO/2799v2rt3r2nPnj2mQ4cOmY4ePWo6ceJEqm3btpk2btxoWrVqlWnhwoWm6dOnm0aMGGEaNGhQpXv37lXq1atXlfr27WsaOnSoafjw4ZXZs2eb5s+fb5o3b56Jz2fOnGmaO3euaenSpalWr15tWrt2bZXWrFljWrdunWnz5s2VLVu2mNavX2/iGbyX98yZM8eke2AIk6tXr5qS48ePV9CxY8dM/C2Q+/btM+3cudOUB5p7Dh48WNm+fbtp06ZNpsWLF5tmzZplArjUu3dvU8+ePaukhvDAaRzuUwPGigMWAW3r1q0mdQKgAXzlypVV4vMNGzaYVI8/GvjgwYNNAO4o4CtWrDAtWbLERBkBD/AdO3akIzaO1GnTppmmTp1aJV0Lt7Nnz5ouXbpkSk6dOpV+eOHCBdP58+fTFsJKsJ3Dhw+bdJ+XYMd7AL9o0SLTuHHjTLIDJNh9+vRJFYGPGjXKRCUFiR5Hb4wWQk8DtGyS8nNvEeCCXRY417YscKlZwHkfoJHKTc9esGCBidGG8H/ACzj30AjUg/dNnDixSvqMUQUvmBAfk3PnzlUuX77cRhEewlKwEAVWSY0gO5G4l5dTyYEDB5oGDBhQpWHDhqWaMmVKlUaPHm2iRwqAPFICCPoHvABwgOYBJ2tgqAog2Q/f0eMpI6ABJhsAZLQB3jNp0qRUuo77fKygUbhu7NixJkaj4sayZctMsfEpY3L69OnKlStXTAKtvyXAAhGwBEsPGnEPwAECcKV/Uv/+/atgS9jNjBkzTPiwB07hywAHWhngvqGKAidQtyxwyVdy8uTJJoY4ebG3BYIYhadSzA1iRqX7ly9fbgIWWQ6gx48fbyJWSJQXYSETJkxI5YGrTjQkHYXGggEZW6JJjQeG10QriZZCQ3jgfMbDeSlQhwwZYvrdgXtLpBOVBg7IRoGrt3Mvw5dh26NHD5OfSQo43idRGSzFB0mJPNpbC38DOA5jZTDYASAQMGlwH7gjcKCiLOAxCDM3oLFgkqinqFecPHnSxOQHkECk93CdbOfMmTMm7KUMcA+7HvAIuAhwUsaOBO4/i8D95M6rJnDl13g6+bZydEmeD3B8HRugUqRLZCdSHnD8kHsouNIxgiSgSQsRw5XZocqQNwulbID2VtIs4HmwDbhSO4Fm8gNEADcbONBbETjQGwEO9CzIPldPyLsBqbxcYip67do10/Xr1023bt0y3bx5M/VzAgQpXr9+/Ux+McoDj7AlCkTwpLAKOADHOrA3LAzwlENTdQDzPIKjJODeQqLqAZfVEYR9aivVgt0wcMEuC9xDz4MdgRPhywBnbaQW8Fqw84ATU1Ae8JEjR7YRdRszZkwlYdkQAZiefO/ePdPdu3dN+gwLIU1jmAKwW7dupq5du7ZRBO5h+5U+BV1A0rDYGlbG+2kYLVwRqAFNMPbAy0CPoBsBLtCo04EjRkVUWeDybHJtZTfkw+0FLuXBbhdwAN+4ccN0+/Zt0/37902PHj0y8ffFixdTEPgUK32Alrp06VJTecCRgiNpXlyjwVIIklnA/aQqC3g96LVgdxpwwW4GcDKHPNgxC2Ek8Tf3e6hZayKIbbksNQqcLIj3MqIAH7OfFDgBMFrIs2fPTE+ePDHpGirOS/1aAqoHHFgo7imqx5YFLu9uFLik+8rAbhngPvUrAlvSypoHrokYeTWVYvOBQEbleL+AM9oagR3lN5vzxGwaa2HeQJmkmsAfPHhgevHihenDhw+mV69epXbDego+iUcxEYgTnHqwAS753fHOAJ4FutnAI3TVs13AAcC2WJwgoKLAEcGQSpA2+tw4blCTf0cAZWEXBc7uD+VhDYjOQDkk6mjAgYmVAPzNmzepmPQwGaJyRYHnwW4UeNaJgKLA68GWVK4ywP2imAddCPj79+9Nnz9/TqVeLpG5sP5MpQHOYZ14aEe+nwccYLIStuTYfQc0IPxZEaCzfMu9xALiDJkP9xSBXRQ475TygHvYBpxgSTrYLOC+VxFss8C3EvBYtjLAqQfWwoy3DfCHDx9WJLKRr1+/mn79+lX59OmT6fXr16Y7d+6YWFkEIgFL9hGHsIctCQZr1QgrYX07HsYhmPIMLR0z20XskscNZJ6t+7GdIrDrAfewI3B/TuWvBK7Pve9LtUDXgx5he3nYZFc0hgGPWcm3b99MP3/+rHz8+NEULQXwHPQEjKbQPoBo3UNwInCJAuozYMU9S/72h0slrc+zUcK9gIyNyPceOGWoBztLtWA3DFyw2wuchSa26poNnJ7Mu1nCZZSw/sIGit9EKQpc9lC018f0UOI91EPrPwlQv3//bhJoejlB8+XLlyZsBxti7VyVoEcBIg84YDTcyetZjOIZXMP32IZWCIHfkcA9tKJWE2H/VsDxX2yDzRGdewRePBrBKAE0J8MUa3w5JTYuCKKsNrLiWAQ4doGUCTEy887yiEPy48ePiuT/ffnyxfT27VvT06dPTQB//PixiRze74VyIBQwBFh/fkXCjvxh//8beIRdC7h6dEsC51QX4DlPzikBKZ4nxxYAzhKynsPzaTgshmt5vrcR3hPPwcdsqSxwiRmnOgqgYx15n5RgJUoFmfQQJAENYIIm6y+a7gOAaX9cfcwC7mEXBc5zygCnR3c0cEZmBM7+ASNMKg1cBWYzV8OSwgOcXo/icPbn0IHIvRQwa0Ob/3srkth+oxzAZThnicaiHsQdAfSBkHPxBGaJ75gr6H6eS7lRJnACowInFoJYvHr+/LmJHvI3Aye+MBp1Dc9lkY9RTv2oj+pvwNlsoEezrgJgggHBxa/ukX3g4UDj5bwUaZQwYti+Y+EMr2ZD2/cMGsn/JsnbQ/RNBSmGeN6eKMEz/hzFH33mu6YBF+DOAk4M6GjgZAVFgWeNlHhvPPWlDIlRBnDKwGjmWf4AbCIotACRPR7S5GW0qv8dI5/xcJ5Brg5coErRdgiwceWSoxlStB9/pt3L5/6kjoBFWAmd5R/wJgBXeX2skagPoMiGZAt8RvkZUUCM9/hz9XFpIv4SA+AKsol/KBfGwOSPmMXf+MSfyVFg7AF4bOFJ8UgGazPA9Uc0aDgai+/4G+uik/iKtgpwTa60Nt9pwD3sosDxdw9cmRIpK9cSoPROyhR/yo61+J+nSxqZWATlxw4YQfxNvf2PEriH5QqpFnAp0YtpeVoI32KYAkj/JwugMeLvNpngEDx1H0EQ6+J5NAC2Q/r57t27FCyzW7773YCzfNAywFkYk4CpZWIW1ZgTsAlCJoWlMDo0mvB34gm2QKPE7MevucS5ASOUevjAF0H79XC/gaLMjhkp6aQB52X+V8USBVQl4g88SRV5Od/zDJ4psPguvZJe62H/qcDVKTscuAoo3wWgLEICGhBZJEO6Ji6Y4eEIwNH/PXBGZjxz4y0LsATdOAfgGTSOTxW5hkZn5MYtQoD71PE/7eyQCLd/dVcAAAAASUVORK5CYII=&quot;,
            [[43.60734015888762, 3.767730367838], [43.81775157192726, 4.3423021816889245]],
            {&quot;crossOrigin&quot;: false, &quot;interactive&quot;: false, &quot;opacity&quot;: 0.7, &quot;zindex&quot;: 1}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);


        var layer_control_16772d432ce6954add8c65224c2554ae = {
            base_layers : {
                &quot;openstreetmap&quot; : tile_layer_8ebb7bdf73674872bc467e44a609f99a,
            },
            overlays :  {
                &quot;January&quot; : image_overlay_06ec3d8b1ee2472e0dee5d2551cc2776,
                &quot;February&quot; : image_overlay_a89a77ffe85cc350b6d5e694476e9cb5,
                &quot;March&quot; : image_overlay_9be1b876fc220c2b22bc75a44a970ca5,
                &quot;April&quot; : image_overlay_d1707f042f199bc480911a410925d0ed,
                &quot;May&quot; : image_overlay_047cf7a4389ddf495e7be954fc800a97,
                &quot;June&quot; : image_overlay_35a70477cdb1aec56ce558f46abbb54b,
                &quot;July&quot; : image_overlay_aa3f48d500ceccff387b5b78900cb83f,
                &quot;August&quot; : image_overlay_6610ffc3843fc286b5f922d0ff999b64,
                &quot;September&quot; : image_overlay_ba0f0d14f5a326f6d0321ceb4226f921,
            },
        };
        L.control.layers(
            layer_control_16772d432ce6954add8c65224c2554ae.base_layers,
            layer_control_16772d432ce6954add8c65224c2554ae.overlays,
            {&quot;autoZIndex&quot;: true, &quot;collapsed&quot;: true, &quot;position&quot;: &quot;topright&quot;}
        ).addTo(map_c411caddeadbebf70f7cebfd52330cd5);

</script> </html>” width=”700” height=”500”style=”border:none !important;” “allowfullscreen” “webkitallowfullscreen” “mozallowfullscreen”></iframe>

Plot temporal series for the data

We calculate the mean temperature per month

data_interpolated
<xarray.DataArray 'LST_Day_1km' (time: 12, y: 48, x: 92)>
array([[[   9.69,    9.69,    9.69, ...,   16.55,   16.71,   16.71],
        [   9.69,    9.69,    9.69, ...,   16.55,   16.71,   16.71],
        [   9.69,    9.69,    9.69, ...,   16.73,   16.69,   16.69],
        ...,
        [  15.51,   15.21,   15.21, ...,   11.71,   11.35,   11.35],
        [  15.47,   15.55,   15.55, ...,   11.71,   11.35,   11.35],
        [  15.47,   15.55,   15.55, ...,   13.47,   14.21,   14.21]],

       [[  12.93,   12.93,   12.93, ...,   17.77,   20.53,   20.53],
        [  12.93,   12.93,   12.93, ...,   20.53,   20.53,   20.53],
        [  12.93,   12.93,   12.93, ...,   20.53,   20.53,   20.53],
        ...,
        [  16.67,   17.83,   17.83, ...,   13.59,   14.23,   14.23],
        [  17.83,   17.83,   17.83, ...,   13.59,   14.23,   14.23],
        [  17.83,   17.83,   17.83, ...,   15.63,   15.97,   15.97]],

       [[  23.83,   23.27,   23.27, ...,   31.33,   31.65,   31.65],
        [  24.17,   23.93,   23.93, ...,   31.33,   31.65,   31.65],
        [  24.17,   23.93,   23.93, ...,   31.13,   30.55,   30.55],
        ...,
...
        ...,
        [-273.15, -273.15, -273.15, ..., -273.15, -273.15, -273.15],
        [-273.15, -273.15, -273.15, ..., -273.15, -273.15, -273.15],
        [-273.15, -273.15, -273.15, ..., -273.15, -273.15, -273.15]],

       [[   8.65,    8.65,    8.65, ...,    9.33,    9.35,    9.35],
        [   8.63,    8.63,    8.63, ...,    9.33,    9.35,    9.35],
        [   8.63,    8.63,    8.63, ...,    9.55,    9.47,    9.47],
        ...,
        [  10.17,   10.23,   10.23, ...,    6.01,    6.09,    6.09],
        [  10.21,   10.17,   10.17, ...,    6.01,    6.09,    6.09],
        [  10.21,   10.17,   10.17, ...,    6.91,    7.67,    7.67]],

       [[  16.45,   16.07,   16.07, ...,   13.49,   13.49,   13.49],
        [  16.59,   16.17,   16.17, ...,   13.49,   13.49,   13.49],
        [  16.59,   16.17,   16.17, ...,   13.49,   13.49,   13.49],
        ...,
        [  16.65,   16.99,   16.99, ...,   12.91,   12.91,   12.91],
        [  16.63,   16.79,   16.79, ...,   12.91,   12.91,   12.91],
        [  16.63,   16.79,   16.79, ...,   12.91,   12.91,   12.91]]])
Coordinates:
  * y            (y) float64 6.303e+06 6.302e+06 ... 6.28e+06 6.279e+06
  * x            (x) float64 7.622e+05 7.628e+05 ... 8.072e+05 8.078e+05
  * time         (time) datetime64[ns] 2021-01-31 2021-02-28 ... 2021-12-31
    spatial_ref  int64 0
Attributes:
    _FillValue:  -273.15

To plot temporal series of the temperature, we must calculate the averaage of temperatures per month in each raster. However, we need to deal with the nodata values first, since a high negative value can influence the mean value.

# Remove data where we have absolute kelvin -273.15. This will add np.nan where temperatures are -273.15
data_with_nan = [data_month.data for data_month in data_interpolated.where(data != -273.15)]

# Remove np.nan values
data_without_nan = [data_month[~np.isnan(data_month)] for data_month in data_with_nan]

# Calculate the mean with the present values
data_temp_mean = [np.nanmean(month_data) for month_data in data_without_nan ]

data_temp_mean
[11.963057065217424,
 16.11516304347829,
 26.789483695652205,
 19.04933423913047,
 29.070221920289885,
 33.015901268115975,
 39.149814311594234,
 35.1516077898551,
 19.180919384058,
 -273.1500000000001,
 9.299877717391333,
 14.698691123188437]

The month of October has nodata values only. Lets replace it by the average of the temperatures of September and November

data_temp_mean[9] = np.mean( np.array([ data_temp_mean[8], data_temp_mean[10] ]), axis=0 )
data_temp_mean
[11.963057065217424,
 16.11516304347829,
 26.789483695652205,
 19.04933423913047,
 29.070221920289885,
 33.015901268115975,
 39.149814311594234,
 35.1516077898551,
 19.180919384058,
 14.240398550724667,
 9.299877717391333,
 14.698691123188437]

We create a dataframe with the columns temperature and month

import pandas as pd
df_temp = pd.DataFrame({'month': list(months.keys()), 'temperature':data_temp_mean})
df_temp
month temperature
0 January 11.963057
1 February 16.115163
2 March 26.789484
3 April 19.049334
4 May 29.070222
5 June 33.015901
6 July 39.149814
7 August 35.151608
8 September 19.180919
9 October 14.240399
10 november 9.299878
11 December 14.698691

Plot data with seaborn

import seaborn as sns
sns.set_theme(rc={'figure.figsize':(8.7,7.27)})
sns.lineplot(x="month", y="temperature", data=df_temp)
plt.xticks(rotation=15)
plt.title(f'Montpellier verage temp in {year}')
plt.show()

png

Make a animation from the output raster

Create folders

print("Current Working Directory " , os.getcwd())
import matplotlib.pyplot as pyplot
in_dir = os.getcwd() + '/temp'
os.chdir(in_dir)
print("Input Directory " , in_dir)

Current Working Directory  /content
Input Directory  /content/temp
out_dir = os.getcwd() + '/animation/'
print("Animation out dir",out_dir)
if not os.path.exists(out_dir): os.makedirs(out_dir)
file_list = [file for file in os.listdir() if file.endswith('.tiff')]
print(os.listdir())

Animation out dir /content/temp/animation/
['modis_tile_10_.tiff', 'modis_tile_02_.tiff', 'modis_tile_03_.tiff', 'modis_tile_12_.tiff', 'modis_tile_05_.tiff', 'modis_tile_06_.tiff', 'modis_tile_11_.tiff', 'modis_tile_08_.tiff', 'animation', 'modis_tile_09_.tiff', 'modis_tile_04_.tiff', 'modis_tile_01_.tiff', 'modis_tile_07_.tiff']
for f in file_list:
    raw_file = rasterio.open(f)
    file = raw_file.read(1)

    cmap='turbo'
    plt.figure(figsize=(12,15))
    plt.imshow(file,cmap=cmap,vmin=0,vmax=60)
    parts = f.split("_")
    file_name = parts[2]
    print(f'Processing: {file_name}')
    plt.title(f"Temperature {file_name}", fontsize=28,fontweight='bold')
    cur_axes = plt.gca()
    cur_axes.axes.get_xaxis().set_visible(False)
    cur_axes.axes.get_yaxis().set_visible(False)
    cb = plt.colorbar(fraction=0.046, pad=0.04,shrink=0.3)
    cb.set_label('Temperature time series', size=24)
    cb.ax.tick_params(labelsize=18)
    plt.tight_layout()
    print(out_dir)
    plt.savefig('{}{}_temperature.png'.format(out_dir,file_name,dpi=150))
    plt.close()

Processing: 10
/content/temp/animation/
Processing: 02
/content/temp/animation/
Processing: 03
/content/temp/animation/
Processing: 12
/content/temp/animation/
Processing: 05
/content/temp/animation/
Processing: 06
/content/temp/animation/
Processing: 11
/content/temp/animation/
Processing: 08
/content/temp/animation/
Processing: 09
/content/temp/animation/
Processing: 04
/content/temp/animation/
Processing: 01
/content/temp/animation/
Processing: 07
/content/temp/animation/

Create animation

import imageio.v2 as imageio

def make_gif(input_folder, save_file_path):
    episode_frames = []
    time_per_step = 4
    for root, _, files in os.walk(input_folder):
        file_paths = [os.path.join(root,file) for file in files]
        file_paths = sorted (file_paths,key=lambda x:os.path.getmtime(x))
        episode_frames = [imageio.imread(file_path) for file_path in file_paths if file_path.endswith('.png')]
    episode_frames = np.array(episode_frames)
    imageio.mimsave(save_file_path,episode_frames,duration=6)

make_gif(out_dir,os.path.join(out_dir,"Animation_Temp.gif"))

Visualize the animation

from IPython.display import Image, display
with open(out_dir + "Animation_Temp.gif", 'rb') as f:
    display(Image(data=f.read(),format='png',width = 600, height = 600))

png

Read More

Unleashing the power of graph analysis with Networkx, NzViz and PyViz

Introduction

Graph netowrks have got a lot of attention in the past decade. Almost every type of data can be represented as a graph. In that sense, newtwork analysis can be a powerfull and usesuful approach as an exloratory data analysis approach to graph neural networks (GNN). This post will focus on Network Analysis of a real data network presented by (Valentin et al., 2024). In the next section, I will make a brief introduction about the dataset and the context on which it was produced before we dive into the graph analysis part

Disease outbreak

In the domain of epidemiology, an outbreak is defined as a sudden increase in the disease frequency, related to time, place, and observed population (Reintjes and Zanuzdana, 2009) Deasease outbreaks are usually reported by agencies that are spread around the globe. Thus, these agencies play a key role to minimize effects of any epidemics, , since if outbreak events are reported early, govermantal organizations may act in time to avoid the rapid spread rate of a possible infectious deasese.An agency must declares and event when an alert is triggered and confirmed by local authorities. So the presence of agencies in all geographical levels (local, regional and global) is essential to stop an eventual epidemic. In veterinary epidemiology, the suveillance of this events is more complex, since animal farms are usually on remote areas and latelly reported by certain type of agencies, but not by others.

Thus, the use of graph network visualization tool could give us some insights about why and how agencies declare a certain type of event and how the declaration of events are connected.

The dataset

The official dataset can be found here. This dataset contains sveral csv files, but will focus on 3 main files :

  • The sources.csv : represents the data about the source of outbreak declaration,
  • The event.csv : represents the data about each outbreak event,
  • The event_source.csv : represents all outbreak declared events

If we transform our data into a graph, source and event will be considered as node types, and event_source are edges of our graph. But before transforming the dataset into a graph, it is good to know which type the of graph we are working with

Graph types

In graph theory, graphs can be classed by different ways (e.g. interconnectivity, number of edges, etc). An exaustive list (with a visual representation) with several types of graphs can be found here. In our case, every source of outbreak declares (edge) many events. With means that for every source of outbreak, there is at least one event associated to the that source. And one event can be associated to many sources. This type of graph is call a bipartite graph :


Examples of bipartite graphs
Source https://mathworld.wolfram.com/BipartiteGraph.html


Bipartite graphs are considered to be complete (or biclique) when every graph vertex of the first set is connected to every vertex of the second set. This is not our case, since not every source declares all events

Buiding the graph with Networkx

Networkx is a python library for network analysis, aloow creation, manipulation and many other features when working with simple and complex network. In order to create a graph using Networkx, we must map the the source dataframe into a tuple. Given the formal bipartite graph definition :

\[G = (U,V,E)\]

where \(U\) and \(V\) represent bi parts (nodes) of the graph, and \(E\) denoting the edges of graph. In that sense, must first build some mapping to feed the graph creation using Networkx. This is done separetelly for each node : the mappings are just list of tuples (for each type of node), where each node can contain features, represented by a dictionary of attribute name and value. For the nodes (source and event) creation, the code is the same :

source_nodes_mappings = [(s,{"type_source":ts}) for s, ts in zip(list(source['source']),list(source['type_source']))]
event_nodes_mappings = [(e,{"region":r}) for e, r in zip(list(events['id_event']),list(events['region']))]

and the output with the top four nodes of event would be :

[(323, {'region': 'EUR'}), (309, {'region': 'EUR'}), (353, {'region': 'EUR'}), (325, {'region': 'EUR'})]

Notice that for each node, we added a dictionary containing the attribute name and its value. We could as many as we want. But for the sake of simplicity, we added the type of source (for source nodes) and region (for event nodes). And it’s done for nodes.

To build the edges, we must build a triplet containing \((S,T,E)\) where \(S\) is the source node (in our case, the source dataframe), \(T\) for target node (event dataframe) and \(E\) for every edge, represented by the event_source data frame. Since the edges are all maped in the event_source dataframe, we can proceed applying the following code :

source_event_mapping_list = [(s,e,{"num":n}) for s, e, n in zip(list(events_source_early['source']),
                                                                list(events_source_early['id_event']),
                                                                list(events_source_early['n']))]

and the output with the top four edges of edges would be :

[('Channel NewsAsia', 293, {'num': 1}), ('DesMoinesRegister.com', 293, {'num': 1}),
('Devdiscourse', 293, {'num': 2}), ('Free Malaysia Today', 293, {'num': 1})]

As for nodes, we also decided to add the only one attribute (n) of each edge in the graph. Finally, we use the defined mapings to initialize a network using networkx :

import networkx as nx

# Initialize the graph
G = nx.Graph()

# Add nodes with the node attribute "bipartite"
G.add_nodes_from(source_nodes, bipartite="Source Type",color="blue" )
G.add_nodes_from(event_nodes, bipartite="Region", color="red")

# Add edges only between nodes of opposite node sets
G.add_edges_from(source_event_mapping_list)

print(list(Ge.nodes)[:10])
print(list(Ge.edges)[:4])

The output :

['OIE', 'Malaysia Vet Auth', 'Gulf Business News', 'Vietnam Vet Auth', 'Reuters', 'The Star Online', 
'Sabbah Vet Auth', 'Devdiscourse', 'Deccan Chronicle', 'The Borneo Post']
[('OIE', 293), ('OIE', 294), ('OIE', 295), ('OIE', 296)]

NetworkX will use the exacly names of nodes and edges as it is in when we output the network. In our case, we want source nodes to be named by type_source and event nodes to be named by region. So we must relabel them

Relabel nodes and edges

To relabel nodes and edges, we will use the relabel function provided by Networkx. To apply the function, we must provide the Graph we want to change and a mapping (as a list of tuples) containg the old and the new name of nodes and edges. In our case, we already have a list containing the type of sources : source_nodes_mappings. To create our mapping we do:

# Crete source mapping from a list
mapping_source_to_type_source = dict([(source[0], source[1].get("type_source")) for source in source_nodes_mappings])

# Relabel source nodes
G = nx.relabel_nodes(Gai, mapping_source_to_type_source_ai)

# Show graph nodes after relabel
print(G.nodes[:15])

which gives us now the nodes by type of source :

['int. org.', 'vet. auth.', 'online news', 'press agency', 'off. auth.', 'radio/TV', 'research org.', 
'social platform', 'laboratory', 323, 309, 353, 325, 304, 326]

Notice howver that our event type nodes are still numbers. So we repeat the process for events, adding a mapping for region of events :

# Crete events mapping from a list
mapping_events_to_region = dict([(event[0], event[1].get("region")) for event in event_nodes_mappings])

# Relabel event nodes
G = nx.relabel_nodes(G, mapping_events_to_region)

# Show graph nodes after relabel
print(G.nodes)

And the result :

['int. org.', 'vet. auth.', 'online news', 'press agency', 'off. auth.', 'radio/TV', 'research org.', 
'social platform', 'laboratory','EUR', 'WPR', 'SEAR', 'EMR', 'AFR', 'AMR']

Indicators of centrality (centrality mesure)

In network analysis, indicators of centrality are really useful when we want to determine the importance of distinct nodes in a graph. In essence, these indicators are based on calculations of some sort of distance, to later assign numbers or rakings to nodes within a graph corresponding to their network postion. They reflect a unit’s prominence; in different substantive settings for instance, this may be its structural power, status, prestige, or visibility (Marsden, 2005). There are several centrality mesures in the literature of graph theory, but for the sake of brevity, I will (try to) explain three among them : Degree Centrality, Betweeness centrality, and Closeness Centrality

1. Degree Centrality

The degree centrality algorithm is the easiest centrality algorithm to calculate in graph data science, as it just considers the total number of relationships connected to a particular node. This centrality can be used to find nodes in a social network with many direct friendships and is often considered a measure of activity, e.g. Facebook. In social networks, for instance, the degree centrality can be used to find the the most “popular” individuals (nodes), as the ones with the largest number of followers or friends. The algorithm definition is given by :

\[C_{D}(i) = {\sum_{j=1\\ j\neq i}^{N}x_{ij}}\]

2. Betweeness centrality

Betweenness centrality is a way of detecting the amount of influence a node has over the flow of information in a graph. The algorithm calculates shortest paths between all pairs of nodes in a graph. Each node receives a score, based on the number of shortest paths that pass through the node. The algorithm definition is given by :

\[C_{B}(i) = \sum_{s\neq v\neq t} \frac{\sigma_{st}(v)} {\sigma_{st}}\]

where \(\sigma_{st}\) is the total number of shortest paths from node \(s\) to node \(t\) and \(\sigma_{st}(v)\) is the number of those paths that pass through \(v\) as long as \(v\) is not the end point

3. Closeness Centrality

A closeness measure conceives of a unit as central to the extent that it is related to other units via short geodesics:

\[C_{C}(i) = \frac{N - 1} {\sum_{j=1}^{N}d_{ij}}\]

where \(d_{ij}\) is the distance (length of the shortest path) of the shortest path between nodes i and j. Since the closiness centrality captures the average between each node and every other node in the graph, we use the normalized form of the average by multiplying \(1 / \sum_{j=1}^{N}d_{ij}\) by \(N - 1\).

Low closeness centrality means that a node is directly connected (or “just one connection away”) from most others in the network. In contrast, nodes in very peripheral locations may have high closeness centrality scores, indicating the high number of connections they need to take to connect to distant others in the network (Hansen et al., 2020)

Calculating indicators of centrality

Hopefullly, we won’t have to implement by hand the algorithms shown above. The great python package Networkx has ready-to-use implemetations of these and many other algorithms.

For each centrality measure discussed, we will compute and them in a plot. So we first define a helper function that will take two arguments : the desired centrality measure and the name of centrality label we want to show in the plot:

def generate_centrality_df(centrality_type,graph,plot_col_name):
  centrality = None
  if centrality_type == 'degree' :
    centrality = nx.degree_centrality(graph)
  elif centrality_type == 'betweenness':
    centrality = nx.betweenness_centrality(graph)
  elif centrality_type == 'closeness':
    centrality = nx.closeness_centrality(graph)
  assert centrality is not None, "Invalid/unsupported centrality type"
  df = pd.DataFrame.from_dict(centrality,orient='index',columns=[plot_col_name])
  return df

And we plot the three indicators of centrality:

# Define plot area
fig, axes = plt.subplots(nrows=1, ncols=3,figsize=(17, 5))

# Define dataframe based on graph data
degree_df_ai = generate_centrality_df('degree',Gai,'Degree centrality')
betweenness_df_ppa = generate_centrality_df('betweenness',Gai,'Betweenness')
closeness_df_ppa = generate_centrality_df('closeness',Gai,'Closeness')

# Plot centrality of all graph deaseases
degree_df_ai.sort_values('Degree centrality',ascending=False).plot(kind="bar",color='green',ax=axes[0],rot=75)
betweenness_df_ppa.sort_values('Betweenness',ascending=False).plot(kind="bar",color='blue',ax=axes[1],rot=80)
closeness_df_ppa.sort_values('Closeness',ascending=False).plot(kind="bar",color='red',ax=axes[2],rot=80)


Indicators of centrality


As we can see..

Visualize graph with NzViz

fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(18, 10))

# Define some legend parms
font = {'family': 'serif',
        'color':  'darkred',
        'weight': 'normal',
        'size': 16,
        }

# Define the circos plot
nv.circos(Gai, 
          group_by="label",node_color_by="label",node_size_by="degree",
          node_enc_kwargs={"size_scale": 0.35},edge_lw_by='count')
annotate.circos_labels(Gai, group_by="label", layout="rotate")
legend_dict = dict(ncol=1, loc='upper right', fancybox=True, shadow=True,bbox_to_anchor=(1, 1))
annotate.node_colormapping(Gai, color_by="label",legend_kwargs=legend_dict)

and the plot output :


Circos plot of Avian Influenza's report events by type of reporting source and event region


Plot an interactive graphs with Pyviz

Comunity detection

Plot an interactive graph with Pyviz


References

  1. Marsden, P.V., 2005. Network Analysis. Encyclopedia of Social Measurement 819–825. https://doi.org/https://doi.org/10.1016/B0-12-369398-5/00409-6

  2. Reintjes, R., Zanuzdana, A., 2009. “Outbreak Investigations.” Modern Infectious Disease Epidemiology: Concepts, Methods, Mathematical Models, and Public Health. Concepts, Methods, Mathematical Models, and Public Health 159–176. https://doi.org/https://doi.org/10.1007/978-0-387-93835-6_9

  3. Hansen, D.L., Shneiderman, B., Smith, M.A., Himelboim, I., 2020. Chapter 3 - Social network analysis: Measuring, mapping, and modeling collections of connections. Analyzing Social Media Networks with NodeXL (Second Edition) 31–51. https://doi.org/https://doi.org/10.1016/B978-0-12-817756-3.00003-0

  4. Valentin, S., Boudoua, B., Sewalk, K., Arınık, N., Roche, M., Lancelot, R., Arsevska, E., 2024. Dissemination of information in event-based surveillance, a case study of Avian Influenza 18, e0285341. https://doi.org/https://doi.org/10.1371/journal.pone.0285341

Read More