スレッドの同期

突然ですが、こんなソースがあったとしたら、
結果に何が表示されるか分かりますか？

volatile int g = 0; // メモリから読むことを強制

WORD WINAPI tp( LPVOID v ) {
    for( int i = 0; i < 10000; i ++ ) g++;
    return 0;
}

int WINAPI WinMain(HINSTANCE hi, HINSTANCE hpi, LPSTR command, int show) {
    const int tmax = 64;
    HANDLE th[ tmax ];
    for( int i = 0; i < tmax; i ++ ) {
        DWORD tid;
        th[ i ] = CreateThread( 0, 0, (LPTHREAD_START_ROUTINE)tp, 0, 0, &tid );
    }
    // 終了待ち
    WaitForMultipleObjects( tmax, th, TRUE, INFINITE );

    // 止めてgを確認
    __asm int 3;
    return 0;
}

tpはgを10000回インクリメントするだけ。
WinMainはgを実行するスレッドを64個作るだけです。

やってみれば分かりますが、答えは「環境によって不定」です。
マルチスレッドを知らない人は「640000」だと思ったんじゃないでしょうか。

VC++2008付属コンパイラだと、これを

mov eax, 0x2710 // 10進で10000
mov ecx, 1
loop:
add [g], ecx
sub eax, ecx
jne loop
ret

こんな感じでコンパイルしてくれますが、
マルチコアCPUだとgへのaddが同時に起こることがあります。
その場合1回のaddと同じ結果になるので、
競合したaddの結果が失われることになり、結果としてインクリメントが失敗します。
つまり正確な答えは、「マルチコア/CPUの場合640000以下、それ以外は640000」
となります。

さて、これをアトミック（重複せずに、1つ1つ、等という意味）に処理しなければ
正常な結果は得られません。これを同期処理と言います。

同期にはミューテックス、セマフォなどなどありまして、
MicrosoftがWin32APIとしていろいろ提供しています。
その中でもとても単純で応用力のあるAPIが
InterlockedExchangeという関数。以下MSDNより。

InterlockedExchange
指定された 1 個の変数の内容ともう 1 つの値の交換を一括して行います。
この関数は、複数のスレッドが同じ変数を同時に使うことを防止します。

LONG InterlockedExchange(
LPLONG Target, // 交換に使われる変数
LONG Value // 新しい値
);
戻り値：Target パラメータが指す変数の、交換前の値が返ります。

この関数は同時に実行されないことが保証されているので、
こんな感じで同期処理が取れます。

long flag = 0;
WORD WINAPI tp( LPVOID v ) {
    for( ;; ) {
        if( InterlockedExchange( &flag, 1 ) == 0 ) break;
        Sleep( 1 );
    }

    for( int i = 0; i < 10000; i ++ ) g++;
    InterlockedExchange( &flag, 0 );
    return 0;
}

例としてスレッド1,2がtpを実行するとし、
スレッド2が先にInterlockedExchangeに入ったとします。

1.スレッド1はInterlockedExchangeで待機状態になります。
2.スレッド2がInterlockedExchangeから戻り値0を取得し、処理続行します。
3.スレッド1がInterlockedExchangeから戻り値1を受け取ります。
4.0でなかったので、一瞬眠って再度InterlockedExchangeを試行します。
　※）スレッド2がflagを0にするまで3,4を繰り返します。
5.スレッド2がflagを0にします
6.スレッド1がInterlockedExchangeから戻り値0を取得し、処理続行します。

とまあこんな流れになります。ややこしいですけど。

ミューテックスにしろセマフォにしろ、結局は
同時に処理してはいけない状況への対策なので、
上記処理を応用すれば実現できます。たぶん。

余談）
類似関数としてInterlockedExchangeAddってAPIがあるんですが、
MSDNでの説明が酷いんですよね。以下MSDNより。

InterlockedExchangeAdd
加数変数への増分値の原子加算を実行します。

これで意味が分かった人はエスパーだと思います。

拡大縮小

拡大縮小といえば画像処理関係では避けて通れない話題。
きれいに計算する方法は昔からいろいろ議論されてます。

拡大と縮小では計算方法…というか、考え方が変わるので、
（拡大は元データからの穴をいかに補完するかが重要ですし、
　縮小は元データをいかに余さず反映するかが重要です）
これを両方できて、かつきれいなアルゴリズムというのは
結構少なかったりします。

で、そのアルゴリズムの中でも、かなり単純で、実装も簡単で、
ほぼ最速で、縮小も拡大もできる便利なアルゴリズムを紹介します。

考え方としては「転送先画像」を主体として、
転送先画像のピクセルpd(dx,dy)に対して、
最も近い転送元画像のピクセルps(sx,sy)を求め、
それを転送先画像のピクセルとする、となります。

日本語下手ですいませんが、要するにソースで書けばこうなります。

vbpicture vp( L"c:\\a.jpg" );
int sw = vp.W, sh = vp.H;
long *ps = vp.GetpData();

int dw = 600, dh = 800;
vbpicture vd( w, h );
long *pd = vd.GetpData();

