2009年からJAXAとNASA/JPLのチームは、RRVにて合同で代替校正キャンペーンを実施している。

代替校正とは、衛星がRRV上空を通過する時刻に同期して、地上にてフィールド携帯型分光光度計（ASD Field Spec）を使用して地表面反射率の測定を行うキャンペーンである。両者を比較することで衛星センサの輝度校正を行うことを目的とする。

RRVはネバダ州中東部、標高約1400mに位置する長さ30㎞、幅20㎞の大きな乾燥湖であり、地表面反射率がほぼ一様であること、GOSATやGOSAT-2の観測可能な地表範囲（以下、フットプリント）が収まることから、代替校正に適した場所である。地上の分光計で測定した地表面反射率と気象条件等を入力として放射伝達モデルを用いて計算した理論値と、衛星が観測したスペクトルとの比較を行うことで、センサの輝度校正が可能となり、地上測定に同期して複数の衛星がRRVを観測することで、衛星をまたがったクロスキャリブレーションも可能となる。

図2-1 RRVでの代替校正イメージ図
（実線はGOSAT path36、点線はpath37からの観測時の衛星、太陽の幾何情報を示す）

図2-2 GOSAT軌道とRRVとの位置関係

図2-3 フォワード計算と衛星観測スペクトルの比較のイメージ図

(1)MODIS BRDFプロダクト

RRVでのキャンペーンでは、図3-1に示す測定点につき、分光計を用いて地表面反射率を測定する。測定波長は400～2500nmである。反射率の測定例を図3-2に示す。

現地の測定点（図3-1）で直下方向に測定した地表面反射率を使用し、衛星の方角（天頂角及び方位角）方向の反射率を求めるには、BRDFモデルを使用する必要がある。RRVは乾燥した湖であり、森林のような反射率のホットスポットが存在しないことから、MODISのBRDFプロダクトの半経験値モデルであるLambertian（等方反射成分）+Ross Thick（体積散乱成分）+Li Sparse（幾何散乱成分）の3つのカーネルにて表現することとした。例を図3-3に示す。

Ross Thickカーネルは、植生などの複雑な地表における体積散乱を表し、Li Sparseカーネルは高さのある樹木等が独立して散乱する状況を表す。

MODISプロダクトとしては、MCD43A1(MODIS/Terra+Aqua BRDF/Albedo Model Parameters Daily L3 Global - 500m)を使用した。

体積散乱成分のRoss Thickカーネルの計算式は以下のようになる。

ここで、cosξは以下の式で表される。

幾何散乱成分のLi Sparseカーネルは、太陽に照射される樹幹による散乱と影とならなかった地表面の散乱の和として、以下のようにあらわされる。

hは樹木の高さ、bは樹幹部の鉛直方向の半径、rは水平方向の半径であり、MODISでは以下のように定義される。

図3-1 RRVにおける地表面反射率測定点とGOSAT path36からの観測時のフットプリント例
（ナディア観測の場合衛星/センサに指定した観測点の緯度経度を中心に直径10㎞程度の円となる）

図3-2 2024/06/16 M03サイトでの地表面反射率測定結果

図3-3 GOSAT-2 2024/06/16 20:57:04(UTC) path72からのRRV観測時のRoss Thick及びLi Sparseカーネル
（赤が太陽位置、水色が衛星位置）

(2)衛星のフットプリントへの拡張

RRVの地表面反射率測定サイトにて「ナディア（直上）方向から測定した“単一地点での地表面反射率”」から“衛星のフットプリントに対応する地表面反射率”を求めるために、フットプリント内に含まれるMODISグリッドを使用し、式(3-13)により算出した。

図3-4 地表面反射率測定サイト(M03)から衛星フットプリントに対応する地表面反射率への拡張方法

図3-5 GOSAT-2 2024/06/16 20:57:04(UTC) path72からのRRV観測時のBRDF補正結果
（地上では測定サイト（M03）にてナディア方向（真上）からASD Field Spec（JAXAFS4-1）にて地表面反射率を測定）

(1)気温・気圧・相対湿度プロファイル

RRVにおける気圧、気温、及び相対湿度などの情報は、現地で観測されたラジオゾンデデータを用いて取得している。ラジオゾンデとは、気象観測用の装置であり、気球にセンサを取り付け、温度、気圧、湿度等の高度分布を測定する装置である。

