diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
index 62e972c4059f4856eaea9acb38db0c3a7bb95185..a9882bc13cd80aeb6e707af51ea78bbe032b4b2b 100644
--- a/.gitlab-ci.yml
+++ b/.gitlab-ci.yml
@@ -11,6 +11,7 @@ variables:
   TEST_HBA_DATASET_NAME: "test_data.tar.gz"
   CALIBRATOR_HBA_RESULTS_NAME: "results_calibrator.tar.gz"
   TARGET_HBA_RESULTS_NAME: "results_target.tar.gz"
+  TARGET_HBA_SELFCAL_RESULTS_NAME: "results_target_selfcal.tar.gz"
   PIP_CACHE_DIR: "$CI_PROJECT_DIR/.cache/pip"
   BUILD_DOCKER_IMAGE: "1"
 
@@ -105,6 +106,7 @@ download_data:
     - wget -nv https://support.astron.nl/software/ci_data/linc/$TEST_HBA_DATASET_NAME -O $TEST_HBA_DATASET_NAME && tar xfz $TEST_HBA_DATASET_NAME && rm -f $TEST_HBA_DATASET_NAME
     - wget -nv https://support.astron.nl/software/ci_data/linc/$CALIBRATOR_HBA_RESULTS_NAME -O $CALIBRATOR_HBA_RESULTS_NAME && tar xfz $CALIBRATOR_HBA_RESULTS_NAME && rm -f $CALIBRATOR_HBA_RESULTS_NAME
     - wget -nv https://support.astron.nl/software/ci_data/linc/$TARGET_HBA_RESULTS_NAME -O $TARGET_HBA_RESULTS_NAME && tar xfz $TARGET_HBA_RESULTS_NAME && rm -f $TARGET_HBA_RESULTS_NAME
+    - wget -nv https://support.astron.nl/software/ci_data/linc/$TARGET_HBA_SELFCAL_RESULTS_NAME -O $TARGET_HBA_SELFCAL_RESULTS_NAME && tar xfz $TARGET_HBA_SELFCAL_RESULTS_NAME && rm -f $TARGET_HBA_SELFCAL_RESULTS_NAME
   artifacts:
     paths:
     - data
@@ -167,6 +169,20 @@ run_hba_target:
       - hba_target_logs.tar.gz
     when: on_failure
 
+run_hba_target_selfcal:
+  stage: tests
+  needs: ["versioning", "download_data"]
+  image: $INTEGRATION_IMAGE
+  script:
+    - cwltool --no-container --preserve-environment PATH --preserve-environment LINC_DATA_ROOT --preserve-environment PYTHONPATH --outdir results --leave-tmpdir --tmpdir-prefix /tmp/run_hba_target/ workflows/HBA_target.cwl test_jobs/HBA_target_selfcal.json
+    - test_jobs/check_workflow_results.py results /builds/RD/LINC/data/results_target_selfcal
+  after_script:
+    - find /tmp/run_hba_target -name "*.log" -print0 | tar czf hba_target_selfcal_logs.tar.gz --null -T -
+  artifacts:
+    paths:
+      - hba_target_selfcal_logs.tar.gz
+    when: on_failure
+
 build_doc:
   stage: docs
   needs: ["versioning", "download_data"]
diff --git a/docs/source/target.rst b/docs/source/target.rst
index 451a0486eef01c69201fa3c703382f2a97ed637d..c7fc87212d5b9cbcb5818ea27a6e4897cce9aa49 100644
--- a/docs/source/target.rst
+++ b/docs/source/target.rst
@@ -283,8 +283,13 @@ User-defined parameter configuration
 - ``use_target``: enable downloading of a target skymodel (default: ``true``)
 - ``skymodel_source``: choose the target skymodel from `TGSS ADR`_ or the new `Global Sky Model`_ (GSM) (default: ``TGSS``)
 - ``skymodel_fluxlimit``: limits the input skymodel to sources that exceed the given flux density limit in Jy (default: None for **HBA**, i.e. all sources of the catalogue will be kept, and 1.0 for **LBA**)
-- ``selfcal``: perform extensive self-calibration according to the `LiLF`_ scheme (recommended for **LBA** observations) (default: ``false``)
-- ``selfcal_region``: ds9-compatible region file to select the image regions used for the self-calibration
+- ``selfcal``: perform self-calibration (default: ``false``)
+- ``selfcal_strategy``: sets the strategy for selfcal. If set to ``HBA``. If set to ``LBA``, selfcal will perform extensive self-calibration according to the `LiLF`_ scheme (recommended for **LBA** observations). (default: ``HBA``)
+- ``selfcal_region``: ds9-compatible region file to select the image regions used for the self-calibration in case of LBA self-calibration.
+- ``selfcal_hba_uvlambdamin``: specify the minimum uv-distance in units of wavelength to be used when performing selfcal with HBA (default: 200)
+- ``selfcal_hba_imsize``: specifies the image size in pixels, as a list, to use during HBA self-calibration (default: ``[20000, 20000]``).
+- ``output_channels_per_chunk``: HBA only. Sets the number of frequency channels to chunk data in after self-calibration (default: 20).
+- ``calib_nchan``: number of channels to be combined when calibration (default: ``0`` (one solution per group) if `selfcal = false`, otherwise `1` (one solution per channel))
 
 A comprehensive explanation of the baseline selection syntax can be found `here`_.
 
@@ -324,11 +329,10 @@ A comprehensive explanation of the baseline selection syntax can be found `here`
 - ``avg_freqresolution`` : intermediate frequency resolution of the data after averaging (default: 48.82kHz, which translates to 4 channels per subband)
 - ``avg_timeresolution_concat``: final time resolution of the data in seconds after averaging and concatenation (default: 8)
 - ``avg_freqresolution_concat``: final frequency resolution of the data after averaging and concatenation (default: 97.64kHz, which translates to 2 channels per subband)
-- ``num_SBs_per_group``: make concatenated measurement-sets with that many subbands, choose a high number if running LBA (default: 10)
 
 *Concatenating of the target data*
 
-- ``num_SBs_per_group``: make concatenated measurement-sets with that many subbands (default: 10)
+- ``num_SBs_per_group``: make concatenated measurement-sets with that many subbands (default: 10 normally, -1 for HBA selfcal)
 - ``reference_stationSB``: station-subband number to use as reference for grouping (default: ``None`` -> use lowest frequency input data as reference)
 - ``chunkduration``: Duration (in seconds) after which to start writing a next measurement set while concatenating (default: 0.0, no chunking in time)
 
