расчет NDVI.

[Kolpa4ok]

всем добрый день:)пишу дипломную работу для расчета NDVI.
насколько я понял наса сейчас сохраняет файлы в hdf5. у меня есть hdf-файл, вида MOD03.A2014158.0715.006.2014158133339.hdf, прибора MODIS, и мне надо каким-то магическим способом вытянуть из него только значения красной и инфракрасной области излучения.
проблема в том, что я еще плохо понимаю структуру hdf-файла:(
на данный момент у меня есть небольшой кусок кода, который не может открыть файл(
может кто-нибудь пытался реализовать свое по и сможет меня направить на правильное направление?
вот мой код:
using HDF5DotNet;

namespace hdf5
{
class Program
{
public class H5Dataset
{
public H5GroupId id;
public String Name;
public ulong lIndex;
public long[] dims;
public H5Dataset(H5GroupId _idg, String strN, ulong _idx)
{
id = _idg;
Name = strN;
lIndex = _idx;
dims = null;
}

}
public class H5Group
{
public H5GroupId id;
public H5FileId idFile;
public String Name;
public ulong lIndex;
public List<H5Dataset> datasets;
public H5Group(H5GroupId _id, String strN, ulong _idx)
{
id = _id;
Name = strN;
lIndex = _idx;
datasets = null;
}

}
static void Main(string[] args)
{
H5FileId idFile = H5F.open("MOD03.A2014158.0715.006.2014158133339.hdf", H5F.OpenMode.ACC_RDONLY);
H5GroupId idGroup = H5G.open(idFile, "/");
List<H5Group> groups = new List<H5Group>();
ulong lIndex = 0;

}
}
}

[nadiopt]

код лучше прятать за кнопу code
поиск по форуму рулит
http://gis-lab.info/docs/hdf.html
http://gis-lab.info/qa/hdfimport-sa.html
ну и вообще http://gis-lab.info/scripts/search.cgi? ... 1%81%D0%BA!

[bolotoved]

HDF можно читать GDAL'ом
Там главное выяснить, как называются нужные вам слои (subdatasets).

Вот кусок кода на питоне с использованием GDAL, как раз для экстракции каналов, в т.ч. красного и инфракрасного:

Код: Выделить всё

# -*- coding: utf-8 -*-
import os, re
modis_dir = '/data/pool/MOD13/'

os.chdir(modis_dir)

# читаем метаданные HDF
scene_gdalinfo = os.popen("gdalinfo " + 'MOD13Q1.A2012225.h21v02.005.2012242083001.hdf').read()

# Смотрим, как называется слой NDVI (SUBDATASET_1)
ndvi_layer = ''.join(re.findall('SUBDATASET_1_NAME.*?\n', scene_gdalinfo))
ndvi_layer = ndvi_layer.replace('\n', '')
ndvi_layer = ndvi_layer.replace('SUBDATASET_1_NAME=', '')
# Смотрим, как называется слой red (SUBDATASET_4)
red_layer = ''.join(re.findall('SUBDATASET_4_NAME.*?\n', scene_gdalinfo))
red_layer = red_layer.replace('\n', '')
red_layer = red_layer.replace('SUBDATASET_4_NAME=', '')
# Смотрим, как называется слой nir (SUBDATASET_5)
nir_layer = ''.join(re.findall('SUBDATASET_5_NAME.*?\n', scene_gdalinfo))
nir_layer = nir_layer.replace('\n', '')
nir_layer = nir_layer.replace('SUBDATASET_5_NAME=', '')
# Смотрим, как называется слой mir (SUBDATASET_7)
mir_layer = ''.join(re.findall('SUBDATASET_7_NAME.*?\n', scene_gdalinfo))
mir_layer = mir_layer.replace('\n', '')
mir_layer = mir_layer.replace('SUBDATASET_7_NAME=', '')

# Преобразуем в TIF
os.system('gdal_translate -of GTiff '  + '"' + ndvi_layer + '"' +' ndvi.tif')
os.system('gdal_translate -of GTiff '  + '"' + red_layer + '"' +' 01_red.tif')
os.system('gdal_translate -of GTiff '  + '"' + nir_layer + '"' +' 02_nir.tif')
os.system('gdal_translate -of GTiff '  + '"' + mir_layer + '"' +' 07_mir.tif')