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 fd70be043a60cbcca4c2b9949cdee74d288f7a51..32f7308c33b09d97ce8d304fcff67dfa734539aa 100644
--- a/build/source/actors/file_access_actor/fortran_code/writeOutputFromOutputStructure.f90
+++ b/build/source/actors/file_access_actor/fortran_code/writeOutputFromOutputStructure.f90
@@ -125,20 +125,24 @@ subroutine writeOutput(handle_ncid, num_steps, start_gru, max_gru, err) bind(C,
   do iGRU=start_gru, max_gru
     stepCounter(:) = outputTimeStep(iGRU)%dat(:) ! We want to avoid updating outputTimeStep
     do iStep=1, num_steps
-      call writeBasin(ncid,iGRU,stepCounter(:),iStep,bvar_meta, &
-              outputStructure(1)%bvarStat(1)%gru(iGRU)%hru(indxHRU)%var, &
-              outputStructure(1)%bvarStruct(1)%gru(iGRU)%hru(indxHRU)%var, bvarChild_map, err, cmessage)
-      
       call writeTime(ncid,outputTimeStep(iGRU)%dat(:),iStep,time_meta, &
               outputStructure(1)%timeStruct(1)%gru(iGRU)%hru(indxHRU)%var,err,cmessage)
     end do ! istep
   end do ! iGRU
 
+
+  numGRU = max_gru-start_gru + 1
   ! ****************************************************************************
-  ! *** write data
+  ! *** write basin data
   ! ****************************************************************************
-  numGRU = max_gru-start_gru + 1 
+  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), &
+                  bvarChild_map,err,cmessage)
 
+  ! ****************************************************************************
+  ! *** write data
+  ! ****************************************************************************
   do iStruct=1,size(structInfo)
     select case(trim(structInfo(iStruct)%structName))
       case('forc')
