diff --git a/build.gradle.kts b/build.gradle.kts index 6cb9e1d..ef99846 100644 --- a/build.gradle.kts +++ b/build.gradle.kts @@ -8,6 +8,7 @@ plugins { `maven-publish` alias(libs.plugins.nyx) + alias(libs.plugins.jmh) alias(libs.plugins.axion.release) } @@ -96,7 +97,7 @@ repositories { dependencies { api(libs.jetbrains.annotations) implementation(libs.slf4j.api) - + jmh(libs.slf4j.simple) testImplementation(libs.bundles.junit) } diff --git a/gradle/libs.versions.toml b/gradle/libs.versions.toml index 3fed5e6..4c82a76 100644 --- a/gradle/libs.versions.toml +++ b/gradle/libs.versions.toml @@ -1,6 +1,6 @@ [versions] nyx = "0.2.3" -jmh = "0.7.2" +jmh = "0.7.3" axion-release = "1.18.12" jetbrains-annotations = "26.0.1" @@ -22,6 +22,7 @@ axion-release = { id = "pl.allegro.tech.build.axion-release", version.ref = "axi jetbrains-annotations = { group = "org.jetbrains", name = "annotations", version.ref = "jetbrains-annotations" } slf4j-api = { group = "org.slf4j", name = "slf4j-api", version.ref = "slf4j-api" } +slf4j-simple = { group = "org.slf4j", name = "slf4j-simple", version.ref = "slf4j-api" } junit-jupiter-api = { group = "org.junit.jupiter", name = "junit-jupiter-api", version.ref = "junit-jupiter" } junit-jupiter-engine = { group = "org.junit.jupiter", name = "junit-jupiter-engine", version.ref = "junit-jupiter" } diff --git a/settings.gradle.kts b/settings.gradle.kts index b8f7892..8a006eb 100644 --- a/settings.gradle.kts +++ b/settings.gradle.kts @@ -1 +1,8 @@ -rootProject.name = "seismic" \ No newline at end of file +rootProject.name = "seismic" + +pluginManagement { + repositories { + gradlePluginPortal() + mavenCentral() + } +} diff --git a/src/jmh/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSamplerBenchmark.java b/src/jmh/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSamplerBenchmark.java new file mode 100644 index 0000000..09b8b35 --- /dev/null +++ b/src/jmh/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSamplerBenchmark.java @@ -0,0 +1,32 @@ +package com.dfsek.seismic.algorithms.sampler.noise; + +import org.openjdk.jmh.annotations.Benchmark; +import org.openjdk.jmh.annotations.Fork; +import org.openjdk.jmh.annotations.Group; +import org.openjdk.jmh.annotations.Measurement; +import org.openjdk.jmh.annotations.Scope; +import org.openjdk.jmh.annotations.State; +import org.openjdk.jmh.annotations.Warmup; + +@State(Scope.Benchmark) +public class CellularSamplerBenchmark { + private final NoiseFunction cellular = new CellularSampler(); + + @Benchmark + @Group("cellular") + @Fork(1) + @Warmup(iterations = 1, time = 1) + @Measurement(iterations = 2, time = 5) + public void benchmarkCellular3D() { + cellular.getNoiseRaw(0, 0, 0, 0); + } + + @Benchmark + @Group("cellular") + @Fork(1) + @Warmup(iterations = 1, time = 1) + @Measurement(iterations = 2, time = 5) + public void benchmarkCellular2D() { + cellular.getNoiseRaw(0, 0, 0); + } +} diff --git a/src/jmh/java/com/dfsek/seismic/math/trigonometry/TrigonometricFunctionsBenchmark.java b/src/jmh/java/com/dfsek/seismic/math/trigonometry/TrigonometricFunctionsBenchmark.java new file mode 100644 index 0000000..489d46c --- /dev/null +++ b/src/jmh/java/com/dfsek/seismic/math/trigonometry/TrigonometricFunctionsBenchmark.java @@ -0,0 +1,28 @@ +package com.dfsek.seismic.math.trigonometry; + +import org.openjdk.jmh.annotations.Benchmark; +import org.openjdk.jmh.annotations.Fork; +import org.openjdk.jmh.annotations.Measurement; +import org.openjdk.jmh.annotations.Scope; +import org.openjdk.jmh.annotations.State; +import org.openjdk.jmh.annotations.Warmup; + + +@State(Scope.Benchmark) +public class TrigonometricFunctionsBenchmark { + @Benchmark + @Fork(1) + @Warmup(iterations = 1, time = 1) + @Measurement(iterations = 2, time = 5) + public void benchmarkAtan2Java() { + double _x = Math.atan2(1.0, 1.0); + } + + @Benchmark + @Fork(1) + @Warmup(iterations = 1, time = 1) + @Measurement(iterations = 2, time = 5) + public void benchmarkAtan2Fast() { + double _x = TrigonometryFunctions.fastAtan2(1, 1); + } +} diff --git a/src/main/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSampler.java b/src/main/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSampler.java index a364db2..3bfd61e 100644 --- a/src/main/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSampler.java +++ b/src/main/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSampler.java @@ -11,6 +11,7 @@ import com.dfsek.seismic.algorithms.sampler.noise.simplex.OpenSimplex2Sampler; import com.dfsek.seismic.math.arithmetic.ArithmeticFunctions; import com.dfsek.seismic.math.floatingpoint.FloatingPointFunctions; +import com.dfsek.seismic.math.trigonometry.TrigonometryFunctions; import com.dfsek.seismic.type.DistanceFunction; import com.dfsek.seismic.type.sampler.Sampler; @@ -199,207 +200,236 @@ public class CellularSampler extends NoiseFunction { private DistanceFunction distanceFunction = DistanceFunction.EuclideanSq; private ReturnType returnType = ReturnType.Distance; private double jitterModifier = 1.0; - private Sampler noiseLookup; - private boolean saltLookup; + private FinalSampler inner; + public CellularSampler() { noiseLookup = new OpenSimplex2Sampler(); } + private FinalSampler getInner() { + if(inner == null) { + inner = new FinalSampler(distanceFunction, returnType, jitterModifier, noiseLookup, saltLookup, frequency, salt); + } + + return inner; + } + + @Override + public double getNoiseRaw(long seed, double x, double y) { + return getInner().getNoiseRaw(seed, x, y); + } + + @Override + public double getNoiseRaw(long seed, double x, double y, double z) { + return getInner().getNoiseRaw(seed, x, y, z); + } + public void setDistanceFunction(DistanceFunction distanceFunction) { this.distanceFunction = distanceFunction; + this.inner = null; } public void setJitterModifier(double jitterModifier) { this.jitterModifier = jitterModifier; + this.inner = null; } public void setNoiseLookup(Sampler noiseLookup) { this.noiseLookup = noiseLookup; + this.inner = null; } public void setReturnType(ReturnType returnType) { this.returnType = returnType; + this.inner = null; } public void setSaltLookup(boolean saltLookup) { this.saltLookup = saltLookup; + this.inner = null; } - @Override - public double getNoiseRaw(long sl, double x, double y) { - int seed = (int) sl; - int xr = FloatingPointFunctions.round(x); - int yr = FloatingPointFunctions.round(y); - - double distance0 = Double.MAX_VALUE; - double distance1 = Double.MAX_VALUE; - double distance2 = Double.MAX_VALUE; - - int closestHash = 0; - - double cellularJitter = 0.43701595 * jitterModifier; - - int xPrimed = (xr - 1) * NoiseFunction.PRIME_X; - int yPrimedBase = (yr - 1) * NoiseFunction.PRIME_Y; - - double centerX = x; - double centerY = y; - - for(int xi = xr - 1; xi <= xr + 1; xi++) { - int yPrimed = yPrimedBase; - - for(int yi = yr - 1; yi <= yr + 1; yi++) { - int hash = HashingFunctions.hashPrimeCoords(seed, xPrimed, yPrimed); - int idx = hash & (255 << 1); - - double vecX = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_2D[idx], cellularJitter, xi - x); - double vecY = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_2D[idx | 1], cellularJitter, yi - y); - - double newDistance = switch(distanceFunction) { - case Euclidean, EuclideanSq -> ArithmeticFunctions.fma(vecX, vecX, vecY * vecY); - case Manhattan -> Math.abs(vecX) + Math.abs(vecY); - case Hybrid -> (Math.abs(vecX) + Math.abs(vecY)) + ArithmeticFunctions.fma(vecX, vecX, vecY * vecY); - }; - - distance1 = Math.max(Math.min(distance1, newDistance), distance0); - if(newDistance < distance0) { - distance0 = newDistance; - closestHash = hash; - centerX = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_2D[idx], cellularJitter, xi) / frequency; - centerY = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_2D[idx | 1], cellularJitter, yi) / frequency; - } else if(newDistance < distance1) { - distance2 = distance1; - distance1 = newDistance; - } else if(newDistance < distance2) { - distance2 = newDistance; - } - yPrimed += NoiseFunction.PRIME_Y; - } - xPrimed += NoiseFunction.PRIME_X; - } + record FinalSampler(DistanceFunction distanceFunction, ReturnType returnType, double jitterModifier, Sampler noiseLookup, + boolean saltLookup, double frequency, long salt) { + public double getNoiseRaw(long sl, double x, double y) { + int seed = (int) sl; + int xr = FloatingPointFunctions.round(x); + int yr = FloatingPointFunctions.round(y); - if(distanceFunction == DistanceFunction.Euclidean && returnType != ReturnType.CellValue) { - distance0 = Math.sqrt(distance0); + double distance0 = Double.MAX_VALUE; + double distance1 = Double.MAX_VALUE; + double distance2 = Double.MAX_VALUE; - if(returnType != ReturnType.Distance) { - distance1 = Math.sqrt(distance1); - } - } + int closestHash = 0; - return switch(returnType) { - case CellValue -> closestHash * (1 / 2147483648.0); - case Distance -> distance0 - 1; - case Distance2 -> distance1 - 1; - case Distance2Add -> (distance1 + distance0) * 0.5 - 1; - case Distance2Sub -> distance1 - distance0 - 1; - case Distance2Mul -> distance1 * distance0 * 0.5 - 1; - case Distance2Div -> distance0 / distance1 - 1; - case NoiseLookup -> noiseLookup.getSample(sl - (saltLookup ? 0 : salt), centerX, centerY); - case LocalNoiseLookup -> noiseLookup.getSample(sl - (saltLookup ? 0 : salt), x / frequency - centerX, y / frequency - centerY); - case Distance3 -> distance2 - 1; - case Distance3Add -> (distance2 + distance0) * 0.5 - 1; - case Distance3Sub -> distance2 - distance0 - 1; - case Distance3Mul -> distance2 * distance0 - 1; - case Distance3Div -> distance0 / distance2 - 1; - case Angle -> Math.atan2(y / frequency - centerY, x / frequency - centerX); - }; - } + double cellularJitter = 0.43701595 * jitterModifier; - @Override - public double getNoiseRaw(long sl, double x, double y, double z) { - int seed = (int) sl; - int xr = FloatingPointFunctions.round(x); - int yr = FloatingPointFunctions.round(y); - int zr = FloatingPointFunctions.round(z); - - double distance0 = Double.MAX_VALUE; - double distance1 = Double.MAX_VALUE; - double distance2 = Double.MAX_VALUE; - int closestHash = 0; - - double cellularJitter = 0.39614353 * jitterModifier; - - int xPrimed = (xr - 1) * NoiseFunction.PRIME_X; - int yPrimedBase = (yr - 1) * NoiseFunction.PRIME_Y; - int zPrimedBase = (zr - 1) * NoiseFunction.PRIME_Z; - - double centerX = x; - double centerY = y; - double centerZ = z; + int xPrimed = (xr - 1) * NoiseFunction.PRIME_X; + int yPrimedBase = (yr - 1) * NoiseFunction.PRIME_Y; - for(int xi = xr - 1; xi <= xr + 1; xi++) { - int yPrimed = yPrimedBase; + double centerX = x; + double centerY = y; - for(int yi = yr - 1; yi <= yr + 1; yi++) { - int zPrimed = zPrimedBase; + for(int xi = xr - 1; xi <= xr + 1; xi++) { + int yPrimed = yPrimedBase; - for(int zi = zr - 1; zi <= zr + 1; zi++) { - int hash = HashingFunctions.hashPrimeCoords(seed, xPrimed, yPrimed, zPrimed); - int idx = hash & (255 << 2); + for(int yi = yr - 1; yi <= yr + 1; yi++) { + int hash = HashingFunctions.hashPrimeCoords(seed, xPrimed, yPrimed); + int idx = hash & (255 << 1); - double vecX = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx], cellularJitter, xi - x); - double vecY = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx | 1], cellularJitter, yi - y); - double vecZ = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx | 2], cellularJitter, zi - z); + double vecX = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_2D[idx], cellularJitter, xi - x); + double vecY = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_2D[idx | 1], cellularJitter, yi - y); double newDistance = switch(distanceFunction) { - case Euclidean, EuclideanSq -> ArithmeticFunctions.fma(vecX, vecX, - ArithmeticFunctions.fma(vecY, vecY, vecZ * vecZ)); - case Manhattan -> Math.abs(vecX) + Math.abs(vecY) + Math.abs(vecZ); - case Hybrid -> (Math.abs(vecX) + Math.abs(vecY) + Math.abs(vecZ)) + ArithmeticFunctions.fma(vecX, vecX, - ArithmeticFunctions.fma(vecY, vecY, vecZ * vecZ)); + case Euclidean, EuclideanSq -> ArithmeticFunctions.fma(vecX, vecX, vecY * vecY); + case Manhattan -> Math.abs(vecX) + Math.abs(vecY); + case Hybrid -> (Math.abs(vecX) + Math.abs(vecY)) + ArithmeticFunctions.fma(vecX, vecX, vecY * vecY); }; - distance1 = Math.max(Math.min(distance1, newDistance), distance0); - if(newDistance < distance0) { + if (newDistance < distance0) { + distance2 = distance1; + distance1 = distance0; distance0 = newDistance; closestHash = hash; - centerX = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx], cellularJitter, xi) / frequency; - centerY = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx | 1], cellularJitter, yi) / frequency; - centerZ = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx | 2], cellularJitter, zi) / frequency; - } else if(newDistance < distance1) { + centerX = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_2D[idx], cellularJitter, xi) / frequency; + centerY = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_2D[idx | 1], cellularJitter, yi) / frequency; + } else if (newDistance < distance1) { distance2 = distance1; distance1 = newDistance; - } else if(newDistance < distance2) { + } else if (newDistance < distance2) { distance2 = newDistance; } - zPrimed += NoiseFunction.PRIME_Z; + yPrimed += NoiseFunction.PRIME_Y; } - yPrimed += NoiseFunction.PRIME_Y; + xPrimed += NoiseFunction.PRIME_X; } - xPrimed += NoiseFunction.PRIME_X; - } - if(distanceFunction == DistanceFunction.Euclidean && returnType != ReturnType.CellValue) { - distance0 = Math.sqrt(distance0); + if(distanceFunction == DistanceFunction.Euclidean && returnType != ReturnType.CellValue) { + distance0 = Math.sqrt(distance0); - if(returnType != ReturnType.Distance) { - distance1 = Math.sqrt(distance1); + if(returnType != ReturnType.Distance) { + distance1 = Math.sqrt(distance1); + } } + + return switch(returnType) { + case CellValue -> closestHash * (1 / 2147483648.0); + case Distance -> distance0 - 1; + case Distance2 -> distance1 - 1; + case Distance2Add -> (distance1 + distance0) * 0.5 - 1; + case Distance2Sub -> distance1 - distance0 - 1; + case Distance2Mul -> distance1 * distance0 * 0.5 - 1; + case Distance2Div -> distance0 / distance1 - 1; + case NoiseLookup -> noiseLookup.getSample(sl - (saltLookup ? 0 : salt), centerX, centerY); + case LocalNoiseLookup -> noiseLookup.getSample(sl - (saltLookup ? 0 : salt), x / frequency - centerX, + y / frequency - centerY); + case Distance3 -> distance2 - 1; + case Distance3Add -> (distance2 + distance0) * 0.5 - 1; + case Distance3Sub -> distance2 - distance0 - 1; + case Distance3Mul -> distance2 * distance0 - 1; + case Distance3Div -> distance0 / distance2 - 1; + case Angle -> TrigonometryFunctions.fastAtan2(y / frequency - centerY, x / frequency - centerX); + }; } - return switch(returnType) { - case CellValue -> closestHash * (1 / 2147483648.0); - case Distance -> distance0 - 1; - case Distance2 -> distance1 - 1; - case Distance2Add -> (distance1 + distance0) * 0.5 - 1; - case Distance2Sub -> distance1 - distance0 - 1; - case Distance2Mul -> distance1 * distance0 * 0.5 - 1; - case Distance2Div -> distance0 / distance1 - 1; - case NoiseLookup -> noiseLookup.getSample(sl - (saltLookup ? 0 : salt), centerX, centerY, centerZ); - case LocalNoiseLookup -> noiseLookup.getSample(sl - (saltLookup ? 0 : salt), x / frequency - centerX, y / frequency - centerY, - z / frequency - centerZ); - case Distance3 -> distance2 - 1; - case Distance3Add -> (distance2 + distance0) * 0.5 - 1; - case Distance3Sub -> distance2 - distance0 - 1; - case Distance3Mul -> distance2 * distance0 - 1; - case Distance3Div -> distance0 / distance2 - 1; - case Angle -> Math.atan2(y / frequency - centerY, x / frequency - centerX); - }; + public double getNoiseRaw(long sl, double x, double y, double z) { + int seed = (int) sl; + int xr = FloatingPointFunctions.round(x); + int yr = FloatingPointFunctions.round(y); + int zr = FloatingPointFunctions.round(z); + + double distance0 = Double.MAX_VALUE; + double distance1 = Double.MAX_VALUE; + double distance2 = Double.MAX_VALUE; + int closestHash = 0; + + double cellularJitter = 0.39614353 * jitterModifier; + + int xPrimed = (xr - 1) * NoiseFunction.PRIME_X; + int yPrimedBase = (yr - 1) * NoiseFunction.PRIME_Y; + int zPrimedBase = (zr - 1) * NoiseFunction.PRIME_Z; + + double centerX = x; + double centerY = y; + double centerZ = z; + + for(int xi = xr - 1; xi <= xr + 1; xi++) { + int yPrimed = yPrimedBase; + + for(int yi = yr - 1; yi <= yr + 1; yi++) { + int zPrimed = zPrimedBase; + + for(int zi = zr - 1; zi <= zr + 1; zi++) { + int hash = HashingFunctions.hashPrimeCoords(seed, xPrimed, yPrimed, zPrimed); + int idx = hash & (255 << 2); + + double vecX = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx], cellularJitter, xi - x); + double vecY = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx | 1], cellularJitter, yi - y); + double vecZ = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx | 2], cellularJitter, zi - z); + + double newDistance = switch(distanceFunction) { + case Euclidean, EuclideanSq -> ArithmeticFunctions.fma(vecX, vecX, + ArithmeticFunctions.fma(vecY, vecY, vecZ * vecZ)); + case Manhattan -> Math.abs(vecX) + Math.abs(vecY) + Math.abs(vecZ); + case Hybrid -> (Math.abs(vecX) + Math.abs(vecY) + Math.abs(vecZ)) + ArithmeticFunctions.fma(vecX, vecX, + ArithmeticFunctions.fma(vecY, vecY, vecZ * vecZ)); + }; + + if (newDistance < distance0) { + distance2 = distance1; + distance1 = distance0; + distance0 = newDistance; + + closestHash = hash; + centerX = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx], cellularJitter, xi) / frequency; + centerY = ArithmeticFunctions.fma(CellularSampler.RAND_VECS_3D[idx | 1], cellularJitter, yi) / frequency; + } else if (newDistance < distance1) { + distance2 = distance1; + distance1 = newDistance; + } else if (newDistance < distance2) { + distance2 = newDistance; + } + zPrimed += NoiseFunction.PRIME_Z; + } + yPrimed += NoiseFunction.PRIME_Y; + } + xPrimed += NoiseFunction.PRIME_X; + } + + if(distanceFunction == DistanceFunction.Euclidean && returnType != ReturnType.CellValue) { + distance0 = Math.sqrt(distance0); + + if(returnType != ReturnType.Distance) { + distance1 = Math.sqrt(distance1); + } + } + + return switch(returnType) { + case CellValue -> closestHash * (1 / 2147483648.0); + case Distance -> distance0 - 1; + case Distance2 -> distance1 - 1; + case Distance2Add -> (distance1 + distance0) * 0.5 - 1; + case Distance2Sub -> distance1 - distance0 - 1; + case Distance2Mul -> distance1 * distance0 * 0.5 - 1; + case Distance2Div -> distance0 / distance1 - 1; + case NoiseLookup -> noiseLookup.getSample(sl - (saltLookup ? 0 : salt), centerX, centerY, centerZ); + case LocalNoiseLookup -> noiseLookup.getSample(sl - (saltLookup ? 0 : salt), x / frequency - centerX, + y / frequency - centerY, + z / frequency - centerZ); + case Distance3 -> distance2 - 1; + case Distance3Add -> (distance2 + distance0) * 0.5 - 1; + case Distance3Sub -> distance2 - distance0 - 1; + case Distance3Mul -> distance2 * distance0 - 1; + case Distance3Div -> distance0 / distance2 - 1; + case Angle -> TrigonometryFunctions.fastAtan2(y / frequency - centerY, x / frequency - centerX); + }; + } } + public enum ReturnType { CellValue, Distance, diff --git a/src/main/java/com/dfsek/seismic/math/trigonometry/TrigonometryConstants.java b/src/main/java/com/dfsek/seismic/math/trigonometry/TrigonometryConstants.java index 5cc3612..4fcb0e6 100644 --- a/src/main/java/com/dfsek/seismic/math/trigonometry/TrigonometryConstants.java +++ b/src/main/java/com/dfsek/seismic/math/trigonometry/TrigonometryConstants.java @@ -7,6 +7,7 @@ public class TrigonometryConstants { * its diameter. */ public static final double PI = Math.PI; + public static final double HALF_PI = PI / 2; /** * The {@code double} value that is closer than any other superior value to @@ -35,4 +36,15 @@ public class TrigonometryConstants { * words, tau is double pi . */ public static final double TAU_SUP = Double.longBitsToDouble(Double.doubleToRawLongBits(TrigonometryConstants.TAU) + 1); + + /** + * Polynomial defs for atan2 taken from https://mazzo.li/ + */ + public static final double a1 = 0.99997726; + public static final double a3 = -0.33262347; + public static final double a5 = 0.19354346; + public static final double a7 = -0.11643287; + public static final double a9 = 0.05265332; + public static final double a11 = -0.01172120; + } diff --git a/src/main/java/com/dfsek/seismic/math/trigonometry/TrigonometryFunctions.java b/src/main/java/com/dfsek/seismic/math/trigonometry/TrigonometryFunctions.java index 7d63ca3..125fd6d 100644 --- a/src/main/java/com/dfsek/seismic/math/trigonometry/TrigonometryFunctions.java +++ b/src/main/java/com/dfsek/seismic/math/trigonometry/TrigonometryFunctions.java @@ -1,5 +1,12 @@ package com.dfsek.seismic.math.trigonometry; +import com.dfsek.seismic.math.arithmetic.ArithmeticFunctions; +import com.dfsek.seismic.util.VMConstants; + +import static com.dfsek.seismic.math.trigonometry.TrigonometryConstants.HALF_PI; +import static com.dfsek.seismic.math.trigonometry.TrigonometryConstants.PI; + + public class TrigonometryFunctions { /** * Returns the trigonometric sine of an angle. @@ -67,4 +74,23 @@ public static double csc(double angle) { public static double cot(double angle) { return TrigonometryFunctions.cos(angle) / TrigonometryFunctions.sin(angle); } + + /** + * A fast but less precise way to implement Math.atan2 + * Taken and adapted from https://mazzo.li/ + */ + public static double fastAtan2(double y, double x) { + boolean swap = Math.abs(x) < Math.abs(y); + double atanInput = (swap ? x : y) / (swap ? y : x); + + // Approximate atan + double zSq = atanInput * atanInput; + double res = atanInput * ArithmeticFunctions.fma(zSq, ArithmeticFunctions.fma(zSq, ArithmeticFunctions.fma(zSq, ArithmeticFunctions.fma(zSq, ArithmeticFunctions.fma(zSq, TrigonometryConstants.a11, TrigonometryConstants.a9), TrigonometryConstants.a7), TrigonometryConstants.a5), TrigonometryConstants.a3), TrigonometryConstants.a1); + + res = swap ? (HALF_PI - res) : res; + res = x < 0.0 ? (PI - res) : res; + + res = Math.abs(res) * Math.signum(y); + return res; + } } diff --git a/src/test/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSamplerTest.java b/src/test/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSamplerTest.java new file mode 100644 index 0000000..5fedff9 --- /dev/null +++ b/src/test/java/com/dfsek/seismic/algorithms/sampler/noise/CellularSamplerTest.java @@ -0,0 +1,21 @@ +package com.dfsek.seismic.algorithms.sampler.noise; + +import org.junit.jupiter.api.Test; + +import static org.junit.jupiter.api.Assertions.*; + + +class CellularSamplerTest { + + @Test + void getNoiseRaw() { + NoiseFunction sampler = new CellularSampler(); + assertEquals(-0.8090170594460182, sampler.getNoiseRaw(12, 12, 456), 0.0001); + } + + @Test + void getNoiseRaw3D() { + NoiseFunction sampler = new CellularSampler(); + assertEquals(-0.8430703036518714, sampler.getNoiseRaw(0, 5674, 43, 423), 0.0001); + } +} \ No newline at end of file