2011年5月15日日曜日

2直線の交点の計算

y=ax+b で定義される2直線をランダムに生成し、その2直線が交わる交点 (x,y) を求める例題です。

反復により一番誤差の少ない点を求める方法と、セカント法による方法を比較しています。
(逆行列を計算してと言う方法もありますが、ここでは単純な計算例ということで。。)

セカント法は、非線型方程式(f(x)=0)の解を求める手法ですが、f(x) の定義部に当たる nw_func() の定義を修正すれば直線方程式の解を求めることもできます。

例題では、両手法の反復回数、探索範囲を同じにしていますが、セカント法の方が誤差を小さく計算できていることが分かります。

実行結果

基礎統計量の計算

統計に必要となる基礎的な統計量の計算と標準化されたデータ列に対する統計量の計算例です。

実行結果

2項分布による確率計算

2項分布を参考に例題を用意しています。

2項分布
次の1~3を満たす試行をベルヌーイ試行という。

1. 各試行において、その事象が発生するか否かのみを問題にする。
2. 各試行は統計的に独立する。
3. 対象とする事象が発生する確率は、各試行を通じて一定である。

1 回の試行において、ある事象 X が発生する確率を p とする。
n 回のベルヌーイ試行列において、ちょうど i 回事象 X が発生する確率は、

で表され、このとき X の確率分布を2項分布といい、X 〜 B(n,p) と表す。

例題
1. 2項分布 B(8,0.4), B(8,0.2) の確率分布を求めよ
2. プラモデルが25セットあり、そのうち2セットは部品がかけているとする。客が任意に3セット選ぶとき、いずれも完全なセットである確率を求めよ
3. ある試薬をネズミに注射すると一定期間のうちに 40% が死亡するという。8匹のネズミにその試薬を注射した場合、
(a) 8匹全部が一定期間以上生存する確率を求めよ
(b) 95% の確率で,8匹のうち何匹以上が生存するということができるか
実行結果

相関行列(係数)の計算

2つの変量から相関係数を計算し、相関行列を計算する例題です。

2つの確率変数(変量)x と y の関連性を調べるために、n 組のデータ (x1,y1),...,(xn,yn) から次式によって(標本)相関係数 rxy を計算します。



例題では、変量 x1 と x3 の相関係数が 0.616 となり、x2 との相関より高いことが分かります。

実行結果

2011年5月14日土曜日

構造体データのバイナリファイル出力、読み取り

構造体データをバイナリファイルとして出力し、バイナリファイルから構造体データを読み取る例題です。

固定サイズの構造体の場合
読み取りの際には、構造体データのサイズ毎に読み取り、データが無くなるまで処理を行っています。
実行結果
サイズが変動する構造体の場合

上記の例題では、読み書きのサイズには sizeof(entry) にてサイズを指定していますが、文字列の長さが可変の場合には注意が必要で、文字列の長さに応じて構造体のサイズを毎回計算する必要があります。

以下の例題では、構造体のサイズが毎回異なる場合のファイル出力、読み込み、表示を行なっています。

ファイルフォーマット:
[格納する名前の数][リストのタイトル][格納数に応じた名前の長さ、文字列の繰り返し]
[格納する名前の数][リストのタイトル][格納数に応じた名前の長さ、文字列の繰り返し]
[格納する名前の数][リストのタイトル][格納数に応じた名前の長さ、文字列の繰り返し]
→ 格納する名前の数が読み取れなくなるまで構造体を読み取る
実行結果
既存のフォーマットファイルの場合(.bmp/.wav ファイル)

RFC 等により既定されている既存のファイルフォーマットは、構造体の形式を予め調べておくことにより入出力することができます(画像:ビットマップ形式PNM形式、音声:AAC/ADTS形式)。

以下では、Windows PCM フォーマット(拡張子 .wav)のファイルを読み込み、チャンネル数やサンプリング周波数等の情報の表示を行なっています(WAVE/PCM 音声の生成)。
実行結果

ファイルサイズの取得

ファイルサイズの取得例です。

stat 関数を利用し、ファイルやディレクトリの状態を取得していますが、処理系によっては利用できません(sys/stat.h が必要)。
ファイルポインタの位置から取得する方法との計算時間も比較していますが、結果としては stat 関数を利用する方が高速であることが言えます。

複数のファイルサイズのチェックが必要な場合や、ファイルサイズが非常に大きな場合などには stat 関数を用いる方が良さそうです。


実行結果

2分木による番号データのソート

入力したデータ(ファイル入力)から2分木を作成し、番号順にソートした結果を表示する例題です。


実行結果

btree.dat
2169 michael
1922 mike
2374 arthor
343 ken
1080 john
1080 carl
1278 david
2231 dianne
1175 alan
1584 andreas
1430 ben