しかしながら、近年のRRVキャンペーン日には、ラジオゾンデによる観測を行わない年もあるため、その際には代替としてNational Centers for Environmental Prediction（以下、NCEP）による再解析データ（https://psl.noaa.gov/data/gridded/data.ncep.reanalysis2.html）を利用している。

ラジオゾンデデータは、地表面から高度約22,280mまで、2秒間隔でデータを取得しており、NCEPの気象データは17層（1000, 925, 850, 700, 600, 500, 400, 300, 250, 200, 150, 100, 70, 50, 30, 20, 10hPa）の気圧層で定義されている。

そのため、LIDORTに入力する前処理として、RRVの地表面気圧と各層の下端の気圧の比が一定間隔になるように20～24層程度にリサンプリングしてから処理を行う。その結果、ガス濃度や濃度変化が大きい地表面付近（下層）はより細かく、濃度が薄い上層は粗く分割することで、処理精度の向上を図っている（図4-1参照）。

図4-1 気圧比間隔一定の層へのリサンプリング

ラジオゾンデおよびNCEP再解析データにおいて、水蒸気量は通常、相対湿度または比湿の形式で提供される。CO₂やCH₄と同様に、水蒸気量をH₂Oとしてのガスプロファイルとして利用するためには、まずこれらの測定値を体積混合比（Volume Mixing Ratio, VMR）に変換する必要がある。この変換には以下の式を用いる。また、飽和水蒸気圧の計算にはTetensの式を使用した。

(2)エアロソル

(a)エアロソルによる散乱

通常複雑な形状のエアロソルの散乱は厳密にはT-Matrix法を使用して求めるが、本研究ではエアロソルの粒形は球形であると近似し、事前に粒形分布（Particle Size Distribution）と複素屈折率（実部が屈折、虚部が吸収を表す）を入力に種類ごとに散乱特性を計算する。エアロソルを球形として光の散乱を扱う理論をミー散乱と呼ぶ。

エアロソル種ごとに消散係数、単一散乱アルベド、位相関数の算出を行う。

本研究では以下の種類のエアロソルを扱っている。

・雲（Water Cloud, Ice Cloud）
・成層圏エアロソル（Stratospheric aerosol）
・黒色炭素（ブラックカーボン）（Black Carbon）
・ダスト（Dust）
・有機炭素（Organic Carbon）
・硫酸塩（Sulfate）
・海塩（Sea Salt）

消散係数は吸収係数と散乱係数の和であり、(4-5)式のようにあらわされ、媒質の内部を伝搬する放射は、複素屈折率の実部、虚部が0ではない一様な媒質の中で、吸収及び散乱される。

• β_e：消散係数
• β_s：散乱係数
• β_a：吸収係数

消散係数β_eを使用すると、点s₁とs₂ (s₂≥s₁)の間の消散は以下のようにあらわされる。

(4-6)式の右辺カッコ内の積分量は点s₁とs₂の間の光路長(Optical Path)と呼ばれ、大気を鉛直方向に積分した場合、光学的厚さ(Aerosol Optical Thickness : AOT)または光学的深さ(Aerosol Optical Depth : AOD)と呼ばれる。

単一散乱アルベドは(4-8)式のようにあらわされ、エアロソルによる吸収＋散乱における散乱の割合を示す。

位相関数とは、光がエアロソル粒子によって散乱されるとき、前方や後方、横方向にどれくらいの強度で放射されるかを示すもので、角度依存性を示す関数である。

位相関数は散乱角と入射角の間の散乱角Θにのみ依存し、ルジャンドル多項式に展開することができる。

Mは位相関数の精度に依存するため、ルジャンドル展開係数のβ_lは十分小さくなるところまで計算し、それ以降は打ち切る。

エアロソルの位相関数の例を図4-2に示す。

図4-2 粒子の半径を10μm、複素屈折率を1.33257+1.67E08i、波長550nmとした場合の位相関数例

(b)エアロソルの地域分布と高度分布

ア)地域分布とRRV処理での扱い

LIDORTを使用したフォワード計算では、エアロソルの種類は4章(2)(a)で述べた7種類を扱うことができる。

