import numpy as np
import cv2
import matplotlib.pyplot as plt
from sklearn.metrics import mean_absolute_error, mean_squared_error, precision_score, recall_score, roc_curve, auc
from skimage.transform import resize  # For resizing the mask
from keras.models import load_model
from keras_preprocessing.image import img_to_array, load_img
from sklearn.metrics import accuracy_score, f1_score

# 1. 加载u-net 模型
model = load_model(r"C:\Users\jack8\OneDrive\LSTM and ARIMA\LSTM\lstm_model.h5")

# 2. 定义辅助函数
def image_to_patches(image, patch_size):
    patches = []
    positions = []
    img_height, img_width = image.shape[:2]
    patch_height, patch_width = patch_size
    for y in range(0, img_height, patch_height):
        for x in range(0, img_width, patch_width):
            patch = image[y:y+patch_height, x:x+patch_width]
            # 如果最后的块尺寸不足，进行填充
            if patch.shape[0] != patch_height or patch.shape[1] != patch_width:
                pad_height = patch_height - patch.shape[0]
                pad_width = patch_width - patch.shape[1]
                patch = np.pad(patch, ((0, pad_height), (0, pad_width), (0, 0)), mode='constant')
            patches.append(patch)
            positions.append((y, x))
    return np.array(patches), positions

def patches_to_image(patches, positions, image_shape):
    reconstructed_image = np.zeros(image_shape, dtype=patches.dtype)
    patch_height, patch_width = patches.shape[1:3]
    for patch, (y, x) in zip(patches, positions):
        y_end = min(y + patch_height, image_shape[0])
        x_end = min(x + patch_width, image_shape[1])
        h = y_end - y
        w = x_end - x
        reconstructed_image[y:y_end, x:x_end] = patch[:h, :w]
    return reconstructed_image

# 3. 加载并预处理输入图像
# 加载大图
input_img_path = r"C:\Users\jack8\OneDrive\LSTM and ARIMA\LSTM\data\south\Finish\segmentation-vidas\2024.jpg"
new_img = load_img(input_img_path)
new_img_array = img_to_array(new_img) / 255.0

# 获取原始图像的尺寸
original_height, original_width = new_img_array.shape[:2]

# 定义模型输入尺寸
patch_size = (640, 640)  # 根据您的模型的输入尺寸调整

# 将大图切割成小块
patches, positions = image_to_patches(new_img_array, patch_size)

# 4. 对每个小块进行预测
predicted_patches = []
for patch in patches:
    patch_input = np.expand_dims(patch, axis=0)
    prediction = model.predict(patch_input)
    predicted_patches.append(prediction[0])
predicted_patches = np.array(predicted_patches)  # Shape: (num_patches, patch_height, patch_width, num_classes)

# 将预测结果二值化
prediction_binary_patches = (predicted_patches > 0.5).astype(np.uint8)

# 重建预测的掩膜图
prediction_mask = patches_to_image(prediction_binary_patches, positions, (original_height, original_width, 1))

# 5. 加载并预处理实际的掩膜图像（mask）
mask_img_path = r"C:\Users\jack8\OneDrive\LSTM and ARIMA\LSTM\data\south\Finish\cover\2024.png"
mask_img = load_img(mask_img_path, color_mode='grayscale')
mask_img_array = img_to_array(mask_img) / 255.0
mask_img_array = mask_img_array[:, :, 0]  # 只保留单通道

# 确保实际掩膜的尺寸与预测掩膜一致
if mask_img_array.shape != prediction_mask.shape[:2]:
    # 如果尺寸不一致，调整实际掩膜的尺寸
    mask_img_array_resized = resize(mask_img_array, prediction_mask.shape[:2], preserve_range=True)
else:
    mask_img_array_resized = mask_img_array

# 为了确保掩膜是二值的，应用阈值
mask_img_array_resized = (mask_img_array_resized > 0.5).astype(np.uint8)

# 确保预测概率掩膜是连续的概率值
prediction_prob_mask = patches_to_image(predicted_patches, positions, (original_height, original_width, 1))
prediction_prob_mask = prediction_prob_mask[:, :, 0]  # 只取单通道

# 展平成一维数组
y_true_roc = mask_img_array_resized.flatten()
y_scores = prediction_prob_mask.flatten()

# 检查 y_true 是否为二进制
unique_values = np.unique(y_true_roc)
print("Unique values in y_true_roc:", unique_values)

if not set(unique_values).issubset({0, 1}):
    print("Error: y_true_roc contains non-binary values.")
    # 可以选择强制二值化
    y_true_roc = (y_true_roc > 0.5).astype(int)

# 6. 应用闭运算
# 将预测的掩膜图像（prediction_mask）转换为 0 和 255 的图像，用于 OpenCV 函数
prediction_mask_uint8 = (prediction_mask[:, :, 0] * 255).astype(np.uint8)

# 定义结构元素的大小和形状
kernel_size = 29
kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (kernel_size, kernel_size))

# 应用闭运算
closed_mask = cv2.morphologyEx(prediction_mask_uint8, cv2.MORPH_CLOSE, kernel)

