3D : Individual Curvature

#!/usr/bin/env python3

import os
import vedo

import numpy as np
import pygranite as pg  # trajectory library
import pgt  # pygranite toolkit


def heightmap(csvf):
    return np.genfromtxt(csvf, delimiter=",")


if __name__ == "__main__":
    root = os.path.abspath(os.path.dirname(__file__))

    # simulation properties
    settings = pg.IntegratorSettings()
    settings.Space = pg.Space3D
    settings.MinimumAliveParticles = 16
    settings.WindfieldMode = pg.WindfieldMode.Dynamic
    settings.Integrator = pg.Integrator.ClassicRungeKutta
    settings.DeltaT = 1
    settings.DataTimeDistance = 60
    settings.MaximumSimulationTime = 1000
    settings.GridScale = [25, 25, 25]
    settings.SaveInterval = 1
    settings.CurvatureMode = pg.CurvatureMode.IndividualAndTotalCurvature

    # load heightmap
    hm = heightmap(root + "/../data/pyr_1000.csv")
    settings.Topography = hm

    # generate start positions
    spawn_min = np.array([2000, 2500, 150])
    spawn_max = np.array([2300, 2600, 170])

    inp = np.random.rand(1024, 3)
    for i in range(inp.shape[0]):
        inp[i] = ((spawn_max - spawn_min) * inp[i] + spawn_min)

    start_set = pg.TrajectorySet(inp)

    with pgt.Timer("Individual Curvature"):
        # initialize windfield loader
        directory = root + "fields_scaled/"
        loader = pgt.netCDFLoader3D([
            directory + f"{i}.nc" for i in range(len(os.listdir(directory)))
        ])
        loader.PrintProgress = False

        # computation
        integrator = pg.TrajectoryIntegrator(
            settings, loader, start_set)
        result_set: pg.TrajectorySet = integrator.compute()

        # visualization using vedo
        total_curvatures = result_set.totalCurvatures()
        curvatures = result_set.curvatures()

        print(total_curvatures.shape) 
        print(curvatures.shape)