Added directional blur for final image

SDFMapCreator/Program.cs

@@ -41,7 +41,6 @@ public class Program {
         Images.Add(new(pixels, pixels.GetLength(0), pixels.GetLength(1)));
         */
         
-
         //check if all the images in Images are the same resolution
         if (Images.Select(img => (img.Width, img.Height)).Distinct().Count() > 1) {
             Console.WriteLine("Error: Not all images have the same resolution.");
@@ -122,7 +121,22 @@ public class Program {
             }
             currStep += stepIncrement;
         }
+        
         finalImage.SaveImage("final.png");
+
+        
+        // apply directional blur
+        var iterations = 10;
+        var finalImageMask = new Vector3[width, height];
+        finalImageMask = SelfMask(Images[^1].Pixels, width, height);
+        
+        for (int i = 0; i < iterations; i++){
+            Console.WriteLine($"Applying directional blur {i}...");
+            var radius = 2f;
+            finalImage = DirectionalBlur(finalImage, finalImageMask, radius, 0.5f, 1f);
+        }
+        
+        finalImage.SaveImage("finalBlur.png");
         Console.WriteLine("Done!");
     }
 
@@ -247,7 +261,6 @@ public class Program {
     private static float EuclideanDistance(Vector2 a, Vector2 b) =>
         MathF.Sqrt(MathF.Pow(a.X - b.X, 2) + MathF.Pow(a.Y - b.Y, 2));
     
-
     static Vector3[,] GetABMask(Vector3[,] A, Vector3[,] B, uint resX, uint resY) {
         var temp = new Vector3[resX, resY];
         Parallel.For(0, resX * resY, parallelOptions, (i) =>
@@ -296,10 +309,89 @@ public class Program {
         return false;
     }
 
+    static Vector3 SampleBilinear(Vector3[,] image, Vector2 pixelCoord) {
+        var width = (uint)image.GetLength(0);
+        var height = (uint)image.GetLength(1);
+
+        var x = pixelCoord.X;
+        var y = pixelCoord.Y;
+        
+        var x0 = (int)x;
+        var y0 = (int)y;
+        
+        var x1 = Math.Min(x0 + 1, width - 1);
+        var y1 = Math.Min(y0 + 1, height - 1);
+        
+        var dx = x - x0;
+        var dy = y - y0;
+        
+        return Lerp(
+            Lerp(image[x0, y0], image[x1, y0], dx),
+            Lerp(image[x0, y1], image[x1, y1], dx),
+            dy
+        );
+    }
+    
+    static Vector3[,] DirectionalBlur(Vector3[,] image, Vector3[,] mask, float radius = 3f, float step = .5f, float sigma = 1f)
+    {
+        var width = (uint)image.GetLength(0);
+        var height = (uint)image.GetLength(1);
+        var output = new Vector3[width, height];
+        
+        Parallel.For(0, (int)width * height, parallelOptions, (i) =>
+        {
+            var x = i % width;
+            var y = i / width;
+            var gradient = new Vector2(0, 0);
+            
+            if(mask[x, y].X == 0)
+            {
+                output[x, y] = image[x, y];
+                return;
+            }
+            
+            for (int dx = -1; dx <= 1; dx++)
+            {
+                for (int dy = -1; dy <= 1; dy++)
+                {
+                    if (x + dx < 0 || x + dx >= width || y + dy < 0 || y + dy >= height) continue;
+                    gradient += new Vector2(dx, dy) * image[x + dx, y + dy].X;
+                }
+            }
+
+            if (gradient == Vector2.Zero)
+            {
+                output[x, y] = image[x, y];
+                return;
+            }
+
+            gradient = Vector2.Normalize(gradient);
+            float sum = 0;
+            
+            // now we follow the direction line and sample the image for length;
+            for(float l = -radius; l <= radius; l += step)
+            {
+                var pixel = new Vector2(x, y) + gradient * l;
+                // check if the pixel is out of bounds
+                if (pixel.X < 0 || pixel.X >= width || pixel.Y < 0 || pixel.Y >= height) continue;
+                var sample = SampleBilinear(image, pixel);
+                float weight = (float)Math.Exp(-(l * l) / (2 * sigma * sigma));
+                output[x, y] += sample * weight;
+                sum += weight;
+            }
+            
+            output[x, y] /= sum;
+        });
+        
+        return output;
+    }
+    
     static T Lerp<T>(T a, T b, float t)
         where T : INumber<T>, IMultiplyOperators<T, float, T>, IAdditionOperators<T, T, T>
         => a * (1 - t) + b * t;
     
+    static Vector3 Lerp(Vector3 a, Vector3 b, float t) => a * (1 - t) + b * t;
+
     static T Remap<T>(T value, T min, T max, T newMin, T newMax)
         where T : INumber<T>, ISubtractionOperators<T, T, T>, IMultiplyOperators<T, T, T>, IAdditionOperators<T, T, T>
         => (value - min) / (max - min) * (newMax - newMin) + newMin;