diff --git a/configs/anvil/config.20170926.FCT2.A_WCYCL1850S.ne30_oECv3.anvil b/configs/anvil/config.20170926.FCT2.A_WCYCL1850S.ne30_oECv3.anvil index f3f32244b..7892f47c2 100644 --- a/configs/anvil/config.20170926.FCT2.A_WCYCL1850S.ne30_oECv3.anvil +++ b/configs/anvil/config.20170926.FCT2.A_WCYCL1850S.ne30_oECv3.anvil @@ -59,9 +59,9 @@ htmlSubdirectory = html # option: # ./run_mpas_analysis config.analysis --generate \ # all,no_ocean,all_timeSeries -# The climatologyMapAntarcticMelt task is disabled because this run did not -# include ice-shelf cavities.` -generate = ['all', 'no_climatologyMapAntarcticMelt'] +# All tasks with tag "landIceCavities" are disabled because this run did not +# include land-ice cavities. +generate = ['all', 'no_landIceCavities'] [climatology] ## options related to producing climatologies, typically to compare against diff --git a/configs/cori/config.20171201.default.GMPAS-IAF.T62_oRRS30to10v3wLI.cori-knl b/configs/cori/config.20171201.default.GMPAS-IAF.T62_oRRS30to10v3wLI.cori-knl index 0043afc35..d21913142 100644 --- a/configs/cori/config.20171201.default.GMPAS-IAF.T62_oRRS30to10v3wLI.cori-knl +++ b/configs/cori/config.20171201.default.GMPAS-IAF.T62_oRRS30to10v3wLI.cori-knl @@ -164,6 +164,14 @@ depths = ['top', -200, -400, -600, -800, 'bot'] # plot on polar plot polarPlot = False +[timeSeriesAntarcticMelt] +## options related to plotting time series of melt below Antarctic ice shelves + +# list of ice shelves to plot or ['all'] for all 106 ice shelves and regions. +# See "regionNames" in the ice shelf masks file in regionMaskDirectory for +# details. +iceShelvesToPlot = ['Peninsula', 'West Antarctica', 'East Antarctica', 'Larsen_C', 'Filchner', 'Ronne', 'Brunt_Stancomb', 'Fimbul', 'Amery', 'Totten', 'Ross_West', 'Ross_East', 'Getz', 'Thwaites', 'Pine_Island', 'Abbot', 'George_VI'] + [regions] ## options related to ocean regions used in several analysis modules diff --git a/configs/edison/config.20170807.beta1.G_oQU240.edison b/configs/edison/config.20170807.beta1.G_oQU240.edison index a03932a01..161c18819 100644 --- a/configs/edison/config.20170807.beta1.G_oQU240.edison +++ b/configs/edison/config.20170807.beta1.G_oQU240.edison @@ -60,9 +60,9 @@ htmlSubdirectory = html # option: # ./run_mpas_analysis config.analysis --generate \ # all,no_ocean,all_timeSeries -# The climatologyMapAntarcticMelt task is disabled because this run did not -# include ice-shelf cavities.` -generate = ['all', 'no_climatologyMapAntarcticMelt'] +# All tasks with tag "landIceCavities" are disabled because this run did not +# include land-ice cavities. +generate = ['all', 'no_landIceCavities'] [climatology] ## options related to producing climatologies, typically to compare against diff --git a/configs/edison/config.20170915.beta2.A_WCYCL1850S.ne30_oECv3_ICG.edison b/configs/edison/config.20170915.beta2.A_WCYCL1850S.ne30_oECv3_ICG.edison index a4d2d124d..fb89b8531 100644 --- a/configs/edison/config.20170915.beta2.A_WCYCL1850S.ne30_oECv3_ICG.edison +++ b/configs/edison/config.20170915.beta2.A_WCYCL1850S.ne30_oECv3_ICG.edison @@ -75,9 +75,9 @@ htmlSubdirectory = html # option: # ./run_mpas_analysis config.analysis --generate \ # all,no_ocean,all_timeSeries -# The climatologyMapAntarcticMelt task is disabled because this run did not -# include ice-shelf cavities.` -generate = ['all', 'no_climatologyMapAntarcticMelt'] +# All tasks with tag "landIceCavities" are disabled because this run did not +# include land-ice cavities. +generate = ['all', 'no_landIceCavities'] [climatology] ## options related to producing climatologies, typically to compare against diff --git a/configs/edison/config.20171102.beta3rc02_1850.ne30_oECv3_ICG.edison b/configs/edison/config.20171102.beta3rc02_1850.ne30_oECv3_ICG.edison index 910214884..18fbffb70 100644 --- a/configs/edison/config.20171102.beta3rc02_1850.ne30_oECv3_ICG.edison +++ b/configs/edison/config.20171102.beta3rc02_1850.ne30_oECv3_ICG.edison @@ -61,9 +61,9 @@ htmlSubdirectory = html # option: # ./run_mpas_analysis config.analysis --generate \ # all,no_ocean,all_timeSeries -# The climatologyMapAntarcticMelt task is disabled because this run did not -# include ice-shelf cavities.` -generate = ['all', 'no_climatologyMapAntarcticMelt'] +# All tasks with tag "landIceCavities" are disabled because this run did not +# include land-ice cavities. +generate = ['all', 'no_landIceCavities'] [climatology] ## options related to producing climatologies, typically to compare against diff --git a/configs/edison/config.B_low_res_ice_shelves_1696_JWolfe_layout_Edison b/configs/edison/config.B_low_res_ice_shelves_1696_JWolfe_layout_Edison index 31fb86e4b..7d4bbeaed 100644 --- a/configs/edison/config.B_low_res_ice_shelves_1696_JWolfe_layout_Edison +++ b/configs/edison/config.B_low_res_ice_shelves_1696_JWolfe_layout_Edison @@ -168,6 +168,14 @@ depths = ['top', -200, -400, -600, -800, 'bot'] # plot on polar plot polarPlot = False +[timeSeriesAntarcticMelt] +## options related to plotting time series of melt below Antarctic ice shelves + +# list of ice shelves to plot or ['all'] for all 106 ice shelves and regions. +# See "regionNames" in the ice shelf masks file in regionMaskDirectory for +# details. +iceShelvesToPlot = ['Peninsula', 'West Antarctica', 'East Antarctica', 'Larsen_C', 'Filchner', 'Ronne', 'Brunt_Stancomb', 'Fimbul', 'Amery', 'Totten', 'Ross_West', 'Ross_East', 'Getz', 'Thwaites', 'Pine_Island', 'Abbot', 'George_VI'] + [regions] ## options related to ocean regions used in several analysis modules diff --git a/configs/lanl/config.MatchBoth_orig b/configs/lanl/config.MatchBoth_orig index 92ec3e87a..7ec34e297 100644 --- a/configs/lanl/config.MatchBoth_orig +++ b/configs/lanl/config.MatchBoth_orig @@ -55,9 +55,9 @@ baseDirectory = /dir/to/analysis/output # option: # ./run_mpas_analysis config.analysis --generate \ # all,no_ocean,all_timeSeries -# The climatologyMapAntarcticMelt task is disabled because this run did not -# include ice-shelf cavities.` -generate = ['all', 'no_climatologyMapAntarcticMelt'] +# All tasks with tag "landIceCavities" are disabled because this run did not +# include land-ice cavities. +generate = ['all', 'no_landIceCavities'] [climatology] ## options related to producing climatologies, typically to compare against diff --git a/configs/olcf/config.20170313.beta1.A_WCYCL1850S.ne30_oECv3_ICG.edison b/configs/olcf/config.20170313.beta1.A_WCYCL1850S.ne30_oECv3_ICG.edison index 76a6fc05a..b57a1eef4 100644 --- a/configs/olcf/config.20170313.beta1.A_WCYCL1850S.ne30_oECv3_ICG.edison +++ b/configs/olcf/config.20170313.beta1.A_WCYCL1850S.ne30_oECv3_ICG.edison @@ -60,9 +60,9 @@ htmlSubdirectory = html # option: # ./run_mpas_analysis config.analysis --generate \ # all,no_ocean,all_timeSeries -# The climatologyMapAntarcticMelt task is disabled because this run did not -# include ice-shelf cavities.` -generate = ['all', 'no_climatologyMapAntarcticMelt'] +# All tasks with tag "landIceCavities" are disabled because this run did not +# include land-ice cavities. +generate = ['all', 'no_landIceCavities'] [climatology] ## options related to producing climatologies, typically to compare against diff --git a/configs/olcf/config.20170915.beta2.A_WCYCL1850S.ne30_oECv3_ICG.edison b/configs/olcf/config.20170915.beta2.A_WCYCL1850S.ne30_oECv3_ICG.edison index 8ede7c504..49c4b749c 100644 --- a/configs/olcf/config.20170915.beta2.A_WCYCL1850S.ne30_oECv3_ICG.edison +++ b/configs/olcf/config.20170915.beta2.A_WCYCL1850S.ne30_oECv3_ICG.edison @@ -77,9 +77,9 @@ htmlSubdirectory = html # to list all task names and their tags # an equivalent syntax can be used on the command line to override this # option: -# ./run_mpas_analysis config.analysis --generate \ -# all,no_ocean,all_timeSeries -generate = ['all', 'no_climatologyMapAntarcticMelt'] +# All tasks with tag "landIceCavities" are disabled because this run did not +# include land-ice cavities. +generate = ['all', 'no_landIceCavities'] [climatology] ## options related to producing climatologies, typically to compare against diff --git a/configs/olcf/config.GMPAS-IAF_oRRS18to6v3.titan b/configs/olcf/config.GMPAS-IAF_oRRS18to6v3.titan index b903d1294..582c94e86 100644 --- a/configs/olcf/config.GMPAS-IAF_oRRS18to6v3.titan +++ b/configs/olcf/config.GMPAS-IAF_oRRS18to6v3.titan @@ -64,9 +64,9 @@ htmlSubdirectory = html # option: # ./run_mpas_analysis config.analysis --generate \ # all,no_ocean,all_timeSeries -# The climatologyMapAntarcticMelt task is disabled because this run did not -# include ice-shelf cavities.` -generate = ['all', 'no_climatologyMapAntarcticMelt'] +# All tasks with tag "landIceCavities" are disabled because this run did not +# include land-ice cavities. +generate = ['all', 'no_landIceCavities'] [climatology] ## options related to producing climatologies, typically to compare against diff --git a/configs/theta/config.20171031.tenYearTest.GMPAS-IAF.T62_oEC60to30v3wLI.60layer.theta b/configs/theta/config.20171031.tenYearTest.GMPAS-IAF.T62_oEC60to30v3wLI.60layer.theta index e295f0475..328e09307 100644 --- a/configs/theta/config.20171031.tenYearTest.GMPAS-IAF.T62_oEC60to30v3wLI.60layer.theta +++ b/configs/theta/config.20171031.tenYearTest.GMPAS-IAF.T62_oEC60to30v3wLI.60layer.theta @@ -172,4 +172,4 @@ depths = ['top', -200, -400, -600, -800, 'bot'] [timeSeriesAntarcticMelt] # a list of ice shelves to plot -iceShelvesToPlot = ['Peninsula', 'West Antarctica', 'East Antarctica', 'Larsen_C', 'Larsen_D', 'Ronne', 'Filchner', 'Brunt_Stancomb', 'Fimbul', 'Riiser-Larsen', 'Baudouin', 'Amery', 'West', 'Shackleton', 'Totten', 'Mertz', 'Ross_East', 'Ross_West', 'Getz', 'Dotson', 'Crosson', 'Thwaites', 'Pine_Island', 'Abbot', 'Stange', 'George_VI', 'Wilkins'] +iceShelvesToPlot = ['Antarctica', 'Peninsula', 'West Antarctica', 'East Antarctica', 'Larsen_C', 'Filchner-Ronne', 'Brunt_Stancomb', 'Fimbul', 'Amery', 'Totten', 'Ross', 'Getz', 'Thwaites', 'Pine_Island', 'Abbot', 'George_VI'] diff --git a/mpas_analysis/config.default b/mpas_analysis/config.default index e9ca8708e..0a39892b4 100644 --- a/mpas_analysis/config.default +++ b/mpas_analysis/config.default @@ -605,6 +605,17 @@ normArgsDifference = {'linthresh': 1., 'linscale': 0.5, 'vmin': -100., colorbarTicksDifference = [-100., -50., -20., -10., -5., -2., -1., 0., 1., 2., 5., 10., 20., 50., 100.] +[timeSeriesAntarcticMelt] +## options related to plotting time series of melt below Antarctic ice shelves + +# list of ice shelves to plot or ['all'] for all 106 ice shelves and regions. +# See "regionNames" in the ice shelf masks file in regionMaskDirectory for +# details. +iceShelvesToPlot = ['Antarctica'] + +# Number of months over which to compute moving average +movingAverageMonths = 1 + [climatologyMapSoseTemperature] ## options related to plotting climatology maps of Antarctic ## potential temperature at various levels, including the sea floor against diff --git a/mpas_analysis/ocean/__init__.py b/mpas_analysis/ocean/__init__.py index a9465cc98..06c4b722b 100644 --- a/mpas_analysis/ocean/__init__.py +++ b/mpas_analysis/ocean/__init__.py @@ -10,3 +10,4 @@ from .index_nino34 import IndexNino34 from .streamfunction_moc import StreamfunctionMOC from .meridional_heat_transport import MeridionalHeatTransport +from .time_series_antarctic_melt import TimeSeriesAntarcticMelt diff --git a/mpas_analysis/ocean/time_series_antarctic_melt.py b/mpas_analysis/ocean/time_series_antarctic_melt.py new file mode 100644 index 000000000..a50751eb1 --- /dev/null +++ b/mpas_analysis/ocean/time_series_antarctic_melt.py @@ -0,0 +1,356 @@ +from __future__ import absolute_import, division, print_function, \ + unicode_literals + +import os +import xarray +import numpy + +from ..shared.analysis_task import AnalysisTask + +from ..shared.constants import constants + +from ..shared.plot.plotting import timeseries_analysis_plot + +from ..shared.io import open_mpas_dataset + +from ..shared.io.utility import build_config_full_path, make_directories + +from ..shared.html import write_image_xml + +import csv + + +class TimeSeriesAntarcticMelt(AnalysisTask): + """ + Performs analysis of the time-series output of Antarctic sub-ice-shelf + melt rates. + + Attributes + ---------- + + mpasTimeSeriesTask : ``MpasTimeSeriesTask`` + The task that extracts the time series from MPAS monthly output + + Authors + ------- + Xylar Asay-Davis, Stephen Price + """ + + def __init__(self, config, mpasTimeSeriesTask): # {{{ + """ + Construct the analysis task. + + Parameters + ---------- + config : instance of MpasAnalysisConfigParser + Contains configuration options + + mpasTimeSeriesTask : ``MpasTimeSeriesTask`` + The task that extracts the time series from MPAS monthly output + + Authors + ------- + Xylar Asay-Davis + """ + # first, call the constructor from the base class (AnalysisTask) + super(TimeSeriesAntarcticMelt, self).__init__( + config=config, + taskName='timeSeriesAntarcticMelt', + componentName='ocean', + tags=['timeSeries', 'melt', 'landIceCavities']) + + self.mpasTimeSeriesTask = mpasTimeSeriesTask + + self.run_after(mpasTimeSeriesTask) + + # }}} + + def setup_and_check(self): # {{{ + """ + Perform steps to set up the analysis and check for errors in the setup. + + Raises + ------ + IOError + If files are not present + + Authors + ------- + Xylar Asay-Davis + """ + # first, call setup_and_check from the base class (AnalysisTask), + # which will perform some common setup, including storing: + # self.inDirectory, self.plotsDirectory, self.namelist, self.streams + # self.calendar + super(TimeSeriesAntarcticMelt, self).setup_and_check() + + self.check_analysis_enabled( + analysisOptionName='config_am_timeseriesstatsmonthly_enable', + raiseException=True) + + config = self.config + + landIceFluxMode = self.namelist.get('config_land_ice_flux_mode') + if landIceFluxMode not in ['standalone', 'coupled']: + raise ValueError('*** timeSeriesAntarcticMelt requires ' + 'config_land_ice_flux_mode \n' + ' to be standalone or coupled. Otherwise, no ' + 'melt rates are available \n' + ' for plotting.') + + mpasMeshName = config.get('input', 'mpasMeshName') + regionMaskDirectory = config.get('regions', 'regionMaskDirectory') + + self.regionMaskFileName = '{}/{}_iceShelfMasks.nc'.format( + regionMaskDirectory, mpasMeshName) + + if not os.path.exists(self.regionMaskFileName): + raise IOError('Regional masking file {} for Antarctica melt-rate ' + 'calculation does not exist'.format( + self.regionMaskFileName)) + + # Load mesh related variables + try: + self.restartFileName = self.runStreams.readpath('restart')[0] + except ValueError: + raise IOError('No MPAS-O restart file found: need at least one ' + 'restart file for Antarctic melt calculations') + + # get a list of timeSeriesStats output files from the streams file, + # reading only those that are between the start and end dates + self.startDate = config.get('timeSeries', 'startDate') + self.endDate = config.get('timeSeries', 'endDate') + + self.inputFile = self.mpasTimeSeriesTask.outputFile + + self.variableList = \ + ['timeMonthly_avg_landIceFreshwaterFlux'] + self.mpasTimeSeriesTask.add_variables(variableList=self.variableList) + + iceShelvesToPlot = config.getExpression('timeSeriesAntarcticMelt', + 'iceShelvesToPlot') + + with xarray.open_dataset(self.regionMaskFileName) as dsRegionMask: + regionNames = [bytes.decode(name) for name in + dsRegionMask.regionNames.values] + nRegions = dsRegionMask.dims['nRegions'] + + if 'all' in iceShelvesToPlot: + iceShelvesToPlot = regionNames + regionIndices = [iRegion for iRegion in range(nRegions)] + + else: + regionIndices = [] + for regionName in iceShelvesToPlot: + if regionName not in regionNames: + raise ValueError('Unknown ice shelf name {}'.format( + regionName)) + + iRegion = regionNames.index(regionName) + regionIndices.append(iRegion) + + self.regionIndices = regionIndices + self.iceShelvesToPlot = iceShelvesToPlot + self.xmlFileNames = [] + + for prefix in ['melt_flux', 'melt_rate']: + for regionName in iceShelvesToPlot: + regionName = regionName.replace(' ', '_') + self.xmlFileNames.append( + '{}/{}_{}.xml'.format(self.plotsDirectory, prefix, + regionName)) + return # }}} + + def run_task(self): # {{{ + """ + Performs analysis of the time-series output of Antarctic sub-ice-shelf + melt rates. + + Authors + ------- + Xylar Asay-Davis + """ + + self.logger.info("\nPlotting Antarctic melt rate time series...") + + self.logger.info(' Load melt rate data...') + + config = self.config + calendar = self.calendar + + # Load data: + ds = open_mpas_dataset(fileName=self.inputFile, + calendar=calendar, + variableList=self.variableList, + startDate=self.startDate, + endDate=self.endDate) + + # Load observations from multiple files and put in dictionary based + # on shelf keyname + observationsDirectory = build_config_full_path(config, + 'oceanObservations', + 'meltSubdirectory') + obsFileNameDict = {'Rignot et al. (2013)': + 'Rignot_2013_melt_rates.csv', + 'Rignot et al. (2013) SS': + 'Rignot_2013_melt_rates_SS.csv'} + + obsDict = {} # dict for storing dict of obs data + for obsName in obsFileNameDict: + obsFileName = '{}/{}'.format(observationsDirectory, + obsFileNameDict[obsName]) + obsDict[obsName] = {} + obsFile = csv.reader(open(obsFileName, 'rU')) + next(obsFile, None) # skip the header line + for line in obsFile: # some later useful values commented out + shelfName = line[0] + # surveyArea = line[1] + meltFlux = float(line[2]) + meltFluxUncertainty = float(line[3]) + meltRate = float(line[4]) + meltRateUncertainty = float(line[5]) + # actualArea = float( line[6] ) # actual area here is in sq km + + # build dict of obs. keyed to filename description + # (which will be used for plotting) + obsDict[obsName][shelfName] = { + 'meltFlux': meltFlux, + 'meltFluxUncertainty': meltFluxUncertainty, + 'meltRate': meltRate, + 'meltRateUncertainty': meltRateUncertainty} + + # If areas from obs file used need to be converted from sq km to sq m + + # work on data from simulations + freshwaterFlux = ds.timeMonthly_avg_landIceFreshwaterFlux + + mainRunName = config.get('runs', 'mainRunName') + movingAverageMonths = config.getint('timeSeriesAntarcticMelt', + 'movingAverageMonths') + + iceShelvesToPlot = self.iceShelvesToPlot + + dsRestart = xarray.open_dataset(self.restartFileName) + areaCell = dsRestart.landIceFraction.isel(Time=0)*dsRestart.areaCell + + dsRegionMask = xarray.open_dataset(self.regionMaskFileName) + + # select only those regions we want to plot + dsRegionMask = dsRegionMask.isel(nRegions=self.regionIndices) + cellMasks = dsRegionMask.regionCellMasks + nRegions = dsRegionMask.dims['nRegions'] + + # convert from kg/s to kg/yr + totalMeltFlux = constants.sec_per_year * \ + (cellMasks*areaCell*freshwaterFlux).sum(dim='nCells') + + totalArea = (cellMasks*areaCell).sum(dim='nCells') + + # from kg/m^2/yr to m/yr + meltRates = (1./constants.rho_fw) * (totalMeltFlux/totalArea) + + # convert from kg/yr to GT/yr + totalMeltFlux /= constants.kg_per_GT + + outputDirectory = build_config_full_path(config, 'output', + 'timeseriesSubdirectory') + + make_directories(outputDirectory) + + self.logger.info(' Make plots...') + for iRegion in range(nRegions): + + regionName = iceShelvesToPlot[iRegion] + + # get obs melt flux and unc. for shelf (similar for rates) + obsMeltFlux = [] + obsMeltFluxUnc = [] + obsMeltRate = [] + obsMeltRateUnc = [] + for obsName in obsDict: + if regionName in obsDict[obsName]: + obsMeltFlux.append( + obsDict[obsName][regionName]['meltFlux']) + obsMeltFluxUnc.append( + obsDict[obsName][regionName]['meltFluxUncertainty']) + obsMeltRate.append( + obsDict[obsName][regionName]['meltRate']) + obsMeltRateUnc.append( + obsDict[obsName][regionName]['meltRateUncertainty']) + else: + # append NaN so this particular obs won't plot + self.logger.warning('{} observations not available for ' + '{}'.format(obsName, regionName)) + obsMeltFlux.append(None) + obsMeltFluxUnc.append(None) + obsMeltRate.append(None) + obsMeltRateUnc.append(None) + + title = regionName.replace('_', ' ') + + regionName = regionName.replace(' ', '_') + + xLabel = 'Time (yr)' + yLabel = 'Melt Flux (GT/yr)' + + timeSeries = totalMeltFlux.isel(nRegions=iRegion) + + filePrefix = 'melt_flux_{}'.format(regionName) + figureName = '{}/{}.png'.format(self.plotsDirectory, filePrefix) + + timeseries_analysis_plot(config, [timeSeries], movingAverageMonths, + title, xLabel, yLabel, figureName, + lineStyles=['b-'], lineWidths=[1.2], + legendText=[mainRunName], + calendar=calendar, obsMean=obsMeltFlux, + obsUncertainty=obsMeltFluxUnc, + obsLegend=list(obsDict.keys())) + + caption = 'Running Mean of Total Melt Flux under Ice ' \ + 'Shelves in the {} Region'.format(title) + write_image_xml( + config=config, + filePrefix=filePrefix, + componentName='Ocean', + componentSubdirectory='ocean', + galleryGroup='Antarctic Melt Time Series', + groupLink='antmelttime', + gallery='Total Melt Flux', + thumbnailDescription=title, + imageDescription=caption, + imageCaption=caption) + + xLabel = 'Time (yr)' + yLabel = 'Melt Rate (m/yr)' + + timeSeries = meltRates.isel(nRegions=iRegion) + + filePrefix = 'melt_rate_{}'.format(regionName) + figureName = '{}/{}.png'.format(self.plotsDirectory, filePrefix) + + timeseries_analysis_plot(config, [timeSeries], movingAverageMonths, + title, xLabel, yLabel, figureName, + lineStyles=['b-'], lineWidths=[1.2], + legendText=[mainRunName], + calendar=calendar, obsMean=obsMeltRate, + obsUncertainty=obsMeltRateUnc, + obsLegend=list(obsDict.keys())) + + caption = 'Running Mean of Area-averaged Melt Rate under Ice ' \ + 'Shelves in the {} Region'.format(title) + write_image_xml( + config=config, + filePrefix=filePrefix, + componentName='Ocean', + componentSubdirectory='ocean', + galleryGroup='Antarctic Melt Time Series', + groupLink='antmelttime', + gallery='Area-averaged Melt Rate', + thumbnailDescription=title, + imageDescription=caption, + imageCaption=caption) + # }}} + +# }}} + +# vim: foldmethod=marker ai ts=4 sts=4 et sw=4 ft=python diff --git a/mpas_analysis/shared/constants/constants.py b/mpas_analysis/shared/constants/constants.py index d649a51f9..9f901aadf 100644 --- a/mpas_analysis/shared/constants/constants.py +++ b/mpas_analysis/shared/constants/constants.py @@ -61,4 +61,7 @@ # density of freshwater (kg/m^3) rho_fw = 1000. +# kilograms per gigatonne +kg_per_GT = 1e12 + # vim: foldmethod=marker ai ts=4 sts=4 et sw=4 ft=python diff --git a/mpas_analysis/shared/plot/plotting.py b/mpas_analysis/shared/plot/plotting.py index 78d72dbb4..33d5258b9 100644 --- a/mpas_analysis/shared/plot/plotting.py +++ b/mpas_analysis/shared/plot/plotting.py @@ -33,7 +33,9 @@ def timeseries_analysis_plot(config, dsvalues, N, title, xlabel, ylabel, fileout, lineStyles, lineWidths, legendText, calendar, maxPoints=None, titleFontSize=None, - figsize=(15, 6), dpi=None, maxXTicks=20): + figsize=(15, 6), dpi=None, maxXTicks=20, + obsMean=None, obsUncertainty=None, + obsLegend=None, legendLocation='lower left'): """ Plots the list of time series data sets and stores the result in an image @@ -87,9 +89,20 @@ def timeseries_analysis_plot(config, dsvalues, N, title, xlabel, ylabel, This may need to be adjusted depending on the figure size and aspect ratio. + obsMean, obsUncertainty : list of float, optional + Mean values and uncertainties for observations to be plotted as error + bars. The two lists must have the same number of elements. + + obsLegend : list of str, optional + The label in the legend for each element in ``obsMean`` (and + ``obsUncertainty``) + + legendLocation : str, optional + The location of the legend (see ``pyplot.legend()`` for details) + Authors ------- - Xylar Asay-Davis, Milena Veneziani + Xylar Asay-Davis, Milena Veneziani, Stephen Price """ if dpi is None: @@ -117,7 +130,26 @@ def timeseries_analysis_plot(config, dsvalues, N, title, xlabel, ylabel, lineStyles[dsIndex], linewidth=lineWidths[dsIndex], label=legendText[dsIndex]) - plt.legend(loc='lower left') + + if obsMean is not None: + obsCount = len(obsMean) + assert(len(obsUncertainty) == obsCount) + + # space the observations along the time line, leaving gaps at either + # end + start = np.amin(minDays) + end = np.amax(maxDays) + obsTimes = np.linspace(start, end, obsCount+2)[1:-1] + obsSymbols = ['o', '^', 's', 'D', '*'] + for iObs in range(obsCount): + if obsMean[iObs] is not None: + plt.errorbar(obsTimes[iObs], obsMean[iObs], + yerr=obsUncertainty[iObs], + fmt=obsSymbols[np.mod(iObs, len(obsSymbols))], + ecolor='k', + capthick=2, label=obsLegend[iObs]) + + plt.legend(loc=legendLocation) ax = plt.gca() @@ -133,7 +165,7 @@ def timeseries_analysis_plot(config, dsvalues, N, title, xlabel, ylabel, if yaxLimits[0]*yaxLimits[1] < 0: indgood = np.where(np.logical_not(np.isnan(mean))) x = mean['Time'][indgood] - plt.plot(x, np.zeros(np.size(x)), 'k-', linewidth=1.2) + plt.plot(x, np.zeros(np.size(x)), 'k-', linewidth=1.2, zorder=1) plot_xtick_format(plt, calendar, minDays, maxDays, maxXTicks) diff --git a/preprocess_masks/make_ice_shelf_masks.py b/preprocess_masks/make_ice_shelf_masks.py new file mode 100755 index 000000000..4e4136506 --- /dev/null +++ b/preprocess_masks/make_ice_shelf_masks.py @@ -0,0 +1,40 @@ +#!/usr/bin/env python + +''' +Make a mask file for a given mesh from ice-shelf geometric features. + +The -m flag is used to specify the name of the ACME mesh to which the +masks should be applied. + +Requires: + * a local link to the MPAS mask creator MpasMaskCreator.x + * a local link to a mesh file named _mesh.nc describing the + desired mesh + * the region file iceShelves.geojson produced by running + ./driver_scripts/setup_ice_shelves.py in the geometric_features repo + +Produces: + * _iceShelfMasks.nc, the mask file + +Author: Xylar Asay-Davis +''' + +import subprocess +import argparse + + +parser = \ + argparse.ArgumentParser(description=__doc__, + formatter_class=argparse.RawTextHelpFormatter) +parser.add_argument('-m', '--mesh_name', dest='mesh_name', + help='The ACME name of the mesh', metavar='MESH_NAME', + required=True) +args = parser.parse_args() + +meshFileName = '{}_mesh.nc'.format(args.mesh_name) +maskFileName = '{}_iceShelfMasks.nc'.format(args.mesh_name) +regionFileName = 'iceShelves.geojson' + +tempRegionMaskFile = 'tempRegionMasks.nc' +subprocess.check_call(['./MpasMaskCreator.x', meshFileName, maskFileName, + '-f', regionFileName]) diff --git a/run_mpas_analysis b/run_mpas_analysis index b3f048e76..f1b13b590 100755 --- a/run_mpas_analysis +++ b/run_mpas_analysis @@ -77,14 +77,15 @@ def build_analysis_list(config): # {{{ analyses.append(ocean.ClimatologyMapMLD(config, oceanClimatolgyTask)) analyses.append(ocean.ClimatologyMapSST(config, oceanClimatolgyTask)) analyses.append(ocean.ClimatologyMapSSS(config, oceanClimatolgyTask)) - analyses.append(ocean.ClimatologyMapAntarcticMelt(config, - oceanClimatolgyTask)) analyses.append(ocean.ClimatologyMapSoseTemperature(config, oceanClimatolgyTask)) analyses.append(ocean.ClimatologyMapSoseSalinity(config, oceanClimatolgyTask)) + analyses.append(ocean.ClimatologyMapAntarcticMelt(config, + oceanClimatolgyTask)) analyses.append(oceanTimeSeriesTask) + analyses.append(ocean.TimeSeriesAntarcticMelt(config, oceanTimeSeriesTask)) analyses.append(ocean.TimeSeriesOHC(config, oceanTimeSeriesTask)) analyses.append(ocean.TimeSeriesSST(config, oceanTimeSeriesTask)) analyses.append(ocean.MeridionalHeatTransport(config, oceanClimatolgyTask))