<jupyter_start><jupyter_text>221 - PMAp Product - Aerosol Optical Depth - Level 2 >> 20 - DATA DISCOVERY SEE ALSO* **30 - CASE STUDIES**  - [317 - Indonesia Fires 2015](../30_case_studies/317_fire_indonesia_2015.ipynb)  - [318 - Indonesia Fires 2020](../30_case_studies/318_fire_indonesia_2020.ipynb)     - [319 - Portugal Fires 2020](../30_case_studies/319_fire_portugal_2020.ipynb)   2.1.5 AC SAF - GOME-2 - Level 2 - Absorbing Aerosol Height (AAH) data The data of the GOME-2 `Absorbing Aerosol Height (AAH)` are provided by KNMI in the framework of the [EUMETSAT Satellite Application Facility on Atmospheric Composition Monitoring (AC SAF)](https://acsaf.org/index.html).The Absorbing Aerosol Height (AAH) has been released in 2020. It uses the Absorbing Aerosol Index (AAI) product for aerosol detection, in order to derive the actual height of the absorbing aerosol layer.The data can be retrieved from the [AC-SAF data access portal](https://acsaf.org/offline_access.php).Find more information on the GOME-2 Level 2 AAH data product processed by KNMI in the [product user manual](https://acsaf.org/docs/pum/Product_User_Manual_AAH_Aug_2020.pdf).  Indonesian Fires 2015In 2015, devastating wildfires in Indonesia broke out. As reported by [Reuters](https://www.reuters.com/article/us-indonesia-politics-haze-analysis-idUSKCN1UV0ZJ), "Between June and October that year, about 2.6 million hectares (6.4 million acres) of land burned in Indonesia, mainly on Sumatra and Borneo islands, said a 2016 World Bank report. A Harvard University study linked the 2015 haze to more than 100,000 premature deaths in Indonesia, Malaysia and Singapore."  Module outline:* [1 - Load and browse AC-SAF GOME-2 Level 2 Absorbing Aerosol Height data](load_browse_aah)* [2 - Open a HDF5 group content as xarray.Dataset object](open_aah)* [3 - Build a xarray.DataArray object for a Metop-B GOME-2 Level 2 parameter](build_xr_aah)* [4 - Combine multiple HDF5 files in a xarray.DataArray object](combine_aah)* [5 - Filter AAH values based on Absorbing Aerosol Index (AAI) information](filter_aah)* [6 - Visualize filtered AAH values](visualize_aah)   Load required libraries<jupyter_code>import os
from netCDF4 import Dataset
import xarray as xr
import numpy as np
import pandas as pd
from datetime import datetime

from IPython.display import HTML
# Python libraries for visualisation
from matplotlib import pyplot as plt
from matplotlib import animation


import cartopy.crs as ccrs
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER<jupyter_output><empty_output><jupyter_text>Load helper functions<jupyter_code>%run ../functions.ipynb<jupyter_output><empty_output><jupyter_text>Load and browse AC-SAF GOME-2 Level 2 Absorbing Aerosol Height data The AC-SAF Metop-B GOME-2 Level 2 Absorbing Aerosol Height data is disseminated in the `HDF5` format. [xarray](http://xarray.pydata.org/en/stable/index.html) and [netCDF4](https://unidata.github.io/netcdf4-python/netCDF4/index.html) are two Python libraries you can use to access and manipulate datasets in `HDF5` format.  Discover the structure of AC-SAF Metop-B GOME-2 Level 2 product data files The AAH data is in the folder `/eodata/gome2/aah/`. The data in the `HDF5` file is organized in groups. AC-SAF Metop-B GOME-2 Level 2 data contain four groups: * `DATA` * Absorbing Aerosol Height information* `GEOLOCATION` * Date & Time (all retrievals) * 4 Corner Coordinates & Centre of Ground Pixel * Solar, Line-of-Sight Zenith, and realtive Azimuth angles at TOA * Orbit information* `METADATA` * General Metadata* `PRODUCT_SPECIFIC_METADATA` * Product Metadata  You can load a `HDF5` file with the `Dataset` constructor of the [netCDF4](https://unidata.github.io/netcdf4-python/netCDF4/index.html) library. Load the file in `diskless non-persistent` mode by setting the `diskless` argument to `True` and the `persist` argument to `False`.The resulting object is a `netCDF4.Dataset` object, which acts like a Python dictionary. Thus, with the `keys` function you can list the different groups the file contains.<jupyter_code>file_name = '../eodata/gome2/aah/2015/09/09/S-O3M_GOME_ARS_02_M01_20150909005657Z_20150909015657Z_R_O_20201120210202Z.hdf5'

file = Dataset(file_name, diskless=True, persist=False)

groups = file.groups.keys()
groups<jupyter_output><empty_output><jupyter_text>Each group contains multiple variables. You can print out the variables of each group nicely formatted. This gives you a helpful overview of available parameters. You will mostly be interested in parameters of the `GEOLOCATION` and `DATA` groups.The overview shows you, that one Metop-B GOME-2 Level 2 data file contains information on `AAH_AbsorbingAerosolHeight`.<jupyter_code>for group in groups:
    print('{}:'.format(group))
    # variables in a group
    variables = sorted(file[group].variables.keys())
    # pretty print of variables in 3 columns
    for v1, v2, v3 in zip(variables[::3], variables[1::3], variables[2::3]):
        print('   - {:<30}- {:<30}- {:<}'.format(v1, v2, v3))
    if (len(variables) % 3) == 1:
        print('   - {:<30}'.format(variables[-1]))
    elif (len(variables) % 3) == 2:
        print('   - {:<30}- {:<}'.format(*variables[-2:]))
    # white line between one line and another
    print()<jupyter_output>DATA:
   - AAH_AbsorbingAerosolHeight    - AAH_AbsorbingAerosolHeightError- AAH_AbsorbingAerosolPressure
   - AAH_AbsorbingAerosolPressureError- AAH_ErrorFlag                 - AAH_NElements
   - AAH_QualityInput              - AAH_QualityProcessing         - AAH_RegimeFlag
   - AAI                           - DegradationCorrectedResidue   - NElements
   - PMD_CloudFraction             - PMD_SceneHomogeneity          - PMD_SnowIceFrac
   - QualityInput                  - QualityProcessing             - SunGlintFlag
   - UncorrectedResidue            

GEOLOCATION:
   - IndexInScan                   - LatitudeCenter                - LatitudeCorner
   - LineOfSightAzimuthAngle       - LineOfSightZenithAngle        - LongitudeCenter
   - LongitudeCorner               - NElements                     - NrOfPixelsInScan
   - RelAzimuthAngle               - ScanDirection                 - ScannerAngle
   - ScatteringAngle               - SolarAzimuthAngle             - SolarZenithAngle
   - S[...]<jupyter_text>Open a `HDF5` group content as `xarray.Dataset` object In a next step, you can make use of the `xarray.backends.NetCDF4DataStore` feature of the [xarray](http://xarray.pydata.org/en/stable/index.html) library to open the `HDF5` group you are interested in as a xarray `Dataset` object. It is important to know that you can only get hold of one `HDF5` group item at the time.As an example, you can load the `DATA` group and open it as xarray `Dataset`. You see that the data has two dimensions: `phony_dim_0` and `phony_dim_1` and 19 data variables.<jupyter_code>data = file['DATA']
xr_data = xr.open_dataset(xr.backends.NetCDF4DataStore(data))
xr_data<jupyter_output><empty_output><jupyter_text>Select and retrieve one specific parameter from a `xarray.Dataset` object Let's say, you are interested in the Absorbing Aerosol Height parameter. With `xarray`, you can retrieve the parameter using square brackets `[]` and the name of the parameter of interest. The resulting object is a `xarray.DataArray`, which has the shape of a one-dimensional vector.The `xarray.DataArray` holds also data attributes, which help a better understanding of the parameter, e.g. you get information about minimum and maximum values of the parameter and the unit of the data values.<jupyter_code>aah = xr_data['AAH_AbsorbingAerosolHeight']
aah<jupyter_output><empty_output><jupyter_text>Build a `xarray.DataArray` object for a Metop-B GOME-2 Level 2 parameter For easier handling and plotting of the parameters, you can combine geolocation information and data values in a xarray `DataArray` structure, which can be created with the xarray `DataArray` constructor. You can make use of the function `generate_xr_from_1D_vec` which takes the following kwargs:* `file` - the data file you want to load the data information from* `lat_path` - internal paths to the latitude information of the HDF5 file* `lon_path` - internal paths to the longitude information of the HDF5 file* `variable` - xarray.DataArray with the variable information* `parameter_name` - Name of the variable* `longname` - Long name information of the variable* `no_of_dims` - Number of dimensions of the data file* `unit` - unit of the variableWe take the centre of the ground pixel as coordinates of the `DataArray`. The values can be retrieved from the internal paths `GEOLOCATION/LatitudeCentre` and `GEOLOCATION/LongitudeCentre` respectively.The resulting object is a xarray `DataArray` object with latitude and longitude information and the AAH values as data values. The `DataArray` object contains 32 x 588 data pixels.<jupyter_code>aah_xr = generate_xr_from_1D_vec(file = file, 
                                 lat_path = 'GEOLOCATION/LatitudeCenter', 
                                 lon_path = 'GEOLOCATION/LongitudeCenter', 
                                 variable = aah, 
                                 parameter_name = aah.name, 
                                 longname = aah.Title,
                                 no_of_dims = 2, 
                                 unit = aah.Unit)

aah_xr<jupyter_output><empty_output><jupyter_text>You can now continue to visualize the `DataArray` object. You can use the function `visualize_scatter()` to do it. The function combines matplotlib's `scatter` functionality with [Cartopy], a powerful library offering different projections for geographical data.The `visualize_scatter` function takes the following kwargs:* `xr_dataarray` - xarray.DataArray object* `conversion_factor` - a number that can be used to convert data values* `projection` - choose a cartopy projection* `vmin` - Minimum value in the data range* `vmax` - Maximum value in the data range* `point_size` - size of the scatters* `color_scale` - string of matplotlib's color ramp references* `unit` - Unit of the data* `title` - Title of the plot<jupyter_code>visualize_scatter(xr_dataarray=aah_xr, 
             conversion_factor=1, 
             projection=ccrs.PlateCarree(),
             vmin=0,
             vmax=12,
             point_size=10,
             color_scale='viridis',
             unit=aah_xr.units, 
             title=aah_xr.long_name)<jupyter_output><empty_output><jupyter_text>Combine multiple HDF5 files in a `xarray.DataArray` object AC-SAF Metop-B GOME-2 Level 2 data files for one individual day are stored in multiple files. In order to combine all sensed ground tracks for one individual day, you can apply the function `load_l2_data_xr`, which loads multiple `HDF5` files and returns a `xarray.DataArray` object with all ground pixels of all files. The function makes internally use of the function `generate_xr_from_1D_vec`. Let us load all the AC-SAF Metop-B GOME-2 Level 2 data files for 9 September 2015. The function takes the following kwargs:* `directory` - file directory containing the `HDF5` files* `internal_filepath` - internal filepath containing the data information, e.g. `DATA`* `parameter` - name of the parameter* `lat_path` - path to the latitude information in the `HDF5` file* `lon_path` - path to the longitude information in the `HDF5` file* `no_of_dims` - Number of dimensions of the file* `paramname` - Name of variable* `unit` - Unit of the variable* `longname` - Long name of the variable<jupyter_code>directory='../eodata/gome2/aah/2015/09/09'
aah_da_mf = load_l2_data_xr(directory= directory,
                           internal_filepath='DATA',
                           parameter='AAH_AbsorbingAerosolHeight',
                           lat_path='GEOLOCATION/LatitudeCenter',
                           lon_path='GEOLOCATION/LongitudeCenter',
                           no_of_dims=2,
                           paramname=aah.name,
                           unit=aah.Unit,
                           longname=aah.Title)

aah_da_mf<jupyter_output><empty_output><jupyter_text>If you visualize the data array containing ground pixel information from multiple files, you see that for one day, almost all part on Earth are covered by Metop-B GOME-2. You can again use the function `visualize_scatter` to visualize the data as a scatter plot.<jupyter_code>visualize_scatter(xr_dataarray=aah_da_mf, 
             conversion_factor=1,
             projection=ccrs.Mollweide(), 
             vmin=0, 
             vmax=aah_da_mf.max(), 
             point_size=10, 
             color_scale='viridis', 
             unit=aah_da_mf.units, 
             title=str(aah_da_mf.long_name)[0:24])<jupyter_output><empty_output><jupyter_text>Filter AAH values based on Absorbing Aerosol Index (AAI) information Absorbing Aerosol Height information is strongly linked with Absorbing Aerosol Index values. AAI values below 2 index points correspond to scenes with too low amounts of aerosol to result in a reliable Absorbing Aerosol Height (AAH). These situations have already been flagged in the data file. However, even AAI values larger than 2 index points but smaller than 4 index points are often not thick enough for a reliable AAH retrieval. For this reason, we need to flag only the AAH values where the AAI values are larger than 4.As a first step, you have to load the `AAI` values, which are part of the data files and can be found as well in the `DATA` group. You can use the function `load_l2_data_xr` to load the AAI values for all files for the same day.<jupyter_code>aai = load_l2_data_xr(directory=directory,
                                    internal_filepath='DATA', 
                                    parameter='AAI', 
                                    lat_path='GEOLOCATION/LatitudeCenter', 
                                    lon_path='GEOLOCATION/LongitudeCenter', 
                                    no_of_dims=2, 
                                    paramname='AAI', 
                                    unit='-', 
                                    longname='Absorbing Aerosol Index')
aai<jupyter_output><empty_output><jupyter_text>In a next step, you can mask the AAH data array with the help of the previously loaded data array that contains the AAI values. You can use the function `generate_masked_array`. You want to filter out all AAH values, where AAI is less than 4 index points. As a result, only 592 AAH values remain valid.<jupyter_code>aah_da_masked = generate_masked_array(aah_da_mf, aai, 4, '>', drop=True)
aah_da_masked<jupyter_output><empty_output><jupyter_text>Visualize filtered AAH values In a last step, you can visualize again the data array with the filtered AAH information.<jupyter_code>visualize_scatter(xr_dataarray=aah_da_masked, 
             conversion_factor=1,
             projection=ccrs.Orthographic(120), 
             vmin=0, 
             vmax=10, 
             point_size=10, 
             color_scale='viridis', 
             unit=aah_da_masked.units[0:2], 
             title=str(aah_da_masked.long_name)[0:24] + ' on 2015-09-09'
                 )<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.07647717541663922, 'fraction_numerical': 0.020749621237072656, 'mean_word_length': 3.3155201819215465, 'pattern_counts': {'":': 0, '<': 58, '<?xml version=': 0, '>': 55, 'https://': 6, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 16, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '13', 'dir': 'jupyter-structured-clean-dedup', 'id': '377658', 'license': 'MIT', 'n_tokens_mistral': 4729, 'n_tokens_neox': 4426, 'n_words': 1556, 'path': '20_data_discovery/215_Metop-B_GOME-2_AAH_L2_load_browse.ipynb', 'repo_name': 'trivedi-c/atm_Practical3', 'size': '15181'}