0
0

ゲーム数学Pythonコード

Last updated at Posted at 2024-08-02

ボックスミューラ法

import numpy as np
import matplotlib.pyplot as plt

def box_muller(n):
    # Generate uniform random variables U1 and U2
    u1 = np.random.uniform(0, 1, n)
    u2 = np.random.uniform(0, 1, n)
    
    # Apply the Box-Muller transform
    z0 = np.sqrt(-2 * np.log(u1)) * np.cos(2 * np.pi * u2)
    z1 = np.sqrt(-2 * np.log(u1)) * np.sin(2 * np.pi * u2)
    
    return z0, z1

# Number of random variables to generate
n = 1000
z0, z1 = box_muller(n)

# Plotting the results
plt.figure(figsize=(12, 6))

# Plot histogram of Z0
plt.subplot(1, 2, 1)
plt.hist(z0, bins=30, density=True, alpha=0.6, color='g')
plt.title('Histogram of Z0')
plt.xlabel('Value')
plt.ylabel('Density')

# Plot histogram of Z1
plt.subplot(1, 2, 2)
plt.hist(z1, bins=30, density=True, alpha=0.6, color='b')
plt.title('Histogram of Z1')
plt.xlabel('Value')
plt.ylabel('Density')

plt.tight_layout()
plt.show()



LCG

import matplotlib.pyplot as plt

def lcg(seed, a, c, m, n):
    values = []
    x = seed
    for _ in range(n):
        x = (a * x + c) % m
        values.append(x / m)  # [0, 1) の範囲にスケーリング
    return values

# 例として、シード=1, a=1103515245, c=12345, m=2^31 を使用
seed = 1
a = 1103515245
c = 12345
m = 2**31
n = 100  # 生成する乱数の数

random_values = lcg(seed, a, c, m, n)

# 乱数をプロットする
plt.figure(figsize=(10, 5))
plt.plot(random_values, marker='o', linestyle='-', markersize=4)
plt.title('LCG Random Numbers')
plt.xlabel('Index')
plt.ylabel('Random Value')
plt.grid(True)
plt.show()

import matplotlib.pyplot as plt

def pokemon_lcg(seed, A, B, M, n):
    values = []
    x = seed
    for _ in range(n):
        x = (A * x + B) % M
        # 上位4桁を抽出
        high_bits = x >> 16
        values.append(high_bits)
    return values

# ポケモンシリーズのLCGパラメータ
A = 0x41c64e6d
B = 0x6073
M = 0x100000000  # 2^32

# シード値と生成する乱数の数
seed = 1
n = 100  # 生成する乱数の数

# 乱数を生成
random_values = pokemon_lcg(seed, A, B, M, n)

# 乱数をプロットする
plt.figure(figsize=(10, 5))
plt.plot(random_values, marker='o', linestyle='-', markersize=4)
plt.title('Pokémon LCG Random Numbers')
plt.xlabel('Index')
plt.ylabel('Pseudo-Random Value')
plt.grid(True)
plt.show()

クォータニオン、共役クォータニオン、逆クォータニオン

import numpy as np
import matplotlib.pyplot as plt

# Define a quaternion class
class Quaternion:
    def __init__(self, w, x, y, z):
        self.w = w
        self.x = x
        self.y = y
        self.z = z

    def conjugate(self):
        return Quaternion(self.w, -self.x, -self.y, -self.z)
    
    def norm_squared(self):
        return self.w**2 + self.x**2 + self.y**2 + self.z**2
    
    def inverse(self):
        norm_sq = self.norm_squared()
        return Quaternion(self.w / norm_sq, -self.x / norm_sq, -self.y / norm_sq, -self.z / norm_sq)
    
    def __str__(self):
        return f"{self.w} + {self.x}i + {self.y}j + {self.z}k"

# Create a quaternion
q = Quaternion(1, 2, 3, 4)

# Compute conjugate and inverse
q_conj = q.conjugate()
q_inv = q.inverse()

# Print the results
print(f"Quaternion: {q}")
print(f"Conjugate: {q_conj}")
print(f"Inverse: {q_inv}")

