国立防災科学技術センター研究報告第43号

1989年3月

528．235／236

リモートセンシング画像処理基本ソフトウェアの開発
I．広域衛星画像の地理的変換

矢崎忍＊・建石隆太郎＊＊

国立防災科学技術センター

Developme11t
inof

RemoteFoumdamemta1Softwaren

Semsi皿gI．Geographica1Transformation

of

Image

WidePm㏄ssing

View

Sate11ite

Image

By

Shinobll

YAZAKI

Mガo舳ばθ∫θκんCθ肋7〃1）｛∫α∫肪肋閉ガo〃，切舳

Rylltaro
Chiba

TATEISHI

University

AbstractAs

oneof

theseriesof

deve1opment

programoffomdamenta1softwaresofremote sensing data proccesing，correspondence between sate11ite image andgeographical1ocation is dea1t with．The polynomina1formu1a of transformation
basedonthatareview

satellite

GCP

data

imagessuch

Forimages，formula

those

va1id

asfor

LANDSAT

image

MOS−1VTIR

based

on

theor

is

NOAA

notsufficientforwideneccesary．In

this

AVHRR．

map

projectionisformu1a

betweenimage

report．

1）

Composition

oftransformation

geographica1coordinate
2）

Superposion

3）

Transformationdifferent

aregivenfor

And

as

coordinate

ofand

eva1uation

geographicalofone

map

ofgridsatellite

coeficients

inonimage．

image

sate11itetoanotherthecoordinate

andformula．

image

coordinate

ofprojection．

imagesthesystemofunified
isMercator

andprojection

LCC

projection．

coordinatefor

resu1tingimage，thesquare

proposed．

1．はじめ1こ

1987年より日本の人工衛星MOS−1（Marine第4研究部，計測研究室

榊千葉大学

一223一

Observation

Sate1lite一ユ） の受信国立防災科学技術センター研究報告第43号

1989年3月が始まり，当センターにおいても解析が始まった．MO

ME

S

S

R（Multispectra1Electronic（Visib1e
Scanningand

ThermalSelf−Scanning

InfraredRadiometer），VTIRRadiometer）およびMSR（MicrowaveRadiometer）を搭載している．このうちME能，狭域の放射計で仕様はLAND

S−1は3つの地上撮像センサー

SATに近い．V

TISS

Rは可視，近赤外の高分解

Rは可視〜熱赤外の中分解能，広域の放射言十で海面・地表・雲頂温度分布の測定を目的にしている．またM

S

Rはマイクロ波の低分解能，広域の放射言十で大気中の水分量や雪，氷の分布などの測定を目的にしている．

一方，既に運用されているアメリカの気象衛星NOAAの画像は，その観測波長等の仕様がMO

S−／VTI

Rと一部共通しており，VTI

Rと相互補完的に利用できる．広域の気象，気候変動の予測等，防災分野におけるこれらのデータの利用が期待される．

リモートセンンングの分野で，当センターはこれまでLANDSATのMS

SやTM，あるいはSPOTの画像を用いて，おもに陸域を対象に災害状況の検出や災害に関連した環境の調査を行ってきた．これらの衛星画像はおもに可視，近赤外の比較的狭い範囲（数10㎞〜100㎞）の画像であった．これに対してMO

S−1／V

T

I

R，MS

RやNOAA／AVHRR

は数100㎞〜数1000㎞におよぶ広い範囲の画像であり，画像と地上位置の対応において従来のLANDSAT等と異なる手法を要する．
本報告は，人工衛星リモートセンシングにおける基本的ソフトウェアの開発の一環として，
MO

S−1VT

I

R，MS

RおよびNO

AAを対象に実施した広域画像の地上位置との対応および幾何学的変換にっいて述べたものである．

2．広域衛星画像と地図投影

衛星画像を扱う場合画像上の各ピクセルと地上における位置の対応は最も基本である．衛星画像とそれに対応する地球上の位置の関係を示したものが図1である．LAND
MO

S−1

SATやMESSRなどにより陸域の地表を解析対象とする場合，通常良く晴れた日の雲の少ない画像を扱うので画像上にGCP（groundcontro1point）を容易にとることができる．さらにこれらの場合，一度に解析する領域の大きさは高々100㎞×100㎞程度であるので地球の曲率の影響は小さく，経度，緯度と画像座標とはG

