照度分布計算その7 - 補助クラスの実装

回折計算によるレーザ光学系の照度分布の計算の実装。前回は計算の中心となる複素数2次元配列のクラスをやった。まあこんなもんだろ。次はこのクラスのインスタンスを初期化するためのクラス。

まず、SAmplitudeFunction。これは、例えばレーザビームの中の振幅分布を返したり、絞りの、半径Rの内側は透過率1でその外は透過率0というような透過率分布を返すためのクラス。 まず、抽象クラス

@interface SAmplitudeFunction : NSObject {
    double  shiftv;
    double  shifth;
}

- (id)initWithVerticalShift:(double)shiftVertical
         andHorizontalShift:(double)shiftHorizontal;

- (void)setVerticalShift:(double)shiftVertical
      andHorizontalShift:(double)shiftHorizontal;

- (complex)amplitudeAtVertical:(double)v andHorizontal:(double)h;
- (double)intensityAtVertical:(double)v andHorizontal:(double)h;
@end

-amplitudeAtVertical:andHorizontal:メソッドは直交座標上で指定された位置での振幅の値を返す。SAmplitudeFunctionそのものはどの場所も0を返すだけ。

具体的な、例えばGaussian振幅分布を返すSGaussianは、このサブクラスとして実装する。

@interface SGaussian : SAmplitudeFunction {
    double  av;
    double  ah;
}
- (id)initForLaser;
- (id)initGaussianFWHMVertical:(double)fwhmVertical
                 andHorizontal:(double)fwhmHorizontal
             withVerticalShift:(double)shiftVertical
            andHorizontalShift:(double)shiftHorizontal;
@end

Gaussianの縦横のFWHM（Full Width of Half Maximum）の値を与えて初期化する。この実装は

static double   gcoef = 1.1774100225154746910;//sqrt(2.0 * log(2.0));
@implementation SGaussian
- (id)initGaussianFWHMVertical:(double)fwhmVertical
                 andHorizontal:(double)fwhmHorizontal
             withVerticalShift:(double)shiftVertical
            andHorizontalShift:(double)shiftHorizontal
{
    self = [super initWithVerticalShift:shiftVertical
                     andHorizontalShift:shiftHorizontal];
    av = gcoef / fwhmVertical;
    ah = gcoef / fwhmHorizontal;
    return self;
}

中身はただ設定するだけ。親クラスの-amplitudeAtVertical:andHorizontal:を上書きして

- (complex)amplitudeAtVertical:(double)v andHorizontal:(double)h
{
    double  cv = (v + shiftv) * av;
    double  ch = (h + shifth) * ah;
    return complexBy(exp(- cv * cv - ch * ch), 0.0);
}

で具体的な値を計算する。座標の値はmmを使うことにする。

普通はコントローラクラスがユーザからのパラメータ設定を読んで、このinit*を使って設定するのが当たり前のやり方だけど、今回のような急いで計算したいときはまだるこしい。パラメータの設定と読み出しのためのシングルトンを作ったので、自分でかってにそれを読んで自力で初期化するようにした。たとえばこのクラスはレーザのビームの振幅分布を読む

- (id)initLaser;

というのを準備する。これは

- (id)initForLaser
{
    SParameters *par = [SParameters currentParameters];
    double      fwhmv = [[par querryFor:kFwhmv] doubleValue];
    double      fwhmh = [[par querryFor:kFwhmh] doubleValue];
    return [self initGaussianFWHMVertical:fwhmv
                            andHorizontal:fwhmh
                        withVerticalShift:0.0
                       andHorizontalShift:0.0];
}

というように、SParametersのシングルトンをアクセスして自分の必要なパラメータを読み出してそれで初期化する。こうすればコントローラクラスがパラメータの詳細を知る必要がなくなって簡単になる（はず）。

他のアパチャやピンホールの透過率分布も全く同じにする。例えばピンホールは

@interface SCircularAperture : SAmplitudeFunction {
    double  radius2;
}
- (id)initForPinhole;
- (id)initCircularAperture:(double)CircluarRadius
         withVerticalShift:(double)shiftVertical
        andHorizontalShift:(double)shiftHorizontal;
@end

みたいなの。中身は例えば

@implementation SCircularAperture
- (id)initCircularAperture:(double)CircluarRadius
         withVerticalShift:(double)shiftVertical
        andHorizontalShift:(double)shiftHorizontal
{
    self = [super initWithVerticalShift:shiftVertical
                     andHorizontalShift:shiftHorizontal];
    radius2 = CircluarRadius * CircluarRadius;
    return self;
}

