실습 5 — Otsu 로 자동 임계값 정하기

import cv2
import numpy as np

gray = cv2.imread("ch08_parts.png", cv2.IMREAD_GRAYSCALE)
# THRESH_OTSU: 임계값을 자동으로 고른다. thresh 자리의 0 은 무시된다
ret, otsu = cv2.threshold(gray, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)
print("Otsu 가 고른 임계값:", int(ret))
print("검은 칸(부품) 비율:", round(float((otsu == 0).mean()), 3))

def my_otsu(gray):
    """후보 t 마다 '왼쪽 비율 × 오른쪽 비율 × (두 평균의 차)²' 를 계산해 가장 큰 t 를 돌려준다."""
    # np.bincount: 0~255 의 값마다 칸이 몇 개인지 센다(히스토그램과 같다)
    hist = np.bincount(gray.ravel(), minlength=256)
    values = np.arange(256)
    total = hist.sum()
    best_t, best_score = 0, -1.0
    for t in range(256):
        left, right = hist[:t + 1], hist[t + 1:]          # t 이하, t 보다 큰 무리
        if left.sum() == 0 or right.sum() == 0:           # 한쪽이 비면 나눈 것이 아니다
            continue
        w0, w1 = left.sum() / total, right.sum() / total  # 두 무리의 비율
        m0 = (left * values[:t + 1]).sum() / left.sum()   # 왼쪽 무리의 평균 밝기
        m1 = (right * values[t + 1:]).sum() / right.sum() # 오른쪽 무리의 평균 밝기
        score = w0 * w1 * (m0 - m1) ** 2
        if score > best_score:                            # 동점이면 먼저 찾은(작은) t 를 남긴다
            best_t, best_score = t, score
    return best_t

print("직접 계산한 임계값:", my_otsu(gray))
print("OpenCV 와 같은가:", my_otsu(gray) == int(ret))

# 본문의 픽셀 8개로도 견준다
pixels = np.array([[30, 40, 50, 90, 120, 160, 170, 180]], dtype=np.uint8)
r8, _ = cv2.threshold(pixels, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)
print("픽셀 8개 - 직접:", my_otsu(pixels), "/ OpenCV:", int(r8))

import matplotlib.pyplot as plt

plt.rcParams["font.family"] = "Malgun Gothic"
plt.rcParams["axes.unicode_minus"] = False

hist = cv2.calcHist([gray], [0], None, [32], [0, 256]).ravel()
fig, axes = plt.subplots(1, 2, figsize=(11, 3.4))
axes[0].bar(np.arange(32) * 8, hist, width=8, align="edge")
axes[0].axvline(ret, color="red", linestyle="--", label=f"Otsu 임계값 {int(ret)}")
axes[0].set_xlim(0, 256)
axes[0].set_xlabel("밝기")
axes[0].set_ylabel("픽셀 수")
axes[0].legend(loc="upper center")
axes[1].imshow(otsu, cmap="gray", vmin=0, vmax=255)
axes[1].set_title("Otsu 결과")
axes[1].axis("off")
plt.tight_layout()
plt.show()

# 씨앗 8 로 고정한 노이즈(평균 0, 표준편차 35)를 칸마다 더하고 0~255 로 자른다
rng = np.random.default_rng(8)
noisy = np.clip(gray + rng.normal(0, 35, gray.shape), 0, 255).astype(np.uint8)
# 6차시에서 배운 가우시안 블러: 7×7 네모 안에서 가운데에 무게를 더 주어 섞는다
blur = cv2.GaussianBlur(noisy, (7, 7), 0)

r1, direct = cv2.threshold(noisy, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)
r2, blurred = cv2.threshold(blur, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)
# 노이즈 없는 원본의 Otsu 결과(otsu)와 같은 칸의 비율로 깨끗함을 잰다
for name, r, out in [("바로 Otsu", r1, direct), ("블러 뒤 Otsu", r2, blurred)]:
    print(f"{name}: 고른 값 {int(r)}, 원본 결과와 같은 칸 {round(float((out == otsu).mean()), 3)}")

fig, axes = plt.subplots(1, 3, figsize=(11, 3))
for ax, (title, im) in zip(axes, [("노이즈 낀 그림", noisy), (f"바로 Otsu ({int(r1)})", direct),
                                  (f"블러 뒤 Otsu ({int(r2)})", blurred)]):
    ax.imshow(im, cmap="gray", vmin=0, vmax=255)
    ax.set_title(title)
    ax.axis("off")
plt.tight_layout()
plt.show()