#########################################################
#
#   GT2DLIB  2-Dimensional Grid Transformation LIBrary
#            2次元格子の座標変換
#
#########################################################

● 概要

本ライブラリーは2次元の補間、並びにそれを用いた2次元の座標変換を扱う。
座標間の対応は離散的な格子点で指定し、間は補間する。補間アルゴリズムは
双線形補間を用いる。このため逆変換は解析的に計算される。将来的には双3
次補間もサポートする可能性があるが、その場合逆変換は反復法によることに
なる。

補間のサブルーチンは、座標変換のみならず一般の用途に適するよう作ってある。

座標変換は、(ξ,η) -> (x,y) において、離散的な点 ξ_i (i=0,1,2,..),
η_j (j=0,1,2,..) において

      x_i,j = x(ξ_i,η_i)
      y_i,j = y(ξ_i,η_i)

という形で表される格子点での対応関係により定義される。従って「正変換」
においては、変換元の座標のデータは2つの1次元配列、変換先の座標のデータ
は2つの2次元配列となる。逆変換はその逆である。なお、2つの2次元配列から
2つの2次元配列へと対応づけられる座標変換は、中間的に1次元配列ベースの
「矩形」格子を設定すれば、逆変換、正変換の2段階で行うことが出来る。こ
れは、幾何的には、基準となるデカルト座標系を用いることに相当する。

● サブルーチンのリスト

* 双線形補間 (BLI: bilinear interpolation)

  G2FBLI   双線形補間を行う

* 2変数の双線形補間 (BL2)

  G2FBL2  2変数の双線形補間を行う
  G2IBL2  2変数の双線形補間(G2FBL2)の逆演算

  注： 2変数 --> 2変数だから一般に逆演算が存在する

* 座標変換 (CTR: coordinate transform) 
  格子点間の対応に基づく座標変換. BL2を利用.

  G2SCTR  変換設定 (初期化)
  G2FCTR  変換
  G2ICTR  逆変換
  G2QCTR  初期化済か否かの問合わせ. 

* 内外判定

  LG2INQ  四辺形の内外判定

* 下請

  G2SGRD  格子データの蓄積
  G2QGRD  指定した位置のグリッドボックス（四辺形）の問合わせ
          （g2sgrd, g2qgrd は C で書いてある）

----------------------------------------------------------------------
* 双線形補間
----------------------------------------------------------------------

● SUBROUTINE G2FBLI

1. 機能

  双線形補間. 入力座標 (P, Q) はともに [0, 1] の区間に正規化されたものをとる.

2. 呼出し方法

      CALL G2FBLI(P, Q, Z00, Z10, Z01, Z11, Z)

3. パラメターの説明

      P, Q (R)  補間すべき点の座標. ４隅が (0,0), (1,0), (1,1), (0,1) 
         となるよう正規化されてるものとする. (入力)

      Z00, Z10, Z01, Z11 (R)  P,Q の４隅(上述)、即ち、それぞれ (P,Q) =
         (0,0), (1,0), (1,1), (0,1) における変数の値 (入力) 

      Z (R) 補間結果 (出力)

4. 備考

  (a) 中身は

         Z = (1-P)*(1-Q)*Z00 + P*(1-Q)*Z10 + (1-P)*Q*Z01 + P*Q*Z11

      である(一行のみ).

  (b) P, Q とも [0,1] の区間内にある場合は補間, そうでなければ外挿
      となる。

  (c) パラメターの説明から明らかであろうが, (P,Q) =
      (0,0), (1,0), (1,1), (0,1) の場合、出力はそれぞれ
      Z00, Z10, Z01, Z11 に等しくなる。

----------------------------------------------------------------------
* 2変数の双線形補間
----------------------------------------------------------------------

● SUBROUTINE G2FBL2

1. 機能

  2変数同時の双線形補間. 入力座標 (P, Q) はともに [0, 1] の区間に正規化
  されたものをとる.

2. 呼出し方法

      CALL G2FBL2(P,Q, X00,X10,X01,X11, Y00,Y10,Y01,Y11, X,Y)

3. パラメターの説明

      P, Q (R)  補間すべき点の座標. ４隅が (0,0), (1,0), (1,1), (0,1) 
         となるよう正規化されてるものとする. (入力)

      X00, X10, X01, X11 (R)  P,Q の４隅(上述)、即ち、それぞれ (P,Q) =
         (0,0), (1,0), (1,1), (0,1) における一つめの変数の値 (入力) 

      Y00, Y10, Y01, Y11 (R)  P,Q の４隅(上述)、即ち、それぞれ (P,Q) =
         (0,0), (1,0), (1,1), (0,1) における二つめの変数の値 (入力) 

      X (R) X00, X10, X01, X11 に関する補間結果 (出力)

      Y (R) Y00, Y10, Y01, Y11 に関する補間結果 (出力)

4. 備考

  (a) G2FBLI を2回呼んでいるだけである。よって

         X = (1-P)*(1-Q)*X00 + P*(1-Q)*X10 + (1-P)*Q*X01 + P*Q*X11
         Y = (1-P)*(1-Q)*Y00 + P*(1-Q)*Y10 + (1-P)*Q*Y01 + P*Q*Y11

  (b) 2変数から2変数への変換であるから逆変換が存在する。
      逆変換を行うサブルーチンは G2IBL2 である。

