# This software is open source software available under the BSD-3 license.
#
# Copyright (c) 2020 Triad National Security, LLC. All rights reserved.
# Copyright (c) 2020 Lawrence Livermore National Security, LLC. All rights
# reserved.
# Copyright (c) 2020 UT-Battelle, LLC. All rights reserved.
#
# Additional copyright and license information can be found in the LICENSE file
# distributed with this code, or at
# https://raw.githubusercontent.com/MPAS-Dev/MPAS-Analysis/master/LICENSE
"""
An analysis subtasks for plotting comparison of 2D model fields against
observations.
"""
# Authors
# -------
# Luke Van Roekel, Xylar Asay-Davis, Milena Veneziani
from __future__ import absolute_import, division, print_function, \
unicode_literals
import xarray as xr
import numpy as np
from pyremap.descriptor import interp_extrap_corner
from mpas_analysis.shared import AnalysisTask
from mpas_analysis.shared.plot import plot_global_comparison, \
plot_polar_projection_comparison
from mpas_analysis.shared.html import write_image_xml
from mpas_analysis.shared.climatology import \
get_remapped_mpas_climatology_file_name
from mpas_analysis.shared.climatology.comparison_descriptors import \
get_comparison_descriptor
from mpas_analysis.ocean.utility import nans_to_numpy_mask
[docs]class PlotClimatologyMapSubtask(AnalysisTask): # {{{
"""
An analysis task for plotting 2D model fields against observations.
Attributes
----------
season : str
A season (key in ``shared.constants.monthDictionary``) to be
plotted.
comparisonGridName : {'latlon', 'antarctic', 'arctic'}
The name of the comparison grid to plot.
remapMpasClimatologySubtask : ``RemapMpasClimatologySubtask``
The subtask for remapping the MPAS climatology that this subtask
will plot
remapObsClimatologySubtask : ``RemapObservedClimatologySubtask``
The subtask for remapping the observational climatology that this
subtask will plot
secondRemapMpasClimatologySubtask : ``RemapMpasClimatologySubtask``
A second subtask for remapping another MPAS climatology to plot
in the second panel and compare with in the third panel
removeMean : bool, optional
If True, a common mask for the model and reference data sets is
computed (where both are valid) and the mean over that mask is
subtracted from both the model and reference results. This is
useful for data sets where the desire is to compare the spatial
pattern but the mean offset is not meaningful (e.g. SSH)
outFileLabel : str
The prefix on each plot and associated XML file
fieldNameInTitle : str
The name of the field being plotted, as used in the plot title
mpasFieldName : str
The name of the variable in the MPAS timeSeriesStatsMonthly output
obsFieldName : str
The name of the variable to use from the observations file
observationTitleLabel : str
the title of the observations subplot
diffTitleLabel : str, optional
the title of the difference subplot
unitsLabel : str
the units of the plotted field, to be displayed on color bars
imageCaption : str
the caption when mousing over the plot or displaying it full
screen
galleryGroup : str
the name of the group of galleries in which this plot belongs
groupSubtitle : str
the subtitle of the group in which this plot belongs (or blank
if none)
groupLink : str
a short name (with no spaces) for the link to the gallery group
galleryName : str
the name of the gallery in which this plot belongs
depth : {None, float, 'top', 'bot'}
Depth at which to perform the comparison, 'top' for the sea surface
'bot' for the sea floor
configSectionName : str
the name of the section where the color map and range is defined
"""
# Authors
# -------
# Luke Van Roekel, Xylar Asay-Davis, Milena Veneziani
[docs] def __init__(self, parentTask, season, comparisonGridName,
remapMpasClimatologySubtask, remapObsClimatologySubtask=None,
secondRemapMpasClimatologySubtask=None, controlConfig=None,
depth=None, removeMean=False, subtaskName=None): # {{{
'''
Construct one analysis subtask for each plot (i.e. each season and
comparison grid) and a subtask for computing climatologies.
Parameters
----------
parentTask : ``AnalysisTask``
The parent (master) task for this subtask
season : str
A season (key in ``shared.constants.monthDictionary``) to be
plotted.
comparisonGridName : {'latlon', 'antarctic', 'arctic'}
The name of the comparison grid to plot.
remapMpasClimatologySubtask : ``RemapMpasClimatologySubtask``
The subtask for remapping the MPAS climatology that this subtask
will plot
remapObsClimatologySubtask : ``RemapObservedClimatologySubtask``, optional
The subtask for remapping the observational climatology that this
subtask will plot
secondRemapMpasClimatologySubtask : ``RemapMpasClimatologySubtask``, optional
A second subtask for remapping another MPAS climatology to plot
in the second panel and compare with in the third panel
controlConfig : ``MpasAnalysisConfigParser``, optional
Configuration options for a control run (if any)
depth : {float, 'top', 'bot'}, optional
Depth the data is being plotted, 'top' for the sea surface
'bot' for the sea floor
removeMean : bool, optional
If True, a common mask for the model and reference data sets is
computed (where both are valid) and the mean over that mask is
subtracted from both the model and reference results. This is
useful for data sets where the desire is to compare the spatial
pattern but the mean offset is not meaningful (e.g. SSH)
subtaskName : str, optinal
The name of the subtask. If not specified, it is
``plot<season>_<comparisonGridName>`` with a suffix indicating the
depth being sliced (if any)
'''
# Authors
# -------
# Xylar Asay-Davis
self.season = season
self.depth = depth
self.comparisonGridName = comparisonGridName
self.remapMpasClimatologySubtask = remapMpasClimatologySubtask
self.remapObsClimatologySubtask = remapObsClimatologySubtask
self.secondRemapMpasClimatologySubtask = \
secondRemapMpasClimatologySubtask
self.controlConfig = controlConfig
self.removeMean = removeMean
if depth is None:
self.depthSuffix = ''
else:
self.depthSuffix = 'depth_{}'.format(depth)
if subtaskName is None:
subtaskName = 'plot{}_{}'.format(season, comparisonGridName)
if depth is not None:
subtaskName = '{}_{}'.format(subtaskName, self.depthSuffix)
config = parentTask.config
taskName = parentTask.taskName
tags = parentTask.tags
# call the constructor from the base class (AnalysisTask)
super(PlotClimatologyMapSubtask, self).__init__(
config=config, taskName=taskName, subtaskName=subtaskName,
componentName='ocean', tags=tags)
# this task should not run until the remapping subtasks are done, since
# it relies on data from those subtasks
self.run_after(remapMpasClimatologySubtask)
if remapObsClimatologySubtask is not None:
self.run_after(remapObsClimatologySubtask)
if secondRemapMpasClimatologySubtask is not None:
self.run_after(secondRemapMpasClimatologySubtask)
# }}}
[docs] def set_plot_info(self, outFileLabel, fieldNameInTitle, mpasFieldName,
refFieldName, refTitleLabel, unitsLabel,
imageCaption, galleryGroup, groupSubtitle, groupLink,
galleryName, diffTitleLabel='Model - Observations',
configSectionName=None):
# {{{
"""
Store attributes related to plots, plot file names and HTML output.
Parameters
----------
outFileLabel : str
The prefix on each plot and associated XML file
fieldNameInTitle : str
The name of the field being plotted, as used in the plot title
mpasFieldName : str
The name of the variable in the MPAS timeSeriesStatsMonthly output
refFieldName : str
The name of the variable to use from the observations or reference
file
refTitleLabel : str
the title of the observations or reference subplot
unitsLabel : str
the units of the plotted field, to be displayed on color bars
imageCaption : str
the caption when mousing over the plot or displaying it full
screen
galleryGroup : str
the name of the group of galleries in which this plot belongs
groupSubtitle : str
the subtitle of the group in which this plot belongs (or blank
if none)
groupLink : str
a short name (with no spaces) for the link to the gallery group
galleryName : str
the name of the gallery in which this plot belongs
diffTitleLabel : str, optional
the title of the difference subplot
configSectionName : str, optional
the name of the section where the color map and range is defined,
default is the name of the task
"""
# Authors
# -------
# Xylar Asay-Davis
self.outFileLabel = outFileLabel
self.fieldNameInTitle = fieldNameInTitle
self.mpasFieldName = mpasFieldName
self.refFieldName = refFieldName
self.refTitleLabel = refTitleLabel
self.diffTitleLabel = diffTitleLabel
self.unitsLabel = unitsLabel
# xml/html related variables
self.imageCaption = imageCaption
self.galleryGroup = galleryGroup
self.groupSubtitle = groupSubtitle
self.groupLink = groupLink
self.galleryName = galleryName
if configSectionName is None:
self.configSectionName = self.taskName
else:
self.configSectionName = configSectionName
season = self.season
depth = self.depth
if depth is None:
self.fieldNameInTitle = fieldNameInTitle
self.thumbnailDescription = season
elif depth == 'top':
self.fieldNameInTitle = 'Sea Surface {}'.format(fieldNameInTitle)
self.thumbnailDescription = '{} surface'.format(season)
elif depth == 'bot':
self.fieldNameInTitle = 'Sea Floor {}'.format(fieldNameInTitle)
self.thumbnailDescription = '{} floor'.format(season)
else:
self.fieldNameInTitle = '{} at z={} m'.format(fieldNameInTitle,
depth)
self.thumbnailDescription = '{} z={} m'.format(season, depth)
# }}}
def setup_and_check(self): # {{{
"""
Perform steps to set up the analysis and check for errors in the setup.
"""
# Authors
# -------
# Xylar Asay-Davis
# first, call setup_and_check from the base class (AnalysisTask),
# which will perform some common setup, including storing:
# self.runDirectory , self.historyDirectory, self.plotsDirectory,
# self.namelist, self.runStreams, self.historyStreams,
# self.calendar
super(PlotClimatologyMapSubtask, self).setup_and_check()
config = self.config
self.startYear = config.getint('climatology', 'startYear')
self.endYear = config.getint('climatology', 'endYear')
self.startDate = config.get('climatology', 'startDate')
self.endDate = config.get('climatology', 'endDate')
mainRunName = config.get('runs', 'mainRunName')
self.xmlFileNames = []
prefixPieces = [self.outFileLabel]
if self.comparisonGridName != 'latlon':
prefixPieces.append(self.comparisonGridName)
prefixPieces.append(mainRunName)
if self.depth is not None:
prefixPieces.append(self.depthSuffix)
years = 'years{:04d}-{:04d}'.format(self.startYear, self.endYear)
prefixPieces.extend([self.season, years])
self.filePrefix = '_'.join(prefixPieces)
self.xmlFileNames.append('{}/{}.xml'.format(self.plotsDirectory,
self.filePrefix))
# }}}
def run_task(self): # {{{
"""
Plots a comparison of E3SM/MPAS output to SST/TEMP, SSS/SALT or MLD
observations or a control run
"""
# Authors
# -------
# Luke Van Roekel, Xylar Asay-Davis, Milena Veneziani
season = self.season
depth = self.depth
comparisonGridName = self.comparisonGridName
self.logger.info("\nPlotting 2-d maps of {} climatologies for {} on "
"the {} grid...".format(self.fieldNameInTitle,
season, comparisonGridName))
# first read the model climatology
remappedFileName = \
self.remapMpasClimatologySubtask.get_remapped_file_name(
season=season, comparisonGridName=comparisonGridName)
remappedModelClimatology = xr.open_dataset(remappedFileName)
if depth is not None:
if str(depth) not in remappedModelClimatology.depthSlice.values:
raise KeyError('The climatology you are attempting to perform '
'depth slices of was originally created '
'without depth {}. You will need to delete and '
'regenerate the climatology'.format(depth))
remappedModelClimatology = remappedModelClimatology.sel(
depthSlice=str(depth), drop=True)
# now the observations or control run
if self.remapObsClimatologySubtask is not None:
remappedFileName = self.remapObsClimatologySubtask.get_file_name(
stage='remapped', season=season,
comparisonGridName=comparisonGridName)
remappedRefClimatology = xr.open_dataset(remappedFileName)
elif self.secondRemapMpasClimatologySubtask is not None:
remappedFileName = \
self.secondRemapMpasClimatologySubtask.get_remapped_file_name(
season=season, comparisonGridName=comparisonGridName)
remappedRefClimatology = xr.open_dataset(remappedFileName)
elif self.controlConfig is not None:
climatologyName = self.remapMpasClimatologySubtask.climatologyName
remappedFileName = \
get_remapped_mpas_climatology_file_name(
self.controlConfig, season=season,
componentName=self.componentName,
climatologyName=climatologyName,
comparisonGridName=comparisonGridName,
op=self.remapMpasClimatologySubtask.op)
remappedRefClimatology = xr.open_dataset(remappedFileName)
controlStartYear = self.controlConfig.getint('climatology',
'startYear')
controlEndYear = self.controlConfig.getint('climatology',
'endYear')
if controlStartYear != self.startYear or \
controlEndYear != self.endYear:
self.refTitleLabel = '{}\n(years {:04d}-{:04d})'.format(
self.refTitleLabel, controlStartYear, controlEndYear)
else:
remappedRefClimatology = None
if remappedRefClimatology is not None and depth is not None:
depthIndex = -1
for index, depthSlice in enumerate(
remappedRefClimatology.depthSlice.values):
try:
depthSlice = depthSlice.decode("utf-8")
except AttributeError:
pass
if depthSlice == str(depth):
depthIndex = index
if depthIndex == -1:
raise KeyError('The climatology you are attempting to perform '
'depth slices of was originally created '
'without depth {}. You will need to delete and '
'regenerate the climatology'.format(depth))
remappedRefClimatology = remappedRefClimatology.isel(
depthSlice=depthIndex, drop=True)
if self.removeMean:
if remappedRefClimatology is None:
remappedModelClimatology[self.mpasFieldName] = \
remappedModelClimatology[self.mpasFieldName] - \
remappedModelClimatology[self.mpasFieldName].mean()
else:
masked = remappedModelClimatology[self.mpasFieldName].where(
remappedRefClimatology[self.refFieldName].notnull())
remappedModelClimatology[self.mpasFieldName] = \
remappedModelClimatology[self.mpasFieldName] - \
masked.mean()
masked = remappedRefClimatology[self.refFieldName].where(
remappedModelClimatology[self.mpasFieldName].notnull())
remappedRefClimatology[self.refFieldName] = \
remappedRefClimatology[self.refFieldName] - \
masked.mean()
if self.comparisonGridName == 'latlon':
self._plot_latlon(remappedModelClimatology, remappedRefClimatology)
elif self.comparisonGridName == 'antarctic' or \
self.comparisonGridName == 'arctic':
self._plot_polar(remappedModelClimatology,
remappedRefClimatology)
# }}}
def _plot_latlon(self, remappedModelClimatology, remappedRefClimatology):
# {{{
""" plotting a global lat-lon data set """
season = self.season
config = self.config
configSectionName = self.configSectionName
mainRunName = config.get('runs', 'mainRunName')
modelOutput = nans_to_numpy_mask(
remappedModelClimatology[self.mpasFieldName].values)
lon = remappedModelClimatology['lon'].values
lat = remappedModelClimatology['lat'].values
lonTarg, latTarg = np.meshgrid(lon, lat)
if remappedRefClimatology is None:
refOutput = None
bias = None
else:
refOutput = nans_to_numpy_mask(
remappedRefClimatology[self.refFieldName].values)
bias = modelOutput - refOutput
filePrefix = self.filePrefix
outFileName = '{}/{}.png'.format(self.plotsDirectory, filePrefix)
title = '{} ({}, years {:04d}-{:04d})'.format(
self.fieldNameInTitle, season, self.startYear,
self.endYear)
plot_global_comparison(config,
lonTarg,
latTarg,
modelOutput,
refOutput,
bias,
configSectionName,
fileout=outFileName,
title=title,
modelTitle='{}'.format(mainRunName),
refTitle=self.refTitleLabel,
diffTitle=self.diffTitleLabel,
cbarlabel=self.unitsLabel)
caption = '{} {}'.format(season, self.imageCaption)
write_image_xml(
config,
filePrefix,
componentName='Ocean',
componentSubdirectory='ocean',
galleryGroup='Global {}'.format(self.galleryGroup),
groupSubtitle=self.groupSubtitle,
groupLink='global_{}'.format(self.groupLink),
gallery=self.galleryName,
thumbnailDescription=self.thumbnailDescription,
imageDescription=caption,
imageCaption=caption)
# }}}
def _plot_polar(self, remappedModelClimatology,
remappedRefClimatology): # {{{
""" plotting an Arctic or Antarctic data set """
season = self.season
comparisonGridName = self.comparisonGridName
config = self.config
configSectionName = self.configSectionName
mainRunName = config.get('runs', 'mainRunName')
oceanMask = remappedModelClimatology['validMask'].values
self.landMask = np.ma.masked_array(
np.ones(oceanMask.shape),
mask=np.logical_not(np.isnan(oceanMask)))
modelOutput = nans_to_numpy_mask(
remappedModelClimatology[self.mpasFieldName].values)
if remappedRefClimatology is None:
refOutput = None
bias = None
else:
refOutput = nans_to_numpy_mask(
remappedRefClimatology[self.refFieldName].values)
bias = modelOutput - refOutput
comparisonDescriptor = get_comparison_descriptor(
config, comparisonGridName)
x = comparisonDescriptor.xCorner
y = comparisonDescriptor.yCorner
filePrefix = self.filePrefix
outFileName = '{}/{}.png'.format(self.plotsDirectory, filePrefix)
title = '{} ({}, years {:04d}-{:04d})'.format(
self.fieldNameInTitle, season, self.startYear,
self.endYear)
if self.comparisonGridName == 'antarctic':
hemisphere = 'south'
else:
# arctic
hemisphere = 'north'
plot_polar_projection_comparison(
config,
x,
y,
self.landMask,
modelOutput,
refOutput,
bias,
fileout=outFileName,
colorMapSectionName=configSectionName,
title=title,
modelTitle='{}'.format(mainRunName),
refTitle=self.refTitleLabel,
diffTitle=self.diffTitleLabel,
cbarlabel=self.unitsLabel,
hemisphere=hemisphere)
upperGridName = comparisonGridName[0].upper() + comparisonGridName[1:]
caption = '{} {}'.format(season, self.imageCaption)
write_image_xml(
config,
filePrefix,
componentName='Ocean',
componentSubdirectory='ocean',
galleryGroup='{} {}'.format(upperGridName,
self.galleryGroup),
groupSubtitle=self.groupSubtitle,
groupLink='{}_{}'.format(comparisonGridName, self.groupLink),
gallery=self.galleryName,
thumbnailDescription=self.thumbnailDescription,
imageDescription=caption,
imageCaption=caption)
# }}}
# }}}
# vim: foldmethod=marker ai ts=4 sts=4 et sw=4 ft=python