Записки океанолога - обработка и визуализация данных
NAO индекс в netCDF формате
Задача: перевести Индекс Северо-Атлантической Осциляции (NAO) из ASCII в netCDF формат
Решение: используем модули Python - PyNio, numpy, time
Индекском NAO (а также его близким родственником AO) пользуется огромное количество народу, но найти его в netCDF формате мне не удалось. Пришлось делать самому. Кому нужен просто файлик - вот он NAO index in netCDF format (up to 2011.04). Кто хочет посмотреть на очередной пример использования Nio для создания netCDF файла, велкам под кат.
Более подробно про NAO можно почитать в википедии. Если коротко, то NAO это индекс, характеризующий изменчивость атмосферного давления, с которым любят связывать (зачастую успешно) различные явления в океане и атмосфере - от ледовитости СЛО, до частоты тропических циклонов. Его упоминание можно найти чуть не в каждой статье посвященной Арктике или северной Атлантике. Тем более странно, что при такой популярности версии NAO в формате netCDF в интернете не нашлось, обычно он распространяется в виде текстовых файлов. Чтобы исправить эту несправедливость я написал небольшой скрипт, при помощи которого можно в принципе конвертировать любой текстовый файл вида:
Год Месяц Значение
в netCDF.
Значения индекса брались с этого сайта (там же можете найти и AO).
Сам скрипт:
#to netCDF file
#
#Created by Nikolay Koldunov
#koldunovn@gmail.com
#Description in Russian at http://koldunov.net/?p=521
import numpy
import Nio
import os
import time
def timetostep(start_date_time,time_step, present_time):
""" Convert date to amount of timesteps since start date
Usage:
timetostep(start_date_time, time_step, present_time )
Input:
start_date_time - should be string in form of YYYYDDMMhhmmss
time_step - model timestep (deltaT) in seconds
present_time - should be string in form of YYYYDDMMhhmmss
Output: time step
"""
start_date_time_in_python_format = time.strptime(start_date_time,"%Y%m%d%H%M%S")
start_date_time_in_seconds = time.mktime(start_date_time_in_python_format)
present_time_in_python_format = time.strptime(present_time,"%Y%m%d%H%M%S")
present_time_in_seconds = time.mktime(present_time_in_python_format)
seconds_from_start_date = present_time_in_seconds - start_date_time_in_seconds
present_time_step = seconds_from_start_date/time_step
return present_time_step
os.system("rm NAO_conv.nc")
input_file = './norm.nao.monthly.b5001.current.ascii'
ifile = open(input_file, 'r')
lines = ifile.readlines()
nao = numpy.array(())
ttime = numpy.array(())
for line in lines[:]:
ttime = numpy.append(ttime, timetostep("19480101000000", 3600 , line.split()[0]+line.split()[1].zfill(2)+"15000000"))
nao = numpy.append(nao,float(line.split()[2]))
opt = Nio.options()
opt.PreFill = False
opt.HeaderReserveSpace = 4000
f = Nio.open_file("NAO_conv.nc","w",opt)
f.title = "NAO index in netCDF format"
f.source = "http://www.cpc.ncep.noaa.gov/products/precip/CWlink/pna/nao.shtml"
f.author = "Nikolay Koldunov, koldunovn@gmail.com"
f.url = "http://koldunov.net/?p=521"
f.create_dimension('time',ttime.shape[0])
f.create_variable('time','d',('time',))
f.variables['time'].units = "hours since 1948-01-01 00:00:00"
f.variables['time'].calendar = "proleptic_gregorian"
f.variables['time'][:] = ttime[:]
f.create_variable('NAO','d',('time',))
f.variables['NAO'].long_name = "NAO index"
f.variables['NAO'].units = "non dimensional"
f.variables['NAO'][:] = nao
f.close()
Небольшие пояснения. Функция timetostep работает очень похоже на то что описано в этом посте, только наоборот :) У нее страшные имена переменных, но уж так мне захотелось сделать в тот момент когда я ее писал. Нужна она нам для того чтобы рассчитать для каждой даты, на которую у нас есть индекс (считаем что индекс задан на 15 число каждого месяца), количество часов прошедшее с "начального" момента времени, который будет стоять в аттрибуте units нашей переменной time. Напомню что в netCDF файле время может задаваться разными способами, и один из них "hours since" - количество часов, прошедшее с какого-нибудь момента времени. В нашем случае это "1948-01-01 00:00:00".
PyNio будет ругаться если файл, уже существует, так что стираем его, если он у нас вдруг уже есть.
Далее стандартно открываем текстовый файл, считываем его построчно в переменную lines, создаем пустые пока переменные nao и ttime. В цикле обрабатываем каждую строку. К переменной time присоединяем результат работы функции timetostep, которой частично передаем значения полученные из нашего текстового файла, а именно год и месяц. Заметьте, что при передаче месяца используется zfill(2), чтобы номер месяца всегда состоял из двух цифр (например не 3, а 03 для марта). Величина временного шага задана в 3600 - что соответствует одному часу в секундах.
Дальше создается netCDF файл, подобно тому как это описано в этом посте. Измерение (dimension) у нас только одно - время. Единицы измеререния - часы с 1 января 1948 года. Можно, при желании, сделать и "seconds since 1948-01-01 00:00:00", тогда нужно будет поменять 3600 на 1 при вызове функции timetostep.
Если вы сконвертируете при помощи этого скрипта какие ни будь полезные индексы, то с удовольствием размещу их netCDF файлы здесь, либо на Oceanographers.RU, для использования людьми которые не так нежно любят Python как мы с вами :)
- 26 Май 15:25
Задача: отобразить наши данные на Google Earth