ソートアルゴリズムの比較

ソートアルゴリズムを各種準備し(昇順、降順含む)、その計算時間を比較した例です。

単純に値の大小を比較し並べ替えるシンプルなソート、C言語でクイックソートを実装したもの、qsort 関数によりクイックソートを行うもの、の 3 種類準備しています。
100,200,400,800,1600,3200,6400,12800,25600,51200,102400 個(計11回の処理時間比較)の値を並べ替えた際の計算時間結果は以下のグラフになっています(Y軸はログスケール)。

(黒:シンプルソート、赤:クイックソート、青:qsort 関数)

グラフを見てみると、それぞれ指数関数的に計算時間が増えていることが分かりますが、それぞれの傾きが異なり、qsort よりも自前のクイックソートを用いた処理時間が短いことが分かります。
(コンパイルオプションには -ffast-math -O4 を指定)

CPUや処理系にも依存する部分ですが、102400 個程度の配列では 1 秒にも満たない時間なので、どちらでも良さそうですが、非常に大きな配列を扱う場合には参考になるかと思います。


実行結果

クイックソートによるソート

C言語の標準関数である qsort() 関数を用いずにクイックソートを行う例です。

入力データは乱数を用いています。
ソートされた結果は0から乱数の最大数までの昇順の数になります。

実行結果

2011年5月10日火曜日

qsort を用いたソート

C言語の標準関数である qsort() 関数を用いたソートの例です。

qsort の引数:
void qsort(void * base, size_t nmemb, size_t size, int (* compar)(const void *, const void *) );
void* base : データ元となるデータを指し示すポインタ用配列
size_t nmemb: 配列の数
size_t size: データの型のサイズ(sizeof で指定)
int (*compar)(const void*, const void*)): ソートで用いる大小比較用関数
qsort() の関数内では、データ元となるデータもしくは配列をそれぞれポインタで指し示し、ポインタが指し示すデータの順番を入れ替えていくことになります。

例題では、
文字列のデータを指し示すポインタ用配列:char* c_table[3];
数値列のデータを指し示すポインタ用配列: int* i_table[5];
のそれぞれのアドレスを
c_table[0] = a; → 配列aの先頭アドレス
i_table[i] = &idata[i]; → idata[]配列のi番目のデータのアドレス

というように指定し、qsort() 関数が指し示されたアドレスの順番を入れ替えていくことになります。

qsort() で指定する大小規定用の関数 compare() は、文字列、数値列用それぞれに準備しています。
compare() の関数は、2つのデータの大小を正負の値を使って表現するような関数を自分で規定するものですが、今回は、単に文字列の長さと、データの大きさで指定しています。

実行結果

2011年5月7日土曜日

2次方程式の解法

2次方程式の解を得る例題です。

引数により与えられたパラメータにより2次方程式を規定し、その方程式の解を計算します。

実行結果

組合せ数 (nCr) の計算

組合せ数(nCr)の計算例です(漸化式、Pascalの3角形による方法)。

n 個の数から r 個を選択する場合の組合せ数を計算します。
通常、組合わせ数 nCr は以下のような漸化式により計算することができます。

また、重複がある場合の組合せ数は以下のように計算することができます。

例では、選択する数に重複がある場合とない場合(選択する数が同じでも順番のみが異なっている場合)に分けて計算しています。

重複がある場合には値が非常に大きくなってしまうため、どの辺までオーバーフローなく計算できているか、若干疑問です。
漸化式による計算と、Pascal の方法による計算を比較していますが、Pascal の方法による計算の方が、加算により計算しているのでオーバーフローは少ないかと思われます。

全ての順列、組み合わせを配列として取得する例は、こちら

実行結果

相関係数の計算

相関係数の計算例です。



相関係数を計算する数値列が同じ場合、相関係数は 1 となり、
相関係数を計算する数値列の符号が逆の場合、相関係数は -1 となります。

相関行列の計算は こちら

実行結果

ビットシフトによる掛算

ビットシフトによる掛算の例題です。

計算をループさせ、通常の掛け算との計算時間を比較していますが、通常の計算の方が早くなっています。
また、単純なビットシフトにより2の冪乗の計算も行っていますが、pow 関数による計算の方が早くもなっています。

特殊な計算が必要な場合以外には、掛け算として用いるのは得策ではないようです。


実行結果

経過時刻の計算/出力(クラス化)

経過時刻の計算/出力例です。

開始時間・終了時間を設定し、計画時刻を出力します。

時刻は "20041125172211" のように、
2004:年,11:月,25:日,17:時,22:分,11:秒
を記した文字列で出力されます。