@@ -170,6 +174,7 @@ subroutine writeOutput(handle_ncid, num_steps, start_gru, max_gru, err) bind(C,
     if(err/=0)then; message=trim(message)//trim(cmessage)//'['//trim(structInfo(iStruct)%structName)//']'; return; endif
   end do  ! (looping through structures)
 
+
   do iFreq = 1,maxvarFreq
     outputTimeStep(start_gru)%dat(iFreq) = outputTimeStep(start_gru)%dat(iFreq) + outputTimeStepUpdate(iFreq) 
   end do ! ifreq
@@ -300,6 +305,7 @@ subroutine writeData(ncid,outputTimestep,outputTimestepUpdate,maxLayers,nSteps,
     do iVar = 1,size(meta)
       stepCounter = 0
 
+      ! Todo: wrap this in a function
       if (meta(iVar)%varName=='time' .and. structName == 'forc')then
         ! get variable index
         err = nf90_inq_varid(ncid%var(iFreq),trim(meta(iVar)%varName),ncVarID)
@@ -394,6 +400,7 @@ subroutine writeScalar(ncid, outputTimestep, outputTimestepUpdate, nSteps, minGR
           if(.not.outputStructure(1)%finalizeStats(1)%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
         end do ! iStep
       end do ! iGRU  
     
@@ -401,7 +408,9 @@ subroutine writeScalar(ncid, outputTimestep, outputTimestepUpdate, nSteps, minGR
       if (outputTimeStepUpdate(iFreq) /= stepCounter ) then
         print*, "ERROR Missmatch in Steps - stat doubleVec"
         print*, "   outputTimeStepUpdate(iFreq) = ", outputTimeStepUpdate(iFreq)
+        print*, "   outputTimeStepUpdate", outputTimeStepUpdate
         print*, "   stepCounter = ", stepCounter
+        print*, "   iFreq = ", iFreq
         return
       endif
     class default; err=20; message=trim(message)//'stats must be scalarv and of type gru_hru_doubleVec'; return
@@ -530,7 +539,9 @@ end subroutine
 ! **************************************************************************************
 ! public subroutine writeBasin: write basin-average variables
 ! **************************************************************************************
-subroutine writeBasin(ncid,iGRU,outputTimestep,iStep,meta,stat,dat,map,err,message)
+subroutine writeBasin(ncid,outputTimestep,outputTimestepUpdate,nSteps,&
+                      minGRU, maxGRU, numGRU, &
+                      meta,stat,dat,map,err,message)
   USE data_types,only:var_info                       ! metadata type
   USE var_lookup,only:maxVarStat                     ! index into stats structure
   USE var_lookup,only:iLookVarType                   ! index into type structure
@@ -541,12 +552,15 @@ subroutine writeBasin(ncid,iGRU,outputTimestep,iStep,meta,stat,dat,map,err,messa
 
   ! declare dummy variables
   type(var_i)   ,intent(in)     :: ncid              ! file ids
-  integer(i4b)  ,intent(in)     :: iGRU              ! GRU index
   integer(i4b)  ,intent(inout)  :: outputTimestep(:) ! output time step
-  integer(i4b)  ,intent(in)     :: iStep            ! number of steps in forcing file
+  integer(i4b)  ,intent(inout)  :: outputTimestepUpdate(:) ! number of HRUs in the run domain
+  integer(i4b)  ,intent(in)     :: nSteps            ! number of timeSteps
+  integer(i4b)  ,intent(in)     :: minGRU            ! minGRU index to write
+  integer(i4b)  ,intent(in)     :: maxGRU            ! maxGRU index to write - probably not needed
+  integer(i4b)  ,intent(in)     :: numGRU            ! number of GRUs to write
   type(var_info),intent(in)     :: meta(:)           ! meta data
-  type(time_dlength),intent(in) :: stat(:)           ! stats data
-  type(time_dlength),intent(in) :: dat(:)            ! timestep data
+  class(*)      ,intent(in)     :: stat              ! stats data
+  class(*)      ,intent(in)     :: dat               ! timestep data
   integer(i4b)  ,intent(in)     :: map(:)            ! map into stats child struct
   integer(i4b)  ,intent(out)    :: err               ! error code
   character(*)  ,intent(out)    :: message           ! error message
@@ -564,7 +578,7 @@ subroutine writeBasin(ncid,iGRU,outputTimestep,iStep,meta,stat,dat,map,err,messa
     if(.not.outFreq(iFreq)) cycle
 
     ! 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(1)%gru(1)%hru(1)%tim(iStep)%dat(iFreq)) cycle
 
     ! loop through model variables
     do iVar = 1,size(meta)
@@ -579,11 +593,14 @@ subroutine writeBasin(ncid,iGRU,outputTimestep,iStep,meta,stat,dat,map,err,messa
       select case (meta(iVar)%varType)
 
         case (iLookVarType%scalarv)
-          err = nf90_put_var(ncid%var(iFreq),meta(iVar)%ncVarID(iFreq),(/stat(map(iVar))%tim(iStep)%dat(iFreq)/),start=(/iGRU,outputTimestep(iFreq)/),count=(/1,1/))
+          call writeScalar(ncid, outputTimeStep, outputTimeStepUpdate, nSteps, &
+                           minGRU, maxGRU, numGRU, iFreq, iVar, meta, stat, map, &
+                           err, message)
+          ! err = nf90_put_var(ncid%var(iFreq),meta(iVar)%ncVarID(iFreq),(/stat(map(iVar))%tim(iStep)%dat(iFreq)/),start=(/iGRU,outputTimestep(iFreq)/),count=(/1,1/))
 
         case (iLookVarType%routing)
           if (iFreq==1 .and. outputTimestep(iFreq)==1) then
-            err = nf90_put_var(ncid%var(iFreq),meta(iVar)%ncVarID(iFreq),(/dat(iVar)%tim(iStep)%dat/),start=(/1/),count=(/1000/))
+            ! err = nf90_put_var(ncid%var(iFreq),meta(iVar)%ncVarID(iFreq),(/dat(iVar)%tim(iStep)%dat/),start=(/1/),count=(/1000/))
           end if
         case default
           err=40; message=trim(message)//"unknownVariableType[name='"//trim(meta(iVar)%varName)//"';type='"//trim(get_varTypeName(meta(iVar)%varType))//    "']"; return
@@ -593,7 +610,7 @@ subroutine writeBasin(ncid,iGRU,outputTimestep,iStep,meta,stat,dat,map,err,messa
       if (err.ne.0) message=trim(message)//trim(meta(iVar)%varName)//'_'//trim(get_statName(iStat))
       call netcdf_err(err,message); if (err/=0) return
     end do ! iVar
-    outputTimeStep(iFreq) = outputTimeStep(iFreq) + 1
+    ! outputTimeStep(iFreq) = outputTimeStep(iFreq) + 1
   end do ! iFreq
 
 end subroutine writeBasin
@@ -602,48 +619,48 @@ end subroutine writeBasin
     ! public subroutine writeTime: write current time to all files
     ! **************************************************************************************
 subroutine writeTime(ncid,outputTimestep,iStep,meta,dat,err,message)
-USE data_types,only:var_info                       ! metadata type
-USE var_lookup,only:iLookStat                      ! index into stat structure
-implicit none
+  USE data_types,only:var_info                       ! metadata type
+  USE var_lookup,only:iLookStat                      ! index into stat structure
+  implicit none
 
-! declare dummy variables
-type(var_i)   ,intent(in)     :: ncid              ! file ids
-integer(i4b)  ,intent(inout)  :: outputTimestep(:) ! output time step
-integer(i4b)  ,intent(in)     :: iStep
-type(var_info),intent(in)     :: meta(:)           ! meta data
-type(time_i)  ,intent(in)     :: dat(:)            ! timestep data
-integer(i4b)  ,intent(out)    :: err               ! error code
-character(*)  ,intent(out)    :: message           ! error message
-! local variables
-integer(i4b)                  :: iVar              ! variable index
-integer(i4b)                  :: iFreq             ! frequency index
-integer(i4b)                  :: ncVarID           ! used only for time
-! initialize error control
-err=0;message="f-writeTime/"
-! loop through output frequencies
-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
-
-  ! loop through model variables
-  do iVar = 1,size(meta)
+  ! declare dummy variables
+  type(var_i)   ,intent(in)     :: ncid              ! file ids
+  integer(i4b)  ,intent(inout)  :: outputTimestep(:) ! output time step
+  integer(i4b)  ,intent(in)     :: iStep
+  type(var_info),intent(in)     :: meta(:)           ! meta data
+  type(time_i)  ,intent(in)     :: dat(:)            ! timestep data
+  integer(i4b)  ,intent(out)    :: err               ! error code
+  character(*)  ,intent(out)    :: message           ! error message
+  ! local variables
+  integer(i4b)                  :: iVar              ! variable index
+  integer(i4b)                  :: iFreq             ! frequency index
+  integer(i4b)                  :: ncVarID           ! used only for time
+  ! initialize error control
+  err=0;message="f-writeTime/"
+  ! loop through output frequencies
+  do iFreq=1,maxvarFreq
 
-    ! check instantaneous
-  if (meta(iVar)%statIndex(iFreq)/=iLookStat%inst) cycle
-    print*, "Time Data", dat(iVar)%tim(iStep)
-      ! get variable id in file
-    err = nf90_inq_varid(ncid%var(iFreq),trim(meta(iVar)%varName),ncVarID)
-    if (err/=0) message=trim(message)//trim(meta(iVar)%varName); call netcdf_err(err,message)
-    if (err/=0) then; err=20; return; end if
-
-    ! add to file
-    err = nf90_put_var(ncid%var(iFreq),ncVarID,(/dat(iVar)%tim(iStep)/),start=(/outputTimestep(iFreq)/),count=(/1/))
-    if (err/=0) message=trim(message)//trim(meta(iVar)%varName);call netcdf_err(err,message)
-    if (err/=0) then; err=20; return; end if
-
-  end do ! iVar
-end do ! iFreq
+    ! 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
+
+    ! loop through model variables
+    do iVar = 1,size(meta)
+
+      ! check instantaneous
+    if (meta(iVar)%statIndex(iFreq)/=iLookStat%inst) cycle
+      print*, "Time Data", dat(iVar)%tim(iStep)
+        ! get variable id in file
+      err = nf90_inq_varid(ncid%var(iFreq),trim(meta(iVar)%varName),ncVarID)
+      if (err/=0) message=trim(message)//trim(meta(iVar)%varName); call netcdf_err(err,message)
+      if (err/=0) then; err=20; return; end if
+
+      ! add to file
+      err = nf90_put_var(ncid%var(iFreq),ncVarID,(/dat(iVar)%tim(iStep)/),start=(/outputTimestep(iFreq)/),count=(/1/))
+      if (err/=0) message=trim(message)//trim(meta(iVar)%varName);call netcdf_err(err,message)
+      if (err/=0) then; err=20; return; end if
+
+    end do ! iVar
+  end do ! iFreq
 
 
 end subroutine writeTime
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 7ab98ae61c9d4fbef195ee49e7e83bf43d81713d..490d7120c859d43fc1667075b997c6a23fba2a59 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
@@ -3,47 +3,15 @@ USE, intrinsic :: iso_c_binding
 USE nrtype
 USE data_types
 
+!TODO: This module is only used for writeParamToNetCDF the others are not.
+
+
 implicit none
 public::writeParamToNetCDF
 public::writeDataToNetCDF
 public::writeBasinToNetCDF
 public::writeTimeToNetCDF
 
-! type var_dlength_ptr
-!   type(var_dlength), pointer
-! type(var_dlength)     forc_stat_all_grus
-! type(var_dlength)     prog_stat_all_grus
-! type(var_dlength)     diag_stat_all_grus
-! type(var_dlength)     flux_stat_all_grus
-! type(var_dlength)     indx_stat_all_grus
-! type(var_dlength)     bvar_stat_all_grus
-
-! type(var_i)           time_struct
-! type(var_d)           forc_struct
-! type(var_d)           attr_struct
-! type(var_i)           type_struct
-! type(var_i8)          id_struct
-
-! type(var_ilength)     indx_struct
-! type(var_dlength)     mpar_struct
-! type(var_dlength)     prog_struct
-! type(var_dlength)     diag_struct
-! type(var_dlength)     flux_struct
-
-! type(var_d)           bpar_struct
-! type(var_dlength)     bvar_struct
-
-! type(var_d)           dpar_struct
-! type(flagVec)         finalize_stats
-! type(var_i)           output_timesteps_all_grus
-
-
-
-! type(flagVec),     pointer, dimension(:,:), allocatable, public :: finalize_stats_all_grus
-! type(var_i),       pointer, dimension(:,:), allocatable, public :: output_timesteps_all_grus
-! type(var_dlength), pointer, dimension(:,:), allocatable, public :: bvar_stats_all_grus
-! type(var_dlength), pointer, dimension(:,:), allocatable, public :: bvar_struct_all_grus
-! type(var_i),       pointer, dimension(:,:), allocatable, public :: time_struct_all_grus
 
 contains
 
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 b6a82d66ddb6fb9e930cc4500222364c1e43a964..2485e6a09ff7127acfab13a543b73f5f390319b2 100644
--- a/build/source/actors/hru_actor/fortran_code/hru_writeOutput.f90
+++ b/build/source/actors/hru_actor/fortran_code/hru_writeOutput.f90
@@ -48,6 +48,7 @@ USE var_lookup,only:iLookDIAG                 ! named variables for local column
 USE var_lookup,only:iLookPROG                 ! named variables for local column model prognostic variables
 USE var_lookup,only:iLookINDEX                ! named variables for local column index variables
 USE var_lookup,only:iLookFreq                 ! named variables for the frequency structure
+USE var_lookup,only:iLookBVAR                 ! named variables for basin parameters
 USE get_ixname_module,only:get_freqName       ! get name of frequency from frequency index
 
 
@@ -265,6 +266,12 @@ subroutine writeHRUToOutputStructure(&
 
  ! 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(:)
+  print*, ""
+  print*, "indxGRU" , indxGRU
+  print*, "averageRoutedRunoff", bvarStruct%var(iLookBVAR%averageRoutedRunoff)%dat(1)
+  print*, "outputStep", outputStep
+  print*, ""
+
  ! ****************************************************************************
  ! *** calculate output statistics
  ! ****************************************************************************
diff --git a/build/source/actors/hru_actor/fortran_code/outputStrucWrite.f90 b/build/source/actors/hru_actor/fortran_code/outputStrucWrite.f90
index 3b4501e0bfbade74cac318d63a83e552eec08fa8..21707dc5db0dbfba4082ab435d3b3e147f84164d 100755
--- a/build/source/actors/hru_actor/fortran_code/outputStrucWrite.f90
+++ b/build/source/actors/hru_actor/fortran_code/outputStrucWrite.f90
@@ -77,10 +77,10 @@ public::writeRestart
 integer(i4b),parameter      :: maxSpectral=2              ! maximum number of spectral bands
 contains
 
- ! **********************************************************************************************************
- ! public subroutine writeParm: write model parameters
- ! **********************************************************************************************************
-  subroutine writeParm(indxGRU,indxHRU,ispatial,struct,meta,structName,err,message)
+! **********************************************************************************************************
+! public subroutine writeParm: write model parameters
+! **********************************************************************************************************
+subroutine writeParm(indxGRU,indxHRU,ispatial,struct,meta,structName,err,message)
   ! USE globalData,only:ncid                      ! netcdf file ids
   USE data_types,only:var_info                    ! metadata info
   USE var_lookup,only:iLookStat                   ! index in statistics vector
@@ -363,10 +363,12 @@ subroutine writeBasin(indxGRU,indxHRU,iStep,finalizeStats,&
 
     case (iLookVarType%scalarv)
       outputStructure(1)%bvarStat(1)%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%dat(iFreq) = stat(map(iVar))%dat(iFreq)
-
+      print*, "iFreq" , iFreq
+      print*, "outputStructure(1)%bvarStat(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(iFreq)", outputStructure(1)%bvarStat(1)%gru(indxGRU)%hru(indxHRU)%var(map(iVar))%tim(iStep)%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)
+      print*, "outputStructure(1)%bvarStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(iFreq)", outputStructure(1)%bvarStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar)%tim(iStep)%dat(iFreq)
      end if
 
     case default
diff --git a/utils/netcdf/OutputVerification/checkOutput.py b/utils/netcdf/OutputVerification/checkOutput.py
index 5079637703395fc899b0f3a15da75303a1222797..5a85051d278baf2edaf4262e7906be96bf47b91d 100644
--- a/utils/netcdf/OutputVerification/checkOutput.py
+++ b/utils/netcdf/OutputVerification/checkOutput.py
@@ -80,15 +80,15 @@ def get_output_vars(model_output_file):
 
 
 
-num_hru = 10
+num_hru = 2
 print("Checking output for", num_hru, "HRUs")
 # dataset_1 = "/scratch/kck540/Summa_Sundials/non-actors/SummaOriginal_G000001-000010_timestep.nc"
 # dataset_2 = "/scratch/kck540/Summa_Sundials/actors/SummaActorsGRU1-10_timestep.nc"
 
-dataset_1 = "/scratch/kck540/Summa_Sundials/non-actors/SummaOriginal-BE_G000001-000010_timestep.nc"
-dataset_2 = "/scratch/kck540/Summa_Sundials/actors/SummaActors-BEGRU1-10_timestep.nc"
+dataset_1 = "/scratch/kck540/Summa_Sundials/non-actors/SummaOriginal-BE_G000001-000002_timestep.nc"
+dataset_2 = "/scratch/kck540/Summa_Sundials/actors/SummaActors-BEGRU1-2_timestep.nc"
 
-model_output_file = "/gladwell/kck540/Sundials_Settings/outputControl.txt"
+model_output_file = "/gladwell/kck540/Summa_Distributed_Settings/Summa_Settings/outputControl.txt"
 
 output_vars = get_output_vars(model_output_file)
 verify_data(dataset_1, dataset_2, num_hru, output_vars)