read_mesons supports legacy

compatible with old versions mesons

src/juobs_reader.jl

@@ -7,75 +7,87 @@ function read_global_header(path::String)
     return a
 end
 
-function read_CHeader(path::String)
+function read_CHeader(path::String; legacy::Bool=false)
     gh = read_global_header(path)
     
     data = open(path, "r")
     seek(data, gh.hsize)
     
-    aux_f = zeros(Float64, 6)
-    aux_i = zeros(Int32, 4)
-    theta = zeros(Float64, 6)
     a = Vector{CHeader}(undef, gh.ncorr)
-    for k = 1:gh.ncorr
-        read!(data, aux_f)
-        read!(data, theta)
-
-        qs1 = read(data, Int32)
-        if qs1 != 0
-            qn1 = read(data, Int32)
-            qeps1 = read(data, Float64)
-            q1 = Sm(qs1, qn1, qeps1, 1)
-        else
-            q1 = Sm(qs1, 1)
-        end
+    if !legacy
+        aux_f = zeros(Float64, 6)
+        aux_i = zeros(Int32, 4)
+        theta = zeros(Float64, 6)
+
+        for k = 1:gh.ncorr
+            read!(data, aux_f)
+            read!(data, theta)
+
+            qs1 = read(data, Int32)
+            if qs1 != 0
+                qn1 = read(data, Int32)
+                qeps1 = read(data, Float64)
+                q1 = Sm(qs1, qn1, qeps1, 1)
+            else
+                q1 = Sm(qs1, 1)
+            end
 
-        qs2 = read(data, Int32)
-        if qs2 != 0
-            qn2 = read(data, Int32)
-            qeps2 = read(data, Float64)
-            q2 = Sm(qs2, qn2, qeps2, 1)
-        else
-            q2 = Sm(qs2, 1)
-        end
+            qs2 = read(data, Int32)
+            if qs2 != 0
+                qn2 = read(data, Int32)
+                qeps2 = read(data, Float64)
+                q2 = Sm(qs2, qn2, qeps2, 1)
+            else
+                q2 = Sm(qs2, 1)
+            end
 
 
-        gs1 = read(data, Int32)
-        if gs1 != 0 && gs1 != 3 && gs1 != 4
-            gn1 = read(data, Int32)
-            geps1 = read(data, Float64)
-            g1 = Sm(gs1, gn1, geps1, 2)
-        elseif gs1 == 3 || gs1 == 4
-            g1 = Sm(gs1, q1.niter, q1.eps, 2)
-        else
-            g1 = Sm(gs1, 2)
+            gs1 = read(data, Int32)
+            if gs1 != 0 && gs1 != 3 && gs1 != 4
+                gn1 = read(data, Int32)
+                geps1 = read(data, Float64)
+                g1 = Sm(gs1, gn1, geps1, 2)
+            elseif gs1 == 3 || gs1 == 4
+                g1 = Sm(gs1, q1.niter, q1.eps, 2)
+            else
+                g1 = Sm(gs1, 2)
+            end
+            
+            gs2 = read(data, Int32)
+            if gs2 != 0 && gs2 != 3 && gs2 != 4
+                gn2 = read(data, Int32)
+                geps2 = read(data, Float64)
+                g2 = Sm(gs2, gn2, geps2, 2)
+            elseif gs1 == 3 || gs1 == 4
+                g2 = Sm(gs2, q2.niter, q2.eps, 2)
+            else
+                g2 = Sm(gs2, 2)
+            end
+            
+            read!(data, aux_i)
+            a[k] = CHeader(aux_f, aux_i, theta, [q1, q2, g1, g2])
         end
-        
-        gs2 = read(data, Int32)
-        if gs2 != 0 && gs2 != 3 && gs2 != 4
-            gn2 = read(data, Int32)
-            geps2 = read(data, Float64)
-            g2 = Sm(gs2, gn2, geps2, 2)
-        elseif gs1 == 3 || gs1 == 4
-            g2 = Sm(gs2, q2.niter, q2.eps, 2)
-        else
-            g2 = Sm(gs2, 2)
+    else
+        aux_f = zeros(Float64, 4)
+        aux_i = zeros(Int32, 4)
+        for k = 1:gh.ncorr
+            read!(data, aux_f)
+            read!(data, aux_i)
+            a[k] = CHeader(aux_f, aux_i)
         end
