plot_north_american_categories.py 4.64 KB
Newer Older
1
"""
2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
Created on Nov 29, 2012

@author: abuddenberg

The purpose of this program is to plot the categories of statistical
significance for the various emissions scenarios in the North American
area.  This plot is used in conjection with the accompanying contour plots
to spot-check anomalies caused by the smoothing in those plots (as well as to
verify the data itself).

The categories are as follows:

1.0 denotes areas little change
2.0 denotes areas of statistical uncertainty
3.0 denotes areas of statistical significance

Fun stuff
19
"""
20 21 22 23 24 25 26
from scipy.io.netcdf import netcdf_file
from mpl_toolkits.basemap import Basemap
from numpy import meshgrid
from numpy.ma import masked_equal
import numpy as np

import matplotlib.pyplot as plt
27
from config import GLOBAL_PRECIP_FILES, NA_PRECIP_FILES, SEASONS
28

29
for infilename, outfilename in NA_PRECIP_FILES:
30 31 32
    nc = netcdf_file(infilename) 
    
    lat_data = nc.variables['lat'].data 
33
    lon_data = nc.variables['lon'].data
abuddenberg's avatar
abuddenberg committed
34
<<<<<<< HEAD
abuddenberg's avatar
Merging  
abuddenberg committed
35

36 37
    fig = plt.figure(figsize=(25,16), dpi=100, tight_layout=True)
    
abuddenberg's avatar
Merging  
abuddenberg committed
38
    for i, season in enumerate(['Winter', 'Spring', 'Summer', 'Fall']): #['Winter', 'Spring', 'Summer', 'Fall']
abuddenberg's avatar
abuddenberg committed
39
=======
40 41 42
    
    fig = plt.figure(figsize=(25,16), dpi=100, tight_layout=True)
    
43
    for i, season in enumerate(['Winter', 'Spring', 'Summer', 'Fall']):
abuddenberg's avatar
abuddenberg committed
44
>>>>>>> dbd7e5f7b2a0c4d96dda5d2c140efbb6a052e6ad
45 46 47 48
        data_var, signif_var = SEASONS[season]
        
        data = nc.variables[data_var].data
        signif = nc.variables[signif_var].data
abuddenberg's avatar
abuddenberg committed
49
<<<<<<< HEAD
abuddenberg's avatar
Merging  
abuddenberg committed
50

abuddenberg's avatar
abuddenberg committed
51
=======
52
        
abuddenberg's avatar
abuddenberg committed
53
>>>>>>> dbd7e5f7b2a0c4d96dda5d2c140efbb6a052e6ad
54
        ax = fig.add_subplot(231 + i)
55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72
        plt.title(season)
        m = Basemap(
            projection='aea',
                    lon_0=-96,
                    lat_0=37.5,
                    lat_1=29.5,
                    lat_2=45.5,
        #            lat_ts=median(lats),
            llcrnrlat=12,
            urcrnrlat=80,
            llcrnrlon=-135,
            urcrnrlon=-25,
            resolution='l',area_thresh=10000
        )
        
        m.drawcoastlines()
        m.drawstates()
        m.drawcountries()
abuddenberg's avatar
abuddenberg committed
73
<<<<<<< HEAD
abuddenberg's avatar
Merging  
abuddenberg committed
74

75
        lons, lats = meshgrid(lon_data, lat_data)
abuddenberg's avatar
Merging  
abuddenberg committed
76 77 78 79 80 81 82 83 84 85 86 87

        #pcolor coordinates should refer to the lower left corner of the grid box;
        # the data refers to the center of the grid box.  Shift the lon and lat south and west by half the height and
        # width of the grid box to obtain the correct representation.
        grid_width = lons[0][1] - lons[0][0]
        grid_height = lats[1][0] - lats[0][0]

        lons_shifted = (lons - 0.5 * grid_width)
        lats_shifted = (lats - 0.5 * grid_height)

        x, y = m(lons, lats)
        x_shifted, y_shifted = m(lons_shifted, lats_shifted)
abuddenberg's avatar
abuddenberg committed
88
=======
89 90 91
        
        lons, lats = meshgrid(lon_data, lat_data)
        x,y = m(lons, lats)
abuddenberg's avatar
abuddenberg committed
92
>>>>>>> dbd7e5f7b2a0c4d96dda5d2c140efbb6a052e6ad
93 94 95 96 97 98 99 100 101 102 103
    
        #Build boolean masks of the gridpoint for each category
        stipples_mask = np.ma.getmask(np.ma.masked_equal(signif, 1.))
        zeros_mask = np.ma.getmask(masked_equal(data, 0.0))

        both_mask = np.ma.mask_or(stipples_mask, zeros_mask)
        third_cat_mask = ~both_mask
        
        #There's got to be a better way of doing this than copying the array
        data = np.ma.masked_array(data)
        data.mask = stipples_mask
104
        data = np.ma.masked_array(data.filled(3.0)) #3.0 denotes areas of statistical significance; red
105 106
        
        data.mask = zeros_mask
107
        data = np.ma.masked_array(data.filled(1.0)) #1.0 denotes areas little change; blue
abuddenberg's avatar
abuddenberg committed
108
<<<<<<< HEAD
abuddenberg's avatar
Merging  
abuddenberg committed
109

110
        data.mask = third_cat_mask
111
        data = np.ma.masked_array(data.filled(2.0)) #2.0 denotes areas of statistical uncertainty; green
abuddenberg's avatar
Merging  
abuddenberg committed
112 113 114 115 116 117

        weird = m.pcolor(x_shifted, y_shifted, data)
        m.colorbar(weird, location='right',  pad="5%")

        m.scatter(x, y, 3, marker='o')

abuddenberg's avatar
abuddenberg committed
118
=======
119 120
                
        data.mask = third_cat_mask
121
        data = np.ma.masked_array(data.filled(2.0)) #2.0 denotes areas of statistical uncertainty; green
122 123 124 125 126
        
 
        weird = m.pcolor(x,y, data)
        m.colorbar(weird,location='right',pad="5%")
        
abuddenberg's avatar
abuddenberg committed
127
>>>>>>> dbd7e5f7b2a0c4d96dda5d2c140efbb6a052e6ad
128 129 130 131 132 133
        #Tests for overlap (There shouldn't be any)
#        print np.any(np.logical_and(stipples_mask, zeros_mask))
#        print np.any(np.logical_and(third_cat_mask, zeros_mask))
#        print np.any(np.logical_and(third_cat_mask, stipples_mask))
        
    
abuddenberg's avatar
abuddenberg committed
134
<<<<<<< HEAD
abuddenberg's avatar
Merging  
abuddenberg committed
135
    # plt.savefig('../dist/' + outfilename.format('categories'), format='eps', dpi=200)
abuddenberg's avatar
abuddenberg committed
136
=======
137
    # plt.savefig('../dist/' + outfilename.format('north_american_categories'), format='eps', dpi=200)
abuddenberg's avatar
abuddenberg committed
138
>>>>>>> dbd7e5f7b2a0c4d96dda5d2c140efbb6a052e6ad
139
    plt.show()
140