[dennou-ruby:002444] Re: integrate

[Takeshi Horinouchi <horinout@xxxxxxxxxxxxxxxxxx>]

竹広さん

堀之内です．

返事遅くてすみません．
 
> 平均操作でなく積分用の重みを用いるのに 1 票です. 
> spmodel library ではそのように座標補助変数を定義してます. 
> そうすると gtool4 規約を変えた方がいいのかしら. 
> 
> ガウス重みは積分重みの代用にならない, ということのココロは?
> -90 〜 90 度までの積分には重みとして OK ですよね. 
> 積分範囲が中途半端な場合は違うんだけ? 

x という軸にそっての積分は，\int f(x) dx ですよね．もしも南北の
座標軸が x = sinφ なら，ガウス重みはまさに積分重みなんですが，
恐らく φ を使ってるんじゃないですか？ （しかも単位は度で．）と
なると \int f(φ) dφ の形にあいません（\int f(φ) (cosφ) dφ 
ですね）．もちろんこれは緯度円の長さの重みをかけた積分で，こうす
るのが物理的に妥当なケースが多いのはわかってますが，だからいって，
ガウス重みを（sinφでなく）φに関する積分重みにするのは一貫性に
欠けます．（またもや一貫性の話ですみません．が，実際問題として，
一貫性が崩れるといろんなところにしわ寄せが来るので，気をつけざる
を得ません．）

ところで，もしもある物理量を経度に関し積分した上で緯度に関して積
分するということをすれば，最初の段階で緯度円の長さの重みがかかる
ため，緯度に関する積分は sinφ でなく φ に対してすることで，全
球積分になります（半径の2乗のファクターを除き）．つまり，本来全
球平均をとるときには，全球に関し「積分」した上で 4π なり 
4πr^2 で割るべきところを，手っ取り早く経度にそって平均し，つい
で緯度にそって(重みつき)「平均」をすることで安直に実現しようとい
うのが，gtool3からの「伝統」だと言えます．

ちなみに，これは結構根が深くて，現在の Grid クラスには複数の軸に
またがる「座標系」の概念が入ってないということにも繋がります．つ
まり，各軸が独立で，その間の関係が表現されない．この問題は最初の
設計当初から意識してたのですが，うまいやり方が思いつかずペンディ
ングにしてました．下手にやるとデータ処理の過程ですぐ一貫性がおか
しくなる．

ということで，座標軸が φ の場合は，ガウス重みは，積分の重み用の
名札（前回の私の提案では integ_weight）とは別の名札をつけて頂く
必要があります．

> 補助変数は具体的にどのようにしたら使えるのでしょうか. 
> 簡単な例を示していただけるとありがたいです. 

ax を Axis クラスのオブジェクトとすると，

   ax.set_aux('integ_weight', varray)

とします．varray は1次元VArrayです．それを読み出すときには，

   wgt = ax.get_aux('integ_weight')

などとします．なければ nil が返ります．名前は任意に決められます
ので，きちんと解釈するには規約がいります．まだ特別な名前として決
めたものはないので，integ_weight が第一号になる可能性大です．
なお，gtool4 規約など，個別のファイル形式に応じたルールとは
対応関係がはっきりしてればいいだけで，名前が一致する必要はありま
せん．

> NetCDF 出力効率とか難しいことは考えず, 自分の仕事に必要なものを
> そろえるべく現在走っておりますので適宜添削していただけたらと思います. 

じつは効率をあげるには，これまで議論してきたように重みが軸の補助
変数として取り込まれていればいいんです．その場合コードに何も手を
入れなくても，Grid 出力時にちゃんと書出されるという特典もありま
す．（というか，始めからこういうことを想定して Axis や Grid が作っ
てあります．）


> > >  * 現在の Gphys::mean メソッドはもしかして Gphys::Axis の 
> > >    average メソッドを使っているわけではないのでしょうか? 
> > 
> > お察しの通り mean と average は別物です。mean は NArray の
> > mean を使うので、重みも何もありません。average のほうは、
> > サイクリック対応してません。gphys/axis.rb に入っている
> > average のコード以外には何も存在しません --- ということが
> > わかれば分かりますよね。
> 
> 送ってから気づいたですが, サイクリック対応は
> 積分重みをあたえてあげれば単純に対応できますね. 
>  
> > >    Gphys.cut して mean をとるとエラーにならないけど, 
> > >    Gphys.cut して integerate をとると軸とデータの数が違うという
> > >    エラーが出てしまいます. 
> > > 
> > > cut して mean を使うとそのようなエラーが出ないので別物なのかな? 
> > > と思ったわけでして...
> > 
> > ほんとですか。（そのうち調べてみます。）
> 
> これは重みをちゃんと同じ長さであたえてなかったからかもしれません. 
> もう一度確認してみます. 
> --