# Function to plot quaternion components
def plot_quaternion(q, title):
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.bar(['w', 'x', 'y', 'z'], [q.w, q.x, q.y, q.z], color=['r', 'g', 'b', 'y'])
    ax.set_title(title)
    ax.set_xlabel('Components')
    ax.set_ylabel('Value')
    ax.set_zlabel('Magnitude')
    plt.show()

# Plot the quaternion, conjugate, and inverse
plot_quaternion(q, 'Quaternion')
plot_quaternion(q_conj, 'Conjugate Quaternion')
plot_quaternion(q_inv, 'Inverse Quaternion')

ロドリゲスの回転公式


import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

def rodrigues_rotation(v, k, theta):
    k = k / np.linalg.norm(k)  # Normalize the rotation axis
    v_rot = (v * np.cos(theta) + 
             np.cross(k, v) * np.sin(theta) + 
             k * np.dot(k, v) * (1 - np.cos(theta)))
    return v_rot

# Define the original vector and the rotation axis
v = np.array([1, 0, 0])
k = np.array([0, 0, 1])
theta = np.pi / 4  # 45 degrees

# Compute the rotated vector
v_rot = rodrigues_rotation(v, k, theta)

# Plotting
fig = plt.figure(figsize=(12, 6))
ax = fig.add_subplot(111, projection='3d')

# Original vector
ax.quiver(0, 0, 0, v[0], v[1], v[2], color='r', label='Original Vector')

# Rotated vector
ax.quiver(0, 0, 0, v_rot[0], v_rot[1], v_rot[2], color='b', label='Rotated Vector')

# Rotation axis
ax.quiver(0, 0, 0, k[0], k[1], k[2], color='g', linestyle='--', label='Rotation Axis')

ax.set_xlim([-1, 1])
ax.set_ylim([-1, 1])
ax.set_zlim([-1, 1])
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.set_title('Rodrigues Rotation Formula')
ax.legend()
plt.show()





import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d.art3d import Line3DCollection

def plot_transformation(T, origin=[0, 0, 0], label=''):
    # Extract rotation matrix and translation vector
    R = T[:3, :3]
    d = T[:3, 3]
    
    # Define the base coordinate system
    base = np.array([origin, origin + np.array([1, 0, 0]), 
                     origin + np.array([0, 1, 0]), 
                     origin + np.array([0, 0, 1])])
    
    # Define the transformed coordinate system
    transformed = np.array([origin, 
                             origin + R @ np.array([1, 0, 0]) + d, 
                             origin + R @ np.array([0, 1, 0]) + d, 
                             origin + R @ np.array([0, 0, 1]) + d])
    
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    
    # Plot the base coordinate system
    base_lines = [[base[0], base[i]] for i in range(1, 4)]
    ax.add_collection3d(Line3DCollection(base_lines, colors='k', linestyles='--'))
    
    # Plot the transformed coordinate system
    transformed_lines = [[transformed[0], transformed[i]] for i in range(1, 4)]
    ax.add_collection3d(Line3DCollection(transformed_lines, colors='r', linestyles='-'))
    
    # Plot origin and axes
    ax.plot([base[0, 0], base[1, 0]], [base[0, 1], base[1, 1]], [base[0, 2], base[1, 2]], 'k--')
    ax.plot([base[0, 0], base[2, 0]], [base[0, 1], base[2, 1]], [base[0, 2], base[2, 2]], 'k--')
    ax.plot([base[0, 0], base[3, 0]], [base[0, 1], base[3, 1]], [base[0, 2], base[3, 2]], 'k--')
    
    ax.plot([transformed[0, 0], transformed[1, 0]], [transformed[0, 1], transformed[1, 1]], [transformed[0, 2], transformed[1, 2]], 'r-')
    ax.plot([transformed[0, 0], transformed[2, 0]], [transformed[0, 1], transformed[2, 1]], [transformed[0, 2], transformed[2, 2]], 'r-')
    ax.plot([transformed[0, 0], transformed[3, 0]], [transformed[0, 1], transformed[3, 1]], [transformed[0, 2], transformed[3, 2]], 'r-')
    
    ax.set_xlabel('X')
    ax.set_ylabel('Y')
    ax.set_zlabel('Z')
    ax.set_title(f'Transformation Matrix: {label}')
    plt.legend(['Base Coordinate System', 'Transformed Coordinate System'])
    plt.show()

