Визуализация скалярного поля в PyNGL
Задача: визуализировать концентрацию льда
Инструменты: Python с модулем PyNGL
Продолжаем тему визуализации при помощи питоновского модуля PyNGL. Подробно о том что это такое и зачем оно нам можно посмотреть в . Здесь же мы пошагово рассмотрим создание такой вот картинки:
Исходные данные находятся в netCDF файле, который можно .
Начнём работать сразу со скрипта который что то нам нарисует:
import Ngl
import Nio
import numpy
file = Nio.open_file("NCEP_MR109_sico_1980-1999_mean_09.nc")
concentration_nio = file.variables["SICOMO"]
lon = file.variables["lon"]
lat = file.variables["lat"]
concentration = concentration_nio[0,0,:,:]
lon = lon[:]
lat = lat[:]
#
# Open a workstation.
#
rlist = Ngl.Resources()
rlist.wkColorMap = 'posneg_1'
wks_type = "ps"
wks = Ngl.open_wks(wks_type,"contours",rlist)
resources = Ngl.Resources()
resources.sfXArray = lon[:]
resources.sfYArray = lat[:]
map = Ngl.contour_map(wks,concentration[:,:],resources)
Назовём его contours.py. Происходит тут следующее - мы открываем файл NCEP_MR109_sico_1980-1999_mean_09.nc, вытягиваем из него переменные
- SICOMO - концентрация льда
- lon - долгота
- lat - широта
Напомню что посмотреть на переменные находящиеся в открытом файле можно набрав в ipython:
Далее переменные преобразуются в NymPy матрицы. А вот после этого идёт интересное, мы создаём рабочую станцию:
- rlist = Ngl.Resources() - создаём переменную в которой будут жить ресурсы относящиеся к рабочей станции
- rlist.wkColorMap = 'posneg_1' - выбираем цветовую палитру
- wks_type = "ps" - выбираем выходной формат картинки
- wks = Ngl.open_wks(wks_type,"contours",rlist) - создаём рабочую станцию.
Цветовую палитру (color map) можно либо выбрать из . Палитра это мегаважно, потому что вы для своего рисунка сможете использовать только цвета из этой палитры. Все цвета пронумерованы, и именно по номерам (в данном случае) мы будем к ним обращаться. Первые два цвета обычно белые, поскольку по умолчанию цвет 0 - задний план (background) рабочей станции (грубо всего рисунка), цвет 1 - передний план рабочей станции (foreground). Есть ещё цвет -1, который обозначает прозрачность.
Дальше:
- resources = Ngl.Resources() создаём переменную в которой будут все ресурсы относящиеся к карте
- resources.sfXArray = lon[:] сообщаем где у нас содержатся долготы
- resources.sfYArray = lat[:] то же с широтами
Запустив скрипт получим такое вот чудо-юдо:
Что то похожее на изолинии в районе полюсов наблюдается (там где должен быть лёд) но разобрать ничего, конечно, нельзя. Понятно что дефолтные значения нас не устраивают, будем менять. Все последующие куски кода нужно вставлять в ваш скрипт перед командой на отрисовку:
map = Ngl.contour_map(wks,concentration[:,:],resources).
И всё что мы будем делать это изменять эту самую resources.
Давайте сконцентрируемся на Арктике, для этого добавим в скрипт:
resources.mpDataBaseVersion = "MediumRes"
resources.mpLimitMode = "LatLon"
resources.mpMinLonF = 0
resources.mpMaxLonF = 360
resources.mpMinLatF = 60
resources.mpMaxLatF = 90
В принципе названия ресурсов и их значения говорят сами за себя. Мы выбираем стереографическую проекцию из , не самое высокое разрешение береговой линии и сообщаем что выбирать границы карты мы будем по широтам и долготам. Дальше задаём границы карты которая будет отображаться.
Получаем:
Как это не грустно, но уже этот зародыш карты лучше чем многое из того что можно встретить даже в топовых научных журналах :).
Хоть это и не по международным правилам, но давайте повернём нашу карту так чтобы Россия была внизу, так многим привычнее, и разрешим заливку.
resources.mpFillOn = True
Мы задали центральную долготу и сказали что хотим залить нашу карту.
Получаем:
Тут нужно сказать о заливке. Большинство областей на карте принадлежит к так называемым группам. Основные - Океан, Суша, Воды суши. Контролируется она ресурсом mpFillColors которому передаются цвета заливки различных групп. Первые четыре это:
resources.mpFillColors = [Не принадлежащие ни какой группе, Океан, Суша, Вода окружённая землёй]
Мы ресурс не задавали, поэтому для нашей карты были использованы дефолтные значения:
resources.mpFillColors = [16, 10, 8, 10]
и убедитесь что океан и суша залиты именно этими цветами.
Но допустим что у нас есть собственное мнение по поводу цвета заливки, а кое что мы вообще заливать не хотим. Добавляем:
resources.mpFillColors = [0,-1,igray,-1]
resources.cnLineDrawOrder = "Predraw"
Мы хотим залить землю серым цветом, а океан и воду суши оставить прозрачными. Но серого у нас в палитре нет. Это не беда, поскольку существует функция Ngl.new_color которая добавляет в конец существующей палиты ещё один цвет. Цвет определяется в RGB но значения не от 0 до 255 а от 0 до 1. Не очень понятно почему так, но как уж есть. Во второй строчке используем ресурс mpFillColors, в результате заливка:
- Области не принадлежащие никакой группе = 0 = цвет бэкграунда (у нас белый)
- Океан = -1 = прозрачность
- Суша = igray = серый заданный нами
- Воды суши = -1 = прозрачность
С таким же успехом вместо igray можно написать 21 (в нашей палитре 20 цветов и мы добавили ещё один - серый).
Также последней строчкой мы заставляем изолинии рисоваться вначале, чтобы потом всякие каки, возникающие у берегов, спокойно перекрыть сушей 
Получаем:
Давайте теперь раскрасим всё это дело:
resources.cnFillDrawOrder = "Predraw"
resources.cnLineLabelsOn = False
Разрешаем заливку контуров, говорим что контуры заливаются первыми и заодно выключаем не нужные теперь подписи на изолиниях.
Получаем:
В принципе уже живенько ). Белая полоска связана с тем что, насколько я понимаю, тут сходятся края матрицы. Бороться с этим можно, но я пока не знаю как
Давайте теперь чуть поиграем с цветами:
resources.nglSpreadColorEnd = -2
resources.cnLevelSelectionMode = "ExplicitLevels" # Define own levels.
resources.cnLevels = numpy.arange(0.,1.,0.1)
Начнём рисовать с цвета номер 8, а закончим вторым с конца, поскольку первый с конца это добавленный нами серый. Последние две строчки показывают как задавать свои собственные уровни, но у нас они адекватно определились, так что для нас в этом смысле ничего не поменяется.
Получаем:
Ну и, наконец, давайте займёмся легендой, которая сейчас уж больно огромная:
resources.lbTitleFontHeightF = 0.018
resources.lbLabelFontHeightF = 0.012
resources.lbTitleOffsetF = -0.27
resources.lbBoxMinorExtentF = 0.15
resources.pmLabelBarOrthogonalPosF = -0.01
resources.lbOrientation = "Horizontal"
lbTitleString отвечает за заголовок легенды, и мы вместо того чтобы явно задать строчку, используем информацию которая содержится в исходном файле. У переменной concentration_nio, которая связана с SICOMO есть атрибут __dict__. Он содержит такой словарь:
Комбинируя long_name и units мы получаем название нашей легенды.
Остальные параметры интуитивно понятны. lbOrientation делает легенду горизонтальной.
Ну и напоследок вспомним что у нас модельные данные, а значит всё на сетке:
resources.cnLinesOn = False
Включаем другой метод заполнения и выключаем изолинии.
Теперь избавимся от дурацкой белой линии посередине. Для этого надо изменить код тут:
resources.sfYArray = Ngl.add_cyclic(lat[:])
и в самом конце:
Функция Ngl.add_cyclic дополняет наш массив ещё одной колонкой, то есть ещё одним рядом "долготы". Огромное спасибо Маше за наводку на эту функцию!
Получаем:
По моему очень живенько
Полностью скрипт:
import Ngl
import Nio
import numpy
file = Nio.open_file("NCEP_MR109_sico_1980-1999_mean_09.nc")
concentration_nio = file.variables["SICOMO"]
lon = file.variables["lon"]
lat = file.variables["lat"]
concentration = concentration_nio[0,0,:,:]
lon = lon[:]
lat = lat[:]
#
# Open a workstation.
#
rlist = Ngl.Resources()
rlist.wkColorMap = 'posneg_1'
wks_type = "ps"
wks = Ngl.open_wks(wks_type,"contours",rlist)
resources = Ngl.Resources()
resources.sfXArray = Ngl.add_cyclic(lon[:])
resources.sfYArray = Ngl.add_cyclic(lat[:])
resources.mpProjection = "Stereographic"
resources.mpDataBaseVersion = "MediumRes"
resources.mpLimitMode = "LatLon"
resources.mpMinLonF = 0
resources.mpMaxLonF = 360
resources.mpMinLatF = 60
resources.mpMaxLatF = 90
resources.mpCenterLonF = 100.
resources.mpFillOn = True
igray = Ngl.new_color(wks,0.7,0.7,0.7)
resources.mpFillColors = [0,-1,igray,-1]
resources.cnLineDrawOrder = "Predraw"
resources.cnFillOn = True
resources.cnFillDrawOrder = "Predraw"
resources.cnLineLabelsOn = False
resources.nglSpreadColorStart = 8
resources.nglSpreadColorEnd = -2
resources.cnLevelSelectionMode = "ExplicitLevels" # Define own levels.
resources.cnLevels = numpy.arange(0.,1.,0.1)
resources.lbTitleString = concentration_nio.__dict__['long_name']+', '+concentration_nio.__dict__['units']
resources.lbTitleFontHeightF = 0.018
resources.lbLabelFontHeightF = 0.012
resources.lbTitleOffsetF = -0.27
resources.lbBoxMinorExtentF = 0.15
resources.pmLabelBarOrthogonalPosF = -0.01
resources.lbOrientation = "Horizontal"
#resources.cnFillMode = "RasterFill"
#resources.cnLinesOn = False
map = Ngl.contour_map(wks,Ngl.add_cyclic(concentration[:,:]),resources)
tut rasskazano kak mozhno beluy polosochku ubrat’
2masha
Спасибо, исправил пост, теперь финальная версия без дырки )