フォワード計算を行う上での仕組みとして、雲及び成層圏エアロソルの2種類に加え、MERRA-2再解析データを基に事前の解析で地域ごと、月ごとにエアロソルを優勢な順にランク付けを行った結果の上位2種類データの計4種類のエアロソルが存在しているとし、AOTを配分することとする。

ただし、RRVキャンペーンは晴天時に実施されていることと、RRVにおいては、「ダスト」が優勢であることを考慮し、「ダスト」のエアロソルのみが存在するという仮定でフォワード計算を実施した。

図4-3 MERRA-2再解析データよりNASAが作成した地域ごとのエアロソル優位順

イ)高度分布

エアロソルの高度方向の分布は、エアロソルの種類ごとに、ピーク高度、標準偏差を定義し、ガウシアンの形状で分布していると仮定して処理を行う。

• S(x)：エアロソル高度分布形状
• x_a：エアロソルピーク高度
• σ_a：ガウシアン幅(1σ)

x=P/P_surfaceとすると単位圧力あたりの光学的厚さは(4-12)式で表される。

AOT_refはAERONETによる実際の観測値である値であり、解析対象波長ウインドウに最も近い波長のAOTをリファレンスとして使用する。Q_λは波長λにおける消散効率である。

図4-4 RRVでのエアロソル（ダスト）の分布
（AOT=AERONETによる観測値(この例では0.05@755nm), Peak Height=0.9, Width(σ_a)=0.05と仮定）

(3)大気分子

(a)大気分子による散乱（レイリー散乱）

入射する電磁波の波長λに比べ、散乱体の粒子の半径rが非常に小さい場合に生じる散乱をレイリー散乱という。散乱理論では で定義されるサイズパラメータから大気分子によるレイリー散乱なのか、エアロソルによるミー散乱なのか判断が可能である。レイリー散乱はx«1の場合に相当するため、大気分子(r<1nm)による可視光(λ~500nm)の散乱はレイリー散乱とみなすことができる。波長500nm付近にはNO₂の吸収帯があり、NO₂を観測対象とするセンサ（TROPOMI、TEMPO及びTANSO-3）はこの波長域をウインドウとして持つため、フォワード計算をする際にはレイリー散乱の影響を考慮することが必須となる。

図4-5 散乱体の粒子径と波長との関係

レイリー散乱の位相関数のルジャンドル多項式展開係数は となる。δは散乱における非対称因子であるが、レイリー散乱は前方と後方にほぼ対称な散乱パターンとなるため、δの値はほぼ0となり、(4-13)式のようにあらわされる。

図4-6 レイリー散乱の位相関数の3次元表示

レイリー散乱の吸収断面積σ_R (m²/molecules)は、以下のモデルで計算する。

• λ：波長[μm]
• ρ：脱分極因子(0.029)
• N_s：標準状態（273.15K、1013.25 mbar）における大気分子密度（2.6871019 molecules m^-3）
• n_s：標準状態（273.15K、1013.25 mbar）における標準大気の屈折率

(4)対象ガス

(a)ガス濃度プロファイル

CO₂、CH₄などの温室効果ガス濃度のプロファイルについては、NASAの大気研究のための航空機実験プロジェクトであるAJAX（Alpha Jet Atmospheric eXperiment）の観測データを利用している。

本来、時期や場所の異なるガス濃度プロファイルを入力にする場合は、以下の4種類の補正を行う必要がある。ただし、今回の一連の解析作業では、RRV上空で過去の同時期に取得されたデータを利用しているため、長期トレンド補正のみを施している。

• ●高度リサンプリング
  リサンプリングは対流圏界面高度に基づいて実施する。
• ●緯度勾配補正（目標地点が異なる場合）
  体積混合比 (VMR) は、プロファイルの緯度と観測地点の緯度の違いに基づいて補正する。成層圏および中間圏ではガス分布が赤道を中心に対称であると仮定し、地表付近では反対称であると仮定する。
• ●長期トレンド補正（目標年が異なる場合）
  基準日から観測日までの経過時間を考慮して補正する。
• ●季節変動補正（目標季節が異なる場合）
  体積混合比 (VMR) の補正において季節による変動を考慮する。

AJAXの観測データの例を図4-7に示す。

図4-7 2016年7月1日のRRV上空でのAJAXによる観測データ（左：CO₂、右：CH₄）

(b)対象ガスの吸収断面積（HITRANデータベース）

