2013年2月24日日曜日

EPD(フルパッケージ版)におけるCythonの使い方

概要

EPDのフルパッケージ版にはCythonが含まれている。
CythonはPythonのコードをCコードに変換してくれるパッケージ。
ループ構造や条件分岐などで特に処理の遅いPythonコードをCythonでCに変換しコンパイルすることで処理の高速化が可能となったりするらしい。

とりあえず関数などをCythonでコンパイルすると良いらしい。これにより生成される実行ファイルをPython側で使うことでPythonで作った関数よりも高速に動作させることが可能となる。
もっとも、NumpyやScipyなどのパッケージは処理を計算ライブラリに投げることで高速化を実現しており、CPU資源をかなり使い切っているためCythonの恩恵は少ないはず。
別途記事を作るかもしれないが、行列-行列計算などではほとんどといっていいほど違いは出ない(今のところ時間計測をPythonで行なっているため、アプリケーションの実行時間で比べれば、また違った結果が見えるかもしれないが)。
Cythonのイメージのようなもの
とにかく、条件によっては性能が100倍以上に上昇する場合もあるので、覚えておいて損はないと思う。
EPDの場合、なんにも考えなくとも、パス等が通っているため楽に使うことができる。
ただ、Windows 7 64bitの場合、デフォルトでは32bitの実行ファイルが生成されてしまう。
これを64bitにするためにかなり調べてググったのでメモとしてまとめておく。

使い方

はじめに、C言語化したいPythonの関数を.pyx形式で作成する。
ここで、Cythonを使うためにいくつかのコードを関数に追加する必要がある。
高速化の恩恵を受けられるかどうかにも関わってくるらしいので覚えておきたい。
Pyxのコードは、例えば次のようになる。

#行列のべき乗計算
import numpy as np
cimport numpy as np
cimport cython
DOUBLE = np.float64
INT = np.int64
ctypedef np.float64_t DOUBLE_t
ctypedef np.int64_t INT_t
def MyMatPow(np.ndarray[DOUBLE_t,ndim=2] a, INT_t n):
    cdef int i
    cdef np.ndarray[DOUBLE_t,ndim=2] tempa = a
    for i in xrange(n-1):
        tempa = np.dot(tempa, a)
    return tempa
こんなかんじでソースを作っておく。

Windows 7の場合のコンパイル方法(多分XPとかでも同じ)を以下に記す。

コンパイル@32bit

EPDだと何も考えずにそのままコンパイルできる。
と言っても、いろいろめんどくさいので、setup.pyというファイルをPYXソースと同一のディレクトリに配置しておく。setup.pyの中身はこんな感じで動くと思う。

'''$ python setup.py build_ext --inplace'''
from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext
import numpy as np
setup(
    cmdclass = {'build_ext': build_ext},
    ext_modules = [Extension("出力ファイルの名前(ex: example)", ["PYXソースの名前(ex: example.pyx)"], language="c", include_dirs=[np.get_include()])]
)

とりあえずコマンドライン上でPYXソースのあるディレクトリ上に移動して、
python setup.py build_ext --inplace
と打てばコンパイルしたファイルが生成される。ここで出力ファイルの名前とPYXソースの名前は合わせておく。例えば、example.pyxという名前のソースを作った場合、出力ファイルの名前にはexampleと書く。
後は、同一ディレクトリにpythonソースを置き、
import example
とかいう風にしてやれば、example.pyxの中身の関数を使うことができる。

コンパイル @64bit

64bitでコンパイルする場合、windwosだとVisual Studioと言うか、 Windows SDK C/C++ compilerが必要になる。
とりあえず、本家に記載のあった通り、Microsoft Windows SDK for Windows 7 and .NET Framework 3.5 SP1辺りをインストールしておく。
とりあえず、C:\Program Files\Microsoft SDKs\Windows\v7.0というフォルダができていれば良い。
次に、
バッチファイルでも作って、

