wizualizacja baselineow

import math
import numpy as np
import pandas as pd
from typing import Tuple

def find_baseline_carti(mask: np.ndarray) -> Tuple[Tuple[int, int], Tuple[int, int]]:
    """
    Funkcja znajdująca baseline, który ma najwięcej punktów pokrywających się
    z górnym obrysem bony roof. Algorytm iteruje po pikselach wysuniętych
    najbardziej na prawo punktów bony roof, które są oddalone o 30 pikseli
    w górę i 10 pikseli w dół w celu zmniejszenia ilości obliczeń.

    Funkcja przyjmuje:
    - maskę w formacie png z zaznaczonymi strukturami
    Funkcja zwraca:
    - współrzędne prawego i lewego punktu baseline
    """
    height, width = mask.shape

    label_upper_contour_bony = 5
    label_upper_contour_carti = 4
    binary_mask_bony = (mask == label_upper_contour_bony).astype(np.uint8)
    binary_mask_carti = (mask == label_upper_contour_carti).astype(np.uint8)
    upper_contour_bony = []

    # Znalezienie punktów górnego obrysu bony roof
    for x in range(width):
        column_bony = binary_mask_bony[:, x]
        if np.any(column_bony):
            y_bony = np.where(column_bony)[0][0]
            upper_contour_bony.append((x, y_bony))

    # Znalezienie najbardziej na prawo wysuniętego punktu carti roof
    x_max_carti = -1
    for x in range(width):
        column_carti = binary_mask_carti[:, x]
        if np.any(column_carti):
            x_max_carti = x  # Ustalamy x_max_carti na najbardziej wysunięty na prawo punkt
            y_max_carti = np.where(column_carti)[0][0]  # Używamy tylko y, ale nie potrzebujemy

    # Ustalanie prawego punktu bony roof (+1 piksel w prawo)
    rightmost_point_x = x_max_carti + 1 if x_max_carti != -1 else None

    # Znalezienie punktu bony roof dla współrzędnej x = x_max_carti + 1
    bony_point_y = None
    if rightmost_point_x is not None and rightmost_point_x < width:
        column_bony = binary_mask_bony[:, rightmost_point_x]
        if np.any(column_bony):
            bony_point_y = np.where(column_bony)[0][0]  # Znajdź punkt na bony roof w kolumnie x_max_carti + 1

    # Definiowanie punktów startowych (iteracja 30 pikseli w górę i 10 w dół od punktu bony roof)
    start_points = []
    if bony_point_y is not None:
        start_points = [(rightmost_point_x, y) for y in range(max(0, bony_point_y - 5),
                                                              min(height, bony_point_y + 5))]

    best_line = []
    max_overlap = 0
    rightmost_point = (rightmost_point_x, bony_point_y) if rightmost_point_x is not None else (None, None)

    for x0, y0 in start_points:
        for angle in range(360):
            line_points = []

            angle_rad = math.radians(angle)
            sin_angle = math.sin(angle_rad)
            cos_angle = math.cos(angle_rad)

            # Definiujemy maksymalną liczbę kroków
            n_max = 2 * (width + height)
            n = np.arange(n_max)

            x = x0 + n * cos_angle
            y = y0 - n * sin_angle

            x_rounded = np.round(x).astype(int)
            y_rounded = np.round(y).astype(int)

            # Filtrujemy punkty znajdujące się w granicach obrazu
            valid_mask = (x_rounded >= 0) & (x_rounded < width) & \
                         (y_rounded >= 0) & (y_rounded < height)

            x_valid = x_rounded[valid_mask]
            y_valid = y_rounded[valid_mask]

            line_points = list(zip(x_valid, y_valid))

            # Obliczenie overlapu między linią a konturem
            overlap = len(set(line_points) & set((p for p in upper_contour_bony if p[0] <= rightmost_point_x)))
            if overlap > max_overlap:
                max_overlap = overlap
                best_line = line_points

    # Znalezienie najbardziej na lewo wysuniętego punktu linii
    leftmost_point = None
    if best_line:
        for point in best_line:
            x, y = point
            if mask[y, x] == label_upper_contour_bony:
                if leftmost_point is None or x < leftmost_point[0]:
                    leftmost_point = (x, y)

    if leftmost_point is None:
        leftmost_point = (None, None)

    return rightmost_point, leftmost_point