C

Pデータから最小自乗法によって係数を決定した1次または低次の変換式により精度よく変換できる．

これに対してMO

S−1／VTIR，MS

RやNOAAなどの広域の画像の場合，地表ばかりではなく海面や雲も解析の対象にするので必ずしもつねにGCPをとることができるとは限らない．また，よく晴れた画像で陸域にGCPをとることができても，この場合地球の曲率が無視できないので変換式は高次になる．これに対して陸域の範囲は全画像においては片寄っているので，GCPデータのみから画像全体で有効な変換式を構成することは困難＿224一

リモートセンシング画像処理基本ソフトウェアの開発

I

広域衛星画像の地理的変換一矢崎・建石

画像座標（未補正）

衡屋の位置、姿勢←地轟の形状、白転

地球座標（λ，g）

地図座標（エ，ツ）

画像座標（1↓，u）（地図投影済）図1

衛星画像と地球上の位置との関係

Fig・Fig．1

である．またM

STransformation

processfromsatel1ite

Rの画像では分解能が低くもともとGCしたがって，これらの画像に対してはG

C

imagetogeographicalcoordinate

Pをとることかできない。

Pを用いない方法，すなわち衛是と地球の位置関係のデータから原画像の各ピクセルに対応する地上位置を言十算する方法によらなければならない．我々が利用できる衛星画像には地図投影のなされていない原画像と規定の地図上に投影された画像があるが，本稿では地図投影された画像におけるピクセルと地上位置の対応，すなわち図1の枠で囲んだ部分についてのみ述べる．

地図投影には多くの図法があるが，ここではMOS一で用いられているLCC図法（LambertConformalConic，ランベルト正軸円錐図法）とMOS−1およびNOAAで用いられているメルカトール図法をとりあげ，これらに対して

川

地球座標から画像座標への変換および逆変換

12〕衛星画像への等緯度線，等経度線および海岸線の書き入れ

13〕投影法の異なる衛星画像の相互の変換および任意の地図上への変換を実施した．

一225一

国立防災科学技術センター研究報告第43号

1989年3月

3．地球座標と画像座標との変換

地上の位置は経度と緯度で表わせるが，ここではこれを地球座標と呼んで（λ，？）で表わすことにする（λ，？はそれぞれ東，北方向を正とする）．一方，画像上のピクセル位置はピクセル値とライン値で表わせる．ここではこれを画像座標と呼んで（u，V）で表わすことにする．画像座標は画像の左上を原点（u＝V＝ユ）として右方および下方をそれぞれ正にとる．地上位置と画像上の位置との対応を知るということは，地球座標と画像座標との変換式を求めることに他ならない．

地図投影された画像の場合地球座標から画像座標への変換は図1に示すように地図座標を介して行われる．地図座標というのは投影面上にとった直交座標で（X，y）で表わすことにする．なお，ここでは地球をベッセル回転楕円体と仮定し赤道半径をRで，離心率をeで表わすことにする．

メルカトール図法とLtanC

C図法の地図投影変換式では？はつねにπg
一十一42ユーeSin91＋eSin9π

ユ十e

Sin91−eSin？

e−2（ユ）

三2（2）またはtan？

一一一42という形で現れる．そこで簡単のため（1）をf（乎）で表わすことにする（このとき（2）は1／f（甲）で表わされる）．すなわち，

f（・）一㎞［云・二ユーeSin91＋eSin？

一1・／叶一云三2ellllll；「（・）f（？）の逆関数f−1（t）も現れるがこれは解析的には表わせない（f■1（t）については付録を参照）．

以下に，メルカトール図法およびLCC図法についてそれぞれの変換式を示す．

＿226一

リモートセンシング画像処理基本ソフトウェアの開発

ψI

広域衛星画像の地理的変換一矢崎・建石yV

図2

メルカトール図法における各座標問の関係

Fig．2

Relationofon

Mercator

thethreecoordinate
projection

systems
planeX

ψ0

λλ0
3．1

メルカトール図法

経緯度（λo，？o）の点を地図座標原点（x＝y＝0）にとり緯線に沿って東向きにx軸を，経線に沿って北向きにy軸をとる（図2）．このとき（λ，g）と（x，y）の関係はメルカトール図法の地図投影法の公式によって

x＝R（λ一λo）（4）

