СПОСОБ ОПРЕДЕЛЕНИЯ БОНИТЕТА НАСАЖДЕНИЙ Российский патент 2009 года по МПК A01G23/00 

Изобретение относится к лесному хозяйству и может найти применение при дистанционном мониторинге лесов космическими средствами на обширных площадях.

Бонитет является показателем, характеризующим качество условий произрастания леса, прежде всего его продуктивность. В России используют единую для страны бонитировочную шкалу классов бонитета [см., например, Общесоюзные нормативы для таксации лесов, справочник, стр.122, таблица 32 Бонитировочные шкалы для семенных и порослевых насаждений по М.М. Орлову, стр.123, Бонитировочные шкалы по группам древесных пород (ВНИИЛМ), аналог].

Классификационным признаком бонитировочной шкалы является средняя высота Н[м] насаждений в соответствующем возрасте. Дробность отечественной бонитировочной шкалы составляет по высоте 4 м, по запасу ≈20%, в соответствии с чем выделены следующие классы бонитета: Iа, I, II, III, IV, V, Va.

Входными параметрами отечественной бонитировочной шкалы являются высота и возраст насаждений. За столетний период составления бонитировочной шкалы накоплены материалы, позволяющие ее уточнить и расширить.

Английская, американская и немецкая бонитировочные шкалы квалификационным признаком имеют максимальный средний годичный прирост насаждений [см., например, Н.П.Анучин, Лесная таксация, учебник, 5-е издание, Москва: Лесная промышленность, стр.230-232, Бонитирование насаждений в США и Западной Европе - аналог]. В английской (американской) бонитировочной шкале разница в величине прироста, равная 20 куб. футам на 1-ом акре, принята в качестве ступени (интервала) между двумя классами производительности. В соответствии с этим для столетних сосновых насаждений установлены следующие шесть классов производительности: 160, 140, 120, 100, 80, 60 куб. футов на акр, что соответствует метрическому интервалу 1,37 м³/га в год. При этом все шесть классов шкалы производительности находятся в пропорциях:

класс бонитета: Ia I II III IV Va

прирост в год, м³/га: 10,96 9,59 8,22 6,85 5,48 4,11.

Недостатками известных аналогов являются:

- неоперативность, большая трудоемкость перечислительной таксации при установлении классов бонитета;

- трудность расчета и субъективизм определения средней высоты насаждения многоярусных насаждений.

Ближайшим аналогом к заявляемому способу является Патент RU №2277325, 2006 г. «Способ определения прироста запаса насаждений».

Способ ближайшего аналога включает:

- получение изображений лесных массивов в виде зависимости яркости I(х,у) от пространственных координат;

- преобразование функции яркости изображения в цифровую матрицу дискретных отсчетов, размерностью |m×n| элементов;

- программный расчет пространственного спектра изображения анализируемой матрицы и определение прикрепляющей точки D_cp=1/F_cp огивы насаждения анализируемого участка;

- расчет запаса насаждений по числовым характеристикам анализируемой матрицы и массовым таблицам;

- вычисление площади рельефа древесного полога (S_p) анализируемого участка с помощью программного анализа матрицы |m×n|;

- вычисление прироста запаса насаждения (Р_м) в соответствии с формулой

где F_cp1, F_cp2 - средняя частота пространственного спектра матриц изображений, полученных в начале и конце временного интервала;

S_P1>S_P2 - площади рельефа древесных пологов анализируемого участка в том же временном интервале;

D_cp - средний диаметр крон древесного полога;

F_cp - средняя частота пространственного спектра матриц изображений.

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

К недостаткам ближайшего аналога следует отнести:

- невозможность непосредственного использования способа для бонитировки насаждений;

- не все признаки продуктивности леса, измеряемые дистанционным способом (например, фотосинтезирующий индекс), используются при оценках.