エアロソルについては4章(2)で述べたが、ガスの吸収断面積とは、光がガス（CO₂やCH₄など）に当たった場合、どの程度その光が吸収されるかを表す指標である。

また図4-8に示すとおり、ガスの種類により、吸収する光の波長帯が異なる。

この特性を使用し、対象ガスが存在する媒質と存在しない媒質に同様の強さの光を入射させ、通過後の光を観測し、どの波長の光がどれだけ減光しているかを計測することで、対象ガスがどれだけ光路上に存在したかを定量的に推定できる。

このようなインバース計算（リトリーバル）をする際には、どのガスによる影響かを識別するために、複数のガスの吸収が重複しない波長帯を選ぶことが望ましい。

図4-8 主な観測対象ガスの吸収波長帯

本研究のフォワード計算では、ハーバード・スミソニアン天体物理学センターによって研究、データ提供されており、多くの研究分野で利用されている標準的なデータベースであるHITRANを使用している。（https://hitran.org/）

ガスによる吸収度合いは、温度や圧力によって変化し、その結果として吸収線形状も変化する。これには、分子の熱運動によるドップラー効果と、分子間衝突による圧力広がりの2つの要因があり、これらを組み合わせたものがVoigtプロファイルである。Voigtプロファイルを使用することで、高度方向に分けた各層の気温・気圧に対応した吸収断面積データを用いて、各層における吸収特性を高精度に反映させることが可能となる。

図4-9 気圧1019.0hPa 気温286.3Kにおけるメタンの吸収断面積（Strong CH₄帯2150-2430nm）

4章で述べた処理対象の気象条件に基づいて吸収断面積を求め、計測されたガスのVMR（体積混合比）を基に次の式で対象ガスの光学的厚さを計算する。

(c)対象ガスの光学的厚さ

(1)の気温・気圧のプロファイルの各層の上端、下端の圧力を用いて、層内の吸収断面積を一定としようとするとその間のガスの吸収断面積の変化が大きくなるため、誤差が生じる。そのため、各層をさらに10層のサブレイヤーに分割し、その間は線形に吸収断面積が変化すると仮定し、光学的厚さを算出する。

• P：気圧[Pa]
• T：気温[K]
• q：局所比湿
• g：重力加速度[m/s²]
• M_dry：乾燥空気のモル質量[kg/mol]
• K(P_cen,T_cen)：サブレイヤー中央の気体吸収係数

各サブレイヤーについて光学的厚さを求めたのち、シンプソンの公式で積分することで、層内の光学的厚さを求める。

(5)大気分子、ガスとエアロソルの光学特性

LIDORTのフォワード計算では、大気分子によるレイリー散乱、微量ガスによる吸収と散乱、エアロソルによる散乱は気象条件としては与えた各層内では均一であるとの仮定のもと、計算を行う。この場合の光学的な特性の入力は以下のようになる。

本研究では、太陽照度データとして、TSIS-1 HSRS(TSIS-1 Hybrid Solar Reference Spectrum)のデータを使用している。

TSIS-1は国際宇宙ステーションに搭載するために、NASAにより、2017年12月15日に打ち上げられ、地球（地球大気への突入前の）への太陽光の入力、つまり太陽放射照度を正確に測定するものである。

高波長分解能（数種類の波長分解能）でコロラド大学から提供されている。

太陽照度データの波長分解能でフォワード計算を実施し、7章で説明する装置関数にて畳み込みを行うことで、理論スペクトルを生成することができる。

GOSATの通常観測に使用される拡散板表面は常に宇宙環境にさらされているため、時間経過とともに劣化する。その劣化の程度を評価するための基準（リファレンス）として、ひと月に一度、拡散板裏面を使用した太陽照度校正運用を実施している。取得した太陽照度データの本来の使用目的は、通常観測（拡散板表面）の観測精度を保つためのものであるが、宇宙環境にさらされていない拡散板裏面で取得した高波長分解能かつ長期（16年以上）の太陽照度のデータは非常に貴重であり、TSIS-1 HSRSデータの精度向上に貢献する可能性があるため、今後コロラド大学のチームと協力しデータ解析を実施していく。

TSIS-1 HSRSのスペクトル図と、GOSATの波長範囲を図5-1に示す。

図5-1 TSIS-1 HSRS（太陽放射照度のスペクトルとGOSAT-2の観測波長帯）

