import cv2
import numpy as np
import matplotlib.pyplot as plt

def niblack_threshold(image_path, window_size=15, k=-0.2):
    """
    实现Niblack阈值算法
    :param image_path: 输入图像路径
    :param window_size: 邻域窗口大小
    :param k: 阈值参数
    :return: 原图和Niblack二值化后的图像
    """
    # 读取图像
    img = cv2.imread(image_path)
    gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
    gray = gray.astype(np.float64)
    
    # 获取图像尺寸
    height, width = gray.shape
    
    # 计算邻域半径
    radius = window_size // 2
    
    # 初始化输出图像
    binary = np.zeros_like(gray)
    
    # 避免边界问题，只处理内部像素
    for i in range(radius, height - radius):
        for j in range(radius, width - radius):
            # 提取邻域
            neighborhood = gray[i-radius:i+radius+1, j-radius:j+radius+1]
            
            # 计算邻域的均值和标准差
            mean_val = np.mean(neighborhood)
            std_val = np.std(neighborhood)
            
            # 计算Niblack阈值
            threshold = mean_val + k * std_val
            
            # 应用阈值
            if gray[i, j] > threshold:
                binary[i, j] = 255
            else:
                binary[i, j] = 0
    
    return img, cv2.cvtColor(binary.astype(np.uint8), cv2.COLOR_GRAY2BGR)

def manual_niblack_threshold(image_path, window_size=15, k=-0.2):
    """
    手动实现Niblack阈值算法
    :param image_path: 输入图像路径
    :param window_size: 邻域窗口大小
    :param k: 阈值参数
    :return: 二值化后的图像
    """
    # 读取图像
    img = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE)
    img = img.astype(np.float64)
    
    # 获取图像尺寸
    height, width = img.shape
    
    # 计算邻域半径
    radius = window_size // 2
    
    # 初始化输出图像
    binary = np.zeros_like(img)
    
    # 使用向量化操作提高效率
    # 首先计算均值图
    mean_img = np.zeros_like(img)
    std_img = np.zeros_like(img)
    
    for i in range(radius, height - radius):
        for j in range(radius, width - radius):
            neighborhood = img[i-radius:i+radius+1, j-radius:j+radius+1]
            mean_img[i, j] = np.mean(neighborhood)
            std_img[i, j] = np.std(neighborhood)
    
    # 计算阈值图
    threshold_img = mean_img + k * std_img
    
    # 应用阈值
    binary[img > threshold_img] = 255
    
    return cv2.cvtColor(img.astype(np.uint8), cv2.COLOR_GRAY2BGR), cv2.cvtColor(binary.astype(np.uint8), cv2.COLOR_GRAY2BGR)

def optimized_niblack_threshold(image_path, window_size=15, k=-0.2):
    """
    优化的Niblack阈值算法（使用积分图像加速计算）
    :param image_path: 输入图像路径
    :param window_size: 邻域窗口大小
    :param k: 阈值参数
    :return: 原图和优化的Niblack二值化结果
    """
    # 读取图像
    img = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE)
    img = img.astype(np.float64)
    
    height, width = img.shape
    radius = window_size // 2
    
    # 使用积分图像来快速计算局部均值和方差
    integral_img = np.zeros((height + 1, width + 1))
    integral_img_sq = np.zeros((height + 1, width + 1))
    
    # 构建积分图像
    for i in range(1, height + 1):
        for j in range(1, width + 1):
            integral_img[i, j] = img[i-1, j-1] + integral_img[i-1, j] + integral_img[i, j-1] - integral_img[i-1, j-1]
            integral_img_sq[i, j] = img[i-1, j-1]**2 + integral_img_sq[i-1, j] + integral_img_sq[i, j-1] - integral_img_sq[i-1, j-1]
    
    # 初始化输出图像
    binary = np.zeros_like(img)
    
    # 计算局部均值和标准差
    for i in range(radius, height - radius):
        for j in range(radius, width - radius):
            # 使用积分图像计算局部和
            x1, y1 = i - radius, j - radius
            x2, y2 = i + radius + 1, j + radius + 1
            
            region_sum = integral_img[x2, y2] - integral_img[x1, y2] - integral_img[x2, y1] + integral_img[x1, y1]
            region_sum_sq = integral_img_sq[x2, y2] - integral_img_sq[x1, y2] - integral_img_sq[x2, y1] + integral_img_sq[x1, y1]
            
            # 计算均值和方差
            region_area = window_size * window_size
            mean_val = region_sum / region_area
            variance = (region_sum_sq / region_area) - (mean_val ** 2)
            std_val = np.sqrt(max(variance, 0))  # 确保非负
            
            # 计算阈值
            threshold = mean_val + k * std_val
            
            # 应用阈值
            if img[i, j] > threshold:
                binary[i, j] = 255
            else:
                binary[i, j] = 0
    
    return cv2.cvtColor(img.astype(np.uint8), cv2.COLOR_GRAY2BGR), cv2.cvtColor(binary.astype(np.uint8), cv2.COLOR_GRAY2BGR)

