diff --git a/n3fit/src/n3fit/scripts/evolven3fit.py b/n3fit/src/n3fit/scripts/evolven3fit.py
index c9b506b4e5..277de17858 100644
--- a/n3fit/src/n3fit/scripts/evolven3fit.py
+++ b/n3fit/src/n3fit/scripts/evolven3fit.py
@@ -7,11 +7,13 @@
 import pathlib
 import sys
 
+from evolven3fit import cli, eko_utils, evolve, utils
 import numpy as np
 
 from eko.runner.managed import solve
-from evolven3fit import cli, eko_utils, evolve, utils
 from n3fit.io.writer import XGRID
+from nnpdf_data import theory_cards
+from nnpdf_data.theorydbutils import fetch_theory
 from validphys.loader import FallbackLoader, Loader
 
 _logger = logging.getLogger(__name__)
@@ -174,7 +176,7 @@ def main():
     else:
         # If we are in the business of producing an eko, do some checks before starting:
         # 1. load the nnpdf theory early to check for inconsistent options and theory problems
-        nnpdf_theory = loader.check_theoryID(args.theoryID).get_description()
+        nnpdf_theory = fetch_theory(theory_cards, args.theoryID)
         if nnpdf_theory.get("ModEv") == "TRN" and args.ev_op_iterations is not None:
             raise ValueError("ev_op_iterations is not accepted with ModEv=TRN solution")
 
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/data.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/data.yaml
new file mode 100644
index 0000000000..255639204c
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/data.yaml
@@ -0,0 +1,3 @@
+data_central:
+- 3.50000000e+06
+- 4.53000000e+06
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/filter.py b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/filter.py
new file mode 100644
index 0000000000..512d5a1d30
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/filter.py
@@ -0,0 +1,88 @@
+"""
+filter.py module for ATLAS_WPWM_13TEV dataset
+When running `python filter.py` the relevant uncertainties , data and kinematics yaml
+file will be created in the `nnpdf_data/commondata/ATLAS_WPWM_13TEV` directory.
+"""
+
+import yaml
+from filter_utils import get_kinematics, get_data_values, get_systematics
+from nnpdf_data.filter_utils.utils import prettify_float
+
+yaml.add_representer(float, prettify_float)
+
+
+def filter_ATLAS_WPWM_13TEV_TOT_data_kinetic():
+    """
+    This function writes the central values and kinematics to yaml files.
+    """
+
+    kin = get_kinematics()
+    # only keep first 2 as the last bin is for Z observable
+    central_values = list(get_data_values())[:-1]
+
+    data_central_yaml = {"data_central": central_values}
+
+    kinematics_yaml = {"bins": kin}
+
+    # write central values and kinematics to yaml file
+    with open("data.yaml", "w") as file:
+        yaml.dump(data_central_yaml, file, sort_keys=False)
+
+    with open("kinematics.yaml", "w") as file:
+        yaml.dump(kinematics_yaml, file, sort_keys=False)
+
+
+def filter_ATLAS_WPWM_13TEV_TOT_systematics():
+    """
+    This function writes the systematics to a yaml file.
+    """
+
+    with open("metadata.yaml", "r") as file:
+        metadata = yaml.safe_load(file)
+
+    systematics = get_systematics()
+
+    # error definition
+    error_definitions = {}
+    errors = []
+
+    for sys in systematics:
+        if sys[0]['name'] == 'stat':
+            error_definitions[sys[0]['name']] = {
+                "description": f"{sys[0]['name']}",
+                "treatment": "ADD",
+                "type": "UNCORR",
+            }
+
+        elif sys[0]['name'] == 'ATLAS_LUMI':
+            error_definitions["ATLASLUMI13"] = {
+                "description": f"ATLASLUMI13",
+                "treatment": "MULT",
+                "type": "ATLASLUMI13",
+            }
+
+        else:
+            error_definitions[sys[0]['name']] = {
+                "description": f"{sys[0]['name']}",
+                "treatment": "ADD",
+                "type": f"{sys[0]['name']}",
+            }
+
+    for i in range(metadata['implemented_observables'][0]['ndata']):
+        error_value = {}
+
+        for sys in systematics:
+            error_value[sys[0]['name']] = float(sys[0]['values'][i])
+
+        errors.append(error_value)
+
+    uncertainties_yaml = {"definitions": error_definitions, "bins": errors}
+
+    # write uncertainties
+    with open(f"uncertainties.yaml", 'w') as file:
+        yaml.dump(uncertainties_yaml, file, sort_keys=False)
+
+
+if __name__ == "__main__":
+    filter_ATLAS_WPWM_13TEV_TOT_data_kinetic()
+    filter_ATLAS_WPWM_13TEV_TOT_systematics()
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/filter_utils.py b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/filter_utils.py
new file mode 100644
index 0000000000..68a97d7b6c
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/filter_utils.py
@@ -0,0 +1,135 @@
+"""
+This module contains helper functions that are used to extract the uncertainties, kinematics and data values 
+from the rawdata files.
+"""
+
+import yaml
+import numpy as np
+from nnpdf_data.filter_utils.utils import decompose_covmat
+
+MW2 = 80.385**2
+UNIT_CONVERSION = 1000000
+TABLES = [9, 8, 11]  # order is W-, W+, Z
+
+
+def get_kinematics():
+    """
+    returns the kinematics in the form of a list of dictionaries.
+    """
+    kin = []
+
+    for _ in range(2):
+
+        kin_value = {
+            'm_W2': {'min': None, 'mid': MW2, 'max': None},
+            'sqrts': {'min': None, 'mid': 13000.0, 'max': None},
+        }
+
+        kin.append(kin_value)
+
+    return kin
+
+
+def get_data_values():
+    """
+    returns the central data values in the form of a list.
+    """
+    name_data = lambda tab: f"rawdata/HEPData-ins1436497-v1-Table_{tab}.yaml"
+
+    data_central = []
+
+    for tab in TABLES:
+        with open(name_data(tab), 'r') as file:
+            input = yaml.safe_load(file)
+        values = input['dependent_variables'][0]['values']
+        data_central.append(values[0]['value'] * UNIT_CONVERSION)
+
+    return data_central
+
+
+def get_uncertainties():
+    """
+    Returns array of shape (3,3)
+    Each row corresponds to a different observable: (W-, W+, Z)
+    Each column corresponds to a different systematic: (stat, sys, lumi)
+
+    See table 3 of paper: https://arxiv.org/abs/1603.09222
+    """
+
+    name_data = lambda tab: f"rawdata/HEPData-ins1436497-v1-Table_{tab}.yaml"
+
+    uncertainties = []
+
+    for tab in TABLES:
+        with open(name_data(tab), 'r') as file:
+            input = yaml.safe_load(file)
+        errors = input['dependent_variables'][0]['values'][0]['errors']
+        uncertainties.append(
+            np.array([errors[0]['symerror'], errors[1]['symerror'], errors[2]['symerror']])
+        )
+
+    return np.array(uncertainties) * UNIT_CONVERSION
+
+
+def get_correlation_matrix():
+    """
+    See extra material page: https://atlas.web.cern.ch/Atlas/GROUPS/PHYSICS/PAPERS/STDM-2015-03/tabaux_03.pdf
+
+    Note that this does not include the normalisation uncertainty due to the luminosity.
+    """
+
+    correlation_matrix = np.ones((3, 3))
+    correlation_matrix[0, 1] = 0.93
+    correlation_matrix[1, 0] = correlation_matrix[0, 1]
+    correlation_matrix[0, 2] = 0.18
+    correlation_matrix[2, 0] = correlation_matrix[0, 2]
+    correlation_matrix[1, 2] = 0.19
+    correlation_matrix[2, 1] = correlation_matrix[1, 2]
+
+    return correlation_matrix
+
+
+def get_covariance_matrices():
+    """
+    For the systematics see Table 3 of paper: https://arxiv.org/abs/1603.09222
+
+    Returns:
+    --------
+    tuple: (cov_matrix_no_lumi, lumi_cov)
+    cov_matrix_no_lumi: np.array, the sum of stat and syst covmats -> to be decomposed into artificial systematics
+    lumi_cov: np.array, the lumi covmat. This is correlated between experiments so needs to be saved with type: SPECIAL
+    """
+    corr_matrix = get_correlation_matrix()
+    uncertainties = get_uncertainties()
+
+    # build correlated systematics covariance
+    sys = np.array([uncertainties[i, 1] for i in range(3)])
+    cov_sys = corr_matrix * np.outer(sys, sys)
+
+    # array of lumi uncertainties
+    lumi_unc = np.array([uncertainties[i, 2] for i in range(3)])
+
+    # array of stat uncertainties
+    stat = np.array([uncertainties[i, 0] for i in range(3)])
+
+    return stat, cov_sys, lumi_unc
+
+
+def get_systematics():
+    stat, cov_sys, lumi_unc = get_covariance_matrices()
+
+    # decompose sys covmat
+    syst_unc = decompose_covmat(cov_sys)
+
+    uncertainties = []
+
+    # store only systematics for W+ and W-
+    for i in range(3):
+        uncertainties.append(
+            [{"name": f"ATLAS_WZ_TOT_13TEV_{i}", "values": [syst_unc[0, i], syst_unc[1, i]]}]
+        )
+
+    uncertainties.append([{"name": "stat", "values": [stat[0], stat[1]]}])
+    uncertainties.append([{"name": "ATLAS_LUMI", "values": [lumi_unc[0], lumi_unc[1]]}])
+
+    return uncertainties
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/kinematics.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/kinematics.yaml
new file mode 100644
index 0000000000..cc57c2bf4c
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/kinematics.yaml
@@ -0,0 +1,17 @@
+bins:
+- m_W2:
+    min: null
+    mid: 6.46174823e+03
+    max: null
+  sqrts:
+    min: null
+    mid: 13000.0
+    max: null
+- m_W2:
+    min: null
+    mid: 6.46174823e+03
+    max: null
+  sqrts:
+    min: null
+    mid: 13000.0
+    max: null
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/metadata.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/metadata.yaml
index 77f74cc9a5..f2b77796a0 100644
--- a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/metadata.yaml
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/metadata.yaml
@@ -1,6 +1,6 @@
 setname: ATLAS_WPWM_13TEV
 version: 1