LIDORTは以下の基本放射伝達方程式に従っており、次の3つの現象が同時に起きている状態と考える。

① 入射した放射輝度が一次散乱や吸収により減衰する現象
② (μ^',φ^')方向から(μ,φ)方向に多重散乱により放射輝度が増加する現象
③ 入射太陽光の直達成分が一次散乱を受けることによる(μ,φ)方向への放射輝度が増加する現象

式中のI(x,μ,φ)は、入射した放射輝度が一次散乱や吸収を受けることで生じる減衰①を表す。また、②および③は放射輝度の増加要因であり、放射源関数J(x,μ,φ)と呼ぶ。

右辺の第1項は(μ^',φ^')方向から(μ,φ)方向に多重散乱される放射輝度の増加分（多重散乱寄与のソース項）を表し、右辺の第2項は入射太陽光の直達成分が一次散乱を受けることによる(μ,φ)方向への放射輝度の増加分を表す。

さらに、LIDORTは擬似球面大気を模擬できる機能を備えており、球面幾何を反映した精密な放射伝達解析が可能である。これにより、平行平面大気を仮定する場合に比べて、より現実的なシナリオに対応した計算が行え、地球大気における光学的厚さの変動や散乱の影響をより正確に評価できる。擬似球面モデルの導入により、地球の曲率が放射輝度に与える影響を考慮し、精度の高い放射伝達計算が可能となる。

平行平面と仮定した場合、4章で述べた各気圧層への太陽光の入射角は一定であるとみなされるが、実際には球面であるため、各気圧層は同心円状の構造を成し、層ごとに入射角が異なる。このため、太陽や衛星の天頂角が大きくなるほど、球面を模擬する効果が大きくなる。

装置関数（Instrument Line Shape Function）とは、観測装置が持つ特有の光学的特性を数理的に表現したものであり、観測装置が測定するスペクトルデータに与える影響を示す関数である。温室効果ガス観測衛星においては、地表面で反射された光のスペクトルを観測し、吸収される波長に基づいて大気中のガス濃度を算出するが、その際、装置が持つ分解能や検出特性に起因するぼやけや広がりの影響を補正するために装置関数が使用される。

装置関数は、装置の光学的設計や構成要素によって決定されるが、GOSATやGOSAT-2の場合はSinc関数、OCO-2, OCO-3, TROPOMI等はガウス関数やスーパーガウス関数で近似される。

GOSAT-2の装置関数を図7-1に示す。

図7-1 GOSAT-2の装置関数（S偏光、中心波数: 12950 [cm^-1]）

図7-2 GOSAT-2の装置関数の中心波数箇所拡大図

一方、衛星Sentinel-5Pに搭載されているラインセンサ(TROPOMI)であれば、次のようなスーパーガウス関数の形状となる。TROPOMIの装置関数を図7-3に示す。

図7-3 TROPOMIの装置関数(detector: 200, center wavelength: 730.032[nm])

LIDORTで算出した理想的なモデルスペクトルに対して、この装置関数を畳み込み処理することで、衛星が観測するスペクトルを再現することができる。

図7-4 GOSAT-2の装置関数での畳み込みイメージ図

GOSAT-2の例で示すと、O₂A帯スペクトルはFTS-2装置関数で畳み込むことによって次のように変化する。

図7-5 GOSAT-2 装置関数を畳み込み前のO₂A帯スペクトル図

図7-6 GOSAT-2 装置関数を畳み込み後のO₂A帯スペクトル図

(1)理論スペクトルと衛星センサの観測スペクトルの比較

LIDORTを使用しフォワード計算した理論スペクトルと衛星が観測したスペクトル（L1Bプロダクトに含まれるスペクトル）について、GOSAT-2の例を図8-1に示す。

Band1のO₂A帯については、スペクトルの肩の部分で5~6%程度観測スペクトル側の輝度値が高く、吸収線にリンクした残差が残っていることがわかる。Band2のWeak CO₂帯では、残差はほぼなく、衛星の観測スペクトルが理論値とほぼ一致していると言える。Band3のStrong CO₂帯では、Band2と同様に肩の部分では残差は発生していないが、吸収線にリンクした差分がわずかにみられる。