Задача, решаемая заявленным способом, состоит в повышении точности и производительности дистанционной бонитировки насаждений на обширных площадях. Это осуществляется путем совместного анализа измерений коэффициента спектральной яркости (КСЯ) в зеленом и красном участках видимого диапазона спектра электромагнитных волн и полученных синхронно в тех же участках спектра двумерных изображений.

Технический результат достигается тем, что в способе определения бонитета насаждений, включающем получение изображений лесных массивов в виде цифровой матрицы |m×n| дискретных отсчетов зависимости амплитуды сигнала А(х,у) от пространственных координат, расчет пространственного спектра матрицы, определение средней частоты пространственного спектра (F_cp) и диаметра кроны среднего дерева D_cp=1/F_cp, вычисление площади рельефа древесного полога S_P анализируемого участка, расчет запаса V[м³/га] по числовым характеристикам матрицы и массовым таблицам, отслеживание изменения расчетных параметров в те же вегетационные периоды на временном интервале в несколько лет, дополнительно осуществляют измерения коэффициентов спектральной яркости (КСЯ) в зеленом (G) и красном (R) участках видимого диапазона с одновременным получением двумерных изображений в этих же участках, что позволяет рассчитать фотосинтезирующий индекс Ф=G/R как отношение средних значений КСЯ в зеленом и красном участках спектра. С использованием фотосинтезирующего индекса находят прирост продуктивности ΔV на временном интервале наблюдений

По величине прироста продуктивности δV/δt=ΔV/L, высоте среднего дерева насаждения H_cpI=H(D_cp1) и шкале бонитировочной характеристики устанавливают класс бонитета, где

L - интервал наблюдения, [лет];

Ф₁,Ф₂ - расчетные значения фотосинтезирующего индекса в начале и конце интервала наблюдений;

D_cp1, D_cp2 - диаметр кроны среднего дерева в начале и конце интервала наблюдения;

П₁/П₂ - полнота насаждения в начале и конце интервала наблюдений.

Изобретение иллюстрируется чертежами, где:

фиг.1 - зависимость коэффициентов спектральной яркости насаждений от длины волны излучения при различном проективном покрытии;

фиг.2 - сечение рельефа древесного полога вертикальной плоскостью;

фиг.3 - зависимость полноты насаждения от соотношения площади рельефа (S_p) и геометрической площади участка (S_о):

фиг.4 - огибающая пространственного спектра функции сигнала;

фиг.5 - зависимость высоты древостоя (Н_ср) от диаметра крон(D_ср);

фиг.6 - шкала бонитировочных характеристик;

фиг.7 - функциональная схема устройства, реализующая способ.

Техническая сущность способа состоит в следующем.

Как следует из способов-аналогов, основным квалификационным признаком бонитировки насаждений является продуктивность, размерность величины м³/га в год. Простейшая формула исчисления запаса перечислительным способом включает произведение трех параметров: высоты, толщины ствола модельного дерева участка и количества деревьев на участке. Чтобы с математической точки зрения постановка задачи считалась корректной, необходимо при дистанционном зондировании также измерять как минимум три независимых параметра. В заявляемом способе независимыми измеряемыми параметрами являются: фотосинтезирующий индекс Ф, диаметр кроны среднего дерева D_cp, полнота насаждений П. Их произведение (Ф×D_ср×П) отражает продуктивность (объем продуцирующей биомассы). Сама градация ступени бонитета имеет величину: 1,37 м³/га в год. Точно отследить изменение такой величины при дистанционном зондировании за счет изменения высоты и толщины отдельного дерева затруднительно. В то же время, за счет ежегодного роста ветвей, молодых побегов увеличивается объем продуцирующей биомассы, существенно изменяются такие интегральные характеристики насаждений как рельеф древесного полога, пространственный спектр изображения, фитоценометрические параметры. Чтобы вычислить прирост продуктивности, нет необходимости рассчитывать детальные таксационные характеристики. Для статистической устойчивости и точности результатов оценки следует оперировать среднестатистическими характеристиками насаждений. Рассмотрим измерение и расчет каждого из параметров, входящих в произведение (Ф×D_ср×П).

