文章目录
【基础】傅里叶变换
傅里叶变换（Fourier Transform，FT）：是一种数学工具，用于将一个信号（图像）分解为多个正弦波（不同的频率成分），以帮助分析图像的频率特征（强度和位置）。
核心思想：任何信号（如图像）都可以表示为若干不同频率的正弦波（sine waves）叠加。每个正弦波有不同的频率、幅度和相位，它们共同组成了原始图像。傅里叶变换通过数学手段计算出这些频率成分，从而将图像从空间域转换到频域。
- 快速傅氏变换（Fast Fourier Transform，FFT）
- 离散傅里叶变换（Discrete Fourier transform，DFT）
0.1、傅里叶变换的复数结果（空间域 -> 频率域）
空间域（Spatial domain）：通常我们在空间域中处理图像，即通过直接操作图像的像素值来实现图像处理。空间域关注的是图像的“空间”特征，如颜色、亮度、边缘等。频率域（Frequency Domain）：频域关注的是图像的频率成分，即图像中不同频率（或周期）变化的强度。频率成分可以理解为图像的不同细节层次：低频部分代表平滑的区域（如背景），而高频部分代表细节（如边缘、纹理）。
傅里叶变换将图像从空间域转换到频率域，生成的结果是一个复数数组，由实部和虚部组成。
每个复数表示频域中某个频率成分的幅度和相位信息，具体而言：
- 实部：表示该频率分量的余弦成分。
- 虚部：表示该频率分量的正弦成分。
傅里叶变换的复数形式，可以通过以下公式表示： F ( u ， v ) = R （ u ， v ） + j ∗ I （ u ， v ） { F(u，v) = R（u，v） + j * I（u，v） } F(u，v)=R（u，v）+j∗I（u，v）
其中：
- F ( u , v ) { F(u, v) } F(u,v)：是频率点 ( u , v ) { (u, v) } (u,v) 的实部（对应图像的余弦分量）。
- I ( u , v ) { I(u, v) } I(u,v)：是频率点 ( u , v ) { (u, v) } (u,v) 的虚部（对应图像的正弦分量）。
- j { j } j：是虚数单位（即 j 2 = − 1 { j^2 = -1 } j2=−1）
0.2、幅度与相位（幅度图像 + 相位图像）
-
幅度（Magnitude）：表示的是每个频率分量的强度信息，即该频率成分在图像中所占的比重。- 计算公式： M a g n i t u d e = s q r t ( R ( u , v ) 2 + I ( u , v ) 2 ) { Magnitude = sqrt(R(u, v)^2 + I(u, v)^2) } Magnitude=sqrt(R(u,v)2+I(u,v)2)
- 幅度图像（Magnitude Image）：反映了图像的不同频率成分的强度。通过操作幅度图像，可以对图像的低频或高频部分进行滤波（如低通滤波、高通滤波），从而实现去噪、平滑、锐化等效果。
-
相位（Phase）：表示的是每个频率分量的位相信息，即该频率成分的 " 开始位置 "。- 计算公式： Phase = arctan ( I ( u , v ) R ( u , v ) ) = tan − 1 ( y x ) \text{Phase} = \arctan\left(\frac{I(u, v)}{R(u, v)}\right) = \tan^{-1}\left(\frac{y}{x}\right) Phase=arctan(R(u,v)I(u,v))=tan−1(xy)
- 相位图像（Phase Image）：描述了不同频率成分在图像中的相对位置和空间分布。改变图像的相位会直接影响图像的结构和形状。例如，即使幅度相同，改变相位也会导致图像的形状发生显著变化。
- 位相信息：频谱图像的中心区域为低频部分，边缘区域为高频部分。
（1）在傅里叶变换后的频谱图像中，通常我们将实部和虚部进行分离显示，或者通过幅度图像展示图像的频率成分。（2）如果要重建原始图像，必须使用幅度和相位信息。
0.3、幅度图像与频谱图像的区别
幅度图像：是一个经过傅里叶变换后的包含频率信息的数值数组（幅度部分），表示各个频率分量的强度（没有空间信息）。频谱图像：是一个经过傅里叶变换后的图像（幅度部分 + 相位部分），展示了各个频率分量的强度和位置分布。通常以灰度图或伪彩色图的方式呈现，能够直观地展示频率成分。
0.4、什么是频谱图像？？？
频谱图像（Spectrogram 或 Frequency Spectrum）：通过傅里叶变换，将图像从空间域转换到频域后得到的图像。是幅度图像的一种可视化形式，通常显示每个频率分量在频域中的位置（相位）和强度（幅度）。
- （1）在显示频谱图像时，通常会将低频部分移动到图像中心（通过 fftshift 操作），这样可以更直观地观察频率成分的分布。
低频部分：位于频谱图像的中心，代表图像中的大范围平缓过渡的区域，如平滑区域，图像的背景、大致轮廓和结构等。
- 低频成分反映了图像的基本外观和主要构图。例如，图像中的背景、平坦区域和较大的结构（如建筑物、天空）主要由低频成分决定。
- 低频成分是图像的大致 " 框架 " ，决定了图像的基本外观和主要构图。
高频部分：位于频谱图像的边缘，代表图像中的小尺度的结构或急剧变化的区域，这些通常是图像中的细节，如图像的边缘、纹理、噪声或锐利的结构变化等。
- 高频成分越多，图像的细节越丰富，表现得越清晰和锐利。
- （2）在显示频谱图像时，通常使用对数尺度，以便增强对比度，特别是当频谱中的低频成分非常强时，高频部分的信号往往显得不明显。
0.4.1、应用领域
- 去噪：
高频噪声：这些噪声通常出现在频谱图像的边缘部分。通过低通滤波（去除高频成分），可以有效去除图像中的细节噪声，保持图像的平滑过渡和背景。周期性噪声（如条纹、网格等）：这种噪声常常呈现出明显的周期性特征，通过分析频谱图像中的频率分量，可以发现噪声的周期性，进而对其进行抑制。- 滤波处理：通过过滤某些频率成分，可以增强或抑制图像的特定特征。
低通滤波：去除图像中的高频成分，适用于图像的平滑、去噪和模糊处理。低频成分通常对应图像的背景和大体结构，去除高频可以使图像更加平滑。高通滤波：去除图像中的低频成分，通常用于增强图像的细节、锐化图像。高频成分通常对应图像中的细节部分，如边缘、纹理等，保留这些部分可以提高图像的清晰度。- 特征分析：在频谱图像中，某些特定的频率分量或模式可能对应图像中的周期性结构（如纹理、图案等）。通过分析频谱图像中的这些模式，可以识别图像中的重复性或周期性特征，帮助进行纹理分析、模式识别或表面检测等任务。
- 图像增强：
- 提高图像的清晰度：通过增强频域中的高频成分，可以让图像更清晰，突出图像中的细节。
- 改善对比度：通过调节频域中的低频部分，可以增强图像的整体对比度，使得图像更加生动和鲜明。
0.4.2、频谱图像显示：水平+竖直
0.4.3、频谱图像处理：水平+竖直（实测）
一、基于FFT的去条纹算法 —— 用于去除图像中的周期性条纹噪声
1.1、基本介绍（Destripe）
基于FFT（快速傅里叶变换）的去条纹算法（Destripe Algorithm）：是一种频域滤波方法。通过分析图像的频域特征，去除图像中的周期性噪声（如条纹、网格等），从而增强图像的清晰度。
核心思想：通过傅里叶变换将图像从空间域转换到频域，在频域中定位并去除条纹噪声对应的频率分量，然后通过逆傅里叶变换将图像还原到空间域，从而得到去条纹后的图像。
" Destripe " 是去除条纹噪声的常用术语，由 "de - " （表示去除、反向）与 " stripe "（条纹）组合而成的词。
1.2、条纹噪声
条纹噪声：是图像处理中常见的一种噪声类型，通常具有明显的周期性和方向性，在图像中呈现出规则的条状或网格状结构，且这些条纹具有较为均匀的间隔和方向。
1.2.1、噪声类型（竖条纹 + 横条纹）
条纹噪声具有固定的频率、方向和幅度。在频域中表现为明显的峰值或亮线。
根据其方向性分为不同的类型，以下是它们的详细介绍：
竖条纹（Vertical Stripe）：是指沿图像垂直方向（即列方向）出现的亮暗相间的线条状图案。
- 在频域中，竖条纹对应的频率分量表现为水平方向的亮线或峰值。
- 可以通过设计水平方向的带阻滤波器来去除。
横条纹（Horizontal Stripe）：是指沿图像水平方向（即行方向）出现的亮暗相间的线条状图案。
- 在频域中，横条纹对应的频率分量表现为垂直方向的亮线或峰值。
- 可以通过设计垂直方向的带阻滤波器来去除。
详细图解请看：【什么是频谱图像？？？】
1.2.2、产生原因
条纹噪声的产生主要有以下几个原因：
- （1）扫描设备问题：在显微镜扫描过程中，扫描线的周期性不稳定或光源的不稳定性可能导致图像中出现条纹噪声。特别是在扫描频率不匹配或光源不均匀时，条纹噪声往往更加明显。
- （2）设备特性：显微镜的硬件（如光源频率、滤波器、镜头等）限制，尤其在共聚焦显微镜或超分辨显微镜中，可能引起不同方向的周期性条纹噪声。
- （3）传感器特性：成像传感器（如CMOS或CCD）的工作原理和电源稳定性可能引起周期性条纹，尤其在高分辨率图像中，这些条纹与细节相互叠加，难以分辨。
- （4）荧光标记和染料的不均匀性：荧光成像中，标记或染料的浓度不均匀、光源照射不均匀等因素也可能引入周期性纹理，表现为条纹噪声。
- （5）物理或光学干涉：光源的不均匀性或光学系统中的干涉效应常常导致图像中出现规则的条纹，尤其在低对比度或背景区域，这种噪声尤为明显。
1.2.3、影响
条纹噪声会显著影响图像的清晰度和准确性，具体表现在以下方面：
- （1）细节丢失：条纹噪声可能会掩盖图像中的真实细节，尤其是在处理细胞、组织切片或类器官等复杂结构时，噪声可能与细胞结构、分子标记等重要特征混淆。
- （2）误导性结果：在进行定量分析或数据提取时，条纹噪声可能导致错误的结果或偏差。例如，在对细胞形态、组织分布等进行自动分析时，噪声可能影响算法的准确性。
- （3）影响后续分析：图像的噪声如果不被有效去除，后续的图像分析（如分割、配准、特征提取等）可能会变得复杂和不可靠。
1.2.4、需求
随着荧光成像技术的进步，尤其在生物医学研究、临床诊断以及高分辨率显微镜成像中，对高质量图像的需求越来越高。因此，去除条纹噪声成为了图像预处理中的重要任务。通过有效的去条纹算法，可以大大提高图像的质量，使得后续的分析和实验结果更加准确和可靠。
1.3、处理步骤
以下是去除条纹的处理步骤，主要利用傅里叶变换和频域滤波方法进行处理。
## （1）离散傅里叶变换（DFT），将空间域图像转换为频域图像。 ———— cv2.DFT_COMPLEX_OUTPUT：指定输出为复数形式，包含实部和虚部
image_dft = cv2.dft(np.float32(image), flags=cv2.DFT_COMPLEX_OUTPUT)
## （2）将频谱的低频部分移到频谱图像的中心
dst_center = np.fft.fftshift(image_dft)
"""（3）频谱图像处理...(滤波、去噪等)"""
## （4）将频谱的低频部分移回原来位置
dst_center1 = np.fft.ifftshift(dst_center)
## （5）逆离散傅里叶变换（IDFT），将频域图像恢复为空间域图像。
# 输出是一个复数数组，其数值范围通常比原始输入图像大得多，因为傅里叶变换和反变换会引入缩放因子。
image_idft = cv2.idft(dst_center1, flags=cv2.DFT_SCALE)
## （6）计算复数数组的模值 ———— magnitude(x, y) = sqrt(x^2 + y^2)
image_magnitude = cv2.magnitude(image_idft[:, :, 0], image_idft[:, :, 1])
（1）读取图像
使用OpenCV或PIL等工具读取图像。并转换为灰度图像。如果图像是彩色图像，需要先转换为灰度图像。
import cv2
image_path = "path_to_your_image.png" # 输入图像路径
image = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE) # 读取灰度图像
（2）前处理（可选） —— 决定了图像质量
对图像进行预处理操作，例如归一化、对比度调整或边缘填充。边缘填充可以减少FFT的边界效应。
# 归一化图像[0, 255]
image = cv2.normalize(image, None, 0, 255, cv2.NORM_MINMAX)
# 边缘填充（可选）
rows, cols = image.shape
nrows = cv2.getOptimalDFTSize(rows)
ncols = cv2.getOptimalDFTSize(cols)
padded_image = cv2.copyMakeBorder(image, 0, nrows - rows, 0, ncols - cols, cv2.BORDER_CONSTANT, value=0)
（3）傅里叶变换：cv2.dft() + np.fft.fftshift()
傅里叶变换将图像从空间域转换到频域，频域中的低频部分代表图像的整体结构，而高频部分代表细节和噪声。
image_dft = cv2.dft(np.float32(image), flags=cv2.DFT_COMPLEX_OUTPUT) # 计算傅里叶变换
dst_center = np.fft.fftshift(image_dft) # 将零频率分量移到中心
# 计算幅度谱（用于可视化）
magnitude_spectrum = 20 * np.log(cv2.magnitude(dft_shifted[:, :, 0], dft_shifted[:, :, 1]))
（4）频域滤波（过滤频率分量） —— mask决定了图像质量（核心）
在频域中，条纹通常表现为周期性的高频分量。通过设计滤波器（如滤波器、掩模），剔除特定频率分量，进而去除图像中的条纹噪声。
- 滤波器
低通滤波：去除高频成分，用于平滑图像，去除细节噪声。高通滤波：去除低频成分，增强细节和边缘。带通滤波：保留一定频率范围的成分，去除其他无关频率成分。带阻滤波：特定频率范围的成分可以被抑制，用于去除周期性噪声（如条纹、网格等）。- 掩模：通常为一个长矩形（常用）、半圆形、锥形等区域，通过保留低频部分并过滤条纹相关的高频部分。
# 获取图像的中心点
rows, cols = image.shape
crow, ccol = rows // 2, cols // 2
# 创建掩模，设置频谱中心为0（即去除中心区域的频率分量）
mask = np.ones((rows, cols, 2), np.uint8) # 初始化掩模为1
mask[crow-30:crow+30, ccol-30:ccol+30] = 0 # 设置中心区域为0，表示保留低频，去除高频
# 应用掩模
dst_center_filtered = dst_center * mask
（5）逆傅里叶变换：np.fft.ifftshift() + cv2.idft()
通过逆傅里叶变换将频域图像转换回空间域，恢复去除条纹噪声后的图像。使用 np.fft.ifftshift() 将低频部分移回原始位置，然后通过 cv2.idft() 进行逆变换。
dst_center1 = np.fft.ifftshift(dst_center_filtered) # 将零频率分量移回原位置。
image_idft = cv2.idft(dst_center1, flags=cv2.DFT_SCALE) # 将图像还原到空间域。
image_idft = cv2.magnitude(image_idft[:, :, 0], image_idft[:, :, 1]) # 计算模值
（6）后处理（可选） —— 决定了图像质量
对去条纹后的图像进行后处理操作，例如归一化、对比度增强或锐化。
# 归一化图像[0, 255]
image_idft = cv2.normalize(image_idft, None, 0, 255, cv2.NORM_MINMAX)
# 数据类型转换
image_idft = image_idft.astype(image.dtype) # 方法一
# image_idft = np.uint8(image_idft) # 方法二
二、项目实战
2.0、使用说明
实测结论（效果最优）：
（1）在自然图像背景下，需要同时去除频谱图像的横向与竖向条纹（去除中心点）；
（2）在荧光图像背景下，只需要去除频谱图像的横向条纹（去除中心点）；
2.1、在自然图像背景下（2D图像）
方法一（代码复现） —— 基于FFT的去条纹算法（性能最优）
在python中使用opencv进行dft和idft去除图像条纹
import os
import numpy as np
import cv2
import tifffile
import matplotlib.pyplot as plt
def striping_algorithm(image, line_width=20, center_radius=64, is_vertical_stripe=True, is_horizontal_stripe=True, is_multi=True):
"""对输入的灰度图像进行离散傅里叶变换（DFT） ———— 通过掩模过滤特定频率分量，去除图像中的条纹噪声。
参数:
image: 输入的灰度图像，NumPy数组。
line_width: 线宽，用于控制中心区域的半径。
center_radius: 中心区域的半径，用于控制保留的频率范围。
is_horizontal_stripe: 是否去除横向条纹（默认True）。
is_vertical_stripe: 是否去除竖向条纹（默认True）。
is_multi: 是否去除频谱中的N条水平线
"""
# （1）离散傅里叶变换：将空间域图像转换到频域，并将低频分量移动到频谱中心
image_dft = cv2.dft(np.float32(image), flags=cv2.DFT_COMPLEX_OUTPUT) # 进行DFT
dst_center = np.fft.fftshift(image_dft) # 将频谱的低频部分移动到中心
# （2）创建掩模：用于过滤特定频率分量（1表示保留，0表示过滤）
rows, cols = image.shape # 获取图像尺寸
mask = np.ones(shape=(rows, cols, 2), dtype=np.uint8) # 初始化全0掩模（保留所有频率）
row_center, col_center = int(rows / 2), int(cols / 2) # 计算图像中心坐标
if is_vertical_stripe:
"""去除竖向条纹，对应频谱中的中心水平线条 ———— 将中心水平线频率分量置0（过滤低频），保留中心区域的高频分量"""
mask[int(row_center - int(line_width // 2)): int(row_center + int(line_width // 2) + 1), :] = 0 # 清空中心水平线
mask[int(row_center - int(line_width // 2)): int(row_center + int(line_width // 2) + 1),
(col_center - center_radius): (col_center + center_radius)] = 1 # 保留中心区域
if is_horizontal_stripe:
"""去除横向条纹，对应频谱中的中心垂直线条 ———— 将中心垂直线频率分量置0（过滤低频），保留中心区域的高频分量"""
mask[:, int(col_center - int(line_width // 2)): int(col_center + int(line_width // 2) + 1)] = 0 # 清空中心垂直线
mask[(row_center - center_radius): (row_center + center_radius),
int(col_center - int(line_width // 2)): int(col_center + int(line_width // 2) + 1)] = 1 # 保留中心区域
if is_multi:
"""去除竖向条纹，如果频谱有多条水平线 ———— 需要手动指定对称线条中位于上侧的中心位置"""
row_value1 = 55 # 325/55
row_value2 = rows - row_value1
mask[int(row_value1 - int(line_width // 2)): int(row_value1 + int(line_width // 2) + 1), :] = 0 # 清空水平线（上下对称）
mask[int(row_value2 - int(line_width // 2)): int(row_value2 + int(line_width // 2) + 1), :] = 0 # 清空水平线（上下对称）
# （3）应用掩模：对频域图像进行滤波
dst_center_mask = dst_center * mask # 将掩模应用到频域图像
# （4）可视化频域图像
amplitude_image = 20 * np.log(cv2.magnitude(dst_center_mask[:, :, 0],
dst_center_mask[:, :, 1])) # 计算幅度谱并进行对数压缩，以便可视化
# （5）反傅里叶变换：将频域图像转换回空间域
dst_origin = np.fft.ifftshift(dst_center_mask) # 将频率重新移回原始位置
image_idft = cv2.idft(dst_origin, flags=cv2.DFT_SCALE) # 进行反傅里叶变换
image_idft = cv2.magnitude(image_idft[:, :, 0], image_idft[:, :, 1]) # 获取幅度，得到重建的空间域图像
image_idft = image_idft.astype(image.dtype)
print(f"输出图像的最大最小像素值：{[np.max(image_idft), np.min(image_idft)]}")
# （6）保存结果
is_save = True
if is_save:
# tifffile.imwrite(os.path.join(os.path.dirname(image_path),
# f"amplitude_image{line_width}_{center_radius}.tif"), amplitude_image)
tifffile.imwrite(os.path.join(os.path.dirname(image_path),
f"result{line_width}_{center_radius}.tif"), image_idft)
# （7）显示结果：原图、掩模、频域图像和重建图像
plt.subplot(131), plt.imshow(image, cmap='gray'), plt.axis('off'), plt.title('Original Image')
# plt.subplot(132), plt.imshow(mask[:, :, 0], cmap='gray'), plt.axis('off'), plt.title('mask')
plt.subplot(132), plt.imshow(amplitude_image, cmap='gray'), plt.axis('off'), plt.title('DFT Image (Frequency Domain)')
plt.subplot(133), plt.imshow(image_idft, cmap='gray'), plt.axis('off'), plt.title('IDFT Image (Reconstructed)')
plt.show()
return image_idft
"""
注意：只支持8bit图像（RGB）
备注：16位图像的像素值范围是0-65535，而8位图像的像素值范围是0-255
"""
if __name__ == '__main__':
# （1）读取图像并转换为灰度图像
# image_path = r"F:\AI\test.tif"
# image = tifffile.imread(image_path)
image_path = r"F:\py\test2.png"
image = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE)
if image is None:
raise FileNotFoundError(f"无法读取图像文件: {image_path}")
# （2）调用 striping_algorithm 函数
value = 1
if value == 1: # 去除竖向 + 横向条纹
image_idft = striping_algorithm(image,
line_width=1,
center_radius=1,
is_vertical_stripe=True,
is_horizontal_stripe=True,
is_multi=False)
elif value == 2: # 去除竖向条纹
image_idft = striping_algorithm(image,
line_width=5,
center_radius=64,
is_vertical_stripe=True,
is_horizontal_stripe=False,
is_multi=False)
elif value == 3: # 去除横向条纹
image_idft = striping_algorithm(image,
line_width=5,
center_radius=64,
is_vertical_stripe=False,
is_horizontal_stripe=True,
is_multi=False)
方法二（代码复现） —— 基于带通滤波器的去条纹算法 —— （比较模糊）
- 结果显示
- 实现方式
方法三（代码复现） —— 基于引导滤波器的去条纹算法 —— （效果减半）
import numpy as np
import cv2
def Normalize_img(img):
"""该函数用于将图像归一化至[0, 1]的范围。
参数:
img: 输入的图像，可以是任何维度的NumPy数组。
返回:
normalized_image: 归一化后的图像，像素值在[0, 1]范围内。
"""
# 计算图像的最小值和最大值
min_n = np.min(img) # 最小值
max_n = np.max(img) # 最大值
normalized_image = (img - min_n) / (max_n - min_n) # 归一化公式[0, 1]
return normalized_image
def guideFilter(I, p, winSize, eps):
"""引导滤波算法，用于图像平滑 ———— 用于将图像I的局部统计信息应用到图像p上。
参数:
I: 引导图像（通常为灰度图像或者边缘信息图像）
p: 输入图像，通常是待平滑或增强的图像。
winSize: 滤波窗口大小，tuple形式，例如(9, 9)。
eps: 正则化常数，防止分母为0，控制平滑程度。
返回:
q: 平滑处理后的图像。
"""
# 对图像I和p进行均值平滑处理
mean_I = cv2.blur(I, winSize) # I的均值平滑
mean_p = cv2.blur(p, winSize) # p的均值平滑
mean_II = cv2.blur(I * I, winSize) # I*I的均值平滑
mean_Ip = cv2.blur(I * p, winSize) # I*p的均值平滑
# 计算方差和协方差
var_I = mean_II - mean_I * mean_I # 方差
cov_Ip = mean_Ip - mean_I * mean_p # 协方差
# 计算相关因子a和b
a = cov_Ip / (var_I + eps) # 相关因子a
b = mean_p - a * mean_I # 相关因子b
# 对a和b进行均值平滑
mean_a = cv2.blur(a, winSize) # 对a进行均值平滑
mean_b = cv2.blur(b, winSize) # 对b进行均值平滑
# 根据平滑后的a和b生成最终结果q
q = mean_a * I + mean_b # 最终输出
return q
def main(image, is_vertical_stripe=True, is_horizontal_stripe=True):
"""主函数，执行行引导滤波、列引导滤波或两者的操作。
参数:
image: 输入的灰度图像，NumPy数组。
is_row: 布尔值，是否执行行引导滤波。默认为True。
is_column: 布尔值，是否执行列引导滤波。默认为True。
"""
I = image / 255.0 # 将图像归一化，像素值范围为[0, 1]
p = I # 引导图像p与I相同
"""（1）行引导滤波处理"""
if is_horizontal_stripe:
eps = 20 # 引导滤波的平滑常数
winSize = (9, 1) # 行向上的窗口大小，9表示窗口大小为9，1表示沿着行方向平滑
guideFilterImg_row = guideFilter(I, p, winSize, eps) # 进行引导滤波处理
guideFilterImg_row = p - guideFilterImg_row # 计算图像差异
else:
guideFilterImg_row = p # 如果不需要行滤波，直接使用原图像
"""（2）列引导滤波处理"""
if is_vertical_stripe:
eps = 20 # 引导滤波的平滑常数
winSize = (1, 15) # 列向上的窗口大小，15表示窗口大小为15，1表示沿着列方向平滑
guideFilter_img_lie = guideFilter(p, guideFilterImg_row, winSize, eps) # 进行引导滤波处理
guideFilterImg_col = p - guideFilter_img_lie # 计算图像差异
guideFilter_img_lie = Normalize_img(guideFilter_img_lie) # 归一化处理
else:
guideFilterImg_col = guideFilterImg_row # 如果不需要列滤波，直接使用行滤波后的结果
"""（3）还原为8bit图像 ———— 将归一化后的列引导滤波图像恢复到[0, 255]范围，并转换为8bit"""
guideFilterImg = Normalize_img(guideFilterImg_col) * 255 # 归一化后乘以255
guideFilterImg[guideFilterImg > 255] = 255 # 确保像素值不超过255
guideFilterImg = np.round(guideFilterImg) # 四舍五入
guideFilterImg = guideFilterImg.astype(np.uint8) # 转换为8位无符号整数图像
"""（4）显示图像 ———— 显示原图、行方向引导滤波图像、列方向引导滤波图像及最终图像"""
import matplotlib.pyplot as plt
plt.subplot(221), plt.imshow(image, cmap='gray'), plt.axis('off'), plt.title('Original Image') # 显示原图
if is_horizontal_stripe:
plt.subplot(222), plt.imshow(guideFilterImg_row, cmap='gray'), plt.axis('off'), plt.title('Row Guide Filter Image') # 显示行引导滤波结果
if is_vertical_stripe:
plt.subplot(223), plt.imshow(guideFilter_img_lie, cmap='gray'), plt.axis('off'), plt.title('Column Guide Filter Image') # 显示列引导滤波结果
plt.subplot(224), plt.imshow(guideFilterImg, cmap='gray'), plt.axis('off'), plt.title('Guide Filter Image') # 显示最终结果
plt.show()
if __name__ == '__main__':
image_path = r"F:\py\test.png" # 图像路径
image = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE) # 读取图像并转换为灰度图
if image is None:
raise FileNotFoundError(f"无法读取图像文件: {image_path}") # 如果图像读取失败，抛出异常
###################################################################################
# 调用主函数，控制是否执行行和列引导滤波
main(image, is_vertical_stripe=True, is_horizontal_stripe=True) # 行+列引导滤波
# main(image, is_vertical_stripe=True, is_horizontal_stripe=False) # 行引导滤波
# main(image, is_vertical_stripe=False, is_horizontal_stripe=True) # 列引导滤波
# main(image, is_vertical_stripe=False, is_horizontal_stripe=False) # 不进行任何滤波
2.2、在荧光图像背景下（类器官、脑图像等）
2.2.2、【 Fiji 复现】 —— 基于FFT的去条纹算法（性能最优）
- （1）打开图像
- （2）高斯滤波
- （3）傅里叶变换（FFT）
- （4）频率滤波（使用Easing ROI自动标定尺寸）
- （5）反傅里叶变换（Invert FFT）
2.2.3、【代码复现】 —— 基于FFT的去条纹算法（性能最优）
测试结果请看：【参数测试】
import os
import numpy as np
import cv2
import tifffile
import matplotlib.pyplot as plt
def fft_filter(img, line_width=1, center_radius=0,
is_vertical_stripe=True, is_horizontal_stripe=False, is_multi=False, is_visualization=True):
"""对输入的灰度图像进行离散傅里叶变换（DFT） ———— 通过掩模过滤特定频率分量，去除图像中的条纹噪声。
参数:
img: 输入的灰度图像，NumPy数组。
line_width: 线宽，用于控制中心区域的半径。
center_radius: 中心区域的半径，用于控制保留的频率范围。
is_horizontal_stripe: 是否去除横向条纹（默认True）。
is_vertical_stripe: 是否去除竖向条纹（默认True）。
is_multi: 是否去除频谱中的N条水平线
"""
# （1）离散傅里叶变换：将空间域图像转换到频域，并将低频分量移动到频谱中心
dft_image = cv2.dft(np.float32(img), flags=cv2.DFT_COMPLEX_OUTPUT) # 进行DFT
dst_center = np.fft.fftshift(dft_image) # 将频谱的低频部分移动到中心
# （2）创建掩模：用于过滤特定频率分量（1表示保留，0表示过滤）
rows, cols = img.shape # 获取图像尺寸
mask = np.ones(shape=(rows, cols, 2), dtype=np.uint8) # 初始化全0掩模（保留所有频率）
row_center, col_center = int(rows / 2), int(cols / 2) # 计算图像中心坐标
if is_vertical_stripe:
"""去除竖向条纹，对应频谱中的中心水平线条 ———— 将中心水平线频率分量置0（过滤低频），保留中心区域的高频分量"""
mask[int(row_center - int(line_width // 2)): int(row_center + int(line_width // 2) + 1), :] = 0 # 清空中心水平线
mask[int(row_center - int(line_width // 2)): int(row_center + int(line_width // 2) + 1),
(col_center - center_radius): (col_center + center_radius)] = 1 # 保留中心区域
############################################################################################
# """（无效操作）在水平线的头尾分别保留 10 个像素"""
# mask[int(row_center - int(line_width // 2)): int(row_center + int(line_width // 2) + 1), :100] = 1
# mask[int(row_center - int(line_width // 2)): int(row_center + int(line_width // 2) + 1), -100:] = 1
############################################################################################
if is_horizontal_stripe:
"""去除横向条纹，对应频谱中的中心垂直线条 ———— 将中心垂直线频率分量置0（过滤低频），保留中心区域的高频分量"""
mask[:, int(col_center - int(line_width // 2)): int(col_center + int(line_width // 2) + 1)] = 0 # 清空中心垂直线
mask[(row_center - center_radius): (row_center + center_radius),
int(col_center - int(line_width // 2)): int(col_center + int(line_width // 2) + 1)] = 1 # 保留中心区域
############################################################################################
# """（无效操作）在竖直线的头尾分别保留 10 个像素"""
# mask[:100, int(col_center - int(line_width // 2)): int(col_center + int(line_width // 2) + 1)] = 1
# mask[-100:, int(col_center - int(line_width // 2)): int(col_center + int(line_width // 2) + 1)] = 1
############################################################################################
if is_multi:
"""去除竖向条纹，如果频谱有多条水平线 ———— 需要手动指定对称线条中位于上侧的y坐标，用于定位（暂未启用代码识别）"""
row_value1 = 55 # 325/55
row_value2 = rows - row_value1
mask[int(row_value1 - int(line_width // 2)): int(row_value1 + int(line_width // 2) + 1), :] = 0 # 清空水平线（上下对称）
mask[int(row_value2 - int(line_width // 2)): int(row_value2 + int(line_width // 2) + 1), :] = 0 # 清空水平线（上下对称）
# （3）应用掩模：对频域图像进行滤波
dst_center_mask = dst_center * mask # 将掩模应用到频域图像
############################################################################################
# """（无效操作）裁剪中心区域的周边区域，取其像素线条替换黑色线条"""
# # 中心水平线条
# over_area = line_width
# print(int(row_center - int(line_width // 2)), int(row_center + int(line_width // 2) + 1))
# print(int(col_center - int(line_width // 2)), int(col_center + int(line_width // 2) + 1))
#
# dst_center_mask[int(row_center - int(line_width // 2)): int(row_center + int(line_width // 2) + 1), :, 0] = \
# dst_center_mask[int(row_center - int(line_width // 2)) + over_area: int(row_center + int(line_width // 2) + 1) + over_area, :, 0]
# dst_center_mask[int(row_center - int(line_width // 2)): int(row_center + int(line_width // 2) + 1), :, 1] = \
# dst_center_mask[int(row_center - int(line_width // 2)) + over_area: int(row_center + int(line_width // 2) + 1) + over_area, :, 1]
#
# # 中心垂直线条
# dst_center_mask[:, int(col_center - int(line_width // 2)): int(col_center + int(line_width // 2) + 1), 0] = \
# dst_center_mask[:, int(col_center - int(line_width // 2)) + over_area: int(col_center + int(line_width // 2) + 1) + over_area, 0]
# dst_center_mask[:, int(col_center - int(line_width // 2)): int(col_center + int(line_width // 2) + 1), 1] = \
# dst_center_mask[:, int(col_center - int(line_width // 2)) + over_area: int(col_center + int(line_width // 2) + 1) + over_area, 1]
############################################################################################
# （4）可视化频域图像
magnitude_result = cv2.magnitude(dst_center_mask[:, :, 0], dst_center_mask[:, :, 1])
magnitude_result[magnitude_result == 0] = 1e-10 # 将零值替换为一个非常小的正数：1e-10
amplitude_image = 20 * np.log(magnitude_result) # 1e-10取log得到负值
amplitude_image[amplitude_image < 0] = 0 # 将负值替换为0
# cv2.namedWindow('amplitude_image', cv2.WINDOW_NORMAL) # 创建命名窗口
# cv2.imshow('amplitude_image', amplitude_image)
# cv2.waitKey(0) # 延迟一秒后自动关闭图像
# cv2.destroyAllWindows() # （同时）摧毁所有图窗
############################################################################################
# amplitude_image = 20 * np.log(cv2.magnitude(dst_center_mask[:, :, 0],
# dst_center_mask[:, :, 1])) # 计算幅度谱并进行对数压缩，以便可视化
# 警告提示：RuntimeWarning: divide by zero encountered in log
# 原因分析：np.log() 函数的输入值为零或非常接近零，而零的对数是未定义的（负无穷），因此会触发这个警告。
############################################################################################
# （5）反傅里叶变换：将频域图像转换回空间域
dst_origin = np.fft.ifftshift(dst_center_mask) # 将频率重新移回原始位置
idft_image = cv2.idft(dst_origin, flags=cv2.DFT_SCALE) # 进行反傅里叶变换
idft_image = cv2.magnitude(idft_image[:, :, 0], idft_image[:, :, 1]) # 获取幅度，得到重建的空间域图像
idft_image = idft_image.astype(img.dtype)
# （6）显示结果：原图、掩模、频域图像和重建图像
if is_visualization:
plt.subplot(131), plt.imshow(img, cmap='gray'), plt.axis('off'), plt.title('Original Image (Spatial Domain)')
plt.subplot(132), plt.imshow(amplitude_image, cmap='gray'), plt.axis('off'), plt.title('DFT Image (Frequency Domain)')
plt.subplot(133), plt.imshow(idft_image, cmap='gray'), plt.axis('off'), plt.title('IDFT Image (Reconstructed)')
# plt.pause(5) # 显示图像并暂停N秒
# plt.close() # 关闭图像窗口
plt.show() # 显示图像
return amplitude_image, idft_image
def connected_area(image_stack, min_area=50, max_area=50000):
from skimage.measure import label, regionprops
# （1）连通区域算法
# from skimage.filters import threshold_otsu
# print("image_stack.shape = ", image_stack.shape) # [通道, 高度, 宽度] = [depth, height, width]
# image_stack[image_stack <= gray_threshold] = 0 # 灰度强度阈值化
# otsu_threshold = threshold_otsu(image_stack) # 使用大津阈值法自适应计算图像的灰度阈值
# binary_image = image_stack > otsu_threshold # 根据阈值将图像转换为二值图像
binary_image = image_stack > 1 # 根据阈值将图像转换为二值图像
labeled_image = label(binary_image) # 对二值图像进行标记连通区域操作，生成标记图像
regions = regionprops(labeled_image) # 提取标记图像中的连通区域的属性信息
# （2）遍历所有区域
num = 0 # 统计连通区域的数量
for i, region in enumerate(regions):
# print(f"区域 {i + 1} 的大小: {region.area}")
if region.area < min_area or region.area > max_area:
# 获取小区域中每个像素的坐标，然后定位图像并置零
for coord in region.coords:
image_stack[coord[0], coord[1]] = 0
else:
num += 1
print(f"{num}/{len(regions)} 区域 {i + 1} 的大小: {region.area}")
return image_stack
def main(image, index=-1):
# （1）图像前处理: 高斯滤波 ———— （下面函数的设置参数）效果完全等同于（Fiji-Process-Filters-GaussianBlur-sigma=1.50）
GaussianBlur_image = cv2.GaussianBlur(image, (0, 0), sigmaX=1.5, sigmaY=0, borderType=cv2.BORDER_DEFAULT)
######################################################################
# （2）调用 striping_algorithm 函数
line_width = 1
center_radius = 0
is_vertical_stripe = True
is_horizontal_stripe = False
is_multi = False
is_visualization = False
amplitude_image, idft_image = fft_filter(img=GaussianBlur_image,
line_width=line_width,
center_radius=center_radius,
is_vertical_stripe=is_vertical_stripe,
is_horizontal_stripe=is_horizontal_stripe,
is_multi=is_multi,
is_visualization=is_visualization)
print(f"输出图像{index}的最大最小像素值：{[np.max(idft_image), np.min(idft_image)]}")
######################################################################
# # 循环测试不同的 center_radius 值
# center_radius_list = [1, 4, 16, 64, 256] # 可以根据需要调整
# for center_radius in center_radius_list:
# print(f"Testing with center_radius = {center_radius}")
# result = striping_algorithm(img,
# line_width=5,
# center_radius=center_radius,
# is_vertical_stripe=True,
# is_horizontal_stripe=True,
# is_multi=True)
# # 循环测试不同的 line_width 值
# line_width_list = [1, 3, 5, 7, 9, 11, 13, 20, 40] # 可以根据需要调整
# for line_width in line_width_list:
# print(f"Testing with line_width = {line_width}")
# amplitude_image, idft_image = fft_filter(GaussianBlur_image,
# line_width=line_width,
# center_radius=center_radius,
# is_vertical_stripe=is_vertical_stripe,
# is_horizontal_stripe=is_horizontal_stripe,
# is_multi=is_multi,
# is_visualization=is_visualization)
######################################################################
# （3）图像后处理: 连通区域（只保留指定面积范围内的像素值，完成去噪） ———— （对比）阈值分割: （1）需要手动调参（2）可能剔除有效区域的部分低像素值
denosie_image = connected_area(idft_image.copy(), min_area=50, max_area=50000)
# （4）保存结果
is_save = True
if is_save:
extension = f"{line_width}_{center_radius}_{is_vertical_stripe}_{is_horizontal_stripe}_{is_multi}_{index}"
tifffile.imwrite(os.path.join(folder_path, f"GaussianBlur_image{extension}.tif"), GaussianBlur_image)
tifffile.imwrite(os.path.join(folder_path, f"amplitude_image{extension}.tif"), amplitude_image)
tifffile.imwrite(os.path.join(folder_path, f"fft_image{extension}.tif"), idft_image)
tifffile.imwrite(os.path.join(folder_path, f"denosie_image{extension}.tif"), denosie_image)
"""
注意：只支持8bit图像
备注：16位图像的像素值范围是0-65535，而8位图像的像素值范围是0-255
"""
if __name__ == '__main__':
# image_path = r"F:\AI\directional_filter\test\MAX_Layer00_ROI00000_Loop00000_20250106161354700.tif"
image_path = r"F:\AI\directional_filter\test\MAX_Layer00_ROI00000_Loop00005_20250113142749411-1.tif" # 树枝状（min设置11，有大致效果）
# image_path = r"F:\AI\directional_filter\test\MAX_Layer00_ROI00000_Loop00170_20250109082420679.tif"
# image_path = r"F:\AI\directional_filter\test\MAX_Layer00_ROI00001_Loop00005_20250113142921377-2.tif"
image = tifffile.imread(image_path)
# image_path = r"F:\AI\directional_filter\sample\e91f9c67f6737002b97828d5b391b6fa.png"
# image = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE) # 读取图像并转换为灰度图像
if image is None:
raise FileNotFoundError(f"无法读取图像文件: {image_path}")
print("图像尺寸=", image.shape)
# （1）新建文件夹
folder_path = os.path.join(os.path.dirname(image_path), "result")
if os.path.exists(folder_path): # 判断路径是否存在
import shutil
shutil.rmtree(folder_path) # 删除目录树（递归删除目录及其下的所有文件和子目录）
os.makedirs(folder_path) # 重建文件夹
if image.ndim == 3:
for slice_index, slice_image in enumerate(image):
main(image[slice_index], slice_index) # 逐切片处理
elif image.ndim == 2:
main(image)
三、基于FFT的去条纹算法 - 参数测试 —— 【通用去条纹算法的核心因素】
3.0、使用说明
通过调参 + 定参 + 前图像处理（高斯滤波） + 后图像处理（灰度强度阈值分割），最终得到解决方案
3.1、频率滤波：线条类型 + 是否去除中心点 + 中心点半径 + 线宽
3.1.1、mask去除线条：水平 or 竖直 or 水平+竖直（保留中心点） —— 中心点半径
结论：
（1）仅去除水平或垂直，将出现反方向条纹；
（2）同时去除水平+竖直，去条纹效果最好；
（3）同时去除多条水平，有一定变化但很小；
（4）中心点半径越小，去条纹效果越好，但伪影越多；
（5）中心点半径越大，去条纹效果越差，但伪影越少；
3.1.2、mask去除线条：中心水平 + 中心竖直（保留中心点） —— 线宽
结论：（6）线宽越大（保留中心点），去条纹效果越差，伪影略微加重；
3.1.3、mask去除线条：中心水平 + 中心竖直（保留中心点） —— 裁剪中心区域的周边区域，取其像素线条替换黑色线条
结论：效果没有改善，没有结论。
3.1.4、mask去除线条：水平 or 竖直 or 水平+竖直（去除中心点）
结论：效果没有改善，没有结论。
3.2、前处理
结论：对输入图像进行高斯滤波，然后使用mask去除中心水平线条（去除中心点），再通过不同的线宽提取最优值，最终得到最优去条纹和去背景效果。
备注：由于去除中心点，导致图像的低频区域消失，因为具有去背景效果。
详细请看：【什么是频谱图像？？？】。
3.2.1、mask去除线条：水平 or 竖直 or 水平+竖直 —— 高斯滤波 + 是否去除中心点
结论：去除中心水平线条（去除中心点）效果最优
3.2.2、mask去除线条：中心水平（去除中心点） —— 高斯滤波 + 线宽
结论：
（1）线宽控制了 " 去伪影 " 效果
（2）没有最优参数，根据图像不同而有所不同（常用参数：线宽=1）
3.2.3、mask去除线条：中心水平（去除中心点） —— 高斯滤波 + sigma
结论：
（1）当sigma = 1.5 或 2 时（区别不大），效果最优
（2）当sigma > 3 及其以上，由于过度平滑，导致信息丢失，出现失真现象。
备注：sigma越大，高斯函数越平缓，滤波器覆盖范围更广，图像更平滑。
3.3、后处理
可变参数：
（1）高斯滤波的sigma参数（已测试，sigma=1.5或2参数最优）
（2）中心水平线条的线宽（需人工核查） —— 对于大部分图像，参数为1
（3）灰度强度阈值（需人工核查） —— 连通区域面积阈值（经验值：5050000）
分割效果对比：
- 连通区域:
（1）需要测试一个通用面积范围（不同需求，面积不同）
（2）只保留指定面积范围内的像素值，完成去噪- 阈值分割:
（1）需要手动调参
（2）可能剔除有效区域的部分低像素值
3.3.1、灰度阈值分割
结论：
（1）灰度阈值控制了 " 分割 " 效果，差异性较大；
（2）没有最优阈值参数，根据图像不同而有所不同；
3.3.2、连通区域分割
结论：
（1）通过区域面积去除噪声，通过经验值设置
（2）不同类型的图像不同而有所不同；