y＝R｛1nf（？）一1nf（？o）｝
で与えられる．

メルカトール図法の場合通常，画像座標を経緯線に沿ってとるので，ここでもそのようにとる（図2）．このとき（x，y）から（u，v）への変換は

X

u：uo＋一d
V＝VOyd（5）−で表わされる．ここで（uo，vo）は地図座標原点（x：y：0）の画像座標，dは？＝0（赤道）におけるピクセルの大きさである（ピクセルは正方形とする）．

（4）と（5）からx，yを消去して
λu＝U＋一D（6）

1

v＝V一一1nf（乎）D＿227一国立防災科学技術センター研究報告第43号

1989年3月

また，これをλ，gについて解いてλ＝D（u−U）（7）？＝f・1［exp｛一D（v−V）｝］

が得られる．ここに

D＝一dRλ。

U＝uo一一（8）D1

V＝vo＋一1nf（ρo）
Dである．これらの変換式からわかるように，この場合の地球座標と画像座標との変換のパラメーターは3つである．

NOAAAVHR

R，87，6．17モードS画像の例

この場合地図投影のパラメーターとしてλ、＝135．O（℃）：画像原点（u＝v＝1）の経度？1＝44．0（d＝）：

〃緯度

3．0（㎞）：赤道におけるピクセルサイズが与えられている．これから（8）により

D＝0−0004704が得られ，（6）に（λ，g）：（λ1，g1），（u，v）＝（1，1）を代入することにより

U＝一5007．80

V＝1812．74が得られる．

3．2

LC

ここではMO

C図法

S−1で用いられている2基線正角円錐図法について述べる．これは軸が地軸と一致し2っの緯度（これらをア1，？。（？1＜乎。）とする）で地球と交わる円錐面に投影する図法である．

一228一

I

リモートセンシング画像処理基本ソフトウェアの開発

広域衛星画像の地理的変換一矢崎・建石

ψ

y図3
uL

Fig．3

C

C図法における各座標間の関係

RelationofthreeonLCC

Projection

thecoordinatesystems

planeX汽λV
λ0

正軸図法であるから，経線は投影面上で常に直線になる．そこで地図座標として経緯度（λ。，乎。）の点を原点（x＝y：0）とし経線に沿って北向きにy軸を，y軸に直交して東向きにx軸をとることにする（図3）．このとき（λ，g）と（x，y）の関係はLCC
図法の地図投影法の公式によって

x＝κ｛f（？）｝μSinμ（λ一λO）（9）

一μy＝κ［lf（乎。）1−lf（乎）1■μ…μ（λ一λ・）］

ここでμ，κはで与えられる．

1一㎞／

COS？1

1−e2Sin2乎2

COS92

1−e2Sin2乎1

1万／・㎞f（？。）
f（91）
（ユ0）

Rcos甲1｛f（91）｝κ＝μ1−e2Sin2ψ1で定義される定数である．

一229一国立防災科学技術センター研究報告第43号

1989年3月

地図投影画像はその定義から投影面上の像を形を変えずに適当に縮尺したものであるから，この変換は原点移動，縮小および回転によって表わされる．すなわち

1（・C・・δ一ySinδ）
u＝uo＋一d
（11）

1d

V＝VO

（・Si・δ十yCOSδ）

ここで（uo，vo）は地図座標原点（x＝y＝O）の画像座標，dは基準緯度（g＝g1，g＝乎2）におけるピクセルのサイズ，

δはV軸のy軸に対するに傾き（y軸に対して一V軸が時計回りの向きのとき正とする）である（図3）．

（9）と（11）からx，yを消去して

1

一μu＝U＋一｛f（甲）｝sin（μλ十」）

D（12）

1

一μv＝V＋一｛f（9）｝cos（μλ十」）D

また，これをλ，？にっいて解いて1

λ＝＿
μarCtanu−U

一∠v−V
（13）

H…
が得られる．

1D2（（・一U）2＋（・一V）2）

ここで

一230一2μリモートセンシング画像処理基本ソフトウェアの開発

I

広域衛星画像の地理的変換一矢崎・建石

dD＝一■
1

一μU＝uo一一｛f（？o）｝
Dsinδ
（14）

1

一μV＝vo一一｛f（？o）｝
Dcosδ

