From 9c15b0918872235c678a90d0c880941207de3eb4 Mon Sep 17 00:00:00 2001
From: KyleKlenk <kyle.c.klenk@gmail.com>
Date: Fri, 18 Aug 2023 16:26:15 -0600
Subject: [PATCH] fixed issue where output was not written properly. Also
 cleaned up the output strucutre and removed some of the calls to subroutines
 that are uneeded.

---
 .../file_access_actor/file_access_actor.hpp   |   2 -
 .../file_access_actor_subroutine_wrappers.hpp |   6 -
 .../cpp_code/file_access_actor.cpp            |  11 -
 .../fortran_code/cppwrap_fileAccess.f90       |   4 -
 .../fortran_code/output_structure.f90         | 168 ++++-----
 .../fortran_code/read_icondFromStructure.f90  | 280 --------------
 .../writeOutputFromOutputStructure.f90        |  40 +-
 .../fortran_code/write_to_netcdf.f90          |  80 ----
 .../hru_actor/fortran_code/hru_modelwrite.f90 |  32 +-
 .../hru_actor/fortran_code/hru_restart.f90    |  25 --
 .../hru_actor/fortran_code/hru_setup.f90      |  12 +-
 .../fortran_code/hru_writeOutput.f90          |   2 +-
 build/source/netcdf/modelwrite.f90            |  41 +--
 build/source/netcdf/read_icond.f90            | 344 ------------------
 14 files changed, 119 insertions(+), 928 deletions(-)
 delete mode 100644 build/source/actors/file_access_actor/fortran_code/read_icondFromStructure.f90
 delete mode 100644 build/source/netcdf/read_icond.f90