Наиболее тесная зависимость между содержанием хлорофилла в хвое, листьях (ответственного за фотосинтез) и коэффициентом спектральной яркости растительных покровов соответствует зональному отношению

[см., например, Выгодская Н.Н., Горшкова Н.Н. Теория и эксперимент в дистанционных исследованиях растительности. Ленинград: Гидрометеоиздат, 1987, стр.123-130].

Зависимость КСЯ растительных покровов в видимом диапазоне от изменения проективного покрытия иллюстрируется графиками на фиг.1 [см., например, Л.И.Чапурский «Отражательные свойства природных объектов в диапазоне 400-2500 нм», ч. I, МО СССР, 1986 г., стр.44-53, Коэффициенты спектральной яркости растительного покрова. КСЯ крон деревьев].

До 95% лучистой энергии солнца поглощается древесной растительностью в процессе фотосинтеза, и лишь небольшой максимум отражения на длине волны λ=550 нм придает растительности зеленую окраску. Чем больше объем фитомассы, тем больше разница КСЯ на длинах волн λ=550 нм и λ=670 нм (см. фиг.1).

Ежегодный прирост объема фитомассы изменяет проективное покрытие (полноту) древостоя, интегральной характеристикой которого является площадь рельефа древесного полога S_P. Изменение площади рельефа древесного полога в зависимости от полноты П древостоя иллюстрируется графиками на фиг.2. Относительную полноту древостоя рассчитывают через отношение площади рельефа S_r древесного полога к геометрической площади участка П=П(S_р/S_о). Зависимость иллюстрируется графиком на фиг.3 [см., например, Патент RU № 2294622, 2007 г. «Способ определения полноты древостоев»]. Площадь рельефа вычисляют программным методом обработки функции сигнала изображения A(x,y) по аналитическому выражению

где m,n - число строк, столбцов анализируемой матрицы изображения; х,у - текущие координаты функции сигнала A(х,у);

σ² - дисперсия сигнала матрицы |m×n| элементов.

Текст программы вычисления площади рельефа [см. Патент RU №2855357, 2005 г. «Способ определения площади рельефа»] приведен ниже в примере реализации.

Следующей интегральной характеристикой, отслеживаемой при дистанционном зондировании, является распределение крон деревьев по их диаметру D. Для чего вычисляют огибающую пространственного спектра матрицы изображения |m×n| элементов путем двумерного Фурье-преобразования функции сигнала А(х,у). На фиг.4 показаны огибающие пространственных спектров сигналов изображений двух участков насаждений, исходного и прирастающего, полученных в способе ближайшего аналога.

Вычисление спектра Фурье матрицы |m×n| элементов, функции сигнала А (х,у) является стандартной процедурой специализированного программного обеспечения, [см. Пакет программ для обработки изображений в науках о Земле «ER MAPPER Reference», Earth Resource Mapping Pty Lid, Western Australia, 6005, 1998, p.295].

