平成
27
年度
学士学位論文
超解像処理を題材とした
プログラム高速化の効果の検証
Verification of Effectiveness of Program Optimization
Using Super-resolution Processing as a Case Study
1160347
畑中 涼
指導教員
松崎 公紀
2016
年
2
月
5
日
要 旨
超解像処理を題材とした
プログラム高速化の効果の検証
畑中 涼
近年の計算機はキャッシュやパイプライン処理などにより高速化を実現している．そのた
め，計算を高速化するには，そのような機能を有効活用する必要がある．
本研究では，超解像化を行う画像処理プログラムを題材にし，最適化手法をいくつか適用
してプログラム高速化を行う．まず，元のプログラムである
Java
プログラムを
C++
プログ
ラムに書き換えた．その後，書き換えた
C++
プログラムに対して各最適化手法を施すこと
で高速化を狙った．評価実験では，元の
Java
プログラムとそのまま書き換えた
C++
プログ
ラム，最適化を施した高速
C++
プログラムの
3
種のプログラムの処理時間を計測し，比較
した．この結果，高速
C++
プログラムの処理時間は
Java
プログラムの処理時間の
27%
にま
で短縮された．また，各最適化手法別で処理時間を計測し，処理時間と書換えコストの観点
から比較した．その結果，題材にしたプログラムでは，コンパイラに定数伝播を促す書き換
えが最も効率的な最適化手法であることを確認した．
以上の結果から，最適化を行うことで処理時間は短くなり，高速なプログラムにすること
が可能であることが確認できた．また，書換えコストを考慮した場合，コンパイラに定数
伝播を促す書換えを行うことで，最も効率的にプログラムの高速化を実現できることがわ
かった．
キーワード
プログラムの書換え，高速化，最適化，画像処理
Abstract
Verification of Effectiveness of Program Optimization
Using Super-resolution Processing as a Case Study
Hatanaka Ryo
Recent computers realizes high performance with the cache and pipelining mechanisms.
Therefore, in order to perform high-performance computing, it is necessary to use such
mecha-nisms effectively.
In this study, I tried to accelerate super-resolution processing program by applying program
optimization. First, I translated the Java program into a C++ program. Then, I tried to accelerate
the translated C++ program by four program-rewriting techniques. I measured the processing
time of the existing Java program, the translated C++ program, and the accelerated program.
As a result, processing time of the accelerated C++ program reduced to 27% of that of the
existing Java program. Furthermore, I measured and compared processing time and rewriting
cost for each optimization technique. As a result, I found that to promote the compiler to do the
constant propagation is the most efficient optimization.
From these results, We can shorten a processing time by performing program
optimiza-tions. In particular to promote the compiler to do the constant propagation is the most efficient
optimization with respect to the rewriting cost.
目次
第
1
章
はじめに
1
1.1
背景と目的
. . . .
1
1.2
論文の構成
. . . .
2
第
2
章
最適化の対象とする超解像処理
プログラム
3
2.1
超解像処理
. . . .
3
2.2
既存の超解像処理プログラム
. . . .
4
第
3
章
高速化のための最適化手法
5
3.1
定数伝播
. . . .
5
3.2
キャッシュメモリの活用
. . . .
5
3.3
条件分岐の削除
. . . .
6
3.4
ループアンローリング
. . . .
7
第
4
章
超解像処理プログラムの高速化
8
4.1
Java
から
C++
への書き換え
. . . .
8
4.2
定数伝播可能にする書換え
. . . .
9
4.3
キャッシュメモリの活用
. . . 10
4.4
条件分岐の削除
. . . 11
4.4.1 new x
及び
new y
に関する条件分岐
. . . 12
4.4.2 lx
及び
ly
に関する条件分岐
. . . 13
4.4.3
x
と
y
に関する条件分岐
. . . 14
4.4.4
条件分岐が必要のない
x
と
y
の範囲
. . . 14
4.5
ループアンローリング
. . . 15
目次
4.6
C++
プログラムのコンパイル
. . . 16
第
5
章
評価実験
18
5.1
実験環境
. . . 18
5.2
元のプログラムとの比較
. . . 18
5.2.1
実験内容
. . . 18
5.2.2
実験結果と考察
. . . 19
5.3
各最適化手法による比較
. . . 19
5.3.1
実験内容
. . . 19
5.3.2
実験結果と考察
. . . 21
第
6
章
まとめと今後の課題
23
謝辞
25
参考文献
26
付録
A
既存の
Java
プログラム
(mri sr.java)
27
付録
B
高速化
C++
プログラム
(mri sr opt.cpp)
31
図目次
2.1
複数枚の観測画像を用いた超解像
. . . .
3
4.1
各プログラムの計算順序
. . . 11
4.2
条件分岐が不要な範囲
. . . 15
4.3
条件分岐を削除した範囲
. . . 15
5.1
プログラム別の処理時間の差を示したグラフ
. . . 20
5.2
各最適化手法別のプログラムの処理時間の差を示したグラフ
. . . 21
表目次
4.1 Java
プログラムに用意されていた関数
. . . .
8
4.2
定数化した変数
. . . .
9
4.3 g++
に準備されている最適化オプション
. . . 17
5.1
プログラムを実行する計算機
. . . 18
5.2
実行するプログラムについて
. . . 19
5.3
プログラム別の処理時間の計測結果
. . . 19
5.4
実行するプログラムについて
. . . 20
5.5
最適化手法別のプログラムの処理時間の計測結果
. . . 21
第
1
章
はじめに
1.1
背景と目的
近年では，技術の発達により，計算処理能力が高い
CPU
や大容量メモリが開発されてい
る．このような現代の技術動向に伴い，高性能な計算機が普及している．そのため，計算量
が多いアプリケーションでも高速に実行する環境が整っている．計算量が多い処理の一つと
して，画像処理があげられる．画像処理は，膨大なデータを規則的に処理することが多いた
め，一般的に計算コストが非常に高いことが知られている．例えば，縦
1280×
横
960
画素
の画像に対し
1
画素毎に単純な処理を行うだけで，少なくとも
1,228,800
回の処理を行う必
要がある．
画像処理の中に超解像技術というものがある．この技術はデジタルカメラやテレビといっ
た身近なものにも使われており，近年注目されている．
現代の高性能な計算機はキャッシュやパイプライン処理などの機能により，高速で計算す
ることを可能にしている．そのため高速な処理を計算機に行わすためには，それらの機能を
有効活用する必要がある．
本研究では既存の
Java
で書かれた超解像処理プログラムを題材に，様々な最適化手法を
用いて，プログラムの高速化を行う．本研究は，プログラムの処理内容はそのままにして，
書き換えることで処理時間の高速化を行うことを目的とする．
1.2
論文の構成
1.2
論文の構成
本論文の構成は次の通りである．第
2
章では既存の超解像プログラムの説明を行う．第
3
章では高速化のための最適化手法についての説明を行う．第
4
章では既存の超解像プログラ
ムの高速化についての説明を行う．第
5
章では評価実験による書換えコストの観点を踏まえ
た高速化の効果から考察を行う．第
6
章では本論文のまとめを行い，今後の課題を述べる．
第
2
章
最適化の対象とする超解像処理
プログラム
2.1
超解像処理
超解像とは，平たくいうと，
「解像度の低い画像
(
画素数が少ない画像
)
」を，
「解像度の高い
画像
(
画素数が多い画像
)
」にする技術である
[1]
．単に画像を引き伸ばすと，見た目がギザ
ギザの荒い画像になる．そのため，なめらかに画像を拡大するには，元の画像にはなかった
画素の情報が必要になる．周辺の画素の輝度値や色情報から値を推測して補間する
[1]
．近
年，デバイスの表示できる解像度が高くなってきているため，この超解像技術が注目されて
いる．この技術はテレビやデジタルカメラといった身近な製品に搭載されている．
図
2.1
複数枚の観測画像を用いた超解像
2.2
既存の超解像処理プログラム
超解像処理には元となる観測画像が必要となる．
1
枚の観測画像から超解像画像を生成す
る手法はパターンマッチングを利用したものがある．また，複数枚の観測画像から高解像度
画像を生成する手法も存在する
(
図
2.1)
．この手法は画像観測時にサンプルした画素のズレ
を用いて超解像を行う
[2]
．
先述した通り，超解像処理は画像処理の一種であり計算コストが高い．処理に用いる画像
の画素数が多いほど計算量は増加し，処理時間は長くなる．
2.2
既存の超解像処理プログラム
本研究で題材とする
Java
プログラムは
fMRI
画像を対象とした超解像処理プログラムで
ある．また，その手法は複数枚の観測画像を用いるものである．
このプログラムは
CSV
形式に加工された観測画像を読み込み，二次元配列に値を代入す
る．その配列とアフィン変換に用いるパラメータを関数
img_iteration
に渡し，超解像処
理が行われる．関数
img_iteration
では以下の様な処理が行われる．
•
注目する超解像画像の座標
(x, y)
に対応する観測画像の座標
(new_x, new_y)
を求
める
• (x, y)
の
4
近傍画素をダウンサンプリングし，観測画像の画素に対応する画素値を求
める
•
誤差値
dEdp
を元に超解像画像の画素値を更新する
関数
img_iteration
では上記の処理を
200
回繰り返している．これは反復法にって誤差値
dEdp
の最小化問題に取り組んでいるためである．
第
3
章
高速化のための最適化手法
3.1
定数伝播
定数伝播とは，数式内の変数値を既知の定数に変換するコンパイラが行う最適化である．
コンパイラが定数伝播を施した結果，ループの繰り返し回数が静的に定まり，ループアン
ローリングの効果を高めたりする
[3]
などの効果がある．この最適化をコンパイラが施すよ
うにプログラムを書き換えることで，プログラムの実行速度の向上が見込める．
以下の擬似コードを用いて定数伝播による最適化の一例を示す．
int x = 16; int y = 8 - x / 2; return y * (32 / x + 2);このコードに一度定数伝播を適用すると次のようになる．
int x = 16; int y = 8 - 16 / 2; return y * (32 / 16 + 2);定数伝播による定数への変換と同時に定数畳み込みによる定数式の単純化が行われることが
多く，更に単純化される．
3.2
キャッシュメモリの活用
命令処理の際，必要なデータに加えていくつか先までのデータをキャッシュメモリに保存
される．次に実行される命令で必要なデータがキャッシュメモリに存在する場合，それを読
み込むことが可能である．このことを，キャッシュヒットといい，短時間でアクセスするこ
3.3
条件分岐の削除
とが可能になる．これは，キャッシュメモリがメインメモリよりも
CPU
に近い位置にある
ためである．
反対に，必要なデータがキャッシュメモリに存在しない場合もある．これはキャッシュの
ヒットミスといい，必要なデータをメインメモリ等に探しに行く．メインメモリは
CPU
か
らの距離があるためアクセスに時間を要する．そのためヒットミスが生じるほど，無駄なア
クセスが発生し，処理が滞るため，全体の処理時間が遅くなる．
連続しているメモリへのアクセスはキャッシュが効率良く働くが，連続していない離れば
なれなメモリへのアクセスではキャッシュのヒットミスが生じ，アクセスに時間がかかるよ
うになる
[4]
．この現象は，数値計算で頻繁に用いられる行列どうしの計算でよく起こる．
行列の値を格納している配列操作を工夫することで，キャッシュのヒットミスを解消するこ
とができる．
3.3
条件分岐の削除
CPU
には分岐が起こるかどうかを予測して処理を進め，予想がはずれた場合に分岐命令
以下の命令を破棄する
[5]
機能がある．このことを分岐予測という．分岐予測はパイプライ
ン処理における制御ハザードを解消し，性能を高めるための機能である．
条件分岐がある場合，命令と分岐の履歴を持つブランチテーブルバッファにアドレスと結
果を書き込み，次からそれを利用して予測する．この履歴があるため連続して同じ側に分岐
する，または分岐しないとパターン化しているものは，正確に予測することができる．その
ため，ランダムで分岐する場合は履歴が一切推測要素にならないため，分岐予測の成功率は
半々である．
このようにして
CPU
は条件分岐に対して分岐予測をして性能を高めている．しかし，分
岐予測ミスによって逆にオーバーヘッドが生じる場合もある．そのため高速なプログラムに
するためには，不要な分岐命令を削除することが必要になる．
3.4
ループアンローリング
3.4
ループアンローリング
ループを展開することでループ回数を減らし，繰り返しごとに行われる終了条件のテスト
を削除することができる．これをループアンローリングという．
例えば，以下の擬似コード
for ( i = 0; i < 100; i = i + 1 ) { x[i] = y[i]; }を
for ( i = 0; i < 100; i = i + 4 ) { x[i] = y[i]; x[i +1] = y[i +1]; x[i +2] = y[i +2]; x[i +3] = y[i +3]; }と展開すると，判定とジャンプ回数が半分になる．また，
y[i]
と
y[i+1]
などで共通部分
式除去の機会が増える可能性がある
[6]
．
このようなループ展開は
SIMD
の観点からも大きなメリットがある．
128
ビットのデータ
を受け付けるプロセッサで上記のループを計算すると仮定する．配列
x
及び
y
が
int
型の場
合，ループ展開前のループ内の命令は
32
ビットのコピーのみである．つまりプロセッサは
すべての計算を
100
回のクロックで計算する．それに比べ，ループ展開したループ内では
32
ビットの命令が
4
つある．プロセッサは
128
ビットまで受け付けることができるため，
4
つの命令を
1
回のクロックで計算することが可能である．そのため，プロセッサはすべての
命令を
25
回のクロックで計算することができる．
ループアンローリングによってコード自体は長くなるが副次的な効果が大きく，プログラ
ムの高速化が可能である．
第
4
章
超解像処理プログラムの高速化
4.1 Java
から
C++
への書き換え
本研究では
Java
プログラムに最適化を行ったのではなく，一般的に
Java
より高速に実行
が可能と言われる
C++
に書き換えてから最適化を行った．
書き換えを行う上で超解像処理とは関係のない部分でいくつか変更した．元の
Java
プロ
グラムには関数
initialize
と関数
clear
が用意されていた．それぞれの関数を表
4.1
に
まとめた．
関数
initialize
は配列のすべての値を
2047
に初期化する関数である．
C++
では配列を
vector <int > result ( width * height , 2047);
と宣言することで，関数
initialize
と同じ働きをする．このため，関数
initialize
は
不要であると判断し
C++
プログラムには書かなかった．
また関数
clear
は配列
next
に対してのみ利用されているが，
next
は計算過程で使われ
ることがなく，毎回上書きされている．そのため，すべての値を
0
に初期化する必要は無い
と判断し，
C++
プログラムには書かなかった．
表
4.1 Java
プログラムに用意されていた関数
関数名
引数
備考
initialize int
型の配列
配列内すべての値を中央値
(2047)
に初期化
clear
int
型の配列
配列内すべての値を
0
に初期化
4.2
定数伝播可能にする書換え
また，読み込むファイル名を指定する上で
int
型の変数をそのまま
string
型に変換する必
要があった．そのため
C++
プログラムには
int
型の値を
string
型に変換する関数
itos
を追
記した．
string itos ( int i) { ostringstream s; s << i; return s.str (); }
4.2
定数伝播可能にする書換え
本研究では，以下の両方の条件を満たす変数に注目し，定数化を行った．
•
定義された後，値が変化しない変数
•
計算時，頻繁に用いられる変数
これらの条件に適した変数に対して定数化を行い，コンパイラの最適化変換で定数伝播が行
われるようにした．また，今回は関数の引数に対して定数伝播が適用されないことを考慮し
て書換えを行った．
具 体 的 に 変 更 し た 変 数 を 表
4.2
に ま と め た ．元 の プ ロ グ ラ ム は 以 下 の よ う に し て
file_n
は引数の二次元配列の大きさ
(img.size())
から求め，
scale
は
file_n
の平
方根
(sqrt(file_n))
で得ていた．
int file_n = img.size (); // = 4 int scale = (int)sqrt( file_n ); // = 2
表
4.2
定数化した変数
変数名
備考
file_n
観測画像の枚数
scale
超解像画像の倍率
4.3
キャッシュメモリの活用
それらの変数は処理中に値は変化せず，処理中の計算で多く用いられているため，それぞれ
の変数を以下のように定数化した．
const int file_n = FILE_N ; // define FILE_N 4 const int scale = SCALE ; // define SCALE 2
また，配列
palam
は
main
関数で宣言し，関数
img_iteration
に引数として渡してい
た．これを以下のように
img_iteration
内で宣言することで，定数伝播が適用されるよう
にした．
4.3
キャッシュメモリの活用
題材にしたプログラムでは，配列に値を格納し，その値を用いて計算を行っている．この
プログラムは，毎回離れた領域を参照しているため，キャッシュのヒットミスが頻繁に起き
ていることが予想できる．この問題は，ループの順序を逆にすることで解決することができ
る．このプログラムでは，変数
x
と
y
の値の間に依存関係がないため，ループの順序を入れ
替えることが可能である．既存のプログラムでは
for ( int x = 0; x < width ; x++ ) { for ( int y = 0; y < height ; y++ ) {
と書かれていたが，このループの順序を
for ( int y = 0; y < height ; y++ ) { for ( int x = 0; x < width ; x++ ) {
と書き換えて入れ替えた．
それぞれのプログラムは，図
4.1
のような順序で二重のループ文によって繰り返し計算し
ている．この繰り返しの入れ替えを行うだけで，計算中に発生する参照はメモリ上の連続し
た領域をアクセスするようになり，プログラム実行速度の向上が見込める．
ソースコード
4.1
プログラム中の
if
文
(mri sr.cpp
，
77
−
97
行目
)
4.4
条件分岐の削除
82 if ( lx < 0 || lx >= width || ly < 0 || ly >= width ) pl += 0;
83 else pl += result [lx + ( width *ly )];
90 if ( 0 < x && x < width -1 && 0 < y && y < height -1 ) { …
95 // update by deDp
96 if ( dEdp > 0 ) next[x + ( width *y)] = result [x + ( width *y)] - ( dEdp +2) / 4;
97 else next [x + ( width *y)] = result [x + ( width *y)] - (dEdp -2) / 4;
4.4
条件分岐の削除
既存のプログラムで出現する条件分岐はソースコード
4.1
の通りである．この中で，条件
式に
x
と
y
が使われていない最後の
if
文に関しては書き換えて省くことはできない．これ
は，この
if
文の判定項目が
1
ループ中に変化する変数
dEdp
であるためである．しかし残
りの
if
文の判定項目は
1
ループ中に変化することがない変数
x
と
y
であるため，
x
と
y
の取
りうる値によっては
if
文が必要なくなる．以下から，各
if
文の条件式に注目し，
if
文が必要
ない特定の
x
と
y
の範囲を求めていく．ただし，
x
と
y
は同じ範囲を示すことが明らかなの
で，
x
についてだけ範囲を求めることにする．
[1]
元のプログラム
(
非連続的な領域
)
[2]
最適化を施したプログラム
(
連続的な領域
)
図
4.1
各プログラムの計算順序
4.4
条件分岐の削除
4.4.1 new x
及び
new y
に関する条件分岐
まずは最初の，変数
new x
及び
new y
の値で判別される
4
つの
if
文に注目した．これら
の
if
文は
new x
及び
new y
が
0
未満または
384
以上の場合に
contune
し，
0
以上
384
未満
の場合に続けて処理を行うという条件分岐である．それらの変数は
x
及び
y
，
dx
，
dy
の値に
よって値が変化する．
dx
と
dy
については直前で定義されている
(
ソースコード
4.2)
．また，
new x
及び
new y
については
if
文の直前で定義されている
(
ソースコード
4.1)
．これらの式
から，この
if
文についての
x
と
y
の取りうる範囲を求める．
ソースコード
4.2 dx
と
dy
の定義
(mri sr.cpp
，
76
行目
)
1 int dx = param [l][4] , dy = param [l ][5];
width = height = 768
，
scale = 2
，
palam[l][4] = {0, 1}
，
palam[l][5] = {0, 1}
，
dx = palam[l][4] = {0, 1}
，
dy = palam[l][5] = {0, 1}
x
と
y
の範囲を求める上で必要になる式をまとめた．上記の式を利用して
new x
について
展開する．
new x =
x + dx
scale
=
x + dx
2
(4.1)
new x < 0
(4.2)
384 ! new x
(4.3)
式
4.1
と
if
文の判定式
(
式
4.2
，
4.3)
より，
x
の取りうる範囲を求める．式
4.1
と式
4.2
から
x + dx
2
<
0
x + dx < 0
が得られる．ここで，
x
のとりうる値を多く考える必要がある．つまり，
dx
が最小となる値
を考える必要があるため，ここでは
0
として考える．よって式
4.4
が得られる．
4.4
条件分岐の削除
次に，式
4.1
と式
4.3
から
384 !
x + dx
2
768 ! x + dx
となるが，先程と同様に
x
の取りうる値を多く考える必用があるため，ここでは
dx
が最大
となる値として考え計算する．
dx
の最大の値は
1
であるため，
768 ! x + 1
width − 1 ! x
(4.5)
となる．
以上から，この
if
文については以下のようにまとめられる．
• x < 0
または
width − 1 ! x
のとき，処理を行わず
continue
• 0 ! x < width − 1
のとき，処理を継続
4.4.2 lx
及び
ly
に関する条件分岐
続いて，
lx
と
ly
の値が条件式に関係している
if
文について，
x
の取りうる範囲を求める．
以下に示す式が，この
if
文の
x
に関係する条件式である．
lx < 0
(4.6)
width ! lx
(4.7)
まずは式
4.6
の場合の
x
の取りうる範囲について考える．この場合，
lx
が最小である場合
を考える必要があるので，
lx
は
lx = (new x × scale) − dx
(4.8)
となる．よって式
4.6
と式
4.8
から式
4.9
が得らる．
(new x × scale) − dx < 0
x < 0
(4.9)
4.4
条件分岐の削除
次に，式
4.7
の場合の
x
の取りうる範囲について考える．この場合，
lx
が最大である場合
を考えるため
lx
は
lx = ((new x + 1) × scale) − dx − 1
(4.10)
となる．よって式
4.7
と式
4.10
から式
4.11
が得られる．
width ! ((new x + 1) × scale) − dx − 1
width − 3 ! x
(4.11)
以上から，この
if
文は以下の様にまとめられる．
• x < 0
または
width − 3 ! x
のとき，
pl += 0
の処理
• 0 ! x < widtn − 3
のとき，
pl += result[lx + (width*ly)]
の処理
4.4.3 x
と
y
に関する条件分岐
最後に，
x
と
y
が条件式に含まれる
if
文について，
x
の取りうる範囲を求める．この
if
文
では条件式がそのまま
x
の取りうる範囲なので，その範囲は
0 < x
(4.12)
x < width − 1
(4.13)
となる．
よって，この
if
文は以下のようにまとめられる．
• x ! 0
または
width − 1 ! x
のとき，処理なし
• 0 < x < width − 1
のとき，処理あり
4.4.4
条件分岐が必要のない
x
と
y
の範囲
最後に，第
4.4.1
項から第
4.4.3
項で考えた全ての条件分岐を削除できる
x
の取りうる範
4.5
ループアンローリング
! g Ɔõ´üdj Ɔõ´ü ! ! üú†ûü´ üú†ûü´dj図
4.2
条件分岐が不要な範囲
! g Ɔõ´üdk Ɔõ´ü ! ! üú†ûü´ üú†ûü´dk図
4.3
条件分岐を削除した範囲
x
と同様に計算して求められる．
0 < x < width − 3
(4.14)
0 < y < height − 3
(4.15)
図
4.2
は計算して得られた条件分岐が不要な範囲を示している．図から，濃い影で示す超解
像画像の縁付近の画素に対しては条件分岐が必要であり，それ以外の影のない部分の画素に
ついては必要ないことがわかる．そのため，条件分岐を削除できる場所は超解像画像の縁付
近の画素以外である．本研究では図
4.3
の影のない範囲について，不要な条件分岐を削除し
高速化を図った．
4.5
ループアンローリング
プログラム内の一部の
for
文
(
ソースコード
4.3)
をループアンローリングを行った．
for
文
を手書きで展開し書き換えたプログラムの一部をソースコード
4.4
に示す．
ソースコード
4.3
元のソースコード
(mri sr.cpp
，
75-89
行目
)
75 for ( int l = 0; l < file_n ; l++ ) {
int dx = param [l][4] , dy = param [l ][5];
…
dEdp += pl - ql; }
4.6 C++
プログラムのコンパイル
ソースコード
4.4
書き換えたソースコード
(mri sr opt.cpp
，
167
−
206
行目
)
167 int l = 0;
int dx = param [l][4] , dy = param [l ][5]; int new_x = (x+dx) / scale ;
int new_y = (y+dy) / scale ; int pl = 0;
for ( int lx = ( new_x * scale ) - dx; lx < (( new_x +1)* scale ) - dx; lx ++ ) { for ( int ly = ( new_y * scale ) - dy; ly < (( new_y +1)* scale ) - dy; ly ++ ) {
pl += result [lx + ( width *ly )]; }
}
pl /= scale * scale ;
int ql = img[l][ new_x + (384 * new_y )]; dEdp += pl - ql;
l = 1; dx = param [l ][4]; dy = param [l ][5]; … …
ql = img [l][ new_x + (384 * new_y )]; dEdp += pl - ql; // l = 3
条件分岐の削除をした後にこの書き換えを行ったため，ループアンローリングを行う対象
の
for
文
(
ソースコード
4.3
と同様の
for
文
)
は複数箇所に現れる．本研究では，実行速度に
大きな変化が現れることが予想される計算量が最も多い，縁画素以外の画素の計算部分の
for
文のみ書き換えた．
4.6 C++
プログラムのコンパイル
トランスレートした
C++
プログラムでは
C++11
で追加されたクラスを用いている部分が
ある．そのクラスは
std::initializer list
であり，初期化リストを扱うクラスである．プログラ
ム中ではパラメータの値を
vector
で宣言する際，
vector < vector <int > > P {{0 ,0 ,0 ,0 ,1 ,1} , {0 ,0 ,0 ,0 ,0 ,1} , {0 ,0 ,0 ,0 ,1 ,0} , {0 ,0 ,0 ,0 ,0 ,0}};
というように初期化リストを用いている．そのため，プログラムのコンパイル時に
C++11
の機能を有効化させるオプション
-std=c++11
を付ける必要がある．
また，他のオプションとして最適化オプションをつけた．
C++
コンパイラの
g++
には最適
化オプションが準備されている．表
4.3
に最適化オプションの一部をまとめた．本研究では
4.6 C++
プログラムのコンパイル
表
4.3 g++
に準備されている最適化オプション
オプション
概要
-O0
最適化を行わない
-O, -O1
コードサイズと実行時間を小さくする
-O2
-O(-O1)
より強力な最適化
-O3
-O2
より強力な最適化
-Os
サイズの最適化
のコンパイル時には以下の様にしてコンパイルを行っている．
✓
✏
$ g++ -std=c++11 -O3 mri_sr_opt.cpp
✒
✑
第
5
章
評価実験
5.1
実験環境
本研究の評価実験では異なる
2
台の計算機を使用した．それぞれの計算機については表
5.1
にまとめた．
OS
は
2
台とも同じ
ubuntu14.04 LTS
である．
5.2
元のプログラムとの比較
5.2.1
実験内容
本実験では，題材にした既存の
Java
プログラム，
C++
に書き換えたプログラム，最適化
を施したプログラムについて行った．画像処理の本処理部
(
関数
img iteration)
における処理
時間を計測し，比較した．本実験では同じ処理を
20
回繰り返し，
1
回あたりの平均処理時
間を求めたものを評価対象とする．各プログラムについては表
5.2
にまとめた．表中の行数
は超解像処理部の関数
img_iteration
の行数であり，プログラム全体の行数ではない．
表
5.1
プログラムを実行する計算機
計算機
CPU (
クロック周波数
)
メモリ
matsu-lab98 IntelXeon Processor E5-2620 v3 (2.4GHz)
32GB
matsu-lab99 IntelCore i7-4770 (3.4GHz)
32GB
5.3
各最適化手法による比較
表
5.2
実行するプログラムについて
ファイル名
行数
備考
mri sr.java
47
既存の
Java
コード
mri sr.cpp
39 C++
に書き換えただけの
C++
コード
mri sr opt.cpp
169
各最適化手法を施した
C++
コード
表
5.3
プログラム別の処理時間の計測結果
計算機
matsu-lab98
matsu-lab99
プログラム
処理時間
(ms)
向上率
処理時間
(ms)
向上率
mri sr.java
16155
-
12334
-mri sr.cpp
11041.8
32 %
8880.09
28 %
mri sr opt.cpp
4245.18
74 %
3491.66
72 %
5.2.2
実験結果と考察
実行環境別の各プログラムの処理時間を表
5.3
にまとめた．求められた数値からグラフ化
したものが図
5.1
である． 計測結果から高速化が確認でき，本研究で行った最適化が処理時
間の短縮につながったことがわかる．最適化を施した
C++
プログラムの処理時間は
Java
プ
ログラムの
70%
以上の向上が確認できた．以上より，最適化手法をプログラムに施すことに
よって高速化の効果を得られたことがわかった．
5.3
各最適化手法による比較
5.3.1
実験内容
続いて，各最適化手法別に処理時間の比較を行った．本研究で行った各最適化手法を一つ
づつ施したプログラムを用意し，それぞれの処理時間を計測した．比較対象は
Java
を
C++
に書き換えたのみのプログラム
(mri sr.cpp)
である．表
5.4
は，用意したプログラムについ
5.3
各最適化手法による比較
図
5.1
プログラム別の処理時間の差を示したグラフ
表
5.4
実行するプログラムについて
ファイル名
行数
施してある最適化手法
mri sr.cpp
39
最適化なし
mri sr opt cache.cpp
39
キャッシュメモリの活用
mri sr opt const.cpp
40
変数の定数化
(
定数伝播
)
mri sr opt jump.cpp
144
条件分岐の削除
mri se opt jump unroll.cpp
168
条件分岐の削除，ループアンローリング
てまとめたものである．実装の都合上ループアンローリングについては，条件分岐の削除と
同時に施す必要があるため，ループアンローリング単独での比較はできなかった．そのため
ループアンローリングを施したプログラムの比較対象は，条件分岐の削除を行ったプログラ
ム
(mri sr opt jump.cpp)
である．
5.3
各最適化手法による比較
表
5.5
最適化手法別のプログラムの処理時間の計測結果
計算機
matsu-lab98
matsu-lab99
プログラム
処理時間
(ms)
向上率
処理時間
(ms)
向上率
mri sr.cpp
11041.8
-
8880.90
-mri sr opt cache.cpp
10391.0
6 %
8601.92
3 %
mri sr opt const.cpp
7372.64
33 %
5760.69
35 %
mri sr opt jump.cpp
9895.32
10 %
7994.95
10 %
mri sr opt jump unroll.cpp
9689.97
12 %
7772.12
12 %
図
5.2
各最適化手法別のプログラムの処理時間の差を示したグラフ
5.3.2
実験結果と考察
最適化手法別のプログラムの処理時間を表
5.5
にまとめた．また，求められた数値からグ
ラフ化したものが図
5.2
である． 計測結果から，本研究で題材にしたプログラムでは，定数
伝播を狙った変数の定数化が最も高速化の効果を得られ，処理時間は約
35%
向上した．書き
換える行数も少なく，向上率が高いことから高速化においては効率的な最適化であると考え
5.3
各最適化手法による比較
る．次に効果が得られた最適化手法は条件分岐の削除であり，処理時間は約
10%
という結
果が得られた．この最適化は，コードの書き換えるコストを考慮すると，あまり効率的では
ないと考えられる．また，キャッシュメモリの活用に関しては高速化の効果が小さく，向上
率は約
6%
であった．向上率は低いが，プログラムの書き換えとしては最も単純であるため，
手頃な最適化であることがわかる．しかし，生成する超解像画像の解像度をより大きい物に
した場合，扱う配列も大きくなる．このような場合においては，処理時間の向上率にも変化
があるのではないかと考える．
また，ループアンローリングを施したプログラムに関しても向上率は低く，約
3%
という
結果が得られた．本研究で取り扱ったプログラムではループアンローリングによる高速化
は，極端に現れなかった．これは
SIMD
最適化が適用できない限り，大きな影響は出なかっ
たと考える．
第
6
章
まとめと今後の課題
本研究では計算量が多いとされる画像処理プログラムの高速化を行った．題材とした超解
像処理
Java
プログラムを元に，最適化を行い処理時間を計測し比較した．これにより，い
くつかの最適化による高速化の効果を検証した．
まずは
Java
プログラムを
C++
プログラムにトランスレートを行った．その作成した
C++
プログラムを元に各最適化手法を施して高速化
C++
プログラムを作成した．
本研究では
4
種の最適化手法を用いて，プログラムを書き換えた．まずは定数伝播を利用
した最適化を施した．プログラム内で値が変化しない変数を定数化することでコンパイラ
に定数伝播をするように書き換え，高速化を図った．次にキャッシュメモリの活用を行った．
書き換える内容は単純で，二重
for
文の外側と内側のループを入れ替えるだけである．また，
条件分岐の削除を行うことで更に高速化を狙った．
if
文の条件式をすべて満たすような変数
の組み合わせを計算し，算出した．プログラムが複雑にはなるが，部分的に条件分岐が必要
のないプログラムにした．同時に，ループアンローリングによる最適化を施した．プログラ
ム中の繰り返し文の一部を展開することで高速化を図った．
評価実験では，
Java
プログラム，トランスレートしただけの
C++
プログラム，最適化を
行った高速化
C++
プログラムの処理時間を計測し，比較した．実験結果から，高速化
C++
プログラムの処理時間は
70%
ほど向上が確認できた．また，各最適化別に高速化の効果を比
較するため，各最適化をそれぞれ施したプログラムを用意し，処理時間を計測し比較した．
実験結果から，確認できた定数伝播をうながす変数を定数化したプログラムが最も高速化の
効果を得られた．その処理時間はトランスレートした
C++
プログラムと比較して
35%
の向
上していることを確認した．以上のことから，本研究で施した最適化によってプログラムの
高速化を実現した．また，本研究で行った最適化の中では値が変化しない変数を定数化し，
定数伝播をさせるものが最も効果的であることがわかった．
今後は，本研究で行っていない最適化できる部分を検討し，更に高速化を図ることが考え
られる．また，並列化を実装することによる高速化も可能であると考える．これらを検討し
て，更に高速化が可能であるかを考える必要がある．
謝辞
本研究を行う際に，ご指導いただいた松崎公紀准教授に深く感謝いたします．また，副査
をしていただいた鵜川始陽准教授，高田喜朗准教授にも御礼申し上げます．最後に，研究中
に助言をいただいた松崎研究室の修士の先輩方や学友に感謝いたします．
参考文献
[1]
延原肇
,
ムハマドハリス
,
渡邉拓也
,
鋤先星太，
“
超解像処理アルゴリズム入門
”
，
「
Inter-face
」
2015
年
6
月号，
pp.28–38
，
CQ
出版，
2015
．
[2]
宮崎玲奈，
“fMRI
画像を対象とした超解像技術に関する研究
”
，高知工科大学 情報学群
卒業論文，
2015
．
[3]
今城哲二
,
布広永示
,
岩澤京子
,
千葉雄司，
“
コンパイラとバーチャルマシン
”
，オーム社，
2004
．
[4]
片山善夫，
“C
プログラム高速化研究班 コードを高速化する
20
の実験と達人の技
”
，ユ
ニバーサル・シェル・プログラミング研究所，
2012
．
[5]
坂井修一，
“
コンピュータアーキテクチャ
”
，コロナ社，
2004
．
[6]
中田育男
,
渡邊坦
,
佐々政孝
,
瀧本宗宏，
“
コンパイラの基盤技術と実践 コンパイラ・イ
ンフラストラクチャ
COINS
を用いて
”
，朝倉書店，
2008
．
付録
A
既存の
Java
プログラム
(mri sr.java)
1 import javax . imageio .*;
2 import java . lang .*;
3 import java .io .*;
4 import java . awt . image . BufferedImage ;
5 import java . util . ArrayList ;
6 import java . util . StringTokenizer ;
7 /* 複 数 枚 のf m r i画 像(1つ のc s v f i l eに ま と ま っ て い る)を 入 力 と し て 超 解 像 を 行 う */
8 /* csv 行intput_nx列64 */
9 /* java mri_sr */
10 public class mri_sr {
11 static int alpha = 4;
12
13 public static void main ( String args []) throws Exception {
14 if ( args. length != 1) {
15 System . out . println (" error : shape ?0:1 ");
16 return ; 17 } 18 19 int [][] P = {{0 ,0 ,0 ,0 ,1 ,1} ,{0 ,0 ,0 ,0 ,0 ,1} , {0 ,0 ,0 ,0 ,1 ,0} , {0 ,0 ,0 ,0 ,0 ,0}}; // ア フ ィ ン 変 換 の パ ラ メ ー タ 指 定 20 // (1 ,1)移 動 (0 ,1)移 動 (1 ,0)移 動 (0 ,0)移 動
21 int input_n = P. length ; // L Rの フ ァ イ ル 数
22 int shape = Integer . valueOf ( args [0]);
23
24 // c s vフ ァ イ ル(lr)の 読 み 込 み
25 int [][] img = new int[ input_n ][384 * 384];
26 int boxcel_n = 384 * 384; // 200704
27
28 BufferedReader br = new BufferedReader ( new FileReader ("./ inputcsv /lr" + input_n + " _shape " + shape + ". csv " ));
29 for (int i = 0; i < input_n ; i++) { // L R画 像 の 枚 数 分
30 //読み込んだファイルを1行 ず つ 処 理 す る
31 StringTokenizer token = new StringTokenizer (br. readLine (), ",");
34 img [i][j] = data ; // pxcel value 35 } 36 } 37 br. close (); 38 long total = 0; 39 int [] sr = {};
40 for ( int loop = 0; loop < 20; loop ++ ) {
41 long start = System . currentTimeMillis ();
42 sr = img_iteration (img , P);
43 long end = System . currentTimeMillis ();
44 System . out . println (" count " + ( loop +1) + ": Elapsed time = " + ( end - start ) + "ms .");
45 total = total + ( end - start );
46 }
47 System . out . println (" Average Elapsed time = " + total / 20 + "ms.");
48 // S R i m a g eフ ァ イ ル の 書 き 出 し
49 analys ana = new analys (sr );
50 ana . to_img ("./ result / il_shape " + shape + "_a" + alpha );
51 ana . to_csv ("./ result / il_shape " + shape + "_a" + alpha );
52 }
53
54 public static int [] img_iteration (int [][] img , int [][] param ) throws Exception {
55 int file_n = img. length ; // 画 像 数 = 移 動 パ ラ メ ー タ 数
56 int scale = (int)Math .sqrt ( file_n );
57 int width = (int)Math .sqrt (img [0]. length ) * scale ; // 超 解 像 画 像 の 幅 と 高 さ
58 int height = width ;
59 int [] result = new int[ width * height ]; // 反 復 法 用 の デ ー タ
60 initialize ( result ); // 中 央 値:2047に 初 期 化0
61 for (int test = 0; test < 200; test ++) { // と り あ え ず1 0 0回 ま わ す
62 int [] next = new int[ width * height ];
63 clear ( next );
64 for (int x = 0; x < width ; x++) {
65 for (int y = 0; y < height ; y++) {
66 int dEdp = 0; // E = Sum (p_i - q_i )ˆ2 と し た と き の dE/dp
67 for (int l = 0; l < file_n ; l++) {
68 int dx = param [l][4] , dy = param [l ][5]; // ↓ 平 行 移 動
69 int newx = (x + dx) / scale ; if (newx < 0) continue ; if ( newx >= 384) continue ;
70 int newy = (y + dy) / scale ; if (newy < 0) continue ; if ( newy >= 384) continue ;
71 int pl = 0; // ↓ ダ ウ ン サ ン プ リ ン グ
72 for (int lx = newx * scale - dx; lx < (newx + 1) * scale - dx; lx ++) {
73 for (int ly = newy * scale - dy; ly < (newy + 1) * scale - dy; ly ++) {
74 if (lx < 0 || lx >= width || ly < 0 || ly >= width ) {
75 pl += 0;
76 } else {
80 }
81 pl /= scale * scale ;
82 int ql = img[l][ newx + 384 * newy ];
83 dEdp += pl - ql; // 差 分 を 合 計
84 }
85 if (x > 0 && y > 0 && x < width - 1 && y < height - 1) {
86 int pn = ( result [x -1 + width * y] + result [x+1 + width * y]
87 + result [x + width * (y -1)] + result [x + width * (y +1)]) /4;
88 dEdp += alpha * ( result [x + width * y] - pn );
89 }
90 if (dEdp > 0) { // 最 低1は 変 化 さ せ る
91 next [x + width * y] = result [x + width * y] - ( dEdp + 2) / 4; //
誤 差 分 画 素 値 修 正(四 捨 五 入)
92 } else {
93 next [x + width * y] = result [x + width * y] - ( dEdp - 2) / 4;
94 } 95 } 96 } // 200回 97 result = next ; 98 } 99 return result ; 100 } 101
102 public static int [] md (int [] img , int [] param , int n) {
103 int width = (int)Math .sqrt (img. length );
104 int height = width ;
105 int [] result_m = new int[ width * height ];
106 int dx = param [4]; // xの 移 動 分
107 int dy = param [5]; // y "
108 // m o t i o n付 加
109 for (int x = 0; x < width ; x++) {
110 for (int y = 0; y < height ; y++) {
111 if (x - dx >= 0 && y - dy >= 0 && x - dx < width && y - dy < height ) {
112 int px = img [(x - dx) + (y - dy) * width ];
113 result_m [x + y * width ] = px; // 移 動
114 }
115 }
116 }
117 // Downsamp
118 int [] result_md = new int [384 * 384];
119 int pass = (int)Math. sqrt(n);
120
121 for (int y = 0; y < height ; y = y + pass ) {
122 for (int x = 0; x < width ; x = x + pass) {
123 int sum = 0;
124 for (int i = 0; i < pass ; i++) {
125 for (int j = 0; j < pass; j++) {
126 sum += result_m [(x + i) + width * (y + j )];
127 }
129 sum /= n; // 平 均 値
130 result_md [(x / pass ) + 384 * (y / pass )] = sum ;
131 }
132 }
133 return result_md ;
134 }
135
136 static int diff (int [] a, int [] b) {
137 int ret = 0;
138 // System .out. println ("a. length = " + a. length );
139 int width = (int)Math .sqrt (a. length );
140 for (int x = 0; x < width ; x++) {
141 for (int y = 0; y < width ; y++) {
142 int av = a[x + y * width ];
143 int bv = b[x + y * width ];
144 ret += (av - bv) * (av - bv );
145 }
146 }
147 return ret;
148 }
149
150 public static void initialize (int [] data ) {l
151 for (int i = 0; i < data . length ; i++) {
152 data [i] = 2047;
153 }
154 return ;
155 }
156 /* 1次 元 配 列 初 気 化 */
157 public static void clear (int [] data) {
158 for (int i = 0; i < data . length ; i++) {
159 data [i] = 0;
160 }
161 return ;
162 }
付録
B
高速化
C++
プログラム
(mri sr opt.cpp)
1 /* 2 * 20151218 Java -> C p pに 書 き 換 え た だ け ． 3 * xとyの ル ー プ を 逆 に 4 * i m g _ i t e r a t i o n内 の file_n , scale を 定 数 化 5 * 画 像 の 縁 の 画 素(４ 画 素)分( i f文 あ り)と 中 間 の 画 素 分( i f文 な し)に 分 け た 6 * 最 も 処 理 が 多 い 中 間 画 素 部 分 に お い て ，f i l e _ nの 数 for ( l < file_n )を す べ て 展 開 7 * p a r a mを 引 数 で 渡 す の で は な く ， ロ ー カ ル 変 数 に し て ，lを 定 数 化 8 */ 9 # include <stdio .h> 10 # include <iostream > 11 # include <math .h> 12 # include <stdlib .h> 13 # include <fstream > 14 # include <string > 15 # include <sstream > 16 # include <vector > 17 # include <time .h> 18 19 # define SIZE 384 20 # define SHAPE "0" 21 # define LOOP 20 22 # define FILE_N 4 23 # define SCALE 2 24 int alpha = 4; 25 using namespace std ;26 string itos ( int i);
27 vector <int > img_iteration ( vector < vector <int > > img );
28
29 int main () {
30 // vector <vector <int > > P {{0 ,0 ,0 ,0 ,1 ,1} , {0 ,0 ,0 ,0 ,0 ,1} , {0 ,0 ,0 ,0 ,1 ,0} , {0 ,0 ,0 ,0 ,0 ,0}};
31 int input_n = 4;
32 string shape = SHAPE ;
35 vector < vector <int > > img ( input_n , vector <int >( SIZE * SIZE ));
36 int boxcel_n = SIZE * SIZE;
37
38 /* Input CSV */
39 string num = itos ( input_n );
40 ifstream ifs ("./ inputcsv /lr" + num + " _shape " + shape + ". csv ");
41 if( !ifs ) {
42 cout << " Error : Input data file not found " << endl ; return 0;
43 }
44 string str ;
45 int x = 0, y = 0;
46 while ( getline (ifs , str) ) {
47 string tmp ;
48 istringstream stream ( str );
49 y = 0;
50 while ( getline (stream , tmp , ’,’) ) {
51 img [x][y] = ( stoi ( tmp ));
52 y ++;
53 }
54 x ++;
55 } /* Input CSV end */
56
57 int lop = LOOP;
58 clock_t total = 0;
59 while (lop --) {
60 clock_t start = clock ();
61 vector <int > sr = img_iteration ( img );
62 clock_t finish = clock ();
63 // for ( int x = 0; x < sr. size (); x++ ) {
64 // if (x != 0) cout << " ,";
65 // cout << sr[x];
66 // }
67 cout << " case " << LOOP - lop << ":: duration = " << (( double )( finish - start ) / CLOCKS_PER_SEC ) * 1000 << "ms .\n";
68 total = total + ( finish - start );
69 }
70 cout << " average duration = " << (( double )( total / LOOP ) / CLOCKS_PER_SEC ) * 1000 << "ms .\n";
71 return 0;
72 }
73
74 vector <int > img_iteration ( vector < vector <int > > img ) {
75 const int file_n = FILE_N ;
76 const int scale = SCALE ;
77 const int width = (int)sqrt (img [0]. size ()) * scale ;
78 const int height = width ;
79 vector <int > result ( width * height , 2047);
80 vector <int > next ( width * height , 0);
83 for ( int y = 0; y < 4; y++ ) {
84 for ( int x = 0; x < width ; x++ ) {
85 int dEdp = 0;
86 for ( int l = 0; l < file_n ; l++ ) {
87 int dx = param [l][4] , dy = param [l ][5];
88 int new_x = (x+dx) / scale ; if ( new_x < 0 ) continue ; if ( 384 <= new_x ) continue ;
89 int new_y = (y+dy) / scale ; if ( new_y < 0 ) continue ; if ( 384 <= new_y ) continue ;
90 int pl = 0;
91 for ( int lx = ( new_x * scale ) - dx; lx < (( new_x +1)* scale ) - dx; lx ++ ) {
92 for ( int ly = ( new_y * scale ) - dy; ly < (( new_y +1)* scale ) - dy; ly ++ ) {
93 // asm (" nop ");
94 if ( lx < 0 || lx >= width || ly < 0 || ly >= width ) pl += 0;
95 else pl += result [lx + ( width *ly )];
96 }
97 }
98 pl /= scale * scale ;
99 int ql = img[l][ new_x + (384 * new_y )];
100 dEdp += pl - ql;
101 }
102 if ( 0 < x && x < width -1 && 0 < y && y < height -1 ) {
103 int pn = ( result [(x -1) + ( width *y)] + result [(x+1) + ( width *y)]
104 + result [x + ( width *(y -1))] + result [x + ( width *(y +1))] ) / 4;
105 dEdp += alpha * ( result [x+( width *y)] - pn );
106 }
107 if ( dEdp > 0 ) next [x + ( width *y)] = result [x + ( width *y)] - (dEdp +2) / 4;
108 else next [x + ( width *y)] = result [x + ( width *y)] - (dEdp -2) / 4;
109 }
110 }
111 for ( int y = height -4; y < height ; y++ ) {
112 for ( int x = 0; x < width ; x++ ) {
113 int dEdp = 0;
114 for ( int l = 0; l < file_n ; l++ ) {
115 int dx = param [l][4] , dy = param [l ][5];
116 int new_x = (x+dx) / scale ; if ( new_x < 0 ) continue ; if ( 384 <= new_x ) continue ;
117 int new_y = (y+dy) / scale ; if ( new_y < 0 ) continue ; if ( 384 <= new_y ) continue ;
118 int pl = 0;
119 for ( int lx = ( new_x * scale ) - dx; lx < (( new_x +1)* scale ) - dx; lx ++ ) {
120 for ( int ly = ( new_y * scale ) - dy; ly < (( new_y +1)* scale ) - dy; ly ++ ) {
121 if ( lx < 0 || lx >= width || ly < 0 || ly >= width ) pl += 0;
123 }
124 }
125 pl /= scale * scale ;
126 int ql = img[l][ new_x + (384 * new_y )];
127 dEdp += pl - ql;
128 }
129 if ( 0 < x && x < width -1 && 0 < y && y < height -1 ) {
130 int pn = ( result [(x -1) + ( width *y)] + result [(x+1) + ( width *y)]
131 + result [x + ( width *(y -1))] + result [x + ( width *(y +1))] ) / 4;
132 dEdp += alpha * ( result [x+( width *y)] - pn );
133 }
134 if ( dEdp > 0 ) next [x + ( width *y)] = result [x + ( width *y)] - (dEdp +2) / 4;
135 else next [x + ( width *y)] = result [x + ( width *y)] - (dEdp -2) / 4;
136 }
137 }
138 for ( int y = 4; y < height -4; y++ ) {
139 for ( int x = 0; x < 4; x++ ) {
140 int dEdp = 0;
141 for ( int l = 0; l < file_n ; l++ ) {
142 int dx = param [l][4] , dy = param [l ][5];
143 int new_x = (x+dx) / scale ; if ( new_x < 0 ) continue ; if ( 384 <= new_x ) continue ;
144 int new_y = (y+dy) / scale ; if ( new_y < 0 ) continue ; if ( 384 <= new_y ) continue ;
145 int pl = 0;
146 for ( int lx = ( new_x * scale ) - dx; lx < (( new_x +1)* scale ) - dx; lx ++ ) {
147 for ( int ly = ( new_y * scale ) - dy; ly < (( new_y +1)* scale ) - dy; ly ++ ) {
148 if ( lx < 0 || lx >= width || ly < 0 || ly >= width ) pl += 0;
149 else pl += result [lx + ( width *ly )];
150 }
151 }
152 pl /= scale * scale ;
153 int ql = img[l][ new_x + (384 * new_y )];
154 dEdp += pl - ql;
155 }
156 if ( 0 < x && x < width -1 && 0 < y && y < height -1 ) {
157 int pn = ( result [(x -1) + ( width *y)] + result [(x+1) + ( width *y)]
158 + result [x + ( width *(y -1))] + result [x + ( width *(y +1))] ) / 4;
159 dEdp += alpha * ( result [x+( width *y)] - pn );
160 }
161 if ( dEdp > 0 ) next [x + ( width *y)] = result [x + ( width *y)] - (dEdp +2) / 4;
164
165 for ( int x = 4; x < width -4; x++ ) {
166 int dEdp = 0;
167 int l = 0;
168 int dx = param [l][4] , dy = param [l ][5];
169 int new_x = (x+dx) / scale ;
170 int new_y = (y+dy) / scale ;
171 int pl = 0;
172 for ( int lx = ( new_x * scale ) - dx; lx < (( new_x +1)* scale ) - dx; lx ++ ) {
173 for ( int ly = ( new_y * scale ) - dy; ly < (( new_y +1)* scale ) - dy; ly ++ ) {
174 pl += result [lx + ( width *ly )];
175 }
176 }
177 pl /= scale * scale ;
178 int ql = img[l][ new_x + (384 * new_y )];
179 dEdp += pl - ql;
180 l = 1; dx = param [l ][4]; dy = param [l ][5]; new_x = (x+dx) / scale ; new_y = (y+dy) / scale ; pl = 0;
181 for ( int lx = ( new_x * scale ) - dx; lx < (( new_x +1)* scale ) - dx; lx ++ ) {
182 for ( int ly = ( new_y * scale ) - dy; ly < (( new_y +1)* scale ) - dy; ly ++ ) {
183 pl += result [lx + ( width *ly )];
184 }
185 }
186 pl /= scale * scale ;
187 ql = img [l][ new_x + (384 * new_y )];
188 dEdp += pl - ql;
189 l = 2; dx = param [l ][4]; dy = param [l ][5]; new_x = (x+dx) / scale ; new_y = (y+dy) / scale ; pl = 0;
190 for ( int lx = ( new_x * scale ) - dx; lx < (( new_x +1)* scale ) - dx; lx ++ ) {
191 for ( int ly = ( new_y * scale ) - dy; ly < (( new_y +1)* scale ) - dy; ly ++ ) {
192 pl += result [lx + ( width *ly )];
193 }
194 }
195 pl /= scale * scale ;
196 ql = img [l][ new_x + (384 * new_y )];
197 dEdp += pl - ql;
198 l = 3; dx = param [l ][4]; dy = param [l ][5]; new_x = (x+dx) / scale ; new_y = (y+dy) / scale ; pl = 0;
199 for ( int lx = ( new_x * scale ) - dx; lx < (( new_x +1)* scale ) - dx; lx ++ ) {
200 for ( int ly = ( new_y * scale ) - dy; ly < (( new_y +1)* scale ) - dy; ly ++ ) {
201 pl += result [lx + ( width *ly )];
202 }
203 }
204 pl /= scale * scale ;
205 ql = img [l][ new_x + (384 * new_y )];
206 dEdp += pl - ql;
