Изменение проекций карты при использовании базовой карты

Я изменил проекции карты в базовой карте раньше, так что я знаю, что это должно быть легко исправить, но мне кажется, что ничего не работает. Я использовал сетку сетки, и отображал свои значения x,y и т. Д., И я просто получаю искаженные или сумасшедшие графики. Я думаю, что это как-то связано с тем фактом, что данные, которые я использую, автоматически устанавливаются для построения на конформной системе Ламберта (что мне не нужно), также размеры указаны в км, а не в широтах и ​​долях. Я не уверен, куда идти отсюда...

Источник данных: http://thredds.ucar.edu/thredds/dodsC/grib/NCEP/RAP/CONUS_13km/RR_CONUS_13km_20150415_0600.grib2/GC.html

Вот мой рабочий код. У меня есть куча комментариев, которые я безуспешно пробовал.

import numpy as np
import math as m
import urllib2
import time
import datetime as dt
import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1 import make_axes_locatable
from mpl_toolkits.basemap import Basemap, shiftgrid
from matplotlib.colors import LinearSegmentedColormap
from pydap.client import open_url
from pydap.proxy import ArrayProxy
import scipy

data_url = 'http://thredds.ucar.edu/thredds/dodsC/grib/NCEP/RAP/CONUS_13km/RR_CONUS_13km_' + '20150415' + '_' + '01' + '00.grib2/GC'

print('Getting Data from URL:\n\n    "{0}"\n'.format(data_url))

# Create Array of all data from URL
dataset = open_url(data_url)

# Map Projection Info

proj_attributes = dataset['LambertConformal_Projection'].attributes
rsphere = proj_attributes['earth_radius']

lat_0 = proj_attributes['latitude_of_projection_origin']
lon_0 = proj_attributes['longitude_of_central_meridian'] 
lat_1 = proj_attributes['standard_parallel']

llcrnrlat = 16.28100013732909 # (1,1)
llcrnrlon = 360-126.18 # (1,1)

urcrnrlat = 55.552133975329625 # (614,428)
urcrnrlon = 360-59.15590040502627 # (614,248)

x = np.array(dataset['x'][:])
y = np.array(dataset['y'][:])

def xy_converter(var):
    """
    Downloads entered variable (x or y) coordinates 
    and converts from m to km.  Inputs for var 
    should be 'x' or 'y'.
    """
    values = dataset[var][:]
    data_array = values * 1000
    newarray = data_array + abs(data_array.min())
    return newarray

# Download x & y coord. and convert m to km
x = xy_converter('x')
y = xy_converter('y')

# Temp Contour
temp_2m = dataset['Temperature_height_above_ground'].array[1,:,:,:]-273.
temp_2m = temp_2m * (9./5.) + 32.
temp_2m    = temp_2m.squeeze() 

#plot
fig = plt.figure(figsize=(11,11))
ax = fig.add_subplot(1,1,1)

map = Basemap(projection='lcc', lat_0 = lat_0, lon_0 = lon_0,
                                     llcrnrlon = llcrnrlon, llcrnrlat = llcrnrlat,
                                     urcrnrlat = urcrnrlat, urcrnrlon = urcrnrlon,
                                     area_thresh = 1000., rsphere = rsphere, resolution='i')
map.drawcoastlines(linewidth=0.3)
map.drawcountries(linewidth=0.3)
map.drawcounties(linewidth=0.1)
map.drawstates(linewidth=0.3)
map.drawmapboundary(linewidth=0.5)

#lons,lats = basemap_parameters.map(x,y)
#lon,lat = basemap_parameters.map(lons,lats,inverse=True)

#ny = range(len(y)); nx = range(len(x))
#ny = temp_2m.shape[0]; nx = temp_2m.shape[1]
#lons, lats = map(x, y) # get lat/lons of ny by nx evenly space grid.
#xx, yy = map(lons, lats)

levels = np.linspace(-42,122,320)
ticks = [-60,-50,-40,-30,-20,-10,0,10,20,30,40,50,60,70,80,90,100,110,120]  

plot = plt.contourf(x,y,temp_2m,levels,cmap='jet',extend='both')

# Set Colorbar Text Color
color_bar = map.colorbar(plot)
color_bar.set_ticks(ticks)

# CONUS
plt.xlim(x[30],x[440])
plt.ylim(y[30],y[290])    
plt.savefig('/home/public_html/conus_temp.png', dpi=100, bbox_inches='tight', pad_inches = .05)
Источник

1 ответ

  • Прежде всего, вы не должны определять свой экземпляр базовой карты как map - это слежка за встроенным в Python map() функция, и просто вызовет путаницу.
  • Во-вторых, ваш x а также y загруженные значения имеют разную длину (451 и 337 значений соответственно). Это нужно исправить, прежде чем мы попробуем что-нибудь еще.
  • Я не уверен, почему вы конвертируете свои данные x и y (проецируемые координаты в LCC?) Путем деления на 1000.

В любом случае разумный подход заключается в следующем:

Загрузите ваши данные и преобразуйте координаты из LCC в lon/lat с помощью Pyproj:

import pyproj
lc = pyproj.Proj("+proj=lcc +lat_1=33 +lat_2=45 +lat_0=40 +lon_0=-97 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs")
lons, lats = lc(x, y, inverse=True)
# lons & lats are now unprojected WGS84

Преобразуйте ваши координаты в желаемые координаты проекции карты:

# now you can get ll and ur lons and lats
# set up your basemap with whatever projection you'd like, e.g.
m = Basemap(
    projection = 'merc',
    llcrnrlon = lllon, llcrnrlat = lllat, urcrnrlon = urlon, urcrnrlat = urlat,
    resolution='h')
# project lons & lats into map coordinates
projected_lons, projected_lats = m(x, y)

Теперь мы можем сделать другие вещи:

Я не уверен, что ваши данные температуры уже интерполированы в сетку. Если это не так:

from matplotlib.mlab import griddata
# set up a square grid with the same extents as our measured data
numcols, numrows = 1000, 1000
xi = np.linspace(projected_lons.min(), projected_lons.max(), numcols)
yi = np.linspace(projected_lats.min(), projected_lats.max(), numrows)
# get lon and lat coords of our grid points
xi, yi = np.meshgrid(xi, yi)
# interpolate
x, y, z = projected_lons, projected_lats, temp2m
zi = griddata(x, y, z, xi, yi)

Настройте другие детали вашей карты, как обычно, и нанесите контур, используя контур (не используйте карту цветов струи!)

conf = m.contourf(xi, yi, zi, cmap='coolwarm', ax=ax)
Источник
Другие вопросы по тегам