From ed42ce9348a603eb230701de88f610af25e9caa7 Mon Sep 17 00:00:00 2001
From: Kyle <kyle.c.klenk@gmail.com>
Date: Wed, 26 Oct 2022 23:22:57 +0000
Subject: [PATCH] got fortran routines for file_access_actor to assemble
parameters from nc file
---
.../file_access_actor/file_access_actor.hpp | 12 +-
.../file_access_actor_subroutine_wrappers.hpp | 19 +-
build/includes/global/message_atoms.hpp | 1 +
build/makefile_sundials | 2 +-
.../cpp_code/file_access_actor.cpp | 69 ++--
.../fortran_code/initOutputStruc.f90 | 10 -
.../fortran_code/read_attribute.f90 | 10 +-
.../fortran_code/read_param.f90 | 355 ++++++++++++++++++
.../fortran_code/read_param_all_hru.f90 | 26 +-
build/source/actors/hru_actor/hru_actor.cpp | 63 ++--
build/source/engine/pOverwrite.f90 | 120 +++---
11 files changed, 534 insertions(+), 153 deletions(-)
create mode 100644 build/source/actors/file_access_actor/fortran_code/read_param.f90
diff --git a/build/includes/file_access_actor/file_access_actor.hpp b/build/includes/file_access_actor/file_access_actor.hpp
index cb34d17..5deecd1 100644
--- a/build/includes/file_access_actor/file_access_actor.hpp
+++ b/build/includes/file_access_actor/file_access_actor.hpp
@@ -34,9 +34,13 @@ struct file_access_state {
std::vector<std::vector<int>> type_arrays_for_hrus;
std::vector<std::vector<long int>> id_arrays_for_hrus;
- // vector of handles for parameters
- std::vector<void *> mpar_struct_handles;
- std::vector<void *> bpar_struct_handles;
+ // Variables for handling parameters file
+ std::vector<void*> mpar_structs_for_hrus;
+ std::vector<std::vector<double>> bpar_arrays_for_hrus;
+ std::vector<std::vector<double>> dpar_arrays_for_hrus;
+ int dpar_array_size;
+ int bpar_array_size;
+ bool param_file_exists;
// Timing Variables
@@ -53,5 +57,7 @@ int write(stateful_actor<file_access_state>* self, int listIndex);
// Read in the attributes for all HRUs that are in the run-domain
void readAttributes(stateful_actor<file_access_state>* self);
+// read in the parameters for all HRUs that are in the run-domain
+void readParameters(stateful_actor<file_access_state>* self);
} // end namespace
\ No newline at end of file
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 c6a8d4e..1df3e82 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
@@ -30,11 +30,11 @@ extern "C" {
void def_output(void* handle_ncid, int* startGRU, int* numGRU, int* numHRU, int* err);
- void Write_HRU_Param(void* handle_ncid, int* indxGRU, int* indxHRU, int* err);
+ // void Write_HRU_Param(void* handle_ncid, int* indxGRU, int* indxHRU, int* err);
// void readAttributeFileAccessActor(int* num_gru, int* err);
- void overwriteParam(int* num_gru, int* err);
+ // void overwriteParam(int* num_gru, int* err);
// void readParamFileAccessActor(int* start_gru, int* num_gru, int* err);
@@ -60,12 +60,25 @@ extern "C" {
// Attributes Files
void openAttributeFile(int* att_ncid, int* err);
- void getNumVar(int* attr_ncid, int* num_var_attr, int* err);
+ void getNumVarAttr(int* attr_ncid, int* num_var_attr, int* err);
void closeAttributeFile(int* attr_ncid, int* err);
void readAttributeFromNetCDF(int* attr_ncid, int* index_gru, int* index_hru, int* num_var,
void* attr_array, void* type_array, void* id_array, int* err);
+
+ // Parameters File
+ void openParamFile(int* param_ncid, bool* param_file_exists, int* err);
+ void getNumVarParam(int* param_ncid, int* num_var_param, int* err);
+ void closeParamFile(int* param_ncid, int* err);
+ void getParamSizes(int* dpar_array_size, int* bpar_array_size);
+
+ void overwriteParam(int* index_gru, int* index_hru, int* num_var_attr,
+ void* type_array, void* dpar_array, void* handle_mpar_struct, void* bpar_array,
+ int* err);
+
+ void readParamFromNetCDF(int* param_ncid, int* index_gru, int* index_hru, int* start_index_gru,
+ int* num_var_param, int* bpar_array_size, void* handle_mpar_struct, void* bpar_array, int* err);
}
diff --git a/build/includes/global/message_atoms.hpp b/build/includes/global/message_atoms.hpp
index e820fad..d9c4d1d 100644
--- a/build/includes/global/message_atoms.hpp
+++ b/build/includes/global/message_atoms.hpp
@@ -39,6 +39,7 @@ CAF_BEGIN_TYPE_ID_BLOCK(summa, first_custom_type_id)
CAF_ADD_ATOM(summa, write_param)
CAF_ADD_ATOM(summa, restart_failures)
CAF_ADD_ATOM(summa, serialized_hru_data)
+ CAF_ADD_ATOM(summa, get_attributes)
// HRU Actor
CAF_ADD_ATOM(summa, run_hru)
CAF_ADD_ATOM(summa, start_hru)
diff --git a/build/makefile_sundials b/build/makefile_sundials
index 8d544da..60f38af 100644
--- a/build/makefile_sundials
+++ b/build/makefile_sundials
@@ -165,7 +165,7 @@ SUMMA_FILEACCESS_INTERFACE = \
deallocateOutputStruc.f90 \
cppwrap_fileAccess.f90 \
read_attribute.f90 \
- read_param_all_hru.f90 \
+ read_param.f90 \
write_to_netcdf.f90
FILEACCESS_INTERFACE = $(patsubst %, $(FILE_ACCESS_DIR)/%, $(SUMMA_FILEACCESS_INTERFACE))
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 c0696ba..e8bfca7 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
@@ -128,6 +128,17 @@ behavior file_access_actor(stateful_actor<file_access_state>* self, int start_gr
}
},
+ [=] (get_attributes, int ref_gru, caf::actor actor_to_respond) {
+ // ref_gru will always be 1 index too high (FORTRAN arrays start at 1)
+ std::vector<double> attr_array_to_send = self->state.attr_arrays_for_hrus[ref_gru-1];
+ std::vector<int> type_array_to_send = self->state.type_arrays_for_hrus[ref_gru-1];
+ std::vector<long int> id_array_to_send = self->state.id_arrays_for_hrus[ref_gru-1];
+
+ self->send(actor_to_respond, get_attributes_v, attr_array_to_send,
+ type_array_to_send, id_array_to_send);
+
+ },
+
[=](write_output, int index_gru, int index_hru,
// statistic structures
std::vector<std::vector<double>> forc_stat, std::vector<std::vector<double>> prog_stat, std::vector<std::vector<double>> diag_stat,
@@ -326,28 +337,22 @@ void initalizeFileAccessActor(stateful_actor<file_access_state>* self) {
Init_OutputStruct(self->state.handle_forcing_file_info, &self->state.outputStrucSize,
&self->state.num_gru, &self->state.err);
- // // Read In all of the attribres for the number of GRUs in the run Domian
- // readAttributeFileAccessActor(&self->state.num_gru, &err);
+ // Read in the attribute information for the HRUs to request
+ readAttributes(self);
+
+ readParameters(self);
+
+ // Noah-MP table information
+ // overwriteParam(&self->state.num_gru, &err);
// if (err != 0) {
- // aout(self) << "ERROR: FILE_ACCESS_ACTOR readAttributeFilAccessActor() \n";
- // std::string function = "readAttributeFileAccessActor";
+ // aout(self) << "ERROR: FILE_ACCESS_ACTOR overwriteParam() \n";
+ // std::string function = "overwriteParam";
// self->send(self->state.parent, file_access_actor_err_v, function);
// self->quit();
// return;
// }
- // Noah-MP table information
- overwriteParam(&self->state.num_gru, &err);
- if (err != 0) {
- aout(self) << "ERROR: FILE_ACCESS_ACTOR overwriteParam() \n";
- std::string function = "overwriteParam";
- self->send(self->state.parent, file_access_actor_err_v, function);
- self->quit();
- return;
- }
-
-
- // // Read in all of the parmeters for the number of GRUs in the run Domain
+ // Read in all of the parmeters for the number of GRUs in the run Domain
// readParamFileAccessActor(&self->state.start_gru, &self->state.num_gru, &err);
// if (err != 0) {
// aout(self) << "ERROR: FILE_ACCESS_ACTOR readParamFileAccessActor() \n";
@@ -405,9 +410,8 @@ int readForcing(stateful_actor<file_access_state>* self, int currentFile) {
void readAttributes(stateful_actor<file_access_state>* self) {
int err = 0;
openAttributeFile(&self->state.attribute_ncid, &err);
- getNumVar(&self->state.attribute_ncid, &self->state.num_var_in_attributes_file, &err);
- for (int index_gru = self->state.start_gru;
- index_gru < self->state.num_gru + self->state.start_gru; index_gru++) {
+ getNumVarAttr(&self->state.attribute_ncid, &self->state.num_var_in_attributes_file, &err);
+ for (int index_gru = 1; index_gru < self->state.num_gru + 1; index_gru++) {
std::vector<double> attr_array(self->state.num_var_in_attributes_file);
std::vector<int> type_array(self->state.num_var_in_attributes_file);
@@ -427,8 +431,31 @@ void readAttributes(stateful_actor<file_access_state>* self) {
}
-// void readParameters(stateful_actor<file_access_state>* self) {
+void readParameters(stateful_actor<file_access_state>* self) {
+ int err = 0;
+ int index_hru = 1;
+ getParamSizes(&self->state.dpar_array_size, &self->state.bpar_array_size);
+ // openParamFile();
+ for (int index_gru = 1; index_gru < self->state.num_gru + 1; index_gru++) {
+
+ std::vector<double> dpar_array(self->state.dpar_array_size);
+ void* handle_mpar_struct = new_handle_var_dlength();
+ std::vector<double> bpar_array(self->state.dpar_array_size);
+
+
+ overwriteParam(&index_gru, &index_hru, &self->state.num_var_in_attributes_file,
+ &self->state.attr_arrays_for_hrus[index_gru-1][0],
+ &dpar_array[0], handle_mpar_struct, &bpar_array[0], &err);
-// }
+
+ // readParamFromNetCDF();
+
+ self->state.dpar_arrays_for_hrus.push_back(dpar_array);
+ self->state.mpar_structs_for_hrus.push_back(handle_mpar_struct);
+ self->state.bpar_arrays_for_hrus.push_back(bpar_array);
+ }
+ // closeParamFile();
+
+}
} // end namespace
\ No newline at end of file
diff --git a/build/source/actors/file_access_actor/fortran_code/initOutputStruc.f90 b/build/source/actors/file_access_actor/fortran_code/initOutputStruc.f90
index f109f8d..9bbaa93 100644
--- a/build/source/actors/file_access_actor/fortran_code/initOutputStruc.f90
+++ b/build/source/actors/file_access_actor/fortran_code/initOutputStruc.f90
@@ -47,16 +47,9 @@ subroutine initalizeOutput(forcFileInfo, maxSteps, num_gru, err)
! Primary Data Structures (scalars)
- allocate(outputStructure(1)%attrStruct(1))
- allocate(outputStructure(1)%typeStruct(1))
- allocate(outputStructure(1)%idStruct(1))
allocate(outputStructure(1)%mparStruct(1))
allocate(outputStructure(1)%bparStruct(1))
allocate(outputStructure(1)%dparStruct(1))
-
- 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)%mparStruct(1)%gru(num_gru))
allocate(outputStructure(1)%bparStruct(1)%gru(num_gru))
allocate(outputStructure(1)%dparStruct(1)%gru(num_gru))
@@ -67,9 +60,6 @@ subroutine initalizeOutput(forcFileInfo, maxSteps, num_gru, err)
num_hru = gru_struc(iGRU)%hruCount
! Primary Data Structures (scalars)
- 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)%mparStruct(1)%gru(iGRU)%hru(num_hru))
allocate(outputStructure(1)%dparStruct(1)%gru(iGRU)%hru(num_hru))
diff --git a/build/source/actors/file_access_actor/fortran_code/read_attribute.f90 b/build/source/actors/file_access_actor/fortran_code/read_attribute.f90
index 3fd01e0..7a0aba8 100644
--- a/build/source/actors/file_access_actor/fortran_code/read_attribute.f90
+++ b/build/source/actors/file_access_actor/fortran_code/read_attribute.f90
@@ -5,7 +5,7 @@ USE nrtype
implicit none
private
public::openAttributeFile
-public::getNumVar
+public::getNumVarAttr
public::closeAttributeFile
public::readAttributeFromNetCDF
@@ -14,8 +14,6 @@ contains
subroutine openAttributeFile(attr_ncid, err) bind(C, name="openAttributeFile")
USE netcdf
USE netcdf_util_module,only:nc_file_open ! open netcdf file
- USE netcdf_util_module,only:nc_file_close ! close netcdf file
- USE netcdf_util_module,only:netcdf_err ! netcdf error handling function
! Attribute File
USE summaActors_FileManager,only:SETTINGS_PATH ! define path to settings files (e.g., parameters, soil and veg. tables)
USE summaActors_FileManager,only:LOCAL_ATTRIBUTES ! name of model initial attributes file
@@ -39,12 +37,12 @@ subroutine openAttributeFile(attr_ncid, err) bind(C, name="openAttributeFile")
end subroutine
-subroutine getNumVar(attr_ncid, num_var, err) bind(C, name="getNumVar")
+subroutine getNumVarAttr(attr_ncid, num_var, err) bind(C, name="getNumVarAttr")
USE netcdf
USE netcdf_util_module,only:netcdf_err ! netcdf error handling function
implicit none
integer(c_int),intent(in) :: attr_ncid
- integer(c_int),intent(out) :: num_var
+ integer(c_int),intent(out) :: num_var
integer(c_int),intent(out) :: err
! local variables
@@ -58,7 +56,7 @@ subroutine getNumVar(attr_ncid, num_var, err) bind(C, name="getNumVar")
return
endif
-end subroutine getNumVar
+end subroutine getNumVarAttr
subroutine closeAttributeFile(attr_ncid, err) bind(C, name="closeAttributeFile")
USE netcdf_util_module,only:nc_file_close
diff --git a/build/source/actors/file_access_actor/fortran_code/read_param.f90 b/build/source/actors/file_access_actor/fortran_code/read_param.f90
new file mode 100644
index 0000000..a56efe3
--- /dev/null
+++ b/build/source/actors/file_access_actor/fortran_code/read_param.f90
@@ -0,0 +1,355 @@
+module read_param_module
+ USE, intrinsic :: iso_c_binding
+ USE nrtype
+ implicit none
+ private
+ public::openParamFile
+ public::getNumVarParam
+ public::closeParamFile
+ public::getParamSizes
+ public::overwriteParam
+ public::readParamFromNetCDF
+ contains
+
+subroutine openParamFile(param_ncid, param_file_exists, err) bind(C, name="openParamFile")
+ USE netcdf
+ USE netcdf_util_module,only:nc_file_open
+ ! Parameters File
+ USE summaActors_FileManager,only:SETTINGS_PATH ! path for metadata files
+ USE summaActors_FileManager,only:PARAMETER_TRIAL ! file with parameter trial values
+
+ implicit none
+ integer(c_int),intent(out) :: param_ncid
+ logical(c_bool),intent(out) :: param_file_exists
+ integer(c_int),intent(out) :: err
+
+ ! local variables
+ character(LEN=1024) :: infile ! input filename
+ character(len=256) :: message ! error message
+ character(len=1024) :: cmessage ! error message for downwind routine
+ err=0; message="read_param.f90 - openParamFile"
+ ! **********************************************************************************************
+ ! * open files, etc.
+ ! **********************************************************************************************
+
+ infile = trim(SETTINGS_PATH)//trim(PARAMETER_TRIAL) ! build filename
+
+ ! check whether the user-specified file exists and warn if it does not
+ inquire(file=trim(infile),exist=param_file_exists)
+ if (.not.param_file_exists) then
+ write(*,'(A)') NEW_LINE('A')//'!! WARNING: trial parameter file not found; proceeding instead with other default parameters; check path in file manager input if this was not the desired behavior'//NEW_LINE('A')
+ return
+ endif
+
+ ! open trial parameters file if it exists
+ call nc_file_open(trim(infile),nf90_nowrite,param_ncid,err,cmessage)
+ if(err/=0)then
+ message=trim(message)//trim(cmessage)
+ print*, message
+ return
+ end if
+
+end subroutine openParamFile
+
+subroutine getNumVarParam(param_ncid, num_var, err) bind(C, name="getNumVarParam")
+ USE netcdf
+ USE netcdf_util_module,only:netcdf_err ! netcdf error handling function
+ implicit none
+ integer(c_int),intent(in) :: param_ncid
+ integer(c_int),intent(out) :: num_var
+ integer(c_int),intent(out) :: err
+
+ ! local variables
+ character(len=256) :: message ! error message
+
+ err=0; message="read_param.f90 - getNumVarAndDims/"
+ ! get the number of variables in the parameter file
+ err=nf90_inquire(param_ncid, nVariables=num_var)
+ call netcdf_err(err,message)
+ if (err/=0) then
+ err=20
+ print*, message
+ return
+ end if
+
+end subroutine getNumVarParam
+
+subroutine closeParamFile(param_ncid, err) bind(C, name="closeParamFile")
+ USE netcdf_util_module,only:nc_file_close
+ implicit none
+ integer(c_int),intent(in) :: param_ncid
+ integer(c_int),intent(out) :: err
+ ! local variables
+ character(len=256) :: message
+
+ err=0; message="read_param.f90 - closeParamFile/"
+
+ call nc_file_close(param_ncid,err,message)
+ if(err/=0)then
+ message=trim(message)
+ print*, message
+ return
+ end if
+
+end subroutine closeParamFile
+
+! get the sizes of the arrays for dpar_array bpar_array
+subroutine getParamSizes(dpar_array_size, bpar_array_size) bind(C, name="getParamSizes")
+ USE var_lookup,only:maxvarMpar ! model parameters: maximum number variables
+ USE var_lookup,only:maxvarBpar ! model parameters: maximum number variables
+
+ implicit none
+ integer(c_int),intent(out) :: dpar_array_size
+ integer(c_int),intent(out) :: bpar_array_size
+
+
+ dpar_array_size = maxvarMpar
+ bpar_array_size = maxvarBpar
+
+
+end subroutine getParamSizes
+
+subroutine overwriteParam(index_gru, index_hru, num_var_attr, type_array, &
+ dpar_array, handle_mpar_struct, bpar_array, err) bind(C, name="overwriteParam")
+ USE var_lookup,only:maxvarMpar ! model parameters: maximum number variables
+ USE var_lookup,only:maxvarBpar ! model parameters: maximum number variables
+ USE var_lookup,only:iLookTYPE ! named variables to index elements of the data vectors
+
+ ! global data
+ USE globalData,only:gru_struc
+ USE globalData,only:localParFallback ! local column default parameters
+ USE globalData,only:basinParFallback ! basin-average default parameter
+ USE globalData,only:mpar_meta
+ USE data_types,only:var_dlength
+
+ USE pOverwrite_module,only:pOverwrite ! module to overwrite default parameter values with info from the Noah tables
+ USE allocspace_module,only:allocLocal
+
+
+ implicit none
+ integer(c_int),intent(in) :: index_gru
+ integer(c_int),intent(in) :: index_hru
+ integer(c_int),intent(in) :: num_var_attr ! size of type array
+ ! arrays
+ integer(c_int),intent(in) :: type_array(num_var_attr)
+ real(c_double),intent(out) :: dpar_array(maxvarMpar)
+ type(c_ptr),intent(in),value :: handle_mpar_struct
+ real(c_double),intent(out) :: bpar_array(maxvarBpar)
+ ! error control
+ integer(c_int), intent(out) :: err
+
+ ! local variables
+ type(var_dlength),pointer :: mpar_struct ! model parameters
+ integer(i4b) :: iVar
+ integer(i4b) :: iDat
+
+ character(len=256) :: message
+ ! ---------------------------------------------------------------------------------------
+ ! * Convert From C++ to Fortran
+ ! ---------------------------------------------------------------------------------------
+ call c_f_pointer(handle_mpar_struct, mpar_struct)
+ ! Start subroutine
+ err=0; message="read_param.f90 - overwriteParam"
+
+ ! initalize the structure with allocatable components
+ call allocLocal(mpar_meta,mpar_struct, &
+ gru_struc(index_gru)%hruInfo(index_hru)%nSnow,&
+ gru_struc(index_gru)%hruInfo(index_hru)%nSoil,&
+ err,message); if(err/=0)then; message=trim(message); print*, message; return; endif;
+
+ ! Set the basin parameters with the default values
+ do ivar=1, size(localParFallback)
+ dpar_array(iVar) = localParFallback(iVar)%default_val
+ end do
+
+ call pOverwrite(type_array(iLookTYPE%vegTypeIndex), & ! vegetation category
+ type_array(iLookTYPE%soilTypeIndex),& ! soil category
+ dpar_array(:),err,message) ! default model parameters
+
+ do ivar=1, size(localParFallback)
+ do iDat=1, size(mpar_struct%var(iVar)%dat)
+ mpar_struct%var(iVar)%dat(iDat) = dpar_array(iVar)
+ end do
+ end do
+
+ do iVar=1, size(basinParFallback)
+ bpar_array(iVar) = basinParFallback(iVar)%default_val
+ end do
+
+end subroutine overwriteParam
+
+
+subroutine readParamFromNetCDF(param_ncid, index_gru, index_hru, start_index_gru, &
+ num_vars, bpar_array_size, handle_mpar_struct, bpar_array, err) bind(C, name="readParamFromNetCDF")
+ USE netcdf
+ USE netcdf_util_module,only:netcdf_err ! netcdf error handling function
+
+ USE data_types,only:var_dlength
+ USE get_ixname_module,only:get_ixparam,get_ixbpar ! access function to find index of elements in structure
+
+ USE globalData,only:index_map,gru_struc ! mapping from global HRUs to the elements in the data structures
+ USE globalData,only:integerMissing ! missing integer
+
+ implicit none
+ ! dummy variables
+ integer(c_int),intent(in) :: param_ncid
+ integer(c_int),intent(in) :: index_gru
+ integer(c_int),intent(in) :: index_hru
+ integer(c_int),intent(in) :: start_index_gru
+ integer(c_int),intent(in) :: num_vars
+ integer(c_int),intent(in) :: bpar_array_size
+ type(c_ptr), intent(in),value :: handle_mpar_struct
+ real(c_double), intent(out) :: bpar_array(bpar_array_size)
+ integer(c_int), intent(out) :: err
+ ! define local variables
+ type(var_dlength),pointer :: mpar_struct ! model parameters
+
+ character(len=256) :: message ! error message
+ character(len=1024) :: cmessage ! error message for downwind routine
+ integer(i4b) :: localHRU_ix ! index of HRU within data structure
+ integer(i4b) :: ixParam ! index of the model parameter in the data structure
+
+ ! indices/metadata in the NetCDF file
+ integer(i4b) :: num_dims ! number of dimensions
+
+ integer(i4b) :: ivarid ! variable index
+ character(LEN=64) :: dimName ! dimension name
+
+ character(LEN=64) :: parName ! parameter name
+ integer(i4b) :: nSoil_file ! number of soil layers in the file
+ integer(i4b) :: idim_list(2) ! list of dimension ids
+ ! data in the netcdf file
+ integer(i4b) :: parLength ! length of the parameter data
+ real(dp),allocatable :: parVector(:) ! model parameter vector
+ integer(i4b) :: fHRU ! index of HRU in input file
+ integer(i4b) :: netcdf_index
+
+ ! ---------------------------------------------------------------------------------------
+ ! * Convert From C++ to Fortran
+ ! ---------------------------------------------------------------------------------------
+ call c_f_pointer(handle_mpar_struct, mpar_struct)
+ err=0; message="read_param.f90 - readParamFromNetCDF/"
+
+
+ ! **********************************************************************************************
+ ! * read the local parameters and the basin parameters
+ ! **********************************************************************************************
+ do ivarid=1,num_vars
+ ! get the parameter name
+ err=nf90_inquire_variable(param_ncid, ivarid, name=parName)
+ call netcdf_err(err,message)
+ if (err/=0) then
+ err=20
+ print*,message
+ return
+ end if
+
+ ! get the local parameters
+ ixParam = get_ixparam( trim(parName) )
+ if(ixParam/=integerMissing)then
+ ! **********************************************************************************************
+ ! * read the local parameters
+ ! **********************************************************************************************
+ ! get the variable shape
+ err=nf90_inquire_variable(param_ncid, ivarid, nDims=num_dims, dimids=idim_list)
+ if(err/=0)then
+ message=trim(message)//trim(cmessage)
+ print*, message
+ return
+ end if
+
+ ! get the length of the depth dimension (if it exists)
+ if(num_dims==2)then
+ ! get the information on the 2nd dimension for 2-d variables
+ err=nf90_inquire_dimension(param_ncid, idim_list(2), dimName, nSoil_file)
+ if(err/=0)then
+ message=trim(message)//trim(cmessage)
+ print*, message
+ return
+ end if
+
+ ! check that it is the depth dimension
+ if(trim(dimName)/='depth')then
+ message=trim(message)//'expect 2nd dimension of 2-d variable to be depth (dimension name = '//trim(dimName)//')'
+ err=20; return
+ endif
+
+ ! ! check that the dimension length is correct
+ if(size(mpar_struct%var(ixParam)%dat) /= nSoil_file)then
+ message=trim(message)//'unexpected number of soil layers in parameter file'
+ err=20; return
+ endif
+
+ ! define parameter length
+ parLength = nSoil_file
+ else
+ parLength = 1
+ endif ! if two dimensions
+
+ ! allocate space for model parameters
+ allocate(parVector(parLength),stat=err)
+ if(err/=0)then
+ message=trim(message)//'problem allocating space for parameter vector'
+ err=20; return
+ endif
+
+
+ localHRU_ix=index_map(index_hru)%localHRU_ix
+ fHRU = gru_struc(index_gru)%hruInfo(localHRU_ix)%hru_nc
+ ! read parameter data
+ select case(num_dims)
+ case(1); err=nf90_get_var(param_ncid, ivarid, parVector, start=(/fHRU/), count=(/1/) )
+ case(2); err=nf90_get_var(param_ncid, ivarid, parVector, start=(/fHRU,1/), count=(/1,nSoil_file/) )
+ case default; err=20; message=trim(message)//'unexpected number of dimensions for parameter '//trim(parName)
+ end select
+
+ ! error check for the parameter read
+ if(err/=0)then
+ message=trim(message)//trim(cmessage)
+ print*, message
+ return
+ end if
+
+ ! populate parameter structures
+ select case(num_dims)
+ case(1); mpar_struct%var(ixParam)%dat(:) = parVector(1) ! also distributes scalar across depth dimension
+ case(2); mpar_struct%var(ixParam)%dat(:) = parVector(:)
+ case default; err=20; message=trim(message)//'unexpected number of dimensions for parameter '//trim(parName)
+ end select
+
+ ! deallocate space for model parameters
+ deallocate(parVector,stat=err)
+ if(err/=0)then
+ message=trim(message)//'problem deallocating space for parameter vector'
+ print*, message
+ err=20; return
+ endif
+
+ ! **********************************************************************************************
+ ! * read the basin parameters
+ ! **********************************************************************************************
+
+ ! get the basin parameters
+ else
+ ! get the parameter index
+ ixParam = get_ixbpar( trim(parName) )
+
+ ! allow extra variables in the file that are not used
+ if(ixParam==integerMissing) cycle
+
+ ! read parameter data
+ netcdf_index = start_index_gru + index_gru - 1
+ err=nf90_get_var(param_ncid, ivarid, bpar_array(ixParam), start=(/netcdf_index/))
+ if(err/=0)then
+ message=trim(message)//trim(cmessage)
+ print*, message
+ return
+ end if
+ endif
+
+ end do ! (looping through the parameters in the NetCDF file)
+
+end subroutine readParamFromNetCDF
+
+
+end module read_param_module
\ No newline at end of file
diff --git a/build/source/actors/file_access_actor/fortran_code/read_param_all_hru.f90 b/build/source/actors/file_access_actor/fortran_code/read_param_all_hru.f90
index 0e49f7c..56d8353 100644
--- a/build/source/actors/file_access_actor/fortran_code/read_param_all_hru.f90
+++ b/build/source/actors/file_access_actor/fortran_code/read_param_all_hru.f90
@@ -112,22 +112,22 @@ subroutine read_param_file_access_actor(startGRU,num_gru,err) bind(C, name="read
! * open files, etc.
! **********************************************************************************************
- infile = trim(SETTINGS_PATH)//trim(PARAMETER_TRIAL) ! build filename
+ ! infile = trim(SETTINGS_PATH)//trim(PARAMETER_TRIAL) ! build filename
- ! check whether the user-specified file exists and warn if it does not
- inquire(file=trim(infile),exist=fexist)
- if (.not.fexist) then
- write(*,'(A)') NEW_LINE('A')//'!! WARNING: trial parameter file not found; proceeding instead with other default parameters; check path in file manager input if this was not the desired behavior'//NEW_LINE('A')
- return
- endif
+ ! ! check whether the user-specified file exists and warn if it does not
+ ! inquire(file=trim(infile),exist=fexist)
+ ! if (.not.fexist) then
+ ! write(*,'(A)') NEW_LINE('A')//'!! WARNING: trial parameter file not found; proceeding instead with other default parameters; check path in file manager input if this was not the desired behavior'//NEW_LINE('A')
+ ! return
+ ! endif
- ! open trial parameters file if it exists
- call nc_file_open(trim(infile),nf90_nowrite,ncid,err,cmessage)
- if(err/=0)then; message=trim(message)//trim(cmessage); return; end if
+ ! ! open trial parameters file if it exists
+ ! call nc_file_open(trim(infile),nf90_nowrite,ncid,err,cmessage)
+ ! if(err/=0)then; message=trim(message)//trim(cmessage); return; end if
- ! get the number of variables in the parameter file
- err=nf90_inquire(ncid, nDimensions=nDims, nVariables=nVars)
- call netcdf_err(err,message); if (err/=0) then; err=20; return; end if
+ ! ! get the number of variables in the parameter file
+ ! err=nf90_inquire(ncid, nDimensions=nDims, nVariables=nVars)
+ ! call netcdf_err(err,message); if (err/=0) then; err=20; return; end if
! initialize the number of HRUs
nHRU_file=integerMissing
diff --git a/build/source/actors/hru_actor/hru_actor.cpp b/build/source/actors/hru_actor/hru_actor.cpp
index 0be9c65..879daf5 100644
--- a/build/source/actors/hru_actor/hru_actor.cpp
+++ b/build/source/actors/hru_actor/hru_actor.cpp
@@ -52,12 +52,38 @@ behavior hru_actor(stateful_actor<hru_state>* self, int refGRU, int indxGRU,
}
- Initialize_HRU(self);
+ initHRU(&self->state.indxGRU, &self->state.num_steps, self->state.handle_lookupStruct, self->state.handle_forcStat,
+ self->state.handle_progStat, self->state.handle_diagStat, self->state.handle_fluxStat, self->state.handle_indxStat,
+ self->state.handle_bvarStat, self->state.handle_timeStruct, self->state.handle_forcStruct, self->state.handle_attrStruct,
+ self->state.handle_typeStruct, self->state.handle_idStruct,self->state.handle_indxStruct, self->state.handle_mparStruct,
+ self->state.handle_progStruct, self->state.handle_diagStruct, self->state.handle_fluxStruct,self->state.handle_bparStruct,
+ self->state.handle_bvarStruct, self->state.handle_dparStruct, self->state.handle_startTime, self->state.handle_finshTime,
+ self->state.handle_refTime,self->state.handle_oldTime, &self->state.err);
+ if (self->state.err != 0) {
+ aout(self) << "Error: HRU_Actor - Initialize - HRU = " << self->state.indxHRU <<
+ " - indxGRU = " << self->state.indxGRU << " - refGRU = "<< self->state.refGRU << std::endl;
+ aout(self) << "Error = " << self->state.err << "\n";
+ self->quit();
+ }
+
+
+ // Get attributes
+ self->send(self->state.file_access_actor, get_attributes_v, self->state.refGRU, self);
+
+
self->state.hru_timing.updateEndPoint("total_duration");
self->send(self, start_hru_v);
return {
// Starts the HRU and tells it to ask for data from the file_access_actor
+ [=](get_attributes, std::vector<double> attr_array, std::vector<int> type_array,
+ std::vector<long int> id_array) {
+
+ aout(self) << "Received Attribute Information \n";
+
+ },
+
+
[=](start_hru) {
self->state.hru_timing.updateStartPoint("total_duration");
@@ -218,43 +244,8 @@ behavior hru_actor(stateful_actor<hru_state>* self, int refGRU, int indxGRU,
void Initialize_HRU(stateful_actor<hru_state>* self) {
self->state.hru_timing.updateStartPoint("init_duration");
-
- initHRU(&self->state.indxGRU,
- &self->state.num_steps,
- self->state.handle_lookupStruct,
- self->state.handle_forcStat,
- self->state.handle_progStat,
- self->state.handle_diagStat,
- self->state.handle_fluxStat,
- self->state.handle_indxStat,
- self->state.handle_bvarStat,
- self->state.handle_timeStruct,
- self->state.handle_forcStruct,
- self->state.handle_attrStruct,
- self->state.handle_typeStruct,
- self->state.handle_idStruct,
- self->state.handle_indxStruct,
- self->state.handle_mparStruct,
- self->state.handle_progStruct,
- self->state.handle_diagStruct,
- self->state.handle_fluxStruct,
- self->state.handle_bparStruct,
- self->state.handle_bvarStruct,
- self->state.handle_dparStruct,
- self->state.handle_startTime,
- self->state.handle_finshTime,
- self->state.handle_refTime,
- self->state.handle_oldTime,
- &self->state.err);
- if (self->state.err != 0) {
- aout(self) << "Error: HRU_Actor - Initialize - HRU = " << self->state.indxHRU <<
- " - indxGRU = " << self->state.indxGRU << " - refGRU = "<< self->state.refGRU << std::endl;
- aout(self) << "Error = " << self->state.err << "\n";
- self->quit();
- return;
- }
// Need to send a message to the file_access_actor for the data
diff --git a/build/source/engine/pOverwrite.f90 b/build/source/engine/pOverwrite.f90
index d90b5bd..af0670c 100755
--- a/build/source/engine/pOverwrite.f90
+++ b/build/source/engine/pOverwrite.f90
@@ -26,71 +26,71 @@ public::pOverwrite
contains
- ! ************************************************************************************************
- ! public subroutine pOverwrite: use Noah tables to overwrite default model parameters
- ! ************************************************************************************************
- subroutine pOverwrite(ixVeg,ixSoil,defaultParam,err,message)
- ! SUMMA dictionary
- USE var_lookup,only:iLookPARAM ! named variables for elements of the data structures
- ! Noah table dimensions
- USE module_sf_noahlsm, only: LUCATS ! dimension of the vegetation tables (number of land use catagories)
- USE module_sf_noahlsm, only: NSLTYPE ! dimension of the soil tables
- ! Noah vegetation tables
- USE NOAHMP_VEG_PARAMETERS, only: Z0MVT ! Noah-MP: momentum roughness length (m)
- USE NOAHMP_VEG_PARAMETERS, only: HVT ! Noah-MP: height at top of canopy (m)
- USE NOAHMP_VEG_PARAMETERS, only: HVB ! Noah-MP: height at bottom of canopy (m)
- USE NOAHMP_VEG_PARAMETERS, only: DLEAF ! Noah-MP: characteristic leaf dimension (m)
- USE NOAHMP_VEG_PARAMETERS, only: VCMX25 ! Noah-MP: maximum Rubisco carboxylation rate (umol m-2 s-1)
- USE NOAHMP_VEG_PARAMETERS, only: MP ! Noah-MP: slope of conductance-photosynthesis relationship (-)
- ! Noah soil tables
- USE module_sf_noahlsm, only: theta_res, theta_sat, vGn_alpha, vGn_n, k_soil ! van Genutchen soil parameters
- USE module_sf_noahlsm, only: REFSMC ! Noah-MP: reference volumetric soil moisture content (-)
- USE module_sf_noahlsm, only: WLTSMC ! Noah-MP: volumetric soil moisture content when plants are wilting (-)
- implicit none
- ! define input
- integer(i4b),intent(in) :: ixVeg ! vegetation category
- integer(i4b),intent(in) :: ixSoil ! soil category
- ! define output
- real(dp),intent(inout) :: defaultParam(:) ! default model parameters
- integer(i4b),intent(out) :: err ! error code
- character(*),intent(out) :: message ! error message
- ! Start procedure here
- err=0; message="pOverwrite/"
+! ************************************************************************************************
+! public subroutine pOverwrite: use Noah tables to overwrite default model parameters
+! ************************************************************************************************
+subroutine pOverwrite(ixVeg,ixSoil,defaultParam,err,message)
+ ! SUMMA dictionary
+ USE var_lookup,only:iLookPARAM ! named variables for elements of the data structures
+ ! Noah table dimensions
+ USE module_sf_noahlsm, only: LUCATS ! dimension of the vegetation tables (number of land use catagories)
+ USE module_sf_noahlsm, only: NSLTYPE ! dimension of the soil tables
+ ! Noah vegetation tables
+ USE NOAHMP_VEG_PARAMETERS, only: Z0MVT ! Noah-MP: momentum roughness length (m)
+ USE NOAHMP_VEG_PARAMETERS, only: HVT ! Noah-MP: height at top of canopy (m)
+ USE NOAHMP_VEG_PARAMETERS, only: HVB ! Noah-MP: height at bottom of canopy (m)
+ USE NOAHMP_VEG_PARAMETERS, only: DLEAF ! Noah-MP: characteristic leaf dimension (m)
+ USE NOAHMP_VEG_PARAMETERS, only: VCMX25 ! Noah-MP: maximum Rubisco carboxylation rate (umol m-2 s-1)
+ USE NOAHMP_VEG_PARAMETERS, only: MP ! Noah-MP: slope of conductance-photosynthesis relationship (-)
+ ! Noah soil tables
+ USE module_sf_noahlsm, only: theta_res, theta_sat, vGn_alpha, vGn_n, k_soil ! van Genutchen soil parameters
+ USE module_sf_noahlsm, only: REFSMC ! Noah-MP: reference volumetric soil moisture content (-)
+ USE module_sf_noahlsm, only: WLTSMC ! Noah-MP: volumetric soil moisture content when plants are wilting (-)
+ implicit none
+ ! define input
+ integer(i4b),intent(in) :: ixVeg ! vegetation category
+ integer(i4b),intent(in) :: ixSoil ! soil category
+ ! define output
+ real(dp),intent(out) :: defaultParam(:) ! default model parameters
+ integer(i4b),intent(out) :: err ! error code
+ character(*),intent(out) :: message ! error message
+ ! Start procedure here
+ err=0; message="pOverwrite/"
- ! define vegetation class
- if(ixVeg < 1)then; err=20; message=trim(message)//'index for vegetation type must be > 0'; return; end if
- if(ixVeg > LUCATS)then
- write(message,'(2(a,i0),a)')trim(message)//'index for vegetation type is greater than dimension of vegetation table [ixVeg = ', ixVeg, &
- '; LUCATS = ', LUCATS, ']'
- err=20; return
- end if
+ ! define vegetation class
+ if(ixVeg < 1)then; err=20; message=trim(message)//'index for vegetation type must be > 0'; return; end if
+ if(ixVeg > LUCATS)then
+ write(message,'(2(a,i0),a)')trim(message)//'index for vegetation type is greater than dimension of vegetation table [ixVeg = ', ixVeg, &
+ '; LUCATS = ', LUCATS, ']'
+ err=20; return
+ end if
- ! define soil class
- if(ixSoil < 1)then; err=20; message=trim(message)//'index for soil type must be > 0'; return; end if
- if(ixSoil > NSLTYPE)then
- write(message,'(2(a,i0),a)')trim(message)//'index for soil type is greater than dimension of soil table [ixSoil = ', ixSoil, &
- '; NSLTYPE = ', NSLTYPE, ']'
- err=20; return
- end if
+ ! define soil class
+ if(ixSoil < 1)then; err=20; message=trim(message)//'index for soil type must be > 0'; return; end if
+ if(ixSoil > NSLTYPE)then
+ write(message,'(2(a,i0),a)')trim(message)//'index for soil type is greater than dimension of soil table [ixSoil = ', ixSoil, &
+ '; NSLTYPE = ', NSLTYPE, ']'
+ err=20; return
+ end if
- ! include parameters from the vegetation tables
- defaultParam(iLookPARAM%heightCanopyTop) = HVT(ixVeg) ! Noah-MP: height at top of canopy (m)
- defaultParam(iLookPARAM%heightCanopyBottom) = HVB(ixVeg) ! Noah-MP: height at bottom of canopy (m)
- defaultParam(iLookPARAM%z0Canopy) = Z0MVT(ixVeg) ! Noah-MP: momentum roughness length (m)
- defaultParam(iLookPARAM%leafDimension) = DLEAF(ixVeg) ! Noah-MP: characteristic leaf dimension (m)
- defaultParam(iLookPARAM%vcmax25_canopyTop) = VCMX25(ixVeg) ! Noah-MP: maximum Rubisco carboxylation rate (umol m-2 s-1)
- defaultParam(iLookPARAM%cond2photo_slope) = MP(ixVeg) ! Noah-MP: slope of conductance-photosynthesis relationship (-)
+ ! include parameters from the vegetation tables
+ defaultParam(iLookPARAM%heightCanopyTop) = HVT(ixVeg) ! Noah-MP: height at top of canopy (m)
+ defaultParam(iLookPARAM%heightCanopyBottom) = HVB(ixVeg) ! Noah-MP: height at bottom of canopy (m)
+ defaultParam(iLookPARAM%z0Canopy) = Z0MVT(ixVeg) ! Noah-MP: momentum roughness length (m)
+ defaultParam(iLookPARAM%leafDimension) = DLEAF(ixVeg) ! Noah-MP: characteristic leaf dimension (m)
+ defaultParam(iLookPARAM%vcmax25_canopyTop) = VCMX25(ixVeg) ! Noah-MP: maximum Rubisco carboxylation rate (umol m-2 s-1)
+ defaultParam(iLookPARAM%cond2photo_slope) = MP(ixVeg) ! Noah-MP: slope of conductance-photosynthesis relationship (-)
- ! include parameters from the soil tables
- defaultParam(iLookPARAM%k_soil) = k_soil(ixSoil) ! hydraulic conductivity (m s-1)
- defaultParam(iLookPARAM%theta_res) = theta_res(ixSoil) ! residual volumetric liquid water content (-)
- defaultParam(iLookPARAM%theta_sat) = theta_sat(ixSoil) ! soil porosity (-)
- defaultParam(iLookPARAM%vGn_alpha) = vGn_alpha(ixSoil) ! van Genutchen "alpha" parameter (m-1)
- defaultParam(iLookPARAM%vGn_n) = vGn_n(ixSoil) ! van Genutchen "n" parameter (-)
- defaultParam(iLookPARAM%critSoilTranspire) = REFSMC(ixSoil) ! Noah-MP: reference volumetric soil moisture content (-)
- defaultParam(iLookPARAM%critSoilWilting) = WLTSMC(ixSoil) ! Noah-MP: volumetric soil moisture content when plants are wilting (-)
+ ! include parameters from the soil tables
+ defaultParam(iLookPARAM%k_soil) = k_soil(ixSoil) ! hydraulic conductivity (m s-1)
+ defaultParam(iLookPARAM%theta_res) = theta_res(ixSoil) ! residual volumetric liquid water content (-)
+ defaultParam(iLookPARAM%theta_sat) = theta_sat(ixSoil) ! soil porosity (-)
+ defaultParam(iLookPARAM%vGn_alpha) = vGn_alpha(ixSoil) ! van Genutchen "alpha" parameter (m-1)
+ defaultParam(iLookPARAM%vGn_n) = vGn_n(ixSoil) ! van Genutchen "n" parameter (-)
+ defaultParam(iLookPARAM%critSoilTranspire) = REFSMC(ixSoil) ! Noah-MP: reference volumetric soil moisture content (-)
+ defaultParam(iLookPARAM%critSoilWilting) = WLTSMC(ixSoil) ! Noah-MP: volumetric soil moisture content when plants are wilting (-)
- end subroutine pOverwrite
+end subroutine pOverwrite
end module pOverwrite_module
--
GitLab