開始時間または終了時間が設定されていない場合、開始時間と終了時間の設定が逆の場合、"Illegal WorkTime" と表示されます。


実行結果

アセンブラによる足し算

アセンブラによる足し算の例題です。

計算をループさせ通常の足し算と計算時間を比較していますが、通常の計算の方が速いようです。

実行結果

約数の計算

整数から約数を計算する例題です。

実行結果

2011年5月5日木曜日

libavcodec による動画のサムネイル作成

動画ファイルから動画位置を0から1までの割合で指定し、サムネイルを作成し保存する例題です。

引数により動画位置や出力ファイル(デフォルト JPG、指定により PPM)を指定し、サムネイルを作成します。

AVFormatContext により動画の長さ(TimeStamp)が取得できない際(MKV 等で多い)には動画位置の指定が上手くいかない場合がありますが、それ以外の場合には上手く取得されています。

コンパイルの際には、libavcodec 等のヘッダやライブラリを指定する必要があるため、makefile による指定方法も掲載しています。
libavcodec によりデコードされ取得される画像から、PPM で保存する場合には RGB888 に変換する必要があります。
変換には、swscale を利用しています。

avthumb.cpp

makefile
実行結果

2011年5月4日水曜日

FFMPEG/libavcodec のインストール

C言語から FFMPEG に含まれる libavcodec 等のライブラリを利用するためのインストール手順です。

OSX では MacPortsfink (現在は brew)などから簡単にインストールすることができますが、ここではあえて最新の FFMPEG のソースからビルドし、インストールを行う手順を記載しておきます。

Windows の場合は、インストールには cygwin/mingw に含まれる gcc や ar コマンド等が必要なので、そちらのインストールが必要になります(説明は省略)。
また、最新のソースを取得するために git コマンドをインストールしておくと便利です(GUI版があったり、cygwin/mingw の環境に含まれている場合があります)。
cygwin からのインストールはビルド時に失敗することが多いので、mingw からインストールすることをお勧めしておきます。

その他、必要となるライブラリ等(mingw ベース)
・MinGW/MSYS
・TDM's GCC/mingw32 Builds
・MSYS Base System
・nasm/yasm
・zlib/libbz2/libgsm
・libfaac/libfaad -> あれば便利という程度
・pthreads-win32

インストールコマンド手順(git のレポジトリ先は ffmpeg のサイトから):

configure のオプションには --enable-gpl/--enable-nonfree のチェックを付けていますが、商用のライブラリを作成する際には --disable-gpl/--disable-nonfree としたり GPL が含まれるライブラリを外したりして、LGPL とする工夫が必要になります。

mingw 環境での LGPL ビルド例

ライブラリが GPL か LGPL かどうかは configure 後に表示されますが、ビルド後以下のコマンドから商用には利用できないことを確認できます(GPL/LGPL の際には GPL のライセンス等が表示されます)。

ビルドや configure に失敗する場合には、手元の環境に必要なライブラリが含まれていなかったりする場合がありますので、その際には必要なライブラリを順次インストール後、再度 configure/ビルドを行う必要があります。


mingw 環境でのインストールメモ

乱数発生に関する考察

乱数発生に関する考察です。

Numerical Recipes in C, Second Edition (1992)(7. Random Numbers) に記載されていた、乱数を発生するアルゴリズムに若干の疑念を覚えたので検証するための例題を用意しました。

以下の例題では、上記 random() 関数と、srand/rand 関数による乱数発生のアルゴリズムによる乱数発生の比較を行っています。
乱数発生を繰り返し行い、ヒストグラムを作成し、それに要する時間と平均値、標準偏差を計算しています。

実行結果

例題では、0から99までの乱数を発生する計算を 100×1000 回行っていますが、平均値は双方のアルゴリズム共に 1000 となっています。
その標準偏差を比較すると、random() 関数では srand/rand 関数のアルゴリズムよりも標準偏差が大きく、バラつきが多くなっています。
また、計算時間も srand/rand 関数の方が短くなっています。

ヒストグラムは下記のグラフのようになっています。
グラフでは、私が使っていたアルゴリズムの方(青線)が、最小値/最大値の発生率が低く、そのせいでバラつきが大きい結果になっていることが分かります。

これまで普通に使ってきましたが、標準関数を利用するようにしていかないといけないですね。

以下、同様に Numerical Recipes in C に掲載されていた乱数発生の例です。
上記と同様にヒストグラムを作成した結果、最大最小値の発生率が低くバラつきが大きい結果になりますので、標準関数を利用した方が良い結果になります。

2011年5月1日日曜日

再帰呼び出しによる階乗の計算

再帰呼び出しによる階乗(1*2*3*....*n)の計算例です。


実行結果