-version_comment: Port of old commondata
+version_comment: New implementation, as for the old one, data is taken from Table 3 of the paper
 nnpdf_metadata:
   nnpdf31_process: DY CC
   experiment: ATLAS
@@ -8,10 +8,11 @@ arXiv:
   url: https://arxiv.org/abs/1603.09222
   journal: Phys. Lett. B759 (2016) 601
 iNSPIRE:
-  url: ''
+  url: 'https://inspirehep.net/literature/1436497'
 hepdata:
-  url: ''
-  version: -1
+  url: 'https://www.hepdata.net/record/73640, https://www.hepdata.net/record/73641'
+
+  version: 1
 implemented_observables:
 - observable_name: TOT
   observable:
@@ -19,11 +20,11 @@ implemented_observables:
     label: ATLAS $W$ inclusive 13 TeV
     units: ''
   process_type: INC
-  tables: []
-  npoints: []
+  tables: [9, 8, 11]
+  npoints: [2]
   ndata: 2
   plotting:
-    kinematics_override: inc_sqrt_scale
+    kinematics_override: identity
     dataset_label: ATLAS $W$ inclusive 13 TeV
     y_label: $\sigma^{fid}$ (fb)
     figure_by:
@@ -37,34 +38,30 @@ implemented_observables:
       - W
     plot_x: ' '
   kinematic_coverage:
-  - k1
-  - k2
-  - k3
+  - m_W2
+  - sqrts
   kinematics:
     variables:
-      k1:
-        description: Variable k1
-        label: k1
-        units: ''
-      k2:
-        description: Variable k2
-        label: k2
-        units: ''
-      k3:
-        description: Variable k3
-        label: k3
-        units: ''
-    file: kinematics_TOT.yaml
+      m_W2:
+        description: W boson squared mass
+        label: '$m_W^2$'
+        units: 'GeV$^2$'
+      sqrts:
+        description: center of mass energy
+        label: '$\sqrt{s}$'
+        units: 'GeV'
+    file: kinematics.yaml
   theory:
     conversion_factor: 1.0
     operation: 'null'
     FK_tables:
     - - ATLAS_WZ_TOT_13TEV-ATLASWZTOT13TEV81PB_WM_tot
       - ATLAS_WZ_TOT_13TEV-ATLASWZTOT13TEV81PB_WP_tot
-  data_uncertainties: []
+  data_uncertainties:
+    - uncertainties.yaml
   variants:
     legacy:
       data_uncertainties:
       - uncertainties_legacy_TOT.yaml
-  data_central: data_legacy_TOT.yaml
+  data_central: data.yaml
   ported_from: ATLAS_W_TOT_13TEV
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/rawdata/HEPData-ins1436497-v1-Table_11.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/rawdata/HEPData-ins1436497-v1-Table_11.yaml
new file mode 100644
index 0000000000..6c400ca342
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/rawdata/HEPData-ins1436497-v1-Table_11.yaml
@@ -0,0 +1,20 @@
+dependent_variables:
+- header: {name: SIG, units: NB}
+  qualifiers:
+  - {name: ABS(ETARAP(C=ELECTRON)), value: < 2.5}
+  - {name: ABS(ETARAP(C=MUON)), value: < 2.5}
+  - {name: MZ, units: GEV, value: 66-116}
+  - {name: PT(C=ELECTRON), units: GEV, value: '> 25'}
+  - {name: PT(C=MUON), units: GEV, value: '> 25'}
+  - {name: RE, value: P P --> ( Z0 < E+ E- + MU+ MU- > + GAMMA* < E+ E- + MU+ MU-
+      > ) X}
+  values:
+  - errors:
+    - {label: stat, symerror: 0.003}
+    - {label: sys, symerror: 0.006}
+    - {label: sys, symerror: 0.016}
+    value: 0.779
+independent_variables:
+- header: {name: SQRT(S), units: GEV}
+  values:
+  - {value: 13000.0}
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/rawdata/HEPData-ins1436497-v1-Table_8.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/rawdata/HEPData-ins1436497-v1-Table_8.yaml
new file mode 100644
index 0000000000..68cd381750
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/rawdata/HEPData-ins1436497-v1-Table_8.yaml
@@ -0,0 +1,20 @@
+dependent_variables:
+- header: {name: SIG, units: NB}
+  qualifiers:
+  - {name: ABS(ETARAP(C=ELECTRON)), value: < 2.5}
+  - {name: ABS(ETARAP(C=MUON)), value: < 2.5}
+  - {name: MT, units: GEV, value: '> 50'}
+  - {name: PT(C=ELECTRON), units: GEV, value: '> 25'}
+  - {name: PT(C=MUON), units: GEV, value: '> 25'}
+  - {name: PT(C=NU), units: GEV, value: '> 25'}
+  - {name: RE, value: P P --> W+ < E+ NUE + MU+ NUMU > X}
+  values:
+  - errors:
+    - {label: stat, symerror: 0.01}
+    - {label: sys, symerror: 0.09}
+    - {label: sys, symerror: 0.1}
+    value: 4.53
+independent_variables:
+- header: {name: SQRT(S), units: GEV}
+  values:
+  - {value: 13000.0}
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/rawdata/HEPData-ins1436497-v1-Table_9.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/rawdata/HEPData-ins1436497-v1-Table_9.yaml
new file mode 100644
index 0000000000..efd8411b2a
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/rawdata/HEPData-ins1436497-v1-Table_9.yaml
@@ -0,0 +1,20 @@
+dependent_variables:
+- header: {name: SIG, units: NB}
+  qualifiers:
+  - {name: 'ABS(ETARAP(C=ELECTRON)}', value: < 2.5}
+  - {name: ABS(ETARAP(C=MUON)), value: < 2.5}
+  - {name: MT, units: GEV, value: '> 50'}
+  - {name: PT(C=ELECTRON), units: GEV, value: '> 25'}
+  - {name: PT(C=MUON), units: GEV, value: '> 25'}
+  - {name: PT(C=NU), units: GEV, value: '> 25'}
+  - {name: RE, value: P P --> W- < E- NUEBAR + MU- NUMUBAR > X}
+  values:
+  - errors:
+    - {label: stat, symerror: 0.01}
+    - {label: sys, symerror: 0.07}
+    - {label: sys, symerror: 0.07}
+    value: 3.5
+independent_variables:
+- header: {name: SQRT(S), units: GEV}
+  values:
+  - {value: 13000.0}
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/uncertainties.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/uncertainties.yaml
new file mode 100644
index 0000000000..e57bb9f4bb
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_WPWM_13TEV/uncertainties.yaml
@@ -0,0 +1,32 @@
+definitions:
+  ATLAS_WZ_TOT_13TEV_0:
+    description: ATLAS_WZ_TOT_13TEV_0
+    treatment: ADD
+    type: ATLAS_WZ_TOT_13TEV_0
+  ATLAS_WZ_TOT_13TEV_1:
+    description: ATLAS_WZ_TOT_13TEV_1
+    treatment: ADD
+    type: ATLAS_WZ_TOT_13TEV_1
+  ATLAS_WZ_TOT_13TEV_2:
+    description: ATLAS_WZ_TOT_13TEV_2
+    treatment: ADD
+    type: ATLAS_WZ_TOT_13TEV_2
+  stat:
+    description: stat
+    treatment: ADD
+    type: UNCORR
+  ATLASLUMI13:
+    description: ATLASLUMI13
+    treatment: MULT
+    type: ATLASLUMI13
+bins:
+- ATLAS_WZ_TOT_13TEV_0: 6.80483165e+04
+  ATLAS_WZ_TOT_13TEV_1: 1.64142169e+04
+  ATLAS_WZ_TOT_13TEV_2: -1.03156424e+01
+  stat: 10000.0
+  ATLAS_LUMI: 70000.0
+- ATLAS_WZ_TOT_13TEV_0: 8.91232992e+04
+  ATLAS_WZ_TOT_13TEV_1: -1.25312809e+04
+  ATLAS_WZ_TOT_13TEV_2: -6.74297267e+01
+  stat: 10000.0
+  ATLAS_LUMI: 100000.0
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/data.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/data.yaml
new file mode 100644
index 0000000000..06a6e43bd9
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/data.yaml
@@ -0,0 +1,2 @@
+data_central:
+- 779000.0
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/filter.py b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/filter.py
new file mode 100644
index 0000000000..1ccca1ede4
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/filter.py
@@ -0,0 +1,96 @@
+"""
+filter.py module for ATLAS_WPWM_13TEV dataset
+When running `python filter.py` the relevant uncertainties , data and kinematics yaml
+file will be created in the `nnpdf_data/commondata/ATLAS_WPWM_13TEV` directory.
+"""
+
+import yaml
+from filter_utils import get_data_values, get_systematics
+from nnpdf_data.filter_utils.utils import prettify_float
+
+yaml.add_representer(float, prettify_float)
+
+MZ2 = 91.1876**2
+
+
+def filter_ATLAS_Z0_13TEV_TOT_data_kinetic():
+    """
+    This function writes the central values and kinematics to yaml files.
+    """
+
+    kin = [
+        {
+            'm_Z2': {'min': None, 'mid': MZ2, 'max': None},
+            'sqrts': {'min': None, 'mid': 13000.0, 'max': None},
+        }
+    ]
+
+    # only keep the last entry corresponding to Z observable
+    central_values = [get_data_values()[-1]]
+
+    data_central_yaml = {"data_central": central_values}
+
+    kinematics_yaml = {"bins": kin}
+
+    # write central values and kinematics to yaml file
+    with open("data.yaml", "w") as file:
+        yaml.dump(data_central_yaml, file, sort_keys=False)
+
+    with open("kinematics.yaml", "w") as file:
+        yaml.dump(kinematics_yaml, file, sort_keys=False)
+
+
+def filter_ATLAS_Z0_13TEV_TOT_systematics():
+    """
+    This function writes the systematics to a yaml file.
+    """
+
+    with open("metadata.yaml", "r") as file:
+        metadata = yaml.safe_load(file)
+
+    systematics = get_systematics()
+
+    # error definition
+    error_definitions = {}
+    errors = []
+
+    for sys in systematics:
+        if sys[0]['name'] == 'stat':
+            error_definitions[sys[0]['name']] = {
+                "description": f"{sys[0]['name']}",
+                "treatment": "ADD",
+                "type": "UNCORR",
+            }
+
+        elif sys[0]['name'] == 'ATLAS_LUMI':
+            error_definitions["ATLASLUMI13"] = {
+                "description": f"ATLASLUMI13",
+                "treatment": "MULT",
+                "type": "ATLASLUMI13",
+            }
+
+        else:
+            error_definitions[sys[0]['name']] = {
+                "description": f"{sys[0]['name']}",
+                "treatment": "ADD",
+                "type": f"{sys[0]['name']}",
+            }
+
+    for i in range(metadata['implemented_observables'][0]['ndata']):
+        error_value = {}
+
+        for sys in systematics:
+            error_value[sys[0]['name']] = float(sys[0]['values'][i])
+
+        errors.append(error_value)
+
+    uncertainties_yaml = {"definitions": error_definitions, "bins": errors}
+
+    # write uncertainties
+    with open(f"uncertainties.yaml", 'w') as file:
+        yaml.dump(uncertainties_yaml, file, sort_keys=False)
+
+
+if __name__ == "__main__":
+    filter_ATLAS_Z0_13TEV_TOT_data_kinetic()
+    filter_ATLAS_Z0_13TEV_TOT_systematics()
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/filter_utils.py b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/filter_utils.py
new file mode 100644
index 0000000000..63f31870ba
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/filter_utils.py
@@ -0,0 +1,115 @@
+"""
+This module contains helper functions that are used to extract the uncertainties, kinematics and data values 
+from the rawdata files.
+"""
+
+import yaml
+import numpy as np
+from nnpdf_data.filter_utils.utils import decompose_covmat
+
+
+UNIT_CONVERSION = 1000000
+TABLES = [9, 8, 11]  # order is W-, W+, Z
+
+
+def get_data_values():
+    """
+    returns the central data values in the form of a list.
+    """
+    name_data = lambda tab: f"rawdata/HEPData-ins1436497-v1-Table_{tab}.yaml"
+
+    data_central = []
+
+    for tab in TABLES:
+        with open(name_data(tab), 'r') as file:
+            input = yaml.safe_load(file)
+        values = input['dependent_variables'][0]['values']
+        data_central.append(values[0]['value'] * UNIT_CONVERSION)
+
+    return data_central
+
+
+def get_uncertainties():
+    """
+    Returns array of shape (3,3)
+    Each row corresponds to a different observable: (W-, W+, Z)
+    Each column corresponds to a different systematic: (stat, sys, lumi)
+
+    See table 3 of paper: https://arxiv.org/abs/1603.09222
+    """
+
+    name_data = lambda tab: f"rawdata/HEPData-ins1436497-v1-Table_{tab}.yaml"
+
+    uncertainties = []
+
+    for tab in TABLES:
+        with open(name_data(tab), 'r') as file:
+            input = yaml.safe_load(file)
+        errors = input['dependent_variables'][0]['values'][0]['errors']
+        uncertainties.append(
+            np.array([errors[0]['symerror'], errors[1]['symerror'], errors[2]['symerror']])
+        )
+
+    return np.array(uncertainties) * UNIT_CONVERSION
+
+
+def get_correlation_matrix():
+    """
+    See extra material page: https://atlas.web.cern.ch/Atlas/GROUPS/PHYSICS/PAPERS/STDM-2015-03/tabaux_03.pdf
+
+    Note that this does not include the normalisation uncertainty due to the luminosity.
+    """
+
+    correlation_matrix = np.ones((3, 3))
+    correlation_matrix[0, 1] = 0.93
+    correlation_matrix[1, 0] = correlation_matrix[0, 1]
+    correlation_matrix[0, 2] = 0.18
+    correlation_matrix[2, 0] = correlation_matrix[0, 2]
+    correlation_matrix[1, 2] = 0.19
+    correlation_matrix[2, 1] = correlation_matrix[1, 2]
+
+    return correlation_matrix
+
+
+def get_covariance_matrices():
+    """
+    For the systematics see Table 3 of paper: https://arxiv.org/abs/1603.09222
+
+    Returns:
+    --------
+    tuple: (cov_matrix_no_lumi, lumi_cov)
+    cov_matrix_no_lumi: np.array, the sum of stat and syst covmats -> to be decomposed into artificial systematics
+    lumi_cov: np.array, the lumi covmat. This is correlated between experiments so needs to be saved with type: SPECIAL
+    """
+    corr_matrix = get_correlation_matrix()
+    uncertainties = get_uncertainties()
+
+    # build correlated systematics covariance
+    sys = np.array([uncertainties[i, 1] for i in range(3)])
+    cov_sys = corr_matrix * np.outer(sys, sys)
+
+    # array of lumi uncertainties
+    lumi_unc = np.array([uncertainties[i, 2] for i in range(3)])
+
+    # array of stat uncertainties
+    stat = np.array([uncertainties[i, 0] for i in range(3)])
+
+    return stat, cov_sys, lumi_unc
+
+
+def get_systematics():
+    stat, cov_sys, lumi_unc = get_covariance_matrices()
+
+    # decompose sys covmat
+    syst_unc = decompose_covmat(cov_sys)
+
+    uncertainties = []
+
+    # store only systematics for W+ and W-
+    for i in range(3):
+        uncertainties.append([{"name": f"ATLAS_WZ_TOT_13TEV_{i}", "values": [syst_unc[2, i]]}])
+
+    uncertainties.append([{"name": "stat", "values": [stat[2]]}])
+    uncertainties.append([{"name": "ATLAS_LUMI", "values": [lumi_unc[2]]}])
+
+    return uncertainties
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/kinematics.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/kinematics.yaml
new file mode 100644
index 0000000000..a881171bd6
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/kinematics.yaml
@@ -0,0 +1,9 @@
+bins:
+- m_Z2:
+    min: null
+    mid: 8.31517839e+03
+    max: null
+  sqrts:
+    min: null
+    mid: 13000.0
+    max: null
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/metadata.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/metadata.yaml
index 4ca09d36f7..e5f3d4b91a 100644
--- a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/metadata.yaml
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/metadata.yaml
@@ -1,6 +1,6 @@
 setname: ATLAS_Z0_13TEV
 version: 1