1＝δ一μλoである．したがって，この場合地球座標と画像座標との変換を規定するパラメーターはU，V，D，」およびμの5っである．
VTlR，87．6．17画像の例

この場合変換パラメーターを求める方法はいろいろあるが，ここではそのうちの1つを示す．用いるデータは？1＝20（。）：

基準緯度19。＝50（。）：〃λo＝139．35（。E）：？o＝35．98（。N）

：

d（㎞）

：

ピクセルサイズ＝16．00（o）：

画像軸の傾きxw＝一63．160164（㎞）：WRδ：O．909

yw＝34．636581

（km）

〃

緯度

Sセンターのx座標〃uw＝1787．73である．

地図座標原点経度：

vw二2132．99

2：：y〃〃〃u〃v〃

まず？1，？2から（10）によりμ：0．580483■＝12684．6（km）が得られる．つぎに（ユ1）においてd，δが与えられているので，

一231＿

1っの点の地図座標と画

国立防災科学技術センター研究報告第43号

ユ989年3月

像座標があればu。，voが求まるが，このような点としてWRSセンターの値が利用できる．

すなわち（1ユ）に（X，y）＝（XW，yw），（u，V）：（uw，Vw）を代入してuo；1865．0v〇二2150．5が求まる．最後に，（14）により

D＝7．1662＊10−5
U＝一742．ユ

V：一6941．7」：

一64．89（o）が得られる．uo，vo，d，δはここで述べたやりかたの他に，チックマークデータなどから最小自乗法により求めることもできる（この場合は地図座標原点は任意である）．

4．衛星画像への等緯度線，等経度線および海岸線の書き入れ

VT

I

RやNO

AAの画像範囲のような広域でみると平均的に地表の6〜8割は雲に覆われている．特に陸地付近は雲が出やすく海岸線がはっきり見える画像はむしろ少ない．した
がってこれらの画像では目視によって地上位置との対応をっけるのは容易ではない．このとき画像上に海岸線や経緯線を重ね合わせると，この対応は容易になる．

等経度線，等緯度線および海岸線を経度，緯度で表わしたデータとして用意しておけば，3．で求めたような経度，緯度から画像座標への変換式により，任意の図法により地図投影された画像上に重ねて描くことができる．

図4〜図6にそれぞれVTI

R，MS

RおよびNわせたものを示す、＿232＿OA

Aの画像に経緯線と海岸線を重ね合

I

リモートセンシング画像処理基本ソフトウェアの開発

図4

MOS一ユV

Fig・4

VTI

図5

MOS一ユM

SR画像（LFig．5

MSR

and

RTI

image

image

R画像（LandCC図法）と地図・経緯線との重ね合わせ

geographicalC

広域衛星画像の地理的変換一矢崎・建石

map

and

grid

superposition（theLCC

Projection）

LCC

projection）

C図法）と地図・経緯線との重ね合わせ

geographical

map

and

grid

一233＿

superposition（the国立防災科学技術センター研究報告

図6

NOAA

Fig．6

NOAA

第43号

1989年3月

AVHRR画像（メルカトール図法）と地図・経緯線との重ね合わせ
image

and

geographicalmap

andgridsuperposition（theMercatorprojection）

5．衛星画像の相互の交換

衛星画像を画像座標系の異なる別の衛星画像の上に変換して2つの画像を幾何学的に重ね合わせることはいろいろな場合に必要になる．1っには同種の画像を比較するためである．
例えばVT

I

RとVT

I

R，あるいはVT

I

RとN

O

AAA

V

H

R

Rの観測日時の違う2つの画像をラジオメトリック特性を一致させた後，差をとることにより経時変化を抽出することができる．もう1っは，同一の衛星の異なるセンサーによる2つの画像を重ね合わせる

ことにより1っの多バンドデータとして利用するためである．MOS−1のVTIRとMSRは同時に観測を行い，それぞれ4っのバンドをもっが，これを重ね合わせることにより共通する部分では8バンドのデータとして利用することができる．V

M

S

T

I

Rはおもに温度情報を，Rはおもに大気中の水分量の情報を持っので相互に補正情報として利用できる．

2つの画像を重ね合わせるということは2つの画像座標間の変換を求めるということである．それぞれの画像座標と地球座標の変換が求められていれば，画像座標問の変換は一方の画像座標を地球座標に変換してそれを他方の画像座標に変換することにより得られる．すな