За время существования лесной науки получены устойчивые статистические зависимости между диаметром кроны, высотой и площадью поперечного сечения ствола дерева (см. фиг.5) [см. аналог. Общесоюзные нормативы для таксации лесов, справочник, стр.147, табл.42, Н.П.Анучин [Таксация, учебник, 5-е издание. Лесная промышленность, 1982 г., стр.105, 296-298, таблица 41; Патент RU №2133565, 1999 г., Способ таксации насаждений, таблица 2, лист 15-16].

Главнейшие древесные породы разделены на две группы:

- первая группа: сосна, лиственница, береза, осина, ольха, H_I=7D^1,2;

- вторая группа: ель, пихта, кедр, ясень, бук, дуб, H_II=5D^1,1;

площадь сечения g=120D^0,8 (см²), а запас вычисляют по массовым таблицам, связывающим видовые числа f_i, ступень толщины (d_1,3i), площади сечений g=π d² _{1,3 i}/4 и количество деревьев в насаждении N_i данной ступени толщины из соотношения:

V=Σf_i Н_i g_i N_i или для среднего модельного дерева: V=f_cp H_сp g_cp N_cp.

На временном интервале наблюдений в несколько лет, в те же вегетационные периоды года, проводят измерения основных параметров заявляемого способа Ф, D, П и оценивают прирост: ΔV=(Ф₂ D_cp2 П₂/ Ф₁ D_cp1 П₁ -1)V.

Полученные величины Н - высоты древостоя и ΔV являются входными в калибровочную шкалу бонитета насаждений, иллюстрированную графиками фиг.6. Пересечение прямой (H=const) с функцией прироста ΔV/L однозначно определяет класс бонитета.

Пример реализации способа

Заявляемый способ может быть реализован на базе устройства, функциональная схема которого показана на рисунке 7. Устройство содержит орбитальный комплекс наблюдения 1 типа Международной космической станции (МКС) с установленной на ее борту поворотной платформой 2 (±15°, от надира) и гиперспетрометром 3 (типа Астрогон-1). Съемка запланированных участков лесных массивов 4, включение гиперспектрометром 3 и выставка поворотной платформы 2 осуществляет бортовой комплекс управления 5 (БКУ) по командам, передаваемым из Центра управления полетом 6 (ЦУП) по радиолинии управления 7. Измерительную информацию зондирования лесных массивов записывают на бортовое устройство хранения информации 8 (типа видеомагнитофона «Нива») и в сеансах видимости МКС с наземных пунктов передают по автономным каналам передачи данных 9 на наземные пункты приема информации 10 (ППИ), где осуществляется запись переданных массивов на устройство хранения информации 11 (типа видеомагнитофона «Арктур»). Информацию с ППИ передают по наземным каналам связи в центр тематической обработки 12, где осуществляют выделение кадров по служебным признакам. Скомпонованные массивы изображений лесных участков по запросам потребителей передают в Региональные центры учета лесных ресурсов, где создают их долговременный архив 13 на базе запоминающих устройств (типа FT -120). Программную обработку изображений лесных участков и расчет бонитета насаждений осуществляют на ПЭВМ 14 в стандартном наборе элементов: процессора 15, оперативного запоминающего устройства (ОЗУ) 16, винчестера 17, дисплея 18, принтера 19, клавиатуры 20. Программу расчета площади рельефа записывают на винчестер 17. Доступ к результатам расчетов и долговременному архиву осуществляется через сеть Интернет 21.

По расчетной величине площади рельефа S_p и геометрической площади S_o определяют полноту насаждения П. Из графиков фиг.2,3 исходная полнота насаждения П₁=0,5. Затем стандартной программой Фурье-преобразования вычисляют огибающую пространственного спектра. Амплитудно-частотная характеристика (АЧХ) пространственного спектра изображения анализируемого участка иллюстрируется графиками фиг.4. Средняя частота пространственного спектра исходного участка насаждения, которая делит площадь под огибающей пространственного спектра пополам, составляет F_cp=0,4, диаметр кроны среднего дерева D_cp=1/F_cp=2,5 м.

Из графиков фиг.5 для древостоя I группы этому диаметру кроны соответствует средняя высота насаждений Н_{ср I}(D_{cp 1})≈21 м. При исходной высоте модельного дерева H=21 м, полноте П=0,5, площади сечения g=120D^0,8 (см²) или g=120D^0,810^-4 м², среднем количестве деревьев на одном гектаре N_cp=SП/π(D_cp ²/4)=10000*0,5*4/3,14*

(2,5²/4)=1020 шт., видовом числе модельного дерева 0,4 [см. Анучин Н.П., стр.105] запас составит:

V=f_cpН_срg_cpN_cp=0,4*21*(252*10^-4)*1020=216(м³/га).

Соотношение Ф₁=G/R=1,7. Те же значения на временном интервале пять лет (фиг.1, 2, 3, 4, 5) составили значения: Ф₂=1,8; П₂=0,53; F_cp2=0,37; D_cp2=2,7.

Откуда ΔV=216(1,8*2,7*0,53/1,7*2,5*0,5-1)=48 (м³/га).

ΔV/5 лет=9,6 (м³/га) в год.

Пересечение прямой H₁=21 и δV/δt=9,6 (м³/га) в год однозначно соответствуют (график 6) II классу бонитета.

По сравнению с ближайшим аналогом, за счет введения фотосинтезирующего индекса, способ обладает высокой чувствительностью, что делает возможным проводить оценку бонитета насаждения на более коротком интервале времени (1-3 года).

Эффективность способа характеризуется также такими качественными показателями, как высокая производительность, оперативность, точность, документальность, масштабность.

Способ включает получение изображений лесных массивов в виде цифровой матрицы |m×n| дискретных отсчетов зависимости амплитуды сигнала А(х,у) от пространственных координат. Кроме того, проводят расчет пространственного спектра матрицы, определение средней частоты пространственного спектра (F_cp) и диаметра кроны среднего дерева (D_cp=1/F_cp). Также проводят вычисление площади рельефа древесного полога (S_p) анализируемого участка, расчет запаса V[м³/га] по числовым характеристикам матрицы и массовым таблицам, отслеживание изменения расчетных параметров в те же вегетационные периоды на временном интервале в несколько лет. В способе осуществляют синхронные измерения коэффициентов спектральной яркости (КСЯ) в зеленом (G) и красном (R) участках видимого диапазона спектра с одновременным получением двумерных изображений в этих же участках. Рассчитывают фотосинтезирующий индекс Ф=G/R как отношение средних значений КСЯ в зеленой и красной полосах. Находят прирост продуктивности δV/δt=ΔV/L на временном интервале в несколько лет,

По величине прироста продуктивности δV/δt=ΔV/L, высоте среднего дерева насаждения H_cpI=H(D_cp1) и шкале бонитировочной характеристики устанавливают класс бонитета, где

L - интервал наблюдения, [лет];

Ф₁, Ф₂ - расчетные значения фотосинтезирующего индекса в начале и конце интервала наблюдений;

D_cp1, D_cp2 - диаметр кроны среднего дерева в начале и конце интервала наблюдения;

П₁, П₂ - полнота насаждения в начале и конце интервала наблюдений. Способ позволит повысить точность и производительность определения бонитета насаждений. 7 ил.

Способ определения бонитета насаждений, включающий получение изображений лесных массивов в виде цифровой матрицы |m×n| дискретных отсчетов зависимости амплитуды сигнала А(х,у) от пространственных координат, расчет пространственного спектра матрицы, определение средней частоты пространственного спектра (F_cp) и диаметра кроны среднего дерева (D_cp=1/F_cp), вычисление площади рельефа древесного полога (S_p) анализируемого участка, расчет запаса V(м³/га) по числовым характеристикам матрицы и массовым таблицам, отслеживание изменения расчетных параметров в те же вегетационные периоды на временном интервале в несколько лет, отличающийся тем, что осуществляют синхронные измерения коэффициентов спектральной яркости (КСЯ) в зеленом (G) и красном (R) участках видимого диапазона спектра с одновременным получением двумерных изображений в этих же участках, рассчитывают фотосинтезирующий индекс Ф=G/R как отношение средних значений КСЯ в зеленой и красной полосах, находят прирост продуктивности δV/δt=ΔV/L на временном интервале в несколько лет, по величине прироста продуктивности δV/δt=ΔV/L, высоте среднего дерева насаждения H_cpI=H(D_cp1) и шкале бонитировочной характеристики устанавливают класс бонитета, где L - интервал наблюдения, (лет);
Ф₁, Ф₂- расчетные значения фотосинтезирующего индекса в начале и конце интервала наблюдений;
D_cp1, D_cp2 - диаметр кроны среднего дерева в начале и конце интервала наблюдения;
П₁, П₂ - полнота насаждения в начале и конце интервала наблюдений.