rwf rework

<W> is a primary observable

src/juobs_obs.jl

@@ -487,16 +487,25 @@ function comp_t0(Y::YData, plat::Vector{Int64}; L::Int64, pl::Bool=false,
         end
     end
 
-    Ysl = isnothing(rw) ? Ysl : apply_rw(Ysl, rw)
     xmax = size(Ysl, 2)
     nt0 = t0_guess(t, Ysl, plat, L)
     
     dt0 = iseven(npol) ? Int64(npol / 2) : Int64((npol+1)/ 2)
     Y_aux = Matrix{uwreal}(undef, xmax, 2*dt0+1)
+    
+    if !isnothing(rw)
+        Ysl_r, W = apply_rw(Ysl, rw)
+        W_obs = uwreal(W, id)
+    end
+
     for i = 1:xmax
         k = 1
         for j = nt0-dt0:nt0+dt0
-            Y_aux[i, k] = uwreal(Ysl[:, i, j], id)
+            if isnothing(rw)
+                Y_aux[i, k] = uwreal(Ysl[:, i, j], id)
+            else
+                Y_aux[i, k] = uwreal(Ysl_r[:, i, j], id) / W_obs
+            end
             k = k + 1
         end
     end
@@ -577,18 +586,30 @@ function comp_t0(Y::Vector{YData}, plat::Vector{Int64}; L::Int64, pl::Bool=false
         end
     end
 
-    Ysl = isnothing(rw) ? Ysl : apply_rw(Ysl, rw)
     tmp = Ysl[1]
-    [tmp = cat(tmp, Ysl[k], dims=1) for k = 2:nr]
+    [tmp = cat(tmp, Ysl[k], dims=1) for k=2:nr]
     nt0 = t0_guess(t, tmp, plat, L)
     xmax = size(tmp, 2)
 
     dt0 = iseven(npol) ? Int64(npol / 2) : Int64((npol+1) / 2)
     Y_aux = Matrix{uwreal}(undef, xmax, 2*dt0+1)
+
+    if !isnothing(rw)
+        Ysl_r, W = apply_rw(Ysl, rw)
+        tmp_r = Ysl_r[1]
+        tmp_W = W[1]
+        [tmp_r = cat(tmp_r, Ysl_r[k], dims=1) for k = 2:nr]
+        [tmp_W = cat(tmp_W, W[k], dims=1) for k = 2:nr]
+        W_obs = uwreal(tmp_W, id, replica)
+    end
     for i = 1:xmax
         k = 1
         for j = nt0-dt0:nt0+dt0
-            Y_aux[i, k] = uwreal(tmp[:, i, j], id, replica)
+            if isnothing(rw)
+                Y_aux[i, k] = uwreal(tmp[:, i, j], id, replica)
+            else
+                Y_aux[i, k] = uwreal(tmp_r[:, i, j], id, replica) / W_obs
+            end
             k = k + 1
         end
     end

src/juobs_tools.jl

@@ -4,8 +4,8 @@ function apply_rw(data::Array{Float64}, W::Matrix{Float64})
     W1 = W[1, 1:nc]
     W2 = W[2, 1:nc]
 
-    data_r = data .* W1 .* W2 / mean(W1 .* W2)
-    return data_r
+    data_r = data .* W1 .* W2
+    return (data_r, W1 .* W2)
 end
 
 function apply_rw(data::Vector{<:Array{Float64}}, W::Vector{Matrix{Float64}})
@@ -16,17 +16,9 @@ function apply_rw(data::Vector{<:Array{Float64}}, W::Vector{Matrix{Float64}})
 
     rw1 = [W[k][1, 1:nc[k]] for k=1:length(W)]
     rw2 = [W[k][2, 1:nc[k]] for k=1:length(W)]
-    rw1_cat = rw1[1]
-    rw2_cat = rw2[1]
-
-    for k = 2:length(W)
-        rw1_cat = cat(rw1_cat, rw1[k], dims=1)
-        rw2_cat = cat(rw2_cat, rw2[k], dims=1)
-    end
-    
-    rw_mean = mean(rw1_cat .* rw2_cat)
-    data_r = [data[k] .* rw1[k].* rw2[k] / rw_mean for k=1:length(data)]
-    return data_r
+    rw = [rw1[k] .* rw2[k] for k=1:length(W)]
+    data_r = [data[k] .* rw[k] for k=1:length(data)]
+    return (data_r, rw)
 end
 
 function check_corr_der(obs::Corr, derm::Vector{Corr})   
@@ -92,11 +84,16 @@ corr_pp_r = corr_obs.(data, rw=[rw1, rw2])
 function corr_obs(cdata::CData; real::Bool=true, rw::Union{Array{Float64, 2}, Nothing}=nothing, L::Int64=1)
 
     real ? data = cdata.re_data ./ L^3 : data = cdata.im_data ./ L^3
-    data_r = isnothing(rw) ? data : apply_rw(data, rw)
-
     nt = size(data)[2]
     obs = Vector{uwreal}(undef, nt)
-    [obs[x0] = uwreal(data_r[:, x0], cdata.id) for x0 = 1:nt]
+
+    if isnothing(rw)
+        [obs[x0] = uwreal(data[:, x0], cdata.id) for x0 = 1:nt]
+    else
+        data_r, W = apply_rw(data, rw)
+        W_obs = uwreal(W, cdata.id)
+        [obs[x0] = uwreal(data_r[:, x0], cdata.id) / W_obs for x0 = 1:nt]
+    end
     return Corr(obs, cdata)
 end
 
@@ -113,15 +110,24 @@ function corr_obs(cdata::Array{CData, 1}; real::Bool=true, rw::Union{Array{Array
     end
 
     real ? data = getfield.(cdata, :re_data) ./ L^3 : data = getfield.(cdata, :im_data) ./ L^3
-    data_r = isnothing(rw) ? data : apply_rw(data, rw)
-
-    tmp = data_r[1]
-    [tmp = cat(tmp, data_r[k], dims=1) for k = 2:nr]
 
     nt = size(data[1])[2]
     obs = Vector{uwreal}(undef, nt)
 
-    [obs[x0] = uwreal(tmp[:, x0], id[1], replica) for x0 = 1:nt]
+    if isnothing(rw)
+        tmp = data[1]
+        [tmp = cat(tmp, data[k], dims=1) for k = 2:nr]
+        [obs[x0] = uwreal(tmp[:, x0], id[1], replica) for x0 = 1:nt]
+    else
+        data_r, W = apply_rw(data, rw)
+        tmp = data_r[1]
+        tmp_W = W[1]
+        [tmp = cat(tmp, data_r[k], dims=1) for k = 2:nr]
+        [tmp_W = cat(tmp_W, W[k], dims=1) for k = 2:nr]
+
+        W_obs = uwreal(tmp_W, id[1], replica)
+        [obs[x0] = uwreal(tmp[:, x0], id[1], replica) / W_obs for x0 = 1:nt]
+    end
 
     return Corr(obs, cdata)
 end