わち，2つの画像A，Bの画像座標をそれぞれ（u。，v。），（uA，VA）→（λ，9）一芋（uB，VB）と表わせる．

一234一（uB，vB）とすると（＊）

リモートセンシング画像処理基本ソフトウェアの開発

I

広域衛星画像の地理的変換一矢崎・建石

画像変換を行なうには各画素毎にこれらの言十算をしなければならないので，無駄な計算はできるだけ省いて計算量が少なくなるようにしておくことは有用である．そのためには（・）の2つの変換式からλ，？を消去して直接u。，v。からu。，v。への変換式を求めておけばよい．ここでは特に利用度の高い，メルカトールとLCC，およびLCCどうしにっいて画像座標間の変換式を求めておく．

5．1

メルカトール画像とL

C

C画像との変換

メルカトール画像の画像座標を（uM，vM），LC

C画像の画像座標を（uL，vL）とす

ると，（uM，vM）と（uL，v。）との変換は（7），（ユ2）からλ，アを消去して

ユuL＝U1＋一exp（μ1vM）

sin（μ1uM＋」1）

D1（15）

1vL＝V1＋一D1exp（μ1vM）cos（μ1uM＋∠1）

およびこれをuM，VMについて解いて・・一÷／・一t叶11；1一・1／
（16）

V

M＝ln

［D12｛（uL−U1）2＋（vL−V1）2｝］2μ1と得られる．ここで

μ1＝μL

DM

D1＝DLexp（μ1VM）」1＝∠L一μ1VM

（17）

U1＝UL
V1＝VL
一235一

国立防災科学技術センター研究報告第43号

1989年3月である．

（ユ7）でサフィックスのMおよびLはそれぞれ（8），（ユ4）で与えられたメルカトールとL

CNOAA

Cの変換のパラメータを示す．

AVHRRのメルカトール投影画像は日本付近を約1500㎞（東西）×1400㎞（南北）（画素数512×480）の5つのブロックに分けて提供される．

NOAA，87．6．17とVT

lR，87．6．17の重ね合せての例

NOAAとVTIRのバラメータは3．1および3．2で求めてあるのでこれらを（ユ7）に代入してμI二2．7306＊ユ0−4Dl＝1．1756＊lO

4∠1：13．46（。）

U1＝一742．11
V1＝一6941．70が得られる．

ここではNOAAをVT
図8はNO

IRの上に変換した例を示す．図7はVT

I

Rのバンド3画像を，AAチャネル4画像をそれぞれ温度表示したもので（観測波長は共に10．5〜1ユ．5μm），図8は図7の約4時間後の画像である．図9にNO

AAの温度からV

T

IRの温度を引いた差の画像を示す．この図で明るいところは4時間の問に温度の上昇したところを，暗いところは温度の下降したところを示している．

図7

V

TI Rバンド3から求められた温度分布（関東周辺）（観測日時：1987．6．ユ710：20（JST）頃）

Fig．7

Temperature

distribution

inandband−3at10：20（JST），17．

aromd

theKanto

Jun．ユ987

＿236一

area

due

tobrightness

of

VTIR

広域衛星画像の地理的変換

リモートセンシング画像処理基本ソフトウユアの開発

図﹇8g8

図9Fig．9

NOAA

矢崎・建石

AVHRRチャンネル4から求められた温度分布（VTIR（図7））との重ね合わせ画像（観測日時：1987．6．ユ714：1O（JST）頃）
Temperature

distribution

at14：10（JST），17．Jun．1987，transformed

duetobrightness

ofto

NOAA
the

chame1−4

VTIR

image

coordinate

NOAAの温度（6，17，14：1O）とVTIRの温度（6．ユ7，ユ0：20）との差

Temperature

theandtheVTIR

difference
imageofbetween
FNOAA

ig．7237一

imageof

Fig．8

of

Fig・7

国立防災科学技術センター研究報告

5．2

LC第43号

1989年3月

Cどうしの変換

ここでは2つの基準緯度（ア1，乎2）は共通とする．MO

S一のL

C

C画像はすべてこの場合にあてはめる．このときは地図投影面が共通であるから画像座標どうしは原点移動，回転および拡大・縮小の関係，すなわちヘルマート変換の関係になる．したがって，2つのL

C

