Merge pull request #1145 from cta-observatory/write_images_dl1ab

Write modified images from dl1ab in output DL1 files and fix `add_config_metadata`

lstchain/image/modifier.py

@@ -299,9 +299,6 @@ def calculate_noise_parameters(simtel_filename, data_dl1_filename,
     pedestal_calibrator.image_extractors[ped_config['charge_product']].apply_integration_correction = True
     shower_calibrator.image_extractors[shower_extractor_type].apply_integration_correction = True
 
-
-
-
     # Pulse integration window width of the (biased) extractor for showers:
     shower_extractor_window_width = config[config['image_extractor']]['window_width']
 

lstchain/io/config.py

@@ -126,3 +126,24 @@ def dump_config(config, filename, overwrite=False):
         raise FileExistsError(f"File {filename} exists, use overwrite=True")
     with open(filename, 'w') as file:
         json.dump(config, file, indent=2)
+
+
+def includes_image_modification(config):
+    """
+    Check if the image modifier has been used in the given configuration.
+
+    Parameters
+    ----------
+    config : `dict`
+        The configuration dictionary to check.
+
+    Returns
+    -------
+    `bool`
+        `True` if the image modifier has been used, `False` otherwise.
+    """
+    imconfig = config.get('image_modifier', {})
+    increase_nsb = imconfig.get("increase_nsb", False)
+    increase_psf = imconfig.get("increase_psf", False)
+
+    return increase_nsb or increase_psf

lstchain/io/io.py

@@ -10,6 +10,9 @@
 import tables
 from tables import open_file
 from tqdm import tqdm
+import json
+from traitlets.config.loader import DeferredConfigString, LazyConfigValue
+from pathlib import PosixPath
 
 import astropy.units as u
 from astropy.table import Table, vstack, QTable
@@ -31,6 +34,8 @@
     ThrownEventsHistogram,
 )
 