diff --git a/scripts/MSChunker.py b/scripts/MSChunker.py
new file mode 100755
index 0000000000000000000000000000000000000000..916ebf117a1ee620346611eec6833f7bd1bf1945
--- /dev/null
+++ b/scripts/MSChunker.py
@@ -0,0 +1,629 @@
+#!/usr/bin/env python
+"""
+Portable module to chunk MeasurementSets fractionally by time or frequency.
+
+Some code adapted from Rapthor: https://git.astron.nl/RD/rapthor
+"""
+import argparse
+import logging
+import os
+import subprocess
+
+import casacore.tables as pt
+import numpy as np
+from astropy.time import Time
+
+logging.basicConfig(
+    format="%(levelname)s:%(asctime)s %(name)s ---- %(message)s",
+    datefmt="%m/%d/%Y %H:%M:%S",
+)
+
+
+def normalize_ra(num: float) -> float:
+    """
+    Normalize RA to be in the range [0, 360).
+
+    Based on https://github.com/phn/angles/blob/master/angles.py
+
+    Parameters
+    ----------
+    num : float
+        The RA in degrees to be normalized.
+
+    Returns
+    -------
+    res : float
+        RA in degrees in the range [0, 360).
+    """
+    lower = 0.0
+    upper = 360.0
+    res = num
+    if num > upper or num == lower:
+        num = lower + abs(num + upper) % (abs(lower) + abs(upper))
+    if num < lower or num == upper:
+        num = upper - abs(num - lower) % (abs(lower) + abs(upper))
+    res = lower if num == upper else num
+
+    return res
+
+
+def normalize_dec(num: float) -> float:
+    """
+    Normalize Dec to be in the range [-90, 90].
+
+    Based on https://github.com/phn/angles/blob/master/angles.py
+
+    Parameters
+    ----------
+    num : float
+        The Dec in degrees to be normalized.
+
+    Returns
+    -------
+    res : float
+        Dec in degrees in the range [-90, 90].
+    """
+    lower = -90.0
+    upper = 90.0
+    res = num
+    total_length = abs(lower) + abs(upper)
+    if num < -total_length:
+        num += np.ceil(num / (-2 * total_length)) * 2 * total_length
+    if num > total_length:
+        num -= np.floor(num / (2 * total_length)) * 2 * total_length
+    if num > upper:
+        num = total_length - num
+    if num < lower:
+        num = -total_length - num
+    res = num
+
+    return res
+
+
+def concat_time_command(msfiles: list, output_file: str) -> list:
+    """
+    Construct command to concatenate files in time using TAQL
+
+    Parameters
+    ----------
+    msfiles : list of str
+        List of MS filenames to be concatenated
+    output_file : str
+        Filename of output concatenated MS
+
+    Returns
+    -------
+    cmd : list of str
+        Command to be executed by subprocess.run()
+    """
+    cmd = [
+        "taql",
+        "select",
+        "from",
+        "[{}]".format(",".join(msfiles)),
+        "giving",
+        '"{}"'.format(output_file),
+        "AS",
+        "PLAIN",
+    ]
+    return cmd
+
+
+def convert_mjd(mjd_sec: float) -> str:
+    """
+    Converts MJD to casacore MVTime
+
+    Parameters
+    ----------
+    mjd_sec : float
+        MJD time in seconds
+
+    Returns
+    -------
+    mvtime : str
+        Casacore MVTime string
+    """
+    t = Time(mjd_sec / 3600 / 24, format="mjd", scale="utc")
+    date, hour = t.iso.split(" ")
+    year, month, day = date.split("-")
+    d = t.datetime
+    month = d.ctime().split(" ")[1]
+
+    return "{0}{1}{2}/{3}".format(day, month, year, hour)
+
+
+class Observation(object):
+    """
+    The Observation object contains various MS-related parameters
+
+    Parameters
+    ----------
+    ms_filename : str
+        Filename of the MS file
+    starttime : float, optional
+        The start time of the observation (in MJD seconds). If None, the start time
+        is the start of the MS file
+    endtime : float, optional
+        The end time of the observation (in MJD seconds). If None, the end time
+        is the end of the MS file
+    startfreq : float, optional
+        The start freq of the observation (in Hz). If None, the start freq
+        is the start of the MS file
+    endfreq : float, optional
+        The end freq of the observation (in Hz). If None, the end freq
+        is the end of the MS file
+    """
+
+    def __init__(
+        self,
+        ms_filename: str,
+        starttime: float = None,
+        endtime: float = None,
+        startfreq: float = None,
+        endfreq: float = None,
+    ):
+        self.ms_filename = ms_filename
+        self.name = os.path.basename(self.ms_filename.rstrip("/"))
+        self.log = logging.getLogger("Observation:{}".format(self.name))
+        self.log.setLevel(logging.INFO)
+        self.starttime = starttime
+        self.endtime = endtime
+        self.startfreq = startfreq
+        self.endfreq = endfreq
+        self.parameters = {}
+        self.scan_ms()
+
+        # Define the infix for filenames
+        if self.startsat_startofms and self.goesto_endofms:
+            # Don't include starttime if observation covers full MS
+            self.infix = ""
+        else:
+            # Include starttime to avoid naming conflicts
+            self.infix = ".mjd{}".format(int(self.starttime))
+
+    def scan_ms(self):
+        """
+        Scans input MS and stores info
+        """
+        # Get time info
+        tab = pt.table(self.ms_filename, ack=False)
+        if self.starttime is None:
+            self.starttime = np.min(tab.getcol("TIME"))
+        else:
+            valid_times = np.where(tab.getcol("TIME") >= self.starttime)[0]
+            if len(valid_times) == 0:
+                raise ValueError(
+                    "Start time of {0} is greater than the last time in the "
+                    "MS".format(self.starttime)
+                )
+            self.starttime = tab.getcol("TIME")[valid_times[0]]
+
+        # DPPP takes ceil(startTimeParset - startTimeMS), so ensure that our start time is
+        # slightly less than the true one (if it's slightly larger, DPPP sets the first
+        # time to the next time, skipping the first time slot)
+        self.starttime -= 0.1
+        if self.starttime > np.min(tab.getcol("TIME")):
+            self.startsat_startofms = False
+        else:
+            self.startsat_startofms = True
+        if self.endtime is None:
+            self.endtime = np.max(tab.getcol("TIME"))
+        else:
+            valid_times = np.where(tab.getcol("TIME") <= self.endtime)[0]
+            if len(valid_times) == 0:
+                raise ValueError(
+                    "End time of {0} is less than the first time in the "
+                    "MS".format(self.endtime)
+                )
+            self.endtime = tab.getcol("TIME")[valid_times[-1]]
+        if self.endtime < np.max(tab.getcol("TIME")):
+            self.goesto_endofms = False
+        else:
+            self.goesto_endofms = True
+        self.timepersample = tab.getcell("EXPOSURE", 0)
+        self.numsamples = int(
+            np.ceil((self.endtime - self.starttime) / self.timepersample)
+        )
+        tab.close()
+
+        # Get frequency info
+        sw = pt.table(self.ms_filename + "::SPECTRAL_WINDOW", ack=False)
+        self.referencefreq = sw.col("REF_FREQUENCY")[0]
+        channels = sw.col("CHAN_FREQ")[0]
+        if self.startfreq is None:
+            self.startfreq = np.min(channels)
+        if (self.endfreq is None) or (self.endfreq > channels.max()):
+            self.endfreq = np.max(channels)
+        self.startchan = int(np.argwhere(channels == self.startfreq))
+        self.endchan = int(np.argwhere(channels == self.endfreq))
+        if self.endfreq >= channels.max():
+            self.endchan += 1
+        self.numchannels = sw.col("NUM_CHAN")[0]
+        self.channelwidth = sw.col("CHAN_WIDTH")[0][0]
+        sw.close()
+
+        # Get pointing info
+        obs = pt.table(self.ms_filename + "::FIELD", ack=False)
+        self.ra = normalize_ra(np.degrees(float(obs.col("REFERENCE_DIR")[0][0][0])))
+        self.dec = normalize_dec(np.degrees(float(obs.col("REFERENCE_DIR")[0][0][1])))
+        obs.close()
+
+        # Get station names and diameter
+        ant = pt.table(self.ms_filename + "::ANTENNA", ack=False)
+        self.stations = ant.col("NAME")[:]
+        self.diam = float(ant.col("DISH_DIAMETER")[0])
+        if "HBA" in self.stations[0]:
+            self.antenna = "HBA"
+        elif "LBA" in self.stations[0]:
+            self.antenna = "LBA"
+        else:
+            self.log.warning(
+                "Antenna type not recognized (only LBA and HBA data "
+                "are supported at this time)"
+            )
+        ant.close()
+
+        # Find mean elevation and FOV
+        el_values = pt.taql(
+            "SELECT mscal.azel1()[1] AS el from " + self.ms_filename + " limit ::10000"
+        ).getcol("el")
+        self.mean_el_rad = np.mean(el_values)
+
+
+class MSChunker:
+    """Handles chunking a MeasurementSet in time."""
+
+    def __init__(self, msin: list, fraction: float = 1.0):
+        """Handles chunking a MeasurementSet
+
+        Parameters
+        ----------
+        msin : list
+            List of MS filenames to be chunked.
+        fraction : float
+            Fraction of time to take.
+        """
+        self.time_fraction = fraction
+        if type(msin) is str:
+            mslist = [msin]
+        else:
+            mslist = msin
+        self.full_observations = []
+        self.mschunks = {}
+        for ms in mslist:
+            obs = Observation(ms)
+            self.full_observations.append(obs)
+            self.mschunks[obs.name] = {"chunks": [], "parsets": []}
+        self.log = logging.getLogger("MSChunker")
+        self.log.setLevel(logging.INFO)
+
+    def chunk_observations_in_time(self, mintime: float, data_fraction: float = 1.0):
+        """
+        Break observations into smaller time chunks.
+
+        Chunking is done if the specified data_fraction < 1.
+
+        Parameters
+        ----------
+        data_fraction : float, optional
+            Fraction of data to use during processing
+        """
+        if data_fraction < 1.0:
+            self.log.info("Calculating time chunks")
+            self.observations = []
+            for obs in self.full_observations:
+                tottime = obs.endtime - obs.starttime
+                if data_fraction < min(1.0, mintime / tottime):
+                    obs.log.warning(
+                        "The specified value of data_fraction ({0:0.3f}) results in a "
+                        "total time for this observation that is less than the "
+                        "minimum timestep. The data fraction will be increased "
+                        "to {1:0.3f} to ensure the minimum timestep requirement is "
+                        "met.".format(data_fraction, min(1.0, mintime / tottime))
+                    )
+                nchunks = int(np.ceil(data_fraction / (mintime / tottime)))
+                if nchunks == 1:
+                    # Center the chunk around the midpoint (which is generally the most
+                    # sensitive, near transit)
+                    midpoint = obs.starttime + tottime / 2
+                    chunktime = min(tottime, max(mintime, data_fraction * tottime))
+                    if chunktime < tottime:
+                        sub_obs = Observation(
+                            obs.ms_filename,
+                            starttime=midpoint - chunktime / 2,
+                            endtime=midpoint + chunktime / 2,
+                        )
+                        self.observations.append(sub_obs)
+                    else:
+                        self.observations.append(obs)
+                else:
+                    obs.log.info("Splitting MS in {:d} chunks.".format(nchunks))
+                    steptime = (
+                        mintime * (tottime / mintime - nchunks) / nchunks
+                        + mintime
+                        + 0.1
+                    )
+                    starttimes = np.arange(obs.starttime, obs.endtime, steptime)
+                    endtimes = np.arange(
+                        obs.starttime + mintime, obs.endtime + mintime, steptime
+                    )
+                    for starttime, endtime in zip(starttimes, endtimes):
+                        if endtime > obs.endtime:
+                            starttime = obs.endtime - mintime
+                            endtime = obs.endtime
+                        sub_obs = Observation(
+                            obs.ms_filename, starttime=starttime, endtime=endtime
+                        )
+                        if sub_obs.name.endswith(".ms"):
+                            sub_obs.name = sub_obs.name.replace(
+                                ".ms", "_mjd{:f}.ms".format(sub_obs.starttime)
+                            )
+                        if sub_obs.name.endswith(".MS"):
+                            sub_obs.name = sub_obs.name.replace(
+                                ".MS", "_mjd{:f}.MS".format(sub_obs.starttime)
+                            )
+                        else:
+                            sub_obs.name = sub_obs.name + "_mjd{:f}".format(
+                                sub_obs.starttime
+                            )
+                        self.observations.append(sub_obs)
+                        self.mschunks[obs.name]["chunks"].append(sub_obs)
+        else:
+            self.observations = self.full_observations[:]
+
+    def chunk_observations_in_freq(
+        self, nchan: int = 0, round_freq_to_int: bool = False
+    ):
+        """
+        Break observations into smaller frequency chunks.
+
+        Chunking is done if the specified number of channels > 0.
+
+        Parameters
+        ----------
+        nchan : int, optional
+            Fraction of data to use during processing.
+
+        round_freq_to_int : int, optional
+            Round the frequency to the nearest integer value.
+        """
+        if nchan > 0:
+            self.log.info("Calculating frequency chunks")
+            self.observations = []
+            for obs in self.full_observations:
+                if nchan > obs.numchannels:
+                    obs.log.warning(
+                        "The specified number of channels exceeds the number of channels in the MeasurementSet."
+                    )
+                nchunks = int(obs.numchannels / nchan) if nchan > 0 else 1
+                obs.log.info("Splitting MS in {:d} chunks.".format(nchunks))
+                if nchunks == 1:
+                    self.observations.append(obs)
+                else:
+                    startfreqs = np.arange(
+                        obs.startfreq, obs.endfreq, obs.channelwidth * nchan
+                    )
+                    endfreqs = np.arange(
+                        obs.startfreq + obs.channelwidth * nchan,
+                        obs.endfreq + obs.channelwidth * nchan,
+                        obs.channelwidth * nchan,
+                    )
+                    for startfreq, endfreq in zip(startfreqs, endfreqs):
+                        sub_obs = Observation(
+                            obs.ms_filename, startfreq=startfreq, endfreq=endfreq
+                        )
+                        dfreq = (sub_obs.startfreq + sub_obs.endfreq) / 2e6
+                        if round_freq_to_int:
+                            dfreq = int(dfreq)
+                        if sub_obs.name.endswith(".ms"):
+                            sub_obs.name = sub_obs.name.replace(
+                                ".ms",
+                                "_{:f}MHz.ms".format(dfreq),
+                            )
+                        elif sub_obs.name.endswith(".MS"):
+                            sub_obs.name = sub_obs.name.replace(
+                                ".MS",
+                                "_{:f}MHz.MS".format(dfreq),
+                            )
+                        else:
+                            sub_obs.name = sub_obs.name + "_{:f}MHz".format(dfreq)
+                        self.observations.append(sub_obs)
+                        self.mschunks[obs.name]["chunks"].append(sub_obs)
+        else:
+            self.observations = self.full_observations[:]
+
+    def make_parsets_time(self):
+        """
+        Generate parsets to create the time chunks.
+        """
+        PARSET = """numthreads=12
+
+msin = {name}
+msin.starttime = {stime}
+msin.endtime = {etime}
+
+msout = {name_out}
+msout.storagemanager = dysco
+
+steps=[]
+"""
+        self.log.info("Writing chunk parsets")
+        for fobs in self.full_observations:
+            for i, obs in enumerate(self.mschunks[fobs.name]["chunks"]):
+                pname = "split_" + fobs.name + "_chunk{:02d}.parset".format(i)
+                with open(pname, "w") as f:
+                    parset = PARSET.format(
+                        name=fobs.name,
+                        stime=convert_mjd(obs.starttime),
+                        etime=convert_mjd(obs.endtime),
+                        name_out=obs.name,
+                    )
+                    # self.mschunks[obs.name]['parsets'][obs.name + "_{:f}".format(obs.starttime)] = pname
+                    self.mschunks[fobs.name]["parsets"].append(pname)
+                    f.write(parset)
+
+    def make_parsets_freq(self):
+        """
+        Generate parsets to create the time chunks.
+        """
+        PARSET = """numthreads=12
+
+msin = {name}
+msin.startchan = {sfreq}
+msin.nchan = {nchan}
+
+msout = {name_out}
+msout.storagemanager = dysco
+
+steps=[]
+"""
+        self.log.info("Writing chunk parsets")
+        for fobs in self.full_observations:
+            for i, obs in enumerate(self.mschunks[fobs.name]["chunks"]):
+                pname = "split_" + fobs.name + "_chunk{:02d}.parset".format(i)
+                with open(pname, "w") as f:
+                    parset = PARSET.format(
+                        name=fobs.name,
+                        sfreq=obs.startchan,
+                        nchan=(obs.endchan - obs.startchan),
+                        name_out=obs.name,
+                    )
+                    # self.mschunks[obs.name]['parsets'][obs.name + "_{:f}".format(obs.starttime)] = pname
+                    self.mschunks[fobs.name]["parsets"].append(pname)
+                    f.write(parset)
+
+    def run_parsets(self):
+        """
+        Run all the parsets, generating chunks.
+        """
+        self.log.info("Writing parsets")
+        for mschunk in self.mschunks.keys():
+            print(mschunk)
+            for parset in self.mschunks[mschunk]["parsets"]:
+                self.log.info("Running {:s}".format(parset))
+                subprocess.run(["DP3", parset])
+
+    def concat_chunks(self) -> list:
+        """
+        Concatenate all generated chunks.
+
+        Returns
+        -------
+        msout : list
+            Name of the concatenated MeasurementSet as <input name>.<int(data_fraction*100)>pc.ms'
+        """
+        self.log.info("Concatenating chunks")
+        msouts = []
+        for mschunk in self.mschunks.keys():
+            chunknames = [obs.name for obs in self.mschunks[mschunk]["chunks"]]
+            msout = mschunk + ".{:d}pc.ms".format(int(self.time_fraction * 100))
+            cmd = concat_time_command(chunknames, msout)
+            subprocess.run(cmd)
+            self.log.info("Concatenated MS written to {:s}".format(msout))
+            msouts.append(msout)
+        return msouts
+
+
+def chunk_and_concat(
+    mslist,
+    mode: str,
+    fraction: float = 1.0,
+    mintime: int = 1,
+    nchan: int = 0,
+    round_freq: bool = False,
+) -> list:
+    """
+    Main entry point. Splits a MeasurementSet in time chunks and concatenate all generated chunks.
+    Parameters
+    ----------
+    mslist : list
+        List of MS filenames to be chunked.
+    fraction : float
+        Fraction of time to take.
+    mintime : int
+        Minimum time span per chunk in seconds.
+    mode : str
+        Mode of chunking: time or frequency.
+    round_freq : bool
+        Round the frequency to an integer value.
+
+    Returns
+    -------
+    msout : list
+        List of MS chunks.
+    """
+    if mode == "time":
+        data = MSChunker(mslist, fraction)
+        data.chunk_observations_in_time(
+            data_fraction=data.time_fraction, mintime=mintime
+        )
+        data.make_parsets_time()
+        epochs = []
+        for epoch in data.mschunks:
+            epochs.append([obs.name for obs in data.mschunks[epoch]["chunks"]])
+        data.run_parsets()
+        data.concat_chunks()
+        return epochs
+    elif mode == "frequency":
+        data = MSChunker(mslist)
+        data.chunk_observations_in_freq(nchan, round_freq_to_int=round_freq)
+        data.make_parsets_freq()
+        epochs = []
+        for epoch in data.mschunks:
+            epochs.append([obs.name for obs in data.mschunks[epoch]["chunks"]])
+        data.run_parsets()
+        return epochs
+
+
+if __name__ == "__main__":
+    parser = argparse.ArgumentParser(
+        description="Chunk a MeasurementSet in time or frequency."
+    )
+    parser.add_argument("ms", nargs="+", help="Input MeasurementSet(s).")
+    parser.add_argument(
+        "--timefraction",
+        type=float,
+        default=1.0,
+        help="Fraction of data to split off. Default: 1.0",
+    )
+    parser.add_argument(
+        "--mintime",
+        required=False,
+        type=int,
+        default=-1,
+        help="Minimum time in seconds. Default: -1 (all)",
+    )
+    parser.add_argument(
+        "--chan_per_chunk",
+        required=False,
+        type=int,
+        default=0,
+        help="Number of channels per output frequency chunk. Default: 0 (all)",
+    )
+    parser.add_argument(
+        "--mode",
+        required=True,
+        type=str,
+        choices=["time", "frequency"],
+        help="Chunk in time or frequency.",
+    )
+
+    parser.add_argument(
+        "--round_freq",
+        required=False,
+        action='store_true',
+        help="Round the frequency in the file name to an integer value.",
+    )
+
+    options = parser.parse_args()
+    if (options.timefraction < 1) and (options.chan_per_chunk > 0):
+        print("Splitting time and frequency simultaneously is not supported.")
+    elif options.timefraction < 1:
+        chunk_and_concat(
+            options.ms,
+            fraction=options.timefraction,
+            mintime=options.mintime,
+            mode=options.mode,
+        )
+    elif options.chan_per_chunk > 0:
+        chunk_and_concat(options.ms, mode=options.mode, nchan=options.chan_per_chunk, round_freq=options.round_freq)
diff --git a/scripts/sort_times_into_freqGroups.py b/scripts/sort_times_into_freqGroups.py
index 1d77c941cbc9e82fb4089875593fc8050fed0633..5ba72a841e52f36cd4281a5cce691173166352ed 100755
--- a/scripts/sort_times_into_freqGroups.py
+++ b/scripts/sort_times_into_freqGroups.py
@@ -172,8 +172,8 @@ def main(MSfile, numSB=10, DP3fill=True, stepname=None, mergeLastGroup=False, tr
         logging.warning("Bandwidth of files is smaller than 51% of the minimum frequency step between two files! (More than about half the data is missing.)")
     #the new output map
 
-    # add 1% of the SB badwidth in case maxfreq might be "exactly" on a group-border
-    maxfreq = np.max(freqliste)+freq_width*0.51
+    # add 1% of the SB badwidth in case maxfreq might be "exactly" on a group-border # set back to exactly 50%, otherwise we add dummy data unnecessarily if using numSB = -1
+    maxfreq = np.max(freqliste)+freq_width*0.50
     if firstSB != None:
         if freqliste[0] < 100e6:
             # LBA Data
diff --git a/steps/LoSoTo.Plot.cwl b/steps/LoSoTo.Plot.cwl
index 0a9c63bbcdcebd5a2302a99b391ff7c0dbcd9092..082620476554f8aeb1ace1691aecaa73bb4a9ddc 100644
--- a/steps/LoSoTo.Plot.cwl
+++ b/steps/LoSoTo.Plot.cwl
@@ -18,7 +18,7 @@ requirements:
         entry: $(get_losoto_config('PLOT').join('\n'))
       - entryname: $(inputs.input_h5parm.basename)
         entry: $(inputs.input_h5parm)
-        writable: true
+        writable: false
       - entryname: run_losoto_plot.sh
         entry: |
           #!/bin/bash
@@ -48,6 +48,8 @@ inputs:
   - id: execute
     type: boolean?
     default: true
+  - id: selfcal_strategy
+    type: string?
   - id: soltab
     type: string
     doc: "Tabs to plot"
diff --git a/steps/applybeam.cwl b/steps/applybeam.cwl
index f360dcdd1bbde1d4252bf3b8810b3e3877fb86da..3c66b11f049e27eed28ed80c60edee8df36f624b 100644
--- a/steps/applybeam.cwl
+++ b/steps/applybeam.cwl
@@ -54,22 +54,26 @@ inputs:
        prefix: msout.storagemanager.databitrate=
        separate: false
   - id: updateweights
-    type: string?
+    type: boolean?
     inputBinding:
       position: 0
-      prefix: applybeam.updateweights=
+      prefix: applybeam.updateweights=True
       separate: false
   - id: usechannelfreq
-    type: string?
+    type: boolean?
+    default: true
     inputBinding:
+      valueFrom: $(!self)
       position: 0
-      prefix: applybeam.usechannelfreq=
+      prefix: applybeam.usechannelfreq=False
       separate: false
   - id: invert
-    type: string?
+    type: boolean?
+    default: true
     inputBinding:
+      valueFrom: $(!self)
       position: 0
-      prefix: applybeam.invert=
+      prefix: applybeam.invert=False
       separate: false
   - id: beammode
     type: string?
diff --git a/steps/chunk_ms.cwl b/steps/chunk_ms.cwl
new file mode 100644
index 0000000000000000000000000000000000000000..ad0b33c3b57e5d8bf2fbb57761d3f1a5b38073d6
--- /dev/null
+++ b/steps/chunk_ms.cwl
@@ -0,0 +1,67 @@
+class: CommandLineTool
+cwlVersion: v1.2
+id: chunkms
+label: ChunkMS
+
+baseCommand:
+  - MSChunker.py
+
+inputs:
+  - id: ms
+    type: Directory
+    inputBinding:
+      position: 6
+    doc: Input MeasurementSet.
+  - id: mode
+    type: string
+    inputBinding:
+      position: 1
+      prefix: --mode
+    doc: Mode to split the MeasurementSet in, time or frequency.
+  - id: timefraction
+    type: float?
+    default: 1.0
+    inputBinding:
+      position: 2
+      prefix: --timefraction
+    doc: The fraction of data to split off, evenly distributed in time.
+  - id: mintime
+    type: int?
+    default: -1
+    inputBinding:
+      position: 3
+      prefix: --mintime
+    doc: Minimum number of seconds per split off time chunk.
+  - id: chan_per_chunk
+    type: int?
+    default: 0
+    inputBinding:
+      position: 4
+      prefix: --chan_per_chunk
+    doc: Number of channels per split off frequency chunk.
+  - id: round_freq
+    type: boolean
+    default: true
+    inputBinding:
+      position: 5
+      prefix: --round_freq
+    doc: Rounds the frequency in the directory name to an integer.
+
+outputs:
+  - id: msouts
+    type: Directory[]
+    doc: Output MeasurementSets.
+    outputBinding:
+      glob: '*MHz.ms'
+
+hints:
+ - class: DockerRequirement
+   dockerPull: astronrd/linc
+
+requirements:
+ - class: InitialWorkDirRequirement
+   listing:
+    - entry: $(inputs.ms)
+ - class: InlineJavascriptRequirement
+stderr: chunk_ms.log
+
diff --git a/steps/ddecal.cwl b/steps/ddecal.cwl
index 016ff570ffc73da0d7dd5bab511576098234a75f..d8fd2aead14c2c322331fdaa439653a669c3e7e6 100644
--- a/steps/ddecal.cwl
+++ b/steps/ddecal.cwl
@@ -235,6 +235,35 @@ inputs:
       position: 0
       prefix: ddecal.approxtolerance=
       separate: false
+  - id: sourcedb
+    type:
+      - File?
+      - Directory?
+    inputBinding:
+      position: 0
+      prefix: ddecal.sourcedb=
+      separate: false
+  - id: usebeammodel
+    default: false
+    type: boolean?
+    inputBinding:
+      position: 0
+      prefix: ddecal.usebeammodel=True
+      separate: false
+  - id: usechannelfreq
+    default: true
+    type: boolean?
+    inputBinding:
+      valueFrom: $(!self)
+      position: 0
+      prefix: ddecal.usechannelfreq=False
+      separate: false
+  - id: beammode
+    type: string?
+    inputBinding:
+      position: 0
+      prefix: ddecal.beammode=
+      separate: false
 
 #--------------------
   - id: save2json
diff --git a/steps/filter_predict.cwl b/steps/filter_predict.cwl
index ab3aedde6451a40c16cafd011fec075879719301..ca71e944ab3b2b45ccb404b74af75bf011cf4a17 100644
--- a/steps/filter_predict.cwl
+++ b/steps/filter_predict.cwl
@@ -55,6 +55,13 @@ inputs:
     inputBinding:
       position: 0
       prefix: predict.usebeammodel=True
+  - id: usechannelfreq
+    default: true
+    type: boolean?
+    inputBinding:
+      valueFrom: $(!self)
+      position: 0
+      prefix: predict.usechannelfreq=False
   - default: false
     id: onebeamperpatch
     type: boolean?
@@ -118,7 +125,6 @@ outputs:
       glob: 'filter_predict*.log'
 arguments:
   - steps=[filter,predict,count]
-  - predict.usechannelfreq=False
   - msout=.
 requirements:
   - class: InitialWorkDirRequirement
diff --git a/steps/gaincal.cwl b/steps/gaincal.cwl
index e8061e72ed5d920d6b7b03cd69ffaf182726241f..2582d0801be6d09888210767fd36875b16e12579 100644
--- a/steps/gaincal.cwl
+++ b/steps/gaincal.cwl
@@ -114,9 +114,11 @@ inputs:
   - default: true
     id: usechannelfreq
     type: boolean?
+    default: true
     inputBinding:
+      valueFrom: $(!self)
       position: 0
-      prefix: gaincal.usechannelfreq=True
+      prefix: gaincal.usechannelfreq=False
       separate: false
   - default: array_factor
     id: beammode
diff --git a/steps/predict.cwl b/steps/predict.cwl
index 33114b9c5b725edd093b175e1cb098106cc8a0d3..647f6bc9f25add1c477f6ebc139ca030bd1a582d 100644
--- a/steps/predict.cwl
+++ b/steps/predict.cwl
@@ -56,11 +56,12 @@ inputs:
       position: 0
       prefix: predict.usebeammodel=True
   - id: usechannelfreq
-    default: false
+    default: true
     type: boolean?
     inputBinding:
+      valueFrom: $(!self)
       position: 0
-      prefix: predict.usechannelfreq=True
+      prefix: predict.usechannelfreq=False
   - id: onebeamperpatch
     default: false
     type: boolean?
diff --git a/steps/wsclean.cwl b/steps/wsclean.cwl
index e435d8bb849b0e102bfbeac5f041edc9cdf7da6c..8a759bbaad19691ae1fe28d0418b46330e4f9fe4 100644
--- a/steps/wsclean.cwl
+++ b/steps/wsclean.cwl
@@ -16,7 +16,7 @@ inputs:
       - 2500
       - 2500
     type: 
-      - int[]
+      - int[]?
     inputBinding:
       position: 1
       shellQuote: false
@@ -315,6 +315,10 @@ inputs:
   - id: fits_image
     type: File?
 outputs:
+  - id: sourcelist
+    type: File?
+    outputBinding:
+      glob: [$(inputs.image_name)-sources.txt]
   - id: dirty_image
     type: File?
     outputBinding:
diff --git a/test_jobs/HBA_target_selfcal.json b/test_jobs/HBA_target_selfcal.json
new file mode 100644
index 0000000000000000000000000000000000000000..6661ede0a96cf4ebd4506ede02c80775513ed6b3
--- /dev/null
+++ b/test_jobs/HBA_target_selfcal.json
@@ -0,0 +1,102 @@
+{
+    "msin": [
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB000_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB001_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB002_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB003_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB004_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB005_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB006_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB007_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB008_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB009_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB010_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB011_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB012_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB013_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB014_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB015_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB016_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB017_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB018_uv.MS"
+        },
+        {
+            "class": "Directory",
+            "path": "../data/L667520_SB019_uv.MS"
+        }
+    ],
+    "avg_timeresolution_concat": 4,
+    "ncores": 12,
+    "cal_solutions": {
+        "class": "File",
+        "path": "../data/results_calibrator/cal_solutions.h5"
+    },
+    "A-Team_skymodel": {
+        "class": "File",
+        "path": "/usr/local/share/linc/skymodels/Ateam_LBA_CC.skymodel"
+    },
+    "calibrator_path_skymodel": {
+        "class": "Directory",
+        "path": "/usr/local/share/linc/skymodels"
+    },
+    "selfcal": true,
+    "selfcal_strategy": "HBA",
+    "selfcal_hba_imsize": [2000,2000],
+    "num_SBs_per_group": -1
+}
diff --git a/workflows/HBA_target.cwl b/workflows/HBA_target.cwl
index 33a79200fd1cf22874da50313a2bb6a4eaef2d53..27edfe6b3d3170df9e930fd8cdad731757e98dda 100644
--- a/workflows/HBA_target.cwl
+++ b/workflows/HBA_target.cwl
@@ -108,7 +108,8 @@ inputs:
     default: 97.64kHz
   - id: num_SBs_per_group
     type: int?