-version_comment: Port of old commondata
+version_comment: New implementation, as for the old one, data is taken from Table 3 of the paper
 nnpdf_metadata:
   nnpdf31_process: DY NC
   experiment: ATLAS
@@ -19,11 +19,11 @@ implemented_observables:
     label: ATLAS $Z$ inclusive 13 TeV
     units: ''
   process_type: INC
-  tables: []
-  npoints: []
+  tables: [9, 8, 11]
+  npoints: [1]
   ndata: 1
   plotting:
-    kinematics_override: inc_sqrt_scale
+    kinematics_override: identity
     dataset_label: ATLAS $Z$ inclusive 13 TeV
     y_label: $\sigma^{fid}$ (fb)
     figure_by:
@@ -35,33 +35,29 @@ implemented_observables:
       - Z
     plot_x: ' '
   kinematic_coverage:
-  - k1
-  - k2
-  - k3
+  - m_Z2
+  - sqrts
   kinematics:
     variables:
-      k1:
-        description: Variable k1
-        label: k1
-        units: ''
-      k2:
-        description: Variable k2
-        label: k2
-        units: ''
-      k3:
-        description: Variable k3
-        label: k3
-        units: ''
-    file: kinematics_TOT.yaml
+      m_Z2:
+        description: Z boson mass squared
+        label: '$m_{Z}^2$'
+        units: 'GeV$^2$'
+      sqrts:
+        description: center of mass energy
+        label: '$\sqrt{s}$'
+        units: 'GeV'
+    file: kinematics.yaml
   theory:
     conversion_factor: 1.0
     operation: 'null'
     FK_tables:
     - - ATLAS_WZ_TOT_13TEV-ATLASWZTOT13TEV81PB_Z_tot