-        
-        read!(data, aux_i)
-        a[k] = CHeader(aux_f, aux_i, theta, [q1, q2, g1, g2])
     end
     close(data)
     return a
 end
 
 @doc raw"""
-read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing)
+read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)
 
 This faction read a mesons dat file at a given path and returns a vector of CData structures for different masses and Dirac structures.
 Dirac structures g1 and/or g2 can be passed as string arguments in order to filter correaltors.
 ADerrors id can be specified as argument. If is not specified, the id is fixed according to the ensemble name (example: "H400"-> id = 400)
 
+*For the old version (without smearing, distance preconditioning and theta) set legacy=true
 Examples:
 ```@example
 read_mesons(path)
@@ -83,12 +95,13 @@ read_mesons(path, "G5")
 read_mesons(path, nothing, "G5")
 read_mesons(path, "G5", "G5")
 read_mesons(path, "G5", "G5", id=1)
+read_mesons(path, "G5_d2", "G5_d2", legacy=true)
+
 ```
 """
-function read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing)
-    
-    isnothing(g1) ? t1=nothing : t1 = findfirst(x-> x==g1, gamma_name) - 1
-    isnothing(g2) ? t2=nothing : t2 = findfirst(x-> x==g2, gamma_name) - 1
+function read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union{String, Nothing}=nothing; id::Union{Int64, Nothing}=nothing, legacy::Bool=false)    
+    t1 = isnothing(g1) ? nothing : findfirst(x-> x==g1, gamma_name) - 1
+    t2 = isnothing(g2) ? nothing : findfirst(x-> x==g2, gamma_name) - 1
     if isnothing(id)
         bname = basename(path)
         m = findfirst(r"[A-Z][0-9]{3}r[0-9]{3}", bname)
@@ -97,7 +110,7 @@ function read_mesons(path::String, g1::Union{String, Nothing}=nothing, g2::Union
 
     data = open(path, "r")
     g_header = read_global_header(path)
-    c_header = read_CHeader(path)
+    c_header = read_CHeader(path, legacy=legacy)
 
     ncorr = g_header.ncorr
     tvals = g_header.tvals

src/juobs_types.jl

@@ -12,7 +12,13 @@
 const noise_name=["Z2", "GAUSS", "U1", "Z4"]
 const gamma_name = ["G0", "G1", "G2", "G3",
 "invalid", "G5", "1", "G0G1", "G0G2", "G0G3",
-"G0G5", "G1G2", "G1G3", "G1G5", "G2G3", "G2G5", "G3G5"]
+"G0G5", "G1G2", "G1G3", "G1G5", "G2G3", "G2G5", "G3G5", 
+"G0_d1", "G1_d1", "G2_d1", "G3_d1",
+"invalid", "G5_d1", "1_d1", "G0G1_d1", "G0G2_d1", "G0G3_d1",
+"G0G5_d1", "G1G2_d1", "G1G3_d1", "G1G5_d1", "G2G3_d1", "G2G5_d1", "G3G5_d1",
+"G0_d2", "G1_d2", "G2_d2", "G3_d2",
+"invalid", "G5_d2", "1_d2", "G0G1_d2", "G0G2_d2", "G0G3_d2",
+"G0G5_d2", "G1G2_d2", "G1G3_d2", "G1G5_d2", "G2G3_d2", "G2G5_d2", "G3G5_d2"]
 const qs_name=["Local", "Wuppertal", "3D Gradient Flow", "Gradient Flow"]
 const gs_name=["Local", "APE", "3D Wilson Flow", "Quark 3D Gradient Flow", "Quark Gradient Flow"]
 
@@ -96,7 +102,30 @@ mutable struct CHeader
         
         return a
     end
+    function CHeader(aux_f::Vector{Float64}, aux_i::Vector{Int32})
+        a = new()
+        a.k1 = aux_f[1]
+        a.k2 = aux_f[2]
+        a.mu1 = aux_f[3]
+        a.mu2 = aux_f[4]
+        a.dp1 = 0.0
+        a.dp2 = 0.0
+        a.type1 = aux_i[1]
+        a.type2 = aux_i[2]
+        a.x0 = aux_i[3]
+        a.is_real = aux_i[4]
+
+        a.theta1 = zeros(3)
+        a.theta2 = zeros(3)
+        a.q1 = Sm(0, 1)
+        a.q2 = Sm(0, 1)
+        a.g1 = Sm(0, 2)
+        a.g2 = Sm(0, 2)
+        a.hsize = 8*4 + 4*4
+        a.dsize = 16 - 8* a.is_real 
 
+        return a
+    end
 end
 
 mutable struct CData