-    default: 10
+  - id: calib_nchan
+    type: int?
   - id: reference_stationSB
     type: int?
     default: null
@@ -161,6 +162,20 @@ inputs:
   - id: selfcal
     type: boolean?
     default: false
+  - id: selfcal_strategy
+    type: string?
+    default: 'HBA'
+  - id: selfcal_hba_imsize
+    type: int[]?
+    default: [20000,20000]
+  - id: selfcal_hba_uvlambdamin
+    type: float?
+    default: 200.
+  - id: output_channels_per_chunk
+    type: int?
+    default: 20
+  - id: selfcal_region
+    type: File?
   - id: chunkduration
     type: float?
     default: 0.0
@@ -169,8 +184,6 @@ inputs:
   - id: make_structure_plot
     type: boolean?
     default: true
-  - id: selfcal_region
-    type: File?
   - id: skymodel_fluxlimit
     type: float?
 outputs:
@@ -267,6 +280,8 @@ steps:
         source: avg_freqresolution_concat
       - id: num_SBs_per_group
         source: num_SBs_per_group
+      - id: calib_nchan
+        source: calib_nchan
       - id: reference_stationSB
         source: reference_stationSB
       - id: ionex_server
@@ -301,16 +316,24 @@ steps:
         source: aoflag_freqconcat
       - id: selfcal
         source: selfcal