-  data_uncertainties: []
+  data_uncertainties: 
+    - uncertainties.yaml
   variants:
     legacy:
       data_uncertainties:
       - uncertainties_legacy_TOT.yaml
-  data_central: data_legacy_TOT.yaml
+  data_central: data.yaml
   ported_from: ATLAS_Z_TOT_13TEV
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/rawdata/HEPData-ins1436497-v1-Table_11.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/rawdata/HEPData-ins1436497-v1-Table_11.yaml
new file mode 100644
index 0000000000..6c400ca342
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/rawdata/HEPData-ins1436497-v1-Table_11.yaml
@@ -0,0 +1,20 @@
+dependent_variables:
+- header: {name: SIG, units: NB}
+  qualifiers:
+  - {name: ABS(ETARAP(C=ELECTRON)), value: < 2.5}
+  - {name: ABS(ETARAP(C=MUON)), value: < 2.5}
+  - {name: MZ, units: GEV, value: 66-116}
+  - {name: PT(C=ELECTRON), units: GEV, value: '> 25'}
+  - {name: PT(C=MUON), units: GEV, value: '> 25'}
+  - {name: RE, value: P P --> ( Z0 < E+ E- + MU+ MU- > + GAMMA* < E+ E- + MU+ MU-
+      > ) X}
+  values:
+  - errors:
+    - {label: stat, symerror: 0.003}
+    - {label: sys, symerror: 0.006}
+    - {label: sys, symerror: 0.016}
+    value: 0.779
+independent_variables:
+- header: {name: SQRT(S), units: GEV}
+  values:
+  - {value: 13000.0}
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/rawdata/HEPData-ins1436497-v1-Table_8.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/rawdata/HEPData-ins1436497-v1-Table_8.yaml
new file mode 100644
index 0000000000..68cd381750
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/rawdata/HEPData-ins1436497-v1-Table_8.yaml
@@ -0,0 +1,20 @@
+dependent_variables:
+- header: {name: SIG, units: NB}
+  qualifiers:
+  - {name: ABS(ETARAP(C=ELECTRON)), value: < 2.5}
+  - {name: ABS(ETARAP(C=MUON)), value: < 2.5}
+  - {name: MT, units: GEV, value: '> 50'}
+  - {name: PT(C=ELECTRON), units: GEV, value: '> 25'}
+  - {name: PT(C=MUON), units: GEV, value: '> 25'}
+  - {name: PT(C=NU), units: GEV, value: '> 25'}
+  - {name: RE, value: P P --> W+ < E+ NUE + MU+ NUMU > X}
+  values:
+  - errors:
+    - {label: stat, symerror: 0.01}
+    - {label: sys, symerror: 0.09}
+    - {label: sys, symerror: 0.1}
+    value: 4.53
+independent_variables:
+- header: {name: SQRT(S), units: GEV}
+  values:
+  - {value: 13000.0}
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/rawdata/HEPData-ins1436497-v1-Table_9.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/rawdata/HEPData-ins1436497-v1-Table_9.yaml
new file mode 100644
index 0000000000..efd8411b2a
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/rawdata/HEPData-ins1436497-v1-Table_9.yaml
@@ -0,0 +1,20 @@
+dependent_variables:
+- header: {name: SIG, units: NB}
+  qualifiers:
+  - {name: 'ABS(ETARAP(C=ELECTRON)}', value: < 2.5}
+  - {name: ABS(ETARAP(C=MUON)), value: < 2.5}
+  - {name: MT, units: GEV, value: '> 50'}
+  - {name: PT(C=ELECTRON), units: GEV, value: '> 25'}
+  - {name: PT(C=MUON), units: GEV, value: '> 25'}
+  - {name: PT(C=NU), units: GEV, value: '> 25'}
+  - {name: RE, value: P P --> W- < E- NUEBAR + MU- NUMUBAR > X}
+  values:
+  - errors:
+    - {label: stat, symerror: 0.01}
+    - {label: sys, symerror: 0.07}
+    - {label: sys, symerror: 0.07}
+    value: 3.5
+independent_variables:
+- header: {name: SQRT(S), units: GEV}
+  values:
+  - {value: 13000.0}
diff --git a/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/uncertainties.yaml b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/uncertainties.yaml
new file mode 100644
index 0000000000..485c31e9e7
--- /dev/null
+++ b/nnpdf_data/nnpdf_data/commondata/ATLAS_Z0_13TEV/uncertainties.yaml
@@ -0,0 +1,27 @@
+definitions:
+  ATLAS_WZ_TOT_13TEV_0:
+    description: ATLAS_WZ_TOT_13TEV_0
+    treatment: ADD
+    type: ATLAS_WZ_TOT_13TEV_0
+  ATLAS_WZ_TOT_13TEV_1:
+    description: ATLAS_WZ_TOT_13TEV_1
+    treatment: ADD
+    type: ATLAS_WZ_TOT_13TEV_1
+  ATLAS_WZ_TOT_13TEV_2:
+    description: ATLAS_WZ_TOT_13TEV_2
+    treatment: ADD
+    type: ATLAS_WZ_TOT_13TEV_2
+  stat:
+    description: stat
+    treatment: ADD
+    type: UNCORR
+  ATLASLUMI13:
+    description: ATLASLUMI13
+    treatment: MULT
+    type: ATLASLUMI13
+bins:
+- ATLAS_WZ_TOT_13TEV_0: 1.13953999e+03
+  ATLAS_WZ_TOT_13TEV_1: -1.14718978e+02
+  ATLAS_WZ_TOT_13TEV_2: 5.88967641e+03
+  stat: 3000.0
+  ATLAS_LUMI: 16000.0