lib.plotting_functions ====================== .. py:module:: lib.plotting_functions .. autoapi-nested-parse:: . Generic computation and plotting helper functions Functions --------- load_dataset() generalized load dataset method used for plotting/analysis functions use_this_norm() switches matplotlib color normalization method get_difference_colors(values) Provide a color norm and colormap assuming `values` is a difference field. mask_land_or_ocean(arr, msk, use_nan=False) Apply a land or ocean mask to provided variable. get_central_longitude(*args) Determine central longitude for maps. global_average(fld, wgt, verbose=False) pure numpy global average. spatial_average(indata, weights=None, spatial_dims=None) Compute spatial average wgt_rmse(fld1, fld2, wgt): Calculate the area-weighted RMSE. annual_mean(data, whole_years=False, time_name='time'): Calculate annual averages from time series data. seasonal_mean(data, season=None, is_climo=None): Calculates the time-weighted seasonal average (or average over all time). domain_stats(data, domain): Provides statistics in specified region. make_polar_plot(wks, case_nickname, base_nickname, case_climo_yrs, baseline_climo_yrs, d1:xr.DataArray, d2:xr.DataArray, difference:Optional[xr.DataArray]=None, domain:Optional[list]=None, hemisphere:Optional[str]=None, **kwargs): Make a stereographic polar plot for the given data and hemisphere. plot_map_vect_and_save(wks, case_nickname, base_nickname, case_climo_yrs, baseline_climo_yrs, plev, umdlfld_nowrap, vmdlfld_nowrap, uobsfld_nowrap, vobsfld_nowrap, udiffld_nowrap, vdiffld_nowrap, **kwargs): Plots a vector field on a map. plot_map_and_save(wks, case_nickname, base_nickname, case_climo_yrs, baseline_climo_yrs, mdlfld, obsfld, diffld, **kwargs): Map plots of `mdlfld`, `obsfld`, and their difference, `diffld`. pres_from_hybrid(psfc, hya, hyb, p0=100000.): Converts a hybrid level to a pressure vert_remap(x_mdl, p_mdl, plev) Interpolates to specified pressure levels. lev_to_plev(data, ps, hyam, hybm, P0=100000., new_levels=None, convert_to_mb=False) Interpolate model hybrid levels to specified pressure levels. pmid_to_plev(data, pmid, new_levels=None, convert_to_mb=False) Interpolate `data` from hybrid-sigma levels to isobaric levels using provided mid-level pressures. zonal_mean_xr(fld) Average over all dimensions except `lev` and `lat`. validate_dims(fld, list_of_dims) Checks if specified dimensions are in a DataArray lat_lon_validate_dims(fld) Check if input field has lat and lon. zm_validate_dims(fld) Check for dimensions for zonal average. zonal_plot(lat, data, ax=None, color=None, **kwargs) Make a line plot or pressure-latitude plot of `data`. meridional_plot(lon, data, ax=None, color=None, **kwargs) Make a line plot or pressure-longitude plot of `data`. prep_contour_plot Preparation for making contour plots. plot_zonal_mean_and_save zonal mean plot plot_meridional_mean_and_save meridioanl mean plot square_contour_difference Produce filled contours of fld1, fld2, and their difference with square axes. .. rubric:: Notes This module includes several "private" methods intended for internal use only. _plot_line(axobject, xdata, ydata, color, **kwargs) Create a generic line plot _meridional_plot_line _zonal_plot_line _zonal_plot_preslat _meridional_plot_preslon Attributes ---------- .. autoapisummary:: lib.plotting_functions.empty_message lib.plotting_functions.props lib.plotting_functions.seasons Functions --------- .. autoapisummary:: lib.plotting_functions.my_formatwarning lib.plotting_functions.load_dataset lib.plotting_functions.use_this_norm lib.plotting_functions.get_difference_colors lib.plotting_functions.mask_land_or_ocean lib.plotting_functions.get_central_longitude lib.plotting_functions.global_average lib.plotting_functions.spatial_average lib.plotting_functions.wgt_rmse lib.plotting_functions.annual_mean lib.plotting_functions.seasonal_mean lib.plotting_functions.domain_stats lib.plotting_functions.make_polar_plot lib.plotting_functions.plot_map_vect_and_save lib.plotting_functions.plot_map_and_save lib.plotting_functions.pres_from_hybrid lib.plotting_functions.vert_remap lib.plotting_functions.lev_to_plev lib.plotting_functions.pmid_to_plev lib.plotting_functions.zonal_mean_xr lib.plotting_functions.validate_dims lib.plotting_functions.lat_lon_validate_dims lib.plotting_functions.zm_validate_dims lib.plotting_functions._plot_line lib.plotting_functions._meridional_plot_line lib.plotting_functions._zonal_plot_line lib.plotting_functions._zonal_plot_preslat lib.plotting_functions._meridional_plot_preslon lib.plotting_functions.zonal_plot lib.plotting_functions.meridional_plot lib.plotting_functions.prep_contour_plot lib.plotting_functions.plot_zonal_mean_and_save lib.plotting_functions.plot_meridional_mean_and_save lib.plotting_functions.square_contour_difference Module Contents --------------- .. py:function:: my_formatwarning(msg, *args, **kwargs) Issue `msg` as warning. .. py:data:: empty_message :value: Multiline-String .. raw:: html
Show Value .. code-block:: python """No Valid Data Points""" .. raw:: html
.. py:data:: props .. py:data:: seasons .. py:function:: load_dataset(fils) This method exists to get an xarray Dataset from input file information that can be passed into the plotting methods. :param fils: strings or paths to input file(s) :type fils: list :rtype: xr.Dataset .. rubric:: Notes When just one entry is provided, use `open_dataset`, otherwise `open_mfdatset` .. py:function:: use_this_norm() Just use the right normalization; avoids a deprecation warning. .. py:function:: get_difference_colors(values) Provide a color norm and colormap assuming this is a difference field. :param values: can be either the data field or a set of specified contour levels. :type values: array-like :returns: * *dnorm* -- Matplotlib color nomalization * *cmap* -- Matplotlib colormap .. rubric:: Notes Uses 'OrRd' colormap for positive definite, 'BuPu_r' for negative definite, and 'RdBu_r' centered on zero if there are values of both signs. .. py:function:: mask_land_or_ocean(arr, msk, use_nan=False) Apply a land or ocean mask to provided variable. :param arr: the xarray variable to apply the mask to. :type arr: xarray.DataArray :param msk: the xarray variable that contains the land or ocean mask, assumed to be the same shape as "arr". :type msk: xarray.DataArray :param use_nan: argument for whether to set the missing values to np.nan values instead of the defaul "-999." values. :type use_nan: bool, optional :returns: **arr** -- Same as input `arr` but masked as specified. :rtype: xarray.DataArray .. py:function:: get_central_longitude(*args) Determine central longitude for maps. Allows an arbitrary number of arguments. If any of the arguments is an instance of `AdfDiag`, then check whether it has a `central_longitude` in `diag_basic_info`. _This case takes precedence._ _Else_, if any of the arguments are scalars in [-180, 360], assumes the FIRST ONE is the central longitude. There are no other possible conditions, so if none of those are met, returns the default value of 180. :param \*args: Any number of objects to check for `central_longitude`. After that, looks for the first number between -180 and 360 in the args. :type \*args: tuple .. rubric:: Notes This allows a script to, for example, allow a config file to specify, but also have a preference: `get_central_longitude(AdfObj, 30.0)` .. py:function:: global_average(fld, wgt, verbose=False) A simple, pure numpy global average. :param fld: an input ndarray :type fld: np.ndarray :param wgt: a 1-dimensional array of weights, should be same size as one dimension of `fld` :type wgt: np.ndarray :param verbose: prints information when `True` :type verbose: bool, optional :rtype: weighted average of `fld` .. py:function:: spatial_average(indata, weights=None, spatial_dims=None) Compute spatial average. :param indata: input data :type indata: xr.DataArray :param weights: the weights to apply, see Notes for default behavior :type weights: np.ndarray or xr.DataArray, optional :param spatial_dims: list of dimensions to average, see Notes for default behavior :type spatial_dims: list, optional :returns: weighted average of `indata` :rtype: xr.DataArray .. rubric:: Notes When `weights` is not provided, tries to find sensible values. If there is a 'lat' dimension, use `cos(lat)`. If there is a 'ncol' dimension, looks for `area` in `indata`. Otherwise, set to equal weights. Makes an attempt to identify the spatial variables when `spatial_dims` is None. Will average over `ncol` if present, and then will check for `lat` and `lon`. When none of those three are found, raise an AdfError. .. py:function:: wgt_rmse(fld1, fld2, wgt) Calculate the area-weighted RMSE. :param fld1: 2-dimensional spatial fields with the same shape. They can be xarray DataArray or numpy arrays. :type fld1: array-like :param fld2: 2-dimensional spatial fields with the same shape. They can be xarray DataArray or numpy arrays. :type fld2: array-like :param wgt: the weight vector, expected to be 1-dimensional, matching length of one dimension of the data. :type wgt: array-like :returns: * *float* -- root mean squared error * *Notes* * ```rmse = sqrt( mean( (fld1 - fld2)**2 ) )``` .. py:function:: annual_mean(data, whole_years=False, time_name='time') Calculate annual averages from monthly time series data. :param data: monthly data values with temporal dimension :type data: xr.DataArray or xr.Dataset :param whole_years: whether to restrict endpoints of the average to start at first January and end at last December :type whole_years: bool, optional :param time_name: name of the time dimension, defaults to `time` :type time_name: str, optional :returns: **result** -- `data` reduced to annual averages :rtype: xr.DataArray or xr.Dataset .. rubric:: Notes This function assumes monthly data, and weights the average by the number of days in each month. `result` includes an attribute that reports the date range used for the average. .. py:function:: seasonal_mean(data, season=None, is_climo=None) Calculates the time-weighted seasonal average (or average over all time). :param data: data to be averaged :type data: xarray.DataArray or xarray.Dataset :param season: the season to extract from `data` If season is `ANN` or None, average all available time. :type season: str, optional :param is_climo: If True, expects data to have time or month dimenion of size 12. If False, then 'time' must be a coordinate, and the `time.dt.days_in_month` attribute must be available. :type is_climo: bool, optional :returns: the average of `data` in season `season` :rtype: xarray.DataArray or xarray.Dataset .. rubric:: Notes If the data is a climatology, the code will make an attempt to understand the time or month dimension, but will assume that it is ordered from January to December. If the data is a climatology and is just a numpy array with one dimension that is size 12, it will assume that dimension is time running from January to December. .. py:function:: domain_stats(data, domain) Provides statistics in specified region. :param data: data values :type data: xarray.DataArray :param domain: the domain specification as: [west_longitude, east_longitude, south_latitude, north_latitude] :type domain: list or tuple or numpy.ndarray :returns: * **x_region_mean** (*float*) -- the regional area-weighted average * **x_region_max** (*float*) -- the maximum value in the region * **x_region_min** (*float*) -- the minimum value in the region .. rubric:: Notes Currently assumes 'lat' is a dimension and uses `cos(lat)` as weight. Should use `spatial_average` .. seealso:: :py:obj:`spatial_average` .. py:function:: make_polar_plot(wks, case_nickname, base_nickname, case_climo_yrs, baseline_climo_yrs, d1: xarray.DataArray, d2: xarray.DataArray, difference: Optional[xarray.DataArray] = None, pctchange: Optional[xarray.DataArray] = None, domain: Optional[list] = None, hemisphere: Optional[str] = None, obs=False, **kwargs) Make a stereographic polar plot for the given data and hemisphere. :param wks: output file path :type wks: str or Path :param case_nickname: short case name for `d1` :type case_nickname: str :param base_nickname: short case name for `d2` :type base_nickname: str :param case_climo_yrs: years for case `d1`, used for annotation :type case_climo_yrs: list :param baseline_climo_yrs: years for case `d2`, used for annotation :type baseline_climo_yrs: list :param d1: input data, must contain dimensions `lat` and `lon` :type d1: xr.DataArray :param d2: input data, must contain dimensions `lat` and `lon` :type d2: xr.DataArray :param difference: data to use as the difference, otherwise `d2 - d1` :type difference: xr.DataArray, optional :param pctchange: :type pctchange: xr.DataArray, optional data to use as the percent change :param domain: the domain to plot, specified as [west_longitude, east_longitude, south_latitude, north_latitude] Defaults to pole to 45-degrees, all longitudes :type domain: list, optional :param hemisphere: Hemsiphere to plot :type hemisphere: {'NH', 'SH'}, optional :param kwargs: variable-dependent options for plots, See Notes. :type kwargs: dict, optional .. rubric:: Notes - Uses contourf. No contour lines (yet). - kwargs is checked for: + `colormap` + `contour_levels` + `contour_levels_range` + `diff_contour_levels` + `diff_contour_range` + `diff_colormap` + `units` .. py:function:: plot_map_vect_and_save(wks, case_nickname, base_nickname, case_climo_yrs, baseline_climo_yrs, plev, umdlfld_nowrap, vmdlfld_nowrap, uobsfld_nowrap, vobsfld_nowrap, udiffld_nowrap, vdiffld_nowrap, obs=False, **kwargs) This plots a vector plot. Vector fields constructed from x- and y-components (u, v). :param wks: output file path :type wks: str or Path :param case_nickname: short name for case :type case_nickname: str :param base_nickname: short name for base case :type base_nickname: str :param case_climo_yrs: list of years in case climatology, used for annotation :type case_climo_yrs: list :param baseline_climo_yrs: list of years in base case climatology, used for annotation :type baseline_climo_yrs: list :param plev: if not None, label denoting the pressure level :type plev: str or float or None :param umdlfld_nowrap: input data for case, the x- and y- components of the vectors :type umdlfld_nowrap: xarray.DataArray :param vmdlfld_nowrap: input data for case, the x- and y- components of the vectors :type vmdlfld_nowrap: xarray.DataArray :param uobsfld_nowrap: input data for base case, the x- and y- components of the vectors :type uobsfld_nowrap: xarray.DataArray :param vobsfld_nowrap: input data for base case, the x- and y- components of the vectors :type vobsfld_nowrap: xarray.DataArray :param udiffld_nowrap: input difference data, the x- and y- components of the vectors :type udiffld_nowrap: xarray.DataArray :param vdiffld_nowrap: input difference data, the x- and y- components of the vectors :type vdiffld_nowrap: xarray.DataArray :param kwargs: variable-specific options, See Notes :type kwargs: dict, optional .. rubric:: Notes kwargs expected to be a variable-specific section, possibly provided by an ADF Variable Defaults YAML file. Currently it is inspected for: - `central_longitude` - `var_name` - `case_name` - `baseline` - `tiString` - `tiFontSize` - `units` _Note_ The title string constructed by kwargs appears to not be used. .. py:function:: plot_map_and_save(wks, case_nickname, base_nickname, case_climo_yrs, baseline_climo_yrs, mdlfld, obsfld, diffld, pctld, obs=False, **kwargs) Plots mdlfld, obsfld, diffld in a 3-row panel plot of maps. :param wks: Output file path. :type wks: str or Path :param case_nickname: Short name for case. :type case_nickname: str :param base_nickname: Short name for base case. :type base_nickname: str :param case_climo_yrs: List of years in case climatology, used for annotation. :type case_climo_yrs: list :param baseline_climo_yrs: List of years in base case climatology, used for annotation. :type baseline_climo_yrs: list :param mdlfld: Input data for case. :type mdlfld: xarray.DataArray :param obsfld: Input data for base case. :type obsfld: xarray.DataArray :param diffld: Input difference data. :type diffld: xarray.DataArray :param pctld: Input percent difference data. :type pctld: xarray.DataArray :param obs: If True, use observational data. Default is False. :type obs: bool, optional :param \*\*kwargs: Variable-specific options. See Notes. :type \*\*kwargs: dict, optional .. rubric:: Notes The `kwargs` parameter is expected to be a variable-specific section, possibly provided by an ADF Variable Defaults YAML file. Currently it is inspected for: - colormap (str): Name of matplotlib colormap. - contour_levels (list or tuple): Explicit values or a tuple (min, max, step). - diff_colormap - diff_contour_levels - tiString (str): Title string. - tiFontSize (int): Title font size. - mpl (dict): Any matplotlib kwargs that should be passed along. Example YAML for `mpl` options:: mpl: subplots: figsize: (3, 9) contourf: levels: 15 cmap: Blues colorbar: shrink: 0.4 This is experimental, and if you find yourself doing much with this, you probably should write a new plotting script that does not rely on this module. When these are not provided, colormap is set to 'coolwarm' and limits/levels are set by data range. .. py:function:: pres_from_hybrid(psfc, hya, hyb, p0=100000.0) Calculates pressure field pressure derived with the formula: ```p = a(k)*p0 + b(k)*ps``` :param psfc: surface pressure :param hya: hybrid-sigma A and B coefficients :param hyb: hybrid-sigma A and B coefficients :param p0: reference pressure, defaults to 100000 Pa :type p0: optional :rtype: pressure, size is same as `psfc` with `len(hya)` levels .. py:function:: vert_remap(x_mdl, p_mdl, plev) Apply simple 1-d interpolation to a field :param x_mdl: input data :type x_mdl: xarray.DataArray or numpy.ndarray :param p_mdl: pressure field, same shape as `x_mdl` :type p_mdl: xarray.DataArray or numpy.ndarray :param plev: the new pressures :type plev: xarray.DataArray or numpy.ndarray :returns: `x_mdl` interpolated to `plev` :rtype: output .. rubric:: Notes Interpolation done in log pressure .. py:function:: lev_to_plev(data, ps, hyam, hybm, P0=100000.0, new_levels=None, convert_to_mb=False) Interpolate model hybrid levels to specified pressure levels. :param data: :param ps: surface pressure :param hyam: hybrid-sigma A and B coefficients :param hybm: hybrid-sigma A and B coefficients :param P0: reference pressure, defaults to 100000 Pa :type P0: float, optional :param new_levels: 1-D array containing pressure levels in Pascals (Pa). If not specified, then the levels will be set to the GeoCAT defaults, which are (in hPa): `1000, 925, 850, 700, 500, 400, 300, 250, 200, 150, 100, 70, 50, 30, 20, 10, 7, 5, 3, 2, 1` :type new_levels: numpy.ndarray, optional :param convert_to_mb: If True, then vertical (lev) dimension will have values of mb/hPa, otherwise the units are Pa. :type convert_to_mb: bool, optional :returns: data interpolated to new pressure levels :rtype: data_interp_rename .. rubric:: Notes The function `interp_hybrid_to_pressure` used here is dask-enabled, and so can potentially be sped-up via the use of a DASK cluster. .. py:function:: pmid_to_plev(data, pmid, new_levels=None, convert_to_mb=False) Interpolate data from hybrid-sigma levels to isobaric levels. :param data: field with a 'lev' coordinate :type data: xarray.DataArray :param pmid: the pressure field (Pa), same shape as `data` :type pmid: xarray.DataArray :param new_levels: the output pressure levels (Pa), defaults to standard levels :type new_levels: optional :param convert_to_mb: flag to convert output to mb (i.e., hPa), defaults to False :type convert_to_mb: bool, optional :returns: **output** -- `data` interpolated onto `new_levels` :rtype: xarray.DataArray .. py:function:: zonal_mean_xr(fld) Average over all dimensions except `lev` and `lat`. .. py:function:: validate_dims(fld, list_of_dims) Check if specified dimensions are in a DataArray. :param fld: field to check for named dimensions :type fld: xarray.DataArray :param list_of_dims: list of strings that specifiy the dimensions to check for :type list_of_dims: list :returns: dict with keys that are "has_{x}" where x is the name from `list_of_dims` and values that are boolean :rtype: dict .. py:function:: lat_lon_validate_dims(fld) Check if input field has lat and lon. :param fld: data with named dimensions :type fld: xarray.DataArray :returns: True if lat and lon are both dimensions, False otherwise. :rtype: bool .. seealso:: :py:obj:`validate_dims` .. py:function:: zm_validate_dims(fld) Check for dimensions for zonal average. Looks for dimensions called 'lev' and 'lat'. :param fld: field to check for lat and/or lev dimensions :type fld: xarray.DataArray :returns: * *tuple* -- (has_lat, has_lev) each are bool * *None* -- If 'lat' is not in dimensions, returns None. .. py:function:: _plot_line(axobject, xdata, ydata, color, **kwargs) Create a generic line plot and check for some ways to annotate. .. py:function:: _meridional_plot_line(ax, lon, data, color, **kwargs) Create line plot with longitude as the X-axis. .. py:function:: _zonal_plot_line(ax, lat, data, color, **kwargs) Create line plot with latitude as the X-axis. .. py:function:: _zonal_plot_preslat(ax, lat, lev, data, **kwargs) Create plot with latitude as the X-axis, and pressure as the Y-axis. .. py:function:: _meridional_plot_preslon(ax, lon, lev, data, **kwargs) Create plot with longitude as the X-axis, and pressure as the Y-axis. .. py:function:: zonal_plot(lat, data, ax=None, color=None, **kwargs) Make zonal plot Determine which kind of zonal plot is needed based on the input variable's dimensions. :param lat: latitude :param data: input data :param ax: axes object to use :type ax: Axes, optional :param color: color for the curve :type color: str or mpl color specification :param kwargs: plotting options :type kwargs: dict, optional .. rubric:: Notes Checks if there is a `lev` dimension to determine if it is a lat-pres plot or a line plot. .. py:function:: meridional_plot(lon, data, ax=None, color=None, **kwargs) Make meridional plot Determine which kind of meridional plot is needed based on the input variable's dimensions. :param lon: longitude :param data: input data :param ax: axes object to use :type ax: Axes, optional :param color: color for the curve :type color: str or mpl color specification :param kwargs: plotting options :type kwargs: dict, optional .. rubric:: Notes Checks if there is a `lev` dimension to determine if it is a lon-pres plot or a line plot. .. py:function:: prep_contour_plot(adata, bdata, diffdata, pctdata, **kwargs) Preparation for making contour plots. Prepares for making contour plots of adata, bdata, diffdata, and pctdata, which is presumably the difference between adata and bdata. - set colormap from kwargs or defaults to coolwarm - set contour levels from kwargs or 12 evenly spaced levels to span the data - normalize colors based on specified contour levels or data range - set option for linear or log pressure when applicable - similar settings for difference, defaults to symmetric about zero - separates Matplotlib kwargs into their own dicts :param adata: the data to be plotted :param bdata: the data to be plotted :param diffdata: the data to be plotted :param pctdata: the data to be plotted :param kwargs: plotting options :type kwargs: dict, optional :returns: a dict with the following: - 'subplots_opt': mpl kwargs for subplots - 'contourf_opt': mpl kwargs for contourf - 'colorbar_opt': mpl kwargs for colorbar - 'diff_colorbar_opt' : mpl kwargs for difference colorbar - 'normdiff': color normalization for difference panel - 'cmapdiff': colormap for difference panel - 'levelsdiff': contour levels for difference panel - 'cmap1': color map for a and b panels - 'norm1': color normalization for a and b panels - 'levels1' : contour levels for a and b panels - 'plot_log_p' : true/false whether to plot log(pressure) axis :rtype: dict .. py:function:: plot_zonal_mean_and_save(wks, case_nickname, base_nickname, case_climo_yrs, baseline_climo_yrs, adata, bdata, has_lev, log_p, obs=False, **kwargs) This is the default zonal mean plot :param adata: The vertical coordinate (lev) must be pressure levels. :type adata: data to plot ([lev], lat, [lon]). :param bdata: - For 2-d variables (reduced to (lat,)): + 2 panels: (top) zonal mean, (bottom) difference - For 3-D variables (reduced to (lev,lat)): + 3 panels: (top) zonal mean adata, (middle) zonal mean bdata, (bottom) difference + pcolormesh/contour plot :type bdata: baseline or observations to plot adata against. :param kwargs -> optional dictionary of plotting options: ** Expecting this to be variable-specific section, possibly provided by ADF Variable Defaults YAML file.** :param - colormap -> str: :param name of matplotlib colormap: :param - contour_levels -> list of explict values or a tuple: :type - contour_levels -> list of explict values or a tuple: (min, max, step) :param - diff_colormap: :param - diff_contour_levels: :param - tiString -> str: :param Title String: :param - tiFontSize -> int: :param Title Font Size: :param - mpl -> dict: + Organize these by the mpl function. In this function (`plot_map_and_save`) we will check for an entry called `subplots`, `contourf`, and `colorbar`. So the YAML might looks something like: ``` mpl: subplots: figsize: (3, 9) contourf: levels: 15 cmap: Blues colorbar: shrink: 0.4 ``` :param This should be any matplotlib kwargs that should be passed along. Keep reading: + Organize these by the mpl function. In this function (`plot_map_and_save`) we will check for an entry called `subplots`, `contourf`, and `colorbar`. So the YAML might looks something like: ``` mpl: subplots: figsize: (3, 9) contourf: levels: 15 cmap: Blues colorbar: shrink: 0.4 ``` .. py:function:: plot_meridional_mean_and_save(wks, case_nickname, base_nickname, case_climo_yrs, baseline_climo_yrs, adata, bdata, has_lev, latbounds=None, obs=False, **kwargs) Default meridional mean plot :param wks: the figure object to plot in :param case_nickname: short name of `adata` case, use for annotation :type case_nickname: str :param base_nickname: short name of `bdata` case, use for annotation :type base_nickname: str :param case_climo_yrs: years in the `adata` case, use for annotation :type case_climo_yrs: list :param baseline_climo_yrs: years in the `bdata` case, use for annotation :type baseline_climo_yrs: list: :param adata: data to plot ([lev], [lat], lon). The vertical coordinate (lev) must be pressure levels. :type adata: xarray.DataArray :param bdata: baseline or observations to plot adata against. It must have the same dimensions and vertical levels as adata. :type bdata: xarray.DataArray :param has_lev: whether lev dimension is present :type has_lev: bool :param latbounds: indicates latitude bounds to average over if it is a number, assume symmetric about equator, otherwise expects `slice(south, north)` defaults to `slice(-5,5)` :type latbounds: numbers.Number or slice, optional :param kwargs: optional dictionary of plotting options, See Notes :type kwargs: dict, optional .. rubric:: Notes - For 2-d variables (reduced to (lon,)): + 2 panels: (top) meridional mean, (bottom) difference - For 3-D variables (reduced to (lev,lon)): + 3 panels: (top) meridonal mean adata, (middle) meridional mean bdata, (bottom) difference + pcolormesh/contour plot - kwargs -> optional dictionary of plotting options ** Expecting this to be variable-specific section, possibly provided by ADF Variable Defaults YAML file.** - colormap -> str, name of matplotlib colormap - contour_levels -> list of explicit values or a tuple: (min, max, step) - diff_colormap -> str, name of matplotlib colormap used for different plot - diff_contour_levels -> list of explicit values or a tuple: (min, max, step) - tiString -> str, Title String - tiFontSize -> int, Title Font Size - mpl -> dict, This should be any matplotlib kwargs that should be passed along. Keep reading: + Organize these by the mpl function. In this function (`plot_meridional_mean_and_save`) we will check for an entry called `subplots`, `contourf`, and `colorbar`. So the YAML might looks something like: ``` mpl: subplots: figsize: (3, 9) contourf: levels: 15 cmap: Blues colorbar: shrink: 0.4 ``` .. py:function:: square_contour_difference(fld1, fld2, **kwargs) Produce filled contours of fld1, fld2, and their difference with square axes. Intended use is latitude-by-month to show the annual cycle. Example use case: use climo files to get data, take zonal averages, rename "time" to "month" if desired, and then provide resulting DataArrays to this function. :param fld1: 2-dimensional DataArrays with same shape :type fld1: xarray.DataArray :param fld2: 2-dimensional DataArrays with same shape :type fld2: xarray.DataArray :param \*\*kwargs: optional keyword arguments this function _only checks_ `kwargs` for `case1name`, `case2name` :type \*\*kwargs: dict, optional :returns: figure object :rtype: fig .. rubric:: Notes Assumes `fld1.shape == fld2.shape` and `len(fld1.shape) == 2` Will try to label the cases by looking for `case1name` and `case2name` in `kwargs`, and then `fld1['case']` and `fld2['case']` (i.e., attributes) If no case name is found proceeds with empty strings. **IF THERE IS A BETTER CONVENTION WE SHOULD USE IT.** Each panel also puts the Min/Max values into the title string. Axis labels are upper-cased names of the coordinates of `fld1`. Ticks are automatic with the exception that if the first coordinate is "month" and is length 12, use `np.arange(1,13)`.