Картографическая проекция и принудительная интерполяция

У меня странное поведение при использовании разных проекций базовой карты.

Сетка измерений, которую я хочу нанести на карту мира, имеет форму [181,83]. Это означает, что у меня есть значения для каждой точки 2 ° / 2 ° в диапазоне от -180 ° до 180 ° долготы и от -82 ° до 82 ° широты.

from mpl_toolkits.basemap import Basemap
import numpy as np
measurementgrid = np.random.random_sample((181,83))
m = Basemap(projection='cyl',llcrnrlon=-180, llcrnrlat=-82, urcrnrlon=180, urcrnrlat=82, resolution='l')
m.drawcoastlines()
m.drawparallels(np.arange(-90,90,30),labels=[1,0,0,0])
m.drawmeridians(np.arange(-180,180,45), labels=[0,0,0,1])
data, x, y = m.transform_scalar(measurementgrid.T, lons=np.arange(-180,182,2), lats=np.arange(-82,84,2), nx = 181, ny = 83, returnxy=True, order=0)
m.imshow(data, origin='lower', interpolation='none')

Используя цилиндрическую проекцию, возвращенная сетка данных совпадает с сеткой измерений, и все в порядке. Если я изменю проекцию на «мельницу», результирующие интерполированные данные будут отличаться от их источника.

Есть ли способ построить сетку измерений как есть, но с учетом меняющейся проекции?


python matplotlib matplotlib-basemap
person nit    schedule 12.08.2013    source источник
comment
В вашем примере кода не показано, откуда вы получаете measurementgrid (он также не показывает import numpy as np, но мне удалось это угадать!   -  person Steve Barnes    schedule 12.08.2013
comment
Могу я предложить опубликовать свое исправление, чтобы помочь кому-нибудь еще, у кого есть аналогичная проблема!   -  person Steve Barnes    schedule 12.08.2013
comment
Я не исправил проблему, но отредактировал свой пост согласно вашему предложению. Так что проблема все еще присутствует!   -  person nit    schedule 13.08.2013

Ответы (1)


arrow_upward
3
arrow_downward

Во-первых, я бы посоветовал вам начать использовать pcolormesh, а не imshow. Pcolormesh должен быть вашим инструментом визуализации при попытке отобразить данные с координатной сеткой в ​​прямоугольниках (как contourf при рисовании данных с координатной сеткой в ​​виде областей с одинаковыми значениями).

Чтобы использовать pcolormesh, вы должны передать координаты x и y углов ваших данных, поэтому:

x = np.linspace(-180, 180, 182)
y = np.linspace(-90, 90, 84)
m.pcolormesh(x, y, data)

Но с помощью базовой карты вы всегда должны преобразовывать координаты в систему координат карты - в конечном итоге это может означать, что координаты x и y должны быть двухмерными массивами, поэтому мы делаем это и конвертируем:

x = np.linspace(-180, 180, 182)
y = np.linspace(-90, 90, 84)
x, y = np.meshgrid(x, y)
converted_x, converted_y = m(x, y)
m.pcolormesh(converted_x, converted_y, data)

Это означает, что теперь вы можете изменить прогноз, и ваши данные будут отображены в нужном месте. Например, я изменил проекцию на «Робин» (Робинсон) и получил следующую картинку:

Результат моего исправления

К сожалению, pcolormesh предназначен для смежных блоков данных, которые, если вы выберете проекцию, которая не имеет одинаковой центральной долготы (также известной как «lon_0»), вы получите плохие результаты. Например, я изменил проекцию на projection='robin',lon_0=180 и получил следующее изображение:

Плохая строка даты

Это потому, что шкала дат в настоящее время не обрабатывается Basemap, и, насколько я могу судить, без серьезной перезаписи - никогда не будет.

Хорошая новость заключается в том, что эта область беспокоила меня в течение долгого времени, поэтому я начал писать новый пакет, чтобы справиться с этой и многими другими причудами картографии для научной визуализации. В результате получился новый пакет под названием cartopy, который выполняет гораздо больше работы за вас, так что такие вещи, как строка даты, «просто работают»:

import numpy as np
import cartopy.crs as ccrs
import matplotlib.pyplot as plt

x = np.linspace(-180, 180, 182)
y = np.linspace(-90, 90, 84)
measurement_grid = np.random.random_sample((83, 181)) * y[:-1, np.newaxis] ** 2

plt.axes(projection=ccrs.Robinson(central_longitude=180))
plt.pcolormesh(x, y, measurement_grid, transform=ccrs.PlateCarree())
plt.gca().coastlines()
plt.show()

Картография

Хотя я не предлагаю вам перейти на cartopy сейчас (установка и производительность все еще в процессе), но стоит знать, что пакет существует, и я ожидаю, что в будущем он станет все более и более привлекательным, когда вы столкнетесь с такими видами. вопросов. http://scitools.org.uk/cartopy/docs/latest

Также стоит отметить, что многие проблемы, возникающие при научной визуализации данных с координатной привязкой, связаны с обработкой данных, их координатами и лежащими в их основе системами координат, поэтому был написан другой пакет, который реализует модель данных для инкапсуляции всех данных. из этой сложной информации в единый объект, который затем может быть передан подпрограммам построения графиков для простого взаимодействия. Я снова рекомендую вам взглянуть на него http://scitools.org.uk/iris/docs/latest.

HTH

person pelson    schedule 15.08.2013