# Define a rotation matrix and translation vector
theta = np.pi / 6  # 30 degrees
R = np.array([[np.cos(theta), -np.sin(theta), 0],
              [np.sin(theta),  np.cos(theta), 0],
              [0,             0,            1]])
d = np.array([1, 2, 3])  # Translation vector

# Construct the transformation matrix
T = np.eye(4)
T[:3, :3] = R
T[:3, 3] = d

# Plot the transformation matrix
plot_transformation(T, label='30 Degrees Rotation Around Z Axis with Translation')


import numpy as np
from scipy.spatial.transform import Rotation as R
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

def slerp(q1, q2, t):
    """
    Perform Spherical Linear Interpolation (SLERP) between two quaternions.
    
    :param q1: The starting quaternion (4-element array).
    :param q2: The ending quaternion (4-element array).
    :param t: Interpolation parameter (0 <= t <= 1).
    :return: Interpolated quaternion.
    """
    q1 = R.from_quat(q1)
    q2 = R.from_quat(q2)

    # Compute the cosine of the angle between the two quaternions
    dot_product = np.dot(q1.as_quat(), q2.as_quat())
    
    # If the dot product is negative, the quaternions are more than 90 degrees apart,
    # so we invert one of them to take the shorter path
    if dot_product < 0.0:
        q2 = R.from_quat(-q2.as_quat())
        dot_product = -dot_product

    # Clamp dot product to the range [-1, 1] to avoid numerical errors
    dot_product = np.clip(dot_product, -1.0, 1.0)

    # Compute the angle between the quaternions
    theta_0 = np.arccos(dot_product)
    theta = theta_0 * t

    # Compute the second quaternion in the interpolation
    q2 = q2.as_quat()
    q1 = q1.as_quat()
    q2 -= dot_product * q1
    q2 = q2 / np.linalg.norm(q2)

    # Compute the interpolated quaternion
    return R.from_quat(np.cos(theta) * q1 + np.sin(theta) * q2).as_quat()

def plot_quaternions(q1, q2, num_points=10):
    """
    Plot the quaternions using SLERP interpolation.
    
    :param q1: The starting quaternion (4-element array).
    :param q2: The ending quaternion (4-element array).
    :param num_points: Number of interpolation points.
    """
    # Define interpolation parameter t
    ts = np.linspace(0, 1, num_points)
    
    # Compute interpolated quaternions
    quats = np.array([slerp(q1, q2, t) for t in ts])
    
    # Plot quaternions on a sphere
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    
    # Plot initial and final rotation axes
    for rotation in quats:
        # Convert quaternion to rotation matrix
        rotation_matrix = R.from_quat(rotation).as_matrix()
        
        # Define the base unit vectors
        origin = np.zeros(3)
        x_axis = np.array([1, 0, 0])
        y_axis = np.array([0, 1, 0])
        z_axis = np.array([0, 0, 1])
        
        # Apply the rotation
        x_rot = rotation_matrix @ x_axis
        y_rot = rotation_matrix @ y_axis
        z_rot = rotation_matrix @ z_axis
        
        # Plot the transformed coordinate axes
        ax.quiver(*origin, *x_rot, color='r', length=0.5, normalize=True)
        ax.quiver(*origin, *y_rot, color='g', length=0.5, normalize=True)
        ax.quiver(*origin, *z_rot, color='b', length=0.5, normalize=True)
    
    ax.set_xlabel('X')
    ax.set_ylabel('Y')
    ax.set_zlabel('Z')
    ax.set_title('SLERP Interpolation of Quaternions')
    plt.show()

# Example quaternions (start and end)
q1 = [1, 0, 0, 0]  # Identity quaternion
q2 = [0, 0, 0, 1]  # 180-degree rotation around Z axis

# Plot SLERP interpolation
plot_quaternions(q1, q2, num_points=30)





import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# 初期の3D座標
point = np.array([1, 1, 1])

