#!/usr/bin/env python
# coding: utf-8

# # matplotlibで図を描く，保存する，アニメーションさせる
# 2015, Yuto Tokunaga

# In[1]:


get_ipython().run_line_magic('matplotlib', 'inline')
import numpy as np
import scipy as sp
import scipy.constants as C
import matplotlib.pyplot as plt
import statsmodels.api as sm
from statsmodels.sandbox.regression.predstd import wls_prediction_std


# この記事は[OUCC Advent Calendar 2015 - Adventar](http://www.adventar.org/calendars/730)の12日目の本来の記事です．溜まっているクソみたいな量のタスクを無視して執筆されました．
# 
# おはようございます．綺麗なグラフを作成することにしか実験レポートのモチベーションを感じられないゆんたんです．みなさんmatplotlib使ってますか？てかmatplotlib以外で何使うの？Excel？(笑)
# 
# 冒頭から各方面(特に部長)に火種を撒いてしまった当記事ですが，まずはmatplotlibの紹介からしようと思います．
# 
# 対象となるmatplotlibのバージョンは1.5です．この記事ではjupyter notebookを使っていますが，jupyterを使わなくてもmatplotlibは使うことができます．Python環境やjupyter，matplotlibのインストール手順などについては割愛します．
# 
# ## [matplotlib](http://matplotlib.org/)とは
# Pythonで使えるグラフ(散布図，~~ここにグラフの種類をずらずら書く~~)作成ライブラリです．描けるグラフの種類については[matplotlib: python plotting — Matplotlib 1.5.0 documentation](http://matplotlib.org/#matlab)なんかを見てみてください．
# 
# まあこんな説明だけだと味気ないので色々実際に描いてみることにしましょう．コードも一緒に載せるので雰囲気を掴んでみてください．
# 
# ### 線でつなぐ，軸ラベルを設定する，凡例を入れる
# 方形波のフーリエ級数展開

# In[2]:


def furie_approx(func, nmaxs=[1, 3, 5], xr=(0, 2 * np.pi), points=200):
    x = np.linspace(float(xr[0]), float(xr[1]), points)
    for nmax in nmaxs:
        res = []
        for j in x:
            res.append(func(j, nmax))
        plt.plot(x, res, '-', label=r"nmax={0}".format(nmax))
    plt.xlim(xr)
    plt.legend(loc='best')

def func2(x, nmax):
    sumn = 0
    for n in range(1, nmax + 1):
        sumn += np.sin((2*n-1)*x) / (2*n-1)
    return 1 / 2 + 2 / np.pi * sumn

furie_approx(func2, nmaxs=[2, 10, 100], xr=(-0.5, 2 * np.pi+0.5))
plt.xticks([0,np.pi/2,np.pi,np.pi*3/2,np.pi*2], [r"$0$", r"$\frac{\pi}{2}$", r"$\pi$", r"$\frac{3\pi}{2}$", r"$2\pi$"], fontsize=20)


# ### 線やマーカーのスタイルを変更する，並べる
# `plot`の第三引数でマーカーや線のスタイルを指定できます．`'-'`で普通の線，`'--'`で破線など．`subplot`で図を並べます．

# In[3]:


def marker(xs):
    plt.plot(xs, xs, '.', label=r"$y=x$")
    plt.plot(xs, xs**2, '+', label=r"$y=x^2$")
    plt.plot(xs, xs**3, 'x', label=r"$y=x^3$")
    plt.plot(xs, xs**5, '1', label=r"$y=x^5$")
    plt.plot(xs, np.exp(xs), '2', label=r"$y=\mathrm{e}^x$")
    
def line(xs):
    plt.plot(xs, xs, '-')
    plt.plot(xs, xs**2, '--')
    plt.plot(xs, xs**3, '-.')
    plt.plot(xs, xs**5, ':')
    plt.plot(xs, np.exp(xs), '--')

plt.figure(figsize=(10,4))
plt.subplot(121)
xs = np.linspace(0.1, 1.1, 11)
marker(xs)
plt.legend(loc='best')

xs = np.linspace(0, 1.2, 100)
line(xs)
plt.ylim(-0.5,3)

plt.subplot(122)
xs = np.linspace(1, 9, 9)
marker(xs)
plt.legend(loc='best')

