""" 単位格子を描画するスクリプト。 概要: 結晶学における単位格子を3Dで可視化する機能を提供します。 詳細説明: 任意の格子定数と角度から単位格子の形状を計算し、Matplotlibを使用して描画します。 基本ベクトル（a, b, c）の描画には、パースペクティブを考慮したカスタム3D矢印クラスを使用します。 関連リンク: draw_unit_cell_usage """ import numpy as np import matplotlib.pyplot as plt from mpl_toolkits.mplot3d import Axes3D from matplotlib.patches import FancyArrowPatch from mpl_toolkits.mplot3d import Axes3D, proj3d class Arrow3D(FancyArrowPatch): """ 概要: Matplotlibの3Dプロットでカスタムの3D矢印を描画するためのクラス。 詳細説明: FancyArrowPatch を継承し、3D空間内の座標を2Dスクリーン座標に変換して描画することで、 常に適切なパースペクティブで矢印が表示されるようにします。 """ def __init__(self, xs, ys, zs, *args, **kwargs): """ 概要: Arrow3Dクラスのコンストラクタ。 引数: :param xs: 矢印のX座標のリストまたはタプル。 :type xs: list[float]またはtuple[float] :param ys: 矢印のY座標のリストまたはタプル。 :type ys: list[float]またはtuple[float] :param zs: 矢印のZ座標のリストまたはタプル。 :type zs: list[float]またはtuple[float] :param args: FancyArrowPatch に渡される追加の引数。 :type args: tuple :param kwargs: FancyArrowPatch に渡される追加のキーワード引数。 :type kwargs: dict """ super().__init__((0, 0), (0, 0), *args, **kwargs) self._verts3d = xs, ys, zs def draw(self, renderer): """ 概要: 矢印をレンダリングします。 詳細説明: 3D座標を2Dスクリーン座標に変換し、変換された座標を使用して矢印を描画します。 引数: :param renderer: Matplotlibのレンダラーオブジェクト。 :type renderer: matplotlib.backend_bases.RendererBase 戻り値: :returns: なし :rtype: None """ xs3d, ys3d, zs3d = self._verts3d proj = self.axes.get_proj() xs, ys, zs = proj3d.proj_transform(xs3d, ys3d, zs3d, proj) self.set_positions((xs[0], ys[0]), (xs[1], ys[1])) super().draw(renderer) def do_3d_projection(self, renderer=None): """ 概要: 3Dプロジェクションを実行し、矢印の深度を計算します。 詳細説明: 3D座標を2Dスクリーン座標に変換し、Z軸方向の最小値（深度）を返すことで、 3Dビューにおけるオブジェクトの描画順序を適切に管理します。 引数: :param renderer: Matplotlibのレンダラーオブジェクト（オプション）。 :type renderer: matplotlib.backend_bases.RendererBase 戻り値: :returns: 矢印のZ軸方向の最小値（深度）。 :rtype: numpy.float """ xs3d, ys3d, zs3d = self._verts3d proj = self.axes.get_proj() xs, ys, zs = proj3d.proj_transform(xs3d, ys3d, zs3d, proj) self.set_positions((xs[0], ys[0]), (xs[1], ys[1])) return np.min(zs) def draw_vector(ax, vec, color, label, fontsize): """ 概要: 3D空間にベクトルと対応するラベルを描画します。 詳細説明: Arrow3D クラスを利用して、原点から指定されたベクトルまでの矢印を描画し、 ベクトルの先端にラベルを配置します。 引数: :param ax: 3DプロットのAxesオブジェクト。 :type ax: matplotlib.axes.Axes :param vec: 描画するベクトルの3D座標（[x, y, z]）。 :type vec: numpy.ndarray :param color: 矢印とラベルの色。 :type color: str :param label: ベクトルに付けるラベル。 :type label: str :param fontsize: ラベルのフォントサイズ。 :type fontsize: int 戻り値: :returns: なし :rtype: None """ arrow = Arrow3D([0, vec[0]], [0, vec[1]], [0, vec[2]], mutation_scale=20, lw=2, arrowstyle="-|>", color=color) ax.add_artist(arrow) ax.text(vec[0] * 1.05, vec[1] * 1.05, vec[2] * 1.05, label, fontsize = fontsize, color = color) def lattice_vectors(a, b, c, alpha, beta, gamma): """ 概要: 格子定数と角度から基本格子ベクトルを計算します。 詳細説明: 結晶学で用いられる格子定数 (a, b, c) と軸間角 (alpha, beta, gamma) を用いて、 直交座標系における3つの基本格子ベクトル (va, vb, vc) を計算します。 引数: :param a: 格子定数aの長さ。 :type a: float :param b: 格子定数bの長さ。 :type b: float :param c: 格子定数cの長さ。 :type c: float :param alpha: b軸とc軸の間の角度（度数）。 :type alpha: float :param beta: a軸とc軸の間の角度（度数）。 :type beta: float :param gamma: a軸とb軸の間の角度（度数）。 :type gamma: float 戻り値: :returns: 3つの基本格子ベクトル (va, vb, vc)。各ベクトルは [x, y, z] の形式の numpy.ndarray。 :rtype: tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray] """ alpha_r, beta_r, gamma_r = np.radians([alpha, beta, gamma]) va = np.array([a, 0, 0]) vb = np.array([b * np.cos(gamma_r), b * np.sin(gamma_r), 0]) cx = c * np.cos(beta_r) cy = c * (np.cos(alpha_r) - np.cos(beta_r) * np.cos(gamma_r)) / np.sin(gamma_r) cz = np.sqrt(max(0, c**2 - cx**2 - cy**2)) # 負の平方根を避けるための安全策 vc = np.array([cx, cy, cz]) return va, vb, vc def set_equal_aspect(ax, points): """ 概要: 3Dプロットのアスペクト比を均等に設定し、描画範囲を調整します。 詳細説明: 描画される点の最大・最小座標に基づいて、X, Y, Z軸の表示範囲を同じにします。 引数: :param ax: 3DプロットのAxesオブジェクト。 :type ax: matplotlib.axes.Axes :param points: 描画される全ての点のNumpy配列。形状は (N, 3)。 :type points: numpy.ndarray 戻り値: :returns: なし :rtype: None """ xlim = [np.min(points[:,0]), np.max(points[:,0])] ylim = [np.min(points[:,1]), np.max(points[:,1])] zlim = [np.min(points[:,2]), np.max(points[:,2])] max_range = max( xlim[1] - xlim[0], ylim[1] - ylim[0], zlim[1] - zlim[0] ) / 2 mid_x = np.mean(xlim) mid_y = np.mean(ylim) mid_z = np.mean(zlim) ax.set_xlim(mid_x - max_range, mid_x + max_range) ax.set_ylim(mid_y - max_range, mid_y + max_range) ax.set_zlim(mid_z - max_range, mid_z + max_range) ax.set_box_aspect([1,1,1]) def draw_unit_cell_plot(a, b, c, alpha, beta, gamma, elev=5.90, azim=-67.55, fontsize=24, color='blue'): """ 概要: 指定された格子定数と角度を持つ単位格子を3Dで描画します。 詳細説明: lattice_vectors 関数を使用して基本格子ベクトルを計算し、それらから単位格子の8つの頂点を生成します。 その後、Matplotlibを用いて単位格子の辺、頂点、および基本格子ベクトルを3D空間に描画します。 引数: :param a: 格子定数aの長さ。 :type a: float :param b: 格子定数bの長さ。 :type b: float :param c: 格子定数cの長さ。 :type c: float :param alpha: b軸とc軸の間の角度（度数）。 :type alpha: float :param beta: a軸とc軸の間の角度（度数）。 :type beta: float :param gamma: a軸とb軸の間の角度（度数）。 :type gamma: float :param elev: 3Dビューの仰角（度数）。デフォルトは5.90。 :type elev: float :param azim: 3Dビューの方位角（度数）。デフォルトは-67.55。 :type azim: float :param fontsize: ラベルのフォントサイズ。デフォルトは24。 :type fontsize: int :param color: ベクトルの色。デフォルトは 'blue'。 :type color: str 戻り値: :returns: なし（グラフを表示し、PNGとして保存します） :rtype: None """ va, vb, vc = lattice_vectors(a, b, c, alpha, beta, gamma) origin = np.array([0, 0, 0]) points = [ origin, va, vb, vc, va + vb, va + vc, vb + vc, va + vb + vc ] points = np.array(points) edges = [ (0,1), (0,2), (0,3), (1,4), (1,5), (2,4), (2,6), (3,5), (3,6), (4,7), (5,7), (6,7) ] fig = plt.figure(figsize=(8,6)) ax = fig.add_subplot(111, projection='3d') ax.view_init(elev=elev, azim=azim) # インタラクティブな角度確認用のハンドラ（直接実行時のみ意味を持つ） def on_draw(event): curr_elev = ax.elev curr_azim = ax.azim print(f"elev: {curr_elev:.2f}, azim: {curr_azim:.2f}") fig.canvas.mpl_connect('draw_event', on_draw) # 単位格子の辺（細い黒線） for i,j in edges: ax.plot(*zip(points[i], points[j]), color='black', linewidth=0.5) # 頂点（黒い小さな●） ax.scatter(points[:,0], points[:,1], points[:,2], color='black', s=20) # 基本ベクトル（太い矢印）とラベル for vec, label in zip([va, vb, vc], ['a', 'b', 'c']): draw_vector(ax, vec, color, label, fontsize) # グラフの見た目を整える ax.set_axis_off() ax.grid(False) ax.xaxis.pane.set_visible(False) ax.yaxis.pane.set_visible(False) ax.zaxis.pane.set_visible(False) ax.set_facecolor((1,1,1,0)) # 背景透明 set_equal_aspect(ax, points) plt.tight_layout() plt.savefig("unit_cell.png", dpi=300, bbox_inches='tight', transparent=True) plt.show() def main(): """ 概要: メインの実行ルーチン。特定の格子定数で単位格子を描画します。 """ # 始点の初期値やフォントサイズなどの設定 params = { 'a': 5.0, 'b': 5.5, 'c': 4.5, 'alpha': 80, 'beta': 70, 'gamma': 100, 'elev': 5.90, 'azim': -67.55, 'fontsize': 24, 'color': 'blue' } print(f"Drawing unit cell with a={params['a']}, b={params['b']}, c={params['c']}") print(f"Angles: alpha={params['alpha']}, beta={params['beta']}, gamma={params['gamma']}") draw_unit_cell_plot(**params) if __name__ == "__main__": main()