diff --git a/Content/MetaPointMap.umap b/Content/MetaPointMap.umap
index 7e4733c7f9d52777e794102217da9ff21ea628ba..a90f205a3061d3582938745c1b4c3550588e57fa 100644
Binary files a/Content/MetaPointMap.umap and b/Content/MetaPointMap.umap differ
diff --git a/Source/MetaCastBachelor/DensityField.cpp b/Source/MetaCastBachelor/DensityField.cpp
index ab49740d352689c4683859bce08b5445cde2753b..6c6af6b8a896b8c791f2ee8170ca6e2ceccfb4ba 100644
--- a/Source/MetaCastBachelor/DensityField.cpp
+++ b/Source/MetaCastBachelor/DensityField.cpp
@@ -13,15 +13,13 @@ FDensityField::FDensityField(): XNum(0), YNum(0), ZNum(0), XStep(0), YStep(0), Z
void FDensityField::InitializeDensityField(const APointCloud* PointCloud, const FVector& MinBounds, const FVector& MaxBounds, const FIntVector AxisNumbers)
{
InitializeDataStructures(MinBounds, MaxBounds, AxisNumbers.X, AxisNumbers.Y, AxisNumbers.Z);
-
- constexpr float InfluenceRadius = 1.0f; // Size of a voxel
- constexpr float Sigma = 1.0f; // Standard deviation for the Gaussian kernel
- CalculateVoxelDensities(PointCloud, InfluenceRadius, Sigma);
- CalculateAndStoreGradients();
-
VoxelPointLookupTable = new FVoxelPointLookupTable();
VoxelPointLookupTable->Initialize(PointCloud->PositionVectors, this);
+
+ CalculateVoxelDensities_2(PointCloud, PointCloud->InfluenceRadius);
+ //CalculateVoxelDensitiesByNumber(PointCloud, 4);
+ CalculateAndStoreGradients();
}
void FDensityField::InitializeDataStructures(const FVector& MinBounds, const FVector& MaxBounds, const int32 XAxisNum, const int32 YAxisNum, const int32 ZAxisNum)
@@ -185,6 +183,196 @@ void FDensityField::CalculateVoxelDensitiesByNumber(const APointCloud* PointClou
UE_LOG(LogTemp, Log, TEXT("Density calculation completed in %f seconds."), ElapsedTime);
}
+/*
+void FDensityField::CalculateVoxelDensities_2(const APointCloud* PointCloud, const float InfluenceRadius)
+{
+ if (!PointCloud)
+ return;
+
+ const double StartTime = FPlatformTime::Seconds();
+ UE_LOG(LogTemp, Log, TEXT("Starting density calculation."));
+
+ TArray<float> SmoothingLengths;
+ SmoothingLengths.Init(InfluenceRadius, PointCloud->PositionVectors.Num());
+
+ // Clear existing densities
+ for (auto& Node : VoxelList)
+ {
+ Node.SetVoxelDensity(0.0);
+ Node.SetClosePointsNumber(0.0);
+ }
+
+ // Iterate over each particle to update smoothing lengths based on local densities
+ UpdateSmoothingLengths(PointCloud, SmoothingLengths, InfluenceRadius);
+
+ // Iterate over each particle for density calculation
+ for (int i = 0; i < PointCloud->PositionVectors.Num(); i++)
+ {
+ const FVector& ParticlePosition = PointCloud->PositionVectors[i];
+ const float ParticleSmoothingLength = SmoothingLengths[i];
+ const float ParticleGaussianDenominator = 2 * ParticleSmoothingLength * ParticleSmoothingLength;
+
+ const int32 StartX = FMath::Max(0, FMath::FloorToInt((ParticlePosition.X - ParticleSmoothingLength - PointCloud->MinBounds.X) / XStep));
+ const int32 EndX = FMath::Min(XNum - 1, FMath::FloorToInt((ParticlePosition.X + ParticleSmoothingLength - PointCloud->MinBounds.X) / XStep));
+ const int32 StartY = FMath::Max(0, FMath::FloorToInt((ParticlePosition.Y - ParticleSmoothingLength - PointCloud->MinBounds.Y) / YStep));
+ const int32 EndY = FMath::Min(YNum - 1, FMath::FloorToInt((ParticlePosition.Y + ParticleSmoothingLength - PointCloud->MinBounds.Y) / YStep));
+ const int32 StartZ = FMath::Max(0, FMath::FloorToInt((ParticlePosition.Z - ParticleSmoothingLength - PointCloud->MinBounds.Z) / ZStep));
+ const int32 EndZ = FMath::Min(ZNum - 1, FMath::FloorToInt((ParticlePosition.Z + ParticleSmoothingLength - PointCloud->MinBounds.Z) / ZStep));
+
+ for (int32 Z = StartZ; Z <= EndZ; ++Z)
+ for (int32 Y = StartY; Y <= EndY; ++Y)
+ for (int32 X = StartX; X <= EndX; ++X)
+ {
+ const int32 Index = GridPositionToIndex(X, Y, Z);
+ if (VoxelList.IsValidIndex(Index))
+ {
+ FVector NodePosition = VoxelList[Index].GetVoxelPosition();
+ const float Distance = FVector::Dist(NodePosition, ParticlePosition);
+ if (Distance <= ParticleSmoothingLength)
+ {
+ const float Weight = exp(-(Distance * Distance) / ParticleGaussianDenominator);
+ VoxelList[Index].AddToVoxelDensity(Weight);
+ VoxelList[Index].SetClosePointsNumber(VoxelList[Index].GetClosePointsNumber() + 1.0);
+ }
+ }
+ }
+ }
+
+ double MaxDensity = 0.0;
+ for (const auto& Node : VoxelList)
+ if (Node.GetVoxelDensity() > MaxDensity)
+ MaxDensity = Node.GetVoxelDensity();
+
+ const double EndTime = FPlatformTime::Seconds();
+ const double ElapsedTime = EndTime - StartTime;
+ UE_LOG(LogTemp, Log, TEXT("Density calculation completed in %f seconds."), ElapsedTime);
+}*/
+
+void FDensityField::SetInitialDensityEstimates(const APointCloud* PointCloud, const float InfluenceRadius)
+{
+ // Iterate over each particle
+ for (const FVector& ParticlePosition : PointCloud->PositionVectors)
+ {
+ const int32 StartX = FMath::Max(0, FMath::FloorToInt((ParticlePosition.X - InfluenceRadius - PointCloud->MinBounds.X) / XStep));
+ const int32 EndX = FMath::Min(XNum - 1, FMath::FloorToInt((ParticlePosition.X + InfluenceRadius - PointCloud->MinBounds.X) / XStep));
+ const int32 StartY = FMath::Max(0, FMath::FloorToInt((ParticlePosition.Y - InfluenceRadius - PointCloud->MinBounds.Y) / YStep));
+ const int32 EndY = FMath::Min(YNum - 1, FMath::FloorToInt((ParticlePosition.Y + InfluenceRadius - PointCloud->MinBounds.Y) / YStep));
+ const int32 StartZ = FMath::Max(0, FMath::FloorToInt((ParticlePosition.Z - InfluenceRadius - PointCloud->MinBounds.Z) / ZStep));
+ const int32 EndZ = FMath::Min(ZNum - 1, FMath::FloorToInt((ParticlePosition.Z + InfluenceRadius - PointCloud->MinBounds.Z) / ZStep));
+
+ for (int32 Z = StartZ; Z <= EndZ; ++Z)
+ {
+ for (int32 Y = StartY; Y <= EndY; ++Y)
+ {
+ for (int32 X = StartX; X <= EndX; ++X)
+ {
+ const int32 Index = GridPositionToIndex(X, Y, Z);
+ if (Index >= 0 && Index < VoxelList.Num())
+ {
+ FVector NodePosition = VoxelList[Index].GetVoxelPosition();
+ const float Distance = FVector::Dist(NodePosition, ParticlePosition);
+
+ if (Distance <= InfluenceRadius)
+ {
+ const float NormalizedDistance = Distance / InfluenceRadius;
+
+ double Weight = 1- (NormalizedDistance * NormalizedDistance);
+
+ Weight *= 0.5968310365946 / (InfluenceRadius * InfluenceRadius * InfluenceRadius);
+ VoxelList[Index].AddToVoxelDensity(Weight);
+ }
+ }
+ }
+ }
+ }
+ }
+}
+
+void FDensityField::CalculateVoxelDensities_2(const APointCloud* PointCloud, const float InfluenceRadius)
+{
+ if (!PointCloud)
+ return;
+
+ // Clear existing densities and counts
+ for (auto& Node : VoxelList)
+ {
+ Node.SetVoxelDensity(0.0);
+ Node.SetClosePointsNumber(0.0);
+ }
+
+ SetInitialDensityEstimates(PointCloud, InfluenceRadius);
+
+ TArray<float> SmoothingLengths;
+ SmoothingLengths.Init(InfluenceRadius, PointCloud->PositionVectors.Num());
+
+ // Iterate over each particle
+ int Count = 0;
+ for (const FVector& ParticlePosition : PointCloud->PositionVectors)
+ {
+ const int VoxelIndex = WorldPositionToIndex(ParticlePosition);
+ const float TotalVoxelDensity = IndexToVoxelDensity(VoxelIndex);
+ const int NumberInVoxel = VoxelPointLookupTable->GetPointsInVoxel(VoxelIndex).Num();
+ const float InterpolatedDensity = FUtilities::InterpolateDensityAtPosition(this, ParticlePosition);
+
+ const float SmoothingLength = FMath::Min(FMathf::Pow(TotalVoxelDensity / NumberInVoxel / InterpolatedDensity,(1.0 / 3.0)), 5 * GetStep().X);
+ SmoothingLengths[Count] = SmoothingLength;
+
+ Count++;
+ }
+
+ for (int i = 0; i < 100; i++)
+ {
+ UE_LOG(LogTemp, Log, TEXT("Smoothing Length: %f"), SmoothingLengths[i]);
+ }
+
+ for (auto& Node : VoxelList)
+ {
+ Node.SetVoxelDensity(0.0);
+ Node.SetClosePointsNumber(0.0);
+ }
+
+ // Iterate over each particle
+ Count = 0;
+ for (const FVector& ParticlePosition : PointCloud->PositionVectors)
+ {
+ const float NewInfluenceRadius = SmoothingLengths[Count];
+ const int32 StartX = FMath::Max(0, FMath::FloorToInt((ParticlePosition.X - NewInfluenceRadius - PointCloud->MinBounds.X) / XStep));
+ const int32 EndX = FMath::Min(XNum - 1, FMath::FloorToInt((ParticlePosition.X + NewInfluenceRadius - PointCloud->MinBounds.X) / XStep));
+ const int32 StartY = FMath::Max(0, FMath::FloorToInt((ParticlePosition.Y - NewInfluenceRadius - PointCloud->MinBounds.Y) / YStep));
+ const int32 EndY = FMath::Min(YNum - 1, FMath::FloorToInt((ParticlePosition.Y + NewInfluenceRadius - PointCloud->MinBounds.Y) / YStep));
+ const int32 StartZ = FMath::Max(0, FMath::FloorToInt((ParticlePosition.Z - NewInfluenceRadius - PointCloud->MinBounds.Z) / ZStep));
+ const int32 EndZ = FMath::Min(ZNum - 1, FMath::FloorToInt((ParticlePosition.Z + NewInfluenceRadius - PointCloud->MinBounds.Z) / ZStep));
+
+ for (int32 Z = StartZ; Z <= EndZ; ++Z)
+ {
+ for (int32 Y = StartY; Y <= EndY; ++Y)
+ {
+ for (int32 X = StartX; X <= EndX; ++X)
+ {
+ const int32 Index = GridPositionToIndex(X, Y, Z);
+ if (Index >= 0 && Index < VoxelList.Num())
+ {
+ FVector NodePosition = VoxelList[Index].GetVoxelPosition();
+ const float Distance = FVector::Dist(NodePosition, ParticlePosition);
+
+ if (Distance <= NewInfluenceRadius)
+ {
+ const float NormalizedDistance = Distance / NewInfluenceRadius;
+
+ double Weight = 1 - (NormalizedDistance * NormalizedDistance);
+
+ Weight *= (NewInfluenceRadius * NewInfluenceRadius * NewInfluenceRadius);
+ VoxelList[Index].AddToVoxelDensity(FUtilities::SmoothingKernel(Distance, NewInfluenceRadius));
+ }
+ }
+ }
+ }
+ }
+
+ Count++;
+ }
+}
+
void FDensityField::CalculateAndStoreGradients()
{
diff --git a/Source/MetaCastBachelor/DensityField.h b/Source/MetaCastBachelor/DensityField.h
index c854dcdecc9bbb3e0fbf913733f6f46cbbdee95e..5a3cd269d9704f36933c56edfb8e44e5c564eaf4 100644
--- a/Source/MetaCastBachelor/DensityField.h
+++ b/Source/MetaCastBachelor/DensityField.h
@@ -55,6 +55,8 @@ public:
// INITIALIZATION FUNCTIONS
void InitializeDensityField(const APointCloud* PointCloud, const FVector& MinBounds, const FVector& MaxBounds, const FIntVector AxisNumbers);
+ void CalculateVoxelDensities_2(const APointCloud* PointCloud, const float InfluenceRadius);
+ void SetInitialDensityEstimates(const APointCloud* PointCloud, const float InfluenceRadius);
// DEBUG FUNCTIONS
void DrawDebugVoxelDensity(const UWorld* World, const AActor* Actor, float Duration, float Scale, float DensityThreshold) const;
diff --git a/Source/MetaCastBachelor/MetaPoint.cpp b/Source/MetaCastBachelor/MetaPoint.cpp
index cde81c999f8b30ab9f8e5808b402b6e4983f317c..6b45ec3eb10b7998cb43a397de22dac85c24867f 100644
--- a/Source/MetaCastBachelor/MetaPoint.cpp
+++ b/Source/MetaCastBachelor/MetaPoint.cpp
@@ -55,9 +55,11 @@ void UMetaPoint::HandleMetaSelectReleased(const FInputActionInstance& Instance)
//ProceduralMesh->ClearAllMeshSections();
//ProceduralMesh->SetVisibility(false);
+ MyPointCloud->ColorPointsInVoxels(FloodedIndices);
+
AsyncTask(ENamedThreads::Type::AnyBackgroundHiPriTask, [this]()
{
- MyPointCloud->ColorPointsInVoxels(FloodedIndices);
+
});
}
@@ -101,14 +103,17 @@ void UMetaPoint::InitMetaPointSelection()
World = GetWorld();
LastHandInput = SelectionWorldPosition;
+ AbortMarchingCubes = true;
+ AccumulatedTime = 0.0f;
}
void UMetaPoint::SelectParticles(const FVector& InputPosition)
{
if (AccumulatedTime >= 1 / MetaCastPerSecond)
{
- AccumulatedTime = -1000;
-
+ AccumulatedTime = -10000;
+
+ FScopeLock Lock(&FloodFillingGuard);
AsyncTask(ENamedThreads::AnyBackgroundThreadNormalTask, [this]()
{
const FVector SelectionWorldPosition = SelectionObject->GetComponentLocation();
@@ -187,8 +192,14 @@ void UMetaPoint::GenerateVoxelMeshWithCubes(const TArray<int32> Voxels) const
});
}
-void UMetaPoint::GenerateVoxelMeshSmooth(TArray<int32> Voxels)
+void UMetaPoint::GenerateVoxelMeshSmooth(const TArray<int32> Voxels)
{
+ if(IsMarchingCubesRunning)
+ return;
+
+ IsMarchingCubesRunning = true;
+
+ FScopeLock Lock(&ProceduralMeshGuard);
AsyncTask(ENamedThreads::AnyBackgroundThreadNormalTask, [this, Voxels]()
{
FVector Min(FLT_MAX, FLT_MAX, FLT_MAX);
@@ -202,9 +213,13 @@ void UMetaPoint::GenerateVoxelMeshSmooth(TArray<int32> Voxels)
UE::Geometry::FMarchingCubes MarchingCubes;
MarchingCubes.Bounds = UE::Geometry::TAxisAlignedBox3(Min, Max);
- MarchingCubes.CubeSize = MyDensityField->GetStep().X;
+ MarchingCubes.CubeSize = MarchingCubeSize == 0 ? MyDensityField->GetStep().X : MarchingCubeSize;
+ MarchingCubes.CancelF = [this]() {
+ return AbortMarchingCubes;
+ };
MarchingCubes.IsoValue = 0;
+ AbortMarchingCubes = false;
// Define the implicit function
MarchingCubes.Implicit = [this, Voxels](const FVector3d& Pos) {
const FVector PosConverted = FVector(Pos.X, Pos.Y, Pos.Z);
@@ -254,12 +269,11 @@ void UMetaPoint::GenerateVoxelMeshSmooth(TArray<int32> Voxels)
AsyncTask(ENamedThreads::Type::GameThread, [this, Vertices, Indices, Normals, VertexColors]()
{
+ FScopeLock MeshLock(&ProceduralMeshGuard);
ProceduralMesh->CreateMeshSection(0, Vertices, Indices, Normals, TArray<FVector2D>(), VertexColors, TArray<FProcMeshTangent>(), false);
ProceduralMesh->SetVisibility(true);
AccumulatedTime = 0.0f;
+ IsMarchingCubesRunning = false;
});
});
-
-
-
}
\ No newline at end of file
diff --git a/Source/MetaCastBachelor/MetaPoint.h b/Source/MetaCastBachelor/MetaPoint.h
index ecb85c63abdc3ef3818c920993345c2302df7df8..20cd30ed39349630e349e66a2d68736b6683695b 100644
--- a/Source/MetaCastBachelor/MetaPoint.h
+++ b/Source/MetaCastBachelor/MetaPoint.h
@@ -24,6 +24,12 @@ class UMetaPoint : public UMetaCastBaseline
FVector LastHandInput;
UPROPERTY(EditAnywhere)
float MinimalHandMovement = 0.01;
+ bool IsMarchingCubesRunning = false;
+ bool AbortMarchingCubes = false;
+
+ // Critical section to protect shared variables
+ mutable FCriticalSection FloodFillingGuard;
+ mutable FCriticalSection ProceduralMeshGuard;
UPROPERTY()
UProceduralMeshComponent* ProceduralMesh;
@@ -36,8 +42,11 @@ class UMetaPoint : public UMetaCastBaseline
TArray<int32> FloodedIndices;
TArray<int32> RevisitIndices;
+
UPROPERTY(EditAnywhere)
- int MetaCastPerSecond = 1;
+ float MarchingCubeSize = 0;
+ UPROPERTY(EditAnywhere)
+ int MetaCastPerSecond = 10;
public:
UMetaPoint();
@@ -55,5 +64,4 @@ public:
void GenerateVoxelMeshSmooth(const TArray<int32> Voxels);
virtual void TickComponent(float DeltaTime, ELevelTick TickType, FActorComponentTickFunction* ThisTickFunction) override;
-};
-
+};
\ No newline at end of file
diff --git a/Source/MetaCastBachelor/PointCloud.cpp b/Source/MetaCastBachelor/PointCloud.cpp
index d319ce91f39e29ecdb9443e723229bac0bd796cc..aad4992d3ac917b71487bc944061f228b45301f8 100644
--- a/Source/MetaCastBachelor/PointCloud.cpp
+++ b/Source/MetaCastBachelor/PointCloud.cpp
@@ -134,7 +134,7 @@ void APointCloud::DrawVoxel(const int Index, const float Time) const
}
-void APointCloud::ColorPointsInVoxels(const TArray<int32> VoxelIndices)
+void APointCloud::ColorPointsInVoxels(const TArray<int32>& VoxelIndices)
{
FScopeLock Lock(&DataGuard);
if (!MyDensityField)
diff --git a/Source/MetaCastBachelor/PointCloud.h b/Source/MetaCastBachelor/PointCloud.h
index 6aae1269053c7110a6e69381a2e2bea7bb3600cf..90e89d217240ff7301b4135764d0527b2c1df935 100644
--- a/Source/MetaCastBachelor/PointCloud.h
+++ b/Source/MetaCastBachelor/PointCloud.h
@@ -12,6 +12,9 @@ class APointCloud : public AActor
GENERATED_BODY()
public:
+ UPROPERTY(EditAnywhere)
+ float InfluenceRadius = 3.0f;
+
FDensityField* MyDensityField;
UPROPERTY()
TArray<FVector> PositionVectors;
@@ -64,7 +67,7 @@ public:
void UpdateSelection();
void DrawVoxel(const int Index, float Time) const;
- void ColorPointsInVoxels(const TArray<int32> VoxelIndices);
+ void ColorPointsInVoxels(const TArray<int32>& VoxelIndices);
UFUNCTION(BlueprintCallable)
void ReadPointCloudFromFile(FFilePath FileNamePoints, FFilePath FileNameFlags);
diff --git a/Source/MetaCastBachelor/Utilities.cpp b/Source/MetaCastBachelor/Utilities.cpp
index 91ab449b6ad2f9f62b3a0a4295633ab4c3af2e30..76b42e35cda7a3a9d9cb502a09d6f3188a2160dd 100644
--- a/Source/MetaCastBachelor/Utilities.cpp
+++ b/Source/MetaCastBachelor/Utilities.cpp
@@ -116,7 +116,7 @@ FVector FUtilities::FollowGradientToMaximum(const FDensityField* DensityField, c
int Iterations = 0;
constexpr int MaxIterations = 100; // Prevent infinite loops
- constexpr float StepSize = 0.005f; // Define a reasonable step size
+ constexpr float StepSize = 0.05f; // Define a reasonable step size
while (Iterations < MaxIterations)
{
@@ -250,7 +250,7 @@ TArray<int32> FUtilities::FloodFilling_2(const FDensityField* DensityField, cons
}else
{
- const float Multiplier = 1 + FMath::Pow(DistanceToStart - MaxDistance, 2);
+ const float Multiplier = 1 + FMath::Pow(DistanceToStart - MaxDistance, 1.1);
if(NeighborDensity > Threshold * Multiplier)
{
Stack.Push(NeighborIndex);