Skip to content

Update plotting routines in QC script #355

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
apcraig merged 5 commits into from
Aug 24, 2019
Merged

Update plotting routines in QC script #355

Changes from 2 commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
bd74075
Add an option to plot the ice thickness data on a map
Aug 22, 2019
971dbad
Update the t-test script to include a new option to plot the mean ice…
Aug 22, 2019
52c2076
Modify the QC script to always try to plot the data on a map. Add ne…
Aug 23, 2019
dde0402
Add a brief paragraph about the plotting functionality of the QC script
Aug 23, 2019
c22e32b
Fix typo in documentation
Aug 23, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
74 changes: 56 additions & 18 deletions configuration/scripts/tests/QC/cice.t-test.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -358,48 +358,63 @@ def skill_test(path_a, fname, data_a, data_b, num_files, hemisphere):
logger.info('Quadratic Skill Test Failed for %s Hemisphere', hemisphere)
return False

def plot_data(data, lat, lon, units):
def plot_data(data, lat, lon, units, case):
'''This function plots CICE data and creates a .png file (ice_thickness_map.png).'''

try:
logger.info('Creating map of the data (ice_thickness_map.png)')
# Load the necessary plotting libraries
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
from mpl_toolkits.axes_grid1 import make_axes_locatable
except ImportError:
logger.warning('Error loading necessary Python modules in plot_data function')
return

# Suppress Matplotlib deprecation warnings
import warnings
warnings.filterwarnings("ignore", category=UserWarning)

# Create the figure and axis
fig = plt.figure(figsize=(12, 8))
ax = fig.add_axes([0.05, 0.08, 0.9, 0.9])
fig, axes = plt.subplots(nrows=1, ncols=2,figsize=(14, 8))

# Plot the northern hemisphere data as a scatter plot
# Create the basemap, and draw boundaries
m = Basemap(projection='kav7', lon_0=180., resolution='l')
m.drawmapboundary(fill_color='white')
plt.sca(axes[0])
m = Basemap(projection='npstere', boundinglat=35,lon_0=270, resolution='l')
m.drawcoastlines()
m.fillcontinents()
m.drawcountries()

# Plot the data as a scatter plot
x, y = m(lon, lat)
sc = m.scatter(x, y, c=data, cmap='jet', lw=0)
sc = m.scatter(x, y, c=data, cmap='jet', lw=0, s=4)

m.drawmeridians(np.arange(0, 360, 60), labels=[0, 0, 0, 1], fontsize=10)
m.drawparallels(np.arange(-90, 90, 30), labels=[1, 0, 0, 0], fontsize=10)
m.drawparallels(np.arange(-90.,120.,15.),labels=[1,0,0,0]) # draw parallels
m.drawmeridians(np.arange(0.,420.,30.),labels=[1,1,1,1]) # draw meridians

plt.title('CICE Ice Thickness')
# Plot the southern hemisphere data as a scatter plot
plt.sca(axes[1])
m = Basemap(projection='spstere', boundinglat=-45,lon_0=270, resolution='l')
m.drawcoastlines()
m.fillcontinents()
m.drawcountries()

x, y = m(lon, lat)
sc = m.scatter(x, y, c=data, cmap='jet', lw=0, s=4)

# Create the colorbar and add Pass / Fail labels
divider = make_axes_locatable(ax)
cax = divider.append_axes("bottom", size="5%", pad=0.5)
cb = plt.colorbar(sc, cax=cax, orientation="horizontal", format="%.2f")
m.drawparallels(np.arange(-90.,120.,15.),labels=[1,0,0,0]) # draw parallels
m.drawmeridians(np.arange(0.,420.,30.),labels=[1,1,1,1]) # draw meridians

plt.suptitle('CICE Mean Ice Thickness\n{}'.format(case), y=0.95)

# Make some room at the bottom of the figure, and create a colorbar
fig.subplots_adjust(bottom=0.2)
cbar_ax = fig.add_axes([0.11,0.1,0.8,0.05])
cb = plt.colorbar(sc, cax=cbar_ax, orientation="horizontal", format="%.2f")
cb.set_label(units, x=1.0)

