add metadata to output netcdf file

azspiceEnv

@@ -6,6 +6,7 @@
 module load csc/2025.3.20 gcc/12.2.0 mpich/4.2.3
 module load netcdf-c/4.9.2 netcdf-fortran/4.6.1
 module load xios/2701
+module load nccmp
 
 export FCFLAGS="-g $FFLAGS"
 export LDFLAGS="$LDFLAGS -lxios -lnetcdf -lnetcdff -lstdc++"

xios_examples/gen_netcdf.py

@@ -136,24 +136,24 @@ def create_ncfile_unstructured(ncmeshout, meshin_file, meshin_varname, nlev, fun
 
     for face_coord in meshin_var.face_coordinates.split(" "):
         face_coordvar = ncmeshin.variables[face_coord]
-        if face_coordvar.standard_name == 'longitude':
+        if face_coordvar.standard_name in ['longitude', 'projection_x_coordinate']:
             face_lon = face_coordvar[:]
-        elif face_coordvar.standard_name == 'latitude':
+        elif face_coordvar.standard_name in ['latitude', 'projection_y_coordinate']:
             face_lat = face_coordvar[:]
 
     for node_coord in meshin_var.node_coordinates.split(" "):
         node_coordvar = ncmeshin.variables[node_coord]
-        if node_coordvar.standard_name == 'longitude':
+        if node_coordvar.standard_name in ['longitude', 'projection_x_coordinate']:
             node_lon = node_coordvar[:]
-        elif node_coordvar.standard_name == 'latitude':
+        elif node_coordvar.standard_name in ['latitude', 'projection_y_coordinate']:
             node_lat = node_coordvar[:]
 
     if 'edge_coordinates' in meshin_var.ncattrs():
         for edge_coord in meshin_var.edge_coordinates.split(" "):
             edge_coordvar = ncmeshin.variables[edge_coord]
-            if edge_coordvar.standard_name == 'longitude':
+            if edge_coordvar.standard_name in ['longitude', 'projection_x_coordinate']:
                 edge_lon = edge_coordvar[:]
-            elif edge_coordvar.standard_name == 'latitude':
+            elif edge_coordvar.standard_name in ['latitude', 'projection_y_coordinate']:
                 edge_lat = edge_coordvar[:]
 
     meshout_varname = 'Mesh2d'

xios_examples/setAttrCoord/Makefile

@@ -0,0 +1,34 @@
+# Make file for the write demonstartion XIOS programme
+# Targets provided our detailed below...
+#
+# all: (default) Build the write programme
+# clean: Delete all final products and working files
+# run: run the programme
+#
+# Environment Variables expected by this MakeFile:
+#
+# FC: mpif90
+# FCFLAGS: -g & include files for netcdf & xios
+# LD_FLAGS: for xios, netcdf, netcdff, stfc++
+# LD_LIBRARY_PATH: for netCDF & XIOS libs
+# XIOS_BINDIR: The directory for XIOS binary files
+
+.PHONY: all, clean, run
+
+all: write
+
+# fortran compilation
+%.o: %.F90
+	$(FC) $(FCFLAGS) -c $<
+
+# fortran linking
+write: write.o
+	$(FC) -o write.exe write.o $(LDFLAGS) \
+		&& ln -fs $(XIOS_BINDIR)/xios_server.exe .
+
+run:
+	mpiexec -n 1 ./write.exe : -n 1 ./xios_server.exe
+
+# cleanup
+clean:
+	rm -f *.exe *.o *.mod *.MOD *.out *.err

xios_examples/setAttrCoord/Readme.md

@@ -0,0 +1,3 @@
+Demonstration of user customised metadata from XIOS
+
+This includes spatial and temproal metadata definitions and global attributes.

xios_examples/setAttrCoord/axis_check.xml

@@ -0,0 +1,20 @@
+<context>
+
+  <calendar type="Gregorian"/>
+
+  <grid_definition>
+    <grid id="oax_grid">
+      <axis id="x" />
+      <axis id="y" />
+      <axis id="alt" />
+    </grid>
+
+  </grid_definition>
+
+  <file_definition type="one_file">
+    <file id="din" name="data_input" output_freq="1ts" mode="read" enabled=".true.">
+      <field id="odatax" name="original_data" grid_ref="oax_grid" operation="instant" read_access=".true." />
+    </file>
+  </file_definition>
+
+</context>

xios_examples/setAttrCoord/data_input.cdl

@@ -0,0 +1,37 @@
+netcdf domain_input {
+dimensions:
+        x = 5 ;
+        y = 5 ;
+        alt = 1 ;
+variables:
+        float x(x) ;
+                x:standard_name = "projection_x_coordinate" ;
+                x:units = "m";
+        float y(y) ;
+                y:standard_name = "projection_y_coordinate" ;
+                y:units = "m";
+        float alt(alt) ;
+                alt:standard_name = "altitude" ;
+                alt:units = "metres";                
+        double original_data(alt, y, x) ;
+                original_data:long_name = "data values" ;
+                original_data:units = "1";
+
+// global attributes:
+                :title = "Input data for XIOS output." ;
+
+data:
+
+ x = 303000, 305000, 307000, 309000, 311000 ;
+
+ y = 107000, 109000, 111000, 113000, 115000 ;
+
+ alt = 50 ;
+
+ original_data =  0,  2,  4,  6,  8,
+                  2,  4,  6,  8, 10,
+                  4,  6,  8, 10, 12,
+                  6,  8, 10, 12, 14,
+                  8, 10, 12, 14, 16 ;
+
+}

xios_examples/setAttrCoord/expected_domain_output.cdl

@@ -0,0 +1,110 @@
+netcdf data_output {
+dimensions:
+	axis_nbounds = 2 ;
+	x = 5 ;
+	y = 5 ;
+	alt = 1 ;
+	time_counter = UNLIMITED ; // (2 currently)
+variables:
+	float x(x) ;
+		x:name = "x" ;
+		x:standard_name = "projection_x_coordinate" ;
+		x:long_name = "x coordinate of projection" ;
+		x:units = "m" ;
+	float y(y) ;
+		y:name = "y" ;
+		y:standard_name = "projection_y_coordinate" ;
+		y:long_name = "y coordinate of projection" ;
+		y:units = "m" ;
+	float alt(alt) ;
+		alt:name = "alt" ;
+		alt:standard_name = "altitude" ;
+		alt:units = "metres" ;
+	float forecast_reference_time ;
+		forecast_reference_time:name = "forecast_reference_time" ;
+		forecast_reference_time:standard_name = "forecast_reference_time" ;
+		forecast_reference_time:units = "seconds since 2022-02-02 12:00:00" ;
+	double time_instant(time_counter) ;
+		time_instant:standard_name = "time" ;
+		time_instant:long_name = "Time axis" ;
+		time_instant:calendar = "gregorian" ;
+		time_instant:units = "seconds since 2022-02-02 12:00:00" ;
+		time_instant:time_origin = "2022-02-02 12:00:00" ;
+		time_instant:bounds = "time_instant_bounds" ;
+	double time_instant_bounds(time_counter, axis_nbounds) ;
+	double time_counter(time_counter) ;
+		time_counter:axis = "T" ;
+		time_counter:standard_name = "time" ;
+		time_counter:long_name = "Time axis" ;
+		time_counter:calendar = "gregorian" ;
+		time_counter:units = "seconds since 2022-02-02 12:00:00" ;
+		time_counter:time_origin = "2022-02-02 12:00:00" ;
+		time_counter:bounds = "time_counter_bounds" ;
+	double time_counter_bounds(time_counter, axis_nbounds) ;
+	double original_data(time_counter, alt, y, x) ;
+		original_data:long_name = "Arbitrary data values" ;
+		original_data:units = "1" ;
+		original_data:online_operation = "instant" ;
+		original_data:interval_operation = "1 h" ;
+		original_data:interval_write = "1 h" ;
+		original_data:cell_methods = "time: point" ;
+		original_data:coordinates = "time_instant forecast_reference_time" ;
+		original_data:grid_mapping = "osgb: x y egm2008:alt" ;
+	short osgb ;
+		osgb:online_operation = "once" ;
+		osgb:_FillValue = -32767s ;
+		osgb:missing_value = -32767s ;
+		osgb:coordinates = "" ;
+		osgb:crs_wkt = "PROJCRS[\"OSGB36 / British National Grid\",BASEGEOGCRS[\"OSGB36\",DATUM[\"Ordnance Survey of Great Britain 1936\",ELLIPSOID[\"Airy 1830\",6377563.396,299.3249646,LENGTHUNIT[\"metre\",1,ID[\"EPSG\",9001]],ID[\"EPSG\",7001]],ID[\"EPSG\",6277]],ID[\"EPSG\",4277]],CONVERSION[\"British National Grid\",METHOD[\"Transverse Mercator\",ID[\"EPSG\",9807]],PARAMETER[\"Latitude of natural origin\",49,ANGLEUNIT[\"degree\",0.0174532925199433,ID[\"EPSG\",9102]],ID[\"EPSG\",8801]],PARAMETER[\"Longitude of natural origin\",-2,ANGLEUNIT[\"degree\",0.0174532925199433,ID[\"EPSG\",9102]],ID[\"EPSG\",8802]],PARAMETER[\"Scale factor at natural origin\",0.9996012717,SCALEUNIT[\"unity\",1,ID[\"EPSG\",9201]],ID[\"EPSG\",8805]],PARAMETER[\"False easting\",400000,LENGTHUNIT[\"metre\",1,ID[\"EPSG\",9001]],ID[\"EPSG\",8806]],PARAMETER[\"False northing\",-100000,LENGTHUNIT[\"metre\",1,ID[\"EPSG\",9001]],ID[\"EPSG\",8807]],ID[\"EPSG\",19916]],CS[Cartesian,2,ID[\"EPSG\",4400]],AXIS[\"Easting (E)\",east],AXIS[\"Northing (N)\",north],LENGTHUNIT[\"metre\",1,ID[\"EPSG\",9001]],ID[\"EPSG\",27700]]" ;
+	short egm2008 ;
+		egm2008:online_operation = "once" ;
+		egm2008:_FillValue = -32767s ;
+		egm2008:missing_value = -32767s ;
+		egm2008:coordinates = "" ;
+		egm2008:crs_wkt = "VERTCRS[\"EGM2008 height\",VDATUM[\"EGM2008 geoid\",ID[\"EPSG\",1027]],CS[vertical,1,ID[\"EPSG\",6499]],AXIS[\"Gravity-related height (H)\",up],LENGTHUNIT[\"metre\",1,ID[\"EPSG\",9001]],GEOIDMODEL[\"WGS 84 to EGM2008 height (1)\",ID[\"EPSG\",3858]],GEOIDMODEL[\"WGS 84 to EGM2008 height (2)\",ID[\"EPSG\",3859]],ID[\"EPSG\",3855]]" ;
+
+// global attributes:
+		:description = "LFRic file format v0.2.0" ;
+		:Conventions = "CF-1.6, UGRID" ;
+		:timeStamp = "2025-Jun-19 08:40:45 GMT" ;
+		:uuid = "4166ab59-b790-423d-9fa3-7693aa16a04a" ;
+		:name = "Attribute demonstration" ;
+		:title = "Attribute demonstration" ;
+data:
+
+ x = 303000, 305000, 307000, 309000, 311000 ;
+
+ y = 107000, 109000, 111000, 113000, 115000 ;
+
+ alt = 50 ;
+
+ forecast_reference_time = 2.71296e+07 ;
+
+ time_instant = 27133200, 27136800 ;
+
+ time_instant_bounds =
+  27133200, 27133200,
+  27136800, 27136800 ;
+
+ time_counter = 27133200, 27136800 ;
+
+ time_counter_bounds =
+  27133200, 27133200,
+  27136800, 27136800 ;
+
+ original_data =
+  0, 2, 4, 6, 8,
+  2, 4, 6, 8, 10,
+  4, 6, 8, 10, 12,
+  6, 8, 10, 12, 14,
+  8, 10, 12, 14, 16,
+  0, 2, 4, 6, 8,
+  2, 4, 6, 8, 10,
+  4, 6, 8, 10, 12,
+  6, 8, 10, 12, 14,
+  8, 10, 12, 14, 16 ;
+
+ osgb = _ ;
+
+ egm2008 = _ ;
+}

xios_examples/setAttrCoord/iodef.xml

@@ -0,0 +1,9 @@
+<?xml version="1.0"?>
+<simulation>
+
+<context id="main" src="main.xml"/>
+<context id="axis_check" src="axis_check.xml"/>
+
+<context id="xios" src="xios.xml"/>
+
+</simulation>

xios_examples/setAttrCoord/main.xml

@@ -0,0 +1,65 @@
+<context>
+
+  <calendar type="Gregorian"/>
+
+  <scalar_definition>
+    <scalar id="metadata_scalar"/>
+    <scalar id="frt" name="forecast_reference_time" standard_name="forecast_reference_time" />
+  </scalar_definition>
+
+  <axis_definition>
+    <axis id="x" />
+    <axis id="y" />
+    <axis id="alt" unit="metres" standard_name="altitude" />
+  </axis_definition>
+
+
+  <grid_definition>
+
+    <grid id="original_grid">
+      <axis axis_ref="x" />
+      <axis axis_ref="y" />
+      <axis axis_ref="alt" />
+      <scalar scalar_ref="frt" />
+    </grid>
+
+  </grid_definition>
+
+  <field_definition prec="8">
+    <field id="odata" name="original_data" grid_ref="original_grid" long_name="Arbitrary data values" unit="1"  />
+  </field_definition>
+
+
+  <file_definition type="one_file">
+    <file id="data_input" output_freq="1ts" mode="read" enabled=".true.">
+      <field id="odatain" name="original_data" grid_ref="original_grid" operation="instant" read_access=".true." />
+
+    </file>
+    <file id="data_output" output_freq="1ts">
+      <field_group operation="instant">
+        <field field_ref="odata" >
+	  <variable name="grid_mapping" type="string">osgb: x y egm2008:alt</variable>
+        </field>
+
+      </field_group>
+
+      <field_group operation="once" scalar_ref="metadata_scalar" >
+        <field id="hcrs" name="osgb" prec="2" default_value="-32767">
+          <variable name="crs_wkt" type="string">PROJCRS["OSGB36 / British National Grid",BASEGEOGCRS["OSGB36",DATUM["Ordnance Survey of Great Britain 1936",ELLIPSOID["Airy 1830",6377563.396,299.3249646,LENGTHUNIT["metre",1,ID["EPSG",9001]],ID["EPSG",7001]],ID["EPSG",6277]],ID["EPSG",4277]],CONVERSION["British National Grid",METHOD["Transverse Mercator",ID["EPSG",9807]],PARAMETER["Latitude of natural origin",49,ANGLEUNIT["degree",0.0174532925199433,ID["EPSG",9102]],ID["EPSG",8801]],PARAMETER["Longitude of natural origin",-2,ANGLEUNIT["degree",0.0174532925199433,ID["EPSG",9102]],ID["EPSG",8802]],PARAMETER["Scale factor at natural origin",0.9996012717,SCALEUNIT["unity",1,ID["EPSG",9201]],ID["EPSG",8805]],PARAMETER["False easting",400000,LENGTHUNIT["metre",1,ID["EPSG",9001]],ID["EPSG",8806]],PARAMETER["False northing",-100000,LENGTHUNIT["metre",1,ID["EPSG",9001]],ID["EPSG",8807]],ID["EPSG",19916]],CS[Cartesian,2,ID["EPSG",4400]],AXIS["Easting (E)",east],AXIS["Northing (N)",north],LENGTHUNIT["metre",1,ID["EPSG",9001]],ID["EPSG",27700]]</variable>
+      </field>
+        <field id="vcrs" name="egm2008" prec="2" default_value="-32767">
+          <variable name="crs_wkt" type="string">VERTCRS["EGM2008 height",VDATUM["EGM2008 geoid",ID["EPSG",1027]],CS[vertical,1,ID["EPSG",6499]],AXIS["Gravity-related height (H)",up],LENGTHUNIT["metre",1,ID["EPSG",9001]],GEOIDMODEL["WGS 84 to EGM2008 height (1)",ID["EPSG",3858]],GEOIDMODEL["WGS 84 to EGM2008 height (2)",ID["EPSG",3859]],ID["EPSG",3855]]</variable>
+      </field>
+
+      </field_group>
+
+
+      <variable name="name" type="string">Attribute demonstration</variable>
+
+      <variable name="title" type="string">Attribute demonstration</variable>
+
+
+    </file>
+  </file_definition>
+
+</context>

xios_examples/setAttrCoord/test_write_metadata.py

@@ -0,0 +1,70 @@
+import copy
+import glob
+import netCDF4
+import numpy as np
+import os
+import subprocess
+import unittest
+
+import xios_examples.shared_testing as xshared
+
+this_path = os.path.realpath(__file__)
+this_dir = os.path.dirname(this_path)
+
+class TestWriteMetadata(xshared._TestCase):
+    test_dir = this_dir
+    transient_inputs = ['data_input.nc']
+    transient_outputs = ['data_output.nc']
+    executable = './write.exe'
+    rtol = 5e-03
+
+    @classmethod
+    def make_a_write_test(cls, inf, nc_method='cdl_files',
+                          nclients=1, nservers=1):
+        """
+        this function makes a test case and returns it as a test function,
+        suitable to be dynamically added to a TestCase for running.
+
+        """
+        # always copy for value, don't pass by reference.
+        infcp = copy.copy(inf)
+        # expected by the fortran XIOS program's main.xml
+        inputfile = cls.transient_inputs[0]
+        outputfile = cls.transient_outputs[0]
+        ## This runs in XIOS3 but ignores the FRT coordinate in the output
+        ## ToDo understand and fix.
+        @unittest.skipIf(os.environ.get('MVER', '') == 'XIOS3/trunk',
+                         "skipping for XIOS3 trunk failing to write scalar frt")
+        def test_write_metadata(self):
+            # create a netCDF file using nc_method
+            cls.make_netcdf(infcp, inputfile, nc_method=nc_method)
+            cls.run_mpi_xios(nclients=nclients, nservers=nservers)
+
+            # load the result netCDF file
+            runfile = '{}/{}'.format(cls.test_dir, outputfile)
+            assert(os.path.exists(runfile))
+            run_cmd = ['ncdump', runfile]
+            ncd = subprocess.run(run_cmd, check=True, capture_output=True)
+
+            with open(f'{this_dir}/expected_domain_output.cdl') as fin:
+                exptd = fin.read()
+            emsg = ''
+            for eline, aline in zip(exptd.split('\n'),
+                                    ncd.stdout.decode("utf-8").split('\n')):
+                if eline != aline:
+                    # skip timestamp and uuid
+                    # until we work out how not to set these
+                    if not (':timeStamp' in eline or ':uuid' in eline):
+                        emsg += eline + '\n' + aline + '\n\tmismatch\n\n'
+
+            self.assertFalse(emsg, msg='\n\n' + emsg[0:300] + '\n... ...')
+        return test_write_metadata
+
+
+f = f'{this_dir}/data_input.cdl'
+# unique name for the test
+tname = 'test_{}'.format(os.path.splitext(os.path.basename(f))[0])
+# add the test as an attribute (function) to the test class
+
+setattr(TestWriteMetadata, tname,
+        TestWriteMetadata.make_a_write_test(f))