{"cells":[{"metadata":{},"cell_type":"markdown","source":"## Note\n\nI created this kernel at immediately after the start of this competition.\n\nLater Paul Mooney released the following kernel for the explanation of band of dataset.  \n\nhttps://www.kaggle.com/paultimothymooney/explore-image-metadata-s5p-gfs-gldas\n\nSo please take this kernel as one example of dataset EDA."},{"metadata":{},"cell_type":"markdown","source":"# About This Kernel\nI'll summarize how to understand sample tiffiles in terms of it's bands.\n\nFirst, I try to check which band in the tifffiles corresponds to which band on the Dataset bands.\n\nSecond, I'll get data by Google Earth Engin API with specific band.\n\nI am glad to help you understand data!"},{"metadata":{},"cell_type":"markdown","source":"## Load Required Library"},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport os\nimport rasterio as rio\nimport folium\nimport tifffile as tiff ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from kaggle_secrets import UserSecretsClient\nfrom google.oauth2.credentials import Credentials\nimport ee","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Snippets for visualize data"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"def overlay_image_on_puerto_rico(file_name,band_layer,lat,lon,zoom):\n    band = rio.open(file_name).read(band_layer)\n    m = folium.Map([lat, lon], zoom_start=zoom)\n    folium.raster_layers.ImageOverlay(\n        image=band,\n        bounds = [[18.6,-67.3,],[17.9,-65.2]],\n        colormap=lambda x: (1, 0, 0, x),\n    ).add_to(m)\n    return m","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Checking Shape of tif files"},{"metadata":{"trusted":true},"cell_type":"code","source":"gfs_img = tiff.imread(\"/kaggle/input/ds4g-environmental-insights-explorer/eie_data/gfs/gfs_2018070100.tif\")\ngldas_img = tiff.imread(\"/kaggle/input/ds4g-environmental-insights-explorer/eie_data/gldas/gldas_20180702_0600.tif\")\ns5p_no2_img = tiff.imread(\"/kaggle/input/ds4g-environmental-insights-explorer/eie_data/s5p_no2/s5p_no2_20180701T161259_20180707T175356.tif\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(\"gfs shape: \", gfs_img.shape)\nprint(\"gldas shape: \", gldas_img.shape)\nprint(\"s5p_no2 shape: \", s5p_no2_img.shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Check Bands of sample tiff "},{"metadata":{},"cell_type":"markdown","source":"## **Sentinel-5P OFFL NO2: Offline Nitrogen Dioxide**"},{"metadata":{},"cell_type":"markdown","source":"According to \"Sentinel 5P OFFL NO2 Dataset\" (below URL), there are 12 Bands in this data.\n\nhttps://developers.google.com/earth-engine/datasets/catalog/COPERNICUS_S5P_OFFL_L3_NO2"},{"metadata":{"trusted":true},"cell_type":"code","source":"NO2_Bands = pd.read_csv(\"../input/ds4gbands-explanation/Sentinel-5P_OFFL _NO2_Bands.csv\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":false},"cell_type":"code","source":"band1 = s5p_no2_img[:,:,0].reshape([70300,])\nband2 = s5p_no2_img[:,:,1].reshape([70300,])\nband3 = s5p_no2_img[:,:,2].reshape([70300,])\nband4 = s5p_no2_img[:,:,3].reshape([70300,])\nband5 = s5p_no2_img[:,:,4].reshape([70300,])\nband6 = s5p_no2_img[:,:,5].reshape([70300,])\nband7= s5p_no2_img[:,:,6].reshape([70300,])\nband8 = s5p_no2_img[:,:,7].reshape([70300,])\nband9 = s5p_no2_img[:,:,8].reshape([70300,])\nband10 = s5p_no2_img[:,:,9].reshape([70300,])\nband11 = s5p_no2_img[:,:,10].reshape([70300,])\nband12 = s5p_no2_img[:,:,11].reshape([70300,])\nbands = [band1, band2, band3, band4, band5, band6, band7, band8, band9, band10, band11, band12]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"NO2_Bands","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"for i, band in enumerate(bands):\n    print(\"band: \", i+1)\n    print(\"max: \", max(band))\n    print(\"min: \", min(band))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"We can find that the bands seems same order order of Dataset web site.\n\nWhen we want to retrieve information of NO2, it would be nice to get information on bands 1 to 4."},{"metadata":{"trusted":true},"cell_type":"code","source":"image = '/kaggle/input/ds4g-environmental-insights-explorer/eie_data/s5p_no2/s5p_no2_20180701T161259_20180707T175356.tif'\nlatitude=18.1429005246921; longitude=-65.4440010699994\noverlay_image_on_puerto_rico(image,band_layer=1,lat=latitude,lon=longitude,zoom=8)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## **GFS: Global Forecast System 384-Hour Predicted Atmosphere Data**"},{"metadata":{"_cell_guid":"","_uuid":"","trusted":true},"cell_type":"markdown","source":"According to \"Global Forecast System 384-Hour Predicted Atmosphere Data\" (below URL), there are 9 Bands in this data.\n\nhttps://developers.google.com/earth-engine/datasets/catalog/NOAA_GFS0P25"},{"metadata":{"trusted":true},"cell_type":"code","source":"GFS_Bands = pd.read_csv(\"../input/ds4gbands-explanation/GFS-Global_Forecast_System_384-Hour_Predicted_Atmosphere_Data_Bands.csv\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":false},"cell_type":"code","source":"band1 = gfs_img[:,:,0].reshape([70300,])\nband2 = gfs_img[:,:,1].reshape([70300,])\nband3 = gfs_img[:,:,2].reshape([70300,])\nband4 = gfs_img[:,:,3].reshape([70300,])\nband5 = gfs_img[:,:,4].reshape([70300,])\nband6 = gfs_img[:,:,5].reshape([70300,])\nbands = [band1, band2, band3, band4, band5, band6]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"GFS_Bands","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"for i, band in enumerate(bands):\n    print(\"band: \", i+1)\n    print(\"max: \", max(band))\n    print(\"min: \", min(band))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"This dataset has different bands from website.\n\nWe can find that the bands seems same order order of Dataset web site, except for the last one.\n\nI can't understand what the last band corresponded to..."},{"metadata":{},"cell_type":"markdown","source":"## **Global Land Data Assimilation System**"},{"metadata":{},"cell_type":"markdown","source":"According to \"GLDAS-2.1: Global Land Data Assimilation System\" (below URL), there are 36 Bands in this data.\n\nhttps://developers.google.com/earth-engine/datasets/catalog/NASA_GLDAS_V021_NOAH_G025_T3H"},{"metadata":{"trusted":true},"cell_type":"code","source":"GLDAS_Bands = pd.read_csv(\"../input/ds4gbands-explanation/GLDAS-2_1- Global_Land_Data_Assimilation_System_Bands.csv\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"band1 = gldas_img[:,:,0].reshape([70300,])\nband2 = gldas_img[:,:,1].reshape([70300,])\nband3 = gldas_img[:,:,2].reshape([70300,])\nband4 = gldas_img[:,:,3].reshape([70300,])\nband5 = gldas_img[:,:,4].reshape([70300,])\nband6 = gldas_img[:,:,5].reshape([70300,])\nband7= gldas_img[:,:,6].reshape([70300,])\nband8 = gldas_img[:,:,7].reshape([70300,])\nband9 = gldas_img[:,:,8].reshape([70300,])\nband10 = gldas_img[:,:,9].reshape([70300,])\nband11 = gldas_img[:,:,10].reshape([70300,])\nband12 = gldas_img[:,:,11].reshape([70300,])\nbands = [band1, band2, band3, band4, band5, band6, band7, band8, band9, band10, band11, band12]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"GLDAS_Bands","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"for i, band in enumerate(bands):\n    print(\"band: \", i+1)\n    print(\"max: \", np.nanmax(band))\n    print(\"min: \", np.nanmin(band))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"In this dataset, bands seems be shaked.\n\nWhen using band information, we may need to choose carefully."},{"metadata":{"trusted":true},"cell_type":"markdown","source":"# Get data by Google Earth Engin API"},{"metadata":{},"cell_type":"markdown","source":"While referring to the great kernels, I'll get 'NO2_column_number_density' data of \"COPERNICUS/S5P/NRTI/L3_NO2\" dataset.\n\nYou have to setup for using Google Earth Engin API:\nhttps://www.kaggle.com/paultimothymooney/how-to-get-started-with-the-earth-engine-data\n\nIf you use User Secrets first, ferer:\nhttps://www.kaggle.com/product-feedback/114053"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"user_secret = \"earth_engine\" # Your user secret, defined in the add-on menu of the notebook editor\nrefresh_token = UserSecretsClient().get_secret(user_secret)\ncredentials = Credentials(\n        None,\n        refresh_token=refresh_token,\n        token_uri=ee.oauth.TOKEN_URI,\n        client_id=ee.oauth.CLIENT_ID,\n        client_secret=ee.oauth.CLIENT_SECRET,\n        scopes=ee.oauth.SCOPES)\nee.Initialize(credentials=credentials)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def add_ee_layer(self, ee_image_object, vis_params, name):\n  # https://github.com/google/earthengine-api/blob/master/python/examples/ipynb/ee-api-colab-setup.ipynb\n  map_id_dict = ee.Image(ee_image_object).getMapId(vis_params)\n  folium.raster_layers.TileLayer(\n    tiles = map_id_dict['tile_fetcher'].url_format,\n    attr = 'Map Data &copy; <a href=\"https://earthengine.google.com/\">Google Earth Engine</a>',\n    name = name,\n    overlay = True,\n    control = True\n  ).add_to(self)\n\ndef plot_ee_data_on_map(dataset,column,begin_date,end_date,minimum_value,maximum_value,latitude,longitude,zoom):\n    # https://github.com/google/earthengine-api/blob/master/python/examples/ipynb/ee-api-colab-setup.ipynb\n    folium.Map.add_ee_layer = add_ee_layer\n    vis_params = {\n      'min': minimum_value,\n      'max': maximum_value,\n      'palette': ['006633', 'E5FFCC', '662A00', 'D8D8D8', 'F5F5F5']}\n    my_map = folium.Map(location=[latitude,longitude], zoom_start=zoom, height=500)\n    s5p = ee.ImageCollection(dataset).filterDate(\n        begin_date, end_date)\n    my_map.add_ee_layer(s5p.first().select(column), vis_params, 'Color')\n    my_map.add_child(folium.LayerControl())\n    display(my_map)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"dataset = \"COPERNICUS/S5P/NRTI/L3_NO2\"\ncolumn = 'NO2_column_number_density'\nbegin_date = '2018-07-08'\nend_date = '2018-07-14'\nminimum_value = -0.00051\nmaximum_value = 0.0192\nlatitude = 60.17\nlongitude = 25.94\nzoom = 5\nplot_ee_data_on_map(dataset,column,begin_date,end_date,minimum_value,maximum_value,latitude,longitude,zoom)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Challenge the application of Google Earth Engin API "},{"metadata":{},"cell_type":"markdown","source":"I'll get ImageCollection of NO2_column_number_density cut out by time rangetime series. And I'll reduces the ImageCollection by calculating the median of all values at each pixel. Finally, I'll try to visualize the result. \n\nWe can refer Google Earth Engin API here!\n\nhttps://developers.google.com/earth-engine/api_docs"},{"metadata":{"trusted":true},"cell_type":"code","source":"dataset_name = \"COPERNICUS/S5P/NRTI/L3_NO2\"\ncolumn = 'NO2_column_number_density'\nbegin_date = '2017-03-07'\nend_date = '2018-12-31'\nminimum_value = -0.00051\nmaximum_value = 0.0192\nzoom = 2","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"r = 1\nregion = [longitude-r,latitude-r,longitude+r,latitude+r]\nregion=ee.Geometry.Rectangle(region);","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"L3_NO2 = ee.ImageCollection(dataset_name);\ndataset = L3_NO2.filterBounds(region).filterDate('2017-03-07','2018-12-31')\nprint(dataset)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"image = dataset.median()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"print(image)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"image.getInfo()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def add_ee_layer_median(self, ee_image_object,vis_params,name):\n  map_id_dict = ee.Image(ee_image_object).getMapId(vis_params)\n  folium.raster_layers.TileLayer(\n    tiles = map_id_dict['tile_fetcher'].url_format,\n    attr = 'Map Data &copy; <a href=\"https://earthengine.google.com/\">Google Earth Engine</a>',\n    name = name,\n    overlay = True,\n    control = True,\n  ).add_to(self)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def plot_ee_data_on_map_median(image,minimum_value,maximum_value,latitude,longitude,zoom):\n    folium.Map.add_ee_layer = add_ee_layer_median\n    vis_params = {\n      'min': minimum_value,\n      'max': maximum_value,\n      'bands': ['NO2_column_number_density'],\n      'palette': ['006633', 'E5FFCC', '662A00', 'D8D8D8', 'F5F5F5']\n    }\n\n    my_map = folium.Map(location=[latitude,longitude], zoom_start=zoom, height=500)\n    my_map.add_ee_layer(image, vis_params, 'NO2')\n    my_map.add_child(folium.LayerControl())\n    display(my_map)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plot_ee_data_on_map_median(image,minimum_value,maximum_value,latitude,longitude,zoom)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Future work"},{"metadata":{},"cell_type":"markdown","source":"It looks like I have successfully plotted the results on the map. But it is too wide! 😢\n\nNext, I'll try to control Google Earth Engine API and folium."},{"metadata":{},"cell_type":"markdown","source":"Thak you for reading!"}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}