"""PyVista 3D インタラクティブビューアモジュール． シーンのメッシュ・LED・センサを PyVista で可視化する． 照度ベイク結果 (irradiance) があればメッシュを照度で着色する． ## 背面カリング（backface_cull 引数） show_scene の `backface_cull=True` を指定すると，ベイク着色経路・ニュートラルグレー背景の 各 add_mesh に `culling="back"` を渡し，カメラに背を向けた面を描画しない VTK 標準の レンダリング時カリングを有効化する． 効果: - 不透明表示（bake_opacity=1.0）でも手前の壁が自動的に消えて内部が見える （dollhouse 表示）． - cull_outward（幾何的な静的カリング）と組み合わせると， 外向き面の除去（静的）＋ カメラ向きに応じた動的カリングの 2 段構成になる． - bake_opacity<1 の半透明表示でも手前の面による混色ブレンドが軽減される． 注意（巻き順依存）: 面の法線は STL の巻き順に依存する．cull_outward で内向き face のみが残っている場合， 内向き face の法線はシーン中心向き = 外から覗くとカメラに背を向けることになり， 期待どおり手前の壁が消える． ただし巻き順が乱れた STL では面の一部が欠けて見える可能性がある （SmTIAS シーンの検証 STL では妥当）． マーカー（LED/センサ球・矢印）には適用しない． ## 投影モード（--render 引数） luminance_map（2D 輝度マップ, np.ndarray）が指定された場合、 評価板メッシュ（id="evaluation_board"）の各頂点を orthographic センサ規約で 2D 輝度マップへ射影し、per-vertex scalar で着色する． 投影規約の詳細は visualize/projection.py の docstring を参照． 投影モードでは: - 評価板が主役：viridis カラーマップ + scalar bar 1 本（単位: Irradiance (W/m²)） - 輝度マップは E_board = π * luminance_map / ρ で照度 (W/m²) に換算する （ρ = 評価板材料の拡散反射率、Lambertian 拡散面前提） - 他メッシュは半透明ニュートラル（グレー, opacity=0.2） - 配色は percentile クリップで高ダイナミックレンジ潰れを防止 ## ハイブリッドモード（luminance_map + irradiance 同時指定時） luminance_map と irradiance dict が**両方**指定された場合： - 評価板メッシュ（id 特定は projection.find_eval_board_mesh_id）: 投影輝度を E_board = π * L / ρ で照度 (W/m²) に換算し、 壁ベイク irradiance と**共有スケール・共有カラーバー 1 本**で着色する． - 評価板以外の全メッシュ（壁・筐体等）: irradiance dict の per-face 照度で着色（viridis・半透明 bake_opacity）． 評価板 E_board を含めた全データで算出した共有 clim を適用する． - depth peeling 有効化（半透明メッシュ混在のため）． ## ベイクモード（irradiance dict のみ指定時） irradiance dict が指定された場合は各メッシュの per-face 照度で着色する． **全メッシュ共有カラーマップ（共有 clim）**: 全メッシュ全 face の値を結合して percentile クリップ（p1–p99）で [vmin, vmax] を算出し、全 add_mesh に渡す． これにより全メッシュが同一スケールで着色される（path ベイクで壁面と評価板を比較可能）． スカラバーは 1 本のみ（単位: Irradiance (W/m²)）． mode="direct"（shadow ray）でも mode="path"（多重反射込み）でも共通して使用する． ## log スケール着色（bake_scale="log"、既定） 壁面ベイクの irradiance は高ダイナミックレンジ（LED 至近: ~4400 W/m²、壁一般: 数十〜数百）を 持つため、linear + percentile クリップでは少数の超高輝度 face が clim を支配し 大多数の壁が viridis 最下色（濃紫）に圧縮される． bake_scale="log" 時の処理順序（ベイク着色経路のみ）: 1. hide_below フィルタ: irradiance <= hide_below の cell を除去（linear 値基準）． 2. log 変換: log10(max(irradiance, eps)) を全残 face に適用（eps=1e-3）． 3. log 空間で共有 clim 算出（percentile クリップ） → scalar bar タイトル「Irradiance (W/m², log10)」． 4. bake_denoise メディアンデノイズ: log 空間スカラに対して face 隣接メディアン（N 回）． MC ノイズ（firefly＝一面だけ飛び抜けた値）をエッジを保ったまま除去する． bake_scale="linear" 時は log 変換なしの cell_data スカラに適用． 5. bake_smooth_radius ガウス平滑化: log 空間スカラに対して物理半径指定のガウス重み付き 平均平滑化（ボクセル splat/gather 方式，O(N)，密度ヘテロ耐性あり， メモリ O(N + ユニークボクセル数) で有界）． bake_denoise の後・bake_smooth の前． 三角形スケールの MC ノイズを消し，radius より大きい照明構造は保つ． 6. bake_smooth 平滑化: log 空間スカラに対して Laplacian 拡散（N 回 往復）． 7. cell → point 補間（Gouraud シェーディング用）． 8. add_mesh． bake_scale="linear" 時の処理順序: 1. bake_denoise メディアンデノイズ: cell_data スカラ（linear 値）に適用． 2. bake_smooth_radius ガウス平滑化: cell_data スカラ（linear 値）に適用（bake_denoise の後）． 3. bake_smooth Laplacian 拡散（N 回 往復）． 4. hide_below フィルタ（平滑化後）． 5. add_mesh． bake_scale="linear" のデノイズ適用は平滑化前（cell_data がある段階）で行う（後方互換）． ## 評価板照度換算（E_board） 評価板の輝度マップ L（W/(m²·sr)）は Lambertian 拡散面の関係式により 入射照度 E（W/m²）に換算できる: E = π * L / ρ ρ: 評価板材料の拡散反射率（scene_config の materials から取得、既定 0.85） Lambertian 拡散面では L = ρ * E / π が成立（SPEC_03「筐体（材料）」参照）． この換算により評価板テクスチャと壁ベイク照度が同一物理量（W/m²）になり、 共有スケール・単一カラーバーで直接比較できる． ## 表示スタイルプリセット（style 引数） show_scene の `style` 引数で発表スライド用の体裁に切り替えられる． - "default"（既定）: 現行表示そのまま（後方互換）．座標軸ウィジェット・ bounds グリッドを表示し、背景は PyVista 既定色、scalar bar は既定体裁． センサマーカーを表示する．カラーマップ既定は viridis． - "slide": 発表スライド用クリーン体裁．背景を純白にし、座標軸ウィジェットと bounds グリッドを非表示にし、アンチエイリアスを有効化し、scalar bar を 黒文字でクリーンな体裁にする．センサマーカー（青い球＋視線矢印）を非表示にする （hide_sensor_marker=True と同等）．カラーマップ既定は hot． scalar bar の体裁は `_build_scalar_bar_args` で 1 箇所に集約し、viewer 内の 全 scalar bar 追加箇所（投影テクスチャのダミーバー・per-vertex フォールバック・ 壁ベイク）に同一の scalar_bar_args を適用する． ## カラーマップ選択（cmap 引数） show_scene の `cmap` 引数でカラーマップ名を指定できる． - None（既定）: style に応じて自動決定する． - "default" → "viridis"（後方互換） - "slide" → "hot"（照明均一性評価と同一のカラーマップ） - 文字列指定: 全スタイルでそのカラーマップを使用する． 不正な名前は warning を出して style 既定にフォールバックする． """ from __future__ import annotations import logging from pathlib import Path from typing import TYPE_CHECKING import numpy as np from sm_tias_lightsim.geometry.loader import load_stl from sm_tias_lightsim.geometry.subdivide import subdivide_mesh_to_max_edge from sm_tias_lightsim.geometry.transforms import scale_vertices from sm_tias_lightsim.scene.schema import ( AreaLight, DirectionalLight, PointLight, SceneConfig, ) if TYPE_CHECKING: import pyvista as pv _logger = logging.getLogger(__name__) # YAML 座標 (mm) → m 変換係数 _MM_TO_M = 0.001 # 方向ベクトルの長さがこれ以下なら縮退とみなしてスキップする _DIR_NORM_EPSILON = 1e-9 # 光源・センサを示す球の半径（m 単位） _MARKER_SPHERE_RADIUS = 0.005 # AreaLight の放射方向矢印のスケール（m 単位） _AREA_LIGHT_ARROW_SCALE = 0.02 # DirectionalLight の矢印をシーン対角線長に対して拡大する比率 _DIRECTIONAL_ARROW_SCALE_RATIO = 0.3 # センサの視線矢印をシーン対角線長に対してクランプする比率 _SENSOR_ARROW_SCALE_RATIO = 0.2 # 投影モードで評価板以外のメッシュの不透明度 _BACKGROUND_MESH_OPACITY = 0.2 # 直接ベイクモードで着色メッシュ以外の不透明度 _UNBAKED_MESH_OPACITY = 0.5 # 評価板反射率のフォールバック値（SPEC_03 基準設計値 printer_white） _FALLBACK_RHO = 0.85 # log10 変換の下限：irradiance が 0 や負の face を log10 で -3 に揃える _LOG10_EPS = 1e-3 # サポートする表示スタイル（未知値は "default" にフォールバック） _VALID_STYLES = ("default", "slide") # slide スタイルの scalar bar フォントサイズ（黒文字クリーン体裁） _SLIDE_SCALAR_BAR_TITLE_FONT_SIZE = 20 _SLIDE_SCALAR_BAR_LABEL_FONT_SIZE = 14 def _build_scalar_bar_args(title: str, style: str) -> dict: """scalar bar の表示引数（scalar_bar_args）をスタイルに応じて組み立てる． viewer 内には scalar bar を追加する箇所が複数ある（投影テクスチャのダミーバー・ per-vertex フォールバック・壁ベイク）ため、体裁をこの 1 箇所に集約して全箇所に 同一の引数を適用する． Args: title: scalar bar のタイトル文字列． style: 表示スタイル（"default" または "slide"）． Returns: plotter.add_mesh の scalar_bar_args に渡す辞書． "slide" のとき黒文字・フォントサイズ調整・fmt="%.2f" のクリーン体裁を付与する． "default" のときは title のみ（現行どおり・後方互換）． """ if style == "slide": return { "title": title, "color": "black", "title_font_size": _SLIDE_SCALAR_BAR_TITLE_FONT_SIZE, "label_font_size": _SLIDE_SCALAR_BAR_LABEL_FONT_SIZE, "fmt": "%.2f", } return {"title": title} def _build_pyvista_polydata(vertices: np.ndarray, faces: np.ndarray) -> pv.PolyData: """頂点・面配列から PyVista PolyData を構築する． PyVista の faces 形式は各行が [n_vertices, v0, v1, v2, ...] の 1D フラット配列． Args: vertices: 頂点配列 (N_v, 3)，m 単位． faces: 三角形面配列 (N_f, 3)，頂点インデックス． Returns: PyVista PolyData オブジェクト． """ import pyvista as pv n_faces = len(faces) # PyVista 形式: [3, v0, v1, v2, 3, v0, v1, v2, ...] pv_faces = np.hstack([np.full((n_faces, 1), 3, dtype=np.int64), faces]).flatten() return pv.PolyData(vertices.astype(np.float32), pv_faces) def _get_scene_bbox_center_and_diag( scene_config: SceneConfig, project_root: Path ) -> tuple[np.ndarray, float]: """シーン全メッシュのバウンディングボックス中心と対角線長を返す． DirectionalLight の矢印位置を決めるために使う． メッシュが空なら原点と 1.0 を返す． Args: scene_config: SceneConfig． project_root: STL 相対パスの基点ディレクトリ． Returns: (center (3,), diagonal_length) のタプル． """ all_vertices: list[np.ndarray] = [] for mesh_cfg in scene_config.geometry.meshes: try: raw_path = Path(mesh_cfg.path) stl_path = raw_path if raw_path.is_absolute() else project_root / mesh_cfg.path loaded = load_stl(stl_path, mesh_cfg.unit) scaled = scale_vertices(loaded, target_unit="m") all_vertices.append(scaled.vertices) except Exception as exc: _logger.warning("バウンディングボックス計算でメッシュを読み込めませんでした: %s", exc) if not all_vertices: return np.zeros(3), 1.0 combined = np.vstack(all_vertices) mins = combined.min(axis=0) maxs = combined.max(axis=0) center = (mins + maxs) / 2.0 diagonal = float(np.linalg.norm(maxs - mins)) return center, max(diagonal, 1e-6) def _apply_percentile_clip(values: np.ndarray, clip_low: float, clip_high: float) -> np.ndarray: """percentile クリップを適用して高 DR 潰れを防ぐ． Args: values: スカラ配列． clip_low: 下限 percentile（%）． clip_high: 上限 percentile（%）． Returns: クリップ済みの配列（元と同 dtype）． """ p_low = float(np.percentile(values, clip_low)) p_high = float(np.percentile(values, clip_high)) return np.clip(values, p_low, p_high) def _get_eval_board_reflectance( scene_config: SceneConfig, eval_board_id: str, ) -> float: """評価板メッシュの拡散反射率 ρ をスカラ値で返す． scene_config の materials から評価板の material を特定し、 DiffuseMaterial の reflectance を Rec.709 輝度でスカラ化して返す． 取得できない場合は 0.85 をフォールバック値とする（SPEC_03 基準設計値）． Args: scene_config: 検証済み SceneConfig． eval_board_id: 評価板メッシュの id 文字列． Returns: 拡散反射率 ρ（スカラ, 0 < ρ ≤ 1）． """ from sm_tias_lightsim.scene.schema import DiffuseMaterial try: eval_mesh_cfg = next(m for m in scene_config.geometry.meshes if m.id == eval_board_id) material = scene_config.materials.get(eval_mesh_cfg.material) if material is None or not isinstance(material, DiffuseMaterial): _logger.warning( "評価板材料が DiffuseMaterial でないため ρ=%.2f にフォールバックします" "（material=%s）", _FALLBACK_RHO, eval_mesh_cfg.material, ) return _FALLBACK_RHO ref = material.reflectance if isinstance(ref, (int, float)): rho = float(ref) else: # RGB タプル → Rec.709 輝度でスカラ化 # Rec.709: Y = 0.2126 R + 0.7152 G + 0.0722 B r, g, b = float(ref[0]), float(ref[1]), float(ref[2]) rho = 0.2126 * r + 0.7152 * g + 0.0722 * b if rho <= 0.0 or rho > 1.0: _logger.warning( "評価板の reflectance=%.4f が有効範囲外のため ρ=%.2f にフォールバックします", rho, _FALLBACK_RHO, ) return _FALLBACK_RHO _logger.info( "評価板 reflectance 取得: material='%s' raw_ref=%s ρ=%.4f", eval_mesh_cfg.material, ref, rho, ) return rho except StopIteration: _logger.warning( "評価板メッシュ '%s' が見つかりません．ρ=%.2f にフォールバックします", eval_board_id, _FALLBACK_RHO, ) return _FALLBACK_RHO except Exception as exc: _logger.warning( "評価板反射率の取得中にエラー: %s．ρ=%.2f にフォールバックします", exc, _FALLBACK_RHO, ) return _FALLBACK_RHO def _make_irradiance_texture( irradiance_map: np.ndarray, vmin: float, vmax: float, bake_scale: str, log10_eps: float, cmap_name: str = "viridis", ) -> pv.Texture: """照度マップ（W/m²）を共有スケールで着色した PyVista Texture を生成する． 評価板の E_board（照度マップ, W/m²）を壁ベイクと同一の [vmin, vmax] および log/linear 変換規約で RGB uint8 画像化し pv.Texture として返す． 共有 clim（壁ベイクと評価板の結合 percentile）は呼び出し側で算出済みの vmin/vmax をそのまま受け取る（この関数は再算出しない）． 処理手順: 1. bake_scale=="log" のとき log10(max(irradiance_map, log10_eps)) を適用する． vmin/vmax も log 空間の値（呼び出し側が変換済み）であること． 2. [vmin, vmax] で [0, 1] 正規化する． 3. cmap_name のカラーマップで RGBA → RGB uint8 (H, W, 3) 画像化する． 4. pv.Texture を生成して返す． テクスチャ画像の行列方向は irradiance_map の (row, col) と同一（PyVista が VTK の 左下原点 UV 系でサンプリングするため、UV 側で v 反転を行う）． Args: irradiance_map: 2D 照度マップ (H, W)（W/m²）． vmin: 正規化の下限（log 空間のとき log10 済み値）． vmax: 正規化の上限（log 空間のとき log10 済み値）． bake_scale: "log" または "linear"．壁ベイクと同一の変換規約を適用する． log10_eps: log 変換時のゼロ除算防止下限値（bake_scale=="log" 時のみ使用）． cmap_name: matplotlib カラーマップ名（既定 "viridis"）．壁ベイクと同一の カラーマップを使うことで評価板テクスチャと壁ベイクの配色が一致する． Returns: pv.Texture（cmap_name 着色・RGB uint8）． """ import matplotlib import pyvista as pv irr2d = irradiance_map.astype(np.float64) if bake_scale == "log": # 壁ベイクと同一の log10 変換 values = np.log10(np.maximum(irr2d, log10_eps)) else: values = irr2d # [vmin, vmax] → [0, 1] 正規化（vmin==vmax の場合はゼロ除算防止） if vmax <= vmin: normalized = np.zeros_like(values) else: normalized = np.clip((values - vmin) / (vmax - vmin), 0.0, 1.0) # cmap_name で RGBA (H, W, 4) → RGB uint8 (H, W, 3) # matplotlib 3.7+ では colormaps[] を使う（get_cmap は非推奨） colormap = matplotlib.colormaps[cmap_name] rgba = colormap(normalized) # (H, W, 4) float64 rgb = (rgba[:, :, :3] * 255).astype(np.uint8) texture = pv.Texture(rgb) _logger.info( "E_board テクスチャ生成: shape=%s bake_scale=%s cmap=%s " "E_board 値範囲=[%.4e, %.4e] 正規化 clim=[%.4f, %.4f]", rgb.shape, bake_scale, cmap_name, float(irr2d.min()), float(irr2d.max()), vmin, vmax, ) return texture def _prepare_projection_scalars( scene_config: SceneConfig, project_root: Path, luminance_map: np.ndarray, clip_low: float, clip_high: float, ) -> tuple[str | None, np.ndarray | None, np.ndarray | None, np.ndarray | None]: """投影モード用に評価板 id・per-vertex 輝度スカラ・UV 座標・カメラ向き face マスクを準備する． 評価板メッシュを特定して STL をロードし、orthographic 投影で per-vertex 輝度スカラ配列と UV テクスチャ座標を計算して返す． また、各 face がカメラ（センサ）に面しているかを判定する bool マスクも返す． Args: scene_config: 検証済み SceneConfig． project_root: STL 相対パスの基点ディレクトリ． luminance_map: sm-render 出力の 2D 輝度マップ (H, W) または (H, W, C)． clip_low: percentile クリップの下限（%）． clip_high: percentile クリップの上限（%）． Returns: (eval_board_id, projected_scalars, vertex_uv, camera_face_mask) のタプル． projected_scalars は per-vertex 輝度スカラ (N_v,)（フォールバック用）． vertex_uv は per-vertex UV テクスチャ座標 (N_v, 2) float32． camera_face_mask は per-face bool マスク (N_f,)（True = カメラ向き内側面）． 失敗時は (None, None, None, None) を返す（例外は投げない）． """ from sm_tias_lightsim.visualize.projection import ( compute_face_camera_mask, compute_orthographic_uv, find_eval_board_mesh_id, project_luminance_to_vertices, ) eval_board_id = find_eval_board_mesh_id(scene_config) _logger.info("投影モード: 評価板メッシュ id = '%s'", eval_board_id) eval_mesh_cfg = next(m for m in scene_config.geometry.meshes if m.id == eval_board_id) try: raw_path = Path(eval_mesh_cfg.path) stl_path = raw_path if raw_path.is_absolute() else project_root / eval_mesh_cfg.path eval_loaded = load_stl(stl_path, eval_mesh_cfg.unit) eval_scaled = scale_vertices(eval_loaded, target_unit="m") # 輝度マップの形状を取得（UV 算出に使用） lum_shape = (luminance_map.shape[0], luminance_map.shape[1]) projected_scalars = project_luminance_to_vertices( eval_scaled.vertices, luminance_map, scene_config.sensor, clip_low=clip_low, clip_high=clip_high, ) vertex_uv = compute_orthographic_uv( eval_scaled.vertices, lum_shape, scene_config.sensor, ) camera_face_mask = compute_face_camera_mask( eval_scaled.vertices, eval_scaled.faces, scene_config.sensor, ) _logger.info( "投影スカラ: min=%.4f, max=%.4f, std=%.4f", float(projected_scalars.min()), float(projected_scalars.max()), float(projected_scalars.std()), ) _logger.info( "UV 座標: u=[%.4f, %.4f] v=[%.4f, %.4f]", float(vertex_uv[:, 0].min()), float(vertex_uv[:, 0].max()), float(vertex_uv[:, 1].min()), float(vertex_uv[:, 1].max()), ) n_camera = int(camera_face_mask.sum()) n_total = len(camera_face_mask) _logger.info( "評価板 face 判定: カメラ向き（内側）=%d / 全 face=%d（外側 face=%d は除外）", n_camera, n_total, n_total - n_camera, ) print( f" 評価板 face: カメラ向き（内側・投影着色）={n_camera} / " f"全 face={n_total}（外側={n_total - n_camera} はニュートラル）" ) print( f" UV: u=[{vertex_uv[:, 0].min():.4f}, {vertex_uv[:, 0].max():.4f}]" f" v=[{vertex_uv[:, 1].min():.4f}, {vertex_uv[:, 1].max():.4f}]" ) except Exception as exc: _logger.warning("投影計算に失敗しました: %s．シーンレイアウトのみ表示します．", exc) return None, None, None, None return eval_board_id, projected_scalars, vertex_uv, camera_face_mask def _get_triangle_faces(mesh: pv.DataSet) -> np.ndarray | None: """PolyData または UnstructuredGrid から三角形 face 配列 (N_f, 3) を取得する． cull_outward 有効時の `extract_cells()` は UnstructuredGrid を返すため， `.faces` プロパティが存在しない．両型を統一的に扱うためのヘルパ． Args: mesh: pyvista.PolyData または pyvista.UnstructuredGrid． Returns: 三角形面配列 (N_f, 3) int64，または取得不能・三角形以外混在時は None． 三角形以外のセルが混在する場合は WARNING を出して None を返す （呼び出し側でデノイズをスキップする安全網）． """ import pyvista as pv if isinstance(mesh, pv.PolyData): # PolyData: 全面が三角形なら faces は [3, v0, v1, v2, ...] の 1D 配列 n_cells = mesh.n_cells if n_cells == 0: return None try: faces_flat = mesh.faces # 各面の先頭は頂点数カウント（三角形なら 3） # reshape(-1, 4) で [n_verts, v0, v1, v2] の 2D 配列に変換 if len(faces_flat) == n_cells * 4: faces_2d = faces_flat.reshape(-1, 4) if np.any(faces_2d[:, 0] != 3): _logger.warning( "PolyData に三角形以外の face が混在しています．デノイズをスキップします．" ) return None return faces_2d[:, 1:].astype(np.int64) else: _logger.warning( "PolyData.faces の形状が三角形メッシュと一致しません" "（len=%d, n_cells=%d）．デノイズをスキップします．", len(faces_flat), n_cells, ) return None except Exception as exc: _logger.warning( "PolyData.faces 取得に失敗しました: %s．デノイズをスキップします．", exc ) return None elif isinstance(mesh, pv.UnstructuredGrid): # UnstructuredGrid: cells_dict から VTK_TRIANGLE (celltype=5) を取得 try: cells_dict = mesh.cells_dict except Exception as exc: _logger.warning( "UnstructuredGrid.cells_dict 取得に失敗しました: %s．デノイズをスキップします．", exc, ) return None _vtk_triangle = 5 # VTK_TRIANGLE if _vtk_triangle not in cells_dict: _logger.warning( "UnstructuredGrid に VTK_TRIANGLE (type=5) が見つかりません" "（存在するセル型: %s）．デノイズをスキップします．", list(cells_dict.keys()), ) return None tri_faces = cells_dict[_vtk_triangle] # (N_f, 3) # 三角形以外が混在するセル型があれば警告 non_tri_types = [k for k in cells_dict if k != _vtk_triangle] if non_tri_types: _logger.warning( "UnstructuredGrid に三角形以外のセル型が混在しています" "（types=%s）．デノイズをスキップします．", non_tri_types, ) return None # セル数と faces 行数の整合チェック n_cells = mesh.n_cells if len(tri_faces) != n_cells: _logger.warning( "UnstructuredGrid: cells_dict の三角形数 (%d) と n_cells (%d) が不一致．" "デノイズをスキップします．", len(tri_faces), n_cells, ) return None return tri_faces.astype(np.int64) else: _logger.warning("未対応のメッシュ型 %s．デノイズをスキップします．", type(mesh).__name__) return None def denoise_face_scalars_median( faces: np.ndarray, values: np.ndarray, n_iter: int, ) -> np.ndarray: """per-face スカラに対してメディアンフィルタデノイズを適用する（ベクトル化実装）． 各 face の値を「頂点を共有する隣接 face 群（自身含む）の中央値」で置き換える操作を n_iter 回反復する．firefly（一面だけ飛び抜けた MC ノイズ）をエッジを保ったまま除去する． 実装方式（2段近似）: 1. 頂点→face の incidence を利用し「各頂点に接続する face 群の中央値」を 頂点スカラとして計算する（np.argsort ベクトル化：O(N_f log N_f)）． 2. 各 face の新しいスカラを「その 3 頂点の頂点スカラの平均値」で決定する． 厳密な「face 近傍 face の中央値」ではなく頂点経由の 2 段近似だが、 firefly 除去性能は同等でありフル Python ループ（O(N_f × 隣接数)）と比較して 数十万 face でも数秒以内に完了する． Args: faces: 三角形面配列 (N_f, 3)，頂点インデックス（非負整数）． values: per-face スカラ配列 (N_f,)，float． n_iter: デノイズ反復回数．0 の場合は入力をそのまま返す（恒等）． Returns: デノイズ後の per-face スカラ配列 (N_f,)，入力と同 dtype． Note: face 数 ~40 万・n_iter=3 で概ね 1〜3 秒（CPU，NumPy のみ，scipy 不要）． """ if n_iter <= 0: return values.copy() n_verts = int(faces.max()) + 1 out = values.astype(np.float64) for _ in range(n_iter): # ---- Step 1: 各頂点に接続する face スカラの中央値を計算する ---- # faces を (N_f, 3) → (N_f * 3,) に展開し vertex インデックスと face 値を並べる # これで "vertex → その頂点が属する face のスカラ値" のペアが得られる vert_ids = faces.ravel() # (N_f * 3,) face_vals = np.repeat(out, 3) # (N_f * 3,) — 各 face の値を 3 頂点分に複製 # vert_ids でソートして同じ頂点に属する face 値を連続させる sort_idx = np.argsort(vert_ids, kind="stable") sorted_verts = vert_ids[sort_idx] sorted_vals = face_vals[sort_idx] # 各頂点ブロックの開始インデックスを取得する（np.searchsorted で境界を特定） vert_range = np.arange(n_verts) starts = np.searchsorted(sorted_verts, vert_range, side="left") ends = np.searchsorted(sorted_verts, vert_range, side="right") # 各頂点ブロックの中央値を計算する（ブロック長が異なるためループが必要だが # 頂点数は face 数より小さく、かつ np.median はブロック内のみなので許容範囲） vert_medians = np.empty(n_verts, dtype=np.float64) for vi in range(n_verts): s, e = starts[vi], ends[vi] if s < e: vert_medians[vi] = np.median(sorted_vals[s:e]) else: # 孤立頂点（どの face にも属さない）は 0 にしておく vert_medians[vi] = 0.0 # ---- Step 2: 各 face の新スカラを 3 頂点の中央値の平均で決定する ---- # faces (N_f, 3) で vert_medians をインデックス引き → (N_f, 3) → 行平均 out = vert_medians[faces].mean(axis=1) return out.astype(values.dtype) def smooth_face_scalars_radius( face_centroids: np.ndarray, values: np.ndarray, radius: float, face_areas: np.ndarray | None = None, ) -> np.ndarray: """per-face スカラを物理半径指定のガウス重み付き平均で平滑化する． 各 face の値を，同一メッシュ内で face 重心間距離が radius 以内の face 群の 距離ガウス重み付き平均（面積重みも掛ける）で置き換える 1 パス処理． 三角形スケールの MC ノイズ（パラパラ感）を消し，radius より大きい照明構造は保つ． 実装方式（ボクセル splat/gather，純 NumPy，O(N)）: 1. ボクセルサイズ h = radius / 2 で各 face 重心を整数ボクセル座標に量子化する． 2. splat: np.add.at を用いてボクセルごとに Σ(w_area)，Σ(w_area·v)， Σ(w_area·centroid) をボクセル集計配列（n_voxels 行）に集約する． 3. LUT: ボクセル ID → 集計配列インデックスを int32 配列で O(1) 解決する （unique_ids の [min, max] 区間を配列で覆い，存在しないボクセルを -1 とする）． 4. gather: 各軸 -2..+2 の最大 125 オフセット（ユークリッド 1.5×radius 足切り後 ~117）を Python ループで回し，各回は n_voxels 要素の一括ベクトル演算で ボクセル重心間距離ガウス重み × ボクセル総面積重みの加重平均を算出する． 処理単位が n_voxels（≪ N）なので，高密度クラスタ（N ≫ n_voxels）で高速． 5. face への配布: ボクセル結果を inverse_idx で O(N) に配布する． 精度上の近似: gather 段でボクセル重心間距離をガウス重みに使うため （同一ボクセル内の全 face に同じ重みが適用される）， 厳密な per-face 距離ガウス加重と比べて ~1% 程度の誤差が生じる （mean_rel_err < 1%，max_rel_err < 5%）．照明可視化の平滑化用途では許容範囲． メモリ消費は O(N + ユニークボクセル数 + LUT span) で有界． LUT span はシーン座標範囲 / h² 程度（典型的に ~60 MB @ int32）． サブディビジョン後の高密度クラスタ（N=17 万で K=7378 になる密度）でも メモリ爆発しない（旧 cKDTree 密行列方式は 9.82 GiB でクラッシュしていた）． 計算量: O(~117 × n_voxels + N)（n_voxels ≤ N）． 性能目安: N_f ≈ 20 万・radius = 0.008 m・CPU で 3〜6 秒（密度ヘテロ耐性あり）． N_f ≈ 10 万（実際の SmTIAS シーン規模）で ~2 秒． Args: face_centroids: face 重心座標 (N_f, 3)，m 単位． values: per-face スカラ配列 (N_f,)，float． radius: 平滑化半径，m 単位．sigma = radius / 2 のガウス核を使用する． face_areas: per-face 面積 (N_f,)，m² 単位（省略可能）． 指定時はガウス重みに面積を掛けて大面積 face の寄与を重視する． None のとき面積重みなし（均等ガウス重み）． Returns: 平滑化後の per-face スカラ配列 (N_f,)，入力と同 dtype． radius <= 0 または非有限の場合は入力のコピーを返す（WARNING ログ出力）． Note: face_centroids と values の長さが一致することを前提とする （不一致時は ValueError が伝播する）． """ if not (np.isfinite(radius) and radius > 0.0): _logger.warning( "smooth_face_scalars_radius: radius=%.4e は正の有限値でないため恒等を返します．", radius, ) return values.copy() n_faces = len(face_centroids) if n_faces == 0: return values.copy() vals = values.astype(np.float64) # h = radius / 2 → sigma = h のガウス核を使用する # このサイズにより radius カットオフを覆うオフセットは各軸 -2..+2 の 5 段で済む h = radius / 2.0 neg_inv_2sigma2 = -1.0 / (2.0 * h * h) # sigma = h # face_areas が None のときは均等重み（全 face の面積 = 1.0 として扱う） w_area = face_areas.astype(np.float64) if face_areas is not None else np.ones(n_faces) # ---- Step 1: ボクセル座標への量子化 ---- # 各 face 重心を floor(centroid / h) の整数 3D 座標に変換する coords_min = face_centroids.min(axis=0) # (3,) — 座標原点シフト用 shifted = face_centroids - coords_min # (N, 3) — 非負化 voxel_ijk = np.floor(shifted / h).astype(np.int64) # (N, 3) # ---- Step 2: ボクセル ID の計算とユニーク化 ---- # 3D ボクセル座標を 1D キーに変換する（stride は軸最大値から算出） ijk_max = voxel_ijk.max(axis=0) + 1 # (3,) — 各軸のボクセル数 stride_j = int(ijk_max[2]) stride_i = int(ijk_max[1]) * stride_j voxel_ids = ( voxel_ijk[:, 0] * stride_i + voxel_ijk[:, 1] * stride_j + voxel_ijk[:, 2] ) # (N,) int64 # ユニークなボクセル ID と face→ボクセルの逆引きインデックスを取得する unique_ids, inverse_idx = np.unique(voxel_ids, return_inverse=True) n_voxels = len(unique_ids) # ---- LUT（ルックアップテーブル）の構築 ---- # ボクセル ID → 集計配列インデックスの O(1) 変換テーブル． # unique_ids の [min, max] 区間を int32 配列で覆い，存在しないボクセルを -1 とする． # メモリ: span × 4 bytes（int32）．サブディビジョン後の典型的なシーンで 50〜100 MB 程度． uid_min = int(unique_ids[0]) uid_max = int(unique_ids[-1]) lut_span = uid_max - uid_min + 1 lut = np.full(lut_span, -1, dtype=np.int32) lut[unique_ids - uid_min] = np.arange(n_voxels, dtype=np.int32) # ---- Step 3: splat（ボクセル単位の集約）---- # 各ボクセルに対して Σ(w_area)，Σ(w_area · v)，Σ(w_area · centroid) を集計する vox_w = np.zeros(n_voxels, dtype=np.float64) vox_wv = np.zeros(n_voxels, dtype=np.float64) vox_wc = np.zeros((n_voxels, 3), dtype=np.float64) np.add.at(vox_w, inverse_idx, w_area) np.add.at(vox_wv, inverse_idx, w_area * vals) np.add.at(vox_wc, inverse_idx, w_area[:, np.newaxis] * face_centroids) # ボクセル内加重重心と面積加重平均値を算出する（ゼロ除算防止） safe_w = np.where(vox_w > 0.0, vox_w, 1.0) vox_centroid = vox_wc / safe_w[:, np.newaxis] # (n_voxels, 3) vox_mean_val = vox_wv / safe_w # (n_voxels,) # ---- ボクセル ijk 座標を unique_ids から逆変換する ---- # gather ループ内でのオフセット計算に使用する vox_vi = unique_ids // stride_i # (n_voxels,) vox_vj = (unique_ids % stride_i) // stride_j vox_vk = unique_ids % stride_j # ---- Step 4: gather（ボクセル単位の近傍オフセットループ）---- # ボクセル間演算で加重平均の分子・分母を集計し，最後に face へ配布する． # gather の処理単位を N（face 数）ではなく n_voxels（ユニークボクセル数）にすることで # サブディビジョン後の高密度クラスタ（N ≫ n_voxels）で大幅に高速化できる． # 精度上の近似: face ごとの重心距離ではなくボクセル重心間距離でガウス重みを計算する # （同一ボクセル内の全 face に同じ重みが適用される）． # この近似により per-face 厳密実装と比べて ~数% の誤差が生じるが， # 照明可視化の平滑化用途では許容範囲（mean_rel_err < 1%，max_rel_err < 5% 程度）． # # 近傍オフセット集合: 各軸 -2..+2（ユークリッド距離 1.5*radius 以内に限定） # h = radius/2 なので ±2h = ±radius の範囲を覆い，対角遠方（sqrt(3)*2h≈1.73r）を足切り _cutoff_sq = (1.5 * radius) ** 2 valid_offsets: list[tuple[int, int, int]] = [] for di in range(-2, 3): for dj in range(-2, 3): for dk in range(-2, 3): if (di * h) ** 2 + (dj * h) ** 2 + (dk * h) ** 2 <= _cutoff_sq: valid_offsets.append((di, dj, dk)) vox_num = np.zeros(n_voxels, dtype=np.float64) vox_den = np.zeros(n_voxels, dtype=np.float64) for di, dj, dk in valid_offsets: # 各ボクセルに対して近傍ボクセルの ID を計算する n_ids = (vox_vi + di) * stride_i + (vox_vj + dj) * stride_j + (vox_vk + dk) # LUT で近傍ボクセルの集計インデックスを O(1) で取得する shifted_ids = n_ids - uid_min in_range = (shifted_ids >= 0) & (shifted_ids < lut_span) safe_shifted = np.where(in_range, shifted_ids, 0) idx = lut[safe_shifted] # (n_voxels,) int32，存在しないボクセルは -1 valid = in_range & (idx >= 0) if not valid.any(): continue # 近傍ボクセル重心と自ボクセル重心の距離を算出してガウス重みを計算する safe_idx = np.where(valid, idx, 0).astype(np.int64) diff = vox_centroid - vox_centroid[safe_idx] # (n_voxels, 3) dist_sq = (diff * diff).sum(axis=1) # (n_voxels,) gauss_w = np.exp(dist_sq * neg_inv_2sigma2) # (n_voxels,) # ガウス重み × 近傍ボクセル総面積重み（存在しない近傍は 0） combined_w = np.where(valid, gauss_w * vox_w[safe_idx], 0.0) # (n_voxels,) # 加重和の分子・分母に加算する vox_num += combined_w * np.where(valid, vox_mean_val[safe_idx], 0.0) vox_den += combined_w # ボクセル単位の加重平均（ゼロ除算防止: 重み総和が 0 なら face 元値を使う） vox_zero = vox_den <= 0.0 vox_out = np.where(vox_zero, vox_mean_val, vox_num / np.where(vox_zero, 1.0, vox_den)) # ボクセル結果を face に配布する（O(N)） out = vox_out[inverse_idx] _logger.debug( "smooth_face_scalars_radius (voxel): N=%d radius=%.4e h=%.4e " "n_voxels=%d lut_span=%d n_offsets=%d", n_faces, radius, h, n_voxels, lut_span, len(valid_offsets), ) return out.astype(values.dtype) def _smooth_scalar_on_polydata( polydata: pv.PolyData, scalar_name: str, n_iter: int, ) -> pv.PolyData: """メッシュ上のスカラ値を Laplacian 平滑化する（頂点座標は変えない）． cell_data → point_data 変換（隣接セル平均）を n_iter 回繰り返すことで スカラの Laplacian 平滑化に相当する効果を得る． 毎回 point_data → cell_data → point_data と往復することで拡散が進む． Args: polydata: 入力 PolyData．cell_data[scalar_name] が存在すること． scalar_name: 平滑化対象のスカラ名．cell_data に格納されている前提． n_iter: 平滑化イテレーション数．0 の場合は point_data 化のみ行う（生）． Returns: point_data[scalar_name] を持つ新しい PolyData（座標は入力と同一）． cell_data のスカラは削除される（point_data のみ残す）． """ # 入力は cell_data にスカラが格納されていると仮定する # まず cell_data → point_data に変換する（隣接セル平均） pd = polydata.cell_data_to_point_data() for _ in range(n_iter): # point_data → cell_data（隣接頂点平均）→ point_data（隣接セル平均） # この往復でスカラが隣接方向へ 1 ステップ拡散する pd = pd.point_data_to_cell_data() pd = pd.cell_data_to_point_data() return pd def show_scene( scene_config: SceneConfig, project_root: Path, irradiance: dict[str, np.ndarray] | None = None, luminance_map: np.ndarray | None = None, *, off_screen: bool = False, screenshot: Path | None = None, clip_low: float = 1.0, clip_high: float = 99.0, bake_opacity: float = 0.5, hide_below: float = 0.0, bake_denoise: int = 0, bake_smooth: int = 4, bake_scale: str = "log", cull_outward: bool = True, hide_sensor_marker: bool = False, hide_led_markers: bool = False, style: str = "default", cmap: str | None = None, subdivide_max_edge: float | None = None, bake_smooth_radius: float | None = None, backface_cull: bool = False, ) -> None: """PyVista でシーンとベイク結果を表示する． メッシュ・LED・センサ・座標軸を含む 3D インタラクティブビューアを開く． 表示モードの優先順位: 1. irradiance のみ指定（luminance_map=None）の場合（ベイクモード・既定）: **全メッシュ（評価板含む）**を per-face 照度で統一着色（viridis）． 全メッシュに bake_opacity を一律適用して内部の光分布が透けて見える． 全メッシュ共有 clim（percentile クリップ）・単一 scalar bar． 2. luminance_map と irradiance が**両方**指定された場合（ハイブリッドモード）: 評価板メッシュ: 実レンダ輝度を orthographic 投影・不透明 評価板以外の全メッシュ: irradiance per-face 照度で着色（viridis・bake_opacity 半透明） 壁ベイクと評価板 E_board の共有 clim（percentile クリップ）・単一 scalar bar depth peeling 有効化（半透明メッシュ混在のため）． 3. luminance_map のみ指定（irradiance=None）の場合（投影モード）: 評価板メッシュに 2D 輝度マップを orthographic 投影して per-vertex 着色． 他メッシュは半透明ニュートラル（opacity=0.2）． scalar bar は評価板の 1 本のみ（単位: Luminance）． bake_opacity は影響しない（投影モード専用の透明度制御がある）． 4. 両方 None の場合: シーンレイアウトのみ（--no-bake 相当）． off_screen=True のときは GUI ウィンドウを開かず，screenshot が指定されていれば保存する． screenshot=None かつ off_screen=True のときは何も出力しない（no-op）． Args: scene_config: 検証済み SceneConfig． project_root: STL 相対パスの基点ディレクトリ． irradiance: {mesh_id: irradiance_per_face (N_f,)} の辞書． None かつ luminance_map も None なら着色なし（後方互換）． luminance_map: sm-render 出力の 2D 輝度マップ (H, W) または (H, W, C)． 指定時は投影モードで評価板に投影する（irradiance より優先）． off_screen: True のとき GUI ウィンドウを開かない（screenshot 専用）． screenshot: 保存先 PNG パス．指定時は off_screen で screenshot を保存する． clip_low: percentile クリップの下限（%）．既定 1． clip_high: percentile クリップの上限（%）．既定 99． bake_opacity: ベイクモードでのメッシュ不透明度（0=透明〜1=不透明）．既定 0.5． 投影モード（luminance_map 指定時）には影響しない． hide_below: ベイク着色で irradiance がこの値以下の face を描画から除外する閾値． 既定 0.0（irradiance <= 0.0 の face を非表示＝無照度の外殻を消す）． 負値指定で実質全 face を表示（後方互換的に全面表示したい場合）． 投影モード・ベイク計算ロジックには影響しない． bake_denoise: ベイク着色のメディアンデノイズ反復回数（既定 0=なし）． N>0 のとき `denoise_face_scalars_median` を N 回反復し、MC ノイズ（firefly= 一面だけ飛び抜けた値）をエッジを保ったまま除去する． bake_scale="log" 時は log 変換後のスカラに対して適用する（log 空間でのメディアンは 単調変換に対して順序保存なので線形空間と同等の効果を持つ）． bake_scale="linear" 時は cell_data スカラ（linear 値）に直接適用する． bake_smooth（Laplacian 拡散）の前に適用される（処理順: denoise → smooth）． 負値は 0 にクランプする． 投影モード・ベイク計算ロジックには影響しない． bake_smooth: ベイク着色メッシュのスカラ Laplacian 平滑化イテレーション数． N>0 のとき cell_data_to_point_data / point_data_to_cell_data の往復を N 回繰り返し スカラを拡散させて三角形モザイクノイズを抑制する． 0 のとき平滑化なし（cell_data → point_data 変換のみ行い Gouraud 補間は有効）． 負値は 0 にクランプする． bake_scale="log" 時は log 変換後のスカラに対して平滑化する（bake_denoise の後）． 投影モード・ベイク計算ロジックには影響しない（評価板投影は常に per-vertex）． bake_scale: ベイク着色の輝度正規化スケール． "log"（既定）: log10 変換後に共有 clim / percentile を算出する． 高ダイナミックレンジ潰れを防ぎ、壁面にコントラストが付く． eps=1e-3 を下限として log10(max(irradiance, eps)) を適用． scalar bar タイトルは「Irradiance (W/m², log10)」． "linear": 従来どおり（後方互換）． scalar bar タイトルは「Irradiance (W/m²)」． 評価板の実レンダ投影（luminance_map）には影響しない． cull_outward: True（既定）のとき、ベイク着色経路の各メッシュについて 「シーン中心側を向く内向き face のみ」を描画し、外向き face はカリングする． 外向き（≈0 の外殻）を除くことで暗殻の濁り/見かけ段差を解消する （PoC poc_cull.py で検証済み）． 共有 clim もカリング後の内向き lit 値で percentile を算出する（PoC と同方式）． False のとき全面描画（後方互換）． --render の評価板テクスチャ投影には影響しない（壁ベイク側にのみ適用）． scene_center は描画対象メッシュ頂点全体の平均を 1 回算出して全メッシュ共通に使う． 巻き順依存の注意: 法線は STL の巻き順から算出するため、巻き順が 揃っていない STL では内向き/外向きが反転する可能性がある （poc_cull.py 検証 STL では妥当）． hide_sensor_marker: True のとき、センサマーカー（青い球＋青い視線矢印）の add_mesh を スキップする（デモ撮影用）．メッシュ着色・scalar bar・座標軸・bounds はそのまま残す． False（既定）のとき従来どおりセンサマーカーを表示する（後方互換）． style="slide" のとき、このフラグに関わらずセンサマーカーは非表示になる （実装: style 正規化後に `hide_sensor_marker = hide_sensor_marker or (style == "slide")` で既存フラグへ畳み込む）． hide_led_markers: True のとき、LED マーカー（PointLight 球 / AreaLight plane＋矢印 / DirectionalLight 矢印）の add_mesh をスキップする（デモ撮影用）． False（既定）のとき従来どおり LED マーカーを表示する（後方互換）． slide スタイルでも LED マーカーは表示される（--no-led-markers で個別に制御可能）． style: 表示スタイルプリセット．keyword-only． "default"（既定）: 現行表示そのまま（後方互換）．座標軸ウィジェット・ bounds グリッドを表示し、背景は PyVista 既定色、scalar bar は既定体裁． センサマーカーを表示する．cmap 未指定時は "viridis"． "slide": 発表スライド用クリーン体裁．背景を純白（white）にし、 座標軸ウィジェット（add_axes）と bounds グリッド（show_bounds）を非表示にし、 アンチエイリアスを有効化し（例外が出る環境では warning フォールバック）、 scalar bar を黒文字・フォントサイズ調整・fmt="%.2f" のクリーン体裁にする． センサマーカー（青い球＋視線矢印）を自動的に非表示にする． cmap 未指定時は "hot"（照明均一性評価と同一のカラーマップ）． 未知の style 値は "default" にフォールバックして warning する （bake_scale と同じパターン）． cmap: カラーマップ名（matplotlib 互換）．keyword-only． None（既定）: style に応じて自動決定する． "default" → "viridis"（後方互換） "slide" → "hot"（照明均一性評価と同一のカラーマップ） 文字列指定: 全スタイルでそのカラーマップを使用する．壁ベイク・評価板テクスチャ・ ダミースカラバーのすべてに同一のカラーマップが適用される（配色一致）． 不正な名前（matplotlib に存在しない）は warning を出して style 既定に フォールバックする（bake_scale / style と同じパターン）． subdivide_max_edge: ベイク着色経路の各メッシュを表示前に三角形の最大辺長が この値以下になるまで再分割する（m 単位）．keyword-only． None（既定）のとき分割なし（後方互換）． baker.py に渡した値と同じ値を指定することで，irradiance 配列の長さと PyVista PolyData の n_cells が一致する （一致しない場合は WARNING + ニュートラルグレーにフォールバック）． 投影モード（--render の評価板テクスチャ）には影響しない． bake_smooth_radius: ベイク着色を物理半径 radius（m 単位）のガウス重み付き平均で 平滑化する．keyword-only．None（既定）のとき適用しない（後方互換）． 三角形スケールの MC ノイズを消し，radius より大きい照明構造は保つ． bake_scale="log" 時は log 変換後のスカラ（log 空間）に適用する． bake_scale="linear" 時は cell_data スカラ（linear 値）に適用する． 適用タイミング: bake_denoise の後・bake_smooth（Laplacian 拡散）の前． 投影モード・ベイク計算ロジックには影響しない． backface_cull: True のとき，ベイク着色経路・ニュートラルグレー背景の add_mesh に `culling="back"` を渡し，カメラに背を向けた面を描画しない（背面カリング）． keyword-only．False（既定）のとき従来どおり全面描画（後方互換）． bake_opacity=1.0 と組み合わせると，手前の壁が自動的に隠れて内部が混色なしで 見える（dollhouse 表示・スライド/デモ向き）． マーカー（LED/センサ球・矢印）には適用しない． 注意: 面の法線は STL の巻き順依存．cull_outward で内向き face のみが残っている 場合，内向き face はシーン中心向き = 外から覗くと背面になるため期待どおりに 動く．巻き順が乱れた STL では面が部分的に欠けて見える可能性がある． """ import pyvista as pv # bake_opacity を [0, 1] にクランプ bake_opacity = max(0.0, min(1.0, bake_opacity)) # bake_denoise を 0 以上にクランプ bake_denoise = max(0, bake_denoise) # bake_smooth を 0 以上にクランプ bake_smooth = max(0, bake_smooth) # bake_scale を正規化（未知値は "log" にフォールバック） if bake_scale not in ("linear", "log"): _logger.warning("未知の bake_scale='%s'、'log' にフォールバックします．", bake_scale) bake_scale = "log" # style を正規化（未知値は "default" にフォールバック） if style not in _VALID_STYLES: _logger.warning("未知の style='%s'、'default' にフォールバックします．", style) style = "default" # cmap を解決する（None のとき style に応じて既定値を選択） _style_default_cmap = "hot" if style == "slide" else "viridis" if cmap is None: _resolved_cmap = _style_default_cmap else: # 指定された cmap 名が matplotlib に存在するか確認する import matplotlib try: matplotlib.colormaps[cmap] # 存在チェック（KeyError が来たら不正） _resolved_cmap = cmap except KeyError: _logger.warning( "未知の cmap='%s'、style 既定 '%s' にフォールバックします．", cmap, _style_default_cmap, ) _resolved_cmap = _style_default_cmap _logger.info("カラーマップ: cmap=%s（style=%s）", _resolved_cmap, style) # slide スタイルのときセンサマーカーを自動的に非表示にする # 既存フラグへ畳み込むことで hide_sensor_marker の他の制御パスと干渉しない hide_sensor_marker = hide_sensor_marker or (style == "slide") # screenshot 指定があれば off_screen を強制する if screenshot is not None: off_screen = True plotter = pv.Plotter(off_screen=off_screen) # slide スタイル: 背景純白 + アンチエイリアス（off_screen の screenshot にも反映される） if style == "slide": plotter.set_background("white") try: plotter.enable_anti_aliasing() _logger.info("アンチエイリアスを有効化しました（style=slide）．") except Exception as exc: _logger.warning("アンチエイリアスの有効化に失敗しました（フォールバック）: %s", exc) # ベイクモード / ハイブリッドモード時に depth peeling を有効化 # （半透明メッシュの深度ブレンド崩れ対策） # 例外が出る環境もあるため try/except でフォールバックする if irradiance is not None: try: plotter.enable_depth_peeling() _logger.info("depth peeling を有効化しました．") except Exception as exc: _logger.warning("depth peeling の有効化に失敗しました（フォールバック）: %s", exc) # モード判定 # use_projection: luminance_map を評価板に投影するか # use_bake_walls: irradiance で評価板以外の壁メッシュを着色するか use_projection = luminance_map is not None # ハイブリッドモード（luminance_map + irradiance 同時指定）では # 投影モードの評価板着色＋ベイクの壁着色を同時に行う use_bake_walls = irradiance is not None # 投影モードの場合：評価板メッシュ id・per-vertex スカラ・UV・カメラ向き face マスクを準備する projected_scalars: np.ndarray | None = None vertex_uv: np.ndarray | None = None camera_face_mask: np.ndarray | None = None eval_board_id: str | None = None proj_texture: pv.Texture | None = None # 評価板照度マップ E_board (W/m²)：輝度 L を E = π*L/ρ で変換したもの e_board_map: np.ndarray | None = None eval_board_rho: float | None = None if use_projection: eval_board_id, projected_scalars, vertex_uv, camera_face_mask = _prepare_projection_scalars( scene_config, project_root, luminance_map, clip_low, clip_high, # type: ignore[arg-type] ) if eval_board_id is None: # 投影計算失敗 → レイアウト表示のみにフォールバック use_projection = False else: # 評価板の反射率 ρ を scene_config から取得する # ρ が取得できない場合は 0.85（SPEC_03 基準設計値）にフォールバックする eval_board_rho = _get_eval_board_reflectance(scene_config, eval_board_id) _eval_material = next( m.material for m in scene_config.geometry.meshes if m.id == eval_board_id ) print(f" 評価板 reflectance: ρ={eval_board_rho:.4f} （material='{_eval_material}'）") # 輝度マップ L (W/(m²·sr)) → 照度マップ E_board (W/m²) # Lambertian 拡散面: L = ρ * E / π → E = π * L / ρ # （SPEC_03「筐体（材料）」の E ≈ pixel × π / 0.85 に準拠） if luminance_map is not None: if luminance_map.ndim == 3: lum2d = luminance_map[:, :, 0].astype(np.float64) else: lum2d = luminance_map.astype(np.float64) e_board_map = (np.pi * lum2d / eval_board_rho).astype(np.float64) _logger.info( "E_board 換算: ρ=%.4f L=[%.4e, %.4e] E_board=[%.4e, %.4e] mean=%.4e", eval_board_rho, float(lum2d.min()), float(lum2d.max()), float(e_board_map.min()), float(e_board_map.max()), float(e_board_map.mean()), ) print( f" E_board 換算 (E=πL/ρ): " f"min={e_board_map.min():.4e} max={e_board_map.max():.4e} " f"mean={e_board_map.mean():.4e} (W/m^2)" ) # ---- cull_outward プレパス: scene_center と内向き face マスクを事前計算 ---- # cull_outward=True のとき、ベイク着色の各メッシュについて内向き face マスクを計算し # inward_masks に保持する（clim 算出とメッシュ描画の両方で再利用する）． # scene_center は全描画対象メッシュ頂点の平均を 1 回だけ計算して共有する（PoC と同方式）． inward_masks: dict[str, np.ndarray] = {} # {mesh_id: bool (N_f,)} — cull_outward 用 # {mesh_id: (vertices_m, faces)} — プレパスで読み込んだ頂点・面のキャッシュ _prepass_vf: dict[str, tuple[np.ndarray, np.ndarray]] = {} if cull_outward and use_bake_walls and irradiance is not None: from sm_tias_lightsim.visualize.projection import compute_inward_face_mask # Step 1: 全ベイク対象メッシュの頂点を集めて scene_center を計算する _all_verts_for_center: list[np.ndarray] = [] for mesh_cfg in scene_config.geometry.meshes: if mesh_cfg.id not in irradiance: continue try: raw_path = Path(mesh_cfg.path) stl_path = raw_path if raw_path.is_absolute() else project_root / mesh_cfg.path loaded = load_stl(stl_path, mesh_cfg.unit) scaled = scale_vertices(loaded, target_unit="m") verts_for_prepass, faces_for_prepass = scaled.vertices, scaled.faces # subdivide_max_edge 指定時: viewer でも同じ分割を適用して # irradiance 配列長と n_cells を一致させる if subdivide_max_edge is not None: verts_for_prepass, faces_for_prepass = subdivide_mesh_to_max_edge( verts_for_prepass, faces_for_prepass, subdivide_max_edge ) _prepass_vf[mesh_cfg.id] = (verts_for_prepass, faces_for_prepass) _all_verts_for_center.append(verts_for_prepass) except Exception as exc: _logger.warning( "cull_outward プレパス: メッシュ '%s' 読み込みエラー（スキップ）: %s", mesh_cfg.id, exc, ) if _all_verts_for_center: _scene_center = np.vstack(_all_verts_for_center).mean(axis=0) _logger.info( "cull_outward: scene_center (m)=[%.4f, %.4f, %.4f]", float(_scene_center[0]), float(_scene_center[1]), float(_scene_center[2]), ) print( f" cull_outward: scene_center (m)=[" f"{_scene_center[0]:.4f}, {_scene_center[1]:.4f}, {_scene_center[2]:.4f}]" ) else: _scene_center = np.zeros(3) _logger.warning("cull_outward: 全頂点収集失敗、scene_center=原点にフォールバック") # Step 2: 各メッシュの内向き face マスクを計算する for mesh_cfg in scene_config.geometry.meshes: if mesh_cfg.id not in irradiance: continue if mesh_cfg.id not in _prepass_vf: continue verts_m, faces = _prepass_vf[mesh_cfg.id] try: mask = compute_inward_face_mask(verts_m, faces, _scene_center) inward_masks[mesh_cfg.id] = mask n_inward = int(mask.sum()) n_total = len(faces) _logger.info( "cull_outward: mesh='%s' inward(描画)=%d outward(カリング)=%d / total=%d", mesh_cfg.id, n_inward, n_total - n_inward, n_total, ) print( f" cull_outward: {mesh_cfg.id}:" f" inward(描画)={n_inward} outward(カリング)={n_total - n_inward}" f" / total={n_total}" ) except Exception as exc: _logger.warning( "cull_outward: mesh='%s' のマスク計算に失敗しました" "（全面描画にフォールバック）: %s", mesh_cfg.id, exc, ) # ---- 共有 clim の算出（ベイクモード / ハイブリッドモード / 投影モード単独）---- # 算出に使うデータ: # ハイブリッドモード（use_bake_walls + use_projection）: # 壁ベイク全 face irradiance ∪ E_board 全画素 → 共有 [vmin, vmax] # ベイクモード（use_bake_walls のみ）: # 壁ベイク全 face irradiance のみ # 投影モード単独（use_projection のみ、use_bake_walls=False）: # E_board のみ → 評価板用の [vmin, vmax] # cull_outward=True のとき壁ベイクの clim 算出は内向き face のみを使う（PoC と同方式）． # bake_scale="log" 時は log 変換後（log 空間）の値で算出する． shared_vmin: float | None = None shared_vmax: float | None = None bake_scalar_bar_added = False # 共有 scalar bar を 1 本だけ追加するフラグ # clim 算出用の値リスト（壁ベイク + E_board を結合する） _clim_value_chunks: list[np.ndarray] = [] if use_bake_walls and irradiance is not None: # 壁ベイク face irradiance を追加 # cull_outward=True のとき内向き face のみを使う（カリング後の lit 値で締める） for mesh_id, irr in irradiance.items(): if len(irr) == 0: continue if cull_outward and mesh_id in inward_masks: # カリング後の内向き face のみで clim を算出する（PoC と同方式） mask = inward_masks[mesh_id] if mask.sum() > 0: _clim_value_chunks.append(irr[mask].astype(np.float64)) else: _clim_value_chunks.append(irr.astype(np.float64)) if use_projection and e_board_map is not None: # 評価板の E_board（照度マップ全画素）を追加 _clim_value_chunks.append(e_board_map.ravel().astype(np.float64)) if _clim_value_chunks: all_values_combined = np.concatenate(_clim_value_chunks) if bake_scale == "log": # log 変換後の値で clim を算出する（log 空間の percentile） all_values_for_clim = np.log10(np.maximum(all_values_combined, _LOG10_EPS)) else: all_values_for_clim = all_values_combined # ---- lit フィルタ: 無照度（floor 値）の face を clim 算出から除外 ---- # log モード: log 変換後の値が log10(_LOG10_EPS) に張り付いた face は # 実質ゼロ照度（無照度 floor）であり、percentile を歪める原因になる． # log10(_LOG10_EPS) = -3.0 より僅かに大きい値のみを lit とみなし clim 算出に使う． # linear モード: linear 値が _LOG10_EPS（=1e-3 W/m²）より大きい face のみを lit とする． # 一貫性のため log モードと同じ絶対閾値を基準とする． # いずれも _LOG10_EPS を閾値の根拠とすることでハードコードを避ける． if bake_scale == "log": _lit_threshold_log = np.log10(_LOG10_EPS) _lit_mask = all_values_for_clim > _lit_threshold_log else: _lit_mask = all_values_for_clim > _LOG10_EPS _n_lit = int(_lit_mask.sum()) _n_total = len(_lit_mask) _logger.info( "lit フィルタ: lit=%d / total=%d (%.1f%%)", _n_lit, _n_total, 100.0 * _n_lit / max(_n_total, 1), ) if _n_lit >= 2: # lit 値のみで percentile を算出（無照度 floor による clim 過大を防ぐ） _values_for_percentile = all_values_for_clim[_lit_mask] _clim_source = "lit のみ" else: # lit 要素が極端に少ない場合は全結合で算出（ゼロ除算/空配列回避） _values_for_percentile = all_values_for_clim _clim_source = "全結合（lit 不足によるフォールバック）" _logger.warning( "lit face 数が不足（%d 個）のため全結合で clim を算出します（フォールバック）．", _n_lit, ) shared_vmin = float(np.percentile(_values_for_percentile, clip_low)) shared_vmax = float(np.percentile(_values_for_percentile, clip_high)) _logger.info( "共有 clim (%s 空間, %s): vmin=%.4e vmax=%.4e (percentile %.1f–%.1f)", bake_scale, _clim_source, shared_vmin, shared_vmax, clip_low, clip_high, ) print( f" 共有 clim ({bake_scale} 空間, {_clim_source}): " f"vmin={shared_vmin:.4e} vmax={shared_vmax:.4e} " f"(p{clip_low:.0f}-p{clip_high:.0f})" ) # 投影モードのテクスチャは共有 clim が確定した後に生成する if use_projection and e_board_map is not None and shared_vmin is not None: try: proj_texture = _make_irradiance_texture( e_board_map, shared_vmin, shared_vmax, # type: ignore[arg-type] bake_scale, _LOG10_EPS, cmap_name=_resolved_cmap, ) except Exception as exc: _logger.warning( "E_board テクスチャ生成に失敗しました（フォールバック: per-vertex）: %s", exc ) proj_texture = None # ---- メッシュ追加 ---- # scalar bar タイトルを bake_scale に応じて決定（投影・ベイク共通） if bake_scale == "log": _scalar_bar_title = "Irradiance (W/m^2, log10)" else: _scalar_bar_title = "Irradiance (W/m^2)" # scalar_bar_args を style に応じて 1 箇所で組み立て、全 scalar bar 追加箇所に適用する _scalar_bar_args = _build_scalar_bar_args(_scalar_bar_title, style) for mesh_cfg in scene_config.geometry.meshes: # プレパスでキャッシュ済みの場合はそれを再利用し、STL の二度読みを回避する # （プレパスでは subdivide_max_edge 分割済み頂点・面が格納されている） if mesh_cfg.id in _prepass_vf: _cached_verts, _cached_faces = _prepass_vf[mesh_cfg.id] else: try: raw_path = Path(mesh_cfg.path) stl_path = raw_path if raw_path.is_absolute() else project_root / mesh_cfg.path loaded = load_stl(stl_path, mesh_cfg.unit) scaled = scale_vertices(loaded, target_unit="m") _cached_verts, _cached_faces = scaled.vertices, scaled.faces # cull_outward プレパスに通らないメッシュ（投影専用等）にも # subdivide_max_edge が指定されていれば同じ分割を適用する if subdivide_max_edge is not None and mesh_cfg.id in (irradiance or {}): _cached_verts, _cached_faces = subdivide_mesh_to_max_edge( _cached_verts, _cached_faces, subdivide_max_edge ) except Exception as exc: _logger.warning("メッシュ '%s' を読み込めませんでした: %s", mesh_cfg.id, exc) continue polydata = _build_pyvista_polydata(_cached_verts, _cached_faces) if use_projection and mesh_cfg.id == eval_board_id and projected_scalars is not None: # ハイブリッドモード（--render 指定時のみ）: 評価板をカメラ向き内側 face のみ投影着色 # render 未指定（全面統一ベイク）時はここには来ない（use_projection=False） # カメラ向き face（内側面）と非カメラ向き face（外側面）に分割して別々に描画する if camera_face_mask is not None and len(camera_face_mask) == polydata.n_cells: # cell インデックスリストで分割する all_cell_ids = np.arange(polydata.n_cells, dtype=np.int64) inner_ids = all_cell_ids[camera_face_mask] outer_ids = all_cell_ids[~camera_face_mask] # ---- 内側 face（カメラ向き）: UV テクスチャで E_board 照度マップを貼る ---- if len(inner_ids) > 0: inner_poly = polydata.extract_cells(inner_ids) if proj_texture is not None and vertex_uv is not None: # ---- テクスチャ方式（主経路）---- # extract_cells 後の頂点インデックスを vtkOriginalPointIds で取得し、 # 元メッシュ全体の vertex_uv から対応 UV を抽出する orig_ids = inner_poly.point_data["vtkOriginalPointIds"] inner_uv = vertex_uv[orig_ids] # (N_inner_v, 2) inner_poly.point_data["TextureCoordinates"] = inner_uv inner_poly.active_texture_coordinates = inner_uv plotter.add_mesh( inner_poly, texture=proj_texture, show_scalar_bar=False, label=mesh_cfg.id, opacity=1.0, ) # カラーバーを別途追加（テクスチャは自動 colorbar を作らない） # 壁ベイクと共有の単一カラーバー（E_board + 壁 irradiance 同一スケール） # ハイブリッドモードでは壁ベイクが先に追加する場合があるため未追加のみ both_clim = shared_vmin is not None and shared_vmax is not None if not bake_scalar_bar_added and both_clim: # ダミースカラ付き不可視アクターで cmap カラーバーを 1 本追加 _dummy = inner_poly.copy() _dummy.point_data["_irr_bar"] = np.linspace( shared_vmin, shared_vmax, _dummy.n_points, dtype=np.float32 ) plotter.add_mesh( _dummy, scalars="_irr_bar", cmap=_resolved_cmap, clim=[shared_vmin, shared_vmax], show_scalar_bar=True, scalar_bar_args=_scalar_bar_args, opacity=0.0, ) bake_scalar_bar_added = True else: # ---- per-vertex スカラ方式（テクスチャ生成失敗時フォールバック）---- # E_board を per-vertex にサンプリングして着色する _logger.warning( "テクスチャが無効のため per-vertex E_board スカラに" "フォールバックします．" ) orig_ids = inner_poly.point_data["vtkOriginalPointIds"] # E_board_map から projected_scalars と同じ bilinear サンプリングで # per-vertex 値を取得する（近似: projected_scalars に π/ρ をスケール） if e_board_map is not None and eval_board_rho is not None: _e_board_scalars = projected_scalars * (np.pi / eval_board_rho) else: _e_board_scalars = projected_scalars inner_poly.point_data["irradiance"] = _e_board_scalars[orig_ids] show_bar = not bake_scalar_bar_added _kwargs_fb: dict = { "scalars": "irradiance", "cmap": _resolved_cmap, "show_scalar_bar": show_bar, "scalar_bar_args": _scalar_bar_args, "label": mesh_cfg.id, "opacity": 1.0, } if shared_vmin is not None and shared_vmax is not None: _kwargs_fb["clim"] = [shared_vmin, shared_vmax] plotter.add_mesh(inner_poly, **_kwargs_fb) if show_bar: bake_scalar_bar_added = True # ---- 外側 face（非カメラ向き）: ニュートラルグレーで描画 ---- if len(outer_ids) > 0: outer_poly = polydata.extract_cells(outer_ids) _outer_kwargs: dict = { "color": "lightgray", "opacity": _BACKGROUND_MESH_OPACITY, "label": f"{mesh_cfg.id}:outer", } if backface_cull: _outer_kwargs["culling"] = "back" plotter.add_mesh(outer_poly, **_outer_kwargs) else: # マスクが取得できなかった場合: 全面を per-vertex E_board スカラで描画 _logger.warning( "カメラ向き face マスクが無効です（全面を投影着色）． face_mask=%s, n_cells=%d", camera_face_mask is None, polydata.n_cells, ) if e_board_map is not None and eval_board_rho is not None: _e_board_scalars_full = projected_scalars * (np.pi / eval_board_rho) else: _e_board_scalars_full = projected_scalars polydata.point_data["irradiance"] = _e_board_scalars_full show_bar = not bake_scalar_bar_added _kwargs_fb2: dict = { "scalars": "irradiance", "cmap": _resolved_cmap, "show_scalar_bar": show_bar, "scalar_bar_args": _scalar_bar_args, "label": mesh_cfg.id, "opacity": 1.0, } if shared_vmin is not None and shared_vmax is not None: _kwargs_fb2["clim"] = [shared_vmin, shared_vmax] plotter.add_mesh(polydata, **_kwargs_fb2) if show_bar: bake_scalar_bar_added = True elif use_bake_walls and irradiance is not None and mesh_cfg.id in irradiance: # ベイクモード（全面統一）/ ハイブリッドモード: per-face 照度で着色 # 全面統一ベイク: 評価板を含む全メッシュがここに来る # ハイブリッドモード: 評価板以外のメッシュがここに来る # （評価板は上の if ブランチで投影着色済み or irradiance に含まれない） # bake_opacity を一律適用して内部の光分布が透けて見えるようにする irr = irradiance[mesh_cfg.id].astype(np.float32) # ---- 安全網: irradiance 長と n_cells の不一致チェック ---- # bake_cache 読み込み時に subdivide 設定が変わった場合等に発生しうる． # 不一致のメッシュはニュートラルグレーで描画して処理を続行する． if len(irr) != polydata.n_cells: _logger.warning( "メッシュ '%s': irradiance 長 (%d) と n_cells (%d) が不一致．" " subdivide 設定とキャッシュの不整合の可能性があります．" " ニュートラルグレーで描画します（irradiance スキップ）．", mesh_cfg.id, len(irr), polydata.n_cells, ) print( f" 警告: {mesh_cfg.id}: irradiance 長({len(irr)})" f" ≠ n_cells({polydata.n_cells})，" "ニュートラルグレーで描画します（subdivide 設定を確認してください）．" ) plotter.add_mesh( polydata, color="lightgray", opacity=bake_opacity, label=mesh_cfg.id, ) continue # ---- cull_outward: 内向き face のみを extract_cells して描画 ---- # 外向き face（シーン外を向く ≈0 の外殻）をカリングし濁り/暗殻を除去する． # 値閾値（hide_below）でなく幾何的向きで判定するため、 # ノイズで暗い lit 面に穴が開かない（PoC poc_cull.py で確認済み）． if cull_outward and mesh_cfg.id in inward_masks: _cull_mask = inward_masks[mesh_cfg.id] _inner_ids = np.where(_cull_mask)[0].astype(np.int64) if len(_inner_ids) == 0: _logger.info( "cull_outward: メッシュ '%s' は内向き face が 0 のため非表示．", mesh_cfg.id, ) continue polydata = polydata.extract_cells(_inner_ids) irr = irr[_cull_mask] polydata.cell_data["irradiance"] = irr # ---- 処理順序（bake_scale="log" 時）---- # 1. hide_below フィルタ（linear 値基準） # 2. log 変換（log10(max(irr, eps))） # 3. bake_smooth 平滑化（log 空間スカラに対して） # 4. cell → point 補間 → add_mesh # bake_scale="linear" 時は従来どおり（平滑化 → hide_below → cell→point） if bake_scale == "log": # ---- Step 1: hide_below フィルタ（linear 値基準、cell_data 上で実行）---- n_faces_before = polydata.n_cells if hide_below >= 0.0: filtered_cd = polydata.threshold( value=hide_below, scalars="irradiance", invert=False ) else: filtered_cd = polydata n_faces_after = filtered_cd.n_cells _logger.info( "threshold フィルタ: mesh='%s' hide_below=%.4e 残 face=%d / 全 face=%d", mesh_cfg.id, hide_below, n_faces_after, n_faces_before, ) if n_faces_after == 0: _logger.info( "メッシュ '%s' は全 face が閾値以下のため非表示にします．", mesh_cfg.id ) continue # ---- Step 2: log 変換（log 空間スカラを cell_data に上書き）---- irr_after_filter = filtered_cd.cell_data["irradiance"].astype(np.float32) log_irr = np.log10(np.maximum(irr_after_filter, _LOG10_EPS)).astype(np.float32) filtered_cd.cell_data["irradiance"] = log_irr _logger.info( "log 変換: mesh='%s' log10 irradiance: min=%.4f max=%.4f std=%.4f", mesh_cfg.id, float(log_irr.min()), float(log_irr.max()), float(log_irr.std()), ) # ---- Step 3: bake_denoise メディアンデノイズ（log 空間・cell_data）---- # cull_outward 後の polydata の faces を使う（extract_cells 後は UnstructuredGrid） if bake_denoise > 0: _faces_for_denoise = _get_triangle_faces(filtered_cd) _irr_before_denoise = filtered_cd.cell_data["irradiance"].astype(np.float64) if _faces_for_denoise is None: _logger.warning( "メディアンデノイズ (log): mesh='%s' の face 取得に失敗しました．" "デノイズをスキップします．", mesh_cfg.id, ) elif len(_faces_for_denoise) != len(_irr_before_denoise): _logger.warning( "メディアンデノイズ (log): mesh='%s' の face 数 (%d) と" " irradiance 数 (%d) が不一致．デノイズをスキップします．", mesh_cfg.id, len(_faces_for_denoise), len(_irr_before_denoise), ) else: _irr_denoised = denoise_face_scalars_median( _faces_for_denoise, _irr_before_denoise, bake_denoise ) filtered_cd.cell_data["irradiance"] = _irr_denoised.astype(np.float32) _logger.info( "メディアンデノイズ (log): mesh='%s' bake_denoise=%d " "std: %.4f → %.4f", mesh_cfg.id, bake_denoise, float(_irr_before_denoise.std()), float(_irr_denoised.std()), ) # ---- Step 4: bake_smooth_radius ガウス平滑化（log 空間・cell_data）---- # bake_denoise の後・bake_smooth（Laplacian 拡散）の前に適用する if bake_smooth_radius is not None: _irr_before_radius = filtered_cd.cell_data["irradiance"].astype(np.float64) _centroids_log = filtered_cd.cell_centers().points.astype(np.float64) _areas_log: np.ndarray | None = None try: _areas_log = ( filtered_cd.compute_cell_sizes(length=False, area=True, volume=False) .cell_data["Area"] .astype(np.float64) ) except Exception as _e: _logger.warning( "ガウス平滑化 (log): mesh='%s' の face 面積取得に失敗しました" "（面積重みなし）: %s", mesh_cfg.id, _e, ) _irr_radius_smoothed = smooth_face_scalars_radius( _centroids_log, _irr_before_radius, bake_smooth_radius, _areas_log, ) filtered_cd.cell_data["irradiance"] = _irr_radius_smoothed.astype(np.float32) _logger.info( "ガウス平滑化 (log): mesh='%s' radius=%.4f m std: %.4f → %.4f", mesh_cfg.id, bake_smooth_radius, float(_irr_before_radius.std()), float(_irr_radius_smoothed.std()), ) # ---- Step 5: bake_smooth 平滑化（log 空間スカラ・cell_data を入力）---- smoothed_pd = _smooth_scalar_on_polydata(filtered_cd, "irradiance", bake_smooth) _logger.info( "スカラ平滑化 (log): mesh='%s' bake_smooth=%d " "point_data irradiance: min=%.4f max=%.4f std=%.4f", mesh_cfg.id, bake_smooth, float(smoothed_pd.point_data["irradiance"].min()), float(smoothed_pd.point_data["irradiance"].max()), float(smoothed_pd.point_data["irradiance"].std()), ) filtered = smoothed_pd else: # linear モード: bake_denoise → bake_smooth → hide_below → cell→point # デノイズと平滑化はフルメッシュ（threshold 前）のスカラに対して行う # （threshold 後に行うと欠損セルの境界で近傍平均が歪むため） if bake_denoise > 0: _faces_lin = _get_triangle_faces(polydata) _irr_lin_before = polydata.cell_data["irradiance"].astype(np.float64) if _faces_lin is None: _logger.warning( "メディアンデノイズ (linear): mesh='%s' の face 取得に失敗しました．" "デノイズをスキップします．", mesh_cfg.id, ) elif len(_faces_lin) != len(_irr_lin_before): _logger.warning( "メディアンデノイズ (linear): mesh='%s' の face 数 (%d) と" " irradiance 数 (%d) が不一致．デノイズをスキップします．", mesh_cfg.id, len(_faces_lin), len(_irr_lin_before), ) else: _irr_lin_denoised = denoise_face_scalars_median( _faces_lin, _irr_lin_before, bake_denoise ) polydata.cell_data["irradiance"] = _irr_lin_denoised.astype(np.float32) _logger.info( "メディアンデノイズ (linear): mesh='%s' bake_denoise=%d " "std: %.4e → %.4e", mesh_cfg.id, bake_denoise, float(_irr_lin_before.std()), float(_irr_lin_denoised.std()), ) # ---- ガウス平滑化 (linear): bake_denoise の後・bake_smooth の前に適用 ---- if bake_smooth_radius is not None: _irr_before_radius_lin = polydata.cell_data["irradiance"].astype(np.float64) _centroids_lin = polydata.cell_centers().points.astype(np.float64) _areas_lin: np.ndarray | None = None try: _areas_lin = ( polydata.compute_cell_sizes(length=False, area=True, volume=False) .cell_data["Area"] .astype(np.float64) ) except Exception as _e: _logger.warning( "ガウス平滑化 (linear): mesh='%s' の face 面積取得に失敗しました" "（面積重みなし）: %s", mesh_cfg.id, _e, ) _irr_radius_smoothed_lin = smooth_face_scalars_radius( _centroids_lin, _irr_before_radius_lin, bake_smooth_radius, _areas_lin, ) polydata.cell_data["irradiance"] = _irr_radius_smoothed_lin.astype(np.float32) _logger.info( "ガウス平滑化 (linear): mesh='%s' radius=%.4e m std: %.4e → %.4e", mesh_cfg.id, bake_smooth_radius, float(_irr_before_radius_lin.std()), float(_irr_radius_smoothed_lin.std()), ) smoothed_pd = _smooth_scalar_on_polydata(polydata, "irradiance", bake_smooth) _logger.info( "スカラ平滑化: mesh='%s' bake_smooth=%d " "point_data irradiance: min=%.4e max=%.4e std=%.4e", mesh_cfg.id, bake_smooth, float(smoothed_pd.point_data["irradiance"].min()), float(smoothed_pd.point_data["irradiance"].max()), float(smoothed_pd.point_data["irradiance"].std()), ) n_faces_before = smoothed_pd.n_cells if hide_below >= 0.0: pd_for_threshold = smoothed_pd.point_data_to_cell_data() filtered_cd = pd_for_threshold.threshold( value=hide_below, scalars="irradiance", invert=False ) if filtered_cd.n_cells > 0: filtered = filtered_cd.cell_data_to_point_data() else: filtered = filtered_cd else: filtered = smoothed_pd n_faces_after = filtered.n_cells _logger.info( "threshold フィルタ: mesh='%s' hide_below=%.4e 残 face=%d / 全 face=%d", mesh_cfg.id, hide_below, n_faces_after, n_faces_before, ) if n_faces_after == 0: _logger.info( "メッシュ '%s' は全 face が閾値以下のため非表示にします．", mesh_cfg.id ) continue # _scalar_bar_title はメッシュループ前に bake_scale に応じて定義済み show_bar = not bake_scalar_bar_added add_mesh_kwargs: dict = { "scalars": "irradiance", "cmap": _resolved_cmap, "show_scalar_bar": show_bar, "scalar_bar_args": _scalar_bar_args, "label": mesh_cfg.id, "opacity": bake_opacity, # Gouraud（頂点補間）シェーディングで三角形境界の硬いエッジを除去する "smooth_shading": True, # カラーマップを頂点間で線形補間する（フラット単色塗りを防ぐ） "interpolate_before_map": True, } # 共有 clim がある場合は全壁メッシュで同一スケールを適用する if shared_vmin is not None and shared_vmax is not None: add_mesh_kwargs["clim"] = [shared_vmin, shared_vmax] # backface_cull=True のとき背面カリングを有効化する（dollhouse 表示） if backface_cull: add_mesh_kwargs["culling"] = "back" plotter.add_mesh(filtered, **add_mesh_kwargs) if show_bar: bake_scalar_bar_added = True elif use_projection: # 投影モード（ハイブリッドでない）または # ハイブリッドモードで irradiance に含まれないメッシュ: # 評価板以外は半透明ニュートラルグレー（背景扱い） _proj_bg_kwargs: dict = { "color": "lightgray", "opacity": _BACKGROUND_MESH_OPACITY, "label": mesh_cfg.id, } if backface_cull: _proj_bg_kwargs["culling"] = "back" plotter.add_mesh(polydata, **_proj_bg_kwargs) else: # 着色なし: グレー半透明（--no-bake または bake 未完了メッシュ） plotter.add_mesh( polydata, color="lightgray", opacity=_UNBAKED_MESH_OPACITY, label=mesh_cfg.id, ) # ---- シーン全体のバウンディングボックス（DirectionalLight 用）---- bbox_center, bbox_diag = _get_scene_bbox_center_and_diag(scene_config, project_root) # ---- LED / 光源表示 ---- if not hide_led_markers: for light in scene_config.lights: if isinstance(light, PointLight): pos_m = [c * _MM_TO_M for c in light.position] sphere = pv.Sphere(radius=_MARKER_SPHERE_RADIUS, center=pos_m) plotter.add_mesh(sphere, color="yellow", label=f"LED:{light.id}") elif isinstance(light, AreaLight): center_m = [c * _MM_TO_M for c in light.center] dir_vec = np.array(light.direction, dtype=np.float64) dir_norm = np.linalg.norm(dir_vec) if dir_norm > _DIR_NORM_EPSILON: dir_vec = dir_vec / dir_norm if light.shape == "disk": plane = pv.Disc( center=center_m, inner=0.0, outer=light.radius * _MM_TO_M, normal=dir_vec.tolist(), c_res=64, ) else: size_m = [light.size[0] * _MM_TO_M, light.size[1] * _MM_TO_M] plane = pv.Plane( center=center_m, direction=dir_vec.tolist(), i_size=size_m[0], j_size=size_m[1], ) plotter.add_mesh(plane, color="orange", opacity=0.7, label=f"AreaLight:{light.id}") arrow = pv.Arrow( start=center_m, direction=dir_vec.tolist(), scale=_AREA_LIGHT_ARROW_SCALE, ) plotter.add_mesh(arrow, color="orange", label=f"AreaLight:{light.id}:dir") elif isinstance(light, DirectionalLight): # バウンディングボックス中心から矢印 1 本で進行方向を示す dir_vec = np.array(light.direction, dtype=np.float64) dir_norm = np.linalg.norm(dir_vec) if dir_norm < _DIR_NORM_EPSILON: continue dir_vec = dir_vec / dir_norm arrow_scale = bbox_diag * _DIRECTIONAL_ARROW_SCALE_RATIO arrow = pv.Arrow( start=bbox_center.tolist(), direction=dir_vec.tolist(), scale=arrow_scale, ) plotter.add_mesh(arrow, color="lightyellow", label=f"DirLight:{light.id}") # ---- センサ表示 ---- sensor_pos_m = [c * _MM_TO_M for c in scene_config.sensor.position] sensor_target_m = [c * _MM_TO_M for c in scene_config.sensor.target] if not hide_sensor_marker: sensor_sphere = pv.Sphere(radius=_MARKER_SPHERE_RADIUS, center=sensor_pos_m) plotter.add_mesh(sensor_sphere, color="blue", label="sensor") sensor_pos_arr = np.array(sensor_pos_m) sensor_target_arr = np.array(sensor_target_m) view_dir = sensor_target_arr - sensor_pos_arr view_len = np.linalg.norm(view_dir) if view_len > _DIR_NORM_EPSILON: view_dir_n = view_dir / view_len sensor_arrow = pv.Arrow( start=sensor_pos_m, direction=view_dir_n.tolist(), scale=min(float(view_len), bbox_diag * _SENSOR_ARROW_SCALE_RATIO), ) plotter.add_mesh(sensor_arrow, color="blue", label="sensor:view") # ---- 座標軸と境界表示 ---- # slide スタイルではクリーンな見た目のため座標軸ウィジェットと bounds グリッドを出さない if style != "slide": plotter.add_axes() plotter.show_bounds(grid="back", location="outer") # ---- 投影モード: 評価板正面が見えるカメラ視点を設定 ---- if use_projection: # センサ位置から評価板を見る視点に設定 # センサは -Z 側から +Z 方向を向いているので、 # ビューワは少し引いて斜め上から俯瞰する位置に設定する sensor_pos = np.array(sensor_pos_m) sensor_tgt = np.array(sensor_target_m) # センサ位置からさらに -Z 方向（センサの裏側）に引いた位置 cam_offset = sensor_pos - sensor_tgt cam_dist = np.linalg.norm(cam_offset) if cam_dist < 1e-9: cam_dist = 0.1 # カメラを評価板正面かつ斜め上から cam_pos = sensor_pos + cam_offset / cam_dist * 0.25 cam_pos[1] += 0.05 # 少し上から plotter.camera.position = cam_pos.tolist() plotter.camera.focal_point = sensor_tgt.tolist() plotter.camera.up = (0.0, 1.0, 0.0) plotter.reset_camera() # ---- 表示 / 保存 ---- if screenshot is not None: plotter.show(screenshot=str(screenshot)) elif not off_screen: plotter.show() # off_screen=True かつ screenshot=None の場合は no-op