diff --git a/docs/source/Benchmarking.rst b/docs/source/Benchmarking.rst new file mode 100644 index 00000000..1d6b2f7a --- /dev/null +++ b/docs/source/Benchmarking.rst @@ -0,0 +1,1553 @@ +.. |br| raw:: html + +
+ +.. _bmk: + +############ +Benchmarking +############ + +The GEOS-Chem Support Team uses GCPy to produce comparison plots and +summary tables from GEOS-Chem benchmark simulations. In this chapter +we will describe this capability of GCPy. + +.. _bmk-scripts: + +====================================== +Location of benchmark plotting scripts +====================================== + +The source code for creating benchmark plots is located in the +:file:`gcpy/benchmark` directory tree. + +.. table:: **Contents of the gcpy/benchmark directory** + + +-------------------------+--------------------------------------------+ + | File or folder | Description | + +=========================+============================================+ + | ``run_benchmark.py`` | Benchmark driver script | + +-------------------------+--------------------------------------------+ + | ``benchmark_slurm.sh`` | Bash script to submit ``run_benchmark,py`` | + | | as a SLURM batch job | + +-------------------------+--------------------------------------------+ + | ``cloud/`` | Directory containing template config files | + | | (in YAML format) for 1-hour and 1-month | + | | benchmark plot jobs on the AWS cloud. | + +-------------------------+--------------------------------------------+ + | ``config/`` | Directory containing editable config files | + | | (in YAML format) for 1-month and 1-year | + | | benchmark plot jobs. | + +-------------------------+--------------------------------------------+ + | ``__init__.py`` | Python import script | + +-------------------------+--------------------------------------------+ + | ``modules/`` | Contains Python modules imported into the | + | | ``run_benchmark.py`` script. | + +-------------------------+--------------------------------------------+ + | ``README.md`` | Readme file in Markdown format | + +-------------------------+--------------------------------------------+ + +.. note:: + + As of this writing, the benchmarking scripts still use several + :ref:`plotting ` and :ref:`tabling + ` functions from module + :file:`gcpy.benchmark_funcs`. We are currently in the process of + moving the functions contained in :file:`gcpy.benchmark_funcs` to + the :file:`gcpy/benchmark/modules` directory. + +.. _bmk-steps: + +=============================== +Steps to create benchmark plots +=============================== + +Follow these instructions to create comparison plots and summary +tables from GEOS-Chem benchmark simulations. + +#. Copy a configuration file from the :file:`gcpy/benchmark/config` + directory. + + In this example we will use the configuration file that will create + plots from 1-year full-chemistry benchmark + simulations. (Configuration files for other benchmark types have a + similar layout.) + + .. code-block:: console + + $ cp /path/to/GCPy/gcpy/benchmark/config/1yr_fullchem_benchmark.yml . + + |br| + +#. Edit the :literal:`paths` section of the configuration file to + specify the proper directory paths for your system. + + .. code-block:: yaml + + # Configuration for 1-year FullChemBenchmark + # + # paths: + # main_dir: High-level directory containing ref & dev rundirs + # results_dir: Directory where plots/tables will be created + # weights_dir: Path to regridding weights + # spcdb_dir: Folder in which the species_database.yml file is + # located. If set to "default", then will look for + # species_database.yml in one of the Dev rundirs. + # obs_data_dir: Path to observational data (for models vs obs plots) + # + paths: + main_dir: /path/to/benchmark/main/dir # EDIT AS NEEDED + results_dir: /path/to/BenchmarkResults # EDIT AS NEEDED + weights_dir: /n/holyscratch01/external_repos/GEOS-CHEM/gcgrid/data/ExtData/GCHP/RegriddingWeights + spcdb_dir: default + obs_data_dir: /path/to/observational/data + + |br| + +#. Edit the :literal:`data` section to specify the directories (and + labels) for the Ref and Dev versions for GEOS-Chem Classic and GCHP. + + .. code-block:: yaml + + # + # data: Contains configurations for ref and dev runs + # version: Version string (must not contain spaces) + # dir: Path to run directory + # outputs_subdir: Subdirectory w/ GEOS-Chem diagnostic files + # restarts_subdir: Subdirectory w/ GEOS-Chem restarts + # bmk_start: Simulation start date (YYYY-MM-DDThh:mm:ss) + # bmk_end: Simulation end date (YYYY-MM-DDThh:mm:ss) + # resolution: GCHP resolution string + # + data: + ref: + gcc: + version: GCC_ref + dir: GCC_ref + outputs_subdir: OutputDir + restarts_subdir: Restarts + bmk_start: "2019-01-01T00:00:00" + bmk_end: "2020-01-01T00:00:00" + gchp: + version: GCC_dev + dir: GCC_dev + outputs_subdir: OutputDir + restarts_subdir: Restarts + bmk_start: "2019-01-01T00:00:00" + bmk_end: "2020-01-01T00:00:00" + is_pre_14.0: False + resolution: c24 + dev: + gcc: + version: GCC_dev + dir: GCC_dev + outputs_subdir: OutputDir + restarts_subdir: Restarts + bmk_start: "2019-01-01T00:00:00" + bmk_end: "2020-01-01T00:00:00" + gchp: + version: GCC_dev + dir: GCC_dev + restarts_subdir: Restarts + bmk_start: "2019-01-01T00:00:00" + bmk_end: "2020-01-01T00:00:00" + is_pre_14.0: False + resolution: c24 + + |br| + +#. Edit the :literal:`comparisons` section to specify the types of + comparisons you would like to perform. + + .. code-block:: yaml + + # + # comparisons: Specifies the comparisons to perform. + # + comparisons: + gcc_vs_gcc: + run: True + dir: GCC_version_comparison + tables_subdir: Tables + gchp_vs_gcc: + run: True + dir: GCHP_GCC_comparison + tables_subdir: Tables + gchp_vs_gchp: + run: True + dir: GCHP_version_comparison + tables_subdir: Tables + gchp_vs_gcc_diff_of_diffs: + run: True + dir: GCHP_GCC_diff_of_diffs + +#. Edit the :literal:`outputs` section to select the plots and tables + that you would like to generate. + + .. code-block:: yaml + + # + # outputs: Specifies the plots and tables to generate + # + outputs: + plot_conc: True + plot_emis: True + emis_table: True + plot_jvalues: True + plot_aod: True + mass_table: True + ops_budget_table: False + aer_budget_table: True + Ox_budget_table: True + ste_table: True # GCC only + OH_metrics: True + plot_models_vs_obs: True + plot_options: + by_spc_cat: True + by_hco_cat: True + + |br| + +#. Edit the :literal:`n_cores` setting if you wish to change the + number of computational cores to use. If not, leave + :literal:`n_cores` set to :literal:`-1`, which will use as many + cores as possible. + + .. code-block:: yaml + + # + # n_cores: Specify the number of cores to use. + # -1: Use $OMP_NUM_THREADS cores + # -2: Use $OMP_NUM_THREADS - 1 cores + # -N: Use $OMP_NUM_THREADS - (N-1) cores + # 1: Disable parallelization (use a single core) + # + n_cores: -1 + + |br| + +#. Run the :file:`run.benchmark.py` script. You may do this in 2 + ways: + + #. Direct execution from the command line: + + .. code-block:: console + + (gcpy_env) $ python -m gcpy.benchmark.run_benchmark + 1yr_fullchem_benchmark.yml + + #. Batch execution with the SLURM scheduler. First, copy the + :file:`benchmark_slurm.sh` script to your current directory: + + .. code-block:: console + + (gcpy_env) $ cp /path/to/GCPy/gcpy/benchmark/benchmark_slurm.sh . + + Next, edit your local copy of :file:`benchmark_slurm.sh` to + specify your SLURM partition name, number of cores, the name of + your Python environment and the configuration file to use. + + .. code-block:: bash + + #!/bin/bash + + #SBATCH -c 8 + #SBATCH -N 1 + #SBATCH -t 0-4:00 + #SBATCH -p seas_compute,shared + #SBATCH --mem=100000 + #SBATCH --mail-type=END + + #============================================================================ + # This us a sample SLURM script that you can use to run the GCPy + # benchmark plotting code as a SLURM batch job. + # + # You can modify the SLURM parameters above for your setup. + # + # Tip: Using less cores can reduce the amount of memory required. + #============================================================================ + + # Apply all bash initialization settings + . ~/.bashrc + + # Make sure to set multiple threads; Joblib will use multiple + # cores to parallelize certain plotting operations. + export OMP_NUM_THREADS=$SLURM_CPUS_PER_TASK + export OMP_STACKSIZE=500m + + # Turn on Python environment (edit for your setup) + mamba activate gcpy_env + + # Specify a YAML file with benchmark options + # Uncomment the file that you wish: + #config="1mo_benchmark.yml" + config="1yr_fullchem_benchmark.yml" + #config="1yr_tt_benchmark.yml" + + # Call the run_benchmark script to make the plots + python -m gcpy.benchmark.run_benchmark "${config}" > benchmark.log 2>&1 + + # Turn off python environment + mamba deactivate + + exit 0 + + Lastly, start the SLURM batch execution with this command: + + .. code-block:: console + + $ sbatch benchmark_slurm.sh + +.. _bmk-funcs-plot: + +============================ +Benchmark plotting functions +============================ + +Module :code:`gcpy.benchmark_funcs` contains several functions for +creating plots and tables from GEOS-Chem benchmark simulations. The +specific outputs generated have been requested by the `GEOS-Chem +Steering Committee `_ in +order to facilitate comparing benchmark output from different model +versions. + +In this section, we will describe functions that create comparison +plots from GEOS-Chem benchmark simulation output. The functions to +create summary tables will be described :ref:`in a separate section +`. + +.. note:: + + We are working towards moving all benchmark-related source code to + the :file:`gcpy/benchmark/` directory tree. For the time being, + the :file:`benchmark_funcs.py` script is located in the + :file:`/path/to/GCPy/gcpy/` directory. + +.. table:: **Functions creating comparison plots from benchmark + simulation output** + + +-----------------------------------------------+----------------------------------------------+ + | Function | Type of 6-panel comparison plot created | + +===============================================+==============================================+ + | ``make_benchmark_aod_plots()`` | Comparison plots for aerosol optical depth | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_conc_plots()`` | Species concentration | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_emis_plots()`` | Emissions (by species and catgegory) | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_jvalue_plots()`` | Comparison plots for J-values (photolysis) | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_wetdep_plots()`` | Comparison plots for species wet deposition | + +-----------------------------------------------+----------------------------------------------+ + +The functions listed above create comparison plots of most GEOS-Chem +output variables divided into specific categories, e.g. species +categories such as :literal:`Aerosols` or :literal:`Bromine` for the +:literal:`SpeciesConcVV` diagnostic. In eachcategory, these function +create single level PDFs for the surface and 500hPa and zonal +mean PDFs for the entire atmosphere and only the stratosphere (defined +a 1-100hPa). For :code:`make_benchmark_emis_plots()`, only single +level plots at the surface are produced. All of these plotting +functions include bookmarks within the generated PDFs that point to +the pages containing each plotted quantity. Thus these functions serve +as tools for quickly creating comprehensive plots comparing two +GEOS-Chem runs. These functions are used to create the publicly +available plots for 1-month and 1-year benchmarks of new versions of +GEOS-Chem. + +Many of the plotting functions listed above use pre-defined lists of +variables in YAML files. If one dataset includes a variable but the +other dataset does not, the data for that variable in the latter +dataset will be considered to be NaN and will be plotted as such. + +.. _bmk-funcs-plot-aod: + +make_benchmark_aod_plots +------------------------ + +This function creates column optical depth plots using the Aerosols +diagnostic output. + +.. code-block:: python + + def make_benchmark_aod_plots( + ref, + refstr, + dev, + devstr, + varlist=None, + dst="./benchmark", + subdst=None, + cmpres=None, + overwrite=False, + verbose=False, + log_color_scale=False, + sigdiff_files=None, + weightsdir='.', + n_job=-1, + time_mean=False, + spcdb_dir=os.path.dirname(__file__) + ): + """ + Creates PDF files containing plots of column aerosol optical + depths (AODs) for model benchmarking purposes. + + Args: + ref: str + Path name for the "Ref" (aka "Reference") data set. + refstr: str + A string to describe ref (e.g. version number) + dev: str + Path name for the "Dev" (aka "Development") data set. + This data set will be compared against the "Reference" + data set. + devstr: str + A string to describe dev (e.g. version number) + + Keyword Args (optional): + varlist: list of str + List of AOD variables to plot. If not passed, then all + AOD variables common to both Dev and Ref will be plotted. + Use the varlist argument to restrict the number of + variables plotted to the pdf file when debugging. + Default value: None + dst: str + A string denoting the destination folder where a + PDF file containing plots will be written. + Default value: ./benchmark. + subdst: str + A string denoting the sub-directory of dst where PDF + files containing plots will be written. In practice, + subdst is only needed for the 1-year benchmark output, + and denotes a date string (such as "Jan2016") that + corresponds to the month that is being plotted. + Default value: None + cmpres: string + Grid resolution at which to compare ref and dev data, e.g. '1x1.25' + overwrite: bool + Set this flag to True to overwrite files in the + destination folder (specified by the dst argument). + Default value: False. + verbose: bool + Set this flag to True to print extra informational output. + Default value: False + log_color_scale: bool + Set this flag to True to enable plotting data (not diffs) + on a log color scale. + Default value: False + sigdiff_files: list of str + Filenames that will contain the list of quantities having + having significant differences in the column AOD plots. + These lists are needed in order to fill out the benchmark + approval forms. + Default value: None + weightsdir: str + Directory in which to place (and possibly reuse) xESMF regridder + netCDF files. + Default value: '.' + n_job: int + Defines the number of simultaneous workers for parallel plotting. + Set to 1 to disable parallel plotting. Value of -1 allows the + application to decide. + Default value: -1 + spcdb_dir: str + Directory of species_datbase.yml file + Default value: Directory of GCPy code repository + time_mean : bool + Determines if we should average the datasets over time + Default value: False + """ + +.. _bmk-funcs-plot-conc: + +make_benchmark_conc_plots +------------------------- + +This function creates species concentration plots using the +SpeciesConc diagnostic output by default. In particular: + +- This function is the only benchmark plotting function that supports + diff-of-diffs plotting, in which 4 datasets are passed and the + differences between two groups of :literal:`Ref` datasets vs. two + groups of :literal:`Dev` datasets is plotted (typically used for + comparing changes in GCHP vs. changes in GEOS-Chem Classic across + model versions). |br| + |br| + +- This is also the only benchmark plotting function that sends plots + to separate folders based on category (as denoted by the + plot_by_spc_cat flag). The full list of species categories is + denoted in `benchmark_categories.yml + `_ + (included in GCPy). |br| + |br| + +- In this function, parallelization occurs at the species category + level. In all other functions, parallelization occurs within calls + to :code:`compare_single_level()` and :code:`compare_zonal_mean()`.= + +.. code-block:: python + + def make_benchmark_conc_plots( + ref, + refstr, + dev, + devstr, + dst="./benchmark", + subdst=None, + overwrite=False, + verbose=False, + collection="SpeciesConc", + benchmark_type="FullChemBenchmark", + cmpres=None, + plot_by_spc_cat=True, + restrict_cats=[], + plots=["sfc", "500hpa", "zonalmean"], + use_cmap_RdBu=False, + log_color_scale=False, + sigdiff_files=None, + normalize_by_area=False, + cats_in_ugm3=["Aerosols", "Secondary_Organic_Aerosols"], + areas=None, + refmet=None, + devmet=None, + weightsdir='.', + n_job=-1, + second_ref=None, + second_dev=None, + time_mean=False, + spcdb_dir=os.path.dirname(__file__) + ): + """ + Creates PDF files containing plots of species concentration + for model benchmarking purposes. + + Args: + ref: str + Path name for the "Ref" (aka "Reference") data set. + refstr: str + A string to describe ref (e.g. version number) + dev: str + Path name for the "Dev" (aka "Development") data set. + This data set will be compared against the "Reference" + data set. + devstr: str + A string to describe dev (e.g. version number) + + Keyword Args (optional): + dst: str + A string denoting the destination folder where a PDF + file containing plots will be written. + Default value: ./benchmark + subdst: str + A string denoting the sub-directory of dst where PDF + files containing plots will be written. In practice, + subdst is only needed for the 1-year benchmark output, + and denotes a date string (such as "Jan2016") that + corresponds to the month that is being plotted. + Default value: None + overwrite: bool + Set this flag to True to overwrite files in the + destination folder (specified by the dst argument). + Default value: False + verbose: bool + Set this flag to True to print extra informational output. + Default value: False + collection: str + Name of collection to use for plotting. + Default value: "SpeciesConc" + benchmark_type: str + A string denoting the type of benchmark output to plot, options are + FullChemBenchmark, TransportTracersBenchmark, or CH4Benchmark. + Default value: "FullChemBenchmark" + cmpres: string + Grid resolution at which to compare ref and dev data, e.g. '1x1.25' + plot_by_spc_cat: logical + Set this flag to False to send plots to one file rather + than separate file per category. + Default value: True + restrict_cats: list of strings + List of benchmark categories in benchmark_categories.yml to make + plots for. If empty, plots are made for all categories. + Default value: empty + plots: list of strings + List of plot types to create. + Default value: ['sfc', '500hpa', 'zonalmean'] + log_color_scale: bool + Set this flag to True to enable plotting data (not diffs) + on a log color scale. + Default value: False + normalize_by_area: bool + Set this flag to true to enable normalization of data + by surfacea area (i.e. kg s-1 --> kg s-1 m-2). + Default value: False + cats_in_ugm3: list of str + List of benchmark categories to to convert to ug/m3 + Default value: ["Aerosols", "Secondary_Organic_Aerosols"] + areas: dict of xarray DataArray: + Grid box surface areas in m2 on Ref and Dev grids. + Default value: None + refmet: str + Path name for ref meteorology + Default value: None + devmet: str + Path name for dev meteorology + Default value: None + sigdiff_files: list of str + Filenames that will contain the lists of species having + significant differences in the 'sfc', '500hpa', and + 'zonalmean' plots. These lists are needed in order to + fill out the benchmark approval forms. + Default value: None + weightsdir: str + Directory in which to place (and possibly reuse) xESMF regridder + netCDF files. + Default value: '.' + n_job: int + Defines the number of simultaneous workers for parallel plotting. + Set to 1 to disable parallel plotting. Value of -1 allows the + application to decide. + Default value: -1 + second_ref: str + Path name for a second "Ref" (aka "Reference") data set for + diff-of-diffs plotting. This dataset should have the same model + type and grid as ref. + Default value: None + second_dev: str + Path name for a second "Ref" (aka "Reference") data set for + diff-of-diffs plotting. This dataset should have the same model + type and grid as ref. + Default value: None + spcdb_dir: str + Directory of species_datbase.yml file + Default value: Directory of GCPy code repository + time_mean : bool + Determines if we should average the datasets over time + Default value: False + """ + +.. _bmk-funcs-plot-emis: + +make_benchmark_emis_plots +------------------------- + +This function generates plots of total emissions using output from +:file:`HEMCO_diagnostics.*` (for GEOS-Chem Classic) and/or +:file:`GCHP.Emissions.*` output files. + +.. code-block:: python + + def make_benchmark_emis_plots( + ref, + refstr, + dev, + devstr, + dst="./benchmark", + subdst=None, + plot_by_spc_cat=False, + plot_by_hco_cat=False, + benchmark_type="FullChemBenchmark", + cmpres=None, + overwrite=False, + verbose=False, + flip_ref=False, + flip_dev=False, + log_color_scale=False, + sigdiff_files=None, + weightsdir='.', + n_job=-1, + time_mean=False, + spcdb_dir=os.path.dirname(__file__) + ): + """ + Creates PDF files containing plots of emissions for model + benchmarking purposes. This function is compatible with benchmark + simulation output only. It is not compatible with transport tracers + emissions diagnostics. + + Args: + ref: str + Path name for the "Ref" (aka "Reference") data set. + refstr: str + A string to describe ref (e.g. version number) + dev: str + Path name for the "Dev" (aka "Development") data set. + This data set will be compared against the "Reference" + data set. + devstr: str + A string to describe dev (e.g. version number) + + Keyword Args (optional): + dst: str + A string denoting the destination folder where + PDF files containing plots will be written. + Default value: './benchmark + subdst: str + A string denoting the sub-directory of dst where PDF + files containing plots will be written. In practice, + and denotes a date string (such as "Jan2016") that + corresponds to the month that is being plotted. + Default value: None + plot_by_spc_cat: bool + Set this flag to True to separate plots into PDF files + according to the benchmark species categories (e.g. Oxidants, + Aerosols, Nitrogen, etc.) These categories are specified + in the YAML file benchmark_species.yml. + Default value: False + plot_by_hco_cat: bool + Set this flag to True to separate plots into PDF files + according to HEMCO emissions categories (e.g. Anthro, + Aircraft, Bioburn, etc.) + Default value: False + benchmark_type: str + A string denoting the type of benchmark output to plot, options are + FullChemBenchmark, TransportTracersBenchmark, or CH4Benchmark. + Default value: "FullChemBenchmark" + cmpres: string + Grid resolution at which to compare ref and dev data, e.g. '1x1.25' + overwrite: bool + Set this flag to True to overwrite files in the + destination folder (specified by the dst argument). + Default value: False + verbose: bool + Set this flag to True to print extra informational output. + Default value: False + flip_ref: bool + Set this flag to True to reverse the vertical level + ordering in the "Ref" dataset (in case "Ref" starts + from the top of atmosphere instead of the surface). + Default value: False + flip_dev: bool + Set this flag to True to reverse the vertical level + ordering in the "Dev" dataset (in case "Dev" starts + from the top of atmosphere instead of the surface). + Default value: False + log_color_scale: bool + Set this flag to True to enable plotting data (not diffs) + on a log color scale. + Default value: False + sigdiff_files: list of str + Filenames that will contain the lists of species having + significant differences in the 'sfc', '500hpa', and + 'zonalmean' plots. These lists are needed in order to + fill out the benchmark approval forms. + Default value: None + weightsdir: str + Directory in which to place (and possibly reuse) xESMF regridder + netCDF files. + Default value: '.' + n_job: int + Defines the number of simultaneous workers for parallel plotting. + Set to 1 to disable parallel plotting. + Value of -1 allows the application to decide. + Default value: -1 + spcdb_dir: str + Directory of species_datbase.yml file + Default value: Directory of GCPy code repository + time_mean : bool + Determines if we should average the datasets over time + Default value: False + + Remarks: + (1) If both plot_by_spc_cat and plot_by_hco_cat are + False, then all emission plots will be placed into the + same PDF file. + + (2) Emissions that are 3-dimensional will be plotted as + column sums. + column sums. + """ + +.. _bmk-funcs-plot-jvalue: + +make_benchmark_jvalue_plots +--------------------------- + +This function generates plots of J-values using the :literal:`JValues` +GEOS-Chem output files. + +.. code-block:: python + + def make_benchmark_jvalue_plots( + ref, + refstr, + dev, + devstr, + varlist=None, + dst="./benchmark", + subdst=None, + local_noon_jvalues=False, + cmpres=None, + plots=["sfc", "500hpa", "zonalmean"], + overwrite=False, + verbose=False, + flip_ref=False, + flip_dev=False, + log_color_scale=False, + sigdiff_files=None, + weightsdir='.', + n_job=-1, + time_mean=False, + spcdb_dir=os.path.dirname(__file__) + ): + """ + Creates PDF files containing plots of J-values for model + benchmarking purposes. + + Args: + ref: str + Path name for the "Ref" (aka "Reference") data set. + refstr: str + A string to describe ref (e.g. version number) + dev: str + Path name for the "Dev" (aka "Development") data set. + This data set will be compared against the "Reference" + data set. + devstr: str + A string to describe dev (e.g. version number) + + Keyword Args (optional): + varlist: list of str + List of J-value variables to plot. If not passed, + then all J-value variables common to both dev + and ref will be plotted. The varlist argument can be + a useful way of restricting the number of variables + plotted to the pdf file when debugging. + Default value: None + dst: str + A string denoting the destination folder where a + PDF file containing plots will be written. + Default value: ./benchmark. + subdst: str + A string denoting the sub-directory of dst where PDF + files containing plots will be written. In practice, + subdst is only needed for the 1-year benchmark output, + and denotes a date string (such as "Jan2016") that + corresponds to the month that is being plotted. + Default value: None + local_noon_jvalues: bool + Set this flag to plot local noon J-values. This will + divide all J-value variables by the JNoonFrac counter, + which is the fraction of the time that it was local noon + at each location. + Default value: False + cmpres: string + Grid resolution at which to compare ref and dev data, e.g. '1x1.25' + plots: list of strings + List of plot types to create. + Default value: ['sfc', '500hpa', 'zonalmean'] + overwrite: bool + Set this flag to True to overwrite files in the + destination folder (specified by the dst argument). + Default value: False. + verbose: bool + Set this flag to True to print extra informational output. + Default value: False + flip_ref: bool + Set this flag to True to reverse the vertical level + ordering in the "Ref" dataset (in case "Ref" starts + from the top of atmosphere instead of the surface). + Default value: False + flip_dev: bool + Set this flag to True to reverse the vertical level + ordering in the "Dev" dataset (in case "Dev" starts + from the top of atmosphere instead of the surface). + Default value: False + log_color_scale: bool + Set this flag to True if you wish to enable plotting data + (not diffs) on a log color scale. + Default value: False + sigdiff_files: list of str + Filenames that will contain the lists of J-values having + significant differences in the 'sfc', '500hpa', and + 'zonalmean' plots. These lists are needed in order to + fill out the benchmark approval forms. + Default value: None + weightsdir: str + Directory in which to place (and possibly reuse) xESMF regridder + netCDF files. + Default value: '.' + n_job: int + Defines the number of simultaneous workers for parallel plotting. + Set to 1 to disable parallel plotting. Value of -1 allows the + application to decide. + Default value: -1 + spcdb_dir: str + Directory of species_datbase.yml file + Default value: Directory of GCPy code repository + time_mean : bool + Determines if we should average the datasets over time + Default value: False + + Remarks: + Will create 4 files containing J-value plots: + (1 ) Surface values + (2 ) 500 hPa values + (3a) Full-column zonal mean values. + (3b) Stratospheric zonal mean values + These can be toggled on/off with the plots keyword argument. + + At present, we do not yet have the capability to split the + plots up into separate files per category (e.g. Oxidants, + Aerosols, etc.). This is primarily due to the fact that + we archive J-values from GEOS-Chem for individual species + but not family species. We could attempt to add this + functionality later if there is sufficient demand. + """ + +.. _bmk-funcs-plot-wetdep: + +make_benchmark_wetdep_plots +--------------------------- + +This function generates plots of wet deposition using +:literal:`WetLossConv` and :literal:`WetLossLS` GEOS-Chem output files. +It is currently primarily used for 1-Year Transport Tracer benchmarks, +plotting values for the following species as defined in +`benchmark_categories.yml +`_ +(included in GCPY). + +.. code-block:: python + + def make_benchmark_wetdep_plots( + ref, + refstr, + dev, + devstr, + collection, + dst="./benchmark", + cmpres=None, + datestr=None, + overwrite=False, + verbose=False, + benchmark_type="TransportTracersBenchmark", + plots=["sfc", "500hpa", "zonalmean"], + log_color_scale=False, + normalize_by_area=False, + areas=None, + refmet=None, + devmet=None, + weightsdir='.', + n_job=-1, + time_mean=False, + spcdb_dir=os.path.dirname(__file__) + ): + """ + Creates PDF files containing plots of species concentration + for model benchmarking purposes. + + Args: + ref: str + Path name for the "Ref" (aka "Reference") data set. + refstr: str + A string to describe ref (e.g. version number) + dev: str + Path name for the "Dev" (aka "Development") data set. + This data set will be compared against the "Reference" + data set. + devstr: str + A string to describe dev (e.g. version number) + collection: str + String name of collection to plot comparisons for. + + Keyword Args (optional): + dst: str + A string denoting the destination folder where a PDF + file containing plots will be written. + Default value: ./benchmark + datestr: str + A string with date information to be included in both the + plot pdf filename and as a destination folder subdirectory + for writing plots + Default value: None + benchmark_type: str + A string denoting the type of benchmark output to plot, options are + FullChemBenchmark, TransportTracersBenchmark, or CH4Benchmark. + Default value: "FullChemBenchmark" + overwrite: bool + Set this flag to True to overwrite files in the + destination folder (specified by the dst argument). + Default value: False. + verbose: bool + Set this flag to True to print extra informational output. + Default value: False. + plots: list of strings + List of plot types to create. + Default value: ['sfc', '500hpa', 'zonalmean'] + normalize_by_area: bool + Set this flag to true to enable normalization of data + by surfacea area (i.e. kg s-1 --> kg s-1 m-2). + Default value: False + areas: dict of xarray DataArray: + Grid box surface areas in m2 on Ref and Dev grids. + Default value: None + refmet: str + Path name for ref meteorology + Default value: None + devmet: str + Path name for dev meteorology + Default value: None + n_job: int + Defines the number of simultaneous workers for parallel plotting. + Set to 1 to disable parallel plotting. Value of -1 allows the + application to decide. + Default value: -1 + spcdb_dir: str + Directory of species_datbase.yml file + Default value: Directory of GCPy code repository + time_mean : bool + Determines if we should average the datasets over time + Default value: False + """ + +.. _bmk-funcs-table: + +=========================== +Benchmark tabling functions +=========================== + +.. table:: **Functions creating summary tables from benchmark + simulation output** + + +-----------------------------------------------+----------------------------------------------+ + | Function | Type of summary table created | + +===============================================+==============================================+ + | ``make_benchmark_aerosol_tables()`` | Global aerosol burdens (1yr benchmarks only) | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_emis_tables()`` | Emissions (by species & inventory) | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_mass_tables()`` | Total mass of each species | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_mass_accumulation_tables()`` | Mass accumulation for each species | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_mass_conservation_table()`` | Total mass of a single species at hourly | + | | intervals (to check mass conservation) | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_oh_metrics()`` | Global OH metrics (mean OH, CH4 lifetime, | + | | methylchloroform lifetime) | + +-----------------------------------------------+----------------------------------------------+ + | ``make_benchmark_operations_budget()`` | Total mass of each species after each | + | | operation (transport, mixing, etc.) | + +-----------------------------------------------+----------------------------------------------+ + +The functions listed above create summary tables for quantities such as +total mass of species, total mass of emissions, and OH metrics. + +Many of these functions use pre-defined lists of variables in YAML +files. If one dataset includes a variable but the other dataset does +not, the data for that variable in the latter dataset will be +considered to be NaN and will be plotted as such. + +.. _bmk-funcs-table-aer: + +make_benchmark_aerosol_tables +----------------------------- + +This function creates tables of global aerosol budgets and burdens from GEOS-Chem +1-year full-chemistry benchmark simulation output. + +.. code-block:: python + + def make_benchmark_aerosol_tables( + devdir, + devlist_aero, + devlist_spc, + devlist_met, + devstr, + year, + days_per_mon, + dst='./benchmark', + overwrite=False, + is_gchp=False, + spcdb_dir=os.path.dirname(__file__) + ): + """ + Compute FullChemBenchmark aerosol budgets & burdens + + Args: + devdir: str + Path to development ("Dev") data directory + devlist_aero: list of str + List of Aerosols collection files (different months) + devlist_spc: list of str + List of SpeciesConc collection files (different months) + devlist_met: list of str + List of meteorology collection files (different months) + devstr: str + Descriptive string for datasets (e.g. version number) + year: str + The year of the benchmark simulation (e.g. '2016'). + days_per_month: list of int + List of number of days per month for all months + + Keyword Args (optional): + dst: str + Directory where budget tables will be created. + Default value: './benchmark' + overwrite: bool + Overwrite burden & budget tables? (default=True) + Default value: False + is_gchp: bool + Whether datasets are for GCHP + Default value: False + spcdb_dir: str + Directory of species_datbase.yml file + Default value: Directory of GCPy code repository + + """ + +.. _bmk-funcs-table-emis: + +make_benchmark_emis_tables +-------------------------- + +This function creates tables of emissions (by species and by +inventory) from the output of GEOS-Chem benchmark simulations. + +.. code-block:: python + + def make_benchmark_emis_tables( + reflist, + refstr, + devlist, + devstr, + dst="./benchmark", + benchmark_type="FullChemBenchmark", + refmet=None, + devmet=None, + overwrite=False, + ref_interval=[2678400.0], + dev_interval=[2678400.0], + spcdb_dir=os.path.dirname(__file__) + ): + """ + Creates a text file containing emission totals by species and + category for benchmarking purposes. + + Args: + reflist: list of str + List with the path names of the emissions file or files + (multiple months) that will constitute the "Ref" + (aka "Reference") data set. + refstr: str + A string to describe ref (e.g. version number) + devlist: list of str + List with the path names of the emissions file or files + (multiple months) that will constitute the "Dev" + (aka "Development") data set + devstr: str + A string to describe dev (e.g. version number) + + Keyword Args (optional): + dst: str + A string denoting the destination folder where the file + containing emissions totals will be written. + Default value: ./benchmark + benchmark_type: str + A string denoting the type of benchmark output to plot, options are + FullChemBenchmark, TransportTracersBenchmark or CH4Benchmark. + Default value: "FullChemBenchmark" + refmet: str + Path name for ref meteorology + Default value: None + devmet: str + Path name for dev meteorology + Default value: None + overwrite: bool + Set this flag to True to overwrite files in the + destination folder (specified by the dst argument). + Default value: False + ref_interval: list of float + The length of the ref data interval in seconds. By default, interval + is set to [2678400.0], which is the number of seconds in July + (our 1-month benchmarking month). + Default value: [2678400.0] + dev_interval: list of float + The length of the dev data interval in seconds. By default, interval + is set to [2678400.0], which is the number of seconds in July + (our 1-month benchmarking month). + Default value: [2678400.0] + spcdb_dir: str + Directory of species_datbase.yml file + Default value: Directory of GCPy code repository + + """ + +.. _bmk-funcs-table-mass: + +make_benchmark_mass_tables +-------------------------- + +This function creates tables of total mass for species in two +different GEOS-Chem benchmark simulations. + +.. code-block:: python + + def make_benchmark_mass_tables( + ref, + refstr, + dev, + devstr, + varlist=None, + dst="./benchmark", + subdst=None, + overwrite=False, + verbose=False, + label="at end of simulation", + spcdb_dir=os.path.dirname(__file__), + ref_met_extra=None, + dev_met_extra=None + ): + """ + Creates a text file containing global mass totals by species and + category for benchmarking purposes. + + Args: + reflist: str + Pathname that will constitute + the "Ref" (aka "Reference") data set. + refstr: str + A string to describe ref (e.g. version number) + dev: list of str + Pathname that will constitute + the "Dev" (aka "Development") data set. The "Dev" + data set will be compared against the "Ref" data set. + devstr: str + A string to describe dev (e.g. version number) + + Keyword Args (optional): + varlist: list of str + List of variables to include in the list of totals. + If omitted, then all variables that are found in either + "Ref" or "Dev" will be included. The varlist argument + can be a useful way of reducing the number of + variables during debugging and testing. + Default value: None + dst: str + A string denoting the destination folder where the file + containing emissions totals will be written. + Default value: ./benchmark + subdst: str + A string denoting the sub-directory of dst where PDF + files containing plots will be written. In practice, + subdst is only needed for the 1-year benchmark output, + and denotes a date string (such as "Jan2016") that + corresponds to the month that is being plotted. + Default value: None + overwrite: bool + Set this flag to True to overwrite files in the + destination folder (specified by the dst argument). + Default value: False + verbose: bool + Set this flag to True to print extra informational output. + Default value: False. + spcdb_dir: str + Directory of species_datbase.yml file + Default value: Directory of GCPy code repository + ref_met_extra: str + Path to ref Met file containing area data for use with restart files + which do not contain the Area variable. + Default value: '' + dev_met_extra: str + Path to dev Met file containing area data for use with restart files + which do not contain the Area variable. + Default value: '' + """ + +.. _bmk-funcs-table-accum: + +make_benchmark_mass_accumulation_tables +--------------------------------------- + +This function creates tables of mass accumulation over time for species in two +different GEOS-Chem benchmark simulations. + +.. code-block:: python + + def create_mass_accumulation_table( + refdatastart, + refdataend, + refstr, + refperiodstr, + devdatastart, + devdataend, + devstr, + devperiodstr, + varlist, + met_and_masks, + label, + trop_only=False, + outfilename="GlobalMassAccum_TropStrat.txt", + verbose=False, + spcdb_dir=os.path.dirname(__file__) + ): + """ + Creates a table of global mass accumulation for a list of species in + two data sets. The data sets, which typically represent output from two + different model versions, are usually contained in netCDF data files. + + Args: + refdatastart: xarray Dataset + The first data set to be compared (aka "Reference"). + refdataend: xarray Dataset + The first data set to be compared (aka "Reference"). + refstr: str + A string that can be used to identify refdata + (e.g. a model version number or other identifier). + refperiodstr: str + Ref simulation period start and end + devdatastart: xarray Dataset + The second data set to be compared (aka "Development"). + devdataend: xarray Dataset + The second data set to be compared (aka "Development"). + devstr: str + A string that can be used to identify the data set specified + by devfile (e.g. a model version number or other identifier). + devperiodstr: str + Ref simulation period start and end + varlist: list of strings + List of species concentation variable names to include + in the list of global totals. + met_and_masks: dict of xarray DataArray + Dictionary containing the meterological variables and + masks for the Ref and Dev datasets. + label: str + Label to go in the header string. Can be used to + pass the month & year. + + Keyword Args (optional): + trop_only: bool + Set this switch to True if you wish to print totals + only for the troposphere. + Default value: False (i.e. print whole-atmosphere totals). + outfilename: str + Name of the text file which will contain the table of + emissions totals. + Default value: "GlobalMass_TropStrat.txt" + verbose: bool + Set this switch to True if you wish to print out extra + informational messages. + Default value: False + spcdb_dir: str + Directory of species_datbase.yml file + Default value: Directory of GCPy code repository + + Remarks: + This method is mainly intended for model benchmarking purposes, + rather than as a general-purpose tool. + + Species properties (such as molecular weights) are read from a + YAML file called "species_database.yml". + """ + + +.. _bmk-funcs-table-cons: + +make_benchmark_mass_conservation_table +-------------------------------------- + +This function creates a timeseries table of the global mass of the +:literal:`PassiveTracer` species. Usually used with output from +1-year TransportTracers benchmark simulations. + +.. code-block:: python + + def make_benchmark_mass_conservation_table( + datafiles, + runstr, + dst="./benchmark", + overwrite=False, + areapath=None, + spcdb_dir=os.path.dirname(__file__) + ): + """ + Creates a text file containing global mass of the PassiveTracer + from Transport Tracer simulations across a series of restart files. + + Args: + datafiles: list of str + Path names of restart files. + runstr: str + Name to put in the filename and header of the output file + refstr: str + A string to describe ref (e.g. version number) + dev: str + Path name of "Dev" (aka "Development") data set file. + The "Dev" data set will be compared against the "Ref" data set. + devmet: list of str + Path name of dev meteorology data set. + devstr: str + A string to describe dev (e.g. version number) + + Keyword Args (optional): + dst: str + A string denoting the destination folder where the file + containing emissions totals will be written. + Default value: "./benchmark" + overwrite: bool + Set this flag to True to overwrite files in the + destination folder (specified by the dst argument). + Default value: False + areapath: str + Path to a restart file containing surface area data. + Default value: None + spcdb_dir: str + Path to the species_database.yml + Default value: points to gcpy/gcpy folder + """ + +.. _bmk-funcs-table-oh: + +make_benchmark_oh_metrics +------------------------- + +This function generates a table of OH metrics from GEOS-Chem benchmark +simulation output. + +.. code-block:: python + + def make_benchmark_oh_metrics( + ref, + refmet, + refstr, + dev, + devmet, + devstr, + dst="./benchmark", + overwrite=False, + ): + """ + Creates a text file containing metrics of global mean OH, MCF lifetime, + and CH4 lifetime for benchmarking purposes. + + Args: + ref: str + Path name of "Ref" (aka "Reference") data set file. + refmet: str + Path name of ref meteorology data set. + refstr: str + A string to describe ref (e.g. version number) + dev: str + Path name of "Dev" (aka "Development") data set file. + The "Dev" data set will be compared against the "Ref" data set. + devmet: list of str + Path name of dev meteorology data set. + devstr: str + A string to describe dev (e.g. version number) + + Keyword Args (optional): + dst: str + A string denoting the destination folder where the file + containing emissions totals will be written. + Default value: "./benchmark" + overwrite: bool + Set this flag to True to overwrite files in the + destination folder (specified by the dst argument). + Default value: False + """ + +.. _bmk-funcs-table-ops: + +make_benchmark_operations_budget +-------------------------------- + +Creates a table with the change in species mass after each GEOS-Chem +operation, using output from GEOS-Chem benchmark simulations. + +.. code-block:: python + + def make_benchmark_operations_budget( + refstr, + reffiles, + devstr, + devfiles, + ref_interval, + dev_interval, + benchmark_type=None, + label=None, + col_sections=["Full", "Trop", "PBL", "Strat"], + operations=[ + "Chemistry", "Convection", "EmisDryDep", + "Mixing", "Transport", "WetDep" + ], + compute_accum=True, + compute_restart=False, + require_overlap=False, + dst='.', + species=None, + overwrite=True, + verbose=False, + spcdb_dir=os.path.dirname(__file__) + ): + """ + Prints the "operations budget" (i.e. change in mass after + each operation) from a GEOS-Chem benchmark simulation. + + Args: + refstr: str + Labels denoting the "Ref" versions + reffiles: list of str + Lists of files to read from the "Ref" version. + devstr: str + Labels denoting the "Dev" versions + devfiles: list of str + Lists of files to read from "Dev" version. + interval: float + Number of seconds in the diagnostic interval. + + Keyword Args (optional): + benchmark_type: str + A string denoting the type of benchmark output to plot, options are + FullChemBenchmark, TransportTracersBenchmark, or CH4Benchmark. + Default value: None + label: str + Contains the date or date range for each dataframe title. + Default value: None + col_sections: list of str + List of column sections to calculate global budgets for. May + include Strat eventhough not calculated in GEOS-Chem, but Full + and Trop must also be present to calculate Strat. + Default value: ["Full", "Trop", "PBL", "Strat"] + operations: list of str + List of operations to calculate global budgets for. Accumulation + should not be included. It will automatically be calculated if + all GEOS-Chem budget operations are passed and optional arg + compute_accum is True. + Default value: ["Chemistry","Convection","EmisDryDep", + "Mixing","Transport","WetDep"] + compute_accum: bool + Optionally turn on/off accumulation calculation. If True, will + only compute accumulation if all six GEOS-Chem operations budgets + are computed. Otherwise a message will be printed warning that + accumulation will not be calculated. + Default value: True + compute_accum: bool + Optionally turn on/off accumulation calculation. If True, will + only compute accumulation if all six GEOS-Chem operations budgets + are computed. Otherwise a message will be printed warning that + accumulation will not be calculated. + Default value: True + compute_restart: bool + Optionally turn on/off calculation of mass change based on restart + file. Only functional for "Full" column section. + Default value: False + require_overlap: bool + Whether to calculate budgets for only species that are present in + both Ref or Dev. + Default value: False + dst: str + Directory where plots & tables will be created. + Default value: '.' (directory in which function is called) + species: list of str + List of species for which budgets will be created. + Default value: None (all species) + overwrite: bool + Denotes whether to overwrite existing budget file. + Default value: True + verbose: bool + Set this switch to True if you wish to print out extra + informational messages. + Default value: False + """ + "" diff --git a/docs/source/Getting-Started-with-GCPy.rst b/docs/source/Getting-Started-with-GCPy.rst index ad93963e..6fe74133 100644 --- a/docs/source/Getting-Started-with-GCPy.rst +++ b/docs/source/Getting-Started-with-GCPy.rst @@ -20,12 +20,29 @@ Requirements #. Windows Subsystem for Linux (running in Microsoft Windows 11) #. MacOS -To install GCPy, you will need a distribution of the :program:`Mamba` -package manager. :program:`Mamba` is a drop-in replacement for the -widely-used :program:`Conda` package manager. You will use -:program:`Mamba` to create a Python environment for GCPy, which -contains the Python interpreter plus other Python packages that GCPY -uses. +To install GCPy, you will need: + +- **EITHER** a distribution of the :program:`Mamba` package manager +- **OR** a distribution of the :program:`Conda` package manager. + +:program:`Mamba` is a fast drop-in replacement for the +widely-used :program:`Conda` package manager. We recommend using +:program:`Mamba` to create a Python environment for GCPy. This +environment will contain a version of the Python interpreter +(in this case, Python 3.9) plus packages upon which GCPy depends. + +.. note:: + + If your system has an existing :program:`Conda` installation, and/or + you do not wish to upgrade from :program:`Conda` to + :program:`Mamba`, you may create the Python environment for GCPy + with :program:`Conda`. See the following sections for detailed + instructions. + +.. _requirements-mamba: + +Check if Mamba is installed +--------------------------- Check if you already have :program:`Mamba` on your system: @@ -35,16 +52,63 @@ Check if you already have :program:`Mamba` on your system: If :program:`Mamba` has been installed, you will see output similar to this: -.. code-block:: +.. code-block:: console mamba version X.Y.Z conda version A.B.C -If you see this output, you may skip ahead to the :ref:`gcpy_install` -section. Otherwise, read further for instructions on how to install -:program:`Mamba` on your system. +If you see this output, you may skip ahead to the :ref:`gcpy-install` +section. + +.. _requirements-conda: + +Check if Conda is installed +--------------------------- + +If your system does not have :program:`Mamba` installed, check if +:program:`Conda` is already present on your system: + +.. code-block:: console + + $ conda --version + +If a :program:`Conda` version exists, you will see its version number +printed to the screen: + +.. code-block:: console + + conda version A.B.C + +If neither :program:`Conda` or :program:`Mamba` are installed, we +recommend installing the :program:`Mamba` package manager yourself. +Please proceed to the :ref:`mamba-install` section for instructions. + +.. _requirements-conda-older: + +Additional setup for older Conda versions +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +If your :program:`Conda` version is earlier than 23.7, you will need +to do the following additional steps. + +.. code-block:: console -.. _mamba_install: + $ conda install -n base conda-libmamba-solver + $ conda config --set solver libmamba + +This will install the fast :program:`Mamba` environment solver into +your :program:`Conda` base environment. Using the :program:`Mamba` +solver within :program:`Conda` will speed up the Python environment +creation considerably. + +.. note:: + + The :program:`Mamba` environment solver is used by default in + :program:`Conda` 23.7 and later. + +You may now skip ahead to the :ref:`gcpy-install` section. + +.. _mamba-install: ================== Install MambaForge @@ -85,7 +149,7 @@ MacOS #. Exit your current terminal session and open a new terminal session. This will apply the changes. -You may now skip ahead to the :ref:`gcpy_install` section. +You may now skip ahead to the :ref:`gcpy-install` section. Linux and Windows Subsystem for Linux @@ -162,7 +226,7 @@ Linux and Windows Subsystem for Linux - Press ENTER to confirm the location - Press CTRL-C to abort the installation - Or specify a different location below - [/home/bob/mambaforge] >>> + [/home/YOUR-USER-NAME/mambaforge] >>> In most cases, it should be OK to accept the default installation location. But on some systems, users may be encouraged to install @@ -179,7 +243,7 @@ Linux and Windows Subsystem for Linux :program:`MambaForge` will downlad and install Python software packages into the :file:`pkgs` subfolder of the root installation path. Similarly, when you :ref:`create Python - environments `, these will be installed to the + environments `, these will be installed to the :file:`envs` subfolder of the root installation path. |br| @@ -198,7 +262,7 @@ Linux and Windows Subsystem for Linux As long as your :envvar:`PYTHONPATH` environment variable only contains the path to the root-level GCPy folder, you may safely ignore this. (More on :envvar:`PYTHONPATH` in the :ref:`next - section `.) |br| + section `.) |br| |br| #. Tell the installer to initialize :program:`MambaForge`. @@ -219,15 +283,15 @@ Linux and Windows Subsystem for Linux #. Exit your current terminal session. Start a new terminal session to apply the updates. You are now ready to install GCPy. -.. _gcpy_install: +.. _gcpy-install: ================================= Install GCPy and its dependencies ================================= -Once :program:`Mamba` has been installed, use it to ccreate a Python -environment for GCPy. (Please return to :ref:`mamba_install` if you -have not yet installed :program:`Mamba`.) +Once you have made sure that :program:`Mamba` (or :program:`Conda`) is +present on your system, you may create a Python environment for GCPy. +Follow these steps: #. **Download the GCPy source code.** @@ -253,21 +317,24 @@ have not yet installed :program:`Mamba`.) maintain a set of package dependencies compatible with GCPy without interfering with Python packages you use for other work. You can create a Python virtual environment from anywhere on your - system. It will be stored in your :program:`Mamba` installation - rather than the directory from which you create it. + system. It will be stored in your :program:`Mamba` (or + :program:`Conda` installation rather than the directory from which + you create it). You can create a Python virtual environment using a file that lists all packages and their versions to be included in the environment. GCPy includes such as file, :file:`environment.yml`, located in the top-level directory of the package. - Run the following command at the command prompt to create a virtual + Run one of the following commands at the command prompt to create a virtual environment for use with GCPy. You can name environment whatever you wish. This example names it :file:`gcpy_env`. .. code-block:: console - $ mamba env create -n gcpy_env --file=environment.yml + $ mamba env create -n gcpy_env --file=environment.yml # If using Mamba + + $ conda env create -n gcpy_env --file=environment.yml # If using Conda A list of packages to be downloaded will be displayed. A confirmation message will ask you if you really wish to install all @@ -275,17 +342,21 @@ have not yet installed :program:`Mamba`.) :command:`n` to abort. Once successfully created you can activate the environment with - this command: + one of these commands: .. code-block:: console - $ mamba activate gcpy_env + $ mamba activate gcpy_env # If using Mamba + + $ conda activate gcpy_env # If using Conda - To exit the environment, use this command: + To exit the environment, use one of these commands: .. code-block:: console - $ mamba deactivate + $ mamba deactivate # If using Mamba + + $ conda deactivate # If using Conda |br| @@ -295,7 +366,7 @@ have not yet installed :program:`Mamba`.) locations of Python libraries on your system that were not installed by :program:`Mamba`. - Add the path to your GCPy source code folder :file:`~/.bashrc` file: + Add the path to your GCPy source code folder :file:`~/.bashrc` file: .. code-block:: bash @@ -310,6 +381,28 @@ have not yet installed :program:`Mamba`.) to apply the change. |br| |br| +#. **Set the** :envvar:`MPLBACKEND` **environment variable** + + The environment variable :envvar:`MPLBACKEND` specifies the X11 + backend that the Matplotlib package will use to render plots to the + screen. + + Add this line to your :file:`~/.bashrc` file on your local PC/Mac + and on any remote computer systems where you will use GCPy: + + .. code-block:: bash + + export MPLBACKEND=tkagg + + And then use: + + .. code-block:: console + + $ source ~/.bashrc + + to apply the change. |br| + |br| + #. **Perform a simple test:** Run the following commands in your terminal to check if the diff --git a/docs/source/Guide-to-Useful-Capabilities.rst b/docs/source/Guide-to-Useful-Capabilities.rst index 09b095c7..473adb5c 100644 --- a/docs/source/Guide-to-Useful-Capabilities.rst +++ b/docs/source/Guide-to-Useful-Capabilities.rst @@ -1,4 +1,3 @@ - .. _capabilities: ######################## @@ -11,7 +10,7 @@ function documentation. .. _capabilities-spatial: ================ -Spatial Plotting +Spatial plotting ================ One hallmark of GCPy is easy-to-use spatial plotting of GEOS-Chem @@ -23,28 +22,29 @@ mean for all layers of the atmosphere. .. _capabilities-spatial-single: -Single Panel Plots +Single panel plots ------------------ -Single panel plots are generated through the -:file:`plot.single_panel()` function. :file:`plot.single_panel()` uses -Matplotlib and Cartopy plotting capabilities while handling certain -behind the scenes operations that are necessary for plotting GEOS-Chem -data, particularly for cubed-sphere and/or zonal mean data. +Single panel plots are generated through the :file:`single_panel()` +function (located in module :file:`gcpy.plot.single_panel`). This +function uses Matplotlib and Cartopy plotting capabilities while +handling certain behind the scenes operations that are necessary for +plotting GEOS-Chem data, particularly for cubed-sphere and/or zonal +mean data. .. code:: python import xarray as xr - import gcpy.plot as gcplot import matplotlib.pyplot as plt + from gcpy.plot.single_panel import single_panel # Read data ds = xr.open_dataset( 'GEOSChem.Restart.20160701_0000z.nc4' ) - # plot surface Ozone over the North Pacific - gcplot.single_panel( + # Plot surface Ozone over the North Pacific + single_panel( ds['SpeciesRst_O3'].isel(lev=0), title='Surface Ozone over the North Pacific', extent=[80, -90, -10, 60] @@ -56,8 +56,8 @@ data, particularly for cubed-sphere and/or zonal mean data. .. code:: python - #plot global zonal mean of Ozone - gcplot.single_panel( + # Plot global zonal mean of Ozone + single_panel( ds['SpeciesRst_O3'], plot_type='zonal_mean', title='Global Zonal Mean of Ozone' @@ -68,15 +68,16 @@ data, particularly for cubed-sphere and/or zonal mean data. .. image:: _static/images/single\_panel\_zonal\_mean.png :align: center -`Click here `__ for an example single panel plotting script. -`Click here `__ for detailed documentation for ``single_panel()``. +:ref:`Click here ` for an example single panel plotting script. +:ref:`Click here ` for detailed documentation for +:code:`single_panel()`. .. _capabilities-spatial-sixpanel: -Six Panel Comparison Plots +Six-panel comparison plots -------------------------- -Six panel plots are used to compare results across two different model +Six-panel plots are used to compare results across two different model runs. Single level and zonal mean plotting options are both available. The two model runs do not need to be the same resolution or even the same grid type (GEOS-Chem Classic and GCHP output can be mixed at will). @@ -84,8 +85,9 @@ same grid type (GEOS-Chem Classic and GCHP output can be mixed at will). .. code:: python import xarray as xr - import gcpy.plot as gcplot import matplotlib.pyplot as plt + from gcpy.plot.compare_single_level import compare_single_level + from gcpy.plot.compare_zonal_mean import compare_zonal_mean # Read data gcc_ds = xr.open_dataset( @@ -95,8 +97,8 @@ same grid type (GEOS-Chem Classic and GCHP output can be mixed at will). 'GCHP.SpeciesConc.20160716_1200z.nc4' ) - #Plot comparison of surface ozone over the North Pacific - gcplot.compare_single_level( + # Plot comparison of surface ozone over the North Pacific + compare_single_level( gcc_ds, 'GEOS-Chem Classic', gchp_ds, @@ -109,11 +111,12 @@ same grid type (GEOS-Chem Classic and GCHP output can be mixed at will). .. image:: _static/images/six\_panel\_single\_level.png :align: center + :width: 80% .. code:: python - #Plot comparison of global zonal mean ozone - gcplot.compare_zonal_mean( + # Plot comparison of global zonal mean ozone + compare_zonal_mean( gcc_ds, 'GEOS-Chem Classic', gchp_ds, @@ -123,34 +126,34 @@ same grid type (GEOS-Chem Classic and GCHP output can be mixed at will). plt.show() .. image:: _static/images/six\_panel\_zonal\_mean.png + :align: center + :width: 80% -`Click here `__ for an example six panel plotting script. -`Click here `__ -for complete documentation for ``compare_single_level()`` and ``compare_zonal_mean()``. +:ref:`Click here ` for an example six panel plotting +script. :ref:`Click here ` for complete documentation +for :code:`compare_single_level()` and :code:`compare_zonal_mean()`. .. _capabilities-spatial-benchmark: -Comprehensive Benchmark Plotting +Comprehensive benchmark plotting -------------------------------- The GEOS-Chem Support Team uses comprehensive plotting functions from -:file:`benchmark.py` to generate full plots of benchmark +module :file:`gcpy.benchmark_funcs` to generate full plots of benchmark diagnostics. Functions like -:file:`benchmark.make_benchmark_conc_plots` by default create plots -for every variable in a given collection -(e.g. :literal:`SpeciesConc`) at multiple vertical levels (surface, -500hPa, zonal mean) and divide plots into separate folders based on -category (e.g. Chlorine, Aerosols). The GCST uses full benchmark -plotting / table scripts similar to `this example -`__ to produce plots and tables for official -model benchmarks. Full documentation for the benchmark plotting -functions can be found -`here `__. +:ref:`gcpy.benchmark_funcs.make_benchmark_conc_plots() +` by default create plots for every variable +in a given collection (e.g. :literal:`SpeciesConc`) at multiple +vertical levels (surface, 500hPa, zonal mean) and divide plots into +separate folders based on category (e.g. Chlorine, Aerosols). The +GEOS-Chem Support Team uses benchmark plotting and tabling table +scripts (described in our :ref:`Benchmarking ` chapter) to +produce plots and tables for official model benchmarks. .. _capabilities-table: ============== -Table Creation +Table creation ============== GCPy has several dedicated functions for tabling GEOS-Chem output data @@ -159,7 +162,7 @@ used for model benchmarking purposes. .. _capabilities-table-budget: -Budget Tables +Budget tables ------------- Currently, budget tables can be created for "operations" (table shows @@ -168,39 +171,35 @@ the GEOS-Chem :literal:`Budget` diagnostics) or in overall averages for different aerosols or the Transport Tracers simulation. Operations budget tables are created using the -:file:`benchmark.make_benchmark_operations_budget` function and appear as -follows: +:ref:`gcpy.benchmark_funcs.make_benchmark_operations_budget() +` function and appear as follows: .. image:: _static/images/budget\_table.png :align: center -Full documentation for operations budget table creation can be found -`here `__. - .. _capabilities-tables-mass: -Mass Tables +Mass tables ----------- -The :file:`benchmark.make_benchmark_mass_tables` function uses species -concentrations and info from meteorology files to generate the total -mass of species in certain segments of the atmosphere (currently -global or only the troposphere). An example table is shown below: +The :ref:`gcpy.benchmark_funcs.make_benchmark_mass_tables() +` function uses species concentrations and info +from meteorology files to generate the total mass of species in +certain segments of the atmosphere (currently global or only the +troposphere). An example table is shown below: .. image:: _static/images/mass\_table.png :align: center -Full documentation for mass table creation can be found -`here `__. - .. _capabilities-tables-emissions: -Emissions Tables +Emissions tables ---------------- -The :file:`benchmark.make_benchmark_emis_tables` function creates -tables of total emissions categorized by species or by -inventory. Examples of both emissions table types are shown below: +The :ref:`gcpy.benchmark_funcs.make_benchmark_emis_tables() +` function creates tables of total emissions +categorized by species or by inventory. Examples of both emissions +table types are shown below: .. image:: _static/images/emissions\_totals.png :align: center @@ -208,9 +207,6 @@ inventory. Examples of both emissions table types are shown below: .. image:: _static/images/inventory\_totals.png :align: center -Full documentation for emissions table creation can be found `here -`__. - .. _capabilities-regridding: ========== @@ -219,24 +215,25 @@ Regridding .. _capabilities-regridding-rules: -General Regridding Rules +General regridding rules ------------------------ -GCPy supports regridding between all horizontal GEOS-Chem grid types, including -latitude/longitude grids (the grid format of GEOS-Chem Classic), +GCPy supports regridding between all horizontal GEOS-Chem grid types, +including latitude/longitude grids (the grid format of GEOS-Chem Classic), standard cubed-sphere (the standard grid format of GCHP), and stretched-grid (an optional grid format in GCHP). GCPy contains several horizontal regridding functions built off of xESMF. GCPy automatically handles most regridding needs when plotting GEOS-Chem data. -:file:`gcpy.file_regrid` allows you to regrid GEOS-Chem Classic files between -different grid resolutions and can be called from the command line or as a -function. +:ref:`gcpy.file_regrid() ` allows you to regrid +GEOS-Chem Classic and GCHP files between different grid resolutions +and can be called from the command line or as a function. -:file:`gcpy.regrid_restart_file` allows you to regrid GCHP files between -between different grid resolutions and grid types (standard and stretched -cubed-sphere grids), and can be called from the command line. +:ref:`gcpy.regrid_restart_file ` allows you to regrid +GCHP files between between different grid resolutions and grid +types (standard and stretched cubed-sphere grids), and can be +called from the command line. The 72-level and 47-level vertical grids are pre-defined in GCPy. Other vertical grids can also be defined if you provide `the A @@ -244,12 +241,13 @@ and B coefficients of the hybrid vertical grid `__. When plotting data of differing grid types or horizontal resolutions -using :file:`compare_single_level` or :file:`compare_zonal_mean`, you +using :ref:`compare_single_level() ` +or :ref:`compare_zonal_mean() `, you can specify a comparison resolution using the :literal:`cmpres` argument. This resolution will be used for the difference panels in each plot (the bottom four panels rather than the top two raw data panels). If you do not specify a comparison resolution, GCPy will automatically choose one. -For more extensive regridding information, visit the `detailed -regridding documentation `__. +For more extensive regridding information, visit the :ref:`detailed +regridding documentation `. diff --git a/docs/source/Plotting.rst b/docs/source/Plotting.rst index 1a22c844..2496a0e9 100644 --- a/docs/source/Plotting.rst +++ b/docs/source/Plotting.rst @@ -4,54 +4,423 @@ Plotting ######## -This page describes in depth the plotting capabilities of GCPy, +This page describes in depth the general plotting capabilities of GCPy, including possible argument values for every plotting function. -.. _plot-csl-and-czm: +For information about GCPy functions that are specific to the +GEOS-Chem benchmark workflow, please see our :ref:`Benchmarking ` +chapter. -=========================================== -compare_single_level and compare_zonal_mean -=========================================== +.. _plot-six-panel: -:code:`gcpy.plot.compare_single_level()` and -:code:`gcpy.plot.compare_zonal_mean()` both generate six panel -plots comparing variables between two datasets. They share significant -structural overlap both in output appearance and code -implementation. This section gives an overview of the components of -the plots generated by these functions, their shared arguments, and -features unique to each function. +========================== +Six-panel comparison plots +========================== -.. _plot-csl-and-czm-shared: +The functions listed below generate six-panel plots comparing +variables between two datasets: -Shared structure ----------------- ++-----------------------------+------------------------------------+ +| Plotting function | Located in GCPy module | ++=============================+====================================+ +| ``compare_single_level()`` | ``gcpy.plot.compare_single_level`` | ++-----------------------------+------------------------------------+ +| ``compare_zonal_mean()`` | ``gcpy.plot.compare_zonal_mean`` | ++-----------------------------+------------------------------------+ Both :code:`compare_single_level()` and :code:`compare_zonal_mean()` generate a six panel plot for each variable passed. These plots can either be saved to PDFs or generated sequentially for visualization in the Matplotlib GUI using :code:`matplotlib.pyplot.show()`. Each plot uses data passed from a reference (:literal:`Ref`) dataset -and a development (:literal:`Dev`) dataset. +and a development (:literal:`Dev`) dataset. Both functions share +significant structural overlap both in output appearance and code +implementation. + +You can import these routines into your code with these statements: + +.. code-block:: python + + from gcpy.plot.compare_single_level import compare_single_level + from gcpy.plot.compare_zonal_mean import compare_zonal_mean Each panel has a title describing the type of panel, a colorbar for the values plotted in that panel, and the units of the data plotted in that panel. The upper two panels of each plot show actual values from -the Ref (left) and Dev (right) datasets for a given variable. The -middle two panels show the difference (:literal:`Dev - Ref`) between -the values in the Dev dataset and the values in the Ref dataset. The -left middle panel uses a full dynamic color map, while the right -middle panel caps the color map at the 5th and 95th percentiles. The -bottom two panels show the ratio (:literal:`Dev/Ref`) between the -values in the Dev dataset and the values in the Ref Dataset. -The left bottom panel uses a full dynamic color map, while the right -bottom panel caps the color map at 0.5 and 2.0. +the :literal:`Ref` (left) and :literal:`Dev` (right) datasets for a +given variable. The middle two panels show the difference +(:literal:`Dev - Ref`) between the values in the :literal:`Dev` +dataset and the values in the :literal:`Ref` dataset. The left middle +panel uses a full dynamic color map, while the right middle panel caps +the color map at the 5th and 95th percentiles. The bottom two panels +show the ratio (:literal:`Dev/Ref`) between the values in the Dev +dataset and the values in the Ref Dataset. The left bottom panel uses +a full dynamic color map, while the right bottom panel caps the color +map at 0.5 and 2.0. + +.. _plot-csl: + +Function :code:`compare_single_level` +------------------------------------- + +The :code:`compare_single_level` function accepts takes the following +arguments: + +.. code-block:: python + + def compare_single_level( + refdata, + refstr, + devdata, + devstr, + varlist=None, + ilev=0, + itime=0, + refmet=None, + devmet=None, + weightsdir='.', + pdfname="", + cmpres=None, + match_cbar=True, + normalize_by_area=False, + enforce_units=True, + convert_to_ugm3=False, + flip_ref=False, + flip_dev=False, + use_cmap_RdBu=False, + verbose=False, + log_color_scale=False, + extra_title_txt=None, + extent=None, + n_job=-1, + sigdiff_list=None, + second_ref=None, + second_dev=None, + spcdb_dir=os.path.dirname(__file__), + sg_ref_path='', + sg_dev_path='', + ll_plot_func='imshow', + **extra_plot_args + ): + """ + Create single-level 3x2 comparison map plots for variables common + in two xarray Datasets. Optionally save to PDF. + + Args: + refdata: xarray dataset + Dataset used as reference in comparison + refstr: str + String description for reference data to be used in plots + devdata: xarray dataset + Dataset used as development in comparison + devstr: str + String description for development data to be used in plots + + Keyword Args (optional): + varlist: list of strings + List of xarray dataset variable names to make plots for + Default value: None (will compare all common variables) + ilev: integer + Dataset level dimension index using 0-based system. + Indexing is ambiguous when plotting differing vertical grids + Default value: 0 + itime: integer + Dataset time dimension index using 0-based system + Default value: 0 + refmet: xarray dataset + Dataset containing ref meteorology + Default value: None + devmet: xarray dataset + Dataset containing dev meteorology + Default value: None + weightsdir: str + Directory path for storing regridding weights + Default value: None (will create/store weights in + current directory) + pdfname: str + File path to save plots as PDF + Default value: Empty string (will not create PDF) + cmpres: str + String description of grid resolution at which + to compare datasets + Default value: None (will compare at highest resolution + of ref and dev) + match_cbar: bool + Set this flag to True if you wish to use the same colorbar + bounds for the Ref and Dev plots. + Default value: True + normalize_by_area: bool + Set this flag to True if you wish to normalize the Ref + and Dev raw data by grid area. Input ref and dev datasets + must include AREA variable in m2 if normalizing by area. + Default value: False + enforce_units: bool + Set this flag to True to force an error if Ref and Dev + variables have different units. + Default value: True + convert_to_ugm3: bool + Whether to convert data units to ug/m3 for plotting. + Default value: False + flip_ref: bool + Set this flag to True to flip the vertical dimension of + 3D variables in the Ref dataset. + Default value: False + flip_dev: bool + Set this flag to True to flip the vertical dimension of + 3D variables in the Dev dataset. + Default value: False + use_cmap_RdBu: bool + Set this flag to True to use a blue-white-red colormap + for plotting the raw data in both the Ref and Dev datasets. + Default value: False + verbose: bool + Set this flag to True to enable informative printout. + Default value: False + log_color_scale: bool + Set this flag to True to plot data (not diffs) + on a log color scale. + Default value: False + extra_title_txt: str + Specifies extra text (e.g. a date string such as "Jan2016") + for the top-of-plot title. + Default value: None + extent: list + Defines the extent of the region to be plotted in form + [minlon, maxlon, minlat, maxlat]. + Default value plots extent of input grids. + Default value: [-1000, -1000, -1000, -1000] + n_job: int + Defines the number of simultaneous workers for parallel + plotting. Set to 1 to disable parallel plotting. + Value of -1 allows the application to decide. + Default value: -1 + sigdiff_list: list of str + Returns a list of all quantities having significant + differences (where |max(fractional difference)| > 0.1). + Default value: None + second_ref: xarray Dataset + A dataset of the same model type / grid as refdata, + to be used in diff-of-diffs plotting. + Default value: None + second_dev: xarray Dataset + A dataset of the same model type / grid as devdata, + to be used in diff-of-diffs plotting. + Default value: None + spcdb_dir: str + Directory containing species_database.yml file. + Default value: Path of GCPy code repository + sg_ref_path: str + Path to NetCDF file containing stretched-grid info + (in attributes) for the ref dataset + Default value: '' (will not be read in) + sg_dev_path: str + Path to NetCDF file containing stretched-grid info + (in attributes) for the dev dataset + Default value: '' (will not be read in) + ll_plot_func: str + Function to use for lat/lon single level plotting with + possible values 'imshow' and 'pcolormesh'. imshow is much + faster but is slightly displaced when plotting from + dateline to dateline and/or pole to pole. + Default value: 'imshow' + extra_plot_args: various + Any extra keyword arguments are passed through the + plotting functions to be used in calls to pcolormesh() (CS) + or imshow() (Lat/Lon). + """ + +and generates a comparison plot such as: + +.. image:: _static/images/six\_panel\_single\_level.png + :align: center + :width: 80% + +.. _plot-czm: + +Function :code:`compare_zonal_mean` +----------------------------------- + +.. code-block:: python + + def compare_zonal_mean( + refdata, + refstr, + devdata, + devstr, + varlist=None, + itime=0, + refmet=None, + devmet=None, + weightsdir='.', + pdfname="", + cmpres=None, + match_cbar=True, + pres_range=None, + normalize_by_area=False, + enforce_units=True, + convert_to_ugm3=False, + flip_ref=False, + flip_dev=False, + use_cmap_RdBu=False, + verbose=False, + log_color_scale=False, + log_yaxis=False, + extra_title_txt=None, + n_job=-1, + sigdiff_list=None, + second_ref=None, + second_dev=None, + spcdb_dir=os.path.dirname(__file__), + sg_ref_path='', + sg_dev_path='', + ref_vert_params=None, + dev_vert_params=None, + **extra_plot_args + ): + """ + Creates 3x2 comparison zonal-mean plots for variables + common in two xarray Datasets. Optionally save to PDF. + + Args: + refdata: xarray dataset + Dataset used as reference in comparison + refstr: str + String description for reference data to be used in plots + devdata: xarray dataset + Dataset used as development in comparison + devstr: str + String description for development data to be used in plots + + Keyword Args (optional): + varlist: list of strings + List of xarray dataset variable names to make plots for + Default value: None (will compare all common 3D variables) + itime: integer + Dataset time dimension index using 0-based system + Default value: 0 + refmet: xarray dataset + Dataset containing ref meteorology + Default value: None + devmet: xarray dataset + Dataset containing dev meteorology + Default value: None + weightsdir: str + Directory path for storing regridding weights + Default value: None (will create/store weights in + current directory) + pdfname: str + File path to save plots as PDF + Default value: Empty string (will not create PDF) + cmpres: str + String description of grid resolution at which + to compare datasets + Default value: None (will compare at highest resolution + of Ref and Dev) + match_cbar: bool + Set this flag to True to use same the colorbar bounds + for both Ref and Dev plots. + Default value: True + pres_range: list of two integers + Pressure range of levels to plot [hPa]. The vertical axis + will span the outer pressure edges of levels that contain + pres_range endpoints. + Default value: [0, 2000] + normalize_by_area: bool + Set this flag to True to to normalize raw data in both + Ref and Dev datasets by grid area. Input ref and dev + datasets must include AREA variable in m2 if normalizing + by area. + Default value: False + enforce_units: bool + Set this flag to True force an error if the variables in + the Ref and Dev datasets have different units. + Default value: True + convert_to_ugm3: str + Whether to convert data units to ug/m3 for plotting. + Default value: False + flip_ref: bool + Set this flag to True to flip the vertical dimension of + 3D variables in the Ref dataset. + Default value: False + flip_dev: bool + Set this flag to True to flip the vertical dimension of + 3D variables in the Dev dataset. + Default value: False + use_cmap_RdBu: bool + Set this flag to True to use a blue-white-red colormap for + plotting raw reference and development datasets. + Default value: False + verbose: logical + Set this flag to True to enable informative printout. + Default value: False + log_color_scale: bool + Set this flag to True to enable plotting data (not diffs) + on a log color scale. + Default value: False + log_yaxis: bool + Set this flag to True if you wish to create zonal mean + plots with a log-pressure Y-axis. + Default value: False + extra_title_txt: str + Specifies extra text (e.g. a date string such as "Jan2016") + for the top-of-plot title. + Default value: None + n_job: int + Defines the number of simultaneous workers for parallel + plotting. Set to 1 to disable parallel plotting. + Value of -1 allows the application to decide. + Default value: -1 + sigdiff_list: list of str + Returns a list of all quantities having significant + differences (where |max(fractional difference)| > 0.1). + Default value: None + second_ref: xarray Dataset + A dataset of the same model type / grid as refdata, + to be used in diff-of-diffs plotting. + Default value: None + second_dev: xarray Dataset + A dataset of the same model type / grid as devdata, + to be used in diff-of-diffs plotting. + Default value: None + spcdb_dir: str + Directory containing species_database.yml file. + Default value: Path of GCPy code repository + sg_ref_path: str + Path to NetCDF file containing stretched-grid info + (in attributes) for the ref dataset + Default value: '' (will not be read in) + sg_dev_path: str + Path to NetCDF file containing stretched-grid info + (in attributes) for the dev dataset + Default value: '' (will not be read in) + ref_vert_params: list(AP, BP) of list-like types + Hybrid grid parameter A in hPa and B (unitless). + Needed if ref grid is not 47 or 72 levels. + Default value: None + dev_vert_params: list(AP, BP) of list-like types + Hybrid grid parameter A in hPa and B (unitless). + Needed if dev grid is not 47 or 72 levels. + Default value: None + extra_plot_args: various + Any extra keyword arguments are passed through the + plotting functions to be used in calls to pcolormesh() + (CS) or imshow() (Lat/Lon). + """ + +and generates a comparison plot such as: + +.. image:: _static/images/six\_panel\_zonal\_mean.png + :align: center + :width: 80% + +.. _plot-shared: + +Shared structure +---------------- Both :code:`compare_single_level()` and :code:`compare_zonal_mean()` have four positional (required) arguments. -Arguments: -~~~~~~~~~~ - .. option:: refdata : xarray.Dataset Dataset used as reference in comparison @@ -70,33 +439,16 @@ Arguments: String description for development data to be used in plots OR list containing [dev1str, dev2str] for diff-of-diffs plots -:option:`refstr` and :option:`devstr` title the top two panels of +:option:`refstr` and :option:`devstr` title the top two panels of each six panel plot. -A basic script that calls :code:`compare_zonal_mean()` or -:code:`compare_single_level()` looks like: +Functions :code:`compare_single_level()` and +:code:`compare_zonal_mean()` share many arguments. Some of these +arguments are plotting options that change the format of the plots: -.. code-block:: python - - #!/usr/bin/env python - - import xarray as xr - import gcpy.plot as gcplot - import matplotlib.pyplot as plt - - file1 = '/path/to/ref' - file2 = '/path/to/dev' - ds1 = xr.open_dataset(file1) - ds2 = xr.open_dataset(file2) - gcplot.compare_zonal_mean(ds1, 'Ref run', ds2, 'Dev run') - #gcplot.compare_single_level(ds1, 'Ref run', ds2, 'Dev run') - plt.show() - -:code`compare_single_level()` and :code:`compare_zonal_mean()` also -share many keyword arguments. Some of these arguments are plotting -options that change the format of the plots, e.g. choosing to convert -units to ug/m\ :sup:`3`, which are automatically handled if you do not -specify a value for that argument. +For example, you may wish to convert units to ug/m\ :sup:`3` when +generating comparison plots of aerosol species. Activate this +option by setting the keyword argument :literal:`convert_to_ugm3=True`. Other arguments are necessary to achieve a correct plot depending on the format of :literal:`refdata` and :literal:`devdata` and require @@ -105,1096 +457,258 @@ specify if one of the datasets should be flipped vertically if Z coordinates in that dataset do not denote decreasing pressure as Z index increases, otherwise the vertical coordinates between your two datasets may be misaligned and result in an undesired plotting -outcome. +outcome. This may be done with by setting the boolean options +:literal:`flip_ref=True` and/or :literal:`flip_dev=True`. The :literal:`n_job` argument governs the parallel plotting settings of :code:`compare_single_level()` and :code:`compare_zonal_mean()` . -GCPy uses the joblib library to create plots in parallel. Due to +GCPy uses the JobLib library to create plots in parallel. Due to limitations with matplotlib, this parallelization creates plots (pages) in parallel rather than individual panels on a single page. Parallel plot creation is not enabled when you do not save to a -PDF. The default value of :literal:`n_job=-1` allows the function call to -automatically scale up to, at most, the number of cores available on -your system. On systems with higher (12+) core counts, the max number of cores is not typically reached because of the process handling mechanics of joblib. -However, on lower-end systems with lower core counts or less available memory, it is advantageous to use :literal:`n_job` to limit the max number of processes. - -Shared keyword arguments: -~~~~~~~~~~~~~~~~~~~~~~~~~ - -.. option:: varlist : list of str - - List of xarray dataset variable names to make plots for - - Default value: None (will compare all common variables) - -.. option:: itime : int - - Dataset time dimension index using 0-based system. Can only plot - values from one time index in a single function call. - - Default value: 0 - -.. option:: refmet : xarray.Dataset - - Dataset containing ref meteorology. Needed for area-based - normalizations / ug/m3 unit conversions. - - Default value: None - -.. option:: devmet : xarray.Dataset - - Dataset containing dev meteorology. Needed for area-based - normalizations and/or ug/m3 unit conversions. - - Default value: None - -.. option:: weightsdir : str - - Directory path for storing regridding weight files generated by - xESMF. - - Default value: None (will create/store weights in current directory) - -.. option:: pdfname : str - - File path to save plots as PDF. - - Default value: Empty string (will not create PDF) - -.. option:: cmpres : str - - String description of grid resolution at which to compare - datasets. The possible formats are 'int' (e.g. '48' for c48) for - a cubed-sphere resolution or 'latxlon' (e.g. '4x5') for a - lat/lon resolution. - - Default value: None (will compare at highest resolution of Ref - and Dev) - -.. option:: match_cbar : bool - - Set this flag to True to use same the colorbar bounds for both - Ref and Dev plots. This only applies to the top two panels of - each plot. - - Default value: True - -.. option:: normalize_by_area : bool - - Set this flag to True to to normalize raw data in both Ref and - Dev datasets by grid area. Either input ref and dev datasets - must include AREA variable in m2 if normalizing by area, or - refmet and devmet datasets must include Met_AREAM2 variable. - - Default value: False - -.. option:: enforce_units : bool - - Set this flag to True force an error if the variables in the Ref - and Dev datasets have different units. - - Default value: True - -.. option:: convert_to_ugm3 : bool - - Whether to convert data units to ug/m3 for plotting. refmet and - devmet cannot be None if converting to ug/m3. - - Default value: False - -.. option:: flip_ref : bool - - Set this flag to True to flip the vertical dimension of 3D - variables in the Ref dataset. - - Default value: False - -.. option:: flip_dev : bool - - Set this flag to True to flip the vertical dimension of 3D - variables in the Dev dataset. - - Default value: False - -.. option:: use_cmap_RdBu : bool - - Set this flag to True to use a blue-white-red colormap for - plotting raw ref and dev data (the top two panels). - - Default value: False - -.. option:: verbose : bool - - Set this flag to True to enable informative printout. - - Default value: False - -.. option:: log_color_scale : bool - - Set this flag to True to enable plotting data (only the top two - panels, not diffs) on a log color scale. - - Default value: False - -.. option:: extra_title_txt : str - - Specifies extra text (e.g. a date string such as "Jan2016") - for the top-of-plot title. - - Default value: None - -.. option:: n_job : int - - Defines the number of simultaneous workers for parallel - plotting. Only applicable when saving to PDF. - Set to 1 to disable parallel plotting. Value of -1 allows the - application to decide. - - Default value: -1 - -.. option:: sigdiff_list : list of str - - Returns a list of all quantities having significant - differences (where \|max(fractional difference)\| > 0.1). - - Default value: [] - -.. option:: second_ref : xarray.Dataset - - A dataset of the same model type / grid as refdata, to be used - in diff-of-diffs plotting. - - Default value: None - -.. option:: second_dev : xarray.Dataset - - A dataset of the same model type / grid as devdata, to be used - in diff-of-diffs plotting. - - Default value: None - -.. option:: spcdb_dir : str - - Directory containing species_database.yml file. This file is - used for unit conversions to ug/m3. GEOS-Chem run directories - include a copy of this file which may be more up-to-date than - the version included with GCPy. - - Default value: Path of GCPy code repository - -.. option:: sg_ref_path : str - - Path to NetCDF file containing stretched-grid info (in - attributes) for the ref dataset. - - Default value: '' (will not be read in) - -.. option:: sg_dev_path : str +PDF. The default value of :literal:`n_job=-1` allows the function call +to automatically scale up to, at most, the number of cores available +on your system. - Path to NetCDF file containing stretched-grid info (in - attributes) for the dev dataset. +.. note:: - Default value: '' (will not be read in) + On systems with higher (12+) core counts, the maximum number of + cores is not typically reached because of the process handling + mechanics of JobLib. However, on lower-end systems with lower core + counts or less available memory, it is advantageous to use + :literal:`n_job` to limit the max number of processes. -compare_single_level --------------------- + Due to how Python handles memory management on Linux systems, using + more cores may result in memory not returned to the system after + the plots are created. Requesting fewer cores with + :literal:`n_job` may help to avoid this situation. -.. code-block:: python - - def compare_single_level(refdata, refstr, devdata, devstr, - varlist=None, ilev=0, itime=0, - refmet=None, devmet=None, weightsdir='.', - pdfname="", cmpres=None, match_cbar=True, - normalize_by_area=False, enforce_units=True, - convert_to_ugm3=False, flip_ref=False, flip_dev=False, - use_cmap_RdBu=False, verbose=False, log_color_scale=False, - extra_title_txt=None, extent = [-1000, -1000, -1000, -1000], - n_job=-1, sigdiff_list=[], second_ref=None, second_dev=None, - spcdb_dir=os.path.dirname(__file__), sg_ref_path='', sg_dev_path='', - ll_plot_func='imshow', **extra_plot_args - ): - - -:code:`compare_single_level()` features several keyword arguments that -are not relevant to :code:`compare_zonal_mean()`, including specifying -which level to plot, the lat/lon extent of the plots, and which -underlying :code:`matplotlib.plot` function to use for plotting. - -Function-specific keyword arguments: -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -.. option:: ilev : int - - Dataset level dimension index using 0-based system - - Default value: 0 - -.. option:: extent : list of float - - Defines the extent of the region to be plotted in form - [minlon, maxlon, minlat, maxlat]. Default value plots extent of input grids. - - Default value: [-1000, -1000, -1000, -1000] - -.. option:: ll_plot_func : str - - Function to use for lat/lon single level plotting with possible - values 'imshow' and 'pcolormesh'. imshow is much faster but is - slightly displaced when plotting from dateline to dateline - and/or pole to pole. +.. _plot-six-panel-example: - Default value: 'imshow' +Example script +-------------- -.. option:: **extra_plot_args - - Any extra keyword arguments are passed through the plotting - functions to be used in calls to :code:`pcolormesh()` (CS) or - :code:`imshow()` (Lat/Lon). - -compare_zonal_mean ------------------- - -.. code-block:: python - - def compare_zonal_mean(refdata, refstr, devdata, devstr, - varlist=None, itime=0, refmet=None, devmet=None, - weightsdir='.', pdfname="", cmpres=None, - match_cbar=True, pres_range=[0, 2000], - normalize_by_area=False, enforce_units=True, - convert_to_ugm3=False, flip_ref=False, flip_dev=False, - use_cmap_RdBu=False, verbose=False, log_color_scale=False, - log_yaxis=False, extra_title_txt=None, n_job=-1, sigdiff_list=[], - second_ref=None, second_dev=None, spcdb_dir=os.path.dirname(__file__), - sg_ref_path='', sg_dev_path='', ref_vert_params=[[],[]], - dev_vert_params=[[],[]], **extra_plot_args - ): - -:code:`compare_zonal_mean()` features several keyword arguments that -are not relevant to :code:`compare_single_level()`, including -specifying the pressure range to plot (defaulting to the complete -atmosphere), whether the y-axis of the plots (pressure) should be in -log format, and hybrid vertical grid parameters to pass if one or more -of Ref and Dev do not use the typical 72-level or 47-level grids. - -Function-specific keyword arguments: -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -.. option:: pres_range : list of ints - - Pressure range of levels to plot [hPa]. The vertical axis will - span the outer pressure edges of levels that contain pres_range - endpoints. - - Default value: [0,2000] - -.. option:: log_yaxis : bool - - Set this flag to True if you wish to create zonal mean - plots with a log-pressure Y-axis. - - Default value: False - -.. option:: ref_vert_params : list of list-like types - - Hybrid grid parameter A in hPa and B (unitless). Needed if ref - grid is not 47 or 72 levels. - - Default value: [[], []] - -.. option:: dev_vert_params : list of list-like types - - Hybrid grid parameter A in hPa and B (unitless). Needed if dev - grid is not 47 or 72 levels. - - Default value: [[], []] - -.. option:: **extra_plot_args - - Any extra keyword arguments are passed through the plotting - functions to be used in calls to :code:`pcolormesh()`. - -============ -Single_panel -============ - -.. code-block:: python - - def single_panel(plot_vals, ax=None, plot_type="single_level", - grid={}, gridtype="", title="fill",comap=WhGrYlRd, - norm=[],unit="",extent=(None, None, None, None), - masked_data=None,use_cmap_RdBu=False, - log_color_scale=False, add_cb=True, - pres_range=[0, 2000], pedge=np.full((1, 1), -1), - pedge_ind=np.full((1,1), -1), log_yaxis=False, - xtick_positions=[], xticklabels=[], proj=ccrs.PlateCarree(), - sg_path='', ll_plot_func="imshow", vert_params=[[],[]], - pdfname="", return_list_of_plots=False **extra_plot_args - ): - - -:code:`gcpy.plot.single_panel()` is used to create plots containing -only one panel of GEOS-Chem data. This function is used within -:code:`compare_single_level()` and :code:`compare_zonal_mean()` to -generate each panel plot. It can also be called directly on its own to -quickly plot GEOS-Chem data in zonal mean or single level format. +Here is a basic script that calls both :code:`compare_zonal_mean()` and +:code:`compare_single_level()`: .. code-block:: python #!/usr/bin/env python import xarray as xr - import gcpy.plot as gcplot import matplotlib.pyplot as plt + from gcpy.plot.compare_single_level import compare_single_level + from gcpy.plot.compare_zonal_mean import compare_zonal_mean - ds = xr.open_dataset('GEOSChem.SpeciesConc.20160701_0000z.nc4') - #get surface ozone - plot_data = ds['SpeciesConc_O3'].isel(lev=0) - - gcplot.single_panel(plot_data) + file1 = '/path/to/ref' + file2 = '/path/to/dev' + ds1 = xr.open_dataset(file1) + ds2 = xr.open_dataset(file2) + compare_zonal_mean(ds1, 'Ref run', ds2, 'Dev run') + plt.show() + compare_single_level(ds1, 'Ref run', ds2, 'Dev run') plt.show() -Currently :code:`single_panel()` expects data with a 1-length ( or -non-existent) time dimension, as well as a 1-length or non-existent Z -dimension for single level plotting, so you'll need to do some -pre-processing of your input data as shown in the above code snippet. - -:code:`single_panel()` contains a few amenities to help with plotting -GEOS-Chem data, including automatic grid detection for lat/lon or -standard cubed-sphere xarray :code:`DataArray`-s. You can also pass NumPy -arrays to plot, though you'll need to manually pass grid info in this -case. - -Arguments: ----------- - -In addition to the specific arguments listed below, any other keyword -arguments will be forwarded to :code:`matplotlib.pyplot.imshow()` / -:code:`matplotlib.pyplot.pcolormesh()`. - -.. option:: plot_vals : xarray.DataArray or numpy array - - Single data variable GEOS-Chem output to plot - -.. option:: ax : matplotlib axes - - Axes object to plot information - - Default value: None (Will create a new axes) - -.. option:: plot_type : str - - Either "single_level" or "zonal_mean" - - Default value: "single_level" - -.. option:: grid : dict - - Dictionary mapping plot_vals to plottable coordinates - - Default value: {} (will attempt to read grid from plot_vals) - -.. option:: gridtype : str - - "ll" for lat/lon or "cs" for cubed-sphere - - Default value: "" (will automatically determine from grid) - -.. option:: title : str - - Title to put at top of plot - - Default value: "fill" (will use name attribute of plot_vals - if available) - -.. option:: comap : matplotlib Colormap - - Colormap for plotting data values - - Default value: WhGrYlRd - -.. option:: norm : list - - List with range [0..1] normalizing color range for matplotlib methods - - Default value: [] (will determine from plot_vals) - -.. option:: unit : str - - Units of plotted data - - Default value: "" (will use units attribute of plot_vals if available) - -.. option:: extent : tuple (minlon, maxlon, minlat, maxlat) - - Describes minimum and maximum latitude and longitude of input data - - Default value: (None, None, None, None) (Will use full extent - of plot_vals if plot is single level. - -.. option:: masked_data : numpy array - - Masked area for avoiding near-dateline cubed-sphere plotting issues - - Default value: None (will attempt to determine from plot_vals) - -.. option:: use_cmap_RdBu : bool - - Set this flag to True to use a blue-white-red colormap - - Default value: False - -.. option:: log_color_scale : bool - - Set this flag to True to use a log-scale colormap - - Default value: False - -.. option:: add_cb : bool - - Set this flag to True to add a colorbar to the plot - - Default value: True - -.. option:: pres_range : list of int - - Range from minimum to maximum pressure for zonal mean plotting - - Default value: [0, 2000] (will plot entire atmosphere) - -.. option:: pedge : numpy array - - Edge pressures of vertical grid cells in plot_vals for zonal mean plotting - - Default value: np.full((1, 1), -1) (will determine automatically) - -.. option:: pedge_ind : numpy array - - Index of edge pressure values within pressure range in - plot_vals for zonal mean plotting - - Default value: np.full((1, 1), -1) (will determine automatically) - -.. option:: log_yaxis : bool - - Set this flag to True to enable log scaling of pressure in - zonal mean plots - - Default value: False - -.. option:: xtick_positions : list of float - - Locations of lat/lon or lon ticks on plot - - Default value: [] (will place automatically for zonal mean plots) - -.. option:: xticklabels : list of str - - Labels for lat/lon ticks - - Default value: [] (will determine automatically from xtick_positions) - -.. option:: sg_path : str - - Path to NetCDF file containing stretched-grid info (in - attributes) for plot_vals - - Default value: '' (will not be read in) - -.. option:: ll_plot_func : str - - Function to use for lat/lon single level plotting with - possible values 'imshow' and 'pcolormesh'. - imshow is much faster but is slightly displaced when plotting - from dateline to dateline and/or pole to pole. - - Default value: 'imshow' - -.. option:: vert_params : list(AP, BP) of list-like types - - Hybrid grid parameter A in hPa and B (unitless). Needed if - grid is not 47 or 72 levels. - - Default value: [[], []] - -.. option:: pdfname : str - - File path to save plots as PDF - - Default value: "" (will not create PDF) - -.. option:: extra_plot_args : various - - Any extra keyword arguments are passed to calls to - pcolormesh() (CS) or imshow() (Lat/Lon). - - -Function-specific return value: -------------------------------- - -:code:`single_panel()` returns the following object: - -.. option:: plot : matplotlib plot - - Plot object created from input - -.. _plot-bmk: - -============================ -Benchmark Plotting Functions -============================ - -:code:`gcpy.benchmark` contains several functions for plotting -GEOS-Chem output in formats requested by the GEOS-Chem Steering -Committee. The primary use of these functions is to create plots of -most GEOS-Chem output variables divided into specific categories, -e.g. species categories such as Aerosols or Bromine for the -SpeciesConc diagnostic. In each category, these functions create -single level PDFs for the surface and 500hPa and zonal mean PDFs for -the entire a tmosphere and only the stratosphere (defined a 1-100hPa). -For :code:`make_benchmark_emis_plots()`, only single level plots at -the surface are produced. -All of these plotting functions include bookmarks within the generated PDFs that point to the pages containing each plotted quantity. -Thus these functions serve as tools for quickly creating comprehensive plots comparing two GEOS-Chem runs. These functions are used to create -the publicly available plots for 1-month and 1-year benchmarks of new versions of GEOS-Chem. - -Many of these functions use pre-defined (via YAML files included in GCPy) lists of variables. If one dataset includes a variable but the other dataset does not, -the data for that variable in the latter dataset will be considered to be NaN and will be plotted as such. - -Shared structure of benchmark functions ---------------------------------------- - -Each of the :code:`gcpy.benchmark.make_benchmark_*_plots()` functions -requires 4 arguments to specify the ref and dev datasets. - -Shared arguments: -~~~~~~~~~~~~~~~~~ - -.. option:: ref: str - - Path name for the "Ref" (aka "Reference") data set. - -.. option:: refstr : str - - A string to describe ref (e.g. version number) - -.. option:: dev : str - - Path name for the "Dev" (aka "Development") data set. - This data set will be compared against the "Reference" data set. - -.. option:: devstr : str - - A string to describe dev (e.g. version number) - -Note that the :literal:`ref` and :literal:`dev` arguments in -:code:`make_benchmark_*_plots()` are the paths to NetCDF files, rather -than xarray Datasets as in :code:`compare_single_level()` and -:code:`compare_zonal_mean()`. The :code:`make_benchmark_*_plots()` -functions internally open these files as xarray Datasets and pass -those datasets to :code:`compare_single_level()` and -:code:`compare_zonal_mean()`. - -The benchmark plotting functions share several keyword -arguments. Keyword arguments that do not share the same purpose across -benchmark plotting functions have :literal:`NOTE:` in the description. - -Shared keyword arguments: -~~~~~~~~~~~~~~~~~~~~~~~~~ - -.. option:: dst : str - - A string denoting the destination folder where a - PDF file containing plots will be written. - - Default value: ./benchmark. - -.. option:: subdst : str - - A string denoting the sub-directory of dst where PDF - files containing plots will be written. In practice, - subdst is only needed for the 1-year benchmark output, - and denotes a date string (such as "Jan2016") that - corresponds to the month that is being plotted. - NOTE: Not available in wetdep_plots - - Default value: None - -.. option:: overwrite : bool - - Set this flag to True to overwrite previously created files in the - destination folder (specified by the dst argument). - - Default value: False. - -.. option:: verbose : bool - - Set this flag to True to print extra informational output. - - Default value: False. - -.. option:: log_color_scale: bool - - Set this flag to True to enable plotting data (the top two panels - of each plot, not diffs) on a log color scale. - - Default value: False - -.. option:: sigdiff_files : list of str - - Filenames that will contain the list of quantities having - significant differences between datasets. Three files are used: - one for surface, one for 500hPa, and one for zonal mean. - These lists are needed in order to fill out the benchmark - approval forms. - - .. note:: Not available in wetdep_plots - - Default value: None - -.. option:: spcdb_dir : str - - Directory containing species_database.yml file. This file is - used for unit conversions to ug/m3. GEOS-Chem run directories - include a copy of this file which may be more up-to-date than - the version included with GCPy. - - Default value: Path of GCPy code repository - -.. option:: weightsdir : str - - Directory in which to place (and possibly reuse) xESMF regridder - netCDF files. - - Default value: '.' - -.. option:: n_job : int - - Defines the number of simultaneous workers for parallel plotting. - Set to 1 to disable parallel plotting. Value of -1 allows the - application to decide. - - .. note:: - - In :code:`make_benchmark_conc_plots()`, parallelization occurs - at the species category level. In all other functions, - parallelization occurs within calls to :code:`compare_single_level()` - and :code:`compare_zonal_mean()`. - - Default value: -1 in :code:`make_benchmark_conc_plots`, 1 in all others - -make_benchmark_aod_plots ------------------------- - -.. code-block:: python - - def make_benchmark_aod_plots(ref, refstr, dev, devstr, varlist=None, - dst="./benchmark", subdst=None, overwrite=False, verbose=False, - log_color_scale=False, sigdiff_files=None, weightsdir='.', n_job=-1, - spcdb_dir=os.path.dirname(__file__) - ): - - """ - Creates PDF files containing plots of column aerosol optical - depths (AODs) for model benchmarking purposes. - """ - -Function-specific keyword args: - -.. option:: varlist : list of str - - List of AOD variables to plot. If not passed, then all - AOD variables common to both Dev and Ref will be plotted. - Use the varlist argument to restrict the number of - variables plotted to the pdf file when debugging. - - Default value: None - - -This function creates column optical depth plots using the Aerosols diagnostic output. - -make_benchmark_conc_plots -------------------------- - -.. code-block:: python - - def make_benchmark_conc_plots(ref, refstr, dev, devstr, dst="./benchmark", - subdst=None, overwrite=False, verbose=False, collection="SpeciesConc", - benchmark_type="FullChemBenchmark", plot_by_spc_cat=True, restrict_cats=[], - plots=["sfc", "500hpa", "zonalmean"], use_cmap_RdBu=False, log_color_scale=False, - sigdiff_files=None, normalize_by_area=False, cats_in_ugm3=["Aerosols", "Secondary_Organic_Aerosols"], - areas=None, refmet=None, devmet=None, weightsdir='.', n_job=-1, second_ref=None - second_dev=None, spcdb_dir=os.path.dirname(__file__) - ): - - """ - Creates PDF files containing plots of species concentration - for model benchmarking purposes. - """ - -Function-specific keyword arguments: -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -.. option:: collection : str - - Name of collection to use for plotting. - - Default value: "SpeciesConc" - -.. option:: benchmark_type: str - - A string denoting the type of benchmark output to plot, - either FullChemBenchmark or TransportTracersBenchmark. - - Default value: "FullChemBenchmark" - -.. option:: plot_by_spc_cat: logical - - Set this flag to False to send plots to one file rather - than separate file per category. - - Default value: True - -.. option:: restrict_cats : list of str - - List of benchmark categories in benchmark_categories.yml to make - plots for. If empty, plots are made for all categories. - - Default value: empty - -.. option:: plots : list of str - - List of plot types to create. - - Default value: ['sfc', '500hpa', 'zonalmean'] - -.. option:: normalize_by_area: bool - - Set this flag to true to enable normalization of data - by surfacea area (i.e. kg s-1 --> kg s-1 m-2). - - Default value: False - -.. option:: cats_in_ugm3: list of str - - List of benchmark categories to to convert to ug/m3 - - Default value: ["Aerosols", "Secondary_Organic_Aerosols"] - -.. option:: areas : dict of xarray DataArray: - - Grid box surface areas in m2 on Ref and Dev grids. - - Default value: None - -.. option:: refmet : str - - Path name for ref meteorology - - Default value: None - -.. option:: devmet : str - - Path name for dev meteorology - - Default value: None - -.. option:: second_ref: str - - Path name for a second "Ref" (aka "Reference") data set for - diff-of-diffs plotting. This dataset should have the same model - type and grid as ref. - - Default value: None - -.. option:: second_dev: str - - Path name for a second "Ref" (aka "Reference") data set for - diff-of-diffs plotting. This dataset should have the same model - type and grid as ref. - - Default value: None - - -This function creates species concentration plots using the -:literal:`SpeciesConc` diagnostic output by default. This function is the -only benchmark plotting function that supports diff-of-diffs plotting, -in which 4 datasets are passed and the differences between two groups -of Ref datasets vs. two groups of Dev datasets is plotted (typically -used for comparing changes in GCHP vs. changes in GEOS-Chem Classic -across model versions). This is also the only benchmark plotting -function that sends plots to separate folders based on category -(as denoted by the :literal:`plot_by_spc_cat` flag). The full list of -species categories is denoted in :file:`benchmark_categories.yml` -(included in GCPy) as follows: - -.. code-block:: python +.. _plot-single-panel: - """ - FullChemBenchmark: - Aerosols: - Dust: DST1, DST2, DST3, DST4 - Inorganic: NH4, NIT, SO4 - OC_BC: BCPI, BCPO, OCPI, OCPO - SOA: Complex_SOA, Simple_SOA - Sea_Salt: AERI, BrSALA, BrSALC, ISALA, ISALC, NITs, - SALA, SALAAL, SALACL, SALC, SALCAL, SALCCL, SO4s - Bromine: Bry, BrOx, Br, Br2, BrCl, BrNO2, BrNO3, BrO, - CH3Br, CH2Br2, CHBr3, HOBr, HBr - Chlorine: Cly, ClOx, Cl, ClO, Cl2, Cl2O2, ClOO, ClNO2, ClNO3, - CCl4, CFCs, CH3Cl, CH2Cl2, CH3CCl3, CHCl3, HOCl, HCl, Halons, HCFCs, OClO - Iodine: Iy, IxOy, I, I2, IBr, ICl, IO, ION, IONO2, CH3I, CH2I2, - CH2ICl, CH2IBr, HI, HOI, OIO - Nitrogen: NOy, NOx, HNO2, HNO3, HNO4, MPAN, NIT, 'NO', NO2, NO3, - N2O5, MPN, PAN, PPN, N2O, NHx, NH3, NH4, MENO3, ETNO3, IPRNO3, NPRNO3 - Oxidants: O3, CO, OH, NOx - Primary_Organics: - Alcohols: EOH, MOH - Biogenics: ISOP, MTPA, MTPO, LIMO - HCs: ALK4, BENZ, CH4, C2H6, C3H8, PRPE, TOLU, XYLE - ROy: H2O2, H, H2, H2O, HO2, O1D, OH, RO2 - Secondary_Organic_Aerosols: - Complex_SOA: TSOA0, TSOA1, TSOA2, TSOA3, ASOA1, ASOA2, ASOA3, - ASOAN, TSOG0, TSOG1, TSOG2, TSOG3, ASOG1, ASOG2, ASOG3 - Isoprene_SOA: INDIOL, LVOCOA, SOAIE, SOAGX - Simple_SOA: SOAP, SOAS - Secondary_Organics: - Acids: ACTA - Aldehydes: ALD2, CH2O, HPALDs, MACR - Epoxides: IEPOX - Ketones: ACET, MEK, MVK - Nitrates: ISOPN - Other: GLYX, HCOOH, MAP, RCHO - Peroxides: MP - Sulfur: SOx, DMS, OCS, SO2, SO4 - TransportTracersBenchmark: - RnPbBeTracers: Rn222, Pb210, Pb210s, Be7, Be7s, Be10, Be10s - TransportTracers: PassiveTracer, SF6, CH3I, aoa, aoa_bl, aoa_nh, - CO_25, CO_50, e90, e90_s, e90_n, st80_25, stOX +================== +Single panel plots +================== - """ +Function :code:`single_panel()` (contained in GCPy module +:code:`gcpy.plot.single_panel`) is used to create plots containing +only one panel of GEOS-Chem data. This function is used within +:code:`compare_single_level()` and :code:`compare_zonal_mean()` to +generate each panel plot. It can also be called directly on its +own to quickly plot GEOS-Chem data in zonal mean or single level +format. +.. _plot-single-panel-func: -make_benchmark_emis_plots +Function: :code:`single_panel` +------------------------------ +Function :code:`single_panel()` accepts the following arguments: .. code-block:: python - def make_benchmark_emis_plots(ref, refstr, dev, devstr, dst="./benchmark", - subdst=None, plot_by_spc_cat=False, plot_by_hco_cat=False, overwrite=False, - verbose=False, flip_ref=False, flip_dev=False, log_color_scale=False, - sigdiff_files=None, weightsdir='.', n_job=-1, spcdb_dir=os.path.dirname(__file__) + def single_panel( + plot_vals, + ax=None, + plot_type="single_level", + grid=None, + gridtype="", + title="fill", + comap=WhGrYlRd, + norm=None, + unit="", + extent=None, + masked_data=None, + use_cmap_RdBu=False, + log_color_scale=False, + add_cb=True, + pres_range=None, + pedge=np.full((1, 1), -1), + pedge_ind=np.full((1, 1), -1), + log_yaxis=False, + xtick_positions=None, + xticklabels=None, + proj=ccrs.PlateCarree(), + sg_path='', + ll_plot_func="imshow", + vert_params=None, + pdfname="", + weightsdir='.', + vmin=None, + vmax=None, + return_list_of_plots=False, + **extra_plot_args ): - """ - Creates PDF files containing plots of emissions for model - benchmarking purposes. This function is compatible with benchmark - simulation output only. It is not compatible with transport tracers - emissions diagnostics. - - Remarks: - -------- - (1) If both plot_by_spc_cat and plot_by_hco_cat are - False, then all emission plots will be placed into the - same PDF file. - - (2) Emissions that are 3-dimensional will be plotted as - column sums. - """ - -Function-specific keyword args: -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -.. option:: plot_by_spc_cat : bool - - Set this flag to True to separate plots into PDF files - according to the benchmark species categories (e.g. Oxidants, - Aerosols, Nitrogen, etc.) These categories are specified - in the YAML file benchmark_species.yml. - - Default value: False - -.. option:: plot_by_hco_cat : bool - - Set this flag to True to separate plots into PDF files - according to HEMCO emissions categories (e.g. Anthro, - Aircraft, Bioburn, etc.) - - Default value: False - -.. option:: flip_ref : bool - - Set this flag to True to reverse the vertical level - ordering in the "Ref" dataset (in case "Ref" starts - from the top of atmosphere instead of the surface). + """ + Core plotting routine -- creates a single plot panel. + + Args: + plot_vals: xarray.DataArray, numpy.ndarray, or dask.array.Array + Single data variable GEOS-Chem output to plot + + Keyword Args (Optional): + ax: matplotlib axes + Axes object to plot information + Default value: None (Will create a new axes) + plot_type: str + Either "single_level" or "zonal_mean" + Default value: "single_level" + grid: dict + Dictionary mapping plot_vals to plottable coordinates + Default value: {} (will attempt to read grid from plot_vals) + gridtype: str + "ll" for lat/lon or "cs" for cubed-sphere + Default value: "" (will automatically determine from grid) + title: str + Title to put at top of plot + Default value: "fill" (will use name attribute of plot_vals + if available) + comap: matplotlib Colormap + Colormap for plotting data values + Default value: WhGrYlRd + norm: list + List with range [0..1] normalizing color range for matplotlib + methods. Default value: None (will determine from plot_vals) + unit: str + Units of plotted data + Default value: "" (will use units attribute of plot_vals + if available) + extent: tuple (minlon, maxlon, minlat, maxlat) + Describes minimum and maximum latitude and longitude of input + data. Default value: None (Will use full extent of plot_vals + if plot is single level). + masked_data: numpy array + Masked area for avoiding near-dateline cubed-sphere plotting + issues Default value: None (will attempt to determine from + plot_vals) + use_cmap_RdBu: bool + Set this flag to True to use a blue-white-red colormap + Default value: False + log_color_scale: bool + Set this flag to True to use a log-scale colormap + Default value: False + add_cb: bool + Set this flag to True to add a colorbar to the plot + Default value: True + pres_range: list(int) + Range from minimum to maximum pressure for zonal mean + plotting. Default value: [0, 2000] (will plot entire + atmosphere) + pedge: numpy array + Edge pressures of vertical grid cells in plot_vals + for zonal mean plotting. Default value: np.full((1, 1), -1) + (will determine automatically) + pedge_ind: numpy array + Index of edge pressure values within pressure range in + plot_vals for zonal mean plotting. + Default value: np.full((1, 1), -1) (will determine + automatically) + log_yaxis: bool + Set this flag to True to enable log scaling of pressure in + zonal mean plots. Default value: False + xtick_positions: list(float) + Locations of lat/lon or lon ticks on plot + Default value: None (will place automatically for + zonal mean plots) + xticklabels: list(str) + Labels for lat/lon ticks + Default value: None (will determine automatically from + xtick_positions) + proj: cartopy projection + Projection for plotting data + Default value: ccrs.PlateCarree() + sg_path: str + Path to NetCDF file containing stretched-grid info + (in attributes) for plot_vals. + Default value: '' (will not be read in) + ll_plot_func: str + Function to use for lat/lon single level plotting with + possible values 'imshow' and 'pcolormesh'. imshow is much + faster but is slightly displaced when plotting from dateline + to dateline and/or pole to pole. Default value: 'imshow' + vert_params: list(AP, BP) of list-like types + Hybrid grid parameter A in hPa and B (unitless). Needed if + grid is not 47 or 72 levels. Default value: None + pdfname: str + File path to save plots as PDF + Default value: "" (will not create PDF) + weightsdir: str + Directory path for storing regridding weights + Default value: "." (will store regridding files in + current directory) + vmin: float + minimum for colorbars + Default value: None (will use plot value minimum) + vmax: float + maximum for colorbars + Default value: None (will use plot value maximum) + return_list_of_plots: bool + Return plots as a list. This is helpful if you are using + a cubedsphere grid and would like access to all 6 plots + Default value: False + extra_plot_args: various + Any extra keyword arguments are passed to calls to + pcolormesh() (CS) or imshow() (Lat/Lon). + + Returns: + plot: matplotlib plot + Plot object created from input + """ + +Function :code:`single_panel()` expects data with a 1-length (or +non-existent) :literal:`T` (time) dimension, as well as a +1-length or non-existent :literal:`Z` (vertical level) dimension. - Default value: False - -.. option:: flip_dev : bool - - Set this flag to True to reverse the vertical level - ordering in the "Dev" dataset (in case "Dev" starts - from the top of atmosphere instead of the surface). - - Default value: False - - -This function generates plots of total emissions using output from ``HEMCO_diagnostics`` (for GEOS-Chem Classic) and/or ``GCHP.Emissions`` output files. - - -make_benchmark_jvalue_plots ---------------------------- - -.. code-block:: python - - def make_benchmark_jvalue_plots(ref, refstr, dev, devstr, varlist=None, - dst="./benchmark", subdst=None, local_noon_jvalues=False, - plots=["sfc", "500hpa", "zonalmean"],overwrite=False, verbose=False, - flip_ref=False, flip_dev=False, log_color_scale=False, sigdiff_files=None, - weightsdir='.', n_job=-1, spcdb_dir=os.path.dirname(__file__) - ): - """ - Creates PDF files containing plots of J-values for model - benchmarking purposes. - - Remarks: - -------- - Will create 4 files containing J-value plots: - (1 ) Surface values - (2 ) 500 hPa values - (3a) Full-column zonal mean values. - (3b) Stratospheric zonal mean values - These can be toggled on/off with the plots keyword argument. - - At present, we do not yet have the capability to split the - plots up into separate files per category (e.g. Oxidants, - Aerosols, etc.). This is primarily due to the fact that - we archive J-values from GEOS-Chem for individual species - but not family species. We could attempt to add this - functionality later if there is sufficient demand. - """ - - -Function-specific keyword args: -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -.. option:: varlist : list of str - - List of J-value variables to plot. If not passed, - then all J-value variables common to both dev - and ref will be plotted. The varlist argument can be - a useful way of restricting the number of variables - plotted to the pdf file when debugging. - - Default value: None - -.. option:: local_noon_jvalues : bool - - Set this flag to plot local noon J-values. This will - divide all J-value variables by the JNoonFrac counter, - which is the fraction of the time that it was local noon - at each location. - - Default value: False - -.. option:: plots : list of strings - - List of plot types to create. - - Default value: ['sfc', '500hpa', 'zonalmean'] - -.. option:: flip_ref : bool - - Set this flag to True to reverse the vertical level - ordering in the "Ref" dataset (in case "Ref" starts - from the top of atmosphere instead of the surface). - - Default value: False - -.. option:: flip_dev : bool - - Set this flag to True to reverse the vertical level - ordering in the "Dev" dataset (in case "Dev" starts - from the top of atmosphere instead of the surface). - - Default value: False - - -This function generates plots of J-values using the -:literal:`JValues` GEOS-Chem output files. +:code:`single_panel()` contains a few amenities to help with plotting +GEOS-Chem data, including automatic grid detection for lat/lon or +standard cubed-sphere xarray :code:`DataArray`-s. You can also pass NumPy +arrays to plot, though you'll need to manually pass grid info in this +case (with the :literal:`gridtype`, :literal:`pedge`, and +:literal:`pedge_ind` keyword arguments). -make_benchmark_wetdep_plots ---------------------------- +The sample script shown below shows how you can data at a single level and +timestep from an :literal:`xarray.DataArray` object. .. code-block:: python - def make_benchmark_wetdep_plots(ref, refstr, dev, devstr, collection, - dst="./benchmark", datestr=None, overwrite=False, verbose=False, - benchmark_type="TransportTracersBenchmark", plots=["sfc", "500hpa", "zonalmean"], - log_color_scale=False, normalize_by_area=False, areas=None, refmet=None, - devmet=None, weightsdir='.', n_job=-1, spcdb_dir=os.path.dirname(__file__) - ): - """ - Creates PDF files containing plots of species concentration - for model benchmarking purposes. - """ - -Function-specific keyword args: -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -.. option:: datestr : str - - A string with date information to be included in both the - plot pdf filename and as a destination folder subdirectory - for writing plots - - Default value: None - -.. option:: benchmark_type: str - - A string denoting the type of benchmark output to plot, - either FullChemBenchmark or TransportTracersBenchmark. - - Default value: "FullChemBenchmark" - -.. option:: plots : list of strings - - List of plot types to create. - - Default value: ['sfc', '500hpa', 'zonalmean'] - -.. option:: normalize_by_area: bool - - Set this flag to true to enable normalization of data - by surfacea area (i.e. kg s-1 --> kg s-1 m-2). - - Default value: False - -.. option:: areas : dict of xarray DataArray: - - Grid box surface areas in m2 on Ref and Dev grids. - - Default value: None - -.. option:: refmet : str - - Path name for ref meteorology - - Default value: None - -.. option:: devmet : str - - Path name for dev meteorology + #!/usr/bin/env python - Default value: None + import xarray as xr + import matplotlib.pyplot as plt + from gcpy.plot.single_panel import single_panel -This function generates plots of wet deposition using -:literal:`WetLossConv` and :literal:`WetLossLS` GEOS-Chem output files. -It is currently primarily used for 1-Year Transport Tracer benchmarks, -plotting values for the following species as defined in -:file:`benchmark_categories.yml`: + # Read data from a file into an xr.Dataset object + dset = xr.open_dataset('GEOSChem.SpeciesConc.20160701_0000z.nc4') -.. code-block:: python + # Extract ozone (v/v) from the xr.Dataset object, + # for time=0 (aka first timestep) and lev=0 (aka surface) + sfc_o3 = dset['SpeciesConcVV_O3'].isel(time=0).isel(lev=0) - """ - WetLossConv: Pb210, Pb210Strat, Be7, Be7Strat, Be10, Be10Strat - WetLossLS: Pb210, Pb210Strat, Be7, Be7Strat, Be10, Be10Strat - """ + # Plot the data! + single_panel(sfc_o3) + plt.show() diff --git a/docs/source/Regridding.rst b/docs/source/Regridding.rst index 439977d7..4a4eac4e 100644 --- a/docs/source/Regridding.rst +++ b/docs/source/Regridding.rst @@ -51,6 +51,12 @@ directly from the command line (:ref:`see the examples below `) or from a Python script or interpreter (:code:`gcpy.file_regrid.file_regrid()`) +.. note:: + + For regridding to or from GCHP stretched-grid restart files, we + recommend using the :ref:`offline regridding weights method + `. + The syntax of :code:`file_regrid` is as follows: .. code-block:: python @@ -70,7 +76,7 @@ The syntax of :code:`file_regrid` is as follows: Regrids an input file to a new horizontal grid specification and saves it as a new file. """ -s + gcpy.file_regrid required arguments: ------------------------------------ diff --git a/docs/source/Single_panel b/docs/source/Single_panel deleted file mode 100644 index e69de29b..00000000 diff --git a/docs/source/Single_panel.rst b/docs/source/Single_panel.rst index 3a7d45d1..37e1d213 100644 --- a/docs/source/Single_panel.rst +++ b/docs/source/Single_panel.rst @@ -1,77 +1,123 @@ +.. _single-panel: + +This example script may also be found at `gcpy/examples/plotting/plot_single_panel.py `_. + +##################### Single Panel Plotting -===================== +##################### .. code-block:: python - #!/usr/bin/env python """ Global and Regional Single Panel Plots -------------------------------------- - This example script demonstrates the core single panel plotting capabilities of GCPy, - including global and regional single level plots as well as global zonal mean plots. - The example data described here is in lat/lon format, but the same code works equally - well for cubed-sphere (GCHP) data. + This example script demonstrates the core single panel plotting + capabilities of GCPy, including global and regional single level plots + as well as global zonal mean plots. + + The example data described here is in lat/lon format, but the same code + works equally well for cubed-sphere (GCHP) data. + + For full documentation on the plotting capabilities of GCPy + (including full argument lists), please see the GCPy documentation + at https://gcpy.readthedocs.io. """ - - #xarray allows us to read in any NetCDF file, the format of most GEOS-Chem diagnostics, - #as an xarray Dataset import xarray as xr - ds = xr.open_dataset('GEOSChem.Restart.20160701_0000z.nc4') - - #You can easily view the variables available for plotting using xarray. - #Each of these variables has its own xarray DataArray within the larger Dataset container. - print(ds.data_vars) - - #Most variables have some sort of prefix; in this example all variables are - #prefixed with 'SpeciesRst_'. We'll select the DataArray for ozone. - da = ds.SpeciesRst_O3 - - #Printing a DataArray gives a summary of the dimensions and attributes of the data. - print(da) - #This Restart file has a time dimension of size 1, with 72 vertical levels, - #46 latitude indicies, and 72 longitude indices. - import gcpy.plot as gcplot - - """ - Single level plots - ------------------ - """ - #gcpy.single_panel is the core plotting function of GCPy, able to create a one panel zonal mean or - #single level plot. Here we will create a single level plot of ozone at ~500 hPa. - #We must manually index into the level that we want to plot (index 22 in the standard 72-layer - #and 47-layer GMAO vertical grids). - slice_500 = da.isel(lev=22) - - #single_panel has many arguments which can be optionally specified. The only argument you must always - #pass to a call to single_panel is the DataArray that you want to plot. - #By default, the created plot includes a colorbar with units read from the DataArray, an automatic title - #(the data variable name in the DataArray), and an extent equivalent to the full lat/lon extent of the DataArray import matplotlib.pyplot as plt - gcplot.single_panel(slice_500) - plt.show() - - #You can specify a specific area of the globe you would like plotted using the 'extent' argument, - #which uses the format [min_longitude, max_longitude, min_latitude, max_latitude] with bounds [-180, 180, -90, 90] - gcplot.single_panel(slice_500, extent=[50, -90, -10, 60]) - plt.show() - - #Other commonly used arguments include specifying a title and a colormap (defaulting to a White-Green-Yellow-Red colormap) - #You can find more colormaps at https://matplotlib.org/tutorials/colors/colormaps.html - gcplot.single_panel(slice_500, title='500mb Ozone over the North Pacific', comap = plt.cm.viridis, - log_color_scale=True, extent=[80, -90, -10, 60]) - plt.show() - - """ - Zonal Mean Plotting - ------------------- - """ - - #Use the plot_type argument to specify zonal_mean plotting - gcplot.single_panel(da, plot_type="zonal_mean") - plt.show() - - #You can specify pressure ranges in hPa for zonal mean plot (by default every vertical level is plotted) - gcplot.single_panel(da, pres_range=[0, 100], log_yaxis=True, log_color_scale=True) - plt.show() - + from gcpy.plot.single_panel import single_panel + + + def main(): + """ + Example routine to create single panel plots. + """ + + # xarray allows us to read in any NetCDF file, the format of + # GEOS-Chem diagnostics as an xarray Dataset + dset = xr.open_dataset('GEOSChem.Restart.20160701_0000z.nc4') + + # You can easily view the variables available for plotting + # using xarray. Each of these variables has its own xarray + # DataArray within the larger Dataset container. + print(dset.data_vars) + + # Most variables have some sort of prefix; in this example all + # variables are prefixed with 'SpeciesRst_'. We'll select the + # DataArray for ozone. + darr = dset.SpeciesRst_O3 + + # Printing a DataArray gives a summary of the dimensions and attributes + # of the data. + print(darr) + + # This Restart file has a time dimension of size 1, with 72 vertical levels, + #46 latitude indicies, and 72 longitude indices. + + + # ================== + # Single-level Plots + # ================== + + # gcpy.single_panel is the core plotting function of GCPy, able to + # create a one panel zonal mean or single level plot. Here we will + # create a single level plot of ozone at ~500 hPa. We must manually + # index into the level that we want to plot (index 22 in the standard + # 72-layer and 47-layer GMAO vertical grids). + slice_500 = darr.isel(lev=22) + + # single_panel has many arguments which can be optionally specified. + # The only argument you must always pass to a call to single_panel is + # the DataArray that you want to plot. By default, the created plot + # includes a colorbar with units read from the DataArray, an + # automatic title (the data variable name in the DataArray), and + # an extent equivalent to the full lat/lon extent of the DataArray + single_panel(slice_500) + plt.show() + + #You can specify a specific area of the globe you would like plotted + # using the 'extent' argument, which uses the format [min_longitude, + # max_longitude, min_latitude, max_latitude] with bounds + # [-180, 180, -90, 90] + single_panel(slice_500, extent=[50, -90, -10, 60]) + plt.show() + + # Other commonly used arguments include specifying a title and a + # colormap (defaulting to a White-Green-Yellow-Red colormap) + #You can find more colormaps at + # https://matplotlib.org/tutorials/colors/colormaps.html + single_panel( + slice_500, + title='500mb Ozone over the North Pacific', + comap=plt.get_cmap("viridis"), + log_color_scale=True, + extent=[80, -90, -10, 60] + ) + plt.show() + + # =================== + # Zonal Mean Plotting + # =================== + + # Use the plot_type argument to specify zonal_mean plotting + single_panel( + darr, + plot_type="zonal_mean" + ) + plt.show() + + #You can specify pressure ranges in hPa for zonal mean plot + # (by default every vertical level is plotted) + single_panel( + darr, + pres_range=[0, 100], + log_yaxis=True, + log_color_scale=True + ) + plt.show() + + + + # Only execute when we run as a standalone script + if __name__ == '__main__': + main() diff --git a/docs/source/Six_panel.rst b/docs/source/Six_panel.rst index 240a787c..553cc0f9 100644 --- a/docs/source/Six_panel.rst +++ b/docs/source/Six_panel.rst @@ -1,61 +1,119 @@ +.. _six-panel: + +################## Six Panel Plotting -================== +################## -.. code-block:: python +This example script may also be found at `gcpy/exampls/plotting/plot_comparisons.py `_. +.. code-block:: python #!/usr/bin/env python """ Six Panel Comparison Plots -------------------------------------- - This example script demonstrates the comparitive plotting capabilities of GCPy, - including single level plots as well as global zonal mean plots. - These comparison plots are frequently used to evaluate results from different runs / versions - of GEOS-Chem, but can also be used to compare results from different points in one run that - are stored in separate xarray datasets. - The example data described here is in lat/lon format, but the same code works equally - well for cubed-sphere (GCHP) data. + This example script demonstrates the comparitive plotting + capabilities of GCPy, including single level plots as well as + global zonal mean plots. These comparison plots are frequently + used to evaluate results from different runs / versions of + GEOS-Chem, but can also be used to compare results from different + points in one run that are stored in separate xarray datasets. + + The example data described here is in lat/lon format, but the same + code works equally well for cubed-sphere (GCHP) data. """ - - #xarray allows us to read in any NetCDF file, the format of most GEOS-Chem diagnostics, - #as an xarray Dataset import xarray as xr - ref_ds = xr.open_dataset('first_run/GEOSChem.Restart.20160801_0000z.nc4') - dev_ds = xr.open_dataset('second_run/GEOSChem.Restart.20160801_0000z.nc4') - - import gcpy.plot as gcplot - - """ - Single level plots - ------------------ - """ - - #compare_single_level generates sets of six panel plots for data at a specified level in your datasets. - #By default, the level at index 0 (likely the surface) is plotted. Here we will plot data at ~500 hPa, - #which is located at index 21 in the standard 72-level and 47-level GMAO vertical grids. - ilev=21 - - #You likely want to look at the same variables across both of your datasets. If a variable is in - #one dataset but not the other, the plots will show NaN values for the latter. - #You can pass variable names in a list to these comparison plotting functions (otherwise all variables will plot). - varlist = ['SpeciesRst_O3', 'SpeciesRst_CO2'] - - #compare_single_level has many arguments which can be optionally specified. The first four arguments are required. - #They specify your first xarray Dataset, the name of your first dataset, your second xarray Dataset, and the name of - #your second dataset. Here we will also pass a specific level and the names of the variables you want to plot. import matplotlib.pyplot as plt - gcplot.compare_single_level(ref_ds, 'Dataset 1', dev_ds, 'Dataset 2', ilev=ilev, varlist=varlist) - plt.show() - - #Using plt.show(), you can view the plots interactively. You can also save out the plots to a PDF. - gcplot.compare_single_level(ref_ds, 'Dataset 1', dev_ds, 'Dataset 2', ilev=ilev, varlist=varlist, pdfname='single_level.pdf') - - """ - Zonal Mean Plotting - ------------------- - """ - #compare_zonal_mean generates sets of six panel plots containing zonal mean data across your dataset. - #compare_zonal_mean shares many of the same arguments as compare_single_level. - #You can specify pressure ranges in hPa for zonal mean plotting (by default every vertical level is plotted) - gcplot.compare_zonal_mean(ref_ds, 'Dataset 1', dev_ds, 'Dataset 2', pres_range=[0, 100], varlist=varlist, pdfname='zonal_mean.pdf') - + from gcpy.constants import skip_these_vars + from gcpy.plot.compare_single_level import compare_single_level + from gcpy.plot.compare_zonal_mean import compare_zonal_mean + + + def main(): + """ + Example function to create six-panel comparison plots. + """ + + # xarray allows us to read in any NetCDF file, the format of + # GEOS-Chem diagnostics, #as an xarray Dataset + # + # The skip_these_vars list avoids trying to read certain + # GCHP variables that cause data read issues. + ref_ds = xr.open_dataset( + 'first_run/GEOSChem.Restart.20160801_0000z.nc4', + drop_variables=skip_these_vars + ) + dev_ds = xr.open_dataset( + 'second_run/GEOSChem.Restart.20160801_0000z.nc4', + drop_variables=skip_these_vars + ) + + # ================== + # Single level plots + # ================== + + # compare_single_level generates sets of six panel plots for + # data at a specified level in your datasets. By default, the + # level at index 0 (likely the surface) is plotted. Here we will + # plot data at ~500 hPa, which is located at index 21 in the + # standard 72-level and 47-level GMAO vertical grids. + ilev=21 + + # You likely want to look at the same variables across both of + # your datasets. If a variable is in one dataset but not the other, + # the plots will show NaN values for the latter. You can pass + # variable names in a list to these comparison plotting functions + # (otherwise all variables will plot). + varlist = ['SpeciesRst_O3', 'SpeciesRst_CO2'] + + # compare_single_level has many arguments which can be optionally + # specified. The first four arguments are required. They specify + # your first xarray Dataset, the name of your first dataset, + # your second xarray Dataset, and the name of your second dataset. + # Here we will also pass a specific level and the names of the + # variables you want to plot. + compare_single_level( + ref_ds, + 'Dataset 1', + dev_ds, + 'Dataset 2', + ilev=ilev, + varlist=varlist + ) + plt.show() + + # Using plt.show(), you can view the plots interactively. + # You can also save out the plots to a PDF. + compare_single_level( + ref_ds, + 'Dataset 1', + dev_ds, + 'Dataset 2', + ilev=ilev, + varlist=varlist, + pdfname='single_level.pdf' + ) + + # ================== + # Zonal Mean Plots + # ================== + + # compare_zonal_mean generates sets of six panel plots containing + # zonal mean data across your dataset. compare_zonal_mean shares + # many of the same arguments as compare_single_level. You can + # specify pressure ranges in hPa for zonal mean plotting (by + # default every vertical level is plotted) + compare_zonal_mean( + ref_ds, + 'Dataset 1', + dev_ds, + 'Dataset 2', + pres_range=[0, 100], + varlist=varlist, + pdfname='zonal_mean.pdf' + ) + + + # Only execute when we run as a standalone script + if __name__ == '__main__': + main() diff --git a/docs/source/Tabling.rst b/docs/source/Tabling.rst deleted file mode 100644 index d09d940f..00000000 --- a/docs/source/Tabling.rst +++ /dev/null @@ -1,436 +0,0 @@ -.. _table: - -####### -Tabling -####### - -This page describes the tabling capabilities of GCPy, including -possible argument values for every tabling function. These functions are primarily used for model benchmarking purposes. All tables are printed to text files. - - -================ -Emissions tables -================ - -.. code-block:: python - - - def make_benchmark_emis_tables(reflist, refstr, devlist, - devstr, dst="./benchmark", refmet=None, devmet=None, - overwrite=False, ref_interval=[2678400.0], dev_interval=[2678400.0], - spcdb_dir=os.path.dirname(__file__) - ): - """ - Creates a text file containing emission totals by species and - category for benchmarking purposes. - """ - -Arguments: ----------- - -.. option:: reflist: list of str - - List with the path names of the emissions file or files - (multiple months) that will constitute the "Ref" - (aka "Reference") data set. - -.. option:: refstr : str - - A string to describe ref (e.g. version number) - -.. option:: devlist : list of str - - List with the path names of the emissions file or files - (multiple months) that will constitute the "Dev" - (aka "Development") data set - -.. option:: devstr : str - - A string to describe dev (e.g. version number) - - -Keyword arguments: ------------------- - -.. option:: dst : str - - A string denoting the destination folder where the file - containing emissions totals will be written. - - Default value: ./benchmark - -.. option:: refmet : str - - Path name for ref meteorology - - Default value: None - -.. option:: devmet : str - - Path name for dev meteorology - - Default value: None - -.. option:: overwrite : bool - - Set this flag to True to overwrite files in the - destination folder (specified by the dst argument). - - Default value: False - -.. option:: ref_interval : list of float - - The length of the ref data interval in seconds. By default, interval - is set to [2678400.0], which is the number of seconds in July - (our 1-month benchmarking month). - - Default value: [2678400.0] - -.. option:: dev_interval : list of float - - The length of the dev data interval in seconds. By default, interval - is set to [2678400.0], which is the number of seconds in July - (our 1-month benchmarking month). - - Default value: [2678400.0] - -.. option:: spcdb_dir : str - - Directory of species_datbase.yml file - - Default value: Directory of GCPy code repository - - -:code:`gcpy.benchmark.make_benchmark_emis_tables()` generates tables -of total emissions categorized by species or by inventory. These -tables contain total global emissions over the lengths of the Ref and -Dev datasets, as well as the differences between totals across the two -datasets. Passing a list of datasets as Ref or Dev (e.g. multiple -months of emissions files) will result in printing totals emissions -summed across all files in the list. Make sure to update the -literal:`ref_interval` and/or :literal:`dev_interval` arguments if you -pass input that does not correspond with 1 31 day month. - -=========== -Mass Tables -=========== - -.. code-block:: python - - def make_benchmark_mass_tables(ref, refstr, dev, devstr, - varlist=None, dst="./benchmark", subdst=None, overwrite=False, - verbose=False, label="at end of simulation", - spcdb_dir=os.path.dirname(__file__), - ref_met_extra='', dev_met_extra='' - ): - """ - Creates a text file containing global mass totals by species and - category for benchmarking purposes. - """ - -Arguments: ----------- - -.. option:: reflist : str - - Pathname that will constitute - the "Ref" (aka "Reference") data set. - -.. option:: refstr : str - - A string to describe ref (e.g. version number) - -.. option:: dev : list of str - - Pathname that will constitute - the "Dev" (aka "Development") data set. The "Dev" - data set will be compared against the "Ref" data set. - -.. option:: devstr : str - - A string to describe dev (e.g. version number) - - -Keyword arguments: ------------------- - -.. option:: varlist : list of str - - List of variables to include in the list of totals. - If omitted, then all variables that are found in either - "Ref" or "Dev" will be included. The varlist argument - can be a useful way of reducing the number of - variables during debugging and testing. - - Default value: None - -.. option:: dst : str - - A string denoting the destination folder where the file - containing emissions totals will be written. - - Default value: ./benchmark - -.. option:: subdst : str - - A string denoting the sub-directory of dst where PDF - files containing plots will be written. In practice, - subdst is only needed for the 1-year benchmark output, - and denotes a date string (such as "Jan2016") that - corresponds to the month that is being plotted. - - Default value: None - -.. option:: overwrite : bool - - Set this flag to True to overwrite files in the - destination folder (specified by the dst argument). - - Default value: False - -.. option:: verbose : bool - - Set this flag to True to print extra informational output. - - Default value: False. - -.. option:: spcdb_dir : str - - Directory of species_datbase.yml file - - Default value: Directory of GCPy code repository - -.. option:: ref_met_extra : str - - Path to ref Met file containing area data for use with restart files - which do not contain the Area variable. - Default value : '' - -.. option:: dev_met_extra : str - - Path to dev Met file containing area data for use with restart files - which do not contain the Area variable. - - Default value: '' - - -:code:`gcpy.benchmark.make_benchmark_mass_tables` is used to create -global mass tables of GEOS-Chem species from a -:literal:`Restart` file. This function will create one table of total -mass by species from the earth's surface to the top of the -stratosphere and one table for only the troposphere. -The tables contain total mass for each of the ref and dev datasets in -Gg, as well as absolute and percentage difference between the two -datasets. If your restart files do not contain an Area variable -(:literal:`AREA` for GEOS-Chem Classic or :literal:`Met_AREAM2` for -GCHP) then you will need to use the :literal:`ref_met_extra` and/or -:literal:`dev_met_extra` arguments to pass the paths of NetCDF files -containing the corresponding area variables (usually contained in -meteorology diagnostic output). - -======================== -Operations Budget Tables -======================== - -.. code-block:: python - - def make_benchmark_operations_budget(refstr, reffiles, devstr, - devfiles, ref_interval, dev_interval, benchmark_type=None, - label=None, col_sections=["Full", "Trop", "PBL", "Strat"], - operations=["Chemistry","Convection","EmisDryDep","Mixing", - "Transport","WetDep"], compute_accum=True, - require_overlap=False, dst='.', species=None, overwrite=True - ): - """ - Prints the "operations budget" (i.e. change in mass after - each operation) from a GEOS-Chem benchmark simulation. - """ - - -Arguments: ----------- - -.. option:: refstr : str - - Labels denoting the "Ref" versions - -.. option:: reffiles : list of str - - Lists of files to read from the "Ref" version. - -.. option:: devstr : str - - Labels denoting the "Dev" versions - -.. option:: devfiles : list of str - - Lists of files to read from "Dev" version. - -.. option:: interval : float - - Number of seconds in the diagnostic interval. - - -Keyword arguments: ------------------- - -.. option:: benchmark_type : str - - "TransportTracersBenchmark" or "FullChemBenchmark". - - Default value: None - -.. option:: label : str - - Contains the date or date range for each dataframe title. - - Default value: None - -.. option:: col_sections : list of str - - List of column sections to calculate global budgets for. May - include Strat eventhough not calculated in GEOS-Chem, but Full - and Trop must also be present to calculate Strat. - - Default value: ["Full", "Trop", "PBL", "Strat"] - -.. option:: operations : list of str - - List of operations to calculate global budgets for. Accumulation - should not be included. It will automatically be calculated if - all GEOS-Chem budget operations are passed and optional arg - compute_accum is True. - - Default value: ["Chemistry","Convection","EmisDryDep", - "Mixing","Transport","WetDep"] - -.. option:: compute_accum : bool - - Optionally turn on/off accumulation calculation. If True, will - only compute accumulation if all six GEOS-Chem operations budgets - are computed. Otherwise a message will be printed warning that - accumulation will not be calculated. - - Default value: True - -.. option:: require_overlap : bool - - Whether to calculate budgets for only species that are present in - both Ref or Dev. - - Default value: False - -.. option:: dst : str - - Directory where plots & tables will be created. - - Default value: '.' (directory in which function is called) - -.. option:: species : list of str - - List of species for which budgets will be created. - - Default value: None (all species) - -.. option:: overwrite : bool - - Denotes whether to overwrite existing budget file. - - Default value: True - -:code:`gcpy.benchmark.make_benchmark_operations_budget()` creates -tables of budgets for species separated by model operation. The tables -show budgets for each of the ref and dev datasets in Gg, as well as -absolute and percentage difference between the two datasets. -Note that total accumulation across all operations will only be -printed if you set :code:`compute_accum==True` and -all operations are included in :literal:`operations`. Note also that -when using the non-local mixing scheme (default), :literal:`'Mixing'` -includes emissions and dry deposition applied below the -PBL. :literal:`'EmisDryDep'` therefore only captures fluxes above the -PBL. When using full mixing, :literal:`'Mixing'` and -:literal:`'EmisDryDep'` are fully separated. - -=========================== -Aerosol Budgets and Burdens -=========================== - -.. code-block:: python - - def make_benchmark_aerosol_tables(devdir, devlist_aero, devlist_spc, - devlist_met, devstr, year, days_per_mon, dst='./benchmark', - overwrite=False, is_gchp=False, spcdb_dir=os.path.dirname(__file__) - ): - """ - Compute FullChemBenchmark aerosol budgets & burdens - """ - - -Arguments: ----------- - -.. option:: devdir: str - - Path to development ("Dev") data directory - -.. option:: devlist_aero : list of str - - List of Aerosols collection files (different months) - -.. option:: devlist_spc : list of str - - List of SpeciesConc collection files (different months) - -.. option:: devlist_met : list of str - - List of meteorology collection files (different months) - -.. option:: devstr : str - - Descriptive string for datasets (e.g. version number) - -.. option:: year : str - - The year of the benchmark simulation (e.g. '2016'). - -.. option:: days_per_mon : list of int - - List of number of days per month for all months - - -Keyword arguments: ------------------- - -.. option:: dst : str - - Directory where budget tables will be created. - - Default value: './benchmark' - -.. option:: overwrite : bool - - Overwrite burden & budget tables? (default=True) - - Default value: False - -.. option:: is_gchp : bool - - Whether datasets are for GCHP - - Default value: False - -.. option:: spcdb_dir : str - - Directory of species_datbase.yml file - - Default value: Directory of GCPy code repository - - -:code:`gcpy.benchmark.make_benchmark_aerosol_tables()` generates two -different tables using output from a single dataset. One -contains annual mean aerosol burdens in Tg in the stratosphere, -troposphere, and combined stratosphere and troposphere. The other -table shows annual global mean AOD in the stratosphere, troposphere, -and combined stratosphere and troposphere. Aerosol species used are -pre-defined in :code:`aod_species.yml`: BCPI, OCPI, SO4, DST1, SALA, -and SALC. diff --git a/docs/source/benchmark_plotting.rst b/docs/source/benchmark_plotting.rst deleted file mode 100644 index e2fd80e4..00000000 --- a/docs/source/benchmark_plotting.rst +++ /dev/null @@ -1,146 +0,0 @@ -############################ -Benchmark Plotting / Tabling -############################ - -Below is an example configuration file used to input the desired -options for the comprehensive benchmark comparison script -:code:`run_benchmark.py`. Additional configuration file examples can -be found in the :file:`benchmarks` directory of GCpy. - -The :file:`run_benchmark.py` script allows one to perform benchmark -comparisons between any simulation duration supplied in the -configuration file provided the ref and dev simulations time periods -match. Additionally, if the durations specified are exactly one year, -then the corresponding :literal:`bmk_type` specialty comparison script -will be run (either :file:`run_1yr_fullchem_benchmark.py` or -:file:`run_1yr_tt_benchmark.py`). Any other duration will run the standard -suite of benchmark comparisons. - -To generate plots from a 1-month benchmark simulation, you would call -:file:`run_benchmark.py` as follows: - -.. code-block:: console - - (gcpy_env) $ run_benchmark.py 1mo_benchmark.yml - -Where :file:`1mo_benchmark.yml` contains the following inputs: - -.. code-block:: yaml - - &--- - # ===================================================================== - # Benchmark configuration file (**EDIT AS NEEDED**) - # customize in the following manner: - # (1) Edit the path variables so that they point to folders w/ model data - # (2) Edit the version strings for each benchmark simulation - # (3) Edit the switches that turn on/off creating of plots and tables - # (4) If necessary, edit labels for the dev and ref versions - # Note: When doing GCHP vs GCC comparisions gchp_dev will be compared - # to gcc_dev (not gcc_ref!). This ensures consistency in version names - # when doing GCHP vs GCC diff-of-diffs (mps, 6/27/19) - # ===================================================================== - # - # Configuration for 1 month FullChemBenchmark - # - # paths: - # main_dir: High-level directory containing ref & dev rundirs - # results_dir: Directory where plots/tables will be created - # weights_dir: Path to regridding weights - # spcdb_dir: Folder in which the species_database.yml file is - # located. If set to "default", then will look for - # species_database.yml in one of the Dev rundirs. - # - paths: - main_dir: /n/holyscratch01/external_repos/GEOS-CHEM/gcgrid/geos-chem/validation/gcpy_test_data/1mon - results_dir: /path/to/BenchmarkResults - weights_dir: /n/holyscratch01/external_repos/GEOS-CHEM/gcgrid/gcdata/ExtData/GCHP/RegriddingWeights - spcdb_dir: default - # - # data: Contains configurations for ref and dev runs - # version: Version string (must not contain spaces) - # dir: Path to run directory - # outputs_subdir: Subdirectory w/ GEOS-Chem diagnostic files - # restarts_subdir: Subdirectory w/ GEOS-Chem restarts - # bmk_start: Simulation start date (YYYY-MM-DDThh:mm:ss) - # bmk_end: Simulation end date (YYYY-MM-DDThh:mm:ss) - # resolution: GCHP resolution string - # - data: - ref: - gcc: - version: GCC_ref - dir: GCC_ref - outputs_subdir: OutputDir - restarts_subdir: Restarts - bmk_start: "2019-07-01T00:00:00" - bmk_end: "2019-08-01T00:00:00" - gchp: - version: GCHP_ref - dir: GCHP_ref - outputs_subdir: OutputDir - restarts_subdir: Restarts - bmk_start: "2019-07-01T00:00:00" - bmk_end: "2019-08-01T00:00:00" - is_pre_13.1: False - is_pre_14.0: False - resolution: c24 - dev: - gcc: - version: GCC_dev - dir: GCC_dev - outputs_subdir: OutputDir - restarts_subdir: Restarts - bmk_start: "2019-07-01T00:00:00" - bmk_end: "2019-08-01T00:00:00" - gchp: - version: GCHP_dev - dir: GCHP_dev - outputs_subdir: OutputDir - restarts_subdir: Restarts - bmk_start: "2019-07-01T00:00:00" - bmk_end: "2019-08-01T00:00:00" - is_pre_13.1: False - is_pre_14.0: False - resolution: c24 - # - # options: Specify the types of comparisons to perform - # - options: - bmk_type: FullChemBenchmark - gcpy_test: True # Specify if this is a gcpy test validation run - comparisons: - gcc_vs_gcc: - run: True # True to run this comparison - dir: GCC_version_comparison - tables_subdir: Tables - gchp_vs_gcc: - run: True - dir: GCHP_GCC_comparison - tables_subdir: Tables - gchp_vs_gchp: - run: True - dir: GCHP_version_comparison - tables_subdir: Tables - gchp_vs_gcc_diff_of_diffs: - run: True - dir: GCHP_GCC_diff_of_diffs - # - # outputs: Types of output to generate (plots/tables) - # - outputs: - plot_conc: True - plot_emis: True - emis_table: True - plot_jvalues: True - plot_aod: True - mass_table: True - ops_budget_table: False - OH_metrics: True - ste_table: True # GCC only - plot_options: # Plot concentrations and emissions by category? - by_spc_cat: True - by_hco_cat: True - -YAML configuration files for 1-year benchmarks -(:file:`1yr_fullchem_benchmark.yml`, :file:`1yr_tt_benchmark.yml`) are -also provided in the :file:`benchmarks` folder. diff --git a/docs/source/bpch_to_nc.rst b/docs/source/bpch_to_nc.rst deleted file mode 100644 index d7dd9325..00000000 --- a/docs/source/bpch_to_nc.rst +++ /dev/null @@ -1,188 +0,0 @@ -Convert BPCH to NetCDF -====================== - -.. code-block:: python - - - #!/usr/bin/env python - ''' - Example script that illustrates how to create a netCDF file - from an old GEOS-Chem binary punch ("bpch") file. - ''' - - # Imports - import gcpy - import xarray as xr - import xbpch as xb - import warnings - - # Suppress harmless run-time warnings (mostly about underflow in division) - warnings.filterwarnings('ignore', category=RuntimeWarning) - warnings.filterwarnings('ignore', category=UserWarning) - - # ---------------------------------------------------------------------- - # User configurable settings (EDIT THESE ACCORDINGLY) - # ---------------------------------------------------------------------- - - # Name of Bpch file - bpchfile = '/path/to/bpch/file' - - # tracerinfo.dat and diaginfo,dat fiels - tinfo_file = '/path/to/tracerinfo.dat' - dinfo_file = '/path/to/diaginfo.dat' - - # Name of netCDF file - ncfile = '/path/to/netcdf/file' - - # Date string for the time:units attribute - datestr = 'YYYY-MM-DD' - - # Number of seconds in the diagnostic interval (assume 1-month) - interval = 86400.0 * 31.0 - - # ---------------------------------------------------------------------- - # Open the bpch file and save it into an xarray Dataset object - # NOTE: For best results, also specify the corresponding - # tracerinfo.dat diaginfo.dat metadata files. - # ---------------------------------------------------------------------- - try: - ds = xb.open_bpchdataset(filename=bpchfile, - tracerinfo_file=tinfo_file, - diaginfo_file=dinfo_file) - except FileNotFoundError: - print('Could not find file {}'.format(bpchfile)) - raise - - # ---------------------------------------------------------------------- - # Further manipulate the Dataset - # ---------------------------------------------------------------------- - - # Transpose the order of the xarray Dataset object read by - # xbpch so that its dimensions will be in the same order as - # Dataset objects read from netCDF files. - ds = ds.transpose() - - # Convert the bpch variable names to the same naming - # convention as the netCDF ("History") diagnostics. - ds = gcpy.convert_bpch_names_to_netcdf_names(ds) - - # xbpch does not include a time dimension, so we'll add one here - coords = ds.coords - coords['time'] = 0.0 - - # ------------------------------------------------------------------ - # Further edit variable attributes - # ------------------------------------------------------------------ - for v in ds.data_vars.keys(): - - # Append time to the data array - ds[v] = xr.concat([ds[v]], 'time') - - # Add long_name attribute for COARDS netCDF compliance - ds[v].attrs['long_name'] = ds[v].attrs['full_name'] - - # Remove some extraneous attributes that xbpch sets - del ds[v].attrs['name'] - del ds[v].attrs['full_name'] - del ds[v].attrs['scale_factor'] - del ds[v].attrs['hydrocarbon'] - del ds[v].attrs['tracer'] - del ds[v].attrs['category'] - del ds[v].attrs['chemical'] - del ds[v].attrs['original_shape'] - del ds[v].attrs['origin'] - del ds[v].attrs['number'] - del ds[v].attrs['molwt'] - del ds[v].attrs['C'] - - # Make the units attribute consistent with the units - # attribute from the GEOS-Chem History diagnostics - # NOTE: There probably is a more Pythonic way to code - # this, but this will work for sure. - if 'ug/m3' in ds[v].units: - ds[v].attrs['units'] = 'ug m-3' - if 'ug Celsius/m3' in ds[v].units: - ds[v].attrs['units'] = 'ug C m-3' - if 'count/cm3' in ds[v].units: - ds[v].attrs['units'] = 'molec m-3' - if 'cm/s' in ds[v].units: - ds[v].attrs['units'] = 'cm s-1' - if 'count/cm2/s' in ds[v].units: - ds[v].attrs['units'] = 'molec cm-2 s-1' - if 'kg/m2s' in ds[v].units: - ds[v].attrs['units'] = 'kg m-2 s-1' - if 'kg/m2/s' in ds[v].units: - ds[v].attrs['units'] = 'kg m-2 s-1' - if 'kg/s' in ds[v].units: - ds[v].attrs['units'] = 'kg s-1' - if 'W/m2' in ds[v].units: - ds[v].attrs['units'] = 'W m-2' - if 'm/s' in ds[v].units: - ds[v].attrs['units'] = 'm s-1' - if 'Pa/s' in ds[v].units: - ds[v].attrs['units'] = 'Pa s-1' - if 'g/kg' in ds[v].units: - ds[v].attrs['units'] = 'g kg-1' - if v.strip() == 'TotalOC': - ds[v].attrs['units'] = 'ug m-3' - if v.strip() in [ 'HO2concAfterChem']: - ds[v].attrs['units'] = 'ppb' - if v.strip() in ['O1DconcAfterChem', - 'O3PconcAfterChem', - 'OHconcAfterChem']: - ds[v].attrs['units'] = 'molec cm-3' - if v.strip() in ['Loss_CO', 'Prod_CO', - 'Loss_Ox', 'Prod_Ox', 'Prod_SO4']: - ds[v].attrs['units'] = 'molec/cm3/s' - if v.strip() in 'Met_CLDTOPS': - ds[v].attrs['units'] = 'level' - if v.strip() in 'Met_PHIS': - ds[v].attrs['units'] = 'm2 s-1' - if v.strip() in ['Met_PRECCON', 'Met_PRECTOT']: - ds[v].attrs['units'] = 'kg m-2 s-1' - if v.strip() in 'Met_AVGW': - ds[v].attrs['units'] = 'vol vol-1' - if v.strip() in 'Met_AIRNUMDEN': - ds[v].attrs['units'] = 'molec cm-3' - if v.strip() in ['ProdCOfromCH4', 'ProdCOfromNMVOC']: - ds[v].attrs['units'] = 'molec cm-3 s-1' - - # Convert these prodloss diagnostics from kg (bpch) to kg/s - # to be consistent with the GEOS-Chem History diagnostics - # NOTE: Assume a 1-month interval ( - if v.strip() in ['ProdSO4fromH2O2inCloud', 'ProdSO4fromO3inCloud', - 'ProdSO4fromO2inCloudMetal', 'ProdSO4fromO3inSeaSalt', - 'ProdSO4fromHOBrInCloud', 'ProdSO4fromSRO3', - 'ProdSO4fromSRHObr', 'ProdSO4fromO3s']: - ds[v].attrs['units'] = 'kg S s-1' - ds[v] = ds[v] / interval - if v.strip() in ['LossHNO3onSeaSalt']: - ds[v].attrs['units'] = 'kg s-1' - ds[v] = ds[v] / interval - - # ------------------------------------------------------------------ - # Edit attributes for coordinate dimensions - # ------------------------------------------------------------------ - - # Time - ds['time'].attrs['long_name'] = 'time' - ds['time'].attrs['units'] = \ - 'hours since {} 00:00:00.00 UTC'.format(datestr) - ds['time'].attrs['calendar'] = 'standard' - ds['time'].attrs['axis'] = 'T' - - # "lon", "lat", "lev" - ds['lon'].attrs['axis'] = 'X' - ds['lat'].attrs['axis'] = 'Y' - ds['lev'].attrs['axis'] = 'Z' - ds['lev'].attrs['units'] = 'level' - - # Global title - ds.attrs['title'] = 'Created by bpch2nc.py' - ds.attrs['conventions'] = 'COARDS' - ds.attrs['references'] = 'www.geos-chem.org; wiki.geos-chem.org' - - # ------------------------------------------------------------------ - # Create the netCDF file - # ------------------------------------------------------------------ - ds.to_netcdf(ncfile) diff --git a/docs/source/index.rst b/docs/source/index.rst index b280e341..1c38079b 100644 --- a/docs/source/index.rst +++ b/docs/source/index.rst @@ -25,25 +25,23 @@ For documentation on setting up and running GEOS-Chem please see our .. toctree:: :maxdepth: 4 :caption: Usage Details: - + Plotting - Tabling Regridding + Benchmarking .. toctree:: :maxdepth: 1 :caption: Example Scripts: - + Six_panel Single_panel - benchmark_plotting plot_timeseries - bpch_to_nc .. toctree:: :maxdepth: 1 :caption: Help & Reference: - + reference/CONTRIBUTING.md reference/SUPPORT.md editing_these_docs diff --git a/docs/source/plot_timeseries.rst b/docs/source/plot_timeseries.rst index be60ce53..f228e533 100644 --- a/docs/source/plot_timeseries.rst +++ b/docs/source/plot_timeseries.rst @@ -1,8 +1,12 @@ +.. _plot-timeseries: + +############### Plot Timeseries -=============== +############### -.. code-block:: python +This example script may also be found at `gcpy/examples/plotting/plot_single_panel.py `_. +.. code-block:: python #!/usr/bin/env python ''' @@ -14,19 +18,19 @@ Plot Timeseries Page 1: ------- - O3 from the first model layer (from the "SpeciesConc" - diagnostic collection is) plotted in blue. + O3 from the first model layer (from the "SpeciesConc" + diagnostic collection is) plotted in blue. - O3 at 10 meter height (from the "SpeciesConc_10m" - diagnostic collection) is plotted in red. + O3 at 10 meter height (from the "SpeciesConc_10m" + diagnostic collection) is plotted in red. Page 2: ------- - HNO3 from the first model layer (from the SpeciesConc - diagnostic collection is) plotted in blue. + HNO3 from the first model layer (from the SpeciesConc + diagnostic collection is) plotted in blue. - HNO3 at 10 meter height (from the SpeciesConc_10m - diagnostic collection) is plotted in red. + HNO3 at 10 meter height (from the SpeciesConc_10m + diagnostic collection) is plotted in red. You can of course modify this for your own particular applications. @@ -38,15 +42,16 @@ Plot Timeseries ''' # Imports - import gcpy.constants as gcon import os + import warnings import numpy as np import matplotlib.dates as mdates import matplotlib.ticker as mticker import matplotlib.pyplot as plt from matplotlib.backends.backend_pdf import PdfPages import xarray as xr - import warnings + from gcpy import constants + # Tell matplotlib not to look for an X-window, as we are plotting to # a file and not to the screen. This will avoid some warning messages. @@ -58,292 +63,292 @@ Plot Timeseries def find_files_in_dir(path, substrs): - ''' - Returns a list of all files in a directory that match one or more - substrings. - - Args: - ----- - path : str - Path to the directory in which to search for files. - - substrs : list of str - List of substrings used in the search for files. - - Returns: - -------- - file_list : list of str - List of files in the directory (specified by path) - that match all substrings (specified in substrs). - ''' - - # Initialize - file_list = [] - - # Walk through the given data directory. Then for each file found, - # add it to file_list if it matches text in search_list. - for root, directory, files in os.walk(path): - for f in files: - for s in substrs: - if s in f: - file_list.append(os.path.join(root, f)) - - # Return an alphabetically sorted list of files - file_list.sort() - return file_list - - - def find_value_index(seq, val): - ''' - Finds the index of a numpy array that is close to a value. - - Args: - ----- - seq : numpy ndarray - An array of numeric values. - - val : number - The value to search for in seq. - - Returns: - -------- - result : integer - The index of seq that has a value closest to val. - - Remarks: - -------- - This algorithm was found on this page: - https://stackoverflow.com/questions/48900977/find-all-indexes-of-a-numpy-array-closest-to-a-value - ''' - r = np.where(np.diff(np.sign(seq - val)) != 0) - idx = r + (val - seq[r]) / (seq[r + np.ones_like(r)] - seq[r]) - idx = np.append(idx, np.where(seq == val)) - idx = np.sort(idx) - result = np.round(idx) - - # NOTE: xarray needs integer values, so convert here! - return int(result[0]) - - - def read_geoschem_data(path, collections): - ''' - Returns an xarray Dataset containing timeseries data. - - Args: - ----- - path : str - Directory path where GEOS-Chem diagnostic output - files may be found. - - collections: list of str - List of GEOS-Chem collections. Files for these - collections will be read into the xarray Dataset. - - Returns: - -------- - ds : xarray Dataset - A Dataset object containing the GEOS-Chem diagnostic - output corresponding to the collections that were - specified. - ''' - - # Get a list of variables that GCPy should not read. - # These are mostly variables introduced into GCHP with the MAPL v1.0.0 - # update. These variables contain either repeated or non-standard - # dimensions that can cause problems in xarray when combining datasets. - skip_vars = gcon.skip_these_vars - - # Find all files in the given - file_list = find_files_in_dir(path, collections) - - # Return a single xarray Dataset containing data from all files - # NOTE: Need to add combine="nested" for xarray 0.15 and higher - v = xr.__version__.split(".") - if int(v[0]) == 0 and int(v[1]) >= 15: - return xr.open_mfdataset(file_list, - drop_variables=skip_vars, - combine="nested", - concat_dim=None) - else: - return xr.open_mfdataset(file_list, - drop_variables=skip_vars) - - - def plot_timeseries_data(ds, site_coords): - ''' - Plots a timseries of data at a given (lat,lon) location. - - Args: - ----- - ds : xarray Dataset - Dataset containing GEOS-Chem timeseries data. - - site_coords : tuple - Contains the coordinate (lat, lon) of a site location - at which the timeseries data will be plotted. - ''' - - # ---------------------------------------------------------------------- - # Get the GEOS-Chem data for O3 and HNO3 corresponding to the - # location of the observational station. We will save these into - # xarray DataArray objects, which we'll need for plotting. - # - # YOU CAN EDIT THIS FOR YOUR OWN PARTICULAR APPLICATION! - # ---------------------------------------------------------------------- - - # Find the indices corresponding to the site lon and lat - lat_idx = find_value_index(ds.lat.values, site_coords[0]) - lon_idx = find_value_index(ds.lon.values, site_coords[1]) - - # Save O3 from the first level (~60m height) (ppb) into a DataArray - O3_L1 = ds['SpeciesConc_O3'].isel(lon=lon_idx, lat=lat_idx, lev=0) - O3_L1 *= 1.0e9 - O3_L1.attrs['units'] = 'ppbv' - - # Save O3 @ 10m height into a DataArray - O3_10m = ds['SpeciesConc10m_O3'].isel(lon=lon_idx, lat=lat_idx) - O3_10m *= 1.0e9 - O3_10m.attrs['units'] = 'ppbv' - - # Save HNO3 from the first level (~60m height) into a DataArray - HNO3_L1 = ds['SpeciesConc_HNO3'].isel(lon=lon_idx, lat=lat_idx, lev=0) - HNO3_L1 *= 1.0e9 - HNO3_L1.attrs['units'] = 'ppbv' - - # Save HNO3 @ 10m height into a DataArray - HNO3_10m = ds['SpeciesConc10m_HNO3'].isel(lon=lon_idx, lat=lat_idx) - HNO3_10m *= 1.0e9 - HNO3_10m.attrs['units'] = 'ppbv' - - # ---------------------------------------------------------------------- - # Create a PDF file of the plots - # ---------------------------------------------------------------------- - - # Get min & max days of the plot span (for setting the X-axis range). - # To better center the plot, add a cushion of 12 hours on either end. - time = ds['time'].values - datemin = np.datetime64(time[0]) - np.timedelta64(12, 'h') - datemax = np.datetime64(time[-1]) + np.timedelta64(12, 'h') - - # Define a PDF object so that we can save the plots to PDF - pdf = PdfPages('O3_and_HNO3.pdf') - - # Loop over number of desired pages (in this case, 2) - for i in range(0, 2): - - # Create a new figure: 1 plot per page, 2x as wide as high - figs, ax0 = plt.subplots(1, 1, figsize=[12, 6]) - - # ----------------------------- - # Plot O3 on the first page - # ----------------------------- - if i == 0: - - # 1st model level - O3_L1.plot.line(ax=ax0, x='time', color='blue', - marker='o', label='O3 from 1st model level', - linestyle='-') - - # 10 mheight - O3_10m.plot.line(ax=ax0, x='time', color='red', - marker='x', label='O3 at 10m height', - linestyle='-') - - # Set title (has to be after the line plots are drawn) - ax0.set_title('O3 from the 1st model level and at 10m height') - - # Set Y-axis minor tick marks at every 2 ppb (5 intervals) - ax0.yaxis.set_minor_locator(mticker.AutoMinorLocator(5)) - - # Set y-axis title - ax0.set_ylabel('O3 (ppbv)') - - # ----------------------------- - # Plot HNO3 on the second page - # ----------------------------- - if i == 1: - - # 1st model level - HNO3_L1.plot.line(ax=ax0, x='time', color='blue', - marker='o', label='HNO3 from 1st model level', - linestyle='-') - - # 10m height - HNO3_10m.plot.line(ax=ax0, x='time', color='red', - marker='x', label='HNO3 at 10m height', - linestyle='-') - - # Set title (has to be after the line plots are drawn - ax0.set_title('HNO3 from the 1st model level and at 10m height') - - # Set Y-axis minor tick marks at every 0.05 ppb (4 intervals) - ax0.yaxis.set_minor_locator(mticker.AutoMinorLocator(4)) - - # Set y-axis title - ax0.set_ylabel('HNO3 (ppbv)') - - # ----------------------------- - # Set general plot parameters - # ----------------------------- - - # Add the plot legend - ax0.legend() - - # Set the X-axis range - ax0.set_xlim(datemin, datemax) - - # Set the X-axis major tickmarks - locator = mdates.DayLocator() - formatter = mdates.DateFormatter('%d') - ax0.xaxis.set_major_locator(locator) - ax0.xaxis.set_major_formatter(formatter) - - # Set X-axis minor tick marks at noon of each day - # (i.e. split up the major interval into 2 bins) - ax0.xaxis.set_minor_locator(mticker.AutoMinorLocator(2)) - - # Don't rotate the X-axis jtick labels - ax0.xaxis.set_tick_params(rotation=0) + ''' + Returns a list of all files in a directory that match one or more + substrings. + + Args: + ----- + path : str + Path to the directory in which to search for files. + + substrs : list of str + List of substrings used in the search for files. + + Returns: + -------- + file_list : list of str + List of files in the directory (specified by path) + that match all substrings (specified in substrs). + ''' + + # Initialize + file_list = [] + + # Walk through the given data directory. Then for each file found, + # add it to file_list if it matches text in search_list. + for root, directory, files in os.walk(path): + for f in files: + for s in substrs: + if s in f: + file_list.append(os.path.join(root, f)) + + # Return an alphabetically sorted list of files + file_list.sort() + return file_list + + + def find_value_index(seq, val): + ''' + Finds the index of a numpy array that is close to a value. + + Args: + ----- + seq : numpy ndarray + An array of numeric values. + + val : number + The value to search for in seq. + + Returns: + -------- + result : integer + The index of seq that has a value closest to val. + + Remarks: + -------- + This algorithm was found on this page: + https://stackoverflow.com/questions/48900977/find-all-indexes-of-a-numpy-array-closest-to-a-value + ''' + r = np.where(np.diff(np.sign(seq - val)) != 0) + idx = r + (val - seq[r]) / (seq[r + np.ones_like(r)] - seq[r]) + idx = np.append(idx, np.where(seq == val)) + idx = np.sort(idx) + result = np.round(idx) + + # NOTE: xarray needs integer values, so convert here! + return int(result[0]) + + + def read_geoschem_data(path, collections): + ''' + Returns an xarray Dataset containing timeseries data. + + Args: + ----- + path : str + Directory path where GEOS-Chem diagnostic output + files may be found. + + collections: list of str + List of GEOS-Chem collections. Files for these + collections will be read into the xarray Dataset. + + Returns: + -------- + ds : xarray Dataset + A Dataset object containing the GEOS-Chem diagnostic + output corresponding to the collections that were + specified. + ''' + + # Get a list of variables that GCPy should not read. + # These are mostly variables introduced into GCHP with the MAPL v1.0.0 + # update. These variables contain either repeated or non-standard + # dimensions that can cause problems in xarray when combining datasets. + skip_vars = constants.skip_these_vars + + # Find all files in the given + file_list = find_files_in_dir(path, collections) + + # Return a single xarray Dataset containing data from all files + # NOTE: Need to add combine="nested" for xarray 0.15 and higher + v = xr.__version__.split(".") + if int(v[0]) == 0 and int(v[1]) >= 15: + return xr.open_mfdataset(file_list, + drop_variables=skip_vars, + combine="nested", + concat_dim=None) + else: + return xr.open_mfdataset(file_list, + drop_variables=skip_vars) + + + def plot_timeseries_data(ds, site_coords): + ''' + Plots a timseries of data at a given (lat,lon) location. + + Args: + ----- + ds : xarray Dataset + Dataset containing GEOS-Chem timeseries data. + + site_coords : tuple + Contains the coordinate (lat, lon) of a site location + at which the timeseries data will be plotted. + ''' + + # ---------------------------------------------------------------------- + # Get the GEOS-Chem data for O3 and HNO3 corresponding to the + # location of the observational station. We will save these into + # xarray DataArray objects, which we'll need for plotting. + # + # YOU CAN EDIT THIS FOR YOUR OWN PARTICULAR APPLICATION! + # ---------------------------------------------------------------------- + + # Find the indices corresponding to the site lon and lat + lat_idx = find_value_index(ds.lat.values, site_coords[0]) + lon_idx = find_value_index(ds.lon.values, site_coords[1]) + + # Save O3 from the first level (~60m height) (ppb) into a DataArray + O3_L1 = ds['SpeciesConc_O3'].isel(lon=lon_idx, lat=lat_idx, lev=0) + O3_L1 *= 1.0e9 + O3_L1.attrs['units'] = 'ppbv' + + # Save O3 @ 10m height into a DataArray + O3_10m = ds['SpeciesConc10m_O3'].isel(lon=lon_idx, lat=lat_idx) + O3_10m *= 1.0e9 + O3_10m.attrs['units'] = 'ppbv' + + # Save HNO3 from the first level (~60m height) into a DataArray + HNO3_L1 = ds['SpeciesConc_HNO3'].isel(lon=lon_idx, lat=lat_idx, lev=0) + HNO3_L1 *= 1.0e9 + HNO3_L1.attrs['units'] = 'ppbv' + + # Save HNO3 @ 10m height into a DataArray + HNO3_10m = ds['SpeciesConc10m_HNO3'].isel(lon=lon_idx, lat=lat_idx) + HNO3_10m *= 1.0e9 + HNO3_10m.attrs['units'] = 'ppbv' + + # ---------------------------------------------------------------------- + # Create a PDF file of the plots + # ---------------------------------------------------------------------- + + # Get min & max days of the plot span (for setting the X-axis range). + # To better center the plot, add a cushion of 12 hours on either end. + time = ds['time'].values + datemin = np.datetime64(time[0]) - np.timedelta64(12, 'h') + datemax = np.datetime64(time[-1]) + np.timedelta64(12, 'h') + + # Define a PDF object so that we can save the plots to PDF + pdf = PdfPages('O3_and_HNO3.pdf') + + # Loop over number of desired pages (in this case, 2) + for i in range(0, 2): + + # Create a new figure: 1 plot per page, 2x as wide as high + figs, ax0 = plt.subplots(1, 1, figsize=[12, 6]) + + # ----------------------------- + # Plot O3 on the first page + # ----------------------------- + if i == 0: + + # 1st model level + O3_L1.plot.line(ax=ax0, x='time', color='blue', + marker='o', label='O3 from 1st model level', + linestyle='-') + + # 10 mheight + O3_10m.plot.line(ax=ax0, x='time', color='red', + marker='x', label='O3 at 10m height', + linestyle='-') + + # Set title (has to be after the line plots are drawn) + ax0.set_title('O3 from the 1st model level and at 10m height') + + # Set Y-axis minor tick marks at every 2 ppb (5 intervals) + ax0.yaxis.set_minor_locator(mticker.AutoMinorLocator(5)) + + # Set y-axis title + ax0.set_ylabel('O3 (ppbv)') + + # ----------------------------- + # Plot HNO3 on the second page + # ----------------------------- + if i == 1: + + # 1st model level + HNO3_L1.plot.line(ax=ax0, x='time', color='blue', + marker='o', label='HNO3 from 1st model level', + linestyle='-') + + # 10m height + HNO3_10m.plot.line(ax=ax0, x='time', color='red', + marker='x', label='HNO3 at 10m height', + linestyle='-') + + # Set title (has to be after the line plots are drawn + ax0.set_title('HNO3 from the 1st model level and at 10m height') + + # Set Y-axis minor tick marks at every 0.05 ppb (4 intervals) + ax0.yaxis.set_minor_locator(mticker.AutoMinorLocator(4)) + + # Set y-axis title + ax0.set_ylabel('HNO3 (ppbv)') + + # ----------------------------- + # Set general plot parameters + # ----------------------------- + + # Add the plot legend + ax0.legend() + + # Set the X-axis range + ax0.set_xlim(datemin, datemax) + + # Set the X-axis major tickmarks + locator = mdates.DayLocator() + formatter = mdates.DateFormatter('%d') + ax0.xaxis.set_major_locator(locator) + ax0.xaxis.set_major_formatter(formatter) + + # Set X-axis minor tick marks at noon of each day + # (i.e. split up the major interval into 2 bins) + ax0.xaxis.set_minor_locator(mticker.AutoMinorLocator(2)) + + # Don't rotate the X-axis jtick labels + ax0.xaxis.set_tick_params(rotation=0) + + # Center the X-axis tick labels + for tick in ax0.xaxis.get_major_ticks(): + tick.label1.set_horizontalalignment('center') + + # Set X-axis and Y-axis labels + ax0.set_xlabel('Day of July (and August) 2016') + + # ----------------------------- + # Save this page to PDF + # ----------------------------- + pdf.savefig(figs) + plt.close(figs) - # Center the X-axis tick labels - for tick in ax0.xaxis.get_major_ticks(): - tick.label1.set_horizontalalignment('center') + # ---------------------------------------------------------------------- + # Save the PDF file to disk + # ---------------------------------------------------------------------- + pdf.close() - # Set X-axis and Y-axis labels - ax0.set_xlabel('Day of July (and August) 2016') - - # ----------------------------- - # Save this page to PDF - # ----------------------------- - pdf.savefig(figs) - plt.close(figs) - # ---------------------------------------------------------------------- - # Save the PDF file to disk - # ---------------------------------------------------------------------- - pdf.close() + def main(): + ''' + Main program. + ''' + # Path where the data files live + # (YOU MUST EDIT THIS FOR YUR OWN PARTICULAR APPLICATION!) + path_to_data = '/path/to/GEOS-Chem/diagnostic/data/files' + # Get a list of files in the ConcAboveSfc and SpeciesConc collections + # (YOU CAN EDIT THIS FOR YOUR OWN PARTICULAR APPLICATION!) + collections = ['ConcAboveSfc', 'SpeciesConc'] - def main(): - ''' - Main program. - ''' - # Path where the data files live - # (YOU MUST EDIT THIS FOR YUR OWN PARTICULAR APPLICATION!) - path_to_data = '/path/to/GEOS-Chem/diagnostic/data/files' - - # Get a list of files in the ConcAboveSfc and SpeciesConc collections - # (YOU CAN EDIT THIS FOR YOUR OWN PARTICULAR APPLICATION!) - collections = ['ConcAboveSfc', 'SpeciesConc'] - - # Read GEOS-Chem data into an xarray Dataset - ds = read_geoschem_data(path_to_data, collections) + # Read GEOS-Chem data into an xarray Dataset + ds = read_geoschem_data(path_to_data, collections) - # Plot timeseries data at Centerville, AL (32.94N, 87.18W) - # (YOU CAN EDIT THIS FOR YOUR OWN PARTICULAR APPLICATION!) - site_coords = (32.94, -87.18) - plot_timeseries_data(ds, site_coords) + # Plot timeseries data at Centerville, AL (32.94N, 87.18W) + # (YOU CAN EDIT THIS FOR YOUR OWN PARTICULAR APPLICATION!) + site_coords = (32.94, -87.18) + plot_timeseries_data(ds, site_coords) - if __name__ == "__main__": - main() + if __name__ == "__main__": + main()