def find_baseline(mask: np.ndarray) -> Tuple[Tuple[int, int], Tuple[int, int]]:
    """
    Funkcja znajdująca baseline, który ma najwięcej punktów pokrywających się
    z górnym obrysem bony roof. Algorytm iteruje po pikselach wysuniętych
    najbardziej na prawo punktów bony roof, które są oddalone o 10 pikseli
    w górę i w dół w celu zmniejszenia ilości obliczeń.

    Funkcja przyjmuje:
    - maskę w formacie png z zaznaczonymi strukturami
    Funkcja zwraca:
    - współrzędne prawego i lewego punktu baseline
    """
    height, width = mask.shape

    label_upper_contour = 5
    binary_mask = (mask == label_upper_contour).astype(np.uint8)
    upper_contour = []

    # Znalezienie punktów górnego obrysu (konturu)
    for x in range(width):
        column = binary_mask[:, x]
        if np.any(column):
            y = np.where(column)[0][0]
            upper_contour.append((x, y))

    # Znalezienie najbardziej na prawo wysuniętego punktu
    x_max = max(p[0] for p in upper_contour)
    y_max = max(p[1] for p in upper_contour if p[0] == x_max)

    # Iterowanie po pikselach znajdujących się 30 pikseli ponad i 10
    # pikseli podnajbardziej wysuniętymi na prawo punktami bony roof

    start_points = [(x, y) for x, y in upper_contour if x == x_max] +\
        [(x_max, y) for y in range(max(0, y_max - 30),
                                   min(height, y_max + 10))]

    best_line = []
    max_overlap = 0
    rightmost_point = None

    for x0, y0 in start_points:
        for angle in range(360):
            line_points = []

            angle_rad = math.radians(angle)
            sin_angle = math.sin(angle_rad)
            cos_angle = math.cos(angle_rad)

            # Definiujemy maksymalną liczbę kroków
            n_max = 2 * (width + height)
            n = np.arange(n_max)

            x = x0 + n * cos_angle
            y = y0 - n * sin_angle

            x_rounded = np.round(x).astype(int)
            y_rounded = np.round(y).astype(int)

            # Filtrujemy punkty znajdujące się w granicach obrazu
            valid_mask = (x_rounded >= 0) & (x_rounded < width) & \
                (y_rounded >= 0) & (y_rounded < height)

            x_valid = x_rounded[valid_mask]
            y_valid = y_rounded[valid_mask]

            line_points = list(zip(x_valid, y_valid))

            # Obliczenie overlapu między linią a konturem
            overlap = len(set(line_points) & set(upper_contour))
            if overlap > max_overlap:
                max_overlap = overlap
                best_line = line_points
                rightmost_point = (x0, y0)

    # Znalezienie najbardziej na lewo wysuniętego punktu linii
    leftmost_point = None
    if best_line:
        for point in best_line:
            x, y = point
            if mask[y, x] == label_upper_contour:
                if leftmost_point is None or x < leftmost_point[0]:
                    leftmost_point = (x, y)

    if leftmost_point is None:
        leftmost_point = (None, None)

    return rightmost_point, leftmost_point