CD /d C:\Program Files\Microsoft SDKs\Windows\v7.0
set DISTUTILS_USE_SDK=1
C:\Windows\System32\cmd.exe /E:ON /V:ON /T:0E /K "C:\Program Files\Microsoft SDKs\Windows\v7.0\Bin\Setenv.cmd" /x64 /release
と書き、実行すると、文字が緑色のコマンドラインが出てくる。
そうしたら、コマンドライン上でPYXソースのあるディレクトリに移動して、
python setup.py build_ext --inplace
と打つと64bitの実行ライブラリが生成される。
Python上での使い方は32bitの場合と変わらない。


2013年2月23日土曜日

EPDによるパッケージのアップデート方法


EPDはパッケージ管理を一括して行なってくれる。
パッケージの管理はコマンドプロント(またはターミナル)上で
enpkg
を使って行う。

アカデミック版(恐らく有償版も)ではアップデートをする前に認証が必要になる。以下の手順で認証すれば良い。

>enpkg enstaller

>enpkg --userpass
Please enter the email address (or username) and password for your
EPD or EPD Free subscription.  If you are not subscribed to EPD,
just press Enter.

Email (or username): メールアドレス@を入力する
Password:パスワードを打ち込む
Wrote configuration file: C:\どっか
You are logged in as メールアドレス (ユーザー名).
Subscription level: EPD Basic or above

認証が正常に終了すると次のコマンドでパッケージの最新版がインストールできる。
>enpkg packagename

任意のバージョンをインストールしたい場合は、パッケージ名の後ろにバージョン番号を入れる。
>enpkg packagename 1.1.1(version number)

バージョン番号は、
>enpkg -i packagename
で調べる。

アップデート可能なパッケージを調べるには、
>enpkg --whats-new
で一覧が表示される。

分からないことがあれば、
>enpkg -h
でコマンドの説明が表示される。

2013年2月22日金曜日

Python,Matplotlibでギリシャ文字出力方法2

概要

matplotlibのギリシャ文字出力方法については以前書いたが、Tex形式で出力する方法('$\mu$'など)とフォントの保有しているギリシャ文字をUnicodeから指定する方法(u'\u03bc'等)がある。
使用するフォントにも依るが、両者でグラフの見た目が結構変わる。
特にArialやそれ系のフォント(Helvetica欲しい……)等、サンセリフ体のものを使う場合、Tセリフ体の文字が出力されるためTex方式では、ギリシャ文字に対してそれ以外の英数字が浮いてしまう。
気にならなければ構わないのだが、個人的にはグラフ全体に同一のフォントを使用したほうが統一感が出て良いと思う。
まあ、Times New Romanとかそれ系のセリフフォント(Times欲しい……)を使うならそれほど違和感はないかもしれないが、これらのフォントにもギリシャ文字が同梱されているわけで、それを使ったほうがやはり統一感が出て美しいのではないかと思う。

学会によってはギリシャ文字はSymbolを使えというところもあるようなので、その場合は仕方ないし、ギリシャ文字はSymbolじゃないとという人もいるのかもしれない。
ただ、グラフのフォントをIllustratorなどを使って差し替える場合、Symbolを使用していると一括変換ができないが、Unicodeを使えば変換先のフォントがコードに対応したフォントを持っていれば問題なく変換できてしまうというメリットがある。

比較

ギリシャ文字をTex形式とUnicode形式で出力したグラフを作成した。英数字のフォントをArialとした場合とTimes New Romanとした場合の結果を以下に示す。

Arial

Arial、ギリシャ文字はUnicode形式

Arial、ギリシャ文字はTex形式

Times New Roman

Times New Roman、ギリシャ文字はUnicode形式

Times New Roman、ギリシャ文字はTex形式

考察

