<jupyter_start><jupyter_text>31 - AC SAF GOME-2 Level 2 case study >>  2.5 Sentinel-3 OLCI Level-1 data Sentinel-3 OLCI Level-1 data products are available as Full (`OL_1_EFR`) and Reduced Resolution (`OL_1_ERR`) data files. The following notebook shows you how `OL_1_EFR` data are structured and how you can generate a RGB True- and False-Color image.Sentinel-3 OLCI data products are disseminated as `.zip archives`, containing data files in `NetCDF` format.Find more information on the OLCI data products in the [Sentinel-3 OLCI User Guide](https://sentinel.esa.int/web/sentinel/user-guides/sentinel-3-olci).  This module has the following outline:* [1 - Load and browse Sentinel-3 OLCI Level-1 data](load_browse)* [2 - Load OLCI channel information as xarray](channel_info)* [3 - Load OLCI geographic coordinates as xarray](load_geo_coords)* [4 - Select OLCI channels for a RGB image](channel_select)* [5 - Advanced image processing - Normalization and histogram equalization](image_processing)* [6 - Visualize a georeferenced RGB True Color image with pcolormesh](visualize_tc)* [7 - Create a False Color RGB image and plot it as georeferenced image ](false_color)  Load required libraries<jupyter_code>%matplotlib inline

import xarray as xr
import numpy as np
import matplotlib.pyplot as plt

import cartopy.crs as ccrs
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
import cartopy.feature as cfeature
from skimage import exposure

import glob

import warnings
warnings.filterwarnings('ignore')
warnings.simplefilter(action = "ignore", category = RuntimeWarning)<jupyter_output><empty_output><jupyter_text>Load and browse Sentinel-3 OLCI Level-1 data Sentinel-3 Level-1 data are dissiminated as `.zip archives` when downloaded. The unzipped folder contains 31 data files in `NetCDF` format. Data for each channel is stored in a single `NetCDF` file. Additionally, you get information on `qualityFlags`, `time_coordinates` or `geo_coordinates`.You can see the names of the 31 data files by looping through the data directory. You see that the channel information follow the same naming and all end with `_radiance.nc`.<jupyter_code>olci_dir = './eodata/sentinel3/olci/2019/07/25/S3A_OL_1_EFR____20190725T022702_20190725T023002_20190726T055327_0179_047_203_1800_LN1_O_NT_002.SEN3/'
for i in glob.glob(olci_dir+'*.nc'):
    tmp = i.split('/')
    print(tmp[8])<jupyter_output><empty_output><jupyter_text>Load OLCI channel information as `xarray`  Load one single channel As a first step, you can load one channel with xarray's function `open_dataset`. This will help you to understand how the data is structured.You see that the data of each channel is a two dimensional data array, with `rows` and `columns` as dimensions.<jupyter_code>olci_xr = xr.open_dataset(olci_dir+'Oa01_radiance.nc')
olci_xr<jupyter_output><empty_output><jupyter_text>Load all channel information into one xarray object The next step is to bring the individual data files together into one `xarray` data object. You can do this with xarray's function `open_mfdataset`. Since all channel files have the same data structure, you can combine all channels by coordinates using the keyword argument `by_coords`. The result is a `xarray` object with 21 data variables. Each channel is represented as a data variable.<jupyter_code>olci_xr_mf = xr.open_mfdataset(olci_dir+'*_radiance.nc', combine='by_coords')
olci_xr_mf<jupyter_output><empty_output><jupyter_text>Example to plot one channel To get an impression of the image, you can simply plot one channel of the xarray object with xarray's function  `.plot.imshow()`. Let's plot channel 8 `Oa08_radiance`. The example shows that you can easily visualize the image itself, but without geographic information.<jupyter_code>fig = plt.figure(figsize=(10,6))
olci_xr_mf['Oa08_radiance'].plot.imshow()<jupyter_output><empty_output><jupyter_text>Load OLCI geographic coordinates as `xarray` If you want to georeference your image, you have to load the geographic coordinates file. You can load it as xarray with `open_dataset`. The file is called `geo_coordinates.nc`. You see that the file contains three variables: `latitude`, `longitude` and `altitude`.Let's store the latitude and longitude data as `lat` and `lon` variables respectively.<jupyter_code>olci_geo_coords = xr.open_dataset(olci_dir+'geo_coordinates.nc')
olci_geo_coords
lat = olci_geo_coords.latitude.data
lon = olci_geo_coords.longitude.data<jupyter_output><empty_output><jupyter_text>Select OLCI channels for a RGB image Depending on the combination of different OLCI channels, your RGB image might highlight specific phenomena. The channel combination for a True Color RGB image could be:- Red: `Oa08_radiance`- Green: `Oa06_radiance`- Blue: `Oa04_radiance`Let's define a function called [select_channels_for_rgb](./ltpy_functions.ipynbrgb_channels), which makes the channel selection more flexible. The function returns the three bands individually.<jupyter_code>def select_channels_for_rgb(xarray, red_channel, green_channel, blue_channel):
    """ 
    Selects the channels / bands of a multi-dimensional xarray for red, green and blue composites.
    
    Parameters:
        xarray(xarray Dataset): xarray Dataset object that stores the different channels / bands.
        red_channel(str): Name of red channel to be selected
        green_channel(str): Name of green channel to be selected
        blue_channel(str): Name of blue channel to be selected

    Returns:
        Three xarray DataArray objects with selected channels / bands
    """  
    return xarray[red_channel], xarray[green_channel], xarray[blue_channel]
red, green, blue = select_channels_for_rgb(olci_xr_mf, 'Oa08_radiance', 'Oa06_radiance', 'Oa04_radiance')
red<jupyter_output><empty_output><jupyter_text>Advanced image processing - Normalization and histogram equalization  Normalization A common operation in image processing is the normalization of data values. Normalization changes the range of pixel intesity and can improve the constrast.Let's define a function called [normalize](./ltpy_functions.ipynbnormalize), which normalizes a numpy array into a scale between 0.0 and 1.0.<jupyter_code>def normalize(array):
    """ 
    Normalizes a numpy array / xarray DataArray object value to values between 0 and 1.
    
    Parameters:
        xarray(numpy array or xarray DataArray): xarray DataArray or numpy array object.

    Returns:
        Normalized array
    """ 
    array_min, array_max = array.min(), array.max()
    return ((array - array_min)/(array_max - array_min))<jupyter_output><empty_output><jupyter_text>You can now apply the function [normalize](./ltpy_functions.ipynbnormalize) to each RGB channel. At the end, you can bring the three channels together into one `rgb` array with `numpy.dstack`. By verifying the shape of the resulting array, you see that the `rgb` array has now three dimensions.<jupyter_code>redn = normalize(red)
greenn = normalize(green)
bluen = normalize(blue)
rgb = np.dstack((redn, greenn, bluen))
rgb.shape<jupyter_output><empty_output><jupyter_text>You can now plot the rgb image to see if maybe further image processing methods could be applied. If you visualize the image with `plt.imshow()`, you see that the constrast is not very strong. A further step is to sharpen the contrast with the help of a **Histogram equalization**.<jupyter_code>fig = plt.figure(figsize=(10,6))
plt.imshow(rgb)<jupyter_output><empty_output><jupyter_text>Histogram equalization Histogram equalization is a method in image processing that adjusts the contrast using the image's histogram. Python's `skikit-learn` library has useful tools to make a histogram equalization quite straighforward. The `skimage` library provides a function `exposure.equalize_adaphtist()` which can be applied to the rgb data array.<jupyter_code>rgb = exposure.equalize_adapthist(rgb)<jupyter_output><empty_output><jupyter_text>You can again plot the rgb image with `plt.imshow()`. You can see that the histogram equalization improved the contrast of the image.<jupyter_code>fig = plt.figure(figsize=(10,6))
plt.imshow(rgb)<jupyter_output><empty_output><jupyter_text>Visualize a georeferenced RGB True Color image with `pcolormesh` The final step is to georeference the `RGB` image. Therefore, you need the extracted `lat` and `lon` information extracted above.Let's plot the image with matplotlib's `pcolormesh` function. To be able to visualize the array as RGB, you have to map first the rgb array to a colour array.<jupyter_code>mesh_rgb = rgb[:, :-1, :]
colorTuple = mesh_rgb.reshape((mesh_rgb.shape[0] * mesh_rgb.shape[1]), 3)
colorTuple = np.insert(colorTuple, 3, 1.0, axis=1)<jupyter_output><empty_output><jupyter_text>The last step is to visualize the color array and add additional information, e.g. gridlines and coastlines. As you want to reuse the plotting code again, you can define a function called [visualize_s3_pcolormesh](./ltpy_functions.ipynbvisualize_s3_pcolormesh).<jupyter_code>def visualize_s3_pcolormesh(color_array, array, latitude, longitude, title):
    """ 
    Visualizes a numpy array (Sentinel-3 data) with matplotlib's 'pcolormesh' function as RGB image.
    
    Parameters:
        color_array (numpy MaskedArray): any numpy MaskedArray, e.g. loaded with the NetCDF library and the Dataset function
        longitude (numpy Array): array with longitude values
        latitude (numpy Array) : array with latitude values
        title (str): title of the resulting plot
    """
    fig=plt.figure(figsize=(20, 12))

    ax=plt.axes(projection=ccrs.Mercator())
    ax.coastlines()

    gl = ax.gridlines(draw_labels=True, linestyle='--')
    gl.xlabels_top=False
    gl.ylabels_right=False
    gl.xformatter=LONGITUDE_FORMATTER
    gl.yformatter=LATITUDE_FORMATTER
    gl.xlabel_style={'size':14}
    gl.ylabel_style={'size':14}

    img1 = plt.pcolormesh(longitude, latitude, array*np.nan, color=colorTuple,
                          clip_on = True,
                          edgecolors=None,
                          zorder=0,
                          transform=ccrs.PlateCarree())
    ax.set_title(title, fontsize=20, pad=20.0)
    plt.show()
visualize_s3_pcolormesh(colorTuple,red, lat, lon, 'Sentinel-3 OLCI Level-1 True Color RGB - "25 Jul 2019"')<jupyter_output><empty_output><jupyter_text>Create a False Color RGB image and plot it as georeferenced image  A false color RGB image can be applied if specific phenomena should be highlighted. A channel combination for Sentinel-3 OLCI data could be:- Red: `Oa17_radiance`- Green: `Oa05_radiance`- Blue: `Oa02_radiance`This selection is advantageous to identify burnt areas and wildfires. The combination highlights healthy vegeation as red and burnt areas as black.Let's repeat the steps from above to visualize the image as False Color RGB composite.  Select RGB channels The first step is to select the channels `Oa17_radiance`, `Oa05_radiance` and `Oa02_radiance` as red, green and blue channels respectively.<jupyter_code>red_fc, green_fc, blue_fc = select_channels_for_rgb(olci_xr_mf, 'Oa17_radiance', 'Oa05_radiance', 'Oa02_radiance')<jupyter_output><empty_output><jupyter_text>Normalize You also want to normalize these channels and stack them afterwards into a three-dimensional array.<jupyter_code>redn_fc = normalize(red_fc)
greenn_fc = normalize(green_fc)
bluen_fc = normalize(blue_fc)
rgb_fc = np.dstack((redn_fc, greenn_fc, bluen_fc))<jupyter_output><empty_output><jupyter_text>Histogram equalization Apply histogram equalization.<jupyter_code>rgb_fc = exposure.equalize_adapthist(rgb_fc)<jupyter_output><empty_output><jupyter_text>Let's have a sneak-peak to the un-georefrenced image with `plt.imshow()`<jupyter_code>fig = plt.figure(figsize=(10,6))
plt.imshow(rgb_fc)<jupyter_output><empty_output><jupyter_text>Visualize with function `visualize_s3_pcolormesh`<jupyter_code>mesh_rgb = rgb_fc[:, :-1, :]
colorTuple = mesh_rgb.reshape((mesh_rgb.shape[0] * mesh_rgb.shape[1]), 3)
colorTuple = np.insert(colorTuple, 3, 1.0, axis=1)
visualize_s3_pcolormesh(colorTuple, red_fc, lat, lon, 'Sentinel-3 OLCI Level-1 False Color RGB - "25 Jul 2019"')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.08277221218590902, 'fraction_numerical': 0.017901133191000163, 'mean_word_length': 5.349843587069865, 'pattern_counts': {'":': 0, '<': 77, '<?xml version=': 0, '>': 79, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '19', 'dir': 'jupyter-structured-clean-dedup', 'id': '468379', 'license': 'MIT', 'n_tokens_mistral': 3951, 'n_tokens_neox': 3685, 'n_words': 1264, 'path': '25_ltpy_Sentinel3_OLCI_L1.ipynb', 'repo_name': 'federicofierli1/LTPy_Notebooks_satellite_copernicus', 'size': '12178'}