update test

Examples/Tests/diff_lumi_diag/analysis.py

@@ -5,15 +5,20 @@
 # In that case, the differential luminosity can be calculated analytically.
 
 import os
+import re
 
 import numpy as np
-from read_raw_data import read_reduced_diags_histogram
+from openpmd_viewer import OpenPMDTimeSeries
 
-# Extract the differential luminosity from the file
-_, _, E_bin, bin_data = read_reduced_diags_histogram(
-    "./diags/reducedfiles/DifferentialLuminosity_beam1_beam2.txt"
-)
-dL_dE_sim = bin_data[-1]  # Differential luminosity at the end of the simulation
+# Extract the 1D differential luminosity from the file
+filename = "./diags/reducedfiles/DifferentialLuminosity_beam1_beam2.txt"
+with open(filename) as f:
+    # First line: header, contains the energies
+    line = f.readline()
+    E_bin = np.array(list(map(float, re.findall("=(.*?)\(", line))))
+data = np.loadtxt(filename)
+dE_bin = E_bin[1] - E_bin[0]
+dL_dE_sim = data[-1, 2:]  # Differential luminosity at the end of the simulation
 
 # Beam parameters
 N = 1.2e10
@@ -33,21 +38,38 @@
     * np.exp(-((E_bin - 2 * E_beam) ** 2) / (2 * sigma_E**2))
 )
 
+# Extract the 2D differential luminosity from the file
+series = OpenPMDTimeSeries("./diags/reducedfiles/DifferentialLuminosity2d_beam1_beam2/")
+d2L_dE1_dE2, info = series.get_field("d2L_dE1_dE2", iteration=80)
+E1, E2 = info.E1, info.E2
+dE1, dE2 = info.dE1, info.dE2
+
 # Extract test name from path
 test_name = os.path.split(os.getcwd())[1]
 print("test_name", test_name)
 
 # Pick tolerance
 if "leptons" in test_name:
-    tol = 1e-2
+    tol1 = 1e-2
+    tol2 = 1e-5
 elif "photons" in test_name:
     # In the photons case, the particles are
     # initialized from a density distribution ;
     # tolerance is larger due to lower particle statistics
-    tol = 6e-2
+    tol1 = 6e-2
+    tol2 = 1e-5
+
+# Check that the 1D diagnostic and analytical result match
+error1 = abs(dL_dE_sim - dL_dE_th).max() / abs(dL_dE_th).max()
+print("Relative error: ", error1)
+print("Tolerance: ", tol1)
+
+# Check that the 2D and 1D diagnostics match
+error2 = abs(np.sum(d2L_dE1_dE2) * dE2 * dE1 - np.sum(dL_dE_sim) * dE_bin) / abs(
+    np.sum(dL_dE_sim) * dE_bin
+)
+print("Relative error: ", error2)
+print("Tolerance: ", tol2)
 
-# Check that the simulation result and analytical result match
-error = abs(dL_dE_sim - dL_dE_th).max() / abs(dL_dE_th).max()
-print("Relative error: ", error)
-print("Tolerance: ", tol)
-assert error < tol
+assert error1 < tol1
+assert error2 < tol2

Examples/Tests/diff_lumi_diag/inputs_base_3d

@@ -107,8 +107,8 @@ DifferentialLuminosity2d_beam1_beam2.type = DifferentialLuminosity2D
 DifferentialLuminosity2d_beam1_beam2.intervals = 80
 DifferentialLuminosity2d_beam1_beam2.species = beam1 beam2
 DifferentialLuminosity2d_beam1_beam2.bin_number_1 = 128
-DifferentialLuminosity2d_beam1_beam2.bin_max_1 = 1.05*beam_energy_eV
+DifferentialLuminosity2d_beam1_beam2.bin_max_1 = 1.45*beam_energy_eV
 DifferentialLuminosity2d_beam1_beam2.bin_min_1 = 0
 DifferentialLuminosity2d_beam1_beam2.bin_number_2 = 128
-DifferentialLuminosity2d_beam1_beam2.bin_max_2 = 1.05*beam_energy_eV
+DifferentialLuminosity2d_beam1_beam2.bin_max_2 = 1.45*beam_energy_eV
 DifferentialLuminosity2d_beam1_beam2.bin_min_2 = 0