def old_find_baseline(class_mask: np.ndarray) -> Tuple[
        Tuple[int, int] | None, Tuple[int, int] | None]:
    """Bazując na położeniu bony i carti roof określa punkt początkowy
    oraz końcowy baseline dla maksymalnej odległości do 5 pikseli.
    Algorytm przeszukuje obszar jednego piksela dookoła (czyli 3x3),
    a jeżeli nie znajdzie punktów, to zwiększa obszar poszukiwań aż
    do progu jakim jest 5 pikseli. Dodatkowo na początku algorytm
    zawęża obszar poszukiwań tylko do obszarów bony roof i carti
    roof oraz pikselom im sąsiadującym

    Args:
        class_mask (np.ndarray): Maska wygenerowana z txt

    Returns:
        Współrzędne punktów baseline
    """
    mask_4_indices = np.argwhere(class_mask == 4)
    mask_5_indices = np.argwhere(class_mask == 5)
    if mask_4_indices.size == 0 or mask_5_indices.size == 0:
        return None, None

    min_y = max(0, min(np.min(mask_4_indices[:, 0]), np.min(mask_5_indices[
        :, 0])))
    max_y = min(class_mask.shape[0], max(np.max(mask_4_indices[:, 0]), np.max(
        mask_5_indices[:, 0])))
    min_x = max(0, min(np.min(mask_4_indices[:, 1]), np.min(
        mask_5_indices[:, 1])))
    max_x = min(class_mask.shape[1], max(np.max(mask_4_indices[:, 1]), np.max(
        mask_5_indices[:, 1])))

    for distance in range(1, 6):
        right_point = None
        left_point = None

        for y in range(min_y, max_y + 1):
            for x in range(min_x, max_x + 1):
                surroundings = class_mask[max(0, y - distance):min(
                    class_mask.shape[0], y + distance + 1),
                                          max(0, x - distance):min(
                                        class_mask.shape[1], x + distance + 1)]
                if (surroundings == 5).any() and (surroundings == 4).any():
                    if right_point is None or (x > right_point[0] or (
                            x == right_point[0] and y > right_point[1])):
                        right_point = (x, y)
                    if left_point is None or (x < left_point[0] or (
                            x == left_point[0] and y > left_point[1])):
                        left_point = (x, y)

        if right_point is not None and left_point is not None:
            if distance == 5 and right_point == left_point:
                left_point = (left_point[0] - 1, left_point[1])

            if right_point == left_point:
                continue

            return right_point, left_point

    return None, None


import os
import cv2
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image

# Słownik kolorów dla poszczególnych klas
CLASS_COLORS = {
    0: {'color_rgb': (0, 0, 0), 'label': 'background'},
    1: {'color_rgb': (255, 0, 0), 'label': 'chondroosseous border'},
    2: {'color_rgb': (0, 255, 0), 'label': 'femoral head'},
    3: {'color_rgb': (0, 0, 255), 'label': 'labrum'},
    4: {'color_rgb': (255, 255, 0), 'label': 'cartilagineous roof'},
    5: {'color_rgb': (0, 255, 255), 'label': 'bony roof'},
    6: {'color_rgb': (159, 2, 250), 'label': 'bony rim'},
    7: {'color_rgb': (255, 132, 0), 'label': 'lower limb'},
    8: {'color_rgb': (255, 0, 255), 'label': 'baseline'},
    9: {'color_rgb': (66, 135, 245), 'label': 'lower limb template'},
    10: {'color_rgb': (255, 69, 0), 'label': 'lower limb - alg. v3'}
}

