Repository: esmvalgroup/esmvalcore
Issue: concatenation of cubes with overlapping data

**Describe the bug**
Thanks for all the efforts so far of concatenating overlapping data (as e.g. #280)! 
For certain datasets however this does not work satisfactorily. E.g. for tas
```
 - {project: CMIP5, dataset: GISS-E2-H-CC, ensemble: r1i1p1, exp: [historical, rcp85], start_year: 1950, end_year: 2100}
```
and uses the files:

```
historical/mon/atmos/Amon/r1i1p1/v20120723/tas/tas_Amon_GISS-E2-H-CC_historical_r1i1p1_195101-201012.nc
historical/mon/atmos/Amon/r1i1p1/v20120723/tas/tas_Amon_GISS-E2-H-CC_historical_r1i1p1_190101-195012.nc
rcp85/mon/atmos/Amon/r1i1p1/v20160512/tas/tas_Amon_GISS-E2-H-CC_rcp85_r1i1p1_200601-205012.nc
rcp85/mon/atmos/Amon/r1i1p1/v20160512/tas/tas_Amon_GISS-E2-H-CC_rcp85_r1i1p1_205101-210012.nc
```

The initial concatenation concatenates the first 2 and the last cube to one (from 1950-2010 and 2051-2100) and leaves the 3rd as it is. _concatenate_overlapping_cubes subsequently ignores this 3rd one.

This is caught later on and throws:
```
esmvalcore.cmor.check.CMORCheckError: There were errors in variable tas:
time: Frequency mon does not match input data
```

Any idea how to approach this issue?

**Please attach**
- [x] The recipe that you are trying to run, you can find a copy in the `run` directory in the output directory  [concatenate_bug_yml.txt](https://github.com/ESMValGroup/ESMValCore/files/4343526/concatenate_bug_yml.txt)

- [x] The `main_log_debug.txt` file, this can also be found in the `run` directory in the output directory
[main_log_debug.txt](https://github.com/ESMValGroup/ESMValCore/files/4343523/main_log_debug.txt)

Changed files:
tests/integration/preprocessor/_io/test_concatenate.py, esmvalcore/preprocessor/_io.py

diff --git a/esmvalcore/preprocessor/_io.py b/esmvalcore/preprocessor/_io.py
index 51fa0c201..203fddb5d 100644
--- a/esmvalcore/preprocessor/_io.py
+++ b/esmvalcore/preprocessor/_io.py
@@ -142,25 +142,21 @@ def _fix_cube_attributes(cubes):
         cube.attributes = attributes
 
 
-def concatenate(cubes):
-    """Concatenate all cubes after fixing metadata."""
-    _fix_cube_attributes(cubes)
-
+def _by_two_concatenation(cubes):
+    """Perform a by-2 concatenation to avoid gaps."""
     concatenated = iris.cube.CubeList(cubes).concatenate()
-    if len(concatenated) == 2:
-        try:
-            concatenated[0].coord('time')
-            concatenated[1].coord('time')
-        except iris.exceptions.CoordinateNotFoundError:
-            pass
+    if len(concatenated) == 1:
+        return concatenated[0]
+    else:
+        concatenated = _concatenate_overlapping_cubes(concatenated)
+        if len(concatenated) == 2:
+            _get_concatenation_error(concatenated)
         else:
-            concatenated = _concatenate_overlapping_cubes(concatenated)
+            return concatenated[0]
 
-    if len(concatenated) == 1:
-        cube = concatenated[0]
-        _fix_aux_factories(cube)
-        return cube
 
+def _get_concatenation_error(cubes):
+    """Raise an error for concatenation."""
     # Concatenation not successful -> retrieve exact error message
     try:
         iris.cube.CubeList(cubes).concatenate_cube()
@@ -168,17 +164,37 @@ def concatenate(cubes):
         msg = str(exc)
     logger.error('Can not concatenate cubes into a single one: %s', msg)
     logger.error('Resulting cubes:')
-    for cube in concatenated:
+    for cube in cubes:
         logger.error(cube)
-        try:
-            time = cube.coord('time')
-        except iris.exceptions.CoordinateNotFoundError:
-            pass
-        else:
-            logger.error('From %s to %s', time.cell(0), time.cell(-1))
+        time = cube.coord("time")
+        logger.error('From %s to %s', time.cell(0), time.cell(-1))
+
     raise ValueError(f'Can not concatenate cubes: {msg}')
 
 
+def concatenate(cubes):
+    """Concatenate all cubes after fixing metadata."""
+    _fix_cube_attributes(cubes)
+
+    if len(cubes) > 1:
+        # order cubes by first time point
+        try:
+            cubes = sorted(cubes, key=lambda c: c.coord("time").cell(0).point)
+        except iris.exceptions.CoordinateNotFoundError as exc:
+            msg = "One or more cubes {} are missing".format(cubes) + \\
+                  " time coordinate: {}".format(str(exc))
+            raise ValueError(msg)
+
+        # iteratively concatenate starting with first cube
+        result = cubes[0]
+        for cube in cubes[1:]:
+            result = _by_two_concatenation([result, cube])
+
+    _fix_aux_factories(result)
+
+    return result
+
+
 def save(cubes, filename, optimize_access='', compress=False, **kwargs):
     """
     Save iris cubes to file.
@@ -364,6 +380,11 @@ def _concatenate_overlapping_cubes(cubes):
     # get time end points
     time_1 = cubes[0].coord('time')
     time_2 = cubes[1].coord('time')
+    if time_1.units != time_2.units:
+        raise ValueError(
+            f"Cubes\\n{cubes[0]}\\nand\\n{cubes[1]}\\ncan not be concatenated: "
+            f"time units {time_1.units}, calendar {time_1.units.calendar} "
+            f"and {time_2.units}, calendar {time_2.units.calendar} differ")
     data_start_1 = time_1.cell(0).point
     data_start_2 = time_2.cell(0).point
     data_end_1 = time_1.cell(-1).point
diff --git a/tests/integration/preprocessor/_io/test_concatenate.py b/tests/integration/preprocessor/_io/test_concatenate.py
index 1fbace4eb..e9d2fbfd1 100644
--- a/tests/integration/preprocessor/_io/test_concatenate.py
+++ b/tests/integration/preprocessor/_io/test_concatenate.py
@@ -253,12 +253,22 @@ class TestConcatenate(unittest.TestCase):
 
     def test_concatenate_with_overlap_2(self):
         """Test a more generic case."""
-        self._add_cube([65., 75.], [3., 200.])
-        self._add_cube([65., 75.], [1000., 7000.])
+        self._add_cube([65., 75., 100.], [9., 10., 11.])
+        self._add_cube([65., 75., 100.], [7., 8., 9.])
         concatenated = _io.concatenate(self.raw_cubes)
         np.testing.assert_array_equal(
             concatenated.coord('time').points,
-            np.array([1., 2., 3., 4., 5., 6., 1000., 7000.]))
+            np.array([1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 11.]))
+
+    def test_concatenate_with_overlap_3(self):
+        """Test a more generic case."""
+        self._add_cube([65., 75., 100.], [9., 10., 11.])
+        self._add_cube([65., 75., 100., 100., 100., 112.],
+                       [7., 8., 9., 10., 11., 12.])
+        concatenated = _io.concatenate(self.raw_cubes)
+        np.testing.assert_array_equal(
+            concatenated.coord('time').points,
+            np.array([1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12.]))
 
     def test_concatenate_with_overlap_same_start(self):
         """Test a more generic case."""
@@ -302,6 +312,26 @@ class TestConcatenate(unittest.TestCase):
         with self.assertRaises(ConcatenateError):
             _io.concatenate(cubes_single_ovlp)
 
+    def test_concatenate_no_time_coords(self):
+        """Test a more generic case."""
+        time_coord_1 = DimCoord([1.5, 5., 7.],
+                                var_name='time',
+                                standard_name='time',
+                                units='days since 1950-01-01')
+        cube1 = Cube([33., 55., 77.],
+                     var_name='sample',
+                     dim_coords_and_dims=((time_coord_1, 0), ))
+        ap_coord_2 = DimCoord([1., 5., 7.],
+                              var_name='air_pressure',
+                              standard_name='air_pressure',
+                              units='m',
+                              attributes={'positive': 'down'})
+        cube2 = Cube([33., 55., 77.],
+                     var_name='sample',
+                     dim_coords_and_dims=((ap_coord_2, 0), ))
+        with self.assertRaises(ValueError):
+            _io.concatenate([cube1, cube2])
+
     def test_concatenate_with_order(self):
         """Test a more generic case."""
         time_coord_1 = DimCoord([1.5, 2., 5., 7.],
@@ -379,6 +409,7 @@ class TestConcatenate(unittest.TestCase):
     def test_fail_metadata_differs(self):
         """Test exception raised if two cubes have different metadata."""
         self.raw_cubes[0].units = 'm'
+        self.raw_cubes[1].units = 'K'
         with self.assertRaises(ValueError):
             _io.concatenate(self.raw_cubes)