Description
This example shows: how to create a new GRIB message by cloning an existing message.
Source code
grib_clone.f90
! (C) Copyright 2005- ECMWF. ! ! This software is licensed under the terms of the Apache Licence Version 2.0 ! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. ! ! In applying this licence, ECMWF does not waive the privileges and immunities granted to it by ! virtue of its status as an intergovernmental organisation nor does it submit to any jurisdiction. ! ! ! FORTRAN 90 Implementation: grib_clone ! ! Description: how to create a new GRIB message by cloning ! an existing message. ! ! program clone use eccodes implicit none integer :: err,i integer :: nx, ny integer :: infile,outfile integer :: igrib_in integer :: igrib_out character(len=2) :: step real(kind=8), dimension(:,:), allocatable :: field2D call codes_open_file(infile,'../../data/constant_field.grib1','r') call codes_open_file(outfile,'out.clone.grib1','w') ! A new GRIB message is loaded from file. ! igrib is the GRIB id to be used in subsequent calls call codes_grib_new_from_file(infile,igrib_in) call codes_get(igrib_in,'Ni', nx) call codes_get(igrib_in,'Nj', ny) allocate(field2D(nx,ny),stat=err) if (err .ne. 0) then print*, 'Failed to allocate ', nx*ny, ' values' STOP end if ! clone the constant field to create 4 new GRIB messages do i=0,18,6 call codes_clone(igrib_in, igrib_out) write(step,'(i2)') i ! Careful: stepRange is a string (could be 0-6, 12-24, etc.) ! use adjustl to remove blank from the left. call codes_set(igrib_out,'stepRange',adjustl(step)) call generate_field(field2D) ! use pack to create 1D values call codes_set(igrib_out,'values',pack(field2D, mask=.true.)) ! write cloned messages to a file call codes_write(igrib_out,outfile) call codes_release(igrib_out) end do call codes_release(igrib_in) call codes_close_file(infile) call codes_close_file(outfile) deallocate(field2D) contains !====================================== subroutine generate_field(gfield2D) real(kind=8), dimension(:,:) :: gfield2D call random_number(gfield2D) end subroutine generate_field !====================================== end program clone
grib_clone.py
# # (C) Copyright 2005- ECMWF. # # This software is licensed under the terms of the Apache Licence Version 2.0 # which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. # # In applying this licence, ECMWF does not waive the privileges and immunities # granted to it by virtue of its status as an intergovernmental organisation # nor does it submit to any jurisdiction. # import traceback import sys import random from eccodes import * INPUT = '../../data/constant_field.grib1' OUTPUT = 'out.clone.grib' VERBOSE = 1 # verbose error reporting def example(): fin = open(INPUT, 'rb') fout = open(OUTPUT, 'wb') gid = codes_grib_new_from_file(fin) assert codes_is_missing(gid, 'Ni') == False assert codes_is_missing(gid, 'Nj') == False nx = codes_get(gid, 'Ni') ny = codes_get(gid, 'Nj') for step in range(0, 24, 6): clone_id = codes_clone(gid) codes_set(clone_id, 'step', step) values = [random.random() for i in range(nx * ny)] codes_set_values(clone_id, values) codes_write(clone_id, fout) codes_release(clone_id) codes_release(gid) fin.close() fout.close() def main(): try: example() except CodesInternalError as err: if VERBOSE: traceback.print_exc(file=sys.stderr) else: sys.stderr.write(err.msg + '\n') return 1 if __name__ == "__main__": sys.exit(main())
grib_clone.c
/* * (C) Copyright 2005- ECMWF. * * This software is licensed under the terms of the Apache Licence Version 2.0 * which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. * * In applying this licence, ECMWF does not waive the privileges and immunities granted to it by * virtue of its status as an intergovernmental organisation nor does it submit to any jurisdiction. */ /* * C Implementation: grib_clone * * Description: How to create a new GRIB message by cloning * an existing message. * */ #include <stdio.h> #include "eccodes.h" static void usage(const char *app) { fprintf(stderr,"Usage is: %s input_file ouput_file\n", app); } int main(int argc, char *argv[]) { FILE *in = NULL; FILE *out = NULL; codes_handle *source_handle = NULL; const void *buffer = NULL; size_t size = 0; int err = 0; if (argc != 3) { usage(argv[0]); return 1; } in = fopen(argv[1],"r"); out = fopen(argv[2],"w"); if (!in || !out) { perror("ERROR: unable to open files"); fclose(out); fclose(in); return 1; } /* loop over the messages in the source grib and clone them */ while ((source_handle = codes_handle_new_from_file(0, in, PRODUCT_GRIB, &err))!=NULL) { codes_handle *clone_handle = codes_handle_clone(source_handle); if (clone_handle == NULL) { perror("ERROR: could not clone field"); return 1; } /* This is the place where you may wish to modify the clone */ /* E.g. CODES_CHECK(codes_set_long(clone_handle, "centre", 250),0); etc... */ /* get the coded message in a buffer */ CODES_CHECK(codes_get_message(clone_handle,&buffer,&size),0); /* write the buffer to a file */ if(fwrite(buffer,1,size,out) != size) { perror(argv[1]); return 1; } codes_handle_delete(clone_handle); codes_handle_delete(source_handle); } fclose(out); fclose(in); return 0; }