+      - id: selfcal_strategy
+        source: selfcal_strategy
+      - id: selfcal_hba_imsize
+        source: selfcal_hba_imsize
+      - id: selfcal_region
+        source: selfcal_region
       - id: chunkduration
         source: chunkduration
       - id: wsclean_tmpdir
         source: wsclean_tmpdir
       - id: make_structure_plot
         source: make_structure_plot
-      - id: selfcal_region
-        source: selfcal_region
       - id: skymodel_fluxlimit
         source: skymodel_fluxlimit
+      - id: output_channels_per_chunk
+        source: output_channels_per_chunk
+      - id: selfcal_hba_uvlambdamin
+        source: selfcal_hba_uvlambdamin
     out:
       - id: logfiles
       - id: msout
diff --git a/workflows/LBA_target.cwl b/workflows/LBA_target.cwl
index 646da039d6339046ba56ea6a1804a36f2e4bc0ec..49f9d45b0b8358620922caeebb292e75611d71ff 100644
--- a/workflows/LBA_target.cwl
+++ b/workflows/LBA_target.cwl
@@ -108,7 +108,8 @@ inputs:
     default: 48.82kHz
   - id: num_SBs_per_group
     type: int?
-    default: -1
+  - id: calib_nchan
+    type: int?
   - id: reference_stationSB
     type: int?
     default: null
@@ -161,6 +162,14 @@ inputs:
   - id: selfcal
     type: boolean?
     default: true
+  - id: selfcal_strategy
+    type: string?
+    default: 'LBA'
+  - id: selfcal_region
+    type: File?
+  - id: selfcal_hba_imsize
+    type: int[]?
+    default: [20000,20000]
   - id: chunkduration
     type: float?
     default: 3600.0
@@ -169,8 +178,6 @@ inputs:
   - id: make_structure_plot
     type: boolean?
     default: true
-  - id: selfcal_region
-    type: File?
   - id: skymodel_fluxlimit
     type: float?
     default: 1.0
@@ -268,6 +275,8 @@ steps:
         source: avg_freqresolution_concat
       - id: num_SBs_per_group
         source: num_SBs_per_group
+      - id: calib_nchan
+        source: calib_nchan
       - id: reference_stationSB
         source: reference_stationSB
       - id: ionex_server
@@ -302,14 +311,18 @@ steps:
         source: aoflag_freqconcat
       - id: selfcal
         source: selfcal
+      - id: selfcal_strategy
+        source: selfcal_strategy
+      - id: selfcal_hba_imsize
+        source: selfcal_hba_imsize
+      - id: selfcal_region
+        source: selfcal_region
       - id: chunkduration
         source: chunkduration
       - id: wsclean_tmpdir
         source: wsclean_tmpdir
       - id: make_structure_plot
         source: make_structure_plot
-      - id: selfcal_region
-        source: selfcal_region
       - id: skymodel_fluxlimit
         source: skymodel_fluxlimit
     out:
diff --git a/workflows/linc_calibrator/apply_calibrate_bp.cwl b/workflows/linc_calibrator/apply_calibrate_bp.cwl
index 69cfef3a1068eb392a7845486e4e2d5962ff3b1f..7b9001f1769ea760bea8923c31991620b8cd5022 100644
--- a/workflows/linc_calibrator/apply_calibrate_bp.cwl
+++ b/workflows/linc_calibrator/apply_calibrate_bp.cwl
@@ -121,11 +121,11 @@ steps:
       - id: databitrate
         default: 0
       - id: updateweights
-        default: 'true'
+        default: true
       - id: usechannelfreq
-        default: 'false'
+        default: false
       - id: invert
-        default: 'true'
+        default: true
       - id: beammode
         default: element
     out:
diff --git a/workflows/linc_calibrator/apply_calibrate_pa.cwl b/workflows/linc_calibrator/apply_calibrate_pa.cwl
index f14b9d4b0913f1995207bf1a800df0a75281bcc2..6289863a120521ab20c53f4ab6183c794ce0ee41 100644
--- a/workflows/linc_calibrator/apply_calibrate_pa.cwl
+++ b/workflows/linc_calibrator/apply_calibrate_pa.cwl
@@ -82,13 +82,13 @@ steps:
       - id: databitrate
         default: 0
       - id: updateweights
-        default: 'true'
+        default: true
       - id: invert
-        default: 'true'
+        default: true
       - id: beammode
         default: element
       - id: usechannelfreq
-        default: 'false'
+        default: false
       - id: msin
         source: applyPA/msout
       - id: type
diff --git a/workflows/linc_target.cwl b/workflows/linc_target.cwl
index d776ce1f9a6d61814bbf789784e6e67c6ae528ab..6174b8e867f74fc44aa1535fdf44efa862c6989d 100644
--- a/workflows/linc_target.cwl
+++ b/workflows/linc_target.cwl
@@ -108,7 +108,8 @@ inputs:
     default: 97.64kHz
   - id: num_SBs_per_group
     type: int?
-    default: 10
+  - id: calib_nchan
+    type: int?
   - id: reference_stationSB
     type: int?
     default: null
@@ -161,6 +162,18 @@ inputs:
   - id: selfcal
     type: boolean?
     default: false
+  - id: selfcal_strategy
+    type: string?
+    default: 'HBA'
+  - id: selfcal_hba_imsize
+    type: int[]?
+    default: [20000,20000]
+  - id: selfcal_hba_uvlambdamin
+    type: float?
+    default: 200
+  - id: output_channels_per_chunk
+    type: int?
+    default: 20
   - id: chunkduration
     type: float?
     default: 0.0
@@ -225,7 +238,10 @@ steps:
       - id: process_baselines_target
         source: process_baselines_target
       - id: num_SBs_per_group
-        source: num_SBs_per_group
+        source: 
+          - num_SBs_per_group
+          - selfcal
+        valueFrom: '$(self[0] == null ? (self[-1] ? -1 : null) : self[0])'
       - id: reference_stationSB
         source: reference_stationSB
       - id: filter_baselines
@@ -358,6 +374,10 @@ steps:
         source: prep/outh5parm
       - id: selfcal
         source: selfcal
+      - id: selfcal_strategy
+        source: selfcal_strategy
+      - id: selfcal_hba_imsize
+        source: selfcal_hba_imsize
       - id: chunkduration
         source: chunkduration
       - id: wsclean_tmpdir
@@ -372,6 +392,13 @@ steps:
         source: prep/filenames
       - id: selfcal_region
         source: selfcal_region
+      - id: calib_nchan
+        source: 
+          - calib_nchan
+          - selfcal
+        valueFrom: '$(self[0] == null ? (self[1] ? null : 0) : self[0])'
+      - id: selfcal_hba_uvlambdamin
+        source: selfcal_hba_uvlambdamin
     out:
       - id: msout
       - id: outh5parm
@@ -438,6 +465,12 @@ steps:
         source: wsclean_tmpdir
       - id: make_structure_plot
         source: make_structure_plot
+      - id: selfcal
+        source: selfcal
+      - id: selfcal_strategy
+        source: selfcal_strategy
+      - id: output_channels_per_chunk
+        source: output_channels_per_chunk
     out:
       - id: msout
       - id: solutions
@@ -449,4 +482,6 @@ steps:
 requirements:
   - class: SubworkflowFeatureRequirement
   - class: ScatterFeatureRequirement
-  - class: MultipleInputFeatureRequirement
\ No newline at end of file
+  - class: MultipleInputFeatureRequirement
+  - class: StepInputExpressionRequirement
+  - class: InlineJavascriptRequirement
diff --git a/workflows/linc_target/calib_targ.cwl b/workflows/linc_target/calib_targ.cwl
index 945be934c128f08b19eb45a316365990e6aed9b2..d55d0f8cc5455f8353b145a65489fc67b08c709b 100644
--- a/workflows/linc_target/calib_targ.cwl
+++ b/workflows/linc_target/calib_targ.cwl
@@ -9,8 +9,8 @@ inputs:
     type: Directory
   - id: skymodel
     type:
-      - File
-      - Directory
+      - File?
+      - Directory?
   - id: do_smooth
     type: boolean?
     default: false
@@ -20,9 +20,22 @@ inputs:
   - id: gsmcal_step
     type: string?
     default: 'phase'
-  - id: execute
-    type: boolean?
-    default: true
+  - id: smoothnessconstraint
+    type: float?
+    default: 0.
+  - id: smoothnessreffrequency
+    type: float?
+    default: 0.
+  - id: beammode
+    type: string?
+  - id: nchan
+    type: int?
+    default: 1
+  - id: model_column
+    type: string[]?
+    default: []
+  - id: uvlambdamin
+    type: float?
 outputs:
   - id: msout
     outputSource:
@@ -64,13 +77,9 @@ steps:
         source: BLsmooth/msout
       - id: msin_datacolumn
         default: SMOOTHED_DATA
-      - id: msout_name
-        default: '.'
-      - id: blrange
-        default:
-          - 150
-          - 9999999
-      - id: caltype
+      - id: uvlambdamin
+        source: uvlambdamin
+      - id: mode
         source: gsmcal_step
         valueFrom: '$(self == "phase" ? "phaseonly" : self)'
       - id: sourcedb
@@ -80,22 +89,28 @@ steps:
       - id: solint
         default: 1
       - id: nchan
-        default: 0
+        source: nchan
       - id: tolerance
         default: 1e-3
-      - id: propagatesolutions
+      - id: propagate_solutions
         source: propagatesolutions
       - id: usebeammodel
         default: true
       - id: usechannelfreq
         default: true
       - id: beammode
-        default: array_factor
+        source: beammode
+      - id: smoothnessconstraint
+        source: smoothnessconstraint
+      - id: smoothnessreffrequency
+        source: smoothnessreffrequency
+      - id: modeldatacolumns
+        source: model_column
     out:
       - id: msout
       - id: h5parm
       - id: logfile
-    run: ../../steps/gaincal.cwl
+    run: ../../steps/ddecal.cwl
   - id: concat_logfiles_gaincal
     in:
       - id: file_list
diff --git a/workflows/linc_target/dp3_prep_targ.cwl b/workflows/linc_target/dp3_prep_targ.cwl
index 37ddeac9f9d6dcbabfec073411ea47d4fe98bc77..778520f7421a6bf8cf645485227fe2ef3b287490 100644
--- a/workflows/linc_target/dp3_prep_targ.cwl
+++ b/workflows/linc_target/dp3_prep_targ.cwl
@@ -281,6 +281,8 @@ steps:
         default: 0
       - id: filter_baselines
         source: process_baselines_target
+      - id: usechannelfreq
+        default: false
       - id: execute
         source: clipAteam
     out:
diff --git a/workflows/linc_target/finalize.cwl b/workflows/linc_target/finalize.cwl
index f52ae1143bbcc8286fe3cf9b3135de99e49afe52..6dd849a62a363ea8e6960adcf60201139254cd83 100644
--- a/workflows/linc_target/finalize.cwl
+++ b/workflows/linc_target/finalize.cwl
@@ -66,18 +66,29 @@ inputs:
   - id: make_structure_plot
     type: boolean?
     default: true
+  - id: selfcal
+    type: boolean?
+    default: false
+  - id: selfcal_strategy
+    type: string?
+    default: 'HBA'
+  - id: output_channels_per_chunk
+    type: int?
+    default: 20
 outputs:
   - id: msout
     outputSource:
-      - apply_gsmcal/msout
+      - split_into_chunks/msouts
+      - apply_gsmcal_final/msout
     type: Directory[]
+    pickValue: first_non_null
   - id: solutions
     outputSource:
       - h5parm_pointingname/outh5parm
     type: File
   - id: logfiles
     outputSource:
-      - concat_logfiles_applygsm/output
+      - concat_logfiles_applygsm_final/output
       - concat_logfiles_solutions/output
       - concat_logfiles_structure/output
       - concat_logfiles_wsclean/output
@@ -111,12 +122,12 @@ steps:
         default: target
       - id: soltab_in
         source: gsmcal_step
-        valueFrom: $(self+'000')
+        valueFrom: $(self == "phase" || self == "scalarphase" ? 'phase000':self+'000')
       - id: soltab_out
         source:
           - skymodel_source
           - gsmcal_step
-        valueFrom: $(self.join(''))
+        valueFrom: $(self.join('')+"_final")
       - id: filter
         source: process_baselines_target
       - id: bad_antennas
@@ -126,7 +137,7 @@ steps:
       - id: log
     run: ../../steps/add_missing_stations.cwl
     label: add_missing_stations
-  - id: apply_gsmcal
+  - id: apply_gsmcal_final
     in:
       - id: msin
         source: msin
@@ -145,7 +156,7 @@ steps:
         source:
           - skymodel_source
           - gsmcal_step
-        valueFrom: $(self.join(''))
+        valueFrom: $(self.join('')+"_final")
       - id: solset
         default: target
       - id: storagemanager
@@ -157,7 +168,7 @@ steps:
       - id: flagged_fraction_dict
       - id: logfile
     run: ../../steps/applytarget.cwl
-    label: apply_gsmcal
+    label: apply_gsmcal_final
     scatter:
       - msin
       - msout_name
@@ -166,7 +177,7 @@ steps:
     in:
       - id: flagged_fraction_dict
         source:
-          - apply_gsmcal/flagged_fraction_dict
+          - apply_gsmcal_final/flagged_fraction_dict
       - id: filter_station
         default: ''
       - id: state
@@ -175,14 +186,29 @@ steps:
       - id: flagged_fraction_antenna
     run: ./../../steps/findRefAnt_join.cwl
     label: final_flags_join