C画像の画像座標をそれぞれ（u。，v。），（u。，v。）とするとこの変換は

uB＝

auA＋bvA＋c（18）vB＝一buA＋avA＋dおよび

uA＝PuB＋qvB＋r（19）vA＝一qvB＋s

ここで係数a，b，・で与えられる．

auB＋P＝はDA■ COS（∠。一1。）
DBDASin（4パ」。）
b：一DB

c＝

UB−aUA−bVA

d＝VB＋bUA−aVA
（20）

P

DB
＝ 一 COS（1。一1。）
DA

q＝DB
一 Sin （4A一∠B）
DA

r＝

UA−PUB−qVB

s＝

VA＋q

UB−PVB

一238＿リモートセンシング画像処理基本ソフトウェアの開発

I

広域衛星画像の地理的変換一矢崎・建石で与えられる．（20）でサフィックスのAおよびBはそれぞれの画像の（14）で与えられるパラメータを示す．

VT

l

R，87．8．8とMS

R

87．8．8の重ね合せの例

これらの画像の場合g1＝20（。），g2＝50（。）は共通である．3．2と同様にV（画像Aとする）およびM

S

R（画像Bとする）のパラメータを求めると

DA＝7．1662＊lO−5
UA＝一541．75

V＾二一7365．58∠。＝一64．89

（。）

および

DB＝7．8836＊1O−4

UB二一143．64
VB二一357，2ユ∠。＝一64．7ユ（。）が得られる．これらから（20）により

a＝．090900

b＝．O00279
c＝

一92，339d＝312，170

p＝ユ1．00105

q＝一0．03381
r＝s＝1026，377

一3431．075

一239一

TI

R

国立防災科学技術センター研究報告第43号

ユ989年3月が得られる．

図10に88年1月24日のVTI測されたMSR画像をVT

I

R画像（一部）を，図11に図10のVTIR画像の上に変換したものを示す．

図10

VTIRバンド3画像（1988．1．24，日本全土および周辺）

Fig．lO

VTIR

図11

M

Fig．1lMSR

S

R画像と同時に観band−3imageof24．Jan．ユ988，inRバンド1画像（1988．1．24）とVTI

band−1image

andaromd

theJapan

is1ands

R（図10）との重ね合わせ

of24．Jan．1988transformed

一240一

to

theVTIR

coordinateofFig．10

リモートセンシング画像処理基本ソフトウェアの開発

5．3

I

広域衛星画像の地理的変換一矢崎・建石

正方形図法への変換

幾何学的に未補正の画像も含めて，いろいろな地図投影法による広域の衛星画像を処理した結果得られた画像（例えば海面温度分布等）をデータファイルとして蓄積していくようなとき，共通の図法に変換しておいたほうが後の利用のために都合がよい．このような共通の図法の要件としては，画像と地上位置の対応が簡単で，広い範囲で均一な表示がなされ，しかも形や面積の歪のあまり大きくない図法が好ましい．まず，画像と地上位置の対応が簡単で広い範囲で均一な表示がなされているという点では円筒図法がよい．円筒図法のうちで例

えばメルカトール図法は方位を正しく表わすことを目的にした図法であり，したがって方位は正しく表わされるが高緯度で面積が著しく拡大される欠点がある．また面積を正しく表わす図法である等積図法では高緯度で形が著しく歪む．このようにすべてを正しく表わすことができない以上それぞれ一長一短があるが，ここでは正方形図法をとりあげる．

正方形図法は一定の経緯度差を地図上の正方形に対応させるもので，一定の経緯度差の領域が画像の一つのピクセルに対応するので対応関係がきわめて簡単である．また形や面積の歪みもそれほど大きくはない．

この図法では地球座標と画像座標の変換は簡単に

λλOu：十ユD
（21）

V＝90

9十 1Dおよび

λ＝λo＋D（u−1）（22）

甲二ψo−D（v−1）
で与えられる．ここでλo，乎oは画像座標原点の経度，緯度，Dは経緯度で言十ったピクセルの大きさである．

図12に正方形図法による表示例としてV

TI

Rの画像を変換したものを示す．この図では

λo：11O．E，go＝60oN，D＝0．1。である．

＿241一

国立防災科学技術センター研究報告

図12

正方形図法によるVTI

Fig．12

Exampleofthe第43号

1989年3月

Rバンド3画像（1987．6．17）
squaremap