GOSAT-2/FTS-2を含むセンサは、宇宙環境下では徐々に感度が劣化するため、経年劣化補正が行われている。しかし、補正量が過大であると、理論値との乖離が大きくなる可能性がある。一方、Band1に見られるような数%の誤差は、理論スペクトル計算する際に使用する地表面反射率の見積りが低い可能性や、エアロソルによる散乱の影響が過小評価されている可能性もある。

そのため、今後は各国のセンサの輝度値との比較を通じて、理論スペクトルの精度向上を目指すとともに、センサ間でのスペクトルレベルでのバイアス削減を図り、各センサの校正状況を引き続き確認していく。

図8-1 GOSAT-2/FTS-2 2024/6/16 path72からの観測とLIDORTによるフォワード計算との比較
（左からLIDORTフォワード計算（黒）&衛星観測スペクトル（赤）、残差（②-①）、比（②/①）
上からO₂A帯、Weak CO₂帯、Strong CO₂帯）
（JAXA所有Field Spec(JAXAFS4-1)を使用してM03地点で測定した地表面反射率を使用）

(2)RRVキャンペーン2024の測定データを使用した理論スペクトルと衛星観測スペクトルの誤差

RRV2024に参加した6センサについて、フォワード計算を実行し、各センサの観測スペクトルとの誤差を箱髭図にしたものを図8-2に示す。

図8-2 RRVキャンペーン2024の観測スペクトルと理論スペクトルとの誤差

図8-2からわかる様に、GOSAT、 GOSAT-2及びOCO-2については、輝度値の誤差がほぼ10%以内に収まっている。経年劣化補正は各衛星チームがそれぞれ独自の判断で実施しているが、唯一OCO-3については理論スペクトルの方が観測スペクトルより全体的に輝度値が高い傾向があり、経年劣化補正の程度が少なすぎるのでは？という可能性が考えられる。またCH₄帯を持つGOSAT-2、TROPOMIを見るとO₂AやCO₂帯に比べ差が大きく、CH₄帯におけるフォワード計算で使用する反射率が低い可能性も考えられる。また、各衛星とも観測時の幾何条件（衛星天頂角、太陽天頂角、相対方位角など）が異なるため、それも考慮し、理論スペクトルと観測スペクトルとの差、各センサ間の差を削減していく必要がある。

RRVでの代替校正キャンペーンで取得した地表面反射率及びMODISプロダクト、ラジオゾンデまたはNCEP気象データ、AERONETのエアロソルAOT、AJAXのガス濃度を入力として、放射伝達モデルのLIDORTを使用したフォワード計算により、6センサについて理論スペクトルを生成することに成功し、衛星センサが観測するスペクトルの校正・検証を実施するための環境が整った。

現在も課題として残っている吸収線にリンクしたスペクトル残差や地表面反射率の補正方法の改善に取り組み、今後は各国の衛星チームと連携し、このスペクトルを基に輝度校正を行っていくことで、マルチプラットフォーム間のバイアスを削減し、データ統一を進めていく。

今後打ち上げ予定のGOSAT-GWや海外の温室効果ガス観測衛星についても、速やかに理論スペクトルの情報を提供し、校正に役立てていく必要がある。

また、フォワード計算が正確に行えるようになったため、今後は観測スペクトルから対象ガスの濃度を推定する処理（リトリーバル処理）も行い、衛星が観測したスペクトルから正確なガス濃度の算出とデータの提供を行っていきたいと考えている。

本研究で紹介したフォワード計算の手法は、宇宙航空研究開発機構（JAXA）からの委託業務「GOSAT-2スペクトル品質確認」の検討結果を複数衛星のセンサ（マルチプラットフォーム）に発展させたものである。本論文中に掲載された一部の図（図1-1、図2-1、図2-2、図2-3及び図3-1）はJAXA/EORCのサイトにて公開されているものを使用している。

本論文の執筆にあたり、関連する業務において有益な技術的アドバイスをいただき、データの利用と発表について快くご承諾くださったJAXAの須藤洋志氏、塩見慶氏、久世曉彦氏に深く感謝する。

また、放射伝達モデルLIDORTの作者であるRT Solutions, IncのRobert J. D. Spurr氏には、他では得られない貴重な知見の提供と多大なご支援をいただいた。Spurr氏のご協力と専門的なご指導なくして、本論文の完成は不可能であったと感じている。深く感謝の意を表する。