CC BY
14
Региональные проблемы. 2013. Том 16, № 2. С. 84-90.
УДК630*43(571.621)
ИНТЕРПОЛЯЦИЯ КОМПЛЕКСНОГО ПОКАЗАТЕЛЯ ПОЖАРНОЙ ОПАСНОСТИ НА ТЕРРИТОРИИ ЕВРЕЙСКОЙ АВТОНОМНОЙ ОБЛАСТИ И ХАБАРОВСКОГО КРАЯ
В.А. Глаголев, P.M. Коган Институт комплексного анализа региональных проблем ДВО РАН, ул. Шолом-Алейхема 4, г. Биробиджан, 679016, e-mail: [email protected], [email protected]
Исследовано применение методов интерполяции для восстановления значений лесопожарного показателя пожарной опасности по условиям погоды на территории Хабаровского края и Еврейской автономной области. Ключевые слова: комплексный показатель пожарной опасности, пожары, интерполяция.
Погодные условия являются одним из основных факторов, влияющих на вероятность возникновения и скорость распространения пожаров растительности. В настоящее время густота сетей метеостанций на территории Российской Федерации не позволяет в полной мере проводить мониторинг пожароопасной обстановки, поскольку расчет показателей по условиям погоды осуществляется в основном по данным единичной гидрометеостанции (ГМС), зона наблюдений которой составляет 30 км [5, 10]. Применение современных средств дистанционного зондирования позволяет вести территориально распределенное наблюдение за погодными условиями, однако дискретный характер некоторых метеоданных, таких как суточный объем осадков и скорость ветра, а также пропуски информации на спутниковых снимках в случае облачности накладывают ограничения на их использование. Поэтому для нахождения показателей пожарной опасности на растительной территории, особенно в кварталах участковых лесничеств, находящихся на значительном удалении от метеостанций и пунктов охраны, необходимо применять различные методы пространственного восстановления данных.
Целью работы является выбор метода восстановления показателей пожарной опасности по условиям погоды и верификация его на территории Дальнего Востока России для расчета показателей пожарной опасности и вероятности возникновения пожаров растительности. Для этого необходимо выполнить следующие задачи: анализ методов интерполяции и возможности применения их для оценки пожарной опасности по условиям погоды; определение параметров проверки методов; разработка модулей реализации интерполяции и проведения ее кросс-валидации; верифицикация восстановленных значений лесопожарного показателя.
Известны две группы интерполяционных методов с географической привязкой данных [3]: детерминистические и геостатистические. Общим принципом является нахождение весовой функции, отражающей динамическое изменение показателей между узлами интерполирования, а также проверка сравнением восстановленных и фактических данных по выбранным узлам интерполирования. Детерминистические широко используются в раз-
личных областях прикладной и научной деятельности, они базируются на технике создания поверхностей по известным значениям в отдельных точках путем продления закономерностей распределения значений по мере удаления от точек, или путем сглаживания различий. К этой группе относят методы: обратно-взвешенные расстояния [16]; Шепарда [18]; минимизация кривизны [20]; полиномиальная регрессия [15]; базисные радиальные функции [19]; триангуляция Делоне с линейной интерполяцией [12]; естественные соседи [12]; локальная полиноминальная интерполяция [16]. Геостатистические методы базируются на статистических закономерностях и используются для более точного и сложного моделирования поверхностей, включая оценку ошибок и вычисление вероятностных параметров поверхностей, например «Кригинг» [14].
Для осуществления интерполяции необходимо выбрать показатель, имеющий монотонный и непрерывный характер. При количественной оценке пожарной опасности на территории РФ используются два показателя, учитывающих ежедневное высыхание и увлажнение растительных горючих материалов [6]. Первый из них - ле-сопожарный показатель Ь на текущий г-ый день рассчитывается на основе дневной температуры воздуха и точки росы (ф.), измеренных в 13-15 часов дня местного
времени: Ц = — Г,) (1).
Вторым является комплексный показатель Р., для расчета которого используется значение Ь. и сведения об объеме выпавших осадков (х) с 9 часов утра предыдущего дня до 9 часов утра текущего дня. Для «сухих» дней (при х< 3 мм/сут.) Р вычисляется путем ежедневного суммирования Ь данного и предшествующих дней, а при наступлении «мокрых» дней (при х. ? 3 мм/сут.) Р равен Ь. текущего дня:
{Ь(т-), при х- > 3мм/сут.
) + Р- р при хг < 3 мм/сут
или:
Р^Щ^ + Р^Щ), (3),
где коэффициент К(х) задается на основе предельного
Рис. 1. Распределение значений лесопожарного Ь. и комплексного Р. показателей в весенний период 1989 г., по данным метеостанции «Биробиджан»
количества осадков, при котором устраняется возможность загорания напочвенного покрова; он может принимать два значения: при х. более 3 мм/сут. К(х) равен 0; при х. менее 3 мм/сут. К(х} равен 1.
Для интерполирования можно выбрать только показатель L., поскольку функция № 1, по которой он определяется, в отличие от показателя Р. (функции № 2, 3) является монотонной и носит непрерывный характер без ограничений (рис. 1) [7].
Впервые восстановление показателей пожарной опасности описано Н.П. Курбатским [5]. Им предложено находить значения показателя Р. в каждом квартале участкового лесничества по распределению фактических метеоданных соседних метеостанций методом обратно взвешенные расстояния, при этом весовые функции принимают значение аналитической функции комплексного аргумента (полином Лагранжа). Затем было показано, что весовая функция комплексного аргумента в виде полинома третьей степени достаточно достоверно отражает распределение метеоданных по сравнению с полиномом Лагранжа первой и второй степени [2].
Три интерполяционных метода (обратно взвешенные расстояния, локальная полиномиальная интерполяция второй степени и полиномиальная регрессия) использованы для восстановления метеорологического индекса пожарной опасности Fire Weather Index (FWI) на примере территории Онтарио в Канаде [16]. Обработка осуществлялась по данным 26 ГМС, одна из которых взята в качестве контрольной (проверочной). Оценка качества восстановления производилась на основе среднеквадра-тического отклонения фактического значения индекса FWI для контрольной метеостанции. Показано, что эти методы удовлетворительно осуществляют восстановление FWI, при этом лучшим оказался метод - полиномиальная регрессия.
С начала 90-х гг. XX в. для восстановления показателей пожарной опасности начинают пользоваться информацией, поступающей от спутников серии NOAA [9, 11]. В любой точке снимка радиационная температура видимой поверхности и температура точки росы находились с помощью метода триангуляция Делоне с линейной
интерполяцией, а суточный объем осадков восстанавливался методом обратно взвешенные расстояния по фактическим данным смежных ГМС.
В Финляндии с 1996 г. метеорологический индекс пожарной опасности FFFI (Finnish Forest Fire Index) рассчитывается в каждой ячейке размером 10x10 км по восстановленным метеорологическим данным (влажность воздуха, скорость ветра, солнечная радиация и суточный объем осадков) 26 ГМС с помощью интерполяционного метода Кригинг [13].
В перечисленных выше примерах не обоснован выбор метода интерполяции, поэтому нами они были проверены для восстановления показателя L. на территории Хабаровского края и ЕАО. Оценивание качества каждого метода в контрольной точке (т) проводилось по двум параметрам:
среднеквадратическая ошибка интерполяции RMSEm (Root Mean Square Error):
RMSEm=^{LrLinf ;
где N- количество ежедневных расчетов; L. и L.n - значения фактического и интерполяционного значения показателя в i-ый день [3,16];
коэффициент Нэша-Сатклиффа или эффективность метода [17]:
Uh-hnf
Е = 1-—-
т N / _ у ,
где L - среднее значение показателя.
При выборе метода интерполяции учитывались следующие требования: среднеквадратическая ошибка метода RMSEm должна иметь наименьшее значение; коэффициент эффективности метода Ет должен быть больше 0,5 (1 метод является идеальным; > 0,75 - хорошим; > 0,5 удовлетворительным; >0,5 - неудовлетворительным).
Для выполнения поставленных задач разработан и реализован алгоритм интерполяции показателя L. и про-
Контрольная точка, список ГМС, сезон, метод
гС
Выполнять для каждого дня сезона
Расчет значения показателя
Задание координат узлов интерполирования
Выполнение метода интерполирования
Создание массива расчетных и восстановленных значений показателя
X
Г)
3
ü~
Экспорт расчетных значений показателя в универсальный формат базы пространственных данных
Создание виртуального объекта картографического приложения
Сохранение результата интерполирования в универсальном формате обмена
метеор ологических данных (MySQL)
пр о стр анств енных данных (csv)
Расчет среднеквадратической ошибки RMSBm и коэффициента эффективности Вт
Таблица "Интерполяция" (MySQL)
Рис. 2. Блок-схема проверки метода интерполяции показателя пожарной опасности в контрольной точке
верки метода в контрольной точке (рис. 2). Вначале задаются координаты контрольной точки, номер метода интерполяции, год пожароопасного сезона. Контрольная точка представлена отдельной ГМС, и при интерполяции она не используется, в то время как остальные ГМС применяются в качестве узлов интерполяции на каждый день пожароопасного сезона. Затем значения показате-
ля L. рассчитываются в информационной системе оценки и прогноза показателей пожарной опасности [1 ] и экспортируются в текстовый файл с расширением «csv».
Авторский скрипт загружает текстовый файл и передает их виртуальному объекту картографического приложения по технологии Microsoft ActiveX (пример приведен на листинге 1 ).
Sub Main 'Начало скрипта
Dim SurferАрр As Object 1 задать переменную виртуального приложения Set SurferApp = CreateObject("Surfer.Application") 'создать виртуальный объект For YearO = 1960 То 1991 'цикл обработки пожароопасных сезонов с 1960 по 1991 г.г. For DayO = 1 То 214 'цикл обработки дней выделенного пожароопасного сезона (в сезоне 214 дней)
Data =YearO+"\"+DayO+".csv" 'Формирование пути к файлу данных для загрузки значений узлов интерполяции 'обратно-взвешенные расстояния
SurferApp.GridData(DataFile:=Data, xcol:=l, ycol:=2, zcol:=3, ShowReport:=False, Algorithm: =1, OutGrid:=Data+"_InverseDistance.grd") 'Восстановление значений ЛПЗ в ячейке ре1улярной сети (грид) GridFile =Data+"_InverseDistance.grd" 'переменная для хранения имени файла грида OutFile = Data+"_InverseDistance.csv" 'переменная для хранения имени файла данных Surfe г Арр. Gr idC onvert(InGrid: =GridFile, OutGrid: =OutFile, OutFmt:=srfGridFmtXYZ)'конвертация файла грида в файл данных Next 'следующий день Next 'следующий пожароопасный сезон End Sub 'Окончание скрипта
Листинг 1. Программный код авторского скрипта картографического пакета Surfer 10 для интерполяции методом обратно взвешенные расстояния
Таблица 1
Общее количество пожароопасных сезонов судовлетворительным коэффициентом эффективности интерполяции значений показателя!,, с 1960 по 1991 гг. по отношению к контрольным точкам
№ п/п Контрольные точки (ГМС) Метод интерполяции
Обратно взвешенные расстояния I Минимизация кривизны Метод Шепарда Метод естественных соседей Полиномиальная регрессия Базисные радиальные функции Триангуляция с линейной интерполяцией Локальная полиномиальная интерполяция
Количество сезонов с хорошим / удовлетворительным коэффициентом эффективности
Южная часть территории
1. Биробиджан 21/4 2/2 0/0 3/0 0/0 27/1 23/3 0/0 25/3
2. Победа 15/16 11/18 20/7 20/7 12/17 25/6 11/19 14/15 17/14
3. Кукан 14/14 9/19 10/18 14/13 10/18 23/3 9/19 10/18 13/16
4. Хабаровск 10/21 16/15 20/0 16/14 15/16 29/1 12/19 15/16 11/20
5. Сутур 12/20 11/21 17/0 18/8 0/0 26/5 10/22 0/0 12/20
6. Хорское 28/0 28/0 21/0 27/0 0/0 18/0 29/0 0/0 29/0
7. Чегдомын 28/2 27/3 1/0 0/0 30/1 29/0 30/1 29/1 30/0
8. Троицкое 24/3 23/4 24/0 12/0 18/10 25/11 20/7 20/8 19/9
Итого с коэффициентом > 0,5,из 256 сезонов, % 152/80 90/10 127/82 82/18 93/25 46/54 110/4 2 60/40 95/62 61/39 202/2 7 90/10 144/90 91/9 88/58 57/43 156/82 93/7
Центральная часть территории
9. Совгавань 7/0 5/0 0/0 0/0 0/0 4/0 5/0 0/0 0/0
10. Дуки 30/0 29/0 26/0 24/0 28/1 25/0 29/0 27/2 29/0
11. Полина Осипенко 19/10 21/8 25/2 21/0 21/8 29/0 15/14 24/5 17/13
12. Нижне-Тамбовское 10/21 13/18 22/6 21/0 11/20 25/6 11/20 25/6 14/17
13. Тумнин 20/0 17/0 17/0 16/0 17/0 24/0 20/0 17/0 15/0
14. Комсомольск-на-Амуре 26/1 26/0 26/0 25/0 26/0 28/0 26/0 26/0 27/0
Итого с коэффициентом >0,5 из 192 сезонов, % 112/32 75/25 111/26 71/29 91/8 52/48 107/0 56 103/28 68 135/6 73 106/34 73 119/13 69 102/30 69
Северная часть территории
15. Богородское 14/14 14/14 20/7 17/11 23/4 18/0 11/16 22/5 16/11
16. Херпучи 25/2 21/6 15/1 0/0 22/3 24/0 20/7 20/4 20/5
17. Чумикан 3/0 3/0 3/0 1/0 3/0 0/0 3/0 2/0 3/0
18. Николаевск-на-Амуре 16/3 14/3 6/0 0/0 16/1 11/0 16/3 20/1 13/0
19. Мареинское 4/28 11/20 28/1 26/0 9/22 28/2 7/25 9/22 9/23
20. Удское 26/3 21/3 6/0 0/0 0/0 16/0 25/3 0/0 21/2
Итого с коэффициентом >0,5 из 192 сезонов,% 88/50 72 82/46 67 76/8 44 42/11 28 70/30 46 98/2 52 75/54 67 71/32 54 79/41 57
В картографическом приложении при выполнении метода интерполяции создается грид файл *.§гс1, который конвертируется в файл данных *.сз\г, содержащий координаты центров ячеек грида и восстановленные значения лесопожарного показателя Ь..
Следующие шаги алгоритма выполняются в модуле
«Проверка метода интерполяции»: выбирается контрольная ГМС, по ее данным производиться расчет фактических значений показателя Ь. и сравнение с восстановленными значениями, которые находятся в центре ячейки грида. Затем производится расчет среднеквадра-тической ошибки ЯМБЕ ^ и коэффициента эффективнос-
Обратно-взвешенные расстояния
Кригинг
Минимизация кривизны
Полиномиальная регрессия
Естественные соседи
Метод Шепарда
Локальная полиномиальная интерполяция
Триангуляция с линейной интерполяцией
Базисные радиальные функции
Рис. 3. Электронные карты восстановления лесопожарного показателя в контрольной точке «Биробиджан», 10.05.1989 г. Примечание: • метеостанция;---административные границы; — изолинии
Хабаровский
Троицкое
Амурская область
Кукан
Облучье
Биробиджан
гМИДОВИЧИ
Условные обозначения
-граница области
•—« >• государственная граница
• гидрометеостанция Участки растительности | не лесного фонда I | лесного фонда Зоны гидрометеостанции | репрезентативности I ответственности
Ленинское
Екатерин&'Никольское,
масштаб 1:50000«
Рис. 4. Зоны влияния и ответственности метеостанций на территории Еврейской автономной области
Таблица 2
Достоверность прогноза пожаров растительности по погодным и лесорастительным условиям
№ п/п Степень пирологической пожароопасности основных растительных формация I II III IV I-IV
1. Фактическое количество пожаров 24 40 47 14 125
2. Расчетное количество пожаров 21 33 44 10 108
3. Достоверность, % при вероятности 0,95 Текущий день 0,88±0,14 0,83±0,12 0,94±0,07 0,71±0,26 0,86±0,06
1 -ый день прогноза 0,88±0,14 0,85±0,11 0,85±0,1 0,57±0,29 0,82±0,07
2-ый день прогноза 0,88±0,14 0,88±0,11 0,77±0,12 0,57±0,29 0,79±0,07
3-ый день прогноза 0,88±0,14 0,85±0,11 0,77±0,12 0,50±0,29 0,78±0,07
ти Ет в контрольной точке за каждый день пожароопасного сезона и сохранение результатов в таблице «Интерполяция» в базе метеорологических данных, состоящей из 6 числовых атрибутов: код расчета, код метеостанции, год, метод интерполяции, среднеквадратическая ошибка, коэффициент эффективности.
Определение среднемноголетнего значения средне-квадратической ошибки RMSEm, коэффициента эффективности Ет и количества пожароопасных сезонов с удовлетворительным коэффициентом эффективности Ет производится агрегирующими SQL-запросам из таблицы «Интерполяция».
В качестве инструмента интерполяции использован картографический пакет Surfer 10, в котором реализованы перечисленные выше детерминистические и геостатистические методы. Восстановление показателей пожарной опасности осуществлялось в ячейках прямоугольной регулярной сети, наложенной на территорию Хабаровского края и ЕАО размером 0,56°х0,56° в проекции Пужово 1942. Все узлы интерполяции представлены 22 ГМС, из них равнинных - 16 и горных - 6 (высота над уровнем моря свыше 200 м). По двум метеостанциям «Арка» и «Нелькан» (северная часть Хабаровского края) восстановленные значения не проверялись на достоверность ввиду их большой удаленности от оставшихся узлов интерполирования (около 340 км) и значительных отличий климатических условий на этих ГМС от всех остальных.
В каждой точке оценка методов осуществлялась ежедневно с 1960 по 1991 гг., для этого был построен ансамбль из 554658 электронных карт. На рис. 3. представлен пример восстановления значений показателя L. для контрольной точки «Биробиджан» 10 мая 1989 г.
Определение оптимального интерполяционного метода осуществлялось подсчетом случаев, когда коэффициент эффективности больше 0,5 и эффективность являлась удовлетворительной и хорошей (табл. 1).
Для территории южной части Хабаровского края и ЕАО оптимальным оказались несколько методов (обрат-но-взвешенные расстояния, полиномиальная регрессия, локальная полиномиальная интерполяция), поскольку в 90-93% сезонов коэффициенты эффективности превышали 0,5. В других частях Дальнего Востока эти же методы являются лучшими, но количество сезонов с удовлетворительным и хорошим коэффициентом эффективности не превышало 75% , что в большинстве случаев связа-
но с пропусками метеоданных в контрольных и проблемой восстановления значений в крайних точках (№ 9, 17, 18 табл. 1).
Для дальнейших исследований выбран метод обрат-но-взвешенные расстояния, который является приемлемым для всего ДВ региона.
Расчет показателя Р. осуществляется на основе интерполированного значения показателя Ь. и суточного объема осадков х, проявляющего в основном дискретный фронтальный характер. Распределение осадков на территории может быть определено по зонам влияния ГМС, которые представлены окружностями с радиусом 30 км [13,18], либо по зонам ответственности - полигонам сети Тиссена [8]. На рис. 4 приведены зоны влияния и ответственности ГМС на территории ЕАО и южной части Хабаровского края, построенные по данным 10-ти ГМС.
Проверка значений расчетных интерполированных значений комплексного показателя Р. на территории ЕАО выполнена в течение пожароопасного сезона 2010 г. В этот временной период в области возникло 126 лесных пожаров, при этом распределение пожаров по участкам с различными классами опасности растительности выглядит следующим образом: 24 пожара (I класс), 40 пожаров (П класс), 47 пожаров (Ш класс), 14 пожаров (IV класс) (строка № 1 табл. 2). На каждый день сезона построен векторный слой квартальной сети лесхозов, в центрах кварталов которых восстановлены значения Ь ,, определены суточные объемы осадков (х) и рассчитаны показатели Р.. Многомерный набор слоев распределения показателя Р. использован для определения вероятности возникновения пожара по погодным условиям по шкале Т.В. Костыриной [4]. Установлена высокая достоверность прогноза (0,86±0,06%) при вероятности 0,95. В случае использования краткосрочного прогноза на 3 суток достоверность незначительно уменьшалась от 0,82±0,07 до 0,78±0,07, но оставалась достаточно высокой (табл. 2).
Вывод
Восстановление лесопожарного показателя методом обратных взвешенных расстояний совместно с нахождением распределения суточного объема осадков по полигонам Тиссена позволяет получить хорошие результаты прогноза появления пожаров растительности на территории Дальнего Востока России.
Работа выполнена при поддержке гранта ДВО РАН № 12-1-0-09-013 и РФФИ № 12-07-31070.
ЛИТЕРАТУРА:
1. Глаголев В.А., Коган P.M. Информационная система оценки и прогноза пожарной опасности по условиям погоды (на примере Среднего Приамурья) // Вестник ШУ 2009. Т. 314, №5. С. 180-184.
2. Егармин П.А. Система детальной оценки пожарной опасности лесной территории: автореф. дис. ... канд. техн. наук / Сибирский государственный технологический университет. Красноярск, 2005. 21 с.
3. Ковин РВ. Алгоритмическое и программное обеспечение геоинформационной системы для анализа двумерных геополей: автореф. дис. ... канд. техн. наук / Томский политехнический университет. Томск, 2004.19 с.
4. Костырина Т.В. Прогнозирование пожарной опасности в лесах юга Хабаровского края: автореф. дис. ... канд. с.-х. наук / Институт леса им. Сукачева СО РАН. Красноярск, 1978. 23 с.
5. Курбатский Н.П. Прогнозирование лесных пожаров с помощью ЭВМ // Лесное хозяйство. 1976. № 7. С. 51-55.
6. Нестеров В.Г. Горимость лесов и методы ее определения. М.-Л.: Гослесбумиздат, 1949. 74 с.
7. Пискунов Н.С. Дифференциальное и интегральное исчисление. Учебник в 2-х т. СПб.: Мифрил. гл. ред. физ.-мат. лит., 1996. Т. 1. 416 с.
8. Подольская A.C., Ершов Д.В., Шуляк П.П. Применение метода оценки вероятности возникновения лесных пожаров в ИСДМ - Рослесхоз // Современные проблемы дистанционного зондирования Земли из космоса. 2011. Т. 8, № 1. С. 118-126.
9. Понамарев Е.И. Сухинин А.И. Методика картирования и среднесрочного прогнозирования пожарной опасности лесов по условиям погоды // География и природные ресурсы. 2002 . № 4. С. 112-117.
It is investigated the use of interpolation methods for the recovery
autonomous region and in the Khabarovsk territory.
Keywords: fire risk index, fire risk, interpolation.
10. Софронов M.A. Пирологическое районирование в таежной зоне. Новосибирск: Наука, 1990. 204 с.
11. Сухинин А.И. Космические методы в мониторинге лесных пожаров и их последствий // Сибирский экологический журнал. 1996. № 1. С. 85-92.
12. Чернова Л.И. Обработка геопространственной информации при цифровом моделировании топографических задач: дис. ... канд. тех. наук:Иркутск,2006. 123 с.
13. Bailey Т.С. & Gatrell А.С. Interactive spatial data analysis. Longman. 1995. 413 p.
14. Davis J.C. Statistics and data analysis in geology (2nd ed.). New York, 1986.
15. Duchon J. Splines minimizing rotation invariant seminorms in sobolev spaces, constructive theory of functions of several variables. 1976. Vol. 1. P. 85-100.
16. Flannigan M.D., Wotton B.M. A study of the interpolation methods for forest fire danger rating in Canada // Canadian Journal of Forest Research. 1989. №19(8). P. 1059-1066.
17. Nash J.E., Sutcliffe J.V River flow forecasting through conceptual models: 1 Adiscussion ofprinciples // Journal of Hydrology. 1970. Vol.10, № 3. P. 282-290.
18. Robert J. Renka. Multivariate Interpolation of Large Sets of Scattered Data // ACM Transactions on Mathematical Software. 1988. Vol. 14,N.2. P. 139-148.
19. Sirayanone S., Hardy, R.L. The Multiquadric-biharmonic Method as Used for Mineral Resources, Meteorological, and Other Applications // Journal of Applied Sciences and Computations. 1995. Vol. 1. P. 437-475
20. Smith W. H. F., Wessel P. Gridding with continuous curvature splines inTension // Geophysics. 1990. № 55. P. 293-305.
of forest fire risk index values on weather conditions in the Jewish