# 平行移動
translation = np.array([2, 3, 4])
T = np.array([
    [1, 0, 0, translation[0]],
    [0, 1, 0, translation[1]],
    [0, 0, 1, translation[2]],
    [0, 0, 0, 1]
])
homogeneous_point = np.append(point, 1)
translated_point = np.dot(T, homogeneous_point)[:3]

# 回転 (X軸回り)
theta = np.radians(45)
R_x = np.array([
    [1, 0, 0],
    [0, np.cos(theta), -np.sin(theta)],
    [0, np.sin(theta), np.cos(theta)]
])
rotated_point = np.dot(R_x, translated_point)

# スケーリング
scaling_factors = np.array([2, 3, 1])
S = np.diag(np.append(scaling_factors, 1))
scaled_point = np.dot(S[:3, :3], rotated_point)

# 透視投影
f = 1.0  # 遠くのクリッピング面
n = 0.1  # 近くのクリッピング面
aspect = 1.0
P = np.array([
    [1/(f*aspect), 0, 0, 0],
    [0, 1/f, 0, 0],
    [0, 0, -(f+n)/(f-n), -2*f*n/(f-n)],
    [0, 0, -1, 0]
])
projected_point = np.dot(P, np.append(scaled_point, 1))
projected_point /= projected_point[3]  # ホモジニアス座標の正規化

# プロット
fig = plt.figure(figsize=(12, 6))

# 3D座標変換のプロット
ax = fig.add_subplot(121, projection='3d')
ax.scatter(*point, color='b', label='Original Point', s=100)
ax.scatter(*translated_point, color='r', label='Translated Point', s=100)
ax.scatter(*rotated_point, color='g', label='Rotated Point', s=100)
ax.scatter(*scaled_point, color='purple', label='Scaled Point', s=100)
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.set_title('3D Coordinates Transformation')
ax.legend()

# 透視投影後の2Dプロット
ax2 = fig.add_subplot(122)
ax2.scatter(projected_point[0], projected_point[1], color='orange', label='Projected Point', s=100)
ax2.set_xlabel('X')
ax2.set_ylabel('Y')
ax2.set_title('2D Perspective Projection')
ax2.legend()

plt.show()

# 結果の出力
print(f'Original point: {point}')
print(f'Translated point: {translated_point}')
print(f'Rotated point: {rotated_point}')
print(f'Scaled point: {scaled_point}')
print(f'Projected point: {projected_point[:2]}')





import numpy as np
import matplotlib.pyplot as plt

class Ball:
    def __init__(self, position, velocity, radius):
        self.position = np.array(position)
        self.velocity = np.array(velocity)
        self.radius = radius

class Wall:
    def __init__(self, point, normal):
        self.point = np.array(point)
        self.normal = np.array(normal) / np.linalg.norm(normal)  # Normalization to ensure unit normal

def reflect(velocity, normal):
    return velocity - 2 * np.dot(velocity, normal) * normal

def detect_collision(ball, wall):
    # Calculate distance from ball center to wall
    distance = np.dot(ball.position - wall.point, wall.normal)
    return abs(distance) <= ball.radius

def update_ball(ball, wall):
    if detect_collision(ball, wall):
        ball.velocity = reflect(ball.velocity, wall.normal)
        # Move ball to a position just outside the wall to prevent sticking
        ball.position = ball.position + (ball.radius - np.dot(ball.position - wall.point, wall.normal)) * wall.normal

def plot(ball, wall):
    fig, ax = plt.subplots()
    ax.set_xlim(-10, 10)
    ax.set_ylim(-10, 10)
    
    # Plot ball
    circle = plt.Circle(ball.position, ball.radius, color='blue', fill=True)
    ax.add_artist(circle)
    
    # Plot wall
    wall_start = wall.point - 10 * wall.normal
    wall_end = wall.point + 10 * wall.normal
    ax.plot([wall_start[0], wall_end[0]], [wall_start[1], wall_end[1]], 'r-')
    
    plt.gca().set_aspect('equal', adjustable='box')
    plt.grid(True)
    plt.show()

# Example usage
ball = Ball(position=[0, 0], velocity=[1, 1], radius=1)
wall = Wall(point=[2, 2], normal=[-1, 1])

# Update ball position and velocity
update_ball(ball, wall)

# Plot the results
plot(ball, wall)

import numpy as np
import matplotlib.pyplot as plt

