image_to_image.utils.physgen_benchmark

  1"""
  2From: https://github.com/physicsgen/physicsgen/blob/main/eval_scripts/sound_metrics.py
  3
  4Slightly adjusted by Tobia Ippolito
  5
  6Original programmed by Martin Spitznagel
  7"""
  8import argparse
  9import os
 10
 11from PIL import Image
 12import numpy as np
 13import matplotlib.pyplot as plt
 14import pandas as pd
 15
 16import numba
 17from numba import jit
 18
 19from tqdm import tqdm 
 20
 21# for black-white image -> osm conversion
 22import itertools
 23import cv2
 24import shapely.geometry
 25import xml.etree.ElementTree as ET
 26
 27
 28def MAE(y_true, y_pred):
 29    return np.mean(np.abs(y_true - y_pred))
 30
 31def calc_mae(true_path, pred_path):
 32    pred_noisemap = (1 - np.array(
 33        Image.open(pred_path).convert("L"),
 34        dtype=np.float32
 35    ) / 255) * 100
 36
 37    true_noisemap = (1 - np.array(
 38        Image.open(true_path).convert("L"),
 39        dtype=np.float32
 40    ) / 255) * 100
 41    return MAE(true_noisemap, pred_noisemap)
 42
 43def calc_mape(true_path, pred_path):
 44    # Load and process the predicted and true noise maps
 45    pred_noisemap = (1 - np.array(
 46        Image.open(pred_path).convert("L"),
 47        dtype=np.float32
 48    ) / 255) * 100
 49    true_noisemap = (1 - np.array(
 50        Image.open(true_path).convert("L"),
 51        dtype=np.float32
 52    ) / 255) * 100
 53
 54    # Initialize an error map with zeros
 55    error_map = np.zeros_like(true_noisemap, dtype=np.float32)
 56
 57    # Find indices where true noisemap is not 0
 58    nonzero_indices = true_noisemap != 0
 59
 60    # Calculate percentage error where true noisemap is not 0
 61    error_map[nonzero_indices] = np.abs((true_noisemap[nonzero_indices] - pred_noisemap[nonzero_indices]) / true_noisemap[nonzero_indices]) * 100
 62
 63    # For positions where true noisemap is 0 but pred noisemap is not, set error to 100%
 64    zero_true_indices = (true_noisemap == 0) & (pred_noisemap != 0)
 65    error_map[zero_true_indices] = 100
 66
 67    # Calculate the MAPE over the whole image, ignoring positions where both are 0
 68    return np.mean(error_map)
 69
 70
 71@jit(nopython=True)
 72def ray_tracing(image_size, image_map):
 73    visibility_map = np.zeros((image_size, image_size))
 74    source = (image_size // 2, image_size // 2)
 75    for x in range(image_size):
 76        for y in range(image_size):
 77            dx = x - source[0]
 78            dy = y - source[1]
 79            dist = np.sqrt(dx*dx + dy*dy)
 80            steps = int(dist)
 81            if steps == 0:
 82                continue  # Skip the source point itself
 83            step_dx = dx / steps
 84            step_dy = dy / steps
 85
 86            visible = True  # Assume this point is visible unless proven otherwise
 87            ray_x, ray_y = source
 88            for _ in range(steps):
 89                ray_x += step_dx
 90                ray_y += step_dy
 91                int_x, int_y = int(ray_x), int(ray_y)
 92                if 0 <= int_x < image_size and 0 <= int_y < image_size:
 93                    if image_map[int_y, int_x] == 0:
 94                        visible = False
 95                        break
 96            visibility_map[y, x] = visible
 97    return visibility_map
 98
 99def compute_visibility(osm_path, image_size=256):
100    image_map = np.array(Image.open(osm_path).convert('L').resize((image_size, image_size)))
101    image_map = np.where(image_map > 0, 1, 0)
102    visibility_map = ray_tracing(image_size, image_map)
103    pixels_in_sight = np.logical_and(visibility_map == 1, image_map == 1)
104    pixels_not_in_sight = np.logical_and(visibility_map == 0, image_map == 1)
105    pixels_not_in_sight = np.where(image_map == 0, 0, pixels_not_in_sight)
106    pixels_in_sight = np.where(image_map == 0, 0, pixels_in_sight)
107    return pixels_in_sight, pixels_not_in_sight
108
109def masked_mae(true_labels, predictions):
110    # Convert to numpy arrays
111    true_labels = np.array(true_labels)
112    predictions = np.array(predictions)
113    
114    # Create a mask where true_labels are not equal to -1
115    mask = true_labels != -1
116    
117    # Filter arrays with the mask
118    true_labels = true_labels[mask]
119    predictions = predictions[mask]
120    
121    # Compute the MAE and return
122    return MAE(true_labels, predictions)
123
124def masked_mape(true_labels, predictions):
125    # Convert to numpy arrays
126    true_labels = np.array(true_labels)
127    predictions = np.array(predictions)
128
129    # Create a mask to exclude -1
130    mask = true_labels != -1
131    
132    # Apply the mask to filter arrays
133    true_labels_filtered = true_labels[mask]
134    predictions_filtered = predictions[mask]
135    
136    # Initialize an error map with zeros
137    error_map = np.zeros_like(true_labels_filtered, dtype=np.float32)
138
139    # Find indices where true noisemap is not 0
140    nonzero_indices = true_labels_filtered != 0
141
142    # Calculate percentage error where true noisemap is not 0
143    error_map[nonzero_indices] = np.abs((true_labels_filtered[nonzero_indices] - predictions_filtered[nonzero_indices]) / true_labels_filtered[nonzero_indices]) * 100
144
145    # For positions where true noisemap is 0 but pred noisemap is not, set error to 100%
146    zero_true_indices = (true_labels_filtered == 0) & (predictions_filtered != 0)
147    error_map[zero_true_indices] = 100
148
149    # Calculate the MAPE over the whole image, ignoring positions where both are 0
150    return np.mean(error_map)
151
152def calculate_sight_error(true_path, pred_path, osm_path):
153    pred_soundmap = (1 - np.array(
154        Image.open(pred_path).convert("L"),
155        dtype=np.float32
156    ) / 255) * 100
157    true_soundmap = (1 - np.array(
158        Image.open(true_path).convert("L"),
159        dtype=np.float32
160    ) / 255) * 100
161    _, true_pixels_not_in_sight = compute_visibility(osm_path)
162
163    in_sight_soundmap = true_soundmap.copy()
164    not_in_sight_soundmap = true_soundmap.copy()
165    
166    in_sight_pred_soundmap = pred_soundmap.copy()
167    not_in_sight_pred_soundmap = pred_soundmap.copy()
168    
169    #only get the pixels in sight
170    for x in range(256):
171        for y in range(256):
172            if true_pixels_not_in_sight[y, x] == 0:
173                not_in_sight_soundmap[y, x] = -1
174                not_in_sight_pred_soundmap[y, x] = -1
175            else:
176                in_sight_soundmap[y, x] = -1
177                in_sight_pred_soundmap[y, x] = -1
178
179    return masked_mae(in_sight_soundmap, in_sight_pred_soundmap), masked_mae(not_in_sight_soundmap, not_in_sight_pred_soundmap), masked_mape(in_sight_soundmap, in_sight_pred_soundmap), masked_mape(not_in_sight_soundmap, not_in_sight_pred_soundmap)
180
181def evaluate_sample(true_path, pred_path, osm_path=None) -> (float, float, float, float):
182    mae = calc_mae(true_path, pred_path)
183    mape = calc_mape(true_path, pred_path)
184
185    mae_in_sight, mae_not_in_sight, mape_in_sight, mape_not_in_sight = None, None, None, None
186    if osm_path is not None:
187        mae_in_sight, mae_not_in_sight, mape_in_sight, mape_not_in_sight = calculate_sight_error(true_path, pred_path, osm_path)
188    return mae, mape, mae_in_sight, mae_not_in_sight, mape_in_sight, mape_not_in_sight
189
190# main function for evaluation
191if __name__ == "__main__":
192
193    parser = argparse.ArgumentParser()
194    parser.add_argument("--data_dir", type=str, default="data/true")
195    parser.add_argument("--pred_dir", type=str, default="data/pred")
196    parser.add_argument("--osm_dir", type=str, default="none")
197    parser.add_argument("--output", type=str, default="evaluation.csv")
198    args = parser.parse_args()
199
200    data_dir = args.data_dir
201    pred_dir = args.pred_dir
202    osm_dir = args.osm_dir
203    osm_dir = None if osm_dir.lower() == "none" else osm_dir
204    output = args.output
205
206    output_path, _ = os.path.split(output)
207    os.makedirs(output_path, exist_ok=True)
208
209    results = []
210
211    # get files
212    file_names = os.listdir(data_dir)
213
214    # Use tqdm to create a progress bar for the loop
215    for sample_name in tqdm(file_names, total=len(file_names), desc="Evaluating samples"):
216        pred_ = os.path.join(pred_dir, sample_name)
217        real_ = os.path.join(data_dir, sample_name)
218        osm_ = os.path.join(osm_dir, sample_name)
219
220        # Check if prediction is available
221        if not os.path.exists(f"{pred_dir}/{sample_name}"):
222            print(f"Prediction for sample {sample_name} not found.")
223            print(f"{pred_dir}/{sample_name}")
224            continue
225        mae, mape, mae_in_sight, mae_not_in_sight, mape_in_sight, mape_not_in_sight = evaluate_sample(real_, pred_, osm_path=osm_) # adjust prediction naming if needed
226        results.append([sample_name, mae, mape, mae_in_sight, mae_not_in_sight, mape_in_sight, mape_not_in_sight])
227
228    results_df = pd.DataFrame(results, columns=["sample_id", "MAE", "MAPE", "LoS_MAE", "NLoS_MAE", "LoS_wMAPE", "NLoS_wMAPE"])
229    results_df.to_csv(output, index=False)
230    print(results_df[["MAE", "MAPE", "LoS_MAE", "NLoS_MAE", "LoS_wMAPE", "NLoS_wMAPE"]].describe())