Arial等のサンセリフフォントを利用する場合、ギリシャ文字をTex形式で出力するとセリフ体とサンセリフ体が1つのグラフに混在することになる。また、Tex形式で出力したデータのギリシャ文字のフォントはCmmi10及びCmr10となる。これらは通常Illustratorで編集することができない。そのため、Illustratorで加工する場合、自分で適宜Symbolにするなりしなければならない事になる。(ただし、これらのTexフォントをIllustratorで扱えるようにする方法もあるらしい) matplotlibはグラフの描画機能が優秀であるため、Illustratorを介さずに十分な品質のグラフを直接出力することが可能だが、フォントを埋め込む場合はpdf形式しか出力できない(epsだとフォントがアウトライン化されてしまう。何とかepsにフォントを埋め込む方法はないかは調査中)。そのため、面倒なくTexに埋め込むにはIllustratorでeps形式に出力し直さなければならない。Unicode方式ではそのまま保存できるが、Tex方式ではTexのフォントをIllustratorで扱えるようにしないと、フォント調整に一手間かかることになる。 一方で、Unicode形式の場合、出力したい文字のUnicodeを一々調べねばならず面倒である。また、グラフ中に数式を描きたい場合はやはりTex形式の方が優れているように思う。 どちらの方法も一長一短であるという毒にも薬にもならない戯言が結論になるだろうか。

付録

グラフ作成に使用したデータとPythonコードを以下に示す。なお、Pythonコードはあちこちのホームページやブログを参考にさせていただいているが、どこがどこを参考にしたのか全く把握できないため、特に引用の記載はない。誠に申し訳ない。

データ:


sample.txt
# This is a sample data
0 0 5
0.314159265 2.938926261 4.484011233
0.628318531 4.755282581 0.71019761
0.942477796 4.755282581 0.71019761
1.256637061 2.938926261 4.484011233
1.570796327 6.12574E-16 5
1.884955592 -2.938926261 -1.39384129
2.199114858 -4.755282581 -8.800367553
2.513274123 -4.755282581 -8.800367553
2.827433388 -2.938926261 -1.39384129
3.141592654 -1.22515E-15 5
3.455751919 2.938926261 4.484011233
3.769911184 4.755282581 0.71019761
4.08407045 4.755282581 0.71019761
4.398229715 2.938926261 4.484011233
4.71238898 1.83772E-15 5
5.026548246 -2.938926261 -1.39384129
5.340707511 -4.755282581 -8.800367553
5.654866776 -4.755282581 -8.800367553
5.969026042 -2.938926261 -1.39384129
6.283185307 -2.4503E-15 5

Pythonコード

#コメントは#マークを使うか"""hogehoge"""とする.

#ライブラリをインポート
import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
#import datetime #ここらへんのパッケージは昔必要になったことがあったようななかったような
#import matplotlib.colors as colors
#import matplotlib.finance as finance
#import matplotlib.dates as mdates
#import matplotlib.ticker as mticker
#import matplotlib.mlab as mlab
#import matplotlib.font_manager as font_manager


pi = np.pi # 円周率が必要な場合はnumpy(インポート設定からnp.~とする),このように出力できる
e = np.e # 自然対数とかも出力可能
###############################################################################
#                               Graph Setting                                 #
###############################################################################
###############################################################################
mpl.rcParams['font.family'] = 'Arial' #使用するフォント名
mpl.rcParams['pdf.fonttype'] = 42 #このおまじないでフォントが埋め込まれるようになる
params = {'backend': 'ps', # バックエンド設定
          'axes.labelsize': 8, # 軸ラベルのフォントサイズ
          'text.fontsize': 8, # テキストサイズだろう      
          'legend.fontsize': 8, # 凡例の文字の大きさ
          'xtick.labelsize': 8, # x軸の数値の文字の大きさ
          'ytick.labelsize': 8, # y軸の数値の文字の大きさ
          'text.usetex': False, # 使用するフォントをtex用(Type1)に変更
          'figure.figsize': [10/2.54, 6/2.54]} # 出力画像のサイズ(インチなので2.54で割る)
mpl.rcParams.update(params)

###############################################################################
#                                 Data Input                                  #
###############################################################################
###############################################################################
#データ読み込みにはnumpyのloadtxtを使用すると便利
filename = "sample.txt"
data = np.loadtxt(filename, #ファイルの名前
                  skiprows=1, #この場合ファイルの一行目をスキップする
                  dtype={'names':('data1',
                                  'data2',
                                  'data3'),
                        'formats':('f8',
                                   'f8',
                                   'f8'
                               )})