# 将处理后的掩膜图像转换回 0 和 1 的二值图像
closed_mask_binary = (closed_mask > 127).astype(np.uint8)

# 展平成一维数组
true_values = mask_img_array_resized.flatten()
predicted_values = closed_mask_binary.flatten()

# 检查值范围
print("True values min:", np.min(true_values), "max:", np.max(true_values))
print("Predicted values min:", np.min(predicted_values), "max:", np.max(predicted_values))

# 二值化处理
true_values = (true_values > 0.5).astype(np.float32)
predicted_values = (predicted_values > 0.5).astype(np.float32)

# 再次检查值范围
print("After thresholding:")
print("True values min:", np.min(true_values), "max:", np.max(true_values))
print("Predicted values min:", np.min(predicted_values), "max:", np.max(predicted_values))

# 检查数组形状是否一致
print("True values shape:", true_values.shape)
print("Predicted values shape:", predicted_values.shape)

# 计算 MAPE
def mean_absolute_percentage_error(y_true, y_pred):
    # 获取 y_true 不为零的索引
    non_zero_indices = y_true != 0
    # 避免除以零，忽略 y_true 为零的元素
    y_true_non_zero = y_true[non_zero_indices]
    y_pred_non_zero = y_pred[non_zero_indices]
    
    if len(y_true_non_zero) == 0:
        print("Warning: No non-zero elements in y_true. MAPE is undefined.")
        mape = np.nan
    else:
        mape = np.mean(np.abs((y_true_non_zero - y_pred_non_zero) / y_true_non_zero)) * 100
    
    return mape

mape = mean_absolute_percentage_error(true_values, predicted_values)
print(f"MAPE: {mape}%")

# 计算 MAE
mae = mean_absolute_error(true_values, predicted_values)
print(f"MAE: {mae}")

# 计算 RMSE
rmse = np.sqrt(mean_squared_error(true_values, predicted_values))
print(f"RMSE: {rmse}")

# 计算 IoU
intersection = np.logical_and(true_values, predicted_values).sum()
union = np.logical_or(true_values, predicted_values).sum()
iou_score = intersection / union if union != 0 else 1  # 避免除以零
print("IoU (Jaccard Index):", iou_score)

# 计算准确率
accuracy = accuracy_score(true_values, predicted_values)
print("Accuracy:", accuracy)

# 计算 F1 Score
f1 = f1_score(true_values, predicted_values, zero_division=0)
print("F1 Score:", f1)

# 计算 Dice 系数
dice_score = 2 * intersection / (true_values.sum() + predicted_values.sum()) if (true_values.sum() + predicted_values.sum()) != 0 else 1
print("Dice Coefficient:", dice_score)

# 计算精确率和召回率
precision = precision_score(true_values, predicted_values, zero_division=0)
recall = recall_score(true_values, predicted_values, zero_division=0)
print("Precision:", precision)
print("Recall:", recall)

# 可视化结果
plt.figure(figsize=(18, 6))
plt.subplot(1, 3, 1)
plt.title("Ground Truth Mask")
plt.imshow(mask_img_array_resized, cmap='gray')
plt.axis('off')

plt.subplot(1, 3, 2)
plt.title("Prediction Mask")
plt.imshow(prediction_mask[:, :, 0], cmap='gray')
plt.axis('off')

plt.subplot(1, 3, 3)
plt.title("Closed Mask Binary")
plt.imshow(closed_mask_binary, cmap='gray')
plt.axis('off')

plt.show()

# ----------------------- 添加 ROC 曲线的绘制部分 -----------------------

# 为了绘制 ROC 曲线，使用模型预测的概率值和二值化后的真实标签
# y_true_roc 已确保是二进制

# 检查 y_true_roc 是否为二进制
unique_values = np.unique(y_true_roc)
print("Unique values in y_true_roc after binarization:", unique_values)

if not set(unique_values).issubset({0, 1}):
    raise ValueError("y_true_roc 必须是二进制标签 (0 和 1)。请检查掩膜图像的二值化步骤。")

# 计算 ROC 曲线和 AUC
fpr, tpr, thresholds = roc_curve(y_true_roc, y_scores)
roc_auc = auc(fpr, tpr)
print(f"ROC AUC: {roc_auc}")

# 绘制 ROC 曲线
plt.figure(figsize=(8, 6))
plt.plot(fpr, tpr, color='blue', lw=2, label=f'ROC curve (AUC = {roc_auc:.2f})')
plt.plot([0, 1], [0, 1], color='red', lw=2, linestyle='--', label='Random Guess')
plt.xlim([0.0, 1.0])
plt.ylim([0.0, 1.05])
plt.xlabel('False Positive Rate', fontsize=14)
plt.ylabel('True Positive Rate', fontsize=14)
plt.title('Receiver Operating Characteristic (ROC) Curve', fontsize=16)
plt.legend(loc="lower right", fontsize=12)
plt.grid(True)
plt.show()


##請改成其他的都不要，只需要Closed Mask Binary這張圖，並且Closed Mask Binary需和Ground Truth Mask尺寸大小一樣