# Pokémon's level, Individual Values (IVs), and Effort Values (EVs)
level = 50
iv = {'HP': 31, 'Attack': 31, 'Defense': 31, 'SpecialAttack': 31, 'SpecialDefense': 31, 'Speed': 31}
ev = {'HP': 252, 'Attack': 252, 'Defense': 252, 'SpecialAttack': 252, 'SpecialDefense': 252, 'Speed': 252}

# Function to calculate stats
def calculate_stat(base_stat, iv, ev, level, nature_modifier=1.0):
    return ((2 * base_stat + iv + (ev // 4)) * level) // 100 + 5

# Garchomp's base stats
base_stats = {
    'HP': 108,
    'Attack': 130,
    'Defense': 95,
    'SpecialAttack': 80,
    'SpecialDefense': 85,
    'Speed': 102
}

# Calculate stats
stats = {stat: calculate_stat(base_stats[stat], iv[stat], ev[stat], level) for stat in base_stats}

# Plot
fig, ax = plt.subplots()
bars = ax.bar(stats.keys(), stats.values(), color='skyblue')
ax.set_xlabel('Stat')
ax.set_ylabel('Value')
ax.set_title('Garchomp Stats')
plt.xticks(rotation=45)
plt.tight_layout()

# Add value labels
for bar in bars:
    height = bar.get_height()
    ax.text(bar.get_x() + bar.get_width() / 2, height, f'{height}', ha='center', va='bottom')

plt.show()


import numpy as np
import matplotlib.pyplot as plt

# 色違いのポケモンが出る確率と試行回数に基づく少なくとも1回出る確率を計算する関数
def probability_at_least_one(n):
    return 1 - ((n - 1) / n) ** n

# 横軸の分母の回数
n_values = np.arange(1, 10001)  # 1から10000までの分母の回数

# 確率計算
probabilities = [probability_at_least_one(n) for n in n_values]

# 色違い確率の値
n_values_special = [8192, 4096]
probabilities_special = [probability_at_least_one(n) for n in n_values_special]

# プロット
plt.figure(figsize=(12, 8))

# 全ての分母の回数に対する確率
plt.plot(n_values, probabilities, color='blue', label='General Case')

# 特定の分母の回数(1/8192 と 1/4096)に対する確率を点でプロット
plt.scatter(n_values_special, probabilities_special, color='red', zorder=5, label='Special Cases (1/8192 and 1/4096)')

plt.xlabel('n (Number of Trials)')
plt.ylabel('Probability of Getting at Least One Shiny Pokémon')
plt.title('Probability of Getting at Least One Shiny Pokémon vs Number of Trials')
plt.legend()
plt.grid(True)
plt.show()


import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# クォータニオンから回転行列への変換
def quaternion_to_rotation_matrix(q):
    w, x, y, z = q
    return np.array([
        [1 - 2*y**2 - 2*z**2, 2*x*y - 2*z*w, 2*x*z + 2*y*w],
        [2*x*y + 2*z*w, 1 - 2*x**2 - 2*z**2, 2*y*z - 2*x*w],
        [2*x*z - 2*y*w, 2*y*z + 2*x*w, 1 - 2*x**2 - 2*y**2]
    ])

# クォータニオンを正規化する
def normalize_quaternion(q):
    return q / np.linalg.norm(q)

# 回転を適用する
def rotate_vector(v, R):
    return np.dot(R, v)

# 初期位置ベクトル
v = np.array([1, 0, 0])

# クォータニオン (w, x, y, z)
q = np.array([0.7071, 0, 0.7071, 0])
q = normalize_quaternion(q)

# 回転行列を計算する
R = quaternion_to_rotation_matrix(q)

# 回転を適用する
v_rot = rotate_vector(v, R)

# プロット
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

# 元のベクトルをプロット
ax.quiver(0, 0, 0, v[0], v[1], v[2], color='r', label='Original vector')

# 回転後のベクトルをプロット
ax.quiver(0, 0, 0, v_rot[0], v_rot[1], v_rot[2], color='b', label='Rotated vector')

# 軸の設定
ax.set_xlim([-1, 1])
ax.set_ylim([-1, 1])
ax.set_zlim([-1, 1])
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.legend()
ax.set_title('Quaternion Rotation')

plt.show()



import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# クォータニオンから回転行列への変換
def quaternion_to_rotation_matrix(q):
    w, x, y, z = q
    return np.array([
        [1 - 2*y**2 - 2*z**2, 2*x*y - 2*z*w, 2*x*z + 2*y*w],
        [2*x*y + 2*z*w, 1 - 2*x**2 - 2*z**2, 2*y*z - 2*x*w],
        [2*x*z - 2*y*w, 2*y*z + 2*x*w, 1 - 2*x**2 - 2*y**2]
    ])

# クォータニオンを正規化する
def normalize_quaternion(q):
    return q / np.linalg.norm(q)

# クォータニオンの内積
def dot_product(q1, q2):
    return np.dot(q1, q2)

# クォータニオンのSlerp
def slerp(q1, q2, t):
    dot = dot_product(q1, q2)
    if dot < 0.0:
        q1 = -q1
        dot = -dot

    if dot > 0.9995:
        return normalize_quaternion(q1 + t * (q2 - q1))

    theta_0 = np.arccos(dot)
    theta = theta_0 * t

    q3 = normalize_quaternion(q2 - q1 * dot)
    return q1 * np.cos(theta) + q3 * np.sin(theta)

# 初期位置ベクトル
v = np.array([1, 0, 0])

# クォータニオン (w, x, y, z)
q1 = np.array([0.7071, 0, 0.7071, 0])
q2 = np.array([0.7071, 0.7071, 0, 0])
q1 = normalize_quaternion(q1)
q2 = normalize_quaternion(q2)

# プロット
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

# 元のベクトルをプロット
ax.quiver(0, 0, 0, v[0], v[1], v[2], color='r', label='Original vector')

# Slerpで生成されたベクトルをプロット
t_values = np.linspace(0, 1, 100)
for t in t_values:
    qt = slerp(q1, q2, t)
    R = quaternion_to_rotation_matrix(qt)
    v_rot = np.dot(R, v)
    ax.quiver(0, 0, 0, v_rot[0], v_rot[1], v_rot[2], color='b', alpha=0.1)

# 軸の設定
ax.set_xlim([-1, 1])
ax.set_ylim([-1, 1])
ax.set_zlim([-1, 1])
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.legend()
ax.set_title('Spherical Linear Interpolation (Slerp)')

plt.show()



import numpy as np
import matplotlib.pyplot as plt

def point_in_triangle(pt, v1, v2, v3):
    """
    点が三角形の内側にあるかを判定する
    """
    def sign(p1, p2, p3):
        return (p1[0] - p3[0]) * (p2[1] - p3[1]) - (p2[0] - p3[0]) * (p1[1] - p3[1])
    
    d1 = sign(pt, v1, v2)
    d2 = sign(pt, v2, v3)
    d3 = sign(pt, v3, v1)
    
    has_neg = (d1 < 0) or (d2 < 0) or (d3 < 0)
    has_pos = (d1 > 0) or (d2 > 0) or (d3 > 0)
    
    return not (has_neg and has_pos)

# 三角形の頂点
triangle = np.array([[0, 0], [2, 0], [1, 2]])

# テストする点
test_points = np.array([[1, 1], [2, 1], [0.5, 0.5], [1.5, 1.5]])

# 有向線分の起点と方向
line_start = np.array([0.5, 0.5])
line_end = np.array([1.5, 1.5])

# プロット
fig, ax = plt.subplots()

# 三角形のプロット
triangle_path = plt.Polygon(triangle, edgecolor='blue', fill=None, linestyle='-', linewidth=2)
ax.add_patch(triangle_path)

# テストする点のプロット
for pt in test_points:
    if point_in_triangle(pt, *triangle):
        ax.plot(pt[0], pt[1], 'go')  # 内側の点は緑
    else:
        ax.plot(pt[0], pt[1], 'ro')  # 外側の点は赤

# 有向線分のプロット
ax.arrow(line_start[0], line_start[1], line_end[0] - line_start[0], line_end[1] - line_start[1], 
         head_width=0.1, head_length=0.2, fc='k', ec='k')

# 軸の設定
ax.set_xlim([-1, 3])
ax.set_ylim([-1, 3])
ax.set_aspect('equal', 'box')
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_title('Point in Triangle and Directed Line Segment')

plt.grid(True)
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# 環状グラデーションの生成
def ring_gradient(size, inner_radius, outer_radius):
    y, x = np.ogrid[-size[0] // 2:size[0] // 2, -size[1] // 2:size[1] // 2]
    distance = np.sqrt(x*x + y*y)
    
    # 内側と外側の半径に基づいて距離を正規化
    ring = np.clip((distance - inner_radius) / (outer_radius - inner_radius), 0, 1)
    return ring

# 傾きを考慮したグラデーションの生成
def gradient_with_slope(size, slope):
    y, x = np.ogrid[0:size[0], 0:size[1]]
    gradient = x * np.cos(slope) + y * np.sin(slope)
    gradient = gradient - gradient.min()
    gradient = gradient / gradient.max()
    return gradient

# グラデーションの合成
def composite_gradient(ring_grad, linear_grad):
    return ring_grad * linear_grad

# サイズ、内側半径、外側半径、傾きを定義
size = (500, 500)
inner_radius = 50
outer_radius = 200
slope = np.pi / 4  # 45度の傾き

# グラデーションを生成
ring_grad = ring_gradient(size, inner_radius, outer_radius)
linear_grad = gradient_with_slope(size, slope)
combined_grad = composite_gradient(ring_grad, linear_grad)

# プロット
fig, ax = plt.subplots(figsize=(6, 6))
cax = ax.imshow(combined_grad, cmap='viridis', origin='lower')
fig.colorbar(cax)

ax.set_title('Ring Gradient with Linear Gradient and Slope')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

# シミュレーションのパラメータ
mass = 1.0  # 質量
force = np.array([0.0, -9.8])  # 力(重力)
time_step = 0.01  # 時間刻み
total_time = 2.0  # 総シミュレーション時間

# 初期条件
initial_position = np.array([0.0, 10.0])  # 初期位置
initial_velocity = np.array([2.0, 0.0])   # 初期速度

# 時間配列の生成
time_array = np.arange(0, total_time, time_step)

# 位置と速度の配列
position = np.zeros((len(time_array), 2))
velocity = np.zeros((len(time_array), 2))

# 初期条件の設定
position[0] = initial_position
velocity[0] = initial_velocity

# シミュレーションの実行
for i in range(1, len(time_array)):
    acceleration = force / mass  # 加速度の計算
    velocity[i] = velocity[i-1] + acceleration * time_step  # 速度の更新
    position[i] = position[i-1] + velocity[i-1] * time_step  # 位置の更新

# プロット
plt.figure(figsize=(10, 5))
plt.plot(position[:, 0], position[:, 1], label='Trajectory')
plt.xlabel('X Position')
plt.ylabel('Y Position')
plt.title('Projectile Motion')
plt.legend()
plt.grid(True)
plt.show()


import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# 平行移動行列
def translation_matrix(tx, ty, tz):
    return np.array([
        [1, 0, 0, tx],
        [0, 1, 0, ty],
        [0, 0, 1, tz],
        [0, 0, 0, 1]
    ])

# x軸周りの回転行列
def rotation_matrix_x(theta):
    return np.array([
        [1, 0, 0, 0],
        [0, np.cos(theta), -np.sin(theta), 0],
        [0, np.sin(theta), np.cos(theta), 0],
        [0, 0, 0, 1]
    ])

# スケーリング行列
def scaling_matrix(sx, sy, sz):
    return np.array([
        [sx, 0, 0, 0],
        [0, sy, 0, 0],
        [0, 0, sz, 0],
        [0, 0, 0, 1]
    ])

# オブジェクトの頂点
cube_vertices = np.array([
    [-1, -1, -1],
    [1, -1, -1],
    [1, 1, -1],
    [-1, 1, -1],
    [-1, -1, 1],
    [1, -1, 1],
    [1, 1, 1],
    [-1, 1, 1]
])

# 4次元ホモジニアス座標に変換
cube_vertices_hom = np.hstack((cube_vertices, np.ones((cube_vertices.shape[0], 1))))

# トランスフォーメーションのパラメータ
tx, ty, tz = 2, 3, 4
theta = np.pi / 4
sx, sy, sz = 1.5, 1.5, 1.5

# トランスフォーメーション行列の生成
T = translation_matrix(tx, ty, tz)
R = rotation_matrix_x(theta)
S = scaling_matrix(sx, sy, sz)

# 合成トランスフォーメーション行列
TRS = T @ R @ S

# トランスフォーメーションの適用
transformed_vertices_hom = (TRS @ cube_vertices_hom.T).T
transformed_vertices = transformed_vertices_hom[:, :3]

# プロット
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

# 元のキューブのプロット
ax.scatter(cube_vertices[:, 0], cube_vertices[:, 1], cube_vertices[:, 2], color='blue', label='Original Cube')

# 変換後のキューブのプロット
ax.scatter(transformed_vertices[:, 0], transformed_vertices[:, 1], transformed_vertices[:, 2], color='red', label='Transformed Cube')

# ラベルとタイトルの設定
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')
ax.set_title('3D Transformation of a Cube')
ax.legend()

plt.show()



import numpy as np
import matplotlib.pyplot as plt

# 点と三角形の衝突判定
def is_point_in_triangle(P, A, B, C):
    v0 = C - A
    v1 = B - A
    v2 = P - A
    
    dot00 = np.dot(v0, v0)
    dot01 = np.dot(v0, v1)
    dot02 = np.dot(v0, v2)
    dot11 = np.dot(v1, v1)
    dot12 = np.dot(v1, v2)
    
    invDenom = 1 / (dot00 * dot11 - dot01 * dot01)
    u = (dot11 * dot02 - dot01 * dot12) * invDenom
    v = (dot00 * dot12 - dot01 * dot02) * invDenom
    
    return (u >= 0) and (v >= 0) and (u + v < 1)

# 線分と線分の交差判定
def is_intersecting(A, B, C, D):
    def ccw(A, B, C):
        return (C[1] - A[1]) * (B[0] - A[0]) > (B[1] - A[1]) * (C[0] - A[0])
    
    return ccw(A, C, D) != ccw(B, C, D) and ccw(A, B, C) != ccw(A, B, D)

# テストデータ
triangle = np.array([[0, 0], [1, 0], [0, 1]])
point_inside = np.array([0.25, 0.25])
point_outside = np.array([1.5, 1.5])

line1 = np.array([[0, 0], [1, 1]])
line2 = np.array([[0, 1], [1, 0]])
line3 = np.array([[2, 2], [3, 3]])

# 衝突判定の実行
inside = is_point_in_triangle(point_inside, *triangle)
outside = is_point_in_triangle(point_outside, *triangle)
intersecting = is_intersecting(*line1, *line2)
non_intersecting = is_intersecting(*line1, *line3)

# プロット
fig, ax = plt.subplots(1, 2, figsize=(14, 6))

# 点と三角形の衝突判定のプロット
ax[0].plot(*zip(*triangle, triangle[0]), 'bo-', label='Triangle')
ax[0].plot(point_inside[0], point_inside[1], 'go' if inside else 'ro', label='Point Inside' if inside else 'Point Outside')
ax[0].plot(point_outside[0], point_outside[1], 'go' if outside else 'ro', label='Point Inside' if outside else 'Point Outside')
ax[0].set_title('Point-Triangle Collision Detection')
ax[0].legend()
ax[0].set_xlim(-1, 2)
ax[0].set_ylim(-1, 2)
ax[0].grid(True)

# 線分と線分の交差判定のプロット
ax[1].plot(*zip(*line1), 'bo-', label='Line 1')
ax[1].plot(*zip(*line2), 'go-' if intersecting else 'ro-', label='Line 2 (Intersecting)' if intersecting else 'Line 2 (Non-Intersecting)')
ax[1].plot(*zip(*line3), 'ro-', label='Line 3 (Non-Intersecting)')
ax[1].set_title('Line-Line Intersection Detection')
ax[1].legend()
ax[1].set_xlim(-1, 4)
ax[1].set_ylim(-1, 4)
ax[1].grid(True)

plt.show()

0
0
1

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0