● SUBROUTINE G2IBL2

1. 機能

  G2FBL2 (2変数同時の双線形補間) の逆変換. 解析解である.

2. 呼出し方法

      CALL G2IBL2(X,Y, X00,X10,X01,X11, Y00,Y10,Y01,Y11, P,Q)

3. パラメターの説明

      X,Y, X00,X10,X01,X11, Y00,Y10,Y01,Y11 (R) (入力)

      P, Q (R) (出力)

      それぞれの説明は G2FBL2 を参照のこと.

4. 備考

   (a) 連立２次方程式を解析的に解く. (実質的に１次方程式になる場合等、
       多くの場合分けを含む)

----------------------------------------------------------------------
* 座標変換
----------------------------------------------------------------------

● SUBROUTINE G2SCTR

1. 機能

  2次元の座標変換の設定(初期化). 格子点間の対応で定義する.

2. 呼出し方法

      CALL G2SCTR(NX, NY, UXA,UYA, TXA,TYA)

3. パラメターの説明

      NX, NY (R) 配列の長さ (入力)

      UXA (R) 変換元の座標 (独立変数その1). 長さ NX の1次元配列 (入力).

      UYA (R) 変換元の座標 (独立変数その2). 長さ NY の1次元配列 (入力).

      TXA (R) 変換先の座標 (独立変数その1). 長さ NX*NYの2次元配列 (入力).
              第1要素(TX(1,1))を RUNDEF とすることで省略できる. 省略の場合, 
              TX(I,J)=UX(I) が用いられる.
      TYA (R) 変換先の座標 (独立変数その2). 長さ NX*NYの2次元配列 (入力).
              第1要素(TY(1,1))を RUNDEF とすることで省略できる. 省略の場合, 
              TY(I,J)=UY(J) が用いられる.

4. 備考

  (a) 2次元配列 (TXA,TYA) のセーブは、C 言語で書かれた下請の関数が行う。
      動的割り当てによるため、サイズに制限がない。一方、一次元の UXA,
      UYA は静的に確保されたサブルーチン内の配列にセーブされる。最大の
      長さは GRPH2 の UWPACK と同様、コンパイル時に変更可能な定数 
      MAXNGRID である。デフォルト長は 4000。

● SUBROUTINE G2FCTR  

1. 機能

  2次元の座標変換を行う.


2. 呼出し方法

      CALL G2FCTR(UX, UY, TX, TY)

3. パラメターの説明

      UX (R) 変換前の座標（その1） (入力)
      UY (R) 変換前の座標（その2） (入力)
      TX (R) 変換後の座標（その1） (出力)
      TY (R) 変換後の座標（その2） (出力)

4. 備考

  (a) 予め G2SCTR で初期化しておかないとエラーを発生する。

  (b) G2SCTR で設定した UXA, UYA の範囲外の入力に対しては、
      警告を出しつつ、外挿を行う。

● SUBROUTINE G2ICTR  

1. 機能

  G2FCTR (2次元の座標変換) の逆変換を行う.

2. 呼出し方法

      CALL G2ICTR(TX, TY, UX, UY)

3. パラメターの説明

      TX (R) 逆変換前の座標（その1） (入力)
      TY (R) 逆変換前の座標（その2） (入力)
      UX (R) 逆変換後の座標（その1） (出力)
      UY (R) 逆変換後の座標（その2） (出力)

4. 備考

  (a) 予め G2SCTR で初期化しておかないとエラーを発生する。

  (b) G2SCTR で設定した TXA, TYA の範囲外の入力に対しては、
      エラーを発生する。

● SUBROUTINE G2QCTR  

1. 機能

  初期化済か否かの問合わせ.
  なお、初期化せずに G2FCTR, G2ICTR を呼べばエラーがでるようになっているの
  で、問合わせは必須ではない。

2. 呼出し方法

      CALL G2QCTR(LINIT)

3. パラメターの説明

      LINIT (L) 初期化済なら .TRUE. まだなら .FALSE. を返す。(出力)

4. 備考

  (a) G2SCTR が既に呼ばれていれば .TRUE. 

----------------------------------------------------------------------
* 内外判定
----------------------------------------------------------------------

● FUNCTION LG2INQ

1. 機能

  ある点が四辺形の内側にあるか判定する。凸四辺形限定。

2. 呼出し方法

      LIN = LG2INQ(TX,TY, TX00,TX10,TX01,TX11, TY00,TY10,TY01,TY11)

3. パラメターの説明

      TX,TY (R) 判定すべき点の座標 (入力)

      TX00,TX10,TX01,TX11, TY00,TY10,TY01,TY11 (R) 四辺形の4隅の座標 (入力)

      戻り値 (L) 内側または境界上にあれば .TRUE. 外側なら .FALSE. を返す.

4. 備考

  (a) G2ICTR が下請として用いる。

  (b) 凹四辺形には対応していない。座標変換においては、
      凹四辺形は多価性のある変換にしか現われない。多価の場合は、
      逆変換はそもそも一意に決まらないので、利用価値が低い。
      将来的には、凹四辺形に対応してもよい。

----------------------------------------------------------------------
* 下請
----------------------------------------------------------------------

● FUNCTION G2SGRD, FUNCTION G2QGRD

G2SCTR の下請であり、動的な配列確保のため C で書かれている。それ以外に
独立な関数とするべき理由はないので説明は割愛する。