+
+
 log = logging.getLogger(__name__)
 
 __all__ = [
@@ -657,6 +662,50 @@ def add_global_metadata(container, metadata):
         container.meta[k] = item
 
 
+
+
+def serialize_config(obj):
+    """
+    Serialize an object to a JSON-serializable format.
+
+    Parameters
+    ----------
+    obj : object
+        The object to serialize.
+
+    Returns
+    -------
+    object
+        The serialized object.
+
+    Raises
+    ------
+    TypeError
+        If the object is not serializable.
+
+    Notes
+    -----
+    This function serializes an object to a JSON-serializable format. It supports the following types:
+    - LazyConfigValue
+    - DeferredConfigString
+    - PosixPath
+    - numpy.ndarray
+
+    If the object is not one of the above types, a TypeError is raised.
+
+    """
+    if isinstance(obj, LazyConfigValue):
+        return obj.to_dict()
+    elif isinstance(obj, DeferredConfigString):
+        return str(obj)
+    elif isinstance(obj, PosixPath):
+        return obj.as_posix()
+    elif isinstance(obj, np.ndarray):
+        return obj.tolist()
+    else:
+        raise TypeError(f"Type {type(obj).__name__} not serializable")
+
+
 def add_config_metadata(container, configuration):
     """
     Add configuration parameters to a container in container.meta.config
@@ -666,18 +715,7 @@ def add_config_metadata(container, configuration):
     container: `ctapipe.containers.Container`
     configuration: config dict
     """
-    linted_config = str(configuration)
-    linted_config = linted_config.replace("<LazyConfigValue {}>", "None")
-    linted_config = re.sub(r"<LazyConfigValue\svalue=(.*?)>", "\\1", linted_config)
-    linted_config = re.sub(r"DeferredConfigString\((.*?)\)", "\\1", linted_config)
-    linted_config = re.sub(r"PosixPath\((.*?)\)", "\\1", linted_config)
-    linted_config = linted_config.replace("\'", "\"")
-    linted_config = linted_config.replace("None", "\"None\"")
-    linted_config = linted_config.replace("inf", "\"inf\"")
-    linted_config = linted_config.replace("True", "true")
-    linted_config = linted_config.replace("False", "false")
-
-    container.meta["config"] = linted_config
+    container.meta["config"] = json.dumps(configuration, default=serialize_config)
 
 
 def write_subarray_tables(writer, event, metadata=None):
@@ -1261,3 +1299,4 @@ def get_mc_fov_offset(filename):
     mean_offset = min_viewcone + 0.5 * (max_viewcone - min_viewcone)
 
     return mean_offset
+

lstchain/io/tests/test_config.py

@@ -1,6 +1,6 @@
-from lstchain.io import config
-from lstchain.io.config import get_cleaning_parameters
 import tempfile
+from lstchain.io import config
+from lstchain.io.config import get_cleaning_parameters, includes_image_modification
 
 
 def test_get_standard_config():
@@ -52,3 +52,21 @@ def test_dump_config():
         config.dump_config(cfg, file.name, overwrite=True)
         read_cfg = config.read_configuration_file(file.name)
         assert read_cfg['myconf'] == 1
+
+
+def test_includes_image_modification_no_modif():
+    cfg = {}
+    assert not includes_image_modification(cfg)
+    cfg = {"image_modifier": {}}
+    assert not includes_image_modification(cfg)
+
+
+def test_includes_image_modification_with_modif():
+    cfg = {"image_modifier": {"increase_psf": True, "increase_nsb": False}}
+    assert includes_image_modification(cfg)
+    cfg = {"image_modifier": {"increase_nsb": True, "increase_psf": False}}
+    assert includes_image_modification(cfg)
+    cfg = {"image_modifier": {"increase_psf": True, "increase_nsb": True}}
+    assert includes_image_modification(cfg)
+    cfg = {"image_modifier": {"increase_psf": False, "increase_nsb": False}}
+    assert not includes_image_modification(cfg)

lstchain/io/tests/test_io.py

@@ -1,11 +1,15 @@
 import tempfile
-
+import json
+import math
 import numpy as np
 import pandas as pd
 import pytest
 import tables
 from astropy.table import Table, QTable
 from ctapipe.instrument import SubarrayDescription
+from lstchain.io import add_config_metadata
+from pathlib import PosixPath
+from traitlets.config.loader import DeferredConfigString, LazyConfigValue
 
 
 @pytest.fixture
@@ -149,39 +153,86 @@ def test_trigger_type_in_dl1_params(simulated_dl1_file):
 
 def test_extract_simulation_nsb(mc_gamma_testfile):
     from lstchain.io.io import extract_simulation_nsb
+
     nsb = extract_simulation_nsb(mc_gamma_testfile)
     assert np.isclose(nsb[0], 0.246, rtol=0.1)
     assert np.isclose(nsb[1], 0.217, rtol=0.1)
 
 
 def test_remove_duplicated_events():
     from lstchain.io.io import remove_duplicated_events
-
-    d = {'event_id': [1, 2, 3,
-                      1, 2, 4,
-                      1, 2, 3], 
-         'gh_score': [0.1, 0.5, 0.7,
-                      0.5, 0.8, 0.1,
-                      0.9, 0.1, 0.5],
-         'alpha': range(9) 
-     }
+
+    d = {
+        "event_id": [1, 2, 3, 1, 2, 4, 1, 2, 3],
+        "gh_score": [0.1, 0.5, 0.7, 0.5, 0.8, 0.1, 0.9, 0.1, 0.5],
+        "alpha": range(9),
+    }
     df = pd.DataFrame(data=d)
     data1 = QTable.from_pandas(df)
     remove_duplicated_events(data1)
 
-    d2 = {'event_id': [3, 2, 4, 1], 
-          'gh_score': [0.7, 0.8, 0.1, 0.9],
-          'alpha': [2, 4, 5, 6]
-      }
+    d2 = {
+        "event_id": [3, 2, 4, 1],
+        "gh_score": [0.7, 0.8, 0.1, 0.9],
+        "alpha": [2, 4, 5, 6],
+    }
     df2 = pd.DataFrame(data=d2)
-    data2= QTable.from_pandas(df2)
+    data2 = QTable.from_pandas(df2)
 
-    assert np.all(data1==data2)
+    assert np.all(data1 == data2)
 
 
 def test_check_mc_type(simulated_dl1_file):
     from lstchain.io.io import check_mc_type
 
     mc_type = check_mc_type(simulated_dl1_file)
-    assert mc_type == 'diffuse'
-
+    assert mc_type == "diffuse"
+
+
+def test_add_config_metadata():
+    class Container:
+        meta = {}
+
+    lazy_value = LazyConfigValue()
+    lazy_value.update({"key": "new_value"})
+
+    config = {
+        "param1": 1,
+        "param2": "value2",
+        "param3": [1, 2, 3],
+        "param4": {"a": 1, "b": 2},
+        "param5": None,
+        "param6": lazy_value,
+        "param7": DeferredConfigString("some_string"),
+        "param8": PosixPath("/path/to/file"),
+        "param9": np.inf,
+        "param10": True,
+        "param11": False,
+        "param12": np.array([1, 2, 3]),
+    }
+
+    expected_config = {
+        "param1": 1,
+        "param2": "value2",
+        "param3": [1, 2, 3],
+        "param4": {"a": 1, "b": 2},
+        "param5": None,
+        "param6": {"update": {"key": "new_value"}},
+        "param7": "some_string",
+        "param8": "/path/to/file",
+        "param9": math.inf,
+        "param10": True,
+        "param11": False,
+        "param12": [1, 2, 3],
+    }
+
+    container = Container()
+    add_config_metadata(container, config)
+    assert json.loads(container.meta["config"]) == expected_config
+
+    # test also with standard config in case of future changes 
+    from lstchain.io.config import get_standard_config
+    config = get_standard_config()
+    container = Container()
+    add_config_metadata(container, config)
+    assert json.loads(container.meta["config"]) == config

lstchain/scripts/lstchain_dl1ab.py

@@ -15,6 +15,7 @@
 import argparse
 import logging
 from pathlib import Path
+import json
 
 import astropy.units as u
 import numpy as np
@@ -31,13 +32,14 @@
 from lstchain.calib.camera.pixel_threshold_estimation import get_threshold_from_dl1_file
 from lstchain.image.cleaning import apply_dynamic_cleaning
 from lstchain.image.modifier import random_psf_smearer, set_numba_seed, add_noise_in_pixels
-from lstchain.io import get_dataset_keys, copy_h5_nodes, HDF5_ZSTD_FILTERS, add_source_filenames
+from lstchain.io import get_dataset_keys, copy_h5_nodes, HDF5_ZSTD_FILTERS, add_source_filenames, add_config_metadata
 
 from lstchain.io.config import (
     get_cleaning_parameters,
     get_standard_config,
     read_configuration_file,
     replace_config,
+    includes_image_modification,
 )
 from lstchain.io.io import (
     dl1_images_lstcam_key,
@@ -82,7 +84,7 @@
 
 parser.add_argument(
     '--no-image', action='store_true',
-    help='Pass this argument to avoid writing the images in the new DL1 files. Beware, if `increase_nsb` or `increase_psf` are True in the config, the images will not be written.',
+    help='Pass this argument to avoid writing the images in the new DL1 files.',
 )
 
 parser.add_argument(
@@ -148,9 +150,9 @@ def main():
 
     if increase_nsb or increase_psf:
         log.info(f"image_modifier configuration: {imconfig}")
-        log.info("NOTE: Using the image_modifier options means images will "
-                 "not be saved.")
-        args.no_image = True
+        if not args.no_image:
+            log.info("Modified images are saved in the output file.")
+
     if increase_nsb:
         extra_noise_in_dim_pixels = imconfig["extra_noise_in_dim_pixels"]
         extra_bias_in_dim_pixels = imconfig["extra_bias_in_dim_pixels"]
@@ -243,6 +245,14 @@ def main():
 
     with tables.open_file(args.input_file, mode='r') as infile:
         image_table = read_table(infile, dl1_images_lstcam_key)
+        # if the image modifier has been used to produce these images, stop here
+        config_from_image_table = json.loads(image_table.meta['config'])
+        if includes_image_modification(config_from_image_table) and includes_image_modification(config):
+            log.critical(f"\nThe image modifier has already been used to produce the images in file {args.input_file}.\n"
+                        "Re-applying the image modifier is not a good practice, start again from unmodified images please.")
+            sys.exit(1)
+
+        images = image_table['image']
         params = read_table(infile, dl1_params_lstcam_key)
         dl1_params_input = params.colnames
 
@@ -270,10 +280,6 @@ def main():
         if increase_psf:
             set_numba_seed(infile.root.dl1.event.subarray.trigger.col('obs_id')[0])
 
-        image_mask_save = not args.no_image and 'image_mask' in infile.root[dl1_images_lstcam_key].colnames
-        if image_mask_save:
-            image_mask = image_table['image_mask']
-
         new_params = set(parameters_to_update.keys()) - set(params.colnames)
         for p in new_params:
             params[p] = np.empty(len(params), dtype=parameters_to_update[p])
@@ -282,7 +288,6 @@ def main():
             copy_h5_nodes(infile, outfile, nodes=nodes_keys)
             add_source_filenames(outfile, [args.input_file])
 
-
             # need container to use lstchain.io.add_global_metadata and lstchain.io.add_config_metadata
             for k, item in metadata.as_dict().items():
                 outfile.root[dl1_params_lstcam_key].attrs[k] = item
@@ -420,13 +425,17 @@ def main():
                 for p in parameters_to_update:
                     params[ii][p] = u.Quantity(dl1_container[p]).value
 
-                if image_mask_save:
-                    image_mask[ii] = signal_pixels
+                images[ii] = image
+
+                if 'image_mask' in image_table.colnames:
+                    image_table['image_mask'][ii] = signal_pixels
+
 
-            if image_mask_save or catB_calib:
-                # the image table has been modified and needs to be saved
+            add_config_metadata(image_table, config)
+            if not args.no_image:
                 write_table(image_table, outfile, dl1_images_lstcam_key, overwrite=True, filters=HDF5_ZSTD_FILTERS)
 
+            add_config_metadata(params, config)
             write_table(params, outfile, dl1_params_lstcam_key, overwrite=True, filters=HDF5_ZSTD_FILTERS)
 
             # write a cat-B calibrations in DL1b