# Application - Introduction

The `readDiag` package was designed to be a tool that enables easy access to diagnostic files generated by the [Gridpoint Statistical Interpolation (GSI) system](https://github.com/NOAA-EMC/GSI). The package is mainly composed of 2 classes, the first one for reading data (`read_diag`) and the other for creating different types of figures with selected data (`plot_diag`).

## The `read_diag` Class

The `read_diag` class consists of 5 functions:

1. `__init__(self, diagFile, diagFileAnl=None, isisList=None, zlevs=None)`: where `diagFile` is the diagnostic file of the first outer loop with OmF information, and `diagFileAnl` is the diagnostic file of the last outer loop with OmA information. Note that it is not necessary to inform both files; `diagFileAnl` is optional. In this case, the read information (Omf or OmA) will depend on the outer loop (diagnostic file) provided.
2. `overview(self)`: this function creates and returns a dictionary with the information in the file.
3. `pfileinfo(self)`: this function elegantly prints a list with the information in the file.
4. `close(self)`: this function closes the last opened file.
5. `tocsv(self, varName=None, varType=None, dateIni=None, dateFin=None, nHour="06", Level=None, Lay=None, SingleL=None)`: this function generates a CSV file of the OmF and OmA parameters with the following information: date, mean, standard deviation, and total data for the chosen variable and type.

The input data for the function and functionalities are identical to those of the `time_series` function presented in the next section in the `plot_diag` class (item 7 - see also details of functionalities [here](#time_series)). What distinguishes the two functions is that this one writes a CSV file, and the other generates time series figures.

## The `plot_diag` Class

The `plot_diag` class consists of 7 functions:

1. `plot(self, varName, varType, param, mask=None, **kwargs)`: the `plot` function generates a figure for the variable `varName` (e.g., `uv`), `varType` (e.g., `220` (dropssonda)), and `param`, which can be various options such as `param="obs"` for the observation value, `param="omf"` for observation minus background, or `param="oma"` for observation minus analysis. It is also possible to mask the data with the `iuse` variable, which indicates whether the data was (`iuse=1`) or was not (`iuse=-1` - monitored data) used in assimilation. Below is an example executed for `varName="uv"`, `varType=220`, `param="obs"`, and mask `iuse==1"`:

    ![png](../imgs/uv_obs_plot.png)

2. `ptmap(self, varName, varType=None, mask=None, **kwargs)`: the `ptmap` function generates a figure with the location of all observations defined by `varName` (e.g., `uv`) and `varType` (can be a single type or a list, e.g., `[200]` or `[220,221,257]`). If `varType` is not informed, then all types will be included in the figure. It is also possible to mask the data with the `iuse` variable, which indicates whether the data was (`iuse=1`) or was not (`iuse=-1` - monitored data) used in assimilation. Below is an example executed for `varName="uv"`, `varType=[254,242,221,220,257,258,281,280]`, and `mask=None`;

    ![png](../imgs/uv_254_242_221_220_257_258_281_280_ptmap.png)

3. `pvmap(self, varName=None, mask=None, **kwargs)`: the `pvmap` function is similar to the `ptmap` function, with the difference of not specifying the type (`varType`), so you can choose a list of variables, e.g., `["uv","ps","t"]`. Below is an example executed for `varName=['uv','ps','t']`, and `mask=None`;

    ![png](../imgs/uv_ps_t_pvmap.png)

4. `pcount(self,varName,**kwargs)`: the `pcount` function generates a histogram of the data quantity for a specific variable `varName` (e.g., `uv`) and all available types (`varType`);

    ![png](../imgs/hist_uv_kx.png)

5. `kxcount(self,**kwargs)`: the `kxcount` function is similar to `pcount` but does not specify the variable (`varName`). This function generates a histogram with the total data (all variables summed) for all available types (`varType`);

    ![png](../imgs/hist_tot_obs_kx.png)

6. `vcount(self,**kwargs)`: the `vcount` function generates a histogram with the total quantity of data for each variable (`ps`, `t`, `q`, `uv`);

    ![png](../imgs/hist_tot_obs.png)

7. `time_series(self, varName=None, varType=None, dateIni=None, dateFin=None, nHour="06", vminOMA=None, vmaxOMA=None, vminSTD=0.0, vmaxSTD=14.0, Level=None, Lay=None, SingleL=None, Clean=None)`: the `time_series` function can generate 6 different types of figures, depending on the configuration specified in its call. The common feature among the 6 types is the time variation, while the difference between the 6 types is the way to treat the vertical levels. Basically, these ways are distributed between varying vertically (different values for different levels/layers) and fixed vertically (specific level, layer mean, or the entire atmosphere).

The input parameters for the function and each of the figure options will be explained in more detail below.

| Parameter       | Example                                  | Description
| :---            | :---:                                    | :---
| `self`          | `['/home/user/diag_conv_01.2019121000']` | List with all full paths (`path/file_name`) of each time in the time series.
| `varName`       | `uv`                                     | Variable name
| `varType`       | `220`                                    | Variable type
| `dateIni`       | `2019121000`                             | Initial date of the time series
| `dateFin`       | `2019121118`                             | Final date of the time series
| `nHour`         | `6`                                      | Time interval in hours between each file in the `self` list
| `vminOMA`       | `-2.0`                                   | Minimum value of the y-scale (ordinate) for OmF and OmA
| `vmaxOMA`       | `2.0`                                    | Maximum value of the y-scale (ordinate) for OmF and OmA
| `vminSTD`       | `0.0`                                    | Minimum value of the y-scale (ordinate) for the standard deviation
| `vmaxSTD`       | `14.0`                                   | Maximum value of the y-scale (ordinate) for the standard deviation
| `Level`         | `Zlevs`                                  | Value of the level to make the time series; options: numerical value corresponding to the level, e.g., 1000 for 1000 hPa; `Zlevs` to plot by layers (around the standard levels); `None` to plot all levels.
| `Lay`           | 25                                       | Half of the layer size (if `Level="Zlevs"`) in hPa if opting for layer sampling. If `Lay=None`, `Lay` will be internally calculated to fill the entire atmosphere containing the standard levels.
| `SingleL`       | `All`                                    | When `Level` is fixed, e.g., 1000 hPa, it will be considered exactly that level (using the option `SingleL=None`) or at all levels as a single layer (using `SingleL="All"`) or in a layer defined around the `Level` value ranging between `Level-Lay` and `Level+Lay`. If `Lay` is not informed, a default value of 50 hPa will be used.
| `Clean`         | `True` or `False`                        | If `True`, after generating and saving the figure, the figure window is restarted (`plt.clf()`) or closed (`plt.close()`); if `False`, this procedure is eliminated, and the figure remains available for viewing with `plt.show()`.

All figures generated with the option `Level` equal to `None` or `Zlevs` will contain the term `all_levels` in the name; otherwise, it will be `level` or `layer`, depending on the `SingleL` option.

Throughout this notebook, examples with code snippets are shown to illustrate the use of the functions listed above.