+  - id: split_into_chunks
+    in:
+      - id: selfcal
+        source: selfcal
+      - id: selfcal_strategy
+        source: selfcal_strategy
+      - id: msin
+        source: apply_gsmcal_final/msout
+      - id: output_channels_per_chunk
+        source: output_channels_per_chunk
+    out:
+      - id: msouts
+    run: split_into_chunks.cwl
+    when: $(inputs.selfcal && inputs.selfcal_strategy == 'HBA')
+    label: chunk_ms
   - id: average
     in:
       - id: msin
-        source: apply_gsmcal/msout
+        source: apply_gsmcal_final/msout
       - id: msout_name
         linkMerge: merge_flattened
         source:
-          - apply_gsmcal/msout
+          - apply_gsmcal_final/msout
         valueFrom: $(self.nameroot+'_wsclean.ms')
       - id: msin_datacolumn
         default: DATA
@@ -249,7 +275,7 @@ steps:
   - id: uvplot
     in:
       - id: MSfiles
-        source: apply_gsmcal/msout
+        source: apply_gsmcal_final/msout
       - id: output_name
         source: targetname
         valueFrom: $(self)_uv-coverage.png
@@ -324,7 +350,7 @@ steps:
     in:
       - id: input
         source:
-          - apply_gsmcal/logfile
+          - apply_gsmcal_final/logfile
     out:
       - id: output
     run: ../../steps/merge_array_files.cwl
@@ -339,7 +365,7 @@ steps:
         source:
           - skymodel_source
           - gsmcal_step
-        valueFrom: $(self.join(''))
+        valueFrom: $(self.join('')+"_final")
       - id: input_file
         source: insolutions
       - id: squeeze
@@ -386,17 +412,17 @@ steps:
     run: ../../steps/LoSoTo.Structure.cwl
     label: structure_function
     when: $(inputs.soltab.includes("phase") && inputs.execute)
-  - id: concat_logfiles_applygsm
+  - id: concat_logfiles_applygsm_final
     in:
       - id: file_list
         source:
           - merge_array_files/output
       - id: file_prefix
-        default: apply_gsmcal
+        default: apply_gsmcal_final
     out:
       - id: output
     run: ../../steps/concatenate_files.cwl
-    label: concat_logfiles_applygsm
+    label: concat_logfiles_applygsm_final
   - id: concat_logfiles_solutions
     in:
       - id: file_list
@@ -407,7 +433,7 @@ steps:
           - write_solutions/log
           - h5parm_pointingname/log
       - id: file_prefix
-        default: losoto_gsmcal
+        default: losoto_gsmcal_final
     out:
       - id: output
     run: ../../steps/concatenate_files.cwl
@@ -473,3 +499,4 @@ requirements:
   - class: InlineJavascriptRequirement
   - class: ScatterFeatureRequirement
   - class: MultipleInputFeatureRequirement
+  - class: SubworkflowFeatureRequirement
diff --git a/workflows/linc_target/gsmcal.cwl b/workflows/linc_target/gsmcal.cwl
index 38e944dab00a9a505402b7cd56c63b768392a297..e2a9bc5bf61ced766847201cbea5afcfe0d634a9 100644
--- a/workflows/linc_target/gsmcal.cwl
+++ b/workflows/linc_target/gsmcal.cwl
@@ -80,13 +80,25 @@ inputs:
     type: File
   - id: selfcal_region
     type: File?
+  - id: selfcal_strategy
+    type: string?
+    default: 'HBA'
+  - id: selfcal_hba_imsize
+    type: int[]?
+    default: [20000,20000]
+  - id: selfcal_hba_uvlambdamin
+    type: float?
+    default: 200.
+  - id: calib_nchan
+    type: int?
+    default: 1
 outputs:
   - id: msout
     outputSource:
+      - selfcal_target_hba/msout
+      - selfcal_target_lba/msout
       - calibrate_target/msout
-      - selfcal_target/msout
-    pickValue: all_non_null
-    linkMerge: merge_flattened
+    pickValue: first_non_null
     type: Directory[]
   - id: outh5parm
     outputSource:
@@ -94,7 +106,8 @@ outputs:
     type: File
   - id: outsolutions
     outputSource:
-      - selfcal_target/outsolutions
+      - selfcal_target_hba/outsolutions
+      - selfcal_target_lba/outsolutions
       - insolutions
     pickValue: first_non_null
     type: File
@@ -118,7 +131,9 @@ outputs:
       - losoto_plot_Pd2/output_plots
       - losoto_plot_tec/output_plots
       - plot_unflagged/output_imag
-      - selfcal_target/inspection
+      - selfcal_target_lba/inspection
+      - selfcal_target_hba/image_before
+      - selfcal_target_hba/inspection
     type: File[]
     linkMerge: merge_flattened
     pickValue: all_non_null
@@ -135,7 +150,8 @@ outputs:
       - concat_logfiles_blsmooth/output
       - concat_logfiles_unflagged/output
       - aoflag/logfile
-      - selfcal_target/logfiles
+      - selfcal_target_lba/logfiles
+      - selfcal_target_hba/logfiles
     type: File[]
     linkMerge: merge_flattened
     pickValue: all_non_null
@@ -226,7 +242,7 @@ steps:
           - calibrate_target/gaincal.log
         pickValue: all_non_null
       - id: file_prefix
-        default: gaincal
+        default: calibrate_target
       - id: execute
         source: selfcal
     out:
@@ -247,7 +263,7 @@ steps:
           - losoto_plot_Pd2/logfile
         pickValue: all_non_null
       - id: file_prefix
-        default: losoto_gsmcal
+        default: losoto_gsmcal_final
     out:
       - id: output
     run: ../../steps/concatenate_files.cwl
@@ -315,6 +331,14 @@ steps:
     label: concat
     scatter:
       - group_id
+  - id: merge_array_concat
+    in:
+      - id: input
+        source: concat/msout
+    out:
+      - id: output
+    run: ../../steps/merge_array.cwl
+    label: merge_array_concat
   - id: aoflag
     in:
       - id: msin
@@ -334,6 +358,7 @@ steps:
       - id: logfile
     run: ../../steps/aoflag.cwl
     label: aoflag
+
   - id: check_unflagged_fraction
     in:
       - id: msin
@@ -350,6 +375,7 @@ steps:
     label: check_unflagged_fraction
     scatter:
       - msin
+
   - id: merge_array
     in:
       - id: input
@@ -358,14 +384,7 @@ steps:
       - id: output
     run: ../../steps/merge_array.cwl
     label: merge_array
-  - id: merge_array_concat
-    in:
-      - id: input
-        source: concat/msout
-    out:
-      - id: output
-    run: ../../steps/merge_array.cwl
-    label: merge_array_concat
+
   - id: merge_array_files
     in:
       - id: input
@@ -374,6 +393,7 @@ steps:
       - id: output
     run: ../../steps/merge_array_files.cwl
     label: merge_array_files
+
   - id: check_filtered_MS_array
     in:
       - id: input
@@ -382,6 +402,8 @@ steps:
       - id: output
     run: ../../steps/check_filtered_MS_array.cwl
     label: check_filtered_MS_array
+  # END concatenation into a single 48 MHz chunk
+
   - id: calibrate_target
     in:
       - id: max_dp3_threads
@@ -396,8 +418,18 @@ steps:
         source: propagatesolutions
       - id: gsmcal_step
         source: gsmcal_step
-      - id: execute
-        source: selfcal
+      - id: selfcal_strategy
+        source: selfcal_strategy
+      - id: nchan
+        source: calib_nchan
+      - id: smoothnessconstraint
+        default: 5e6
+      - id: smoothnessreffrequency
+        default: 12e7
+      - id: beammode
+        default: array_factor
+      - id: uvlambdamin
+        default: 200
     out:
       - id: msout
       - id: BLsmooth.log
@@ -407,8 +439,60 @@ steps:
     label: calibrate_target
     scatter:
       - msin
-    when: $(!inputs.execute)
-  - id: selfcal_target
+    when: $(inputs.selfcal_strategy == 'HBA')
+
+  - id: selfcal_target_hba
+    in:
+      - id: max_dp3_threads
+        source: max_dp3_threads
+      - id: msin
+        source:
+          - calibrate_target/msout
+        pickValue: all_non_null
+      - id: do_smooth
+        source: do_smooth
+      - id: propagatesolutions
+        source: propagatesolutions
+      - id: preapply_h5parm
+        source:
+          - calibrate_target/outh5parm
+        pickValue: all_non_null
+      - id: refant
+        source: findRefAnt_join/refant
+      - id: skymodel_source
+        source: skymodel_source
+      - id: gsmcal_step
+        source: gsmcal_step
+      - id: process_baselines_target
+        source: process_baselines_target
+      - id: bad_antennas
+        source: identifybadantennas_join/filter_out
+      - id: insolutions
+        source: insolutions
+      - id: wsclean_tmpdir
+        source: wsclean_tmpdir
+      - id: execute
+        source: selfcal
+      - id: selfcal_strategy
+        source: selfcal_strategy
+      - id: selfcal_hba_imsize
+        source: selfcal_hba_imsize
+      - id: calib_nchan
+        source: calib_nchan
+      - id: calib_uvlambdamin
+        source: selfcal_hba_uvlambdamin
+    out:
+      - id: msout
+      - id: outh5parm
+      - id: image_before
+      - id: outsolutions
+      - id: inspection
+      - id: logfiles
+    run: ./selfcal_targ_hba.cwl
+    label: selfcal_target_hba
+    when: $(inputs.execute && inputs.selfcal_strategy == 'HBA')
+
+  - id: selfcal_target_lba
     in:
       - id: max_dp3_threads
         source: max_dp3_threads
@@ -452,36 +536,41 @@ steps:
         source: aoflag_freqconcat
       - id: execute
         source: selfcal
+      - id: selfcal_strategy
+        source: selfcal_strategy
     out:
       - id: msout
       - id: logfiles
       - id: inspection
       - id: outh5parm
       - id: outsolutions
-    run: ./selfcal_targ.cwl
-    label: selfcal_target
-    when: $(inputs.execute)
+    run: ./selfcal_targ_lba.cwl
+    label: selfcal_target_lba
+    when: $(inputs.execute && inputs.selfcal_strategy == 'LBA')
+
   - id: h5parm_collector
     in:
       - id: h5parmFiles
         source:
+          - selfcal_target_hba/outh5parm
+          - selfcal_target_lba/outh5parm
           - calibrate_target/outh5parm
-          - selfcal_target/outh5parm
-        pickValue: all_non_null
-        linkMerge: merge_flattened
+        pickValue: first_non_null
       - id: squeeze
         default: true
       - id: verbose
         default: true
       - id: clobber
         default: true
-      - id: execute
-        source: selfcal
     out:
       - id: outh5parm
       - id: log
     run: ../../steps/H5ParmCollector.cwl
     label: H5parm_collector
+
+  # Collect the selfcal h5parm separately, such that the
+  # plotting steps below can simply plot the original solve
+
   - id: plot_unflagged
     in:
       - id: frequencies
@@ -496,7 +585,8 @@ steps:
   - id: losoto_plot_tec
     in:
       - id: input_h5parm