###############################################################################
#                                 Make Graph                                  #
###############################################################################
###############################################################################
plt.rc('axes', grid=True)#
plt.rc('grid', color='0.8', linestyle='-.', linewidth=0.5)#
bs = plt.figure(facecolor='white')
base_ax = plt.axes([0.15,0.15,0.75,0.75])
ax1 = plt.plot(data['data1'], data['data2'], 'ro-', label='Sample1')#散布図の描画.凡例のラベル設定はここで行う
ax2 = plt.plot(data['data1'], data['data3'], 'bs--', label='Sample2')
lx = plt.xlabel(u'x line [\u03bcm]')#x軸のラベル
ly = plt.ylabel(u'y line [\u03bcm]')#y軸のラベル
tx = plt.xticks(np.arange(0, 2.1*pi, 0.5*pi),#np.arange(0, 2*pi, 0.5*pi)とすると2*piが含まれない
                [u'0',u'0.5\u03C0',u'\u03C0',u'1.5\u03C0',u'2\u03C0']) #xラベルの設定
                #$\pi$表記でTex形式出力も可能
limx = plt.xlim(0, 2*pi) # xの範囲
leg = plt.legend(loc='best', shadow=True, fancybox=True)#凡例出力
leg.get_frame().set_alpha(0.5)#凡例を半透明にする
plt.savefig('sample.pdf')# PDF形式で保存する場合,フォントが埋め込まれる
plt.savefig('sample.eps')# EPS形式で保存する場合,フォントはアウトライン化される
plt.savefig('sample.png', dpi=100)# png形式ではdipの調整ができる.デフォルトでは80だったはず
plt.show()#グラフ出力.最後に持ってくる(savefigがうまく動かなかったりするため)

2013年2月19日火曜日

Matplotlibの凡例設定

昔やったはずなのに忘れていたからメモ。
物忘れが激しくて困る。

import matplotlib.pyplot as plt
plt.legend(loc='best', fancybox=True, shadow=True, ...)
要はmatplotlib.pyplot.legend()関数内で色々設定する。決まった定型文のようなものをコンマで区切って入れていく。
以下、自分が使いそうな関連をまとめておく。

凡例の位置:

基本的には
 loc='~'
と言った感じで作成。本家のページを探したらこんなのがあった。数字を入れてもいいし、loc='best'
のように書いてもいい。
Location StringLocation Code
‘best’0
‘upper right’1
‘upper left’2
‘lower left’3
‘lower right’4
‘right’5
‘center left’6
‘center right’7
‘lower center’8
‘upper center’9
‘center’10
更にbbox_to_anchor=(0.5,0.5)と書くと凡例の位置が指定した箇所からずれる。
どんな風にずれるかはプロットしながら調整すれば良いと思う。

凡例の表示方法:

numpoints=? :
 凡例に表示する点の数を指定する。デフォルトでは2になっているが、点同士を結ばない場合、2つ点が並んで違和感がでるのでよく使うとおもう。
ncol=?:
 凡例の列数。デフォルトでは1。たまに使いたくなる。
mode=“expand”:
 凡例をグラフの横いっぱいに広げるらしい。一度使ったことがある
 fancybox=True:
 凡例の角を丸くする。最初は違和感があったが、なんとなく見映えがいいような気がする
shadow=True:
 凡例に影をつける。これも違和感があったが、かっこいい気がしてきた。
title='~'
 凡例のタイトルを指定する

凡例を半透明にする

関数内でalpha=?としたらエラーが出た。この構文はpyplot.plotやpyplot.xlabelなどに対応しているがpyplot.legend関数には対応していないらしい。
そこで、今のところ、凡例を半透明にするには次のように書いている。
他に方法がないのかはそのうち調べるかも(たぶんしない)
leg = plt.legend('''色々指定する''')
leg.get_frame().set_alpha(0.5)

