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 5bbfe726abd1231020ebac3abd8bc1a26d7a57aa..68c20d046c32ae4f6a3a014ffb74e04e484bf6af 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 @@ -1,42 +1,47 @@ #pragma once extern "C" { - + // initalizeFileAccessActor + void ffile_info(int* indxGRU, void* forcFileInfo, int* numFiles, int* err); + void mDecisions(int* numSteps, int* err); void read_pinit_C(int* err); - void read_vegitationTables(int* err); - void initFailedHRUTracker(int* numGRU); + void def_output(void* handle_ncid, int* startGRU, int* numGRU, int* numHRU, int* err); + + // Attributes Files- called inside initalizeFileAccessActor + void allocateAttributeStructures(int* index_gru, int* index_hru, void* handle_attr_struct, + void* handle_type_struct, void* handle_id_struct, int* err); + void openAttributeFile(int* att_ncid, 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 - called inside initalizeFileAccessActor + void allocateParamStructures(int* index_gru, int* index_hru, void* handle_dpar_struct, + void* handle_mpar_struct, void* handle_bpar_struct, int* err); + 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, int* type_array_size); + void overwriteParam(int* index_gru, int* index_hru, + void* handle_type_struct, void* handle_dpar_struct, void* handle_mpar_struct, + void* handle_bpar_struct, int* err); + void readParamFromNetCDF(int* param_ncid, int* index_gru, int* index_hru, int* start_index_gru, + int* num_var_param, void* handle_mpar_struct, void* _handle_bpar_struct, int* err); void updateFailed(int* indxHRU); void resetFailedArray(); void resetOutputCounter(int* indxGRU); - - void mDecisions_C(int* numSteps, int* err); - - void ffile_info_C(int* indxGRU, void* forcFileInfo, int* numFiles, int* err); - - void Init_OutputStruct(void* forcFileInfo, int* maxSteps, int* nGru, int* err); void FileAccessActor_ReadForcing(void* forcFileInfo, int* currentFile, int* stepsInFile, int* startGRU, int* numGRU, int* err); - void FileAccessActor_WriteOutput(void* handle_ncid, - int* stepsInCurrentFile, int*indxGRU, int*indxHRU, int* err); - void FileAccessActor_DeallocateStructures(void* handle_forcFileInfo, void* handle_ncid); - 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 readAttributeFileAccessActor(int* num_gru, int* err); - - // void overwriteParam(int* num_gru, int* err); - - // void readParamFileAccessActor(int* start_gru, int* num_gru, int* err); // Writing to NETCDF void writeParamToNetCDF(void* handle_ncid, int* index_gru, int* index_hru, @@ -57,33 +62,9 @@ extern "C" { void* handle_time_struct, int* err); - // Attributes Files - void allocateAttributeStructures(int* index_gru, int* index_hru, void* handle_attr_struct, - void* handle_type_struct, void* handle_id_struct, int* err); - - void openAttributeFile(int* att_ncid, 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 allocateParamStructures(int* index_gru, int* index_hru, void* handle_dpar_struct, - void* handle_mpar_struct, void* handle_bpar_struct, int* err); - 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, int* type_array_size); - - void overwriteParam(int* index_gru, int* index_hru, - void* handle_type_struct, void* handle_dpar_struct, void* handle_mpar_struct, - void* handle_bpar_struct, int* err); - - void readParamFromNetCDF(int* param_ncid, int* index_gru, int* index_hru, int* start_index_gru, - int* num_var_param, void* handle_mpar_struct, void* _handle_bpar_struct, int* err); } diff --git a/build/includes/file_access_actor/forcing_file_info.hpp b/build/includes/file_access_actor/forcing_file_info.hpp index 14945aabb5352fa0dad6ec0b9da49f00d0e2910a..3a3ef8c2f2c8565e8f895b9caea2648dfbfe8115 100644 --- a/build/includes/file_access_actor/forcing_file_info.hpp +++ b/build/includes/file_access_actor/forcing_file_info.hpp @@ -1,5 +1,8 @@ #pragma once +#include <vector> + + class Forcing_File_Info { private: int file_ID; @@ -17,4 +20,12 @@ class Forcing_File_Info { void updateNumSteps(int num_steps); +}; + +struct Forcing_Info { + int num_vars; + int num_timesteps; + std::vector<int> index_forc_var; + std::vector<int> ncid_var; + }; \ No newline at end of file diff --git a/build/includes/global/message_atoms.hpp b/build/includes/global/message_atoms.hpp index d9c4d1dd2bc6ea17b1b51b62b0df26ea300bda94..61875e3a1b084077c8d5adeac3124e9e020f0b40 100644 --- a/build/includes/global/message_atoms.hpp +++ b/build/includes/global/message_atoms.hpp @@ -39,7 +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) + CAF_ADD_ATOM(summa, get_attributes_params) // HRU Actor CAF_ADD_ATOM(summa, run_hru) CAF_ADD_ATOM(summa, start_hru) diff --git a/build/includes/job_actor/job_actor_subroutine_wrappers.hpp b/build/includes/job_actor/job_actor_subroutine_wrappers.hpp index d215986be4fb30a84b4c880066e32b848e65c9c8..07622122f95a06bd0d06a57628f922ef82c5d1e6 100644 --- a/build/includes/job_actor/job_actor_subroutine_wrappers.hpp +++ b/build/includes/job_actor/job_actor_subroutine_wrappers.hpp @@ -1,9 +1,6 @@ #pragma once extern "C" { - void initGlobals(char const*str1, int* totalGRUs, int* totalHRUs, - int* numGRUs, int* numHRUs, int* startGRUIndex, int* err); - void setTimesDirsAndFiles(char const* file_manager, int* err); void defineGlobalData(int* start_gru_index, int* err); diff --git a/build/makefile_sundials b/build/makefile_sundials index 4153042a394a23247aabae9998ef62efba54b241..f23974d234895df71e79c59d35ff3836c8a3a0b5 100644 --- a/build/makefile_sundials +++ b/build/makefile_sundials @@ -161,8 +161,6 @@ INTERFACE = $(patsubst %, $(ACTORS_DIR)/global/%, $(SUMMA_INTERFACE)) SUMMA_FILEACCESS_INTERFACE = \ - initOutputStruc.f90 \ - deallocateOutputStruc.f90 \ cppwrap_fileAccess.f90 \ read_attribute.f90 \ read_param.f90 \ 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 c0c452fcd096168259510f17b2215e6a01ebac0a..4e84a05c747248c62f4c2ccabe5b498bb9c9080c 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,13 +128,7 @@ 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]; - // std::vector<double> bpar_array_to_send = self->state.bpar_arrays_for_hrus[ref_gru-1]; - // std::vector<double> dpar_array_to_send = self->state.dpar_arrays_for_hrus[ref_gru-1]; + [=] (get_attributes_params, int ref_gru, caf::actor actor_to_respond) { void* handle_attr_struct = self->state.attr_structs_for_hrus[ref_gru-1]; std::vector<double> attr_struct_to_send = get_var_d(handle_attr_struct); @@ -153,7 +147,7 @@ behavior file_access_actor(stateful_actor<file_access_state>* self, int start_gr void* handle_mpar_struct = self->state.mpar_structs_for_hrus[ref_gru-1]; std::vector<std::vector<double>> mpar_struct_to_send = get_var_dlength(handle_mpar_struct); - self->send(actor_to_respond, get_attributes_v, attr_struct_to_send, + self->send(actor_to_respond, get_attributes_params_v, attr_struct_to_send, type_struct_to_send, id_struct_to_send, bpar_struct_to_send, dpar_struct_to_send, mpar_struct_to_send); @@ -303,8 +297,10 @@ behavior file_access_actor(stateful_actor<file_access_state>* self, int start_gr void initalizeFileAccessActor(stateful_actor<file_access_state>* self) { int indx = 1; int err = 0; - - ffile_info_C(&indx, self->state.handle_forcing_file_info, &self->state.numFiles, &err); + + // read information on model forcing files + ffile_info(&indx, + self->state.handle_forcing_file_info, &self->state.numFiles, &err); if (err != 0) { aout(self) << "Error: ffile_info_C - File_Access_Actor \n"; std::string function = "ffile_info_C"; @@ -313,7 +309,8 @@ void initalizeFileAccessActor(stateful_actor<file_access_state>* self) { return; } - mDecisions_C(&self->state.num_steps, &err); + // save model decisions as named integers + mDecisions(&self->state.num_steps, &err); if (err != 0) { aout(self) << "Error: mDecisions - FileAccess Actor \n"; std::string function = "mDecisions_C"; @@ -366,7 +363,6 @@ void initalizeFileAccessActor(stateful_actor<file_access_state>* self) { } } - int readForcing(stateful_actor<file_access_state>* self, int currentFile) { // Check if we have already loaded this file if(self->state.forcing_file_list[currentFile -1].isFileLoaded()) { @@ -399,7 +395,6 @@ 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); 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 86d95d1ca35ac8cabec9f50cb48c6a46e6d9d46f..18d68a0d4d43516682d1330a9e25a9af52fa5dc0 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 @@ -10,61 +10,11 @@ module cppwrap_fileAccess implicit none - public::ffile_info_C - public::mDecisions_C public::initFailedHRUTracker public::FileAccessActor_ReadForcing contains -subroutine ffile_info_C(indxGRU, handle_forcFileInfo, numFiles, err) bind(C, name='ffile_info_C') - USE ffile_info_module,only:ffile_info - - implicit none - ! dummy variables - - integer(c_int), intent(in) :: indxGRU - type(c_ptr), intent(in), value :: handle_forcFileInfo - integer(c_int), intent(out) :: numFiles - integer(c_int) :: err - ! local variables - type(file_info_array),pointer :: forcFileInfo - character(len=256) :: message - - call c_f_pointer(handle_forcFileInfo, forcFileInfo) - - call ffile_info(indxGRU,forcFileInfo,numFiles,err,message) - if(err/=0)then - message=trim(message) - write(*,*) message - return - endif - -end subroutine ffile_info_C - -! THis subroutine needs to be called after the ffile_info_C subroutine -! It might make more sense to have this in the global data inialization -! subroutine. But here it sits -subroutine mDecisions_C(num_steps,err) bind(C, name='mDecisions_C') - USE mDecisions_module,only:mDecisions ! module to read model decisions - - implicit none - ! dummy variables - integer(c_int),intent(inout) :: num_steps - integer(c_int),intent(inout) :: err ! Error Code - ! local variables - character(len=256) :: message ! error message - character(LEN=256) :: cmessage ! error message of downwind routine - - call mDecisions(num_steps,err,cmessage) - if(err/=0)then - message=trim(message)//trim(cmessage) - print*, message - return - endif - -end subroutine mDecisions_C - ! Read in the inital parameters, from the txt files that are give to summa as input ! LocalParamInfo.txt @@ -107,16 +57,15 @@ subroutine read_vegitationTables(err) bind(C, name="read_vegitationTables") implicit none - integer(c_int),intent(inout) :: err ! Error Code - + integer(c_int),intent(inout) :: err ! Error Code err = 0 ! read Noah soil and vegetation tables - call soil_veg_gen_parm(trim(SETTINGS_PATH)//trim(VEGPARM), & ! filename for vegetation table - trim(SETTINGS_PATH)//trim(SOILPARM), & ! filename for soils table - trim(SETTINGS_PATH)//trim(GENPARM), & ! filename for general table - trim(model_decisions(iLookDECISIONS%vegeParTbl)%cDecision), & ! classification system used for vegetation - trim(model_decisions(iLookDECISIONS%soilCatTbl)%cDecision)) ! classification system used for soils + call soil_veg_gen_parm(trim(SETTINGS_PATH)//trim(VEGPARM), & ! filename for vegetation table + trim(SETTINGS_PATH)//trim(SOILPARM), & ! filename for soils table + trim(SETTINGS_PATH)//trim(GENPARM), & ! filename for general table + trim(model_decisions(iLookDECISIONS%vegeParTbl)%cDecision), & ! classification system used for vegetation + trim(model_decisions(iLookDECISIONS%soilCatTbl)%cDecision)) ! classification system used for soils ! read Noah-MP vegetation tables call read_mp_veg_parameters(trim(SETTINGS_PATH)//trim(MPTABLE), & ! filename for Noah-MP table diff --git a/build/source/actors/file_access_actor/fortran_code/deallocateOutputStruc.f90 b/build/source/actors/file_access_actor/fortran_code/deallocateOutputStruc.f90 deleted file mode 100644 index c0341661c2d659e80e33846a8c853a2ef620494d..0000000000000000000000000000000000000000 --- a/build/source/actors/file_access_actor/fortran_code/deallocateOutputStruc.f90 +++ /dev/null @@ -1,133 +0,0 @@ -module summaActors_deallocateOuptutStruct - USE nrtype - implicit none - - contains - -subroutine deallocateData_output(dataStruct) - USE data_types,only:gru_hru_time_doubleVec, & - gru_hru_time_intVec, & - gru_hru_time_flagVec, & - gru_hru_time_int, & - gru_hru_int, & - gru_hru_time_int8, & - gru_hru_time_double, & - gru_hru_double, & - gru_double - implicit none - class(*),intent(inout) :: dataStruct - ! local variables - integer(i4b) :: iGRU - integer(i4b) :: iHRU - integer(i4b) :: iVar - integer(i4b) :: iTim - - select type(dataStruct) - class is (gru_hru_time_doubleVec) - do iGRU = 1, size(dataStruct%gru(:)) - do iHRU = 1, size(dataStruct%gru(iGRU)%hru(:)) - do iVar = 1, size(dataStruct%gru(iGRU)%hru(iHRU)%var(:)) - do iTim = 1, size(dataStruct%gru(iGRU)%hru(iHRU)%var(iVar)%tim(:)) - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var(iVar)%tim(iTim)%dat) - end do ! Time - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var(iVar)%tim) - end do ! var - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var) - end do ! hru - deallocate(dataStruct%gru(iGRU)%hru) - end do ! gru - deallocate(dataStruct%gru) - - class is (gru_hru_time_intVec) - do iGRU = 1, size(dataStruct%gru(:)) - do iHRU = 1, size(dataStruct%gru(iGRU)%hru(:)) - do iVar = 1, size(dataStruct%gru(iGRU)%hru(iHRU)%var(:)) - do iTim = 1, size(dataStruct%gru(iGRU)%hru(iHRU)%var(iVar)%tim(:)) - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var(iVar)%tim(iTim)%dat) - end do ! Time - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var(iVar)%tim) - end do ! var - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var) - end do ! hru - deallocate(dataStruct%gru(iGRU)%hru) - end do ! gru - deallocate(dataStruct%gru) - - class is (gru_hru_time_flagVec) - do iGRU = 1, size(dataStruct%gru(:)) - do iHRU = 1, size(dataStruct%gru(iGRU)%hru(:)) - do iTim = 1, size(dataStruct%gru(iGRU)%hru(iHRU)%tim(:)) - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%tim(iTim)%dat) - end do ! Time - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%tim) - end do ! hru - deallocate(dataStruct%gru(iGRU)%hru) - end do ! gru - deallocate(dataStruct%gru) - - class is (gru_hru_time_int) - do iGRU = 1, size(dataStruct%gru(:)) - do iHRU = 1, size(dataStruct%gru(iGRU)%hru(:)) - do iVar = 1, size(dataStruct%gru(iGRU)%hru(iHRU)%var(:)) - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var(iVar)%tim) - end do ! var - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var) - end do ! hru - deallocate(dataStruct%gru(iGRU)%hru) - end do ! gru - deallocate(dataStruct%gru) - - class is (gru_hru_int) - do iGRU = 1, size(dataStruct%gru(:)) - do iHRU = 1, size(dataStruct%gru(iGRU)%hru(:)) - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var) - end do ! hru - deallocate(dataStruct%gru(iGRU)%hru) - end do ! gru - deallocate(dataStruct%gru) - - class is (gru_hru_time_int8) - do iGRU = 1, size(dataStruct%gru(:)) - do iHRU = 1, size(dataStruct%gru(iGRU)%hru(:)) - do iVar = 1, size(dataStruct%gru(iGRU)%hru(iHRU)%var(:)) - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var(iVar)%tim) - end do ! var - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var) - end do ! hru - deallocate(dataStruct%gru(iGRU)%hru) - end do ! gru - deallocate(dataStruct%gru) - - class is (gru_hru_time_double) - do iGRU = 1, size(dataStruct%gru(:)) - do iHRU = 1, size(dataStruct%gru(iGRU)%hru(:)) - do iVar = 1, size(dataStruct%gru(iGRU)%hru(iHRU)%var(:)) - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var(iVar)%tim) - end do ! var - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var) - end do ! hru - deallocate(dataStruct%gru(iGRU)%hru) - end do ! gru - deallocate(dataStruct%gru) - - class is (gru_hru_double) - do iGRU = 1, size(dataStruct%gru(:)) - do iHRU = 1, size(dataStruct%gru(iGRU)%hru(:)) - deallocate(dataStruct%gru(iGRU)%hru(iHRU)%var) - end do ! hru - deallocate(dataStruct%gru(iGRU)%hru) - end do ! gru - deallocate(dataStruct%gru) - - class is (gru_double) - do iGRU = 1, size(dataStruct%gru(:)) - deallocate(dataStruct%gru(iGRU)%var) - end do ! gru - deallocate(dataStruct%gru) - - - end select - -end subroutine - -end module \ 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 deleted file mode 100644 index 9bbaa9316140a46a212a3b0cbe5b77054d2f6fa0..0000000000000000000000000000000000000000 --- a/build/source/actors/file_access_actor/fortran_code/initOutputStruc.f90 +++ /dev/null @@ -1,105 +0,0 @@ -module summaActors_initOutputStruct - USE nrtype - implicit none - public::initalizeOutput - contains - -subroutine initalizeOutput(forcFileInfo, maxSteps, num_gru, err) - USE globalData,only:outputStructure - USE globalData,only:time_meta,forc_meta,attr_meta,type_meta ! metadata structures - USE globalData,only:prog_meta,diag_meta,flux_meta,id_meta ! metadata structures - USE globalData,only:mpar_meta,indx_meta ! metadata structures - USE globalData,only:bpar_meta,bvar_meta ! metadata structures - USE globalData,only:statForc_meta ! child metadata for stats - USE globalData,only:statProg_meta ! child metadata for stats - USE globalData,only:statDiag_meta ! child metadata for stats - USE globalData,only:statFlux_meta ! child metadata for stats - USE globalData,only:statIndx_meta ! child metadata for stats - USE globalData,only:statBvar_meta ! child metadata for stats - USE globalData,only:gru_struc - USE globalData,only:structInfo ! information on the data structures - USE alloc_file_access,only:alloc_outputStruc - USE multiconst,only:secprday ! number of seconds in a day - USE data_types,only:file_info_array - USE var_lookup,only:maxvarFreq ! maximum number of output files - - implicit none - type(file_info_array), pointer :: forcFileInfo - integer(i4b), intent(in) :: maxSteps - integer(i4b), intent(in) :: num_gru - integer(i4b), intent(inout) :: err - - ! local variables - integer(i4b) :: nVars - integer(i4b) :: iGRU - integer(i4b) :: iHRU - integer(i4b) :: iStep - integer(i4b) :: nSnow - integer(i4b) :: nSoil - integer(i4b) :: iStruct - character(len=256) :: message - integer(i4b) :: num_hru - - ! Allocate structure to hold output files - if (.not.allocated(outputStructure))then - allocate(outputStructure(1)) - end if - - - ! Primary Data Structures (scalars) - allocate(outputStructure(1)%mparStruct(1)) - allocate(outputStructure(1)%bparStruct(1)) - allocate(outputStructure(1)%dparStruct(1)) - allocate(outputStructure(1)%mparStruct(1)%gru(num_gru)) - allocate(outputStructure(1)%bparStruct(1)%gru(num_gru)) - allocate(outputStructure(1)%dparStruct(1)%gru(num_gru)) - - ! Finalize Stats for writing - - do iGRU = 1, num_gru - num_hru = gru_struc(iGRU)%hruCount - - ! Primary Data Structures (scalars) - allocate(outputStructure(1)%mparStruct(1)%gru(iGRU)%hru(num_hru)) - allocate(outputStructure(1)%dparStruct(1)%gru(iGRU)%hru(num_hru)) - - end do - - do iGRU=1,num_gru - do iHRU=1,gru_struc(iGRU)%hruCount - - ! Get the maximum number of steps needed to initalize the output structure - nVars = maxval(forcFileInfo%ffile_list(:)%nVars) - nSnow = gru_struc(iGRU)%hruInfo(iHRU)%nSnow - nSoil = gru_struc(iGRU)%hruInfo(iHRU)%nSoil - - call alloc_outputStruc(attr_meta,outputStructure(1)%attrStruct(1)%gru(iGRU)%hru(iHRU), & - maxSteps,nSnow,nSoil,err,message); ! local attributes for each HRU - - call alloc_outputStruc(type_meta,outputStructure(1)%typeStruct(1)%gru(iGRU)%hru(iHRU), & - maxSteps,nSnow,nSoil,err,message); ! classification of soil veg etc. - - call alloc_outputStruc(id_meta,outputStructure(1)%idStruct(1)%gru(iGRU)%hru(iHRU), & - maxSteps,nSnow,nSoil,err,message); ! local values of hru gru IDs - - call alloc_outputStruc(mpar_meta,outputStructure(1)%mparStruct(1)%gru(iGRU)%hru(iHRU), & - maxSteps,nSnow,nSoil,err,message); - - call alloc_outputStruc(mpar_meta, outputStructure(1)%dparStruct(1)%gru(iGRU)%hru(iHRU), & - maxSteps,err=err,message=message) - - call alloc_outputStruc(bpar_meta,outputStructure(1)%bparStruct(1)%gru(iGRU), & - maxSteps,nSnow=0,nSoil=0,err=err,message=message); ! basin-average params - ! check errors - if(err/=0)then - message=trim(message)//'initOutputStruc.f90 - [structure = '//trim(structInfo(iStruct)%structName)//']' - print*, "message" - return - endif - end do ! Looping through GRUs - end do - - -end subroutine initalizeOutput - -end module \ No newline at end of file diff --git a/build/source/actors/file_access_actor/fortran_code/read_attribute_all_hru.f90 b/build/source/actors/file_access_actor/fortran_code/read_attribute_all_hru.f90 deleted file mode 100644 index f0912a65236817dddd59c5ee22b09ecc97d2150a..0000000000000000000000000000000000000000 --- a/build/source/actors/file_access_actor/fortran_code/read_attribute_all_hru.f90 +++ /dev/null @@ -1,291 +0,0 @@ -module read_attribute_all_hru - USE, intrinsic :: iso_c_binding - USE nrtype - implicit none - private - public::read_attribute_file_access_actor -contains -subroutine read_attribute_file_access_actor(num_gru,err) bind(C, name="readAttributeFileAccessActor") - USE globalData,only:outputStructure ! Using the output structure as global input for the attribute data. This is so we can hrus can setup params in parallel. - 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 - ! provide access to derived data types - USE data_types,only:var_d ! x%var(:) (i4b) - USE data_types,only:var_i ! x%var(:) integer(8) - USE data_types,only:var_i8 ! x%var(:) (dp) - ! provide access to global data - USE globalData,only:gru_struc ! gru-hru mapping structure - USE globalData,only:attr_meta,type_meta,id_meta ! metadata structures - USE get_ixname_module,only:get_ixAttr,get_ixType,get_ixId ! access function to find index of elements in structure - ! 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 - - - implicit none - - integer(c_int),intent(in) :: num_gru ! id of the HRU - integer(c_int),intent(out) :: err ! error code - - ! Local Variables - character(len=256) :: message ! error message - character(len=256) :: cmessage ! error message for downwind routine - integer(i4b) :: iVar ! loop through varibles in the netcdf file - integer(i4b) :: varType ! type of variable (categorica, numerical, idrelated) - integer(i4b) :: varIndx ! index of variable within its data structure - - ! check structures - integer(i4b) :: iCheck ! index of an attribute name - logical(lgt),allocatable :: checkType(:) ! vector to check if we have all desired categorical values - logical(lgt),allocatable :: checkId(:) ! vector to check if we have all desired IDs - logical(lgt),allocatable :: checkAttr(:) ! vector to check if we have all desired local attributes - - ! netcdf variables - integer(i4b) :: ncID ! netcdf file id - character(LEN=nf90_max_name) :: varName ! character array of netcdf variable name - integer(i4b) :: nVar ! number of variables in netcdf local attribute file - integer(i4b),parameter :: categorical=101 ! named variable to denote categorical data - integer(i4b),parameter :: numerical=102 ! named variable to denote numerical data - integer(i4b),parameter :: idrelated=103 ! named variable to denote ID related data - integer(i4b) :: categorical_var(1) ! temporary categorical variable from local attributes netcdf file - real(dp) :: numeric_var(1) ! temporary numeric variable from local attributes netcdf file - integer(8) :: idrelated_var(1) ! temporary ID related variable from local attributes netcdf file - - integer(i4b) :: iGRU - integer(i4b) :: iHRU - ! attribute file - character(len=256) :: attrFile ! attributes file name - - - ! define mapping variables - - ! Start procedure here - err=0; message="read_attriute_all_hru " - - attrFile = trim(SETTINGS_PATH)//trim(LOCAL_ATTRIBUTES) - - - ! ********************************************************************************************** - ! (1) prepare check vectors - ! ********************************************************************************************** - allocate(checkType(size(type_meta)),checkAttr(size(attr_meta)),checkId(size(id_meta)),stat=err) - if(err/=0)then - err=20 - message=trim(message)//'problem allocating space for variable check vectors' - print*, message - return - endif - checkType(:) = .false. - checkAttr(:) = .false. - checkId(:) = .false. - - ! ********************************************************************************************** - ! (2) open netcdf file - ! ********************************************************************************************** - ! open file - call nc_file_open(trim(attrFile),nf90_noWrite,ncID,err,cmessage) - if(err/=0)then - message=trim(message)//trim(cmessage) - print*, message - return - endif - - ! get number of variables total in netcdf file - err = nf90_inquire(ncID,nvariables=nVar) - call netcdf_err(err,message) - if (err/=0) then - message=trim(message)//'problem with nf90_inquire' - return - endif - ! ********************************************************************************************** - ! (3) read local attributes - ! ********************************************************************************************** - ! loop through variables in netcdf file and pull out local attributes - iCheck = 1 - do iVar = 1,nVar - - ! inqure about current variable name, type, number of dimensions - err = nf90_inquire_variable(ncID,iVar,name=varName) - if(err/=nf90_noerr)then; - message=trim(message)//'problem inquiring variable: '//trim(varName)//'/'//trim(nf90_strerror(err)); - print*, message - return - endif - - ! find attribute name - select case(trim(varName)) - - ! ** categorical data - case('vegTypeIndex','soilTypeIndex','slopeTypeIndex','downHRUindex') - - ! get the index of the variable - varType = categorical - varIndx = get_ixType(varName) - checkType(varIndx) = .true. - - ! check that the variable could be identified in the data structure - if(varIndx < 1)then - err=20; - message=trim(message)//'unable to find variable ['//trim(varName)//'] in data structure'; - print*, message - return; - endif - - - do iGRU=1,num_gru - do iHRU = 1,gru_struc(iGRU)%hruCount - err = nf90_get_var(ncID,iVar,categorical_var,start=(/gru_struc(iGRU)%hruInfo(iHRU)%hru_nc/),count=(/1/)) - if(err/=nf90_noerr)then - message=trim(message)//'problem reading: '//trim(varName) - print*, message - return - end if - outputStructure(1)%typeStruct(1)%gru(iGRU)%hru(iHRU)%var(varIndx) = categorical_var(1) - end do - end do - - ! ** ID related data - case('hruId') - ! get the index of the variable - varType = idrelated - varIndx = get_ixId(varName) - checkId(varIndx) = .true. - - ! check that the variable could be identified in the data structure - if(varIndx < 1)then - err=20 - message=trim(message)//'unable to find variable ['//trim(varName)//'] in data structure' - print*, message - return - endif - - ! get data from netcdf file and store in vector - do iGRU=1,num_gru - do iHRU = 1,gru_struc(iGRU)%hruCount - err = nf90_get_var(ncID,iVar,idrelated_var,start=(/gru_struc(iGRU)%hruInfo(iHRU)%hru_nc/),count=(/1/)) - if(err/=nf90_noerr)then - message=trim(message)//'problem reading: '//trim(varName) - print*, message - return - end if - outputStructure(1)%idStruct(1)%gru(iGRU)%hru(iHRU)%var(varIndx) = idrelated_var(1) - end do - end do - - ! ** numerical data - case('latitude','longitude','elevation','tan_slope','contourLength','HRUarea','mHeight') - - ! get the index of the variable - varType = numerical - varIndx = get_ixAttr(varName) - checkAttr(varIndx) = .true. - - ! check that the variable could be identified in the data structure - if(varIndx < 1)then - err=20; message=trim(message)//'unable to find variable ['//trim(varName)//'] in data structure' - print*, message - return - endif - ! get data from netcdf file and store in vector - - do iGRU=1,num_gru - do iHRU = 1,gru_struc(iGRU)%hruCount - err = nf90_get_var(ncID,iVar,numeric_var,start=(/gru_struc(iGRU)%hruInfo(iHRU)%hru_nc/),count=(/1/)) - if(err/=nf90_noerr)then - message=trim(message)//'problem reading: '//trim(varName) - print*, message - return - end if - outputStructure(1)%attrStruct(1)%gru(iGRU)%hru(iHRU)%var(varIndx) = numeric_var(1) - end do - end do - - ! for mapping varibles, do nothing (information read above) - case('hru2gruId','gruId'); cycle - - ! check that variables are what we expect - case default - message=trim(message)//'unknown variable ['//trim(varName)//'] in local attributes file' - print*,message - err=20 - return - - end select ! select variable - - end do ! (looping through netcdf local attribute file) - - ! ** now handle the optional aspect variable if it's missing - varIndx = get_ixAttr('aspect') - ! check that the variable was not found in the attribute file - if(.not. checkAttr(varIndx)) then - write(*,*) NEW_LINE('A')//'INFO: aspect not found in the input attribute file, continuing ...'//NEW_LINE('A') - do iGRU=1,num_gru - do iHRU = 1, gru_struc(iGRU)%hruCount - outputStructure(1)%attrStruct(1)%gru(iGRU)%hru(iHRU)%var(varIndx) = nr_realMissing ! populate variable with out-of-range value, used later - end do - end do - checkAttr(varIndx) = .true. - endif - - ! TODO: find out why this is here, probably for the lateral flows - varIndx = get_ixTYPE('downkHRU') - checkType(varIndx) = .true. - ! outputStructure(1)%typeStruct(1)%gru(iGRU)%hru(iHRU)%var(varIndx) = 0 - - ! ********************************************************************************************** - ! (4) check that we have all the desired varaibles - ! ********************************************************************************************** - ! check that we have all desired categorical variables - if(any(.not.checkType))then - do iCheck = 1,size(type_meta) - if(.not.checkType(iCheck))then - err=20 - message=trim(message)//'missing variable ['//trim(type_meta(iCheck)%varname)//'] in local attributes file' - print*, message - return - endif - end do - endif - - ! check that we have all desired ID variables - if(any(.not.checkId))then - do iCheck = 1,size(id_meta) - if(.not.checkId(iCheck))then - err=20 - message=trim(message)//'missing variable ['//trim(id_meta(iCheck)%varname)//'] in local attributes file' - print*, message - return - endif - end do - endif - - ! check that we have all desired local attributes - if(any(.not.checkAttr))then - do iCheck = 1,size(attr_meta) - if(.not.checkAttr(iCheck))then; - err=20 - message=trim(message)//'missing variable ['//trim(attr_meta(iCheck)%varname)//'] in local attributes file' - print*, message - return - endif - end do - endif - - ! ********************************************************************************************** - ! (5) close netcdf file - ! ********************************************************************************************** - - call nc_file_close(ncID,err,cmessage) - if (err/=0)then; message=trim(message)//trim(cmessage); return; end if - - ! free memory - deallocate(checkType) - deallocate(checkId) - deallocate(checkAttr) - -end subroutine - - -end module \ No newline at end of file diff --git a/build/source/actors/file_access_actor/fortran_code/read_forcing.f90 b/build/source/actors/file_access_actor/fortran_code/read_forcing.f90 new file mode 100644 index 0000000000000000000000000000000000000000..6c173fd5b17c07f1de5546bdee410fb8c956c2ae --- /dev/null +++ b/build/source/actors/file_access_actor/fortran_code/read_forcing.f90 @@ -0,0 +1,10 @@ +! Module used by the file access actor for reading in the forcing data +module read_forcing_module + USE, intrinsic :: iso_c_binding + USE nrtype + + implicit none + + + +end module read_forcing_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 deleted file mode 100644 index 56d835313013c17070d481232cfe5beec4d7ebbf..0000000000000000000000000000000000000000 --- a/build/source/actors/file_access_actor/fortran_code/read_param_all_hru.f90 +++ /dev/null @@ -1,323 +0,0 @@ -module read_param_all_hru - USE, intrinsic :: iso_c_binding - USE nrtype - implicit none - private - public::read_param_file_access_actor - public::overwriteParam -contains - - -! overwriteParm sets the basin parameters -! It can then overwrite them from information from the Noah-MP tables -subroutine overwriteParam(num_gru, err) bind(C, name="overwriteParam") - USE globalData,only:outputStructure - USE pOverwrite_module,only:pOverwrite ! module to overwrite default parameter values with info from the Noah tables - USE globalData,only:gru_struc - USE globalData,only:localParFallback ! local column default parameters - USE globalData,only:basinParFallback ! basin-average default parameter - USE var_lookup,only:iLookTYPE ! look-up values for classification of veg, soils etc. - - implicit none - integer(c_int),intent(in) :: num_gru ! number of GRUs in the run_domain - integer(c_int),intent(out) :: err ! error code - - ! local - integer(i4b) :: iGRU - integer(i4b) :: iHRU - integer(i4b) :: iVar - integer(i4b) :: iDat - character(len=256) :: message - - err=0; message="overwriteParam" - - ! Need to set the basin parameters with the default values for when we copy - do iGRU=1,num_gru - do iHRU=1,gru_struc(iGRU)%hruCount - do iVar=1, size(localParFallback) - outputStructure(1)%dparStruct(1)%gru(iGRU)%hru(iHRU)%var(iVar) = localParFallback(iVar)%default_val - end do - ! overwrite default model parameters with information from the Noah-MP tables - call pOverwrite(outputStructure(1)%typeStruct(1)%gru(iGRU)%hru(iHRU)%var(iLookTYPE%vegTypeIndex), & ! vegetation category - outputStructure(1)%typeStruct(1)%gru(iGRU)%hru(iHRU)%var(iLookTYPE%soilTypeIndex), & ! soil category - outputStructure(1)%dparStruct(1)%gru(iGRU)%hru(iHRU)%var(:), & ! default model parameters - err,message) - - do iVar=1, size(localParFallback) - do iDat=1, size(outputStructure(1)%mparStruct(1)%gru(iGRU)%hru(iHRU)%var(iVar)%dat) - outputStructure(1)%mparStruct(1)%gru(iGRU)%hru(iHRU)%var(iVar)%dat(iDat) = outputStructure(1)%dparStruct(1)%gru(iGRU)%hru(iHRU)%var(iVar) - end do - end do - end do - do iVar=1,size(basinParFallback) - outputStructure(1)%bparStruct(1)%gru(iGRU)%var(iVar) = basinParFallback(iVar)%default_val - end do - end do - -end subroutine -subroutine read_param_file_access_actor(startGRU,num_gru,err) bind(C, name="readParamFileAccessActor") - ! used to read model initial conditions - USE summaActors_FileManager,only:SETTINGS_PATH ! path for metadata files - USE summaActors_FileManager,only:PARAMETER_TRIAL ! file with parameter trial values - 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 var_lookup,only:iLookPARAM,iLookTYPE,iLookID ! named variables to index elements of the data vectors - USE globalData,only:integerMissing ! missing integer - USE globalData,only:realMissing ! missing real number - - USE netcdf - USE netcdf_util_module,only:nc_file_close ! close netcdf file - USE netcdf_util_module,only:nc_file_open ! open netcdf file - USE netcdf_util_module,only:netcdf_err ! netcdf error handling function - - USE globalData,only:outputStructure - - implicit none - ! define input - - integer(c_int),intent(in) :: startGRU ! starting Index of gru Batch - integer(c_int),intent(in) :: num_gru ! number of GRUs in the run_domain - integer(c_int),intent(out) :: err ! error code - - ! define local variables - character(len=256) :: message ! error message - character(len=1024) :: cmessage ! error message for downwind routine - character(LEN=1024) :: infile ! input filename - 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) :: ncid ! netcdf id - integer(i4b) :: nDims ! number of dimensions - integer(i4b) :: nVars ! number of variables - integer(i4b) :: idimid ! dimension index - integer(i4b) :: ivarid ! variable index - character(LEN=64) :: dimName ! dimension name - character(LEN=64) :: parName ! parameter name - integer(i4b) :: dimLength ! dimension length - integer(i4b) :: nHRU_file ! number of HRUs in the parafile - integer(i4b) :: nGRU_file ! number of GRUs in the parafile - 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 - integer(8),allocatable :: hruId(:) ! HRU identifier in the file - real(dp),allocatable :: parVector(:) ! model parameter vector - logical :: fexist ! inquire whether the parmTrial file exists - integer(i4b) :: fHRU ! index of HRU in input file - integer(i4b) :: iGRU - integer(i4b) :: iHRU - - err=0; message="read_param_all_hru.f90/" - ! ********************************************************************************************** - ! * 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=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 - - ! ! 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 - nGRU_file=integerMissing - - ! get the length of the dimensions - do idimid=1,nDims - ! get the dimension name and length - err=nf90_inquire_dimension(ncid, idimid, name=dimName, len=dimLength) - if(err/=0)then; message=trim(message)//trim(cmessage); return; end if - ! get the number of HRUs - if(trim(dimName)=='hru') nHRU_file=dimLength - if(trim(dimName)=='gru') nGRU_file=dimLength - end do - - ! allocate hruID vector - allocate(hruId(nHRU_file)) - - ! check HRU dimension exists - if(nHRU_file==integerMissing)then - message=trim(message)//'unable to identify HRU dimension in file '//trim(infile) - err=20; return - endif - ! ********************************************************************************************** - ! * read the HRU index - ! ********************************************************************************************** - ! loop through the parameters in the NetCDF file - do ivarid=1,nVars - - ! get the parameter name - err=nf90_inquire_variable(ncid, ivarid, name=parName) - call netcdf_err(err,message) - if (err/=0) then - err=20 - print*, message - return - end if - - ! special case of the HRU id - if(trim(parName)=='hruIndex' .or. trim(parName)=='hruId')then - - ! read HRUs - err=nf90_get_var(ncid, ivarid, hruId) - if(err/=0)then - message=trim(message)//trim(cmessage) - print*, message - return - end if - - endif ! if the HRU id - - end do ! looping through variables in the file - - ! ********************************************************************************************** - ! * read the local parameters and the basin parameters - ! ********************************************************************************************** - do ivarid=1,nVars - ! get the parameter name - err=nf90_inquire_variable(ncid, ivarid, name=parName) - call netcdf_err(err,message); if (err/=0) then; err=20; 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(ncid, ivarid, nDims=nDims, 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(nDims==2)then - ! get the information on the 2nd dimension for 2-d variables - err=nf90_inquire_dimension(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 - - ! TODO: implement this check - ! ! check that the dimension length is correct - ! if(size(outputStructur(1)%mparStruct%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 - - ! loop through all hrus - do iGRU=1, num_gru - do iHRU=1, gru_struc(iGRU)%hruCount - localHRU_ix=index_map(iHRU)%localHRU_ix - fHRU = gru_struc(iGRU)%hruInfo(localHRU_ix)%hru_nc - ! read parameter data - select case(nDims) - case(1); err=nf90_get_var(ncid, ivarid, parVector, start=(/fHRU/), count=(/1/) ) - case(2); err=nf90_get_var(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(nDims) - case(1); outputStructure(1)%mparStruct(1)%gru(iGRU)%hru(localHRU_ix)%var(ixParam)%dat(:) = parVector(1) ! also distributes scalar across depth dimension - case(2); outputStructure(1)%mparStruct(1)%gru(iGRU)%hru(localHRU_ix)%var(ixParam)%dat(:) = parVector(:) - case default; err=20; message=trim(message)//'unexpected number of dimensions for parameter '//trim(parName) - end select - end do - end do - - ! 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 - - ! allocate space for model parameters - allocate(parVector(nGRU_file),stat=err) - if(err/=0)then - message=trim(message)//'problem allocating space for parameter vector' - print*, message - err=20; return - endif - - ! read parameter data - err=nf90_get_var(ncid, ivarid, parVector ) - if(err/=0)then; message=trim(message)//trim(cmessage); return; end if - - ! populate parameter structures - do iGRU=1, num_gru - outputStructure(1)%bparStruct(1)%gru(iGRU)%var(ixParam) = parVector(iGRU+startGRU-1) - end do - - ! 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 - endif - - end do ! (looping through the parameters in the NetCDF file) - - ! Now we must close the netcdf file - call nc_file_close(ncid,err,message) - if(err/=0)then;message=trim(message)//trim(cmessage);return;end if -end subroutine read_param_file_access_actor -end module read_param_all_hru \ No newline at end of file diff --git a/build/source/actors/hru_actor/cpp_code/hru_actor.cpp b/build/source/actors/hru_actor/cpp_code/hru_actor.cpp index cc3bb61708d1cbe2af15596aad9aac34b363041e..42d06d1d4a06113d53e9ad85d51af811e6e6c096 100644 --- a/build/source/actors/hru_actor/cpp_code/hru_actor.cpp +++ b/build/source/actors/hru_actor/cpp_code/hru_actor.cpp @@ -67,14 +67,14 @@ behavior hru_actor(stateful_actor<hru_state>* self, int refGRU, int indxGRU, // Get attributes - self->send(self->state.file_access_actor, get_attributes_v, self->state.refGRU, self); + self->send(self->state.file_access_actor, get_attributes_params_v, self->state.refGRU, self); self->state.hru_timing.updateEndPoint("total_duration"); return { // Starts the HRU and tells it to ask for data from the file_access_actor - [=](get_attributes, std::vector<double> attr_struct, std::vector<int> type_struct, + [=](get_attributes_params, std::vector<double> attr_struct, std::vector<int> type_struct, std::vector<long int> id_struct, std::vector<double> bpar_struct, std::vector<double> dpar_struct, std::vector<std::vector<double>> mpar_struct) { diff --git a/build/source/engine/ffile_info.f90 b/build/source/engine/ffile_info.f90 index 10935780b5651ef654e8b713d4a8bda99fc4243b..3e2bebabe0d54e9720fccd2950b538a550ddd559 100755 --- a/build/source/engine/ffile_info.f90 +++ b/build/source/engine/ffile_info.f90 @@ -19,6 +19,7 @@ ! along with this program. If not, see <http://www.gnu.org/licenses/>. module ffile_info_module +USE, intrinsic :: iso_c_binding USE nrtype USE netcdf USE data_types @@ -26,12 +27,18 @@ USE globalData,only:integerMissing implicit none private public::ffile_info +private::getAndSetVars +private::popFileNames +private::popForcStructVars +private::getAndCheckDataStep contains + + ! ************************************************************************************************ ! public subroutine ffile_info: read information on model forcing files ! ************************************************************************************************ -subroutine ffile_info(indxGRU,forcFileInfo,numFiles,err,message) +subroutine ffile_info(indxGRU,handle_forcFileInfo,num_forcing_files,err) bind(C, name='ffile_info') ! used to read metadata on the forcing data file USE ascii_util_module,only:file_open USE ascii_util_module,only:linewidth @@ -40,43 +47,39 @@ subroutine ffile_info(indxGRU,forcFileInfo,numFiles,err,message) USE summaActors_FileManager,only:SETTINGS_PATH ! path for metadata files USE summaActors_FileManager,only:FORCING_PATH ! path for forcing files USE summaActors_FileManager,only:FORCING_FILELIST ! list of model forcing files - USE summaActors_FileManager,only:FORCING_START - USE globalData,only:data_step USE globalData,only:forc_meta ! forcing metadata USE get_ixname_module,only:get_ixtime,get_ixforce ! identify index of named variable USE ascii_util_module,only:get_vlines ! get a vector of non-comment lines USE ascii_util_module,only:split_line ! split a line into words - USE globalData,only:gru_struc ! gru-hru mapping structure USE time_utils_module,only:extractTime USE globalData,only:time_meta USE allocspace_module,only:allocLocal USE time_utils_module,only:extractTime ! extract time info from units string - USE summaActors_FileManager,only: SIM_START_TM, SIM_END_TM, FORCING_START ! time info from control file module + USE summaActors_FileManager,only: SIM_START_TM, SIM_END_TM! time info from control file module USE var_lookup,only:iLookTIME ! named variables that identify indices in the time structures implicit none ! define input & output - integer(i4b),intent(in) :: indxGRU - type(file_info_array),intent(inout) :: forcFileInfo - integer(i4b),intent(out) :: numFiles - integer(i4b),intent(out) :: err ! error code - character(*),intent(out) :: message ! error message + integer(c_int),intent(in) :: indxGRU + type(c_ptr), intent(in), value :: handle_forcFileInfo + integer(c_int),intent(out) :: num_forcing_files + integer(c_int),intent(out) :: err ! error code ! define local variables + type(file_info_array),pointer :: forcFileInfo ! netcdf file i/o related integer(i4b) :: ncid ! netcdf file id integer(i4b) :: mode ! netCDF file open mode - integer(i4b) :: varid ! netcdf variable id integer(i4b) :: dimId ! netcdf dimension id character(LEN=nf90_max_name) :: varName ! character array of netcdf variable name integer(i4b) :: iNC ! index of a variable in netcdf file integer(i4b) :: nvar ! number of variables in netcdf local attribute file ! the rest character(LEN=linewidth),allocatable :: dataLines(:) ! vector of lines of information (non-comment lines) + character(len=256) :: message ! error message for downwind routine character(len=256) :: cmessage ! error message for downwind routine character(LEN=256) :: infile ! input filename integer(i4b) :: unt ! file unit (free unit output from file_open) - character(LEN=256) :: filenameData ! name of forcing datafile integer(i4b) :: ivar ! index of model variable integer(i4b) :: iFile ! counter for forcing files integer(i4b) :: nFile ! number of forcing files in forcing file list @@ -84,8 +87,6 @@ subroutine ffile_info(indxGRU,forcFileInfo,numFiles,err,message) integer(i4b) :: startIndx ! total number of forcing files defiend in the forcing file list integer(i4b) :: file_nHRU ! number of HRUs in current forcing file integer(i4b) :: nForcing ! number of forcing variables - integer(i4b) :: localHRU_ix ! index of HRU - integer(8) :: ncHruId(1) ! hruID from the forcing files real(dp) :: dataStep_iFile ! data step for a given forcing data file logical(lgt) :: xist ! .TRUE. if the file exists ! Time Variables @@ -94,98 +95,35 @@ subroutine ffile_info(indxGRU,forcFileInfo,numFiles,err,message) type(var_i) :: finishTime real(rkind) :: dsec,dsec_tz integer(i4b) :: ffinfo_index - + + call c_f_pointer(handle_forcFileInfo, forcFileInfo) ! Start procedure here err=0; message="ffile_info/" ! ------------------------------------------------------------------------------------------------------------------ ! (1) read from the list of forcing files ! ------------------------------------------------------------------------------------------------------------------ - + ! build filename for forcing-file list file infile = trim(SETTINGS_PATH)//trim(FORCING_FILELIST) ! open file call file_open(trim(infile),unt,err,cmessage) - if(err/=0)then - message=trim(message)//trim(cmessage) - print*, message - return - end if + if(err/=0)then;message=trim(message)//trim(cmessage);print*,message;return;end if; ! get a list of character strings from non-comment lines call get_vlines(unt,dataLines,err,cmessage) - if(err/=0)then - err=20 - message=trim(message)//trim(cmessage) - print*, message - return - end if - nFile = size(dataLines) + if(err/=0)then; err=20;message=trim(message)//trim(cmessage);print*,message;return;end if; - ! Get the number of forcing files needed - call allocLocal(time_meta, startTime, err=err, message=cmessage) - call allocLocal(time_meta, forcingStart, err=err, message=cmessage) - call allocLocal(time_meta, finishTime, err=err, message=cmessage) - call extractTime(trim(SIM_START_TM), & ! date-time string - startTime%var(iLookTIME%iyyy), & ! year - startTime%var(iLookTIME%im), & ! month - startTime%var(iLookTIME%id), & ! day - startTime%var(iLookTIME%ih), & ! hour - startTime%var(iLookTIME%imin), & ! minute - dsec, & ! second - startTime%var(iLookTIME%ih_tz), & ! time zone hour - startTime%var(iLookTIME%imin_tz), & ! time zone minnute - dsec_tz, & ! time zone seconds - err,cmessage) - if(err/=0)then; err=20; message=trim(message)//trim(cmessage); return; end if - call extractTime(trim(FORCING_START), & ! date-time string - forcingStart%var(iLookTIME%iyyy), & ! year - forcingStart%var(iLookTIME%im), & ! month - forcingStart%var(iLookTIME%id), & ! day - forcingStart%var(iLookTIME%ih), & ! hour - forcingStart%var(iLookTIME%imin), & ! minute - dsec, & ! second - forcingStart%var(iLookTIME%ih_tz), & ! time zone hour - forcingStart%var(iLookTIME%imin_tz), & ! time zone minnute - dsec_tz, & ! time zone seconds - err,cmessage) - if(err/=0)then; err=20; message=trim(message)//trim(cmessage); return; end if - ! put simulation end time information into the time structures - call extractTime(trim(SIM_END_TM), & ! date-time string - finishTime%var(iLookTIME%iyyy), & ! year - finishTime%var(iLookTIME%im), & ! month - finishTime%var(iLookTIME%id), & ! day - finishTime%var(iLookTIME%ih), & ! hour - finishTime%var(iLookTIME%imin), & ! minute - dsec, & ! second - finishTime%var(iLookTIME%ih_tz), & ! time zone hour - finishTime%var(iLookTIME%imin_tz), & ! time zone minnute - dsec_tz, & ! time zone seconds - err,cmessage) ! error control - if(err/=0)then; err=20; message=trim(message)//trim(cmessage); return; end if - - startIndx = (((startTime%var(iLookTIME%iyyy) - forcingStart%var(iLookTIME%iyyy)) * 12) & - + startTime%var(iLookTIME%im) - forcingStart%var(iLookTIME%im) + 1) - - totalFiles = (((finishTime%var(iLookTIME%iyyy) - startTime%var(iLookTIME%iyyy)) * 12) & - + finishTime%var(iLookTIME%im) - startTime%var(iLookTIME%im) + 1) + nFile = size(dataLines) ! allocate space for forcing information if(allocated(forcFileInfo%ffile_list)) deallocate(forcFileInfo%ffile_list) - allocate(forcFileInfo%ffile_list(totalFiles), stat=err) - if(err/=0)then; err=20; message=trim(message)//'problem allocating space for forcFileInfo'; return; end if + allocate(forcFileInfo%ffile_list(nFile), stat=err) + if(err/=0)then; err=20; message=trim(message)//'problem allocating space for forcFileInfo';print*,message;return;endif - ffinfo_index = 1 ! poputate the forcingInfo structure with filenames - do iFile=startIndx,nFile - if (ffinfo_index > totalFiles)then; exit; endif; - ! split the line into "words" (expect one word: the file describing forcing data for that index) - read(dataLines(iFile),*,iostat=err) filenameData - if(err/=0)then; message=trim(message)//'problem reading a line of data from file ['//trim(infile)//']'; return; end if - ! set forcing file name attribute - forcFileInfo%ffile_list(ffinfo_index)%filenmData = trim(filenameData) - ffinfo_index = ffinfo_index + 1 - end do ! (looping through files) + call popFileNames(forcFileInfo,dataLines,nFile,inFile,err,message) + if(err/=0)then; err=20; message=trim(message);print*,message;return;endif ! close ascii file close(unit=unt,iostat=err); if(err/=0)then;message=trim(message)//'problem closing forcing file list'; return; end if @@ -198,7 +136,7 @@ subroutine ffile_info(indxGRU,forcFileInfo,numFiles,err,message) nForcing = size(forc_meta) ! loop through files, and read descriptive information from each file - do iFile=1,totalFiles + do iFile=1,nFile ! ensure allocatable structure components are deallocated if(allocated(forcFileInfo%ffile_list(iFile)%data_id)) deallocate(forcFileInfo%ffile_list(iFile)%data_id) @@ -206,7 +144,7 @@ subroutine ffile_info(indxGRU,forcFileInfo,numFiles,err,message) ! allocate space for structure components allocate(forcFileInfo%ffile_list(iFile)%data_id(nForcing), forcFileInfo%ffile_list(iFile)%varName(nForcing), stat=err) - if(err/=0)then; err=41; message=trim(message)//"problemAllocateStructureElement"; return; end if + if(err/=0)then; err=41; message=trim(message)//"problemAllocateStructureElement"; return;endif; ! initialize variable ids to missing forcFileInfo%ffile_list(iFile)%data_id(:) = integerMissing @@ -215,40 +153,16 @@ subroutine ffile_info(indxGRU,forcFileInfo,numFiles,err,message) infile = trim(FORCING_PATH)//trim(forcFileInfo%ffile_list(iFile)%filenmData) ! check if file exists inquire(file=trim(infile),exist=xist) - if(.not.xist)then - message=trim(message)//"FileNotFound[file='"//trim(infile)//"']" - print*, message - err=10; return - end if + if(.not.xist)then;message=trim(message)//"FileNotFound[file='"//trim(infile)//"']";print*,message;err=10;return;endif; ! open file mode=nf90_NoWrite call nc_file_open(trim(infile), mode, ncid, err, cmessage) - if(err/=0)then - message=trim(message)//trim(cmessage) - print*, message - return - end if - - ! how many variables are there? - err = nf90_inquire(ncid, nvariables=nVar) - call netcdf_err(err,message); if (err/=0) return - - ! set nVar attribute - forcFileInfo%ffile_list(iFile)%nVars = nVar - if(allocated(forcFileInfo%ffile_list(iFile)%var_ix))then - print*, "Already Allocated" - endif - ! allocate space - allocate(forcFileInfo%ffile_list(iFile)%var_ix(nVar), stat=err) - if(err/=0)then - message=trim(message)//'problem allocating space for structure element' - print*, message - err=20; return - endif + if(err/=0)then;message=trim(message)//trim(cmessage);print*,message;return;endif; - ! initialize data structure - forcFileInfo%ffile_list(iFile)%var_ix(:) = integerMissing + ! check how many variables are in the netCDF file + ! populate forcFileInfo with variable information + call getAndSetVars(ncid,forcFileInfo,iFile,nVar,err,cmessage) ! inquire nhru dimension size err = nf90_inq_dimid(ncid,'hru',dimId); if(err/=0)then; message=trim(message)//'cannot find dimension hru'; return; endif @@ -270,64 +184,21 @@ subroutine ffile_info(indxGRU,forcFileInfo,numFiles,err,message) ! if variable is in the forcing vector case('time','pptrate','SWRadAtm','LWRadAtm','airtemp','windspd','airpres','spechum') + + call popForcStructVars(ncid,forcFileInfo,iFile,iNC,varName,err,message) + if(err/=0)then;print*,message;return;endif - ! get variable index - ivar = get_ixforce(trim(varname)) - if(ivar < 0)then; err=40; message=trim(message)//"variableNotFound[var="//trim(varname)//"]"; return; end if - if(ivar>size(forcFileInfo%ffile_list(iFile)%data_id))then; err=35; message=trim(message)//"indexOutOfRange[var="//trim(varname)//"]"; return; end if - - ! put netcdf file variable index in the forcing file metadata structure - err = nf90_inq_varid(ncid, trim(varName), forcFileInfo%ffile_list(iFile)%data_id(ivar)) - if(err/=0)then; message=trim(message)//"problem inquiring forcing variable[var="//trim(varName)//"]"; return; end if - - ! put variable index of the forcing structure in the metadata structure - if(trim(varName)/='time')then - forcFileInfo%ffile_list(iFile)%var_ix(iNC) = ivar - forcFileInfo%ffile_list(iFile)%varName(ivar) = trim(varName) - - ! get first time from file, place into forcFileInfo - else - err = nf90_get_var(ncid,forcFileInfo%ffile_list(iFile)%data_id(ivar),forcFileInfo%ffile_list(iFile)%firstJulDay,start=(/1/)) - if(err/=0)then; message=trim(message)//'problem reading first Julian day'; return; end if - end if ! if the variable name is time - - ! data step case('data_step' ) - ! read data_step from netcdf file - err = nf90_inq_varid(ncid, "data_step", varId); if(err/=0)then; message=trim(message)//'cannot find data_step'; return; end if - err = nf90_get_var(ncid,varid,dataStep_iFile); if(err/=0)then; message=trim(message)//'cannot read data_step'; return; end if - - ! check data_step is the same for all forcing files - if(iFile == 1)then - data_step = dataStep_iFile - else - if(abs(dataStep_iFile - data_step) > epsilon(dataStep_iFile))then - write(message,'(a,i0,a)') trim(message)//'data step for forcing file ',iFile,'differs from the datastep of the first forcing file' - err=20; return - end if - end if - + call getAndCheckDataStep(ncid,dataStep_iFile,iFile,err,message) + if(err/=0)then;print*,message;return;endif + ! HRU id -- required case('hruId') - - ! check to see if hruId exists as a variable, this is a required variable - err = nf90_inq_varid(ncid,trim(varname),varId) - if(err/=0)then; message=trim(message)//'hruID variable not present'; return; endif - - ! check that the hruId is what we expect - ! NOTE: we enforce that the HRU order in the forcing files is the same as in the zLocalAttributes files (too slow otherwise) - do localHRU_ix=1,gru_struc(indxGRU)%hruCount - ! check the HRU is what we expect - err = nf90_get_var(ncid,varId,ncHruId,start=(/gru_struc(indxGRU)%hruInfo(localHRU_ix)%hru_nc/),count=(/1/)) - if(gru_struc(indxGRU)%hruInfo(localHRU_ix)%hru_id /= ncHruId(1))then - write(message,'(a,i0,a,i0,a,i0,a,a)') trim(message)//'hruId for global HRU: ',gru_struc(indxGRU)%hruInfo(localHRU_ix)%hru_nc,' - ', & - ncHruId(1), ' differs from the expected: ',gru_struc(indxGRU)%hruInfo(localHRU_ix)%hru_id, ' in file ', trim(infile) - write(message,'(a)') trim(message)//' order of hruId in forcing file needs to match order in zLocalAttributes.nc' - err=40; return - endif - end do - + + call setHRUID(ncid,indxGRU,varName,inFile,err,message) + if(err/=0)then;print*,message;return;endif + ! OK to have additional variables in the forcing file that are not used case default; cycle end select ! select variable name @@ -351,7 +222,182 @@ subroutine ffile_info(indxGRU,forcFileInfo,numFiles,err,message) end do ! (loop through files) ! Get the number of forcing files we have - numFiles = size(forcFileInfo%ffile_list(:)) + num_forcing_files = size(forcFileInfo%ffile_list(:)) end subroutine ffile_info + + ! check how many variables are in the netCDF file + ! populate forcFileInfo with variable information +subroutine getAndSetVars(ncid,forcFileInfo,iFile,nVar,err,message) + USE netCDF + USE netcdf_util_module,only:netcdf_err + USE globalData,only:integerMissing + + implicit none + integer(i4b),intent(in) :: ncid + type(file_info_array),pointer,intent(in) :: forcFileInfo + integer(i4b),intent(in) :: iFile + integer(i4b),intent(out) :: nVar + integer(i4b),intent(out) :: err + character(*),intent(inout) :: message + + ! how many variables are there? + err = nf90_inquire(ncid, nvariables=nVar) + call netcdf_err(err,message); if (err/=0) return + + ! set nVar attribute + forcFileInfo%ffile_list(iFile)%nVars = nVar + if(allocated(forcFileInfo%ffile_list(iFile)%var_ix))then;print*,"Already Allocated";endif; + ! allocate space + allocate(forcFileInfo%ffile_list(iFile)%var_ix(nVar), stat=err) + if(err/=0)then;message=trim(message)//'problem allocating space for structure element';print*,message;err=20;return;endif; + + ! initialize data structure + forcFileInfo%ffile_list(iFile)%var_ix(:) = integerMissing + + +end subroutine getAndSetVars + +! Private subroutine that filles the array of forcing files with the names +! of the forcing files +subroutine popFileNames(forcFileInfo,dataLines,nFile,inFile,err,message) + USE ascii_util_module,only:linewidth + implicit none + + type(file_info_array),pointer :: forcFileInfo + character(*),intent(out) :: dataLines(:) ! forcingFileList + integer(i4b),intent(in) :: nFile ! number of forcing files + character(*),intent(inout) :: inFile ! file that populated dataLines + integer(i4b),intent(out) :: err + character(*),intent(inout) :: message + ! local variables + character(LEN=256) :: filenameData ! name of forcing datafile + integer(i4b) :: iFile + + do iFile=1,nFile + ! split the line into "words" (expect one word: the file describing forcing data for that index) + read(dataLines(iFile),*,iostat=err) filenameData + if(err/=0)then; message=trim(message)//'problem reading a line of data from file ['//trim(infile)//']'; return; end if + ! set forcing file name attribute + forcFileInfo%ffile_list(iFile)%filenmData = trim(filenameData) + end do ! (looping through files) + + +end subroutine popFileNames + +! populate the variable indexes and names in ffile_list +subroutine popForcStructVars(ncid,forcFileInfo,iFile,iNC,varName,err,message) + USE get_ixname_module,only:get_ixforce ! identify index of named variable + USE netCDF + implicit none + integer(i4b),intent(in) :: ncid + type(file_info_array),pointer :: forcFileInfo + integer(i4b),intent(in) :: iFile + integer(i4b),intent(in) :: iNC ! index of a variable in netcdf file + character(*),intent(in) :: varName + integer(i4b),intent(out) :: err + character(*),intent(inout) :: message + ! local variables + integer(i4b) :: iVar + + + ! get variable index + ivar = get_ixforce(trim(varname)) + if(ivar < 0)then + err=40 + message=trim(message)//"variableNotFound[var="//trim(varname)//"]" + return + end if + if(ivar>size(forcFileInfo%ffile_list(iFile)%data_id))then + err=35 + message=trim(message)//"indexOutOfRange[var="//trim(varname)//"]" + return + end if + + ! put netcdf file variable index in the forcing file metadata structure + err = nf90_inq_varid(ncid, trim(varName), forcFileInfo%ffile_list(iFile)%data_id(ivar)) + if(err/=0)then; message=trim(message)//"problem inquiring forcing variable[var="//trim(varName)//"]"; return; end if + + ! put variable index of the forcing structure in the metadata structure + if(trim(varName)/='time')then + forcFileInfo%ffile_list(iFile)%var_ix(iNC) = ivar + forcFileInfo%ffile_list(iFile)%varName(ivar) = trim(varName) + + ! get first time from file, place into forcFileInfo + else + err = nf90_get_var(ncid,forcFileInfo%ffile_list(iFile)%data_id(ivar),forcFileInfo%ffile_list(iFile)%firstJulDay,start=(/1/)) + if(err/=0)then; message=trim(message)//'problem reading first Julian day'; return; end if + end if ! if the variable name is time + +end subroutine popForcStructVars + +! set the datastep from the first file +! check every file has the same data_step +subroutine getAndCheckDataStep(ncid,dataStep_iFile,iFile,err,message) + USE globalData,only:data_step + USE netCDF + implicit none + integer(i4b),intent(in) :: ncid + real(dp),intent(out) :: dataStep_iFile ! data step for a given forcing data file + integer(i4b),intent(in) :: iFile ! index of current forcing file + integer(i4b),intent(out) :: err + character(*),intent(inout) :: message + ! local variables + integer(i4b) :: varid ! netcdf variable id + + ! read data_step from netcdf file + err = nf90_inq_varid(ncid, "data_step", varId); if(err/=0)then; message=trim(message)//'cannot find data_step'; return; end if + err = nf90_get_var(ncid,varid,dataStep_iFile); if(err/=0)then; message=trim(message)//'cannot read data_step'; return; end if + + ! check data_step is the same for all forcing files + if(iFile == 1)then + data_step = dataStep_iFile + else + if(abs(dataStep_iFile - data_step) > epsilon(dataStep_iFile))then + write(message,'(a,i0,a)') trim(message)//'data step for forcing file ',iFile,'differs from the datastep of the first forcing file' + err=20; return + end if + end if + +end subroutine getAndCheckDataStep + +! get the HRU_ID from the forcing file +! populate the gru_struc with the HRU ID +! Check if the hru is expected +subroutine setHRUID(ncid,indxGRU,varName,inFile,err,message) + USE globalData,only:gru_struc ! gru-hru mapping structure + USE netCDF + implicit none + integer(i4b),intent(in) :: ncid + integer(c_int),intent(in) :: indxGRU + character(*),intent(in) :: varName + character(*),intent(in) :: inFile ! file that populated dataLines + integer(i4b),intent(out) :: err + character(*),intent(inout) :: message + ! local variables + integer(i4b) :: localHRU_ix ! index of HRU + integer(i4b) :: varid ! netcdf variable id + integer(8) :: ncHruId(1) ! hruID from the forcing files + + + ! check to see if hruId exists as a variable, this is a required variable + err = nf90_inq_varid(ncid,trim(varname),varId) + if(err/=0)then; message=trim(message)//'hruID variable not present'; return; endif + + ! check that the hruId is what we expect + ! NOTE: we enforce that the HRU order in the forcing files is the same as in the zLocalAttributes files (too slow otherwise) + do localHRU_ix=1,gru_struc(indxGRU)%hruCount + ! check the HRU is what we expect + err = nf90_get_var(ncid,varId,ncHruId,start=(/gru_struc(indxGRU)%hruInfo(localHRU_ix)%hru_nc/),count=(/1/)) + if(gru_struc(indxGRU)%hruInfo(localHRU_ix)%hru_id /= ncHruId(1))then + write(message,'(a,i0,a,i0,a,i0,a,a)') trim(message)//'hruId for global HRU: ',gru_struc(indxGRU)%hruInfo(localHRU_ix)%hru_nc,' - ', & + ncHruId(1), ' differs from the expected: ',gru_struc(indxGRU)%hruInfo(localHRU_ix)%hru_id, ' in file ', trim(infile) + write(message,'(a)') trim(message)//' order of hruId in forcing file needs to match order in zLocalAttributes.nc' + err=40; return + endif + end do + +end subroutine + + end module ffile_info_module diff --git a/build/source/engine/mDecisions.f90 b/build/source/engine/mDecisions.f90 index 4ead60cf0ca3f420381fb6078b9d9cd7feaa26c3..595f37c5a65379e330e3944d7ece9b1f7f846d28 100755 --- a/build/source/engine/mDecisions.f90 +++ b/build/source/engine/mDecisions.f90 @@ -19,6 +19,7 @@ ! along with this program. If not, see <http://www.gnu.org/licenses/>. module mDecisions_module +USE, intrinsic :: iso_c_binding USE nrtype USE var_lookup, only: maxvarDecisions ! maximum number of decisions implicit none @@ -147,8 +148,8 @@ integer(i4b),parameter,public :: pahaut_76 = 314 ! Pahaut 1976, wi integer(i4b),parameter,public :: meltDripUnload = 321 ! Hedstrom and Pomeroy (1998), Storck et al 2002 (snowUnloadingCoeff & ratioDrip2Unloading) integer(i4b),parameter,public :: windUnload = 322 ! Roesch et al 2001, formulate unloading based on wind and temperature ! look-up values for the choice of energy equation -integer(i4b),parameter,public :: enthalpyFD = 323 ! enthalpyFD -integer(i4b),parameter,public :: closedForm = 324 ! closedForm +integer(i4b),parameter,public :: enthalpyFD = 323 ! enthalpyFD +integer(i4b),parameter,public :: closedForm = 324 ! closedForm ! ----------------------------------------------------------------------------------------------------------- contains @@ -156,7 +157,7 @@ contains ! ************************************************************************************************ ! public subroutine mDecisions: save model decisions as named integers ! ************************************************************************************************ -subroutine mDecisions(num_steps,err,message) +subroutine mDecisions(num_steps,err) bind(C, name='mDecisions') ! model time structures USE multiconst,only:secprday ! number of seconds in a day USE var_lookup,only:iLookTIME ! named variables that identify indices in the time structures @@ -185,12 +186,12 @@ subroutine mDecisions(num_steps,err,message) implicit none ! define output - integer(i4b),intent(out) :: num_steps - integer(i4b),intent(out) :: err ! error code - character(*),intent(out) :: message ! error message + integer(c_int),intent(out) :: num_steps + integer(c_int),intent(out) :: err ! error code ! define local variables + character(len=256) :: message ! error message character(len=256) :: cmessage ! error message for downwind routine - real(rkind) :: dsec,dsec_tz ! second + real(rkind) :: dsec,dsec_tz ! second ! initialize error control err=0; message='mDecisions/' diff --git a/build/source/netcdf/outputStrucWrite.f90 b/build/source/netcdf/outputStrucWrite.f90 deleted file mode 100755 index 2556ed8887205c804c64d08745fe4f6aa6fb3701..0000000000000000000000000000000000000000 --- a/build/source/netcdf/outputStrucWrite.f90 +++ /dev/null @@ -1,703 +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 contains subroutines for writing to the global output strucutre -module outputStrucWrite_module - -! NetCDF types -USE netcdf -USE netcdf_util_module,only:netcdf_err ! netcdf error handling function - -! top-level data types -USE nrtype - -! missing values -USE globalData,only: integerMissing, realMissing - -! provide access to global data -USE globalData,only:gru_struc ! gru->hru mapping structure -USE globalData,only:outputStructure - - -! provide access to the derived types to define the data structures -USE data_types,only:& - ! final data vectors - dlength, & ! var%dat - ilength, & ! var%dat - ! no spatial dimension - var_i, & ! x%var(:) (i4b) - var_i8, & ! x%var(:) integer(8) - var_d, & ! x%var(:) (dp) - var_ilength, & ! x%var(:)%dat (i4b) - var_dlength, & ! x%var(:)%dat (dp) - ! no variable dimension - hru_i, & ! x%hru(:) (i4b) - hru_d, & ! x%hru(:) (dp) - ! gru dimension - gru_int, & ! x%gru(:)%var(:) (i4b) - gru_double, & ! x%gru(:)%var(:) (dp) - gru_intVec, & ! x%gru(:)%var(:)%dat (i4b) - gru_doubleVec, & ! x%gru(:)%var(:)%dat (dp) - ! gru+hru dimension - gru_hru_int, & ! x%gru(:)%hru(:)%var(:) (i4b) - gru_hru_int8, & ! x%gru(:)%hru(:)%var(:) integer(8) - gru_hru_double, & ! x%gru(:)%hru(:)%var(:) (dp) - gru_hru_intVec, & ! x%gru(:)%hru(:)%var(:)%dat (i4b) - gru_hru_doubleVec ! x%gru(:)%hru(:)%var(:)%dat (dp) - -! vector lengths -USE var_lookup, only: maxvarFreq ! number of output frequencies -USE var_lookup, only: maxvarStat ! number of statistics - -implicit none -private -public::writeParm -public::writeData -public::writeBasin -public::writeTime -public::writeRestart -! define dimension lengths -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) - ! USE globalData,only:ncid ! netcdf file ids - USE data_types,only:var_info ! metadata info - USE var_lookup,only:iLookStat ! index in statistics vector - USE var_lookup,only:iLookFreq ! index in vector of model output frequencies - implicit none - - ! declare input variables - integer(i4b) ,intent(in) :: indxGRU ! Index into output Structure - integer(i4b) ,intent(in) :: indxHRU ! Index into output Structure - integer(i4b) ,intent(in) :: iSpatial ! hydrologic response unit - class(*) ,intent(in) :: struct ! data structure - type(var_info),intent(in) :: meta(:) ! metadata structure - character(*) ,intent(in) :: structName ! Name to know which global struct to write to - integer(i4b) ,intent(out) :: err ! error code - character(*) ,intent(out) :: message ! error message - ! local variables - integer(i4b) :: iVar ! loop through variables - - ! initialize error control - err=0;message="outputStrucWrite.f90-writeParm/" - - ! loop through local column model parameters - do iVar = 1,size(meta) - - ! check that the variable is desired - if (meta(iVar)%statIndex(iLookFREQ%timestep)==integerMissing) cycle - - ! initialize message - message=trim(message)//trim(meta(iVar)%varName)//'/' - - ! HRU data - if (iSpatial/=integerMissing) then - select type (struct) - class is (var_i) - if (structName == "type")then - outputStructure(1)%typeStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar) = struct%var(iVar) - end if - class is (var_i8) - - class is (var_d) - if (structName == "attr")then - outputStructure(1)%attrStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar) = struct%var(iVar) - end if - class is (var_dlength) - if (structName == "mpar")then - outputStructure(1)%mparStruct(1)%gru(indxGRU)%hru(indxHRU)%var(iVar) = struct%var(iVar) - end if - - class default; err=20; message=trim(message)//'unknown variable type (with HRU)'; return - end select - call netcdf_err(err,message); if (err/=0) return - - ! GRU data - else - select type (struct) - class is (var_d) - if (structName == "bpar")then - outputStructure(1)%bparStruct(1)%gru(indxGRU)%var(iVar) = struct%var(iVar) ! this will overwrite data - print*, "bpar" - end if - class is (var_i8) - class default; err=20; message=trim(message)//'unknown variable type (no HRU)'; return - end select - end if - call netcdf_err(err,message); if (err/=0) return - - ! re-initialize message - message="writeParm/" - end do ! looping through local column model parameters - -end subroutine writeParm - - ! ************************************************************************************** - ! public subroutine writeData: write model time-dependent data - ! ************************************************************************************** -subroutine writeData(indxGRU,indxHRU,iStep,structName,finalizeStats, & - maxLayers,meta,stat,dat,map,indx,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 - USE var_lookup,only:iLookIndex ! index into index structure - USE var_lookup,only:iLookStat ! index into stat structure - USE globalData,only:outFreq ! output file information - USE get_ixName_module,only:get_varTypeName ! to access type strings for error messages - USE get_ixName_module,only:get_statName ! to access type strings for error messages - implicit none - ! declare dummy variables - integer(i4b) ,intent(in) :: indxGRU - integer(i4b) ,intent(in) :: indxHRU - integer(i4b) ,intent(in) :: iStep - character(*) ,intent(in) :: structName - logical(lgt) ,intent(in) :: finalizeStats(:) ! flags to finalize statistics - integer(i4b) ,intent(in) :: maxLayers ! maximum number of layers - type(var_info),intent(in) :: meta(:) ! meta data - class(*) ,intent(in) :: stat ! stats data - class(*) ,intent(in) :: dat ! timestep data - integer(i4b) ,intent(in) :: map(:) ! map into stats child struct - type(var_ilength) ,intent(in) :: indx ! index data - integer(i4b) ,intent(out) :: err ! error code - character(*) ,intent(out) :: message ! error message - ! local variables - integer(i4b) :: iVar ! variable index - integer(i4b) :: iStat ! statistics index - integer(i4b) :: iFreq ! frequency index - integer(i4b) :: nSnow ! number of snow layers - integer(i4b) :: nSoil ! number of soil layers - integer(i4b) :: nLayers ! total number of layers - ! output arrays - integer(i4b) :: datLength ! length of each data vector - integer(i4b) :: maxLength ! maximum length of each data vector - integer(i4b),parameter :: ixInteger=1001 ! named variable for integer - integer(i4b),parameter :: ixReal=1002 ! named variable for real - ! initialize error control - err=0;message="outputStrucWrite.f90-writeData/" - - ! loop through output frequencies - do iFreq=1,maxvarFreq - - ! skip frequencies that are not needed - if(.not.outFreq(iFreq)) cycle - - ! check that we have finalized statistics for a given frequency - if(.not.finalizeStats(iFreq)) cycle - - ! loop through model variables - do iVar = 1,size(meta) - - ! handle time first - if (meta(iVar)%varName=='time')then - ! 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) - 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 - - ! define the statistics index - iStat = meta(iVar)%statIndex(iFreq) - - ! check that the variable is desired - if (iStat==integerMissing.or.trim(meta(iVar)%varName)=='unknown') cycle - - ! stats output: only scalar variable type - if(meta(iVar)%varType==iLookVarType%scalarv) then - select type(stat) - 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) - case('prog') - outputStructure(1)%progStat(1)%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) - case('flux') - outputStructure(1)%fluxStat(1)%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) - case default - err=21; message=trim(message)//"Stats structure not found"; return - end select - class default; err=20; message=trim(message)//'stats must be scalarv and of type gru_hru_doubleVec'; return - end select ! stat - - ! non-scalar variables: regular data structures - else - - ! 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 - nSnow = indx%var(iLookIndex%nSnow)%dat(1) - outputStructure(1)%indxStruct(1)%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 - - ! get the length of each data vector - select case (meta(iVar)%varType) - case(iLookVarType%wLength); datLength = maxSpectral - case(iLookVarType%midToto); datLength = nLayers - case(iLookVarType%midSnow); datLength = nSnow - case(iLookVarType%midSoil); datLength = nSoil - case(iLookVarType%ifcToto); datLength = nLayers+1 - case(iLookVarType%ifcSnow); datLength = nSnow+1 - case(iLookVarType%ifcSoil); datLength = nSoil+1 - case default; cycle - end select ! vartype - - ! get the data vectors - select type (dat) - 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(:) - case('diag') - outputStructure(1)%diagStruct(1)%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(:) - 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(:) - class default; err=20; message=trim(message)//'data must not be scalarv and either of type var_dlength or var_ilength'; return - end select - - ! get the maximum length of each data vector - select case (meta(iVar)%varType) - case(iLookVarType%wLength); maxLength = maxSpectral - case(iLookVarType%midToto); maxLength = maxLayers - case(iLookVarType%midSnow); maxLength = maxLayers-nSoil - case(iLookVarType%midSoil); maxLength = nSoil - case(iLookVarType%ifcToto); maxLength = maxLayers+1 - case(iLookVarType%ifcSnow); maxLength = (maxLayers-nSoil)+1 - case(iLookVarType%ifcSoil); maxLength = nSoil+1 - case default; cycle - end select ! vartype - - end if ! not scalarv - - ! process error code - if (err/=0) message=trim(message)//trim(meta(iVar)%varName)//'_'//trim(get_statName(iStat)) - call netcdf_err(err,message); if (err/=0) return - - end do ! iVar - end do ! iFreq -end subroutine writeData - -! ************************************************************************************** -! public subroutine writeBasin: write basin-average variables -! ************************************************************************************** -subroutine writeBasin(indxGRU,indxHRU,iStep,finalizeStats,& - outputTimestep,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 - USE globalData,only:outFreq ! output file information - USE get_ixName_module,only:get_varTypeName ! to access type strings for error messages - USE get_ixName_module,only:get_statName ! to access type strings for error messages - implicit none - - ! declare dummy variables - integer(i4b) ,intent(in) :: indxGRU - integer(i4b) ,intent(in) :: indxHRU - integer(i4b) ,intent(in) :: iStep - logical(lgt) ,intent(in) :: finalizeStats(:) ! flags to finalize statistics - integer(i4b) ,intent(in) :: outputTimestep(:) ! output time step - type(var_info),intent(in) :: meta(:) ! meta data - type(dlength) ,intent(in) :: stat(:) ! stats data - type(dlength) ,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 - ! local variables - integer(i4b) :: iVar ! variable index - integer(i4b) :: iStat ! statistics index - integer(i4b) :: iFreq ! frequency index - ! initialize error control - err=0;message="outputStrucWrite.f90-writeBasin/" - - ! loop through output frequencies - do iFreq=1,maxvarFreq - - ! skip frequencies that are not needed - if(.not.outFreq(iFreq)) cycle - - ! check that we have finalized statistics for a given frequency - if(.not.finalizeStats(iFreq)) cycle - - ! loop through model variables - do iVar = 1,size(meta) - - ! define the statistics index - iStat = meta(iVar)%statIndex(iFreq) - - ! check that the variable is desired - if (iStat==integerMissing.or.trim(meta(iVar)%varName)=='unknown') cycle - - ! stats/data output - select data type - 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) - - 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) - end if - - case default - err=40; message=trim(message)//"unknownVariableType[name='"//trim(meta(iVar)%varName)//"';type='"//trim(get_varTypeName(meta(iVar)%varType))// "']"; return - end select ! variable type - - ! process error code - if (err.ne.0) message=trim(message)//trim(meta(iVar)%varName)//'_'//trim(get_statName(iStat)) - if (err/=0) return - - end do ! iVar - end do ! iFreq - -end subroutine writeBasin - - ! ************************************************************************************** - ! public subroutine writeTime: write current time to all files - ! ************************************************************************************** -subroutine writeTime(indxGRU,indxHRU,iStep,finalizeStats,meta,dat,err,message) - USE data_types,only:var_info ! metadata type - USE var_lookup,only:iLookStat ! index into stat structure - implicit none - - ! declare dummy variables - integer(i4b) ,intent(in) :: indxGRU - integer(i4b) ,intent(in) :: indxHRU - integer(i4b) ,intent(in) :: iStep - logical(lgt) ,intent(in) :: finalizeStats(:) ! flags to finalize statistics - type(var_info),intent(in) :: meta(:) ! meta data - integer ,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 - ! initialize error control - err=0;message="outputStrucWrite.f90-writeTime/" - - ! loop through output frequencies - do iFreq=1,maxvarFreq - - ! check that we have finalized statistics for a given frequency - if(.not.finalizeStats(iFreq)) cycle - - ! loop through model variables - do iVar = 1,size(meta) - - ! check instantaneous - 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) - if (err/=0) message=trim(message)//trim(meta(iVar)%varName) - if (err/=0) then; err=20; return; end if - - end do ! iVar - end do ! iFreq - -end subroutine writeTime - -! ********************************************************************************************************* -! public subroutine printRestartFile: print a re-start file -! ********************************************************************************************************* -subroutine writeRestart(filename, & ! intent(in): name of restart file - nGRU, & ! intent(in): number of GRUs - nHRU, & ! intent(in): number of HRUs - prog_meta, & ! intent(in): prognostics metadata - prog_data, & ! intent(in): prognostics data - bvar_meta, & ! intent(in): basin (gru) variable metadata - bvar_data, & ! intent(in): basin (gru) variable data - maxLayers, & ! intent(in): maximum number of layers - maxSnowLayers, & ! intent(in): maximum number of snow layers - indx_meta, & ! intent(in): index metadata - indx_data, & ! intent(in): index data - err,message) ! intent(out): error control - ! -------------------------------------------------------------------------------------------------------- - ! -------------------------------------------------------------------------------------------------------- - ! access the derived types to define the data structures - USE data_types,only:var_info ! metadata - ! access named variables defining elements in the data structures - USE var_lookup,only:iLookINDEX ! named variables for structure elements - USE var_lookup,only:iLookVarType ! named variables for structure elements - USE var_lookup,only:iLookBVAR ! named variables for structure elements - ! constants - USE globalData,only:gru_struc ! gru-hru mapping structures - ! external routines - USE netcdf_util_module,only:nc_file_close ! close netcdf file - USE netcdf_util_module,only:nc_file_open ! open netcdf file - USE globalData,only:nTimeDelay ! number of timesteps in the time delay histogram - - implicit none - ! -------------------------------------------------------------------------------------------------------- - ! input - character(len=256),intent(in) :: filename ! name of the restart file - integer(i4b),intent(in) :: nGRU ! number of GRUs - integer(i4b),intent(in) :: nHRU ! number of HRUs - type(var_info),intent(in) :: prog_meta(:) ! prognostic variable metadata - type(var_dlength),intent(in) :: prog_data ! prognostic vars - type(var_info),intent(in) :: bvar_meta(:) ! basin variable metadata - type(var_dlength),intent(in) :: bvar_data ! basin variables - type(var_info),intent(in) :: indx_meta(:) ! metadata - type(var_ilength),intent(in) :: indx_data ! indexing vars - ! output: error control - integer(i4b),intent(out) :: err ! error code - character(*),intent(out) :: message ! error message - ! -------------------------------------------------------------------------------------------------------- - ! dummy variables - integer(i4b), intent(in) :: maxLayers ! maximum number of total layers - integer(i4b), intent(in) :: maxSnowLayers ! maximum number of snow layers - - ! local variables - integer(i4b) :: ncid ! netcdf file id - integer(i4b),allocatable :: ncVarID(:) ! netcdf variable id - integer(i4b) :: ncSnowID ! index variable id - integer(i4b) :: ncSoilID ! index variable id - - integer(i4b) :: nSoil ! number of soil layers - integer(i4b) :: nSnow ! number of snow layers - integer(i4b) :: maxSnow ! maximum number of snow layers - integer(i4b) :: maxSoil ! maximum number of soil layers - integer(i4b) :: nLayers ! number of total layers - integer(i4b),parameter :: nSpectral=2 ! number of spectal bands - integer(i4b),parameter :: nScalar=1 ! size of a scalar - integer(i4b) :: nProgVars ! number of prognostic variables written to state file - - integer(i4b) :: hruDimID ! variable dimension ID - integer(i4b) :: gruDimID ! variable dimension ID - integer(i4b) :: tdhDimID ! variable dimension ID - integer(i4b) :: scalDimID ! variable dimension ID - integer(i4b) :: specDimID ! variable dimension ID - integer(i4b) :: midSnowDimID ! variable dimension ID - integer(i4b) :: midSoilDimID ! variable dimension ID - integer(i4b) :: midTotoDimID ! variable dimension ID - integer(i4b) :: ifcSnowDimID ! variable dimension ID - integer(i4b) :: ifcSoilDimID ! variable dimension ID - integer(i4b) :: ifcTotoDimID ! variable dimension ID - - character(len=32),parameter :: hruDimName ='hru' ! dimension name for HRUs - character(len=32),parameter :: gruDimName ='gru' ! dimension name for GRUs - character(len=32),parameter :: tdhDimName ='tdh' ! dimension name for time-delay basin variables - character(len=32),parameter :: scalDimName ='scalarv' ! dimension name for scalar data - character(len=32),parameter :: specDimName ='spectral' ! dimension name for spectral bands - character(len=32),parameter :: midSnowDimName='midSnow' ! dimension name for snow-only layers - 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 :: ifcSnowDimName='ifcSnow' ! dimension name for snow-only layers - character(len=32),parameter :: ifcSoilDimName='ifcSoil' ! dimension name for soil-only layers - character(len=32),parameter :: ifcTotoDimName='ifcToto' ! dimension name for layered variables - - integer(i4b) :: cHRU ! count of HRUs - integer(i4b) :: iHRU ! index of HRUs - integer(i4b) :: iGRU ! index of GRUs - integer(i4b) :: iVar ! variable index - logical(lgt) :: okLength ! flag to check if the vector length is OK - character(len=256) :: cmessage ! downstream error message - ! -------------------------------------------------------------------------------------------------------- - - ! initialize error control - err=0; message='writeRestart/' - - ! size of prognostic variable vector - nProgVars = size(prog_meta) - allocate(ncVarID(nProgVars+1)) ! include 1 additional basin variable in ID array (possibly more later) - - ! maximum number of soil layers - maxSoil = gru_struc(1)%hruInfo(1)%nSoil - - ! maximum number of snow layers - maxSnow = maxSnowLayers - - ! create file - !TODO: Verify if this is needed - err = nf90_create(trim(filename),nf90_classic_model,ncid) - message='iCreate[create]'; call netcdf_err(err,message); if(err/=0)return - - ! define dimensions - !TODO: Verify if this is needed - err = nf90_def_dim(ncid,trim(hruDimName) ,nHRU , hruDimID); message='iCreate[hru]' ; call netcdf_err(err,message); if(err/=0)return - err = nf90_def_dim(ncid,trim(gruDimName) ,nGRU , gruDimID); message='iCreate[gru]' ; call netcdf_err(err,message); if(err/=0)return - err = nf90_def_dim(ncid,trim(tdhDimName) ,nTimeDelay , tdhDimID); message='iCreate[tdh]' ; call netcdf_err(err,message); if(err/=0)return - err = nf90_def_dim(ncid,trim(scalDimName) ,nScalar , scalDimID); message='iCreate[scalar]' ; call netcdf_err(err,message); if(err/=0)return - err = nf90_def_dim(ncid,trim(specDimName) ,nSpectral , specDimID); message='iCreate[spectral]'; call netcdf_err(err,message); if(err/=0)return - err = nf90_def_dim(ncid,trim(midSoilDimName),maxSoil ,midSoilDimID); message='iCreate[ifcSoil]' ; call netcdf_err(err,message); if(err/=0)return - err = nf90_def_dim(ncid,trim(midTotoDimName),maxLayers ,midTotoDimID); message='iCreate[midToto]' ; call netcdf_err(err,message); if(err/=0)return - err = nf90_def_dim(ncid,trim(ifcSoilDimName),maxSoil+1 ,ifcSoilDimID); message='iCreate[ifcSoil]' ; call netcdf_err(err,message); if(err/=0)return - err = nf90_def_dim(ncid,trim(ifcTotoDimName),maxLayers+1,ifcTotoDimID); message='iCreate[ifcToto]' ; call netcdf_err(err,message); if(err/=0)return - if (maxSnow>0) err = nf90_def_dim(ncid,trim(midSnowDimName),maxSnow ,midSnowDimID); message='iCreate[ifcSnow]' ; call netcdf_err(err,message); if(err/=0)return - if (maxSnow>0) err = nf90_def_dim(ncid,trim(ifcSnowDimName),maxSnow+1 ,ifcSnowDimID); message='iCreate[ifcSnow]' ; call netcdf_err(err,message); if(err/=0)return - ! re-initialize error control - err=0; message='writeRestart/' - - ! define prognostic variables - do iVar = 1,nProgVars - if (prog_meta(iVar)%varType==iLookvarType%unknown) cycle - - ! define variable - select case(prog_meta(iVar)%varType) - !TODO: Verify if this is needed - case(iLookvarType%scalarv); err = nf90_def_var(ncid,trim(prog_meta(iVar)%varname),nf90_double,(/hruDimID, scalDimID /),ncVarID(iVar)) - case(iLookvarType%wLength); err = nf90_def_var(ncid,trim(prog_meta(iVar)%varname),nf90_double,(/hruDimID, specDimID /),ncVarID(iVar)) - case(iLookvarType%midSoil); err = nf90_def_var(ncid,trim(prog_meta(iVar)%varname),nf90_double,(/hruDimID,midSoilDimID/),ncVarID(iVar)) - case(iLookvarType%midToto); err = nf90_def_var(ncid,trim(prog_meta(iVar)%varname),nf90_double,(/hruDimID,midTotoDimID/),ncVarID(iVar)) - case(iLookvarType%ifcSoil); err = nf90_def_var(ncid,trim(prog_meta(iVar)%varname),nf90_double,(/hruDimID,ifcSoilDimID/),ncVarID(iVar)) - case(iLookvarType%ifcToto); err = nf90_def_var(ncid,trim(prog_meta(iVar)%varname),nf90_double,(/hruDimID,ifcTotoDimID/),ncVarID(iVar)) - case(iLookvarType%midSnow); if (maxSnow>0) err = nf90_def_var(ncid,trim(prog_meta(iVar)%varname),nf90_double,(/hruDimID,midSnowDimID/),ncVarID(iVar)) - case(iLookvarType%ifcSnow); if (maxSnow>0) err = nf90_def_var(ncid,trim(prog_meta(iVar)%varname),nf90_double,(/hruDimID,ifcSnowDimID/),ncVarID(iVar)) - end select - - ! check errors - if(err/=0)then - message=trim(message)//trim(cmessage)//' [variable '//trim(prog_meta(iVar)%varName)//']' - return - end if - - ! add parameter description - !TODO: Verify if this is needed - err = nf90_put_att(ncid,ncVarID(iVar),'long_name',trim(prog_meta(iVar)%vardesc)) - call netcdf_err(err,message) - - ! add parameter units - !TODO: Verify if this is needed - err = nf90_put_att(ncid,ncVarID(iVar),'units',trim(prog_meta(iVar)%varunit)) - call netcdf_err(err,message) - - end do ! iVar - - ! define selected basin variables (derived) -- e.g., hillslope routing - !TODO: Verify if this is needed - err = nf90_def_var(ncid, trim(bvar_meta(iLookBVAR%routingRunoffFuture)%varName), nf90_double, (/gruDimID, tdhDimID /), ncVarID(nProgVars+1)) - err = nf90_put_att(ncid,ncVarID(nProgVars+1),'long_name',trim(bvar_meta(iLookBVAR%routingRunoffFuture)%vardesc)); call netcdf_err(err,message) - err = nf90_put_att(ncid,ncVarID(nProgVars+1),'units' ,trim(bvar_meta(iLookBVAR%routingRunoffFuture)%varunit)); call netcdf_err(err,message) - - ! define index variables - snow - !TODO: Verify if this is needed - err = nf90_def_var(ncid,trim(indx_meta(iLookIndex%nSnow)%varName),nf90_int,(/hruDimID/),ncSnowID); call netcdf_err(err,message) - err = nf90_put_att(ncid,ncSnowID,'long_name',trim(indx_meta(iLookIndex%nSnow)%vardesc)); call netcdf_err(err,message) - err = nf90_put_att(ncid,ncSnowID,'units' ,trim(indx_meta(iLookIndex%nSnow)%varunit)); call netcdf_err(err,message) - - ! define index variables - soil - !TODO: Verify if this is needed - err = nf90_def_var(ncid,trim(indx_meta(iLookIndex%nSoil)%varName),nf90_int,(/hruDimID/),ncSoilID); call netcdf_err(err,message) - err = nf90_put_att(ncid,ncSoilID,'long_name',trim(indx_meta(iLookIndex%nSoil)%vardesc)); call netcdf_err(err,message) - err = nf90_put_att(ncid,ncSoilID,'units' ,trim(indx_meta(iLookIndex%nSoil)%varunit)); call netcdf_err(err,message) - - ! end definition phase - !TODO: Verify if this is needed - err = nf90_enddef(ncid); call netcdf_err(err,message); if (err/=0) return - - ! write variables - do iGRU = 1,nGRU - do iHRU = 1,gru_struc(iGRU)%hruCount - cHRU = gru_struc(iGRU)%hruInfo(iHRU)%hru_ix - do iVar = 1,size(prog_meta) - - ! excape if this variable is not used - if (prog_meta(iVar)%varType==iLookvarType%unknown) cycle - - ! actual number of layers - nSnow = gru_struc(iGRU)%hruInfo(iHRU)%nSnow - nSoil = gru_struc(iGRU)%hruInfo(iHRU)%nSoil - nLayers = nSoil + nSnow - - ! check size - ! NOTE: this may take time that we do not wish to use - okLength=.true. - select case (prog_meta(iVar)%varType) - ! case(iLookVarType%scalarv); okLength = (size(prog_data%gru(iGRU)%hru(iHRU)%var(iVar)%dat) == nScalar ) - ! case(iLookVarType%wlength); okLength = (size(prog_data%gru(iGRU)%hru(iHRU)%var(iVar)%dat) == nSpectral) - ! case(iLookVarType%midSoil); okLength = (size(prog_data%gru(iGRU)%hru(iHRU)%var(iVar)%dat) == nSoil ) - ! case(iLookVarType%midToto); okLength = (size(prog_data%gru(iGRU)%hru(iHRU)%var(iVar)%dat) == nLayers ) - ! case(iLookVarType%ifcSoil); okLength = (size(prog_data%gru(iGRU)%hru(iHRU)%var(iVar)%dat) == nSoil+1 ) - ! case(iLookVarType%ifcToto); okLength = (size(prog_data%gru(iGRU)%hru(iHRU)%var(iVar)%dat) == nLayers+1) - ! case(iLookVarType%midSnow); if (nSnow>0) okLength = (size(prog_data%gru(iGRU)%hru(iHRU)%var(iVar)%dat) == nSnow ) - ! case(iLookVarType%ifcSnow); if (nSnow>0) okLength = (size(prog_data%gru(iGRU)%hru(iHRU)%var(iVar)%dat) == nSnow+1 ) - case(iLookVarType%scalarv); okLength = (size(prog_data%var(iVar)%dat) == nScalar ) - case(iLookVarType%wlength); okLength = (size(prog_data%var(iVar)%dat) == nSpectral) - case(iLookVarType%midSoil); okLength = (size(prog_data%var(iVar)%dat) == nSoil ) - case(iLookVarType%midToto); okLength = (size(prog_data%var(iVar)%dat) == nLayers ) - case(iLookVarType%ifcSoil); okLength = (size(prog_data%var(iVar)%dat) == nSoil+1 ) - case(iLookVarType%ifcToto); okLength = (size(prog_data%var(iVar)%dat) == nLayers+1) - case(iLookVarType%midSnow); if (nSnow>0) okLength = (size(prog_data%var(iVar)%dat) == nSnow ) - case(iLookVarType%ifcSnow); if (nSnow>0) okLength = (size(prog_data%var(iVar)%dat) == nSnow+1 ) - case default; err=20; message=trim(message)//'unknown var type'; return - end select - - ! error check - if(.not.okLength)then - message=trim(message)//'bad vector length for variable '//trim(prog_meta(iVar)%varname) - err=20; return - endif - - ! write data - select case (prog_meta(iVar)%varType) - case(iLookVarType%scalarv); err=nf90_put_var(ncid,ncVarID(iVar),(/prog_data%var(iVar)%dat/),start=(/cHRU,1/),count=(/1,nScalar /)) - case(iLookVarType%wlength); err=nf90_put_var(ncid,ncVarID(iVar),(/prog_data%var(iVar)%dat/),start=(/cHRU,1/),count=(/1,nSpectral/)) - case(iLookVarType%midSoil); err=nf90_put_var(ncid,ncVarID(iVar),(/prog_data%var(iVar)%dat/),start=(/cHRU,1/),count=(/1,nSoil /)) - case(iLookVarType%midToto); err=nf90_put_var(ncid,ncVarID(iVar),(/prog_data%var(iVar)%dat/),start=(/cHRU,1/),count=(/1,nLayers /)) - case(iLookVarType%ifcSoil); err=nf90_put_var(ncid,ncVarID(iVar),(/prog_data%var(iVar)%dat/),start=(/cHRU,1/),count=(/1,nSoil+1 /)) - case(iLookVarType%ifcToto); err=nf90_put_var(ncid,ncVarID(iVar),(/prog_data%var(iVar)%dat/),start=(/cHRU,1/),count=(/1,nLayers+1/)) - case(iLookVarType%midSnow); if (nSnow>0) err=nf90_put_var(ncid,ncVarID(iVar),(/prog_data%var(iVar)%dat/),start=(/cHRU,1/),count=(/1,nSnow /)) - case(iLookVarType%ifcSnow); if (nSnow>0) err=nf90_put_var(ncid,ncVarID(iVar),(/prog_data%var(iVar)%dat/),start=(/cHRU,1/),count=(/1,nSnow+1 /)) - case default; err=20; message=trim(message)//'unknown var type'; return - end select - - ! error check - if (err.ne.0) message=trim(message)//'writing variable:'//trim(prog_meta(iVar)%varName) - call netcdf_err(err,message); if (err/=0) return - err=0; message='writeRestart/' - - end do ! iVar loop - - ! write index variables - !TODO: Verify if this is needed - err=nf90_put_var(ncid,ncSnowID,(/indx_data%var(iLookIndex%nSnow)%dat/),start=(/cHRU/),count=(/1/)) - err=nf90_put_var(ncid,ncSoilID,(/indx_data%var(iLookIndex%nSoil)%dat/),start=(/cHRU/),count=(/1/)) - - end do ! iHRU loop - - ! write selected basin variables - !TODO: Verify if this is needed - err=nf90_put_var(ncid,ncVarID(nProgVars+1),(/bvar_data%var(iLookBVAR%routingRunoffFuture)%dat/), start=(/iGRU/),count=(/1,nTimeDelay/)) - - end do ! iGRU loop - - ! close file - call nc_file_close(ncid,err,cmessage) - if(err/=0)then;message=trim(message)//trim(cmessage);return;end if - - ! cleanup - deallocate(ncVarID) - -end subroutine writeRestart - -end module outputStrucWrite_module