- (complex)amplitudeAtVertical:(double)v andHorizontal:(double)h
{
    double	r2 = (v - shiftv) * (v - shiftv) + (h - shifth) * (h - shifth);
    if (r2 <= radius2)
        return complexBy(1.0, 0.0);
    else
        return complexBy(0.0, 0.0);
}

みたいなの。あともいっしょ。

それから、SCoordinateConverter。さっきのSAmplitudeFunctionが自分でかってに決めた単位で位置を指定しなければいけない。例えばレーザの振幅分布は光軸を原点にしたミリメータ単位の直交座標で指定する。またピンホールは同じだけど単位はミクロンにするのが便利。そこでこのクラスは前回の配列の中身を初期化するときに、横h番目、縦v番目の配列のインデクス値から、SAmplitudeFunctionが使う具体的な長さの座標に変換するための座標変換用のクラス。

@interface SCoordinateConverter : NSObject {
    int     iRange;
    double  indexUnitv;
	double  indexUnith;
	double  shiftv;
	double  shifth;
}
- (id)initWithVerticalUnit:(double)vUnit
            horizontalUnit:(double)hUnit
             verticalShift:(double)vShift
           horizontalShift:(double)hShift
             andIndexRange:(int)indexRange;

- (double)horizontalFromIndex:(int)index;
- (double)verticalFromIndex:(int)index;
- (int)horizontalIndexFrom:(double)hor;   // it's for debugging
- (int)verticalIndexFrom:(double)ver;   // it's also for debugging
- (double)horizontalFractionalIndexFrom:(double)hor;
- (double)verticalFractionalIndexFrom:(double)ver;
@end

これもinit*は定数を設定するだけ。具体的な計算はメソッドでやる。たとえば

- (double)horizontalFromIndex:(int)index
{
    if (index <= iRange / 2)
        return index * indexUnith - shifth;
    else
        return (index - iRange) * indexUnith - shifth;
}

あとも同じ。

このクラスにも自分でかってにやるinit*を作っておく。

- (id)initForLaser;
- (id)initForAperture;
- (id)initForPinhole;
- (id)initForSubstrate;

実装はたとえば

- (id)initForLaser
{
    SParameters	*par = [SParameters currentParameters];
    double      mag = [[par querryFor:kMagnificationOfRelay] doubleValue];
    double      rshifth = [[par querryFor:kRelayShifth] doubleValue];
    double      rshiftv = [[par querryFor:kRelayShiftv] doubleValue];
    double      psiz = [[par querryFor:kPupilAreaSize] doubleValue];
    int         frange = [[par querryFor:kFftSize] intValue];
    double      ipitch = psiz / frange / mag;
	
    return [self initWithVerticalUnit:ipitch
                       horizontalUnit:ipitch
                        verticalShift:rshiftv
                      horizontalShift:rshifth
                        andIndexRange:frange];
}

というようなもの。後も同じ。キー入力が面倒なだけ。

SAmplitudeFunctionとSCoordinateConverterの両方に縦横シフトのインスタンス変数がある。これは片方だけでいい。でも場合によって片方が便利でもう一方は面倒になったりする。それはケースバイケースなんだけど、いくつも使って、それぞれが違うやりかただと間違いの元なので、とりあえず両方作っといて、最終的には見通しのいい方に統一する。一見無駄な作業だけど、書いてる時点でどっちがいいか判断できない場合はこうしといて、都合のいい方だけ使って悪い方は0.0を設定する。この時点で労力を惜しむより、最終的に見通しよくなるようにしたほうが結局は全体の労力を減らすことができる。この例はそれほど大きな違いではないけど、何でも一事が万事。

これで例えばレーザの振幅分布のインスタンス生成は

    SGaussian   *lfunc = [[[SGaussian alloc] initForLaser] autorelease];
    SCoordinateConverter *lconv = [[[SCoordinateConverter alloc] initForLaser] autorelease];
    SAmplitudePlane *laser = [[[SAmplitudePlane alloc]
                               initWithAmplitudeFunction:f
                                     coordinateConverter: lconv
                                                  ofSize:sSize] autorelease];

でいい。一見、文字が多くてめんどくさそうだけど、決まりきったことを書いているだけ。こうしておけばここでバグが発生することは少ない、と見なせる。デバグのときは、それぞれの中身を個別にチェックすればいい。ファイルをあっち行ったりこっち行ったりする必要が少なくなる。これ、重要。

くわー。こんなに丁寧にメモ残してたらコードだけ書くのにくらべて何倍も労力がいる。こんな分かり切ったところで時間かけてるのもったいないがな。さっさといこ。でも、もうおねむ。また明日。