2013年2月17日日曜日

Matplotlibの高解像度画像保存


概要

PNG形式でグラフを保存する場合解像度はデフォルトで80dpi位になっている。
Texを使う場合はベクターデータとしてPDF形式(フォント埋め込み可)かEPS(フォントのアウトライン化)で保存して適宜加工すれば良いが、MS Wordなどを使う場合はPNGで直接出力したほうが早い。
プリント出力する場合、グラフの解像度を上げたいことがある。
しかし、普通にfigure内でdpiを宣言しても上手く行かなかった。
ちょっと調べれば方法は分かったが、メモとして残しておく。
デフォルトの解像度


高解像度(300dpi)


方法

次のように書けば良い。要はsavefig関数の中でdpiを宣言する。
import matplotlib.pyplot as plt
plt.savefig("sample.png",format = 'png', dpi=300)

2013年2月16日土曜日

このBlogについて

このブログはチラシの裏の落書き的メモ置き場です.
Pythonとかそこら辺のチップスのようなものをつらつら書いていく予定です.
飽きるまで.

EPDのMKLをマルチスレッドで実行する方法

概要

EPDでNumpyの行列-行列積をレッツノートJ10(i3 2350m)とデスクトップ(i5 2500k)で実行したが,性能が周波数分程度しか変わっていない.
どうやら,デフォルトではマルチスレッドでプログラムが走らないらしい.それでも,下手にC言語で行列演算をするよりもMKL連携のNumpyの方が高速なのだが,まさかMKLがマルチスレッドに対応していないわけがないと思うので,方法を探した.
方法発見までそこそこ苦労したのでまとめる.

方法

import mkl
mkl.set_num_threads(MyNum)

MyNumにスレッド数を入れる.これだけ.
あとは全部NumpyとMKLが上手くやってくれる.

ベンチ

デスクトップ環境(定格i5 2500k)で超簡単な行列-行列積ベンチマークを行ったので,ソースと実行結果を示す.なお,EPDは64bit版でOSはWindows7 64bit,メモリは16GB.

ソース:

import time import mkl import numpy as np mkl.set_num_threads(4) num2 = 10000 #行列サイズ #計算用の行列作成 starttime = time.clock() A = np.ones((num2,num2),'d') B = np.ones((num2,num2),'d') #行列作成の時間出力 print time.clock() - starttime #行列-行列積の実行 starttime = time.clock() C = np.dot(A, B) # 実行時間と演算性能の出力 print time.clock() - starttime print str(2*num2**3/(time.clock() - starttime)*1e-6) + u" [MFROPS]"

結果:

スレッド数1: 0.465379485104 77.7279291321 25730.7639873 [MFROPS] 
スレッド数2:
0.45703352782 39.7218425794 50350.0858396 [MFROPS]
スレッド数4: 0.358474042039 22.0854201077 90557.3095905 [MFROPS]
スレッド数8: 0.459133079081 22.1900988011 90130.1277687 [MFROPS]
インテルによるとi5 2500kの理論性のはターボ無しで105.6 [GFLOPS] = 105600 [MFLOPS]らしいので,4スレッドで実行すると85.8 %の実行効率が得られたことになる. これほどの性能が得られるとちょっと信じがたい値だ. ループのない簡単な演算を行うコードならCythonは必要無さそう.理論性能をい引き出したとしても,性能向上の余地は17 %しかないということなのだから. そして,Visual Studio C++で行列演算プログラムを作成してもせいぜい10 %しか実行効率を得られない自分は行列演算では大人しく計算ライブラリやPythonのNumpyを使うべきなのだろう.

EPDについて

Enthought Python Distribution (EPD)

Pythonは単体では殆ど使えない.
なので,必要に応じてパッケージをインストールしていくわけだが,パッケージ間の連携やバージョン管理などが面倒くさすぎる.
幸いにして,複数のパッケージをまとめて管理してくれるものがいくつかあるので,それを活用する.
今のところ自分はPythonをMatlabの代わりとして使いたい(実験データの処理とグラフ作成,正直Matlabのするべきことではない)ので,Numpy,Scipy,Matplotlibが必要となる.