xs = np.linspace(0, 10, 100)
line(xs)
plt.ylim(0,1000)


# ### 注釈を入れる，数式を入れる
# `text`と`annotate`で位置を指定して注釈やテキストを入れることができます．ドルマークで囲んだテキストは$\TeX$として解釈されます．

# In[4]:


# 惑星,orbital semimajor axis[AU],orbital period[day]
data = np.loadtxt('sample_orbital.csv', delimiter=',', skiprows=1, dtype=[('obj', 'S8'), ('axis', np.float), ('period', np.float)])

approx1 = lambda x: ((x*149600000000)**3 * 4*np.pi**2 / (C.gravitational_constant*1.989*1e30))**0.5 / (60*60*24)
approx2 = lambda x: ((x*149600000000)**3 * 4*np.pi**2 / (C.gravitational_constant*1.898*1e27))**0.5 / (60*60*24)

plt.figure(figsize=(8,6))
plt.xlabel('Mean distance [AU] $a$')
plt.ylabel('Orbital period [day] $T$')
plt.xscale('log')
plt.xlim((1e-3, 1e2))
plt.yscale('log')

for x,y,s in zip(data['axis'], data['period'], data['obj']):
    plt.plot((x,), (y,), 'o')
    plt.text(x*10**0.15, y*10**-0.05, s.decode('utf-8'))
    
xs=np.logspace(-1,2,4)
plt.plot(xs, approx1(xs), '-.')
xs=np.logspace(-3,-1,3)
plt.plot(xs, approx2(xs), '-.')

plt.title("タイトル")
plt.text(1e-2, 1e5, r"$T^2=ka^3$", bbox={'facecolor':'red', 'alpha':0.5, 'pad':10}, fontsize=20)
plt.annotate('月', xy=(3e-3,30), xytext=(1e-2,1e3), arrowprops=dict(facecolor='black', shrink=0.05))
plt.text(2e-2, 2e3, "Moons of Jupiter", style='italic', rotation=41, fontsize=14)
plt.text(3, 2e4, "Planets", style='italic', rotation=41, fontsize=14)


# 他にも3Dプロットやヒストグラム，ベクトル場，行列，極座標プロットなどとにかく色々できます．[Matplotlib Examples — Matplotlib 1.5.0 documentation](http://matplotlib.org/examples/index.html)にコード例が色々あるのでぜひ見てみてください．
# 
# で，ここからが本題なんですがもう結構尺を使ってしまった．
# 
# ## matplotlibrc
# matplotlibには設定ファイルを読みこませることができます．いつも使う設定は設定ファイルに書いておくと，コードから設定を変更しなくても良いので便利です．設定ファイルは`matplotlibrc`という名前で，windows,macならユーザーディレクトリの`.matplotlib/matplotlibrc`に，linuxなら`$HOME/.config/matplotlib/matplotlibrc`に保存しておきます．設定ファイルのサンプルがmatplotlibをインストールしたディレクトリの`/matplotlib/mpl-data/matplotlibrc`にあるのでそれをコピーして使うと良いです．詳しいことは[Customizing matplotlib — Matplotlib 1.5.0 documentation](http://matplotlib.org/users/customizing.html)にあります．
# 
# ### 日本語の設定
# matplotlibは日本語のフォントを指定しないと，日本語のテキストを描画する際豆腐になってしまいます．`FontProperties`オブジェクトを使う方法もあるのですが，設定ファイルに書いてしまったほうが楽です．
# 
# ```
# font.serif: IPAPMincho
# font.sans-serif: IPAPGothic
# #font.sans-serif: Source Han Sans JP, IPAPGothic
# 
# pdf.fonttype       : 42         # Output Type 3 (Type3) or Type 42 (TrueType)
# ```
# 
# フォントの設定を変更した際は一度`$HOME/.cache/matplotlib/`以下を消してフォントのキャッシュを再作成させる必要があるようです．
# 
# ### スタイルの設定
# 図のデフォルトのスタイルの設定もできます．例えば以下のような設定をするとこの記事の図のスタイルになります．
# 
# ```
# lines.linestyle   : None       # solid line
# lines.marker      : . # the default marker
# lines.markeredgewidth  : 1     # the line width around the marker symbol
# lines.markersize  : 8            # markersize, in points
# lines.antialiased : True         # render lines in antialised (no jaggies)
# 
# patch.linewidth        : 0.5     # edge width in points
# patch.facecolor        : blue
# patch.edgecolor        : eeeeee
# patch.antialiased      : True  
# 
# axes.edgecolor      : bcbcbc   # axes edge color
# axes.grid           : True   # display grid or not
# axes.titlesize      : x-large   # fontsize of the axes title
# axes.labelsize      : large  # fontsize of the x any y labels
# #axes.color_cycle    : 348ABD, A60628, 7A68A6, 467821,D55E00,  CC79A7, 56B4E9, 009E73, F0E442, 0072B2  # color cycle for plot lines
#                                             # as list of string colorspecs:
#                                             # single letter, long name, or
#                                             # web-style hex
# axes.prop_cycle     : (cycler('color', ['348ABD', 'A60628', '7A68A6', '467821', 'D55E00', 'CC79A7', '56B4E9', '009E73', 'F0E442', '0072B2']) + cycler('marker', ['.', '+', 'x', '1', '2', '3', '4', 4, 5, 6]))
# 
# polaraxes.grid      : True    # display grid on polar axes
# axes3d.grid         : True    # display grid on 3d axes
# 
# legend.fancybox      : True  # if True, use a rounded box for the
# legend.framealpha     : 0.8  # Default transparency of the legend frame box
# ```
# 
# matplotlib 1.5から入った`axes.prop_cycle`を使うと，複数の系列をプロットした時に，色やマーカーをどの順番で使うかを指定できます．スタイルのサンプルがいくつか[daler/matplotlibrc](https://github.com/daler/matplotlibrc)にあるので見ながら色々調整すると良いでしょう．
# 
# ## スタイリング
# matplotlibrcを使う以外にも図のスタイルを変更する方法があります．matplotlib 1.4からは`style`パッケージが追加され，以下のようにしてスタイルを変更することができます．