plt.savefig('ice_thickness_map.png', dpi=300)
outfile = 'ice_thickness_{}.png'.format(case.replace('\n- ','_minus_'))
logger.info('Creating map of the data ({})'.format(outfile))
plt.savefig(outfile, dpi=300, bbox_inches='tight')

def plot_two_stage_failures(data, lat, lon):
'''This function plots each grid cell and whether or not it Passed or Failed
Expand Down Expand Up @@ -428,7 +443,7 @@ def plot_two_stage_failures(data, lat, lon):
ax = fig.add_axes([0.05, 0.08, 0.9, 0.9])

# Create the basemap, and draw boundaries
m = Basemap(projection='kav7', lon_0=180., resolution='l')
m = Basemap(projection='moll', lon_0=0., resolution='l')
m.drawmapboundary(fill_color='white')
m.drawcoastlines()
m.drawcountries()
Expand All @@ -440,7 +455,7 @@ def plot_two_stage_failures(data, lat, lon):

# Plot the data as a scatter plot
x, y = m(lon, lat)
sc = m.scatter(x, y, c=int_data, cmap=cm, lw=0, vmin=0, vmax=1)
sc = m.scatter(x, y, c=int_data, cmap=cm, lw=0, vmin=0, vmax=1, s=4)

m.drawmeridians(np.arange(0, 360, 60), labels=[0, 0, 0, 1], fontsize=10)
m.drawparallels(np.arange(-90, 90, 30), labels=[1, 0, 0, 0], fontsize=10)
Expand Down Expand Up @@ -482,8 +497,11 @@ def main():
help='Path to the test history (iceh_inst*) files. REQUIRED')
parser.add_argument('-v', '--verbose', dest='verbose', help='Print debug output?', \
action='store_true')
parser.add_argument('--plot', dest='plot', help='Create maps of ice thickness?', \
action='store_true')

parser.set_defaults(verbose=False)
parser.set_defaults(plot=False)

# If no arguments are provided, print the help message
if len(sys.argv) == 1:
Expand Down Expand Up @@ -528,6 +546,15 @@ def main():
# If test failed, attempt to create a plot of the failure locations
if not PASSED:
plot_two_stage_failures(H1_array, t_lat, t_lon)
if args.plot:
baseDir = os.path.abspath(args.base_dir).rstrip('history/').rstrip(\
'history').split('/')[-1]
testDir = os.path.abspath(args.test_dir).rstrip('history/').rstrip( \
'history').split('/')[-1]
plot_data(np.mean(data_base,axis=0), t_lat, t_lon, 'm', baseDir)
plot_data(np.mean(data_test,axis=0), t_lat, t_lon, 'm', testDir)
plot_data(np.mean(data_base-data_test,axis=0), t_lat, t_lon, 'm', '{}\n- {}'.\
format(baseDir,testDir))
logger.error('Quality Control Test FAILED')
sys.exit(-1)

Expand Down Expand Up @@ -559,6 +586,17 @@ def main():

PASSED_SKILL = PASSED_NH and PASSED_SH

# Plot the ice thickness data for the base and test cases
if args.plot:
baseDir = os.path.abspath(args.base_dir).rstrip('history/').rstrip( \
'history').split('/')[-1]
testDir = os.path.abspath(args.test_dir).rstrip('history/').rstrip( \
'history').split('/')[-1]
plot_data(np.mean(data_base,axis=0), t_lat, t_lon, 'm', baseDir)
plot_data(np.mean(data_test,axis=0), t_lat, t_lon, 'm', testDir)
plot_data(np.mean(data_base-data_test,axis=0), t_lat, t_lon, 'm', '{}\n- {}'.\
format(baseDir,testDir))

logger.info('')
if not PASSED_SKILL:
logger.error('Quality Control Test FAILED')
Expand Down