def compare_niblack_parameters(image_path):
    """
    比较不同参数的Niblack算法效果
    :param image_path: 输入图像路径
    :return: 原图和不同参数的二值化结果
    """
    # 读取图像
    img = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE)
    img = img.astype(np.float64)
    
    # 不同参数组合
    params = [
        {'window_size': 15, 'k': -0.2},
        {'window_size': 15, 'k': -0.5},
        {'window_size': 21, 'k': -0.2},
        {'window_size': 21, 'k': -0.5}
    ]
    
    results = {}
    
    for i, param in enumerate(params):
        window_size = param['window_size']
        k = param['k']
        
        height, width = img.shape
        radius = window_size // 2
        binary = np.zeros_like(img)
        
        for x in range(radius, height - radius):
            for y in range(radius, width - radius):
                neighborhood = img[x-radius:x+radius+1, y-radius:y+radius+1]
                mean_val = np.mean(neighborhood)
                std_val = np.std(neighborhood)
                threshold = mean_val + k * std_val
                
                if img[x, y] > threshold:
                    binary[x, y] = 255
                else:
                    binary[x, y] = 0
        
        results[f'Window:{window_size}, k:{k}'] = cv2.cvtColor(binary.astype(np.uint8), cv2.COLOR_GRAY2BGR)
    
    return cv2.cvtColor(img.astype(np.uint8), cv2.COLOR_GRAY2BGR), results

def niblack_vs_other_methods(image_path):
    """
    比较Niblack算法与其他阈值方法
    :param image_path: 输入图像路径
    :return: 原图和不同方法的二值化结果
    """
    # 读取图像
    img = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE)
    
    # Niblack算法
    niblack_img, niblack_result = manual_niblack_threshold(image_path, 15, -0.2)
    
    # 全局阈值
    _, global_result = cv2.threshold(img, 127, 255, cv2.THRESH_BINARY)
    
    # Otsu阈值
    _, otsu_result = cv2.threshold(img, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)
    
    # 自适应阈值
    adaptive_mean = cv2.adaptiveThreshold(img, 255, cv2.ADAPTIVE_THRESH_MEAN_C, cv2.THRESH_BINARY, 11, 2)
    adaptive_gaussian = cv2.adaptiveThreshold(img, 255, cv2.ADAPTIVE_THRESH_GAUSSIAN_C, cv2.THRESH_BINARY, 11, 2)
    
    results = {
        'Niblack': niblack_result,
        'Global (127)': cv2.cvtColor(global_result, cv2.COLOR_GRAY2BGR),
        'Otsu': cv2.cvtColor(otsu_result, cv2.COLOR_GRAY2BGR),
        'Adaptive Mean': cv2.cvtColor(adaptive_mean, cv2.COLOR_GRAY2BGR),
        'Adaptive Gaussian': cv2.cvtColor(adaptive_gaussian, cv2.COLOR_GRAY2BGR)
    }
    
    return cv2.cvtColor(img, cv2.COLOR_GRAY2BGR), results

def niblack_with_noise_analysis(image_path):
    """
    分析Niblack算法对噪声的鲁棒性
    :param image_path: 输入图像路径
    :return: 原图、含噪声图像和不同噪声水平下的Niblack结果
    """
    # 读取图像
    img = cv2.imread(image_path, cv2.IMREAD_GRAYSCALE)
    
    # 添加不同水平的噪声
    noise_levels = [0.01, 0.05, 0.1]
    results = {}
    
    for level in noise_levels:
        # 添加高斯噪声
        noise = np.random.normal(0, level * 255, img.shape)
        noisy_img = np.clip(img + noise, 0, 255).astype(np.uint8)
        
        # 应用Niblack算法
        _, binary = manual_niblack_threshold_from_array(noisy_img, 15, -0.2)
        results[f'Noise Level: {level}'] = binary
    
    return cv2.cvtColor(img, cv2.COLOR_GRAY2BGR), results

def manual_niblack_threshold_from_array(img_array, window_size=15, k=-0.2):
    """
    从数组输入实现Niblack阈值算法
    :param img_array: 输入图像数组
    :param window_size: 邻域窗口大小
    :param k: 阈值参数
    :return: 二值化后的图像
    """
    img = img_array.astype(np.float64)
    height, width = img.shape
    radius = window_size // 2
    
    # 初始化输出图像
    binary = np.zeros_like(img)
    
    for i in range(radius, height - radius):
        for j in range(radius, width - radius):
            neighborhood = img[i-radius:i+radius+1, j-radius:j+radius+1]
            mean_val = np.mean(neighborhood)
            std_val = np.std(neighborhood)
            threshold = mean_val + k * std_val
            
            if img[i, j] > threshold:
                binary[i, j] = 255
            else:
                binary[i, j] = 0
    
    return cv2.cvtColor(binary.astype(np.uint8), cv2.COLOR_GRAY2BGR)

# 使用示例
if __name__ == "__main__":
    # 注意：需要提供实际的图像路径
    # img, result = niblack_threshold('image.jpg', 15, -0.2)
    # manual_img, manual_result = manual_niblack_threshold('image.jpg', 15, -0.2)
    # opt_img, opt_result = optimized_niblack_threshold('image.jpg', 15, -0.2)
    # comp_img, comp_results = compare_niblack_parameters('image.jpg')
    # vs_img, vs_results = niblack_vs_other_methods('image.jpg')
    # noise_img, noise_results = niblack_with_noise_analysis('image.jpg')
    pass