From eec647615c21e31b169cb144fcd0d8929a62efc8 Mon Sep 17 00:00:00 2001
From: thopri <thopri@noc.ac.uk>
Date: Wed, 25 Mar 2020 13:25:48 +0000
Subject: [PATCH] added U and V grid and constant current tests

---
 inputs/namelist_remote.bdy                |  2 +-
 pynemo/nemo_bdy_extr_tm3.py               | 12 ++--
 pynemo/nemo_bdy_ncgen.py                  |  9 +--
 pynemo/nemo_bdy_ncpop.py                  |  2 +-
 pynemo/profile.py                         |  4 +-
 unit_tests/namelist_unit_test_offset.bdy  | 16 ++---
 unit_tests/namelist_unit_test_orth.bdy    | 18 +++---
 unit_tests/namelist_unit_test_rotated.bdy | 16 ++---
 unit_tests/src_data_unit_tests.ncml       | 10 ----
 unit_tests/test_gen.py                    | 26 ++++-----
 unit_tests/unit_test.py                   | 71 ++++++++++++++++++-----
 unit_tests/unit_tests.ncml                | 28 +++++++++
 12 files changed, 138 insertions(+), 76 deletions(-)
 delete mode 100644 unit_tests/src_data_unit_tests.ncml
 create mode 100644 unit_tests/unit_tests.ncml

diff --git a/inputs/namelist_remote.bdy b/inputs/namelist_remote.bdy
index 72b58c8..bb2a1a7 100644
--- a/inputs/namelist_remote.bdy
+++ b/inputs/namelist_remote.bdy
@@ -56,7 +56,7 @@
     ln_mask_file   = .false.              !  =T : read mask from file
     cn_mask_file   = 'mask.nc'            !  name of mask file 
                                           !  (if ln_mask_file=.TRUE.)
-    ln_dyn2d       = .false.              !  boundary conditions for 
+    ln_dyn2d       = .true.              !  boundary conditions for
                                           !  barotropic fields
     ln_dyn3d       = .false.              !  boundary conditions for 
                                           !  baroclinic velocities
diff --git a/pynemo/nemo_bdy_extr_tm3.py b/pynemo/nemo_bdy_extr_tm3.py
index 632e035..584ff48 100644
--- a/pynemo/nemo_bdy_extr_tm3.py
+++ b/pynemo/nemo_bdy_extr_tm3.py
@@ -894,12 +894,12 @@ class Extract:
             
 #        for v in self.variables:
         for v in self.var_nam:
-            if self.settings['dyn2d']: # Calculate depth averaged velocity
-                tile_dz = np.tile(self.bdy_dz, [len(self.time_counter), 1, 1, 1])
-                tmp_var = np.reshape(self.d_bdy[v][year]['data'][:,:,:], tile_dz.shape)
-                tmp_var = np.nansum(tmp_var * tile_dz, 2) /np.nansum(tile_dz, 2)
-            else: # Replace NaNs with specified fill value
-                tmp_var = np.where(np.isnan(self.d_bdy[v][year]['data'][:,:,:]),
+            #if self.settings['dyn2d']: # Calculate depth averaged velocity
+            #    tile_dz = np.tile(self.bdy_dz, [len(self.time_counter), 1, 1, 1])
+            #    tmp_var = np.reshape(self.d_bdy[v][year]['data'][:,:,:], tile_dz.shape)
+            #    tmp_var = np.nansum(tmp_var * tile_dz, 2) /np.nansum(tile_dz, 2)
+            #else: # Replace NaNs with specified fill value
+            tmp_var = np.where(np.isnan(self.d_bdy[v][year]['data'][:,:,:]),
                                             self.settings['fv'], 
                                             self.d_bdy[v][year]['data'][:,:,:])
                
diff --git a/pynemo/nemo_bdy_ncgen.py b/pynemo/nemo_bdy_ncgen.py
index 96e1e95..c0ccccc 100644
--- a/pynemo/nemo_bdy_ncgen.py
+++ b/pynemo/nemo_bdy_ncgen.py
@@ -240,7 +240,7 @@ def CreateBDYNetcdfFile(filename, N, I, J, K, rw, h, orig, fv, calendar, grd, va
 
         ncid.close()
 
-    if var_nam[0] == 'thetao' or var_nam[0] == 'so':
+    if var_nam[0] == 'thetao' or var_nam[0] == 'so' or var_nam[0] == 'uo' or var_nam[0] == 'vo':
         logging.info('CMEMS variables identified, using default CMEMS variables.....')
         gridNames = ['T', 'I', 'U', 'V', 'E', 'Z'] # All possible grids
 
@@ -286,6 +286,7 @@ def CreateBDYNetcdfFile(filename, N, I, J, K, rw, h, orig, fv, calendar, grd, va
             vartmpID = ncid.createVariable('thetao', 'f4', ('time_counter', 'z', 'yb', 'xb', ),
                                           fill_value=fv)
             varsalID = ncid.createVariable('so', 'f4', ('time_counter', 'z', 'yb', 'xb', ), fill_value=fv)
+            varsshID = ncid.createVariable('zos', 'f4', ('time_counter', 'z', 'yb', 'xb',), fill_value=fv)
             if grd == 'I':
                 varildID = ncid.createVariable('ileadfra', 'f4', ('time_counter', 'yb', 'xb',),
                                                fill_value=fv)
@@ -297,16 +298,16 @@ def CreateBDYNetcdfFile(filename, N, I, J, K, rw, h, orig, fv, calendar, grd, va
             varztID = ncid.createVariable('depthu', 'f4', ('z', 'yb', 'xb', ), fill_value=fv)
             varbtuID = ncid.createVariable('vobtcrtx', 'f4', ('time_counter', 'yb', 'xb', ),
                                            fill_value=fv)
-            vartouID = ncid.createVariable('vozocrtx', 'f4', ('time_counter', 'z', 'yb', 'xb', ),
+            vartouID = ncid.createVariable('uo', 'f4', ('time_counter', 'z', 'yb', 'xb', ),
                                            fill_value=fv)
         elif grd == 'V':
             varztID = ncid.createVariable('depthv', 'f4', ('z', 'yb', 'xb', ))
             varbtvID = ncid.createVariable('vobtcrty', 'f4', ('time_counter', 'yb', 'xb', ),
                                            fill_value=fv)
-            vartovID = ncid.createVariable('vomecrty', 'f4', ('time_counter', 'z', 'yb', 'xb',),
+            vartovID = ncid.createVariable('vo', 'f4', ('time_counter', 'z', 'yb', 'xb',),
                                            fill_value=fv)
         elif grd == 'Z':
-            varsshID = ncid.createVariable('sossheig', 'f4', ('time_counter', 'yb', 'xb', ),
+            varsshID = ncid.createVariable('zos', 'f4', ('time_counter', 'yb', 'xb', ),
                                            fill_value=fv)
             varmskID = ncid.createVariable('bdy_msk', 'f4', ('y', 'x', ), fill_value=fv)
         else:
diff --git a/pynemo/nemo_bdy_ncpop.py b/pynemo/nemo_bdy_ncpop.py
index 1f55c5d..7e67dc2 100644
--- a/pynemo/nemo_bdy_ncpop.py
+++ b/pynemo/nemo_bdy_ncpop.py
@@ -18,7 +18,7 @@ def write_data_to_file(filename, variable_name, data):
     ncid = Dataset(filename, 'a', clobber=False, format='NETCDF4')
     count = data.shape
 
-    three_dim_variables = ['votemper', 'vosaline', 'N1p', 'N3n', 'N5s','thetao','so']
+    three_dim_variables = ['votemper', 'vosaline', 'N1p', 'N3n', 'N5s','thetao','so','uo','vo']
     two_dim_variables = ['sossheig', 'vobtcrtx', 'vobtcrty', 'iicethic', 'ileadfra', 'isnowthi','zos']
 
     if variable_name in three_dim_variables:
diff --git a/pynemo/profile.py b/pynemo/profile.py
index 9110828..ab394b1 100644
--- a/pynemo/profile.py
+++ b/pynemo/profile.py
@@ -490,8 +490,8 @@ def process_bdy(setup_filepath=0, mask_gui=False):
     
     emap = {}
     grd  = [  't',  'u',  'v']
-    pair = [ None, 'uv', 'uv'] # TODO: devolve this to the namelist?
-    
+    #pair = [ None, 'uv', 'uv'] # TODO: devolve this to the namelist?
+    pair = [None, None, None]
     # TODO: The following is a temporary stop gap to assign variables for both CMEMS downloads
     #  and existing variable names. In future we need a slicker way of determining the variables to extract.
     # Perhaps by scraping the .ncml file - this way biogeochemical tracers
diff --git a/unit_tests/namelist_unit_test_offset.bdy b/unit_tests/namelist_unit_test_offset.bdy
index b11b850..10d164b 100644
--- a/unit_tests/namelist_unit_test_offset.bdy
+++ b/unit_tests/namelist_unit_test_offset.bdy
@@ -30,17 +30,17 @@
 !------------------------------------------------------------------------------
 !  grid information 
 !------------------------------------------------------------------------------
-   sn_src_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_src_hgr_zps.nc'
-   sn_src_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_src_zgr_zps.nc'
-   sn_dst_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_offset_dst_hgr_zps.nc'
-   sn_dst_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_offset_dst_zgr_zps.nc'
-   sn_src_msk = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/mask.nc'
-   sn_bathy   = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_offset_dst_bathy.nc'
+   sn_src_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_src_hgr_zps.nc'
+   sn_src_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_src_zgr_zps.nc'
+   sn_dst_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_offset_dst_hgr_zps.nc'
+   sn_dst_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_offset_dst_zgr_zps.nc'
+   sn_src_msk = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/mask.nc'
+   sn_bathy   = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_offset_dst_bathy.nc'
 
 !------------------------------------------------------------------------------
 !  I/O 
 !------------------------------------------------------------------------------
-   sn_src_dir = '/Users/thopri/Projects/PyNEMO/unit_tests/src_data_unit_tests.ncml' ! src_files/'
+   sn_src_dir = '/Users/thopri/Projects/PyNEMO/unit_tests/unit_tests.ncml' ! src_files/'
    sn_dst_dir = '/Users/thopri/Projects/PyNEMO/unit_tests/test_outputs'
    sn_fn      = 'unit_test_offset'             ! prefix for output files
    nn_fv      = -1e20                 !  set fill value for output files
@@ -61,7 +61,7 @@
     ln_mask_file   = .false.              !  =T : read mask from file
     cn_mask_file   = 'mask.nc'            !  name of mask file 
                                           !  (if ln_mask_file=.TRUE.)
-    ln_dyn2d       = .false.              !  boundary conditions for 
+    ln_dyn2d       = .true.              !  boundary conditions for
                                           !  barotropic fields
     ln_dyn3d       = .false.              !  boundary conditions for 
                                           !  baroclinic velocities
diff --git a/unit_tests/namelist_unit_test_orth.bdy b/unit_tests/namelist_unit_test_orth.bdy
index 35e9d9d..8542652 100644
--- a/unit_tests/namelist_unit_test_orth.bdy
+++ b/unit_tests/namelist_unit_test_orth.bdy
@@ -30,17 +30,17 @@
 !------------------------------------------------------------------------------
 !  grid information 
 !------------------------------------------------------------------------------
-   sn_src_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_src_hgr_zps.nc'
-   sn_src_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_src_zgr_zps.nc'
-   sn_dst_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_dst_hgr_zps.nc'
-   sn_dst_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_dst_zgr_zps.nc'
-   sn_src_msk = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/mask.nc'
-   sn_bathy   = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_dst_bathy.nc'
+   sn_src_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_src_hgr_zps.nc'
+   sn_src_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_src_zgr_zps.nc'
+   sn_dst_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_dst_hgr_zps.nc'
+   sn_dst_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_dst_zgr_zps.nc'
+   sn_src_msk = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/mask.nc'
+   sn_bathy   = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_dst_bathy.nc'
 
 !------------------------------------------------------------------------------
 !  I/O 
 !------------------------------------------------------------------------------
-   sn_src_dir = '/Users/thopri/Projects/PyNEMO/unit_tests/src_data_unit_tests.ncml' ! src_files/'
+   sn_src_dir = '/Users/thopri/Projects/PyNEMO/unit_tests/unit_tests.ncml' ! src_files/'
    sn_dst_dir = '/Users/thopri/Projects/PyNEMO/unit_tests/test_outputs'
    sn_fn      = 'unit_test_orth'             ! prefix for output files
    nn_fv      = -1e20                 !  set fill value for output files
@@ -61,9 +61,9 @@
     ln_mask_file   = .false.              !  =T : read mask from file
     cn_mask_file   = 'mask.nc'            !  name of mask file 
                                           !  (if ln_mask_file=.TRUE.)
-    ln_dyn2d       = .false.              !  boundary conditions for
+    ln_dyn2d       = .true.              !  boundary conditions for
                                           !  barotropic fields
-    ln_dyn3d       = .false.              !  boundary conditions for 
+    ln_dyn3d       = .false.              !  boundary conditions for
                                           !  baroclinic velocities
     ln_tra         = .true.               !  boundary conditions for T and S
     ln_ice         = .false.              !  ice boundary condition   
diff --git a/unit_tests/namelist_unit_test_rotated.bdy b/unit_tests/namelist_unit_test_rotated.bdy
index 063affb..f970101 100644
--- a/unit_tests/namelist_unit_test_rotated.bdy
+++ b/unit_tests/namelist_unit_test_rotated.bdy
@@ -30,17 +30,17 @@
 !------------------------------------------------------------------------------
 !  grid information 
 !------------------------------------------------------------------------------
-   sn_src_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_src_hgr_zps.nc'
-   sn_src_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_src_zgr_zps.nc'
-   sn_dst_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_rot_dst_hgr_zps.nc'
-   sn_dst_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_rot_dst_zgr_zps.nc'
-   sn_src_msk = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/mask.nc'
-   sn_bathy   = '/Users/thopri/Projects/PyNEMO/unit_tests/test_data/test_rot_dst_bathy.nc'
+   sn_src_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_src_hgr_zps.nc'
+   sn_src_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_src_zgr_zps.nc'
+   sn_dst_hgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_rot_dst_hgr_zps.nc'
+   sn_dst_zgr = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_rot_dst_zgr_zps.nc'
+   sn_src_msk = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/mask.nc'
+   sn_bathy   = '/Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/test_rot_dst_bathy.nc'
 
 !------------------------------------------------------------------------------
 !  I/O 
 !------------------------------------------------------------------------------
-   sn_src_dir = '/Users/thopri/Projects/PyNEMO/unit_tests/src_data_unit_tests.ncml' ! src_files/'
+   sn_src_dir = '/Users/thopri/Projects/PyNEMO/unit_tests/unit_tests.ncml' ! src_files/'
    sn_dst_dir = '/Users/thopri/Projects/PyNEMO/unit_tests/test_outputs'
    sn_fn      = 'unit_test_rotated'             ! prefix for output files
    nn_fv      = -1e20                 !  set fill value for output files
@@ -61,7 +61,7 @@
     ln_mask_file   = .false.              !  =T : read mask from file
     cn_mask_file   = 'mask.nc'            !  name of mask file 
                                           !  (if ln_mask_file=.TRUE.)
-    ln_dyn2d       = .false.              !  boundary conditions for 
+    ln_dyn2d       = .true.              !  boundary conditions for
                                           !  barotropic fields
     ln_dyn3d       = .false.              !  boundary conditions for 
                                           !  baroclinic velocities
diff --git a/unit_tests/src_data_unit_tests.ncml b/unit_tests/src_data_unit_tests.ncml
deleted file mode 100644
index 3016c52..0000000
--- a/unit_tests/src_data_unit_tests.ncml
+++ /dev/null
@@ -1,10 +0,0 @@
-<ns0:netcdf xmlns:ns0="http://www.unidata.ucar.edu/namespaces/netcdf/ncml-2.2" title="NEMO aggregation">
-  <ns0:aggregation type="union">
-    <ns0:netcdf>
-      <ns0:aggregation dimName="time" name="temperature" type="joinExisting">
-        <ns0:scan location="file://Users/thopri/Projects/PyNEMO/unit_tests/test_data/" regExp=".*T\.nc$" />
-      </ns0:aggregation>
-    </ns0:netcdf>
-  </ns0:aggregation>
-  <ns0:variable name="thetao" orgName="thetao" />
-</ns0:netcdf>
diff --git a/unit_tests/test_gen.py b/unit_tests/test_gen.py
index 0d6f982..9c27327 100644
--- a/unit_tests/test_gen.py
+++ b/unit_tests/test_gen.py
@@ -21,8 +21,8 @@ def _main():
     max_dep = 100
     min_dep = 10
     z_end_dim = 1
-    h_fname = 'unit_tests/test_data/test_src_hgr_zps.nc'
-    z_fname = 'unit_tests/test_data/test_src_zgr_zps.nc'
+    h_fname = 'unit_tests/test_inputs/test_src_hgr_zps.nc'
+    z_fname = 'unit_tests/test_inputs/test_src_zgr_zps.nc'
     grid_h1 = gt.set_hgrid(dx,dy,jpi,jpj)
     grid_z1 = gt.set_zgrid(grid_h1,jpk,max_dep,min_dep,z_end_dim)
     write_coord_H = gt.write_coord_H(h_fname,grid_h1)
@@ -42,14 +42,14 @@ def _main():
     min_dep = 10
     z_end_dim = 1
     sf = 10
-    h_fname = 'unit_tests/test_data/test_dst_hgr_zps.nc'
-    z_fname = 'unit_tests/test_data/test_dst_zgr_zps.nc'
+    h_fname = 'unit_tests/test_inputs/test_dst_hgr_zps.nc'
+    z_fname = 'unit_tests/test_inputs/test_dst_zgr_zps.nc'
     grid_h2 = gt.set_hgrid(dx,dy,jpi,jpj,zoffx,zoffy,sf)
     grid_z2 = gt.set_zgrid(grid_h2,jpk,max_dep,min_dep,z_end_dim)
     write_coord_H = gt.write_coord_H(h_fname,grid_h2)
     write_coord_Z = gt.write_coord_Z(z_fname,grid_h2,grid_z2)
     # write bathy files (constant bathy)
-    bathy_fname = 'unit_tests/test_data/test_dst_bathy.nc'
+    bathy_fname = 'unit_tests/test_inputs/test_dst_bathy.nc'
     bathy = gt.write_bathy(bathy_fname,grid_h2,grid_z2)
     if write_coord_H + write_coord_Z + bathy == 0:
         print("Org child grid generation successful!")
@@ -60,8 +60,8 @@ def _main():
     origin = (8,8)
 
     # rotate grid
-    rot_h_fname = 'unit_tests/test_data/test_rot_dst_hgr_zps.nc'
-    rot_z_fname = 'unit_tests/test_data/test_rot_dst_zgr_zps.nc'
+    rot_h_fname = 'unit_tests/test_inputs/test_rot_dst_hgr_zps.nc'
+    rot_z_fname = 'unit_tests/test_inputs/test_rot_dst_zgr_zps.nc'
     grid_rot = grid_h2.copy()
     grid_rot['latt'], grid_rot['lont'] = gt.rotate_around_point(grid_h2['latt'],grid_h2['lont'],theta,origin)
     grid_rot['latu'], grid_rot['lonu'] = gt.rotate_around_point(grid_h2['latu'], grid_h2['lonu'], theta, origin)
@@ -70,7 +70,7 @@ def _main():
     write_coord_H = gt.write_coord_H(rot_h_fname,grid_rot)
     write_coord_Z = gt.write_coord_Z(rot_z_fname,grid_rot,grid_z2)
     # write bathy files (constant bathy)
-    bathy_fname = 'unit_tests/test_data/test_rot_dst_bathy.nc'
+    bathy_fname = 'unit_tests/test_inputs/test_rot_dst_bathy.nc'
     bathy = gt.write_bathy(bathy_fname,grid_rot,grid_z2)
     if write_coord_H + write_coord_Z + bathy == 0:
         print("Rotated child grid generation Successful!")
@@ -87,14 +87,14 @@ def _main():
     min_dep = 10
     z_end_dim = 1
     sf = 10
-    h_fname = 'unit_tests/test_data/test_offset_dst_hgr_zps.nc'
-    z_fname = 'unit_tests/test_data/test_offset_dst_zgr_zps.nc'
+    h_fname = 'unit_tests/test_inputs/test_offset_dst_hgr_zps.nc'
+    z_fname = 'unit_tests/test_inputs/test_offset_dst_zgr_zps.nc'
     grid_h3 = gt.set_hgrid(dx,dy,jpi,jpj,zoffx,zoffy,sf)
     grid_z3 = gt.set_zgrid(grid_h2,jpk,max_dep,min_dep,z_end_dim)
     write_coord_H = gt.write_coord_H(h_fname,grid_h3)
     write_coord_Z = gt.write_coord_Z(z_fname,grid_h3,grid_z3)
     # write bathy files (constant bathy)
-    bathy_fname = 'unit_tests/test_data/test_offset_dst_bathy.nc'
+    bathy_fname = 'unit_tests/test_inputs/test_offset_dst_bathy.nc'
     bathy = gt.write_bathy(bathy_fname,grid_h3,grid_z3)
     if write_coord_H + write_coord_Z + bathy == 0:
         print("Offset child grid gneration successful!")
@@ -104,7 +104,7 @@ def _main():
     #           grid_h3['lont'],grid_h1['latt'],grid_h1['lont'])
 
     # write boundary files (constant parameters)
-    out_fname = 'unit_tests/test_data/output_boundary' #drop file extension
+    out_fname = 'unit_tests/test_inputs/output_boundary' #drop file extension
     params_t = {'param1': {'name':'thetao','const_value':15.0,'longname':'temperature','units':'degreesC'},
               'param2': {'name':'so','const_value':35.0,'longname':'salinity','units':'PSU'},
               'param3': {'name': 'zos', 'const_value': 1.0, 'longname': 'sea surface height', 'units': 'metres'}
@@ -121,7 +121,7 @@ def _main():
         print('Boundary file generation successful!')
 
     #write_mask
-    mask_fname = 'unit_tests/test_data/mask.nc'
+    mask_fname = 'unit_tests/test_inputs/mask.nc'
     mask = gt.write_mask(mask_fname,grid_h1,grid_z1)
     if mask == 0:
         print('Mask file generation successful!')
diff --git a/unit_tests/unit_test.py b/unit_tests/unit_test.py
index 2f5a5c8..2f041ad 100644
--- a/unit_tests/unit_test.py
+++ b/unit_tests/unit_test.py
@@ -34,7 +34,9 @@ for n in namelist_files:
 
 # perform tests
 def test_temp():
-    test_files = glob.glob('unit_tests/test_outputs/unit_test*')
+    test_files = glob.glob('unit_tests/test_outputs/*bdyT*')
+    if len(test_files) == 0:
+        raise Exception('DONT PANIC: no temperature test files found')
     for t in test_files:
         results = Dataset(t) # open results
         temp = results['thetao'][:]
@@ -45,7 +47,9 @@ def test_temp():
         assert abs(temp_[temp_ != 0.0].min() - 15) <= 0.001
 
 def test_salinty():
-    test_files = glob.glob('unit_tests/test_outputs/unit_test*')
+    test_files = glob.glob('unit_tests/test_outputs/*bdyT*')
+    if len(test_files) == 0:
+        raise Exception('DONT PANIC: no salinity test files found')
     for t in test_files:
         results = Dataset(t)  # open results
         sal = results['so'][:]
@@ -55,18 +59,57 @@ def test_salinty():
         assert abs(sal_[sal_ != 0.0].max() - 35) <= 0.001
         assert abs(sal_[sal_ != 0.0].min() - 35) <= 0.001
 
+def test_ssh():
+    test_files = glob.glob('unit_tests/test_outputs/*bdyT*')
+    if len(test_files) == 0:
+        raise Exception('DONT PANIC: no SSH test files found')
+    for t in test_files:
+        results = Dataset(t)  # open results
+        ssh = results['zos'][:]
+        results.close()
+        ssh_ = np.ma.masked_array(ssh,ssh == -32767.0)
+        assert abs(ssh_[ssh_!=0.0].mean() - 1.0) <= 0.001
+        assert abs(ssh_[ssh_ != 0.0].max() - 1.0) <= 0.001
+        assert abs(ssh_[ssh_ != 0.0].min() - 1.0) <= 0.001
+
+def test_U():
+    test_files = glob.glob('unit_tests/test_outputs/*bdyU*')
+    if len(test_files) == 0:
+        raise Exception('DONT PANIC: no U current test files found')
+    for t in test_files:
+        results = Dataset(t)  # open results
+        U = results['uo'][:]
+        results.close()
+        U_ = np.ma.masked_array(U,U == -32767.0)
+        assert abs(U_[U_!=0.0].mean() - 0.5) <= 0.001
+        assert abs(U_[U_ != 0.0].max() - 0.5) <= 0.001
+        assert abs(U_[U_ != 0.0].min() - 0.5) <= 0.001
+
+def test_V():
+    test_files = glob.glob('unit_tests/test_outputs/*bdyV*')
+    if len(test_files) == 0:
+        raise Exception('DONT PANIC: no V current test files found')
+    for t in test_files:
+        results = Dataset(t)  # open results
+        V = results['vo'][:]
+        results.close()
+        V_ = np.ma.masked_array(V,V == -32767.0)
+        assert abs(V_[V_!=0.0].mean() - 0.5) <= 0.001
+        assert abs(V_[V_ != 0.0].max() - 0.5) <= 0.001
+        assert abs(V_[V_ != 0.0].min() - 0.5) <= 0.001
+
 # clean up test I/O
-files = glob.glob('unit_tests/test_outputs/*')
-for f in files:
-    os.remove(f)
-files = glob.glob('unit_tests/test_outputs/*')
-if len(files) != 0:
-    raise Exception('DONT PANIC: output folder not cleaned')
+def test_rm_out():
+    files = glob.glob('unit_tests/test_outputs/*')
+    for f in files:
+        os.remove(f)
+    files = glob.glob('unit_tests/test_outputs/*')
+    assert len(files) == 0
 
-files = glob.glob('unit_tests/test_data/*')
-for f in files:
-    os.remove(f)
-files = glob.glob('unit_tests/test_data/*')
-if len(files) != 0:
-    raise Exception('DONT PANIC: input folder not cleaned')
 
+def test_rm_in():
+    files = glob.glob('unit_tests/test_inputs/*')
+    for f in files:
+        os.remove(f)
+    files = glob.glob('unit_tests/test_inputs/*')
+    assert len(files) == 0
diff --git a/unit_tests/unit_tests.ncml b/unit_tests/unit_tests.ncml
new file mode 100644
index 0000000..d376a10
--- /dev/null
+++ b/unit_tests/unit_tests.ncml
@@ -0,0 +1,28 @@
+<ns0:netcdf xmlns:ns0="http://www.unidata.ucar.edu/namespaces/netcdf/ncml-2.2" title="NEMO aggregation">
+  <ns0:aggregation type="union">
+    <ns0:netcdf>
+      <ns0:aggregation dimName="time" name="temperature" type="joinExisting">
+        <ns0:scan location="file://Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/" regExp=".*T\.nc$" />
+      </ns0:aggregation>
+    </ns0:netcdf>
+    <ns0:netcdf>
+      <ns0:aggregation dimName="time" name="zonal_velocity" type="joinExisting">
+        <ns0:scan location="file://Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/" regExp=".*U\.nc$" />
+      </ns0:aggregation>
+    </ns0:netcdf>
+    <ns0:netcdf>
+      <ns0:aggregation dimName="time" name="meridian_velocity" type="joinExisting">
+        <ns0:scan location="file://Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/" regExp=".*V\.nc$" />
+      </ns0:aggregation>
+    </ns0:netcdf>
+    <ns0:netcdf>
+      <ns0:aggregation dimName="time" name="sea_surface_height" type="joinExisting">
+        <ns0:scan location="file://Users/thopri/Projects/PyNEMO/unit_tests/test_inputs/" regExp=".*T\.nc$" />
+      </ns0:aggregation>
+    </ns0:netcdf>
+  </ns0:aggregation>
+  <ns0:variable name="thetao" orgName="thetao" />
+  <ns0:variable name="uo" orgName="uo" />
+  <ns0:variable name="vo" orgName="vo" />
+  <ns0:variable name="zos" orgName="zos" />
+</ns0:netcdf>
-- 
GitLab