Merge pull request #22 from warisa-r/ina_test_branch

Noise analysis from Ina

analysis/NoiseAnalysis.jl

@@ -0,0 +1,13 @@
+module NoiseAnalysis
+
+using Distributions
+using LinearAlgebra
+using ShockwaveDetection #: FlowData
+
+export NoiseData
+export apply_noise, apply_noise_to_flow
+export compare_shock_clusters_over_time, compare_shock_positions_over_time, compare_shock_fits_over_time
+
+include("noise.jl")
+
+end

analysis/examples/sod_shock_tube_1d.jl

@@ -0,0 +1,24 @@
+using Euler2D
+using LinearAlgebra
+using ShockwaveProperties
+using Unitful
+using ShockwaveDetection
+include("../NoiseAnalysis.jl")
+using .NoiseAnalysis
+using Distributions
+
+# Create a NoiseData instance
+noise_data = NoiseAnalysis.NoiseData(0.01, Normal(0, 0.2))  #Gaussian distribution
+
+flow_data = FlowData("examples/data/sod_shock_left_1d.tape", false)
+
+# Apply noise to flow data
+noisy_flow_data = NoiseAnalysis.apply_noise_to_flow(flow_data, noise_data)
+
+point_detect_algo = GradientShockDetectionAlgo(0.8)
+
+original_shock_positions_over_time = detect(flow_data, point_detect_algo)
+noisy_shock_positions_over_time = detect(noisy_flow_data, point_detect_algo)
+
+# Detect and compare shock positions with and without noise
+shock_diff = NoiseAnalysis.compare_shock_positions_over_time_1d(original_shock_positions_over_time, noisy_shock_positions_over_time; threshold=10.0)

analysis/examples/sod_shock_tube_2d.jl

@@ -0,0 +1,23 @@
+using ShockwaveDetection
+include("../NoiseAnalysis.jl")
+using .NoiseAnalysis
+using Distributions
+ENV["JULIA_NUM_THREADS"] = "4"
+
+flow_data = FlowData("examples/data/sod_shock_right_2d.tape", false)
+noise_data = NoiseAnalysis.NoiseData(0.01, Normal(0, 0.1))
+noisy_flow_data = NoiseAnalysis.apply_noise_to_flow(flow_data, noise_data)
+
+point_detect_algo = ImageProcessingShockDetectionAlgo(0.8, :prewitt)
+fitting_algo = FittingAlgo(0.1, true)
+dbscan_algo = DBSCANAlgo(5.95, 3, 10)
+
+# Detect shock points, clusters, and fits for original and noisy data
+original_result = detect(flow_data, point_detect_algo, dbscan_algo, fitting_algo)
+noisy_result = detect(noisy_flow_data, point_detect_algo, dbscan_algo, fitting_algo)
+
+rmse_value = NoiseAnalysis.compare_shock_fits_over_time(original_result, noisy_result, flow_data; threshold=10.0)
+println("RMSE between fits, cluster_fit_mismatch_count, empty_residuals_count, shock_data_missing_count, residual_length_mismatch_count: ", rmse_value)
+
+#plot_shock_fits_over_time(flow_data, original_detection, true)
+#plot_shock_fits_over_time(noisy_flow_data, noisy_detection, true)