/*
 * SUMMARY:      ReadRadMap.c - Read radiation images
 * USAGE:        Part of DHSVM
 *
 * AUTHOR:       Bart Nijssen
 * ORG:          University of Washington, Department of Civil Engineering
 * E-MAIL:       nijssen@u.washington.edu
 * ORIG-DATE:    Apr-96
 * Last Change: Wed Jan 27 23:49:56 1999 by Bart Nijssen <nijssen@u.washington.edu>
 * DESCRIPTION:  Read radiation images
 * DESCRIP-END.
 * FUNCTIONS:    ReadRadMap()
 * COMMENTS:     
 */

#ifndef lint
static char vcid[] = "$Id: ReadRadMap.c,v 1.2 2001/03/31 00:10:32 waichler Exp $";
#endif /* lint */

#include <stdio.h>
#include <stdlib.h>
#include "settings.h"
#include "data.h"
#include "DHSVMerror.h"
#include "fileio.h"
#include "functions.h"
#include "constants.h"
#include "rad.h"
#include "sizeofnt.h"

/*****************************************************************************
  Function name: ReadRadMap()

  Purpose      : Read monthly radiation images

  Required     :
    Map      - DHSVM map structure
    RadMap   - address of array pointers to structure in which to store the 
               radiation values for each pixel
    FileName - name of input file

  Returns      : void

  Modifies     : 
    RadMap.BeamCode
    RadMap.DiffuseCode

  Comments     :
*****************************************************************************/

void ReadRadMap(MAPSIZE *Map, int EndStep, int NDaySteps, int NStats, 
                char *FileName, TOPOPIX **TopoMap, RADCLASSPIX **RadMap, 
                METLOCATION *Stat)
{
  const char *Routine = "ReadRadMap";
  int i;			/* counter */
  int x;			/* counter */
  int y;			/* counter */
  int NumberType;		/* number type */
  int NSet;			/* Number of dataset to be read */
  void *Array;

  if (DEBUG)
    printf("Reading radiation map from file: %s\n", FileName);

  NumberType = NC_INT;

  if (!(Array = (unsigned int *) calloc(Map->NY * Map->NX, 
                                        SizeOfNumberType(NumberType))))
    ReportError((char *) Routine, 1);    

  /* Read the encoded direct beam radiation */

  NSet = 0;

  Read2DMatrix(FileName, Array, NumberType, Map->NY, Map->NX, NSet++);

  for (y = 0, i = 0; y < Map->NY; y++) {
    for (x = 0; x < Map->NX; x++, i++) {
      if (InBasin(TopoMap[y][x].Mask)) {
        DecodeIPW(EndStep, NDaySteps, ((unsigned int *)Array)[i], 
                  RadMap[y][x].BeamCode);
      }
    }
  }

  for (i = 0; i < NStats; i++) {
    y = Stat[i].Loc.N;
    x = Stat[i].Loc.E;
    DecodeIPW(EndStep, NDaySteps, 
              ((unsigned int *)Array)[y*Map->NX+x], Stat[i].Data.BeamCode);
  }
  
  /* Read the encoded diffuse radiation */

  Read2DMatrix(FileName, Array, NumberType, Map->NY, Map->NX, NSet);

  for (y = 0, i = 0; y < Map->NY; y++) {
    for (x = 0; x < Map->NX; x++, i++) {
      if (InBasin(TopoMap[y][x].Mask)) {
        DecodeIPW(EndStep, NDaySteps, ((unsigned int *)Array)[i], 
                  RadMap[y][x].DiffuseCode);
      }
    }
  }

  for (i = 0; i < NStats; i++) {
    y = Stat[i].Loc.N;
    x = Stat[i].Loc.E;
    DecodeIPW(EndStep, NDaySteps, 
              ((unsigned int *)Array)[y*Map->NX+x], Stat[i].Data.DiffuseCode);
  }

  free(Array);
}