-        source: h5parm_collector/outh5parm
+        source:
+          - h5parm_collector/outh5parm
       - id: soltab
         default: sol000/tec000
       - id: axesInPlot
@@ -514,6 +604,8 @@ steps:
           - 0.5
       - id: prefix
         default: tec_
+      - id: selfcal_strategy
+        source: selfcal_strategy
       - id: execute
         source: selfcal
     out:
@@ -521,12 +613,13 @@ steps:
       - id: logfile
       - id: parset
     run: ../../steps/LoSoTo.Plot.cwl
-    when: $(inputs.execute)
+    when: $(inputs.execute && inputs.selfcal_strategy == "LBA")
     label: losoto_plot_tec
   - id: losoto_plot_P
     in:
       - id: input_h5parm
-        source: h5parm_collector/outh5parm
+        source:
+          - h5parm_collector/outh5parm
       - id: soltab
         default: sol000/phase000
       - id: axesInPlot
@@ -545,9 +638,11 @@ steps:
         source: findRefAnt_join/refant
       - id: prefix
         default: ph_
+      - id: selfcal_strategy
+        source: selfcal_strategy
       - id: execute
         source: selfcal
-        valueFrom: '$(self ? false : true)'
+        valueFrom: '$(self ? (false || inputs.selfcal_strategy == "HBA") : true)'
     out:
       - id: output_plots
       - id: logfile
@@ -558,7 +653,8 @@ steps:
   - id: losoto_plot_P2
     in:
       - id: input_h5parm
-        source: h5parm_collector/outh5parm
+        source:
+          - h5parm_collector/outh5parm
       - id: soltab
         default: sol000/phase000
       - id: axesInPlot
@@ -578,9 +674,11 @@ steps:
         source: findRefAnt_join/refant
       - id: prefix
         default: ph_
+      - id: selfcal_strategy
+        source: selfcal_strategy
       - id: execute
         source: selfcal
-        valueFrom: '$(self ? false : true)'
+        valueFrom: '$(self ? (false || inputs.selfcal_strategy == "HBA") : true)'
     out:
       - id: output_plots
       - id: logfile
@@ -591,7 +689,8 @@ steps:
   - id: losoto_plot_Pd
     in:
       - id: input_h5parm
-        source: h5parm_collector/outh5parm
+        source:
+          - h5parm_collector/outh5parm
       - id: soltab
         default: sol000/phase000
       - id: axesInPlot
@@ -612,9 +711,11 @@ steps:
         source: findRefAnt_join/refant
       - id: prefix
         default: ph_poldif
+      - id: selfcal_strategy
+        source: selfcal_strategy
       - id: execute
         source: selfcal
-        valueFrom: '$(self ? false : true)'
+        valueFrom: '$(self ? (false || inputs.selfcal_strategy == "HBA") : true)'
     out:
       - id: output_plots
       - id: logfile
@@ -625,7 +726,8 @@ steps:
   - id: losoto_plot_Pd2
     in:
       - id: input_h5parm
-        source: h5parm_collector/outh5parm
+        source:
+          - h5parm_collector/outh5parm
       - id: soltab
         default: sol000/phase000
       - id: axesInPlot
@@ -645,9 +747,11 @@ steps:
         source: findRefAnt_join/refant
       - id: prefix
         default: ph_poldif_
+      - id: selfcal_strategy
+        source: selfcal_strategy
       - id: execute
         source: selfcal
-        valueFrom: '$(self ? false : true)'
+        valueFrom: '$(self ? (false || inputs.selfcal_strategy == "HBA") : true)'
     out:
       - id: output_plots
       - id: logfile
@@ -660,4 +764,4 @@ requirements:
   - class: ScatterFeatureRequirement
   - class: MultipleInputFeatureRequirement
   - class: StepInputExpressionRequirement
