import matplotlib.pyplot as plt
import numpy as np
# %matplotlib widget

radius = 3
theta = np.linspace(0, 2*np.pi, 100)
alpha = np.linspace(0, np.pi, 100)

x = radius*np.outer(np.cos(theta), np.sin(alpha))
y = radius*np.outer(np.sin(theta), np.sin(alpha))
z = radius*np.outer(np.ones(theta.shape), np.cos(alpha))

fig = plt.figure(figsize = (6, 6))
ax = plt.axes(projection = '3d', box_aspect = (1, 1, 1))

ax.plot_surface(x, y, z, color = 'blue')

ax.set_xlim([-3.2, 3.2])
ax.set_ylim([-3.2, 3.2])
ax.set_zlim([-3.2, 3.2])

(-3.2, 3.2)

radius = 3
sides = 8
theta = np.linspace(0, 2*np.pi, 50)
center = np.array([0, 0, 0])
alpha = np.linspace(0, np.pi, sides + 1)

def rotate_circles(x1, y1, z1, rotation_matrix):
    for j in range(y1.shape[0]):
        for k in range(y1.shape[1]):
            vec = np.array([x1[j, k], y1[j, k], z1[j, k]])
            vec_r = rotation_matrix.dot(vec)
            x1[j, k], y1[j, k], z1[j, k] = vec_r[0], vec_r[1], vec_r[2]
    return x1, y1, z1

def plot_circle(center, theta, rotation_matrix):
    radius_m = np.linspace(0, radius, 2)
    x1 = np.outer(radius_m, np.cos(theta))
    y1 = np.outer(radius_m, np.sin(theta))
    z1 = np.zeros(y1.shape)
        
    x1, y1, z1 = rotate_circles(x1, y1, z1, rotation_matrix)
    
        
    ax.plot_surface(x1, y1, z1, color = 'red', alpha = 0.5)

fig = plt.figure(figsize = (6, 6))
ax = plt.axes(projection='3d', box_aspect = (1, 1, 1))

ax.plot(*center, color = 'red', marker = 'o')
ax.text(center[0] + 0.25, center[1] + 0.25, center[2] + 0.25, s = 'O', fontsize = 12)

for i in range(len(alpha) - 1):
    rotation_matrix_x = np.array([
        [1, 0, 0],
        [0, np.cos(alpha[i]), -np.sin(alpha[i])],
        [0, np.sin(alpha[i]), np.cos(alpha[i])],
    ])
    plot_circle(center, theta, rotation_matrix_x)
    
# ax.plot_surface(Xc, Yc, Zc, color = 'blue', alpha = 0.5)
    
ax.set_xlim(-5, 5)
ax.set_ylim(-5, 5)
ax.set_zlim(-5, 5)

(-5.0, 5.0)

radius = 3
shapes = 15
zs = np.linspace(-radius, radius, shapes + 2)
theta = np.linspace(0, 2*np.pi, 50)
center = np.array([0, 0, 0])

fig = plt.figure(figsize = (6, 6))
ax = plt.axes(projection='3d', box_aspect = (1, 1, 1))

ax.plot(*center, color = 'red', marker = 'o')
ax.text(center[0] + 0.25, center[1] + 0.25, center[2] + 0.25, s = 'O', fontsize = 12)

for z in zs:
    radius_m = np.linspace(0, np.sqrt(radius**2 - z**2), 2)
    x1 = np.outer(radius_m, np.cos(theta))
    y1 = np.outer(radius_m, np.sin(theta))
    z1 = np.zeros(y1.shape) + z

    ax.plot_surface(x1, y1, z1, color = 'red', alpha = 0.5)
    
ax.set_xlim(-5, 5)
ax.set_ylim(-5, 5)
ax.set_zlim(-5, 5)

(-5.0, 5.0)

fig = plt.figure(figsize = (6, 6))
ax = plt.axes(projection='3d', box_aspect = (1, 1, 1))
ax.plot(*center, color = 'red', marker = 'o')
ax.text(center[0] + 0.25, center[1] + 0.25, center[2] + 0.25, s = 'O', fontsize = 12)
for i, z in enumerate(zs):
    radius_m = np.linspace(0, np.sqrt(radius**2 - z**2), 2)
    x1 = np.outer(radius_m, np.cos(theta))
    y1 = np.outer(radius_m, np.sin(theta))
    z1 = np.zeros(y1.shape) + z
    if i == 10:
        color = 'blue'
        alpha = 1
    else:
        color = 'red'
        alpha = 0.2
    ax.plot_surface(x1, y1, z1, color = color, alpha = alpha)
ax.set_xlim(-5, 5)
ax.set_ylim(-5, 5)
ax.set_zlim(-5, 5)

(-5.0, 5.0)