diff --git a/configure.ac b/configure.ac index 0c0508a3..affc8c67 100644 --- a/configure.ac +++ b/configure.ac @@ -61,6 +61,7 @@ AC_ARG_ENABLE([quad-precision], AC_PROG_CC([icc gcc nvc]) AM_PROG_CC_C_O +AC_PROG_LN_S ## When autoconf v2.70 is more available, this can be replaced with: ##AC_C__GENERIC GX_C__GENERIC diff --git a/src/Makefile.am b/src/Makefile.am index da17811b..b60e295c 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -39,6 +39,7 @@ bin_PROGRAMS = \ decompress-ncc \ check_mask \ mppncscatter \ + mppdisttool \ ncexists \ $(combine_programs) \ $(grid_programs) @@ -446,3 +447,30 @@ runoff_regrid_parallel_LDADD = \ # ********************************************************* # transfer_to_mosaic_grid_SOURCES = \ transfer-mosaic-grid/transfer_to_mosaic.c + +# ********************************************************* # +# mppdisttool # +# ********************************************************* # +mppdisttool_SOURCES = \ + mpp-disttool/main.c \ + mpp-disttool/cmd_combine.c mpp-disttool/cmd_combine.h \ + mpp-disttool/cmd_scatter.c mpp-disttool/cmd_scatter.h \ + mpp-disttool/cmd_decompress.c mpp-disttool/cmd_decompress.h \ + mpp-disttool/cmd_check.c mpp-disttool/cmd_check.h \ + mpp-disttool/cmd_auto.c mpp-disttool/cmd_auto.h \ + mpp-disttool/compress.c mpp-disttool/compress.h \ + mpp-disttool/nc_utils.c mpp-disttool/nc_utils.h \ + mpp-disttool/domain.c mpp-disttool/domain.h \ + mpp-disttool/strlist.c mpp-disttool/strlist.h \ + mpp-disttool/xmalloc.c mpp-disttool/xmalloc.h + +mppdisttool_CFLAGS = $(NETCDF_CFLAGS) -I$(top_srcdir)/lib/libfrencutils \ + -std=c11 -D_POSIX_C_SOURCE=200809L +mppdisttool_LDADD = libver.a \ + $(top_builddir)/lib/libfrencutils/libfrencutils.a \ + $(NETCDF_LDFLAGS) $(NETCDF_LIBS) -lm $(RPATH_FLAGS) + +install-exec-hook: + for name in mppncscatter mppnccombine combine-ncc scatter-ncc decompress-ncc is-compressed; do \ + $(LN_S) -f mppdisttool $(DESTDIR)$(bindir)/$$name; \ + done diff --git a/src/mpp-disttool/cmd_auto.c b/src/mpp-disttool/cmd_auto.c new file mode 100644 index 00000000..01c868eb --- /dev/null +++ b/src/mpp-disttool/cmd_auto.c @@ -0,0 +1,278 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +/* + * Port of src/combine-restarts/combine_restarts.in. + * Scans CWD for files matching *{res,nc}.[0-9]{4}, groups by base name, + * then routes each group through the appropriate combine path. + */ +#define _POSIX_C_SOURCE 200809L +#include +#include +#include +#include +#include +#include +#include +#include +#include "cmd_auto.h" +#include "cmd_combine.h" +#include "cmd_check.h" +#include "compress.h" +#include "xmalloc.h" + +/* ------------------------------------------------------------------ */ +/* Simple dynamic string array */ +/* ------------------------------------------------------------------ */ +typedef struct { char **items; int n, cap; } strarray_t; + +static void sa_init(strarray_t *a) { a->items=NULL; a->n=0; a->cap=0; } +static void sa_push(strarray_t *a, const char *s) { + if (a->n >= a->cap) { + a->cap = a->cap ? a->cap*2 : 16; + a->items = (char**)realloc(a->items, (size_t)a->cap * sizeof(char*)); + } + a->items[a->n++] = strdup(s); +} +static void sa_free(strarray_t *a) { + for (int i=0;in;i++) free(a->items[i]); + free(a->items); a->n=0; a->cap=0; +} +static int sa_contains(strarray_t *a, const char *s) { + for (int i=0;in;i++) if (strcmp(a->items[i],s)==0) return 1; + return 0; +} +static int sa_cmp(const void *a, const void *b) { + return strcmp(*(const char**)a, *(const char**)b); +} + +/* ------------------------------------------------------------------ */ +/* File grouping */ +/* ------------------------------------------------------------------ */ +/* Match files in CWD whose names contain "res" or "nc" as a word and + end in .[0-9]{4}. Returns base names (suffix stripped) de-duplicated. */ +static strarray_t +find_base_names(int verbose) +{ + strarray_t bases; + sa_init(&bases); + + regex_t re_tail, re_word; + /* .[0-9]{4} at end of name */ + if (regcomp(&re_tail, "\\.[0-9][0-9][0-9][0-9]$", REG_EXTENDED|REG_NOSUB) != 0) return bases; + /* word "res" or "nc" (bounded by non-alphanumeric or start/end) */ + if (regcomp(&re_word, "(^|[^a-zA-Z0-9_])(res|nc)([^a-zA-Z0-9_]|$)", REG_EXTENDED|REG_NOSUB) != 0) { + regfree(&re_tail); return bases; + } + + DIR *dir = opendir("."); + if (!dir) { regfree(&re_tail); regfree(&re_word); return bases; } + + struct dirent *ent; + while ((ent = readdir(dir)) != NULL) { + const char *name = ent->d_name; + if (regexec(&re_tail, name, 0, NULL, 0) != 0) continue; + if (regexec(&re_word, name, 0, NULL, 0) != 0) continue; + + /* Strip trailing .[0-9]{4} to get base name. */ + char base[4096]; + size_t blen = strlen(name); + if (blen < 5) continue; /* Need at least ".NNNN" */ + if (blen - 5 >= sizeof(base)) { + fprintf(stderr, "Warning: filename too long, skipping: %s\n", name); + continue; + } + memcpy(base, name, blen - 5); + base[blen - 5] = '\0'; /* remove ".NNNN" */ + + if (!sa_contains(&bases, base)) { + if (verbose > 1) fprintf(stderr,"debug2: found base '%s'\n", base); + sa_push(&bases, base); + } + } + closedir(dir); + regfree(&re_tail); regfree(&re_word); + + /* Sort for deterministic order. */ + if (bases.n > 0) + qsort(bases.items, (size_t)bases.n, sizeof(char*), sa_cmp); + return bases; +} + +/* Collect input files for a given base name. */ +static strarray_t +collect_inputs(const char *base) +{ + strarray_t files; + sa_init(&files); + + char pattern[4096]; + snprintf(pattern, sizeof(pattern), "\\.[0-9][0-9][0-9][0-9]$"); + regex_t re; + if (regcomp(&re, pattern, REG_EXTENDED|REG_NOSUB) != 0) return files; + + DIR *dir = opendir("."); + if (!dir) { regfree(&re); return files; } + + struct dirent *ent; + size_t blen = strlen(base); + while ((ent = readdir(dir)) != NULL) { + const char *name = ent->d_name; + if (strncmp(name, base, blen) != 0) continue; + if (strlen(name) != blen + 5) continue; /* base + .NNNN */ + if (regexec(&re, name, 0, NULL, 0) != 0) continue; + sa_push(&files, name); + } + closedir(dir); + regfree(&re); + + if (files.n > 0) + qsort(files.items, (size_t)files.n, sizeof(char*), sa_cmp); + return files; +} + +/* ------------------------------------------------------------------ */ +/* cmd_auto */ +/* ------------------------------------------------------------------ */ +int +cmd_auto(int argc, char **argv) +{ + int verbose = 0; + int dryrun = 0; + int removein = 0; + const char *combineNccOpts = NULL; + const char *mppCombOpts = "-h 16384 -m"; + + for (int i = 1; i < argc; i++) { + if (strcmp(argv[i],"-v")==0) verbose++; + else if (strcmp(argv[i],"-n")==0) dryrun = 1; + else if (strcmp(argv[i],"-r")==0) removein = 1; + else if (strcmp(argv[i],"-C")==0 && i+1= 0) { + is_compressed_file = file_is_compressed(probeid) ? 1 : 0; + nc_close(probeid); + } + if (!is_compressed_file) + is_iceberg_file = check_is_iceberg(inputs.items[0]) ? 1 : 0; + + if (verbose > 1) { + fprintf(stderr,"debug2: '%s': compressed=%d iceberg=%d\n", + base, is_compressed_file, is_iceberg_file); + } + + int r = 0; + if (dryrun) { + /* Print what would be executed. */ + if (is_compressed_file) { + printf("combine-ncc"); + if (combineNccOpts) printf(" %s", combineNccOpts); + } else if (is_iceberg_file) { + printf("iceberg_comb.sh"); + } else { + printf("mppnccombine"); + if (mppCombOpts) printf(" %s", mppCombOpts); + printf(" %s", base); + } + for (int f = 0; f < inputs.n; f++) printf(" %s", inputs.items[f]); + if (is_compressed_file) printf(" %s", base); + printf("\n"); + if (removein) { + printf("rm -f"); + for (int f = 0; f < inputs.n; f++) printf(" %s", inputs.items[f]); + printf("\n"); + } + } else { + if (is_compressed_file) { + opts.land_mode = 1; + r = combine_land(&opts); + } else if (is_iceberg_file) { + opts.iceberg_mode = 1; + r = combine_iceberg(&opts); + } else { + opts.mpp_mode = 1; + r = combine_mpp(&opts); + } + if (r != 0) { + fprintf(stderr,"auto: error combining '%s'\n", base); + ret = 1; + } + } + sa_free(&inputs); + } + sa_free(&bases); + return ret; +} diff --git a/src/mpp-disttool/cmd_auto.h b/src/mpp-disttool/cmd_auto.h new file mode 100644 index 00000000..611e93d9 --- /dev/null +++ b/src/mpp-disttool/cmd_auto.h @@ -0,0 +1,28 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_CMD_AUTO_H +#define MPPDISTTOOL_CMD_AUTO_H + +/* Entry point: "mppdisttool auto [-C opts] [-M opts] [-n] [-r] [-v]" + Replicates combine_restarts: scans CWD for *{res,nc}.####, groups by + base name, and dispatches to combine_land / combine_iceberg / combine_mpp. */ +int cmd_auto(int argc, char **argv); + +#endif /* MPPDISTTOOL_CMD_AUTO_H */ diff --git a/src/mpp-disttool/cmd_check.c b/src/mpp-disttool/cmd_check.c new file mode 100644 index 00000000..e3788f60 --- /dev/null +++ b/src/mpp-disttool/cmd_check.c @@ -0,0 +1,117 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#include +#include +#include +#include +#include +#include "cmd_check.h" +#include "compress.h" + +/* Icebergs file format version thresholds (matches iceberg_comb.sh). */ +#define ICEBERG_MAJOR_MIN 0 +#define ICEBERG_MINOR_MIN 1 + +bool +check_is_compressed(const char *filename) +{ + int ncid; + if (nc_open(filename, NC_NOWRITE, &ncid) != NC_NOERR) return false; + bool result = file_is_compressed(ncid); + nc_close(ncid); + return result; +} + +bool +check_is_iceberg(const char *filename) +{ + int ncid; + if (nc_open(filename, NC_NOWRITE, &ncid) != NC_NOERR) return false; + + int major = -1, minor = -1, nfiles = -1; + nc_get_att_int(ncid, NC_GLOBAL, "file_format_major_version", &major); + nc_get_att_int(ncid, NC_GLOBAL, "file_format_minor_version", &minor); + nc_get_att_int(ncid, NC_GLOBAL, "NumFilesInSet", &nfiles); + nc_close(ncid); + + if (major < 0 || minor < 0 || nfiles <= 0) return false; + if (major > ICEBERG_MAJOR_MIN) return true; + if (major == ICEBERG_MAJOR_MIN && minor >= ICEBERG_MINOR_MIN) return true; + return false; +} + +bool +check_iceberg_has_data(const char *filename) +{ + int ncid, unlimdimid; + if (nc_open(filename, NC_NOWRITE, &ncid) != NC_NOERR) return false; + if (nc_inq_unlimdim(ncid, &unlimdimid) != NC_NOERR || unlimdimid < 0) { + nc_close(ncid); + return false; + } + size_t dimlen = 0; + nc_inq_dimlen(ncid, unlimdimid, &dimlen); + nc_close(ncid); + return dimlen > 0; +} + +int +cmd_check(int argc, char **argv) +{ + bool do_compressed = false; + bool do_iceberg = false; + bool quiet = false; + const char *filename = NULL; + + for (int i = 1; i < argc; i++) { + if (strcmp(argv[i], "--compressed") == 0) { do_compressed = true; } + else if (strcmp(argv[i], "--iceberg") == 0) { do_iceberg = true; } + else if (strcmp(argv[i], "--quiet") == 0) { quiet = true; } + else if (strcmp(argv[i], "-h") == 0 || strcmp(argv[i], "--help") == 0) { + printf("Usage: mppdisttool check [--compressed|--iceberg] [--quiet] FILE\n"); + return 0; + } else if (argv[i][0] != '-') { + filename = argv[i]; + } + } + + if (!filename) { + fprintf(stderr, "check: missing file operand\n"); + return 1; + } + if (!do_compressed && !do_iceberg) { + /* Default: behave like is-compressed */ + do_compressed = true; + } + + if (do_compressed) { + bool r = check_is_compressed(filename); + if (!quiet && r) + printf("%s is compressed\n", filename); + return r ? 0 : -1; + } + if (do_iceberg) { + bool r = check_is_iceberg(filename); + if (!quiet && r) + printf("%s is an iceberg restart file\n", filename); + return r ? 0 : 255; + } + return 1; +} diff --git a/src/mpp-disttool/cmd_check.h b/src/mpp-disttool/cmd_check.h new file mode 100644 index 00000000..a2267a26 --- /dev/null +++ b/src/mpp-disttool/cmd_check.h @@ -0,0 +1,41 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_CMD_CHECK_H +#define MPPDISTTOOL_CMD_CHECK_H + +#include + +/* Returns true (exit 0) if any dimension variable in the file has a + "compress" attribute. Matches is-compressed.c logic exactly. */ +bool check_is_compressed(const char *filename); + +/* Returns true if the file is a valid iceberg restart (has + file_format_major_version >= 0 and file_format_minor_version >= 1, + and has NumFilesInSet > 0). */ +bool check_is_iceberg(const char *filename); + +/* Returns true if the iceberg file has icebergs (UNLIMITED dim > 0). */ +bool check_iceberg_has_data(const char *filename); + +/* Entry point for "mppdisttool check [--compressed|--iceberg] [--quiet] FILE" + Returns 0 (compressed/iceberg), -1 (not compressed), 255 (not iceberg). */ +int cmd_check(int argc, char **argv); + +#endif /* MPPDISTTOOL_CMD_CHECK_H */ diff --git a/src/mpp-disttool/cmd_combine.c b/src/mpp-disttool/cmd_combine.c new file mode 100644 index 00000000..33dc928b --- /dev/null +++ b/src/mpp-disttool/cmd_combine.c @@ -0,0 +1,1184 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +/* + * Three combine paths: + * combine_mpp – port of mppnccombine.c (modernised API) + * combine_land – port of combine-ncc.F90 + * combine_iceberg – port of iceberg_comb.sh (no ncrcat/ncatted) + */ +#define _POSIX_C_SOURCE 200809L +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include "cmd_combine.h" +#include "cmd_check.h" +#include "nc_utils.h" +#include "compress.h" +#include "xmalloc.h" + +#define NC_CHECK(call) do { \ + int _err = (call); \ + if (_err != NC_NOERR) { \ + fprintf(stderr, "NetCDF error: %s\n", nc_strerror(_err)); \ + return 1; \ + } \ +} while (0) + +#ifndef MAX_BF +#define MAX_BF 100 +#endif +#ifndef DEFAULT_BF +#define DEFAULT_BF 1 +#endif +#ifndef NC_BLKSZ +#define NC_BLKSZ 65536 +#endif + +/* ================================================================== */ +/* combine_mpp – modern port of mppnccombine.c */ +/* ================================================================== */ + +/* Per-file metadata for mppnccombine algorithm. */ +typedef struct { + int ncid; + int ndims, nvars, ngatts, recdim; + char varname[NC_MAX_VARS][NC_MAX_NAME + 1]; + nc_type datatype[NC_MAX_VARS]; + int varndims[NC_MAX_VARS]; + int vardim[NC_MAX_VARS][NC_MAX_DIMS]; + int natts[NC_MAX_VARS]; + unsigned char vardecomp[NC_MAX_VARS]; + char dimname[NC_MAX_DIMS][NC_MAX_NAME + 1]; + size_t dimsize[NC_MAX_DIMS]; + size_t dimfullsize[NC_MAX_DIMS]; + size_t dimstart[NC_MAX_DIMS]; /* 0-based */ + size_t dimend[NC_MAX_DIMS]; /* 0-based; SIZE_MAX = not decomposed */ + unsigned char varmiss[NC_MAX_VARS]; + unsigned char varmissval[NC_MAX_VARS][8]; +} fileinfo_t; + +static void ***varbuf = NULL; +static size_t mem_allocated = 0; +static size_t estimated_maxrss = 0; +static int g_mem_dry_run = 0; +static int g_print_mem = 0; + +static int g_min(int a, int b) { return a < b ? a : b; } + +static size_t nc_type_sz(nc_type t) { + switch (t) { + case NC_BYTE: case NC_CHAR: return 1; + case NC_SHORT: return 2; + case NC_INT: return 4; + case NC_FLOAT: return 4; + case NC_DOUBLE: return 8; + default: return 0; + } +} + +/* Generic vara read/write by nc_type. */ +static int get_vara_typed(int ncid, int varid, const size_t *st, + const size_t *cnt, void *buf, nc_type t) { + switch (t) { + case NC_BYTE: case NC_CHAR: return nc_get_vara_text(ncid,varid,st,cnt,(char*)buf); + case NC_SHORT: return nc_get_vara_short(ncid,varid,st,cnt,(short*)buf); + case NC_INT: return nc_get_vara_int(ncid,varid,st,cnt,(int*)buf); + case NC_FLOAT: return nc_get_vara_float(ncid,varid,st,cnt,(float*)buf); + case NC_DOUBLE: return nc_get_vara_double(ncid,varid,st,cnt,(double*)buf); + default: return NC_EBADTYPE; + } +} + +static int put_vara_typed(int ncid, int varid, const size_t *st, + const size_t *cnt, const void *buf, nc_type t) { + switch (t) { + case NC_BYTE: case NC_CHAR: return nc_put_vara_text(ncid,varid,st,cnt,(const char*)buf); + case NC_SHORT: return nc_put_vara_short(ncid,varid,st,cnt,(const short*)buf); + case NC_INT: return nc_put_vara_int(ncid,varid,st,cnt,(const int*)buf); + case NC_FLOAT: return nc_put_vara_float(ncid,varid,st,cnt,(const float*)buf); + case NC_DOUBLE: return nc_put_vara_double(ncid,varid,st,cnt,(const double*)buf); + default: return NC_EBADTYPE; + } +} + +static int mpp_flush_decomp(fileinfo_t *outf, int nfiles, int r, int bf, int verbose) { + for (int v = 0; v < outf->nvars; v++) { + if (!outf->vardecomp[v]) continue; + size_t outstart[NC_MAX_DIMS] = {0}; + size_t count[NC_MAX_DIMS]; + int varrecdim = -1; + for (int d = 0; d < outf->varndims[v]; d++) { + outstart[d] = 0; + if (outf->vardim[v][d] == outf->recdim) { + count[d] = 1; varrecdim = d; + } else { + count[d] = outf->dimfullsize[outf->vardim[v][d]]; + } + } + if (varrecdim >= 0) outstart[varrecdim] = (size_t)r; + if (varrecdim < 0 && r > 0) continue; + if (!varbuf || !varbuf[r % bf] || !varbuf[r % bf][v]) continue; + if (put_vara_typed(outf->ncid, v, outstart, count, + varbuf[r % bf][v], outf->datatype[v]) != NC_NOERR) { + fprintf(stderr, "Error: cannot write variable %d.\n", v); + return 1; + } + } + return 0; +} + +static int mpp_process_vars(fileinfo_t *inf, fileinfo_t *outf, + int appendnc, int *nrecs, int *nblocks, int *bf, + int r, int nfiles, int f, int verbose, int missing) { + int recdimsize = (inf->recdim < 0) ? 1 : (int)inf->dimsize[inf->recdim]; + + if (*nrecs == 1) { + *nrecs = recdimsize; + if (*bf >= 1) { + if (*bf > *nrecs) *bf = *nrecs; + } else { + *bf = g_min(MAX_BF, *nrecs); + } + *nblocks = ((*nrecs % *bf) != 0) ? (*nrecs / *bf) + 1 : (*nrecs / *bf); + } else if (recdimsize != *nrecs) { + fprintf(stderr, "Error: different number of records than first input.\n"); + return 1; + } + + /* Allocate varbuf on first call. */ + if (!varbuf) { + /* Validate buffer factor to prevent overflow */ + if (*bf <= 0 || *bf > MAX_BF) { + fprintf(stderr, "Invalid buffer factor: %d\n", *bf); + return 1; + } + /* Check for overflow: bf * sizeof(void**) */ + if ((size_t)(*bf) > SIZE_MAX / sizeof(void**)) { + fprintf(stderr, "Buffer allocation size overflow\n"); + return 1; + } + size_t ptr_array_size = (size_t)(*bf) * sizeof(void**); + /* Check for overflow: bf * NC_MAX_VARS * sizeof(void*) */ + if ((size_t)(*bf) > SIZE_MAX / NC_MAX_VARS || + (size_t)(*bf) * NC_MAX_VARS > SIZE_MAX / sizeof(void*)) { + fprintf(stderr, "Buffer allocation size overflow\n"); + return 1; + } + size_t data_array_size = (size_t)(*bf) * NC_MAX_VARS * sizeof(void*); + /* Check for overflow when adding the two sizes */ + if (ptr_array_size > SIZE_MAX - data_array_size) { + fprintf(stderr, "Buffer allocation size overflow\n"); + return 1; + } + size_t nbytes = ptr_array_size + data_array_size; + if (g_mem_dry_run) estimated_maxrss += nbytes; + varbuf = (void***)calloc(nbytes, 1); + if (!varbuf) { fprintf(stderr,"Cannot allocate varbuf.\n"); exit(1); } + for (int z = 0; z < *bf; z++) + varbuf[z] = (void**)((size_t)varbuf + (*bf)*sizeof(void**) + z*NC_MAX_VARS*sizeof(void*)); + } + + for (int v = 0; v < inf->nvars; v++) { + size_t recsize = 1, recfullsize = 1; + int varrecdim = -1; + size_t instart[NC_MAX_DIMS] = {0}; + size_t outstart[NC_MAX_DIMS] = {0}; + size_t count[NC_MAX_DIMS]; + + for (int d = 0; d < inf->varndims[v]; d++) { + if (inf->vardim[v][d] == inf->recdim) { + count[d] = 1; varrecdim = d; + } else { + count[d] = inf->dimsize[inf->vardim[v][d]]; + recsize *= count[d]; + instart[d] = 0; + outstart[d] = (inf->dimend[inf->vardim[v][d]] == (size_t)(-1)) ? 0 : inf->dimstart[inf->vardim[v][d]]; + recfullsize *= inf->dimfullsize[inf->vardim[v][d]]; + } + } + + /* Skip if not the right record / already written. */ + if (r > 0) { + int dimid_check; + bool is_dim = (nc_inq_dimid(inf->ncid, inf->varname[v], &dimid_check) == NC_NOERR); + if (is_dim) { + if (dimid_check == inf->recdim) { if (f != 0) continue; } + else continue; + } else { + if (varrecdim < 0) continue; + if (!inf->vardecomp[v] && f > 0) continue; + } + } else { + if (!inf->vardecomp[v] && appendnc) continue; + } + + void *values = malloc(nc_type_sz(inf->datatype[v]) * recsize); + if (!values) { + fprintf(stderr, "Cannot allocate values buffer.\n"); + return 1; + } + + if (varrecdim >= 0) instart[varrecdim] = outstart[varrecdim] = (size_t)r; + + if (get_vara_typed(inf->ncid, v, instart, count, values, inf->datatype[v]) != NC_NOERR) { + fprintf(stderr,"Error reading variable '%s'.\n", inf->varname[v]); + free(values); return 1; + } + + if (!inf->vardecomp[v] && !g_mem_dry_run) { + if (put_vara_typed(outf->ncid, v, outstart, count, values, inf->datatype[v]) != NC_NOERR) { + fprintf(stderr,"Error writing variable '%s'.\n", inf->varname[v]); + free(values); return 1; + } + } else { + /* Buffer decomposed var. */ + if (!varbuf[r % (*bf)][v]) { + size_t vsz = nc_type_sz(inf->datatype[v]) * recfullsize; + if (g_mem_dry_run) { estimated_maxrss += vsz; varbuf[r%(*bf)][v] = (void*)"x"; free(values); continue; } + varbuf[r % (*bf)][v] = calloc(vsz, 1); + if (!varbuf[r % (*bf)][v]) { fprintf(stderr,"Cannot allocate varbuf entry.\n"); free(values); return 1; } + mem_allocated += vsz; + if (missing && outf->varmiss[v]) { + /* Fill with missing_value. */ + size_t esz = nc_type_sz(inf->datatype[v]); + for (size_t s = 0; s < recfullsize; s++) + memcpy((char*)varbuf[r%(*bf)][v] + s*esz, outf->varmissval[v], esz); + } + } + + /* Scatter decomposed data into full buffer. */ + int ndv = inf->varndims[v]; + if (ndv == 0) { + memcpy(varbuf[r%(*bf)][v], values, nc_type_sz(inf->datatype[v])); + free(values); continue; + } + + /* Build stride/offset for scattering. */ + size_t imax, jmax=1, kmax=1, lmax=1; + size_t imaxfull, jmaxfull=1, kmaxfull=1; + + imax = inf->dimsize[inf->vardim[v][ndv-1]]; + imaxfull = inf->dimfullsize[inf->vardim[v][ndv-1]]; + if (ndv > 1 && inf->vardim[v][ndv-2] != inf->recdim) { + jmax = inf->dimsize[inf->vardim[v][ndv-2]]; + jmaxfull = inf->dimfullsize[inf->vardim[v][ndv-2]]; + } + if (ndv > 2 && inf->vardim[v][ndv-3] != inf->recdim) { + kmax = inf->dimsize[inf->vardim[v][ndv-3]]; + kmaxfull = inf->dimfullsize[inf->vardim[v][ndv-3]]; + } + if (ndv > 3 && inf->vardim[v][ndv-4] != inf->recdim) { + lmax = inf->dimsize[inf->vardim[v][ndv-4]]; + (void)inf->dimfullsize[inf->vardim[v][ndv-4]]; /* lmaxfull not used in scatter */ + } + + size_t ioff = (outstart[ndv-1]); + size_t joff = (ndv > 1) ? outstart[ndv-2] : 0; + size_t koff = (ndv > 2) ? outstart[ndv-3] : 0; + size_t loff = (ndv > 3) ? outstart[ndv-4] : 0; + if (varrecdim >= 0) { + if (ndv==1) ioff=0; else if (ndv==2) joff=0; + else if (ndv==3) koff=0; else loff=0; + } + + size_t ijfull = imaxfull * jmaxfull; + size_t ijkfull = ijfull * kmaxfull; + size_t esz = nc_type_sz(inf->datatype[v]); + size_t b = 0; + for (size_t l = 0; l < lmax; l++) + for (size_t k = 0; k < kmax; k++) + for (size_t j = 0; j < jmax; j++) + for (size_t i = 0; i < imax; i++) { + size_t off = (i+ioff) + (j+joff)*imaxfull + (k+koff)*ijfull + (l+loff)*ijkfull; + memcpy((char*)varbuf[r%(*bf)][v] + off*esz, (char*)values + b*esz, esz); + b++; + } + } + free(values); + } + return 0; +} + +static int mpp_process_file(const char *ncname, int appendnc, + fileinfo_t *outf, const char *outncname, + int *nfiles, int *nrecs, int *nblocks, int *bf, + int block, int f, int headerpad, int verbose, + int missing, int deflate, int deflation, int shuffle) { + fileinfo_t *inf = XMALLOC(fileinfo_t, 1); + size_t blksz = NC_BLKSZ; + + if (nc_open(ncname, NC_NOWRITE, &inf->ncid) != NC_NOERR) { + fprintf(stderr,"Error: cannot open '%s'\n", ncname); + free(inf); return 1; + } + + int nfiles2 = -1; + nc_get_att_int(inf->ncid, NC_GLOBAL, "NumFilesInSet", &nfiles2); + if (nfiles2 > 0) *nfiles = nfiles2; + + NC_CHECK(nc_inq(inf->ncid, &inf->ndims, &inf->nvars, &inf->ngatts, &inf->recdim)); + + for (int d = 0; d < inf->ndims; d++) { + NC_CHECK(nc_inq_dim(inf->ncid, d, inf->dimname[d], &inf->dimsize[d])); + inf->dimfullsize[d] = inf->dimsize[d]; + inf->dimstart[d] = 0; + inf->dimend[d] = (size_t)(-1); /* not decomposed */ + } + + if (block == 0 && !g_mem_dry_run) { + outf->nvars = inf->nvars; + outf->recdim = inf->recdim; + } + + for (int v = 0; v < inf->nvars; v++) { + NC_CHECK(nc_inq_var(inf->ncid, v, inf->varname[v], &inf->datatype[v], + &inf->varndims[v], inf->vardim[v], &inf->natts[v])); + + int dimid_check; + if (nc_inq_dimid(inf->ncid, inf->varname[v], &dimid_check) == NC_NOERR) { + int decomp[4]; + if (nc_get_att_int(inf->ncid, v, "domain_decomposition", decomp) == NC_NOERR) { + inf->dimfullsize[dimid_check] = (size_t)(decomp[1] - decomp[0] + 1); + inf->dimstart[dimid_check] = (size_t)(decomp[2] - decomp[0]); /* 0-based */ + inf->dimend[dimid_check] = (size_t)(decomp[3] - decomp[0]); /* 0-based */ + } + } + } + + for (int v = 0; v < inf->nvars; v++) { + inf->vardecomp[v] = 0; + for (int d = 0; d < inf->varndims[v]; d++) { + if (inf->dimend[inf->vardim[v][d]] != (size_t)(-1)) { + inf->vardecomp[v] = 1; break; + } + } + if (!appendnc && !g_mem_dry_run) { + if (block == 0) { + outf->varndims[v] = inf->varndims[v]; + for (int d = 0; d < inf->ndims; d++) + outf->dimfullsize[d] = inf->dimfullsize[d]; + for (int d = 0; d < inf->varndims[v]; d++) + outf->vardim[v][d] = inf->vardim[v][d]; + outf->vardecomp[v] = inf->vardecomp[v]; + outf->datatype[v] = inf->datatype[v]; + strcpy(outf->varname[v], inf->varname[v]); + outf->varmiss[v] = 0; + } + } + } + + if (!appendnc && !g_mem_dry_run) { + /* Determine output format. */ + int ncinformat, ncoutformat; + nc_inq_format(inf->ncid, &ncinformat); + nc_inq_format(outf->ncid, &ncoutformat); + + if (ncoutformat == NC_FORMAT_CLASSIC) { + int newmode; + switch (ncinformat) { + case NC_FORMAT_CLASSIC: case NC_FORMAT_64BIT_OFFSET: + newmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + case NC_FORMAT_NETCDF4: + newmode = NC_CLOBBER | NC_NETCDF4; break; + case NC_FORMAT_NETCDF4_CLASSIC: + newmode = NC_CLOBBER | NC_NETCDF4 | NC_CLASSIC_MODEL; break; + default: + newmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + } + nc_close(outf->ncid); + if (nc__create(outncname, newmode, 0, &blksz, &outf->ncid) != NC_NOERR) { + fprintf(stderr,"Error: cannot recreate output '%s'.\n", outncname); + free(inf); return 1; + } + nc_set_fill(outf->ncid, NC_NOFILL, NULL); + } + + /* Define dimensions in output. */ + for (int d = 0; d < inf->ndims; d++) { + int dummy; + size_t dlen = (d == inf->recdim) ? NC_UNLIMITED : inf->dimfullsize[d]; + nc_def_dim(outf->ncid, inf->dimname[d], dlen, &dummy); + } + + /* Define variables + copy attributes. */ + for (int v = 0; v < inf->nvars; v++) { + int newv; + nc_def_var(outf->ncid, inf->varname[v], inf->datatype[v], + inf->varndims[v], inf->vardim[v], &newv); + char attname[NC_MAX_NAME + 1]; + for (int n = 0; n < inf->natts[v]; n++) { + nc_inq_attname(inf->ncid, v, n, attname); + if (missing && strcmp(attname, "missing_value") == 0) { + outf->varmiss[v] = 1; + size_t esz = nc_type_sz(inf->datatype[v]); + nc_get_att(inf->ncid, v, "missing_value", outf->varmissval[v]); + (void)esz; + } + if (strcmp(attname, "domain_decomposition") == 0) continue; + nc_copy_att(inf->ncid, v, attname, outf->ncid, newv); + } + } + + /* Copy global attributes (skip NumFilesInSet). */ + char attname[NC_MAX_NAME + 1]; + for (int n = 0; n < inf->ngatts; n++) { + nc_inq_attname(inf->ncid, NC_GLOBAL, n, attname); + if (strcmp(attname, "NumFilesInSet") == 0) continue; + if (strcmp(attname, "filename") == 0) { + nc_put_att_text(outf->ncid, NC_GLOBAL, attname, + strlen(outncname), outncname); + continue; + } + nc_copy_att(inf->ncid, NC_GLOBAL, attname, outf->ncid, NC_GLOBAL); + } + + if (deflate && deflation > 0) + for (int v = 0; v < inf->nvars; v++) + nc_def_var_deflate(outf->ncid, v, shuffle, deflate, deflation); + + nc__enddef(outf->ncid, (size_t)headerpad, 4, 0, 4); + } + + /* Read / buffer variables. */ + int r = block * (*bf); + int err = 0; + do { + err += mpp_process_vars(inf, outf, appendnc, nrecs, nblocks, bf, + r, *nfiles, f, verbose, missing); + r++; + appendnc = 1; + } while (r < g_min((block + 1) * (*bf), *nrecs)); + + nc_close(inf->ncid); + free(inf); + return err ? 1 : 0; +} + +int +combine_mpp(combine_opts_t *opts) +{ + g_mem_dry_run = opts->mem_dry_run; + g_print_mem = opts->print_mem; + mem_allocated = 0; + estimated_maxrss = 0; + varbuf = NULL; + + int bf = opts->blocking > 0 ? opts->blocking : DEFAULT_BF; + int nstart = opts->nstart; + int nend = opts->nend; + int force = opts->force; + int verbose = opts->verbose; + int missing = opts->missing; + int removein = opts->removein; + int headerpad = opts->headerpad > 0 ? opts->headerpad : 16384; + int deflate = opts->deflate; + int deflation = opts->deflation; + int shuffle = opts->shuffle; + + int nfiles = -1, nrecs = 1, nblocks = 1; + + const char *outncname = opts->outfile; + size_t blksz = NC_BLKSZ; + + fileinfo_t *outf = XMALLOC(fileinfo_t, 1); + memset(outf, 0, sizeof(*outf)); + + struct stat statbuf; + int create_mode = NC_NOCLOBBER; + if (opts->format_64) create_mode |= NC_64BIT_OFFSET; + if (opts->format_n4) create_mode = NC_NOCLOBBER | NC_NETCDF4 | NC_CLASSIC_MODEL; + + if (!opts->appendnc) { + if (stat(outncname, &statbuf) == 0) { + fprintf(stderr,"Error: output file '%s' already exists.\n", outncname); + free(outf); return 1; + } + if (nc__create(outncname, create_mode, 0, &blksz, &outf->ncid) != NC_NOERR) { + fprintf(stderr,"Error: cannot create output '%s'.\n", outncname); + free(outf); return 1; + } + nc_set_fill(outf->ncid, NC_NOFILL, NULL); + } else { + if (nc_open(outncname, NC_WRITE, &outf->ncid) != NC_NOERR) { + fprintf(stderr,"Error: cannot open '%s' for appending.\n", outncname); + free(outf); return 1; + } + } + + int appendnc = opts->appendnc; + int infile_errors = 0, outfile_errors = 0; + int peWidth = -1; + char infilename[2048]; + (void)0; + + /* --- Explicit input files on command line --- */ + if (opts->ninfiles > 0) { + for (int block = 0; block < nblocks; block++) { + int f = 0; + for (int a = 0; a < opts->ninfiles; a++) { + if (stat(opts->infiles[a], &statbuf) != 0) { + if (!force) { + fprintf(stderr,"ERROR: missing file '%s'.\n", opts->infiles[a]); + unlink(outncname); free(outf); return 9; + } + infile_errors = 1; + } + if (verbose) printf(" processing '%s'\n", opts->infiles[a]); + if (mpp_process_file(opts->infiles[a], appendnc, outf, outncname, + &nfiles, &nrecs, &nblocks, &bf, block, f, + headerpad, verbose, missing, deflate, deflation, shuffle)) + infile_errors = 1; + appendnc = 1; f++; + if (f == nfiles || a == opts->ninfiles - 1) { + if (g_mem_dry_run) { + printf("%.0f\n", ceil((float)estimated_maxrss/(1024*1024))); + goto done; + } + for (int r = block * bf; r < g_min((block+1)*bf, nrecs); r++) + outfile_errors += mpp_flush_decomp(outf, nfiles, r, bf, verbose); + f = 0; continue; + } + } + for (int k = 0; k < bf; k++) + for (int v = 0; v < NC_MAX_VARS; v++) + if (varbuf && varbuf[k][v]) { free(varbuf[k][v]); varbuf[k][v] = NULL; } + } + } else { + /* --- Auto-discover files by numeric extension --- */ + if (nend < 0) nend = nstart + 1; + for (int block = 0; block < nblocks; block++) { + int f = 0; + for (int a = nstart; a <= nend; a++) { + if (peWidth < 0) { + sprintf(infilename, "%s.%04d", outncname, a); + if (stat(infilename, &statbuf) == 0) peWidth = 4; + else { + sprintf(infilename, "%s.%06d", outncname, a); + if (stat(infilename, &statbuf) == 0) peWidth = 6; + else continue; + } + } + sprintf(infilename, "%s.%0*d", outncname, peWidth, a); + if (stat(infilename, &statbuf) != 0) { + if (!force) { + fprintf(stderr,"ERROR: missing '%s'.\n", infilename); + unlink(outncname); free(outf); return 9; + } + infile_errors = 1; + } + if (mpp_process_file(infilename, appendnc, outf, outncname, + &nfiles, &nrecs, &nblocks, &bf, block, f, + headerpad, verbose, missing, deflate, deflation, shuffle)) + infile_errors = 1; + if (a == nstart && nfiles > 0 && opts->nend < 0) + nend = nstart + nfiles; + appendnc = 1; f++; + if (f == nfiles || a == nend) { + if (g_mem_dry_run) { + printf("%.0f\n", ceil((float)estimated_maxrss/(1024*1024))); + goto done; + } + for (int r = block*bf; r < g_min((block+1)*bf, nrecs); r++) + outfile_errors += mpp_flush_decomp(outf, nfiles, r, bf, verbose); + f = 0; break; + } + } + for (int k = 0; k < bf; k++) + for (int v = 0; v < NC_MAX_VARS; v++) + if (varbuf && varbuf[k][v]) { free(varbuf[k][v]); varbuf[k][v] = NULL; } + } + } + +done: + if (nc_sync(outf->ncid) != NC_NOERR) outfile_errors++; + if (nc_close(outf->ncid) != NC_NOERR) outfile_errors++; + free(outf); + + if (!infile_errors && !outfile_errors) { + if (removein) { + if (opts->ninfiles > 0) { + for (int a = 0; a < opts->ninfiles; a++) { + if (stat(opts->infiles[a], &statbuf) == 0) unlink(opts->infiles[a]); + } + } else { + for (int a = nstart; a <= nend; a++) { + sprintf(infilename, "%s.%0*d", outncname, peWidth > 0 ? peWidth : 4, a); + if (stat(infilename, &statbuf) == 0) unlink(infilename); + } + } + } + return 0; + } + fprintf(stderr,"Warning: output file may be incomplete.\n"); + return 1; +} + +/* ================================================================== */ +/* combine_land – port of combine-ncc.F90 */ +/* ================================================================== */ +int +combine_land(combine_opts_t *opts) +{ + if (opts->ninfiles < 1) { + fprintf(stderr,"combine --land: no input files specified.\n"); + return 1; + } + + const char *outfile = opts->outfile; + int nfiles = opts->ninfiles; + int *input = XMALLOC(int, nfiles); + int ncid; + + /* Open all input files. */ + int in_format = NC_FORMAT_CLASSIC; + for (int i = 0; i < nfiles; i++) { + if (nc_open(opts->infiles[i], NC_NOWRITE, &input[i]) != NC_NOERR) { + fprintf(stderr,"combine --land: cannot open '%s'.\n", opts->infiles[i]); + free(input); return 1; + } + nc_inq_format(input[i], &in_format); + } + + /* Create output file with matching format. */ + int cmode; + switch (in_format) { + case NC_FORMAT_NETCDF4: cmode = NC_NETCDF4; break; + case NC_FORMAT_NETCDF4_CLASSIC: cmode = NC_NETCDF4 | NC_CLASSIC_MODEL; break; + case NC_FORMAT_64BIT_OFFSET: cmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + default: cmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + } + + size_t blksz = 65536; + if (nc__create(outfile, cmode, 0, &blksz, &ncid) != NC_NOERR) { + fprintf(stderr,"combine --land: cannot create '%s'.\n", outfile); + free(input); return 1; + } + + int ref = input[nfiles - 1]; /* Use last file as template. */ + + /* --- Define dimensions --- */ + int ndims; + nc_inq_ndims(ref, &ndims); + + /* For each dimension determine output length. */ + int *dim_is_compressed = XMALLOC(int, ndims); + int *dim_out_len = XMALLOC(int, ndims); + + for (int dimid = 0; dimid < ndims; dimid++) { + char dname[NC_MAX_NAME + 1]; + size_t dlen; + bool is_unlim; + ncu_inq_dim(ref, dimid, dname, &dlen, &is_unlim); + + int varid, attid; + dim_is_compressed[dimid] = 0; + dim_out_len[dimid] = (int)dlen; + + if (nc_inq_varid(ref, dname, &varid) == NC_NOERR && + nc_inq_attid(ref, varid, "compress", &attid) == NC_NOERR) { + dim_is_compressed[dimid] = 1; + /* Count valid (non-negative) indices across all files. */ + int total_len = 0; + for (int fi = 0; fi < nfiles; fi++) { + size_t slen; + if (ncu_inq_dim(input[fi], dimid, NULL, &slen, NULL) != NC_NOERR) continue; + int *buf = XMALLOC(int, slen); + if (ncu_get_var_int(input[fi], dname, buf) == NC_NOERR) + for (size_t k = 0; k < slen; k++) + if (buf[k] >= 0) total_len++; + free(buf); + } + dim_out_len[dimid] = total_len > 0 ? total_len : 1; + } + + int dummy; + nc_def_dim(ncid, dname, is_unlim ? NC_UNLIMITED : (size_t)dim_out_len[dimid], &dummy); + } + + /* --- Define variables --- */ + int nvars; + nc_inq_nvars(ref, &nvars); + for (int v = 0; v < nvars; v++) + ncu_clone_var(ref, v, ncid); + + /* --- Copy global attributes --- */ + ncu_copy_global_atts(ref, ncid); + + nc__enddef(ncid, 16384, 4, 0, 4); + + /* --- Copy uncompressed variables --- */ + for (int v = 0; v < nvars; v++) { + char vname[NC_MAX_NAME + 1]; + nc_type xtype; + int vndims, dimids[NC_MAX_DIMS]; + size_t dimlens[NC_MAX_DIMS]; + int nrec, recsize; + bool has_records; + + ncu_inq_var(ref, v, vname, &xtype, &vndims, dimids, dimlens, NULL, + NULL, &has_records, NULL, (size_t*)&recsize, &nrec); + + /* Does this variable have any compressed dimension? */ + int has_compressed = 0; + for (int d = 0; d < vndims; d++) + if (dim_is_compressed[dimids[d]]) { has_compressed++; break; } + + if (has_compressed) continue; /* handled later */ + + /* Copy uncompressed variable. */ + if (xtype == NC_CHAR) { + size_t total = (size_t)recsize * (size_t)nrec; + char *text = XMALLOC(char, total); + nc_get_var_text(ref, v, text); + int ov; nc_inq_varid(ncid, vname, &ov); + nc_put_var_text(ncid, ov, text); + free(text); + } else { + double *buf = XMALLOC(double, (size_t)recsize); + for (int rec = 0; rec < nrec; rec++) { + ncu_get_rec(ref, v, rec, buf); + int ov; nc_inq_varid(ncid, vname, &ov); + ncu_put_rec(ncid, ov, rec, buf); + } + free(buf); + } + } + + /* --- Copy compressed variables --- */ + for (int dimid = 0; dimid < ndims; dimid++) { + if (!dim_is_compressed[dimid]) continue; + + char dname[NC_MAX_NAME + 1]; + nc_inq_dimname(ref, dimid, dname); + + /* Get per-file compressed dimension sizes. */ + int *sizes = XMALLOC(int, nfiles); + int buflen = 0; + for (int fi = 0; fi < nfiles; fi++) { + size_t slen; + if (ncu_inq_dim(input[fi], dimid, NULL, &slen, NULL) == NC_NOERR) + sizes[fi] = (int)slen; + else + sizes[fi] = 0; + buflen += sizes[fi]; + } + + if (dim_out_len[dimid] == 0) { free(sizes); continue; } + + /* Build reordering index. */ + int *rank_buf = XMALLOC(int, buflen); + int k = 0; + for (int fi = 0; fi < nfiles; fi++) { + if (sizes[fi] == 0) continue; + if (ncu_get_var_int(input[fi], dname, rank_buf + k) != NC_NOERR) + memset(rank_buf + k, -1, sizes[fi] * sizeof(int)); + k += sizes[fi]; + } + + int *rank = XMALLOC(int, buflen); + rank_ascending(rank_buf, rank, buflen); + + /* Skip leading negatives. */ + int skip = 0; + while (skip < buflen && rank_buf[rank[skip]] < 0) skip++; + for (int i = 0; i < buflen - skip; i++) rank[i] = rank[i + skip]; + for (int i = buflen - skip; i < buflen; i++) rank[i] = -1; + + /* Process variables depending on this compressed dimension. */ + for (int v = 0; v < nvars; v++) { + char vname[NC_MAX_NAME + 1]; + nc_type xtype; + int vndims, dimids[NC_MAX_DIMS]; + size_t dimlens[NC_MAX_DIMS]; + int nrec; + bool has_records; + ncu_inq_var(ref, v, vname, &xtype, &vndims, dimids, dimlens, + NULL, NULL, &has_records, NULL, NULL, &nrec); + + int cdim = -1; + for (int d = 0; d < vndims; d++) + if (dimids[d] == dimid) { cdim = d; break; } + if (cdim < 0) continue; + + int ovarid; + if (nc_inq_varid(ncid, vname, &ovarid) != NC_NOERR) continue; + + /* Product of non-compressed dimensions. */ + size_t slice_count = 1; + for (int d = 0; d < vndims; d++) + if (d != cdim) slice_count *= dimlens[d]; + + double *ibuf = XMALLOC(double, (size_t)buflen); + double *obuf = XMALLOC(double, (size_t)dim_out_len[dimid]); + size_t start[NC_MAX_DIMS], cnt[NC_MAX_DIMS]; + + for (size_t si = 0; si < slice_count; si++) { + /* Determine start indices for this slice. */ + size_t ii = si; + for (int d = 0; d < vndims; d++) { + if (d == cdim) { start[d] = 0; } + else { start[d] = ii % dimlens[d]; ii /= dimlens[d]; } + } + for (int d = 0; d < vndims; d++) cnt[d] = 1; + + /* Read from all input files. */ + k = 0; + for (int fi = 0; fi < nfiles; fi++) { + if (sizes[fi] == 0) continue; + int ivarid; + if (nc_inq_varid(input[fi], vname, &ivarid) != NC_NOERR) { + k += sizes[fi]; continue; + } + cnt[cdim] = (size_t)sizes[fi]; + nc_get_vara_double(input[fi], ivarid, start, cnt, ibuf + k); + k += sizes[fi]; + } + + /* Reorder. */ + for (int i = 0; i < dim_out_len[dimid]; i++) + obuf[i] = ibuf[rank[i]]; + + /* Write to output. */ + cnt[cdim] = (size_t)dim_out_len[dimid]; + nc_put_vara_double(ncid, ovarid, start, cnt, obuf); + } + free(ibuf); free(obuf); + } + free(rank_buf); free(rank); free(sizes); + } + + nc_sync(ncid); nc_close(ncid); + for (int i = 0; i < nfiles; i++) nc_close(input[i]); + free(input); free(dim_is_compressed); free(dim_out_len); + + if (opts->removein) + for (int i = 0; i < nfiles; i++) unlink(opts->infiles[i]); + return 0; +} + +/* ================================================================== */ +/* combine_iceberg – port of iceberg_comb.sh (no ncrcat/ncatted) */ +/* ================================================================== */ +int +combine_iceberg(combine_opts_t *opts) +{ + if (opts->ninfiles < 2) { + fprintf(stderr,"combine --iceberg: need at least 2 input files.\n"); + return 1; + } + const char *outfile = opts->outfile; + + /* Validate: all files must be iceberg restart files. */ + for (int i = 0; i < opts->ninfiles; i++) { + if (!check_is_iceberg(opts->infiles[i])) { + fprintf(stderr,"combine --iceberg: '%s' is not an iceberg restart.\n", + opts->infiles[i]); + return 1; + } + } + + /* Filter to files that actually have icebergs (non-empty UNLIMITED dim). */ + char **combine_files = XMALLOC(char *, opts->ninfiles); + int ncombine = 0; + for (int i = 0; i < opts->ninfiles; i++) + if (check_iceberg_has_data(opts->infiles[i])) + combine_files[ncombine++] = opts->infiles[i]; + + if (ncombine == 0) { + /* No icebergs: copy first input as output. */ + FILE *src = fopen(opts->infiles[0], "rb"); + FILE *dst = fopen(outfile, "wb"); + if (!src || !dst) { + fprintf(stderr,"combine --iceberg: cannot copy file.\n"); + if (src) fclose(src); if (dst) fclose(dst); + free(combine_files); return 1; + } + char buf[65536]; size_t n; + while ((n = fread(buf, 1, sizeof(buf), src)) > 0) + fwrite(buf, 1, n, dst); + fclose(src); fclose(dst); + /* Remove NumFilesInSet attribute. */ + int ncid; + if (nc_open(outfile, NC_WRITE, &ncid) == NC_NOERR) { + nc_del_att(ncid, NC_GLOBAL, "NumFilesInSet"); + nc_close(ncid); + } + free(combine_files); + return 0; + } + + /* Open first non-empty file as template. */ + int ref_ncid; + if (nc_open(combine_files[0], NC_NOWRITE, &ref_ncid) != NC_NOERR) { + fprintf(stderr,"combine --iceberg: cannot open '%s'.\n", combine_files[0]); + free(combine_files); return 1; + } + + /* Determine output format. */ + int in_fmt = NC_FORMAT_CLASSIC; + nc_inq_format(ref_ncid, &in_fmt); + int cmode; + switch (in_fmt) { + case NC_FORMAT_NETCDF4: cmode = NC_NETCDF4; break; + case NC_FORMAT_NETCDF4_CLASSIC: cmode = NC_NETCDF4 | NC_CLASSIC_MODEL; break; + case NC_FORMAT_64BIT_OFFSET: cmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + default: cmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + } + + int out_ncid; + size_t blksz = 65536; + if (nc__create(outfile, cmode, 0, &blksz, &out_ncid) != NC_NOERR) { + fprintf(stderr,"combine --iceberg: cannot create '%s'.\n", outfile); + nc_close(ref_ncid); free(combine_files); return 1; + } + nc_set_fill(out_ncid, NC_NOFILL, NULL); + + /* Query unlimited dimension name and total combined size. */ + int unlimdimid; + char unlimited_name[NC_MAX_NAME + 1] = ""; + /* total_recs not needed; NC_UNLIMITED handles growth automatically. */ + + if (nc_inq_unlimdim(ref_ncid, &unlimdimid) == NC_NOERR && unlimdimid >= 0) + nc_inq_dimname(ref_ncid, unlimdimid, unlimited_name); + + for (int i = 0; i < ncombine; i++) { + int tmpid; + if (nc_open(combine_files[i], NC_NOWRITE, &tmpid) == NC_NOERR) { + int uid; + if (nc_inq_unlimdim(tmpid, &uid) == NC_NOERR && uid >= 0) { + size_t dlen; nc_inq_dimlen(tmpid, uid, &dlen); (void)dlen; + } + nc_close(tmpid); + } + } + + /* Define dimensions in output. */ + int ndims; + nc_inq_ndims(ref_ncid, &ndims); + for (int d = 0; d < ndims; d++) { + char dname[NC_MAX_NAME + 1]; size_t dlen; + nc_inq_dim(ref_ncid, d, dname, &dlen); + int dummy; + if (d == unlimdimid) + nc_def_dim(out_ncid, dname, NC_UNLIMITED, &dummy); + else + nc_def_dim(out_ncid, dname, dlen, &dummy); + } + + /* Define variables. */ + int nvars; + nc_inq_nvars(ref_ncid, &nvars); + for (int v = 0; v < nvars; v++) ncu_clone_var(ref_ncid, v, out_ncid); + + /* Copy global attributes, skip NumFilesInSet. */ + int ngatts; + nc_inq_natts(ref_ncid, &ngatts); + for (int i = 0; i < ngatts; i++) { + char aname[NC_MAX_NAME + 1]; + nc_inq_attname(ref_ncid, NC_GLOBAL, i, aname); + if (strcmp(aname, "NumFilesInSet") == 0) continue; + nc_copy_att(ref_ncid, NC_GLOBAL, aname, out_ncid, NC_GLOBAL); + } + + nc_enddef(out_ncid); + + /* Copy non-record variables from first file. */ + for (int v = 0; v < nvars; v++) { + char vname[NC_MAX_NAME + 1]; + nc_type xtype; + int vndims, dimids[NC_MAX_DIMS], natts; + nc_inq_var(ref_ncid, v, vname, &xtype, &vndims, dimids, &natts); + + bool is_rec = false; + for (int d = 0; d < vndims; d++) + if (dimids[d] == unlimdimid) { is_rec = true; break; } + if (is_rec) continue; + + size_t total = 1; + for (int d = 0; d < vndims; d++) { size_t dl; nc_inq_dimlen(ref_ncid,dimids[d],&dl); total*=dl; } + size_t esz = ncu_type_size(xtype); + void *buf = malloc(total * esz); + ncu_get_vara(ref_ncid, v, NULL, NULL, buf); + int ov; nc_inq_varid(out_ncid, vname, &ov); + /* Build start/count for put. */ + size_t st[NC_MAX_DIMS], cnt[NC_MAX_DIMS]; + for (int d = 0; d < vndims; d++) { st[d]=0; nc_inq_dimlen(ref_ncid,dimids[d],&cnt[d]); } + ncu_put_vara(out_ncid, ov, st, cnt, buf); + free(buf); + } + + /* Concatenate record variables along UNLIMITED dim. */ + size_t out_rec_offset = 0; + for (int fi = 0; fi < ncombine; fi++) { + int in_ncid; + if (nc_open(combine_files[fi], NC_NOWRITE, &in_ncid) != NC_NOERR) continue; + int uid; size_t nrecs_in = 0; + if (nc_inq_unlimdim(in_ncid, &uid) == NC_NOERR && uid >= 0) + nc_inq_dimlen(in_ncid, uid, &nrecs_in); + + for (size_t rec = 0; rec < nrecs_in; rec++) { + for (int v = 0; v < nvars; v++) { + char vname[NC_MAX_NAME + 1]; + nc_type xtype; + int vndims, dimids[NC_MAX_DIMS], natts; + nc_inq_var(ref_ncid, v, vname, &xtype, &vndims, dimids, &natts); + + int rec_dim = -1; + for (int d = 0; d < vndims; d++) + if (dimids[d] == unlimdimid) { rec_dim = d; break; } + if (rec_dim < 0) continue; + + /* Get var from input file. */ + int iv; + if (nc_inq_varid(in_ncid, vname, &iv) != NC_NOERR) continue; + + size_t st_in[NC_MAX_DIMS], cnt[NC_MAX_DIMS]; + size_t st_out[NC_MAX_DIMS]; + for (int d = 0; d < vndims; d++) { + if (d == rec_dim) { + st_in[d] = rec; cnt[d] = 1; + st_out[d] = out_rec_offset + rec; + } else { + st_in[d] = 0; st_out[d] = 0; + nc_inq_dimlen(in_ncid, dimids[d], &cnt[d]); + } + } + size_t total = 1; + for (int d = 0; d < vndims; d++) total *= cnt[d]; + size_t esz = ncu_type_size(xtype); + void *buf = malloc(total * esz); + ncu_get_vara(in_ncid, iv, st_in, cnt, buf); + + int ov; nc_inq_varid(out_ncid, vname, &ov); + ncu_put_vara(out_ncid, ov, st_out, cnt, buf); + free(buf); + } + } + out_rec_offset += nrecs_in; + nc_close(in_ncid); + } + nc_close(ref_ncid); + + /* Remove NumFilesInSet (already excluded above, but ensure it's absent). */ + nc_redef(out_ncid); + nc_del_att(out_ncid, NC_GLOBAL, "NumFilesInSet"); + nc_enddef(out_ncid); + + nc_sync(out_ncid); + nc_close(out_ncid); + free(combine_files); + + if (opts->removein) + for (int i = 0; i < opts->ninfiles; i++) unlink(opts->infiles[i]); + return 0; +} + +/* ================================================================== */ +/* cmd_combine – argument parsing + dispatch */ +/* ================================================================== */ +static void usage_combine(void) { + printf("Usage: mppdisttool combine [-o outfile] [-v] [-r] [--land|--iceberg|--mpp]\n" + " [-h N] [-m] [-k N] [-64|-n4] [-a] [-f] infile...\n"); +} + +int +cmd_combine(int argc, char **argv) +{ + combine_opts_t opts; + memset(&opts, 0, sizeof(opts)); + opts.headerpad = 16384; + opts.blocking = DEFAULT_BF; + opts.nend = -1; + + char **infiles = XMALLOC(char *, argc); + int ninfiles = 0; + + for (int i = 1; i < argc; i++) { + if (strcmp(argv[i],"--land") == 0) opts.land_mode = 1; + else if (strcmp(argv[i],"--iceberg") == 0) opts.iceberg_mode = 1; + else if (strcmp(argv[i],"--mpp") == 0) opts.mpp_mode = 1; + else if (strcmp(argv[i],"-v") == 0) opts.verbose++; + else if (strcmp(argv[i],"-r") == 0) opts.removein = 1; + else if (strcmp(argv[i],"-a") == 0) opts.appendnc = 1; + else if (strcmp(argv[i],"-f") == 0) opts.force = 1; + else if (strcmp(argv[i],"-m") == 0) opts.missing = 1; + else if (strcmp(argv[i],"-64") == 0) opts.format_64 = 1; + else if (strcmp(argv[i],"-n4") == 0) opts.format_n4 = 1; + else if (strcmp(argv[i],"-M") == 0) opts.print_mem = 1; + else if (strcmp(argv[i],"-x") == 0) opts.mem_dry_run = 1; + else if ((strcmp(argv[i],"-o") == 0) && i+1 < argc) { + opts.outfile = argv[++i]; + } else if ((strcmp(argv[i],"-h") == 0) && i+1 < argc) { + opts.headerpad = atoi(argv[++i]); + } else if ((strcmp(argv[i],"-k") == 0) && i+1 < argc) { + opts.blocking = atoi(argv[++i]); + } else if ((strcmp(argv[i],"-n") == 0) && i+1 < argc) { + opts.nstart = atoi(argv[++i]); + } else if ((strcmp(argv[i],"-e") == 0) && i+1 < argc) { + opts.nend = atoi(argv[++i]); + } else if ((strcmp(argv[i],"-d") == 0) && i+1 < argc) { + opts.deflate = 1; opts.deflation = atoi(argv[++i]); + } else if (strcmp(argv[i],"-s") == 0) { + opts.shuffle = 1; + } else if (strcmp(argv[i],"--help") == 0 || strcmp(argv[i],"-?") == 0) { + usage_combine(); free(infiles); return 0; + } else if (argv[i][0] != '-') { + infiles[ninfiles++] = argv[i]; + } + } + + /* First non-option arg is output if --mpp mode with no -o. */ + if (!opts.outfile && ninfiles > 0 && (opts.mpp_mode || (!opts.land_mode && !opts.iceberg_mode))) { + opts.outfile = infiles[0]; + opts.infiles = infiles + 1; + opts.ninfiles = ninfiles - 1; + } else { + opts.infiles = infiles; + opts.ninfiles = ninfiles; + } + + /* If no mode explicitly given, auto-detect from first input file. */ + if (!opts.land_mode && !opts.iceberg_mode && !opts.mpp_mode) { + if (opts.ninfiles > 0) { + const char *probe = opts.infiles[0]; + int probeid; + if (nc_open(probe, NC_NOWRITE, &probeid) == NC_NOERR) { + if (file_is_compressed(probeid)) opts.land_mode = 1; + else if (check_is_iceberg(probe)) opts.iceberg_mode = 1; + else opts.mpp_mode = 1; + nc_close(probeid); + } else opts.mpp_mode = 1; + } else opts.mpp_mode = 1; + } + + int ret; + if (opts.land_mode) ret = combine_land(&opts); + else if (opts.iceberg_mode) ret = combine_iceberg(&opts); + else ret = combine_mpp(&opts); + + free(infiles); + return ret; +} diff --git a/src/mpp-disttool/cmd_combine.h b/src/mpp-disttool/cmd_combine.h new file mode 100644 index 00000000..045a94dc --- /dev/null +++ b/src/mpp-disttool/cmd_combine.h @@ -0,0 +1,60 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_CMD_COMBINE_H +#define MPPDISTTOOL_CMD_COMBINE_H + +/* Options shared by all combine paths. */ +typedef struct { + char *outfile; /* output file name */ + char **infiles; /* input file names */ + int ninfiles; + int verbose; + int removein; /* remove input files on success */ + int appendnc; /* append to existing output */ + int headerpad; /* extra header padding bytes */ + int missing; /* init decomposed vars with missing_value */ + int force; /* combine even if files missing */ + int format_64; /* force 64-bit offset output */ + int format_n4; /* force NETCDF4_CLASSIC output */ + int deflate; /* deflate (NETCDF4 only) */ + int deflation; /* deflation level */ + int shuffle; /* shuffle filter */ + int nstart; /* first PE extension (no explicit infiles) */ + int nend; /* last PE extension (-1 = auto) */ + int blocking; /* record blocking factor */ + int mem_dry_run; /* print estimated memory and exit */ + int print_mem; /* print memory usage statistics */ + /* land path */ + int land_mode; + /* iceberg path */ + int iceberg_mode; + /* MPP (default) */ + int mpp_mode; +} combine_opts_t; + +/* Entry point: "mppdisttool combine [opts] [-o out] infiles..." */ +int cmd_combine(int argc, char **argv); + +/* Individual path entry points (used by cmd_auto). */ +int combine_land(combine_opts_t *opts); +int combine_iceberg(combine_opts_t *opts); +int combine_mpp(combine_opts_t *opts); + +#endif /* MPPDISTTOOL_CMD_COMBINE_H */ diff --git a/src/mpp-disttool/cmd_decompress.c b/src/mpp-disttool/cmd_decompress.c new file mode 100644 index 00000000..becb8a32 --- /dev/null +++ b/src/mpp-disttool/cmd_decompress.c @@ -0,0 +1,268 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +/* + * Port of src/land_utils/decompress-ncc.F90. + * Reads one or more compressed-by-gathering files and writes a single + * decompressed (full-grid) NetCDF file. + */ +#define _POSIX_C_SOURCE 200809L +#include +#include +#include +#include +#include +#include "cmd_decompress.h" +#include "nc_utils.h" +#include "compress.h" +#include "xmalloc.h" + +int +cmd_decompress(int argc, char **argv) +{ + int debug = 0; + bool add_missing = false; + const char **infiles = NULL; + int ninfiles = 0; + const char *outfile = NULL; + + /* Parse arguments. */ + for (int i = 1; i < argc; i++) { + if (strcmp(argv[i],"-v")==0 || strcmp(argv[i],"-D")==0 || + strcmp(argv[i],"--debug-level")==0 || strcmp(argv[i],"--verbosity-level")==0) { + if (i+1 < argc) debug = atoi(argv[++i]); else debug++; + } else if (strcmp(argv[i],"-m")==0 || strcmp(argv[i],"--add-missing-value")==0) { + add_missing = true; + } else if (strcmp(argv[i],"-h")==0 || strcmp(argv[i],"--help")==0) { + printf("Usage: mppdisttool decompress [-v level] [-m] infile... outfile\n"); + return 0; + } else if (argv[i][0] != '-') { + /* Collect positional args; last one is outfile. */ + if (!infiles) infiles = XMALLOC(const char *, argc); + infiles[ninfiles++] = argv[i]; + } + } + + if (ninfiles < 2) { + fprintf(stderr,"decompress: at least one input and one output file required.\n"); + return 1; + } + outfile = infiles[--ninfiles]; /* Last arg is output. */ + + /* Open input files. */ + int *input = XMALLOC(int, ninfiles); + int in_format = NC_FORMAT_CLASSIC; + for (int i = 0; i < ninfiles; i++) { + if (nc_open(infiles[i], NC_NOWRITE, &input[i]) != NC_NOERR) { + fprintf(stderr,"decompress: cannot open '%s'.\n", infiles[i]); + free(input); return 1; + } + nc_inq_format(input[i], &in_format); + } + + int ref = input[ninfiles - 1]; /* Template. */ + + /* Determine output format. */ + int cmode; + switch (in_format) { + case NC_FORMAT_NETCDF4: cmode = NC_NETCDF4; break; + case NC_FORMAT_NETCDF4_CLASSIC: cmode = NC_NETCDF4 | NC_CLASSIC_MODEL; break; + case NC_FORMAT_64BIT_OFFSET: cmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + default: cmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + } + + size_t blksz = 65536; + int ncid; + if (nc__create(outfile, cmode, 0, &blksz, &ncid) != NC_NOERR) { + fprintf(stderr,"decompress: cannot create '%s'.\n", outfile); + free(input); return 1; + } + + /* ---- Define dimensions (skip compressed, include expanded ones) ---- */ + int ndims; + nc_inq_ndims(ref, &ndims); + for (int d = 0; d < ndims; d++) { + char dname[NC_MAX_NAME+1]; size_t dlen; bool is_unlim; + ncu_inq_dim(ref, d, dname, &dlen, &is_unlim); + + /* Skip compressed dimensions. */ + int varid, attid; + if (nc_inq_varid(ref, dname, &varid) == NC_NOERR && + nc_inq_attid(ref, varid, "compress", &attid) == NC_NOERR) + continue; + + int dummy; + nc_def_dim(ncid, dname, is_unlim ? NC_UNLIMITED : dlen, &dummy); + } + + /* ---- Define variables, replacing compressed dims with expanded dims ---- */ + int nvars; + nc_inq_nvars(ref, &nvars); + for (int v = 0; v < nvars; v++) { + char vname[NC_MAX_NAME+1]; + nc_type xtype; + int vndims, dimids[NC_MAX_VARS], natts; + nc_inq_var(ref, v, vname, &xtype, &vndims, dimids, &natts); + if (debug) printf("decompress: defining '%s'\n", vname); + + /* Is this a compressed dim variable? */ + bool is_dim = false, is_comp = false; + int dummy_did; + if (nc_inq_dimid(ref, vname, &dummy_did) == NC_NOERR) is_dim = true; + int attid; + if (is_dim && nc_inq_varid(ref,vname,&dummy_did)==NC_NOERR && + nc_inq_attid(ref,dummy_did,"compress",&attid)==NC_NOERR) is_comp = true; + if (is_dim && is_comp) continue; /* Don't define compressed dim vars. */ + + /* Expand any compressed dims in the variable's dim list. */ + int new_ndims = 0; + int new_dimids[NC_MAX_DIMS]; + for (int d = 0; d < vndims; d++) { + char dname[NC_MAX_NAME+1]; + nc_inq_dimname(ref, dimids[d], dname); + /* Is this dim compressed? */ + int cdims; int cdimids[NC_MAX_DIMS]; size_t cdimlens[NC_MAX_DIMS]; + if (compress_inq_dim(ref, dimids[d], &cdims, cdimids, cdimlens) == NC_NOERR) { + /* Replace with expanded dims. */ + for (int k = 0; k < cdims; k++) { + char cdname[NC_MAX_NAME+1]; + nc_inq_dimname(ref, cdimids[k], cdname); + int new_did; + if (nc_inq_dimid(ncid, cdname, &new_did) == NC_NOERR) + new_dimids[new_ndims++] = new_did; + } + } else { + int new_did; + if (nc_inq_dimid(ncid, dname, &new_did) == NC_NOERR) + new_dimids[new_ndims++] = new_did; + } + } + + int new_varid; + nc_def_var(ncid, vname, xtype, new_ndims, new_dimids, &new_varid); + for (int n = 0; n < natts; n++) { + char attname[NC_MAX_NAME+1]; + nc_inq_attname(ref, v, n, attname); + nc_copy_att(ref, v, attname, ncid, new_varid); + } + + /* Optionally add missing_value if compressed and not present. */ + if (add_missing && is_comp) { + int attid2; + if (nc_inq_attid(ref, v, "missing_value", &attid2) != NC_NOERR && + nc_inq_attid(ref, v, "_FillValue", &attid2) != NC_NOERR) { + double mv; + switch (xtype) { + case NC_DOUBLE: mv = NC_FILL_DOUBLE; break; + case NC_FLOAT: mv = NC_FILL_FLOAT; break; + default: mv = NC_FILL_INT; break; + } + nc_put_att_double(ncid, new_varid, "missing_value", xtype, 1, &mv); + } + } + } + + /* Copy global attributes. */ + ncu_copy_global_atts(ref, ncid); + nc__enddef(ncid, 16384, 4, 0, 4); + + /* Extend record dimension if needed. */ + for (int v = 0; v < nvars; v++) { + size_t cdimlens[NC_MAX_DIMS]; + char vname[NC_MAX_NAME+1]; + nc_type xtype; + int vndims, dimids[NC_MAX_DIMS], natts; + ncu_inq_var(ref, v, vname, &xtype, &vndims, dimids, cdimlens, &natts, NULL, NULL, NULL, NULL, NULL); + /* If has record dim, write a dummy int to extend. */ + int unlimdim_out; + nc_inq_unlimdim(ncid, &unlimdim_out); + if (unlimdim_out < 0) break; + char unlim_name[NC_MAX_NAME+1]; + nc_inq_dimname(ncid, unlimdim_out, unlim_name); + int ov; + if (nc_inq_varid(ncid, vname, &ov) != NC_NOERR) continue; + int ov_ndims, ov_dimids[NC_MAX_DIMS]; + nc_inq_varndims(ncid, ov, &ov_ndims); + nc_inq_vardimid(ncid, ov, ov_dimids); + bool ov_has_rec = false; + for (int d = 0; d < ov_ndims; d++) + if (ov_dimids[d] == unlimdim_out) { ov_has_rec = true; break; } + if (!ov_has_rec) continue; + + /* Get max dimlens of this var. */ + size_t end_idx[NC_MAX_DIMS]; + for (int d = 0; d < ov_ndims; d++) { + size_t dlen; nc_inq_dimlen(ref, dimids[d < vndims ? d : 0], &dlen); + end_idx[d] = dlen > 0 ? dlen - 1 : 0; + } + nc_put_var1_int(ncid, ov, end_idx, (const int[]){0}); + break; + } + + /* ---- Gather compressed data and write to output ---- */ + int out_nvars; + nc_inq_nvars(ncid, &out_nvars); + for (int v = 0; v < out_nvars; v++) { + char vname[NC_MAX_NAME+1]; + nc_type xtype; + size_t vsize; + ncu_inq_var(ncid, v, vname, &xtype, NULL, NULL, NULL, NULL, NULL, NULL, &vsize, NULL, NULL); + if (debug) printf("decompress: processing '%s'\n", vname); + + double *buf = XMALLOC(double, vsize > 0 ? vsize : 1); + bool *mask = XMALLOC(bool, vsize > 0 ? vsize : 1); + + /* Determine missing value. */ + double missing = NC_FILL_DOUBLE; + ncu_inq_var(ncid, v, NULL, &xtype, NULL, NULL, NULL, NULL, NULL, NULL, &vsize, NULL, NULL); + double ocean_val = 0.0; + bool do_ocean = (nc_get_att_double(ncid, v, "ocean_fillvalue", &ocean_val) == NC_NOERR); + if (nc_get_att_double(ncid, v, "missing_value", &missing) != NC_NOERR) + if (nc_get_att_double(ncid, v, "_FillValue", &missing) != NC_NOERR) { + switch (xtype) { + case NC_DOUBLE: missing = NC_FILL_DOUBLE; break; + case NC_FLOAT: missing = NC_FILL_FLOAT; break; + default: missing = NC_FILL_INT; break; + } + } + for (size_t k = 0; k < vsize; k++) buf[k] = missing; + memset(mask, 0, vsize * sizeof(bool)); + + /* Read from all input files. */ + for (int fi = 0; fi < ninfiles; fi++) { + int iv; + if (nc_inq_varid(input[fi], vname, &iv) != NC_NOERR) continue; + compress_get_var_double(input[fi], iv, buf, mask, vsize); + } + + if (do_ocean) + for (size_t k = 0; k < vsize; k++) + if (!mask[k]) buf[k] = ocean_val; + + /* Write to output. */ + nc_put_var_double(ncid, v, buf); + free(buf); free(mask); + } + + nc_sync(ncid); nc_close(ncid); + for (int i = 0; i < ninfiles; i++) nc_close(input[i]); + free(input); + if (infiles) free((void*)infiles); + return 0; +} diff --git a/src/mpp-disttool/cmd_decompress.h b/src/mpp-disttool/cmd_decompress.h new file mode 100644 index 00000000..db29cec7 --- /dev/null +++ b/src/mpp-disttool/cmd_decompress.h @@ -0,0 +1,26 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_CMD_DECOMPRESS_H +#define MPPDISTTOOL_CMD_DECOMPRESS_H + +/* Entry point: "mppdisttool decompress infile... outfile" */ +int cmd_decompress(int argc, char **argv); + +#endif /* MPPDISTTOOL_CMD_DECOMPRESS_H */ diff --git a/src/mpp-disttool/cmd_scatter.c b/src/mpp-disttool/cmd_scatter.c new file mode 100644 index 00000000..069ad071 --- /dev/null +++ b/src/mpp-disttool/cmd_scatter.c @@ -0,0 +1,493 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +/* + * Two scatter paths: + * scatter_mpp – domain decomposition (from mppncscatter.c + domain.c) + * scatter_land – CF compressed-by-gathering scatter (from scatter-ncc.F90) + */ +#define _POSIX_C_SOURCE 200809L +#include +#include +#include +#include +#include +#include "cmd_scatter.h" +#include "domain.h" +#include "nc_utils.h" +#include "compress.h" +#include "xmalloc.h" +#include "strlist.h" + +/* ================================================================== */ +/* scatter_mpp – wraps domain.c logic (ported from mppncscatter) */ +/* ================================================================== */ +static int +scatter_mpp(scatter_opts_t *opts) +{ + int nc; + if (nc_open(opts->filein, NC_NOWRITE, &nc) != NC_NOERR) { + fprintf(stderr,"scatter: cannot open '%s'\n", opts->filein); + return -1; + } + + int ndims, nvars, ngatts, unlimdimid; + nc_inq(nc, &ndims, &nvars, &ngatts, &unlimdimid); + + int format_flag = NC_FORMAT_CLASSIC; + nc_inq_format(nc, &format_flag); + int create_flags; + switch (format_flag) { + case NC_FORMAT_64BIT_OFFSET: create_flags = NC_64BIT_OFFSET; break; + case NC_FORMAT_NETCDF4: create_flags = NC_NETCDF4; break; + case NC_FORMAT_NETCDF4_CLASSIC: create_flags = NC_NETCDF4 | NC_CLASSIC_MODEL; break; + default: create_flags = 0; break; + } + + /* Strip directory prefix from filein for output names. */ + const char *prefix = opts->filein; + const char *p = strrchr(prefix, '/'); + if (p) prefix = p + 1; + + /* Build output filename format string. */ + char outnameformat[512]; + if (opts->prefix && strlen(opts->prefix) > 0) + snprintf(outnameformat, sizeof(outnameformat), "%s/%%s.%%0%dd", + opts->prefix, opts->width); + else + snprintf(outnameformat, sizeof(outnameformat), "%%s.%%0%dd", opts->width); + + ScatterDim *scatterdims[NC_MAX_DIMS]; + scatter_dims(nc, ndims, nvars, scatterdims, opts); + + int nfiles = scatter_get_num_files(opts); + int *ncids = XMALLOC(int, nfiles); + + for (int i = 0; i < nfiles; i++) { + char output[512]; + snprintf(output, sizeof(output), outnameformat, prefix, i + opts->start); + + if (opts->verbose) printf("Info: Creating '%s'\n", output); + + if (!opts->dryrun) { + int nc_err = nc_create(output, NC_CLOBBER | create_flags, &ncids[i]); + if (nc_err != NC_NOERR) { + fprintf(stderr, "scatter: cannot create '%s': %s\n", + output, nc_strerror(nc_err)); + /* Clean up already-created files */ + for (int j = 0; j < i; j++) { + nc_close(ncids[j]); + } + free(ncids); + nc_close(nc); + return -1; + } + int dummy; + nc_set_fill(ncids[i], NC_NOFILL, &dummy); + nc_put_att_int(ncids[i], NC_GLOBAL, "NumFilesInSet", NC_INT, 1, &nfiles); + } + } + + /* Copy global attributes. */ + if (!opts->dryrun) { + char attname[NC_MAX_NAME + 1]; + for (int j = 0; j < ngatts; j++) { + nc_inq_attname(nc, NC_GLOBAL, j, attname); + for (int i = 0; i < nfiles; i++) + nc_copy_att(nc, NC_GLOBAL, attname, ncids[i], NC_GLOBAL); + } + } + + scatter_def_dim(nc, ncids, ndims, scatterdims, opts); + scatter_def_var(nc, ncids, nvars, ndims, scatterdims, opts); + + if (!opts->dryrun) + for (int i = 0; i < nfiles; i++) nc_enddef(ncids[i]); + + scatter_put_var(nc, ncids, ndims, nvars, scatterdims, opts); + + nc_close(nc); + if (!opts->dryrun) + for (int i = 0; i < nfiles; i++) { nc_sync(ncids[i]); nc_close(ncids[i]); } + + scatter_dims_free(scatterdims, ndims); + free(ncids); + return 0; +} + +/* ================================================================== */ +/* scatter_land – port of scatter-ncc.F90 */ +/* ================================================================== */ +static int +scatter_land(const char *infile, int npex, int npey) +{ + int in_ncid; + if (nc_open(infile, NC_NOWRITE, &in_ncid) != NC_NOERR) { + fprintf(stderr,"scatter --land: cannot open '%s'\n", infile); + return 1; + } + + int in_fmt; + nc_inq_format(in_ncid, &in_fmt); + int cmode; + switch (in_fmt) { + case NC_FORMAT_NETCDF4: cmode = NC_NETCDF4; break; + case NC_FORMAT_NETCDF4_CLASSIC: cmode = NC_NETCDF4 | NC_CLASSIC_MODEL; break; + case NC_FORMAT_64BIT_OFFSET: cmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + default: cmode = NC_CLOBBER | NC_64BIT_OFFSET; break; + } + + int nfiles_out = npex * npey; + + /* Find lon/lat/zfull dimensions in input. */ + int ndims; + nc_inq_ndims(in_ncid, &ndims); + + int nlon = 0, nlat = 0, nz = 1; + for (int d = 0; d < ndims; d++) { + char dname[NC_MAX_NAME + 1]; size_t dlen; + nc_inq_dim(in_ncid, d, dname, &dlen); + if (strcmp(dname,"lon") == 0) nlon = (int)dlen; + if (strcmp(dname,"lat") == 0) nlat = (int)dlen; + if (strcmp(dname,"zfull") == 0) nz = (int)dlen; + } + + if (nlon == 0 || nlat == 0) { + fprintf(stderr,"scatter --land: input lacks 'lon' or 'lat' dimension.\n"); + nc_close(in_ncid); return 1; + } + if (nlon % npex != 0 || nlat % npey != 0) { + fprintf(stderr,"scatter --land: domain not evenly divisible by npex/npey.\n"); + nc_close(in_ncid); return 1; + } + int nlon_local = nlon / npex; + int nlat_local = nlat / npey; + + /* Create output files. */ + int *out_ncids = XMALLOC(int, nfiles_out); + int *dimlen_list = XMALLOC(int, nfiles_out); + size_t blksz = 65536; + for (int n = 0; n < nfiles_out; n++) { + char outname[2048]; + snprintf(outname, sizeof(outname), "%s.%04d", infile, n); + int nc_err = nc__create(outname, cmode, 0, &blksz, &out_ncids[n]); + if (nc_err != NC_NOERR) { + fprintf(stderr, "scatter --land: cannot create '%s': %s\n", + outname, nc_strerror(nc_err)); + /* Clean up already-created files */ + for (int j = 0; j < n; j++) { + nc_close(out_ncids[j]); + } + free(out_ncids); + free(dimlen_list); + nc_close(in_ncid); + return 1; + } + } + + /* Count compressed dimension sizes per output tile. */ + + for (int d = 0; d < ndims; d++) { + char dname[NC_MAX_NAME + 1]; size_t dlen; bool is_unlim; + ncu_inq_dim(in_ncid, d, dname, &dlen, &is_unlim); + + int is_compressed = 0; + int varid, attid; + if (nc_inq_varid(in_ncid, dname, &varid) == NC_NOERR && + nc_inq_attid(in_ncid, varid, "compress", &attid) == NC_NOERR) + is_compressed = 1; + + memset(dimlen_list, 0, nfiles_out * sizeof(int)); + + if (is_compressed) { + /* Use compress_get_var_double to get full uncompressed data with mask. */ + int cndims; size_t cdimlens[NC_MAX_DIMS]; + compress_inq_dim(in_ncid, d, &cndims, NULL, cdimlens); + size_t full_size = 1; + for (int k = 0; k < cndims; k++) full_size *= cdimlens[k]; + + double *buf = XMALLOC(double, full_size); + bool *mask = XMALLOC(bool, full_size); + memset(mask, 0, full_size * sizeof(bool)); + compress_get_var_double(in_ncid, varid, buf, mask, full_size); + + int npts = nlon * nlat; + for (size_t l = 0; l < full_size; l++) { + if (!mask[l]) continue; + int ij = (int)(l % (size_t)npts); + int i = ij % nlon; + int j = ij / nlon; + int ii = i / nlon_local; + int jj = j / nlat_local; + int nn = jj * npex + ii; + if (nn < nfiles_out) dimlen_list[nn]++; + } + free(buf); free(mask); + } + + size_t out_dlen = (size_t)dlen; + for (int n = 0; n < nfiles_out; n++) { + int dummy; + if (is_compressed) out_dlen = (size_t)(dimlen_list[n] > 0 ? dimlen_list[n] : 1); + nc_def_dim(out_ncids[n], dname, is_unlim ? NC_UNLIMITED : out_dlen, &dummy); + } + } + + /* Clone variable definitions and global attributes. */ + int nvars; + nc_inq_nvars(in_ncid, &nvars); + for (int v = 0; v < nvars; v++) + for (int n = 0; n < nfiles_out; n++) + ncu_clone_var(in_ncid, v, out_ncids[n]); + + int ngatts; + nc_inq_natts(in_ncid, &ngatts); + for (int i = 0; i < ngatts; i++) { + char aname[NC_MAX_NAME + 1]; + nc_inq_attname(in_ncid, NC_GLOBAL, i, aname); + for (int n = 0; n < nfiles_out; n++) + nc_copy_att(in_ncid, NC_GLOBAL, aname, out_ncids[n], NC_GLOBAL); + } + + for (int n = 0; n < nfiles_out; n++) + nc__enddef(out_ncids[n], 16384, 4, 0, 4); + + /* Query number of records. */ + int unlimdimid, nrec = 1; + if (nc_inq_unlimdim(in_ncid, &unlimdimid) == NC_NOERR && unlimdimid >= 0) { + size_t nrec_sz; + nc_inq_dimlen(in_ncid, unlimdimid, &nrec_sz); + nrec = (int)nrec_sz; + } + + int npts = nlon * nlat; + int npts_local = nlon_local * nlat_local; + + /* For each variable and each record, scatter data. */ + for (int tlev = 0; tlev < nrec; tlev++) { + for (int v = 0; v < nvars; v++) { + char vname[NC_MAX_NAME + 1]; + nc_type xtype; + int vndims, dimids[NC_MAX_DIMS]; + size_t dimlens[NC_MAX_DIMS]; + bool has_records, is_compressed_var; + ncu_inq_var(in_ncid, v, vname, &xtype, &vndims, dimids, dimlens, + NULL, NULL, &has_records, NULL, NULL, NULL); + + if (!has_records && tlev > 0) continue; + + /* Check if first dim is compressed. */ + is_compressed_var = false; + if (vndims > 0) { + char d0name[NC_MAX_NAME+1]; + int varid0, attid0; + nc_inq_dimname(in_ncid, dimids[0], d0name); + if (nc_inq_varid(in_ncid,d0name,&varid0)==NC_NOERR && + nc_inq_attid(in_ncid,varid0,"compress",&attid0)==NC_NOERR) + is_compressed_var = true; + } + + /* Get uncompressed variable size. */ + int cndims_v; size_t cdimlens_v[NC_MAX_DIMS]; + compress_inq_dim(in_ncid, dimids[0], &cndims_v, NULL, cdimlens_v); + + size_t vsize; + if (is_compressed_var) { + /* Use first compressed dim's uncompressed dims + rest. */ + vsize = 1; + for (int k = 0; k < cndims_v; k++) vsize *= cdimlens_v[k]; + for (int d = 1; d < vndims; d++) if (dimids[d] != unlimdimid) vsize *= dimlens[d]; + } else { + vsize = 1; + for (int d = 0; d < vndims; d++) if (dimids[d] != unlimdimid) vsize *= dimlens[d]; + } + + double *buf = XMALLOC(double, vsize > 0 ? vsize : 1); + bool *mask = XMALLOC(bool, vsize > 0 ? vsize : 1); + memset(mask, 0, vsize * sizeof(bool)); + + /* Read variable from input. */ + if (is_compressed_var) { + /* Determine level dimension. */ + int k_dim = -1; + for (int d = 0; d < vndims; d++) { + char ddn[NC_MAX_NAME+1]; + nc_inq_dimname(in_ncid, dimids[d], ddn); + if (strcmp(ddn,"zfull")==0) k_dim = d; + } + int nz_var = (k_dim >= 0) ? nz : 1; + (void)nz_var; /* recsize not used in this path */ + + size_t st[NC_MAX_DIMS], cnt[NC_MAX_DIMS]; + st[0] = (size_t)tlev; cnt[0] = 1; + for (int d = 1; d < vndims; d++) { st[d]=1; cnt[d]=1; } + /* Read all at once. */ + compress_get_var_double(in_ncid, v, buf, mask, vsize); + } else { + if (has_records) { + /* 1D or scalar record variable. */ + size_t st1[1] = {(size_t)tlev}, cnt1[1] = {1}; + nc_get_vara_double(in_ncid, v, st1, cnt1, buf); + } else { + nc_get_var_double(in_ncid, v, buf); + } + } + + /* Scatter to output tiles. */ + for (int n = 0; n < nfiles_out; n++) { + int ii_tile = n % npex; + int jj_tile = n / npex; + + if (is_compressed_var) { + /* Collect points belonging to this tile. */ + double *tbuf = XMALLOC(double, (size_t)npts_local); + bool *tmask = XMALLOC(bool, (size_t)npts_local); + memset(tbuf, 0, npts_local * sizeof(double)); + memset(tmask, 0, npts_local * sizeof(bool)); + int tcount = 0; + + for (int l = 0; l < (int)vsize; l++) { + if (!mask[l]) continue; + int ij = l % npts; + int gi = ij % nlon; + int gj = ij / nlon; + int li_tile = gi / nlon_local; + int lj_tile = gj / nlat_local; + if (li_tile != ii_tile || lj_tile != jj_tile) continue; + int li_local = gi % nlon_local; + int lj_local = gj % nlat_local; + int ll = lj_local * nlon_local + li_local; + if (ll < npts_local) { tbuf[ll] = buf[l]; tmask[ll] = true; tcount++; } + } + if (tcount > 0) { + /* Write packed data to output tile. */ + double *packed = XMALLOC(double, tcount); + int pk = 0; + for (int l = 0; l < npts_local; l++) + if (tmask[l]) packed[pk++] = tbuf[l]; + + int ov; + if (nc_inq_varid(out_ncids[n], vname, &ov) == NC_NOERR) { + size_t st[4]={0,0,0,0}, cnt[4]={0,1,1,1}; + if (has_records) st[unlimdimid >= 0 ? 0 : 0] = (size_t)tlev; + cnt[0] = (size_t)tcount; + nc_put_vara_double(out_ncids[n], ov, st, cnt, packed); + } + free(packed); + } + free(tbuf); free(tmask); + } else { + /* Non-compressed variable: write same value to all tiles. */ + int ov; + if (nc_inq_varid(out_ncids[n], vname, &ov) == NC_NOERR) { + size_t st[1] = {(size_t)tlev}; + size_t cnt[1] = {1}; + if (has_records) + nc_put_vara_double(out_ncids[n], ov, st, cnt, buf); + else + nc_put_var_double(out_ncids[n], ov, buf); + } + } + } + free(buf); free(mask); + } + } + + nc_close(in_ncid); + for (int n = 0; n < nfiles_out; n++) { nc_sync(out_ncids[n]); nc_close(out_ncids[n]); } + free(out_ncids); free(dimlen_list); + return 0; +} + +/* ================================================================== */ +/* cmd_scatter – argument parsing + dispatch */ +/* ================================================================== */ +static void usage_scatter(void) { + printf("Usage: mppdisttool scatter -x N -y N [-i N] [-j N] [-p prefix]\n" + " [-s N] [-w N] [-X dims] [-Y dims] [-n] [-v] infile\n" + " mppdisttool scatter --land -i ndiv_x -j ndiv_y infile\n"); +} + +int +cmd_scatter(int argc, char **argv) +{ + scatter_opts_t opts; + memset(&opts, 0, sizeof(opts)); + opts.width = 4; + + int land_mode = 0; + int npex_land = 0, npey_land = 0; + + int dummy_n; + newstringlist(&opts.xdims, &dummy_n, NC_MAX_DIMS); + newstringlist(&opts.ydims, &dummy_n, NC_MAX_DIMS); + + for (int i = 1; i < argc; i++) { + if (strcmp(argv[i],"--land") == 0) land_mode = 1; + else if (strcmp(argv[i],"-n") == 0) opts.dryrun = 1; + else if (strcmp(argv[i],"-v") == 0 || strcmp(argv[i],"-V")==0) opts.verbose = 1; + else if ((strcmp(argv[i],"-x") == 0) && i+1 < argc) opts.nx = atoi(argv[++i]); + else if ((strcmp(argv[i],"-y") == 0) && i+1 < argc) opts.ny = atoi(argv[++i]); + else if ((strcmp(argv[i],"-i") == 0) && i+1 < argc) { + int v = atoi(argv[++i]); + if (land_mode) npex_land = v; else opts.nxio = v; + } else if ((strcmp(argv[i],"-j") == 0) && i+1 < argc) { + int v = atoi(argv[++i]); + if (land_mode) npey_land = v; else opts.nyio = v; + } else if ((strcmp(argv[i],"-p") == 0) && i+1 < argc) { + opts.prefix = argv[++i]; + } else if ((strcmp(argv[i],"-s") == 0) && i+1 < argc) { + opts.start = atoi(argv[++i]); + } else if ((strcmp(argv[i],"-w") == 0) && i+1 < argc) { + opts.width = atoi(argv[++i]); + } else if ((strcmp(argv[i],"-X") == 0) && i+1 < argc) { + getstringlist(argv[++i], &opts.xdims, &opts.xdims_len); + } else if ((strcmp(argv[i],"-Y") == 0) && i+1 < argc) { + getstringlist(argv[++i], &opts.ydims, &opts.ydims_len); + } else if (strcmp(argv[i],"--help") == 0 || strcmp(argv[i],"-h") == 0) { + usage_scatter(); return 0; + } else if (argv[i][0] != '-') { + opts.filein = argv[i]; + } + } + + if (!opts.filein) { + fprintf(stderr,"scatter: missing input file.\n"); + usage_scatter(); return 1; + } + + int ret; + if (land_mode) { + if (npex_land <= 0 || npey_land <= 0) { + fprintf(stderr,"scatter --land: -i ndiv_x and -j ndiv_y required.\n"); + return 1; + } + ret = scatter_land(opts.filein, npex_land, npey_land); + } else { + if (opts.nx <= 0 || opts.ny <= 0) { + fprintf(stderr,"scatter: -x N and -y N required.\n"); + usage_scatter(); return 1; + } + ret = scatter_mpp(&opts); + } + + freestringlist(&opts.xdims, NC_MAX_DIMS); + freestringlist(&opts.ydims, NC_MAX_DIMS); + return ret; +} diff --git a/src/mpp-disttool/cmd_scatter.h b/src/mpp-disttool/cmd_scatter.h new file mode 100644 index 00000000..3ee1c554 --- /dev/null +++ b/src/mpp-disttool/cmd_scatter.h @@ -0,0 +1,26 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_CMD_SCATTER_H +#define MPPDISTTOOL_CMD_SCATTER_H + +/* Entry point: "mppdisttool scatter [opts] infile" */ +int cmd_scatter(int argc, char **argv); + +#endif /* MPPDISTTOOL_CMD_SCATTER_H */ diff --git a/src/mpp-disttool/compress.c b/src/mpp-disttool/compress.c new file mode 100644 index 00000000..822eacf0 --- /dev/null +++ b/src/mpp-disttool/compress.c @@ -0,0 +1,374 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +/* + * C port of lib/libnfu/nfu_compress.F90 and the rank_ascending / + * mergerank / merge routines from src/land_utils/combine-ncc.F90. + */ +#define _POSIX_C_SOURCE 200809L +#include +#include +#include +#include +#include +#include "compress.h" +#include "nc_utils.h" +#include "xmalloc.h" + +#define NC_CHECK(call) do { \ + int _err = (call); \ + if (_err != NC_NOERR) return _err; \ +} while (0) + +/* ------------------------------------------------------------------ */ +/* compress_inq_dim – check if a dimension is compressed by gathering */ +/* ------------------------------------------------------------------ */ +int +compress_inq_dim(int ncid, int dimid, + int *ndims_out, + int *dimids_out, + size_t *dimlens_out) +{ + char dimname[NC_MAX_NAME + 1]; + char *compress_att = NULL; + int varid; + size_t att_len; + nc_type att_type; + + NC_CHECK(nc_inq_dimname(ncid, dimid, dimname)); + + /* Does the dimension variable exist and have a "compress" attribute? */ + if (nc_inq_varid(ncid, dimname, &varid) != NC_NOERR) + return NC_ENOTATT; + if (nc_inq_att(ncid, varid, "compress", &att_type, &att_len) != NC_NOERR) + return NC_ENOTATT; + if (att_type != NC_CHAR || att_len == 0) + return NC_ENOTATT; + + /* Dynamically allocate buffer based on actual attribute length */ + compress_att = (char *)malloc(att_len + 1); + if (!compress_att) { + fprintf(stderr, "Cannot allocate memory for compress attribute\n"); + return NC_ENOMEM; + } + + int err = nc_get_att_text(ncid, varid, "compress", compress_att); + if (err != NC_NOERR) { + free(compress_att); + return err; + } + compress_att[att_len] = '\0'; + + /* Parse space-separated dimension names in compress_att. + The Fortran code scans from the end, producing them in reverse order. + We do forward scanning, which gives the same set. */ + int n = 0; + char *saveptr = NULL; + char *tok = strtok_r(compress_att, " \t", &saveptr); + while (tok) { + if (ndims_out || dimids_out || dimlens_out) { + if (dimids_out || dimlens_out) { + int did; + int nc_err = nc_inq_dimid(ncid, tok, &did); + if (nc_err != NC_NOERR) { + free(compress_att); + return nc_err; + } + if (dimids_out) dimids_out[n] = did; + if (dimlens_out) { + size_t dlen; + nc_err = nc_inq_dimlen(ncid, did, &dlen); + if (nc_err != NC_NOERR) { + free(compress_att); + return nc_err; + } + dimlens_out[n] = dlen; + } + } + n++; + } + tok = strtok_r(NULL, " \t", &saveptr); + } + if (ndims_out) *ndims_out = n; + free(compress_att); + return NC_NOERR; +} + +/* ------------------------------------------------------------------ */ +/* compress_build_diminfo */ +/* ------------------------------------------------------------------ */ +int +compress_build_diminfo(int ncid, int varid, + diminfo_t *diminfo, int *ndims_out, + size_t *varsize_out) +{ + int ndims, dimids[NC_MAX_VAR_DIMS]; + NC_CHECK(nc_inq_varndims(ncid, varid, &ndims)); + NC_CHECK(nc_inq_vardimid(ncid, varid, dimids)); + + if (ndims_out) *ndims_out = ndims; + + int stride = 1; + for (int i = 0; i < ndims; i++) { + size_t dimlen; + int cdims, cdimids[NC_MAX_DIMS]; + size_t cdimlens[NC_MAX_DIMS]; + + NC_CHECK(nc_inq_dimlen(ncid, dimids[i], &dimlen)); + diminfo[i].length = (int)dimlen; + diminfo[i].idx = NULL; + + int length_for_stride; + if (compress_inq_dim(ncid, dimids[i], + &cdims, cdimids, cdimlens) == NC_NOERR) { + /* Compressed dimension: load index array. */ + diminfo[i].idx = XMALLOC(int, dimlen); + char cname[NC_MAX_NAME + 1]; + NC_CHECK(nc_inq_dimname(ncid, dimids[i], cname)); + NC_CHECK(ncu_get_var_int(ncid, cname, diminfo[i].idx)); + /* Stride contribution = product of uncompressed dim sizes. */ + length_for_stride = 1; + for (int k = 0; k < cdims; k++) + length_for_stride *= (int)cdimlens[k]; + } else { + length_for_stride = diminfo[i].length; + } + diminfo[i].stride = stride; + stride *= length_for_stride; + } + if (varsize_out) { + /* Total compressed size = product of compressed lengths. */ + size_t vsz = 1; + for (int i = 0; i < ndims; i++) + vsz *= (size_t)diminfo[i].length; + *varsize_out = vsz; + } + return NC_NOERR; +} + +void +compress_free_diminfo(diminfo_t *diminfo, int ndims) +{ + for (int i = 0; i < ndims; i++) { + if (diminfo[i].idx) { + free(diminfo[i].idx); + diminfo[i].idx = NULL; + } + } +} + +/* ------------------------------------------------------------------ */ +/* Internal: scatter compressed buffer into uncompressed output array */ +/* ------------------------------------------------------------------ */ +static void +scatter_to_output(const double *buffer, size_t bufsize, + double *data, bool *mask, + const diminfo_t *diminfo, int ndims) +{ + int local_idx[NC_MAX_VAR_DIMS]; + memset(local_idx, 0, ndims * sizeof(int)); + + for (size_t i = 0; i < bufsize; i++) { + /* Compute 0-based flat output index. */ + int ii = 0; + bool valid = true; + for (int n = 0; n < ndims; n++) { + if (diminfo[n].idx != NULL) { + int v = diminfo[n].idx[local_idx[n]]; + if (v < 0) { valid = false; break; } + ii += v * diminfo[n].stride; + } else { + ii += local_idx[n] * diminfo[n].stride; + } + } + if (valid) { + data[ii] = buffer[i]; + if (mask) mask[ii] = true; + } + /* Increment multi-dimensional counter. */ + for (int n = 0; n < ndims; n++) { + local_idx[n]++; + if (local_idx[n] < diminfo[n].length) break; + local_idx[n] = 0; + } + } +} + +/* ------------------------------------------------------------------ */ +/* compress_get_var_double */ +/* ------------------------------------------------------------------ */ +int +compress_get_var_double(int ncid, int varid, + double *data, bool *mask, + size_t uncompressed_size) +{ + return compress_get_var_double_range(ncid, varid, NULL, NULL, + data, mask, uncompressed_size); +} + +int +compress_get_var_double_range(int ncid, int varid, + const size_t *start, const size_t *count, + double *data, bool *mask, + size_t uncompressed_size) +{ + diminfo_t diminfo[NC_MAX_VAR_DIMS]; + int ndims; + size_t varsize; + + NC_CHECK(compress_build_diminfo(ncid, varid, diminfo, &ndims, &varsize)); + + size_t bufsize = (start && count) ? count[0] : varsize; + double *buffer = XMALLOC(double, bufsize); + + int err; + if (start && count) { + err = nc_get_vara_double(ncid, varid, start, count, buffer); + /* For sub-range reads, treat as 1D. */ + diminfo_t di1[1]; + di1[0].idx = diminfo[0].idx; + di1[0].length = (int)count[0]; + di1[0].stride = diminfo[0].stride; + if (err == NC_NOERR) + scatter_to_output(buffer, bufsize, data, mask, di1, 1); + /* Don't free di1[0].idx since it's borrowed from diminfo[0]. */ + } else { + err = nc_get_var_double(ncid, varid, buffer); + if (err == NC_NOERR) + scatter_to_output(buffer, bufsize, data, mask, diminfo, ndims); + } + + free(buffer); + compress_free_diminfo(diminfo, ndims); + return err; +} + +/* ------------------------------------------------------------------ */ +/* compress_put_var_double */ +/* ------------------------------------------------------------------ */ +int +compress_put_var_double(int ncid, int varid, + const double *src, size_t uncompressed_size) +{ + diminfo_t diminfo[NC_MAX_VAR_DIMS]; + int ndims; + size_t varsize; + + NC_CHECK(compress_build_diminfo(ncid, varid, diminfo, &ndims, &varsize)); + + double *buffer = XMALLOC(double, varsize); + int local_idx[NC_MAX_VAR_DIMS]; + memset(local_idx, 0, ndims * sizeof(int)); + + for (size_t i = 0; i < varsize; i++) { + int ii = 0; + for (int n = 0; n < ndims; n++) { + if (diminfo[n].idx != NULL) + ii += diminfo[n].idx[local_idx[n]] * diminfo[n].stride; + else + ii += local_idx[n] * diminfo[n].stride; + } + buffer[i] = src[ii]; + for (int n = 0; n < ndims; n++) { + local_idx[n]++; + if (local_idx[n] < diminfo[n].length) break; + local_idx[n] = 0; + } + } + + int err = nc_put_var_double(ncid, varid, buffer); + free(buffer); + compress_free_diminfo(diminfo, ndims); + return err; +} + +int +compress_put_vara_double(int ncid, int varid, + const size_t *start, const size_t *count, + const double *buf) +{ + return nc_put_vara_double(ncid, varid, start, count, buf); +} + +/* ------------------------------------------------------------------ */ +/* Merge sort – ported from combine-ncc.F90:358-430 */ +/* ------------------------------------------------------------------ */ +static void +merge_arrays(const int *x, + const int *a, int na, + const int *b, int nb, + int *c, int nc_size) +{ + int i = 0, j = 0, k = 0; + while (i < na && j < nb) { + if (x[a[i]] <= x[b[j]]) + c[k++] = a[i++]; + else + c[k++] = b[j++]; + } + while (i < na) + c[k++] = a[i++]; + (void)nc_size; +} + +static void +mergerank(const int *x, int *a, int n, int *t) +{ + if (n < 2) return; + if (n == 2) { + if (x[a[0]] > x[a[1]]) { int v = a[0]; a[0] = a[1]; a[1] = v; } + return; + } + int na = (n + 1) / 2; + int nb = n - na; + mergerank(x, a, na, t); + mergerank(x, a + na, nb, t); + if (x[a[na - 1]] > x[a[na]]) { + memcpy(t, a, na * sizeof(int)); + merge_arrays(x, t, na, a + na, nb, a, n); + } +} + +void +rank_ascending(const int *x, int *idx, int n) +{ + for (int i = 0; i < n; i++) idx[i] = i; + int *t = XMALLOC(int, (n + 1) / 2); + mergerank(x, idx, n, t); + free(t); +} + +/* ------------------------------------------------------------------ */ +/* file_is_compressed */ +/* ------------------------------------------------------------------ */ +bool +file_is_compressed(int ncid) +{ + int ndims; + if (nc_inq_ndims(ncid, &ndims) != NC_NOERR) return false; + for (int dimid = 0; dimid < ndims; dimid++) { + char name[NC_MAX_NAME + 1]; + int varid, attid; + if (nc_inq_dimname(ncid, dimid, name) != NC_NOERR) continue; + if (nc_inq_varid(ncid, name, &varid) == NC_NOERR && + nc_inq_attid(ncid, varid, "compress", &attid) == NC_NOERR) + return true; + } + return false; +} diff --git a/src/mpp-disttool/compress.h b/src/mpp-disttool/compress.h new file mode 100644 index 00000000..73fa2a1e --- /dev/null +++ b/src/mpp-disttool/compress.h @@ -0,0 +1,95 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_COMPRESS_H +#define MPPDISTTOOL_COMPRESS_H + +#include +#include +#include + +/* Mirrors the Fortran diminfo_type from nfu_compress.F90. + For each dimension of a compressed variable: + idx – mapping from compressed position to flat uncompressed index + (NULL for non-compressed dims; values are 0-based; -1 = invalid) + length – size of this dimension in the compressed (input) array + stride – stride of this dimension in the uncompressed output array */ +typedef struct { + int *idx; + int length; + int stride; +} diminfo_t; + +/* Check whether a dimension is CF "compressed by gathering". + Returns NC_NOERR if the dimension variable has a "compress" attribute. + On success, fills ndims_out, dimids_out[0..ndims-1], dimlens_out[0..ndims-1] + with the expanded (uncompressed) dimension info. + Pass NULL for outputs you don't need. */ +int compress_inq_dim(int ncid, int dimid, + int *ndims_out, + int *dimids_out, /* [NC_MAX_DIMS] */ + size_t *dimlens_out); /* [NC_MAX_DIMS] */ + +/* Initialise a diminfo_t array for a compressed variable. + Fills diminfo[0..ndims-1] based on the variable's dimension list. + Allocates diminfo[i].idx for each compressed dimension. + Call compress_free_diminfo() when done. */ +int compress_build_diminfo(int ncid, int varid, + diminfo_t *diminfo, int *ndims_out, + size_t *varsize_out); + +/* Free allocated idx arrays from compress_build_diminfo(). */ +void compress_free_diminfo(diminfo_t *diminfo, int ndims); + +/* Read a compressed variable into a pre-allocated uncompressed array. + 'data' must hold at least uncompressed_size doubles (initialised by caller). + 'mask' (optional) is set to true at positions that were filled. + Pass NULL for mask if not needed. */ +int compress_get_var_double(int ncid, int varid, + double *data, bool *mask, + size_t uncompressed_size); + +/* Variant that reads only a sub-range along the first dimension. + start[0..ndims-1] and count[0..ndims-1] follow nc_get_vara conventions + but must have count[i]==1 for i>0. Pass NULL to read the whole variable. */ +int compress_get_var_double_range(int ncid, int varid, + const size_t *start, const size_t *count, + double *data, bool *mask, + size_t uncompressed_size); + +/* Write an uncompressed array into a compressed variable. */ +int compress_put_var_double(int ncid, int varid, + const double *src, size_t uncompressed_size); + +/* Partial write: nc_put_vara_double for the compressed variable. */ +int compress_put_vara_double(int ncid, int varid, + const size_t *start, const size_t *count, + const double *buf); + +/* ---- Merge-sort helpers (ported from combine-ncc.F90) ---- */ + +/* Rank array x in ascending order: idx[0..n-1] receives indices of x + in non-decreasing order of x values (0-based). */ +void rank_ascending(const int *x, int *idx, int n); + +/* Check if any variable in the file is compressed by gathering. + Returns true if at least one dimension variable has a "compress" attribute. */ +bool file_is_compressed(int ncid); + +#endif /* MPPDISTTOOL_COMPRESS_H */ diff --git a/src/mpp-disttool/domain.c b/src/mpp-disttool/domain.c new file mode 100644 index 00000000..0796e03e --- /dev/null +++ b/src/mpp-disttool/domain.c @@ -0,0 +1,596 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +/* + * Domain decomposition helpers extracted from src/mpp-ncscatter/mppncscatter.c + * and src/mpp-ncscatter/scatterdim.c. + */ +#include +#include +#include +#include +#include +#include +#include "domain.h" +#include "xmalloc.h" + +/* ------------------------------------------------------------------ */ +/* ScatterDim */ +/* ------------------------------------------------------------------ */ +ScatterDim * +ScatterDim_new(int id, size_t len, const char *name, + scatter_t scatter_type, int ndiv) +{ + ScatterDim *p = XMALLOC(ScatterDim, 1); + p->id = id; + p->len = len; + /* Safe string copy with explicit null termination */ + size_t name_len = strlen(name); + if (name_len > NC_MAX_NAME) name_len = NC_MAX_NAME; + memcpy(p->name, name, name_len); + p->name[name_len] = '\0'; + p->scatter_type = scatter_type; + p->scatter_ndiv = (size_t)ndiv; + p->scatter_start = XMALLOC(size_t, ndiv > 0 ? ndiv : 1); + p->scatter_end = XMALLOC(size_t, ndiv > 0 ? ndiv : 1); + p->scatter_len = XMALLOC(size_t, ndiv > 0 ? ndiv : 1); + return p; +} + +void +ScatterDim_free(ScatterDim *p) +{ + if (!p) return; + XFREE(p->scatter_start); + XFREE(p->scatter_end); + XFREE(p->scatter_len); + free(p); +} + +void +scatter_dims_free(ScatterDim *dims[], int ndims) +{ + for (int i = 0; i < ndims; i++) { + ScatterDim_free(dims[i]); + dims[i] = NULL; + } +} + +/* ------------------------------------------------------------------ */ +/* mpp_compute_extent (verbatim port from mppncscatter.c) */ +/* ------------------------------------------------------------------ */ +#define EVEN(n) (!((n) & 1)) +#define ODD(n) (((n) & 1)) +#define MAX(a,b) ((a)>(b)?(a):(b)) + +void +mpp_compute_extent(size_t isg, size_t ieg, size_t ndivs, + size_t *start, size_t *end) +{ + size_t n = ieg - isg + 1; + size_t iss = isg; + size_t ndiv, imax = ieg, ndmax = ndivs, ie, ndmirror; + char symmetrize; + + for (ndiv = 0; ndiv < ndivs; ++ndiv) { + symmetrize = (EVEN(ndivs) && EVEN(n)) || + (ODD(ndivs) && ODD(n)) || + (ODD(ndivs) && EVEN(n) && (ndivs < (n / 2))); + + if (ndiv == 0) { imax = ieg; ndmax = ndivs; } + + if (ndiv < ((ndivs - 1) / 2 + 1)) { + ie = iss + (size_t)ceil(((float)(imax - iss + 1.0)) / (ndmax - ndiv)) - 1; + ndmirror = (ndivs - 1) - ndiv; + if ((ndmirror > ndiv) && symmetrize) { + start[ndmirror] = MAX(isg + ieg - ie, ie + 1); + end[ndmirror] = MAX(isg + ieg - iss, ie + 1); + imax = start[ndmirror] - 1; + ndmax = ndmax - 1; + } + } else { + if (symmetrize) { + iss = start[ndiv]; + ie = end[ndiv]; + } else { + ie = iss + (size_t)ceil(((float)(imax - iss + 1.0)) / (ndmax - ndiv)) - 1; + } + } + start[ndiv] = iss; + end[ndiv] = ie; + iss = ie + 1; + } +} + +/* ------------------------------------------------------------------ */ +/* hyperslabcopy (verbatim port from mppncscatter.c) */ +/* ------------------------------------------------------------------ */ +void +hyperslabcopy(nc_type type, size_t *dimlen, int *dimids, + size_t *start, size_t *count, int ndim, + char *ti, short *si, int *ii, float *fi, double *di, + char *t, short *s, int *iv, float *f, double *d) +{ + size_t k, j, offset = 0, i0, stridek, stridej; + switch (ndim) { + case 0: + switch (type) { + case NC_BYTE: case NC_CHAR: *t = *ti; break; + case NC_SHORT: *s = *si; break; + case NC_INT: *iv = *ii; break; + case NC_FLOAT: *f = *fi; break; + case NC_DOUBLE:*d = *di; break; + default: break; + } + break; + case 1: + switch (type) { + case NC_BYTE: case NC_CHAR: memcpy(t, ti + start[0], count[0] * sizeof(char)); break; + case NC_SHORT: memcpy(s, si + start[0], count[0] * sizeof(short)); break; + case NC_INT: memcpy(iv, ii + start[0], count[0] * sizeof(int)); break; + case NC_FLOAT: memcpy(f, fi + start[0], count[0] * sizeof(float)); break; + case NC_DOUBLE:memcpy(d, di + start[0], count[0] * sizeof(double)); break; + default: break; + } + break; + case 2: + i0 = start[0] * dimlen[dimids[1]]; + stridek = dimlen[dimids[1]]; + switch (type) { + case NC_BYTE: case NC_CHAR: + for (k = 0; k < count[0]; ++k) { memcpy(t+offset, ti+i0+k*stridek+start[1], count[1]*sizeof(char)); offset += count[1]; } break; + case NC_SHORT: + for (k = 0; k < count[0]; ++k) { memcpy(s+offset, si+i0+k*stridek+start[1], count[1]*sizeof(short)); offset += count[1]; } break; + case NC_INT: + for (k = 0; k < count[0]; ++k) { memcpy(iv+offset, ii+i0+k*stridek+start[1], count[1]*sizeof(int)); offset += count[1]; } break; + case NC_FLOAT: + for (k = 0; k < count[0]; ++k) { memcpy(f+offset, fi+i0+k*stridek+start[1], count[1]*sizeof(float)); offset += count[1]; } break; + case NC_DOUBLE: + for (k = 0; k < count[0]; ++k) { memcpy(d+offset, di+i0+k*stridek+start[1], count[1]*sizeof(double)); offset += count[1]; } break; + default: break; + } + break; + case 3: + i0 = start[0]*dimlen[dimids[1]]*dimlen[dimids[2]] + start[1]*dimlen[dimids[2]] + start[2]; + stridek = dimlen[dimids[1]]*dimlen[dimids[2]]; + stridej = dimlen[dimids[2]]; + switch (type) { + case NC_BYTE: case NC_CHAR: + for (k = 0; k < count[0]; ++k) { for (j = 0; j < count[1]; ++j) { memcpy(t+offset, ti+i0+k*stridek+j*stridej, count[2]*sizeof(char)); offset += count[2]; } } break; + case NC_SHORT: + for (k = 0; k < count[0]; ++k) { for (j = 0; j < count[1]; ++j) { memcpy(s+offset, si+i0+k*stridek+j*stridej, count[2]*sizeof(short)); offset += count[2]; } } break; + case NC_INT: + for (k = 0; k < count[0]; ++k) { for (j = 0; j < count[1]; ++j) { memcpy(iv+offset, ii+i0+k*stridek+j*stridej, count[2]*sizeof(int)); offset += count[2]; } } break; + case NC_FLOAT: + for (k = 0; k < count[0]; ++k) { for (j = 0; j < count[1]; ++j) { memcpy(f+offset, fi+i0+k*stridek+j*stridej, count[2]*sizeof(float)); offset += count[2]; } } break; + case NC_DOUBLE: + for (k = 0; k < count[0]; ++k) { for (j = 0; j < count[1]; ++j) { memcpy(d+offset, di+i0+k*stridek+j*stridej, count[2]*sizeof(double)); offset += count[2]; } } break; + default: break; + } + break; + case 4: + i0 = start[1]*dimlen[dimids[2]]*dimlen[dimids[3]] + start[2]*dimlen[dimids[3]] + start[3]; + stridek = dimlen[dimids[2]]*dimlen[dimids[3]]; + stridej = dimlen[dimids[3]]; + switch (type) { + case NC_BYTE: case NC_CHAR: + for (k = 0; k < count[1]; ++k) { for (j = 0; j < count[2]; ++j) { memcpy(t+offset, ti+i0+k*stridek+j*stridej, count[3]*sizeof(char)); offset += count[3]; } } break; + case NC_SHORT: + for (k = 0; k < count[0]; ++k) { for (j = 0; j < count[1]; ++j) { memcpy(s+offset, si+i0+k*stridek+j*stridej, count[2]*sizeof(short)); offset += count[2]; } } break; + case NC_INT: + for (k = 0; k < count[1]; ++k) { for (j = 0; j < count[2]; ++j) { memcpy(iv+offset, ii+i0+k*stridek+j*stridej, count[3]*sizeof(int)); offset += count[3]; } } break; + case NC_FLOAT: + for (k = 0; k < count[1]; ++k) { for (j = 0; j < count[2]; ++j) { memcpy(f+offset, fi+i0+k*stridek+j*stridej, count[3]*sizeof(float)); offset += count[3]; } } break; + case NC_DOUBLE: + for (k = 0; k < count[1]; ++k) { for (j = 0; j < count[2]; ++j) { memcpy(d+offset, di+i0+k*stridek+j*stridej, count[3]*sizeof(double)); offset += count[3]; } } break; + default: break; + } + break; + } +} + +/* ------------------------------------------------------------------ */ +/* Internal scatter helpers */ +/* ------------------------------------------------------------------ */ +static void +get_num_divs(scatter_opts_t *opts, int *nx, int *ny) +{ + if (!opts) { *nx = *ny = 0; return; } + if (opts->nxio && opts->nyio) { *nx = opts->nxio; *ny = opts->nyio; } + else { *nx = opts->nx; *ny = opts->ny; } +} + +int +scatter_get_num_files(scatter_opts_t *opts) +{ + int nx, ny; + if (!opts) return 0; + get_num_divs(opts, &nx, &ny); + return nx * ny; +} + +static void +get_scatter_dims_from_file(int nc, int ndims, int nvars, + ScatterDim *scatterdims[], scatter_opts_t *opt) +{ + int recid; + if (nc_inq_unlimdim(nc, &recid) != NC_NOERR) recid = -1; + + for (int dimid = 0; dimid < ndims; dimid++) { + scatter_t stype = NOSCATTER; + size_t dimlen = 0; + int ndiv = 0; + char name[NC_MAX_NAME + 1]; + int varid; + nc_type att_type; + char att[256]; + + if (nc_inq_dimname(nc, dimid, name) != NC_NOERR) goto done; + + if (dimid == recid) { dimlen = NC_UNLIMITED; goto done; } + if (nc_inq_dimlen(nc, dimid, &dimlen) != NC_NOERR) goto done; + + if (instringlist(opt->xdims, name, opt->xdims_len)) { stype = SCATTERX; goto done; } + if (instringlist(opt->ydims, name, opt->ydims_len)) { stype = SCATTERY; goto done; } + + if (nc_inq_varid(nc, name, &varid) != NC_NOERR) goto done; + + /* Check units attribute for degrees_east / degrees_north */ + int found = 0; + if (nc_inq_atttype(nc, varid, "units", &att_type) == NC_NOERR && + att_type == NC_CHAR) { + memset(att, 0, sizeof(att)); + if (nc_get_att_text(nc, varid, "units", att) == NC_NOERR) { + if (!strncasecmp(att, "degrees_east",12) || !strncasecmp(att,"degree_east",11) || + !strncasecmp(att, "degrees_e",9) || !strncasecmp(att,"degree_e",8) || + !strncasecmp(att, "degreee",7) || !strncasecmp(att,"degreese",8)) { + stype = SCATTERX; found = 1; + } else if (!strncasecmp(att,"degrees_north",13) || !strncasecmp(att,"degree_north",12) || + !strncasecmp(att,"degrees_n",9) || !strncasecmp(att,"degree_n",8) || + !strncasecmp(att,"degreesn",8) || !strncasecmp(att,"degreen",7)) { + stype = SCATTERY; found = 1; + } + } + } + if (!found) { + if (nc_inq_atttype(nc, varid, "cartesian_axis", &att_type) == NC_NOERR && + att_type == NC_CHAR) { + memset(att, 0, sizeof(att)); + if (nc_get_att_text(nc, varid, "cartesian_axis", att) == NC_NOERR) { + if (!strncasecmp(att, "x", 1)) stype = SCATTERX; + else if (!strncasecmp(att, "y", 1)) stype = SCATTERY; + } + } + } + + done: + ndiv = (stype == SCATTERX) ? opt->nx : (stype == SCATTERY) ? opt->ny : 0; + scatterdims[dimid] = ScatterDim_new(dimid, dimlen, name, stype, ndiv); + } +} + +static void +get_scatter_extents(ScatterDim *scatterdims[], int ndims) +{ + for (int i = 0; i < ndims; i++) { + ScatterDim *pd = scatterdims[i]; + if (!pd || pd->scatter_type == NOSCATTER) continue; + mpp_compute_extent(0, pd->len - 1, pd->scatter_ndiv, + pd->scatter_start, pd->scatter_end); + for (size_t j = 0; j < pd->scatter_ndiv; j++) + pd->scatter_len[j] = pd->scatter_end[j] - pd->scatter_start[j] + 1; + } +} + +static void +get_scatter_extents_iolayout(ScatterDim *scatterdims[], int ndims, + int nxio, int nyio) +{ + for (int d = 0; d < ndims; d++) { + ScatterDim *pd = scatterdims[d]; + if (!pd || pd->scatter_type == NOSCATTER) continue; + int ndivio = (pd->scatter_type == SCATTERX) ? nxio : nyio; + int step = (int)pd->scatter_ndiv / ndivio; + size_t *nstart = XMALLOC(size_t, ndivio); + size_t *nend = XMALLOC(size_t, ndivio); + size_t *nlen = XMALLOC(size_t, ndivio); + int i = 0, k = 0; + while (i < (int)pd->scatter_ndiv) { + nstart[k] = pd->scatter_start[i]; + nend[k] = pd->scatter_end[i + step - 1]; + nlen[k] = nend[k] - nstart[k] + 1; + k++; i += step; + } + pd->scatter_ndiv = (size_t)ndivio; + XFREE(pd->scatter_start); pd->scatter_start = nstart; + XFREE(pd->scatter_end); pd->scatter_end = nend; + XFREE(pd->scatter_len); pd->scatter_len = nlen; + } +} + +void +scatter_dims(int nc, int ndims, int nvars, + ScatterDim *scatterdims[], scatter_opts_t *opt) +{ + get_scatter_dims_from_file(nc, ndims, nvars, scatterdims, opt); + get_scatter_extents(scatterdims, ndims); + if (opt->nxio && opt->nyio) { + if (opt->nx % opt->nxio) { + fprintf(stderr, "Error: x divisions not wholly divisible by io-layout x.\n"); + exit(1); + } + if (opt->ny % opt->nyio) { + fprintf(stderr, "Error: y divisions not wholly divisible by io-layout y.\n"); + exit(1); + } + get_scatter_extents_iolayout(scatterdims, ndims, opt->nxio, opt->nyio); + } +} + +/* ------------------------------------------------------------------ */ +/* scatter_def_dim */ +/* ------------------------------------------------------------------ */ +void +scatter_def_dim(int nc, int *ncids, int ndims, + ScatterDim *scatterdims[], scatter_opts_t *opts) +{ + int nx, ny, dummy; + get_num_divs(opts, &nx, &ny); + + for (int dimid = 0; dimid < ndims; dimid++) { + ScatterDim *sd = scatterdims[dimid]; + if (!sd) continue; + for (int yi = 0; yi < ny; yi++) { + for (int xi = 0; xi < nx; xi++) { + int i = yi * nx + xi; + if (opts->dryrun) continue; + switch (sd->scatter_type) { + case NOSCATTER: + nc_def_dim(ncids[i], sd->name, sd->len, &dummy); break; + case SCATTERX: + nc_def_dim(ncids[i], sd->name, sd->scatter_len[xi], &dummy); break; + case SCATTERY: + nc_def_dim(ncids[i], sd->name, sd->scatter_len[yi], &dummy); break; + } + } + } + } +} + +/* ------------------------------------------------------------------ */ +/* scatter_def_var */ +/* ------------------------------------------------------------------ */ +void +scatter_def_var(int nc, int *ncids, int nvars, int ndims, + ScatterDim *scatterdims[], scatter_opts_t *opts) +{ + int nx, ny; + get_num_divs(opts, &nx, &ny); + int scatatt[4]; + scatatt[0] = 1; + + for (int varid = 0; varid < nvars; varid++) { + char varname[NC_MAX_NAME + 1], attname[NC_MAX_NAME + 1]; + nc_type type; + int ndimvar, dimids[NC_MAX_DIMS], natt, newvarid; + + if (nc_inq_var(nc, varid, varname, &type, &ndimvar, dimids, &natt) != NC_NOERR) continue; + + for (int yi = 0; yi < ny; yi++) { + for (int xi = 0; xi < nx; xi++) { + int ifile = yi * nx + xi; + if (!opts->dryrun) { + if (nc_def_var(ncids[ifile], varname, type, ndimvar, dimids, &newvarid) != NC_NOERR) { + fprintf(stderr, "Error: failed to define var '%s' in file %d.\n", varname, ifile); + exit(-1); + } + for (int j = 0; j < natt; j++) { + nc_inq_attname(nc, varid, j, attname); + nc_copy_att(nc, varid, attname, ncids[ifile], newvarid); + } + } + /* Add domain_decomposition attribute for 1D dimension vars. */ + if (ndimvar == 1) { + ScatterDim *sd = scatterdims[dimids[0]]; + if (!sd) continue; + if (sd->scatter_type == SCATTERX && strcmp(varname, sd->name) == 0) { + scatatt[1] = (int)sd->len; + scatatt[2] = (int)sd->scatter_start[xi] + 1; + scatatt[3] = (int)(sd->scatter_start[xi] + sd->scatter_len[xi]); + if (!opts->dryrun) + nc_put_att_int(ncids[ifile], varid, "domain_decomposition", NC_INT, 4, scatatt); + } else if (sd->scatter_type == SCATTERY && strcmp(varname, sd->name) == 0) { + scatatt[1] = (int)sd->len; + scatatt[2] = (int)sd->scatter_start[yi] + 1; + scatatt[3] = (int)(sd->scatter_start[yi] + sd->scatter_len[yi]); + if (!opts->dryrun) + nc_put_att_int(ncids[ifile], varid, "domain_decomposition", NC_INT, 4, scatatt); + } + } + } + } + } +} + +/* ------------------------------------------------------------------ */ +/* scatter_put_var */ +/* ------------------------------------------------------------------ */ +void +scatter_put_var(int nc, int *ncids, int ndims, int nvars, + ScatterDim *scatterdims[], scatter_opts_t *opt) +{ + if (opt->dryrun) return; + + int nx, ny; + get_num_divs(opt, &nx, &ny); + + int recid; + size_t nrec = 0; + nc_inq_unlimdim(nc, &recid); + if (recid != -1) nc_inq_dimlen(nc, recid, &nrec); + + size_t dimlen[NC_MAX_DIMS]; + for (int i = 0; i < ndims; i++) nc_inq_dimlen(nc, i, &dimlen[i]); + + /* Find maximum buffer sizes needed. */ + size_t maxchar = 0, maxshort = 0, maxint = 0, maxfloat = 0, maxdouble = 0; + for (int varid = 0; varid < nvars; varid++) { + nc_type type; + int ndimvar, dimids[NC_MAX_DIMS], natt; + nc_inq_var(nc, varid, NULL, &type, &ndimvar, dimids, &natt); + size_t sz = 1; + for (int i = 0; i < ndimvar; i++) + if (dimids[i] != recid) sz *= dimlen[dimids[i]]; + switch (type) { + case NC_BYTE: case NC_CHAR: if (sz > maxchar) maxchar = sz; break; + case NC_SHORT: if (sz > maxshort) maxshort = sz; break; + case NC_INT: if (sz > maxint) maxint = sz; break; + case NC_FLOAT: if (sz > maxfloat) maxfloat = sz; break; + case NC_DOUBLE: if (sz > maxdouble) maxdouble = sz; break; + default: break; + } + } + + char *tp = maxchar ? (char *)malloc(maxchar * sizeof(char)) : NULL; + char *otp= maxchar ? (char *)malloc(maxchar * sizeof(char)) : NULL; + short *sp = maxshort ? (short *)malloc(maxshort * sizeof(short)) : NULL; + short *osp= maxshort ? (short *)malloc(maxshort * sizeof(short)) : NULL; + int *ip = maxint ? (int *)malloc(maxint * sizeof(int)) : NULL; + int *oip= maxint ? (int *)malloc(maxint * sizeof(int)) : NULL; + float *fp = maxfloat ? (float *)malloc(maxfloat * sizeof(float)) : NULL; + float *ofp= maxfloat ? (float *)malloc(maxfloat * sizeof(float)) : NULL; + double *dp = maxdouble ? (double *)malloc(maxdouble * sizeof(double)) : NULL; + double *odp= maxdouble ? (double *)malloc(maxdouble * sizeof(double)) : NULL; + + size_t outstart[4] = {0, 0, 0, 0}; + size_t instart[4], count[4]; + + /* Process non-record variables. */ + for (int varid = 0; varid < nvars; varid++) { + char varname[NC_MAX_NAME+1]; + nc_type type; + int ndimvar, dimids[NC_MAX_DIMS], natt; + nc_inq_var(nc, varid, varname, &type, &ndimvar, dimids, &natt); + if (recid != -1 && ndimvar > 0 && dimids[0] == recid) continue; + + switch (type) { + case NC_BYTE: case NC_CHAR: nc_get_var_text(nc, varid, tp); break; + case NC_SHORT: nc_get_var_short(nc, varid, sp); break; + case NC_INT: nc_get_var_int(nc, varid, ip); break; + case NC_FLOAT: nc_get_var_float(nc, varid, fp); break; + case NC_DOUBLE: nc_get_var_double(nc, varid, dp); break; + default: break; + } + + for (int yi = 0; yi < ny; yi++) { + for (int xi = 0; xi < nx; xi++) { + int i = yi * nx + xi; + if (ndimvar > 0) { + for (int dimi = 0; dimi < ndimvar; dimi++) { + ScatterDim *sd = scatterdims[dimids[dimi]]; + if (!sd) continue; + switch (sd->scatter_type) { + case NOSCATTER: instart[dimi]=0; count[dimi]=sd->len; break; + case SCATTERX: instart[dimi]=sd->scatter_start[xi]; count[dimi]=sd->scatter_len[xi]; break; + case SCATTERY: instart[dimi]=sd->scatter_start[yi]; count[dimi]=sd->scatter_len[yi]; break; + } + } + hyperslabcopy(type, dimlen, dimids, instart, count, ndimvar, + tp, sp, ip, fp, dp, otp, osp, oip, ofp, odp); + switch (type) { + case NC_BYTE: case NC_CHAR: nc_put_vara_text(ncids[i], varid, outstart, count, otp); break; + case NC_SHORT: nc_put_vara_short(ncids[i], varid, outstart, count, osp); break; + case NC_INT: nc_put_vara_int(ncids[i], varid, outstart, count, oip); break; + case NC_FLOAT: nc_put_vara_float(ncids[i], varid, outstart, count, ofp); break; + case NC_DOUBLE: nc_put_vara_double(ncids[i], varid, outstart, count, odp); break; + default: break; + } + } else { + switch (type) { + case NC_BYTE: case NC_CHAR: nc_put_var_text(ncids[i], varid, tp); break; + case NC_SHORT: nc_put_var_short(ncids[i], varid, sp); break; + case NC_INT: nc_put_var_int(ncids[i], varid, ip); break; + case NC_FLOAT: nc_put_var_float(ncids[i], varid, fp); break; + case NC_DOUBLE: nc_put_var_double(ncids[i], varid, dp); break; + default: break; + } + } + } + } + } + + /* Process record variables. */ + for (size_t reci = 0; reci < nrec; reci++) { + instart[0] = reci; outstart[0] = reci; count[0] = 1; + for (int varid = 0; varid < nvars; varid++) { + char varname[NC_MAX_NAME+1]; + nc_type type; + int ndimvar, dimids[NC_MAX_DIMS], natt, dimids2[NC_MAX_DIMS]; + nc_inq_var(nc, varid, varname, &type, &ndimvar, dimids, &natt); + if (dimids[0] != recid) continue; + int ndimvar2 = ndimvar - 1; + + for (int dimi = 1; dimi < ndimvar; dimi++) { + count[dimi] = dimlen[dimids[dimi]]; + dimids2[dimi-1] = dimids[dimi]; + } + + switch (type) { + case NC_BYTE: case NC_CHAR: nc_get_vara_text(nc,varid,instart,count,tp); break; + case NC_SHORT: nc_get_vara_short(nc,varid,instart,count,sp); break; + case NC_INT: nc_get_vara_int(nc,varid,instart,count,ip); break; + case NC_FLOAT: nc_get_vara_float(nc,varid,instart,count,fp); break; + case NC_DOUBLE: nc_get_vara_double(nc,varid,instart,count,dp); break; + default: break; + } + + for (int yi = 0; yi < ny; yi++) { + for (int xi = 0; xi < nx; xi++) { + int i = yi * nx + xi; + size_t instart2[4]={0,0,0,0}, count2[4]={1,1,1,1}; + for (int dimi = 1; dimi < ndimvar; dimi++) { + ScatterDim *sd = scatterdims[dimids[dimi]]; + if (!sd) continue; + switch (sd->scatter_type) { + case NOSCATTER: instart[dimi]=0; count[dimi]=sd->len; break; + case SCATTERX: instart[dimi]=sd->scatter_start[xi]; count[dimi]=sd->scatter_len[xi]; break; + case SCATTERY: instart[dimi]=sd->scatter_start[yi]; count[dimi]=sd->scatter_len[yi]; break; + } + count2[dimi-1] = count[dimi]; + instart2[dimi-1] = instart[dimi]; + } + hyperslabcopy(type, dimlen, dimids2, instart2, count2, ndimvar2, + tp, sp, ip, fp, dp, otp, osp, oip, ofp, odp); + switch (type) { + case NC_BYTE: case NC_CHAR: nc_put_vara_text(ncids[i],varid,outstart,count,otp); break; + case NC_SHORT: nc_put_vara_short(ncids[i],varid,outstart,count,osp); break; + case NC_INT: nc_put_vara_int(ncids[i],varid,outstart,count,oip); break; + case NC_FLOAT: nc_put_vara_float(ncids[i],varid,outstart,count,ofp); break; + case NC_DOUBLE: nc_put_vara_double(ncids[i],varid,outstart,count,odp); break; + default: break; + } + } + } + } + } + + free(tp); free(otp); free(sp); free(osp); + free(ip); free(oip); free(fp); free(ofp); + free(dp); free(odp); +} diff --git a/src/mpp-disttool/domain.h b/src/mpp-disttool/domain.h new file mode 100644 index 00000000..ad3c19b7 --- /dev/null +++ b/src/mpp-disttool/domain.h @@ -0,0 +1,93 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_DOMAIN_H +#define MPPDISTTOOL_DOMAIN_H + +#include +#include +#include "strlist.h" + +typedef enum { NOSCATTER = 0, SCATTERX, SCATTERY } scatter_t; + +typedef struct { + int id; + size_t len; + char name[NC_MAX_NAME + 1]; + scatter_t scatter_type; + size_t *scatter_start; + size_t *scatter_end; + size_t *scatter_len; + size_t scatter_ndiv; +} ScatterDim; + +/* Options for the MPP scatter operation. */ +typedef struct { + char dryrun; + char *filein; + int nx; + int ny; + int nxio; + int nyio; + char *prefix; + int start; + int verbose; + int width; + char **xdims; + int xdims_len; + char **ydims; + int ydims_len; +} scatter_opts_t; + +/* ScatterDim lifecycle. */ +ScatterDim *ScatterDim_new(int id, size_t len, const char *name, + scatter_t scatter_type, int ndiv); +void ScatterDim_free(ScatterDim *p); +void scatter_dims_free(ScatterDim *dims[], int ndims); + +/* Symmetric MPP tiling (from mppncscatter). + Fills start[0..ndivs-1] and end[0..ndivs-1] with 0-based indices. */ +void mpp_compute_extent(size_t isg, size_t ieg, size_t ndivs, + size_t *start, size_t *end); + +/* Copy a hyperslab of multi-dimensional data (any nc_type). */ +void hyperslabcopy(nc_type type, size_t *dimlen, int *dimids, + size_t *start, size_t *count, int ndim, + char *ti, short *si, int *ii, float *fi, double *di, + char *t, short *s, int *i, float *f, double *d); + +/* Populate scatterdims[] from the open NetCDF file. + Also applies io_layout compression if nxio/nyio are set. */ +void scatter_dims(int nc, int ndims, int nvars, + ScatterDim *scatterdims[], scatter_opts_t *opt); + +/* Output file count from scatter_opts_t. */ +int scatter_get_num_files(scatter_opts_t *opts); + +/* Define dimensions and variables in output file set. */ +void scatter_def_dim(int nc, int *ncids, int ndims, + ScatterDim *scatterdims[], scatter_opts_t *opts); +void scatter_def_var(int nc, int *ncids, int nvars, int ndims, + ScatterDim *scatterdims[], scatter_opts_t *opts); + +/* Write variable data into scattered output files. */ +void scatter_put_var(int nc, int *ncids, int ndims, int nvars, + ScatterDim *scatterdims[], scatter_opts_t *opt); + +#endif /* MPPDISTTOOL_DOMAIN_H */ diff --git a/src/mpp-disttool/main.c b/src/mpp-disttool/main.c new file mode 100644 index 00000000..d3f3a816 --- /dev/null +++ b/src/mpp-disttool/main.c @@ -0,0 +1,164 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +/* + * mppdisttool – single binary replacing seven MPP distribution tools. + * Supports subcommand dispatch and argv[0] backward-compatibility symlinks. + */ +#define _POSIX_C_SOURCE 200809L +#include +#include +#include +#include /* basename() */ +#include "cmd_check.h" +#include "cmd_combine.h" +#include "cmd_scatter.h" +#include "cmd_decompress.h" +#include "cmd_auto.h" + +/* print_version is defined in src/print_version.c and linked via libver.a. */ +extern void print_version(const char *command_name); + +/* ------------------------------------------------------------------ */ +/* Subcommand dispatch table */ +/* ------------------------------------------------------------------ */ +typedef struct { const char *name; int (*fn)(int, char **); } subcmd_t; + +static const subcmd_t commands[] = { + { "check", cmd_check }, + { "combine", cmd_combine }, + { "scatter", cmd_scatter }, + { "decompress", cmd_decompress }, + { "auto", cmd_auto }, + { NULL, NULL } +}; + +static void usage(void) { + printf("mppdisttool – unified MPP distribution tool\n\n" + "Usage: mppdisttool [options]\n\n" + "Subcommands:\n" + " check Test if a file is compressed or an iceberg restart\n" + " combine Combine decomposed NetCDF files into one\n" + " scatter Decompose a NetCDF file into many\n" + " decompress Expand compressed-by-gathering NetCDF files\n" + " auto Auto-detect and combine restart files in CWD\n" + "\nBackward-compatibility symlinks (argv[0] dispatch):\n" + " mppncscatter -> scatter\n" + " mppnccombine -> combine --mpp\n" + " combine-ncc -> combine --land\n" + " scatter-ncc -> scatter --land\n" + " decompress-ncc -> decompress\n" + " is-compressed -> check --compressed\n" + "\nRun 'mppdisttool --help' for subcommand options.\n"); +} + +/* ------------------------------------------------------------------ */ +/* argv[0] backward-compatibility mapping */ +/* ------------------------------------------------------------------ */ +/* Returns extra argv words to prepend, or NULL if no match. */ +static const char *compat_argv0_extra(const char *progname) { + if (strcmp(progname, "mppncscatter") == 0) return "scatter"; + if (strcmp(progname, "mppnccombine") == 0) return "combine\0--mpp"; + if (strcmp(progname, "combine-ncc") == 0) return "combine\0--land"; + if (strcmp(progname, "scatter-ncc") == 0) return "scatter\0--land"; + if (strcmp(progname, "decompress-ncc")== 0) return "decompress"; + if (strcmp(progname, "is-compressed") == 0) return "check\0--compressed"; + return NULL; +} + +/* Build a new argv with extra words prepended after argv[0]. */ +static char **build_compat_argv(int argc, char **argv, + const char *extra, int *new_argc_out) { + /* Count extra words (NUL-terminated, double-NUL ends sequence). */ + int nextra = 0; + const char *p = extra; + while (*p) { nextra++; p += strlen(p) + 1; } + + int new_argc = argc + nextra; + char **new_argv = (char **)malloc((size_t)(new_argc + 1) * sizeof(char *)); + new_argv[0] = argv[0]; + p = extra; + for (int i = 1; i <= nextra; i++) { + new_argv[i] = (char *)p; + p += strlen(p) + 1; + } + for (int i = 1; i < argc; i++) + new_argv[nextra + i] = argv[i]; + new_argv[new_argc] = NULL; + *new_argc_out = new_argc; + return new_argv; +} + +/* ------------------------------------------------------------------ */ +/* main */ +/* ------------------------------------------------------------------ */ +int +main(int argc, char *argv[]) +{ + /* Determine invocation name for argv[0] compat. */ + char *prog_copy = strdup(argv[0]); + const char *progname = basename(prog_copy); + + /* Check for argv[0] backward-compatibility mapping. */ + const char *extra = compat_argv0_extra(progname); + if (extra) { + int new_argc; + char **new_argv = build_compat_argv(argc, argv, extra, &new_argc); + free(prog_copy); + /* Dispatch through subcommand. */ + const char *subcmd = new_argv[1]; + for (const subcmd_t *c = commands; c->name; c++) { + if (strcmp(c->name, subcmd) == 0) { + int ret = c->fn(new_argc - 1, new_argv + 1); + free(new_argv); + return ret; + } + } + free(new_argv); + fprintf(stderr,"mppdisttool: unknown subcommand '%s'\n", subcmd); + return 1; + } + free(prog_copy); + + /* Normal subcommand dispatch. */ + if (argc < 2) { + usage(); + return 1; + } + + const char *subcmd = argv[1]; + + if (strcmp(subcmd,"-V")==0 || strcmp(subcmd,"--version")==0) { + print_version("mppdisttool"); + return 0; + } + if (strcmp(subcmd,"-h")==0 || strcmp(subcmd,"--help")==0) { + usage(); + return 0; + } + + for (const subcmd_t *c = commands; c->name; c++) { + if (strcmp(c->name, subcmd) == 0) + return c->fn(argc - 1, argv + 1); + } + + fprintf(stderr,"mppdisttool: unknown subcommand '%s'\n" + "Run 'mppdisttool --help' for usage.\n", subcmd); + return 1; +} diff --git a/src/mpp-disttool/nc_utils.c b/src/mpp-disttool/nc_utils.c new file mode 100644 index 00000000..3af297bd --- /dev/null +++ b/src/mpp-disttool/nc_utils.c @@ -0,0 +1,305 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#include +#include +#include +#include +#include "nc_utils.h" + +#define NC_CHECK(call) do { \ + int _err = (call); \ + if (_err != NC_NOERR) return _err; \ +} while (0) + +size_t +ncu_type_size(nc_type xtype) +{ + switch (xtype) { + case NC_BYTE: case NC_CHAR: return 1; + case NC_SHORT: return 2; + case NC_INT: return 4; + case NC_FLOAT: return 4; + case NC_DOUBLE: return 8; + default: return 0; + } +} + +int +ncu_clone_dim(int incid, int dimid, int oncid) +{ + char name[NC_MAX_NAME + 1]; + size_t len; + int unlimdimid, newdimid; + + NC_CHECK(nc_inq_dim(incid, dimid, name, &len)); + NC_CHECK(nc_inq_unlimdim(incid, &unlimdimid)); + if (dimid == unlimdimid) + len = NC_UNLIMITED; + + /* If already present in oncid, accept as-is. */ + if (nc_inq_dimid(oncid, name, &newdimid) == NC_NOERR) + return NC_NOERR; + + NC_CHECK(nc_def_dim(oncid, name, len, &newdimid)); + return NC_NOERR; +} + +int +ncu_clone_var(int incid, int varid, int oncid) +{ + char name[NC_MAX_NAME + 1]; + char dname[NC_MAX_NAME + 1]; + char attname[NC_MAX_NAME + 1]; + nc_type xtype; + int ndims, dimids[NC_MAX_VAR_DIMS], natts, ovarid; + + NC_CHECK(nc_inq_var(incid, varid, name, &xtype, &ndims, dimids, &natts)); + + /* Remap dimension IDs to the output file's dimension IDs. */ + for (int i = 0; i < ndims; i++) { + NC_CHECK(nc_inq_dimname(incid, dimids[i], dname)); + NC_CHECK(nc_inq_dimid(oncid, dname, &dimids[i])); + } + + NC_CHECK(nc_def_var(oncid, name, xtype, ndims, dimids, &ovarid)); + + for (int i = 0; i < natts; i++) { + NC_CHECK(nc_inq_attname(incid, varid, i, attname)); + NC_CHECK(nc_copy_att(incid, varid, attname, oncid, ovarid)); + } + return NC_NOERR; +} + +int +ncu_copy_global_atts(int incid, int oncid) +{ + int ngatts; + char attname[NC_MAX_NAME + 1]; + + NC_CHECK(nc_inq_natts(incid, &ngatts)); + for (int i = 0; i < ngatts; i++) { + NC_CHECK(nc_inq_attname(incid, NC_GLOBAL, i, attname)); + NC_CHECK(nc_copy_att(incid, NC_GLOBAL, attname, oncid, NC_GLOBAL)); + } + return NC_NOERR; +} + +int +ncu_get_var_int(int ncid, const char *name, int *data) +{ + int varid; + NC_CHECK(nc_inq_varid(ncid, name, &varid)); + return nc_get_var_int(ncid, varid, data); +} + +int +ncu_inq_var(int ncid, int varid, + char *name_out, + nc_type *xtype_out, + int *ndims_out, + int *dimids_out, + size_t *dimlens_out, + int *natts_out, + bool *is_dim_out, + bool *has_records_out, + size_t *varsize_out, + size_t *recsize_out, + int *nrec_out) +{ + char name[NC_MAX_NAME + 1]; + nc_type xtype; + int ndims, dimids[NC_MAX_VAR_DIMS], natts; + int unlimdim; + + NC_CHECK(nc_inq_var(ncid, varid, name, &xtype, &ndims, dimids, &natts)); + NC_CHECK(nc_inq_unlimdim(ncid, &unlimdim)); + + if (name_out) strcpy(name_out, name); + if (xtype_out) *xtype_out = xtype; + if (ndims_out) *ndims_out = ndims; + if (natts_out) *natts_out = natts; + if (dimids_out) + memcpy(dimids_out, dimids, ndims * sizeof(int)); + + if (is_dim_out) { + int dummy; + *is_dim_out = (nc_inq_dimid(ncid, name, &dummy) == NC_NOERR); + } + if (has_records_out) { + *has_records_out = false; + for (int i = 0; i < ndims; i++) + if (dimids[i] == unlimdim) { *has_records_out = true; break; } + } + + if (dimlens_out || varsize_out || recsize_out || nrec_out) { + size_t vsize = 1, rsize = 1; + for (int i = 0; i < ndims; i++) { + size_t dlen; + NC_CHECK(nc_inq_dimlen(ncid, dimids[i], &dlen)); + if (dimlens_out) dimlens_out[i] = dlen; + vsize *= dlen; + if (dimids[i] != unlimdim) rsize *= dlen; + } + if (varsize_out) *varsize_out = vsize; + if (recsize_out) *recsize_out = rsize; + if (nrec_out) { + *nrec_out = 1; + if (unlimdim != -1) { + bool has_rec = false; + for (int i = 0; i < ndims; i++) + if (dimids[i] == unlimdim) { has_rec = true; break; } + if (has_rec) { + size_t nrec_sz; + NC_CHECK(nc_inq_dimlen(ncid, unlimdim, &nrec_sz)); + *nrec_out = (int)nrec_sz; + } + } + } + } + return NC_NOERR; +} + +int +ncu_inq_dim(int ncid, int dimid, + char *name_out, + size_t *len_out, + bool *is_unlim_out) +{ + char name[NC_MAX_NAME + 1]; + size_t len; + int unlimdim; + + NC_CHECK(nc_inq_dim(ncid, dimid, name, &len)); + if (name_out) strcpy(name_out, name); + if (len_out) *len_out = len; + if (is_unlim_out) { + NC_CHECK(nc_inq_unlimdim(ncid, &unlimdim)); + *is_unlim_out = (dimid == unlimdim); + } + return NC_NOERR; +} + +int +ncu_format_cmode(int incid) +{ + int fmt; + if (nc_inq_format(incid, &fmt) != NC_NOERR) + return NC_64BIT_OFFSET; + switch (fmt) { + case NC_FORMAT_NETCDF4: + return NC_NETCDF4; + case NC_FORMAT_NETCDF4_CLASSIC: + return NC_NETCDF4 | NC_CLASSIC_MODEL; + case NC_FORMAT_64BIT_DATA: + return NC_64BIT_DATA; + /* classic and 64-bit offset both produce 64-bit offset output, + matching mppnccombine behaviour */ + default: + return NC_64BIT_OFFSET; + } +} + +int +ncu_get_vara(int ncid, int varid, const size_t *start, + const size_t *count, void *buf) +{ + nc_type xtype; + NC_CHECK(nc_inq_vartype(ncid, varid, &xtype)); + switch (xtype) { + case NC_BYTE: case NC_CHAR: + return nc_get_vara_text(ncid, varid, start, count, (char *)buf); + case NC_SHORT: + return nc_get_vara_short(ncid, varid, start, count, (short *)buf); + case NC_INT: + return nc_get_vara_int(ncid, varid, start, count, (int *)buf); + case NC_FLOAT: + return nc_get_vara_float(ncid, varid, start, count, (float *)buf); + case NC_DOUBLE: + return nc_get_vara_double(ncid, varid, start, count, (double *)buf); + default: + return NC_EBADTYPE; + } +} + +int +ncu_put_vara(int ncid, int varid, const size_t *start, + const size_t *count, const void *buf) +{ + nc_type xtype; + NC_CHECK(nc_inq_vartype(ncid, varid, &xtype)); + switch (xtype) { + case NC_BYTE: case NC_CHAR: + return nc_put_vara_text(ncid, varid, start, count, (const char *)buf); + case NC_SHORT: + return nc_put_vara_short(ncid, varid, start, count, (const short *)buf); + case NC_INT: + return nc_put_vara_int(ncid, varid, start, count, (const int *)buf); + case NC_FLOAT: + return nc_put_vara_float(ncid, varid, start, count, (const float *)buf); + case NC_DOUBLE: + return nc_put_vara_double(ncid, varid, start, count, (const double *)buf); + default: + return NC_EBADTYPE; + } +} + +int +ncu_get_rec(int ncid, int varid, int rec, void *buf) +{ + int ndims, dimids[NC_MAX_VAR_DIMS], unlimdim; + size_t start[NC_MAX_VAR_DIMS], count[NC_MAX_VAR_DIMS]; + + NC_CHECK(nc_inq_unlimdim(ncid, &unlimdim)); + NC_CHECK(nc_inq_varndims(ncid, varid, &ndims)); + NC_CHECK(nc_inq_vardimid(ncid, varid, dimids)); + + for (int i = 0; i < ndims; i++) { + if (dimids[i] == unlimdim) { + start[i] = (size_t)rec; + count[i] = 1; + } else { + start[i] = 0; + NC_CHECK(nc_inq_dimlen(ncid, dimids[i], &count[i])); + } + } + return ncu_get_vara(ncid, varid, start, count, buf); +} + +int +ncu_put_rec(int ncid, int varid, int rec, const void *buf) +{ + int ndims, dimids[NC_MAX_VAR_DIMS], unlimdim; + size_t start[NC_MAX_VAR_DIMS], count[NC_MAX_VAR_DIMS]; + + NC_CHECK(nc_inq_unlimdim(ncid, &unlimdim)); + NC_CHECK(nc_inq_varndims(ncid, varid, &ndims)); + NC_CHECK(nc_inq_vardimid(ncid, varid, dimids)); + + for (int i = 0; i < ndims; i++) { + if (dimids[i] == unlimdim) { + start[i] = (size_t)rec; + count[i] = 1; + } else { + start[i] = 0; + NC_CHECK(nc_inq_dimlen(ncid, dimids[i], &count[i])); + } + } + return ncu_put_vara(ncid, varid, start, count, buf); +} diff --git a/src/mpp-disttool/nc_utils.h b/src/mpp-disttool/nc_utils.h new file mode 100644 index 00000000..c6f3af11 --- /dev/null +++ b/src/mpp-disttool/nc_utils.h @@ -0,0 +1,86 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_NC_UTILS_H +#define MPPDISTTOOL_NC_UTILS_H + +#include +#include +#include + +/* Clone a dimension from incid to oncid (creates if not present). + If dimension already exists in oncid, verifies size matches. */ +int ncu_clone_dim(int incid, int dimid, int oncid); + +/* Clone a variable definition (and all its attributes) from incid to oncid. + Dimensions in the variable are looked up by name in oncid. */ +int ncu_clone_var(int incid, int varid, int oncid); + +/* Copy all global attributes from incid to oncid. */ +int ncu_copy_global_atts(int incid, int oncid); + +/* Get integer variable data by name. 'data' must be pre-allocated. */ +int ncu_get_var_int(int ncid, const char *name, int *data); + +/* Query variable by id: optionally return name, xtype, ndims, dimids, + dimlens (per-dim lengths), natts, is_dim flag, has_records flag, + varsize (total elements), recsize (elements per record), nrec. + Pass NULL for outputs you don't need. */ +int ncu_inq_var(int ncid, int varid, + char *name_out, /* NC_MAX_NAME+1 or NULL */ + nc_type *xtype_out, + int *ndims_out, + int *dimids_out, /* [NC_MAX_VAR_DIMS] or NULL */ + size_t *dimlens_out, /* [NC_MAX_VAR_DIMS] or NULL */ + int *natts_out, + bool *is_dim_out, + bool *has_records_out, + size_t *varsize_out, + size_t *recsize_out, + int *nrec_out); + +/* Query dimension by id. Pass NULL for outputs not needed. */ +int ncu_inq_dim(int ncid, int dimid, + char *name_out, /* NC_MAX_NAME+1 or NULL */ + size_t *len_out, + bool *is_unlim_out); + +/* Return the nc_create flags (format) appropriate for creating an output + file that matches the format of incid. If incid is classic, returns + NC_64BIT_OFFSET (preserving mppnccombine behaviour). */ +int ncu_format_cmode(int incid); + +/* Bytes per netCDF element type. */ +size_t ncu_type_size(nc_type xtype); + +/* Read a single record of a variable (record dimension = rec, 0-based). + 'buf' must be pre-allocated to recsize elements of the variable type. + Works for any nc_type. */ +int ncu_get_rec(int ncid, int varid, int rec, void *buf); + +/* Write a single record (rec, 0-based) to a variable. */ +int ncu_put_rec(int ncid, int varid, int rec, const void *buf); + +/* Generic vara read/write for any nc_type. */ +int ncu_get_vara(int ncid, int varid, const size_t *start, + const size_t *count, void *buf); +int ncu_put_vara(int ncid, int varid, const size_t *start, + const size_t *count, const void *buf); + +#endif /* MPPDISTTOOL_NC_UTILS_H */ diff --git a/src/mpp-disttool/strlist.c b/src/mpp-disttool/strlist.c new file mode 100644 index 00000000..b4a5f5c0 --- /dev/null +++ b/src/mpp-disttool/strlist.c @@ -0,0 +1,196 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#include "strlist.h" + +int newstringlist(char ***list, int *n, int size) +{ + int i; + if (*list != NULL) return EXIT_FAILURE; + *list = XMALLOC(char *, size); + for (i = 0; i < size; ++i) + (*list)[i] = NULL; + *n = size; + return EXIT_SUCCESS; +} + +int clearstringlist(char **list, int n) +{ + if ((list == NULL) || (n <= 0)) + return 1; + for (; --n >= 0;) { + if (list[n] == NULL) continue; + XFREE(list[n]); + list[n] = NULL; + } + return EXIT_SUCCESS; +} + +void getstringlist(char *optarg, char ***list, int *nitems) +{ + char *p; + int i = 0; + char *saveptr = NULL; + p = optarg; + *nitems = 1; + while (*p++) { if (*p == ',') ++(*nitems); } + if (*list == NULL) return; + for (p = strtok_r(optarg, ",", &saveptr); + p != NULL; + p = strtok_r(NULL, ",", &saveptr)) { + size_t plen = strlen(p); + (*list)[i] = XMALLOC(char, plen + 1); + if ((*list)[i] == NULL) { + fprintf(stderr, "ERROR: Failed to allocate memory for string.\n"); + exit(1); + } + memcpy((*list)[i], p, plen); + (*list)[i][plen] = '\0'; + ++i; + } +} + +void freestringlist(char ***list, int nitems) +{ + if (*list == NULL) return; + for (; nitems > 0; --nitems) + XFREE((*list)[nitems - 1]); + XFREE(*list); + *list = NULL; +} + +void printstrlist(char **list, int n, FILE *f) +{ + int i; + for (i = 0; i < n; ++i) + fprintf(f, "%s ", list[i]); + fprintf(f, "\n"); +} + +int instringlist(char **list, char *str, int nitems) +{ + int i; + if ((list == NULL) || (str == NULL)) return 0; + for (i = 0; i < nitems; ++i) { + if (list[i] == NULL) continue; + if (strcmp(list[i], str) == 0) return 1; + } + return 0; +} + +int addstringtolist(char **list, char *string, int nitems) +{ + int i; + if (string == NULL) return EXIT_FAILED; + if (list == NULL) return EXIT_FAILED; + for (i = 0; i < nitems; ++i) { + if (list[i] == NULL) { + size_t slen = strlen(string); + list[i] = XMALLOC(char, slen + 1); + memcpy(list[i], string, slen); + list[i][slen] = '\0'; + break; + } + } + if (i >= nitems) return EXIT_FAILED; + return EXIT_SUCCESS; +} + +int appendstringtolist(char ***list, char *string, int *nitems) +{ + int newsize; + char **newlist = NULL; + if (string == NULL) return EXIT_FAILED; + if (*list == NULL) return EXIT_FAILED; + newstringlist(&newlist, &newsize, *nitems + 1); + copystrlist(*list, newlist, *nitems, newsize); + size_t slen = strlen(string); + newlist[newsize - 1] = (char *)malloc(sizeof(char) * (slen + 1)); + if (newlist[newsize - 1] == NULL) { + fprintf(stderr, "ERROR: Failed to allocate memory for string.\n"); + exit(1); + } + memcpy(newlist[newsize - 1], string, slen); + newlist[newsize - 1][slen] = '\0'; + freestringlist(list, *nitems); + *list = newlist; + *nitems = newsize; + return EXIT_SUCCESS; +} + +int strlistu(char **list1, char **list2, char **listunion, int n1, int n2, int nu) +{ + int i; + if (listunion == NULL) return EXIT_FAILED; + for (i = 0; i < n1; ++i) { + if (list1[i] == NULL) continue; + if (!instringlist(listunion, list1[i], nu)) + addstringtolist(listunion, list1[i], nu); + } + for (i = 0; i < n2; ++i) { + if (list2[i] == NULL) continue; + if (!instringlist(listunion, list2[i], nu)) + addstringtolist(listunion, list2[i], nu); + } + return EXIT_SUCCESS; +} + +int strlistsd(char **list1, char **list2, char **listdiff, int n1, int n2, int nsd) +{ + int i; + if (listdiff == NULL) return EXIT_FAILED; + if (n1 == 0) return EXIT_SUCCESS; + if (n2 == 0) { + for (i = 0; i < n1; ++i) + addstringtolist(listdiff, list1[i], nsd); + } else { + for (i = 0; i < n1; ++i) { + if (!instringlist(list2, list1[i], n2)) + addstringtolist(listdiff, list1[i], nsd); + } + } + return EXIT_SUCCESS; +} + +int copystrlist(char **listsrc, char **listdst, int nsrc, int ndst) +{ + int i; + if (nsrc > ndst) return EXIT_FAILED; + for (i = 0; i < ndst; ++i) { + XFREE(listdst[i]); + listdst[i] = NULL; + } + for (i = 0; i < nsrc; ++i) { + if (listsrc[i] == NULL) continue; + size_t slen = strlen(listsrc[i]); + listdst[i] = XMALLOC(char, slen + 1); + memcpy(listdst[i], listsrc[i], slen); + listdst[i][slen] = '\0'; + } + return EXIT_SUCCESS; +} + +int getnumstrlist(char **list, int nlist) +{ + int n = 0, i; + for (i = 0; i < nlist; ++i) + if ((list[i] != NULL) && (strlen(list[i]) > 0)) + ++n; + return n; +} diff --git a/src/mpp-disttool/strlist.h b/src/mpp-disttool/strlist.h new file mode 100644 index 00000000..ab10927c --- /dev/null +++ b/src/mpp-disttool/strlist.h @@ -0,0 +1,41 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_STRLIST_H +#define MPPDISTTOOL_STRLIST_H 1 + +#include +#include +#include +#include "xmalloc.h" + +int newstringlist(char ***list, int *n, int size); +int clearstringlist(char **list, int n); +void getstringlist(char *optarg, char ***list, int *nitems); +void freestringlist(char ***list, int nitems); +int instringlist(char **list, char *str, int nitems); +int addstringtolist(char **list, char *string, int nitems); +int appendstringtolist(char ***list, char *string, int *nitems); +void printstrlist(char **list, int n, FILE *f); +int strlistu(char **list1, char **list2, char **listunion, int n1, int n2, int nu); +int strlistsd(char **list1, char **list2, char **listdiff, int n1, int n2, int nsd); +int copystrlist(char **listsrc, char **listdst, int nsrc, int ndst); +int getnumstrlist(char **list, int nlist); + +#endif /* MPPDISTTOOL_STRLIST_H */ diff --git a/src/mpp-disttool/xmalloc.c b/src/mpp-disttool/xmalloc.c new file mode 100644 index 00000000..73d3cba1 --- /dev/null +++ b/src/mpp-disttool/xmalloc.c @@ -0,0 +1,59 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#include "xmalloc.h" +#include + +void * +xmalloc(size_t num) +{ + void *p = malloc(num); + if (!p) { + fprintf(stderr, "Memory exhausted\n"); + exit(EXIT_FATAL); + } + return p; +} + +void * +xrealloc(void *p, size_t num) +{ + void *new; + if (!p) + return xmalloc(num); + new = realloc(p, num); + if (!new) { + fprintf(stderr, "Memory exhausted\n"); + exit(EXIT_FATAL); + } + return new; +} + +void * +xcalloc(size_t num, size_t size) +{ + /* Check for multiplication overflow before allocating */ + if (num > 0 && size > 0 && size > SIZE_MAX / num) { + fprintf(stderr, "Memory allocation overflow prevented\n"); + exit(EXIT_FATAL); + } + void *p = xmalloc(num * size); + memset(p, 0, num * size); + return p; +} diff --git a/src/mpp-disttool/xmalloc.h b/src/mpp-disttool/xmalloc.h new file mode 100644 index 00000000..0081c8d3 --- /dev/null +++ b/src/mpp-disttool/xmalloc.h @@ -0,0 +1,46 @@ +/*********************************************************************** + * GNU Lesser General Public License + * + * This file is part of the GFDL FRE NetCDF tools package (FRE-NCTools). + * + * FRE-NCtools is free software: you can redistribute it and/or modify it under + * the terms of the GNU Lesser General Public License as published by + * the Free Software Foundation, either version 3 of the License, or (at + * your option) any later version. + * + * FRE-NCtools is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with FRE-NCTools. If not, see + * . + **********************************************************************/ +#ifndef MPPDISTTOOL_XMALLOC_H +#define MPPDISTTOOL_XMALLOC_H 1 + +#include +#include +#include + +#define EXIT_FATAL 2 +#define EXIT_FAILED 3 + +#define XCALLOC(type, num) \ + ((type*) xcalloc((size_t)(num), (size_t)sizeof(type))) + +#define XMALLOC(type, num) \ + ((type*) xmalloc((size_t)((num) * sizeof(type)))) + +#define XREALLOC(type, p, num) \ + ((type*) xrealloc((p), (size_t)((num) * sizeof(type)))) + +#define XFREE(stale) \ + do { if (stale) { free(stale); stale = NULL; } } while (0) + +extern void *xcalloc(size_t num, size_t size); +extern void *xmalloc(size_t num); +extern void *xrealloc(void *p, size_t num); + +#endif /* MPPDISTTOOL_XMALLOC_H */