本発明は、森林のCO2フラックス測定法に関する。さらに詳述すると、本発明は、渦相関法(EC)により評価されていた森林全体のCO2吸収量の評価技術の改良に関する。
温暖化ガス削減策として京都メカニズムにおける森林吸収源を利用するには、森林吸収源におけるCO2の吸収量の評価方法を示し、実行する必要がある。温暖化抑制や気候変動における森林の影響を定量化する上では、森林全体のCO2収支の評価(フルカーボンアカウンティング)が重要である。このような森林全体のCO2収支は、CO2の交換量(フラックス)を積算することで評価できるものであり、現在、地球規模のCO2収支の評価・予測が盛んに行われている。このように評価や予測を行うにあたっては、地表面で実測されたエネルギー・物質交換量(フラックス)の検証データとしての利用価値が高いことから、フラックスサイトをネットワーク化し、広域的かつ長期間にわたってフラックスデータをグラウンドトゥルースとして利用できる環境の整備が世界的に進められている。また、現在、IPCC(Intergovernmental Panel on Climatic Change、1988年に設立された国連の組織で、地球温暖化問題についての科学的な議論を行う場)が推奨するCO2吸収量の評価手法があるが、この評価手法は森林のCO2収支のうち土壌からのCO2放出を正確に把握できないことから、森林のCO2収支を正しく評価できる新しい手法が今後は求められるであろうという面もある。
例えば、小規模な公園程度の大きさであればCO2の出入りの量に差が生じるが、数十haあるいは数百haといった大規模な森林であれば出入りの量はほぼ等しいと考えることができる。したがって、森林におけるCO2吸収量を評価するにあたり、森林が横方向に均一だと考えれば、森林に対して横方向から流入する風と横方向に流出していく風のCO2の流量(交換量)はほぼ等しいとして扱うことが可能となる(図1参照)。したがって、この考えに基づけば、横方向の流量は考慮せず、森林の上部空間において上下方向に出入りする空気量などを例えば1年間測定することにより、1年の間に森林がどの程度CO2を吸収したか、あるいは放出したかがわかる(図1参照)。より具体的には、森林の上下方向におけるCO2の交換量速度を把握し続け、最後に積算すればある期間(例えば1年)の吸収量(放出量)がわかる。この吸収ないしは放出の量、つまりはCO2の交換量を単位時間あたりの量として考えたものが「フラックス」、つまりいわば交換速度ということになる(図1参照)。
ここで、森林上におけるCO2の収支を表すNEE(Net Ecosystem Exchange:生態系純交換量)は、CO2フラックスFCを用いて次式で表される。
[数7]
NEE=FC+FS
FSは森林内の空間に一時的に存在するCO2量(林冠内貯留量)で、1日の中で貯留と放出を繰り返し、1日を通した収支は非常に小さい。したがって、1年間のNEEは、CO2フラックスの時系列1年分を積算することで求めることができる。1年間のNEE(mgCO2m-2yr-1)は光合成の同化作用によって森林生態系が1年間に吸収したCO2量を示す(図5参照)。
従来、このCO2フラックスを測定するための手法の一つとして、渦相関法(EC:Eddy-covariance method)が広く用いられている。この渦相関法は、垂直風速の変動(w’)と気温や湿度、CO2密度変動の測定値(c’)の共分散の平均値(w’c’(文中ではアッパーラインを付していないが、算術平均を表しているものとする))から乱流輸送量であるフラックスを求める、つまり、図1に示してあるように、
CO2フラックス=風の乱れ×CO2密度変動=w’c’
という式に従い、CO2フラックスを求めるというものである。この渦相関法による測定手法は、各物理量を直接測定するものであることから最も仮定が少ない測定法とされている。また、この渦相関法を利用する場合の施設としては森林上空の風の流れとCO2の密度の変動を直接的に測定するための施設、例えば森林内に建てられたタワーや当該タワーに設けられたセンサなどが利用されている(例えば、非特許文献1参照)。センサとしては、風の乱れを測定するための超音波風速計や、CO2密度変動を測定するためのCO2/H2O分析計などが用いられる。
しかしながら、従来の測定法の場合、測定精度が十分ではないとの見方がある。すなわち、フラックスは大気下層の地表面付近では大気乱流によって駆動されることから、晴天日の夜間など大気が安定した状態や、大気の対流が大きい不安定な条件では、渦相関法による熱フラックスと放射計による放射フラックスが釣り合わない現象(インバランス)が多発している。そして、このような現象(インバランス)は、主に渦相関法に起因する問題と認識されている。
そこで、本発明は、森林全体におけるCO2の吸収量を高精度に評価するための森林のCO2フラックス測定法を提供することを目的とする。
かかる目的を達成するため、本発明者は種々の検討を行った。森林全体におけるCO2の吸収量を高精度に評価するため、森林全体における時間平均ないしは空間平均のCO2フラックスを測定することを考えた。ただし、上述したように、渦相関法だと測定精度が十分ではないとの見方があり、本発明者はこの問題について検討した結果、森林上空の1点にて測定を行ういわば「点測定」型であることに理由があるのではないかと考え、さらに検討を進めた結果、パス(トランスミッタとレシーバとの設置間隔)を長くとれるシンチロメータを使ってCO2フラックスをいわば面的に広い領域にて測定するという着想を得るに至った。シンチロメータならば、トランスミッタとレシーバとを数十mあるいは100m以上離して設置し、それらの間におけるレーザ光の屈折率変化を測定することにより、少なくとも渦相関法の場合よりも広いエリアについてCO2フラックスを測定することが可能となる。
ただし、シンチロメータの場合だと、温度の乱流強度しか測定することができないために、上述したCO2フラックスを求めるための要素である「CO2密度変動」に関するデータを直接求めることができないという欠点もある。つまり、より広い領域にてフラックス測定することにより測定精度を高めることが期待できても、CO2フラックスを求めるための要素を直接求めることができなければそもそもCO2フラックス測定が不可能となってしまう。
ここで、本発明者は、パスを長くとれるシンチロメータを測定手段として利用しつつCO2フラックスを求めるためさらなる検討を加えた。そこで、測定手法の差異に対するソースエリアの影響を調べるために、ソースエリアおよびシンチロメータと渦相関センサによる消散率を比較することに想到し、さらなる検討の結果、ある知見を得るに至った。
本発明はかかる知見に基づくものであり、請求項1に記載の森林のCO2フラックス測定法は、森林上空のシンチロメータのトランスミッタから照射したレーザ光をシンチロメータのレシーバで受光する際、森林上における乱流あるいは温度変動に起因するシンチレーションに基づく光強度の変動を計測して温度変動の消散率εTを求め、さらに、森林上空の渦相関センサによって、温度に関する変動の分散またはその平方根σT、水蒸気に関する変動の分散またはその平方根σq、CO2に関する変動の分散またはその平方根σcを測定し、その後、以下の関係式
に基づき、水蒸気の消散率εq、CO2の消散率εcをそれぞれ算出し、スカラー量である各消散率εT、εq、εcからフラックスを求めることを特徴とするものである。
例えば森林の上空で空気乱流や熱変動が生じると、光の散乱が生じ、光強度が変動する。そこで、このように相関的に変動した光強度に基づき、乱流運動エネルギーの消散率、あるいは温度変動の消散率を求める(図2参照)。また、空気の乱れとCO2の密度とを掛け合わせる、あるいは空気の乱れと温度変動と掛け合わせることによって熱のフラックスを求めることができる。ここで、フラックスの生成と消散には一対の関係がある。そこで、消散率を求めることによって間接的にフラックスを算出することが可能である。
ただし、この手法においてはCO2のフラックス自体を直接的に求めることはできず、顕熱フラックス(伝熱)と風速シアー(乱流の強さ)という2つが直接的に求められるに過ぎない。そこで、この方法を応用してCO2のフラックスあるいは水蒸気のフラックスを求める手法を新たに導入する。すなわち、上述したようにフラックスの生成と消散とには一対の関係があることから、本発明においては、消散率を求めることによってフラックスを間接的に算出することとしている。この場合、渦相関センサによる計測値からシンチロメータの空間スケールに対応するフラックスを算出する点が特徴的である。
より具体的に説明すると、シンチロメータを用いればフラックスを間接的に測定することは可能だが、このシンチロメータの場合、温度と乱流の強さしかわからない。そこで、温度と水蒸気の消散率の比が、温度と水蒸気の分散の比の二乗と等しいという関係式を導入し、渦相関センサによって計測したデータから、温度に関する変動の分散またはその平方根σTと水蒸気に関する変動の分散またはその平方根σqの比を求め、さらには温度に関する変動の分散またはその平方根σTとCO2に関する変動の分散またはその平方根σcの比を求め、ここから、シンチロメータによっては直接的には測りえない水蒸気の消散率εq、CO2の消散率εcをそれぞれ算出することとしている(図3参照)。こうしてスカラー量である各消散率εT、εq、εcを求めることができれば、これらの結果からフラックス(例えばCO2フラックス)を求めることが可能となる。
なお、渦相関法において熱フラックスの過小評価が見られたことがあった。このような測定手法に対し、本発明の場合、シンチロメータによる熱フラックスは熱収支のインバランスを緩和させ、より妥当な評価結果が得られる。また、CO2フラックスは顕熱および潜熱フラックスと同じ測定原理に基づくことから、シンチロメータの適用によって、より信頼できる値が得られることになる。
ここで、本発明者は検討を重ねた結果、森林上に適用したシンチロメータから得られた乱流運動エネルギーおよびスカラー分散の消散率(フラックス算出のための支配的パラメータ)が、測定パスに直交する風向時に渦相関センサによる消散率を特徴的に上回る現象を知見するとともに、この現象の原因のひとつに、シンチロメータと渦相関センサの間でソースエリアが大きく異なることがあると考えた。請求項2に記載の発明はかかる知見に基づくもので、上記のようなCO2フラックス測定を実施する場合には、請求項2に記載のように、シンチロメータにおけるレーザ光の照射方向たるパスの方向が森林における風向と直交するようにトランスミッタおよびレシーバを配置しまたは配置を変換することが好ましい。シンチロメータによる測定領域を広げるという観点からは、上記のようにパス方向を風向に直交させれば、シンチロメータを利用してCO2フラックスを測定する場合により大きな効果を得ることが可能となる(図4参照)。例えば、風向が一定ないしは安定となっている場合にパスを直交させることとすれば測定精度が高くなる。また、特に測定パスをより長距離化した場合に、シンチロメータは広い領域においてフラックスが導出できるため、より広い面的な測定が可能となるという利点がある。本発明の場合、トランスミッタとレシーバの少なくとも一方を移動可能としておけば風向に応じてパスの方向を適宜変化させることが可能となる。
以上がシンチロメータを利用した森林のCO2フラックス測定法であるが、本発明者はさらに検討を重ねた。シンチロメータを利用すれば上記のような効果が得られる反面、レーザ光の光路が遮られると測定不能となるため、雨天や霧など悪天候時のデータが得られないという面もある。また、シンチロメータの測定パスに風向が直交する際には領域を代表する値が得られる一方で、パスに平行な風向時には渦相関法と同等の評価値となる。フラックスの積算によってCO2吸収量を得るためには、シンチロメータの空間スケールを持つフラックスを長期に安定して算出する必要もあると考えた。そこでさらに検討を重ねた結果、シンチロメータと渦相関法による温度変動の消散率の比εTSAS/εTECが、(i)センサの測定範囲の比(換言すればソースエリアの比、RA=ASAS/AEC )および (ii)大気の乱流強度(σw/u) と相関が高いことから、これら2つの要素をパラメータとして、渦相関法の測定値をシンチロメータ相当のフラックスに変換する手法に想到するに至った。
請求項3に記載の発明はこのような考えに基づき想到したものであり、森林全体のCO2吸収量を測定するためのCO2フラックス測定法において、シンチロメータを用いた場合と渦相関法による場合の温度変動の消散率の比εTSAS/εTECが、シンチロメータを用いた場合の測定範囲ASASと渦相関センサによる場合の測定範囲AEC の比RA(=ASAS/AEC )との相関、および大気の乱流強度σw/u との相関が高いことに基づき、これら2つの要素RA(=ASAS/AEC )およびσw/u をパラメータとする関係式
を、シンチロメータおよび渦相関センサを使って実際に得られた一定期間における測定結果に基づきあらかじめ用意しておくとともに、数式10中における測定範囲の比RAには一定値を与えることによりすべての風向に対するシンチロメータの広い測定範囲を適用することとし、その後、森林の上空に実際に設置した渦相関センサによって大気の乱流強度σw/uを測定し、以下の関係式
に基づき、渦相関法によって求めた顕熱フラックスHECをシンチロメータによって求めた場合に相当する顕熱フラックスHSASに補正し、尚かつ、渦相関法によって求めた潜熱フラックスλEECおよびCO2フラックスFcECに関しても、以下の関係式
に基づき、シンチロメータによって求めた場合に相当する潜熱フラックスλESASおよびCO2フラックスFcSASに補正し、シンチロメータにて測定したのと同様の測定結果を得ることを特徴としている。この場合のRAは風向によって変化するシンチロメータの空間スケール、σw/uは大気の乱流状態を表している。
本発明の場合、上述の(i)と(ii)を変数とする関数(数式10)を、実際にシンチロメータおよび渦相関センサを使った過去の測定結果に基づきいわば経験値としてあらかじめ用意しておく。その後の計測は渦相関法のみを実施し、シンチロメータを使った実測はもはや行わないが、あらかじめ用意された関係式を利用し、フラックスを補正することによってシンチロメータを使った場合と同様の測定結果を得ることができる。つまりは、当該測定対象森林にシンチロメータを実際に設置する必要はなく、あくまでも従来と同様の手順により渦相関センサのみにより実際の測定を実施して測定データを得ることができる。例えば渦相関センサのみによる測定結果は実際の値よりも過小評価される傾向があるが、本発明においては特定の関係式に基づきこの結果を補正ないしは補間し、あたかもシンチロメータを使って得たのと同様、真の値により近い測定結果を算出することができる。
請求項1のCO2フラックス測定法によれば、測定パスを長距離化できるシンチロメータを利用することにより、広い領域においてフラックスを計測することが可能となる。例えば従来の渦相関法の場合、いわば点で測定するため測定領域が極めて狭いのに対し、本発明にかかるフラックス測定法の場合だと大きい領域を測定できるため、狭小領域における変動には依存しない平均的かつ精度の高い測定データを得ることができる。これにより、森林全体のCO2吸収量を高精度に評価することが可能となる。
また、請求項2の発明によると、風向にパスを直交させるか、あるいはこれに近似した状態とすることによってより広い測定領域を確保することができる(図4参照)。これによれば、森林全体のCO2吸収量をより高精度に評価することが可能となる。
請求項3の発明によると、シンチロメータは使わずに、シンチロメータを使ったのと同様の高精度の測定結果を得ることができる。こうした場合、トランスミッタとレシーバのそれぞれを設置するためのタワーを設けなくても済むため、コストを抑えつつ、CO2吸収量を高精度に評価することが可能となる。要は、測定領域の異なる計測値から導いた補正を施すことによって、あたかもシンチロメータにて計測を実施したのと同様の結果を得ることができ、渦相関法による熱フラックスのインバランスを解消に向かわせることが可能となる。
しかも、シンチロメータは光学的な計測機器であるために雨や霧の場合に稼働率が下がるといったように天候に左右されやすいが、実際にシンチロメータを使わずに高精度な結果が得られるようにした本発明にかかる測定法によれば天候の影響を受けずに随時測定することが可能になるという利点もある。
以下、本発明の構成を図面に示す実施の形態に基づいて詳細に説明する。
図1〜図8に本発明の実施の一形態を示す。本発明にかかるCO2測定法は、森林全体におけるCO2フラックスを計測し評価するというもので、森林上空のシンチロメータ1のトランスミッタ2から照射したレーザ光をシンチロメータ1のレシーバ3で受光する際、森林上における乱流あるいは温度変動に起因するシンチレーションに基づく光強度の変動を計測して温度変動の消散率εTを求めるか、あるいは経験値としてあらかじめ用意しておき、さらに、森林上空の渦相関センサによって、温度に関する変動の分散(またはその平方根σT)、水蒸気に関する変動の分散(またはその平方根σq)、CO2に関する変動の分散(またはその平方根σc)を測定し、その後、以下の関係式
に基づき、水蒸気の消散率εq、CO2の消散率εcをそれぞれ算出し、スカラー量である各消散率εT、εq、εcからフラックスを求めるというものである。
まず、ここで本明細書中における各種用語について説明しておく。まず、一般に「フラックス」といった場合、気流の乱れによって鉛直方向に運ばれる乱れのエネルギー(運動量フラックス)や、熱・物質の単位面積あたりの移動速度(顕熱フラックス、潜熱フラックス、CO2フラックス)のことを意味する。本明細書では、単位時間当たりのCO2などの交換量のことをフラックスと呼ぶ。「渦相関法」は、超音波風速計とCO2/H2O変動計を用いて計測した水平および鉛直方向の風速と物質量の変動成分からフラックスを算定する手法である。「消散率」は、大気の乱れのエネルギーや、熱・物質の密度変動の消散を表す値である。接地境界層の乱流エネルギーについての生成と消散に関する収支式から乱流フラックスが算定される。
本実施形態では、まずシンチロメータ1によるフラックスの算出方法の概要について説明する。次に、その手法を拡張し、シンチロメータ1の空間スケールから決定した関数を用いて渦相関法の測定値を補正するというフラックス算出方法について説明する。
シンチロメータ1(本実施形態においてはSAS:Small Aperture Scintillometerを用いている)は、可視レーザーを射出するトランスミッタ2とレーザー強度を検出するレシーバ3の一対で構成されており、乱流によって引き起こされる大気の屈折率の変動を光強度変化として捉え、顕熱フラックス算出のパラメータである温度変動の消散率(εT)を算出する(図2参照)。温度変動の消散率は、渦相関センサによって計測した風速と温度の時系列からも計算することができる。フラックスの収支式にモニンオブコフの相似則を適用することで、消散率は次式によってスカラーフラックスである顕熱H(Wm−2)、潜熱λE(Wm−2)およびCO2フラックスFC(mgCO2m−2s−1)と関連付けられる。
ここで、u*は摩擦速度(ms−1)、ρは空気密度(gm−3)、kkarはカルマン定数(=0.4)、zは測定高さ(m)、cpは空気の定圧比熱(Jg−1K−1)およびλは水の蒸発潜熱(Jg−1)を表す。εq、εcはそれぞれ、比湿(q kgkg−1)変動およびCO2密度(mgCO2m−3)変動の消散率である。φε(ζ)は無次元化された消散率で大気安定度を表すパラメータ
[数19]
ζ=z/Lmo(Lmo:モニンオブコフ長)
の関数として与えられる。
続いて、シンチロメータ1によって決定したu*およびεTと、同時に測定したqおよびCO2密度変動の分散から、εq,εcを求めて潜熱フラックスおよびCO2フラックスを算出する手法の内容について説明する。
シンチロメータ1によるフラックス算定は消散法に基づいて行われる。乱れは生成消散によって成り立っており、渦相関法がスカラー変動や乱流運動エネルギー(気流の乱れのエネルギー:Turbulent Kinetic Energy=TKE)生成の直接測定に相当するのに対し、消散法はスカラー変動やTKEの消散率(エネルギーの消散を表す正の値)からそれらの生成を計算する手法である。消散法は船の動揺など低周波ノイズの影響を受けにくいことから海上でのフラックス測定などに応用されてきた。消散法によるフラックス算出では、まずTKEとスカラーの変動の消散率を求め、次いでフラックスと消散率との収支式からフラックスを算定する。シンチロメータ1はこのうちTKEとスカラー変動の消散率を測定する装置である。
消散法の原理について以下に説明する。水平一様で定常な接地境界層における乱流運動エネルギーTKE(et)およびスカラー変動の収支方程式は以下のように書かれる。
ここでsは気温T、比湿q、やCO2などのスカラーを表し、pは気圧で
[数22]
θ’v=θ(1+0.61q)
は仮温位(ただしθは温位)、etはTKE(et=(u2+v2+w2)/2)を表す。ここで、u,v,wは平均風向、横風方向、鉛直方向の風速成分である。また、この文中ではアッパーラインを付していないが、数式中のバーは平均値、プライムは平均値からの変動を表す。数式20の各項はそれぞれ生成項、浮力項、輸送項、圧力項、そして消散率である。また、数式21の各項は、生成項、輸送項、そして消散率である。両式第1項のu’w’およびw’s’(ただし本文中ではアッパーラインを省略)が求めるべき運動量およびスカラーのフラックスに相当する。接地境界層においては両式とも輸送項および圧力項を無視できると仮定して、数式23で定義される摩擦速度u*、
Monin-Obukov長(Lmo)およびKarman定数(kkar 0.4)を用いて、それぞれを
u* 3/kkar(z−d0) と kkar(z−d0)u*/(w’s’)2 によって無次元化すると、数式20、数式21の各項は測定高さz(m)において、
と表される。φεとφεsはTKEおよびスカラーの無次元消散率を、φmとφsはTKEおよびスカラーの無次元した鉛直勾配を表す。d0 はゼロ面変位で、植生地上での風速プロファイルの対数分布を仮定した際の、高さを測る基準面を上方に修正するためのパラメータである。群落高×0.7という値がよく使われる。さらに、ここでTKEおよびスカラー変動の消散率は次のように定義される(Brutsaert,1982)。
水平一様な接地境界層ではこれらの無次元消散率は、大気安定度を示すパラメータζ=(z−d0)/Lmo(安定時ζ>0−、不安定時ζ<0)の関数として表されることを数式24、数式25は示している。ここで用いられるMonin-Obukov長Lmoは、以下の式で定義される。
ρは空気密度、Tvは仮温度であり、H、λEは地表面における顕熱および潜熱フラックスである。数式26と数式27のs=θおよびs=qに関する式をそれぞれu*、H、λEについて解くと、
で、数式30〜数式32は、無次元消散率の関数形が既知であれば、消散率を与えることで連立方程式として解くことができ、フラックスu*、H、λEが求められる。大気安定度パラメータζによる無次元関数表現は、接地境界層におけるMonin-Obukov similarity(MOS)が成立することに基づいており、物理量の鉛直プロファイルがζで一意に表現できることを示している。例えば数式24、数式25中のφm≡kkarz/u*∂u/∂z、
φs≡kkarzu*/w’s’∂s/∂z をζの関数で表すのにBusinger-Deyer式がよく用いられる。また、無次元消散率の関数形も観測結果を元にいくつか提案されている。例えば、Kader(1992)による経験式は以下である。
温位と比湿の大気安定度に対する応答が相似であれば φεc=φε0=φεq が成立し、数式30、数式31、数式32は繰り返し計算によって解くことが出来る。更にスカラーの無次元消散率の相似性の成立を仮定した場合は、CO2の消散率から、数式31、数式32と同様にしてフラックスを算出することが出来る。
次に、渦相関センサによる消散率の算定について説明する。シンチロメータ1との比較に用いるために、渦相関センサによる計測値から消散率を算定する方法について述べる。消散率の算定には、スペクトル(エネルギーや物質量変動の周波数表現)の慣性小領域における乱流特性を利用して間接的に求める慣性消散法(inertial dissipation method)を用いる。ただし、フーリエ変換の結果は一般に乱れが大きく、消散率の推定誤差要因となりうるため、ここではスペクトルの実空間表現である構造関数にもとづいた消散率の推定法を用いる。TKEおよびスカラー変動のスペクトルに対応する2次の構造関数Dab(r)は次式で定義される。
ここでa,bはそれぞれu,θ,qなどの風速やスカラーを表し、rは実空間における距離である。時系列のデータ間の時間差をτとすると平均風速を用いて、τu=rの関係にある。Kolmogorovの理論によればrが慣性小領域にあるとき2次の構造関数と消散率は以下の関係にある。
ここでεとεsはそれぞれTKEおよびスカラー変動の消散率を表す。αuuはKolmogorov定数(およそ0.55)で、αssはObukhov-Corrsin定数(およそ0.8)と呼ばれる普遍的な定数である。比例定数Cs 2は構造関数定数で、次のaとbの組み合わせ[T:T](気温)、[q:q](比湿)および[T:q]はそれぞれ構造関数定数CT 2、Cq 2およびCTqに対応する。数式36、数式37に示すように2次の構造関数がr2/3に比例していることから、慣性小領域に相当するrの範囲で構造関数にr2/3の曲線をフィットさせることで数式36、数式37から消散率を求め、消散法によるフラックス計算へと進む(図11参照)。超音波風速計のパス長dsより小さい乱流変動は測定されず、測定高さよりも大きい渦はその一部分しか測られないことから、十分に乱流変動計測が保障されるフィット範囲として 2ds<r<(z−d0)/2 を選択する。
次に、シンチロメータ1を用いた消散率の算定方法について述べる。SLS40A(Scintec AG,Germany)は、市販のdisplaced-beam small aperture scintillometer で、波長0.67μmのHeNeレーザーを偏光面が直交する2本のビームに分離し、2光路のシンチレーションからTKEおよび気温変動の消散率を求め、フラックスを算出する装置である。SLSシリーズは単一光源を用いることで動作および計測の安定性と経済性を両立させた製品である。乱流による屈折率変動は大気中を伝搬する単色光の光強度のlog分散として計測することが出来る。独立した2光源と受光機がある場合の光強度のlog共分散B12は次式で表される(Thiemann,1992)。
ここで、xは、パス長Lの光軸上での位置;K=2π/λは光の波数;k乱流の波数;Cn 2は屈折率変動の構造関数定数、l0内部スケール(乱流渦が徐々に小さな渦に崩壊していく際、これ以上小さくなると分子粘性により渦が維持できなくなり熱として分解していくスケールの目安)である。J0およびJ1は0次と1次の第1種ベッセル関数で、ビーム間の距離dおよび受光機の直径Dによって定まる感度分布を表す。検出器単体による分散を表すB1とB2は数式38でd=0とすることで得られる。数式38中のφnは屈折率変動の次元スペクトルで、
であり、f(kl0)は、Hill and Clifford(1978)による屈折率変動の3次元スペクトルモデルである。数式38は、σx 2<0.3の低散乱条件(Rytov approximation)が前提である。式中の被積分関数はシンチロメータ1の感度分布を表す荷重関数で、パスの中央で最大値となる釣鐘型の関数形を取る。数式39を数式38に代入し、xおよびkについて積分することで次式が得られる。
fBはl0,dおよびDの増加に対するB12の減衰を表す関数である。相関係数r12=B12/(B1B2)1/2 =B12/B1=B12/B2 はl0,dおよびDの関数となり、dとDを定数と考えると、相関係数r12から内部スケールl0が求められる。l0がわかれば、B1(あるいはB2)からCn 2が求められる。Scintecのソフトウェアでは、l0に対するrおよびfBをパス長Lに対する数表として用意し、B12とB1およびB2 からリアルタイムでl0とCn 2を計算する。得られたl0から乱流運動エネルギー(TKE)の消散率は次の関係によって求められる。
ここでν(m2s-1)は空気の動粘性係数である。屈折率の変動に最も寄与するのは気温T(K)であり、ついで比湿q(kgkg-1)の変動である。CO2などのガスは大気中に微量しか存在しないため屈折率変動には影響しない。屈折率変動の構造関数定数Cn 2は気温変動の構造関数定数CT 2(K2m−2/3)、比湿変動Cq 2(m−2/3)および気温と比湿T−qの結合構造関数定数CTq(Km−2/3)によって次のように表される(Hill et al.,1980)。
定数ATとAqは光の波長、気圧、気温および比湿の影響を受ける。本実施形態で使用するシンチロメータ1で利用する波長λ=0.67μmにおいては、Andreas(1989)による以下の式を用いる。
ここでP(hPa)は気圧である。光源に単一波長を用いるシンチロメータ1ではCn 2の値が一つしか得られないため数式42からCT 2およびCq 2を一意に求めることが出来ない。この問題に対し、Thiermannは比湿がCn 2に与える影響が気温に比べて1オーダー小さいためこれを無視できるとし、数式42の右辺第1項のみからCT 2を求めた。本実施形態で用いるシンチロメータ1に付属のフラックス計算ソフトウェアはこの計算方法を採用しており、運動量フラックスおよび顕熱フラックスのみが計算される。
また、本実施形態では水蒸気やCO2密度などのスカラーについてもシンチロメータ1に対応する値を得るための手法を適用する。
まず水蒸気補正のための既存手法であるが、数式42から明らかなように、シンチロメータ1は水蒸気の変動に対しても感度を持ち、水蒸気フラックスが顕熱フラックスに比べて非常に大きい条件ではCT 2が過大評価されてしまうことから水蒸気に対する補正の必要性が指摘されている(Green and Hayashi,1998)。水蒸気影響の古典的な補正法はWesely(1976)によって提唱されたボーエン比(顕熱と潜熱の比:H/λE)補正である。ボーエン比補正にはw’T’とw’q’が必要なため、渦相関システムなどのセンサによる測定値を用いる。Moene(2003)はボーエン比補正について精査し、気温と比湿の標準偏差σT,σqの比を用いる手法が、フラックスの比を用いた手法に対して、測定項目が少ないため誤差が小さいことを示している。すなわち
とすることで、数式42は、Cn 2とCT 2に関して次のように変形できる。
Tとqの平均値と分散、そして相関係数RTqを他のセンサで計測することで、水蒸気の影響をキャンセルしたCT 2を計算することが出来る。
次にスカラーフラックスへの拡張であるが、水蒸気補正の際に導入した水蒸気と気温の構造関数定数比と分散比の相似性が成立するなら数式44にしたがって、シンチロメータ1によるCT 2と、渦相関センサによるσT,σqからCq 2を算出できると考えられる。同様にこの相似性をCO2などのガスフラックスまで拡張すると以下が成り立つ。
これにより、渦相関センサによって得られた気温、水蒸気密度およびCO2密度の分散とシンチロメータ1によって得られた気温変動の構造関数定数から水蒸気およびCO2変動の構造関数定数CTq 2,Cc 2が求められる。この算出法は純粋なシンチロメータ1による導出ではなく、渦相関センサとの組み合わせによることを留意する必要がある。しかし、現時点でシンチロメータ1と同様のスケールでの変動を測定する手法は存在しない。数式37にqおよびcに関する構造関数定数を代入して消散率εqとεcが求められ、ε,εTを加えて消散法によるフラックス計算へと進む。
続いて、渦相関法の補正法の誘導について説明する。
以下では、シンチロメータ1によるフラックスおよび消散率にはSASの添え字を、ECに相当する値にはECの添え字を添えて区別することとする。ここで、数式16〜数式18において各フラックスの比をとると数式47〜数式49となり、シンチロメータ1と渦相関センサによる消散率の比(εiSAS/εiEC)と渦相関法によるフラックスより、シンチロメータ1の空間スケールに相当するフラックスが得られる。なお、εiにおける添字iは、何の消散率を表すかによってT,q,cのいずれかの添字が該当することを意味している。
ここで、εTSAS/εTECは、乱流強度(水平風速で標準化した垂直風偏差:σw/u)および風向に応答して変化することが分かっている。また、シンチロメータ1のパスに対する風向の変化はシンチロメータ1の空間スケールに影響する。つまり、風向の変化、すなわちパスに対する風の迎え角が変化することで、空間スケール(測定する範囲)が変化する。ただし、このほかに風速や大気の安定度なども影響するが、空間スケールの変化に最も大きく寄与するのが風向と考えられる。そこで、εiSAS/εiEC を乱流強度とシンチロメータ1の空間スケールを変数とする関数で表すことで、渦相関センサによる計測値からシンチロメータ1の空間スケールに対応するフラックスを算出することが可能となる。
シンチロメータ1の空間スケールに対応するフラックス算出に関しては以下のとおりである。まず、フラックス測定センサによって検出される乱流信号は風上の領域を起源とし、その広さ(ソースエリア)は風速や大気安定度、そしてセンサの測定パス長に対して変化することになる。シンチロメータ1と渦相関の空間スケールの違いを定量的に解析するために、Kormann and Meixner(2001) のフットプリントモデルを応用してシンチロメータ1のソースエリアを算出すると以下のとおりとなる。すなわち、フットプリントF(x)は風上方向にx(m)離れた位置のフラックスへの寄与を表した関数で、風速u(ms−1)、摩擦速度u*(ms−1)測定高さz(m)、大気安定度パラメータζを変数とする。さらに、De Bruin et al.(2002) の計算例に従い、シンチロメータ1のソースエリア(ASAS)を、F(x)とセンサのパス長で表される横風方向y(m)の感度関数D(x,y)との積が描く等値線に囲まれた面積から求める(図6参照)。ここで、D(x,y)はシンチロメータ1が信号を検出する感度がパスの中央部で最も高いことを表す関数で、本実施形態では式形の簡単さから、Wang et al.(1978) が示した感度関数D(x,y)
[数50]
D(x,y)=(y/Lp(1−y/Lp))5/6
を標準化して用いる。ただし、Lpはシンチロメータ1のパス長である。また、風向によって変化するシンチロメータ1のソースエリアは、パスに直交する方角からの風向角(θ)を用いて次のように表すこととする。
シンチロメータ1および渦相関センサのソースエリアの比RA(=ASAS/AEC )は、風向がシンチロメータ1のパスに平行な場合(θ=90°、270°)に1となり、直交風向時に増加する(図7参照)。
次に、領域平均フラックスを算出する。ここでは、まず上記で求めたソースエリア比RAおよび乱流強度σw/uを変数とする関数f(σw/u,RA)の形を定め、一定期間に得られたεTSAS/εTEC の実測値に対する回帰によって、関数の係数を決定する。
次に、妥当と考えられるフラックス評価値を導くRAを求め、その関係が全ての風向に対して成立すると仮定して、数式47より顕熱フラックスHSASを求める。こうして求めたフラックスを領域平均フラックスと定義する。領域平均フラックスは全風向に対してシンチロメータ1の空間スケールを持ち、渦相関センサによる計測値のみでフラックス評価を行うことができるため、シンチロメータ1が測定不能な気象条件のフラックスも評価できる。さらに、水蒸気とCO2が熱と同じ原理で輸送されると考えられることから、
と仮定すると、潜熱とCO2についても同様に領域平均フラックスが計算できる(図8参照)。
なお、上述の実施形態は本発明の好適な実施の一例ではあるがこれに限定されるものではなく、本発明の要旨を逸脱しない範囲において種々変形実施可能である。
中部地方の浅間東麓の54年生の落葉広葉樹林においてCO2を含むフラックス測定を、2002年から2004年を通して行った。ここでは、シンチロメータ1の空間スケールから決定した関数を用いて渦相関法の測定値を補正することにより、領域平均フラックスを算出するという「領域平均フラックスの算出」を実施し、算出した領域平均フラックスによって熱収支が正しく評価されていることを確認した。また、光合成量の計測結果とCO2フラックスから求めた群落光合成量とを比較した。優先種であるダケカンバ個葉の光合成特性と葉の分布をもとに算出した森林の光合成ポテンシャルに対し、渦相関法が約CO2程度過小評価する一方で、領域平均フラックスによる算出値は良好に一致したことから、CO2交換量においても領域平均フラックスによって評価値の信頼性が向上することが分かった。さらに、森林のCO2吸収量の評価を実施した結果、気象条件や樹種が類似した場所で得られた既往の渦相関法による評価結果とほぼ同等であった。以上につき、以下に実施例1として説明する。
ここでは、高さ28mの鋼製足場材組上げ式のタワーを東西に86m離して配置し、東側の塔頂で渦相関法を実施するとともに領域平均フラックス測定のために、2基のタワー間をシンチロメータ1の測定パスとした。渦相関法では顕熱、潜熱、CO2フラックスを測定した。フラックス算出には10Hzサンプリング時系列の30分ブロック平均を用いた。NEEの算出には降雨などの影響を受けにくく、測定値のドリフトが少ない閉光路(closed-path)式CO2分析計(LI-7000, Licor USA)によって計測したCO2密度変動を用いた。内径4mmのデカボンチューブ40mを経由してタワー基部の観測舎に塔頂の大気サンプルを9.0Lmin−1で吸引し、その一部を2.0Lmin−1でCO2分析計に導入した。サンプルエアが塔頂から分析器までに到達する時間遅れは、9秒程度であった。サンプルエアーの吸引によって信号の高周波成分が劣化してしまう問題は、乱流のスケールが大きな森林上ではほとんど考慮しなくて良いことが報告されている。当試験林においても、応答の速い開光路(open9-path)式のCO2分析計と比較して、closed-path式のシステムが問題ない精度でCO2フラックスを計測することを確認した。領域平均フラックスを計算するための関数f(σw/u,RA)を定めるために、シンチロメータ1の測定値が得られた2002年6月から11月までの晴天日の測定値を用いて回帰計算を行った。
次に、欠測データの補間について説明する。フラックスからNEE(Net Ecosystem Exchange:生態系純交換量)を算出するためには、計測器の故障、メインテナンス、停電などさまざまな理由により発生する欠測を補い、連続した時系列を得る必要がある。そのために、1.5時間以内の短期欠測は線形補完を行い、それ以外はCO2フラックスと相関の高い気象要素を用い、あらかじめ求めておいた気温‐呼吸量および光‐光合成の関係から補間を行った(図9参照)。測定値が得られない場合はAMeDASデータから当試験林の日平均日射量および気温を推定した。
また、フラックス測定値との比較を行うために、葉面積、個葉の光合成量の測定を行った。タワーを中心とした60m(N−S)×140m(E−W)の範囲を10m×10mの区画に分け、すべての生態データ取得の基準とした。ここでは、葉面積の季節変動を明らかにするために、開口部1m×1mの落葉採取装置(リタートラップ)を林床から1.3mの高さに千鳥状に19個配置し、9月中旬から11月上旬の落葉期に落葉・落枝を定期的に採取した。また、強風などによる落葉期以外の落葉・落枝を考慮して着葉期の6月から9月上旬まで、平均的な落葉量を示した3個所のリタートラップを用いて落葉量を調査した。あらかじめ求めておいた落葉の乾燥重量と面積の1次回帰式から、採取した落葉サンプルの葉面積を求め、全葉面積から調査区域のLAI(Leaf Area Index:波面積指数)を求めた。光合成量の測定は、渦相関法を実施したタワーの北側に隣接したダケカンバ(Betula ermanii Cham.)の高度9m付近にある任意の個葉を対象に、携帯型光合成蒸散量測定装置(LI-6400, Licor USA)を用いて行った。ダケカンバは、事前の毎木調査によって、当試験林を構成する樹木の地上部(枝+幹)炭素量の44.3%を占める優先種であることが分かっている。当装置は小型のリーフチャンバを備え、チャンバ内の気温、湿度およびCO2濃度を外部環境に追従させて制御できる。ダケカンバ葉が十分に光飽和する光強度から7段階(2000, 1500, 1000, 500, 200, 100, 0μmolm-2s-1)に変化させ、最も光合成能力が高いと思われる7月から9月期間で光強度に対する光合成量の反応(光−光合成反応)を調べた。
以上の結果、以下のような結果が得られた。まず、2002年から2004年にかけての年毎のフラックスデータ欠測率はそれぞれ18%、5%、3%であった。2002年に欠測率が高いのは落雷による停電の影響で、そのほかにはメインテナンスや悪天候、機器の故障などの原因があった。ここで、FLUXNET(Falge et al.,2002)や、Asia Flux(Asia Flux運営委員会、2003) などでは信頼性の高いフラックス測定のために、乱流データのクオリティコントロール(QC)が提言されている。QCでは乱流の成立性が乏しい状況での渦相関法は信頼性が低いとして、摩擦速度u*や時系列データの定常性に基準を設け、基準に満たないデータを棄却して補間の対象とする。このため、標準的な森林フラックスサイトの欠測率は30%程度と言われている。しかし、本実施例では、渦相関法の信頼性が低いとされる乱流の乏しい状態あるいは対流の大きな状態における計測値も含めてシンチロメータ1のデータと比較し、その結果に基づいて領域平均フラックスを算出するため、QCは適用せずすべてのデータを用いた。
続いて、補正関数の決定を行った。ここでは、シンチロメータ1と渦相関センサによる温度変動の消散率比εTSAS/εTEC は、機械的な乱流強度(σw/u)が大きな中立時では、ECセンサによる消散率がSASを上回り、(σw/u)が小さな安定時および不安定時ではSASの消散率がECより大きくなった(図10(a)、図10(b)参照)。消散率とフラックスの絶対値は数式47のように比例関係にあるため、この結果は、対流が卓越する条件や成層した大気条件では渦相関センサがフラックスを過小評価するというKanda et al.(2004)の結果と一致した。また、シンチロメータ1のソースエリアが渦相関センサに比べて大きくなる(RAが大きくなる)につれてεTSAS/εTEC の分散の上端が大きくなる傾向が見られた(図10(c)、図10(d)参照)。センサの空間スケールが大きくなることで、より広い領域を起源とする信号の平均値を検出できるため、空間的スケールの大きな大気の乱流現象を捉えることができ、結果的に上記のσw/uに対するεTSAS/εTEC の関係が強調されると考えられた。
そこで、図10(a)および図10(b)のf(σw/u,RA)が大きい領域の漸近値をオフセットと考え、指数減衰関数(数式54、数式55)によって上記の関係を表現することとした(図10の右欄参照)。ここでは、数式中の係数a1,a2,a3を昼間のおよび夜間について最小二乗法によって決定した。夜間と昼間で係数が異なるのは、ソースエリアの算出に用いたフットプリントモデルが夜間と昼間で異なる関数で表されるためである。
[数54]
f(σw/u,RA)=y0+a1exp[−(σw/u−x0)/a2]
[数55]
x0=(RA−1)/a3
関数f(σw/u,RA)と消散率比εTSAS/εTEC の比較を図12に示す。昼間(図12(a))に比べて夜間(図12(b))の分散が大きいのは、渦相関センサによる温度変動の消散率εTEC の算出精度が乱流状態の乏しい状況で低下したことが原因と考えられた。
次に、領域平均フラックスが得られると考えられる最適なRAを計算するために、年間の熱フラックスの収支とRAの関係を調べた(図13)。昼間の不安定時および夜間安定時共に、RAの増加にしたがってインバランス率I
[数56]
I=(H+λE−Rn)/Rn
が0に近づき、熱収支が釣り合う効果が見られた。特に、夜間安定時はその効果が強かった。潜熱と顕熱の他に林冠内の大気や植物体に蓄えられる熱や地中伝熱( 林冠内貯熱)が熱収支に関与するが、その値は大きくても10数Wm−2であることと、RAが20を超える測定値がほとんど無く、外挿による誤差の拡大を防ぐ観点から、インバランス率=−0.3を示すRAを求めた(夜間10.8、昼間23.1)。これによってf(σw/u,RA)がσw/uの変数として表され、渦相関センサの測定値に基づいて領域平均フラックスが得られた。
ここで、領域平均フラックスの特徴について記す。領域平均フラックスは、σw/uが小さな領域では渦相関法によるフラックスより大きな値を示し、σw/uが大きな領域では渦相関法より小さな値を示す。その結果、夜間の放出および昼間の吸収量ピークが大きく見積もられ、日変動の振幅が拡大する。図14に熱およびCO2フラックスの測定例を示す。熱フラックスでは領域平均フラックスによる熱収支の改善が見られる一方で、大気が中立な条件では、領域平均フラックスが渦相関法を下回る状況も見られた(図14のA参照)。また、熱収支の改善効果が小さい状況も見られた。林冠内貯熱の寄与が大きい時間帯が存在する可能性もあると考えられた。FCについても同様に放出・吸収ピークが拡大する傾向が見られ(図14のB参照)、FCに基づいて計算される群落光合成量や、群落呼吸量も変化することになる。
さらに、領域平均フラックスと土壌呼吸量との比較を行った。土壌中の動物および根の呼吸と微生物による有機物の分解によってCO2が放出される土壌呼吸と、地上部の植物体による呼吸が合算して森林からのCO2放出となり、その放出量を群落呼吸量と呼ぶ。森林上で観測したCO2からは、群落呼吸量に相当する情報が得られる。一方、土壌呼吸量は、チャンバ法などによって実測が可能である。そこで、CO2フラックスの妥当性を検証するために、群落呼吸量と土壌呼吸量を比較した。
まず、群落光合成量Ag、群落呼吸量ReとNEEの間には
[数57]
NEE=−Ag+Re
の関係がある。夜間はAg=0なので、NEEが群落呼吸量となる。群落呼吸量は温度の相関が高く、温度に対する群落呼吸量のパラメタリゼーションには数式58が用いられる(Lloyd and Taylor,1969)。
ここで、T(K)は温度、T
0(K)は任意に決定する基準温度で、R
0はT
0(K)におけるCO
2フラックスを表す。また、Q
10は温度係数で、基準温度におけるCO
2フラックスに対する、基準温度から10K上昇したときのCO
2フラックスの比を表す。気温と群落呼吸量の関係(2004年)によると(図15参照)、見易さのため σ
w/u>0.3以上のデータを示した。そこで、当該森林の群落呼吸量を解析するために、SASおよびECから求められる夜間のCO
2フラックス(F
C)のデータと温度の関係から、パラメータQ
10,R
0を決定した(表1参照)。ここでは、冬季間に積雪することを考慮して、Tには地温ではなく地上高さ3mの位置の気温を用いた。なお、群落呼吸量R
eは、土壌呼吸量(土壌中の生物による呼吸と根の呼吸)および植物地上部の呼吸量を含む。
渦相関法にもとづく群落呼吸量のパラメータQ
10と領域平均フラックスによるQ
10は年間の差異において異なる傾向を示したが、いずれも平均気温との間には明確な傾向は見られなかった。次に、夜間の温度−呼吸量の関係が昼間でも同様に成り立つと仮定して、昼間の気温と数式58から群落呼吸量を推定した。こうして求めた温度と群落呼吸量の関係と後述する群落光合成量とを用いて、数式57から欠測の補間を行った。以前、土壌呼吸チャンバを用いて評価した当試験林の2002年における年間の土壌呼吸量は6.6tCha
−1yr
−1であったのに対し、渦相関法と領域平均フラックスにもとづく群落呼吸量はそれぞれ、7.4tCha
−1yr
−1と7.3tCha
−1yr
−1で(表2参照)、フラックス測定手法間に明確な相違は見られなかった。一方、2003年および2004年では、渦相関法の3.9tCha
−1yr
−1および7.2tCha
−1yr
−1に対して、領域平均フラックスは6.7tCha
−1yr
−1および8.9tCha
−1yr
−1と大きな値を示した。群落呼吸量は、土壌呼吸量に植物地上部の呼吸量が加わるため、フラックスによる評価値がチャンバー法を上回るのは妥当な結果といえる。ただし、温度に対する群落呼吸量の分散が大きく(図15参照)、カーブフィッティングの際の誤差が大きいこととなった。領域平均フラックスから算出した群落呼吸量は土壌呼吸チャンバの測定結果と一致し、領域平均フラックスによって正しく群落呼吸量が評価されることが確かめられた。
続いて、領域平均フラックスと群落光合成量との比較を行った。本実施例では、CO2フラックスの妥当性を検証するために、群落光合成量との比較を行った。CO2フラックスによって評価される群落光合成量は、森林全体の光合成量を表す。一方、葉の光合成量はリーフチャンバによって実測が可能である。しかし、森林は複雑な立体構造を持ち、林床面積あたりにその数倍の面積の葉が存在する。このため、個葉において計測された光合成量を群落光合成量と比較するためには個葉の計測結果を群落に拡張する必要がある。群落光合成量は、群落最上端に位置する葉の最大光合成量と葉面積から評価することができる。
まず第一に、群落光合成量のラメタリゼーションを行った。NEEと群落呼吸量の差は同化作用による群落光合成量(数式59中のAg)を示す。群落光合成量(μmolm-2s-1)Pmは最大光合成量(μmolm-2s-1)を表し、αは初期光合成効率を表す。光合成の単位にはCO2の分子量を乗じた密度表現(mgCO2m-2s-1)も用いられる。数式57に従って入射光量と数式58に昼間の気温を代入して求めた呼吸量を差し引いたFCのプロットに数式59をフィットしてパラメータを決定した。こうして求めたパラメータと入射光量および気温からAgとReを算出し、欠測の補間を行った。
次に、個葉の光合成量を群落光合成量Pcanに換算するための葉面積の算出結果を示す。落葉量から積算した葉面積指数(LAI(Leaf Area Index:波面積指数))は、2002年から2004年までそれぞれ5.6,6.3,6.8で、増加傾向を示した。
また、図16に示した葉面積指数の変化を見ると、展葉が完了した6月から、落葉が開始する9月までの期間も落葉が継続している。当試験林の優先種であるダケカンバは、春の展葉の後に2度目の伸張をすることが認められた。6月から9月の期間中も2次伸張による展葉が継続している場合、その間の落葉量をLAIに算入するとLAIを過大評価することになる。林床面積あたりに存在する実際の葉面積に近い値として、落葉開始直前の9月上旬の値を用いると、それぞれ、5.4,5.6,5.8となり、年間の格差は小さくなったが大きさの順位は変わらなかった。全落葉重量に占めるダケカンバの落葉重量は38.2%,42.1%,42.7%とほぼ41%程度であり、地上部炭素量の割合44.3%と近い値を示した。葉面積換算では29.1%,31.9%,32.5%とほぼ31%程度であった。
続いて、群落光合成量の算出を行った。領域平均フラックスによるPmは、7月および8月に最大値を示し、その値は1.7,1.9,1.8(mgCO2m-2s-1)であった。一方、リーフチャンバによって計測したダケカンバ個葉の最大光合成量の平均値はそれぞれ0.8±0.2,0.9±0.1,0.9±0.1(mgCO2m-2s-1)であった。リーフチャンバで得られたダケカンバ個葉の最大光合成量を群落最上端の最大光合成量とみなすと、Beerの法則に従い、積算葉面積Fの位置で葉群に吸収された光量Iinkは以下で表される
ここで、I0は群落上端の入射光、Kは葉の角度などで決まる群落内の吸光定数で、0.5とする。aは葉の透過率’(=0.14)である。次に、入射光Iinkに対する光合成量Pは数式59のPARをIinkとすることで算出できる。図17に、積算葉面積に対する葉群の吸光量Iinkおよび数式59にIinkを代入して求めた光合成量Pを示す。そして、PをFで積分することにより、群落全体の光合成量Pcanを計算することができる。
群落下層の植物は群落上端に比べて光合成能力が低い(最大光合成量が小さい)のが一般的で、実際にはさらに群落光合成量が小さくなることが考えられた。また、このモデルでは水分不足による気孔閉塞がもたらす昼寝現象などは考慮していない。したがって、本計算結果は群落光合成のポテンシャルと考えるのが適当である。植物体内の水分状態が気孔の開閉に影響したり、乱流の状態によってCO2の交換量が変化するなど、フラックスの瞬時値はほぼ群落光合成のポテンシャル以下に分布すると考えられる。図18に示すように、領域平均フラックスによるAgの分布の上端とPcanがほぼ重なる。このことは、領域平均フラックスが当該森林の光合成特性をよく表しているといえる。
また、図19に当試験林で観測される最大光合成有効放射量2000μmolm-2s-1における群落光合成のポテンシャルとフラックスから導いた群落光合成量を示した。渦相関法の過小評価はCO2フラックスにおいても明らかで、領域平均フラックスによる最大光合成量は、ECによる最大光合成量の約1.5倍と、理論値である群落光合成のポテンシャルに近づいた。この結果、領域平均フラックスによって森林の光合成が正しく評価されることが示された。
続いて、領域平均フラックスの妥当性評価を行った。本発明者は、以前、シンチロメータ1の適用による熱収支式の整合をフラックス測定の妥当性の根拠とし、同じ原理で測定されるCO2についても、フラックスが妥当に測定されると判断したが、本実施例においては更に領域平均フラックスから導いた群落呼吸量および群落光合成量と、群落呼吸量および個葉の光合成特性に基づく群落光合成量を比較した結果、領域平均フラックスによる算出値は、土壌呼吸量および光合成量の実測値とよく一致した。これにより、領域平均フラックスによるCO2評価の妥当性が示された。また、熱フラックスの収支を判断材料として、フラックス測定の妥当性を評価する方法論が支持される結果となった。以上のような本実施例における結果は、当試験地特有の現象を捉えたものであるが、既存のフラックス測定サイトおよび新規測定サイトにおいても、渦相関法とシンチロメータ1と組み合わせた観測によってサイト固有の補正関数を求め、領域平均フラックスを求めることが可能と考えられた。
さらに、領域平均フラックスに基づく森林CO2吸収量の評価を行った。表2に示したように、NEEは大きな年変動を示した。これには、降水量、日射量、気温および水蒸気の飽差などの気象条件に加えて、窒素負荷や病害虫なども複合で影響すると考えられた。例えば、2004年は他の年より落葉時期が早く、9月中に葉量を減らしたにもかかわらず生産量は高かった。領域平均フラックスによるNEEは渦相関法を上回り(図20のA参照)、その結果、群落呼吸量および群落最大光合成量も大きく見積もられた(図20のB参照)。3年間のNEEの平均値は渦相関法では180gCm-2yr-1、領域平均フラックスでは266gCm-2yr-1と約1.5倍になった。以前、当試験林と類似した立地条件の40年生のダケカンバ林で1998年から1999年にかけて実施された渦相関法によるフラックス測定では163gCm-2yr-1で、本実施例の結果に近い値であった。また、熱収支のインバランスを解消させるようにフラックスを補正したところ、214gCm-2yr-1となり約1.3倍となった。生態系呼吸量(Re)と生態系純生産量(NEE)の和で表されるGPP(Gross primary production)について比較すると、従来法で875gCm-2yr-1、補正後では1146gCm-2yr-1で、当試験林の3年間平均値811gCm-2yr-1および領域平均フラックスによる1156gCm-2yr-1と同等の結果であった。これは本実施例とは異なった補正法を用いているが、本実施例の結果と同様に渦相関法による過小評価が起こっていると考えられた。
以上のとおり、本実施例においては、森林におけるCO2吸収量を高い精度で評価することを目的として、シンチロメータ1の空間スケールから決定した関数を用いて渦相関法の測定値を補正することにより、領域平均フラックスを算出した。ここでは、熱収支が整合するように調整した領域平均フラックスと、森林における群落呼吸量や群落光合成量との比較を通して領域平均フラックスの確かさを検証した。その結果、領域平均フラックスを基に算出した群落呼吸量は、以前に開発されたチャンバ法による土壌呼吸量の評価結果との一致を見た。また、従来法である渦相関法から算出した群落光合成量がダケカンバ個葉の計測結果から導いた群落光合成量を大きく下回るのに対し、領域平均フラックスによる群落光合成量は個葉の計測結果とよく対応することが明らかになった。以上より、森林における熱収支および群落のCO2収支において、本発明者が領域平均フラックス評価法として開発した手法が従来の渦相関法の信頼性を高めることが分かった。例えば、従来法によって評価したCO2吸収量の3年間の平均値は180gCm-2yr-1で、類似の立地条件にある落葉広葉樹林の値と同等であった。それに対し、本発明にかかる測定手法を用いて推定した森林によるCO2吸収量の3年間の平均値は266gCm-2yr-1で、渦相関法による評価値の1.5倍であった。
風速、乱流強度、TKE(乱流運動エネルギー)、温度、比湿、CO2の各要素に関し、シンチロメータ1(SCとも表記している)および渦相関センサ(SC)間の消散率比と風向の関係について実験を行った(図21参照)。図中の角度は方位を表し、北を0°としている。グラフ中の値が大きければ大きいほどシンチロメータ1による消散率の比が大きく、例えば風速でいえば、北風ほどシンチロメータ1による消散率の比が大きくなる傾向が強い、という結果が得られた。
摩擦速度u*、顕熱フラックスH、潜熱のフラックスλE、そしてCO2フラックスfcのそれぞれに関し、シンチロメータ1および渦相関センサによるフラックス算出値の比較を行った(図22参照)。図中、四角い記号で示すシンチロメータ1のパスに平行な西風、丸い記号でシンチロメータ1のパスに直交する北風のデータを示した。各期につき、シンチロメータ1による測定結果と渦相関センサによる測定結果の比は図に示すように、シンチロメータ1による測定結果が大きくなる傾向にあった。
シンチロメータ1と渦相関センサのそれぞれの場合の非着葉期におけるインバランス率をグラフに表示することにより、熱収支インバランスが改善されたことの確認を行った(図23参照)。インバランス率RI、熱収支式はそれぞれ図中に示す式で求めることができる。この結果、本発明にかかるCO2フラックス測定法によれば従来法よりもインバランス率が改善されることがわかった(図23参照)。
本発明のフラックス測定法におけるフラックスの定義と測定法を説明するための概略図である。
シンチロメータの測定原理を示す図である。
本発明におけるフラックス算出のデータ処理の流れを示す図である。
シンチロメータのパス長による測定領域の違いを渦相関センサの場合も併せて説明するための図である。
森林生態系のCO2交換の概念を示す図である。
風向がパスに直交するときのフットプリントおよびソースエリアを示す図である。
ソースエリア比RAと風向θの関係を示す図である。
領域平均フラックス算出の手順を示す図である。
欠測データの補間のフローを示す図である。
εTSAS/εTEC に対する乱流強度σw/uおよびソースエリア比RAの関係、ならびに補正関係f(σw/u,RA)と乱流強度σw/uおよびソースエリア比RAの関係を示す図である。
構造関数の例と曲線フィットによる消散率決定の内容(2002年6月10日)を示す図である。
εTSAS/εTEC とf(σw/u,RA)の比較を示す図で、a)は昼間、b)は夜間を表している。
ソースエリア比RAとインバランス率I(ただし、Iは3年間の平均値)の関係を示す図である。
熱フラックスおよびCO2フラックスの日変化における渦相関法ECと領域平均フラックスSASの関係を示す図である。
気温と群落呼吸量の関係(2004年)を示す図であって、見やすさのためσw/u>0.3 の条件部分についてのみ表しているものである。
リター積算に基づいた落葉期のLAIの変化を示す図である。
葉群の吸光量I0および光合成量(P)の積算葉面積(F)に対するプロファイルを示す図である。
群落光合成のポテンシャルおよび群落光合成量の入射光量(PAR)との関係(2004年7月)を示す図である。
群落光合成ポテンシャルと群落光合成量の比較(PAR=2000μmolm-2s-1)を示す図である。
領域平均フラックスと渦相関フラックスによる生態系純交換量(NEE)および群落最大光合成量(Ag)の比較を示す図である。
シンチロメータ(SC)および渦相関センサ(EC)間の消散率比と風向の関係を各要素について示す図である。
摩擦速度u*、顕熱フラックスH、潜熱のフラックスλE、CO2フラックスfcのそれぞれに関し、シンチロメータおよび渦相関センサによるフラックス算出値の比較結果を示す図である。
シンチロメータと渦相関センサのそれぞれの場合の非着葉期におけるインバランス率が改善されたことを示す図である。
符号の説明
1 シンチロメータ
2 トランスミッタ
3 レシーバ