projection

appliedtotheVTIRAAで用いられているメルカトール図法とLCband−3image

6．まとめ
MO

S−1およびNO

C図法に対して，画像座標（ピクセル）と地球座標（経度，緯度）との変換式およびそれに含まれるパラメーターの値の決定のアルゴリズムを求めた．また，その変換式を用いて上記の図法に投影された衛星画像に等経度，等緯度線および海岸線の重ね合わせ，異なる地図投影法による衛星画像間の変換および共通地図座標系としての正方形図法への変換を行った．

これらの手法の開発によりMO

S−1およびN参考OA文

Aの画像を統一的に扱えるようにした．

献

地図投影法の公式は宇宙開発事業団の内部資料より引用した．

（1988年12月19日

＿242一

原稿受理）

リモートセンシング画像処理基本ソフトウェアの開発

I

広域衛星画像の地理的変換一矢崎・建石

付録f−1の計算
f（g）の逆関数f■1（t）は解析的にexplicitな形には表わせないので逐次計算により値を求めなければならない．

（Aユ）

甲：f−1（t）とおくと，定義から

t−f（・）一刈÷・ξ

1−eSin？

1＋eSin甲且2（A2）

これを変形してH舳／t1＋eSin？

1−eSin乎eτ／一÷（A3）となる．この形のままで逐次計算してもよいが，三角関数の計算を多数回行なうことになるので次のようにしたほうが計算時問が少なくてすむ．

a＝e（A4）Sin9

1一…t−／t1＋a1−ae㍉（A5）でa， θを定義すると（A3）はπ（A6）ア＝2θ一一2

一243＿国立防災科学技術センター研究報告となる．

＝

1989年3月（A6）から

（A4），

a第43号

一ecos2（A7）

θとなるが，ここで三角関数の公式

1−tan2θcos2θ＝1＋tan2θを使うと（A7）はまた

tan2θ一1a＝（A8）

etan2θ十1と書ける．さらに（A9）b＝tan2θでbを定義すると， （A8）はb一ユa＝（A1O）

eb＋1と書ける．一方，（A5），（A9）から・一t・［ 1＋ae（Al1）1−aが得られる．
（A10），（Aユ1）はa， bに対する連立方程式になっているから，適当な初期値から始＿244一

リモートセンシング画像処理基本ソフトウェアの開発

I

広域衛星画像の地理的変換一矢崎・建石

めて（A1O）， （A11）を交互に求めていけばよい．漸化式の形に書けば例えば

ao＝O1＋ai−1

bi＝t2
ユーai−1（A12）

bi−1
ai−ebi＋1となる．これをn回繰り返して得られるa。から， （A4）によりψ（t）の第n近似として

a？n（t）＝n（A13）arCSineが求まる．

ここで（A12）の収束性を調べておく必要がある．収束残差9。（t）一ア（t）のtの値のすべて定義域（0≦t＜oo）に亘る最大値をd。，すなわちd、（t）＝maxt1甲、（t）一g（t）1（A14）として，各nにっいてd、を求めると，表A1のようになった．例えば2回の繰り返しで収束残差は地上距離にして80m位である．したがって，MOS−1やNOAAの場合は地上分解能からみてn二2で十分であることが分かる．また，どのような場合にもn＝3で十分であることも分かる．

最後に図A1に（A12），（A13）を用いてf−1（t）を求める関数サブルーチンの例を示す．

一245一

国立防災科学技術センター研究報告

表A第43号

1989年3月

繰り返し回数と収束残差lTab1eA1． Convergence

residuals

in

terms

of

numberof

iterationnd。（単位：度）

地上距離換算（概数）1234

2x10−120km

8×10i4804×10−640c㎜

3×10−83n藺聰

FWCTI㎝門W（τ）
P＾R＾㎜EτεR

E＝0．0816968

P服㎜ET醐ト2（or3）π；丁球τ＾＝E＊（π一1）／（π十1）D01I＝2，＝a1N

B＝π＊（（1＋＾）／（1−A））＊＊E

1トE＊（ト1）／（Bり）＝aF㎜V・＾SIN（＾／E）

肥Tu則

酬D
r1（t）を求めるFORTRANサブルーチン

図A1
Fig．Al．

A

FORTRANprogramtocompute

一246一function＝bi

r1（t）

i