Serendipity/cDNA-image-processing
1
1
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

feat: 简化版可视化拆分为6张独立图片

Browse Source 
01_grid_overlay.png - 网格线叠加
02_col_projection.png - 列投影曲线
03_row_projection.png - 行投影曲线
04_histogram.png - 直方图+Otsu阈值
05_segmentation_raw.png - 逐格分割(后处理前)
06_post_processed.png - 最终二值图
This commit is contained in:
刘航宇
2026-05-08 10:08:00 +08:00
parent 70790c1d3e
commit 99c9d0263a
8 changed files with 95 additions and 45 deletions
Download Patch File Download Diff File Expand all files Collapse all files
src/cDNA_gridding_simple.py
+95 -45
View File
@@ -109,31 +109,25 @@ def draw_grid_lines(gray: np.ndarray):
流程：
Otsu 求自适应百分比 → 列/行投影 → 减阈值 → 过零点配对 → 空隙中点
返回 (纵线x列表, 横线y列表, T, 自适应百分比)。
返回 (纵线, 横线, T, pct, 列投影, 行投影, 减阈值后的列投影, 减阈值后的行投影)
"""
T = otsu_threshold_pixels(gray) # 像素级最佳阈值
pct = T / 255.0 # 归一化为百分比
H, W = gray.shape
# ---- 1. 横轴投影 ----
# 对每一列上所有行的灰度求和 → 长度=宽度的数组
# 斑点所在列亮像素多 → 和较大（曲线凸起）
# 空隙所在列暗像素多 → 和较小（曲线凹陷）
col_profile = np.sum(gray, axis=0).astype(float)
# ---- 2. 纵轴投影 ----
# 对每一行上所有列的灰度求和 → 长度=高度的数组
row_profile = np.sum(gray, axis=1).astype(float)
# ---- 3. 投影阈值 ----
# 振幅 × 自适应百分比 = 动态阈值（不写死，每张图自调整）
col_T = (np.max(col_profile) - np.min(col_profile)) * pct
row_T = (np.max(row_profile) - np.min(row_profile)) * pct
col_T_val = (np.max(col_profile) - np.min(col_profile)) * pct
row_T_val = (np.max(row_profile) - np.min(row_profile)) * pct
# ---- 4. 曲线减去阈值 ----
# 减完后：正 = 斑点区域，负 = 空隙，0 = 分界线（过零点）
col_shifted = col_profile - col_T
row_shifted = row_profile - row_T
col_shifted = col_profile - col_T_val
row_shifted = row_profile - row_T_val
# ---- 5. 过零点配对 → 空隙中线 ----
def find_gap_lines(prof_shifted: np.ndarray) -> np.ndarray:
@@ -177,7 +171,7 @@ def draw_grid_lines(gray: np.ndarray):
x_lines = find_gap_lines(col_shifted)
y_lines = find_gap_lines(row_shifted)
return x_lines, y_lines, T, pct
return x_lines, y_lines, T, pct, col_profile, row_profile, col_shifted, row_shifted, col_T_val, row_T_val
# ================================================================
@@ -237,61 +231,117 @@ def main():
gray = (color.rgb2gray(img[:, :, :3]) * 255).astype(np.uint8)
# ---- 1. 网格划线 ----
x_lines, y_lines, T_otsu, pct = draw_grid_lines(gray)
(x_lines, y_lines, T_otsu, pct,
col_prof, row_prof, col_shifted, row_shifted,
col_T_val, row_T_val) = draw_grid_lines(gray)
print(f"检测到 {len(x_lines)} 条纵线, {len(y_lines)} 条横线")
print(f"Otsu 阈值: T={T_otsu}, 自适应百分比: {pct*100:.1f}%")
# ---- 2. 逐格分割 + 后处理 ----
bw_full = np.zeros_like(gray) # 初始化全黑底图
for i in range(len(y_lines) - 1): # 遍历每一行格子
for j in range(len(x_lines) - 1): # 遍历每一列格子
r1, r2 = y_lines[i], y_lines[i + 1] # 格子行范围
c1, c2 = x_lines[j], x_lines[j + 1] # 格子列范围
blk = gray[r1:r2, c1:c2] # 提取该格子
bw_full = np.zeros_like(gray)
for i in range(len(y_lines) - 1):
for j in range(len(x_lines) - 1):
r1, r2 = y_lines[i], y_lines[i + 1]
c1, c2 = x_lines[j], x_lines[j + 1]
blk = gray[r1:r2, c1:c2]
if blk.size == 0:
continue
# 对单个格子做 Otsu 分割
T = otsu_threshold_pixels(blk)
bw_blk = (blk > T).astype(np.uint8) # 大于 T 为前景
# 每格只保留最大连通域（去噪点）
bw_blk = (blk > T).astype(np.uint8)
bw_blk = keep_largest_object(bw_blk)
bw_full[r1:r2, c1:c2] = bw_blk # 拼回全局图
bw_full[r1:r2, c1:c2] = bw_blk
# 全局去除小连通域（自动判别噪声阈值）
bw_clean = remove_small_objects(bw_full)
# ---- 3. 统计斑点 ----
labeled, num = ndimage.label(bw_clean)
spot_sizes = [int(np.sum(labeled == i)) for i in range(1, num + 1)]
valid = [s for s in spot_sizes if s >= 10] # 过滤极小噪声
valid = [s for s in spot_sizes if s >= 10]
print(f"检测到 {len(valid)} 个斑点")
# ---- 4. 输出三栏图 ----
fig, axes = plt.subplots(1, 3, figsize=(20, 7))
# ---- 4. 可视化输出（每张图独立保存）----
# 左：灰度图 + 网格线
axes[0].imshow(gray, cmap='gray')
# 图1：网格线叠加原图
fig1, ax1 = plt.subplots(figsize=(8, 8))
ax1.imshow(gray, cmap='gray')
for x in x_lines:
axes[0].axvline(x=x, color='lime', linewidth=0.5)
ax1.axvline(x=x, color='lime', linewidth=0.5)
for y in y_lines:
axes[0].axhline(y=y, color='lime', linewidth=0.5)
axes[0].set_title(f'网格划分 ({len(x_lines)}x{len(y_lines)})', fontsize=13)
axes[0].axis('off')
ax1.axhline(y=y, color='lime', linewidth=0.5)
ax1.set_title(f'网格划分 ({len(x_lines)}x{len(y_lines)})', fontsize=13)
ax1.axis('off')
fig1.savefig(os.path.join(OUTPUT_DIR, '01_grid_overlay.png'), dpi=150, bbox_inches='tight')
plt.close(fig1)
# 中：逐格 Otsu 分割结果（后处理前）
axes[1].imshow(bw_full, cmap='gray')
axes[1].set_title('逐格 Otsu 分割（T=' + str(T_otsu) + '）', fontsize=13)
axes[1].axis('off')
# 图2：列投影曲线（带阈值线和过零点标记）
fig2, ax2 = plt.subplots(figsize=(10, 4))
xs = np.arange(len(col_prof))
ax2.plot(xs, col_prof, 'b-', linewidth=0.6, label='col profile')
ax2.axhline(y=col_T_val, color='orange', linestyle='--', linewidth=1,
label=f'threshold X={col_T_val:.0f}')
ax2.plot(xs, col_shifted, 'g-', linewidth=0.6, alpha=0.5, label='after -X')
ax2.fill_between(xs, 0, col_shifted, where=(col_shifted > 0), color='green', alpha=0.1)
ax2.fill_between(xs, 0, col_shifted, where=(col_shifted < 0), color='red', alpha=0.1)
zero_idx = np.where(np.diff(col_shifted > 0) != 0)[0]
for zi in zero_idx[:50]:
ax2.axvline(x=zi, color='purple', linewidth=0.3, alpha=0.5)
for xl in x_lines:
ax2.axvline(x=xl, color='red', linewidth=0.8, alpha=0.7)
ax2.set_title('col projection', fontsize=12)
ax2.set_xlabel('col')
ax2.legend(fontsize=8)
fig2.savefig(os.path.join(OUTPUT_DIR, '02_col_projection.png'), dpi=120, bbox_inches='tight')
plt.close(fig2)
# 右：后处理之后的最终二值图
axes[2].imshow(bw_clean, cmap='gray')
axes[2].set_title(f'后处理结果 ({len(valid)}个斑点)', fontsize=13)
axes[2].axis('off')
# 图3：行投影曲线
fig3, ax3 = plt.subplots(figsize=(10, 4))
ys = np.arange(len(row_prof))
ax3.plot(row_prof, ys, 'b-', linewidth=0.6, label='row profile')
ax3.axvline(x=row_T_val, color='orange', linestyle='--', linewidth=1,
label=f'threshold X={row_T_val:.0f}')
ax3.plot(row_shifted, ys, 'g-', linewidth=0.6, alpha=0.5, label='after -X')
ax3.fill_betweenx(ys, 0, row_shifted, where=(row_shifted > 0), color='green', alpha=0.1)
ax3.fill_betweenx(ys, 0, row_shifted, where=(row_shifted < 0), color='red', alpha=0.1)
zero_idx_r = np.where(np.diff(row_shifted > 0) != 0)[0]
for zi in zero_idx_r[:50]:
ax3.axhline(y=zi, color='purple', linewidth=0.3, alpha=0.5)
for yl in y_lines:
ax3.axhline(y=yl, color='red', linewidth=0.8, alpha=0.7)
ax3.set_title('row projection', fontsize=12)
ax3.set_ylabel('row')
ax3.legend(fontsize=8)
fig3.savefig(os.path.join(OUTPUT_DIR, '03_row_projection.png'), dpi=120, bbox_inches='tight')
plt.close(fig3)
out_path = os.path.join(OUTPUT_DIR, 'gridding_simple.png')
fig.savefig(out_path, dpi=150, bbox_inches='tight')
plt.close(fig)
print(f"保存: {out_path}")
# 图4：灰度直方图 + Otsu 阈值
fig4, ax4 = plt.subplots(figsize=(8, 5))
ax4.hist(gray.ravel(), bins=50, color='gray', edgecolor='black', linewidth=0.3)
ax4.axvline(x=T_otsu, color='red', linestyle='--', linewidth=2,
label=f'Otsu T={T_otsu} (pct={pct*100:.1f}%)')
ax4.set_title('histogram + Otsu threshold', fontsize=12)
ax4.set_xlabel('gray value')
ax4.set_ylabel('pixel count')
ax4.legend()
fig4.savefig(os.path.join(OUTPUT_DIR, '04_histogram.png'), dpi=120, bbox_inches='tight')
plt.close(fig4)
# 图5：逐格 Otsu 分割（后处理前）
fig5, ax5 = plt.subplots(figsize=(8, 8))
ax5.imshow(bw_full, cmap='gray')
ax5.set_title('per-cell Otsu (before post-processing)', fontsize=13)
ax5.axis('off')
fig5.savefig(os.path.join(OUTPUT_DIR, '05_segmentation_raw.png'), dpi=150, bbox_inches='tight')
plt.close(fig5)
# 图6：后处理结果（最终二值图）
fig6, ax6 = plt.subplots(figsize=(8, 8))
ax6.imshow(bw_clean, cmap='gray')
ax6.set_title(f'post-processed ({len(valid)} spots)', fontsize=13)
ax6.axis('off')
fig6.savefig(os.path.join(OUTPUT_DIR, '06_post_processed.png'), dpi=150, bbox_inches='tight')
plt.close(fig6)
print(f"共保存 6 张图到: {OUTPUT_DIR}")
if __name__ == '__main__':