diff --git a/azspiceEnv b/azspiceEnv index c0bc22b..7d82211 100644 --- a/azspiceEnv +++ b/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++" diff --git a/xios_examples/gen_netcdf.py b/xios_examples/gen_netcdf.py index 1aab521..1e847d6 100644 --- a/xios_examples/gen_netcdf.py +++ b/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' diff --git a/xios_examples/setAttrCoord/Makefile b/xios_examples/setAttrCoord/Makefile new file mode 100644 index 0000000..b5e888c --- /dev/null +++ b/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 diff --git a/xios_examples/setAttrCoord/Readme.md b/xios_examples/setAttrCoord/Readme.md new file mode 100644 index 0000000..9f77c70 --- /dev/null +++ b/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. diff --git a/xios_examples/setAttrCoord/__init__.py b/xios_examples/setAttrCoord/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/xios_examples/setAttrCoord/axis_check.xml b/xios_examples/setAttrCoord/axis_check.xml new file mode 100644 index 0000000..31def66 --- /dev/null +++ b/xios_examples/setAttrCoord/axis_check.xml @@ -0,0 +1,20 @@ + + + + + + + + + + + + + + + + + + + + diff --git a/xios_examples/setAttrCoord/data_input.cdl b/xios_examples/setAttrCoord/data_input.cdl new file mode 100644 index 0000000..efe1dad --- /dev/null +++ b/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 ; + +} diff --git a/xios_examples/setAttrCoord/expected_domain_output.cdl b/xios_examples/setAttrCoord/expected_domain_output.cdl new file mode 100644 index 0000000..6bd4bf3 --- /dev/null +++ b/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 = _ ; +} diff --git a/xios_examples/setAttrCoord/iodef.xml b/xios_examples/setAttrCoord/iodef.xml new file mode 100644 index 0000000..37acd4d --- /dev/null +++ b/xios_examples/setAttrCoord/iodef.xml @@ -0,0 +1,9 @@ + + + + + + + + + diff --git a/xios_examples/setAttrCoord/main.xml b/xios_examples/setAttrCoord/main.xml new file mode 100644 index 0000000..483e65c --- /dev/null +++ b/xios_examples/setAttrCoord/main.xml @@ -0,0 +1,65 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + osgb: x y egm2008:alt + + + + + + + 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]] + + + 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]] + + + + + + Attribute demonstration + + Attribute demonstration + + + + + + diff --git a/xios_examples/setAttrCoord/test_write_metadata.py b/xios_examples/setAttrCoord/test_write_metadata.py new file mode 100644 index 0000000..c6e216b --- /dev/null +++ b/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)) diff --git a/xios_examples/setAttrCoord/write.F90 b/xios_examples/setAttrCoord/write.F90 new file mode 100644 index 0000000..0c17128 --- /dev/null +++ b/xios_examples/setAttrCoord/write.F90 @@ -0,0 +1,154 @@ +!----------------------------------------------------------------------------- +! (C) Crown copyright 2025 Met Office. All rights reserved. +! The file LICENCE, distributed with this code, contains details of the terms +! under which the code may be used. +!----------------------------------------------------------------------------- +!> Programme to just write data with coordinate system definition metadata. +!> +program write + use xios + use mpi + USE, INTRINSIC :: ISO_C_BINDING + + implicit none + + integer :: comm = -1 + integer :: rank = -1 + integer :: npar = 0 + + call initialise() + call simulate() + call finalise() +contains + + subroutine initialise() + + type(xios_date) :: origin + type(xios_date) :: start + type(xios_duration) :: tstep + integer :: mpi_error + integer :: lenx + integer :: leny + integer :: lenz + double precision :: frtv + double precision, dimension (:), allocatable :: y_vals, x_vals, altvals + character(len=64) :: t_origin + + origin = xios_date(2022, 2, 2, 12, 0, 0) + start = xios_date(2022, 12, 13, 12, 0, 0) + tstep = xios_hour + + ! Initialise MPI and XIOS + call MPI_INIT(mpi_error) + call xios_initialize('client', return_comm=comm) + + call MPI_Comm_rank(comm, rank, mpi_error) + call MPI_Comm_size(comm, npar, mpi_error) + + ! use the axis_check context to obtain sizing information on all arrays + ! for use in defining the main context interpretation + + call xios_context_initialize('axis_check', comm) + call xios_set_time_origin(origin) + call xios_set_start_date(start) + call xios_set_timestep(tstep) + + call xios_close_context_definition() + + ! fetch sizes of axes from the input file for allocate + call xios_get_axis_attr('x', n_glo=lenx) + call xios_get_axis_attr('y', n_glo=leny) + call xios_get_axis_attr('alt', n_glo=lenz) + + allocate ( x_vals(lenx) ) + allocate ( y_vals(leny) ) + allocate ( altvals(lenz) ) + + ! fetch coordinate value arrays from the input file + call xios_get_axis_attr('x', value=x_vals) + call xios_get_axis_attr('y', value=y_vals) + call xios_get_axis_attr('alt', value=altvals) + + ! finalise axis_check context, no longer in use + call xios_set_current_context('axis_check') + call xios_context_finalize() + + ! initialize the main context for interacting with the data. + call xios_context_initialize('main', comm) + + call xios_set_time_origin(origin) + call xios_set_start_date(start) + call xios_set_timestep(tstep) + + call xios_set_axis_attr("x", n_glo=lenx, n=lenx, begin=0) + call xios_set_axis_attr("x", value=x_vals) + call xios_set_axis_attr("y", n_glo=leny, n=leny, begin=0) + call xios_set_axis_attr("y", value=y_vals) + call xios_set_axis_attr("alt", n_glo=lenz, n=lenz, begin=0) + call xios_set_axis_attr("alt", value=altvals) + call xios_set_file_attr("data_output", convention_str="CF-1.6, UGRID") + call xios_set_file_attr("data_output", description="LFRic file format v0.2.0") + if (.true.) then + call xios_set_axis_attr("x", standard_name="projection_x_coordinate", & + unit="m", long_name="x coordinate of projection") + call xios_set_axis_attr("y", standard_name="projection_y_coordinate", & + unit="m", long_name="y coordinate of projection") + end if + call xios_date_convert_to_string(origin, t_origin) + call xios_set_scalar_attr("frt", unit="seconds since "//t_origin) + frtv = dble(xios_date_convert_to_seconds(start)) + call xios_set_scalar_attr("frt", value=frtv) + + call xios_close_context_definition() + + end subroutine initialise + + subroutine finalise() + + integer :: mpi_error + + ! Finalise all XIOS contexts and MPI + call xios_set_current_context('main') + call xios_context_finalize() + + call xios_finalize() + call MPI_Finalize(mpi_error) + + end subroutine finalise + + subroutine simulate() + + type(xios_date) :: start + type(xios_date) :: current + integer :: ts + integer :: lenx + integer :: leny + integer :: lenz + + ! Allocatable arrays, size is taken from input file + double precision, dimension (:,:,:), allocatable :: inodata + + ! obtain sizing of the grid for the array allocation + + call xios_get_axis_attr('x', n_glo=lenx) + call xios_get_axis_attr('y', n_glo=leny) + call xios_get_axis_attr('alt', n_glo=lenz) + + allocate ( inodata(lenz, leny, lenx) ) + + ! Load data from the input file + call xios_recv_field('odatain', inodata) + + do ts=1, 2 + call xios_update_calendar(ts) + call xios_get_current_date(current) + + ! Send (copy) the original data to the output file. + call xios_send_field('odata', inodata) + enddo + + deallocate (inodata) + + end subroutine simulate + +end program write diff --git a/xios_examples/setAttrCoord/xios.xml b/xios_examples/setAttrCoord/xios.xml new file mode 100644 index 0000000..9ec1df0 --- /dev/null +++ b/xios_examples/setAttrCoord/xios.xml @@ -0,0 +1,22 @@ + + + + + performance + + + 1.0 + + + + + true + + 100 + + + true + + + +