# In[5]:


from matplotlib import style
style.use('fivethirtyeight')


# In[6]:


xs = np.linspace(0,2*np.pi,100)
plt.plot(xs, np.sin(xs), '-')


# 使えるスタイルは`style.available`で調べられます

# In[7]:


style.available


# ## 図の出力
# ### 画面に表示
# `matplotlib.use('backend')`で出力に使用するバックエンド(表示するためのインタフェース)を指定できます．利用できるバックエンドは環境により異なります．詳しくは[What is a backend?](http://matplotlib.org/faq/usage_faq.html#what-is-a-backend)．また，図を表示するには`plt.show()`する必要があります．
# 
# ### jupyter notebookで表示
# notebook上で`%matplotlib inline`を実行するとmatplotlibで描画した図が埋め込まれます．また，`%matplotlib notebook`を実行すると図をマウスで掴んで動かしたり，軸のスケールを変えることができます．しかし移動したりスケールを変えるごとにサーバーで図を再描画して表示する仕様のため，結構重く感じると思います．
# 
# ### ファイルに出力
# `plt.save('filename.ext')`で図をファイルに出力することができます．拡張子としてsvg，png，pdf等指定すればその形式で出力されます．図をあとから拡大・縮小する場合はSVG，PDF等のベクタグラフィックス形式を使いましょう．
# 
# また，$\LaTeX$のTikzというパッケージを用いて，matplotlibで描画した図を$\LaTeX$の文章に埋め込むこともできます．Tikzで扱える形式に変換するには拡張子を`pgf`にしてください．詳しくはこちらに書いてあります．
# 
# [matplotlibはTikz(PGF)出力ができる - 気まぐれLinux](http://yuntan-t.hateblo.jp/entry/2015/08/22/115505)
# 
# この方法を用いると非常に高品質な文章が得られるのですが，複雑な図を入れると$\LaTeX$文章のタイプセットにかなり時間がかかってしまいます．個人的にはPDFで出力して`\includegraphics`で埋め込む方法をおすすめします．PDFで出力しておくと拡大・縮小してもPNGのようにピクセル感が出てくることはないです．
# 
# ## アニメーション
# `matplotlib.animation`パッケージを使えば，パラパラマンガを作ることができます．アニメーションを作成するには2通りの方法があります．
# 
# * 前に用意したデータ(Artist)を用いる `ArtistAnimation`
# * データを更新する関数を用意する `FuncAnimation`
# 
# どちらの方法がやりやすいかは場合によると思います．ちょっと難しめですがWebでサンプルコードがいろいろ見つかるのでチャレンジしてみてください．
# 
# ### jupyter notebookで表示
# 以前まではnotebookに動画を表示しようとするとひと手間必要だった(アニメーションを保存し，base64でエンコードして`video`タグに埋め込むコードを書く必要があった)のですが，matplotlib 1.5から`Animation.to_html5_video`メソッドが追加され，簡単になりました．notebook上で動画を再生するにはmatplotlibrcの`animation.html`を`'html5'`にするか，又はnotebook上で以下の設定をします．