def visualize_masks_with_baseline(masks_dir: str, output_directory: str):
    # Upewnij się, że katalog wyjściowy istnieje
    if not os.path.exists(output_directory):
        os.makedirs(output_directory)

    # Iterowanie po katalogu z maskami
    for mask_filename in os.listdir(masks_dir):
        if mask_filename.endswith('.png'):
            mask_path = os.path.join(masks_dir, mask_filename)
            mask = cv2.imread(mask_path, cv2.IMREAD_GRAYSCALE)  # Wczytaj maskę w skali szarości

            # Znajdź punkty baseline za pomocą nowego algorytmu
            rightmost_point, leftmost_point = find_baseline_carti(mask)

            # Debug: Sprawdzenie punktów
            print(f"Maska: {mask_filename}, Punkt prawy: {rightmost_point}, Punkt lewy: {leftmost_point}")

            # Stwórz obraz RGB, na który będziemy rysować
            mask_rgb = np.zeros((mask.shape[0], mask.shape[1], 3), dtype=np.uint8)

            # Iteracja po klasach w masce i kolorowanie
            for class_value, class_info in CLASS_COLORS.items():
                mask_rgb[mask == class_value] = class_info['color_rgb']

            # Rysowanie białej linii między punktami zwróconymi przez find_baseline
            if rightmost_point and leftmost_point and rightmost_point != (None, None) and leftmost_point != (None, None):
                cv2.line(mask_rgb, rightmost_point, leftmost_point, (255, 255, 255), 1)  # Biała linia
                cv2.circle(mask_rgb, rightmost_point, 1, (255, 255, 255), -1)  # Prawy punkt
                cv2.circle(mask_rgb, leftmost_point, 1, (255, 255, 255), -1)  # Lewy punkt

            # Zapisz wizualizację do katalogu
            output_path = os.path.join(output_directory, f'visualization_{mask_filename}')
            cv2.imwrite(output_path, mask_rgb)

            # Wyświetlenie obrazu
            plt.figure(figsize=(10, 10))
            plt.imshow(cv2.cvtColor(mask_rgb, cv2.COLOR_BGR2RGB))  # Zamiana na RGB dla wyświetlenia
            plt.title(f'Visualization for {mask_filename}')
            plt.axis('off')

# Przykładowe wywołanie funkcji
masks_directory = './Angles/dane/masks_baseline'
output_directory = './Angles/dane/visualizations_carti'
#visualize_masks_with_baseline(masks_directory, output_directory)


import os
import numpy as np
import matplotlib.pyplot as plt
from typing import Tuple, List
import math

def slope(x1: int, y1: int, x2: int, y2: int) -> float:
    """Calculates the slope of the line through two points."""
    # Sprawdzanie czy punkty są None
    if None in (x1, y1, x2, y2):
        print(f"Invalid slope points: ({x1}, {y1}), ({x2}, {y2}). Skipping.")
        return None  # Zwracamy None, aby pominąć obliczenia

    if x2 - x1 == 0:  # Unikamy dzielenia przez 0
        return float('inf')  # Linia pionowa
    return (y2 - y1) / (x2 - x1)

def find_angle(line1_point1: Tuple[int, int],
               line1_point2: Tuple[int, int],
               line2_point1: Tuple[int, int],
               line2_point2: Tuple[int, int],
               return_absolute_value: bool = True) -> float:
    """Finds the angle between two lines."""
    if len(line1_point1) != 2 or len(line1_point2) != 2 or \
            len(line2_point1) != 2 or len(line2_point2) != 2:
        print('Point with wrong number of coordinates')
        return np.NaN  # Zwracamy NaN zamiast zgłaszać wyjątek
        
    M1 = slope(line1_point1[0], line1_point1[1], line1_point2[0], line1_point2[1])
    M2 = slope(line2_point1[0], line2_point1[1], line2_point2[0], line2_point2[1])

    PI = 3.14159265

    # Sprawdzamy, czy M1 i M2 nie są None
    if M1 is None or M2 is None:
        print(f"Invalid slopes: M1={M1}, M2={M2}. Skipping angle calculation.")
        return None

    if M2 != float('inf') and M1 * M2 != -1:
        if return_absolute_value:
            angle = abs((M2 - M1) / (1 + M1 * M2))
        else:
            angle = (M2 - M1) / (1 + M1 * M2)

        ret = math.atan(angle)
        val = (ret * 180) / PI
        return round(val, 1)

    elif M1 * M2 == -1:
        return 90.0
    else:
        return np.NaN

