How to plot contours from a polar stereographic grib2 file in Python
I am trying to plot a CMC grib2 pressure forecast file using matplotlib to plot the pressure contours. The description of the grib2 grid can be found here: https://weather.gc.ca/grib/grib2_reg_10km_e.html. The grib2 file is found in this directory: http://dd.weather.gc.ca/model_gem_regional/10km/grib2/00/000/ and starts with CMC_reg_PRMSL_MSL_0_ps10km followed by the date. It is a grib file containing pressure at mean sea level.
My problem is that I end up having some straight line contours that follow the lines of latitude on top of the actual pressure contours. I thought it might be because I am plotting in PlateCarree as opposed to Geodetic but the contour plot will not allow using Geodetic. The result of my plot is:
Code is as follows:
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import datetime as dt
import cartopy
import cartopy.crs as ccrs
import Nio
gr = Nio.open_file('./data/CMC_reg_PRMSL_MSL_0_ps10km_2018111800_P000.grib2', 'r')
print(gr)
names = gr.variables.keys()
print("Variable Names:", names)
dims = gr.dimensions
print("Dimensions: ", dims)
attr = gr.attributes.keys()
print("Attributes: ", attr)
obs = gr.variables['PRMSL_P0_L101_GST0'][:]
lats = gr.variables["gridlat_0"][:]
lons = gr.variables["gridlon_0"][:]
fig = plt.figure(figsize=(15, 2))
intervals = range(95000, 105000, 400)
ax=plt.axes([0.,0.,1.,1.],projection=ccrs.PlateCarree())
obsobj = plt.contour(lons, lats, obs, intervals, cmap='jet',transform=ccrs.PlateCarree())
states_provinces = cartopy.feature.NaturalEarthFeature(
category='cultural',
name='admin_1_states_provinces_lines',
scale='50m',
facecolor='none')
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
ax.add_feature(states_provinces,edgecolor='gray')
obsobj.clabel()
colbar =plt.colorbar(obsobj)
Any suggestions would be appreciated.
UPDATE
For anyone without PyNIO the following can be used to reproduce using the dump files in the comments section.
Just remove all the references to NIO and replace the lats, lons, obs assignment with the following.
lats = np.load('lats.dump')
lons = np.load('lons.dump')
obs = np.load('obs.dump')
python matplotlib cartopy grib
add a comment |
I am trying to plot a CMC grib2 pressure forecast file using matplotlib to plot the pressure contours. The description of the grib2 grid can be found here: https://weather.gc.ca/grib/grib2_reg_10km_e.html. The grib2 file is found in this directory: http://dd.weather.gc.ca/model_gem_regional/10km/grib2/00/000/ and starts with CMC_reg_PRMSL_MSL_0_ps10km followed by the date. It is a grib file containing pressure at mean sea level.
My problem is that I end up having some straight line contours that follow the lines of latitude on top of the actual pressure contours. I thought it might be because I am plotting in PlateCarree as opposed to Geodetic but the contour plot will not allow using Geodetic. The result of my plot is:
Code is as follows:
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import datetime as dt
import cartopy
import cartopy.crs as ccrs
import Nio
gr = Nio.open_file('./data/CMC_reg_PRMSL_MSL_0_ps10km_2018111800_P000.grib2', 'r')
print(gr)
names = gr.variables.keys()
print("Variable Names:", names)
dims = gr.dimensions
print("Dimensions: ", dims)
attr = gr.attributes.keys()
print("Attributes: ", attr)
obs = gr.variables['PRMSL_P0_L101_GST0'][:]
lats = gr.variables["gridlat_0"][:]
lons = gr.variables["gridlon_0"][:]
fig = plt.figure(figsize=(15, 2))
intervals = range(95000, 105000, 400)
ax=plt.axes([0.,0.,1.,1.],projection=ccrs.PlateCarree())
obsobj = plt.contour(lons, lats, obs, intervals, cmap='jet',transform=ccrs.PlateCarree())
states_provinces = cartopy.feature.NaturalEarthFeature(
category='cultural',
name='admin_1_states_provinces_lines',
scale='50m',
facecolor='none')
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
ax.add_feature(states_provinces,edgecolor='gray')
obsobj.clabel()
colbar =plt.colorbar(obsobj)
Any suggestions would be appreciated.
UPDATE
For anyone without PyNIO the following can be used to reproduce using the dump files in the comments section.
Just remove all the references to NIO and replace the lats, lons, obs assignment with the following.
lats = np.load('lats.dump')
lons = np.load('lons.dump')
obs = np.load('obs.dump')
python matplotlib cartopy grib
Unfortunately, I don't think so (though I'm open to suggestions on how this could be achieved). The varaibles are 824x935 numpy arrays so it is a little difficult to put in the code directly. Since I am not sure why it is behaving the way it is, I don't know how I could create a smaller example to reproduce the problem.
– Campbell
Nov 18 '18 at 18:10
That works. Here are the links to the files (for some reason I was getting errors including them in my question update: drive.google.com/open?id=1A7mZtqFcPCk92WRJ11TxOnHrsniZH_s- drive.google.com/open?id=15RC-IsoojBLIG9VMknTpFP1rPqV2qDN3 drive.google.com/open?id=1An3oXDzev-4hGmapvYykduG293fhC1PH
– Campbell
Nov 18 '18 at 18:46
The problem is the grid, which is completely deformed in cartesian coordinates and winds about the north pole. I plotted the grid index in both dimensions here to visualize the problem. In principle this problem seems similar to this one. However, the one answer does not seem to solve the problem; the other answer (written by me) cannot be applied here due to the gridshape.
– ImportanceOfBeingErnest
Nov 19 '18 at 0:43
I'm still very new to numpy and gridded data... do you think it might be possible to discard the upper 15-20° of latitude? I could get by with just data up to 80°N. I'm just not sure how to find where to make the cut and what the impact would be for the contours.
– Campbell
Nov 19 '18 at 1:12
add a comment |
I am trying to plot a CMC grib2 pressure forecast file using matplotlib to plot the pressure contours. The description of the grib2 grid can be found here: https://weather.gc.ca/grib/grib2_reg_10km_e.html. The grib2 file is found in this directory: http://dd.weather.gc.ca/model_gem_regional/10km/grib2/00/000/ and starts with CMC_reg_PRMSL_MSL_0_ps10km followed by the date. It is a grib file containing pressure at mean sea level.
My problem is that I end up having some straight line contours that follow the lines of latitude on top of the actual pressure contours. I thought it might be because I am plotting in PlateCarree as opposed to Geodetic but the contour plot will not allow using Geodetic. The result of my plot is:
Code is as follows:
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import datetime as dt
import cartopy
import cartopy.crs as ccrs
import Nio
gr = Nio.open_file('./data/CMC_reg_PRMSL_MSL_0_ps10km_2018111800_P000.grib2', 'r')
print(gr)
names = gr.variables.keys()
print("Variable Names:", names)
dims = gr.dimensions
print("Dimensions: ", dims)
attr = gr.attributes.keys()
print("Attributes: ", attr)
obs = gr.variables['PRMSL_P0_L101_GST0'][:]
lats = gr.variables["gridlat_0"][:]
lons = gr.variables["gridlon_0"][:]
fig = plt.figure(figsize=(15, 2))
intervals = range(95000, 105000, 400)
ax=plt.axes([0.,0.,1.,1.],projection=ccrs.PlateCarree())
obsobj = plt.contour(lons, lats, obs, intervals, cmap='jet',transform=ccrs.PlateCarree())
states_provinces = cartopy.feature.NaturalEarthFeature(
category='cultural',
name='admin_1_states_provinces_lines',
scale='50m',
facecolor='none')
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
ax.add_feature(states_provinces,edgecolor='gray')
obsobj.clabel()
colbar =plt.colorbar(obsobj)
Any suggestions would be appreciated.
UPDATE
For anyone without PyNIO the following can be used to reproduce using the dump files in the comments section.
Just remove all the references to NIO and replace the lats, lons, obs assignment with the following.
lats = np.load('lats.dump')
lons = np.load('lons.dump')
obs = np.load('obs.dump')
python matplotlib cartopy grib
I am trying to plot a CMC grib2 pressure forecast file using matplotlib to plot the pressure contours. The description of the grib2 grid can be found here: https://weather.gc.ca/grib/grib2_reg_10km_e.html. The grib2 file is found in this directory: http://dd.weather.gc.ca/model_gem_regional/10km/grib2/00/000/ and starts with CMC_reg_PRMSL_MSL_0_ps10km followed by the date. It is a grib file containing pressure at mean sea level.
My problem is that I end up having some straight line contours that follow the lines of latitude on top of the actual pressure contours. I thought it might be because I am plotting in PlateCarree as opposed to Geodetic but the contour plot will not allow using Geodetic. The result of my plot is:
Code is as follows:
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import datetime as dt
import cartopy
import cartopy.crs as ccrs
import Nio
gr = Nio.open_file('./data/CMC_reg_PRMSL_MSL_0_ps10km_2018111800_P000.grib2', 'r')
print(gr)
names = gr.variables.keys()
print("Variable Names:", names)
dims = gr.dimensions
print("Dimensions: ", dims)
attr = gr.attributes.keys()
print("Attributes: ", attr)
obs = gr.variables['PRMSL_P0_L101_GST0'][:]
lats = gr.variables["gridlat_0"][:]
lons = gr.variables["gridlon_0"][:]
fig = plt.figure(figsize=(15, 2))
intervals = range(95000, 105000, 400)
ax=plt.axes([0.,0.,1.,1.],projection=ccrs.PlateCarree())
obsobj = plt.contour(lons, lats, obs, intervals, cmap='jet',transform=ccrs.PlateCarree())
states_provinces = cartopy.feature.NaturalEarthFeature(
category='cultural',
name='admin_1_states_provinces_lines',
scale='50m',
facecolor='none')
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
ax.add_feature(states_provinces,edgecolor='gray')
obsobj.clabel()
colbar =plt.colorbar(obsobj)
Any suggestions would be appreciated.
UPDATE
For anyone without PyNIO the following can be used to reproduce using the dump files in the comments section.
Just remove all the references to NIO and replace the lats, lons, obs assignment with the following.
lats = np.load('lats.dump')
lons = np.load('lons.dump')
obs = np.load('obs.dump')
python matplotlib cartopy grib
python matplotlib cartopy grib
edited Nov 18 '18 at 19:45
ImportanceOfBeingErnest
128k12137214
128k12137214
asked Nov 18 '18 at 14:30
CampbellCampbell
886
886
Unfortunately, I don't think so (though I'm open to suggestions on how this could be achieved). The varaibles are 824x935 numpy arrays so it is a little difficult to put in the code directly. Since I am not sure why it is behaving the way it is, I don't know how I could create a smaller example to reproduce the problem.
– Campbell
Nov 18 '18 at 18:10
That works. Here are the links to the files (for some reason I was getting errors including them in my question update: drive.google.com/open?id=1A7mZtqFcPCk92WRJ11TxOnHrsniZH_s- drive.google.com/open?id=15RC-IsoojBLIG9VMknTpFP1rPqV2qDN3 drive.google.com/open?id=1An3oXDzev-4hGmapvYykduG293fhC1PH
– Campbell
Nov 18 '18 at 18:46
The problem is the grid, which is completely deformed in cartesian coordinates and winds about the north pole. I plotted the grid index in both dimensions here to visualize the problem. In principle this problem seems similar to this one. However, the one answer does not seem to solve the problem; the other answer (written by me) cannot be applied here due to the gridshape.
– ImportanceOfBeingErnest
Nov 19 '18 at 0:43
I'm still very new to numpy and gridded data... do you think it might be possible to discard the upper 15-20° of latitude? I could get by with just data up to 80°N. I'm just not sure how to find where to make the cut and what the impact would be for the contours.
– Campbell
Nov 19 '18 at 1:12
add a comment |
Unfortunately, I don't think so (though I'm open to suggestions on how this could be achieved). The varaibles are 824x935 numpy arrays so it is a little difficult to put in the code directly. Since I am not sure why it is behaving the way it is, I don't know how I could create a smaller example to reproduce the problem.
– Campbell
Nov 18 '18 at 18:10
That works. Here are the links to the files (for some reason I was getting errors including them in my question update: drive.google.com/open?id=1A7mZtqFcPCk92WRJ11TxOnHrsniZH_s- drive.google.com/open?id=15RC-IsoojBLIG9VMknTpFP1rPqV2qDN3 drive.google.com/open?id=1An3oXDzev-4hGmapvYykduG293fhC1PH
– Campbell
Nov 18 '18 at 18:46
The problem is the grid, which is completely deformed in cartesian coordinates and winds about the north pole. I plotted the grid index in both dimensions here to visualize the problem. In principle this problem seems similar to this one. However, the one answer does not seem to solve the problem; the other answer (written by me) cannot be applied here due to the gridshape.
– ImportanceOfBeingErnest
Nov 19 '18 at 0:43
I'm still very new to numpy and gridded data... do you think it might be possible to discard the upper 15-20° of latitude? I could get by with just data up to 80°N. I'm just not sure how to find where to make the cut and what the impact would be for the contours.
– Campbell
Nov 19 '18 at 1:12
Unfortunately, I don't think so (though I'm open to suggestions on how this could be achieved). The varaibles are 824x935 numpy arrays so it is a little difficult to put in the code directly. Since I am not sure why it is behaving the way it is, I don't know how I could create a smaller example to reproduce the problem.
– Campbell
Nov 18 '18 at 18:10
Unfortunately, I don't think so (though I'm open to suggestions on how this could be achieved). The varaibles are 824x935 numpy arrays so it is a little difficult to put in the code directly. Since I am not sure why it is behaving the way it is, I don't know how I could create a smaller example to reproduce the problem.
– Campbell
Nov 18 '18 at 18:10
That works. Here are the links to the files (for some reason I was getting errors including them in my question update: drive.google.com/open?id=1A7mZtqFcPCk92WRJ11TxOnHrsniZH_s- drive.google.com/open?id=15RC-IsoojBLIG9VMknTpFP1rPqV2qDN3 drive.google.com/open?id=1An3oXDzev-4hGmapvYykduG293fhC1PH
– Campbell
Nov 18 '18 at 18:46
That works. Here are the links to the files (for some reason I was getting errors including them in my question update: drive.google.com/open?id=1A7mZtqFcPCk92WRJ11TxOnHrsniZH_s- drive.google.com/open?id=15RC-IsoojBLIG9VMknTpFP1rPqV2qDN3 drive.google.com/open?id=1An3oXDzev-4hGmapvYykduG293fhC1PH
– Campbell
Nov 18 '18 at 18:46
The problem is the grid, which is completely deformed in cartesian coordinates and winds about the north pole. I plotted the grid index in both dimensions here to visualize the problem. In principle this problem seems similar to this one. However, the one answer does not seem to solve the problem; the other answer (written by me) cannot be applied here due to the gridshape.
– ImportanceOfBeingErnest
Nov 19 '18 at 0:43
The problem is the grid, which is completely deformed in cartesian coordinates and winds about the north pole. I plotted the grid index in both dimensions here to visualize the problem. In principle this problem seems similar to this one. However, the one answer does not seem to solve the problem; the other answer (written by me) cannot be applied here due to the gridshape.
– ImportanceOfBeingErnest
Nov 19 '18 at 0:43
I'm still very new to numpy and gridded data... do you think it might be possible to discard the upper 15-20° of latitude? I could get by with just data up to 80°N. I'm just not sure how to find where to make the cut and what the impact would be for the contours.
– Campbell
Nov 19 '18 at 1:12
I'm still very new to numpy and gridded data... do you think it might be possible to discard the upper 15-20° of latitude? I could get by with just data up to 80°N. I'm just not sure how to find where to make the cut and what the impact would be for the contours.
– Campbell
Nov 19 '18 at 1:12
add a comment |
1 Answer
1
active
oldest
votes
The problem
The problem is that the grid winds around the earth. Hence there will be points on the grid at -180° whose nearst neighbor sits at +180°, i.e. the grid wraps around the antimeridian. The following plots the grid index along both directions. One can see that the first grid row (black) appears on both sides of the plot.
Hence a contour line following the pacific westwards needs to then cross straight through the plot to continue towards japan on the other side of the plot. This will lead to the undesired lines
A solution
A solution is to mask the outer points of the PlateCarree out. Those occur in the middle of the grid. Cutting the grid at coordinates of longitude larger than 179° or smaller than -179°, as well as leaving the north pole out would look like
where the blue denotes the cut out points.
Applying this to the contour plot gives:
import matplotlib.pyplot as plt
import numpy as np
import cartopy
import cartopy.crs as ccrs
lats = np.load('data/lats.dump')
lons = np.load('data/lons.dump')
obs = np.load('data/obs.dump')
intervals = range(95000, 105000, 400)
fig, ax = plt.subplots(figsize=(15,4), subplot_kw=dict(projection=ccrs.PlateCarree()))
fig.subplots_adjust(left=0.03, right=0.97, top=0.8, bottom=0.2)
mask = (lons > 179) | (lons < -179) | (lats > 89)
maskedobs = np.ma.array(obs, mask=mask)
pc = ax.contour(lons, lats, maskedobs, intervals, cmap='jet', transform=ccrs.PlateCarree())
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
colbar =plt.colorbar(pc)
plt.show()
add a comment |
Your Answer
StackExchange.ifUsing("editor", function () {
StackExchange.using("externalEditor", function () {
StackExchange.using("snippets", function () {
StackExchange.snippets.init();
});
});
}, "code-snippets");
StackExchange.ready(function() {
var channelOptions = {
tags: "".split(" "),
id: "1"
};
initTagRenderer("".split(" "), "".split(" "), channelOptions);
StackExchange.using("externalEditor", function() {
// Have to fire editor after snippets, if snippets enabled
if (StackExchange.settings.snippets.snippetsEnabled) {
StackExchange.using("snippets", function() {
createEditor();
});
}
else {
createEditor();
}
});
function createEditor() {
StackExchange.prepareEditor({
heartbeatType: 'answer',
autoActivateHeartbeat: false,
convertImagesToLinks: true,
noModals: true,
showLowRepImageUploadWarning: true,
reputationToPostImages: 10,
bindNavPrevention: true,
postfix: "",
imageUploader: {
brandingHtml: "Powered by u003ca class="icon-imgur-white" href="https://imgur.com/"u003eu003c/au003e",
contentPolicyHtml: "User contributions licensed under u003ca href="https://creativecommons.org/licenses/by-sa/3.0/"u003ecc by-sa 3.0 with attribution requiredu003c/au003e u003ca href="https://stackoverflow.com/legal/content-policy"u003e(content policy)u003c/au003e",
allowUrls: true
},
onDemand: true,
discardSelector: ".discard-answer"
,immediatelyShowMarkdownHelp:true
});
}
});
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
StackExchange.ready(
function () {
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fstackoverflow.com%2fquestions%2f53361962%2fhow-to-plot-contours-from-a-polar-stereographic-grib2-file-in-python%23new-answer', 'question_page');
}
);
Post as a guest
Required, but never shown
1 Answer
1
active
oldest
votes
1 Answer
1
active
oldest
votes
active
oldest
votes
active
oldest
votes
The problem
The problem is that the grid winds around the earth. Hence there will be points on the grid at -180° whose nearst neighbor sits at +180°, i.e. the grid wraps around the antimeridian. The following plots the grid index along both directions. One can see that the first grid row (black) appears on both sides of the plot.
Hence a contour line following the pacific westwards needs to then cross straight through the plot to continue towards japan on the other side of the plot. This will lead to the undesired lines
A solution
A solution is to mask the outer points of the PlateCarree out. Those occur in the middle of the grid. Cutting the grid at coordinates of longitude larger than 179° or smaller than -179°, as well as leaving the north pole out would look like
where the blue denotes the cut out points.
Applying this to the contour plot gives:
import matplotlib.pyplot as plt
import numpy as np
import cartopy
import cartopy.crs as ccrs
lats = np.load('data/lats.dump')
lons = np.load('data/lons.dump')
obs = np.load('data/obs.dump')
intervals = range(95000, 105000, 400)
fig, ax = plt.subplots(figsize=(15,4), subplot_kw=dict(projection=ccrs.PlateCarree()))
fig.subplots_adjust(left=0.03, right=0.97, top=0.8, bottom=0.2)
mask = (lons > 179) | (lons < -179) | (lats > 89)
maskedobs = np.ma.array(obs, mask=mask)
pc = ax.contour(lons, lats, maskedobs, intervals, cmap='jet', transform=ccrs.PlateCarree())
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
colbar =plt.colorbar(pc)
plt.show()
add a comment |
The problem
The problem is that the grid winds around the earth. Hence there will be points on the grid at -180° whose nearst neighbor sits at +180°, i.e. the grid wraps around the antimeridian. The following plots the grid index along both directions. One can see that the first grid row (black) appears on both sides of the plot.
Hence a contour line following the pacific westwards needs to then cross straight through the plot to continue towards japan on the other side of the plot. This will lead to the undesired lines
A solution
A solution is to mask the outer points of the PlateCarree out. Those occur in the middle of the grid. Cutting the grid at coordinates of longitude larger than 179° or smaller than -179°, as well as leaving the north pole out would look like
where the blue denotes the cut out points.
Applying this to the contour plot gives:
import matplotlib.pyplot as plt
import numpy as np
import cartopy
import cartopy.crs as ccrs
lats = np.load('data/lats.dump')
lons = np.load('data/lons.dump')
obs = np.load('data/obs.dump')
intervals = range(95000, 105000, 400)
fig, ax = plt.subplots(figsize=(15,4), subplot_kw=dict(projection=ccrs.PlateCarree()))
fig.subplots_adjust(left=0.03, right=0.97, top=0.8, bottom=0.2)
mask = (lons > 179) | (lons < -179) | (lats > 89)
maskedobs = np.ma.array(obs, mask=mask)
pc = ax.contour(lons, lats, maskedobs, intervals, cmap='jet', transform=ccrs.PlateCarree())
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
colbar =plt.colorbar(pc)
plt.show()
add a comment |
The problem
The problem is that the grid winds around the earth. Hence there will be points on the grid at -180° whose nearst neighbor sits at +180°, i.e. the grid wraps around the antimeridian. The following plots the grid index along both directions. One can see that the first grid row (black) appears on both sides of the plot.
Hence a contour line following the pacific westwards needs to then cross straight through the plot to continue towards japan on the other side of the plot. This will lead to the undesired lines
A solution
A solution is to mask the outer points of the PlateCarree out. Those occur in the middle of the grid. Cutting the grid at coordinates of longitude larger than 179° or smaller than -179°, as well as leaving the north pole out would look like
where the blue denotes the cut out points.
Applying this to the contour plot gives:
import matplotlib.pyplot as plt
import numpy as np
import cartopy
import cartopy.crs as ccrs
lats = np.load('data/lats.dump')
lons = np.load('data/lons.dump')
obs = np.load('data/obs.dump')
intervals = range(95000, 105000, 400)
fig, ax = plt.subplots(figsize=(15,4), subplot_kw=dict(projection=ccrs.PlateCarree()))
fig.subplots_adjust(left=0.03, right=0.97, top=0.8, bottom=0.2)
mask = (lons > 179) | (lons < -179) | (lats > 89)
maskedobs = np.ma.array(obs, mask=mask)
pc = ax.contour(lons, lats, maskedobs, intervals, cmap='jet', transform=ccrs.PlateCarree())
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
colbar =plt.colorbar(pc)
plt.show()
The problem
The problem is that the grid winds around the earth. Hence there will be points on the grid at -180° whose nearst neighbor sits at +180°, i.e. the grid wraps around the antimeridian. The following plots the grid index along both directions. One can see that the first grid row (black) appears on both sides of the plot.
Hence a contour line following the pacific westwards needs to then cross straight through the plot to continue towards japan on the other side of the plot. This will lead to the undesired lines
A solution
A solution is to mask the outer points of the PlateCarree out. Those occur in the middle of the grid. Cutting the grid at coordinates of longitude larger than 179° or smaller than -179°, as well as leaving the north pole out would look like
where the blue denotes the cut out points.
Applying this to the contour plot gives:
import matplotlib.pyplot as plt
import numpy as np
import cartopy
import cartopy.crs as ccrs
lats = np.load('data/lats.dump')
lons = np.load('data/lons.dump')
obs = np.load('data/obs.dump')
intervals = range(95000, 105000, 400)
fig, ax = plt.subplots(figsize=(15,4), subplot_kw=dict(projection=ccrs.PlateCarree()))
fig.subplots_adjust(left=0.03, right=0.97, top=0.8, bottom=0.2)
mask = (lons > 179) | (lons < -179) | (lats > 89)
maskedobs = np.ma.array(obs, mask=mask)
pc = ax.contour(lons, lats, maskedobs, intervals, cmap='jet', transform=ccrs.PlateCarree())
ax.add_feature(cartopy.feature.BORDERS)
ax.coastlines(resolution='10m')
colbar =plt.colorbar(pc)
plt.show()
answered Nov 19 '18 at 3:12
ImportanceOfBeingErnestImportanceOfBeingErnest
128k12137214
128k12137214
add a comment |
add a comment |
Thanks for contributing an answer to Stack Overflow!
- Please be sure to answer the question. Provide details and share your research!
But avoid …
- Asking for help, clarification, or responding to other answers.
- Making statements based on opinion; back them up with references or personal experience.
To learn more, see our tips on writing great answers.
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
StackExchange.ready(
function () {
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fstackoverflow.com%2fquestions%2f53361962%2fhow-to-plot-contours-from-a-polar-stereographic-grib2-file-in-python%23new-answer', 'question_page');
}
);
Post as a guest
Required, but never shown
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
Sign up or log in
StackExchange.ready(function () {
StackExchange.helpers.onClickDraftSave('#login-link');
});
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Sign up using Google
Sign up using Facebook
Sign up using Email and Password
Post as a guest
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Required, but never shown
Unfortunately, I don't think so (though I'm open to suggestions on how this could be achieved). The varaibles are 824x935 numpy arrays so it is a little difficult to put in the code directly. Since I am not sure why it is behaving the way it is, I don't know how I could create a smaller example to reproduce the problem.
– Campbell
Nov 18 '18 at 18:10
That works. Here are the links to the files (for some reason I was getting errors including them in my question update: drive.google.com/open?id=1A7mZtqFcPCk92WRJ11TxOnHrsniZH_s- drive.google.com/open?id=15RC-IsoojBLIG9VMknTpFP1rPqV2qDN3 drive.google.com/open?id=1An3oXDzev-4hGmapvYykduG293fhC1PH
– Campbell
Nov 18 '18 at 18:46
The problem is the grid, which is completely deformed in cartesian coordinates and winds about the north pole. I plotted the grid index in both dimensions here to visualize the problem. In principle this problem seems similar to this one. However, the one answer does not seem to solve the problem; the other answer (written by me) cannot be applied here due to the gridshape.
– ImportanceOfBeingErnest
Nov 19 '18 at 0:43
I'm still very new to numpy and gridded data... do you think it might be possible to discard the upper 15-20° of latitude? I could get by with just data up to 80°N. I'm just not sure how to find where to make the cut and what the impact would be for the contours.
– Campbell
Nov 19 '18 at 1:12