自分は学術計算用Pythonパッケージの管理を一括でしてくれるEPDを使用している.PythonxyのほうがSpyderまで一括でインストールされて楽だが,後述の利用でEPDを選択した.

有償だが機能が少ない(ただし,Numpy,Scipy,Matplotlibと自分に必要なパッケージの揃った)Free版もある.
アカデミックは無償でフルパッケージの利用ができるため,自分はこれを使用している.
フルパッケージ版の特徴としては,Numpy,ScipyにおけるMKL連携,Windows版では64bitバージョンがあること(FreeではWindowsは32bitのみ),豊富なパッケージ(自分には使い道が分からない.Cythonはそのうち勉強しておきたいけど)がある.
当初はFree版でも良いと思っていたが,フルパッケージ版のMKL連携が気に入ったため,アカデミック版用に登録をした.
登録には”.ac”の入ったメールアドレスが必要となる.
具体的な手順はHPに(英語で)書いてあるし,分からないことがあればググれば良いと思う.

2013年2月2日土曜日

Python,Matplotlibでギリシャ文字出力方法


# --- グラフにギリシャ文字を書きたい場合 --- #
方法1:TeXを使う
$ドルマークで囲まれた範囲はtex形式になる$
例:
label=('$\mu$')#
μが出力される
問題点:
フォントがTeXのものになるため,Arialなどとデザインが合わない

方法2:文字コードを指定する
u'\uhoge'
''の前にuをつける(u'hoge')ことで文字コードがUTF-16に指定される模様
\uhogeでhogeに文字コードを指定する
例:
label=u'\u03bc'#使用しているフォントでμが出力される
#03bcがギリシャ文字μのUTF-16進のコード
問題点:
使用しているフォントにその文字がない場合文字が表示されない
Arialで日本語のコードを入力すると□出力を確認
ギリシャ文字を使用するならギリシャ文字を含むフォントが必要
……ギリシャ文字は大体のフォントに含まれている気がする
論文投稿でギリシャ文字はSymbolと規定されていたりする

ギリシャ文字と文字コード対応表
U+
0
1
2
3
4
5
6
7
8
9
A
B
C
D
E
F
0370
Ͱ
ͱ
Ͳ
ͳ
ʹ
͵
Ͷ
ͷ


ͺ
ͻ
ͼ
ͽ
;

0380




΄
΅
Ά
·
Έ
Ή
Ί

Ό

Ύ
Ώ
0390
ΐ
Α
Β
Γ
Δ
Ε
Ζ
Η
Θ
Ι
Κ
Λ
Μ
Ν
Ξ
Ο
03A0
Π
Ρ

Σ
Τ
Υ
Φ
Χ
Ψ
Ω
Ϊ
Ϋ
ά
έ
ή
ί
03B0
ΰ
α
β
γ
δ
ε
ζ
η
θ
ι
κ
λ
μ
ν
ξ
ο
03C0
π
ρ
ς
σ
τ
υ
φ
χ
ψ
ω
ϊ
ϋ
ό
ύ
ώ
Ϗ
03D0
ϐ
ϑ
ϒ
ϓ
ϔ
ϕ
ϖ
ϗ
Ϙ
ϙ
Ϛ
ϛ
Ϝ
ϝ
Ϟ
ϟ
03E0
Ϡ
ϡ














03F0
ϰ
ϱ
ϲ
ϳ
ϴ
ϵ
϶
Ϸ
ϸ
Ϲ
Ϻ
ϻ
ϼ
Ͻ
Ͼ
Ͽ
1F00
1F10




1F20
1F30
Ἷ
1F40




1F50




1F60
1F70
ά
έ
ή
ί
ό
ύ
ώ


1F80
1F90
1FA0
1FB0

Ά
ι
᾿
1FC0

Έ
Ή
1FD0
ΐ


Ί

1FE0
ΰ
Ύ
΅
`
1FF0



Ό
Ώ
´


貼り付け元  <