write_data.c

超强的大尺度水文模拟工具

#include <stdio.h>
#include <stdlib.h>
#include <vicNl.h>

static char vcid[] = "$Id: write_data.c,v 4.1.2.4 2004/05/06 20:26:41 tbohn Exp $";

void write_data(out_data_struct *out_data,
		outfiles_struct *outfiles,
		dmy_struct      *dmy,
		int              dt)
/**********************************************************************
	write_data	Dag Lohmann		Janurary 1996

  This subroutine writes all energy and moisture balance parameters to
  output files.

  OUTPUT:
	evaporation and vapor fluxes in mm/time step
	layer moisture in mm/time step
	runoff in mm/time step
	baseflow in mm/time step
	freezing and thawing depths in cm
	snow depth in cm
	snow water equivlence in mm
	all energy fluxes are in W/m^2

  Modifications:
  5/20/96	Program was modified to account for a variable
		number of soil layers.  It was also modified to
		write out frozen soils data per time step.	KAC
  1/15/97	Program modified to output daily sums, or values
		independant of selected time step.  This aids in
		comparisons between model versions.		KAC
  3/98          Routine modified to output fluxes in PILPS2c 
                ASCII column format                             Dag
  4/30/98       Routine modified to add binary output options for
                improved file speed, and less disk usage for large
		model basins                                    KAC
  7/19/99       modified to output a single binary file containing
                the data selected for the LDAS project         KAC
  8/3/99        modified again to reduce the storage space needed
                for the LDAS output files.  
  1/4/2000      modified to allow both standard and LDAS formatted
                output using a compiler flag                    KAC
  3-12-03   added energy fluxes to snow band output files   KAC
  04-23-2003    modified LDAS SWQ output, so that it is multiplied by
                10 instead of 100 before being converted to a short
                integer.  This reduces stored value precision to 0.1,
                but increases the maximum storable SWQ, which was
                exceeded in previous LDAS simulations.          KAC

**********************************************************************/
{
  extern option_struct options;
#if LINK_DEBUG
  extern debug_struct debug;
#endif

  char               *tmp_cptr;
  short int          *tmp_siptr;
  unsigned short int *tmp_usiptr;
  int                 band, j;
  int                 fidx;
  int                *tmp_iptr;
  float              *tmp_fptr;

#if OPTIMIZE

  /*****************************************************************
    Create optimization output files
      type: ASCII
      columns: 5
               year
	       month
	       day
	       runoff
	       baseflow
               snow_depth
      comment: runoff and baseflow are output as daily sums for all
               defined model time steps.
  ******************************************************************/

  fprintf(outfiles->fluxes,"%04i\t%02i\t%02i\t%.4f\t%.4f\t%.4f\n",
	  dmy->year, dmy->month, dmy->day, out_data->runoff,
	  out_data->baseflow);
  
#elif LDAS_OUTPUT 

  /*************************************************************
    Write output files using LDAS binary format
    - single file, limited output variables, output compressed

  LDAS output file format:
  unsigned short int             year
  char                           month
  char                           day
 (char)                         (hour)
  unsigned short int             prec * 100
  short int                      evap * 100
  float                          runoff
  float                          baseflow
  unsigned short int * Nlayers   moist * 10
  unsigned short int             swq * 10
  short int                      net_short * 10
  short int                      in_long * 10
  short int                      r_net * 10
  short int                      latent * 10
  short int                      sensible * 10
  short int                      grnd_flux * 10
  unsigned short int             albedo * 10000
  short int                      surf_temp * 100
  unsigned short int             rel_humid * 100
  short int                      air_temp*100
  unsigned short int             wind*100
  if FROZEN_SOIL activated
    unsigned short int         ice[Nlayers] * 10     mean
    for each NUM_FRONTS:
      unsigned short int         fdepth[] * 100     mean
      unsigned short int         tdepth[] * 100     mean

  *************************************************************/

  tmp_cptr = (char *)calloc(1,sizeof(char));
  tmp_siptr = (short int *)calloc(1,sizeof(short int));
  tmp_usiptr = (unsigned short int *)calloc(1,sizeof(unsigned short int));
  tmp_iptr = (int *)calloc(1,sizeof(int));
  tmp_fptr = (float *)calloc(1,sizeof(float));

  /************************************
    Output Standard Energy and Moisture Flux Variables
  ************************************/

  /***** Write Binary Fluxes File *****/
  tmp_usiptr[0] = (unsigned short int)dmy->year;
  fwrite(tmp_usiptr,1,sizeof(unsigned short int),outfiles->fluxes);
  tmp_cptr[0] = (char)dmy->month;
  fwrite(tmp_cptr,1,sizeof(char),outfiles->fluxes);
  tmp_cptr[0] = (char)dmy->day;
  fwrite(tmp_cptr,1,sizeof(char),outfiles->fluxes);
  if(dt<24) {
    tmp_cptr[0] = (char)dmy->hour;
    fwrite(tmp_cptr,1,sizeof(char),outfiles->fluxes);
  }
  
  /** water balance components **/
  tmp_usiptr[0] = (unsigned short int)(out_data->prec * 100.);
  fwrite(tmp_usiptr,1,sizeof(unsigned short int),outfiles->fluxes);
  tmp_siptr[0] = (short int)(out_data->evap * 100.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);
  tmp_fptr[0] = (float)out_data->runoff;
  fwrite(tmp_fptr,1,sizeof(float),outfiles->fluxes);
  tmp_fptr[0] = (float)out_data->baseflow;
  fwrite(tmp_fptr,1,sizeof(float),outfiles->fluxes);
  for(j=0;j<options.Nlayer;j++) {    
    tmp_usiptr[0] = (unsigned short int)(out_data->moist[j]*10.);
    fwrite(tmp_usiptr,1,sizeof(unsigned short int),outfiles->fluxes);
  }
  tmp_usiptr[0] = (unsigned short int)(out_data->swq[0]*10.);
  fwrite(tmp_usiptr,1,sizeof(unsigned short int),outfiles->fluxes);
  
  /** energy balance components **/
  tmp_siptr[0] = (short int)(out_data->net_short*10.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);
  tmp_siptr[0] = (short int)(out_data->in_long*10.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);
  tmp_siptr[0] = (short int)(out_data->r_net*10.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);
  tmp_siptr[0] = (short int)(out_data->latent*10.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);
  tmp_siptr[0] = (short int)(out_data->sensible*10.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);
  tmp_siptr[0] = (short int)(out_data->grnd_flux*10.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);
  tmp_usiptr[0] = (unsigned short int)(out_data->albedo*10000.);
  fwrite(tmp_usiptr,1,sizeof(unsigned short int),outfiles->fluxes);
  tmp_siptr[0] = (short int)((out_data->rad_temp-KELVIN)*100.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);
  tmp_usiptr[0] = (unsigned short int)(out_data->rel_humid*100.);
  fwrite(tmp_usiptr,1,sizeof(unsigned short int),outfiles->fluxes);

  /* write air temperature */
  tmp_siptr[0] = (short int)(out_data->air_temp * 100.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);

  /* write wind speed */
  tmp_siptr[0] = (short int)(out_data->wind * 100.);
  fwrite(tmp_siptr,1,sizeof(short int),outfiles->fluxes);

  /* frozen soil depths */
  if(options.FROZEN_SOIL) {
    for(j=0;j<options.Nlayer;j++) {    
      tmp_usiptr[0] = (unsigned short int)(out_data->ice[j]*10.);
      fwrite(tmp_usiptr,1,sizeof(unsigned short int),outfiles->fluxes);
    }
    for(fidx=0;fidx<MAX_FRONTS;fidx++) {
      tmp_usiptr[0] = (unsigned short int)(out_data->fdepth[fidx]*100.);
      fwrite(tmp_usiptr,1,sizeof(unsigned short int),outfiles->fluxes);
      tmp_usiptr[0] = (unsigned short int)(out_data->tdepth[fidx]*100.);
      fwrite(tmp_usiptr,1,sizeof(unsigned short int),outfiles->fluxes);
    }
  }

  free((char *)tmp_cptr);
  free((char *)tmp_siptr);
  free((char *)tmp_usiptr);
  free((char *)tmp_iptr);
  free((char *)tmp_fptr);

#else

  /***************************************************************
    Write output files using default VIC ASCII or BINARY formats
    - multiple files, all variables, no truncation

    see VIC web page for format details:
      www.hydro.washington.edu/Lettenmaier/Models/VIC/VIChome.html
  ***************************************************************/

  tmp_iptr = (int *)calloc(1,sizeof(int));
  tmp_fptr = (float *)calloc(1,sizeof(float));

  /************************************
    Output Frozen Soil Variables
  ************************************/
  if(options.FROZEN_SOIL && options.BINARY_OUTPUT) {
    /***** Write Binary Frozen Soil Output File *****/
    tmp_iptr[0] = dmy->year;
    fwrite(tmp_iptr,1,sizeof(int),outfiles->fdepth);
    tmp_iptr[0] = dmy->month;
    fwrite(tmp_iptr,1,sizeof(int),outfiles->fdepth);
    tmp_iptr[0] = dmy->day;
    fwrite(tmp_iptr,1,sizeof(int),outfiles->fdepth);
    tmp_iptr[0] = dmy->hour;
    fwrite(tmp_iptr,1,sizeof(int),outfiles->fdepth);
    for(fidx=0;fidx<MAX_FRONTS;fidx++) {
      tmp_fptr[0] = (float)out_data->fdepth[fidx];
      fwrite(tmp_fptr,1,sizeof(float),outfiles->fdepth);
      tmp_fptr[0] = (float)out_data->tdepth[fidx];
      fwrite(tmp_fptr,1,sizeof(float),outfiles->fdepth);
    }
    for(j=0;j<options.Nlayer;j++) {
      tmp_fptr[0] = (float)(out_data->ice[j] + out_data->moist[j]);
      fwrite(tmp_fptr,1,sizeof(float),outfiles->fdepth);
    }
  }
  else if(options.FROZEN_SOIL) {
    /***** Write ASCII Frozen Soil Output File *****/
    fprintf(outfiles->fdepth  ,"%04i\t%02i\t%02i\t%02i",
	    dmy->year, dmy->month, dmy->day, dmy->hour);
    for(fidx=0;fidx<MAX_FRONTS;fidx++) {
      fprintf(outfiles->fdepth  ,"\t%.4f\t%.4f",out_data->fdepth[fidx], 
            out_data->tdepth[fidx]);
    }
    for(j=0;j<options.Nlayer;j++) {
      fprintf(outfiles->fdepth,"\t%f", out_data->ice[j] + out_data->moist[j]);
    }
    fprintf(outfiles->fdepth,"\n"); 
  }

  /**************************************
    Ouput Snow Variables
  **************************************/
  if(options.BINARY_OUTPUT) {
    /***** Write Binary Snow Output File *****/
    tmp_iptr[0] = dmy->year;
    fwrite(tmp_iptr,1,sizeof(int),outfiles->snow);
    tmp_iptr[0] = dmy->month;
    fwrite(tmp_iptr,1,sizeof(int),outfiles->snow);
    tmp_iptr[0] = dmy->day;
    fwrite(tmp_iptr,1,sizeof(int),outfiles->snow);
    if(dt<24) {
      tmp_iptr[0] = dmy->hour;
      fwrite(tmp_iptr,1,sizeof(int),outfiles->snow);
    }
    tmp_fptr[0] = (float)out_data->swq[0];
    fwrite(tmp_fptr,1,sizeof(float),outfiles->snow);
    tmp_fptr[0] = (float)out_data->snow_depth[0];
    fwrite(tmp_fptr,1,sizeof(float),outfiles->snow);
    tmp_fptr[0] = (float)out_data->snow_canopy[0];
    fwrite(tmp_fptr,1,sizeof(float),outfiles->snow);
/*     tmp_fptr[0] = (float)out_data->coverage[0]; */
/*     fwrite(tmp_fptr,1,sizeof(float),outfiles->snow); */
    if(options.FULL_ENERGY) {
      tmp_fptr[0] = (float)out_data->advection[0];
      fwrite(tmp_fptr,1,sizeof(float),outfiles->snow);
      tmp_fptr[0] = (float)out_data->deltaCC[0];
      fwrite(tmp_fptr,1,sizeof(float),outfiles->snow);
      tmp_fptr[0] = (float)out_data->snow_flux[0];
      fwrite(tmp_fptr,1,sizeof(float),outfiles->snow);
      tmp_fptr[0] = (float)out_data->refreeze_energy[0];