Jupyter | Matplotlib > 凹凸の激しい形状の可視化 > v0.4, v0.5: ADDA InitField-Yファイル読込み > Sphere | Prism (N=3) | Chevyshev particles

```txt:動作環境
GeForce GTX 1070 (8GB)
ASRock Z170M Pro4S [Intel Z170chipset]
Ubuntu 16.04 LTS desktop amd64
TensorFlow v1.2.1
cuDNN v5.1 for Linux
CUDA v8.0
Python 3.5.2
IPython 6.0.0 -- An enhanced Interactive Python.
gcc (Ubuntu 5.4.0-6ubuntu1~16.04.4) 5.4.0 20160609
GNU bash, version 4.3.48(1)-release (x86_64-pc-linux-gnu)
```

関連: [Jupyter | Matplotlib > 凹凸の激しい形状の可視化 > ADDA:IntField-Yファイルからの読込み v0.1 > A sphere | Chebyshev particle の表示](http://qiita.com/7of9/items/ec9e4bd04a2449e54099)


[Jupyter | Matplotlib > 凹凸の激しい形状の可視化 > v0.1-0.3: 2つの球を並べる | 8つの球 | 3000球 + 色](http://qiita.com/7of9/items/34a6a200934d40d4de11)
の続き。


## v0.4 > ADDA InitField-Y対応

```py:showChebyshev_170910.ipynb
import matplotlib.pyplot as plt
from matplotlib import cm, colors
from mpl_toolkits.mplot3d import Axes3D
from pylab import rcParams
import numpy as np
import time

"""
v0.4 Sep. 10, 2017
  - read ADDA file [IntField-Y]
v0.3 Sep. 10, 2017
  - set colors in X direction
  - increase number of spheres to 3000
v0.2 Sep. 10, 2017
  - show 8 spheres
  - lower the resolution of the sphere (from 100j to 6j)
v0.1 Sep. 10, 2017
  - show 2 spheres
"""

# coding rule: PEP8

rcParams['figure.figsize'] = 15, 10


# reference
# https://stackoverflow.com/questions/31768031/plotting-points-on-the-surface-of-a-sphere-in-pythons-matplotlib


def plot_spheres(xps, yps, zps):
    for elem in zip(xps, yps, zps):
        axp, ayp, azp = elem
        # print(elem)
        dx = x + axp
        dy = y + ayp
        dz = z + azp
        ax.plot_surface(
            dx, dy, dz,  rstride=1, cstride=1, color='c', 
            alpha=1.0, linewidth=0,
            facecolors=plt.cm.Set2((dx - 0) / (50 - 0)))  # 50: arbitrary chosen to set colors

start_time = time.time()
        
# Create a sphere
r = 1
pi = np.pi
cos = np.cos
sin = np.sin
phi, theta = np.mgrid[0.0:pi:6j, 0.0:2.0*pi:6j]
x = r*sin(phi)*cos(theta)
y = r*sin(phi)*sin(theta)
z = r*cos(phi)


# read from ADDA file
INPFILE = 'IntField-Y'
SHIFT_POS = 10.0
dat = np.genfromtxt(INPFILE, delimiter=' ', skip_header=1)
xpar, ypar, zpar = [], [], []
for elem in dat:
    axp, ayp, azp = elem[:3]
    xpar += [(axp + SHIFT_POS) * 3.0]  # 3.0: arbitrary chosen to adjust sphere positions
    ypar += [(ayp + SHIFT_POS) * 3.0]
    zpar += [(azp + SHIFT_POS) * 3.0]

# Set colours and render
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

SKIP_NUM = 2
xp = xpar[::SKIP_NUM]
yp = ypar[::SKIP_NUM]
zp = zpar[::SKIP_NUM]

plot_spheres(xp, yp, zp)

ax.set_xlim([0, 50])
ax.set_ylim([0, 50])
ax.set_zlim([0, 50])
ax.set_aspect("equal")
plt.tight_layout()
plt.show()

duration = time.time() - start_time
print("%.3f" % duration)

```

### Result > Sphere

![qiita.png](https://qiita-image-store.s3.amazonaws.com/0/32870/87bc30f0-1a83-90c4-b6c6-7cfcb51f69f5.png)

なかなか良い[Borg Sphere](http://memory-alpha.wikia.com/wiki/Borg_sphere)だ。

処理時間: 48.053秒。

### Result > Prism (N=3)

![qiita.png](https://qiita-image-store.s3.amazonaws.com/0/32870/fbf60043-51f0-f5f8-50dd-302f7f38d231.png)

処理時間: 106.506秒。

### Result > Chevyshev particle

![qiita.png](https://qiita-image-store.s3.amazonaws.com/0/32870/4e6e7587-cb25-311c-44dd-6c40dfad6117.png)

処理時間: 21.136秒。

やはりこういう形状だったのか。

Chevyshev particleと言っても、パラメータによって色々な形にはなる。

```txt:log
command: './adda -grid 26 -shape chebyshev 0.6 5 -store_int_field '
lambda: 6.283185307
shape: axisymmetric chebyshev particle; size along x-axis (Dx):10.77580685, amplitude eps=0.6, order n=5, initial radius r0/Dx=0.3269304082
box dimensions: 26x26x24
```

## v0.5 > ファイル保存

### 変更点

elevationを変更した場合の画像ファイル保存を追加。

```py
for idx, elev in enumerate(range(60, -61, -30)):
    ax.view_init(elev, -60)
    plt.show()
    plt.draw()
    filename = 'out%d.png' % idx
    fig.savefig(filename)
```

### code全部

```py:showChebyshev_170910.ipynb
import matplotlib.pyplot as plt
from matplotlib import cm, colors
from mpl_toolkits.mplot3d import Axes3D
from pylab import rcParams
import numpy as np
import time

"""
v0.5 Sep. 10, 2017
   - output to files with different elevations
   - plot_spheres() takes [dstax] arg
v0.4 Sep. 10, 2017
  - read ADDA file [IntField-Y]
v0.3 Sep. 10, 2017
  - set colors in X direction
  - increase number of spheres to 3000
v0.2 Sep. 10, 2017
  - show 8 spheres
  - lower the resolution of the sphere (from 100j to 6j)
v0.1 Sep. 10, 2017
  - show 2 spheres
"""

# coding rule: PEP8

rcParams['figure.figsize'] = 15, 10


# reference
# https://stackoverflow.com/questions/31768031/plotting-points-on-the-surface-of-a-sphere-in-pythons-matplotlib


def plot_spheres(xps, yps, zps, dstax):
    for elem in zip(xps, yps, zps):
        axp, ayp, azp = elem
        # print(elem)
        dx = x + axp
        dy = y + ayp
        dz = z + azp
        dstax.plot_surface(
            dx, dy, dz,  rstride=1, cstride=1, color='c',
            alpha=1.0, linewidth=0,
            facecolors=plt.cm.Set2((dx - 0) / (50 - 0)))  # 50: arbitrary chosen to set colors

start_time = time.time()

# Create a sphere
r = 1
pi = np.pi
cos = np.cos
sin = np.sin
phi, theta = np.mgrid[0.0:pi:6j, 0.0:2.0*pi:6j]
x = r*sin(phi)*cos(theta)
y = r*sin(phi)*sin(theta)
z = r*cos(phi)


# read from ADDA file
INPFILE = 'IntField-Y'
SHIFT_POS = 10.0
TOP_VIEW = True  # False: bottom view

dat = np.genfromtxt(INPFILE, delimiter=' ', skip_header=1)
xpar, ypar, zpar = [], [], []
for elem in dat:
    axp, ayp, azp = elem[:3]
    xpar += [(axp + SHIFT_POS) * 3.0]  # 3.0: arbitrary chosen to adjust sphere positions
    ypar += [(ayp + SHIFT_POS) * 3.0]
    zpar += [(azp + SHIFT_POS) * 3.0]

# Set colours and render
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

SKIP_NUM = 2
xp = xpar[::SKIP_NUM]
yp = ypar[::SKIP_NUM]
zp = zpar[::SKIP_NUM]

plot_spheres(xp, yp, zp, dstax=ax)

if not TOP_VIEW:
    ax.view_init(-30, -60)

ax.set_xlim([0, 50])
ax.set_ylim([0, 50])
ax.set_zlim([0, 50])
ax.set_aspect("equal")
plt.tight_layout()
plt.show()

for idx, elev in enumerate(range(60, -61, -30)):
    ax.view_init(elev, -60)
    plt.show()
    plt.draw()
    filename = 'out%d.png' % idx
    fig.savefig(filename)

duration = time.time() - start_time
print("%.3f" % duration)
```

### 画像ファイル閲覧

実行後、以下のようにする。

```
$ eog out*.png
```

カーソルキーの右左で異なるelevationにて見ることができる。

![movie.gif](https://qiita-image-store.s3.amazonaws.com/0/32870/c2ee8760-7aa8-27c8-c0b1-25ae2d882066.gif)

上記は以下のパラメータで作成した。

```
command: './adda -grid 52 -shape chebyshev 0.7 12 -store_int_field '
lambda: 6.283185307
shape: axisymmetric chebyshev particle; size along x-axis (Dx):21.83684002, amplitude eps=0.7, order n=12, initial radius r0/Dx=0.2941176471
box dimensions: 52x52x52
```