python-service/app/unity7.py - main - BELEGO/GIS_APPLICATION

python-service/app/unity7.py
237 строк10 KB
HACKER_LAPTOP
В python-service создал geometric quad-crop: unity7.py детектит 4 угла повёрнутой сцены через cv2.minAreaRect + convexHull и сохраняет quad в geo_images.crop_meta и create_voting_masks.py строит маску из этого же quad'а через cv2.fillPoly — снаружи quad'а прозрачно, внутри чёрный фон, огонь белый, image и mask теперь геометрически синхронны.
20 май 2026, 05:00
20 май 2026, 05:00e98c6b7
# debug_visualize_only.py
import os
import sys
import numpy as np
import cv2
from osgeo import gdal
import matplotlib.pyplot as plt
from PIL import Image
from fastapi import HTTPException

# Определяем корень проекта
current_file_path = os.path.abspath(__file__)
project_root = os.path.dirname(os.path.dirname(os.path.dirname(current_file_path)))
sys.path.append(os.path.dirname(current_file_path))

from Models.GeoImage import GeoImage
from database import SessionLocal


def _detect_scene_quad(rgb_image: np.ndarray, threshold: float = 0.01,
                        morph_kernel: int = 31) -> np.ndarray | None:
    """
    Возвращает 4 точки повёрнутого прямоугольника сцены (quad) в координатах
    исходного raster'а. None если детекция не удалась.

    Алгоритм:
      1. valid = пиксели где хотя бы один RGB > threshold (no-data → False)
      2. morphological CLOSE — заращивает тёмные дырки ВНУТРИ сцены
         (вода, тени), чтобы они не дробили контур
      3. convexHull от крупнейшего external-контура
      4. minAreaRect → boxPoints → 4 точки
    """
    valid = (np.max(rgb_image, axis=2) > threshold).astype(np.uint8) * 255
    kernel = cv2.getStructuringElement(cv2.MORPH_RECT, (morph_kernel, morph_kernel))
    valid_filled = cv2.morphologyEx(valid, cv2.MORPH_CLOSE, kernel)

    contours, _ = cv2.findContours(valid_filled, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
    if not contours:
        return None
    largest = max(contours, key=cv2.contourArea)
    if cv2.contourArea(largest) < 1000:
        return None
    hull = cv2.convexHull(largest)
    rect = cv2.minAreaRect(hull)
    quad = cv2.boxPoints(rect).astype(np.int32)
    return quad


def _aabb_of_quad(quad: np.ndarray, max_w: int, max_h: int) -> dict:
    """{x,y,width,height} — AABB вокруг quad с клипом по границам raster'а."""
    x, y, w, h = cv2.boundingRect(quad)
    x = max(0, int(x))
    y = max(0, int(y))
    w = min(max_w - x, int(w))
    h = min(max_h - y, int(h))
    return {"x": x, "y": y, "width": w, "height": h}


def _spatial_alpha_from_quad(shape: tuple, quad: np.ndarray) -> np.ndarray:
    """Возвращает uint8 mask формы shape: 255 внутри quad, 0 снаружи."""
    alpha = np.zeros(shape, dtype=np.uint8)
    cv2.fillPoly(alpha, [quad], 255)
    return alpha


def _save_crop_meta(base_name: str, src_width: int, src_height: int,
                    bbox: dict, quad: list) -> None:
    """Сохраняет crop_meta (включая quad) в geo_images.crop_meta по productid."""
    try:
        with SessionLocal() as s:
            img = s.query(GeoImage).filter(
                GeoImage.ProductId.ilike(f"%{base_name}%")
            ).first()
            if img:
                img.CropMeta = {
                    "src_width": src_width,
                    "src_height": src_height,
                    "crop": bbox,
                    "quad": quad,
                }
                s.commit()
                print(f"[INFO] crop_meta сохранена для {base_name}: bbox={bbox}, quad={quad}")
            else:
                print(f"[WARN] GeoImage не найден для {base_name} — crop_meta не сохранена")
    except Exception as e:
        print(f"[WARN] не удалось сохранить crop_meta для {base_name}: {e}")


def combine_and_visualize_tiff_bands_by_path(
        relative_folder_path: str,
        target_bands=("B1", "B2", "B3", "B4", "B5", "B6", "B7"),
        visualize_rgb_bands=("B7", "B6", "B2")
):
    """
    Объединяет каналы TIFF в многослойный TIFF и создает PNG с прозрачным фоном.
    Возвращает путь к созданному TIFF файлу.
    """
    # Определяем абсолютный путь к папке с TIFF файлами
    folder_path = os.path.join(project_root, relative_folder_path)

    if not os.path.exists(folder_path):
        raise FileNotFoundError(f"Папка не найдена: {folder_path}")

    # Ищем все TIFF файлы
    all_files = os.listdir(folder_path)
    tif_files = [f for f in all_files if f.lower().endswith(".tif")]

    if not tif_files:
        raise FileNotFoundError(f"В папке нет TIFF файлов: {folder_path}")

    # Находим файлы для каждого целевого канала
    image_files = []
    for band in target_bands:
        found = [f for f in tif_files if f"_{band}." in f]
        if not found:
            raise FileNotFoundError(f"Канал {band} не найден в папке {folder_path}")
        image_files.append(found[0])

    # Формируем имя и путь для объединенного TIFF
    base_name = image_files[0].rsplit('_', 1)[0]
    merged_path = os.path.join(folder_path, f"{base_name}_merged.tif")

    if os.path.exists(merged_path):
        # Переиспользуем существующий merged.tif — но всегда обновляем PNG + crop_meta
        print(f"[INFO] merged.tif существует — переиспользуем: {merged_path}")
        ds = gdal.Open(merged_path)
        if not ds:
            raise ValueError(f"Не удалось открыть существующий merged: {merged_path}")
        width, height = ds.RasterXSize, ds.RasterYSize
        img = np.zeros((height, width, ds.RasterCount), dtype=np.float32)
        for i in range(ds.RasterCount):
            img[:, :, i] = ds.GetRasterBand(i + 1).ReadAsArray()
        ds = None
    else:
        # Загружаем данные для получения размеров
        first_band_path = os.path.join(folder_path, image_files[0])
        first_band_ds = gdal.Open(first_band_path)
        if not first_band_ds:
            raise ValueError(f"Не удалось открыть файл: {first_band_path}")

        width, height = first_band_ds.RasterXSize, first_band_ds.RasterYSize

        # Загружаем все каналы в массив
        img = np.zeros((height, width, len(image_files)), dtype=np.float32)

        for i, filename in enumerate(image_files):
            band_path = os.path.join(folder_path, filename)
            ds = gdal.Open(band_path)
            if ds:
                img[:, :, i] = ds.ReadAsArray()

        # Создаем объединенный TIFF файл
        driver = gdal.GetDriverByName("GTiff")
        ds_out = driver.Create(merged_path, width, height, len(image_files), gdal.GDT_Float32)
        for i in range(len(image_files)):
            ds_out.GetRasterBand(i + 1).WriteArray(img[:, :, i])
        ds_out.FlushCache()
        ds_out = None

    # Создаем словарь для быстрого доступа к индексам каналов
    band_dict = {band: idx for idx, band in enumerate(target_bands)}

    # Функция для нормализации канала
    def normalize_band(band):
        min_val, max_val = np.min(band), np.max(band)
        if max_val > min_val:
            return (band - min_val) / (max_val - min_val)
        return np.zeros_like(band)

    # Выбираем индексы каналов для RGB визуализации
    try:
        rgb_indices = [band_dict[b] for b in visualize_rgb_bands]
    except KeyError as e:
        # Если не все каналы для визуализации есть, возвращаем путь к TIFF без PNG
        return merged_path

    # Создаем RGB изображение
    rgb_bands = [normalize_band(img[:, :, i]) for i in rgb_indices]
    rgb_image = np.dstack(rgb_bands)

    # === GEOMETRIC CROP по quad'у повёрнутой сцены ===
    # 1) Находим 4 угла сцены (rotated quad) — НЕ по цвету пикселей, а по геометрии
    # 2) Spatial alpha = fillPoly(quad)  → внутри quad alpha=1, снаружи=0
    #    (тёмные пиксели вода/тень ВНУТРИ сцены остаются видимыми)
    # 3) AABB = boundingRect(quad)  → slicing-окно для PNG
    quad = _detect_scene_quad(rgb_image)
    if quad is not None:
        bbox = _aabb_of_quad(quad, width, height)
        spatial_alpha = _spatial_alpha_from_quad((height, width), quad)
        alpha_f = spatial_alpha.astype(np.float32) / 255.0
        rgba_image = np.dstack((rgb_image, alpha_f))

        x, y, w, h = bbox["x"], bbox["y"], bbox["width"], bbox["height"]
        rgba_image = rgba_image[y:y + h, x:x + w, :]

        _save_crop_meta(base_name, width, height, bbox, quad.tolist())
    else:
        # fallback — если quad не нашёлся (странный raster): просто RGB без альфы
        print(f"[WARN] не удалось детектировать quad сцены для {base_name} — сохраняем без crop")
        alpha_f = np.ones(rgb_image.shape[:2], dtype=np.float32)
        rgba_image = np.dstack((rgb_image, alpha_f))

    # Сохраняем PNG только в Front/public
    folder_path_png = os.path.join(project_root, "Front", "public")
    os.makedirs(folder_path_png, exist_ok=True)

    png_filename = f"{base_name}_merged.png"
    png_path = os.path.join(folder_path_png, png_filename)

    # Конвертируем в 8-бит и сохраняем через PIL (корректно сохраняет альфа-канал)
    rgba_8bit = (rgba_image * 255).astype(np.uint8)
    img_pil = Image.fromarray(rgba_8bit, mode='RGBA')
    img_pil.save(png_path, 'PNG')

    # Если скрипт запущен напрямую, показываем результат
    if __name__ == "__main__":
        plt.figure(figsize=(12, 10))
        plt.imshow(rgba_image)
        plt.title(f'RGBA с прозрачностью')
        plt.axis('off')
        plt.tight_layout()
        plt.show()

    return merged_path


if __name__ == "__main__":
    # Тестовый вызов для отладки
    test_path = "geo-images/LC09_L1TP_145021_20220507_20230417_02_T1/"

    try:
        result = combine_and_visualize_tiff_bands_by_path(
            relative_folder_path=test_path
        )
        print(f"Создан TIFF: {result}")
    except Exception as e:
        print(f"Ошибка: {e}")