diff --git a/build/includes/file_access_actor/file_access_actor.hpp b/build/includes/file_access_actor/file_access_actor.hpp
index 76056bd..7a4e21c 100644
--- a/build/includes/file_access_actor/file_access_actor.hpp
+++ b/build/includes/file_access_actor/file_access_actor.hpp
@@ -94,8 +94,6 @@ struct file_access_state {
 behavior file_access_actor(stateful_actor<file_access_state>* self, int startGRU, int numGRU, 
    File_Access_Actor_Settings file_access_actor_settings, actor parent);
 
-// Read in the inital conditions for all the HRUs that are in the run-domain
-void readInitConditions(stateful_actor<file_access_state>* self);
 
 void initalizeOutputHandles(stateful_actor<file_access_state>* self);
 
diff --git a/build/includes/file_access_actor/file_access_actor_subroutine_wrappers.hpp b/build/includes/file_access_actor/file_access_actor_subroutine_wrappers.hpp
index eab0b13..002e4c9 100644
--- a/build/includes/file_access_actor/file_access_actor_subroutine_wrappers.hpp
+++ b/build/includes/file_access_actor/file_access_actor_subroutine_wrappers.hpp
@@ -18,12 +18,6 @@ extern "C" {
   void deallocateOutputStructure(int* err);
   void writeOutput_fortran(void* handle_ncid, int* num_steps, int* start_gru, int* max_gru, int* err);
 
-  // Set up global initial conditions
-  void openInitCondFile(int* init_cond_ncid, int* err);
-  void closeInitCondFile(int* init_cond_ncid, int* err);
-  void readInitCond_prog(int* init_cond_ncid, int* start_gru, int* num_gru, int* err);
-  void readInitCond_bvar(int* init_cond_ncid, int* start_gru, int* num_gru, int* err);
-
   void updateFailed(int* indxHRU);
 
   void resetFailedArray();
diff --git a/build/source/actors/file_access_actor/cpp_code/file_access_actor.cpp b/build/source/actors/file_access_actor/cpp_code/file_access_actor.cpp
index 4f91d8c..3b0b55a 100644
--- a/build/source/actors/file_access_actor/cpp_code/file_access_actor.cpp
+++ b/build/source/actors/file_access_actor/cpp_code/file_access_actor.cpp
@@ -44,9 +44,6 @@ behavior file_access_actor(stateful_actor<file_access_state>* self, int start_gr
 
     aout(self) << "Simluations Steps: " << self->state.num_steps << "\n";
 
-
-    // read in the inital conditions for the grus/hrus
-    readInitConditions(self);
     
     // Inital Files Have Been Loaded - Send Message to Job_Actor to Start Simulation
     self->send(self->state.parent, init_gru_v);
@@ -267,12 +264,4 @@ void writeOutput(stateful_actor<file_access_state>* self, Output_Partition* part
     partition->resetReadyToWriteList();
 }
 
-void readInitConditions(stateful_actor<file_access_state>* self) {
-    int err;
-    openInitCondFile(&self->state.init_cond_ncid, &err);
-    readInitCond_prog(&self->state.init_cond_ncid, &self->state.start_gru, &self->state.num_gru, &err);
-    readInitCond_bvar(&self->state.init_cond_ncid, &self->state.start_gru, &self->state.num_gru, &err);
-    closeInitCondFile(&self->state.init_cond_ncid, &err); 
-}
-
 } // end namespace
\ No newline at end of file
diff --git a/build/source/actors/file_access_actor/fortran_code/cppwrap_fileAccess.f90 b/build/source/actors/file_access_actor/fortran_code/cppwrap_fileAccess.f90
index 8937e37..f4b881f 100644
--- a/build/source/actors/file_access_actor/fortran_code/cppwrap_fileAccess.f90
+++ b/build/source/actors/file_access_actor/fortran_code/cppwrap_fileAccess.f90
@@ -401,8 +401,6 @@ subroutine FileAccessActor_DeallocateStructures(handle_forcFileInfo, handle_ncid
   USE globalData,only:forcingDataStruct
   USE globalData,only:vectime
   USE globalData,only:outputTimeStep
-  USE globalData,only:init_cond_prog
-  USE globalData,only:init_cond_bvar
   implicit none
   type(c_ptr),intent(in), value        :: handle_forcFileInfo
   type(c_ptr),intent(in), value        :: handle_ncid
@@ -430,8 +428,6 @@ subroutine FileAccessActor_DeallocateStructures(handle_forcFileInfo, handle_ncid
   deallocate(ncid)
   deallocate(failedHRUs)
   deallocate(outputTimeStep)
-  deallocate(init_cond_prog)
-  if (allocated(init_cond_bvar))then; deallocate(init_cond_bvar); endif;
 end subroutine FileAccessActor_DeallocateStructures
 
 
diff --git a/build/source/actors/file_access_actor/fortran_code/output_structure.f90 b/build/source/actors/file_access_actor/fortran_code/output_structure.f90
index a6757fe..c1e85b1 100644
--- a/build/source/actors/file_access_actor/fortran_code/output_structure.f90
+++ b/build/source/actors/file_access_actor/fortran_code/output_structure.f90
@@ -60,41 +60,41 @@ module output_structure_module
   private::is_var_desired
 
   type, public :: summa_output_type
-    type(gru_hru_z_vLookup)                          :: lookupStruct                  ! x%gru(:)%hru(:)%z(:)%var(:)%lookup(:) -- lookup tables
+    type(gru_hru_z_vLookup)                          :: lookupStruct                   ! x%gru(:)%hru(:)%z(:)%var(:)%lookup(:) -- lookup tables
 
     ! define the statistics structures
-    type(gru_hru_time_doubleVec),allocatable          :: forcStat(:)                   ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model forcing data
-    type(gru_hru_time_doubleVec),allocatable          :: progStat(:)                   ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model prognostic (state) variables
-    type(gru_hru_time_doubleVec),allocatable          :: diagStat(:)                   ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model diagnostic variables
-    type(gru_hru_time_doubleVec),allocatable          :: fluxStat(:)                   ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model fluxes
-    type(gru_hru_time_doubleVec),allocatable          :: indxStat(:)                   ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model indices
-    type(gru_hru_time_doubleVec),allocatable          :: bvarStat(:)                   ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- basin-average variabl
+    type(gru_hru_time_doubleVec)                      :: forcStat                      ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model forcing data
+    type(gru_hru_time_doubleVec)                      :: progStat                      ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model prognostic (state) variables
+    type(gru_hru_time_doubleVec)                      :: diagStat                      ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model diagnostic variables
+    type(gru_hru_time_doubleVec)                      :: fluxStat                      ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model fluxes
+    type(gru_hru_time_doubleVec)                      :: indxStat                      ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model indices
+    type(gru_hru_time_doubleVec)                      :: bvarStat                      ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- basin-average variabl
 
     ! define the primary data structures (scalars)
-    type(gru_hru_time_int),allocatable                :: timeStruct(:)                 ! x%gru(:)%hru(:)%var(:)%tim(:)     -- model time data
-    type(gru_hru_time_double),allocatable             :: forcStruct(:)                 ! x%gru(:)%hru(:)%var(:)%tim(:)     -- model forcing data
+    type(gru_hru_time_int)                            :: timeStruct                    ! x%gru(:)%hru(:)%var(:)%tim(:)     -- model time data
+    type(gru_hru_time_double)                         :: forcStruct                    ! x%gru(:)%hru(:)%var(:)%tim(:)     -- model forcing data
     type(gru_hru_double)                              :: attrStruct                    ! x%gru(:)%hru(:)%var(:)            -- local attributes for each HRU, DOES NOT CHANGE OVER TIMESTEPS
     type(gru_hru_int)                                 :: typeStruct                    ! x%gru(:)%hru(:)%var(:)            -- local classification of soil veg etc. for each HRU, DOES NOT CHANGE OVER TIMESTEPS
     type(gru_hru_int8)                                :: idStruct                      ! x%gru(:)%hru(:)%var(:)
 
     ! define the primary data structures (variable length vectors)
-    type(gru_hru_time_intVec),allocatable             :: indxStruct(:)                 ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model indices
+    type(gru_hru_time_intVec)                         :: indxStruct                    ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model indices
     type(gru_hru_intVec)                              :: indxStruct_init               ! x%gru(:)%hru(:)%var(:)%dat        -- model indices
     type(gru_hru_doubleVec)                           :: mparStruct                    ! x%gru(:)%hru(:)%var(:)%dat        -- model parameters, DOES NOT CHANGE OVER TIMESTEPS TODO: MAYBE
-    type(gru_hru_time_doubleVec),allocatable          :: progStruct(:)                 ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model prognostic (state) variables
+    type(gru_hru_time_doubleVec)                      :: progStruct                    ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model prognostic (state) variables
     type(gru_hru_doubleVec)                           :: progStruct_init               ! x%gru(:)%hru(:)%var(:)%dat        -- model prognostic (state) variables
-    type(gru_hru_time_doubleVec),allocatable          :: diagStruct(:)                 ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model diagnostic variables
-    type(gru_hru_time_doubleVec),allocatable          :: fluxStruct(:)                 ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model fluxes
+    type(gru_hru_time_doubleVec)                      :: diagStruct                    ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model diagnostic variables
+    type(gru_hru_time_doubleVec)                      :: fluxStruct                    ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- model fluxes
 
     ! define the basin-average structures
     type(gru_double)                                  :: bparStruct                    ! x%gru(:)%var(:)                   -- basin-average parameters, DOES NOT CHANGE OVER TIMESTEPS
-    type(gru_hru_time_doubleVec),allocatable          :: bvarStruct(:)                 ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- basin-average variables
+    type(gru_hru_time_doubleVec)                      :: bvarStruct                    ! x%gru(:)%hru(:)%var(:)%tim(:)%dat -- basin-average variables
     type(gru_doubleVec)                               :: bvarStruct_init               ! x%gru(:)%hru(:)%var(:)%dat        -- basin-average variables
     ! define the ancillary data structures
     type(gru_hru_double)                              :: dparStruct                    ! x%gru(:)%hru(:)%var(:)
 
     ! finalize stats structure
-    type(gru_hru_time_flagVec),allocatable            :: finalizeStats(:)              ! x%gru(:)%hru(:)%tim(:)%dat -- flags on when to write to file
+    type(gru_hru_time_flagVec)                        :: finalizeStats                 ! x%gru(:)%hru(:)%tim(:)%dat -- flags on when to write to file
 
 
     type(gru_d)                                       :: upArea
@@ -170,108 +170,62 @@ subroutine initOutputStructure(forcFileInfo, maxSteps, num_gru, err)
   if (.not.allocated(outputStructure))then
     allocate(outputStructure(1))
   else
-    print*, "Already Allocated"
-    return;
+    print*, "Already Allocated"; return;
   end if
 
   ! LookupStructure
-  ! allocate(outputStructure(1)%lookupStruct(1))
-  ! allocate(outputStructure(1)%lookupStruct(1)%gru(num_gru))
 
   ! Statistics Structures
-  allocate(outputStructure(1)%forcStat(1))
-  allocate(outputStructure(1)%progStat(1))
-  allocate(outputStructure(1)%diagStat(1))
-  allocate(outputStructure(1)%fluxStat(1))
-  allocate(outputStructure(1)%indxStat(1))
-  allocate(outputStructure(1)%bvarStat(1))
-  allocate(outputStructure(1)%forcStat(1)%gru(num_gru))
-  allocate(outputStructure(1)%progStat(1)%gru(num_gru))
-  allocate(outputStructure(1)%diagStat(1)%gru(num_gru))
-  allocate(outputStructure(1)%fluxStat(1)%gru(num_gru))
-  allocate(outputStructure(1)%indxStat(1)%gru(num_gru))
-  allocate(outputStructure(1)%bvarStat(1)%gru(num_gru))
-
+  allocate(outputStructure(1)%forcStat%gru(num_gru))
+  allocate(outputStructure(1)%progStat%gru(num_gru))
+  allocate(outputStructure(1)%diagStat%gru(num_gru))
+  allocate(outputStructure(1)%fluxStat%gru(num_gru))
+  allocate(outputStructure(1)%indxStat%gru(num_gru))
+  allocate(outputStructure(1)%bvarStat%gru(num_gru))
   ! Primary Data Structures (scalars)
-  allocate(outputStructure(1)%timeStruct(1))
-  allocate(outputStructure(1)%forcStruct(1))
-  ! allocate(outputStructure(1)%attrStruct)
-  ! allocate(outputStructure(1)%typeStruct(1))
-  ! allocate(outputStructure(1)%idStruct(1))
-  allocate(outputStructure(1)%timeStruct(1)%gru(num_gru))
-  allocate(outputStructure(1)%forcStruct(1)%gru(num_gru))
-  ! allocate(outputStructure(1)%attrStruct(1)%gru(num_gru))
-  ! allocate(outputStructure(1)%typeStruct(1)%gru(num_gru))
-  ! allocate(outputStructure(1)%idStruct(1)%gru(num_gru))
-  
+  allocate(outputStructure(1)%timeStruct%gru(num_gru))
+  allocate(outputStructure(1)%forcStruct%gru(num_gru))
   ! Primary Data Structures (variable length vectors)
-  allocate(outputStructure(1)%indxStruct(1))
-  ! allocate(outputStructure(1)%mparStruct(1))
-  allocate(outputStructure(1)%progStruct(1))
-  allocate(outputStructure(1)%diagStruct(1))
-  allocate(outputStructure(1)%fluxStruct(1))
-  allocate(outputStructure(1)%indxStruct(1)%gru(num_gru))
-  ! allocate(outputStructure(1)%mparStruct(1)%gru(num_gru))
-  allocate(outputStructure(1)%progStruct(1)%gru(num_gru))
-  allocate(outputStructure(1)%diagStruct(1)%gru(num_gru))
-  allocate(outputStructure(1)%fluxStruct(1)%gru(num_gru))
-
+  allocate(outputStructure(1)%indxStruct%gru(num_gru))
+  allocate(outputStructure(1)%progStruct%gru(num_gru))
+  allocate(outputStructure(1)%diagStruct%gru(num_gru))
+  allocate(outputStructure(1)%fluxStruct%gru(num_gru))
   ! Basin-Average structures
-  ! allocate(outputStructure(1)%bparStruct(1))
-  allocate(outputStructure(1)%bvarStruct(1))
-  ! allocate(outputStructure(1)%bvarStruct_init(1))
-  ! allocate(outputStructure(1)%bparStruct(1)%gru(num_gru))
-  allocate(outputStructure(1)%bvarStruct(1)%gru(num_gru))
-  ! allocate(outputStructure(1)%bvarStruct_init(1)%gru(num_gru))
-
-
-  ! define the ancillary data structures
-  ! allocate(outputStructure(1)%dparStruct(1))
-  ! allocate(outputStructure(1)%dparStruct(1)%gru(num_gru))
-
+  allocate(outputStructure(1)%bvarStruct%gru(num_gru))
   ! Finalize Stats for writing
-  allocate(outputStructure(1)%finalizeStats(1))
-  allocate(outputStructure(1)%finalizeStats(1)%gru(num_gru))
-
+  allocate(outputStructure(1)%finalizeStats%gru(num_gru))
+  ! Extras
   allocate(outputStructure(1)%upArea%gru(num_gru))
   
   
   do iGRU = 1, num_gru
     num_hru = gru_struc(iGRU)%hruCount
-    ! LookupStructure
-    ! allocate(outputStructure(1)%lookupStruct(1)%gru(iGRU)%hru(num_hru))
-
     ! Statistics Structures
-    allocate(outputStructure(1)%forcStat(1)%gru(iGRU)%hru(num_hru))
-    allocate(outputStructure(1)%progStat(1)%gru(iGRU)%hru(num_hru))
-    allocate(outputStructure(1)%diagStat(1)%gru(iGRU)%hru(num_hru))
-    allocate(outputStructure(1)%fluxStat(1)%gru(iGRU)%hru(num_hru))
-    allocate(outputStructure(1)%indxStat(1)%gru(iGRU)%hru(num_hru))
-    allocate(outputStructure(1)%bvarStat(1)%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%forcStat%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%progStat%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%diagStat%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%fluxStat%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%indxStat%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%bvarStat%gru(iGRU)%hru(num_hru))
 
     ! Primary Data Structures (scalars)
-    allocate(outputStructure(1)%timeStruct(1)%gru(iGRU)%hru(num_hru))
-    allocate(outputStructure(1)%forcStruct(1)%gru(iGRU)%hru(num_hru))
-    ! allocate(outputStructure(1)%attrStruct(1)%gru(iGRU)%hru(num_hru))
-    ! allocate(outputStructure(1)%typeStruct(1)%gru(iGRU)%hru(num_hru))
-    ! allocate(outputStructure(1)%idStruct(1)%gru(iGRU)%hru(num_hru))
-  
+    allocate(outputStructure(1)%timeStruct%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%forcStruct%gru(iGRU)%hru(num_hru))
+
     ! Primary Data Structures (variable length vectors)
-    allocate(outputStructure(1)%indxStruct(1)%gru(iGRU)%hru(num_hru))
-    ! allocate(outputStructure(1)%mparStruct(1)%gru(iGRU)%hru(num_hru))
-    allocate(outputStructure(1)%progStruct(1)%gru(iGRU)%hru(num_hru))
-    allocate(outputStructure(1)%diagStruct(1)%gru(iGRU)%hru(num_hru))
-    allocate(outputStructure(1)%fluxStruct(1)%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%indxStruct%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%progStruct%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%diagStruct%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%fluxStruct%gru(iGRU)%hru(num_hru))
   
     ! Basin-Average structures
-    allocate(outputStructure(1)%bvarStruct(1)%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%bvarStruct%gru(iGRU)%hru(num_hru))
 
 
    ! define the ancillary data structures
-    ! allocate(outputStructure(1)%dparStruct(1)%gru(iGRU)%hru(num_hru))
 
     ! Finalize Stats for writing
-    allocate(outputStructure(1)%finalizeStats(1)%gru(iGRU)%hru(num_hru))
+    allocate(outputStructure(1)%finalizeStats%gru(iGRU)%hru(num_hru))
 
     allocate(outputStructure(1)%upArea%gru(iGRU)%hru(num_hru))
 
@@ -315,14 +269,14 @@ subroutine initOutputStructure(forcFileInfo, maxSteps, num_gru, err)
         ! allocate space structures
           select case(trim(structInfo(iStruct)%structName))    
             case('time')
-              call alloc_outputStruc(time_meta,outputStructure(1)%timeStruct(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(time_meta,outputStructure(1)%timeStruct%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,err=err,message=message)     ! model forcing data
             case('forc')
               ! Structure
-              call alloc_outputStruc(forc_meta,outputStructure(1)%forcStruct(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(forc_meta,outputStructure(1)%forcStruct%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,message=message);    ! model forcing data
               ! Statistics
-              call alloc_outputStruc(statForc_meta(:)%var_info,outputStructure(1)%forcStat(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(statForc_meta(:)%var_info,outputStructure(1)%forcStat%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,message=message);    ! model forcing data
             case('attr'); cycle;
               ! call allocGlobal(attr_meta, outputStructure(1)%attrStruct(1)%gru(iGRU)%hru(iHRU), err, message)
@@ -331,39 +285,39 @@ subroutine initOutputStructure(forcFileInfo, maxSteps, num_gru, err)
             case('mpar'); cycle;
             case('indx')
               ! Structure
-              call alloc_outputStruc(indx_meta,outputStructure(1)%indxStruct(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(indx_meta,outputStructure(1)%indxStruct%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,str_name='indx',message=message);    ! model variables
               ! Statistics
-              call alloc_outputStruc(statIndx_meta(:)%var_info,outputStructure(1)%indxStat(1)%gru(iGRU)%hru(1), &
+              call alloc_outputStruc(statIndx_meta(:)%var_info,outputStructure(1)%indxStat%gru(iGRU)%hru(1), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,message=message);    ! index vars
             case('prog')
               ! Structure
-              call alloc_outputStruc(prog_meta,outputStructure(1)%progStruct(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(prog_meta,outputStructure(1)%progStruct%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,message=message);    ! model prognostic (state) variables
               ! Statistics
-              call alloc_outputStruc(statProg_meta(:)%var_info,outputStructure(1)%progStat(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(statProg_meta(:)%var_info,outputStructure(1)%progStat%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,message=message);    ! model prognostic 
             case('diag')
               ! Structure
-              call alloc_outputStruc(diag_meta,outputStructure(1)%diagStruct(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(diag_meta,outputStructure(1)%diagStruct%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,message=message);    ! model diagnostic variables
               ! Statistics
-              call alloc_outputStruc(statDiag_meta(:)%var_info,outputStructure(1)%diagStat(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(statDiag_meta(:)%var_info,outputStructure(1)%diagStat%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,message=message);    ! model diagnostic
             case('flux')
               ! Structure
-              call alloc_outputStruc(flux_meta,outputStructure(1)%fluxStruct(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(flux_meta,outputStructure(1)%fluxStruct%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,message=message);    ! model fluxes
               ! Statistics
-              call alloc_outputStruc(statFlux_meta(:)%var_info,outputStructure(1)%fluxStat(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(statFlux_meta(:)%var_info,outputStructure(1)%fluxStat%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=nSnow,nSoil=nSoil,err=err,message=message);    ! model fluxes
             case('bpar'); cycle;
             case('bvar')
               ! Structure
-              call alloc_outputStruc(bvar_meta,outputStructure(1)%bvarStruct(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(bvar_meta,outputStructure(1)%bvarStruct%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=0,nSoil=0,err=err,message=message);  ! basin-average variables
               ! Statistics
-              call alloc_outputStruc(statBvar_meta(:)%var_info,outputStructure(1)%bvarStat(1)%gru(iGRU)%hru(iHRU), &
+              call alloc_outputStruc(statBvar_meta(:)%var_info,outputStructure(1)%bvarStat%gru(iGRU)%hru(iHRU), &
                                       nSteps=maxSteps,nSnow=0,nSoil=0,err=err,message=message);  ! basin-average variables
             case('deriv');  cycle
             case('lookup'); cycle  ! lookup tables
@@ -379,9 +333,9 @@ subroutine initOutputStructure(forcFileInfo, maxSteps, num_gru, err)
       end do  ! looping through data structures
     
       ! Finalize stats structure for writing to output file
-      allocate(outputStructure(1)%finalizeStats(1)%gru(iGRU)%hru(iHRU)%tim(maxSteps))
+      allocate(outputStructure(1)%finalizeStats%gru(iGRU)%hru(iHRU)%tim(maxSteps))
       do iStep = 1, maxSteps
-        allocate(outputStructure(1)%finalizeStats(1)%gru(iGRU)%hru(iHRU)%tim(iStep)%dat(1:maxVarFreq))
+        allocate(outputStructure(1)%finalizeStats%gru(iGRU)%hru(iHRU)%tim(iStep)%dat(1:maxVarFreq))
       end do ! timeSteps
     end do ! Looping through GRUs
   end do
diff --git a/build/source/actors/file_access_actor/fortran_code/read_icondFromStructure.f90 b/build/source/actors/file_access_actor/fortran_code/read_icondFromStructure.f90
deleted file mode 100644
index 6eabadd..0000000
--- a/build/source/actors/file_access_actor/fortran_code/read_icondFromStructure.f90
+++ /dev/null
@@ -1,280 +0,0 @@
-module read_icondFromStructure_module
-USE, intrinsic :: iso_c_binding
-USE nrtype
-
-
-implicit none
-private
-public :: openInitCondFile
-public :: closeInitCondFile
-public :: readInitCond_prog
-public :: readInitCond_bvar
-
-! define single HRU restart file
-integer(i4b), parameter :: singleHRU=1001
-integer(i4b), parameter :: multiHRU=1002
-integer(i4b), parameter :: restartFileType=multiHRU
-contains
-
-subroutine openInitCondFile(init_cond_ncid, err) bind(C, name="openInitCondFile")
-  USE netcdf_util_module,only:nc_file_open               ! open netcdf file
-  USE netcdf
-  ! file paths
-  USE summaFileManager,only:SETTINGS_PATH           ! path to settings files (e.g., Noah vegetation tables)
-  USE summaFileManager,only:STATE_PATH              ! optional path to state/init. condition files (defaults to SETTINGS_PATH)
-  USE summaFileManager,only:MODEL_INITCOND          ! name of model initial conditions file
-  implicit none
-
-  integer(c_int),intent(out)          :: init_cond_ncid
-  integer(c_int),intent(out)          :: err
-
-  character(len=256)                  :: init_cond_fileName
-  character(len=256)                  :: message          ! error message
-  character(len=1024)                 :: cmessage         ! error message for downwind routine
-  ! --------------------------------------------------------------------------------------------------------
-
-  ! define restart file path/name
-  if(STATE_PATH == '') then
-    init_cond_fileName = trim(SETTINGS_PATH)//trim(MODEL_INITCOND)
-  else
-    init_cond_fileName = trim(STATE_PATH)//trim(MODEL_INITCOND)
-  endif
-
-  call nc_file_open(init_cond_fileName,nf90_nowrite,init_cond_ncid,err,cmessage)
-  if (err/=0) then
-    message=trim(message)//trim(cmessage)
-    print(message)
-    return
-  end if
-
-end subroutine openInitCondFile
-
-subroutine closeInitCondFile(init_cond_ncid,err) bind(C, name="closeInitCondFile")
-  USE netcdf_util_module,only:nc_file_close
-  implicit none
-  
-  integer(c_int),intent(out)          :: init_cond_ncid
-  integer(c_int),intent(out)          :: err
-  
-  ! local variables
-  character(len=256)                  :: message
-  ! --------------------------------------------------------------------------------------------------------
-  err=0; message="read_icondFromStructure.f90 - closeInitCondFile/"
-
-  call nc_file_close(init_cond_ncid,err,message)
-  if (err/=0) then
-    message=trim(message)
-    print(message)
-    return
-  end if
-
-end subroutine closeInitCondFile
-
-subroutine readInitCond_prog(init_cond_ncid, start_gru, num_gru, err) bind(C, name="readInitCond_prog")
-  ! Structures to populate
-  USE globalData,only:init_cond_prog
-  ! Netcdf
-  USE netcdf
-  USE netcdf_util_module,only:nc_file_close              ! close netcdf file
-  USE netcdf_util_module,only:netcdf_err                 ! netcdf error handling
-
-  ! metadata
-  USE globalData,only:prog_meta                     ! metadata for prognostic variables
-  USE globalData,only:bvar_meta                     ! metadata for basin (GRU) variables
-  ! var_lookup
-  USE var_lookup,only:iLookVarType                  ! variable lookup structure
-  USE var_lookup,only:iLookPROG                     ! variable lookup structure
-  USE var_lookup,only:iLookPARAM                    ! variable lookup structure
-  USE var_lookup,only:iLookINDEX                    ! variable lookup structure
-  ! access type string
-  USE get_ixName_module,only:get_varTypeName        ! to access type strings for error messages
-  implicit none
-
-  integer(c_int), intent(in)             :: init_cond_ncid
-  integer(c_int), intent(in)             :: start_gru
-  integer(c_int), intent(in)             :: num_gru
-  integer(c_int), intent(out)            :: err
- 
-  ! local variables
-  integer(i4b)                           :: iVar
-  integer(i4b)                           :: dimID        ! varible dimension ids
-  integer(i4b)                           :: ncVarID      ! variable ID in netcdf file
-  integer(i4b)                           :: dimLen       ! data dimensions
-  character(256)                         :: dimName      ! not used except as a placeholder in call to inq_dim function
-  integer(i4b)                           :: fileHRU      ! number of HRUs in file
-
-  character(LEN=256)                     :: icond_file  ! restart/icond_file file name
-  character(len=256)                     :: message
-  character(len=256)                     :: cmessage
-
-  character(len=32),parameter            :: scalDimName   ='scalarv'  ! dimension name for scalar data
-  character(len=32),parameter            :: midSoilDimName='midSoil'  ! dimension name for soil-only layers
-  character(len=32),parameter            :: midTotoDimName='midToto'  ! dimension name for layered varaiables
-  character(len=32),parameter            :: ifcTotoDimName='ifcToto'  ! dimension name for layered varaiables
-  ! --------------------------------------------------------------------------------------------------------
-
-  err=0; message="read_icondFromStructure.f90 - readInitCond_prog"
-  
-  ! get number of HRUs in file
-  err = nf90_inq_dimid(init_cond_ncid,"hru",dimID);
-  if(err/=nf90_noerr)then; message=trim(message)//'problem finding hru dimension/'//trim(nf90_strerror(err)); return; end if
-  err = nf90_inquire_dimension(init_cond_ncid,dimID,len=fileHRU);
-  if(err/=nf90_noerr)then; message=trim(message)//'problem reading hru dimension/'//trim(nf90_strerror(err)); return; end if
-
-  allocate(init_cond_prog(size(prog_meta)))
-
-  ! loop through prognostic variables
-  do iVar = 1,size(prog_meta)
-    ! skip variables that are computed later
-    if(prog_meta(iVar)%varName=='scalarCanopyWat'           .or. &
-       prog_meta(iVar)%varName=='spectralSnowAlbedoDiffuse' .or. &
-       prog_meta(iVar)%varName=='scalarSurfaceTemp'         .or. &
-       prog_meta(iVar)%varName=='mLayerVolFracWat'          .or. &
-       prog_meta(iVar)%varName=='mLayerHeight'                   )then 
-       cycle
-    endif
-
-    ! get variable id
-    err = nf90_inq_varid(init_cond_ncid,trim(prog_meta(iVar)%varName),ncVarID); call netcdf_err(err,message)
-    if(err/=0)then
-      message=trim(message)//': problem with getting variable id, var='//trim(prog_meta(iVar)%varName)
-      print*,message
-      return
-    endif
-
-    ! get variable dimension IDs
-    select case (prog_meta(iVar)%varType)
-      case (iLookVarType%scalarv); err = nf90_inq_dimid(init_cond_ncid,trim(scalDimName)   ,dimID)
-        call netcdf_err(err,message)
-      case (iLookVarType%midSoil); err = nf90_inq_dimid(init_cond_ncid,trim(midSoilDimName),dimID)
-        call netcdf_err(err,message)
-      case (iLookVarType%midToto); err = nf90_inq_dimid(init_cond_ncid,trim(midTotoDimName),dimID)
-        call netcdf_err(err,message)
-      case (iLookVarType%ifcToto); err = nf90_inq_dimid(init_cond_ncid,trim(ifcTotoDimName),dimID)
-        call netcdf_err(err,message)
-    case default
-      message=trim(message)//"unexpectedVariableType[name='"//trim(prog_meta(iVar)%varName)//"';type='"//trim(get_varTypeName(prog_meta(iVar)%varType))//"']"
-      print*, message
-      err=20; return
-    end select
-
-    ! check errors
-    if(err/=0)then
-      message=trim(message)//': problem with dimension ids, var='//trim(prog_meta(iVar)%varName)
-      print*, message
-      return
-    endif
-
-    ! get the dimension length
-    err = nf90_inquire_dimension(init_cond_ncid,dimID,dimName,dimLen); call netcdf_err(err,message)
-    if(err/=0)then;message=trim(message)//': problem getting the dimension length';print*, message;return;endif
-
-    ! initialize the variable data
-    allocate(init_cond_prog(iVar)%var_data(num_gru,dimLen),stat=err)
-    if(err/=0)then;message=trim(message)//'problem allocating HRU variable data';print*, message;return;endif
-
-     ! get data
-    err = nf90_get_var(init_cond_ncid,ncVarID,init_cond_prog(iVar)%var_data, start=(/start_gru,1/),count=(/num_gru,dimLen/)) 
-    call netcdf_err(err,message)
-    if(err/=0)then; message=trim(message)//': problem getting the data for variable '//trim(prog_meta(iVar)%varName); return; endif
-
-  end do
-
-end subroutine readInitCond_prog
-
-subroutine readInitCond_bvar(init_cond_ncid, start_gru, num_gru, err) bind(C, name="readInitCond_bvar")
-  USE globalData,only:init_cond_bvar
-  USE globalData,only: nTimeDelay   ! number of hours in the time delay histogram
-  ! var_lookup
-  USE var_lookup,only:iLookBVAR                     ! variable lookup structure
-  ! metadata structures
-  USE globalData,only:bvar_meta                     ! metadata for basin (GRU) variables
-  ! netcdf
-  USE netcdf
-  USE netcdf_util_module,only:netcdf_err                 ! netcdf error handling
-  implicit none
-  
-  integer(c_int), intent(in)             :: init_cond_ncid
-  integer(c_int), intent(in)             :: start_gru
-  integer(c_int), intent(in)             :: num_gru
-  integer(c_int), intent(out)            :: err
-
-  ! local variables
-  integer(i4b)                           :: nTDH          ! number of points in time-delay histogram
-  integer(i4b)                           :: dimID        ! varible dimension ids
-  integer(i4b)                           :: fileGRU      ! number of GRUs in file
-  integer(i4b)                           :: i          
-  integer(i4b)                           :: iVar     
-  integer(i4b)                           :: dimLen       ! data dimensions
-  character(256)                         :: dimName      ! not used except as a placeholder in call to inq_dim function
-  integer(i4b)                           :: ncVarID      ! variable ID in netcdf file
-  integer(i4b),dimension(1)              :: ndx          ! intermediate array of loop indices
-   
-
-  character(len=256)                     :: message
-  
-  character(len=32),parameter            :: tdhDimName    ='tdh'      ! dimension name for time-delay basin variables
-
-  ! --------------------------------------------------------------------------------------------------------
-  err = 0; message="read_icondFromStructure.f90 - readInitCond_bvar/"
-  if(restartFileType/=singleHRU)then
-    ! get dimension of time delay histogram (TDH) from initial conditions file
-    err = nf90_inq_dimid(init_cond_ncid,"tdh",dimID);
-    
-    if(err/=nf90_noerr)then
-      write(*,*) 'WARNING: routingRunoffFuture is not in the initial conditions file ... using zeros'  ! previously created in var_derive.f90
-      err=nf90_noerr    ! reset this err
-    
-    else 
-      ! the state file *does* have the basin variable(s), so process them
-      err = nf90_inquire_dimension(init_cond_ncid,dimID,len=nTDH);
-      if(err/=nf90_noerr)then
-        message=trim(message)//'problem reading tdh dimension from initial condition file/'//trim(nf90_strerror(err))
-        print*, message
-        return
-      end if
-      
-      ! get number of GRUs in file
-      err = nf90_inq_dimid(init_cond_ncid,"gru",dimID);               if(err/=nf90_noerr)then; message=trim(message)//'problem finding gru dimension/'//trim(nf90_strerror(err)); return; end if
-      err = nf90_inquire_dimension(init_cond_ncid,dimID,len=fileGRU); if(err/=nf90_noerr)then; message=trim(message)//'problem reading gru dimension/'//trim(nf90_strerror(err)); return; end if
-
-      ! check vs hardwired value set in globalData.f90
-      if(nTDH /= nTimeDelay)then
-        write(*,*) 'tdh=',nTDH,' nTimeDelay=',nTimeDelay
-        message=trim(message)//': state file time delay dimension tdh does not match summa expectation of nTimeDelay set in globalData()'
-        return
-      endif
-
-      ndx = (/iLookBVAR%routingRunoffFuture/)   ! array of desired variable indices
-      allocate(init_cond_bvar(size(ndx)))
-      do i = 1, size(ndx)
-        iVar = ndx(i)
-        ! get tdh dimension Id in file (should be 'tdh')
-        err = nf90_inq_dimid(init_cond_ncid,trim(tdhDimName), dimID);
-        if(err/=0)then;message=trim(message)//': problem with dimension ids for tdh vars';print*,message;return;endif
-
-        ! get the tdh dimension length (dimName and dimLen are outputs of this call)
-        err = nf90_inquire_dimension(init_cond_ncid,dimID,dimName,dimLen); call netcdf_err(err,message)
-        if(err/=0)then;message=trim(message)//': problem getting the dimension length for tdh vars';print*,message;return;endif
-
-        ! get tdh-based variable id
-        err = nf90_inq_varid(init_cond_ncid,trim(bvar_meta(iVar)%varName),ncVarID); call netcdf_err(err,message)
-        if(err/=0)then; message=trim(message)//': problem with getting basin variable id, var='//trim(bvar_meta(iVar)%varName); return; endif
-
-        allocate(init_cond_bvar(i)%var_data(num_gru,dimLen),stat=err)
-        if(err/=0)then; print*, 'err= ',err; message=trim(message)//'problem allocating GRU variable data'; return; endif
-
-        ! get data
-        err = nf90_get_var(init_cond_ncid,ncVarID,init_cond_bvar(i)%var_data, start=(/start_gru,1/),count=(/num_gru,dimLen/)); call netcdf_err(err,message)
-        if(err/=0)then; message=trim(message)//': problem getting the data'; return; endif
-      end do
-    endif
-  endif
-
-
-end subroutine readInitCond_bvar
-
-
-
-
-end module read_icondFromStructure_module
\ No newline at end of file
diff --git a/build/source/actors/file_access_actor/fortran_code/writeOutputFromOutputStructure.f90 b/build/source/actors/file_access_actor/fortran_code/writeOutputFromOutputStructure.f90
index 3e159d5..364d145 100644
--- a/build/source/actors/file_access_actor/fortran_code/writeOutputFromOutputStructure.f90
+++ b/build/source/actors/file_access_actor/fortran_code/writeOutputFromOutputStructure.f90
@@ -126,7 +126,7 @@ subroutine writeOutput_fortran(handle_ncid, num_steps, start_gru, max_gru, err)
     stepCounter(:) = outputTimeStep(iGRU)%dat(:) ! We want to avoid updating outputTimeStep
     do iStep=1, num_steps
       call writeTime(ncid,outputTimeStep(iGRU)%dat(:),iStep,time_meta, &
-              outputStructure(1)%timeStruct(1)%gru(iGRU)%hru(indxHRU)%var,err,cmessage)
+              outputStructure(1)%timeStruct%gru(iGRU)%hru(indxHRU)%var,err,cmessage)
     end do ! istep
   end do ! iGRU
 
@@ -137,7 +137,7 @@ subroutine writeOutput_fortran(handle_ncid, num_steps, start_gru, max_gru, err)
   ! ****************************************************************************
   call writeBasin(ncid,outputTimeStep(start_gru)%dat(:),outputTimeStepUpdate,num_steps,&
                   start_gru, max_gru, numGRU, &
-                  bvar_meta,outputStructure(1)%bvarStat(1),outputStructure(1)%bvarStruct(1), &
+                  bvar_meta,outputStructure(1)%bvarStat,outputStructure(1)%bvarStruct, &
                   bvarChild_map,err,cmessage)
 
   ! ****************************************************************************
@@ -148,28 +148,28 @@ subroutine writeOutput_fortran(handle_ncid, num_steps, start_gru, max_gru, err)
       case('forc')
         call writeData(ncid,outputTimeStep(start_gru)%dat(:),outputTimestepUpdate,maxLayers,num_steps,&
                         start_gru, max_gru, numGRU, & 
-                        forc_meta,outputStructure(1)%forcStat(1),outputStructure(1)%forcStruct(1),'forc', &
-                        forcChild_map,outputStructure(1)%indxStruct(1),err,cmessage)
+                        forc_meta,outputStructure(1)%forcStat,outputStructure(1)%forcStruct,'forc', &
+                        forcChild_map,outputStructure(1)%indxStruct,err,cmessage)
       case('prog')
         call writeData(ncid,outputTimeStep(start_gru)%dat(:),outputTimestepUpdate,maxLayers,num_steps,&
                         start_gru, max_gru, numGRU, &
-                        prog_meta,outputStructure(1)%progStat(1),outputStructure(1)%progStruct(1),'prog', &
-                        progChild_map,outputStructure(1)%indxStruct(1),err,cmessage)
+                        prog_meta,outputStructure(1)%progStat,outputStructure(1)%progStruct,'prog', &
+                        progChild_map,outputStructure(1)%indxStruct,err,cmessage)
       case('diag')
         call writeData(ncid,outputTimeStep(start_gru)%dat(:),outputTimestepUpdate,maxLayers,num_steps,&
                         start_gru, max_gru, numGRU, &
-                        diag_meta,outputStructure(1)%diagStat(1),outputStructure(1)%diagStruct(1),'diag', &
-                        diagChild_map,outputStructure(1)%indxStruct(1),err,cmessage)
+                        diag_meta,outputStructure(1)%diagStat,outputStructure(1)%diagStruct,'diag', &
+                        diagChild_map,outputStructure(1)%indxStruct,err,cmessage)
       case('flux')
         call writeData(ncid,outputTimeStep(start_gru)%dat(:),outputTimestepUpdate,maxLayers,num_steps,&
                         start_gru, max_gru, numGRU, &
-                        flux_meta,outputStructure(1)%fluxStat(1),outputStructure(1)%fluxStruct(1),'flux', &
-                        fluxChild_map,outputStructure(1)%indxStruct(1),err,cmessage)
+                        flux_meta,outputStructure(1)%fluxStat,outputStructure(1)%fluxStruct,'flux', &
+                        fluxChild_map,outputStructure(1)%indxStruct,err,cmessage)
       case('indx')
         call writeData(ncid,outputTimeStep(start_gru)%dat(:),outputTimestepUpdate,maxLayers,num_steps,&
                         start_gru, max_gru, numGRU, &
-                        indx_meta,outputStructure(1)%indxStat(1),outputStructure(1)%indxStruct(1),'indx', &
-                        indxChild_map,outputStructure(1)%indxStruct(1),err,cmessage)
+                        indx_meta,outputStructure(1)%indxStat,outputStructure(1)%indxStruct,'indx', &
+                        indxChild_map,outputStructure(1)%indxStruct,err,cmessage)
     end select
     if(err/=0)then; message=trim(message)//trim(cmessage)//'['//trim(structInfo(iStruct)%structName)//']'; return; endif
   end do  ! (looping through structures)
@@ -325,9 +325,9 @@ subroutine writeData(ncid,outputTimestep,outputTimestepUpdate,maxLayers,nSteps,
 
         do iStep = 1, nSteps
           ! check if we want this timestep
-          if(.not.outputStructure(1)%finalizeStats(1)%gru(verifiedGRUIndex)%hru(1)%tim(iStep)%dat(iFreq)) cycle
+          if(.not.outputStructure(1)%finalizeStats%gru(verifiedGRUIndex)%hru(1)%tim(iStep)%dat(iFreq)) cycle
           stepCounter = stepCounter+1
-          timeVec(stepCounter) = outputStructure(1)%forcStruct(1)%gru(verifiedGRUIndex)%hru(1)%var(iVar)%tim(iStep)
+          timeVec(stepCounter) = outputStructure(1)%forcStruct%gru(verifiedGRUIndex)%hru(1)%var(iVar)%tim(iStep)
         end do ! iStep
         err = nf90_put_var(ncid%var(iFreq),ncVarID,timeVec(1:stepCounter),start=(/outputTimestep(iFreq)/),count=(/stepCounter/))
         call netcdf_err(err,message); if (err/=0)then; print*, "err"; return; endif
@@ -396,7 +396,7 @@ subroutine writeScalar(ncid, outputTimestep, outputTimestepUpdate, nSteps, minGR
         stepCounter = 0
         gruCounter = gruCounter + 1
         do iStep = 1, nSteps
-          if(.not.outputStructure(1)%finalizeStats(1)%gru(iGRU)%hru(1)%tim(iStep)%dat(iFreq)) cycle
+          if(.not.outputStructure(1)%finalizeStats%gru(iGRU)%hru(1)%tim(iStep)%dat(iFreq)) cycle
           stepCounter = stepCounter + 1
           realVec(gruCounter, stepCounter) = stat%gru(iGRU)%hru(1)%var(map(iVar))%tim(iStep)%dat(iFreq)
           outputTimeStepUpdate(iFreq) = stepCounter
@@ -412,9 +412,8 @@ subroutine writeScalar(ncid, outputTimestep, outputTimestepUpdate, nSteps, minGR
         print*, "   maxGRU = ", maxGRU
         err = 20
         return
-
-      err = nf90_put_var(ncid%var(iFreq),meta(iVar)%ncVarID(iFreq),realVec(1:gruCounter, 1:stepCounter),start=(/minGRU,outputTimestep(iFreq)/),count=(/numGRU,stepCounter/))
       endif
+      err = nf90_put_var(ncid%var(iFreq),meta(iVar)%ncVarID(iFreq),realVec(1:gruCounter, 1:stepCounter),start=(/minGRU,outputTimestep(iFreq)/),count=(/numGRU,stepCounter/))
     class default; err=20; message=trim(message)//'stats must be scalarv and of type gru_hru_doubleVec'; return
   end select  ! stat
 
@@ -492,11 +491,11 @@ subroutine writeVector(ncid, outputTimestep, maxLayers, nSteps, minGRU, maxGRU,
       ! get the data vectors
       select type (dat)
           class is (gru_hru_time_doubleVec)
-              if(.not.outputStructure(1)%finalizeStats(1)%gru(iGRU)%hru(1)%tim(iStep)%dat(iFreq)) cycle
+              if(.not.outputStructure(1)%finalizeStats%gru(iGRU)%hru(1)%tim(iStep)%dat(iFreq)) cycle
               realArray(gruCounter,1:datLength) = dat%gru(iGRU)%hru(1)%var(iVar)%tim(iStep)%dat(:)
 
           class is (gru_hru_time_intVec)
-              if(.not.outputStructure(1)%finalizeStats(1)%gru(iGRU)%hru(1)%tim(iStep)%dat(iFreq)) cycle
+              if(.not.outputStructure(1)%finalizeStats%gru(iGRU)%hru(1)%tim(iStep)%dat(iFreq)) cycle
               intArray(gruCounter,1:datLength) = dat%gru(iGRU)%hru(1)%var(iVar)%tim(iStep)%dat(:)
           class default; err=20; message=trim(message)//'data must not be scalarv and either of type gru_hru_doubleVec or gru_hru_intVec'; return
       end select
@@ -567,6 +566,7 @@ subroutine writeBasin(ncid,outputTimestep,outputTimestepUpdate,nSteps,&
   ! initialize error control
   err=0;message="f-writeBasin/"
 
+
   ! loop through output frequencies
   do iFreq=1,maxvarFreq
 
@@ -632,7 +632,7 @@ subroutine writeTime(ncid,outputTimestep,iStep,meta,dat,err,message)
   do iFreq=1,maxvarFreq
 
     ! check that we have finalized statistics for a given frequency
-    if(.not.outputStructure(1)%finalizeStats(1)%gru(1)%hru(1)%tim(iStep)%dat(iFreq)) cycle
+    if(.not.outputStructure(1)%finalizeStats%gru(1)%hru(1)%tim(iStep)%dat(iFreq)) cycle
 
     ! loop through model variables
     do iVar = 1,size(meta)
diff --git a/build/source/actors/file_access_actor/fortran_code/write_to_netcdf.f90 b/build/source/actors/file_access_actor/fortran_code/write_to_netcdf.f90
index 959fd3e..4dc4668 100644
--- a/build/source/actors/file_access_actor/fortran_code/write_to_netcdf.f90
+++ b/build/source/actors/file_access_actor/fortran_code/write_to_netcdf.f90
@@ -8,9 +8,6 @@ USE data_types
 
 implicit none
 public::writeParamToNetCDF
-public::writeDataToNetCDF
-public::writeBasinToNetCDF
-public::writeTimeToNetCDF
 public::writeGRUStatistics
 
 contains
@@ -199,83 +196,6 @@ subroutine writeDataToNetCDF(handle_ncid,           &
   end do  ! (looping through structures)
 end subroutine writeDataToNetCDF
 
-subroutine writeBasinToNetCDF(handle_ncid, index_gru, handle_finalize_stats, &
-  handle_output_timestep, handle_bvar_stat, handle_bvar_struct, err) bind(C, name="writeBasinToNetCDF")
-  USE modelwrite_module,only:writeBasin
-  USE globalData,only:bvar_meta                 ! metadata on basin-average variables
-  USE globalData,only:bvarChild_map             ! index of the child data structure: stats bvar
- 
-  implicit none
-  ! dummy variables
-  type(c_ptr),    intent(in), value  :: handle_ncid
-  integer(c_int), intent(in)         :: index_gru
-  type(c_ptr),    intent(in), value  :: handle_finalize_stats
-  type(c_ptr),    intent(in), value  :: handle_output_timestep
-  type(c_ptr),    intent(in), value  :: handle_bvar_stat
-  type(c_ptr),    intent(in), value  :: handle_bvar_struct
-  integer(c_int), intent(out)        :: err
-  ! local pointers for dummy variables
-  type(var_i), pointer               :: ncid
-  type(flagVec), pointer             :: finalize_stats
-  type(var_i),pointer                :: output_timestep
-  type(var_dlength),pointer          :: bvar_stat
-  type(var_dlength),pointer          :: bvar_struct
-  ! local Variables
-  character(len=256)                 :: message
-  ! ---------------------------------------------------------------------------------------
-  ! * Convert From C++ to Fortran
-  ! ---------------------------------------------------------------------------------------
-  call c_f_pointer(handle_ncid, ncid)
-  call c_f_pointer(handle_finalize_stats, finalize_stats)
-  call c_f_pointer(handle_output_timestep, output_timestep)
-  call c_f_pointer(handle_bvar_stat, bvar_stat)
-  call c_f_pointer(handle_bvar_struct, bvar_struct)
-  message="file_access_actor.f90 - writeBasinToNetCDF"
-
-
-  call writeBasin(ncid,index_gru,finalize_stats%dat(:),output_timestep%var(:),&
-      bvar_meta, bvar_stat%var, bvar_struct%var, bvarChild_map, err, message)
-  if(err/=0)then 
-    message=trim(message)//'[bvar]'
-    print*, message
-    return
-  endif 
-
-end subroutine writeBasinToNetCDF
-
-subroutine writeTimeToNetCDF(handle_ncid, handle_finalize_stats, handle_output_timestep, &
-    handle_time_struct, err) bind(C, name="writeTimeToNetCDF")
-  USE modelwrite_module,only:writeTime
-  USE globalData,only:time_meta
-
-  implicit none
-  type(c_ptr), intent(in), value :: handle_ncid
-  type(c_ptr), intent(in), value :: handle_finalize_stats
-  type(c_ptr), intent(in), value :: handle_output_timestep
-  type(c_ptr), intent(in), value :: handle_time_struct
-  integer(c_int), intent(out)    :: err
-
-  type(var_i), pointer           :: ncid
-  type(flagVec), pointer         :: finalize_stats
-  type(var_i),pointer            :: output_timestep
-  type(var_i),pointer            :: time_struct
-  character(len=256)             :: message
-
-  call c_f_pointer(handle_ncid, ncid)
-  call c_f_pointer(handle_finalize_stats, finalize_stats)
-  call c_f_pointer(handle_output_timestep, output_timestep)
-  call c_f_pointer(handle_time_struct, time_struct)
-  message="file_access_actor.f90 - writeTimeToNetCDF"
-
-  call writeTime(ncid, finalize_stats%dat(:),output_timestep%var(:),&
-      time_meta, time_struct%var,err,message)
-  if(err/=0)then 
-    message=trim(message)//'writeTime'
-    print*, message
-    return
-  endif 
-
-end subroutine writeTimeToNetCDF
 
 subroutine writeGRUStatistics(handle_ncid,      &
                               gru_var_ids,      &
diff --git a/build/source/actors/hru_actor/fortran_code/hru_modelwrite.f90 b/build/source/actors/hru_actor/fortran_code/hru_modelwrite.f90
index 49e906e..96c9568 100755
--- a/build/source/actors/hru_actor/fortran_code/hru_modelwrite.f90
+++ b/build/source/actors/hru_actor/fortran_code/hru_modelwrite.f90
@@ -213,7 +213,7 @@ subroutine writeData(indxGRU,indxHRU,iStep,structName,finalizeStats, &
         ! Write the time step values
         select type(dat)      ! forcStruc
           class is (var_d)    ! x%var(:)
-            outputStructure(1)%forcStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep) = dat%var(iVar)
+            outputStructure(1)%forcStruct%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep) = dat%var(iVar)
           class default; err=20; message=trim(message)//'time variable must be of type var_d (forcing data structure)'; return
         end select
       end if  ! id time
@@ -230,15 +230,15 @@ subroutine writeData(indxGRU,indxHRU,iStep,structName,finalizeStats, &
           class is (var_dlength)
             select case(trim(structName))
             case('forc')
-              outputStructure(1)%forcStat(1)%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
+              outputStructure(1)%forcStat%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
             case('prog')
-              outputStructure(1)%progStat(1)%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
+              outputStructure(1)%progStat%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
             case('diag')
-              outputStructure(1)%diagStat(1)%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
+              outputStructure(1)%diagStat%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
             case('flux')
-              outputStructure(1)%fluxStat(1)%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
+              outputStructure(1)%fluxStat%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
             case('indx')
-              outputStructure(1)%indxStat(1)%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
+              outputStructure(1)%indxStat%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat%var(map(iVar))%dat(iFreq)
             case default
               err=21; message=trim(message)//"Stats structure not found"; return
             end select
@@ -250,11 +250,11 @@ subroutine writeData(indxGRU,indxHRU,iStep,structName,finalizeStats, &
 
         ! get the model layers
         nSoil   = indx%var(iLookIndex%nSoil)%dat(1)
-        outputStructure(1)%indxStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iLookIndex%nSoil)%tim(iStep)%dat(1) = nSoil
+        outputStructure(1)%indxStruct%gru(indxGRU)%hru(indxHRU)%var(iLookIndex%nSoil)%tim(iStep)%dat(1) = nSoil
         nSnow   = indx%var(iLookIndex%nSnow)%dat(1)
-        outputStructure(1)%indxStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iLookIndex%nSnow)%tim(iStep)%dat(1) = nSnow
+        outputStructure(1)%indxStruct%gru(indxGRU)%hru(indxHRU)%var(iLookIndex%nSnow)%tim(iStep)%dat(1) = nSnow
         nLayers = indx%var(iLookIndex%nLayers)%dat(1)
-        outputStructure(1)%indxStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iLookIndex%nLayers)%tim(iStep)%dat(1) = nLayers
+        outputStructure(1)%indxStruct%gru(indxGRU)%hru(indxHRU)%var(iLookIndex%nLayers)%tim(iStep)%dat(1) = nLayers
 
         ! get the length of each data vector
         select case (meta(iVar)%varType)
@@ -273,16 +273,16 @@ subroutine writeData(indxGRU,indxHRU,iStep,structName,finalizeStats, &
           class is (var_dlength)
             select case(trim(structName))
               case('prog')
-                outputStructure(1)%progStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(:) = dat%var(iVar)%dat(:)
+                outputStructure(1)%progStruct%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(:) = dat%var(iVar)%dat(:)
               case('diag')
-                outputStructure(1)%diagStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(:) = dat%var(iVar)%dat(:)
+                outputStructure(1)%diagStruct%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(:) = dat%var(iVar)%dat(:)
               case('flux')
-                outputStructure(1)%fluxStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(:) = dat%var(iVar)%dat(:)
+                outputStructure(1)%fluxStruct%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(:) = dat%var(iVar)%dat(:)
               case default
                 err=21; message=trim(message)//'data structure not found for output'
             end select
           class is (var_ilength) 
-            outputStructure(1)%indxStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(:) = dat%var(iVar)%dat(:)
+            outputStructure(1)%indxStruct%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(:) = dat%var(iVar)%dat(:)
           class default; err=20; message=trim(message)//'data must not be scalarv and either of type var_dlength or var_ilength'; return
         end select
 
@@ -362,10 +362,10 @@ subroutine writeBasin(indxGRU,indxHRU,iStep,finalizeStats,&
    select case (meta(iVar)%varType)
 
     case (iLookVarType%scalarv)
-      outputStructure(1)%bvarStat(1)%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat(map(iVar))%dat(iFreq)
+      outputStructure(1)%bvarStat%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat(map(iVar))%dat(iFreq)
     case (iLookVarType%routing)
      if (iFreq==1 .and. outputTimestep(iFreq)==1) then
-      outputStructure(1)%bvarStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(iFreq) = dat(iVar)%dat(iFreq)
+      outputStructure(1)%bvarStruct%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(iFreq) = dat(iVar)%dat(iFreq)
      end if
 
     case default
@@ -417,7 +417,7 @@ subroutine writeTime(indxGRU,indxHRU,iStep,finalizeStats,meta,dat,err,message)
    if (meta(iVar)%statIndex(iFreq)/=iLookStat%inst) cycle
 
    ! add to outputStructure
-   outputStructure(1)%timeStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep) = dat(iVar)
+   outputStructure(1)%timeStruct%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep) = dat(iVar)
    if (err/=0) message=trim(message)//trim(meta(iVar)%varName)
    if (err/=0) then; err=20; return; end if
 
diff --git a/build/source/actors/hru_actor/fortran_code/hru_restart.f90 b/build/source/actors/hru_actor/fortran_code/hru_restart.f90
index a894a0b..cc604d0 100644
--- a/build/source/actors/hru_actor/fortran_code/hru_restart.f90
+++ b/build/source/actors/hru_actor/fortran_code/hru_restart.f90
@@ -131,31 +131,6 @@ subroutine summa_readRestart(&
 
   ! initialize error control
   err=0; message='hru_actor_readRestart/'
-  nGRU = 1
-
-  ! *****************************************************************************
-  ! *** read/check initial conditions
-  ! *****************************************************************************
-  ! read initial conditions
-  ! call read_icond(&
-  !                 indxGRU,                       & ! intent(in):    index of GRU in gru_struc
-  !                 indxHRU,                       & ! intent(in):    index of HRU in gru_struc
-  !                 mparStruct,                    & ! intent(in):    model parameters
-  !                 progStruct,                    & ! intent(inout): model prognostic variables
-  !                 bvarStruct,                    & ! intent(inout): model basin (GRU) variables
-  !                 indxStruct,                    & ! intent(inout): model indices
-  !                 err,cmessage)                    ! intent(out):   error control
-  ! if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
-
-  ! ! check initial conditions
-  ! call check_icond4chm(&
-  !                   indxGRU,                        & ! intent(in):   index of GRU in gru_struc
-  !                   indxHRU,                        & ! intent(in):   index of HRU in gru_struc           
-  !                   progStruct,                     & ! intent(in):   model prognostic (state) variables
-  !                   mparStruct,                     & ! intent(in):   model parameters
-  !                   indxStruct,                     & ! intent(in):   layer indexes
-  !                   err,cmessage)                   ! intent(out):   error control
-  ! if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
 
 
   ! *****************************************************************************
diff --git a/build/source/actors/hru_actor/fortran_code/hru_setup.f90 b/build/source/actors/hru_actor/fortran_code/hru_setup.f90
index 4c4c650..c9edd4c 100644
--- a/build/source/actors/hru_actor/fortran_code/hru_setup.f90
+++ b/build/source/actors/hru_actor/fortran_code/hru_setup.f90
@@ -195,12 +195,14 @@ subroutine setupHRUParam(&
   do ivar=1, size(outputStructure(1)%bvarStruct_init%gru(indxGRU)%var(:))
     bvarStruct%var(ivar)%dat(:) = outputStructure(1)%bvarStruct_init%gru(indxGRU)%var(ivar)%dat(:)
   enddo
-  ! Copy the lookup Struct
-  do i_z=1, size(outputStructure(1)%lookupStruct%gru(indxGRU)%hru(indxHRU)%z(:))
-    do iVar=1, size(outputStructure(1)%lookupStruct%gru(indxGRU)%hru(indxHRU)%z(i_z)%var(:))
-      lookupStruct%z(i_z)%var(ivar)%lookup(:) = outputStructure(1)%lookupStruct%gru(indxGRU)%hru(indxHRU)%z(i_z)%var(iVar)%lookup(:)
+  ! Copy the lookup Struct if its allocated
+  if (allocated(outputStructure(1)%lookupStruct%gru(indxGRU)%hru(indxHRU)%z)) then
+    do i_z=1, size(outputStructure(1)%lookupStruct%gru(indxGRU)%hru(indxHRU)%z(:))
+      do iVar=1, size(outputStructure(1)%lookupStruct%gru(indxGRU)%hru(indxHRU)%z(i_z)%var(:))
+        lookupStruct%z(i_z)%var(ivar)%lookup(:) = outputStructure(1)%lookupStruct%gru(indxGRU)%hru(indxHRU)%z(i_z)%var(iVar)%lookup(:)
+      end do
     end do
-  end do
+  endif
   ! Copy the progStruct_init
   do ivar=1, size(outputStructure(1)%progStruct_init%gru(indxGRU)%hru(indxHRU)%var(:))
     progStruct%var(ivar)%dat(:) = outputStructure(1)%progStruct_init%gru(indxGRU)%hru(indxHRU)%var(ivar)%dat(:)
diff --git a/build/source/actors/hru_actor/fortran_code/hru_writeOutput.f90 b/build/source/actors/hru_actor/fortran_code/hru_writeOutput.f90
index 3fd07da..91f9d67 100644
--- a/build/source/actors/hru_actor/fortran_code/hru_writeOutput.f90
+++ b/build/source/actors/hru_actor/fortran_code/hru_writeOutput.f90
@@ -265,7 +265,7 @@ subroutine writeHRUToOutputStructure(&
   if(err/=0)then; message=trim(message)//trim(cmessage); return; endif
 
  ! If we do not do this looping we segfault - I am not sure why
-  outputStructure(1)%finalizeStats(1)%gru(indxGRU)%hru(indxHRU)%tim(outputStep)%dat(:) = finalizeStats%dat(:)
+  outputStructure(1)%finalizeStats%gru(indxGRU)%hru(indxHRU)%tim(outputStep)%dat(:) = finalizeStats%dat(:)
 
  ! ****************************************************************************
  ! *** calculate output statistics
diff --git a/build/source/netcdf/modelwrite.f90 b/build/source/netcdf/modelwrite.f90
index 36a5900..a4f45e6 100644
--- a/build/source/netcdf/modelwrite.f90
+++ b/build/source/netcdf/modelwrite.f90
@@ -205,32 +205,21 @@ subroutine writeData(ncid, finalize_stats, output_timestep, max_layers, index_gr
         ! get variable index
         err = nf90_inq_varid(ncid%var(iFreq),trim(meta(iVar)%varName),ncVarID)
         call netcdf_err(err,message)
-        if (err/=0) then
-          print*, message
-          return
+        if (err/=0) then; print*, message; return
         endif
 
         select type(dat)
           class is(var_d)
             err = nf90_put_var(ncid%var(iFreq),ncVarID,dat%var(iVar),start=(/output_timestep(iFreq)/))
             call netcdf_err(err,message)
-            if (err/=0) then
-              print*, message
-              return
-            endif
+            if (err/=0) then; print*, message; return; endif
             cycle
-            class default
-              err=20
-              message=trim(message)//'time variable must be of type var_dlength (forcing data structure)'
-              print*, message
-              return
+          class default
+            err=20;message=trim(message)//'time variable must be of type var_dlength (forcing data structure)';print*, message;return
         end select
 
         call netcdf_err(err,message)
-        if (err/=0) then
-          print*, message
-          return
-        endif
+        if (err/=0) then; print*, message;return;endif
       endif
 
       ! define the statistics index
@@ -245,17 +234,10 @@ subroutine writeData(ncid, finalize_stats, output_timestep, max_layers, index_gr
           class is (var_dlength)
             realVec(1) = stat%var(map(iVar))%dat(iFreq)
             err = nf90_put_var(ncid%var(iFreq),meta(iVar)%ncVarID(iFreq),realVec,start=(/index_gru,output_timestep(iFreq)/),count=(/num_gru,1/))
-            if (err/=0) then
-              print*, message
-              return
-            endif
+            if (err/=0) then; print*, message; return; endif
           class default
-            err=20
-            message=trim(message)//'stats must be scalarv and of type var_dlength'
-            print*, message
-            return
+            err=20; message=trim(message)//'stats must be scalarv and of type var_dlength'; print*, message; return
         end select ! stat
-
       else
 
          ! Write the data
@@ -366,15 +348,19 @@ subroutine writeBasin(ncid,iGRU,finalizeStats,outputTimestep,meta,stat,dat,map,e
   ! initialize error control
   err=0;message="f-writeBasin/"
 
+  print*, "WE should see this"
+
   ! loop through output frequencies
   do iFreq=1,maxvarFreq
 
     ! skip frequencies that are not needed
+    print*, "Before outputFreq"
     if(.not.outFreq(iFreq)) cycle
-
+    print*, "After outputFreq"
     ! check that we have finalized statistics for a given frequency
+    print*, "Before stats"
     if(.not.finalizeStats(iFreq)) cycle
-
+    print*, "After stats"
     ! loop through model variables
     do iVar = 1,size(meta)
 
@@ -388,6 +374,7 @@ subroutine writeBasin(ncid,iGRU,finalizeStats,outputTimestep,meta,stat,dat,map,e
       select case (meta(iVar)%varType)
 
         case (iLookVarType%scalarv)
+          print*, "output", stat(map(iVar))%dat(iFreq)
           err = nf90_put_var(ncid%var(iFreq),meta(iVar)%ncVarID(iFreq),(/stat(map(iVar))%dat(iFreq)/),start=(/iGRU,outputTimestep(iFreq)/),count=(/1,1/))
 
         case (iLookVarType%routing)
diff --git a/build/source/netcdf/read_icond.f90 b/build/source/netcdf/read_icond.f90
deleted file mode 100644
index c95af1c..0000000
--- a/build/source/netcdf/read_icond.f90
+++ /dev/null
@@ -1,344 +0,0 @@
-! SUMMA - Structure for Unifying Multiple Modeling Alternatives
-! Copyright (C) 2014-2020 NCAR/RAL; University of Saskatchewan; University of Washington
-!
-! This file is part of SUMMA
-!
-! For more information see: http://www.ral.ucar.edu/projects/summa
-!
-! This program is free software: you can redistribute it and/or modify
-! it under the terms of the GNU General Public License as published by
-! the Free Software Foundation, either version 3 of the License, or
-! (at your option) any later version.
-!
-! This program 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 General Public License
-! along with this program.  If not, see <http://www.gnu.org/licenses/>.
-
-module read_icond_actors_module
-USE, intrinsic :: iso_c_binding
-USE nrtype
-USE netcdf
-USE globalData,only: ixHRUfile_min,ixHRUfile_max
-USE globalData,only: nTimeDelay   ! number of hours in the time delay histogram
-implicit none
-private
-public::read_icond
-public::read_icond_nlayers
-! define single HRU restart file
-integer(i4b), parameter :: singleHRU=1001
-integer(i4b), parameter :: multiHRU=1002
-integer(i4b), parameter :: restartFileType=multiHRU
-contains
-
- ! ************************************************************************************************
- ! public subroutine read_icond_nlayers: read model initial conditions file for number of snow/soil layers
- ! ************************************************************************************************
-subroutine read_icond_nlayers(iconFile,nGRU,indx_meta,err,message)
-  ! --------------------------------------------------------------------------------------------------------
-  ! modules
-  USE nrtype
-  USE var_lookup,only:iLookIndex                        ! variable lookup structure
-  USE globalData,only:gru_struc                         ! gru-hru mapping structures
-  USE netcdf_util_module,only:nc_file_close             ! close netcdf file
-  USE netcdf_util_module,only:nc_file_open              ! close netcdf file
-  USE netcdf_util_module,only:netcdf_err                ! netcdf error handling
-  USE data_types,only:gru_hru_intVec                    ! actual data
-  USE data_types,only:var_info                          ! metadata
-  implicit none
-
-  ! --------------------------------------------------------------------------------------------------------
-  ! variable declarations
-  ! dummies
-  character(*)        ,intent(in)     :: iconFile       ! name of input (restart) file
-  integer(i4b)        ,intent(in)     :: nGRU           ! total # of GRUs in run domain
-  type(var_info)      ,intent(in)     :: indx_meta(:)   ! metadata
-  integer(i4b)        ,intent(out)    :: err            ! error code
-  character(*)        ,intent(out)    :: message        ! returned error message
-
-  ! locals
-  integer(i4b)                        :: ncID                       ! netcdf file id
-  integer(i4b)                        :: dimID                      ! netcdf file dimension id
-  integer(i4b)                        :: fileHRU                    ! number of HRUs in netcdf file
-  integer(i4b)                        :: snowID, soilID             ! netcdf variable ids
-  integer(i4b)                        :: iGRU, iHRU                 ! loop indexes
-  integer(i4b)                        :: iHRU_global                ! index of HRU in the netcdf file
-  integer(i4b),allocatable            :: snowData(:)                ! number of snow layers in all HRUs
-  integer(i4b),allocatable            :: soilData(:)                ! number of soil layers in all HRUs  
-  character(len=256)                  :: cmessage                   ! downstream error message
-
-  ! --------------------------------------------------------------------------------------------------------
-  ! initialize error message
-  err=0
-  message = 'read_icond_nlayers/'
-
-  ! open netcdf file
-  call nc_file_open(iconFile,nf90_nowrite,ncid,err,cmessage);
-  if (err/=0) then; message=trim(message)//trim(cmessage); return; end if
-  
-
-  ! get number of HRUs in file (the GRU variable(s), if present, are processed at the end)
-  err = nf90_inq_dimid(ncID,"hru",dimId);               if(err/=nf90_noerr)then; message=trim(message)//'problem finding hru dimension/'//trim(nf90_strerror(err)); return; end if
-  err = nf90_inquire_dimension(ncID,dimId,len=fileHRU); if(err/=nf90_noerr)then; message=trim(message)//'problem reading hru dimension/'//trim(nf90_strerror(err)); return; end if
-
-  ! allocate storage for reading from file (allocate entire file size, even when doing subdomain run)
-  allocate(snowData(fileHRU))
-  allocate(soilData(fileHRU))
-  snowData = 0
-  soilData = 0
-
-  ! get netcdf ids for the variables holding number of snow and soil layers in each hru
-  err = nf90_inq_varid(ncid,trim(indx_meta(iLookIndex%nSnow)%varName),snowid); call netcdf_err(err,message)
-  err = nf90_inq_varid(ncid,trim(indx_meta(iLookIndex%nSoil)%varName),soilid); call netcdf_err(err,message)
-
-  ! get nSnow and nSoil data (reads entire state file)
-  err = nf90_get_var(ncid,snowid,snowData); call netcdf_err(err,message)
-  err = nf90_get_var(ncid,soilid,soilData); call netcdf_err(err,message)
-
-  ixHRUfile_min=huge(1)
-  ixHRUfile_max=0
-  ! find the min and max hru indices in the state file
-  do iGRU = 1,nGRU
-  do iHRU = 1,gru_struc(iGRU)%hruCount
-  if(gru_struc(iGRU)%hruInfo(iHRU)%hru_nc < ixHRUfile_min) ixHRUfile_min = gru_struc(iGRU)%hruInfo(iHRU)%hru_nc
-  if(gru_struc(iGRU)%hruInfo(iHRU)%hru_nc > ixHRUfile_max) ixHRUfile_max = gru_struc(iGRU)%hruInfo(iHRU)%hru_nc
-  end do
-  end do
-  
-
-  ! loop over grus in current run to update snow/soil layer information
-  do iGRU = 1,nGRU
-  do iHRU = 1,gru_struc(iGRU)%hruCount
-  iHRU_global = gru_struc(iGRU)%hruInfo(iHRU)%hru_nc
-
-  ! single HRU (Note: 'restartFileType' is hardwired above to multiHRU)
-  if(restartFileType==singleHRU) then
-      gru_struc(iGRU)%hruInfo(iHRU)%nSnow = snowData(1)
-      gru_struc(iGRU)%hruInfo(iHRU)%nSoil = soilData(1)
-
-  ! multi HRU
-  else
-      gru_struc(iGRU)%hruInfo(iHRU)%nSnow = snowData(iHRU_global)
-      gru_struc(iGRU)%hruInfo(iHRU)%nSoil = soilData(iHRU_global)
-  endif
-
-  end do
-  end do
-
-
-  ! close file
-  call nc_file_close(ncid,err,cmessage)
-  if(err/=0)then;message=trim(message)//trim(cmessage);return;end if
-
-  ! cleanup
-  deallocate(snowData,soilData)
-
-end subroutine read_icond_nlayers
-
-
-! ************************************************************************************************
-! public subroutine read_icond: read model initial conditions
-! ************************************************************************************************
-subroutine read_icond(&
-                    indxGRU,                       & ! intent(in):    Index of GRU in gru_struc
-                    indxHRU,                       & ! intent(in):    Index of HRU in gru_struc
-                    mparData,                      & ! intent(in):    model parameters
-                    progData,                      & ! intent(inout): model prognostic variables
-                    bvarData,                      & ! intent(inout): model basin (GRU) variables
-                    indxData,                      & ! intent(inout): model indices
-                    err,message)                     ! intent(out):   error control
-  ! --------------------------------------------------------------------------------------------------------
-  ! modules
-  USE nrtype
-  USE var_lookup,only:iLookVarType                       ! variable lookup structure
-  USE var_lookup,only:iLookPROG                          ! variable lookup structure
-  USE var_lookup,only:iLookPARAM                         ! variable lookup structure
-  USE var_lookup,only:iLookBVAR                          ! variable lookup structure
-  USE var_lookup,only:iLookINDEX                         ! variable lookup structure
-  USE globalData,only:prog_meta                          ! metadata for prognostic variables
-  USE globalData,only:bvar_meta                          ! metadata for basin (GRU) variables
-  USE globalData,only:gru_struc                          ! gru-hru mapping structures
-  USE globaldata,only:iname_soil,iname_snow              ! named variables to describe the type of layer
-  USE data_types,only:var_ilength                        ! full integer structure
-  USE data_types,only:var_dlength                        ! double precision structure for a single HRU
-  USE data_types,only:var_info                           ! metadata
-  USE get_ixName_module,only:get_varTypeName             ! to access type strings for error messages
-  USE updatState_module,only:updateSoil                  ! update soil states
-
-  USE netcdf
-
-  USE globalData,only:init_cond_prog
-  USE globalData,only:init_cond_bvar
-
-
-  implicit none
-
-  ! --------------------------------------------------------------------------------------------------------
-  ! variable declarations
-  ! dummies
-  integer(i4b)     ,intent(in)             :: indxGRU      ! index of GRU in gru_struc
-  integer(i4b)     ,intent(in)             :: indxHRU      ! index of HRU in hru_struc
-  type(var_dlength),intent(in)             :: mparData     ! model parameters
-  type(var_dlength),intent(inout)          :: progData     ! model prognostic variables
-  type(var_dlength),intent(inout)          :: bvarData     ! model basin (GRU) variables
-  type(var_ilength),intent(inout)          :: indxData     ! model indices
-  integer(i4b)     ,intent(out)            :: err          ! error code
-  character(*)     ,intent(out)            :: message      ! returned error message
-
-  ! locals
-  character(len=256)                     :: cmessage     ! downstream error message
-  integer(i4b)                           :: fileHRU      ! number of HRUs in file
-  integer(i4b)                           :: fileGRU      ! number of GRUs in file
-  integer(i4b)                           :: iVar, i      ! loop indices
-  integer(i4b),dimension(1)              :: ndx          ! intermediate array of loop indices
-  integer(i4b)                           :: dimID        ! varible dimension ids
-  integer(i4b)                           :: ncVarID      ! variable ID in netcdf file
-  character(256)                         :: dimName      ! not used except as a placeholder in call to inq_dim function
-  integer(i4b)                           :: dimLen       ! data dimensions
-  integer(i4b)                           :: ncID         ! netcdf file ID
-  integer(i4b)                           :: ixFile       ! index in file
-  integer(i4b)                           :: nSoil, nSnow, nToto ! # layers
-  integer(i4b)                           :: nTDH          ! number of points in time-delay histogram
-  integer(i4b)                           :: iLayer,jLayer ! layer indices
-  integer(i4b),parameter                 :: nBand=2       ! number of spectral bands
-
-  character(len=32),parameter            :: scalDimName   ='scalarv'  ! dimension name for scalar data
-  character(len=32),parameter            :: midSoilDimName='midSoil'  ! dimension name for soil-only layers
-  character(len=32),parameter            :: midTotoDimName='midToto'  ! dimension name for layered varaiables
-  character(len=32),parameter            :: ifcTotoDimName='ifcToto'  ! dimension name for layered varaiables
-  character(len=32),parameter            :: tdhDimName    ='tdh'      ! dimension name for time-delay basin variables
-
-  ! --------------------------------------------------------------------------------------------------------
-
-  ! Start procedure here
-  err=0; message="read_icondActors.f90 - read_icond/"
-
-  ! loop through prognostic variables
-  do iVar = 1,size(prog_meta)
-
-    ! skip variables that are computed later
-    if(prog_meta(iVar)%varName=='scalarCanopyWat'          .or. &
-      prog_meta(iVar)%varName=='spectralSnowAlbedoDiffuse' .or. &
-      prog_meta(iVar)%varName=='scalarSurfaceTemp'         .or. &
-      prog_meta(iVar)%varName=='mLayerVolFracWat'          .or. &
-      prog_meta(iVar)%varName=='mLayerHeight'                   ) cycle
-
-
-    ! get the number of layers
-    nSnow = gru_struc(indxGRU)%hruInfo(indxHRU)%nSnow
-    nSoil = gru_struc(indxGRU)%hruInfo(indxHRU)%nSoil
-    nToto = nSnow + nSoil
-
-    ! put the data into data structures and check that none of the values are set to nf90_fill_double
-    select case (prog_meta(iVar)%varType)
-      case (iLookVarType%scalarv)
-        progData%var(iVar)%dat(1)       = init_cond_prog(iVar)%var_data(indxGRU,1)
-        if(abs(progData%var(iVar)%dat(1) - nf90_fill_double) < epsilon(init_cond_prog(iVar)%var_data))then; err=20; endif
-      case (iLookVarType%midSoil)
-        progData%var(iVar)%dat(1:nSoil) = init_cond_prog(iVar)%var_data(indxGRU,1:nSoil)
-        if(any(abs(progData%var(iVar)%dat(1:nSoil) - nf90_fill_double) < epsilon(init_cond_prog(iVar)%var_data)))then; err=20; endif
-      case (iLookVarType%midToto)
-        progData%var(iVar)%dat(1:nToto) = init_cond_prog(iVar)%var_data(indxGRU,1:nToto)
-        if(any(abs(progData%var(iVar)%dat(1:nToto) - nf90_fill_double) < epsilon(init_cond_prog(iVar)%var_data)))then; err=20; endif
-      case (iLookVarType%ifcToto)
-        progData%var(iVar)%dat(0:nToto) = init_cond_prog(iVar)%var_data(indxGRU,1:nToto+1)
-        if(any(abs(progData%var(iVar)%dat(0:nToto) - nf90_fill_double) < epsilon(init_cond_prog(iVar)%var_data)))then; err=20; endif
-      case default
-        message=trim(message)//"unexpectedVariableType[name='"//trim(prog_meta(iVar)%varName)//"';type='"//trim(get_varTypeName(prog_meta(iVar)%varType))//"']"
-        print*,message
-        err=20; return
-    end select
-
-    if(err==20)then; 
-      message=trim(message)//"data set to the fill value (name='"//trim(prog_meta(iVar)%varName)//"')"; 
-      print*, message
-      return;
-    endif
-
-    ! fix the snow albedo
-    if(progData%var(iLookPROG%scalarSnowAlbedo)%dat(1) < 0._dp)then
-      progData%var(iLookPROG%scalarSnowAlbedo)%dat(1) = mparData%var(iLookPARAM%albedoMax)%dat(1)
-    endif
-
-    ! make sure canopy water is positive
-    if(progData%var(iLookPROG%scalarCanopyWat)%dat(1) < 0.0001_rkind)then
-      progData%var(iLookPROG%scalarCanopyliq)%dat(1) = 0.0001_rkind
-    endif
-
-    ! initialize the spectral albedo
-    progData%var(iLookPROG%spectralSnowAlbedoDiffuse)%dat(1:nBand) = progData%var(iLookPROG%scalarSnowAlbedo)%dat(1)
-
-  end do ! end looping through prognostic variables (iVar)
-
-  ! --------------------------------------------------------------------------------------------------------
-  ! (2) set number of layers
-  ! --------------------------------------------------------------------------------------------------------
-
-  ! save the number of layers
-  indxData%var(iLookINDEX%nSnow)%dat(1)   = gru_struc(indxGRU)%hruInfo(indxHRU)%nSnow
-  indxData%var(iLookINDEX%nSoil)%dat(1)   = gru_struc(indxGRU)%hruInfo(indxHRU)%nSoil
-  indxData%var(iLookINDEX%nLayers)%dat(1) = gru_struc(indxGRU)%hruInfo(indxHRU)%nSnow + gru_struc(indxGRU)%hruInfo(indxHRU)%nSoil
-
-  ! set layer type
-  indxData%var(iLookINDEX%layerType)%dat(1:gru_struc(indxGRU)%hruInfo(indxHRU)%nSnow) = iname_snow
-  indxData%var(iLookINDEX%layerType)%dat((gru_struc(indxGRU)%hruInfo(indxHRU)%nSnow+1):(gru_struc(indxGRU)%hruInfo(indxHRU)%nSnow+gru_struc(indxGRU)%hruInfo(indxHRU)%nSoil)) = iname_soil
-
-
-  ! --------------------------------------------------------------------------------------------------------
-  ! (3) update soil layers (diagnostic variables)
-  ! --------------------------------------------------------------------------------------------------------
-  ! loop through soil layers
-  do iLayer = 1,indxData%var(iLookINDEX%nSoil)%dat(1)
-
-    ! get layer in the total vector
-    jLayer = iLayer+indxData%var(iLookINDEX%nSnow)%dat(1)
-
-    ! update soil layers
-    call updateSoil(&
-                    ! input
-                    progData%var(iLookPROG%mLayerTemp          )%dat(jLayer),& ! intent(in): temperature vector (K)
-                    progData%var(iLookPROG%mLayerMatricHead    )%dat(iLayer),& ! intent(in): matric head (m)
-                    mparData%var(iLookPARAM%vGn_alpha          )%dat(iLayer),& ! intent(in): van Genutchen "alpha" parameter
-                    mparData%var(iLookPARAM%vGn_n              )%dat(iLayer),& ! intent(in): van Genutchen "n" parameter
-                    mparData%var(iLookPARAM%theta_sat          )%dat(iLayer),& ! intent(in): soil porosity (-)
-                    mparData%var(iLookPARAM%theta_res          )%dat(iLayer),& ! intent(in): soil residual volumetric water content (-)
-                    1._dp - 1._dp/mparData%var(iLookPARAM%vGn_n)%dat(iLayer),& ! intent(in): van Genutchen "m" parameter (-)
-                    ! output
-                    progData%var(iLookPROG%mLayerVolFracWat    )%dat(jLayer),& ! intent(out): volumetric fraction of total water (-)
-                    progData%var(iLookPROG%mLayerVolFracLiq    )%dat(jLayer),& ! intent(out): volumetric fraction of liquid water (-)
-                    progData%var(iLookPROG%mLayerVolFracIce    )%dat(jLayer),& ! intent(out): volumetric fraction of ice (-)
-                    err,message)                                                                   ! intent(out): error control
-    if (err/=0) then; message=trim(message)//trim(cmessage); return; end if
-
-  end do  ! looping through soil layers
-
-  ! --------------------------------------------------------------------------------------------------------
-  ! (2) now get the basin variable(s)
-  ! --------------------------------------------------------------------------------------------------------
-
-  ! get the index in the file: single HRU
-  if(allocated(init_cond_bvar))then
-
-    ! loop through specific basin variables (currently 1 but loop provided to enable inclusion of others)
-    ndx = (/iLookBVAR%routingRunoffFuture/)   ! array of desired variable indices
-    do i = 1,size(ndx)
-      iVar = ndx(i)
-
-      ! store data in basin var (bvar) structure
-
-      ! put the data into data structures
-      bvarData%var(iVar)%dat(1:nTDH) = init_cond_bvar(i)%var_data(indxGRU,1:nTDH)
-      ! check whether the first values is set to nf90_fill_double
-      if(any(abs(bvarData%var(iVar)%dat(1:nTDH) - nf90_fill_double) < epsilon(init_cond_bvar(i)%var_data)))then; err=20; endif
-      if(err==20)then; message=trim(message)//"data set to the fill value (name='"//trim(bvar_meta(iVar)%varName)//"')"; return; endif
-
-    end do ! end looping through basin variables
-  endif  ! end if case for not being a singleHRU run - gaurded by the structure not being allocated
-
-
-end subroutine read_icond
-
-end module read_icond_actors_module
-- 
GitLab