// 単純拡大縮小アルゴリズム
for( int y = 0; y < h; y ++ ) {
    for( int x = 0; x < w; x ++ ) {
        int sx = x * sw / dw;
        int sy = y * sh / dh;
        pd[ x + y * w ] = ps[ sx + sy * vp.W ];
    }
}

BitBlt( GetDC( 0 ), 0, 0, vd.W, vd.H, vd.hDC, 0, 0, SRCCOPY );

確かニアストレイバーだとかなんとか言うような気がします。
一応自分で考えたアルゴリズムなんで、正式名称とかよく分かりません。

ちなみに実行速度はCeleron2.2GHz 800x600x32 -> 600x800x32で
BitBltまで含めて3.2msecほど。
StretchBltの高速モード（NOT HALFTONE）と同じぐらいかな。
あれよりマシな画質なはずですし、そこそこ使えるのではないでしょうか。

計算コスト的には、転送先へのシーケンシャルアクセスのコストは
どんなアルゴリズムだろうが絶対に避けられませんので、
まあそれなりに優秀だと思います。

自動レベル補正

普通は全画素をYUV色空間に変換して、輝度のminとmaxを求めて、
各画素pに対して(p-min)*最大輝度/(max-min)みたいな計算をして、
RGB色空間に直して再セット、といったアルゴリズムになります。

ただそれだとノイズだとかでうまく補正できなかったりするので、
上下ある程度は無視するような処理を入れたりもします。
厳密にやるなら全輝度のソートとかするといい感じなんですが、
面倒だしそこまでの画質はいらん、ってなことでさくっと作ったのがこちら。

void PictBrightAdjust( vbpicture *pp ) {

    // 中のレベルを大雑把にまとめる
    const int nsmax = 256;
    int ns[ nsmax ];
    memset( ns, 0, nsmax * sizeof( int ) );

    // 初期化
    int r,g,b;
    long *pd = pp->GetpData();

    // 一度輝度を調べる
    for( int y = 0; y < pp->H; y ++ ) {
        for( int x = 0; x < pp->W; x ++ ) {
            long2RGB( pd[ x + y * pp->W ], r, g, b );
            int n = ( r + g + b ) / 3;
            ns[ n ] ++;
        }
    }

    // nsの一覧から、上下5%を除外
    int min, max, pt, co = pp->W * pp->H;
    pt = 0; for( min = 0; min < nsmax; min ++ ) { pt += ns[ min ]; if( pt >= co * 0.05 ) break; }
    pt = 0; for( max = nsmax - 1; max >= 0; max -- ) { pt += ns[ max ]; if( pt >= co * 0.05 ) break; }
    if( max - min <= 0 ) return; // min, maxチェック
    
    // 変換
    for( int y = 0; y < pp->H; y ++ ) {
        for( int x = 0; x < pp->W; x ++ ) {
            long2RGB( pd[ x + y * pp->W ], r, g, b );
            r = ( r - min ) * nsmax / ( max - min ); if( r > 255 ) { r = 255; } else if( r < 0 ) { r = 0; }
            g = ( g - min ) * nsmax / ( max - min ); if( g > 255 ) { g = 255; } else if( g < 0 ) { g = 0; }
            b = ( b - min ) * nsmax / ( max - min ); if( b > 255 ) { b = 255; } else if( b < 0 ) { b = 0; }
            pd[ x + y * pp->W ] = r << 16 | g << 8 | b;
        }
    }
}

輝度の計算に(R+G+B)/3とか使ってます。あほです。
でも適当ならこれでいいです。自然画ならそれなりに見栄えもします。

ただ見栄えはそれなりなんですが、問題は意外（当たり前？）にも遅かったこと。
CPU2GHzのCeleronで800x600x32の画像処理に30msec近くかかってました。

HALFTONEのStretchBlt並に遅いとかどんだけ無駄処理なのよってことで
さくっと高速化してみました。（変換のとこだけ）

    // 変換
    int ptrmax = pp->H * pp->W;
    int nstemp = ( nsmax << 12 ) / ( max - min );
    for( int ptr = 0; ptr < ptrmax; ptr ++ ) {
        long2rgb( pd[ ptr ], r, g, b );
        r = ( ( r - min ) * nstemp ) >> 12; if( r > 255 ) { r = 255; } else if( r < 0 ) { r = 0; }
        g = ( ( g - min ) * nstemp ) >> 12; if( g > 255 ) { g = 255; } else if( g < 0 ) { g = 0; }
        b = ( ( b - min ) * nstemp ) >> 12; if( b > 255 ) { b = 255; } else if( b < 0 ) { b = 0; }
        pd[ ptr ] = r << 16 | g << 8 | b;
    }

まずptrを計算するのをやめ。どうせy,xの要素は使わないんでインクリメントで代用しました。
あと割り算はやっぱり超遅いので、掛けて割る分をまとめてキャッシュ。
割り算をキャッシュする場合、整数演算ですから普通はオーバーフローしちゃうんですけど
このアルゴリズムだと最大で255*255の16ビットまでしか使わないので、
残りの12ビット（64ビットCPU限定ならもっと使えますが、まあそこまでいらないし）を
オーバーフロー対策として使用してます。

それなりに効果はあったみたいで、12〜13msecぐらいまで削れました。
さくっと高速化するならあとはMMX使うのがいいですね。桁あふれ対策もあるし。