def read_doctor_coordinates(file_path: str) -> Tuple[Tuple[int, int], Tuple[int, int]]:
    """Reads doctor coordinates from a text file."""
    with open(file_path, 'r') as file:
        lines = file.readlines()
        
        if len(lines) < 2:
            return None, None
        
        try:
            point1 = tuple(map(int, lines[-2].strip().split()))
            point2 = lines[-1].strip()
            
            if point2 == "0":
                return None, None
            
            point2 = tuple(map(int, point2.split()))
        except ValueError as e:
            print(f"Error reading coordinates from {file_path}: {e}")
            return None, None
        
        return point1, point2

def calculate_angle(pt1: Tuple[int, int], pt2: Tuple[int, int], baseline: Tuple[Tuple[int, int], Tuple[int, int]]) -> float:
    """Calculates the angle between the line formed by pt1 and pt2 and the baseline."""
    # Check if points are None
    if pt1 is None or pt2 is None or baseline[0] is None or baseline[1] is None:
        print(f"Invalid points: pt1={pt1}, pt2={pt2}, baseline={baseline}. Skipping.")
        return None  # Return None to skip the calculation

    angle_doctor = find_angle(pt1, pt2, baseline[0], baseline[1], return_absolute_value=False)
    
    if angle_doctor is None:
        return None
    
    angle_baseline = find_angle(baseline[0], baseline[1], baseline[0], baseline[1], return_absolute_value=True)

    angle_diff = angle_baseline - angle_doctor
    return angle_diff

def compare_baselines_to_doctor(masks_dir: str, coords_dir: str, method: int) -> np.ndarray:
    """Compares baseline methods with doctor's coordinates and generates a histogram."""
    angle_diffs = []

    for filename in os.listdir(masks_dir):
        if not filename.endswith('.png'):
            continue
        
        mask_path = os.path.join(masks_dir, filename)
        coords_path = os.path.join(coords_dir, filename.replace('.png', '.txt'))
        
        doctor_coords = read_doctor_coordinates(coords_path)
        if doctor_coords == (None, None):
            print(f"Skipping file {coords_path}, invalid doctor coordinates.")
            continue
        
        mask = np.array(Image.open(mask_path))
        
        if method == 1:
            rightmost, leftmost = find_baseline_carti(mask)
        elif method == 2:
            rightmost, leftmost = find_baseline(mask)
        elif method == 3:
            rightmost, leftmost = old_find_baseline(mask)
        else:
            raise ValueError("Invalid method specified.")
        
        if rightmost is None or leftmost is None:
            print(f"Skipping file {mask_path}, invalid baseline points.")
            continue

        baseline = (leftmost, rightmost)
        
        angle_diff = calculate_angle(doctor_coords[0], doctor_coords[1], baseline)
        if angle_diff is not None:
            angle_diffs.append(angle_diff)

    return np.array(angle_diffs)

def plot_histogram(angle_diffs: np.ndarray):
    """Plots histogram of angle differences."""
    plt.figure(figsize=(10, 6))
    
    if angle_diffs.size == 0:
        print("No angle differences to plot.")
        return

    mean_angle = np.mean(angle_diffs)
    std_angle = np.std(angle_diffs)
    
    plt.hist(angle_diffs, bins=np.linspace(-15, 15, 30), edgecolor='black')
    
    plt.title(f'Histogram of Angle Differences - Wersja najstarsza\nMean: {mean_angle:.2f}°, Std: {std_angle:.2f}°')
    plt.xlabel('Angle Difference (degrees)')
    plt.ylabel('Frequency')
    plt.xlim(-15, 15)
    plt.grid(axis='y')
    plt.axvline(0, color='red', linestyle='dashed', linewidth=1)
    plt.show()

# Example usage:
masks_directory = './Angles/dane/masks_baseline'
coords_directory = './Angles/dane/templates_baseline'
method_type = 3  # Choose the method (1, 2, or 3)
angle_differences = compare_baselines_to_doctor(masks_directory, coords_directory, method_type)
plot_histogram(angle_differences)