请注意:笔记内容片面粗浅,请读者批判着阅读!
一、图像传感与获取的工程实现
1. 三大传感器类型对比
传感器类型 | 工作原理 | 典型应用场景 | 数学建模 |
---|---|---|---|
单传感器 | 机械扫描逐点采集 | 高精度文档扫描仪 | f ( x , y ) = ∑ s ( t ) δ ( x − x ( t ) ) f(x,y) = \sum s(t)δ(x-x(t)) f(x,y)=∑s(t)δ(x−x(t)) |
条带传感器 | 线阵传感器+运动合成 | 卫星遥感、工业流水线 | g ( t ) = ∫ y = 0 H s ( x , y ) d y g(t) = \int_{y=0}^H s(x,y)dy g(t)=∫y=0Hs(x,y)dy |
阵列传感器 | CMOS/CCD平面光电转换 | 手机摄影、监控摄像头 | I i j = Φ ( E i j ⋅ Δ t ) I_{ij} = \Phi(E_{ij} \cdot Δt) Iij=Φ(Eij⋅Δt) |
注:CMOS传感器量子效率曲线(图2.20)显示其在550nm(绿光)处达到峰值65%
2. 图像形成模型
- 光学系统简化模型:
g ( x , y ) = h ( x , y ) ∗ f ( x , y ) + η ( x , y ) g(x,y) = h(x,y) * f(x,y) + \eta(x,y) g(x,y)=h(x,y)∗f(x,y)+η(x,y)
h ( x , y ) h(x,y) h(x,y)为点扩散函数, η \eta η为传感器噪声(服从泊松分布)。
二、图像数字化:采样与量化
1. 采样:从连续到离散
奈奎斯特临界采样实验
import numpy as np
import matplotlib.pyplot as plt# 生成含高频正弦信号
x = np.linspace(0, 4*np.pi, 1000)
signal = np.sin(20*x) + 0.5*np.sin(50*x) # 最高频率50Hz# 临界采样(100Hz)与欠采样(30Hz)
t_critical = np.arange(0, 4*np.pi, 0.02) # 100Hz > 2*50Hz
t_under = np.arange(0, 4*np.pi, 0.033) # 30Hz < 2*50Hzplt.figure(figsize=(12,4))
plt.plot(x, signal, 'b-', label='原始信号')
plt.stem(t_critical, np.interp(t_critical, x, signal), 'g', label='临界采样')
plt.stem(t_under, np.interp(t_under, x, signal), 'r', label='欠采样')
plt.title("奈奎斯特采样定理验证"), plt.legend()
plt.show()
输出结论:欠采样信号出现低频混叠(红色),验证采样定理重要性
另外说一句,如何中文显示方框或乱码的,加入下面语句即可。
from pylab import mpl
# 设置显示中文字体
mpl.rcParams["font.sans-serif"] = ["SimHei"]
# 设置正常显示符号
mpl.rcParams["axes.unicode_minus"] = False
2. 量化:从连续到分级
医疗影像高位深量化必要性(模拟)
import numpy as np
import matplotlib.pyplot as plt
from pylab import mpl
# 设置显示中文字体
mpl.rcParams["font.sans-serif"] = ["SimHei"]
# 设置正常显示符号
mpl.rcParams["axes.unicode_minus"] = False# 生成16-bit合成医学影像(模拟肺部CT)
def generate_synthetic_ct(size=512):x, y = np.meshgrid(np.linspace(-5, 5, size), np.linspace(-5, 5, size))# 基础组织模拟img = 30000 * np.exp(-(x ** 2 + y ** 2)) # 正态分布背景# 添加血管结构(细微灰度变化)vessels = 5000 * np.abs(np.sin(10 * x) * np.cos(8 * y))# 模拟微小病变(直径10像素,Δ=200 HU)lesion = 200 * (np.sqrt((x - 1) ** 2 + (y - 1) ** 2) < 0.1)return (img + vessels + lesion).astype(np.uint16)original = generate_synthetic_ct()
quantized = (original // 256).astype(np.uint8) # 量化到8-bitdef apply_ct_window(img, window_level=30000, window_width=40000):"""将16-bit DICOM映射到8-bit显示范围"""min_val = window_level - window_width // 2max_val = window_level + window_width // 2img_clipped = np.clip(img, min_val, max_val)return ((img_clipped - min_val) / (max_val - min_val) * 255).astype(np.uint8)# 正确显示16-bit图像(使用肺窗:WL=30000, WW=40000)
display_16bit = apply_ct_window(original)
display_8bit = quantizedplt.subplot(121), plt.imshow(display_16bit, cmap='gray'), plt.title('16-bit(病变可见)')
plt.subplot(122), plt.imshow(display_8bit, cmap='gray'), plt.title('8-bit(病变丢失)')
plt.show()
三、插值算法对比与代码实现
插值算法是一种通过已知数据点预测未知数据的方法。在科学计算、图像处理和工程应用中,插值算法被广泛应用于数据拟合和信号处理。
1. 三种经典插值算法数学原理
三种经典插值算法数学原理对比表
算法 | 数学表达式 | 计算复杂度 | 适用场景 |
---|---|---|---|
最近邻插值 | f ( x , y ) = f ( round ( x ) , round ( y ) ) f(x,y) = f(\text{round}(x), \text{round}(y)) f(x,y)=f(round(x),round(y)) | O ( 1 ) O(1) O(1) | 实时系统、像素艺术 |
双线性插值 | f ( x , y ) = a 0 + a 1 x + a 2 y + a 3 x y f(x,y) = a_0 + a_1x + a_2y + a_3xy f(x,y)=a0+a1x+a2y+a3xy (基于4个邻近点加权平均) | O ( n ) O(n) O(n) | 通用图像缩放(默认方案) |
双三次插值 | f ( x , y ) = ∑ i = 0 3 ∑ j = 0 3 a i j x i y j f(x,y) = \sum_{i=0}^3 \sum_{j=0}^3 a_{ij}x^i y^j f(x,y)=i=0∑3j=0∑3aijxiyj (基于16邻近点) | O ( n 2 ) O(n^2) O(n2) | 高质量放大/医学影像重建 |
核心公式说明
-
最近邻插值
- 直接取最邻近像素值,数学上等价于二维脉冲函数采样
- 计算极快,但会产生锯齿(数学不连续)
-
双线性插值
- 在x/y方向分别线性插值,表达式展开为:
f ( x , y ) = ( x 2 − x ) ( y 2 − y ) f 11 + ( x − x 1 ) ( y 2 − y ) f 21 + ( x 2 − x ) ( y − y 1 ) f 12 + ( x − x 1 ) ( y − y 1 ) f 22 ( x 2 − x 1 ) ( y 2 − y 1 ) f(x,y) = \frac{(x_2 - x)(y_2 - y)f_{11} + (x - x_1)(y_2 - y)f_{21} + (x_2 - x)(y - y_1)f_{12} + (x - x_1)(y - y_1)f_{22}}{(x_2 - x_1)(y_2 - y_1)} f(x,y)=(x2−x1)(y2−y1)(x2−x)(y2−y)f11+(x−x1)(y2−y)f21+(x2−x)(y−y1)f12+(x−x1)(y−y1)f22
- 在x/y方向分别线性插值,表达式展开为:
-
双三次插值
- 基于三次多项式逼近,需满足边界导数连续条件
- 常用核函数:
$$
W(t) = \begin{cases}
(a+2)|t|^3 - (a+3)|t|^2 + 1 & \text{for } |t| \leq 1 \
a|t|^3 - 5a|t|^2 + 8a|t| - 4a & \text{for } 1 < |t| < 2 \
0 & \text{otherwise}
\end{cases}
(通常取 (通常取 (通常取 a = -0.5 $$)
复杂度对比实验
import time
import numpy as np
from skimage import transformdef benchmark_interpolation(method, scale=8.0, num_runs=10):# 预生成测试图像(避免重复生成影响计时)img = np.random.rand(512, 512)# 根据方法调整抗锯齿anti_aliasing = False if method == 'nearest' else True# 预热(首次运行不计时)transform.rescale(img, scale, order={'nearest': 0, 'bilinear': 1, 'bicubic': 3}[method],anti_aliasing=anti_aliasing, anti_aliasing_sigma=0.1)# 正式计时total_time = 0.0for _ in range(num_runs):start = time.time()transform.rescale(img, scale,order={'nearest': 0, 'bilinear': 1, 'bicubic': 3}[method],anti_aliasing=anti_aliasing)total_time += time.time() - startreturn total_time / num_runsprint(f"最近邻耗时: {benchmark_interpolation('nearest'):.4f}s")
print(f"双线性耗时: {benchmark_interpolation('bilinear'):.4f}s")
print(f"双三次耗时: {benchmark_interpolation('bicubic'):.4f}s")
输出结果如下:
最近邻耗时: 0.0929s
双线性耗时: 0.2310s
双三次耗时: 0.7054s
2. Python插值效果对比实验
from skimage import data, transformimg = data.camera()
scale_factor = 8.0 # 放大8倍methods = [('nearest', '最近邻插值'),('bilinear', '双线性插值'), ('bicubic', '双三次插值')
]plt.figure(figsize=(15,5))
for i, (method, name) in enumerate(methods):enlarged = transform.rescale(img, scale_factor, order={'nearest':0, 'bilinear':1, 'bicubic':3}[method], anti_aliasing=False)plt.subplot(1,3,i+1)plt.imshow(enlarged[200:300, 200:300], cmap='gray') # 局部放大plt.title(f'{name}(边缘效果)'), plt.axis('off')
plt.tight_layout()
plt.show()
实验结果:
最近邻:明显锯齿
双线性:平滑但有模糊
双三次:保留细节但计算耗时