From b50e7e16be57c347ce616587df3ace26b4eb26a9 Mon Sep 17 00:00:00 2001 From: Enigmatisms <984041003@qq.com> Date: Sat, 22 Jun 2024 01:42:00 +0800 Subject: [PATCH 01/17] WIP: grid volume and RGB heterogeneous estimator --- bxdf/.gitignore | 3 + bxdf/setup.py | 27 +++++++ {parsers => bxdf}/vol_loader/vol2numpy.cpp | 10 +-- bxdf/volume.py | 89 ++++++++++++++++++++++ parsers/xml_parser.py | 7 +- 5 files changed, 130 insertions(+), 6 deletions(-) create mode 100644 bxdf/.gitignore create mode 100644 bxdf/setup.py rename {parsers => bxdf}/vol_loader/vol2numpy.cpp (91%) create mode 100644 bxdf/volume.py diff --git a/bxdf/.gitignore b/bxdf/.gitignore new file mode 100644 index 0000000..e7f4644 --- /dev/null +++ b/bxdf/.gitignore @@ -0,0 +1,3 @@ +*egg-info/ +dist/ +build/ \ No newline at end of file diff --git a/bxdf/setup.py b/bxdf/setup.py new file mode 100644 index 0000000..6bf7b69 --- /dev/null +++ b/bxdf/setup.py @@ -0,0 +1,27 @@ +from pybind11.setup_helpers import Pybind11Extension +from setuptools import setup + +__version__ = "0.1.0" +cxx_std=11 + +ext_modules = [ + Pybind11Extension("vol_loader", + ["vol_loader/vol2numpy.cpp"], + # Example: passing in the version to the compiled code + define_macros = [('VERSION_INFO', __version__)], + extra_compile_args= ['-g', '-O3'], + ), +] + +setup( + name="vol_loader", + version=__version__, + author="Qianyue He", + description="Volume grid loader", + include_dirs="/usr/include/eigen3/", + long_description="", + ext_modules=ext_modules, + extras_require={"test": "pytest"}, + zip_safe=False, + python_requires=">=3.7", +) diff --git a/parsers/vol_loader/vol2numpy.cpp b/bxdf/vol_loader/vol2numpy.cpp similarity index 91% rename from parsers/vol_loader/vol2numpy.cpp rename to bxdf/vol_loader/vol2numpy.cpp index 41f4631..d873ac6 100644 --- a/parsers/vol_loader/vol2numpy.cpp +++ b/bxdf/vol_loader/vol2numpy.cpp @@ -76,13 +76,14 @@ auto loadVol2Numpy(const std::string& filename) { readVolumeData(filename, volume); py::array_t vol_numpy(volume.size()); + float* const data_ptr = vol_numpy.mutable_data(0); if (volume.channels == 1) { for (int z = 0; z < volume.zres; ++z) { int zy_base = z * volume.yres; for (int y = 0; y < volume.yres; ++y) { int zyx_base = (zy_base + y) * volume.xres; for (int x = 0; x < volume.xres; ++x) { - vol_numpy[zyx_base + x] = volume.data[zyx_base + x]; + data_ptr[zyx_base + x] = volume.data[zyx_base + x]; } } } @@ -93,9 +94,9 @@ auto loadVol2Numpy(const std::string& filename) { int zyx_base = (zy_base + y) * volume.xres; for (int x = 0; x < volume.xres; ++x) { int index = (zyx_base + x) * 3; - vol_numpy[index] = volume.data[index]; - vol_numpy[index + 1] = volume.data[index + 1]; - vol_numpy[index + 2] = volume.data[index + 2]; + data_ptr[index] = volume.data[index]; + data_ptr[index + 1] = volume.data[index + 1]; + data_ptr[index + 2] = volume.data[index + 2]; } } } @@ -107,7 +108,6 @@ auto loadVol2Numpy(const std::string& filename) { } PYBIND11_MODULE(vol_loader, m) { - m.doc() = "Volume grid (.vol / .vdb) loader (to numpy)\n"; m.def("vol_file_to_numpy", &loadVol2Numpy, "Load volume grid from mitsuba3 .vol file (return numpy)\n" diff --git a/bxdf/volume.py b/bxdf/volume.py new file mode 100644 index 0000000..07afc0e --- /dev/null +++ b/bxdf/volume.py @@ -0,0 +1,89 @@ +""" + Grid Volume Scattering Medium (RGB / Float) + @author: Qianyue He + @date: 2024-6-22 +""" + +import numpy as np +import taichi as ti +import xml.etree.ElementTree as xet + +from taichi.math import vec3 + +from bxdf.phase import PhaseFunction +from bxdf.medium import Medium_np +from la.cam_transform import delocalize_rotate +from parsers.general_parser import get, rgb_parse +from sampler.general_sampling import random_rgb +from rich.console import Console +CONSOLE = Console(width = 128) + +try: + from vol_loader import vol_file_to_numpy +except ImportError: + CONSOLE.log("[red]:skeleton: Error [/red]: vol_loader is not found, possible reason:") + CONSOLE.log(f"module not installed. Use 'python ./setup.py install --user' in {'./bxdf/'}") + raise ImportError("vol_loader not found") + +""" TODO: add transform parsing +""" +class GridVolume_np: + __type_mapping = {"none": 0, "mono": 1, "rgb": 2} + def __init__(self, elem: xet.Element, is_world = False): + self.albedo = np.ones(3, np.float32) + self.par = np.zeros(3, np.float32) + self.pdf = np.float32([1., 0., 0.]) + self.phase_type_id = -1 + self.phase_type = "hg" + self.type_id = -1 + self.type_name = "none" + + self.xres = 0 + self.yres = 0 + self.zres = 0 + self.channel = 0 + self.grids = None + + elem_to_query = {"rgb": rgb_parse, "float": lambda el: get(el, "value"), "string": lambda path: vol_file_to_numpy(path)} + if elem is not None: + type_name = elem.get("type") + if type_name in GridVolume_np.__type_mapping: + self.type_id = GridVolume_np.__type_mapping[type_name] + else: + raise NotImplementedError(f"GridVolume type '{type_name}' is not supported.") + self.type_name = type_name + + phase_type = elem.get("phase_type") + if type_name in Medium_np.__type_mapping: + self.type_id = Medium_np.__type_mapping[phase_type] + else: + raise NotImplementedError(f"Phase function type '{phase_type}' is not supported.") + + for tag, query_func in elem_to_query.items(): + tag_elems = elem.findall(tag) + for tag_elem in tag_elems: + name = tag_elem.get("name") + if hasattr(self, name): + self.__setattr__(name, query_func(tag_elem)) + else: + if not is_world: + CONSOLE.log("[yellow]:warning: Warning: default initialization yields , which is a trivial medium.") + self.u_e = self.u_a + self.u_s + + def setup_volume(self, path:str): + self.grids, (self.xres, self.yres, self.zres, self.channel) = vol_file_to_numpy(path) + CONSOLE.log(f"Volume grid loaded, type {self.type_name},\n \ + shape: (x = {self.xres}, y = {self.yres}, z = {self.zres}, c = {self.channel})") + + def export(self): + pass + + def __repr__(self): + return f"" + +@ti.dataclass +class GridVolume: + _type: int + albedo: vec3 # scattering + \ No newline at end of file diff --git a/parsers/xml_parser.py b/parsers/xml_parser.py index 288619a..81205e7 100644 --- a/parsers/xml_parser.py +++ b/parsers/xml_parser.py @@ -86,6 +86,10 @@ def parse_emitters(em_elem: list): raise ValueError(f"Source [{emitter_type}] is not implemented.") return sources, source_id_dict +def parse_volume(directory: str, vol_list: List[xet.Element]): + + pass + def parse_wavefront( directory: str, obj_list: List[xet.Element], bsdf_dict: dict, emitter_dict: dict, texture_dict: List[Texture_np]) -> List[Arr]: @@ -255,6 +259,7 @@ def scene_parsing(directory: str, file: str): shape_nodes = root_node.findall("shape") sensor_node = root_node.find("sensor") world_node = root_node.find("world") + volume_node = root_node.find("volume") assert(sensor_node) emitter_configs, \ emitter_dict = parse_emitters(emitter_nodes) @@ -268,7 +273,7 @@ def scene_parsing(directory: str, file: str): Each mapping might needs a different packing, therefore we need different image packaging and texture info block For normal mapping, non-bidirectional renderers will be simple but not for BDPT roughness is of lower priority - - [ ] Speed up python->taichi conversion + - [x] Speed up python->taichi conversion """ array_info, all_objs, area_lut, has_vertex_normal \ = parse_wavefront(directory, shape_nodes, bsdf_dict, emitter_dict, textures) From d418813428a431a3934539806c3888f49bf93177 Mon Sep 17 00:00:00 2001 From: Enigmatisms <984041003@qq.com> Date: Sat, 22 Jun 2024 12:07:26 +0800 Subject: [PATCH 02/17] Update volume methods --- bxdf/volume.py | 71 +++++++++++++++++++++++++++++++++++++++++++++++--- 1 file changed, 67 insertions(+), 4 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index 07afc0e..9e8c0a9 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -8,7 +8,7 @@ import taichi as ti import xml.etree.ElementTree as xet -from taichi.math import vec3 +from taichi.math import vec2, vec3, mat3 from bxdf.phase import PhaseFunction from bxdf.medium import Medium_np @@ -18,6 +18,8 @@ from rich.console import Console CONSOLE = Console(width = 128) +vec3i = ti.types.vector(3, int) + try: from vol_loader import vol_file_to_numpy except ImportError: @@ -25,6 +27,8 @@ CONSOLE.log(f"module not installed. Use 'python ./setup.py install --user' in {'./bxdf/'}") raise ImportError("vol_loader not found") +__all__ = ["GridVolume_np", "GridVolume"] + """ TODO: add transform parsing """ class GridVolume_np: @@ -81,9 +85,68 @@ def export(self): def __repr__(self): return f"" - + @ti.dataclass class GridVolume: - _type: int - albedo: vec3 # scattering + """ Grid Volume Taichi End definition""" + _type: int + """ grid volume type: 0 (none), 1 (mono-chromatic), 2 (RGB) """ + albedo: vec3 + """ scattering albedo """ + inv_R: mat3 + """ Inverse rotation matrix (world frame) """ + trans: vec3 + """ Translation vector (world frame) """ + mini: vec3 + """ Min bound (AABB) minimum coordinate """ + maxi: vec3 + """ Max bound (AABB) maximum coordinate """ + min_idxs: vec3i # min_indices + """ Minimum index, usually (0, 0, 0) """ + max_idxs: vec3i # max_indices + """ Maximum index for clamping, in case there's out-of-range access """ + majorant: vec3 # max_values + """ Maximum value of sigma_t as majorant, for mono-chromatic """ + ph: PhaseFunction # phase function + """ Phase function for scattering sampling """ + + @ti.func + def intersect_volume(self, ray_o: vec3, ray_d: vec3, max_t: float) -> vec2: + """ Get ray-volume intersection AABB near and far distance """ + near_far = vec2([0, max_t]) + inv_dir = 1.0 / ray_d + + t1s = (self.mini - ray_o) * inv_dir + t2s = (self.maxi - ray_o) * inv_dir + + tmin: vec3 = t1s.min(t2s) + tmax: vec3 = t1s.max(t2s) + + near_far[0] = ti.max(ti.max(0, tmin[0]), ti.max(tmin[1], tmin[2])) + near_far[1] = ti.min(ti.min(max_t, tmax[0]), ti.min(tmax[1], tmax[2])) + return near_far + + @ti.experimental.real_func + def transmittance(self, volume: ti.template, ray_o: vec3, ray_d: vec3, max_t: float) -> vec3: + transm = vec3([1, 1, 1]) + near_far = self.intersect_volume(volume, ray_o, ray_d, max_t) + if near_far[0] >= near_far[1] or near_far[1] <= 0: + return transm + ray_o = self.inv_R @ (ray_o - self.trans) + ray_d = self.inv_R @ ray_d + if self._type: + if self._type == 1: # single channel + return self.evalTrRatioTracking(volume, ray_o, ray_d, near_far) + else: + return self.evalTrRatioTracking3D(volume, ray_o, ray_d, near_far) + return transm + + @ti.func + def evalTrRatioTracking(self, volume: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: + + pass + + @ti.func + def evalTrRatioTracking3D(self, volume: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: + pass \ No newline at end of file From 0cb32facd0531e40b6ce10cf3b78f1d9b1e48cbd Mon Sep 17 00:00:00 2001 From: Enigmatisms <984041003@qq.com> Date: Sun, 23 Jun 2024 01:01:38 +0800 Subject: [PATCH 03/17] WIP: compilation bug. mono-chromatic volume is done. --- bxdf/medium.py | 5 + bxdf/setup.py | 4 +- bxdf/vol_loader/vol2numpy.cpp | 25 +++- bxdf/volume.py | 239 ++++++++++++++++++++++++++++------ parsers/general_parser.py | 14 +- parsers/obj_loader.py | 13 +- parsers/xml_parser.py | 10 +- renderer/bdpt.py | 2 +- renderer/vpt.py | 38 +++++- scenes/cbox/cbox-volgrid.xml | 134 +++++++++++++++++++ tracer/path_tracer.py | 11 +- 11 files changed, 428 insertions(+), 67 deletions(-) create mode 100755 scenes/cbox/cbox-volgrid.xml diff --git a/bxdf/medium.py b/bxdf/medium.py index 6fee29f..d5d76ef 100644 --- a/bxdf/medium.py +++ b/bxdf/medium.py @@ -23,6 +23,11 @@ class Medium_np: __type_mapping = {"hg": 0, "multi-hg": 1, "rayleigh": 2, "mie": 3, "transparent": -1} + + @staticmethod + def is_supported_type(_type: str): + return Medium_np.__type_mapping.get(_type, None) + def __init__(self, elem: xet.Element, is_world = False): """ Without spectral information, Rayleigh scattering here might not be physically-based diff --git a/bxdf/setup.py b/bxdf/setup.py index 6bf7b69..a75d6d8 100644 --- a/bxdf/setup.py +++ b/bxdf/setup.py @@ -1,7 +1,7 @@ from pybind11.setup_helpers import Pybind11Extension from setuptools import setup -__version__ = "0.1.0" +__version__ = "0.2.0" cxx_std=11 ext_modules = [ @@ -9,7 +9,7 @@ ["vol_loader/vol2numpy.cpp"], # Example: passing in the version to the compiled code define_macros = [('VERSION_INFO', __version__)], - extra_compile_args= ['-g', '-O3'], + extra_compile_args= ['-g', '-O3', '-fopenmp'], ), ] diff --git a/bxdf/vol_loader/vol2numpy.cpp b/bxdf/vol_loader/vol2numpy.cpp index d873ac6..95f0645 100644 --- a/bxdf/vol_loader/vol2numpy.cpp +++ b/bxdf/vol_loader/vol2numpy.cpp @@ -6,6 +6,7 @@ * @date 2024-06-21 * @copyright Copyright (c) 2024 */ +#include #include #include #include @@ -71,13 +72,14 @@ bool readVolumeData(const std::string& filename, VolumeData& volume) { } // mitsuba3 vol to numpy -auto loadVol2Numpy(const std::string& filename) { +auto loadVol2Numpy(const std::string& filename, bool force_mono_color = false) { VolumeData volume; readVolumeData(filename, volume); py::array_t vol_numpy(volume.size()); float* const data_ptr = vol_numpy.mutable_data(0); if (volume.channels == 1) { + #pragma omp parallel for num_threads(4) for (int z = 0; z < volume.zres; ++z) { int zy_base = z * volume.yres; for (int y = 0; y < volume.yres; ++y) { @@ -88,15 +90,23 @@ auto loadVol2Numpy(const std::string& filename) { } } } else if (volume.channels == 3) { + #pragma omp parallel for num_threads(4) for (int z = 0; z < volume.zres; ++z) { int zy_base = z * volume.yres; for (int y = 0; y < volume.yres; ++y) { int zyx_base = (zy_base + y) * volume.xres; - for (int x = 0; x < volume.xres; ++x) { - int index = (zyx_base + x) * 3; - data_ptr[index] = volume.data[index]; - data_ptr[index + 1] = volume.data[index + 1]; - data_ptr[index + 2] = volume.data[index + 2]; + if (!force_mono_color) { + for (int x = 0; x < volume.xres; ++x) { + int index = (zyx_base + x) * 3; + data_ptr[index] = volume.data[index]; + data_ptr[index + 1] = volume.data[index + 1]; + data_ptr[index + 2] = volume.data[index + 2]; + } + } else { + for (int x = 0; x < volume.xres; ++x) { + int index = zyx_base + x; + data_ptr[index] = volume.data[index * 3]; + } } } } @@ -111,6 +121,7 @@ PYBIND11_MODULE(vol_loader, m) { m.doc() = "Volume grid (.vol / .vdb) loader (to numpy)\n"; m.def("vol_file_to_numpy", &loadVol2Numpy, "Load volume grid from mitsuba3 .vol file (return numpy)\n" - "Input: filename, input path of the file\n" + "Input: filename, [string] input path of the file\n" + "Input: force_mono_color, [bool] whether to extract only one channel from the volume data, default = False (True only for testing)\n" ); } \ No newline at end of file diff --git a/bxdf/volume.py b/bxdf/volume.py index 9e8c0a9..bfd19ea 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -4,17 +4,19 @@ @date: 2024-6-22 """ +import os import numpy as np import taichi as ti import xml.etree.ElementTree as xet -from taichi.math import vec2, vec3, mat3 - +from typing import Tuple +from taichi.math import vec2, vec3, vec4, mat3 from bxdf.phase import PhaseFunction from bxdf.medium import Medium_np from la.cam_transform import delocalize_rotate -from parsers.general_parser import get, rgb_parse +from parsers.general_parser import get, rgb_parse, transform_parse from sampler.general_sampling import random_rgb +from renderer.constants import ZERO_V3, ONES_V3 from rich.console import Console CONSOLE = Console(width = 128) @@ -29,37 +31,60 @@ __all__ = ["GridVolume_np", "GridVolume"] +FORCE_MONO_COLOR = True + """ TODO: add transform parsing """ class GridVolume_np: - __type_mapping = {"none": 0, "mono": 1, "rgb": 2} + NONE = 0 + MONO = 1 + RGB = 2 + __type_mapping = {"none": NONE, "mono": MONO, "rgb": RGB} def __init__(self, elem: xet.Element, is_world = False): self.albedo = np.ones(3, np.float32) self.par = np.zeros(3, np.float32) self.pdf = np.float32([1., 0., 0.]) self.phase_type_id = -1 self.phase_type = "hg" - self.type_id = -1 self.type_name = "none" + self.type_id = GridVolume_np.NONE self.xres = 0 self.yres = 0 self.zres = 0 self.channel = 0 - self.grids = None - elem_to_query = {"rgb": rgb_parse, "float": lambda el: get(el, "value"), "string": lambda path: vol_file_to_numpy(path)} + self.density_grid = None + + self.scale = None + self.rotation = np.eye(3, dtype = np.float32) + self.offset = np.zeros(3, dtype = np.float32) + self.forward_t = None + self.density_scaling = np.ones(3, dtype = np.float32) + + # phase function + self.par = np.zeros(3, np.float32) + self.pdf = np.float32([1., 0., 0.]) + + elem_to_query = { + "rgb": rgb_parse, + "float": lambda el: get(el, "value"), + "string": lambda el: self.setup_volume(get(el, "path", str)), + "transform": lambda el: transform_parse(el) + } if elem is not None: type_name = elem.get("type") if type_name in GridVolume_np.__type_mapping: self.type_id = GridVolume_np.__type_mapping[type_name] else: + CONSOLE.log(f"[red]Error :skull:[/red]: Volume '{elem.get('name')}' has unsupported type '{type_name}'") raise NotImplementedError(f"GridVolume type '{type_name}' is not supported.") self.type_name = type_name phase_type = elem.get("phase_type") - if type_name in Medium_np.__type_mapping: - self.type_id = Medium_np.__type_mapping[phase_type] + _phase_type_id = Medium_np.is_supported_type(phase_type) + if _phase_type_id is not None: + self.phase_type_id = _phase_type_id else: raise NotImplementedError(f"Phase function type '{phase_type}' is not supported.") @@ -71,16 +96,70 @@ def __init__(self, elem: xet.Element, is_world = False): self.__setattr__(name, query_func(tag_elem)) else: if not is_world: - CONSOLE.log("[yellow]:warning: Warning: default initialization yields , which is a trivial medium.") - self.u_e = self.u_a + self.u_s + raise ValueError("Grid volume can not be initialized by xml node {None}.") + + if self.scale is None: + self.scale = np.eye(3, dtype = np.float32) + else: + self.scale = np.diag(self.scale) + self.forward_t = self.rotation @ self.scale + if self.type_id == GridVolume_np.MONO: + self.density_grid *= self.density_scaling[0] + else: + self.density_grid *= self.density_scaling def setup_volume(self, path:str): - self.grids, (self.xres, self.yres, self.zres, self.channel) = vol_file_to_numpy(path) + if not os.path.exists(path): + CONSOLE.log(f"[red]Error :skull:[/red]: {path} contains no valid volume file.") + raise RuntimeError("Volume file not found.") + density_grid, (self.xres, self.yres, self.zres, self.channel) = vol_file_to_numpy(path, FORCE_MONO_COLOR) + if FORCE_MONO_COLOR: + CONSOLE.log(f"[yellow]Warning[/yellow]: FORCE_MONO_COLOR is True. This only makes sense when we are testing the code.") CONSOLE.log(f"Volume grid loaded, type {self.type_name},\n \ shape: (x = {self.xres}, y = {self.yres}, z = {self.zres}, c = {self.channel})") + return density_grid.reshape((self.zres, self.yres, self.xres)) + def get_shape(self) -> Tuple[int, int, int]: + return (self.zres, self.yres, self.xres) + def export(self): - pass + if self.type_id == GridVolume_np.NONE: + return GridVolume(_type = 0) + aabb_mini, aabb_maxi = self.get_aabb() + maj = self.density_grid.max(axis = (0, 1, 2)) # the shape of density grid: (zres, yres, xres, channels) + return GridVolume( + _type = self.type_id, + albedo = vec3(self.albedo), + inv_T = mat3(np.linalg.inv(self.forward_t)), + trans = vec3(self.offset), + mini = aabb_mini, + maxi = aabb_maxi, + max_idxs = vec3i([self.xres, self.yres, self.zres]) - 1, + majorant = vec3(maj) if self.type_id == GridVolume_np.RGB else maj, + ph = PhaseFunction(_type = self.phase_type_id, par = vec3(self.par), pdf = vec3(self.pdf)) + ) + + def local_to_world(self, point: np.ndarray) -> np.ndarray: + """ Take a point (shape (3) and transform it to the world space) """ + print(point.shape, self.forward_t.shape, self.offset.shape) + return point @ self.forward_t.T + self.offset + + def get_aabb(self): + x, y, z = self.xres, self.yres, self.zres + all_points = np.float32([ + [0, 0, 0], + [x, 0, 0], + [0, y, 0], + [x, y, 0], + [0, 0, z], + [x, 0, z], + [0, y, z], + [x, y, z] + ]) + world_point = self.local_to_world(all_points) + # conservative AABB + return world_point.min(axis = 0) - 1e-4, world_point.max(axis = 0) + 1e-4 + def __repr__(self): return f"= 1 + @ti.func def intersect_volume(self, ray_o: vec3, ray_d: vec3, max_t: float) -> vec2: """ Get ray-volume intersection AABB near and far distance """ @@ -119,34 +200,118 @@ def intersect_volume(self, ray_o: vec3, ray_d: vec3, max_t: float) -> vec2: t1s = (self.mini - ray_o) * inv_dir t2s = (self.maxi - ray_o) * inv_dir - tmin: vec3 = t1s.min(t2s) - tmax: vec3 = t1s.max(t2s) + tmin = ti.min(t1s, t2s) + tmax = ti.max(t1s, t2s) - near_far[0] = ti.max(ti.max(0, tmin[0]), ti.max(tmin[1], tmin[2])) - near_far[1] = ti.min(ti.min(max_t, tmax[0]), ti.min(tmax[1], tmax[2])) + near_far[0] = ti.max(0, tmin.max()) + 1e-4 + near_far[1] = ti.min(max_t, tmax.min()) - 1e-4 return near_far - @ti.experimental.real_func - def transmittance(self, volume: ti.template, ray_o: vec3, ray_d: vec3, max_t: float) -> vec3: + @ti.func + def transmittance(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: float) -> vec3: transm = vec3([1, 1, 1]) - near_far = self.intersect_volume(volume, ray_o, ray_d, max_t) - if near_far[0] >= near_far[1] or near_far[1] <= 0: - return transm - ray_o = self.inv_R @ (ray_o - self.trans) - ray_d = self.inv_R @ ray_d if self._type: - if self._type == 1: # single channel - return self.evalTrRatioTracking(volume, ray_o, ray_d, near_far) - else: - return self.evalTrRatioTracking3D(volume, ray_o, ray_d, near_far) + near_far = self.intersect_volume(ray_o, ray_d, max_t) + if near_far[0] < near_far[1] and near_far[1] > 0: + ray_o = self.inv_T @ (ray_o - self.trans) + ray_d = self.inv_T @ ray_d + if self._type: + if self._type == GridVolume_np.MONO: # single channel + transm = self.eval_tr_ratio_tracking(grid, ray_o, ray_d, near_far), + else: + transm = self.eval_tr_ratio_tracking_3d(grid, ray_o, ray_d, near_far) return transm @ti.func - def evalTrRatioTracking(self, volume: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: + def sample_mfp(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: float) -> vec4: + transm = vec4([1, 1, 1, -1]) + if self._type: + near_far = self.intersect_volume(ray_o, ray_d, max_t) + if near_far[0] < near_far[1] and near_far[1] > 0: + ray_o = self.inv_T @ (ray_o - self.trans) + ray_d = self.inv_T @ ray_d + if self._type: + if self._type == 1: # single channel + transm = self.sample_distance_delta_tracking(grid, ray_o, ray_d, near_far) + else: + transm = self.sample_distance_delta_tracking_3d(grid, ray_o, ray_d, near_far) + return transm + + @ti.func + def density_lookup_3d(self, grid: ti.template(), index: vec3, u_offset: vec3) -> float: + """ Stochastic lookup of density (mono-chromatic volume) """ + index = int(ti.floor(index + (u_offset - 0.5))) + # indexing pattern z, y, x + return ti.select((index >= 0).all() and (index <= self.max_idxs).all(), grid[index[2], index[1], index[0]], ZERO_V3) + + @ti.func + def density_lookup(self, grid: ti.template(), index: vec3, u_offset: vec3) -> float: + """ Stochastic lookup of density (RGB volume) """ + index = int(ti.floor(index + (u_offset - 0.5))) + # indexing pattern z, y, x + return ti.select((index >= 0).all() and (index <= self.max_idxs).all(), grid[index[2], index[1], index[0]], 0) + + @ti.func + def sample_distance_delta_tracking(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec4: + """ Sample distance (mono-chromatic volume) via delta tracking + Note that there is no 'sampling PDF', since we are not going to use it anyway + """ + result = vec4([1, 1, 1, -1]) + if self._type: + inv_maj = 1.0 / self.majorant[0] + + t = near_far[0] + while t < near_far[1]: + t -= ti.log(1.0 - ti.random(float)) * inv_maj + d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ + ti.random(float), ti.random(float), ti.random(float) + ])) - pass + # Scatter upon real collision + if ti.random(float) < d * inv_maj: + result[:3] = self.albedo + result[3] = t + break + return result + + @ti.func + def eval_tr_ratio_tracking(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: + inv_maj = 1.0 / self.majorant[0] + + t = near_far[0] + Tr = 1.0 + while True: + t -= ti.log(1.0 - ti.random(float)) * inv_maj + if t >= near_far[1]: break + d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ + ti.random(float), ti.random(float), ti.random(float) + ])) + + Tr *= ti.max(0, 1.0 - d * inv_maj) + # Russian Roulette + if Tr < 0.1: + if ti.random(float) >= Tr: + Tr = 0.0 + break + Tr = 1.0 + return vec3([Tr, Tr, Tr]) + + @ti.func + def sample_new_rays(self, incid: vec3): + ret_spec = vec3([1, 1, 1]) + ret_dir = incid + ret_pdf = 1.0 + if self.is_scattering(): # medium interaction - evaluate phase function (currently output a float) + local_new_dir, ret_pdf = self.ph.sample_p(incid) # local frame ray_dir should be transformed + ret_dir, _ = delocalize_rotate(incid, local_new_dir) + ret_spec *= ret_pdf + return ret_dir, ret_spec, ret_pdf + + @ti.func + def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec4: + return ZERO_V3, -1 @ti.func - def evalTrRatioTracking3D(self, volume: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: - pass + def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: + return ZERO_V3 \ No newline at end of file diff --git a/parsers/general_parser.py b/parsers/general_parser.py index 42bffc1..0df0f71 100644 --- a/parsers/general_parser.py +++ b/parsers/general_parser.py @@ -53,32 +53,34 @@ def vec3d_parse(elem: xet.Element): # for positions, implicit float -> vec3 conversion is not allowed return parse_str(elem.get("value"), no_else_branch = True) -def transform_parse(transform_elem: xet.Element) -> Tuple[Arr, Arr]: +def transform_parse(transform_elem: xet.Element) -> Tuple[Arr, Arr, Arr]: """ Note that: extrinsic rotation is not supported, meaning that we can only rotate around the object centroid, which is [intrinsic rotation]. Yet, extrinsic rotation can be composed by intrinsic rotation and translation """ - trans_r, trans_t = None, None + trans_r, trans_t, trans_s = None, None, None for child in transform_elem: if child.tag == "translate": trans_t = np.float32([get(child, "x"), get(child, "y"), get(child, "z")]) elif child.tag == "rotate": - rot_type = child.get("type") + rot_type = child.get("type", "euler") if rot_type == "euler": r_angle = get(child, "r") # roll p_angle = get(child, "p") # pitch y_angle = get(child, "y") # yaw trans_r = Rot.from_euler("zxy", (r_angle, p_angle, y_angle), degrees = True).as_matrix() elif rot_type == "quaternion": - trans_r = Rot.from_quat([get(child, "x"), get(child, "y"), get(child, "z"), get(child, "w")]) + trans_r = Rot.from_quat([get(child, "x"), get(child, "y"), get(child, "z"), get(child, "w")]).as_matrix() elif rot_type == "angle-axis": axis: Arr = np.float32([get(child, "x"), get(child, "y"), get(child, "z")]) axis /= np.linalg.norm(axis) * get(child, "angle") / 180. * np.pi - trans_r = Rot.from_rotvec(axis) + trans_r = Rot.from_rotvec(axis).as_matrix() else: raise ValueError(f"Unsupported rotation representation '{rot_type}'") + elif child.tag == "scale": + trans_s = np.float32([get(child, "x"), get(child, "y"), get(child, "z")]) elif child.tag.lower() == "lookat": target_point = parse_str(child.get("target")) origin_point = parse_str(child.get("origin")) @@ -93,7 +95,7 @@ def transform_parse(transform_elem: xet.Element) -> Tuple[Arr, Arr]: raise ValueError(f"Unsupported transformation representation '{child.tag}'") # Note that, trans_r (rotation) is defualt to be intrinsic (apply under the centroid coordinate) # Therefore, do not use trans_r unless you know how to correctly transform objects with it - return trans_r, trans_t # trans_t and trans_r could be None, if is not defined in the object + return trans_r, trans_t, trans_s # trans_t, trans_r and trans_s could be None, if is not defined in the object def parse_sphere_element(elem: xet.Element): sphere_info = np.zeros((1, 2, 3), np.float32) diff --git a/parsers/obj_loader.py b/parsers/obj_loader.py index 0776dc5..54e63bd 100644 --- a/parsers/obj_loader.py +++ b/parsers/obj_loader.py @@ -92,13 +92,24 @@ def calculate_surface_area(meshes: Arr, _type = 0): area_sum = 4. * np.pi * radius ** 2 return area_sum -def apply_transform(meshes: Arr, normals: Arr, trans_r: Arr, trans_t: Arr) -> Arr: +def is_uniform_scaling(scale: Arr) -> bool: + if scale[0] != scale[1] or scale[0] != scale[2]: + return False + return True + +def apply_transform(meshes: Arr, normals: Arr, trans_r: Arr, trans_t: Arr, trans_s: Arr) -> Arr: """ - input normals are of shape (N, 3) - input meshes are of shape (N, 3, 3), and for the last two dims - 3(front): entry index for the vertices of triangles - 3(back): entry index for (x, y, z) """ + if trans_s is not None: + if not is_uniform_scaling(trans_s): + CONSOLE.log("Warning: scaling for meshes should be uniform, otherwise normals should be re-computed.") + CONSOLE.log(f"Scaling factor is adjusted to be: {trans_s[0]}") + trans_s[1] = trans_s[0] + trans_s[2] = trans_s[0] if trans_r is not None: center = meshes.mean(axis = 1).mean(axis = 0) meshes -= center # decentralize diff --git a/parsers/xml_parser.py b/parsers/xml_parser.py index 81205e7..3bb5d47 100644 --- a/parsers/xml_parser.py +++ b/parsers/xml_parser.py @@ -115,8 +115,8 @@ def parse_wavefront( meshes, normals, vns, uvs = extract_obj_info(os.path.join(directory, filepath_child.get("value"))) transform_child = elem.find("transform") if transform_child is not None: - trans_r, trans_t = transform_parse(transform_child) - meshes, normals = apply_transform(meshes, normals, trans_r, trans_t) + trans_r, trans_t, trans_s = transform_parse(transform_child) + meshes, normals = apply_transform(meshes, normals, trans_r, trans_t, trans_s) if vns is not None: has_vertex_normal = True else: # CURRENTLY, only sphere is supported @@ -259,7 +259,10 @@ def scene_parsing(directory: str, file: str): shape_nodes = root_node.findall("shape") sensor_node = root_node.find("sensor") world_node = root_node.find("world") - volume_node = root_node.find("volume") + volume_node = root_node.findall("volume") + if len(volume_node) > 1: + CONSOLE.log(f"There are {len(volume_node)} in total, only the first volume will be kept.") + volume_node = volume_node[:1] assert(sensor_node) emitter_configs, \ emitter_dict = parse_emitters(emitter_nodes) @@ -281,6 +284,7 @@ def scene_parsing(directory: str, file: str): configs['world'] = parse_world(world_node) configs['packed_textures'] = teximgs configs['has_vertex_normal'] = has_vertex_normal + configs['volume'] = volume_node emitter_configs = update_emitter_config(emitter_configs, area_lut) return emitter_configs, array_info, all_objs, configs diff --git a/renderer/bdpt.py b/renderer/bdpt.py index 99f3931..6c1adaf 100644 --- a/renderer/bdpt.py +++ b/renderer/bdpt.py @@ -240,7 +240,7 @@ def random_walk( # Step 2: check for mean free path sampling # Calculate mfp, path_beta = transmittance / PDF - is_mi, it.min_depth, path_beta = self.sample_mfp(it.obj_id, in_free_space, it.min_depth) + is_mi, it.min_depth, path_beta = self.sample_mfp(ray_o, ray_d, it.obj_id, in_free_space, it.min_depth) if it.obj_id < 0 and not is_mi: break # exiting world bound throughput *= path_beta # attenuate first if throughput.max() < 5e-5: break diff --git a/renderer/vpt.py b/renderer/vpt.py index 348347d..c1d5151 100644 --- a/renderer/vpt.py +++ b/renderer/vpt.py @@ -13,6 +13,7 @@ from la.cam_transform import * from tracer.path_tracer import PathTracer from emitters.abtract_source import LightSource +from bxdf.volume import GridVolume_np from parsers.obj_desc import ObjDescriptor from sampler.general_sampling import balance_heuristic @@ -34,6 +35,24 @@ def __init__(self, objects: List[ObjDescriptor], prop: dict, bvh_delay: bool = False ): super().__init__(emitters, array_info, objects, prop, bvh_delay) + self.has_volume = False + self.volume = None + if prop["volume"]: + CONSOLE.log("Loading Grid Volume ...") + volume = GridVolume_np(prop["volume"][0]) + if volume.type_id == GridVolume_np.MONO: + self.density_grid = ti.field(float, volume.get_shape()) + elif volume.type_id == GridVolume_np.RGB: + self.density_grid = ti.Vector.field(3, float, volume.get_shape()) + else: + CONSOLE.log(f"Volume type '{volume.type_name}' is not supported.") + raise RuntimeError("Unsupported volume type.") + self.has_volume = True + self.density_grid.from_numpy(volume.density_grid) + self.volume = volume.export() + CONSOLE.log("Grid volume loaded") + else: + self.density_grid = ti.field(float, (1, 1, 1)) self.world_scattering = self.world.medium._type >= 0 @ti.func @@ -57,11 +76,11 @@ def non_null_surface(self, idx: int): return non_null @ti.func - def sample_mfp(self, idx: int, in_free_space: int, depth: float): + def sample_mfp(self, ray_o: vec3, ray_d: vec3, idx: int, in_free_space: int, depth: float): """ Mean free path sampling, returns: - whether medium is a valid scattering medium / mean free path """ - is_mi = False + is_mi = 0 mfp = depth beta = vec3([1., 1., 1.]) # whether the world is valid for scattering: inside the free space and world has scattering medium @@ -75,10 +94,16 @@ def sample_mfp(self, idx: int, in_free_space: int, depth: float): elif not in_free_space: is_mi, mfp, beta = self.bsdf_field[idx].medium.sample_mfp(depth) # use medium to sample / calculate transmittance + if self.has_volume: # grid volume event might override the world medium event (not physically based, but simple to implement) + result = self.volume.sample_mfp(self.density_grid, ray_o, ray_d, depth) + if result[3] > 0: # in case grid volume is nested with world medium + is_mi = 2 + mfp = result[3] + beta = result[:3] return is_mi, mfp, beta @ti.func - def track_ray(self, ray, start_p, depth): + def track_ray(self, cur_ray, cur_point, depth): """ For medium interaction, check if the path to one point is not blocked (by non-null surface object) And also we need to calculate the attenuation along the path, e.g.: if the ray passes through @@ -86,9 +111,10 @@ def track_ray(self, ray, start_p, depth): FIXME: the speed of this method should be boosted """ tr = vec3([1., 1., 1.]) + if self.has_volume: # grid volume transmittance + tr *= self.volume.transmittance(self.density_grid, cur_point, cur_ray, depth) + in_free_space = True - cur_point = start_p - cur_ray = ray acc_depth = 0.0 for _i in range(7): # maximum tracking depth = 7 (corresponding to at most 2 clouds of smoke) it = self.ray_intersect(cur_ray, cur_point, depth) @@ -154,7 +180,7 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int in_free_space = tm.dot(it.n_g, ray_d) < 0 # Step 3: check for mean free path sampling # Calculate mfp, path_beta = transmittance / PDF - is_mi, it.min_depth, path_beta = self.sample_mfp(it.obj_id, in_free_space, it.min_depth) + is_mi, it.min_depth, path_beta = self.sample_mfp(ray_o, ray_d, it.obj_id, in_free_space, it.min_depth) if it.obj_id < 0 and not is_mi: break # exiting world bound hit_point = ray_d * it.min_depth + ray_o throughput *= path_beta # attenuate first diff --git a/scenes/cbox/cbox-volgrid.xml b/scenes/cbox/cbox-volgrid.xml new file mode 100755 index 0000000..94cb7b8 --- /dev/null +++ b/scenes/cbox/cbox-volgrid.xml @@ -0,0 +1,134 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + \ No newline at end of file diff --git a/tracer/path_tracer.py b/tracer/path_tracer.py index da2168f..db05d63 100644 --- a/tracer/path_tracer.py +++ b/tracer/path_tracer.py @@ -433,10 +433,13 @@ def sample_new_ray(self, ret_pdf = 1.0 is_specular = False if is_mi: - if in_free_space: # sample world medium - ret_dir, ret_spec, ret_pdf = self.world.medium.sample_new_rays(incid) - else: # sample object medium - ret_dir, ret_spec, ret_pdf = self.bsdf_field[it.obj_id].medium.sample_new_rays(incid) + if is_mi > 1: # well this is..., I know it is ugly but... + ret_dir, ret_spec, ret_pdf = self.volume.sample_new_rays(incid) + else: + if in_free_space: # sample world medium + ret_dir, ret_spec, ret_pdf = self.world.medium.sample_new_rays(incid) + else: # sample object medium + ret_dir, ret_spec, ret_pdf = self.bsdf_field[it.obj_id].medium.sample_new_rays(incid) else: # surface sampling if ti.is_active(self.obj_nodes, it.obj_id): # active means the object is attached to BRDF if ti.static(self.brdf_two_sides): From fc12b6da4b3d96d5c93cde17bcb71e54da8135a3 Mon Sep 17 00:00:00 2001 From: Enigmatisms <984041003@qq.com> Date: Sun, 23 Jun 2024 23:37:21 +0800 Subject: [PATCH 04/17] Single channel grid volume rendering completed. --- .vscode/settings.json | 3 +- bxdf/vol_loader/vol2numpy.cpp | 12 ++-- bxdf/volume.py | 80 ++++++++++++---------- renderer/vpt.py | 10 +-- scenes/cbox/cbox-volgrid.xml | 80 +++++++++++++++++----- scenes/meshes/cornell/cbox_blue_light.obj | 15 ++++ scenes/meshes/cornell/cbox_green_light.obj | 15 ++++ scenes/meshes/cornell/cbox_red_light.obj | 15 ++++ 8 files changed, 168 insertions(+), 62 deletions(-) create mode 100644 scenes/meshes/cornell/cbox_blue_light.obj create mode 100644 scenes/meshes/cornell/cbox_green_light.obj create mode 100644 scenes/meshes/cornell/cbox_red_light.obj diff --git a/.vscode/settings.json b/.vscode/settings.json index 3ac7a7e..b8b8cc8 100644 --- a/.vscode/settings.json +++ b/.vscode/settings.json @@ -71,7 +71,8 @@ "xstring": "cpp", "xutility": "cpp", "xlocnum": "cpp", - "xtr1common": "cpp" + "xtr1common": "cpp", + "regex": "cpp" }, "python.analysis.diagnosticSeverityOverrides": { "reportInvalidTypeForm": "none" diff --git a/bxdf/vol_loader/vol2numpy.cpp b/bxdf/vol_loader/vol2numpy.cpp index 95f0645..136f82a 100644 --- a/bxdf/vol_loader/vol2numpy.cpp +++ b/bxdf/vol_loader/vol2numpy.cpp @@ -21,8 +21,9 @@ struct VolumeData { int channels; std::vector data; - inline size_t size() const noexcept { - return xres * yres * zres; + inline size_t size(bool force_mono_color = false) const noexcept { + size_t ch = force_mono_color ? 1 : channels; + return ch * xres * yres * zres; } inline auto shape() const noexcept { @@ -76,7 +77,7 @@ auto loadVol2Numpy(const std::string& filename, bool force_mono_color = false) { VolumeData volume; readVolumeData(filename, volume); - py::array_t vol_numpy(volume.size()); + py::array_t vol_numpy(volume.size(force_mono_color)); float* const data_ptr = vol_numpy.mutable_data(0); if (volume.channels == 1) { #pragma omp parallel for num_threads(4) @@ -105,7 +106,7 @@ auto loadVol2Numpy(const std::string& filename, bool force_mono_color = false) { } else { for (int x = 0; x < volume.xres; ++x) { int index = zyx_base + x; - data_ptr[index] = volume.data[index * 3]; + data_ptr[index] = volume.data[index * 3 + 1]; } } } @@ -113,6 +114,9 @@ auto loadVol2Numpy(const std::string& filename, bool force_mono_color = false) { } else { std::cerr << "Grid channel: <" << volume.channels << "> is not supported, supported channels: [1, 3]" << std::endl; } + + if (force_mono_color) + volume.channels = 1; return std::tuple(vol_numpy, volume.shape()); } diff --git a/bxdf/volume.py b/bxdf/volume.py index bfd19ea..870f77a 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -40,7 +40,7 @@ class GridVolume_np: MONO = 1 RGB = 2 __type_mapping = {"none": NONE, "mono": MONO, "rgb": RGB} - def __init__(self, elem: xet.Element, is_world = False): + def __init__(self, elem: xet.Element): self.albedo = np.ones(3, np.float32) self.par = np.zeros(3, np.float32) self.pdf = np.float32([1., 0., 0.]) @@ -56,6 +56,7 @@ def __init__(self, elem: xet.Element, is_world = False): self.density_grid = None + self.toWorld = None self.scale = None self.rotation = np.eye(3, dtype = np.float32) self.offset = np.zeros(3, dtype = np.float32) @@ -70,7 +71,7 @@ def __init__(self, elem: xet.Element, is_world = False): "rgb": rgb_parse, "float": lambda el: get(el, "value"), "string": lambda el: self.setup_volume(get(el, "path", str)), - "transform": lambda el: transform_parse(el) + "transform": lambda el: self.assign_transform(*transform_parse(el)) } if elem is not None: type_name = elem.get("type") @@ -94,19 +95,21 @@ def __init__(self, elem: xet.Element, is_world = False): name = tag_elem.get("name") if hasattr(self, name): self.__setattr__(name, query_func(tag_elem)) - else: - if not is_world: - raise ValueError("Grid volume can not be initialized by xml node {None}.") - - if self.scale is None: - self.scale = np.eye(3, dtype = np.float32) - else: - self.scale = np.diag(self.scale) - self.forward_t = self.rotation @ self.scale - if self.type_id == GridVolume_np.MONO: - self.density_grid *= self.density_scaling[0] - else: - self.density_grid *= self.density_scaling + + if self.scale is None: + self.scale = np.eye(3, dtype = np.float32) + else: + self.scale = np.diag(self.scale) + self.forward_t = self.rotation @ self.scale + if self.type_id == GridVolume_np.MONO: + self.density_grid *= self.density_scaling[0] + else: + self.density_grid *= self.density_scaling + + def assign_transform(self, trans_r, trans_t, trans_s): + self.rotation = trans_r + self.offset = trans_t + self.scale = trans_s def setup_volume(self, path:str): if not os.path.exists(path): @@ -115,9 +118,7 @@ def setup_volume(self, path:str): density_grid, (self.xres, self.yres, self.zres, self.channel) = vol_file_to_numpy(path, FORCE_MONO_COLOR) if FORCE_MONO_COLOR: CONSOLE.log(f"[yellow]Warning[/yellow]: FORCE_MONO_COLOR is True. This only makes sense when we are testing the code.") - CONSOLE.log(f"Volume grid loaded, type {self.type_name},\n \ - shape: (x = {self.xres}, y = {self.yres}, z = {self.zres}, c = {self.channel})") - return density_grid.reshape((self.zres, self.yres, self.xres)) + return density_grid.reshape((self.zres, self.yres, self.xres, self.channel)) def get_shape(self) -> Tuple[int, int, int]: return (self.zres, self.yres, self.xres) @@ -132,16 +133,15 @@ def export(self): albedo = vec3(self.albedo), inv_T = mat3(np.linalg.inv(self.forward_t)), trans = vec3(self.offset), - mini = aabb_mini, - maxi = aabb_maxi, - max_idxs = vec3i([self.xres, self.yres, self.zres]) - 1, - majorant = vec3(maj) if self.type_id == GridVolume_np.RGB else maj, + mini = vec3(aabb_mini), + maxi = vec3(aabb_maxi), + max_idxs = vec3i([self.xres, self.yres, self.zres]), + majorant = vec3(maj) if self.type_id == GridVolume_np.RGB else vec3([maj, maj, maj]), ph = PhaseFunction(_type = self.phase_type_id, par = vec3(self.par), pdf = vec3(self.pdf)) ) def local_to_world(self, point: np.ndarray) -> np.ndarray: """ Take a point (shape (3) and transform it to the world space) """ - print(point.shape, self.forward_t.shape, self.offset.shape) return point @ self.forward_t.T + self.offset def get_aabb(self): @@ -158,12 +158,14 @@ def get_aabb(self): ]) world_point = self.local_to_world(all_points) # conservative AABB - return world_point.min(axis = 0) - 1e-4, world_point.max(axis = 0) + 1e-4 + return world_point.min(axis = 0) - 0.1, world_point.max(axis = 0) + 0.1 def __repr__(self): + aabb_min, aabb_max = self.get_aabb() return f"" + shape: (x = {self.xres}, y = {self.yres}, z = {self.zres}, c = {self.channel}),\n \ + AABB: Min = {aabb_min}, Max = {aabb_max}>" @ti.dataclass class GridVolume: @@ -181,7 +183,7 @@ class GridVolume: maxi: vec3 """ Max bound (AABB) maximum coordinate """ max_idxs: vec3i # max_indices - """ Maximum index for clamping, in case there's out-of-range access """ + """ Maximum index for clamping, in case there's out-of-range access (xres - 1, yres - 1, zres - 1)""" majorant: vec3 # max_values """ (Scaled) Maximum value of sigma_t as majorant, for mono-chromatic, [0] is the value """ ph: PhaseFunction # phase function @@ -217,7 +219,7 @@ def transmittance(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: fl ray_d = self.inv_T @ ray_d if self._type: if self._type == GridVolume_np.MONO: # single channel - transm = self.eval_tr_ratio_tracking(grid, ray_o, ray_d, near_far), + transm = self.eval_tr_ratio_tracking(grid, ray_o, ray_d, near_far) else: transm = self.eval_tr_ratio_tracking_3d(grid, ray_o, ray_d, near_far) return transm @@ -238,18 +240,24 @@ def sample_mfp(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: float return transm @ti.func - def density_lookup_3d(self, grid: ti.template(), index: vec3, u_offset: vec3) -> float: + def density_lookup_3d(self, grid: ti.template(), index: vec3, u_offset: vec3) -> vec3: """ Stochastic lookup of density (mono-chromatic volume) """ - index = int(ti.floor(index + (u_offset - 0.5))) + idx = ti.cast(ti.floor(index + (u_offset - 0.5)), int) + val = ZERO_V3 + if (idx >= 0).all() and (idx <= self.max_idxs).all(): + val = grid[idx[2], idx[1], idx[0]] # indexing pattern z, y, x - return ti.select((index >= 0).all() and (index <= self.max_idxs).all(), grid[index[2], index[1], index[0]], ZERO_V3) + return val @ti.func - def density_lookup(self, grid: ti.template(), index: vec3, u_offset: vec3) -> float: + def density_lookup(self, grid: ti.template(), index, u_offset: vec3) -> float: """ Stochastic lookup of density (RGB volume) """ - index = int(ti.floor(index + (u_offset - 0.5))) + idx = ti.cast(ti.floor(index + (u_offset - 0.5)), int) + val = 0.0 + if (idx >= 0).all() and (idx <= self.max_idxs).all(): + val = grid[idx[2], idx[1], idx[0]] # indexing pattern z, y, x - return ti.select((index >= 0).all() and (index <= self.max_idxs).all(), grid[index[2], index[1], index[0]], 0) + return val @ti.func def sample_distance_delta_tracking(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec4: @@ -266,7 +274,6 @@ def sample_distance_delta_tracking(self, grid: ti.template(), ray_ol: vec3, ray_ d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) - # Scatter upon real collision if ti.random(float) < d * inv_maj: result[:3] = self.albedo @@ -286,7 +293,6 @@ def eval_tr_ratio_tracking(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3 d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) - Tr *= ti.max(0, 1.0 - d * inv_maj) # Russian Roulette if Tr < 0.1: @@ -309,9 +315,9 @@ def sample_new_rays(self, incid: vec3): @ti.func def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec4: - return ZERO_V3, -1 + return vec4([1, 1, 1, -1]) @ti.func def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: - return ZERO_V3 + return ONES_V3 \ No newline at end of file diff --git a/renderer/vpt.py b/renderer/vpt.py index c1d5151..4483c17 100644 --- a/renderer/vpt.py +++ b/renderer/vpt.py @@ -42,6 +42,7 @@ def __init__(self, volume = GridVolume_np(prop["volume"][0]) if volume.type_id == GridVolume_np.MONO: self.density_grid = ti.field(float, volume.get_shape()) + volume.density_grid = volume.density_grid.squeeze() elif volume.type_id == GridVolume_np.RGB: self.density_grid = ti.Vector.field(3, float, volume.get_shape()) else: @@ -50,8 +51,9 @@ def __init__(self, self.has_volume = True self.density_grid.from_numpy(volume.density_grid) self.volume = volume.export() - CONSOLE.log("Grid volume loaded") + CONSOLE.log(f"Grid volume loaded: {volume}") else: + self.volume = GridVolume_np(None).export() self.density_grid = ti.field(float, (1, 1, 1)) self.world_scattering = self.world.medium._type >= 0 @@ -94,7 +96,7 @@ def sample_mfp(self, ray_o: vec3, ray_d: vec3, idx: int, in_free_space: int, dep elif not in_free_space: is_mi, mfp, beta = self.bsdf_field[idx].medium.sample_mfp(depth) # use medium to sample / calculate transmittance - if self.has_volume: # grid volume event might override the world medium event (not physically based, but simple to implement) + if ti.static(self.has_volume): # grid volume event might override the world medium event (not physically based, but simple to implement) result = self.volume.sample_mfp(self.density_grid, ray_o, ray_d, depth) if result[3] > 0: # in case grid volume is nested with world medium is_mi = 2 @@ -111,7 +113,7 @@ def track_ray(self, cur_ray, cur_point, depth): FIXME: the speed of this method should be boosted """ tr = vec3([1., 1., 1.]) - if self.has_volume: # grid volume transmittance + if ti.static(self.has_volume): # grid volume transmittance tr *= self.volume.transmittance(self.density_grid, cur_point, cur_ray, depth) in_free_space = True @@ -171,7 +173,7 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int # Step 2: ray intersection it = self.ray_intersect(ray_d, ray_o) if it.obj_id < 0: - if not self.world_scattering: break # nothing is hit, break + if ti.static(not self.world_scattering and not self.has_volume): break # nothing is hit, break else: # the world is filled with scattering medium it.min_depth = self.world_bound_time(ray_o, ray_d) in_free_space = True diff --git a/scenes/cbox/cbox-volgrid.xml b/scenes/cbox/cbox-volgrid.xml index 94cb7b8..77b0310 100755 --- a/scenes/cbox/cbox-volgrid.xml +++ b/scenes/cbox/cbox-volgrid.xml @@ -9,7 +9,7 @@ - + @@ -17,6 +17,8 @@ + + @@ -28,8 +30,8 @@ - - + + @@ -46,16 +48,16 @@ - + - - + + - + - - + + @@ -66,16 +68,49 @@ - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + - - @@ -108,16 +143,29 @@ + + - - + + - + - + diff --git a/scenes/meshes/cornell/cbox_blue_light.obj b/scenes/meshes/cornell/cbox_blue_light.obj new file mode 100644 index 0000000..276088e --- /dev/null +++ b/scenes/meshes/cornell/cbox_blue_light.obj @@ -0,0 +1,15 @@ +# Blender v3.0.1 OBJ File: '' +# www.blender.org +o Plane +v 0.042202 1.211291 5.535800 +v 5.470666 1.211291 5.535800 +v 0.042202 0.884762 5.535800 +v 5.470666 0.884762 5.535800 +vt 1.000000 0.000000 +vt 0.000000 1.000000 +vt 0.000000 0.000000 +vt 1.000000 1.000000 +vn 0.0000 0.0000 -1.0000 +s off +f 2/1/1 3/2/1 1/3/1 +f 2/1/1 4/4/1 3/2/1 diff --git a/scenes/meshes/cornell/cbox_green_light.obj b/scenes/meshes/cornell/cbox_green_light.obj new file mode 100644 index 0000000..184111f --- /dev/null +++ b/scenes/meshes/cornell/cbox_green_light.obj @@ -0,0 +1,15 @@ +# Blender v3.0.1 OBJ File: '' +# www.blender.org +o Plane.001 +v 0.042202 2.768821 5.535800 +v 5.470666 2.768821 5.535800 +v 0.042202 2.442292 5.535800 +v 5.470666 2.442292 5.535800 +vt 1.000000 0.000000 +vt 0.000000 1.000000 +vt 0.000000 0.000000 +vt 1.000000 1.000000 +vn 0.0000 0.0000 -1.0000 +s off +f 2/1/1 3/2/1 1/3/1 +f 2/1/1 4/4/1 3/2/1 diff --git a/scenes/meshes/cornell/cbox_red_light.obj b/scenes/meshes/cornell/cbox_red_light.obj new file mode 100644 index 0000000..67c8223 --- /dev/null +++ b/scenes/meshes/cornell/cbox_red_light.obj @@ -0,0 +1,15 @@ +# Blender v3.0.1 OBJ File: '' +# www.blender.org +o Plane.002 +v 0.042202 4.489385 5.535800 +v 5.470666 4.489385 5.535800 +v 0.042202 4.162856 5.535800 +v 5.470666 4.162856 5.535800 +vt 1.000000 0.000000 +vt 0.000000 1.000000 +vt 0.000000 0.000000 +vt 1.000000 1.000000 +vn 0.0000 0.0000 -1.0000 +s off +f 2/1/1 3/2/1 1/3/1 +f 2/1/1 4/4/1 3/2/1 From 2e9430b5e0f775f224d9aa63d332dbf4d65534fe Mon Sep 17 00:00:00 2001 From: Enigmatisms <984041003@qq.com> Date: Mon, 24 Jun 2024 00:28:31 +0800 Subject: [PATCH 05/17] WIP: RGB volume --- bxdf/volume.py | 27 +++++++++++++++++++++++++-- 1 file changed, 25 insertions(+), 2 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index 870f77a..a260c64 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -186,6 +186,8 @@ class GridVolume: """ Maximum index for clamping, in case there's out-of-range access (xres - 1, yres - 1, zres - 1)""" majorant: vec3 # max_values """ (Scaled) Maximum value of sigma_t as majorant, for mono-chromatic, [0] is the value """ + pdf: vec3 # normalized majorant + """ PDF used for MIS (wavelength spectrum) """ ph: PhaseFunction # phase function """ Phase function for scattering sampling """ @@ -246,7 +248,6 @@ def density_lookup_3d(self, grid: ti.template(), index: vec3, u_offset: vec3) -> val = ZERO_V3 if (idx >= 0).all() and (idx <= self.max_idxs).all(): val = grid[idx[2], idx[1], idx[0]] - # indexing pattern z, y, x return val @ti.func @@ -315,8 +316,30 @@ def sample_new_rays(self, incid: vec3): @ti.func def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec4: - return vec4([1, 1, 1, -1]) + result = vec4([1, 1, 1, -1]) + if self._type: + # First, choose one element according to the majorant + inv_maj = 1.0 / ti.max(self.majorant[2], 1e-4) + val = ti.random(float) + if val < self.pdf[0]: + inv_maj = 1.0 / self.majorant[0] + elif val < self.pdf[0] + self.pdf[1]: + inv_maj = 1.0 / self.majorant[1] + + t = near_far[0] + while t < near_far[1]: + t -= ti.log(1.0 - ti.random(float)) * inv_maj + d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ + ti.random(float), ti.random(float), ti.random(float) + ])) + # Scatter upon real collision + if ti.random(float) < d * inv_maj: + result[:3] = self.albedo + result[3] = t + break + return result + @ti.func def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: return ONES_V3 From 865ba112c93d524f2981a6d44216e59d88c24294 Mon Sep 17 00:00:00 2001 From: enigmahe Date: Mon, 24 Jun 2024 16:23:01 +0800 Subject: [PATCH 06/17] Updated windows AdaPT vscode config --- .vscode/settings.json | 211 +++++++++++++++++++++++++++--------------- 1 file changed, 134 insertions(+), 77 deletions(-) diff --git a/.vscode/settings.json b/.vscode/settings.json index 3ac7a7e..db0591b 100644 --- a/.vscode/settings.json +++ b/.vscode/settings.json @@ -1,79 +1,136 @@ { - "files.associations": { - "cctype": "cpp", - "clocale": "cpp", - "cmath": "cpp", - "cstdarg": "cpp", - "cstddef": "cpp", - "cstdio": "cpp", - "cstdlib": "cpp", - "cstring": "cpp", - "ctime": "cpp", - "cwchar": "cpp", - "cwctype": "cpp", - "array": "cpp", - "atomic": "cpp", - "strstream": "cpp", - "*.tcc": "cpp", - "bitset": "cpp", - "chrono": "cpp", - "cinttypes": "cpp", - "codecvt": "cpp", - "complex": "cpp", - "condition_variable": "cpp", - "cstdint": "cpp", - "deque": "cpp", - "forward_list": "cpp", - "list": "cpp", - "unordered_map": "cpp", - "unordered_set": "cpp", - "vector": "cpp", - "exception": "cpp", - "algorithm": "cpp", - "filesystem": "cpp", - "functional": "cpp", - "iterator": "cpp", - "map": "cpp", - "memory": "cpp", - "memory_resource": "cpp", - "numeric": "cpp", - "optional": "cpp", - "random": "cpp", - "ratio": "cpp", - "set": "cpp", - "string": "cpp", - "string_view": "cpp", - "system_error": "cpp", - "tuple": "cpp", - "type_traits": "cpp", - "utility": "cpp", - "fstream": "cpp", - "initializer_list": "cpp", - "iomanip": "cpp", - "iosfwd": "cpp", - "iostream": "cpp", - "istream": "cpp", - "limits": "cpp", - "mutex": "cpp", - "new": "cpp", - "ostream": "cpp", - "sstream": "cpp", - "stdexcept": "cpp", - "streambuf": "cpp", - "thread": "cpp", - "cfenv": "cpp", - "typeindex": "cpp", - "typeinfo": "cpp", - "valarray": "cpp", - "variant": "cpp", - "bit": "cpp", - "stack": "cpp", - "xstring": "cpp", - "xutility": "cpp", - "xlocnum": "cpp", - "xtr1common": "cpp" - }, - "python.analysis.diagnosticSeverityOverrides": { - "reportInvalidTypeForm": "none" - } + "files.associations": { + "cctype": "cpp", + "clocale": "cpp", + "cmath": "cpp", + "cstdarg": "cpp", + "cstddef": "cpp", + "cstdio": "cpp", + "cstdlib": "cpp", + "cstring": "cpp", + "ctime": "cpp", + "cwchar": "cpp", + "cwctype": "cpp", + "array": "cpp", + "atomic": "cpp", + "strstream": "cpp", + "*.tcc": "cpp", + "bitset": "cpp", + "chrono": "cpp", + "cinttypes": "cpp", + "codecvt": "cpp", + "complex": "cpp", + "condition_variable": "cpp", + "cstdint": "cpp", + "deque": "cpp", + "forward_list": "cpp", + "list": "cpp", + "unordered_map": "cpp", + "unordered_set": "cpp", + "vector": "cpp", + "exception": "cpp", + "algorithm": "cpp", + "filesystem": "cpp", + "functional": "cpp", + "iterator": "cpp", + "map": "cpp", + "memory": "cpp", + "memory_resource": "cpp", + "numeric": "cpp", + "optional": "cpp", + "random": "cpp", + "ratio": "cpp", + "set": "cpp", + "string": "cpp", + "string_view": "cpp", + "system_error": "cpp", + "tuple": "cpp", + "type_traits": "cpp", + "utility": "cpp", + "fstream": "cpp", + "initializer_list": "cpp", + "iomanip": "cpp", + "iosfwd": "cpp", + "iostream": "cpp", + "istream": "cpp", + "limits": "cpp", + "mutex": "cpp", + "new": "cpp", + "ostream": "cpp", + "sstream": "cpp", + "stdexcept": "cpp", + "streambuf": "cpp", + "thread": "cpp", + "cfenv": "cpp", + "typeindex": "cpp", + "typeinfo": "cpp", + "valarray": "cpp", + "variant": "cpp", + "bit": "cpp", + "stack": "cpp", + "xstring": "cpp", + "xutility": "cpp", + "xlocnum": "cpp", + "xtr1common": "cpp" + }, + "python.analysis.diagnosticSeverityOverrides": { + "reportInvalidTypeForm": "none" + }, + "C_Cpp_Runner.cCompilerPath": "gcc", + "C_Cpp_Runner.cppCompilerPath": "g++", + "C_Cpp_Runner.debuggerPath": "gdb", + "C_Cpp_Runner.cStandard": "", + "C_Cpp_Runner.cppStandard": "", + "C_Cpp_Runner.msvcBatchPath": "C:/Program Files/Microsoft Visual Studio/VR_NR/Community/VC/Auxiliary/Build/vcvarsall.bat", + "C_Cpp_Runner.useMsvc": false, + "C_Cpp_Runner.warnings": [ + "-Wall", + "-Wextra", + "-Wpedantic", + "-Wshadow", + "-Wformat=2", + "-Wcast-align", + "-Wconversion", + "-Wsign-conversion", + "-Wnull-dereference" + ], + "C_Cpp_Runner.msvcWarnings": [ + "/W4", + "/permissive-", + "/w14242", + "/w14287", + "/w14296", + "/w14311", + "/w14826", + "/w44062", + "/w44242", + "/w14905", + "/w14906", + "/w14263", + "/w44265", + "/w14928" + ], + "C_Cpp_Runner.enableWarnings": true, + "C_Cpp_Runner.warningsAsError": false, + "C_Cpp_Runner.compilerArgs": [], + "C_Cpp_Runner.linkerArgs": [], + "C_Cpp_Runner.includePaths": [], + "C_Cpp_Runner.includeSearch": [ + "*", + "**/*" + ], + "C_Cpp_Runner.excludeSearch": [ + "**/build", + "**/build/**", + "**/.*", + "**/.*/**", + "**/.vscode", + "**/.vscode/**" + ], + "C_Cpp_Runner.useAddressSanitizer": false, + "C_Cpp_Runner.useUndefinedSanitizer": false, + "C_Cpp_Runner.useLeakSanitizer": false, + "C_Cpp_Runner.showCompilationTime": false, + "C_Cpp_Runner.useLinkTimeOptimization": false, + "C_Cpp_Runner.msvcSecureNoWarnings": false } \ No newline at end of file From d85208e94b31bac6dee65757de842c0e703e9299 Mon Sep 17 00:00:00 2001 From: enigmahe Date: Mon, 24 Jun 2024 19:42:04 +0800 Subject: [PATCH 07/17] WIP: vol PT debug (strange shadowing) --- .vscode/c_cpp_properties.json | 30 +++--- .vscode/launch.json | 24 +++++ bxdf/setup.py | 4 +- bxdf/vol_loader/vol2numpy.cpp | 6 +- bxdf/volume.py | 172 +++++++++++++++++++++++++++------- renderer/constants.py | 1 + scenes/cbox/cbox-rgbvol.xml | 136 +++++++++++++++++++++++++++ tracer/setup.py | 4 +- 8 files changed, 319 insertions(+), 58 deletions(-) create mode 100644 .vscode/launch.json create mode 100644 scenes/cbox/cbox-rgbvol.xml diff --git a/.vscode/c_cpp_properties.json b/.vscode/c_cpp_properties.json index 4245cf6..4c50ae3 100644 --- a/.vscode/c_cpp_properties.json +++ b/.vscode/c_cpp_properties.json @@ -1,18 +1,16 @@ { - "configurations": [ - { - "name": "Linux", - "includePath": [ - "${workspaceFolder}/ext/pybind11/include/**", - "/usr/include/python3.8/**", - "/usr/include/eigen3/**" - ], - "defines": [], - "compilerPath": "/usr/bin/gcc", - "cStandard": "c11", - "cppStandard": "c++11", - "intelliSenseMode": "clang-x64" - } - ], - "version": 4 + "configurations": [ + { + "name": "windows-gcc-x64", + "includePath": [ + "${workspaceFolder}/**" + ], + "defines": [], + "compilerPath": "gcc", + "cStandard": "${default}", + "cppStandard": "${default}", + "intelliSenseMode": "windows-gcc-x64" + } + ], + "version": 4 } \ No newline at end of file diff --git a/.vscode/launch.json b/.vscode/launch.json new file mode 100644 index 0000000..09ca61a --- /dev/null +++ b/.vscode/launch.json @@ -0,0 +1,24 @@ +{ + "version": "0.2.0", + "configurations": [ + { + "name": "C/C++ Runner: Debug Session", + "type": "cppdbg", + "request": "launch", + "args": [], + "stopAtEntry": false, + "externalConsole": true, + "cwd": "e:/AdaPT/tracer/bvh", + "program": "e:/AdaPT/tracer/bvh/build/Debug/outDebug", + "MIMode": "gdb", + "miDebuggerPath": "gdb", + "setupCommands": [ + { + "description": "Enable pretty-printing for gdb", + "text": "-enable-pretty-printing", + "ignoreFailures": true + } + ] + } + ] +} \ No newline at end of file diff --git a/bxdf/setup.py b/bxdf/setup.py index a75d6d8..811d065 100644 --- a/bxdf/setup.py +++ b/bxdf/setup.py @@ -9,7 +9,7 @@ ["vol_loader/vol2numpy.cpp"], # Example: passing in the version to the compiled code define_macros = [('VERSION_INFO', __version__)], - extra_compile_args= ['-g', '-O3', '-fopenmp'], + extra_compile_args= ['-g', '-O2', '-fopenmp'], ), ] @@ -18,7 +18,7 @@ version=__version__, author="Qianyue He", description="Volume grid loader", - include_dirs="/usr/include/eigen3/", + include_dirs="E:\\eigen-3.4.0", long_description="", ext_modules=ext_modules, extras_require={"test": "pytest"}, diff --git a/bxdf/vol_loader/vol2numpy.cpp b/bxdf/vol_loader/vol2numpy.cpp index 136f82a..029f84a 100644 --- a/bxdf/vol_loader/vol2numpy.cpp +++ b/bxdf/vol_loader/vol2numpy.cpp @@ -62,6 +62,8 @@ bool readVolumeData(const std::string& filename, VolumeData& volume) { return false; } + printf("Channel: %d\n", volume.channels); + file.seekg(24, std::ios::cur); // Skip the bounding box int numVoxels = volume.xres * volume.yres * volume.zres * volume.channels; @@ -73,7 +75,7 @@ bool readVolumeData(const std::string& filename, VolumeData& volume) { } // mitsuba3 vol to numpy -auto loadVol2Numpy(const std::string& filename, bool force_mono_color = false) { +auto loadVol2Numpy(const std::string& filename, bool force_mono_color) { VolumeData volume; readVolumeData(filename, volume); @@ -118,7 +120,7 @@ auto loadVol2Numpy(const std::string& filename, bool force_mono_color = false) { if (force_mono_color) volume.channels = 1; - return std::tuple(vol_numpy, volume.shape()); + return std::tuple, std::tuple>(vol_numpy, volume.shape()); } PYBIND11_MODULE(vol_loader, m) { diff --git a/bxdf/volume.py b/bxdf/volume.py index a260c64..365cdc8 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -16,7 +16,7 @@ from la.cam_transform import delocalize_rotate from parsers.general_parser import get, rgb_parse, transform_parse from sampler.general_sampling import random_rgb -from renderer.constants import ZERO_V3, ONES_V3 +from renderer.constants import ZERO_V3, ONES_V3, EPSL_V3 from rich.console import Console CONSOLE = Console(width = 128) @@ -31,7 +31,7 @@ __all__ = ["GridVolume_np", "GridVolume"] -FORCE_MONO_COLOR = True +FORCE_MONO_COLOR = False """ TODO: add transform parsing """ @@ -104,6 +104,7 @@ def __init__(self, elem: xet.Element): if self.type_id == GridVolume_np.MONO: self.density_grid *= self.density_scaling[0] else: + print(f"Shape: {self.density_grid.shape}, {self.density_scaling}") self.density_grid *= self.density_scaling def assign_transform(self, trans_r, trans_t, trans_s): @@ -116,10 +117,17 @@ def setup_volume(self, path:str): CONSOLE.log(f"[red]Error :skull:[/red]: {path} contains no valid volume file.") raise RuntimeError("Volume file not found.") density_grid, (self.xres, self.yres, self.zres, self.channel) = vol_file_to_numpy(path, FORCE_MONO_COLOR) + density_grid = GridVolume_np.make_rgb_volume_from_monocolor(density_grid) + self.channel = 3 if FORCE_MONO_COLOR: CONSOLE.log(f"[yellow]Warning[/yellow]: FORCE_MONO_COLOR is True. This only makes sense when we are testing the code.") return density_grid.reshape((self.zres, self.yres, self.xres, self.channel)) + @staticmethod + def make_rgb_volume_from_monocolor(density_grid): + density_grid = np.stack([density_grid, density_grid * 0.2, density_grid * 0.2], axis = -1) + return density_grid + def get_shape(self) -> Tuple[int, int, int]: return (self.zres, self.yres, self.xres) @@ -166,7 +174,27 @@ def __repr__(self): return f"" - + + +@ti.func +def rgb_select(rgb: vec3, channel: int) -> float: + """ Get the value of a specified channel + This implement looks dumb, but I think it has its purpose: + Dynamic indexing will cause the local array get moved to global memory + i.e. the access speed will drop significantly + + I am not sure whether Taichi Lang has the problem or not, yet I know CUDA + has this problem + """ + result = 0.0 + if channel == 0: + result = rgb[0] + elif channel == 1: + result = rgb[1] + else: + result = rgb[2] + return result + @ti.dataclass class GridVolume: """ Grid Volume Taichi End definition""" @@ -220,8 +248,8 @@ def transmittance(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: fl ray_o = self.inv_T @ (ray_o - self.trans) ray_d = self.inv_T @ ray_d if self._type: - if self._type == GridVolume_np.MONO: # single channel - transm = self.eval_tr_ratio_tracking(grid, ray_o, ray_d, near_far) + if self._type == GridVolume_np.RGB: # single channel + transm = self.eval_tr_ratio_tracking_3d(grid, ray_o, ray_d, near_far) else: transm = self.eval_tr_ratio_tracking_3d(grid, ray_o, ray_d, near_far) return transm @@ -235,8 +263,8 @@ def sample_mfp(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: float ray_o = self.inv_T @ (ray_o - self.trans) ray_d = self.inv_T @ ray_d if self._type: - if self._type == 1: # single channel - transm = self.sample_distance_delta_tracking(grid, ray_o, ray_d, near_far) + if self._type == 2: # single channel + transm = self.sample_distance_delta_tracking_3d(grid, ray_o, ray_d, near_far) else: transm = self.sample_distance_delta_tracking_3d(grid, ray_o, ray_d, near_far) return transm @@ -270,8 +298,8 @@ def sample_distance_delta_tracking(self, grid: ti.template(), ray_ol: vec3, ray_ inv_maj = 1.0 / self.majorant[0] t = near_far[0] + t -= ti.log(1.0 - ti.random(float)) * inv_maj while t < near_far[1]: - t -= ti.log(1.0 - ti.random(float)) * inv_maj d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) @@ -280,27 +308,29 @@ def sample_distance_delta_tracking(self, grid: ti.template(), ray_ol: vec3, ray_ result[:3] = self.albedo result[3] = t break + t -= ti.log(1.0 - ti.random(float)) * inv_maj return result @ti.func def eval_tr_ratio_tracking(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: - inv_maj = 1.0 / self.majorant[0] - - t = near_far[0] Tr = 1.0 - while True: - t -= ti.log(1.0 - ti.random(float)) * inv_maj - if t >= near_far[1]: break - d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ - ti.random(float), ti.random(float), ti.random(float) - ])) - Tr *= ti.max(0, 1.0 - d * inv_maj) - # Russian Roulette - if Tr < 0.1: - if ti.random(float) >= Tr: - Tr = 0.0 - break - Tr = 1.0 + if self._type: + inv_maj = 1.0 / self.majorant[0] + + t = near_far[0] + while True: + t -= ti.log(1.0 - ti.random(float)) * inv_maj + if t >= near_far[1]: break + d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ + ti.random(float), ti.random(float), ti.random(float) + ])) + Tr *= ti.max(0, 1.0 - d * inv_maj) + # Russian Roulette + if Tr < 0.1: + if ti.random(float) >= Tr: + Tr = 0.0 + break + Tr = 1.0 return vec3([Tr, Tr, Tr]) @ti.func @@ -318,29 +348,99 @@ def sample_new_rays(self, incid: vec3): def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec4: result = vec4([1, 1, 1, -1]) if self._type: - - # First, choose one element according to the majorant - inv_maj = 1.0 / ti.max(self.majorant[2], 1e-4) - val = ti.random(float) + # mis_pdf = ONES_V3 + + # Step 1: choose one element according to the majorant + channel = 2 + maj = self.majorant[2] + val = ti.random(float) if val < self.pdf[0]: - inv_maj = 1.0 / self.majorant[0] + maj = self.majorant[0] + channel = 0 elif val < self.pdf[0] + self.pdf[1]: - inv_maj = 1.0 / self.majorant[1] + maj = self.majorant[1] + channel = 1 + inv_maj = 1.0 / maj t = near_far[0] - while t < near_far[1]: - t -= ti.log(1.0 - ti.random(float)) * inv_maj - d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ + dist = -ti.log(1.0 - ti.random(float)) * inv_maj + while t + dist < near_far[1]: + t += dist + + d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) + # Scatter upon real collision - if ti.random(float) < d * inv_maj: - result[:3] = self.albedo + n_s = rgb_select(d, channel) + if ti.random(float) < n_s * inv_maj: + result[:3] *= self.albedo result[3] = t + # mis_pdf *= d / ti.max(EPSL_V3, self.majorant) break + + maj_tr = ti.exp(-self.majorant * dist) + result[:3] *= maj_tr / ti.exp(- dist * maj) + # mis_pdf *= maj_tr + dist = -ti.log(1.0 - ti.random(float)) * inv_maj + if t + dist >= near_far[1]: + dist = near_far[1] - t + + maj_tr = ti.exp(-self.majorant * dist) + result[:3] *= maj_tr / ti.exp(- dist * maj) + # mis_pdf *= maj_tr + + # MIS + # mis_numerator = rgb_select(mis_pdf, channel) + # result[:3] *= mis_numerator / mis_pdf.sum() return result @ti.func def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: - return ONES_V3 + Tr = ONES_V3 + if self._type: + # mis_pdf = ONES_V3 + + # Step 1: choose one element according to the majorant + channel = 2 + maj = self.majorant[2] + val = ti.random(float) + if val < self.pdf[0]: + maj = self.majorant[0] + channel = 0 + elif val < self.pdf[0] + self.pdf[1]: + maj = self.majorant[1] + channel = 1 + inv_maj = 1.0 / maj + + t = near_far[0] + while True: + dist = -ti.log(1.0 - ti.random(float)) * inv_maj + + if t + dist >= near_far[1]: + dist = near_far[1] - t # get the actual optical length + Tr *= ti.exp(- self.majorant * dist) / ti.exp(- dist * maj) + break + t += dist + # for mono-chromatic medium, this is 1 + d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ + ti.random(float), ti.random(float), ti.random(float) + ])) + + majorant_tr = ti.exp(- self.majorant * dist) + Tr *= ti.max(ZERO_V3, 1.0 - d * inv_maj) * majorant_tr / ti.exp(- dist * maj) + # mis_pdf *= d * majorant_tr * self.pdf + + # Russian Roulette + if Tr.max() < 0.1: + rr_max_val = Tr.max() + if ti.random(float) >= rr_max_val: + Tr.fill(0) + break + Tr /= rr_max_val + + # evaluate MIS weight + # mis_numerator = rgb_select(mis_pdf, channel) + # Tr *= mis_numerator / mis_pdf.sum() + return Tr \ No newline at end of file diff --git a/renderer/constants.py b/renderer/constants.py index 684c513..9b9f711 100644 --- a/renderer/constants.py +++ b/renderer/constants.py @@ -34,6 +34,7 @@ ZERO_V3 = vec3([0, 0, 0]) ONES_V3 = vec3([1, 1, 1]) +EPSL_V3 = vec3([1e-4, 1e-4, 1e-4]) INVALID = vec3([-1, -1, -1]) AXIS_X = vec3([1, 0, 0]) AXIS_Y = vec3([0, 1, 0]) diff --git a/scenes/cbox/cbox-rgbvol.xml b/scenes/cbox/cbox-rgbvol.xml new file mode 100644 index 0000000..eac48ca --- /dev/null +++ b/scenes/cbox/cbox-rgbvol.xml @@ -0,0 +1,136 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + \ No newline at end of file diff --git a/tracer/setup.py b/tracer/setup.py index ecdae95..4950aca 100644 --- a/tracer/setup.py +++ b/tracer/setup.py @@ -9,7 +9,7 @@ ["bvh/bvh.cpp"], # Example: passing in the version to the compiled code define_macros = [('VERSION_INFO', __version__)], - extra_compile_args= ['-g', '-O3'], + extra_compile_args= ['-g', '-O2'], ), ] @@ -18,7 +18,7 @@ version=__version__, author="Qianyue He", description="BVH constructed via C++ backend", - include_dirs="/usr/include/eigen3/", + include_dirs="E:\\eigen-3.4.0", long_description="", ext_modules=ext_modules, extras_require={"test": "pytest"}, From 1146ba4a50dff87da428ec666fe1da6c16cf280d Mon Sep 17 00:00:00 2001 From: enigmahe Date: Mon, 24 Jun 2024 20:40:23 +0800 Subject: [PATCH 08/17] WIP: result correct. --- bxdf/volume.py | 72 ++++++++++++++++++++------------------------------ 1 file changed, 29 insertions(+), 43 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index 365cdc8..ebf608b 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -42,8 +42,6 @@ class GridVolume_np: __type_mapping = {"none": NONE, "mono": MONO, "rgb": RGB} def __init__(self, elem: xet.Element): self.albedo = np.ones(3, np.float32) - self.par = np.zeros(3, np.float32) - self.pdf = np.float32([1., 0., 0.]) self.phase_type_id = -1 self.phase_type = "hg" self.type_name = "none" @@ -117,15 +115,24 @@ def setup_volume(self, path:str): CONSOLE.log(f"[red]Error :skull:[/red]: {path} contains no valid volume file.") raise RuntimeError("Volume file not found.") density_grid, (self.xres, self.yres, self.zres, self.channel) = vol_file_to_numpy(path, FORCE_MONO_COLOR) - density_grid = GridVolume_np.make_rgb_volume_from_monocolor(density_grid) - self.channel = 3 if FORCE_MONO_COLOR: CONSOLE.log(f"[yellow]Warning[/yellow]: FORCE_MONO_COLOR is True. This only makes sense when we are testing the code.") - return density_grid.reshape((self.zres, self.yres, self.xres, self.channel)) + + density_grid = density_grid.reshape((self.zres, self.yres, self.xres, self.channel)) + density_grid = GridVolume_np.make_colorful_volume(density_grid, self.xres, self.yres, self.zres) + self.channel = 3 + return density_grid @staticmethod - def make_rgb_volume_from_monocolor(density_grid): - density_grid = np.stack([density_grid, density_grid * 0.2, density_grid * 0.2], axis = -1) + def make_colorful_volume(density_grid: np.ndarray, xres: int, yres: int, zres: int): + z_coords = np.linspace(0, 0.9, zres, dtype = np.float32).reshape(-1, 1, 1, 1) + 0.1 + y_coords = np.linspace(0, 0.9, yres, dtype = np.float32).reshape(1, -1, 1, 1) + 0.1 + x_coords = np.linspace(0, 0.9, xres, dtype = np.float32).reshape(1, 1, -1, 1) + 0.1 + density_grid = np.concatenate([ + density_grid * z_coords, + density_grid * y_coords, + density_grid * x_coords + ], axis = -1) return density_grid def get_shape(self) -> Tuple[int, int, int]: @@ -135,7 +142,9 @@ def export(self): if self.type_id == GridVolume_np.NONE: return GridVolume(_type = 0) aabb_mini, aabb_maxi = self.get_aabb() - maj = self.density_grid.max(axis = (0, 1, 2)) # the shape of density grid: (zres, yres, xres, channels) + maj = self.density_grid.max(axis = (0, 1, 2)) + majorant = vec3(maj) if self.type_id == GridVolume_np.RGB else vec3([maj, maj, maj]) + # the shape of density grid: (zres, yres, xres, channels) return GridVolume( _type = self.type_id, albedo = vec3(self.albedo), @@ -144,7 +153,8 @@ def export(self): mini = vec3(aabb_mini), maxi = vec3(aabb_maxi), max_idxs = vec3i([self.xres, self.yres, self.zres]), - majorant = vec3(maj) if self.type_id == GridVolume_np.RGB else vec3([maj, maj, maj]), + majorant = majorant, + pdf = majorant / majorant.sum(), ph = PhaseFunction(_type = self.phase_type_id, par = vec3(self.par), pdf = vec3(self.pdf)) ) @@ -348,11 +358,9 @@ def sample_new_rays(self, incid: vec3): def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec4: result = vec4([1, 1, 1, -1]) if self._type: - # mis_pdf = ONES_V3 - # Step 1: choose one element according to the majorant - channel = 2 - maj = self.majorant[2] + channel = 0 + maj = 1.0 val = ti.random(float) if val < self.pdf[0]: maj = self.majorant[0] @@ -360,6 +368,9 @@ def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, r elif val < self.pdf[0] + self.pdf[1]: maj = self.majorant[1] channel = 1 + else: + channel = 2 + maj = self.majorant[2] inv_maj = 1.0 / maj t = near_far[0] @@ -376,60 +387,39 @@ def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, r if ti.random(float) < n_s * inv_maj: result[:3] *= self.albedo result[3] = t - # mis_pdf *= d / ti.max(EPSL_V3, self.majorant) break - - maj_tr = ti.exp(-self.majorant * dist) - result[:3] *= maj_tr / ti.exp(- dist * maj) - # mis_pdf *= maj_tr dist = -ti.log(1.0 - ti.random(float)) * inv_maj if t + dist >= near_far[1]: dist = near_far[1] - t - - maj_tr = ti.exp(-self.majorant * dist) - result[:3] *= maj_tr / ti.exp(- dist * maj) - # mis_pdf *= maj_tr - - # MIS - # mis_numerator = rgb_select(mis_pdf, channel) - # result[:3] *= mis_numerator / mis_pdf.sum() return result @ti.func def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: Tr = ONES_V3 if self._type: - # mis_pdf = ONES_V3 - # Step 1: choose one element according to the majorant - channel = 2 - maj = self.majorant[2] + maj = 1.0 val = ti.random(float) if val < self.pdf[0]: maj = self.majorant[0] - channel = 0 elif val < self.pdf[0] + self.pdf[1]: maj = self.majorant[1] - channel = 1 + else: + maj = self.majorant[2] inv_maj = 1.0 / maj t = near_far[0] while True: dist = -ti.log(1.0 - ti.random(float)) * inv_maj - if t + dist >= near_far[1]: - dist = near_far[1] - t # get the actual optical length - Tr *= ti.exp(- self.majorant * dist) / ti.exp(- dist * maj) - break + if t + dist >= near_far[1]: break t += dist # for mono-chromatic medium, this is 1 d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) - majorant_tr = ti.exp(- self.majorant * dist) - Tr *= ti.max(ZERO_V3, 1.0 - d * inv_maj) * majorant_tr / ti.exp(- dist * maj) - # mis_pdf *= d * majorant_tr * self.pdf + Tr *= ti.max(ZERO_V3, 1.0 - d * inv_maj) # Russian Roulette if Tr.max() < 0.1: @@ -438,9 +428,5 @@ def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: v Tr.fill(0) break Tr /= rr_max_val - - # evaluate MIS weight - # mis_numerator = rgb_select(mis_pdf, channel) - # Tr *= mis_numerator / mis_pdf.sum() return Tr \ No newline at end of file From 528f395dae6f6c24b441c56f1b6a2e20478120ac Mon Sep 17 00:00:00 2001 From: enigmahe Date: Mon, 24 Jun 2024 21:17:24 +0800 Subject: [PATCH 09/17] Colored volume (yet I think this is strange) --- bxdf/volume.py | 81 +++++++------------------------------ renderer/vpt.py | 8 +--- scenes/cbox/cbox-rgbvol.xml | 1 + 3 files changed, 18 insertions(+), 72 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index ebf608b..40fc562 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -64,12 +64,14 @@ def __init__(self, elem: xet.Element): # phase function self.par = np.zeros(3, np.float32) self.pdf = np.float32([1., 0., 0.]) + self.mono2rgb = False elem_to_query = { "rgb": rgb_parse, "float": lambda el: get(el, "value"), "string": lambda el: self.setup_volume(get(el, "path", str)), - "transform": lambda el: self.assign_transform(*transform_parse(el)) + "transform": lambda el: self.assign_transform(*transform_parse(el)), + "bool": lambda el: get(el, "value", bool) } if elem is not None: type_name = elem.get("type") @@ -93,6 +95,16 @@ def __init__(self, elem: xet.Element): name = tag_elem.get("name") if hasattr(self, name): self.__setattr__(name, query_func(tag_elem)) + + if self.mono2rgb == True: + if self.channel == 3: + CONSOLE.log("Setting 'mono2rgb = True' is meanless for RGB volume. This setting is only meaningful when the volume is mono-chromatic.") + else: + CONSOLE.log("Mono-chromatic volume is converted to RGB volume.") + self.channel = 3 + self.density_grid = GridVolume_np.make_colorful_volume(self.density_grid, self.xres, self.yres, self.zres) + else: + self.density_grid = np.concatenate([self.density_grid] * 3, axis = -1) if self.scale is None: self.scale = np.eye(3, dtype = np.float32) @@ -119,8 +131,6 @@ def setup_volume(self, path:str): CONSOLE.log(f"[yellow]Warning[/yellow]: FORCE_MONO_COLOR is True. This only makes sense when we are testing the code.") density_grid = density_grid.reshape((self.zres, self.yres, self.xres, self.channel)) - density_grid = GridVolume_np.make_colorful_volume(density_grid, self.xres, self.yres, self.zres) - self.channel = 3 return density_grid @staticmethod @@ -258,10 +268,7 @@ def transmittance(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: fl ray_o = self.inv_T @ (ray_o - self.trans) ray_d = self.inv_T @ ray_d if self._type: - if self._type == GridVolume_np.RGB: # single channel - transm = self.eval_tr_ratio_tracking_3d(grid, ray_o, ray_d, near_far) - else: - transm = self.eval_tr_ratio_tracking_3d(grid, ray_o, ray_d, near_far) + transm = self.eval_tr_ratio_tracking_3d(grid, ray_o, ray_d, near_far) return transm @ti.func @@ -273,10 +280,7 @@ def sample_mfp(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: float ray_o = self.inv_T @ (ray_o - self.trans) ray_d = self.inv_T @ ray_d if self._type: - if self._type == 2: # single channel - transm = self.sample_distance_delta_tracking_3d(grid, ray_o, ray_d, near_far) - else: - transm = self.sample_distance_delta_tracking_3d(grid, ray_o, ray_d, near_far) + transm = self.sample_distance_delta_tracking_3d(grid, ray_o, ray_d, near_far) return transm @ti.func @@ -288,61 +292,6 @@ def density_lookup_3d(self, grid: ti.template(), index: vec3, u_offset: vec3) -> val = grid[idx[2], idx[1], idx[0]] return val - @ti.func - def density_lookup(self, grid: ti.template(), index, u_offset: vec3) -> float: - """ Stochastic lookup of density (RGB volume) """ - idx = ti.cast(ti.floor(index + (u_offset - 0.5)), int) - val = 0.0 - if (idx >= 0).all() and (idx <= self.max_idxs).all(): - val = grid[idx[2], idx[1], idx[0]] - # indexing pattern z, y, x - return val - - @ti.func - def sample_distance_delta_tracking(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec4: - """ Sample distance (mono-chromatic volume) via delta tracking - Note that there is no 'sampling PDF', since we are not going to use it anyway - """ - result = vec4([1, 1, 1, -1]) - if self._type: - inv_maj = 1.0 / self.majorant[0] - - t = near_far[0] - t -= ti.log(1.0 - ti.random(float)) * inv_maj - while t < near_far[1]: - d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ - ti.random(float), ti.random(float), ti.random(float) - ])) - # Scatter upon real collision - if ti.random(float) < d * inv_maj: - result[:3] = self.albedo - result[3] = t - break - t -= ti.log(1.0 - ti.random(float)) * inv_maj - return result - - @ti.func - def eval_tr_ratio_tracking(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: - Tr = 1.0 - if self._type: - inv_maj = 1.0 / self.majorant[0] - - t = near_far[0] - while True: - t -= ti.log(1.0 - ti.random(float)) * inv_maj - if t >= near_far[1]: break - d = self.density_lookup(grid, ray_ol + t * ray_dl, vec3([ - ti.random(float), ti.random(float), ti.random(float) - ])) - Tr *= ti.max(0, 1.0 - d * inv_maj) - # Russian Roulette - if Tr < 0.1: - if ti.random(float) >= Tr: - Tr = 0.0 - break - Tr = 1.0 - return vec3([Tr, Tr, Tr]) - @ti.func def sample_new_rays(self, incid: vec3): ret_spec = vec3([1, 1, 1]) diff --git a/renderer/vpt.py b/renderer/vpt.py index 4483c17..e37136e 100644 --- a/renderer/vpt.py +++ b/renderer/vpt.py @@ -40,12 +40,8 @@ def __init__(self, if prop["volume"]: CONSOLE.log("Loading Grid Volume ...") volume = GridVolume_np(prop["volume"][0]) - if volume.type_id == GridVolume_np.MONO: - self.density_grid = ti.field(float, volume.get_shape()) - volume.density_grid = volume.density_grid.squeeze() - elif volume.type_id == GridVolume_np.RGB: - self.density_grid = ti.Vector.field(3, float, volume.get_shape()) - else: + self.density_grid = ti.Vector.field(3, float, volume.get_shape()) + if volume.type_id not in{ GridVolume_np.MONO, GridVolume_np.RGB}: CONSOLE.log(f"Volume type '{volume.type_name}' is not supported.") raise RuntimeError("Unsupported volume type.") self.has_volume = True diff --git a/scenes/cbox/cbox-rgbvol.xml b/scenes/cbox/cbox-rgbvol.xml index eac48ca..e502c90 100644 --- a/scenes/cbox/cbox-rgbvol.xml +++ b/scenes/cbox/cbox-rgbvol.xml @@ -115,6 +115,7 @@ + From 44cfed20de8d97ef5c8f34843eb13e972e6cd7d5 Mon Sep 17 00:00:00 2001 From: enigmahe Date: Tue, 25 Jun 2024 14:44:47 +0800 Subject: [PATCH 10/17] Updated volume rendering. RGB rendering problem persists. --- bxdf/setup.py | 3 ++- bxdf/vol_loader/vol2numpy.cpp | 2 -- bxdf/volume.py | 7 +++---- emitters/abtract_source.py | 2 -- renderer/vpt.py | 2 -- scenes/cbox/cbox-rgbvol.xml | 18 ++++++++++++++++-- scenes/csphere/balls-mono.xml | 2 +- tracer/path_tracer.py | 6 ++++-- tracer/setup.py | 3 ++- 9 files changed, 28 insertions(+), 17 deletions(-) diff --git a/bxdf/setup.py b/bxdf/setup.py index 811d065..8b222fe 100644 --- a/bxdf/setup.py +++ b/bxdf/setup.py @@ -1,5 +1,6 @@ from pybind11.setup_helpers import Pybind11Extension from setuptools import setup +import platform __version__ = "0.2.0" cxx_std=11 @@ -9,7 +10,7 @@ ["vol_loader/vol2numpy.cpp"], # Example: passing in the version to the compiled code define_macros = [('VERSION_INFO', __version__)], - extra_compile_args= ['-g', '-O2', '-fopenmp'], + extra_compile_args= ['-g', '-O3' if platform.system() == "Linux" else '-O2', '-fopenmp'], ), ] diff --git a/bxdf/vol_loader/vol2numpy.cpp b/bxdf/vol_loader/vol2numpy.cpp index 029f84a..7a6747a 100644 --- a/bxdf/vol_loader/vol2numpy.cpp +++ b/bxdf/vol_loader/vol2numpy.cpp @@ -62,8 +62,6 @@ bool readVolumeData(const std::string& filename, VolumeData& volume) { return false; } - printf("Channel: %d\n", volume.channels); - file.seekg(24, std::ios::cur); // Skip the bounding box int numVoxels = volume.xres * volume.yres * volume.zres * volume.channels; diff --git a/bxdf/volume.py b/bxdf/volume.py index 40fc562..d494135 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -114,7 +114,6 @@ def __init__(self, elem: xet.Element): if self.type_id == GridVolume_np.MONO: self.density_grid *= self.density_scaling[0] else: - print(f"Shape: {self.density_grid.shape}, {self.density_scaling}") self.density_grid *= self.density_scaling def assign_transform(self, trans_r, trans_t, trans_s): @@ -153,7 +152,8 @@ def export(self): return GridVolume(_type = 0) aabb_mini, aabb_maxi = self.get_aabb() maj = self.density_grid.max(axis = (0, 1, 2)) - majorant = vec3(maj) if self.type_id == GridVolume_np.RGB else vec3([maj, maj, maj]) + majorant = vec3(maj) if self.type_id == GridVolume_np.RGB else vec3([maj, maj, maj]) + CONSOLE.log(f"Majorant: {majorant}. PDF: {majorant / majorant.sum()}") # the shape of density grid: (zres, yres, xres, channels) return GridVolume( _type = self.type_id, @@ -344,7 +344,7 @@ def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, r @ti.func def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: - Tr = ONES_V3 + Tr = ONES_V3 if self._type: # Step 1: choose one element according to the majorant maj = 1.0 @@ -356,7 +356,6 @@ def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: v else: maj = self.majorant[2] inv_maj = 1.0 / maj - t = near_far[0] while True: dist = -ti.log(1.0 - ti.random(float)) * inv_maj diff --git a/emitters/abtract_source.py b/emitters/abtract_source.py index 761dda0..720733c 100644 --- a/emitters/abtract_source.py +++ b/emitters/abtract_source.py @@ -107,8 +107,6 @@ def sample_hit( if ti.static(HEMISPHERE_SAMPLE_SPHERE): to_hit = (hit_pos - center).normalized() - # the pdf here can be viewed as being both both sa & area measure - # since for a unit sphere, different unit solid angle extends to the same amount of area local_dir, pdf = uniform_sphere() normal, _ = delocalize_rotate(to_hit, local_dir) ret_pos = center + normal * radius diff --git a/renderer/vpt.py b/renderer/vpt.py index e37136e..91f298d 100644 --- a/renderer/vpt.py +++ b/renderer/vpt.py @@ -36,7 +36,6 @@ def __init__(self, ): super().__init__(emitters, array_info, objects, prop, bvh_delay) self.has_volume = False - self.volume = None if prop["volume"]: CONSOLE.log("Loading Grid Volume ...") volume = GridVolume_np(prop["volume"][0]) @@ -49,7 +48,6 @@ def __init__(self, self.volume = volume.export() CONSOLE.log(f"Grid volume loaded: {volume}") else: - self.volume = GridVolume_np(None).export() self.density_grid = ti.field(float, (1, 1, 1)) self.world_scattering = self.world.medium._type >= 0 diff --git a/scenes/cbox/cbox-rgbvol.xml b/scenes/cbox/cbox-rgbvol.xml index e502c90..ac31562 100644 --- a/scenes/cbox/cbox-rgbvol.xml +++ b/scenes/cbox/cbox-rgbvol.xml @@ -117,12 +117,26 @@ + + + + + + + + diff --git a/scenes/csphere/balls-mono.xml b/scenes/csphere/balls-mono.xml index 8a8a410..cdbd0f6 100644 --- a/scenes/csphere/balls-mono.xml +++ b/scenes/csphere/balls-mono.xml @@ -78,7 +78,7 @@ - + diff --git a/tracer/path_tracer.py b/tracer/path_tracer.py index db05d63..4b15140 100644 --- a/tracer/path_tracer.py +++ b/tracer/path_tracer.py @@ -9,7 +9,6 @@ import sys sys.path.append("..") -import tqdm import numpy as np import taichi as ti import taichi.math as tm @@ -22,6 +21,7 @@ from bxdf.brdf import BRDF from bxdf.bsdf import BSDF, BSDF_np +from bxdf.volume import GridVolume_np from bxdf.texture import Texture, Texture_np from parsers.opts import get_options from parsers.obj_desc import ObjDescriptor @@ -31,7 +31,7 @@ from sampler.general_sampling import * from utils.tools import TicToc from tracer.interaction import Interaction -from tracer.ti_bvh import LinearBVH, LinearNode, export_python_bvh +from tracer.ti_bvh import LinearBVH, LinearNode from rich.console import Console CONSOLE = Console(width = 128) @@ -90,6 +90,8 @@ def __init__(self, # maybe we should opt for environment map (env lighting) self.roughness_map = Texture.field() # TODO: this is useless for now + self.volume = GridVolume_np(None).export() + ti.root.dense(ti.i, self.src_num).place(self.src_field) # Light source Taichi storage self.obj_nodes = ti.root.bitmasked(ti.i, self.num_objects) self.obj_nodes.place(self.brdf_field) # BRDF Taichi storage diff --git a/tracer/setup.py b/tracer/setup.py index 4950aca..0646bf6 100644 --- a/tracer/setup.py +++ b/tracer/setup.py @@ -1,5 +1,6 @@ from pybind11.setup_helpers import Pybind11Extension from setuptools import setup +import platform __version__ = "0.1.0" cxx_std=11 @@ -9,7 +10,7 @@ ["bvh/bvh.cpp"], # Example: passing in the version to the compiled code define_macros = [('VERSION_INFO', __version__)], - extra_compile_args= ['-g', '-O2'], + extra_compile_args= ['-g', '-O3' if platform.system() == "Linux" else '-O2'], ), ] From 6ac936f8932ee5ee0d88f87df3fb21cd62a91af5 Mon Sep 17 00:00:00 2001 From: enigmahe Date: Tue, 25 Jun 2024 21:08:57 +0800 Subject: [PATCH 11/17] WIP: throughput sampling (feed throughput) --- bxdf/volume.py | 91 +++++++++++++++++++++++++----------- scenes/cbox/cbox-rgbvol.xml | 6 +-- scenes/cbox/cbox-volgrid.xml | 2 + 3 files changed, 68 insertions(+), 31 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index d494135..25442ce 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -71,7 +71,7 @@ def __init__(self, elem: xet.Element): "float": lambda el: get(el, "value"), "string": lambda el: self.setup_volume(get(el, "path", str)), "transform": lambda el: self.assign_transform(*transform_parse(el)), - "bool": lambda el: get(el, "value", bool) + "bool": lambda el: el.get("value") in {"True", "true"} } if elem is not None: type_name = elem.get("type") @@ -102,9 +102,10 @@ def __init__(self, elem: xet.Element): else: CONSOLE.log("Mono-chromatic volume is converted to RGB volume.") self.channel = 3 + self.type_id = GridVolume_np.RGB self.density_grid = GridVolume_np.make_colorful_volume(self.density_grid, self.xres, self.yres, self.zres) else: - self.density_grid = np.concatenate([self.density_grid] * 3, axis = -1) + self.density_grid = np.concatenate([self.density_grid, self.density_grid, self.density_grid], axis = -1) if self.scale is None: self.scale = np.eye(3, dtype = np.float32) @@ -113,8 +114,10 @@ def __init__(self, elem: xet.Element): self.forward_t = self.rotation @ self.scale if self.type_id == GridVolume_np.MONO: self.density_grid *= self.density_scaling[0] + CONSOLE.log(f"Volume density scaled by {self.density_scaling[0]}.") else: self.density_grid *= self.density_scaling + CONSOLE.log(f"Volume density scaled by {self.density_scaling}.") def assign_transform(self, trans_r, trans_t, trans_s): self.rotation = trans_r @@ -138,9 +141,9 @@ def make_colorful_volume(density_grid: np.ndarray, xres: int, yres: int, zres: i y_coords = np.linspace(0, 0.9, yres, dtype = np.float32).reshape(1, -1, 1, 1) + 0.1 x_coords = np.linspace(0, 0.9, xres, dtype = np.float32).reshape(1, 1, -1, 1) + 0.1 density_grid = np.concatenate([ - density_grid * z_coords, - density_grid * y_coords, - density_grid * x_coords + density_grid, + density_grid, + density_grid, ], axis = -1) return density_grid @@ -308,25 +311,31 @@ def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, r result = vec4([1, 1, 1, -1]) if self._type: # Step 1: choose one element according to the majorant + Tr = 1.0 + pdf = 1.0 + albedo = 0.0 channel = 0 - maj = 1.0 + maj = self.majorant[2] val = ti.random(float) if val < self.pdf[0]: + albedo = self.albedo[0] maj = self.majorant[0] + pdf = self.pdf[0] channel = 0 elif val < self.pdf[0] + self.pdf[1]: + albedo = self.albedo[1] maj = self.majorant[1] + pdf = self.pdf[1] channel = 1 else: - channel = 2 + albedo = self.albedo[2] maj = self.majorant[2] + pdf = self.pdf[2] + channel = 2 inv_maj = 1.0 / maj - t = near_far[0] - dist = -ti.log(1.0 - ti.random(float)) * inv_maj - while t + dist < near_far[1]: - t += dist - + t = near_far[0] - ti.log(1.0 - ti.random(float)) * inv_maj + while t < near_far[1]: d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) @@ -334,47 +343,73 @@ def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, r # Scatter upon real collision n_s = rgb_select(d, channel) if ti.random(float) < n_s * inv_maj: - result[:3] *= self.albedo + Tr *= albedo result[3] = t break - dist = -ti.log(1.0 - ti.random(float)) * inv_maj - if t + dist >= near_far[1]: - dist = near_far[1] - t + + t -= ti.log(1.0 - ti.random(float)) * inv_maj + if self._type == GridVolume_np.RGB: + if channel == 0: + result[:3] = vec3([Tr / pdf, 0, 0]) + elif channel == 1: + result[:3] = vec3([0, Tr / pdf, 0]) + else: + result[:3] = vec3([0, 0, Tr / pdf]) + else: + result[:3] = vec3([Tr, Tr, Tr]) return result @ti.func def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: - Tr = ONES_V3 + transm = ONES_V3 if self._type: # Step 1: choose one element according to the majorant - maj = 1.0 + Tr = 1.0 + pdf = 1.0 + channel = 0 + maj = self.majorant[2] val = ti.random(float) if val < self.pdf[0]: maj = self.majorant[0] + pdf = self.pdf[0] + channel = 0 elif val < self.pdf[0] + self.pdf[1]: maj = self.majorant[1] + pdf = self.pdf[1] + channel = 1 else: maj = self.majorant[2] + pdf = self.pdf[2] + channel = 2 inv_maj = 1.0 / maj + t = near_far[0] while True: - dist = -ti.log(1.0 - ti.random(float)) * inv_maj + t -= ti.log(1.0 - ti.random(float)) * inv_maj - if t + dist >= near_far[1]: break - t += dist + if t >= near_far[1]: break # for mono-chromatic medium, this is 1 d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) - Tr *= ti.max(ZERO_V3, 1.0 - d * inv_maj) + n_s = rgb_select(d, channel) + Tr *= ti.max(0.0, 1.0 - n_s * inv_maj) # Russian Roulette - if Tr.max() < 0.1: - rr_max_val = Tr.max() - if ti.random(float) >= rr_max_val: - Tr.fill(0) + if Tr < 0.1: + if ti.random(float) >= Tr: + Tr = 0.0 break - Tr /= rr_max_val - return Tr + Tr = 1.0 + if self._type == GridVolume_np.RGB: + if channel == 0: + transm = vec3([Tr / pdf, 0, 0]) + elif channel == 1: + transm = vec3([0, Tr / pdf, 0]) + else: + transm = vec3([0, 0, Tr / pdf]) + else: + transm = vec3([Tr, Tr, Tr]) + return transm \ No newline at end of file diff --git a/scenes/cbox/cbox-rgbvol.xml b/scenes/cbox/cbox-rgbvol.xml index ac31562..61bd898 100644 --- a/scenes/cbox/cbox-rgbvol.xml +++ b/scenes/cbox/cbox-rgbvol.xml @@ -68,7 +68,7 @@ - + @@ -110,10 +110,10 @@ - + - + diff --git a/scenes/cbox/cbox-volgrid.xml b/scenes/cbox/cbox-volgrid.xml index 77b0310..5a7911a 100755 --- a/scenes/cbox/cbox-volgrid.xml +++ b/scenes/cbox/cbox-volgrid.xml @@ -148,6 +148,7 @@ + @@ -161,6 +162,7 @@ + From dff234647c46211296d0cbf113bf47450e75f49c Mon Sep 17 00:00:00 2001 From: Enigmatisms <984041003@qq.com> Date: Wed, 26 Jun 2024 00:23:48 +0800 Subject: [PATCH 12/17] Seems to be the problem of early exiting. --- bxdf/volume.py | 47 +++++++++++++++++++++--------------- renderer/bdpt.py | 4 +-- renderer/vpt.py | 12 ++++----- scenes/cbox/cbox-volgrid.xml | 2 +- 4 files changed, 36 insertions(+), 29 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index 25442ce..9e199d7 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -15,8 +15,7 @@ from bxdf.medium import Medium_np from la.cam_transform import delocalize_rotate from parsers.general_parser import get, rgb_parse, transform_parse -from sampler.general_sampling import random_rgb -from renderer.constants import ZERO_V3, ONES_V3, EPSL_V3 +from renderer.constants import ZERO_V3, ONES_V3 from rich.console import Console CONSOLE = Console(width = 128) @@ -263,7 +262,7 @@ def intersect_volume(self, ray_o: vec3, ray_d: vec3, max_t: float) -> vec2: return near_far @ti.func - def transmittance(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: float) -> vec3: + def transmittance(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, thp: vec3, max_t: float) -> vec3: transm = vec3([1, 1, 1]) if self._type: near_far = self.intersect_volume(ray_o, ray_d, max_t) @@ -271,11 +270,11 @@ def transmittance(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: fl ray_o = self.inv_T @ (ray_o - self.trans) ray_d = self.inv_T @ ray_d if self._type: - transm = self.eval_tr_ratio_tracking_3d(grid, ray_o, ray_d, near_far) + transm = self.eval_tr_ratio_tracking_3d(grid, ray_o, ray_d, thp, near_far) return transm @ti.func - def sample_mfp(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: float) -> vec4: + def sample_mfp(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, thp: vec3, max_t: float) -> vec4: transm = vec4([1, 1, 1, -1]) if self._type: near_far = self.intersect_volume(ray_o, ray_d, max_t) @@ -283,7 +282,7 @@ def sample_mfp(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, max_t: float ray_o = self.inv_T @ (ray_o - self.trans) ray_d = self.inv_T @ ray_d if self._type: - transm = self.sample_distance_delta_tracking_3d(grid, ray_o, ray_d, near_far) + transm = self.sample_distance_delta_tracking_3d(grid, ray_o, ray_d, thp, near_far) return transm @ti.func @@ -307,30 +306,34 @@ def sample_new_rays(self, incid: vec3): return ret_dir, ret_spec, ret_pdf @ti.func - def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec4: + def sample_distance_delta_tracking_3d(self, + grid: ti.template(), ray_ol: vec3, ray_dl: vec3, thp: vec3, near_far: vec2 + ) -> vec4: result = vec4([1, 1, 1, -1]) if self._type: # Step 1: choose one element according to the majorant + pdfs = thp * self.pdf + pdfs /= pdfs.sum() Tr = 1.0 pdf = 1.0 albedo = 0.0 channel = 0 maj = self.majorant[2] val = ti.random(float) - if val < self.pdf[0]: + if val < pdfs[0]: albedo = self.albedo[0] maj = self.majorant[0] - pdf = self.pdf[0] + pdf = pdfs[0] channel = 0 - elif val < self.pdf[0] + self.pdf[1]: + elif val < pdfs[0] + pdfs[1]: albedo = self.albedo[1] maj = self.majorant[1] - pdf = self.pdf[1] + pdf = pdfs[1] channel = 1 else: albedo = self.albedo[2] maj = self.majorant[2] - pdf = self.pdf[2] + pdf = pdfs[2] channel = 2 inv_maj = 1.0 / maj @@ -343,8 +346,8 @@ def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, r # Scatter upon real collision n_s = rgb_select(d, channel) if ti.random(float) < n_s * inv_maj: - Tr *= albedo - result[3] = t + Tr *= albedo + result[3] = t break t -= ti.log(1.0 - ti.random(float)) * inv_maj @@ -360,26 +363,30 @@ def sample_distance_delta_tracking_3d(self, grid: ti.template(), ray_ol: vec3, r return result @ti.func - def eval_tr_ratio_tracking_3d(self, grid: ti.template(), ray_ol: vec3, ray_dl: vec3, near_far: vec2) -> vec3: + def eval_tr_ratio_tracking_3d(self, + grid: ti.template(), ray_ol: vec3, ray_dl: vec3, thp: vec3, near_far: vec2 + ) -> vec3: transm = ONES_V3 if self._type: # Step 1: choose one element according to the majorant + pdfs = thp * self.pdf + pdfs /= pdfs.sum() Tr = 1.0 pdf = 1.0 channel = 0 maj = self.majorant[2] val = ti.random(float) - if val < self.pdf[0]: + if val < pdfs[0]: maj = self.majorant[0] - pdf = self.pdf[0] + pdf = pdfs[0] channel = 0 - elif val < self.pdf[0] + self.pdf[1]: + elif val < pdfs[0] + pdfs[1]: maj = self.majorant[1] - pdf = self.pdf[1] + pdf = pdfs[1] channel = 1 else: maj = self.majorant[2] - pdf = self.pdf[2] + pdf = pdfs[2] channel = 2 inv_maj = 1.0 / maj diff --git a/renderer/bdpt.py b/renderer/bdpt.py index 6c1adaf..93006b8 100644 --- a/renderer/bdpt.py +++ b/renderer/bdpt.py @@ -240,7 +240,7 @@ def random_walk( # Step 2: check for mean free path sampling # Calculate mfp, path_beta = transmittance / PDF - is_mi, it.min_depth, path_beta = self.sample_mfp(ray_o, ray_d, it.obj_id, in_free_space, it.min_depth) + is_mi, it.min_depth, path_beta = self.sample_mfp(ray_o, ray_d, throughput, it.obj_id, in_free_space, it.min_depth) if it.obj_id < 0 and not is_mi: break # exiting world bound throughput *= path_beta # attenuate first if throughput.max() < 5e-5: break @@ -392,7 +392,7 @@ def connect_path(self, i: int, j: int, sid: int, tid: int, decomp: int): le = cam_v.beta * fr_cam * fr_lit * lit_v.beta / (depth * depth) ret_time = lit_v.time + cam_v.time if le.max() > 0 and calc_transmittance == True: - tr, track_depth = self.track_ray(connect_dir, track_pos, depth) + tr, track_depth = self.track_ray(connect_dir, track_pos, ONES_V3, depth) le *= tr ret_time += track_depth weight = 0. diff --git a/renderer/vpt.py b/renderer/vpt.py index 91f298d..0cb21f9 100644 --- a/renderer/vpt.py +++ b/renderer/vpt.py @@ -72,7 +72,7 @@ def non_null_surface(self, idx: int): return non_null @ti.func - def sample_mfp(self, ray_o: vec3, ray_d: vec3, idx: int, in_free_space: int, depth: float): + def sample_mfp(self, ray_o: vec3, ray_d: vec3, thp: vec3, idx: int, in_free_space: int, depth: float): """ Mean free path sampling, returns: - whether medium is a valid scattering medium / mean free path """ @@ -91,7 +91,7 @@ def sample_mfp(self, ray_o: vec3, ray_d: vec3, idx: int, in_free_space: int, dep is_mi, mfp, beta = self.bsdf_field[idx].medium.sample_mfp(depth) # use medium to sample / calculate transmittance if ti.static(self.has_volume): # grid volume event might override the world medium event (not physically based, but simple to implement) - result = self.volume.sample_mfp(self.density_grid, ray_o, ray_d, depth) + result = self.volume.sample_mfp(self.density_grid, ray_o, ray_d, thp, depth) if result[3] > 0: # in case grid volume is nested with world medium is_mi = 2 mfp = result[3] @@ -99,7 +99,7 @@ def sample_mfp(self, ray_o: vec3, ray_d: vec3, idx: int, in_free_space: int, dep return is_mi, mfp, beta @ti.func - def track_ray(self, cur_ray, cur_point, depth): + def track_ray(self, cur_ray: vec3, cur_point: vec3, thp: vec3, depth: float): """ For medium interaction, check if the path to one point is not blocked (by non-null surface object) And also we need to calculate the attenuation along the path, e.g.: if the ray passes through @@ -108,7 +108,7 @@ def track_ray(self, cur_ray, cur_point, depth): """ tr = vec3([1., 1., 1.]) if ti.static(self.has_volume): # grid volume transmittance - tr *= self.volume.transmittance(self.density_grid, cur_point, cur_ray, depth) + tr = self.volume.transmittance(self.density_grid, cur_point, cur_ray, thp, depth) in_free_space = True acc_depth = 0.0 @@ -176,7 +176,7 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int in_free_space = tm.dot(it.n_g, ray_d) < 0 # Step 3: check for mean free path sampling # Calculate mfp, path_beta = transmittance / PDF - is_mi, it.min_depth, path_beta = self.sample_mfp(ray_o, ray_d, it.obj_id, in_free_space, it.min_depth) + is_mi, it.min_depth, path_beta = self.sample_mfp(ray_o, ray_d, throughput, it.obj_id, in_free_space, it.min_depth) if it.obj_id < 0 and not is_mi: break # exiting world bound hit_point = ray_d * it.min_depth + ray_o throughput *= path_beta # attenuate first @@ -202,7 +202,7 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int to_emitter = emit_pos - hit_point emitter_d = to_emitter.norm() light_dir = to_emitter / emitter_d - tr, _ = self.track_ray(light_dir, hit_point, emitter_d) + tr, _ = self.track_ray(light_dir, hit_point, throughput, emitter_d) shadow_int *= tr direct_spec = self.eval(it, ray_d, light_dir, is_mi, in_free_space) else: # the only situation for being invalid, is when there is only one source and the ray hit the source diff --git a/scenes/cbox/cbox-volgrid.xml b/scenes/cbox/cbox-volgrid.xml index 5a7911a..17ec7b2 100755 --- a/scenes/cbox/cbox-volgrid.xml +++ b/scenes/cbox/cbox-volgrid.xml @@ -9,7 +9,7 @@ - + From 852d73a696189a874418613e107a0e1cb6cf247e Mon Sep 17 00:00:00 2001 From: enigmahe Date: Wed, 26 Jun 2024 16:22:02 +0800 Subject: [PATCH 13/17] Current state of debugging --- bxdf/volume.py | 34 ++++++++++++++++++++-------------- renderer/vpt.py | 9 +++++++-- scenes/cbox/cbox-rgbvol.xml | 6 +++--- scenes/cbox/cbox-volgrid.xml | 4 ++-- 4 files changed, 32 insertions(+), 21 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index 9e199d7..463f2d6 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -148,13 +148,19 @@ def make_colorful_volume(density_grid: np.ndarray, xres: int, yres: int, zres: i def get_shape(self) -> Tuple[int, int, int]: return (self.zres, self.yres, self.xres) + + def get_majorant(self, guard = 0.2, scale_ratio = 1.05): + maj = self.density_grid.max(axis = (0, 1, 2)) + maj_guard = np.mean(maj) * guard + maj = np.maximum(maj, maj_guard) + maj *= scale_ratio + return vec3(maj) if self.type_id == GridVolume_np.RGB else vec3([maj, maj, maj]) def export(self): if self.type_id == GridVolume_np.NONE: return GridVolume(_type = 0) aabb_mini, aabb_maxi = self.get_aabb() - maj = self.density_grid.max(axis = (0, 1, 2)) - majorant = vec3(maj) if self.type_id == GridVolume_np.RGB else vec3([maj, maj, maj]) + majorant = self.get_majorant() CONSOLE.log(f"Majorant: {majorant}. PDF: {majorant / majorant.sum()}") # the shape of density grid: (zres, yres, xres, channels) return GridVolume( @@ -164,7 +170,7 @@ def export(self): trans = vec3(self.offset), mini = vec3(aabb_mini), maxi = vec3(aabb_maxi), - max_idxs = vec3i([self.xres, self.yres, self.zres]), + max_idxs = vec3i([self.xres - 1, self.yres - 1, self.zres - 1]), majorant = majorant, pdf = majorant / majorant.sum(), ph = PhaseFunction(_type = self.phase_type_id, par = vec3(self.par), pdf = vec3(self.pdf)) @@ -188,7 +194,7 @@ def get_aabb(self): ]) world_point = self.local_to_world(all_points) # conservative AABB - return world_point.min(axis = 0) - 0.1, world_point.max(axis = 0) + 0.1 + return world_point.min(axis = 0) - 0.01, world_point.max(axis = 0) + 0.01 def __repr__(self): @@ -257,8 +263,8 @@ def intersect_volume(self, ray_o: vec3, ray_d: vec3, max_t: float) -> vec2: tmin = ti.min(t1s, t2s) tmax = ti.max(t1s, t2s) - near_far[0] = ti.max(0, tmin.max()) + 1e-4 - near_far[1] = ti.min(max_t, tmax.min()) - 1e-4 + near_far[0] = ti.max(0, tmin.max()) + 1e-5 + near_far[1] = ti.min(max_t, tmax.min()) - 1e-5 return near_far @ti.func @@ -288,10 +294,10 @@ def sample_mfp(self, grid: ti.template(), ray_o: vec3, ray_d: vec3, thp: vec3, m @ti.func def density_lookup_3d(self, grid: ti.template(), index: vec3, u_offset: vec3) -> vec3: """ Stochastic lookup of density (mono-chromatic volume) """ - idx = ti.cast(ti.floor(index + (u_offset - 0.5)), int) - val = ZERO_V3 + idx = ti.cast(ti.floor(index + (u_offset - 0.5)), int) + val = ZERO_V3 if (idx >= 0).all() and (idx <= self.max_idxs).all(): - val = grid[idx[2], idx[1], idx[0]] + val = grid[idx[2], idx[1], idx[0]] return val @ti.func @@ -320,12 +326,12 @@ def sample_distance_delta_tracking_3d(self, channel = 0 maj = self.majorant[2] val = ti.random(float) - if val < pdfs[0]: + if val <= pdfs[0]: albedo = self.albedo[0] maj = self.majorant[0] pdf = pdfs[0] channel = 0 - elif val < pdfs[0] + pdfs[1]: + elif val <= pdfs[0] + pdfs[1]: albedo = self.albedo[1] maj = self.majorant[1] pdf = pdfs[1] @@ -342,7 +348,6 @@ def sample_distance_delta_tracking_3d(self, d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) - # Scatter upon real collision n_s = rgb_select(d, channel) if ti.random(float) < n_s * inv_maj: @@ -376,11 +381,11 @@ def eval_tr_ratio_tracking_3d(self, channel = 0 maj = self.majorant[2] val = ti.random(float) - if val < pdfs[0]: + if val <= pdfs[0]: maj = self.majorant[0] pdf = pdfs[0] channel = 0 - elif val < pdfs[0] + pdfs[1]: + elif val <= pdfs[0] + pdfs[1]: maj = self.majorant[1] pdf = pdfs[1] channel = 1 @@ -392,6 +397,7 @@ def eval_tr_ratio_tracking_3d(self, t = near_far[0] while True: + # problem with coordinates t -= ti.log(1.0 - ti.random(float)) * inv_maj if t >= near_far[1]: break diff --git a/renderer/vpt.py b/renderer/vpt.py index 0cb21f9..7f91fef 100644 --- a/renderer/vpt.py +++ b/renderer/vpt.py @@ -155,10 +155,11 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int color = vec3([0, 0, 0]) throughput = vec3([1, 1, 1]) emission_weight = 1.0 - + in_free_space = True bounce = 0 - while True: + surface_bounce = 0 + while surface_bounce < 1: # for _i in range(self.max_bounce): # Step 1: ray termination test - Only RR termination is allowed max_value = throughput.max() @@ -178,6 +179,10 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int # Calculate mfp, path_beta = transmittance / PDF is_mi, it.min_depth, path_beta = self.sample_mfp(ray_o, ray_d, throughput, it.obj_id, in_free_space, it.min_depth) if it.obj_id < 0 and not is_mi: break # exiting world bound + + if is_mi == 0: + surface_bounce += 1 + hit_point = ray_d * it.min_depth + ray_o throughput *= path_beta # attenuate first if not is_mi and not self.non_null_surface(it.obj_id): diff --git a/scenes/cbox/cbox-rgbvol.xml b/scenes/cbox/cbox-rgbvol.xml index 61bd898..400b400 100644 --- a/scenes/cbox/cbox-rgbvol.xml +++ b/scenes/cbox/cbox-rgbvol.xml @@ -9,7 +9,7 @@ - + @@ -20,7 +20,7 @@ - + @@ -115,7 +115,7 @@ - + diff --git a/scenes/cbox/cbox-volgrid.xml b/scenes/cbox/cbox-volgrid.xml index 17ec7b2..b6d347b 100755 --- a/scenes/cbox/cbox-volgrid.xml +++ b/scenes/cbox/cbox-volgrid.xml @@ -9,7 +9,7 @@ - + @@ -20,7 +20,7 @@ - + From 2267944b706dbc6004bf99e49348bc73dfde6830 Mon Sep 17 00:00:00 2001 From: enigmahe Date: Wed, 26 Jun 2024 18:17:19 +0800 Subject: [PATCH 14/17] Russian Roulette is to blame. --- bxdf/volume.py | 10 ++++++---- renderer/vpt.py | 18 +++++++++--------- scenes/cbox/cbox-volgrid.xml | 4 ++-- 3 files changed, 17 insertions(+), 15 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index 463f2d6..43564ab 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -32,8 +32,6 @@ FORCE_MONO_COLOR = False -""" TODO: add transform parsing -""" class GridVolume_np: NONE = 0 MONO = 1 @@ -323,9 +321,11 @@ def sample_distance_delta_tracking_3d(self, Tr = 1.0 pdf = 1.0 albedo = 0.0 + maj = 0.0 channel = 0 - maj = self.majorant[2] val = ti.random(float) + + # avoid dynamic indexing on GPU (and array might be moved from local registers to global memory) if val <= pdfs[0]: albedo = self.albedo[0] maj = self.majorant[0] @@ -378,9 +378,11 @@ def eval_tr_ratio_tracking_3d(self, pdfs /= pdfs.sum() Tr = 1.0 pdf = 1.0 + maj = 0.0 channel = 0 - maj = self.majorant[2] val = ti.random(float) + + # avoid dynamic indexing on GPU (and array might be moved from local registers to global memory) if val <= pdfs[0]: maj = self.majorant[0] pdf = pdfs[0] diff --git a/renderer/vpt.py b/renderer/vpt.py index 7f91fef..2b94644 100644 --- a/renderer/vpt.py +++ b/renderer/vpt.py @@ -158,13 +158,16 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int in_free_space = True bounce = 0 - surface_bounce = 0 - while surface_bounce < 1: + while True: # for _i in range(self.max_bounce): - # Step 1: ray termination test - Only RR termination is allowed - max_value = throughput.max() - if ti.random(float) > max_value: break - else: throughput *= 1. / ti.max(max_value, 1e-7) # unbiased calculation + # Step 1: ray termination test + if ti.static(self.use_rr): + # Simple Russian Roullete ray termination + max_value = throughput.max() + if ti.random(float) > max_value: break + else: throughput *= 1. / (max_value + 1e-7) # unbiased calculation + else: + if throughput.max() < 1e-5: break # contribution too small, break # Step 2: ray intersection it = self.ray_intersect(ray_d, ray_o) if it.obj_id < 0: @@ -180,9 +183,6 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int is_mi, it.min_depth, path_beta = self.sample_mfp(ray_o, ray_d, throughput, it.obj_id, in_free_space, it.min_depth) if it.obj_id < 0 and not is_mi: break # exiting world bound - if is_mi == 0: - surface_bounce += 1 - hit_point = ray_d * it.min_depth + ray_o throughput *= path_beta # attenuate first if not is_mi and not self.non_null_surface(it.obj_id): diff --git a/scenes/cbox/cbox-volgrid.xml b/scenes/cbox/cbox-volgrid.xml index b6d347b..4fd8d3e 100755 --- a/scenes/cbox/cbox-volgrid.xml +++ b/scenes/cbox/cbox-volgrid.xml @@ -9,7 +9,7 @@ - + @@ -162,7 +162,7 @@ - + From 928c97973798f5b0f85e6793086a9760e7005674 Mon Sep 17 00:00:00 2001 From: Enigmatisms <984041003@qq.com> Date: Thu, 27 Jun 2024 01:17:46 +0800 Subject: [PATCH 15/17] Finished colored volume rendering. --- bxdf/volume.py | 12 ++++++------ renderer/vanilla_renderer.py | 7 ++++--- renderer/vpt.py | 5 +++-- scenes/cbox/cbox-rgbvol.xml | 19 ++++++++++--------- tracer/path_tracer.py | 2 ++ 5 files changed, 25 insertions(+), 20 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index 43564ab..e31cbce 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -134,13 +134,13 @@ def setup_volume(self, path:str): @staticmethod def make_colorful_volume(density_grid: np.ndarray, xres: int, yres: int, zres: int): - z_coords = np.linspace(0, 0.9, zres, dtype = np.float32).reshape(-1, 1, 1, 1) + 0.1 - y_coords = np.linspace(0, 0.9, yres, dtype = np.float32).reshape(1, -1, 1, 1) + 0.1 - x_coords = np.linspace(0, 0.9, xres, dtype = np.float32).reshape(1, 1, -1, 1) + 0.1 + z_coords = np.linspace(0, 4.0, zres, dtype = np.float32).reshape(-1, 1, 1, 1) + 0.1 + y_coords = np.linspace(0, 0.5, yres, dtype = np.float32).reshape(1, -1, 1, 1) + 0.1 + x_coords = np.linspace(0, 0.1, xres, dtype = np.float32).reshape(1, 1, -1, 1) + 0.01 density_grid = np.concatenate([ - density_grid, - density_grid, - density_grid, + density_grid * z_coords, + density_grid * y_coords, + density_grid * x_coords, ], axis = -1) return density_grid diff --git a/renderer/vanilla_renderer.py b/renderer/vanilla_renderer.py index b8a4606..007d854 100644 --- a/renderer/vanilla_renderer.py +++ b/renderer/vanilla_renderer.py @@ -45,13 +45,14 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int color = vec3([0, 0, 0]) contribution = vec3([1, 1, 1]) emission_weight = 1.0 - for _i in range(self.max_bounce): + for bounce in range(self.max_bounce): if it.is_ray_not_hit(): break # nothing is hit, break if ti.static(self.use_rr): # Simple Russian Roullete ray termination max_value = contribution.max() - if ti.random(float) > max_value: break - else: contribution *= 1. / (max_value + 1e-7) # unbiased calculation + if max_value < self.rr_threshold and bounce >= self.rr_bounce_th: + if ti.random(float) > max_value: break + else: contribution *= 1. / (max_value + 1e-7) # unbiased calculation else: if contribution.max() < 1e-4: break # contribution too small, break hit_point = ray_d * it.min_depth + ray_o diff --git a/renderer/vpt.py b/renderer/vpt.py index 2b94644..7527da5 100644 --- a/renderer/vpt.py +++ b/renderer/vpt.py @@ -164,8 +164,9 @@ def render(self, _t_start: int, _t_end: int, _s_start: int, _s_end: int, _a: int if ti.static(self.use_rr): # Simple Russian Roullete ray termination max_value = throughput.max() - if ti.random(float) > max_value: break - else: throughput *= 1. / (max_value + 1e-7) # unbiased calculation + if max_value < self.rr_threshold and bounce >= self.rr_bounce_th: + if ti.random(float) > max_value: break + else: throughput *= 1. / (max_value + 1e-7) # unbiased calculation else: if throughput.max() < 1e-5: break # contribution too small, break # Step 2: ray intersection diff --git a/scenes/cbox/cbox-rgbvol.xml b/scenes/cbox/cbox-rgbvol.xml index 400b400..57c36ae 100644 --- a/scenes/cbox/cbox-rgbvol.xml +++ b/scenes/cbox/cbox-rgbvol.xml @@ -9,7 +9,7 @@ - + @@ -17,10 +17,11 @@ - + + + - @@ -68,7 +69,7 @@ - + @@ -113,14 +114,14 @@ - + - + - - - + + + diff --git a/tracer/path_tracer.py b/tracer/path_tracer.py index 4b15140..5cb70c3 100644 --- a/tracer/path_tracer.py +++ b/tracer/path_tracer.py @@ -64,6 +64,8 @@ def __init__(self, self.stratified_sample = prop['stratified_sampling'] # whether to use stratified sampling self.use_mis = prop['use_mis'] # whether to use multiple importance sampling self.num_shadow_ray = prop['num_shadow_ray'] # number of shadow samples to trace + self.rr_threshold = prop.get('rr_threshold', 0.1) # threshold of employing RR + self.rr_bounce_th = prop.get('rr_bounce_th', 4) # minimum number of bounces to start RR self.brdf_two_sides = prop.get('brdf_two_sides', False) if self.num_shadow_ray > 0: From b25daab838651c8f2a28c0f68407ef8df4ba495a Mon Sep 17 00:00:00 2001 From: Enigmatisms <984041003@qq.com> Date: Fri, 28 Jun 2024 01:14:43 +0800 Subject: [PATCH 16/17] Debug: Erroneous color. --- bxdf/volume.py | 56 ++++++++++++++++++++++++++++++------ scenes/cbox/cbox-rgbvol.xml | 14 ++++----- scenes/cbox/cbox-volgrid.xml | 23 ++++++++------- 3 files changed, 66 insertions(+), 27 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index e31cbce..de7e337 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -7,6 +7,7 @@ import os import numpy as np import taichi as ti +import taichi.math as tm import xml.etree.ElementTree as xet from typing import Tuple @@ -134,15 +135,33 @@ def setup_volume(self, path:str): @staticmethod def make_colorful_volume(density_grid: np.ndarray, xres: int, yres: int, zres: int): - z_coords = np.linspace(0, 4.0, zres, dtype = np.float32).reshape(-1, 1, 1, 1) + 0.1 - y_coords = np.linspace(0, 0.5, yres, dtype = np.float32).reshape(1, -1, 1, 1) + 0.1 - x_coords = np.linspace(0, 0.1, xres, dtype = np.float32).reshape(1, 1, -1, 1) + 0.01 density_grid = np.concatenate([ - density_grid * z_coords, - density_grid * y_coords, - density_grid * x_coords, + density_grid, + density_grid, + density_grid ], axis = -1) - return density_grid + + # 创建一个从0到1的渐变数组，表示x轴从左到右 + half_x = zres // 3 + grad_l = np.linspace(1, 0, half_x, dtype = np.float32) ** 0.25 + grad_r = np.linspace(0, 1, zres - half_x, dtype = np.float32) ** 0.1 + + # 前半部分，从红色（1, 0, 0）到白色（1, 1, 1） + left_half = np.zeros((half_x, 3), dtype = np.float32) + left_half[:, 0] = 1 # 红色通道保持1 + left_half[:, 1] = 1 - grad_l # 绿色通道从0到1 + left_half[:, 2] = 1 - grad_l # 蓝色通道从0到1 + + # 后半部分，从白色（1, 1, 1）到蓝色（0, 0, 1） + right_half = np.zeros((zres - half_x, 3), dtype = np.float32) + right_half[:, 0] = 1 - grad_r # 红色通道从1到0 + right_half[:, 1] = 1 - grad_r # 绿色通道从1到0 + right_half[:, 2] = 1 # 蓝色通道保持1 + + # 将两个部分组合在一起 + color_gradient = np.vstack((left_half, right_half)) + + return density_grid * color_gradient[:, None, None, :] def get_shape(self) -> Tuple[int, int, int]: return (self.zres, self.yres, self.xres) @@ -298,6 +317,25 @@ def density_lookup_3d(self, grid: ti.template(), index: vec3, u_offset: vec3) -> val = grid[idx[2], idx[1], idx[0]] return val + @ti.func + def density_lookup_lerp_3d(self, grid: ti.template(), index: vec3, u_offset: vec3) -> vec3: + """ Stochastic lookup of density (mono-chromatic volume) """ + coord = index + (u_offset - 0.5) + idx = ti.cast(ti.floor(coord), int) + coord -= idx + val = ZERO_V3 + if (idx >= 0).all() and (idx <= self.max_idxs - 1).all(): + v1 = tm.mix(grid[idx[2], idx[1], idx[0]], grid[idx[2], idx[1], idx[0] + 1], coord[0]) + v2 = tm.mix(grid[idx[2] + 1, idx[1], idx[0]], grid[idx[2] + 1, idx[1], idx[0] + 1], coord[0]) + v3 = tm.mix(grid[idx[2], idx[1] + 1, idx[0]], grid[idx[2], idx[1] + 1, idx[0] + 1], coord[0]) + v4 = tm.mix(grid[idx[2] + 1, idx[1] + 1, idx[0]], grid[idx[2] + 1, idx[1] + 1, idx[0] + 1], coord[0]) + + v1 = tm.mix(v1, v2, coord[2]) + v3 = tm.mix(v3, v4, coord[2]) + + val = tm.mix(v1, v3, coord[1]) + return val + @ti.func def sample_new_rays(self, incid: vec3): ret_spec = vec3([1, 1, 1]) @@ -345,7 +383,7 @@ def sample_distance_delta_tracking_3d(self, t = near_far[0] - ti.log(1.0 - ti.random(float)) * inv_maj while t < near_far[1]: - d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ + d = self.density_lookup_lerp_3d(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) # Scatter upon real collision @@ -404,7 +442,7 @@ def eval_tr_ratio_tracking_3d(self, if t >= near_far[1]: break # for mono-chromatic medium, this is 1 - d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ + d = self.density_lookup_lerp_3d(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) diff --git a/scenes/cbox/cbox-rgbvol.xml b/scenes/cbox/cbox-rgbvol.xml index 57c36ae..8791cc8 100644 --- a/scenes/cbox/cbox-rgbvol.xml +++ b/scenes/cbox/cbox-rgbvol.xml @@ -50,13 +50,13 @@ - + - + @@ -69,7 +69,7 @@ - + @@ -114,14 +114,14 @@ - + - - - + + + diff --git a/scenes/cbox/cbox-volgrid.xml b/scenes/cbox/cbox-volgrid.xml index 4fd8d3e..e0e42f9 100755 --- a/scenes/cbox/cbox-volgrid.xml +++ b/scenes/cbox/cbox-volgrid.xml @@ -9,7 +9,7 @@ - + @@ -17,10 +17,11 @@ - + + + - @@ -66,19 +67,19 @@ - + - + - + @@ -104,14 +105,14 @@ - + @@ -160,9 +161,9 @@ - + - + From 9f0bfed73f914bad04d2f935206ff38e54377b53 Mon Sep 17 00:00:00 2001 From: Enigmatisms <984041003@qq.com> Date: Sat, 29 Jun 2024 11:03:16 +0800 Subject: [PATCH 17/17] Volume rendering no MIS completed. --- bxdf/volume.py | 24 +++++++--------- scenes/cbox/cbox-rgbvol.xml | 55 +++++++++++++++++++++++++------------ 2 files changed, 47 insertions(+), 32 deletions(-) diff --git a/bxdf/volume.py b/bxdf/volume.py index de7e337..e637cfa 100644 --- a/bxdf/volume.py +++ b/bxdf/volume.py @@ -141,24 +141,20 @@ def make_colorful_volume(density_grid: np.ndarray, xres: int, yres: int, zres: i density_grid ], axis = -1) - # 创建一个从0到1的渐变数组，表示x轴从左到右 half_x = zres // 3 - grad_l = np.linspace(1, 0, half_x, dtype = np.float32) ** 0.25 - grad_r = np.linspace(0, 1, zres - half_x, dtype = np.float32) ** 0.1 + grad_l = np.linspace(1, 0, half_x, dtype = np.float32) ** 0.65 + grad_r = np.linspace(0, 1, zres - half_x, dtype = np.float32) ** 0.6 - # 前半部分，从红色（1, 0, 0）到白色（1, 1, 1） left_half = np.zeros((half_x, 3), dtype = np.float32) - left_half[:, 0] = 1 # 红色通道保持1 - left_half[:, 1] = 1 - grad_l # 绿色通道从0到1 - left_half[:, 2] = 1 - grad_l # 蓝色通道从0到1 + left_half[:, 0] = 1 - grad_l + left_half[:, 1] = 1 + left_half[:, 2] = 1 - # 后半部分，从白色（1, 1, 1）到蓝色（0, 0, 1） right_half = np.zeros((zres - half_x, 3), dtype = np.float32) - right_half[:, 0] = 1 - grad_r # 红色通道从1到0 - right_half[:, 1] = 1 - grad_r # 绿色通道从1到0 - right_half[:, 2] = 1 # 蓝色通道保持1 + right_half[:, 0] = 1 + right_half[:, 1] = 1 + right_half[:, 2] = 1 - grad_r - # 将两个部分组合在一起 color_gradient = np.vstack((left_half, right_half)) return density_grid * color_gradient[:, None, None, :] @@ -383,7 +379,7 @@ def sample_distance_delta_tracking_3d(self, t = near_far[0] - ti.log(1.0 - ti.random(float)) * inv_maj while t < near_far[1]: - d = self.density_lookup_lerp_3d(grid, ray_ol + t * ray_dl, vec3([ + d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) # Scatter upon real collision @@ -442,7 +438,7 @@ def eval_tr_ratio_tracking_3d(self, if t >= near_far[1]: break # for mono-chromatic medium, this is 1 - d = self.density_lookup_lerp_3d(grid, ray_ol + t * ray_dl, vec3([ + d = self.density_lookup_3d(grid, ray_ol + t * ray_dl, vec3([ ti.random(float), ti.random(float), ti.random(float) ])) diff --git a/scenes/cbox/cbox-rgbvol.xml b/scenes/cbox/cbox-rgbvol.xml index 8791cc8..3f09345 100644 --- a/scenes/cbox/cbox-rgbvol.xml +++ b/scenes/cbox/cbox-rgbvol.xml @@ -9,7 +9,7 @@ - + @@ -18,8 +18,8 @@ - - + + @@ -50,13 +50,13 @@ - + - + @@ -67,18 +67,37 @@ - - - + + + + + + + + + + + + + - + - - - - + + + + + + + + + + + + + @@ -112,16 +131,16 @@ - + - - + + - + - +