Fixed h1_summary and h1_summary2 to correctly construct time_bounds values.
authorForrest Hoffman <forrest@climatemodeling.org>
Mon, 01 Oct 2007 15:12:14 -0400
changeset 12ce4ee911439
parent 0 3c02cce30be8
child 2 978f4510987d
Fixed h1_summary and h1_summary2 to correctly construct time_bounds values.

h1_summary and h1_summary2 previously used the time_bounds from the last time
value from the last input file when summarizing, suggesting that the field
values were appropriate only over that short time range instead of the complete
time period over which the statistics were calculated. C-LAMP Experiment 1
runs used the previous code, so the time_bounds were incorrect in the
statistical summaries produced. C-LAMP Experiment 2 runs will use this new
code for production of statistical summaries.
Makefile
h1_summary.c
h1_summary2.c
     1.1 --- a/Makefile	Wed Sep 26 17:16:40 2007 -0400
     1.2 +++ b/Makefile	Mon Oct 01 15:12:14 2007 -0400
     1.3 @@ -3,12 +3,12 @@
     1.4  # robin1 and phoenix
     1.5  LIBS=-L$(NETCDF)/lib -lnetcdf -lm
     1.6  # Penguins
     1.7 -#LIBS=-L/usr/local/netcdf/netcdf-3.6.1-gcc+pgi/lib -lnetcdf -lm
     1.8 +#LIBS=-L/usr/lib64/netcdf-3 -lnetcdf -lm
     1.9  #
    1.10  # robin1 and phoenix
    1.11  CPPFLAGS=-I$(NETCDF)/include
    1.12  # Penguins
    1.13 -#CPPFLAGS=-I/usr/local/netcdf/netcdf-3.6.1-gcc+pgi/include
    1.14 +#CPPFLAGS=-I/usr/include/netcdf-3
    1.15  # phoenix
    1.16  #CC=cc
    1.17  #CFLAGS=-O -h list=m $(CPPFLAGS)
     2.1 --- a/h1_summary.c	Wed Sep 26 17:16:40 2007 -0400
     2.2 +++ b/h1_summary.c	Mon Oct 01 15:12:14 2007 -0400
     2.3 @@ -82,6 +82,7 @@
     2.4  };
     2.5  
     2.6  static char *time_name = "time";
     2.7 +static char *time_bounds_name = "time_bounds";
     2.8  static char *mcsec_name = "mcsec";
     2.9  static char *history_name = "history";
    2.10  static char *nsamples_name = "nsamples";
    2.11 @@ -783,7 +784,8 @@
    2.12  }
    2.13  
    2.14  void copy_metadata(int in_ncid, struct var *in_var_head,
    2.15 -	struct dim **in_dim_idx, int out_ncid, struct var *out_var_head)
    2.16 +	struct dim **in_dim_idx, int out_ncid, struct var *out_var_head,
    2.17 +	double *tbounds)
    2.18  {
    2.19  	int i, j;
    2.20  	float hr[HOURS_PER_DAY];
    2.21 @@ -819,7 +821,10 @@
    2.22  time-based metadata variable %s\n", in_vnode->name);
    2.23  #endif /* DEBUG */
    2.24  					val = read_timeslice(in_ncid, in_vnode->ncvarid, in_vnode->nctype, in_vnode->ndims, in_vnode->dimids, in_dim_idx, ((input_dim_idx[in_vnode->dimids[0]])->len - 1));
    2.25 -					write_timeslice(out_ncid, out_vnode->ncvarid, out_vnode->nctype, in_vnode->ndims, in_vnode->dimids, in_dim_idx, val, 0);
    2.26 +					if (!strcmp(in_vnode->name, time_bounds_name))
    2.27 +						write_timeslice(out_ncid, out_vnode->ncvarid, out_vnode->nctype, in_vnode->ndims, in_vnode->dimids, in_dim_idx, tbounds, 0);
    2.28 +					else
    2.29 +						write_timeslice(out_ncid, out_vnode->ncvarid, out_vnode->nctype, in_vnode->ndims, in_vnode->dimids, in_dim_idx, val, 0);
    2.30  				}
    2.31  				free(val);
    2.32  				/*
    2.33 @@ -845,11 +850,63 @@
    2.34  	return;
    2.35  }
    2.36  
    2.37 -void open_inout(char *input_fname, char *mean_fname, char *stddev_fname, char *flist, size_t nsamples)
    2.38 +void get_time_bounds(char *first_fname, char *last_fname, double *tbounds)
    2.39 +{
    2.40 +	int time_dimid, time_bounds_varid;
    2.41 +	size_t time_len, time_bounds_index[2];
    2.42 +	double bnd1, bnd2;
    2.43 +
    2.44 +	/* Open first input file */
    2.45 +	wrap_nc(nc_open(first_fname, NC_NOWRITE, &input_ncid));
    2.46 +	/* Get dimension ID of the time dimension */
    2.47 +	wrap_nc(nc_inq_dimid(input_ncid, time_name, &time_dimid));
    2.48 +	/* Get length of time dimension */
    2.49 +	wrap_nc(nc_inq_dimlen(input_ncid, time_dimid, &time_len));
    2.50 +	/* Get variable ID of time_bounds variable */
    2.51 +	wrap_nc(nc_inq_varid(input_ncid, time_bounds_name, &time_bounds_varid));
    2.52 +	/* Set index for reading the first value of time_bounds from the
    2.53 +	 * first timeslice of the first file */
    2.54 +	time_bounds_index[0] = time_bounds_index[1] = 0;
    2.55 +	/* Read the value */
    2.56 +	wrap_nc(nc_get_var1_double(input_ncid, time_bounds_varid,
    2.57 +		time_bounds_index, &bnd1));
    2.58 +	/* If the first and last file are not the same, close the first one and
    2.59 +	 * open the second one */
    2.60 +	if (strcmp(first_fname, last_fname)) {
    2.61 +		/* Close the first input file */
    2.62 +		wrap_nc(nc_close(input_ncid));
    2.63 +		/* Open the last input file */
    2.64 +		wrap_nc(nc_open(last_fname, NC_NOWRITE, &input_ncid));
    2.65 +		/* Get dimension ID of the time dimension */
    2.66 +		wrap_nc(nc_inq_dimid(input_ncid, time_name, &time_dimid));
    2.67 +		/* Get length of time dimension */
    2.68 +		wrap_nc(nc_inq_dimlen(input_ncid, time_dimid, &time_len));
    2.69 +		/* Get variable ID of time_bounds variable */
    2.70 +		wrap_nc(nc_inq_varid(input_ncid, time_bounds_name,
    2.71 +			&time_bounds_varid));
    2.72 +	}
    2.73 +	/* Set index for reading the second value of time_bounds from the
    2.74 +	 * last timeslice of the last file */
    2.75 +	time_bounds_index[0] = time_len - 1; time_bounds_index[1] = 1;
    2.76 +	/* Read the value */
    2.77 +	wrap_nc(nc_get_var1_double(input_ncid, time_bounds_varid,
    2.78 +		time_bounds_index, &bnd2));
    2.79 +
    2.80 +	/* Close the last input file */
    2.81 +	wrap_nc(nc_close(input_ncid));
    2.82 +
    2.83 +	tbounds[0] = bnd1;
    2.84 +	tbounds[1] = bnd2;
    2.85 +
    2.86 +	return;
    2.87 +}
    2.88 +
    2.89 +void open_inout(char *first_fname, char *last_fname, char *mean_fname, char *stddev_fname, char *flist, size_t nsamples)
    2.90  {
    2.91  	char *mean_history_gatt, *stddev_history_gatt,
    2.92  		*mean_prefix="Statistical means from history files:",
    2.93  		*stddev_prefix="Statistical standard deviations from history files:";
    2.94 +	double tbounds[2];
    2.95  
    2.96  	/*
    2.97  	 * Construct strings for history global attributes for the two output
    2.98 @@ -866,8 +923,16 @@
    2.99  	sprintf(mean_history_gatt, "%s%s", mean_prefix, flist);
   2.100  	sprintf(stddev_history_gatt, "%s%s", stddev_prefix, flist);
   2.101  
   2.102 -	/* Open input file */
   2.103 -	wrap_nc(nc_open(input_fname, NC_NOWRITE, &input_ncid));
   2.104 +	/* The two time_bounds values must be handled explicitly by obtaining
   2.105 +	 * the first one from the first file and the last one from the last
   2.106 +	 * file.  These are then passed to copy_metadata() below. */
   2.107 +	get_time_bounds(first_fname, last_fname, tbounds);
   2.108 +#ifdef DEBUG
   2.109 +	fprintf(stderr, "Got back tbounds=%lf,%lf\n", tbounds[0], tbounds[1]);
   2.110 +#endif /* DEBUG */
   2.111 +
   2.112 +	/* Open last input file */
   2.113 +	wrap_nc(nc_open(last_fname, NC_NOWRITE, &input_ncid));
   2.114  	/* Inquire about number of dimensions, variables, global attributes
   2.115  	 * and the ID of the unlimited dimension
   2.116  	 */
   2.117 @@ -909,9 +974,9 @@
   2.118  	wrap_nc(nc_enddef(stddev_ncid));
   2.119  	/* Write out metdata variables that do not depend on "time" */
   2.120  	copy_metadata(input_ncid, input_var_head, input_dim_idx, mean_ncid,
   2.121 -		mean_var_head);
   2.122 +		mean_var_head, tbounds);
   2.123  	copy_metadata(input_ncid, input_var_head, input_dim_idx, stddev_ncid,
   2.124 -		stddev_var_head);
   2.125 +		stddev_var_head, tbounds);
   2.126  
   2.127  	wrap_nc(nc_close(input_ncid));
   2.128  
   2.129 @@ -1513,7 +1578,8 @@
   2.130  	 * only values from the last time slice from the period over which
   2.131  	 * the statistics are computed.
   2.132  	 */
   2.133 -	open_inout(ifnames[(nifnames-1)], mfname, sfname, flist, nsamples);
   2.134 +	open_inout(ifnames[0], ifnames[(nifnames-1)], mfname, sfname, flist,
   2.135 +		nsamples);
   2.136  
   2.137  	compute_stats(nifnames, ifnames, nsamples);
   2.138  
     3.1 --- a/h1_summary2.c	Wed Sep 26 17:16:40 2007 -0400
     3.2 +++ b/h1_summary2.c	Mon Oct 01 15:12:14 2007 -0400
     3.3 @@ -82,6 +82,7 @@
     3.4  };
     3.5  
     3.6  static char *time_name = "time";
     3.7 +static char *time_bounds_name = "time_bounds";
     3.8  static char *mcsec_name = "mcsec";
     3.9  static char *history_name = "history";
    3.10  static char *nsamples_name = "nsamples";
    3.11 @@ -783,7 +784,8 @@
    3.12  }
    3.13  
    3.14  void copy_metadata(int in_ncid, struct var *in_var_head,
    3.15 -	struct dim **in_dim_idx, int out_ncid, struct var *out_var_head)
    3.16 +	struct dim **in_dim_idx, int out_ncid, struct var *out_var_head,
    3.17 +	double *tbounds)
    3.18  {
    3.19  	int i, j;
    3.20  	float hr[HOURS_PER_DAY];
    3.21 @@ -819,7 +821,10 @@
    3.22  time-based metadata variable %s\n", in_vnode->name);
    3.23  #endif /* DEBUG */
    3.24  					val = read_timeslice(in_ncid, in_vnode->ncvarid, in_vnode->nctype, in_vnode->ndims, in_vnode->dimids, in_dim_idx, ((input_dim_idx[in_vnode->dimids[0]])->len - 1));
    3.25 -					write_timeslice(out_ncid, out_vnode->ncvarid, out_vnode->nctype, in_vnode->ndims, in_vnode->dimids, in_dim_idx, val, 0);
    3.26 +					if (!strcmp(in_vnode->name, time_bounds_name))
    3.27 +						write_timeslice(out_ncid, out_vnode->ncvarid, out_vnode->nctype, in_vnode->ndims, in_vnode->dimids, in_dim_idx, tbounds, 0);
    3.28 +					else
    3.29 +						write_timeslice(out_ncid, out_vnode->ncvarid, out_vnode->nctype, in_vnode->ndims, in_vnode->dimids, in_dim_idx, val, 0);
    3.30  				}
    3.31  				free(val);
    3.32  				/*
    3.33 @@ -845,11 +850,63 @@
    3.34  	return;
    3.35  }
    3.36  
    3.37 -void open_inout(char *input_fname, char *mean_fname, char *stddev_fname, char *flist, size_t nsamples)
    3.38 +void get_time_bounds(char *first_fname, char *last_fname, double *tbounds)
    3.39 +{
    3.40 +	int time_dimid, time_bounds_varid;
    3.41 +	size_t time_len, time_bounds_index[2];
    3.42 +	double bnd1, bnd2;
    3.43 +
    3.44 +	/* Open first input file */
    3.45 +	wrap_nc(nc_open(first_fname, NC_NOWRITE, &input_ncid));
    3.46 +	/* Get dimension ID of the time dimension */
    3.47 +	wrap_nc(nc_inq_dimid(input_ncid, time_name, &time_dimid));
    3.48 +	/* Get length of time dimension */
    3.49 +	wrap_nc(nc_inq_dimlen(input_ncid, time_dimid, &time_len));
    3.50 +	/* Get variable ID of time_bounds variable */
    3.51 +	wrap_nc(nc_inq_varid(input_ncid, time_bounds_name, &time_bounds_varid));
    3.52 +	/* Set index for reading the first value of time_bounds from the
    3.53 +	 * first timeslice of the first file */
    3.54 +	time_bounds_index[0] = time_bounds_index[1] = 0;
    3.55 +	/* Read the value */
    3.56 +	wrap_nc(nc_get_var1_double(input_ncid, time_bounds_varid,
    3.57 +		time_bounds_index, &bnd1));
    3.58 +	/* If the first and last file are not the same, close the first one and
    3.59 +	 * open the second one */
    3.60 +	if (strcmp(first_fname, last_fname)) {
    3.61 +		/* Close the first input file */
    3.62 +		wrap_nc(nc_close(input_ncid));
    3.63 +		/* Open the last input file */
    3.64 +		wrap_nc(nc_open(last_fname, NC_NOWRITE, &input_ncid));
    3.65 +		/* Get dimension ID of the time dimension */
    3.66 +		wrap_nc(nc_inq_dimid(input_ncid, time_name, &time_dimid));
    3.67 +		/* Get length of time dimension */
    3.68 +		wrap_nc(nc_inq_dimlen(input_ncid, time_dimid, &time_len));
    3.69 +		/* Get variable ID of time_bounds variable */
    3.70 +		wrap_nc(nc_inq_varid(input_ncid, time_bounds_name,
    3.71 +			&time_bounds_varid));
    3.72 +	}
    3.73 +	/* Set index for reading the second value of time_bounds from the
    3.74 +	 * last timeslice of the last file */
    3.75 +	time_bounds_index[0] = time_len - 1; time_bounds_index[1] = 1;
    3.76 +	/* Read the value */
    3.77 +	wrap_nc(nc_get_var1_double(input_ncid, time_bounds_varid,
    3.78 +		time_bounds_index, &bnd2));
    3.79 +
    3.80 +	/* Close the last input file */
    3.81 +	wrap_nc(nc_close(input_ncid));
    3.82 +
    3.83 +	tbounds[0] = bnd1;
    3.84 +	tbounds[1] = bnd2;
    3.85 +
    3.86 +	return;
    3.87 +}
    3.88 +
    3.89 +void open_inout(char *first_fname, char *last_fname, char *mean_fname, char *stddev_fname, char *flist, size_t nsamples)
    3.90  {
    3.91  	char *mean_history_gatt, *stddev_history_gatt,
    3.92  		*mean_prefix="Statistical means from history files:",
    3.93  		*stddev_prefix="Statistical standard deviations from history files:";
    3.94 +	double tbounds[2];
    3.95  
    3.96  	/*
    3.97  	 * Construct strings for history global attributes for the two output
    3.98 @@ -866,8 +923,16 @@
    3.99  	sprintf(mean_history_gatt, "%s%s", mean_prefix, flist);
   3.100  	sprintf(stddev_history_gatt, "%s%s", stddev_prefix, flist);
   3.101  
   3.102 -	/* Open input file */
   3.103 -	wrap_nc(nc_open(input_fname, NC_NOWRITE, &input_ncid));
   3.104 +	/* The two time_bounds values must be handled explicitly by obtaining
   3.105 +	 * the first one from the first file and the last one from the last
   3.106 +	 * file.  These are then passed to copy_metadata() below. */
   3.107 +	get_time_bounds(first_fname, last_fname, tbounds);
   3.108 +#ifdef DEBUG
   3.109 +	fprintf(stderr, "Got back tbounds=%lf,%lf\n", tbounds[0], tbounds[1]);
   3.110 +#endif /* DEBUG */
   3.111 +
   3.112 +	/* Open last input file */
   3.113 +	wrap_nc(nc_open(last_fname, NC_NOWRITE, &input_ncid));
   3.114  	/* Inquire about number of dimensions, variables, global attributes
   3.115  	 * and the ID of the unlimited dimension
   3.116  	 */
   3.117 @@ -909,9 +974,9 @@
   3.118  	wrap_nc(nc_enddef(stddev_ncid));
   3.119  	/* Write out metdata variables that do not depend on "time" */
   3.120  	copy_metadata(input_ncid, input_var_head, input_dim_idx, mean_ncid,
   3.121 -		mean_var_head);
   3.122 +		mean_var_head, tbounds);
   3.123  	copy_metadata(input_ncid, input_var_head, input_dim_idx, stddev_ncid,
   3.124 -		stddev_var_head);
   3.125 +		stddev_var_head, tbounds);
   3.126  
   3.127  	wrap_nc(nc_close(input_ncid));
   3.128  
   3.129 @@ -1602,7 +1667,8 @@
   3.130  	 * only values from the last time slice from the period over which
   3.131  	 * the statistics are computed.
   3.132  	 */
   3.133 -	open_inout(ifnames[(nifnames-1)], mfname, sfname, flist, nsamples);
   3.134 +	open_inout(ifnames[0], ifnames[(nifnames-1)], mfname, sfname, flist,
   3.135 +		nsamples);
   3.136  
   3.137  	compute_stats(nifnames, ifnames, nsamples);
   3.138