# HG changeset patch
# User Neil Dizon <neil.dizon@helsinki.fi>
# Date 1713696095 -10800
# Node ID befb8d5125cd2094ecfa6c3942978f6ca51a82c9
# Parent  c286925c0f35fcb536ed060dcd1418cb30b26193
added stable interval for PET

diff -r c286925c0f35 -r befb8d5125cd src/PET/ImGenerate.jl
--- a/src/PET/ImGenerate.jl	Sun Apr 21 12:46:04 2024 +0300
+++ b/src/PET/ImGenerate.jl	Sun Apr 21 13:41:35 2024 +0300
@@ -181,7 +181,7 @@
                               params :: NamedTuple)
     
     # Restart the seed to enable comparison across predictors                       
-    Random.seed!(rng,9182737465)   
+    Random.seed!(rng,314159)   
 
     nextv = shake(params)
     v_true = nextv()
@@ -217,20 +217,25 @@
 # PETscan
 ########################################################################
 function generate_radon(im, sz,
-    :: Type{DisplacementConstant},
-    datachannel :: Channel{PetOnlineData{DisplacementConstant}},
-    params :: NamedTuple)
+                        :: Type{DisplacementConstant},
+                        datachannel :: Channel{PetOnlineData{DisplacementConstant}},
+                        params :: NamedTuple)
+
+    # Set up counter and zero factor for stabilisation
+    @assert(params.maxiter ≥ maximum(params.stable_interval))
+    indx = 1
+    zero_factor = indx in params.stable_interval ? 0.0 : 1.0
 
     # Restart the seed to enable comparison across predictors                       
-    Random.seed!(rng,9182737465)  
+    Random.seed!(rng,9182737465)
 
     nextv = shake(params)
-    v_true = nextv()
+    v_true = zero_factor.*nextv()
     v_cumul = copy(v_true)
 
     S_true = generate_radonmask(params)
 
-    theta_true = rotatebytheta(params)
+    theta_true = zero_factor*rotatebytheta(params)
     theta_cumul = copy(theta_true)
 
     while true
@@ -261,13 +266,16 @@
         end
 
         # Next step shake
-        v_true = nextv()
+        v_true = zero_factor.*nextv()
         v_cumul .+= v_true
         # Next theta
-        theta_true = rotatebytheta(params)
+        theta_true = zero_factor*rotatebytheta(params)
         theta_cumul += theta_true
         # Next sinogram mask 
         S_true = generate_radonmask(params)
+        # Update indx and zero factor
+        indx += 1
+        zero_factor = indx in params.stable_interval ? 0.0 : 1.0
     end
 end
 
diff -r c286925c0f35 -r befb8d5125cd src/PET/PET.jl
--- a/src/PET/PET.jl	Sun Apr 21 12:46:04 2024 +0300
+++ b/src/PET/PET.jl	Sun Apr 21 13:41:35 2024 +0300
@@ -102,7 +102,7 @@
     shake_noise_level = 0.1,
     shake = 1.0, 
     rotation_factor = 0.075,                  
-    rotation_noise_level = 0.0075,         
+    rotation_noise_level = 0.00075,         
     α = 1.0,
     ρ̃₀ = 1.0,
     σ̃₀ = 1.0,
@@ -111,13 +111,14 @@
     τ₀ = 0.9,
     λ = 1,
     origsize = 256,             
-    radondims = [128, 64],
+    radondims = [128,64],
     sz = (256,256),
     scale = 1,
     c = 1.0,
     sino_sparsity = 0.5,
     L = 300.0,
     L_experiment = false,
+    stable_interval = union(Set(1500:200),Set(3500:4000)),
 )
 
 const petscan = imgen_shepplogan_radon(p_known₀_pet.sz)