[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[dennou-ruby:002268] Re: Ruby-ISPACK 作成中



北村です.

堀之内さんに言ったら作業変数を隠蔽するような上位のラッパーもほ しい
っていうことですのでそれもおいおいできればいいかなと思います。

自分が利用しているN2PACKの一部だけですが,ラップしてみまし た.
作業変数の隠蔽と,グリッドデータ,スペクトルデータが1次元 配列でなくても
動くことを目標にしています.以下に実装例とテストプログラムを添付 します.


サンプルでは,スペクトル変換のためのオブジェクトを作成したら,そ のオブ
ジェクトにデータを引数として与えれば,諸々の演算を行えるようにし てあり
ます.テストプログラムでは,ガウシアンを与えてスペクトル空間で x微分して
います.それなりにspmodelの雰囲気に近づいてはいると思いま す.
ちなみにx微分,y微分はN2PACKのラップではなく, 自前で用意したものです.


どういう設計にするのがいいのかがよく分からないので,もっと適切な 方法が
ありましたら教えていただければ幸いです.


require "narray"
require "ispack_raw"

class N2Pack
  def initialize(im, jm, nx, ny)
    @im, @jm, @nx, @ny = im, jm, nx, ny
    @itj = NArray.int(5)
    @iti = NArray.int(5)
    @tj = NArray.float(2*jm)
    @ti = NArray.float(2*im)
    @wkfft = NArray.float(im*jm)
    @itj, @tj, @iti, @ti = ISPACK::n2init(@jm, @im)
  end

  def forward(grid)
    shape = grid.shape
    grid.reshape!(grid.total)
    obj = ISPACK::n2g2sa(@ny ,@nx, @jm, @im, grid,
                 @wkfft, @itj, @tj, @iti, @ti)
    grid.reshape!(*shape)
    return obj
  end

  def backward(spectrum)
    shape = spectrum.shape
    spectrum.reshape!(spectrum.total)
    obj = ISPACK::n2s2ga(@ny ,@nx, @jm, @im, spectrum,
             @wkfft, @itj, @tj, @iti, @ti)
    spectrum.reshape!(*shape)
    return obj
  end

def dx(spectrum)
tmp = spectrum.reshape(2*@ny+1, 2*@nx+1)
factor = -(NArray.float(1,tmp.shape[1]).indgen! - (tmp.shape[1] - 1)/2)
tmp = tmp[-1..0, -1..0]
return (tmp*factor).reshape!(tmp.total)
end


def dy(spectrum)
tmp = spectrum.reshape(2*@ny+1, 2*@nx+1)
factor = -(NArray.float(tmp.shape[0],1).indgen! - (tmp.shape[0] - 1)/2)
tmp = tmp[-1..0, -1..0]
return (tmp*factor).reshape!(tmp.total)
end


end

=begin
###
### test program
###

require "numru/dcl"

include NumRu
include NMath

## グラフの描画
def draw(grid)
  DCL.grfrm
  DCL.grswnd(0.0, 1.0, 0.0, 1.0)
  DCL.uspfit
  DCL.grstrf

  DCL.ussttl('X', '', 'Y', '')
  DCL.uelset('ltone', true)
  DCL.uetone(grid)
  DCL.usdaxs
  DCL.udcntr(grid)
end

## 初期設定
i = 64
j = 64
nx = 31
ny = 31
x = NArray.float(1,i).indgen!/i
y = NArray.float(j,1).indgen!/j
grid = exp(-100.0*((x-0.5)*(x-0.5) + (y-0.5)*(y-0.5)))

DCL.gropn(1)

draw(grid)
## Fourier  変換 + 演算
n2pack = N2Pack.new(i, j, nx, ny)
grid = n2pack.backward(n2pack.dx(n2pack.forward(grid))).
       reshape!(j, i).transpose(1,0)

draw(grid)

DCL.grcls
=end

----
Yuji Kitamura, Ph.D. <kitamura@xxxxxxxxxxxxxxxxxx>
Division of Earth and Planetary Sciences,
Graduate School of Science, Kyoto University
Tel: +81-75-753-3933/Fax: +81-75-753-3715