-  - class: InlineJavascriptRequirement
\ No newline at end of file
+  - class: InlineJavascriptRequirement
diff --git a/workflows/linc_target/selfcal_targ_hba.cwl b/workflows/linc_target/selfcal_targ_hba.cwl
new file mode 100644
index 0000000000000000000000000000000000000000..f3f96b34790ef7c9dfd1443e4e189163be784a2e
--- /dev/null
+++ b/workflows/linc_target/selfcal_targ_hba.cwl
@@ -0,0 +1,484 @@
+class: Workflow
+cwlVersion: v1.2
+id: selfcal_target
+label: selfcal_target
+inputs:
+  - id: max_dp3_threads
+    type: int?
+  - id: msin
+    type: Directory[]
+  - id: propagatesolutions
+    type: boolean?
+    default: true
+  - id: preapply_h5parm
+    type: File[]
+  - id: gsmcal_step
+    type: string?
+    default: 'phaseonly'
+  - id: insolutions
+    type: File
+  - id: selfcal_hba_imsize
+    type: int[]?
+    default: [20000,20000]
+  - id: wsclean_tmpdir
+    type: string?
+    default: '/tmp'
+  - id: process_baselines_target
+    type: string?
+    default: '[CR]S*&'
+  - id: bad_antennas
+    type: string?
+    default: '[CR]S*&'
+  - id: skymodel_source
+    type: string?
+    default: 'TGSS'
+  - id: refant
+    type: string?
+    default: 'CS001HBA0'
+    doc: |
+      Reference antenna for phases. By default None.
+  - id: do_smooth
+    type: boolean?
+  - id: calib_nchan
+    type: int?
+    default: 1
+  - id: calib_uvlambdamin
+    type: float?
+    default: 200.
+outputs:
+  - id: msout
+    outputSource:
+      - self_calibrate_target/msout
+    type: Directory[]
+  - id: outh5parm
+    outputSource:
+      - self_calibrate_target/outh5parm
+    type: File[]
+  - id: image_before
+    outputSource:
+      - image_target/image
+    type: File
+  - id: outsolutions
+    outputSource:
+      - write_solutions/outh5parm
+    type: File
+  - id: inspection
+    outputSource:
+      - losoto_plot_gsm_P/output_plots
+      - losoto_plot_gsm_P2/output_plots
+      - losoto_plot_gsm_Pd/output_plots
+      - losoto_plot_gsm_Pd2/output_plots
+    type: File[]
+    linkMerge: merge_flattened
+  - id: logfiles
+    outputSource:
+      - concat_logfiles_losoto/output
+      - concat_logfiles_applygsm/output
+      - concat_logfiles_blsmooth/output
+      - concat_logfiles_image_target/output
+      - concat_logfiles_self_calibrate/output
+    linkMerge: merge_flattened
+    type: File[]
+steps:
+  - id: h5parm_collector_gsm
+    in:
+      - id: h5parmFiles
+        source: preapply_h5parm
+        linkMerge: merge_flattened
+      - id: squeeze
+        default: true
+      - id: verbose
+        default: true
+      - id: clobber
+        default: true
+    out:
+      - id: outh5parm
+      - id: log
+    run: ../../steps/H5ParmCollector.cwl
+    label: h5parm_collector_gsm
+  - id: add_missing_stations
+    in:
+      - id: h5parm
+        source: h5parm_collector_gsm/outh5parm
+      - id: refh5parm
+        source: insolutions
+      - id: solset
+        default: sol000
+      - id: refsolset
+        default: target
+      - id: soltab_in
+        default: phase000
+      - id: soltab_out
+        source:
+          - skymodel_source
+          - gsmcal_step
+        valueFrom: $(self.join(''))
+      - id: filter
+        source: process_baselines_target
+      - id: bad_antennas
+        source: bad_antennas
+    out:
+      - id: outh5parm
+      - id: log
+    run: ../../steps/add_missing_stations.cwl
+    label: add_missing_stations
+  - id: write_solutions
+    in:
+      - id: h5parmFile
+        source: add_missing_stations/outh5parm
+      - id: outsolset
+        default: target
+      - id: insoltab
+        source:
+          - skymodel_source
+          - gsmcal_step
+        valueFrom: $(self.join(''))
+      - id: input_file
+        source: insolutions
+      - id: squeeze
+        default: true
+      - id: verbose
+        default: true
+      - id: history
+        default: true
+    out:
+      - id: outh5parm
+      - id: log
+    run: ../../steps/h5parmcat.cwl
+    label: write_solutions
+
+  - id: apply_gsmcal
+    in:
+      - id: max_dp3_threads
+        source: max_dp3_threads
+      - id: msin
+        source: msin
+      - id: msin_datacolumn
+        default: DATA
+      - id: parmdb
+        source: write_solutions/outh5parm
+      - id: msout_datacolumn
+        default: CORRECTED_DATA
+      - id: storagemanager
+        default: Dysco
+      - id: databitrate
+        default: 0
+      - id: solset
+        default: target
+      - id: correction
+        source:
+          - skymodel_source
+          - gsmcal_step
+        valueFrom: $(self.join(''))
+    out:
+      - id: msout
+      - id: logfile
+    run: ../../steps/applycal.cwl
+    scatter:
+      - msin
+    label: apply_gsmcal
+
+  - id: image_target
+    in:
+      - id: msin
+        source: apply_gsmcal/msout
+      - id: image_scale
+        default: 1.5asec
+      - id: auto_threshold
+        default: 0.5
+      - id: image_size
+        source: selfcal_hba_imsize
+      - id: niter
+        default: 285000
+      - id: nmiter
+        default: 12
+      - id: multiscale
+        default: true
+      - id: use-wgridder
+        default: true
+      - id: mgain
+        default: 0.8
+      - id: parallel-reordering
+        default: 4
+      - id: parallel-deconvolution
+        default: 2048
+      - id: channels-out
+        default: 6
+      - id: join-channels
+        default: True
+      - id: fit-spectral-pol
+        default: 3
+      - id: weighting
+        default: briggs -1.5
+      - id: minuv-l
+        default: 0
+      - id: tempdir
+        source: wsclean_tmpdir
+      - id: image_name
+        default: image_before_selfcal
+      - id: no_model_update
+        default: false
+      # This can be turned on if we do a separate predict after, not sure if that would be faster.
+      # - id: baseline_averaging
+      #   default: 0.4908738521234052
+      - id: save_source_list
+        default: true
+    out:
+      - id: msout
+      - id: dirty_image
+      - id: image
+      - id: logfile
+      - id: sourcelist
+    run: ../../steps/wsclean.cwl
+    label: image_target
+  # END image target field
+
+  # START selfcal
+  - id: self_calibrate_target
+    in:
+      - id: max_dp3_threads
+        source: max_dp3_threads
+      - id: msin
+        source: image_target/msout
+      - id: do_smooth
+        source: do_smooth
+      - id: propagatesolutions
+        source: propagatesolutions
+      - id: gsmcal_step
+        source: gsmcal_step
+      - id: model_column
+        default: ["MODEL_DATA"]
+      - id: nchan
+        source: calib_nchan
+      - id: smoothnessconstraint
+        default: 5e6
+      - id: smoothnessreffrequency
+        default: 12e7
+      - id: beammode
+        default: array_factor
+      - id: uvlambdamin
+        source: calib_uvlambdamin
+    out:
+      - id: msout
+      - id: BLsmooth.log
+      - id: gaincal.log
+      - id: outh5parm
+    run: ./calib_targ.cwl
+    label: self_calibrate_target
+    scatter:
+      - msin
+  # END selfcal
+
+  # START Plot selfcal solutions
+  - id: losoto_plot_gsm_P
+    in:
+      - id: input_h5parm
+        source:
+          - h5parm_collector_gsm/outh5parm
+      - id: soltab
+        default: sol000/phase000
+      - id: axesInPlot
+        default:
+          - time
+          - freq
+      - id: axisInTable
+        default: ant
+      - id: minmax
+        default:
+          - -3.14
+          - 3.14
+      - id: plotFlag
+        default: true
+      - id: refAnt
+        source: refant
+      - id: prefix
+        source: skymodel_source
+        valueFrom: $("ph_"+self+"_")
+    out:
+      - id: output_plots
+      - id: logfile
+      - id: parset
+    run: ../../steps/LoSoTo.Plot.cwl
+    label: losoto_plot_gsm_P
+
+  - id: losoto_plot_gsm_P2
+    in:
+      - id: input_h5parm
+        source:
+          - h5parm_collector_gsm/outh5parm
+      - id: soltab
+        default: sol000/phase000
+      - id: axesInPlot
+        default:
+          - time
+      - id: axisInTable
+        default: ant
+      - id: axisInCol
+        default: pol
+      - id: minmax
+        default:
+          - -3.14
+          - 3.14
+      - id: plotFlag
+        default: true
+      - id: refAnt
+        source: refant
+      - id: prefix
+        source: skymodel_source
+        valueFrom: $("ph_"+self+"_")
+    out:
+      - id: output_plots
+      - id: logfile
+      - id: parset
+    run: ../../steps/LoSoTo.Plot.cwl
+    label: losoto_plot_gsm_P2
+
+  - id: losoto_plot_gsm_Pd
+    in:
+      - id: input_h5parm
+        source:
+          - h5parm_collector_gsm/outh5parm
+      - id: soltab
+        default: sol000/phase000
+      - id: axesInPlot
+        default:
+          - time
+          - freq
+      - id: axisInTable
+        default: ant
+      - id: axisDiff
+        default: pol
+      - id: minmax
+        default:
+          - -3.14
+          - 3.14
+      - id: plotFlag
+        default: true
+      - id: refAnt
+        source: refant
+      - id: prefix
+        source: skymodel_source
+        valueFrom: $("ph_"+self+"_poldif_")
+    out:
+      - id: output_plots
+      - id: logfile
+      - id: parset
+    run: ../../steps/LoSoTo.Plot.cwl
+    label: losoto_plot_gsm_Pd
+
+  - id: losoto_plot_gsm_Pd2
+    in:
+      - id: input_h5parm
+        source:
+          - h5parm_collector_gsm/outh5parm
+      - id: soltab
+        default: sol000/phase000
+      - id: axesInPlot
+        default:
+          - time
+      - id: axisInTable
+        default: ant
+      - id: axisDiff
+        default: pol
+      - id: minmax
+        default:
+          - -3.14
+          - 3.14
+      - id: plotFlag
+        default: true
+      - id: refAnt
+        source: refant
+      - id: prefix
+        source: skymodel_source
+        valueFrom: $("ph_"+self+"_poldif_")
+    out:
+      - id: output_plots
+      - id: logfile
+      - id: parset
+    run: ../../steps/LoSoTo.Plot.cwl
+    label: losoto_plot_gsm_Pd2
+  # END plot selfcal solutions
+
+  - id: concat_logfiles_losoto
+    in:
+      - id: file_list
+        source:
+          - h5parm_collector_gsm/outh5parm
+          - losoto_plot_gsm_P/logfile
+          - losoto_plot_gsm_P2/logfile
+          - losoto_plot_gsm_Pd/logfile
+          - losoto_plot_gsm_Pd2/logfile
+          - add_missing_stations/log
+          - write_solutions/log
+        linkMerge: merge_flattened
+      - id: file_prefix
+        default: losoto_gsmcal
+    out:
+      - id: output
+    run: ../../steps/concatenate_files.cwl
+    label: concat_logfiles_losoto
+
+  - id: merge_array_files
+    in:
+      - id: input
+        source: apply_gsmcal/logfile
+    out:
+      - id: output
+    run: ../../steps/merge_array_files.cwl
+    label: merge_array_files
+
+  - id: concat_logfiles_applygsm
+    in:
+      - id: file_list
+        source:
+          - merge_array_files/output
+        linkMerge: merge_flattened
+      - id: file_prefix
+        default: apply_gsmcal
+    out:
+      - id: output
+    run: ../../steps/concatenate_files.cwl
+    label: concat_logfiles_applygsm
+
+  - id: concat_logfiles_image_target
+    in:
+      - id: file_list
+        source:
+          - image_target/logfile
+        linkMerge: merge_flattened
+      - id: file_prefix
+        default: image_target
+    out:
+      - id: output
+    run: ../../steps/concatenate_files.cwl
+    label: concat_logfiles_image_target
+
+  - id: concat_logfiles_self_calibrate
+    in:
+      - id: file_list
+        source:
+          - self_calibrate_target/gaincal.log
+        linkMerge: merge_flattened
+      - id: file_prefix
+        default: self_calibrate_target
+    out:
+      - id: output
+    run: ../../steps/concatenate_files.cwl
+    label: concat_logfiles_self_calibrate
+
+  - id: concat_logfiles_blsmooth
+    in:
+      - id: file_list
+        source:
+          - self_calibrate_target/BLsmooth.log
+        pickValue: all_non_null
+      - id: file_prefix
+        default: blsmooth_self_target
+    out:
+      - id: output
+    run: ../../steps/concatenate_files.cwl
+    label: concat_logfiles_blsmooth
+
+requirements:
+  - class: ScatterFeatureRequirement
+  - class: SubworkflowFeatureRequirement
diff --git a/workflows/linc_target/selfcal_targ.cwl b/workflows/linc_target/selfcal_targ_lba.cwl
similarity index 100%
rename from workflows/linc_target/selfcal_targ.cwl
rename to workflows/linc_target/selfcal_targ_lba.cwl
diff --git a/workflows/linc_target/split_into_chunks.cwl b/workflows/linc_target/split_into_chunks.cwl
new file mode 100644
index 0000000000000000000000000000000000000000..750651991ab9d50b8c9c71cef72da131cf08e7cb
--- /dev/null
+++ b/workflows/linc_target/split_into_chunks.cwl
@@ -0,0 +1,54 @@
+class: Workflow
+cwlVersion: v1.2
+id: split_into_chunks
+label: split_into_chunks
+
+inputs:
+  - id: msin
+    type: Directory[]
+  - id: output_channels_per_chunk
+    type: int?
+    default: 20
+  - id: min_unflagged_fraction
+    type: float?
+    default: 0.5
+
+outputs:
+  - id: msouts
+    outputSource:
+      - merge_array/output
+    type: Directory[]
+
+steps:
+  - id: split_into_chunks
+    in:
+      - id: ms
+        source: msin
+        linkMerge: merge_flattened
+      - id: mode
+        valueFrom: 'frequency'
+      - id: chan_per_chunk
+        source: output_channels_per_chunk
+      - id: round_freq
+        default: true
+    out:
+      - id: msouts
+    run: ../../steps/chunk_ms.cwl
+    label: chunk_ms
+    scatter:
+      - ms
+
+  - id: merge_array
+    in:
+      - id: input
+        source:
+          - split_into_chunks/msouts
+    out:
+      - id: output
+    run: ../../steps/merge_array.cwl
+    label: merge_array_files
+
+requirements:
+  - class: ScatterFeatureRequirement
+
+
diff --git a/workflows/linc_target/split_into_chunks.cwl.backup b/workflows/linc_target/split_into_chunks.cwl.backup
new file mode 100644
index 0000000000000000000000000000000000000000..3f0a8134704d427796975a2ed94abfb1b49a8e37
--- /dev/null
+++ b/workflows/linc_target/split_into_chunks.cwl.backup
@@ -0,0 +1,99 @@
+class: Workflow
+cwlVersion: v1.2
+id: split_into_chunks
+label: split_into_chunks
+doc: Splits the given MeasurementSets into multiple frequency chunks.
+
+inputs:
+  - id: msin
+    type: Directory[]
+    doc: The input MeasurementSets to split into chunks.
+  - id: freqchans_per_chunk
+    type: int?
+    default: 20
+    doc: The number of frequency channels per output chunk.
+  - id: min_unflagged_fraction
+    type: float?
+    default: 0.5
+    doc: The minimum unflagged fraction required in an _output_ chunk.
+
+outputs:
+  - id: msouts
+    outputSource:
+      - check_unflagged_fraction/msout
+      #- merge_array_unflagged/output
+    type: Directory[]
+    doc: Output chunks as new MeasurementSets that have been split from the original.
+
+steps:
+  - id: split_into_chunks
+    in:
+      - id: ms
+        source: msin
+        linkMerge: merge_flattened
+      - id: mode
+        valueFrom: 'frequency'
+      - id: chan_per_chunk
+        source: freqchans_per_chunk
+      - id: round_freq
+        default: true
+    out:
+      - id: msouts
+    run: ../../steps/chunk_ms.cwl
+    label: chunk_ms
+    scatter:
+      - ms
+
+  - id: merge_array_split
+    in:
+      - id: input
+        source:
+          - split_into_chunks/msouts
+    out:
+      - id: output
+    run: ../../steps/merge_array.cwl
+    label: merge_array_files
+
+  - id: check_unflagged_fraction
+    in:
+      - id: msin
+        source: split_into_chunks/msouts
+        #source: merge_array_split/output
+        linkMerge: merge_flattened
+      - id: min_fraction
+        source: min_unflagged_fraction
+    out:
+      - id: msout
+      - id: frequency
+      - id: unflagged_fraction
+      - id: filenames
+      - id: logfile
+    run: ../../steps/check_unflagged_fraction.cwl
+    label: check_unflagged_fraction
+    scatter:
+      - msin
+
+  # - id: merge_array_unflagged
+  #   in:
+  #     - id: input
+  #       source:
+  #         - check_unflagged_fraction/msout
+  #   out:
+  #     - id: output
+  #   run: ../../steps/merge_array.cwl
+  #   label: merge_array_files
+
+  # - id: check_filtered_MS_array
+  #   in:
+  #     - id: input
+  #       source: merge_array_unflagged/output
+  #   out:
+  #     - id: output
+  #   run: ../../steps/check_filtered_MS_array.cwl
+  #   label: check_filtered_MS_array
+
+
+requirements:
+  - class: ScatterFeatureRequirement
+
+