diff --git a/src/lathe/lathe.py b/src/lathe/lathe.py index ffc0170..3c24cef 100644 --- a/src/lathe/lathe.py +++ b/src/lathe/lathe.py @@ -11,13 +11,15 @@ import numpy as np import typer -from mylogger import log, std_con +from numpy.typing import NDArray from rich.logging import RichHandler from rich.table import Table -from terrain import sample_octaves from typing_extensions import Annotated -from util import Mesh, MeshArray, create_mesh, now, rescale, save_world -from viz import viz + +from .mylogger import log, std_con +from .terrain import sample_octaves +from .util import Mesh, MeshArray, create_mesh, now, rescale, save_world +from .viz import viz # Define globals. @@ -33,7 +35,7 @@ ) # The highest elevation in the world, in meters. The highest peak on Earth is approximately 8700 m. ZRANGE: float = ZMAX - ZMIN ZTILT: float = 23.4 # Used in mapping spherical coordinates to lat-lon coordinates in a Coordinate Reference System (CRS). The present z-axis tilt of the Earth is approximately 23.4 degrees. -ZSCALE: int = 1 # A scaling factor for elevations to make them visible in the plot. +ZSCALE: int = 20 # A scaling factor for elevations to make them visible in the plot. OCEAN_PERCENT: float = 0.55 # Sets the point of elevation 0 as a relative percent of the range from zmin to zmax during rescaling. OCEAN_POINT: float = OCEAN_PERCENT * ZRANGE RECURSION: int = 9 # The number of recursions used in creating the icosphere mesh. 9 yields 2,621,442 points and 5,242,880 cells. Surface area of the Earth is approximately 514 million km2. @@ -262,7 +264,7 @@ def main( std_con.print("Rescaling elevations.\r\n") - rescaled_elevations: np.NDArray = rescale( + rescaled_elevations: NDArray[np.float64] = rescale( elevations=raw_elevations, zmin=zmin, zmax=zmax ) @@ -275,7 +277,9 @@ def main( std_con.print("Applying elevations to mesh.\r\n") log.debug(msg="Generating elevation scalars.") - elevation_scalars: np.NDArray = ((rescaled_elevations * zscale) + radius) / radius + elevation_scalars: NDArray[np.float64] = ( + ((raw_elevations) + radius) / radius + ) * zscale log.debug(msg="") log.debug(msg="Elevation Scalars:") @@ -317,9 +321,9 @@ def main( log.info(msg="Starting vizualization.") viz( world_mesh=world_mesh, - scalars="Elevations", + # scalars="Elevations", radius=radius, - zscale=ZSCALE, + zscale=zscale, zmin=zmin, zmax=zmax, ) diff --git a/src/lathe/terrain.py b/src/lathe/terrain.py index 8a6dcbe..76fa7f7 100644 --- a/src/lathe/terrain.py +++ b/src/lathe/terrain.py @@ -2,10 +2,12 @@ import numpy as np -from numba import prange import opensimplex as osi +from numba import prange from numpy.typing import NDArray +from .mylogger import std_con + def sample_noise( points, roughness, strength, feature_size, radius @@ -15,38 +17,18 @@ def sample_noise( for v in prange(len(rough_verts)): elevations[v] = osi.noise4( - x=rough_verts[v][0] / feature_size, - y=rough_verts[v][1] / feature_size, - z=rough_verts[v][2] / feature_size, + x=rough_verts[v][0], + y=rough_verts[v][1], + z=rough_verts[v][2], w=1 / feature_size, ) + std_con.print(f"Point: {v} ") # ?: Adding +1 to elevation moves negative values in the 0-1 range. Multiplying by 0.5 drags any values > 1 back into the 0-1 range. I'm not sure if multiplying by the radius is the proper thing to do in my next implementation. return (elevations + 1) * 0.5 * strength * radius -def v_sample_noise( - points: NDArray[np.float64], - roughness: float, - strength: float, - feature_size: float, - radius: float, -) -> NDArray[np.float64]: - # Vectorized operations - rough_verts = points * roughness - scaled_verts = rough_verts / feature_size - w_value = 1 / feature_size - - # Vectorized noise computation - elevations = osi.noise4( - x=scaled_verts[:, 0], y=scaled_verts[:, 1], z=scaled_verts[:, 2], w=w_value - ) - - # Adjust elevations - return (elevations + 1) * 0.5 * strength * radius - - def sample_octaves( points, octaves, @@ -74,7 +56,7 @@ def sample_octaves( # Initialize elevations array. - elevations = np.zeros(shape=len(points), dtype=np.float64) + elevations = np.ones(shape=len(points), dtype=np.float64) # NOTE: In my separate-sampling experiment, rough/strength pairs of (1.6, 0.4) (5, 0.2) and (24, 0.02) were good for 3 octaves. The final 3 results were added and then multiplied by 0.4 @@ -88,5 +70,6 @@ def sample_octaves( ) init_roughness *= roughness init_strength *= persistence + std_con.print(f"Octave: {i} ", "\r\n") return elevations diff --git a/src/lathe/viz.py b/src/lathe/viz.py index a13ee02..3237f85 100644 --- a/src/lathe/viz.py +++ b/src/lathe/viz.py @@ -5,11 +5,10 @@ # import external dependencies import pyvista as pv - # import pyvistaqt as pvqt -def viz(world_mesh, scalars, radius, zscale, zmin, zmax): +def viz(world_mesh, radius, zscale, zmin, zmax): logging.debug(msg="Starting viz().") logging.debug(msg="Setup PyVista plotter.") @@ -47,17 +46,17 @@ def viz(world_mesh, scalars, radius, zscale, zmin, zmax): ) annotations = {} - world_mesh.compute_normals(cell_normals=True, point_normals=True, inplace=True) - globe_mesh = world_mesh.warp_by_scalar(factor=-50) + # world_mesh.compute_normals(cell_normals=True, point_normals=True, inplace=True) + # globe_mesh = world_mesh.warp_by_scalar(factor=-50) pl.add_mesh( - globe_mesh, - name="Base Terrain", - scalars=scalars, + world_mesh, + name="Global Topography", + scalars="Elevations", scalar_bar_args=scalar_args, show_scalar_bar=True, annotations=annotations, - style="surface", + # style="surface", smooth_shading=True, show_edges=False, edge_color="red", @@ -65,26 +64,26 @@ def viz(world_mesh, scalars, radius, zscale, zmin, zmax): cmap="cmo.topo", lighting=True, pickable=True, - preference="cell", + # preference="cell", ) - # ocean_shell = pv.ParametricEllipsoid(radius, radius, radius, u_res=300, v_res=300) - - # pl.add_mesh( - # ocean_shell, - # show_edges=False, - # smooth_shading=True, - # color="blue", - # opacity=0.15, - # ) + ocean_shell = pv.ParametricEllipsoid(radius, radius, radius, u_res=300, v_res=300) - pl.show_bounds( - grid=True, - location="back", - axes_ranges=[0, 6400, 0, 6400, 0, 6400], - show_zlabels=True, + pl.add_mesh( + ocean_shell, + show_edges=False, + smooth_shading=True, + color="blue", + opacity=0.15, ) + # pl.show_bounds( + # grid=False, + # location="back", + # axes_ranges=[zmin, zmax, zmin, zmax, zmin, zmax], + # show_zlabels=True, + # ) + pl.camera.zoom(1.25) pl.enable_terrain_style(mouse_wheel_zooms=True)