# In[8]:


import matplotlib
matplotlib.rcParams['animation.html'] = 'html5'


# 引数の`repeat`で繰り返し再生するかを指定できます．以下は`FuncAnimation`を用いて定常波を描画する例です．

# In[9]:


from matplotlib.animation import FuncAnimation

fig = plt.figure()
ax = plt.axes(xlim=(0, 20), ylim=(-2, 2), xlabel=r"$x$", ylabel=r"$y$")
rline, = ax.plot([], [], '-', lw=2)  # 合成波
gline, = ax.plot([], [], '--', lw=2)  # 前進波
bline, = ax.plot([], [], '--', lw=2)  # 後進波
x = np.linspace(0, 20, 1000)

# initialization function: plot the background of each frame
def init():
    rline.set_data([], [])
    gline.set_data([], [])
    bline.set_data([], [])
    return rline, gline, bline

# animation function.  This is called sequentially
def animfunc(t):
    ax.set_title(r"$k=2\pi/10,\omega=2\pi/100,t={0}$".format(t))
    uplus  = np.cos(2 * np.pi / 10 * x - 2 * np.pi / 100 * t)
    uminus = np.cos(2 * np.pi / 10 * x + 2 * np.pi / 100 * t)
    gline.set_data(x, uplus)
    bline.set_data(x, uminus)
    rline.set_data(x, uplus + uminus)
    return rline, gline, bline

FuncAnimation(fig, animfunc, init_func=init, frames=101, interval=100, blit=True, repeat=False)


# アニメーションを使ってライフゲームをやってみた例です． http://nbviewer.ipython.org/gist/yuntan/0ad2431123718ff83ae5
# 
# ## おまけ

# In[10]:


# Based on "Self-Description" from XKCD by Randall Monroe
# http://xkcd.com/688/

style.use('grayscale')

with plt.xkcd():
    plt.figure(figsize=(10, 4))
    plt.subplot(121)
    plt.pie([1,10], colors=['black','white'], center=(-100,0), labels=["FRACTION OF\nTHIS IMAGE\nWHICH IS BLACK","FRACTION OF\nTHIS IMAGE\nWHICH IS WHITE"], startangle=200)
    
    ax = plt.subplot(122)
    plt.bar([1,2], [4,9], align="center")
    plt.xticks([1,2], ['1','2'])
    plt.yticks([],[])
    plt.ylabel("AMOUNT OF BLACK INK\nBY PANEL")
    plt.ylim(0, 10)
    ax.spines['right'].set_color('none')
    ax.spines['top'].set_color('none')
    ax.xaxis.set_ticks_position('bottom')
    plt.suptitle("Self-Description")


# ## 参考
# * [What’s new in matplotlib — Matplotlib 1.5.0 documentation](http://matplotlib.org/users/whats_new.html#display-hook-for-animations-in-the-ipython-notebook)
# * [pythonで散布図アニメーションを試してみた - 株式会社CFlatの明後日スタイルのブログ](http://cflat-inc.hatenablog.com/entry/2014/03/17/214719)
# * [animation — Matplotlib 1.5.0 documentation](http://matplotlib.org/api/animation_api.html)
# * [XKCD Plots in Matplotlib: Going the Whole Way](http://jakevdp.github.io/blog/2013/07/10/XKCD-plots-in-matplotlib/)

# In[ ]: