diff --git a/Apps/MAPL_GridCompSpecs_ACG.py b/Apps/MAPL_GridCompSpecs_ACG.py index b65e51d3a808..990e4c63ee2a 100755 --- a/Apps/MAPL_GridCompSpecs_ACG.py +++ b/Apps/MAPL_GridCompSpecs_ACG.py @@ -73,16 +73,11 @@ def get_rank(self): @staticmethod def internal_name(name): - if name[-1] == '*': - return name[:-1] - else: - return name + return name.replace('*','') + @staticmethod def mangled_name(name): - if name[-1] == '*': - return "'" + name[:-1] + "'//trim(comp_name)" - else: - return "'" + name + "'" + return "'" + name.replace("*","'//trim(comp_name)//'") + "'" # Pointers must be declared regardless of COND status. Deactivated # pointers should not be _referenced_ but such sections should still diff --git a/CHANGELOG.md b/CHANGELOG.md index 3c55472112fd..3423e26c7526 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -17,6 +17,25 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Deprecated +## [2.19.0] - 2022-03-18 + +### Fixed + +- Fixed duration of the clock to be the smaller of the user specified duration and (END_DATE - currTime) +- Fixes for compiling on M1 Macs (remove REAL128) +- Fix for CMake when `esmf` is already a target + +### Added + +- New cmake option USE_EXTDATA2G to enable the next generation of ExtData for development, by default uses 1st generation ExtData +- MAPL_ESMFFieldBundleRead/Write modules are now available in when using MAPL + +### Changed + +- Replaced a wild card "*" in any position of a string in MAPL_GridCompSpecs_ACG.py +- Updated `MAPL_SunGetSolarConstantFromNRLFile` to open NRL Solar Table file only on root and broadcast the tables to all processes. Now all processes do interpolation. +- Add voting interpolation method as optional argument to SimpleBundleRead method + ## [2.18.3] - 2022-03-15 ### Fixed @@ -54,6 +73,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Updates to CircleCI - Added GEOSadas CI ifort build test - Add "like-UFS" build to CI. This is no FLAP and pFlogger, and static build +- Added new `_STAT` and `_IOSTAT` macros a la `_RC` ### Changed @@ -66,6 +86,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Updated `components.yaml`. These changes are to support using Spack to build MAPL - ESMA_cmake v3.10.0 (add `FindESMF.cmake` from NOAA-EMC) - ecbuild geos/v1.2.0 (updat `FindNetCDF.cmake` to that from NOAA-EMC) +- Improved file-not-found error handling in ExtData for non-templated filenames ## [2.17.2] - 2022-02-16 diff --git a/CMakeLists.txt b/CMakeLists.txt index 52ac0ba6c524..5090ffa19888 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -4,7 +4,7 @@ cmake_policy (SET CMP0054 NEW) project ( MAPL - VERSION 2.18.3 + VERSION 2.19.0 LANGUAGES Fortran CXX C) # Note - CXX is required for ESMF # Set the default build type to release @@ -84,8 +84,16 @@ endif() if(NOT TARGET FARGPARSE::fargparse) find_package(FARGPARSE QUIET) endif() -if(NOT TARGET YAFYAML::yafyaml) - find_package(YAFYAML QUIET) + +option(USE_EXTDATA2G "Use ExtData2G" ON) +if(USE_EXTDATA2G) + set (EXTDATA2G_TARGET "MAPL.ExtData2G" CACHE STRING "ExtData2G Target") + if(NOT TARGET YAFYAML::yafyaml) + find_package(YAFYAML REQUIRED) + endif() + message (STATUS "Building with ExtData2G") +else() + set (EXTDATA2G_TARGET "" CACHE STRING "ExtData2G Target") endif() option(BUILD_WITH_PFLOGGER "Build MAPL with pFlogger library support" ON) @@ -114,7 +122,9 @@ if (NOT Baselibs_FOUND) add_definitions(-DH5_HAVE_PARALLEL) endif() - find_package(ESMF MODULE REQUIRED) + if (NOT TARGET esmf) + find_package(ESMF MODULE REQUIRED) + endif () endif () if (BUILD_WITH_PFLOGGER) diff --git a/MAPL/CMakeLists.txt b/MAPL/CMakeLists.txt index d91f0081f870..a7d0f97de35d 100644 --- a/MAPL/CMakeLists.txt +++ b/MAPL/CMakeLists.txt @@ -3,11 +3,13 @@ esma_set_this() esma_add_library (${this} SRCS MAPL.F90 - DEPENDENCIES MAPL.base MAPL.generic MAPL.pfio MAPL_cfio_r4 MAPL.gridcomps MAPL.orbit MAPL.griddedio + DEPENDENCIES MAPL.base MAPL.generic MAPL.pfio MAPL_cfio_r4 MAPL.gridcomps MAPL.orbit MAPL.griddedio ${EXTDATA_TARGET} esmf NetCDF::NetCDF_Fortran MPI::MPI_Fortran $<$:FLAP::FLAP> TYPE ${MAPL_LIBRARY_TYPE} ) +target_compile_definitions (${this} PRIVATE $<$:BUILD_WITH_EXTDATA2G>) + target_include_directories (${this} PUBLIC $) diff --git a/MAPL/MAPL.F90 b/MAPL/MAPL.F90 index 6022c41e7dd3..6751de902c0a 100644 --- a/MAPL/MAPL.F90 +++ b/MAPL/MAPL.F90 @@ -4,11 +4,12 @@ module MAPL use MAPLBase_mod use MAPL_GenericMod use MAPL_VarSpecMiscMod - use MAPL_ExtDataGridCompMod, only: T_EXTDATA_STATE, EXTDATA_WRAP use ESMF_CFIOMod use pFIO use MAPL_GridCompsMod use mapl_StubComponent + use MAPL_ESMFFieldBundleRead + use MAPL_ESMFFieldBundleWrite implicit none end module MAPL diff --git a/Tests/CMakeLists.txt b/Tests/CMakeLists.txt index e3e9ae8c9808..bba1dbb97aa7 100644 --- a/Tests/CMakeLists.txt +++ b/Tests/CMakeLists.txt @@ -17,6 +17,7 @@ if (BUILD_WITH_FLAP) target_link_libraries(ExtDataDriver.x PRIVATE OpenMP::OpenMP_Fortran) endif () set_target_properties(ExtDataDriver.x PROPERTIES Fortran_MODULE_DIRECTORY ${MODULE_DIRECTORY}) + target_compile_definitions (ExtDataDriver.x PRIVATE $<$:BUILD_WITH_EXTDATA2G>) ecbuild_add_executable (TARGET pfio_MAPL_demo.x SOURCES pfio_MAPL_demo.F90) target_link_libraries (pfio_MAPL_demo.x PRIVATE MAPL FLAP::FLAP esmf) diff --git a/Tests/ExtDataDriverGridComp.F90 b/Tests/ExtDataDriverGridComp.F90 index f93681785932..1b394018695d 100644 --- a/Tests/ExtDataDriverGridComp.F90 +++ b/Tests/ExtDataDriverGridComp.F90 @@ -4,7 +4,10 @@ module ExtData_DriverGridCompMod use ESMF use MAPL - use MAPL_ExtDataGridCompMod, only : ExtData_SetServices => SetServices +#if defined(BUILD_WITH_EXTDATA2G) + use MAPL_ExtDataGridComp2G, only : ExtData2G_SetServices => SetServices +#endif + use MAPL_ExtDataGridCompMod, only : ExtData1G_SetServices => SetServices use MAPL_HistoryGridCompMod, only : Hist_SetServices => SetServices use MAPL_Profiler, only : get_global_time_profiler, BaseProfiler @@ -140,6 +143,7 @@ subroutine initialize_gc(gc, import_state, export_state, clock, rc) procedure(), pointer :: root_set_services type(ExtData_DriverGridComp), pointer :: cap class(BaseProfiler), pointer :: t_p + logical :: use_extdata2g _UNUSED_DUMMY(import_state) _UNUSED_DUMMY(export_state) @@ -212,6 +216,7 @@ subroutine initialize_gc(gc, import_state, export_state, clock, rc) !EOR enableTimers = ESMF_UtilStringUpperCase(enableTimers, rc = status) _VERIFY(status) + call MAPL_GetResource(maplobj,use_extdata2g,"USE_EXTDATA2G:",default=.false.,_RC) if (enableTimers /= 'YES') then call MAPL_ProfDisable(rc = status) @@ -319,10 +324,16 @@ subroutine initialize_gc(gc, import_state, export_state, clock, rc) call MAPL_Set(MAPLOBJ, CF=CAP%CF_EXT, RC=STATUS) _VERIFY(STATUS) - - cap%extdata_id = MAPL_AddChild (MAPLOBJ, name = 'EXTDATA', SS = ExtData_SetServices, rc = status) - _VERIFY(status) - + if (use_extdata2g) then +#if defined(BUILD_WITH_EXTDATA2G) + cap%extdata_id = MAPL_AddChild (MAPLOBJ, name = 'EXTDATA', SS = ExtData2G_SetServices, _RC) +#else + _FAIL('ExtData2G requested but not built') +#endif + else + cap%extdata_id = MAPL_AddChild (MAPLOBJ, name = 'EXTDATA', SS = ExtData1G_SetServices, _RC) + end if + end if ! Query MAPL for the the children's for GCS, IMPORTS, EXPORTS diff --git a/base/FileMetadataUtilities.F90 b/base/FileMetadataUtilities.F90 index 13e9e503569d..64428d1a29ad 100644 --- a/base/FileMetadataUtilities.F90 +++ b/base/FileMetadataUtilities.F90 @@ -417,6 +417,7 @@ subroutine get_time_info(this,startTime,startyear,startmonth,startday,starthour, allocate(timeVector,source=tVec,stat=status) _VERIFY(status) end if + _RETURN(_SUCCESS) end subroutine get_time_info @@ -516,6 +517,7 @@ subroutine get_coordinate_info(this,coordinate_name,coordSize,coordUnits,coords, _ASSERT(.false.,"unsupported coordel variable type") end select end if + _RETURN(_SUCCESS) end subroutine get_coordinate_info diff --git a/base/MAPL_SimpleBundleMod.F90 b/base/MAPL_SimpleBundleMod.F90 index 0c27e5ad715c..3e9ab23b9ad2 100644 --- a/base/MAPL_SimpleBundleMod.F90 +++ b/base/MAPL_SimpleBundleMod.F90 @@ -749,7 +749,8 @@ end subroutine MAPL_SimpleBundleDestroy ! Function MAPL_SimpleBundleRead (filename, bundle_name, grid, time, verbose, & - only_vars, expid, rc ) result (self) + only_vars, expid, voting, unusable, rc ) result (self) + use mapl_KeywordEnforcerMod ! !ARGUMENTS: @@ -762,6 +763,8 @@ Function MAPL_SimpleBundleRead (filename, bundle_name, grid, time, verbose, & logical, OPTIONAL, intent(in) :: verbose character(len=*), optional, intent(IN) :: only_vars character(len=*), optional, intent(IN) :: expid + class(KeywordEnforcer), optional, intent(in) :: unusable + logical, optional, intent(in) :: voting integer, OPTIONAL, intent(out) :: rc ! !DESCRIPTION: @@ -781,7 +784,7 @@ Function MAPL_SimpleBundleRead (filename, bundle_name, grid, time, verbose, & Bundle = ESMF_FieldBundleCreate ( name=bundle_name, __RC__ ) call ESMF_FieldBundleSet ( bundle, grid=Grid, __RC__ ) call MAPL_CFIORead ( filename, Time, Bundle, verbose=verbose, & - ONLY_VARS=only_vars, expid=expid, __RC__ ) + ONLY_VARS=only_vars, expid=expid, voting=voting, __RC__ ) self = MAPL_SimpleBundleCreate ( Bundle, __RC__ ) self%bundleAlloc = .true. diff --git a/base/MAPL_sun_uc.F90 b/base/MAPL_sun_uc.F90 index 0773fe2cede3..14b9d7166a6f 100644 --- a/base/MAPL_sun_uc.F90 +++ b/base/MAPL_sun_uc.F90 @@ -2475,54 +2475,45 @@ subroutine MAPL_SunGetSolarConstantFromNRLFile(CLOCK,filename_in,SC,MG,SB,Persis CREATE_TABLE: if (.not. TableCreated) then - ! Open the file - ! ------------- - - filename = trim(filename_in) - - ! Does the file exist? - inquire( FILE=filename, EXIST=found ) - _ASSERT( found ,'Could not find NRL data file '//trim(filename) ) - - UNIT = GETFILE(filename, DO_OPEN=0, form="formatted", rc=status) - _VERIFY(STATUS) - - open(unit=unit, file=filename) + ! First we open the file on root to get the + ! number of lines so we can allocate our tables + ! --------------------------------------------- if (amIRoot) then + ! Open the file + ! ------------- + filename = trim(filename_in) + open(newunit=unit, file=filename, form="formatted", status="old", iostat=status) + _ASSERT(status==0,'Could not find NRL data file '// trim(filename )) + ! Determine length of file ! ------------------------ - call lgr%debug("Scanning the Solar Table to determine number of data points") - numlines = num_lines_in_file(UNIT) - call lgr%debug("Solar Table Data Points: %i0", numlines) - ! Allocate our arrays - ! ------------------- + end if - allocate(yearTable(numlines), source=0, stat=status) - _VERIFY(STATUS) + ! Broadcast the number of lines + ! ----------------------------- + call MAPL_CommsBcast(vm, DATA=numlines, N=1, ROOT=0, _RC) - allocate(doyTable(numlines), source=0, stat=status) - _VERIFY(STATUS) + ! Allocate our arrays on all processes + ! ------------------------------------ - allocate(tsi(numlines), source=0.0, stat=status) - _VERIFY(STATUS) + allocate(yearTable(numlines), source=0, _STAT) + allocate(doyTable(numlines), source=0, _STAT) + allocate(tsi(numlines), source=0.0, _STAT) + allocate(mgindex(numlines), source=0.0, _STAT) + allocate(sbindex(numlines), source=0.0, _STAT) - allocate(mgindex(numlines), source=0.0, stat=status) - _VERIFY(STATUS) + ! Back to root to read in the values + ! ---------------------------------- - allocate(sbindex(numlines), source=0.0, stat=status) - _VERIFY(STATUS) - - ! Read in arrays - ! -------------- + if (amIRoot) then call lgr%debug("Reading the Solar Table") - i = 1 do read(unit,'(A)',iostat=stat) line @@ -2536,51 +2527,75 @@ subroutine MAPL_SunGetSolarConstantFromNRLFile(CLOCK,filename_in,SC,MG,SB,Persis ! Belt and suspenders check that all data was read _ASSERT(size(yearTable) == numlines,"Inconsistency in NRL number of lines") + close(unit, _IOSTAT) + end if - ! Close the file - ! -------------- + ! Broadcast the tables + ! -------------------- - call FREE_FILE(UNIT) + call MAPL_CommsBcast(vm, DATA=yearTable, N=numlines, ROOT=0, _RC) + call MAPL_CommsBcast(vm, DATA=doyTable, N=numlines, ROOT=0, _RC) + call MAPL_CommsBcast(vm, DATA=tsi, N=numlines, ROOT=0, _RC) + call MAPL_CommsBcast(vm, DATA=mgindex, N=numlines, ROOT=0, _RC) + call MAPL_CommsBcast(vm, DATA=sbindex, N=numlines, ROOT=0, _RC) TableCreated = .TRUE. end if CREATE_TABLE - ON_ROOT: if (amIRoot) then + ! Now we need to find the two bracketing days + ! ------------------------------------------- - ! Now we need to find the two bracketing days - ! ------------------------------------------- + ! Get current time + ! ---------------- + call ESMF_ClockGet(CLOCK, CURRTIME=currentTime, _RC) - ! Get current time - ! ---------------- - call ESMF_ClockGet(CLOCK, CURRTIME=currentTime, RC=STATUS) - _VERIFY(STATUS) + call ESMF_TimeGet( currentTime, YY = currentYear, & + MM = currentMon, & + DD = currentDay, & + dayOfYear = currentDOY, _RC) - call ESMF_TimeGet( currentTime, YY = currentYear, & - MM = currentMon, & - DD = currentDay, & - dayOfYear = currentDOY, & - RC = STATUS ) - _VERIFY(STATUS) + ! Test if current time is outside our file + ! ---------------------------------------- - ! Test if current time is outside our file - ! ---------------------------------------- + outOfTable = .FALSE. - outOfTable = .FALSE. + ! First is current year higher than last in file... + if ( currentYear > yearTable(numlines) ) then + outOfTable = .TRUE. + ! ...or if a partial year, are we near the end + else if ( currentYear == yearTable(numlines) .and. currentDOY >= doyTable(numlines)) then + outOfTable = .TRUE. + end if - ! First is current year higher than last in file... - if ( currentYear > yearTable(numlines) ) then - outOfTable = .TRUE. - ! ...or if a partial year, are we near the end - else if ( currentYear == yearTable(numlines) .and. currentDOY >= doyTable(numlines)) then - outOfTable = .TRUE. - end if + ! If we are out of the table... + ! ----------------------------- + + OUT_OF_TABLE: if ( outOfTable ) then + + PERSIST_SOLAR: if ( PersistSolar_ ) then + + ! If we are outOfTable and we have the PersistSolar + ! option we just use the last value in the table... + ! ------------------------------------------------- + + SC = tsi(numlines) + MG = mgindex(numlines) + SB = sbindex(numlines) + + call lgr%debug('Off the end of table, persisting last values') + call lgr%debug(' tsi at end of table: %F8.3', tsi(numlines)) + call lgr%debug(' mgindex at end of table: %F8.6', mgindex(numlines)) + call lgr%debug(' sbindex at end of table: %F9.4', sbindex(numlines)) - ! If we are out of the table and not persisting, we must - ! recenter our day to be based on the last complete Solar Cycle - ! ------------------------------------------------------------- - OUT_OF_TABLE_AND_CYCLE: if ( outOfTable .and. (.not. PersistSolar_) ) then + _RETURN(ESMF_SUCCESS) + + else + + ! If we are out of the table and not persisting, we must + ! recenter our day to be based on the last complete Solar Cycle + ! ------------------------------------------------------------- ! Create an ESMF_Time at start of Cycle 24 ! ---------------------------------------- @@ -2589,9 +2604,7 @@ subroutine MAPL_SunGetSolarConstantFromNRLFile(CLOCK,filename_in,SC,MG,SB,Persis DD = 1, & H = 12, & M = 00, & - S = 00, & - RC = STATUS ) - _VERIFY(STATUS) + S = 00, _RC) ! Create an ESMF_Time at start of Cycle 25 ! ---------------------------------------- @@ -2600,9 +2613,7 @@ subroutine MAPL_SunGetSolarConstantFromNRLFile(CLOCK,filename_in,SC,MG,SB,Persis DD = 1, & H = 12, & M = 00, & - S = 00, & - RC = STATUS ) - _VERIFY(STATUS) + S = 00, _RC) ! Create TimeInterval based on interval ! from start of latest Cycle 25 @@ -2631,126 +2642,86 @@ subroutine MAPL_SunGetSolarConstantFromNRLFile(CLOCK,filename_in,SC,MG,SB,Persis ! Get new currentYear, currentMon, currentDay ! ------------------------------------------- - call ESMF_TimeGet( currentTime, YY = currentYear, & - MM = currentMon, & - DD = currentDay, & - dayOfYear = currentDOY, & - RC = STATUS ) - _VERIFY(STATUS) - + call ESMF_TimeGet( currentTime, YY = currentYear, & + MM = currentMon, & + DD = currentDay, & + dayOfYear = currentDOY, _RC) - ! Debugging Prints - ! ---------------- call lgr%debug('Off the end of table, moving into last complete cycle') call lgr%debug(' Original Year-Mon-Day to Find: %i0.4~-%i0.2~-%i0.2', originalYear,originalMon,originalDay) call lgr%debug(' Original Day of Year: %i0', origDOY) call lgr%debug(' New Year-Mon-Day to Find: %i0.4~-%i0.2~-%i0.2', currentYear,currentMon,currentDay) call lgr%debug(' New Day of Year: %i0', currentDOY) - end if OUT_OF_TABLE_AND_CYCLE - - ! Create an ESMF_Time on noon of current day - ! ------------------------------------------ - call ESMF_TimeSet( noonCurrentDay, YY = currentYear, & - MM = currentMon, & - DD = currentDay, & - H = 12, & - M = 00, & - S = 00, & - RC = STATUS ) - _VERIFY(STATUS) - - ! Figure out bracketing days for interpolation - ! NOTE: nextNoon is mainly for debugging purposes - ! ----------------------------------------------- - call ESMF_TimeIntervalSet(oneDayInterval, D=1, rc=status) - if (currentTime <= noonCurrentDay) then - prevNoon = noonCurrentDay - oneDayInterval - nextNoon = noonCurrentDay - else - prevNoon = noonCurrentDay - nextNoon = noonCurrentDay + oneDayInterval - end if + end if PERSIST_SOLAR - ! Get the DOYs - ! ------------ - call ESMF_TimeGet( prevNoon, YY = prevNoonYear, dayOfYear = prevDOY, rc = status ) - call ESMF_TimeGet( nextNoon, YY = nextNoonYear, dayOfYear = nextDOY, rc = status ) + end if OUT_OF_TABLE - ! Our interpolation factor is based of when we are compared to the next noon - ! -------------------------------------------------------------------------- - intToNextNoon = nextNoon-currentTime - - ! The FAC for interpolating is just the real version - ! of the size of the timeinterval to the next noon - ! -------------------------------------------------- - call ESMF_TimeIntervalGet(intToNextNoon, d_r8=days_r8, rc=STATUS) - _VERIFY(STATUS) - FAC = real(days_r8) - - ! Use our find_file_index function to get the index for previous noon - ! ------------------------------------------------------------------- - INDX1 = find_file_index(numlines, yearTable, prevNoonYear, prevDOY) - INDX2 = INDX1 + 1 + ! Create an ESMF_Time on noon of current day + ! ------------------------------------------ + call ESMF_TimeSet( noonCurrentDay, YY = currentYear, & + MM = currentMon, & + DD = currentDay, & + H = 12, & + M = 00, & + S = 00, _RC) + + ! Figure out bracketing days for interpolation + ! NOTE: nextNoon is mainly for debugging purposes + ! ----------------------------------------------- + call ESMF_TimeIntervalSet(oneDayInterval, D=1, _RC) + if (currentTime <= noonCurrentDay) then + prevNoon = noonCurrentDay - oneDayInterval + nextNoon = noonCurrentDay + else + prevNoon = noonCurrentDay + nextNoon = noonCurrentDay + oneDayInterval + end if - ! If we are outOfTable and we have the PersistSolar - ! option we just use the last value in the table... - ! ------------------------------------------------- - OUT_OF_TABLE_AND_PERSIST: if ( outOfTable .and. PersistSolar_) then + ! Get the DOYs + ! ------------ + call ESMF_TimeGet( prevNoon, YY = prevNoonYear, dayOfYear = prevDOY, _RC) + call ESMF_TimeGet( nextNoon, YY = nextNoonYear, dayOfYear = nextDOY, _RC) - SC = tsi(numlines) - MG = mgindex(numlines) - SB = sbindex(numlines) + ! Our interpolation factor is based of when we are compared to the next noon + ! -------------------------------------------------------------------------- + intToNextNoon = nextNoon-currentTime - ! Debugging Prints - ! ---------------- - call lgr%debug('Off the end of table, persisting last values') - call lgr%debug(' tsi at end of table: %F8.3', tsi(numlines)) - call lgr%debug(' mgindex at end of table: %F8.6', mgindex(numlines)) - call lgr%debug(' sbindex at end of table: %F9.4', sbindex(numlines)) + ! The FAC for interpolating is just the real version + ! of the size of the timeinterval to the next noon + ! -------------------------------------------------- + call ESMF_TimeIntervalGet(intToNextNoon, d_r8=days_r8, _RC) + FAC = real(days_r8) - ! Otherwise we interpolate to the table - ! ------------------------------------- - else + ! Use our find_file_index function to get the index for previous noon + ! ------------------------------------------------------------------- + INDX1 = find_file_index(numlines, yearTable, prevNoonYear, prevDOY) + INDX2 = INDX1 + 1 - ! Linear Interpolation to the given day-of-month - ! ---------------------------------------------- - - SC = tsi(INDX1)*FAC + tsi(INDX2)*(1.0-FAC) - MG = mgindex(INDX1)*FAC + mgindex(INDX2)*(1.0-FAC) - SB = sbindex(INDX1)*FAC + sbindex(INDX2)*(1.0-FAC) - - ! Debugging Prints - ! ---------------- - call lgr%debug(' First DOY to Find: %i3', prevDOY) - call lgr%debug(' file_index for date: %i6', INDX1) - call lgr%debug(' yearTable(date): %i4', yearTable(INDX1)) - call lgr%debug(' doyTable(date): %i3', doyTable(INDX1)) - call lgr%debug(' tsi(date): %f8.3', tsi(INDX1)) - call lgr%debug(' mgindex(date): %f8.6', mgindex(INDX1)) - call lgr%debug(' sbindex(date): %f9.4', sbindex(INDX1)) - - call lgr%debug(' Second DOY to Find: %i3', nextDOY) - call lgr%debug(' file_index for date: %i6', INDX2) - call lgr%debug(' yearTable(date): %i4', yearTable(INDX2)) - call lgr%debug(' doyTable(date): %i3', doyTable(INDX2)) - call lgr%debug(' tsi(date): %f8.3', tsi(INDX2)) - call lgr%debug(' mgindex(date): %f8.6', mgindex(INDX2)) - call lgr%debug(' sbindex(date): %f9.4', sbindex(INDX2)) - - call lgr%debug(' Interpolation Factor: %f8.6', FAC) - end if OUT_OF_TABLE_AND_PERSIST - end if ON_ROOT - - ! Broadcast the values - ! -------------------- + ! Linear Interpolation to the given day-of-month + ! ---------------------------------------------- - call MAPL_CommsBcast(vm, DATA=SC, N=1, ROOT=0, RC=status) - _VERIFY(STATUS) - call MAPL_CommsBcast(vm, DATA=MG, N=1, ROOT=0, RC=status) - _VERIFY(STATUS) - call MAPL_CommsBcast(vm, DATA=SB, N=1, ROOT=0, RC=status) - _VERIFY(STATUS) + SC = tsi(INDX1)*FAC + tsi(INDX2)*(1.0-FAC) + MG = mgindex(INDX1)*FAC + mgindex(INDX2)*(1.0-FAC) + SB = sbindex(INDX1)*FAC + sbindex(INDX2)*(1.0-FAC) + + call lgr%debug(' First DOY to Find: %i3', prevDOY) + call lgr%debug(' file_index for date: %i6', INDX1) + call lgr%debug(' yearTable(date): %i4', yearTable(INDX1)) + call lgr%debug(' doyTable(date): %i3', doyTable(INDX1)) + call lgr%debug(' tsi(date): %f8.3', tsi(INDX1)) + call lgr%debug(' mgindex(date): %f8.6', mgindex(INDX1)) + call lgr%debug(' sbindex(date): %f9.4', sbindex(INDX1)) + + call lgr%debug(' Second DOY to Find: %i3', nextDOY) + call lgr%debug(' file_index for date: %i6', INDX2) + call lgr%debug(' yearTable(date): %i4', yearTable(INDX2)) + call lgr%debug(' doyTable(date): %i3', doyTable(INDX2)) + call lgr%debug(' tsi(date): %f8.3', tsi(INDX2)) + call lgr%debug(' mgindex(date): %f8.6', mgindex(INDX2)) + call lgr%debug(' sbindex(date): %f9.4', sbindex(INDX2)) + + call lgr%debug(' Interpolation Factor: %f8.6', FAC) _RETURN(ESMF_SUCCESS) diff --git a/gridcomps/CMakeLists.txt b/gridcomps/CMakeLists.txt index 3a29088ffec6..6493a3ad2de6 100644 --- a/gridcomps/CMakeLists.txt +++ b/gridcomps/CMakeLists.txt @@ -17,3 +17,6 @@ add_subdirectory(Cap) add_subdirectory(History) add_subdirectory(Orbit) add_subdirectory(ExtData) +if(USE_EXTDATA2G) + add_subdirectory(ExtData2G) +endif() diff --git a/gridcomps/Cap/CMakeLists.txt b/gridcomps/Cap/CMakeLists.txt index c0ab60db408e..07a2fe92b3cb 100644 --- a/gridcomps/Cap/CMakeLists.txt +++ b/gridcomps/Cap/CMakeLists.txt @@ -4,6 +4,7 @@ set (srcs MAPL_CapGridComp.F90 MAPL_NUOPCWrapperMod.F90 CapOptions.F90 + ExternalGCStorage.F90 ) if (BUILD_WITH_FLAP) list (APPEND srcs FlapCLI.F90) @@ -11,10 +12,12 @@ endif() esma_add_library (${this} SRCS ${srcs} DEPENDENCIES MAPL.shared MAPL.constants MAPL.base MAPL.profiler MAPL.history - MAPL.ExtData TYPE ${MAPL_LIBRARY_TYPE}) + MAPL.ExtData ${EXTDATA2G_TARGET} TYPE ${MAPL_LIBRARY_TYPE}) target_link_libraries (${this} PUBLIC GFTL::gftl GFTL_SHARED::gftl-shared esmf NetCDF::NetCDF_Fortran PRIVATE MPI::MPI_Fortran $<$:FLAP::FLAP>) +target_compile_definitions (${this} PRIVATE $<$:BUILD_WITH_EXTDATA2G>) + # CMake has an OpenMP issue with NAG Fortran: https://gitlab.kitware.com/cmake/cmake/-/issues/21280 if (NOT CMAKE_Fortran_COMPILER_ID MATCHES "NAG") target_link_libraries(${this} PRIVATE OpenMP::OpenMP_Fortran) diff --git a/gridcomps/Cap/ExternalGCStorage.F90 b/gridcomps/Cap/ExternalGCStorage.F90 new file mode 100644 index 000000000000..8711c21b8126 --- /dev/null +++ b/gridcomps/Cap/ExternalGCStorage.F90 @@ -0,0 +1,14 @@ +module MAPL_ExternalGCStorage +use esmf +implicit none + +type t_extdata_state + type(ESMF_State) :: expState + type(ESMF_GridComp) :: gc +end type t_extdata_state + +type extdata_wrap + type (t_extdata_state), pointer :: PTR +end type extdata_wrap + +end module MAPL_ExternalGCStorage diff --git a/gridcomps/Cap/MAPL_CapGridComp.F90 b/gridcomps/Cap/MAPL_CapGridComp.F90 index 029356c9c8b6..2b8b92f67530 100644 --- a/gridcomps/Cap/MAPL_CapGridComp.F90 +++ b/gridcomps/Cap/MAPL_CapGridComp.F90 @@ -17,8 +17,10 @@ module MAPL_CapGridCompMod use MAPL_ShmemMod use MAPL_HistoryGridCompMod, only : Hist_SetServices => SetServices use MAPL_HistoryGridCompMod, only : HISTORY_ExchangeListWrap - use MAPL_ExtDataGridCompMod, only : ExtData_SetServices => SetServices - use MAPL_ExtDataGridCompMod, only : T_EXTDATA_STATE, EXTDATA_WRAP +#if defined(BUILD_WITH_EXTDATA2G) + use MAPL_ExtDataGridComp2G, only : ExtData2G_SetServices => SetServices +#endif + use MAPL_ExtDataGridCompMod, only : ExtData1G_SetServices => SetServices use MAPL_ConfigMod use MAPL_DirPathMod use MAPL_KeywordEnforcerMod @@ -28,6 +30,7 @@ module MAPL_CapGridCompMod use gFTL_StringVector use pflogger, only: logging, Logger use MAPL_TimeUtilsMod, only: is_valid_time, is_valid_date + use MAPL_ExternalGCStorage use iso_fortran_env @@ -173,8 +176,8 @@ subroutine initialize_gc(gc, import_state, export_state, clock, rc) integer :: status - type (T_ExtData_STATE), pointer :: ExtData_internal_state => null() - type (ExtData_wrap) :: wrap + type (t_extdata_state), pointer :: ExtData_internal_state => null() + type (extdata_wrap) :: wrap character(len=ESMF_MAXSTR ) :: timerModeStr @@ -210,6 +213,7 @@ subroutine initialize_gc(gc, import_state, export_state, clock, rc) class(BaseProfiler), pointer :: t_p class(Logger), pointer :: lgr type(ESMF_Clock) :: cap_clock + logical :: use_extdata2g _UNUSED_DUMMY(import_state) _UNUSED_DUMMY(export_state) @@ -399,6 +403,7 @@ subroutine initialize_gc(gc, import_state, export_state, clock, rc) !EOR enableTimers = ESMF_UtilStringUpperCase(enableTimers, rc = status) _VERIFY(status) + call MAPL_GetResource(maplobj,use_extdata2g,"USE_EXTDATA2G:",default=.false.,_RC) if (enableTimers /= 'YES') then call MAPL_ProfDisable(rc = status) @@ -566,8 +571,16 @@ subroutine initialize_gc(gc, import_state, export_state, clock, rc) call MAPL_Set(MAPLOBJ, CF=CAP%CF_EXT, RC=STATUS) _VERIFY(STATUS) - cap%extdata_id = MAPL_AddChild (MAPLOBJ, name = 'EXTDATA', SS = ExtData_SetServices, rc = status) - _VERIFY(status) + if (use_extdata2g) then +#if defined(BUILD_WITH_EXTDATA2G) + cap%extdata_id = MAPL_AddChild (MAPLOBJ, name = 'EXTDATA', SS = ExtData2G_SetServices, _RC) +#else + call lgr%error('ExtData2G requested but not built') + _FAIL('ExtData2G requested but not built') +#endif + else + cap%extdata_id = MAPL_AddChild (MAPLOBJ, name = 'EXTDATA', SS = ExtData1G_SetServices, _RC) + end if call t_p%stop('SetService') ! Add NX and NY from AGCM.rc to ExtData.rc as well as name of ExtData rc file @@ -1573,6 +1586,7 @@ subroutine MAPL_ClockInit ( MAPLOBJ, Clock, nsteps, rc) type(ESMF_Time) :: CurrTime ! Current Current Time of Experiment type(ESMF_TimeInterval) :: timeStep ! HEARTBEAT type(ESMF_TimeInterval) :: duration + type(ESMF_TimeInterval) :: maxDuration type(ESMF_Calendar) :: cal character(ESMF_MAXSTR) :: calendar @@ -1824,6 +1838,9 @@ subroutine MAPL_ClockInit ( MAPLOBJ, Clock, nsteps, rc) rc = STATUS ) _VERIFY(STATUS) + maxDuration = EndTime - currTime + if (duration > maxDuration) duration = maxDuration + stopTime = currTime + duration ! initialize model time step diff --git a/gridcomps/ExtData/ExtDataGridCompMod.F90 b/gridcomps/ExtData/ExtDataGridCompMod.F90 index c728b2e77cbc..9fa76a0fa575 100644 --- a/gridcomps/ExtData/ExtDataGridCompMod.F90 +++ b/gridcomps/ExtData/ExtDataGridCompMod.F90 @@ -64,8 +64,6 @@ MODULE MAPL_ExtDataGridCompMod ! !PUBLIC MEMBER FUNCTIONS: PUBLIC SetServices - public T_EXTDATA_STATE - public EXTDATA_WRAP !EOP ! ! !REVISION HISTORY: @@ -212,17 +210,6 @@ MODULE MAPL_ExtDataGridCompMod type (MAPL_ExtData_State), pointer :: PTR => null() end type MAPL_ExtData_WRAP - type T_EXTDATA_STATE - type(ESMF_State) :: expState - type(ESMF_GridComp) :: gc - end type T_EXTDATA_STATE - - ! Wrapper for extracting internal state - ! ------------------------------------- - type EXTDATA_WRAP - type (T_EXTDATA_STATE), pointer :: PTR - end type EXTDATA_WRAP - class(Logger), pointer :: lgr @@ -675,10 +662,6 @@ SUBROUTINE Initialize_ ( GC, IMPORT, EXPORT, CLOCK, rc ) primary%item(totalPrimaryEntries)%cyclic='n' end if - - if ( primary%item(totalPrimaryEntries)%isConst .eqv. .false. ) then - call CreateTimeInterval(primary%item(totalPrimaryEntries),clock,__RC__) - end if end if enddo ! Derived Exports @@ -911,6 +894,10 @@ SUBROUTINE Initialize_ ( GC, IMPORT, EXPORT, CLOCK, rc ) call lgr%debug('ExtData Initialize_(): PrimaryLoop: ') + if ( .not. item%isConst ) then + call CreateTimeInterval(item,clock,__RC__) + end if + item%pfioCollection_id = MAPL_DataAddCollection(item%file,use_file_coords=self%use_file_coords) ! parse refresh template to see if we have a time shift during constant updating @@ -1889,6 +1876,7 @@ subroutine CreateTimeInterval(item,clock,rc) character(len=ESMF_MAXSTR) :: creffTime, ctInt integer :: status + logical :: found creffTime = '' ctInt = '' @@ -1933,6 +1921,10 @@ subroutine CreateTimeInterval(item,clock,rc) else ! couldn't find any tokens so all the data must be on one file call ESMF_TimeIntervalSet(item%frequency,__RC__) + + ! check if non-token file exists + inquire(file=trim(item%file),EXIST=found) + _ASSERT(found,'File ' // trim(item%file) // ' not found') end if else ! Reference time should look like: diff --git a/gridcomps/ExtData2G/CMakeLists.txt b/gridcomps/ExtData2G/CMakeLists.txt new file mode 100644 index 000000000000..6efdc8a2d362 --- /dev/null +++ b/gridcomps/ExtData2G/CMakeLists.txt @@ -0,0 +1,32 @@ +esma_set_this (OVERRIDE MAPL.ExtData2G) + +set (srcs + ExtDataFileStream.F90 + ExtDataRule.F90 + ExtDataDerived.F90 + ExtDataConfig.F90 + ExtDataGridCompNG.F90 + TimeStringConversion.F90 + ExtDataTypeDef.F90 + ExtDataOldTypesCreator.F90 + ExtDataBracket.F90 + ExtDataUpdatePointer.F90 + ExtDataAbstractFileHandler.F90 + ExtDataClimFileHandler.F90 + ExtDataSimpleFileHandler.F90 + ExtDataNode.F90 + ExtDataLgr.F90 + ExtDataConstants.F90 + ExtDataSample.F90 + ExtData_IOBundleMod.F90 + ExtData_IOBundleVectorMod.F90 + ) + + +esma_add_library (${this} SRCS ${srcs} DEPENDENCIES MAPL.shared MAPL.base MAPL.generic MAPL.griddedio TYPE SHARED) +target_link_libraries (${this} PUBLIC GFTL::gftl GFTL_SHARED::gftl-shared YAFYAML::yafyaml esmf NetCDF::NetCDF_Fortran + PRIVATE MPI::MPI_Fortran) +target_include_directories (${this} PUBLIC ${INC_ESMF} ${INC_NETCDF} + $) + +set_target_properties (${this} PROPERTIES Fortran_MODULE_DIRECTORY ${include_${this}}) diff --git a/gridcomps/ExtData2G/ExtDataAbstractFileHandler.F90 b/gridcomps/ExtData2G/ExtDataAbstractFileHandler.F90 new file mode 100644 index 000000000000..ec003f7276a6 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataAbstractFileHandler.F90 @@ -0,0 +1,168 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" +#include "unused_dummy.H" +module MAPL_ExtdataAbstractFileHandler + use ESMF + use yafYaml + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_ExtDataBracket + use MAPL_ExtDataFileStream + use MAPL_ExtDataFileStreamMap + use MAPL_DataCollectionMod + use MAPL_CollectionVectorMod + use MAPL_ExtDataConstants + use MAPL_DataCollectionManagerMod + use MAPL_FileMetadataUtilsMod + use MAPL_TimeStringConversion + use MAPL_StringTemplate + implicit none + private + public :: ExtDataAbstractFileHandler + + type, abstract :: ExtDataAbstractFileHandler + character(:), allocatable :: file_template + type(ESMF_TimeInterval) :: frequency + type(ESMF_Time) :: reff_time + integer :: collection_id + type(ESMF_Time), allocatable :: valid_range(:) + logical :: persist_closest + contains + procedure :: initialize + procedure :: make_metadata + procedure :: get_time_on_file + procedure(get_file_bracket), deferred :: get_file_bracket + end type + + abstract interface + subroutine get_file_bracket(this, input_time, source_time, bracket, rc) + use ESMF + use MAPL_ExtDataBracket + import ExtDataAbstractFileHandler + class(ExtDataAbstractFileHandler), intent(inout) :: this + type(ESMF_Time), intent(in) :: input_time + type(ESMF_Time), intent(in) :: source_time(:) + type(ExtDataBracket), intent(inout) :: bracket + integer, optional, intent(out) :: rc + end subroutine get_file_bracket + + end interface + +contains + + subroutine initialize(this,file_series,persist_closest,unusable,rc) + class(ExtDataAbstractFileHandler), intent(inout) :: this + type(ExtDataFileStream), intent(in) :: file_series + class(KeywordEnforcer), optional, intent(in) :: unusable + logical, optional, intent(in) :: persist_closest + integer, optional, intent(out) :: rc + + integer :: status + + _UNUSED_DUMMY(unusable) + + this%file_template = file_series%file_template + this%frequency = file_series%frequency + this%reff_time = file_series%reff_time + if (allocated(file_series%valid_range)) then + allocate(this%valid_range,source=file_series%valid_range) + end if + this%collection_id = file_series%collection_id + if (present(persist_closest)) then + this%persist_closest = persist_closest + else + this%persist_closest = .false. + end if + + end subroutine initialize + + subroutine get_time_on_file(this,filename,target_time,bracketside,time_index,output_time,unusable,wrap,rc) + class(ExtdataAbstractFileHandler), intent(inout) :: this + character(len=*), intent(inout) :: filename + type(ESMF_Time), intent(in) :: target_time + character(len=*), intent(in) :: bracketside + integer, intent(Out) :: time_index + type(ESMF_Time), intent(out) :: output_time + class (KeywordEnforcer), optional, intent(out) :: unusable + integer, optional, intent(inout) :: wrap + integer, optional, intent(out) :: rc + integer :: status + + type(FileMetadataUtils), pointer :: file_metadata + type(ESMF_Time), allocatable :: time_series(:) + logical :: in_bounds, found_time, wrap_ + integer :: i,num_times + + _UNUSED_DUMMY(unusable) + if (present(wrap)) then + wrap_= .true. + else + wrap_=.false. + end if + time_index=time_not_found + + call this%make_metadata(filename,file_metadata,__RC__) + call file_metadata%get_time_info(timeVector=time_series,__RC__) + num_times = size(time_series) + found_time = .false. + if (bracketside == 'L') then + in_bounds = .not.(target_time < time_series(1)) + if (in_bounds) then + do i=num_times,1,-1 + if (target_time >= time_series(i)) then + output_time = time_series(i) + time_index = i + found_time = .true. + exit + end if + enddo + else + if (wrap_) then + output_time=time_series(num_times) + time_index = num_times + found_time = .true. + wrap = -1 + end if + end if + else if (bracketside == 'R') then + in_bounds = .not.(target_time >= time_series(num_times)) + if (in_bounds) then + do i=1,num_times + if (target_time < time_series(i)) then + output_time = time_series(i) + time_index = i + found_time = .true. + exit + end if + enddo + else + if (wrap_) then + output_time=time_series(1) + time_index = 1 + found_time = .true. + wrap = 1 + end if + end if + else + _ASSERT(.false.,"unknown bracket side") + end if + + _RETURN(_SUCCESS) + + end subroutine get_time_on_file + + subroutine make_metadata(this,file,metadata,rc) + class(ExtdataAbstractFileHandler), intent(inout) :: this + character(len=*), intent(in ) :: file + type(FileMetadataUtils), pointer, intent(inout) :: metadata + integer, optional, intent(out ) :: rc + type(MAPLDataCollection), pointer :: collection => null() + + Collection => DataCollections%at(this%collection_id) + metadata => collection%find(file) + _RETURN(_SUCCESS) + + end subroutine make_metadata + + +end module MAPL_ExtdataAbstractFileHandler diff --git a/gridcomps/ExtData2G/ExtDataBracket.F90 b/gridcomps/ExtData2G/ExtDataBracket.F90 new file mode 100644 index 000000000000..d887b73c8f42 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataBracket.F90 @@ -0,0 +1,270 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" +module MAPL_ExtDataBracket + use ESMF + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_BaseMod, only: MAPL_UNDEF + use MAPL_ExtDataNode + implicit none + private + + public :: ExtDataBracket + + type ExtDataBracket + type(ExtDataNode) :: left_node + type(ExtDataNode) :: right_node + real :: scale_factor = 0.0 + real :: offset = 0.0 + logical :: disable_interpolation = .false. + logical :: intermittent_disable = .false. + logical :: new_file_right + logical :: new_file_left + contains + procedure :: interpolate_to_time + procedure :: time_in_bracket + procedure :: set_parameters + procedure :: get_parameters + procedure :: set_node + procedure :: get_node + procedure :: swap_node_fields + procedure :: reset + end type ExtDataBracket + +contains + + subroutine reset(this) + class(ExtDataBracket), intent(inout) :: this + this%new_file_right=.false. + this%new_file_left =.false. + end subroutine reset + + function time_in_bracket(this,time) result(in_bracket) + class(ExtDataBracket), intent(in) :: this + logical :: in_bracket + type(ESMF_Time), intent(in) :: time + + in_bracket = (this%left_node%time <=time) .and. (time < this%right_node%time) + + end function time_in_bracket + + subroutine set_node(this, bracketside, unusable, field, file, time, time_index, was_set, rc) + class(ExtDataBracket), intent(inout) :: this + character(len=*), intent(in) :: bracketside + class(KeywordEnforcer), optional, intent(in) :: unusable + type(ESMF_Field), optional, intent(in) :: field + character(len=*), optional, intent(in) :: file + integer, optional, intent(in) :: time_index + type(ESMF_Time), optional, intent(in) :: time + logical, optional, intent(in) :: was_set + integer, optional, intent(out) :: rc + + _UNUSED_DUMMY(unusable) + if (bracketside=='L') then + if (present(field)) this%left_node%field=field + if (present(time)) this%left_node%time=time + if (present(time_index)) this%left_node%time_index=time_index + if (present(file)) this%left_node%file=file + if (present(was_set)) this%left_node%was_set=was_set + else if (bracketside=='R') then + if (present(field)) this%right_node%field=field + if (present(time)) this%right_node%time=time + if (present(time_index)) this%right_node%time_index=time_index + if (present(file)) this%right_node%file=file + if (present(was_set)) this%right_node%was_set=was_set + else + _ASSERT(.false.,'wrong bracket side') + end if + _RETURN(_SUCCESS) + + end subroutine set_node + + subroutine get_node(this, bracketside, unusable, field, file, time, time_index, was_set, rc) + class(ExtDataBracket), intent(inout) :: this + character(len=*), intent(in) :: bracketside + class(KeywordEnforcer), optional, intent(in) :: unusable + type(ESMF_Field), optional, intent(out) :: field + character(len=*), optional, intent(out) :: file + integer, optional, intent(out) :: time_index + type(ESMF_Time), optional, intent(out) :: time + logical, optional, intent(out) :: was_set + integer, optional, intent(out) :: rc + + _UNUSED_DUMMY(unusable) + if (bracketside=='L') then + if (present(field)) field=this%left_node%field + if (present(time)) time=this%left_node%time + if (present(time_index)) time_index=this%left_node%time_index + if (present(file)) file=this%left_node%file + if (present(was_set)) was_set=this%left_node%was_set + else if (bracketside=='R') then + if (present(field)) field=this%right_node%field + if (present(time)) time=this%right_node%time + if (present(time_index)) time_index=this%right_node%time_index + if (present(file)) file=this%right_node%file + if (present(was_set)) was_set=this%right_node%was_set + else + _ASSERT(.false.,'wrong bracket side') + end if + _RETURN(_SUCCESS) + + end subroutine get_node + + + subroutine set_parameters(this, unusable, linear_trans, disable_interpolation, left_field, right_field, intermittent_disable, rc) + class(ExtDataBracket), intent(inout) :: this + class(KeywordEnforcer), optional, intent(in) :: unusable + real, optional, intent(in) :: linear_trans(2) + logical, optional, intent(in) :: disable_interpolation + type(ESMF_Field), optional, intent(in) :: left_field + type(ESMF_Field), optional, intent(in) :: right_field + logical, optional, intent(in) :: intermittent_disable + integer, optional, intent(out) :: rc + + _UNUSED_DUMMY(unusable) + if (present(linear_trans)) then + this%offset=linear_trans(1) + this%scale_factor=linear_trans(2) + end if + if (present(disable_interpolation)) this%disable_interpolation = disable_interpolation + if (present(left_field)) this%left_node%field=left_field + if (present(right_field)) this%right_node%field=right_field + if (present(intermittent_disable)) this%intermittent_disable = intermittent_disable + _RETURN(_SUCCESS) + + end subroutine set_parameters + + subroutine get_parameters(this, bracket_side, unusable, field, file, time, time_index, update, rc) + class(ExtDataBracket), intent(inout) :: this + character(len=*), intent(in) :: bracket_side + class(KeywordEnforcer), optional, intent(in) :: unusable + type(ESMF_Field), optional, intent(out) :: field + character(len=*), optional, intent(out) :: file + type(ESMF_Time), optional, intent(out) :: time + integer, optional, intent(out) :: time_index + logical, optional, intent(out) :: update + integer, optional, intent(out) :: rc + + _UNUSED_DUMMY(unusable) + if (bracket_side == 'L') then + if (present(field)) field = this%left_node%field + if (present(file)) file = trim(this%left_node%file) + if (present(time)) time = this%left_node%time + if (present(time_index)) time_index = this%left_node%time_index + if (present(update)) update = this%new_file_left + else if (bracket_side == 'R') then + if (present(field)) field = this%right_node%field + if (present(file)) file = trim(this%right_node%file) + if (present(time)) time = this%right_node%time + if (present(time_index)) time_index = this%right_node%time_index + if (present(update)) update = this%new_file_right + else + _ASSERT(.false.,'invalid bracket side!') + end if + _RETURN(_SUCCESS) + + end subroutine get_parameters + + subroutine interpolate_to_time(this,field,time,rc) + class(ExtDataBracket), intent(inout) :: this + type(ESMF_Field), intent(inout) :: field + type(ESMF_Time), intent(in) :: time + integer, optional, intent(out) :: rc + + type(ESMF_TimeInterval) :: tinv1, tinv2 + real :: alpha + real, pointer :: var2d(:,:) => null() + real, pointer :: var3d(:,:,:) => null() + real, pointer :: var2d_left(:,:) => null() + real, pointer :: var2d_right(:,:) => null() + real, pointer :: var3d_left(:,:,:) => null() + real, pointer :: var3d_right(:,:,:) => null() + integer :: field_rank + integer :: status + + call ESMF_FieldGet(field,dimCount=field_rank,__RC__) + alpha = 0.0 + if ( (.not.this%disable_interpolation) .and. (.not.this%intermittent_disable)) then + tinv1 = time - this%left_node%time + tinv2 = this%right_node%time - this%left_node%time + alpha = tinv1/tinv2 + end if + if (field_rank==2) then + call ESMF_FieldGet(field,localDE=0,farrayPtr=var2d,__RC__) + call ESMF_FieldGet(this%right_node%field,localDE=0,farrayPtr=var2d_right,__RC__) + call ESMF_FieldGet(this%left_node%field,localDE=0,farrayPtr=var2d_left,__RC__) + if (time == this%left_node%time .or. this%disable_interpolation) then + var2d = var2d_left + else if (time == this%right_node%time) then + var2d = var2d_right + else + where( (var2d_left /= MAPL_UNDEF) .and. (var2d_right /= MAPL_UNDEF)) + var2d = var2d_left + alpha*(var2d_right-var2d_left) + elsewhere + var2d = MAPL_UNDEF + endwhere + end if + + if (this%scale_factor == 0.0 .and. this%offset /= 0.0) then + where(var2d /= MAPL_UNDEF) var2d=var2d+this%offset + end if + if (this%scale_factor /= 0.0 .and. this%offset == 0.0) then + where(var2d /= MAPL_UNDEF) var2d=var2d*this%scale_factor + end if + if (this%scale_factor /= 0.0 .and. this%offset /= 0.0) then + where(var2d /= MAPL_UNDEF) var2d=var2d*this%scale_factor+this%offset + end if + + else if (field_rank==3) then + call ESMF_FieldGet(field,localDE=0,farrayPtr=var3d,__RC__) + call ESMF_FieldGet(this%right_node%field,localDE=0,farrayPtr=var3d_right,__RC__) + call ESMF_FieldGet(this%left_node%field,localDE=0,farrayPtr=var3d_left,__RC__) + if (time == this%left_node%time .or. this%disable_interpolation) then + var3d = var3d_left + else if (time == this%right_node%time) then + var3d = var3d_right + else + where( (var3d_left /= MAPL_UNDEF) .and. (var3d_right /= MAPL_UNDEF)) + var3d = var3d_left + alpha*(var3d_right-var3d_left) + elsewhere + var3d = MAPL_UNDEF + endwhere + end if + + if (this%scale_factor == 0.0 .and. this%offset /= 0.0) then + where(var3d /= MAPL_UNDEF) var3d=var3d+this%offset + end if + if (this%scale_factor /= 0.0 .and. this%offset == 0.0) then + where(var3d /= MAPL_UNDEF) var3d=var3d*this%scale_factor + end if + if (this%scale_factor /= 0.0 .and. this%offset /= 0.0) then + where(var3d /= MAPL_UNDEF) var3d=var3d*this%scale_factor+this%offset + end if + + end if + _RETURN(_SUCCESS) + + end subroutine interpolate_to_time + + subroutine swap_node_fields(this,rc) + class(ExtDataBracket), intent(inout) :: this + integer, optional, intent(out) :: rc + integer :: status + integer :: field_rank + real, pointer :: var3d_left(:,:,:),var3d_right(:,:,:) + real, pointer :: var2d_left(:,:),var2d_right(:,:) + + call ESMF_FieldGet(this%left_node%field,dimCount=field_rank,__RC__) + if (field_rank == 2) then + call ESMF_FieldGet(this%right_node%field,localDE=0,farrayPtr=var2d_right,__RC__) + call ESMF_FieldGet(this%left_node%field,localDE=0,farrayPtr=var2d_left,__RC__) + var2d_left = var2d_right + else if (field_rank ==3) then + call ESMF_FieldGet(this%right_node%field,localDE=0,farrayPtr=var3d_right,__RC__) + call ESMF_FieldGet(this%left_node%field,localDE=0,farrayPtr=var3d_left,__RC__) + var3d_left = var3d_right + end if + _RETURN(_SUCCESS) + end subroutine swap_node_fields + +end module MAPL_ExtDataBracket diff --git a/gridcomps/ExtData2G/ExtDataClimFileHandler.F90 b/gridcomps/ExtData2G/ExtDataClimFileHandler.F90 new file mode 100644 index 000000000000..8dc2619aae33 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataClimFileHandler.F90 @@ -0,0 +1,280 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" +module MAPL_ExtdataClimFileHandler + use ESMF + use MAPL_ExtDataAbstractFileHandler + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_DataCollectionMod + use MAPL_CollectionVectorMod + use MAPL_DataCollectionManagerMod + use MAPL_FileMetadataUtilsMod + use MAPL_TimeStringConversion + use MAPL_StringTemplate + use MAPL_ExtDataBracket + use MAPL_ExtDataConstants + implicit none + private + public ExtDataClimFileHandler + + integer, parameter :: CLIM_NULL = -1 + type, extends(ExtDataAbstractFileHandler) :: ExtDataClimFileHandler + integer :: clim_year = CLIM_NULL + contains + procedure :: get_file_bracket + procedure :: get_file + end type + +contains + + subroutine get_file_bracket(this, input_time, source_time, bracket, rc) + class(ExtdataClimFileHandler), intent(inout) :: this + type(ESMF_Time), intent(in) :: input_time + type(ESMF_Time), intent(in) :: source_time(:) + type(ExtDataBracket), intent(inout) :: bracket + integer, optional, intent(out) :: rc + + type(ESMF_Time) :: time + integer :: time_index + character(len=ESMF_MAXPATHLEN) :: current_file + integer :: status + type(ESMF_TimeInterval) :: zero + type(ESMF_Time) :: target_time + + integer :: target_year, original_year,clim_shift,valid_years(2) + integer, allocatable :: source_years(:) + + + if (bracket%time_in_bracket(input_time)) then + _RETURN(_SUCCESS) + end if + + target_time=input_time + _ASSERT(size(this%valid_range) == 2, 'Valid time is not defined so can not do any extrapolation or climatology') + call ESMF_TimeGet(this%valid_range(1),yy=valid_years(1),__RC__) + call ESMF_TimeGet(this%valid_range(2),yy=valid_years(2),__RC__) + if (size(source_time)==2) then + allocate(source_years(2)) + call ESMF_TimeGet(source_time(1),yy=source_years(1),__RC__) + call ESMF_TimeGet(source_time(2),yy=source_years(2),__RC__) + _ASSERT(source_years(1) >= valid_years(1),'source time outide valid range') + _ASSERT(source_years(1) <= valid_years(2),'source time outide valid range') + _ASSERT(source_years(2) >= valid_years(1),'source time outide valid range') + _ASSERT(source_years(2) <= valid_years(2),'source time outide valid range') + end if + + ! shift target year to request source time if specified + ! is TS1 < TM < TS2, if not then extrapolate beyond that + call ESMF_TimeGet(target_time,yy=target_year,__RC__) + original_year=target_year + + !if (size(source_years)>0) then + if (allocated(source_years)) then + if (target_year < source_years(1)) then + target_year = source_years(1) + this%clim_year = target_year + else if (target_year >= source_years(2)) then + target_year = source_years(2) + this%clim_year = target_year + end if + call swap_year(target_time,target_year,__RC__) + else + if (target_year < valid_years(1)) then + target_year = valid_years(1) + this%clim_year = target_year + call swap_year(target_time,target_year,__RC__) + else if (target_year >= valid_years(2)) then + target_year = valid_years(2) + this%clim_year = target_year + call swap_year(target_time,target_year,__RC__) + end if + end if + + ! the target time is contained in the dataset and we are not extrapolating outside of source time selection based on available data + if (this%clim_year == CLIM_NULL) then + + call ESMF_TimeIntervalSet(zero,__RC__) + if (this%frequency == zero) then + current_file = this%file_template + call this%get_time_on_file(current_file,input_time,'L',time_index,time,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found on file") + call bracket%set_node('L',file=current_file,time_index=time_index,time=time,__RC__) + if (bracket%left_node == bracket%right_node) then + call bracket%swap_node_fields(rc=status) + _VERIFY(status) + else + bracket%new_file_left=.true. + end if + call this%get_time_on_file(current_file,input_time,'R',time_index,time,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found on file") + call bracket%set_node('R',file=current_file,time_index=time_index,time=time,__RC__) + bracket%new_file_right=.true. + else + call this%get_file(current_file,target_time,0,__RC__) + call this%get_time_on_file(current_file,target_time,'L',time_index,time,rc=status) + if (time_index == time_not_found) then + call this%get_file(current_file,target_time,-1,__RC__) + call this%get_time_on_file(current_file,target_time,'L',time_index,time,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found on file") + end if + call bracket%set_node('L',file=current_file,time_index=time_index,time=time,__RC__) + if (bracket%left_node == bracket%right_node) then + call bracket%swap_node_fields(rc=status) + _VERIFY(status) + else + bracket%new_file_left=.true. + end if + + call this%get_file(current_file,target_time,0,__RC__) + call this%get_time_on_file(current_file,target_time,'R',time_index,time,rc=status) + if (time_index == time_not_found) then + call this%get_file(current_file,target_time,1,__RC__) + call this%get_time_on_file(current_file,target_time,'R',time_index,time,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found on file") + end if + call bracket%set_node('R',file=current_file,time_index=time_index,time=time,__RC__) + bracket%new_file_right=.true. + end if + + ! the target time has been specified to be a climatology for the year; either we + ! are outside the dataset or we have requested a source time range and are on + ! or outside either end + else + + call ESMF_TimeIntervalSet(zero,__RC__) + if (this%frequency == zero) then + current_file = this%file_template + clim_shift=0 + call this%get_time_on_file(current_file,target_time,'L',time_index,time,wrap=clim_shift,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found on file") + call swap_year(time,original_year+clim_shift,__RC__) + call bracket%set_node('L',file=current_file,time_index=time_index,time=time,__RC__) + if (bracket%left_node == bracket%right_node) then + call bracket%swap_node_fields(rc=status) + _VERIFY(status) + else + bracket%new_file_left=.true. + end if + + clim_shift=0 + call this%get_time_on_file(current_file,target_time,'R',time_index,time,wrap=clim_shift,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found on file") + call swap_year(time,original_year+clim_shift,__RC__) + call bracket%set_node('R',file=current_file,time_index=time_index,time=time,__RC__) + bracket%new_file_right=.true. + + else + + call this%get_file(current_file,target_time,0,__RC__) + call this%get_time_on_file(current_file,target_time,'L',time_index,time,rc=status) + if (time_index == time_not_found) then + call this%get_file(current_file,target_time,-1,__RC__) + call this%get_time_on_file(current_file,target_time,'L',time_index,time,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found on file") + call ESMF_TimeGet(target_time,yy=target_year,__RC__) + if (target_year > this%clim_year) then + call swap_year(time,original_year-1,__RC__) + else + call swap_year(time,original_year,__RC__) + end if + else + call swap_year(time,original_year,__RC__) + end if + if (bracket%left_node == bracket%right_node) then + call bracket%swap_node_fields(rc=status) + _VERIFY(status) + else + bracket%new_file_left=.true. + end if + call bracket%set_node('L',file=current_file,time_index=time_index,time=time,__RC__) + + call this%get_file(current_file,target_time,0,__RC__) + call this%get_time_on_file(current_file,target_time,'R',time_index,time,rc=status) + if (time_index == time_not_found) then + call this%get_file(current_file,target_time,1,__RC__) + call this%get_time_on_file(current_file,target_time,'R',time_index,time,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found on file") + call ESMF_TimeGet(target_time,yy=target_year,__RC__) + if (target_year < this%clim_year) then + call swap_year(time,original_year+1,__RC__) + else + call swap_year(time,original_year,__RC__) + end if + else + call swap_year(time,original_year,__RC__) + end if + call bracket%set_node('R',file=current_file,time_index=time_index,time=time,__RC__) + bracket%new_file_right=.true. + + end if + + end if + + _RETURN(_SUCCESS) + + end subroutine get_file_bracket + + subroutine get_file(this,filename,target_time,shift,rc) + class(ExtdataClimFileHandler), intent(inout) :: this + character(len=*), intent(out) :: filename + type(ESMF_Time) :: target_time + integer, intent(in) :: shift + integer, intent(out), optional :: rc + + type(ESMF_Time) :: ftime + integer :: n,status + logical :: file_found + integer :: new_year + integer(ESMF_KIND_I8) :: interval_seconds + + + call ESMF_TimeIntervalGet(this%frequency,s_i8=interval_seconds) + if (interval_seconds==0) then + ! time is not representable as absolute time interval (month, year etc...) do this + ! brute force way. Not good but ESMF leaves no choice + ftime=this%reff_time + do while (ftime < target_time) + ftime = ftime + this%frequency + enddo + ftime=ftime -this%frequency + shift*this%frequency + else + n = (target_time-this%reff_time)/this%frequency + ftime = this%reff_time+(n+shift)*this%frequency + end if + if (this%clim_year /= CLIM_NULL) then + call ESMF_TimeGet(ftime,yy=new_year,__RC__) + if (new_year/=this%clim_year) then + call swap_year(ftime,this%clim_year,__RC__) + if (shift > 0) then + call swap_year(target_time,this%clim_year-shift) + else if (shift < 0) then + call swap_year(target_time,this%clim_year+shift) + end if + end if + end if + call fill_grads_template(filename,this%file_template,time=ftime,__RC__) + inquire(file=trim(filename),exist=file_found) + _ASSERT(file_found,"get_file did not file a file using: "//trim(this%file_template)) + _RETURN(_SUCCESS) + + end subroutine get_file + + subroutine swap_year(time,year,rc) + type(ESMF_Time), intent(inout) :: time + integer, intent(in) :: year + integer, optional, intent(out) :: rc + logical :: is_leap_year + type(ESMF_Calendar) :: calendar + integer :: status, month, day, hour, minute, second + + is_leap_year=.false. + call ESMF_TimeGet(time,mm=month,dd=day,h=hour,m=minute,s=second,calendar=calendar,__RC__) + if (day==29 .and. month==2) then + is_leap_year = ESMF_CalendarIsLeapYear(calendar,year,__RC__) + if (.not.is_leap_year) day=28 + end if + call ESMF_TimeSet(time,yy=year,mm=month,dd=day,h=hour,m=minute,s=second,__RC__) + _RETURN(_SUCCESS) + end subroutine + +end module MAPL_ExtdataClimFileHandler diff --git a/gridcomps/ExtData2G/ExtDataConfig.F90 b/gridcomps/ExtData2G/ExtDataConfig.F90 new file mode 100644 index 000000000000..4f3d0dcc7212 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataConfig.F90 @@ -0,0 +1,200 @@ +#include "MAPL_ErrLog.h" +module MAPL_ExtDataConfig + use ESMF + use yaFyaml + use gFTL_StringVector + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_ExtDataFileStream + use MAPL_ExtDataFileStreamMap + use MAPL_ExtDataRule + use MAPL_ExtDataRuleMap + use MAPL_ExtDataDerived + use MAPL_ExtDataDerivedMap + use MAPL_ExtDataConstants + use MAPL_ExtDataTimeSample + use MAPL_ExtDataTimeSampleMap + implicit none + private + + type, public :: ExtDataConfig + integer :: debug + type(ExtDataRuleMap) :: rule_map + type(ExtDataDerivedMap) :: derived_map + type(ExtDataFileStreamMap) :: file_stream_map + type(ExtDataTimeSampleMap) :: sample_map + + contains + procedure :: get_item_type + procedure :: get_debug_flag + procedure :: new_ExtDataConfig_from_yaml + end type + +contains + + recursive subroutine new_ExtDataConfig_from_yaml(ext_config,config_file,current_time,unusable,rc) + class(ExtDataConfig), intent(inout), target :: ext_config + character(len=*), intent(in) :: config_file + type(ESMF_Time), intent(in) :: current_time + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(Parser) :: p + type(Configuration) :: config, subcfg, ds_config, rule_config, derived_config, sample_config + type(ConfigurationIterator) :: iter + character(len=:), allocatable :: key + type(ExtDataFileStream) :: ds + type(ExtDataDerived) :: derived + type(ExtDataRule) :: rule,ucomp,vcomp + type(ExtDataTimeSample) :: ts + integer :: status, semi_pos + character(len=:), allocatable :: uname,vname + type(FileStream) :: fstream + + type(ExtDataFileStream), pointer :: temp_ds + type(ExtDataTimeSample), pointer :: temp_ts + type(ExtDataRule), pointer :: temp_rule + type(ExtDataDerived), pointer :: temp_derived + + type(Configuration) :: subconfigs + character(len=:), allocatable :: sub_file + integer :: i + + type(ExtDataTimeSample), pointer :: ts_grr + + _UNUSED_DUMMY(unusable) + + p = Parser('core') + fstream=FileStream(config_file) + config = p%load(fstream) + call fstream%close() + + if (config%has("subconfigs")) then + subconfigs = config%at("subconfigs") + _ASSERT(subconfigs%is_sequence(),'subconfigs is not a sequence') + do i=1,subconfigs%size() + sub_file = subconfigs%of(i) + call new_ExtDataConfig_from_yaml(ext_config,sub_file,current_time,rc=status) + _VERIFY(status) + end do + end if + + if (config%has("Samplings")) then + sample_config = config%of("Samplings") + iter = sample_config%begin() + do while (iter /= sample_config%end()) + call iter%get_key(key) + temp_ts => ext_config%sample_map%at(key) + _ASSERT(.not.associated(temp_ts),"defined duplicate named sample key") + call iter%get_value(subcfg) + ts = ExtDataTimeSample(subcfg,_RC) + _VERIFY(status) + call ext_config%sample_map%insert(trim(key),ts) + call iter%next() + enddo + end if + + if (config%has("Collections")) then + ds_config = config%of("Collections") + iter = ds_config%begin() + do while (iter /= ds_config%end()) + call iter%get_key(key) + temp_ds => ext_config%file_stream_map%at(key) + _ASSERT(.not.associated(temp_ds),"defined duplicate named collection") + call iter%get_value(subcfg) + ds = ExtDataFileStream(subcfg,current_time,_RC) + call ext_config%file_stream_map%insert(trim(key),ds) + call iter%next() + enddo + end if + + if (config%has("Exports")) then + rule_config = config%of("Exports") + iter = rule_config%begin() + do while (iter /= rule_config%end()) + call rule%set_defaults(rc=status) + _VERIFY(status) + call iter%get_key(key) + call iter%get_value(subcfg) + rule = ExtDataRule(subcfg,ext_config%sample_map,key,_RC) + semi_pos = index(key,";") + if (semi_pos > 0) then + call rule%split_vector(key,ucomp,vcomp,rc=status) + uname = key(1:semi_pos-1) + vname = key(semi_pos+1:len_trim(key)) + temp_rule => ext_config%rule_map%at(trim(uname)) + _ASSERT(.not.associated(temp_rule),"duplicated export entry key") + call ext_config%rule_map%insert(trim(uname),ucomp) + temp_rule => ext_config%rule_map%at(trim(vname)) + _ASSERT(.not.associated(temp_rule),"duplicated export entry key") + call ext_config%rule_map%insert(trim(vname),vcomp) + else + temp_rule => ext_config%rule_map%at(trim(key)) + _ASSERT(.not.associated(temp_rule),"duplicated export entry key") + call ext_config%rule_map%insert(trim(key),rule) + end if + call iter%next() + enddo + end if + + if (config%has("Derived")) then + derived_config = config%at("Derived") + iter = derived_config%begin() + do while (iter /= derived_config%end()) + call derived%set_defaults(rc=status) + _VERIFY(status) + call iter%get_key(key) + call iter%get_value(subcfg) + derived = ExtDataDerived(subcfg,_RC) + temp_derived => ext_config%derived_map%at(trim(uname)) + _ASSERT(.not.associated(temp_derived),"duplicated derived entry key") + call ext_config%derived_map%insert(trim(key),derived) + call iter%next() + enddo + end if + + if (config%has("debug")) then + call config%get(ext_config%debug,"debug",rc=status) + _VERIFY(status) + end if + ts_grr =>ext_config%sample_map%at('sample_0') + + _RETURN(_SUCCESS) + end subroutine new_ExtDataConfig_from_yaml + + function get_item_type(this,item_name,unusable,rc) result(item_type) + class(ExtDataConfig), intent(inout) :: this + character(len=*), intent(in) :: item_name + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + integer :: item_type + type(ExtDataRule), pointer :: rule + type(ExtDataDerived), pointer :: derived + + _UNUSED_DUMMY(unusable) + item_type=ExtData_not_found + rule => this%rule_map%at(trim(item_name)) + if (associated(rule)) then + if (allocated(rule%vector_component)) then + if (rule%vector_component=='EW') then + item_type=Primary_Type_Vector_comp2 + else if (rule%vector_component=='NS') then + item_type=Primary_Type_Vector_comp1 + end if + else + item_type=Primary_Type_scalar + end if + end if + derived => this%derived_map%at(trim(item_name)) + if (associated(derived)) then + item_type=derived_type + end if + _RETURN(_SUCCESS) + end function get_item_type + + integer function get_debug_flag(this) + class(ExtDataConfig), intent(inout) :: this + get_debug_flag=this%debug + end function get_debug_flag + +end module MAPL_ExtDataConfig diff --git a/gridcomps/ExtData2G/ExtDataConstants.F90 b/gridcomps/ExtData2G/ExtDataConstants.F90 new file mode 100644 index 000000000000..dd711bf12b74 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataConstants.F90 @@ -0,0 +1,12 @@ +module MAPL_ExtDataConstants +implicit none +private + + integer, parameter, public :: ExtData_Not_Found = 0 + integer, parameter, public :: Primary_Type_Scalar = 1 + integer, parameter, public :: Primary_Type_Vector_comp1 = 2 + integer, parameter, public :: Primary_Type_Vector_comp2 = 3 + integer, parameter, public :: Derived_TYpe = 4 + integer, parameter, public :: time_not_found = -1 + +end module MAPL_ExtDataConstants diff --git a/gridcomps/ExtData2G/ExtDataDerived.F90 b/gridcomps/ExtData2G/ExtDataDerived.F90 new file mode 100644 index 000000000000..86cfbe1d70e1 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataDerived.F90 @@ -0,0 +1,90 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" +module MAPL_ExtDataDerived + use yaFyaml + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + implicit none + private + + type, public :: ExtDataDerived + character(:), allocatable :: expression + character(:), allocatable :: sample_key + contains + procedure :: display + procedure :: set_defaults + end type + + interface ExtDataDerived + module procedure new_ExtDataDerived + end interface + +contains + + function new_ExtDataDerived(config,unusable,rc) result(rule) + type(Configuration), intent(in) :: config + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(ExtDataDerived) :: rule + logical :: is_present + integer :: status + character(len=:), allocatable :: tempc + _UNUSED_DUMMY(unusable) + + + if (allocated(tempc)) deallocate(tempc) + is_present = config%has("function") + _ASSERT(is_present,"no expression found in derived entry") + call config%get(tempc,"function",rc=status) + _VERIFY(status) + rule%expression=tempc + + if (allocated(tempc)) deallocate(tempc) + is_present = config%has("sample") + if (is_present) then + call config%get(tempc,"sample",rc=status) + _VERIFY(status) + rule%sample_key=tempc + end if + + _RETURN(_SUCCESS) + end function new_ExtDataDerived + + + subroutine set_defaults(this,unusable,rc) + class(ExtDataDerived), intent(inout), target :: this + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + this%expression='' + _RETURN(_SUCCESS) + _UNUSED_DUMMY(unusable) + end subroutine set_defaults + + subroutine display(this) + class(ExtDataDerived) :: this + write(*,*)"function: ",trim(this%expression) + end subroutine display + +end module MAPL_ExtDataDerived + +module MAPL_ExtDataDerivedMap + use MAPL_ExtDataDerived + +#include "types/key_deferredLengthString.inc" +#define _value type(ExtDataDerived) +#define _alt + +#define _map ExtDataDerivedMap +#define _iterator ExtDataDerivedMapIterator + +#include "templates/map.inc" + +#undef _iterator +#undef _map + +#undef _alt +#undef _value + +end module MAPL_ExtDataDerivedMap diff --git a/gridcomps/ExtData2G/ExtDataFileStream.F90 b/gridcomps/ExtData2G/ExtDataFileStream.F90 new file mode 100644 index 000000000000..bee7c4208ab5 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataFileStream.F90 @@ -0,0 +1,204 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" +module MAPL_ExtDataFileStream + use ESMF + use yaFyaml + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_TimeStringConversion + use MAPL_DataCollectionMod + use MAPL_CollectionVectorMod + use MAPL_DataCollectionManagerMod + use MAPL_FileMetadataUtilsMod + use MAPL_StringTemplate + implicit none + private + + type, public :: ExtDataFileStream + character(:), allocatable :: file_template + type(ESMF_TimeInterval) :: frequency + type(ESMF_Time) :: reff_time + integer :: collection_id + type(ESMF_Time), allocatable :: valid_range(:) + type(FileMetaData) :: metadata + contains + procedure :: detect_metadata + end type + + interface ExtDataFileStream + module procedure new_ExtDataFileStream + end interface ExtDataFileStream +contains + + function new_ExtDataFileStream(config,current_time,unusable,rc) result(data_set) + type(Configuration), intent(in) :: config + type(ESMF_Time), intent(in) :: current_time + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(ExtDataFileStream) :: data_set + integer :: status + integer :: last_token + integer :: iyy,imm,idd,ihh,imn,isc,idx + character(len=2) :: token + character(len=:), allocatable :: file_frequency, file_reff_time,range_str + logical :: is_present + + _UNUSED_DUMMY(unusable) + + if (config%is_scalar()) then + + else if (config%is_mapping()) then + is_present = config%has("template") + _ASSERT(is_present,"no file template in the collection") + if (is_present) then + call config%get(data_set%file_template,"template",rc=status) + _VERIFY(status) + file_frequency = get_string_with_default(config,"freq") + file_reff_time = get_string_with_default(config,"ref_time") + range_str = get_string_with_default(config,"valid_range") + end if + end if + + if (file_frequency /= '') then + data_set%frequency = string_to_esmf_timeinterval(file_frequency) + else + last_token = index(data_set%file_template,'%',back=.true.) + if (last_token.gt.0) then + token = data_set%file_template(last_token+1:last_token+2) + select case(token) + case("y4") + call ESMF_TimeIntervalSet(data_set%frequency,yy=1,__RC__) + case("m2") + call ESMF_TimeIntervalSet(data_set%frequency,mm=1,__RC__) + case("d2") + call ESMF_TimeIntervalSet(data_set%frequency,d=1,__RC__) + case("h2") + call ESMF_TimeIntervalSet(data_set%frequency,h=1,__RC__) + case("n2") + call ESMF_TimeIntervalSet(data_set%frequency,m=1,__RC__) + end select + else + ! couldn't find any tokens so all the data must be on one file + call ESMF_TimeIntervalSet(data_set%frequency,__RC__) + end if + end if + + if (file_reff_time /= '') then + data_set%reff_time = string_to_esmf_time(file_reff_time) + else + last_token = index(data_set%file_template,'%',back=.true.) + if (last_token.gt.0) then + call ESMF_TimeGet(current_time, yy=iyy, mm=imm, dd=idd,h=ihh, m=imn, s=isc ,__RC__) + token = data_set%file_template(last_token+1:last_token+2) + select case(token) + case("y4") + call ESMF_TimeSet(data_set%reff_time,yy=iyy,mm=1,dd=1,h=0,m=0,s=0,__RC__) + case("m2") + call ESMF_TimeSet(data_set%reff_time,yy=iyy,mm=imm,dd=1,h=0,m=0,s=0,__RC__) + case("d2") + call ESMF_TimeSet(data_set%reff_time,yy=iyy,mm=imm,dd=idd,h=0,m=0,s=0,__RC__) + case("h2") + call ESMF_TimeSet(data_set%reff_time,yy=iyy,mm=imm,dd=idd,h=ihh,m=0,s=0,__RC__) + case("n2") + call ESMF_TimeSet(data_set%reff_time,yy=iyy,mm=imm,dd=idd,h=ihh,m=imn,s=0,__RC__) + end select + else + data_set%reff_time = current_time + end if + end if + + if (range_str /= '') then + idx = index(range_str,',') + _ASSERT(idx/=0,'invalid specification of time range') + if (allocated(data_set%valid_range)) deallocate(data_set%valid_range) + allocate(data_set%valid_range(2)) + data_set%valid_range(1)=string_to_esmf_time(range_str(:idx-1)) + data_set%valid_range(2)=string_to_esmf_time(range_str(idx+1:)) + call ESMF_TimeGet(data_set%reff_time,yy=iyy,mm=imm,dd=idd,h=ihh,m=imn,__RC__) + call ESMF_TimeGet(data_set%valid_range(1),yy=iyy,__RC__) + call ESMF_TimeSet(data_set%reff_time,yy=iyy,mm=imm,dd=idd,h=ihh,m=imn,__RC__) + end if + data_set%collection_id = MAPL_DataAddCollection(data_set%file_template) + + _RETURN(_SUCCESS) + + contains + + function get_string_with_default(config,selector) result(string) + type(Configuration), intent(in) :: config + character(len=*), intent(In) :: selector + character(len=:), allocatable :: string + + if (config%has(selector)) then + string=config%of(selector) + else + string='' + end if + end function + + end function new_ExtDataFileStream + + subroutine detect_metadata(this,metadata_out,time,get_range,rc) + class(ExtDataFileStream), intent(inout) :: this + type(FileMetadataUtils), intent(inout) :: metadata_out + type(ESMF_Time), intent(in) :: time + logical, optional, intent(in) :: get_range + integer, optional, intent(out) :: rc + + logical :: get_range_ + type(MAPLDataCollection), pointer :: collection + type(FileMetadataUtils), pointer :: metadata + type(ESMF_Time), allocatable :: time_series(:) + integer :: status + character(len=ESMF_MAXPATHLEN) :: filename + + if (present(get_range)) then + get_range_ = get_range + else + get_range_ = .false. + end if + + collection => DataCollections%at(this%collection_id) + if (get_range_ .and. (.not.allocated(this%valid_range))) then + if (index('%',this%file_template) == 0) then + metadata => collection%find(this%file_template) + call metadata%get_time_info(timeVector=time_series,__RC__) + allocate(this%valid_range(2)) + this%valid_range(1)=time_series(1) + this%valid_range(2)=time_series(size(time_series)) + end if + end if + + if (get_range_) then + call fill_grads_template(filename,this%file_template,time=this%valid_range(1),__RC__) + else + call fill_grads_template(filename,this%file_template,time=time,__RC__) + end if + metadata => collection%find(filename,__RC__) + metadata_out = metadata + _RETURN(_SUCCESS) + + end subroutine detect_metadata + +end module MAPL_ExtDataFileStream + +module MAPL_ExtDataFileStreamMap + use MAPL_ExtDataFileStream + +#include "types/key_deferredLengthString.inc" +#define _value type(ExtDataFileStream) +#define _alt + +#define _map ExtDataFileStreamMap +#define _iterator ExtDataFileStreamMapIterator + +#include "templates/map.inc" + +#undef _iterator +#undef _map + +#undef _alt +#undef _value + +end module MAPL_ExtDataFileStreamMap diff --git a/gridcomps/ExtData2G/ExtDataGridCompNG.F90 b/gridcomps/ExtData2G/ExtDataGridCompNG.F90 new file mode 100644 index 000000000000..5122856c5a77 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataGridCompNG.F90 @@ -0,0 +1,2288 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_Generic.h" +#include "unused_dummy.H" + +!------------------------------------------------------------------------- +! NASA/GSFC, Global Modeling and Assimilation Office, Code 910.1 ! +!------------------------------------------------------------------------- + + MODULE MAPL_ExtDataGridComp2G + +!BOP +! !MODULE: MAPL_ExtDataGridCompMod - Implements Interface to External Data +! +! !DESCRIPTION: +! +! {\tt MAPL\_ExtDataGridComp} is an ESMF gridded component implementing +! an interface to boundary conditions and other types of external data +! files. +! +! Developed for GEOS-5 release Fortuna 2.0 and later. +! +! !USES: +! + USE ESMF + use gFTL_StringVector + use MAPL_BaseMod + use MAPL_CommsMod + use MAPL_ShmemMod + use ESMFL_Mod + use MAPL_GenericMod + use MAPL_VarSpecMod + use MAPL_CFIOMod + use MAPL_NewArthParserMod + use MAPL_ConstantsMod, only: MAPL_PI,MAPL_PI_R8,MAPL_RADIANS_TO_DEGREES + use MAPL_IOMod, only: MAPL_NCIOParseTimeUnits + use, intrinsic :: iso_fortran_env, only: REAL64 + use linearVerticalInterpolation_mod + use ESMF_CFIOCollectionVectorMod + use ESMF_CFIOCollectionMod + use MAPL_ConfigMod + use MAPL_GridManagerMod + use MAPL_ExtDataNG_IOBundleMod + use MAPL_ExtDataNG_IOBundleVectorMod + use MAPL_ExceptionHandling + use MAPL_DataCollectionMod + use MAPL_CollectionVectorMod + use MAPL_DataCollectionManagerMod + use MAPL_FileMetadataUtilsMod + use pFIO_ClientManagerMod, only : i_Clients + use MAPL_GriddedIOItemMod + use MAPL_GriddedIOItemVectorMod + use MAPL_ExtDataConfig + use MAPL_ExtDataTypeDef + use MAPL_ExtDataOldTypesCreator + use MAPL_StringTemplate + use pflogger, only: logging, Logger + use MAPL_ExtDataLogger + use MAPL_ExtDataConstants + + IMPLICIT NONE + PRIVATE +! +! !PUBLIC MEMBER FUNCTIONS: + + PUBLIC SetServices +!EOP +! +! !REVISION HISTORY: +! +! 12Dec2009 da Silva Design and first implementation. +! +!------------------------------------------------------------------------- + + integer :: Ext_Debug + integer, parameter :: MAPL_ExtDataLeft = 1 + integer, parameter :: MAPL_ExtDataRight = 2 + logical :: hasRun + character(len=ESMF_MAXSTR) :: error_msg_str + + type PrimaryExports + PRIVATE + integer :: nItems = 0 + logical :: have_phis + type(PrimaryExport), pointer :: item(:) => null() + end type PrimaryExports + + type DerivedExports + PRIVATE + integer :: nItems = 0 + type(DerivedExport), pointer :: item(:) => null() + end type DerivedExports + +! Legacy state +! ------------ + type MAPL_ExtData_State + PRIVATE + type(PrimaryExports) :: Primary + type(DerivedExports) :: Derived + ! will add fields from export state to this state + ! will also add new fields that could be mask + ! or primary exports that were not in the export + ! state recieved by ExtData, i.e. fields that are + ! needed by a derived field where the primary fields + ! are not actually required + type(ESMF_State) :: ExtDataState + type(ESMF_Config) :: CF + logical :: active + integer, allocatable :: PrimaryOrder(:) + end type MAPL_ExtData_State + +! Hook for the ESMF +! ----------------- + type MAPL_ExtData_Wrap + type (MAPL_ExtData_State), pointer :: PTR => null() + end type MAPL_ExtData_WRAP + +CONTAINS + + +!------------------------------------------------------------------------- +! NASA/GSFC, Global Modeling and Assimilation Office, Code 910.1 ! +!------------------------------------------------------------------------- +!BOP +! +! !IROUTINE: SetServices --- Sets IRF services for the MAPL_ExtData +! +! !INTERFACE: + + SUBROUTINE SetServices ( GC, RC ) + +! !ARGUMENTS: + + type(ESMF_GridComp), intent(INOUT) :: GC ! gridded component + integer, optional :: RC ! return code + +! !DESCRIPTION: Sets Initialize, Run and Finalize services. +! +! !REVISION HISTORY: +! +! 12Dec2009 da Silva Design and first implementation. +! +!EOP +!------------------------------------------------------------------------- + +! Local derived type aliases +! -------------------------- + type (MAPL_ExtData_State), pointer :: self ! internal, that is + type (MAPL_ExtData_wrap) :: wrap + + character(len=ESMF_MAXSTR) :: comp_name + character(len=ESMF_MAXSTR) :: Iam + integer :: status + +! ------------ + +! Get my name and set-up traceback handle +! --------------------------------------- + Iam = 'SetServices' + call ESMF_GridCompGet( GC, name=comp_name, __RC__ ) + Iam = trim(comp_name) // '::' // trim(Iam) + +! Wrap internal state for storing in GC; rename legacyState +! ------------------------------------- + allocate ( self, stat=STATUS ) + _VERIFY(STATUS) + wrap%ptr => self + +! ------------------------ +! ESMF Functional Services +! ------------------------ + +! Set the Initialize, Run, Finalize entry points +! ---------------------------------------------- + call MAPL_GridCompSetEntryPoint ( GC, ESMF_METHOD_INITIALIZE, Initialize_, __RC__ ) + call MAPL_GridCompSetEntryPoint ( GC, ESMF_METHOD_RUN, Run_, __RC__ ) + call MAPL_GridCompSetEntryPoint ( GC, ESMF_METHOD_FINALIZE, Finalize_, __RC__ ) + +! Store internal state in GC +! -------------------------- + call ESMF_UserCompSetInternalState ( GC, 'MAPL_ExtData_state', wrap, STATUS ) + _VERIFY(STATUS) + + call MAPL_TimerAdd(gc,name="Initialize", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="Run", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="-Read_Loop", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="--CheckUpd", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="--Read", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="--GridCreate", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="--IclientWait", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="--PRead", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="---CreateCFIO", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="---prefetch", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="----add-collection", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="----make-reference", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="----RegridStore", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="----request", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="---IclientDone", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="----RegridApply", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="---read-prefetch", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="--Swap", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="--Bracket", rc=status) + _VERIFY(STATUS) + call MAPL_TimerAdd(gc,name="-Interpolate", rc=status) + _VERIFY(STATUS) +! Generic Set Services +! -------------------- + call MAPL_GenericSetServices ( GC, __RC__ ) + +! All done +! -------- + + _RETURN(ESMF_SUCCESS) + + END SUBROUTINE SetServices + + +!------------------------------------------------------------------------- +! NASA/GSFC, Global Modeling and Assimilation Office, Code 610.1 ! +!------------------------------------------------------------------------- +!BOP +! +! !IROUTINE: Initialize_ --- Initialize MAPL_ExtData +! +! !INTERFACE: +! + + SUBROUTINE Initialize_ ( GC, IMPORT, EXPORT, CLOCK, rc ) + +! !USES: + + implicit NONE + +! !INPUT PARAMETERS: + + type(ESMF_Clock), intent(inout) :: CLOCK ! The clock + +! !OUTPUT PARAMETERS: + + type(ESMF_GridComp), intent(inout) :: GC ! Grid Component + type(ESMF_State), intent(inout) :: IMPORT ! Import State + type(ESMF_State), intent(inout) :: EXPORT ! Export State + integer, intent(out) :: rc ! Error return code: + ! 0 - all is well + ! 1 - + +! !DESCRIPTION: This is a simple ESMF wrapper. +! +! !REVISION HISTORY: +! +! 12Dec2009 da Silva Design and first implementation. +! +!EOP +!------------------------------------------------------------------------- + + type(MAPL_ExtData_state), pointer :: self ! Legacy state + type(ESMF_Grid) :: GRID ! Grid + type(ESMF_Config) :: CF_master ! Universal Config + + character(len=ESMF_MAXSTR) :: comp_name + character(len=ESMF_MAXSTR) :: Iam + integer :: Status + + type(PrimaryExport), pointer :: item + integer :: i + integer :: ItemCount + integer :: PrimaryItemCount, DerivedItemCount + + type(ESMF_Time) :: time + + type (ESMF_Field) :: field,left_field,right_field + integer :: fieldRank, lm + type (ESMF_StateItem_Flag), pointer :: ITEMTYPES(:) + character(len=ESMF_MAXSTR), allocatable :: ITEMNAMES(:) + + real, pointer :: ptr2d(:,:) => null() + real, pointer :: ptr3d(:,:,:) => null() + integer :: idx + type(ESMF_VM) :: vm + type(MAPL_MetaComp),pointer :: MAPLSTATE + + type(ExtDataOldTypesCreator),target :: config_yaml + character(len=:), allocatable :: new_rc_file + logical :: found_in_config + integer :: num_primary,num_derived + integer, allocatable :: item_types(:) + type(StringVector) :: unsatisfied_imports + !class(logger), pointer :: lgr + +! Get my name and set-up traceback handle +! --------------------------------------- + Iam = 'Initialize_' + call ESMF_GridCompGet( GC, name=comp_name, config=CF_master, __RC__ ) + Iam = trim(comp_name) // '::' // trim(Iam) + call MAPL_GetLogger(gc, extdata_lgr, __RC__) + +! Extract relevant runtime information +! ------------------------------------ + call extract_ ( GC, self, CF_master, __RC__) + self%CF = CF_master + +! Start Some Timers +! ----------------- + call MAPL_GetObjectFromGC ( gc, MAPLSTATE, RC=STATUS) + _VERIFY(STATUS) + call MAPL_TimerOn(MAPLSTATE,"TOTAL") + call MAPL_TimerOn(MAPLSTATE,"Initialize") + + call ESMF_ClockGet(CLOCK, currTIME=time, __RC__) + new_rc_file = "extdata.yaml" + config_yaml = ExtDataOldTypesCreator(new_rc_file,time,__RC__) +! Get information from export state +!---------------------------------- + call ESMF_StateGet(EXPORT, ITEMCOUNT=ItemCount, RC=STATUS) + _VERIFY(STATUS) + + ! set ExtData on by default, let user turn it off if they want + call ESMF_ConfigGetAttribute(CF_master,self%active, Label='USE_EXTDATA:',default=.true.,rc=status) + + ! no need to run ExtData if there are no imports to fill + if (ItemCount == 0) then + self%active = .false. + end if + + if (.not.self%active) then + call MAPL_TimerOff(MAPLSTATE,"Initialize") + call MAPL_TimerOff(MAPLSTATE,"TOTAL") + _RETURN(ESMF_SUCCESS) + end if + +! Greetings +! --------- + if (MAPL_am_I_root()) then + print *, TRIM(Iam)//': ACTIVE' + print * + end if + + allocate(ITEMNAMES(ITEMCOUNT), STAT=STATUS) + _VERIFY(STATUS) + allocate(ITEMTYPES(ITEMCOUNT), STAT=STATUS) + _VERIFY(STATUS) + + call ESMF_StateGet(EXPORT, ITEMNAMELIST=ITEMNAMES, & + ITEMTYPELIST=ITEMTYPES, RC=STATUS) + _VERIFY(STATUS) + +! -------- +! Initialize MAPL Generic +! ----------------------- + call MAPL_GenericInitialize ( GC, IMPORT, EXPORT, clock, __RC__ ) + + call extdata_lgr%error("Using ExtData2G, note this is still in BETA stage") + +! --------------------------- +! Parse ExtData Resource File +! --------------------------- + num_primary=0 + num_derived=0 + primaryitemcount=0 + deriveditemcount=0 + allocate(item_types(size(itemnames)),__STAT__) + do i=1,size(itemnames) + item_types(i) = config_yaml%get_item_type(trim(itemnames(i)),rc=status) + _VERIFY(status) + found_in_config = (item_types(i)/= ExtData_not_found) + if (.not.found_in_config) call unsatisfied_imports%push_back(itemnames(i)) + if (item_types(i) == derived_type) then + deriveditemcount=deriveditemcount+1 + else + primaryitemcount=primaryitemcount+1 + end if + enddo + if (unsatisfied_imports%size() > 0) then + do i=1,unsatisfied_imports%size() + call extdata_lgr%error("In ExtData resource file, could not find: "//trim(unsatisfied_imports%at(i))) + enddo + _FAIL("Unsatisfied imports in ExtData") + end if + + ext_debug=config_yaml%get_debug_flag() + allocate(self%primary%item(PrimaryItemCount),__STAT__) + allocate(self%derived%item(DerivedItemCount),__STAT__) + self%primary%nitems = PrimaryItemCount + self%derived%nitems = DerivedItemCount + + self%ExtDataState = ESMF_StateCreate(Name="ExtDataNameSpace",__RC__) + num_primary=0 + num_derived=0 + do i=1,size(itemnames) + if (item_types(i)==Primary_Type_Scalar .or. item_types(i)==Primary_Type_Vector_comp1) then + num_primary=num_primary+1 + call config_yaml%fillin_primary(trim(itemnames(i)),self%primary%item(num_primary),time,clock,__RC__) + else if (item_types(i)==Derived_type) then + num_derived=num_derived+1 + call config_yaml%fillin_derived(trim(itemnames(i)),self%derived%item(num_derived),time,clock,__RC__) + end if + call ESMF_StateGet(Export,trim(itemnames(i)),field,__RC__) + call MAPL_StateAdd(self%ExtDataState,field,__RC__) + enddo +! note: handle case if variables in derived expression need to be allocated! + + PrimaryLoop: do i = 1, self%primary%nItems + + item => self%primary%item(i) + + item%pfioCollection_id = MAPL_DataAddCollection(item%file_template) + +! Read the single step files (read interval equal to zero) +! -------------------------------------------------------- + + if (item%isConst) then + + if (item%vartype == MAPL_FieldItem) then + call ESMF_StateGet(self%ExtDataState,trim(item%name),field,__RC__) + call ESMF_FieldGet(field,dimCount=fieldRank,__RC__) + if (fieldRank == 2) then + call MAPL_GetPointer(self%ExtDataState, ptr2d, trim(item%name),__RC__) + ptr2d = item%const + else if (fieldRank == 3) then + call MAPL_GetPointer(self%ExtDataState, ptr3d, trim(item%name), __RC__) + ptr3d = item%const + endif + else if (item%vartype == MAPL_VectorField) then + call ESMF_StateGet(self%ExtDataState,trim(item%vcomp1),field,__RC__) + call ESMF_FieldGet(field,dimCount=fieldRank,__RC__) + if (fieldRank == 2) then + call MAPL_GetPointer(self%ExtDataState, ptr2d, trim(item%vcomp1),__RC__) + ptr2d = item%const + else if (fieldRank == 3) then + call MAPL_GetPointer(self%ExtDataState, ptr3d, trim(item%vcomp1), __RC__) + ptr3d = item%const + endif + call ESMF_StateGet(self%ExtDataState,trim(item%vcomp2),field,__RC__) + call ESMF_FieldGet(field,dimCount=fieldRank,__RC__) + if (fieldRank == 2) then + call MAPL_GetPointer(self%ExtDataState, ptr2d, trim(item%vcomp2),__RC__) + ptr2d = item%const + else if (fieldRank == 3) then + call MAPL_GetPointer(self%ExtDataState, ptr3d, trim(item%vcomp2), __RC__) + ptr3d = item%const + endif + end if + cycle + end if + + ! get levels, other information + call GetLevs(item,__RC__) + call ESMF_VMBarrier(vm) + ! register collections + item%iclient_collection_id=i_clients%add_ext_collection(trim(item%file_template)) + ! create interpolating fields, check if the vertical levels match the file + if (item%vartype == MAPL_FieldItem) then + + call ESMF_StateGet(self%ExtDataState, trim(item%name), field,__RC__) + call ESMF_FieldGet(field,grid=grid,rank=fieldRank,__RC__) + + lm=0 + if (fieldRank==3) then + call ESMF_FieldGet(field,0,farrayPtr=ptr3d,__RC__) + lm = size(ptr3d,3) + end if + if (item%lm /= lm .and. lm /= 0 .and. item%havePressure) then + item%do_VertInterp = .true. + else if (item%lm /= lm .and. lm /= 0) then + item%do_Fill = .true. + end if + left_field = MAPL_FieldCreate(field,item%var,doCopy=.true.,__RC__) + right_field = MAPL_FieldCreate(field,item%var,doCopy=.true.,__RC__) + call item%modelGridFields%comp1%set_parameters(left_field=left_field,right_field=right_field, __RC__) + if (item%do_fill .or. item%do_vertInterp) then + call createFileLevBracket(item,cf_master,__RC__) + end if + + else if (item%vartype == MAPL_VectorField) then + + ! check that we are not asking for conservative regridding +!!$ if (item%Trans /= MAPL_HorzTransOrderBilinear) then + if (item%Trans /= REGRID_METHOD_BILINEAR) then + _ASSERT(.false.,'No conservative re-gridding with vectors') + end if + + block + integer :: gridRotation1, gridRotation2 + call ESMF_StateGet(self%ExtDataState, trim(item%vcomp1), field,__RC__) + call ESMF_AttributeGet(field, NAME='ROTATION', value=gridRotation1, __RC__) + call ESMF_StateGet(self%ExtDataState, trim(item%vcomp2), field,__RC__) + call ESMF_AttributeGet(field, NAME='ROTATION', value=gridRotation2, __RC__) + _ASSERT(GridRotation1 == gridRotation2,'Grid rotations must match when performing vector re-gridding') + end block + + call ESMF_StateGet(self%ExtDataState, trim(item%vcomp1), field,__RC__) + call ESMF_FieldGet(field,grid=grid,rank=fieldRank,__RC__) + + lm = 0 + if (fieldRank==3) then + call ESMF_FieldGet(field,0,farrayPtr=ptr3d,__RC__) + lm = size(ptr3d,3) + end if + if (item%lm /= lm .and. item%havePressure) then + item%do_VertInterp = .true. + else if (item%lm /= lm .and. lm /= 0) then + item%do_Fill = .true. + end if + + left_field = MAPL_FieldCreate(field,item%fcomp1,doCopy=.true.,__RC__) + right_field = MAPL_FieldCreate(field,item%fcomp1,doCopy=.true.,__RC__) + call item%modelGridFields%comp1%set_parameters(left_field=left_field,right_field=right_field, __RC__) + call ESMF_StateGet(self%ExtDataState, trim(item%vcomp2), field,__RC__) + left_field = MAPL_FieldCreate(field,item%fcomp2,doCopy=.true.,__RC__) + right_field = MAPL_FieldCreate(field,item%fcomp2,doCopy=.true.,__RC__) + call item%modelGridFields%comp2%set_parameters(left_field=left_field,right_field=right_field, __RC__) + + if (item%do_fill .or. item%do_vertInterp) then + call createFileLevBracket(item,cf_master,__RC__) + end if + + end if + + end do PrimaryLoop + +! Check if we have any files that would need to be vertically interpolated +! if so ensure that PS is done first + allocate(self%primaryOrder(size(self%primary%item)),__STAT__) + do i=1,size(self%primary%item) + self%primaryOrder(i)=i + enddo +! check for PS + idx = -1 + if (any(self%primary%item%do_VertInterp .eqv. .true.)) then + do i=1,size(self%primary%item) + if (self%primary%item(i)%name=='PS') then + idx =i + end if + enddo + _ASSERT(idx/=-1,'Surface pressure not present for vertical interpolation') + self%primaryOrder(1)=idx + self%primaryOrder(idx)=1 + self%primary%item(idx)%units = ESMF_UtilStringUppercase(self%primary%item(idx)%units,rc=status) + _ASSERT(trim(self%primary%item(idx)%units)=="PA",'PS must be in units of PA') + end if +! check for PHIS + idx = -1 + if (any(self%primary%item%do_VertInterp .eqv. .true.)) then + do i=1,size(self%primary%item) + if (self%primary%item(i)%name=='PHIS') then + idx =i + end if + enddo + if (idx/=-1) then + self%primaryOrder(2)=idx + self%primaryOrder(idx)=2 + self%primary%have_phis=.true. + end if + end if + + call extdata_lgr%info('*******************************************************') + call extdata_lgr%info('** Variables to be provided by the ExtData Component **') + call extdata_lgr%info('*******************************************************') + do i = 1, ItemCount + call extdata_lgr%info('---- %i0.5~: %a', i, trim(ItemNames(i))) + end do + call extdata_lgr%info('*******************************************************\n') + +! Clean up +! -------- + deallocate(ItemTypes) + deallocate(ItemNames) + +! Set has run to false to we know when we first go to run method it is first call + hasRun = .false. + + call MAPL_TimerOff(MAPLSTATE,"Initialize") + call MAPL_TimerOff(MAPLSTATE,"TOTAL") +! All done +! -------- + + IF ( (Ext_Debug > 0) .AND. MAPL_Am_I_Root() ) THEN + Write(*,*) 'ExtData Initialize_: End' + ENDIF + + _RETURN(ESMF_SUCCESS) + + END SUBROUTINE Initialize_ + + +!------------------------------------------------------------------------- +! NASA/GSFC, Global Modeling and Assimilation Office, Code 610.1 ! +!------------------------------------------------------------------------- +!BOP +! +! !IROUTINE: Run_ --- Runs MAPL_ExtData +! +! !INTERFACE: +! + + SUBROUTINE Run_ ( GC, IMPORT, EXPORT, CLOCK, rc ) + +! !USES: + + implicit NONE + +! !INPUT PARAMETERS: + + type(ESMF_Clock), intent(inout) :: CLOCK ! The clock + +! !OUTPUT PARAMETERS: + + type(ESMF_GridComp), intent(inout) :: GC ! Grid Component + type(ESMF_State), intent(inout) :: IMPORT ! Import State + type(ESMF_State), intent(inout) :: EXPORT ! Export State + integer, intent(out) :: rc ! Error return code: + ! 0 - all is well + ! 1 - + +! !DESCRIPTION: This is a simple ESMF wrapper. +! +! !REVISION HISTORY: +! +! 12Dec2009 da Silva Design and first implementation. +! +!EOP +!------------------------------------------------------------------------- + + type(MAPL_ExtData_state), pointer :: self ! Legacy state + type(ESMF_Config) :: CF ! Universal Config + + character(len=ESMF_MAXSTR) :: comp_name + character(len=ESMF_MAXSTR) :: Iam + integer :: status + + type(PrimaryExport), pointer :: item + type(DerivedExport), pointer :: derivedItem + integer :: i + + type(ESMF_Time) :: time, time0 + type(MAPL_MetaComp), pointer :: MAPLSTATE + + logical :: doUpdate_ + character(len=ESMF_MAXPATHLEN) :: file_processed + logical, allocatable :: doUpdate(:) + type(ESMF_Time), allocatable :: useTime(:) + + integer :: bracket_side + integer :: entry_num + type(IOBundleNGVector), target :: IOBundles + type(IOBundleNGVectorIterator) :: bundle_iter + type(ExtDataNG_IOBundle), pointer :: io_bundle + + _UNUSED_DUMMY(IMPORT) + _UNUSED_DUMMY(EXPORT) + +! Declare pointers to IMPORT/EXPORT/INTERNAL states +! ------------------------------------------------- +! #include "MAPL_ExtData_DeclarePointer___.h" + +! Get my name and set-up traceback handle +! --------------------------------------- + Iam = 'Run_' + call ESMF_GridCompGet( GC, name=comp_name, __RC__ ) + Iam = trim(comp_name) // '::' // trim(Iam) + + +! Call Run for every Child +! ------------------------- +!ALT call MAPL_GenericRunChildren ( GC, IMPORT, EXPORT, CLOCK, __RC__) + + +! Extract relevant runtime information +! ------------------------------------ + call extract_ ( GC, self, CF, __RC__ ) + + if (.not. self%active) then + _RETURN(ESMF_SUCCESS) + end if + + call MAPL_GetObjectFromGC ( gc, MAPLSTATE, RC=STATUS) + _VERIFY(STATUS) + call MAPL_TimerOn(MAPLSTATE,"TOTAL") + call MAPL_TimerOn(MAPLSTATE,"Run") + + call ESMF_ClockGet(CLOCK, currTIME=time0, __RC__) + + +! Fill in the internal state with data from the files +! --------------------------------------------------- + + allocate(doUpdate(self%primary%nitems),stat=status) + _VERIFY(STATUS) + doUpdate = .false. + allocate(useTime(self%primary%nitems),stat=status) + _VERIFY(STATUS) + + call MAPL_TimerOn(MAPLSTATE,"-Read_Loop") + + IF ( (Ext_Debug > 0) .AND. MAPL_Am_I_Root() ) THEN + Write(*,*) 'ExtData Run_: Start' + Write(*,*) 'ExtData Run_: READ_LOOP: Start' + ENDIF + + READ_LOOP: do i = 1, self%primary%nItems + + item => self%primary%item(self%primaryOrder(i)) + + IF ( (Ext_Debug > 0) .AND. MAPL_Am_I_Root() ) THEN + Write(*,*) ' ' + Write(*,'(a,I0.3,a,I0.3,a,a)') 'ExtData Run_: READ_LOOP: variable ', i, ' of ', self%primary%nItems, ': ', trim(item%var) + Write(*,*) ' ==> file: ', trim(item%file_template) + Write(*,*) ' ==> isConst: ', item%isConst + ENDIF + + if (item%isConst) then + IF ( (Ext_Debug > 0) .AND. MAPL_Am_I_Root() ) THEN + Write(*,*) ' ==> Break loop since isConst is true' + ENDIF + cycle + endif + + call MAPL_TimerOn(MAPLSTATE,"--CheckUpd") + + call item%update_freq%check_update(doUpdate(i),time,time0,.not.hasRun,__RC__) + !doUpdate(i) = doUpdate_ .or. (.not.hasRun) + call MAPL_TimerOff(MAPLSTATE,"--CheckUpd") + + DO_UPDATE: if (doUpdate(i)) then + + call item%modelGridFields%comp1%reset() + call item%filestream%get_file_bracket(time,item%source_time, item%modelGridFields%comp1,__RC__) + call IOBundle_Add_Entry(IOBundles,item,self%primaryOrder(i)) + useTime(i)=time + + end if DO_UPDATE + + end do READ_LOOP + + IF ( (Ext_Debug > 0) .AND. MAPL_Am_I_Root() ) THEN + Write(*,*) 'ExtData Run_: READ_LOOP: Done' + ENDIF + + bundle_iter = IOBundles%begin() + do while (bundle_iter /= IoBundles%end()) + io_bundle => bundle_iter%get() + bracket_side = io_bundle%bracket_side + entry_num = io_bundle%entry_index + file_Processed = io_bundle%file_name + item => self%primary%item(entry_num) + + io_bundle%pbundle = ESMF_FieldBundleCreate(rc=status) + _VERIFY(STATUS) + + call MAPL_ExtDataPopulateBundle(item,bracket_side,io_bundle%pbundle,rc=status) + _VERIFY(status) + call bundle_iter%next() + enddo + + call MAPL_TimerOn(MAPLSTATE,"--PRead") + call MAPL_TimerOn(MAPLSTATE,"---CreateCFIO") + call MAPL_ExtDataCreateCFIO(IOBundles, rc=status) + _VERIFY(status) + call MAPL_TimerOff(MAPLSTATE,"---CreateCFIO") + + call MAPL_TimerOn(MAPLSTATE,"---prefetch") + call MAPL_ExtDataPrefetch(IOBundles, rc=status) + _VERIFY(status) + call MAPL_TimerOff(MAPLSTATE,"---prefetch") + _VERIFY(STATUS) + call MAPL_TimerOn(MAPLSTATE,"---IclientDone") + + call i_Clients%done_collective_prefetch() + call i_Clients%wait() + + call MAPL_TimerOff(MAPLSTATE,"---IclientDone") + _VERIFY(STATUS) + + call MAPL_TimerOn(MAPLSTATE,"---read-prefetch") + call MAPL_ExtDataReadPrefetch(IOBundles,rc=status) + _VERIFY(status) + call MAPL_TimerOff(MAPLSTATE,"---read-prefetch") + call MAPL_TimerOff(MAPLSTATE,"--PRead") + + bundle_iter = IOBundles%begin() + do while (bundle_iter /= IOBundles%end()) + io_bundle => bundle_iter%get() + bracket_side = io_bundle%bracket_side + entry_num = io_bundle%entry_index + item => self%primary%item(entry_num) + call MAPL_ExtDataVerticalInterpolate(self,item,bracket_side,rc=status) + _VERIFY(status) + call bundle_iter%next() + enddo + call MAPL_ExtDataDestroyCFIO(IOBundles,rc=status) + _VERIFY(status) + + call MAPL_TimerOff(MAPLSTATE,"-Read_Loop") + + call MAPL_TimerOn(MAPLSTATE,"-Interpolate") + + IF ( (Ext_Debug > 0) .AND. MAPL_Am_I_Root() ) THEN + Write(*,*) 'ExtData Run_: INTERP_LOOP: Start' + ENDIF + + INTERP_LOOP: do i = 1, self%primary%nItems + + item => self%primary%item(self%primaryOrder(i)) + + if (doUpdate(i)) then + + IF ( (Ext_Debug > 0) .AND. MAPL_Am_I_Root() ) THEN + Write(*,*) ' ' + Write(*,'(a)') 'ExtData Run_: INTERP_LOOP: interpolating between bracket times' + Write(*,*) ' ==> variable: ', trim(item%var) + Write(*,*) ' ==> file: ', trim(item%file_template) + ENDIF + + ! finally interpolate between bracketing times + + call MAPL_ExtDataInterpField(item,self%ExtDataState,useTime(i),__RC__) + + endif + + nullify(item) + + end do INTERP_LOOP + + IF ( (Ext_Debug > 0) .AND. MAPL_Am_I_Root() ) THEN + Write(*,*) 'ExtData Run_: INTERP_LOOP: Done' + ENDIF + + call MAPL_TimerOff(MAPLSTATE,"-Interpolate") + + ! now take care of derived fields + do i=1,self%derived%nItems + + derivedItem => self%derived%item(i) + + call derivedItem%update_freq%check_update(doUpdate_,time,time0,.not.hasRun,__RC__) + !doUpdate_ = doUpdate_ .or. (.not.hasRun) + + if (doUpdate_) then + + call CalcDerivedField(self%ExtDataState,derivedItem%name,derivedItem%expression, & + derivedItem%masking,__RC__) + + end if + + end do + + IF ( (Ext_Debug > 0) .AND. MAPL_Am_I_Root() ) THEN + Write(*,*) 'ExtData Run_: End' + ENDIF + +! All done +! -------- + deallocate(doUpdate) + deallocate(useTime) + + if (hasRun .eqv. .false.) hasRun = .true. + call MAPL_TimerOff(MAPLSTATE,"Run") + call MAPL_TimerOff(MAPLSTATE,"TOTAL") + + _RETURN(ESMF_SUCCESS) + + END SUBROUTINE Run_ + +!------------------------------------------------------------------------- +! NASA/GSFC, Global Modeling and Assimilation Office, Code 610.1 ! +!------------------------------------------------------------------------- +!BOP +! +! !IROUTINE: Finalize_ --- Finalize MAPL_ExtData +! +! !INTERFACE: +! + + SUBROUTINE Finalize_ ( GC, IMPORT, EXPORT, CLOCK, rc ) + +! !USES: + + implicit NONE + +! !INPUT PARAMETERS: + + type(ESMF_Clock), intent(inout) :: CLOCK ! The clock + +! !OUTPUT PARAMETERS: + + type(ESMF_GridComp), intent(inout) :: GC ! Grid Component + type(ESMF_State), intent(inout) :: IMPORT ! Import State + type(ESMF_State), intent(inout) :: EXPORT ! Export State + integer, intent(out) :: rc ! Error return code: + ! 0 - all is well + ! 1 - + +! !DESCRIPTION: This is a simple ESMF wrapper. +! +! !REVISION HISTORY: +! +! 12Dec2009 da Silva Design and first implementation. +! +!EOP +!------------------------------------------------------------------------- + + type(MAPL_ExtData_state), pointer :: self ! Legacy state + type(ESMF_Config) :: CF ! Universal Config + + character(len=ESMF_MAXSTR) :: comp_name + character(len=ESMF_MAXSTR) :: Iam + integer :: status + + +! Get my name and set-up traceback handle +! --------------------------------------- + Iam = 'Finalize_' + call ESMF_GridCompGet( GC, name=comp_name, __RC__ ) + Iam = trim(comp_name) // trim(Iam) + +! Finalize MAPL Generic +! --------------------- + call MAPL_GenericFinalize ( GC, IMPORT, EXPORT, CLOCK, __RC__ ) + +! Extract relevant runtime information +! ------------------------------------ + call extract_ ( GC, self, CF, __RC__) + +! Free the memory used to hold the primary export items +! ----------------------------------------------------- + if (associated(self%primary%item)) then + deallocate(self%primary%item) + end if + + +! All done +! -------- + _RETURN(ESMF_SUCCESS) + + end SUBROUTINE Finalize_ + +!....................................................................... + + subroutine extract_ ( GC, self, CF, rc) + + type(ESMF_GridComp), intent(INout) :: GC ! Grid Comp object + + type(MAPL_ExtData_state), pointer :: self ! Legacy state + type(ESMF_Config), intent(out) :: CF ! Universal Config + + integer, intent(out), optional :: rc + +! --- + + character(len=ESMF_MAXSTR) :: comp_name + character(len=ESMF_MAXSTR) :: Iam + integer :: status + + type(MAPL_ExtData_Wrap) :: wrap + +! Get my name and set-up traceback handle +! --------------------------------------- + Iam = 'extract_' + call ESMF_GridCompGet( GC, NAME=comp_name, __RC__ ) + Iam = trim(COMP_NAME) // '::' // trim(Iam) + + If (present(rc)) rc=ESMF_SUCCESS + +! Get my internal state +! --------------------- + call ESMF_UserCompGetInternalState(gc, 'MAPL_ExtData_state', WRAP, STATUS) + _VERIFY(STATUS) + self => wrap%ptr + +! Get the configuration +! --------------------- + call ESMF_GridCompGet ( GC, config=CF, __RC__ ) + + + _RETURN(ESMF_SUCCESS) + + end subroutine extract_ + +! ............................................................................ + + logical function PrimaryExportIsConstant_(item) + + type(PrimaryExport), intent(in) :: item + + if ( item%update_freq%is_single_shot() .or. & + trim(item%file_template) == '/dev/null' ) then + PrimaryExportIsConstant_ = .true. + else + PrimaryExportIsConstant_ = .false. + end if + + end function PrimaryExportIsConstant_ + +! ............................................................................ + + logical function DerivedExportIsConstant_(item) + + type(DerivedExport), intent(in) :: item + + if ( item%update_freq%is_disabled() ) then + DerivedExportIsConstant_ = .true. + else + DerivedExportIsConstant_ = .false. + end if + + end function DerivedExportIsConstant_ + + ! ............................................................................ + + type (ESMF_Time) function timestamp_(time, template, rc) + type(ESMF_Time), intent(inout) :: time + character(len=ESMF_MAXSTR), intent(in) :: template + integer, optional, intent(inout) :: rc + + ! locals + integer, parameter :: DATETIME_MAXSTR_ = 32 + integer :: yy, mm, dd, hs, ms, ss + character(len=DATETIME_MAXSTR_) :: buff, buff_date, buff_time + character(len=DATETIME_MAXSTR_) :: str_yy, str_mm, str_dd + character(len=DATETIME_MAXSTR_) :: str_hs, str_ms, str_ss + + integer :: i, il, ir + integer :: status + + ! test the length of the timestamp template + _ASSERT(len_trim(template) < DATETIME_MAXSTR_,'Timestamp template is greater than Maximum allowed len') + + buff = trim(template) + buff = ESMF_UtilStringLowerCase(buff, __RC__) + + ! test if the template is empty and return the current time as result + if (buff == '-' .or. buff == '--' .or. buff == '---' .or. & + buff == 'na' .or. buff == 'none' .or. buff == 'n/a') then + + timestamp_ = time + else + ! split the time stamp template into a date and time strings + i = scan(buff, 't') + If (.not.(i > 3)) Then + _ASSERT(.False.,'ERROR: Time stamp ' // trim(template) // ' uses the fixed format, and must therefore contain a T') + End If + + buff_date = buff(1:i-1) + buff_time = buff(i+1:) + + ! parse the date string + il = scan(buff_date, '-', back=.false.) + ir = scan(buff_date, '-', back=.true. ) + str_yy = trim(buff_date(1:il-1)) + str_mm = trim(buff_date(il+1:ir-1)) + str_dd = trim(buff_date(ir+1:)) + + ! parse the time string + il = scan(buff_time, ':', back=.false.) + ir = scan(buff_time, ':', back=.true. ) + str_hs = trim(buff_time(1:il-1)) + str_ms = trim(buff_time(il+1:ir-1)) + str_ss = trim(buff_time(ir+1:)) + + ! remove the trailing 'Z' from the seconds string + i = scan(str_ss, 'z') + if (i > 0) then + str_ss = trim(str_ss(1:i-1)) + end if + + ! apply the timestamp template + call ESMF_TimeGet(time, yy=yy, mm=mm, dd=dd, h=hs, m=ms, s=ss, __RC__) + + i = scan(str_yy, '%'); if (i == 0) read (str_yy, '(I4)') yy + i = scan(str_mm, '%'); if (i == 0) read (str_mm, '(I2)') mm + i = scan(str_dd, '%'); if (i == 0) read (str_dd, '(I2)') dd + i = scan(str_hs, '%'); if (i == 0) read (str_hs, '(I2)') hs + i = scan(str_ms, '%'); if (i == 0) read (str_ms, '(I2)') ms + i = scan(str_ss, '%'); if (i == 0) read (str_ss, '(I2)') ss + + call ESMF_TimeSet(timestamp_, yy=yy, mm=mm, dd=dd, h=hs, m=ms, s=ss, __RC__) + end if + + _RETURN(ESMF_SUCCESS) + + end function timestamp_ + + subroutine GetLevs(item, rc) + + type(PrimaryExport) , intent(inout) :: item + integer, optional , intent(out ) :: rc + + integer :: status + + real, allocatable :: levFile(:) + character(len=ESMF_MAXSTR) :: levunits,tlevunits + character(len=:), allocatable :: levname + character(len=:), pointer :: positive + type(Variable), pointer :: var + integer :: i + + positive=>null() + var => null() + if (item%isVector) then + var=>item%file_metadata%get_variable(trim(item%fcomp1)) + _ASSERT(associated(var),"Variable "//TRIM(item%fcomp1)//" not found in file "//TRIM(item%file_template)) + var => null() + var=>item%file_metadata%get_variable(trim(item%fcomp2)) + _ASSERT(associated(var),"Variable "//TRIM(item%fcomp2)//" not found in file "//TRIM(item%file_template)) + else + var=>item%file_metadata%get_variable(trim(item%var)) + _ASSERT(associated(var),"Variable "//TRIM(item%var)//" not found in file "//TRIM(item%file_template)) + end if + + levName = item%file_metadata%get_level_name(rc=status) + _VERIFY(status) + if (trim(levName) /='') then + call item%file_metadata%get_coordinate_info(levName,coordSize=item%lm,coordUnits=tLevUnits,coords=levFile,__RC__) + levUnits=MAPL_TrimString(tlevUnits) + ! check if pressure + item%levUnit = ESMF_UtilStringLowerCase(levUnits) + if (trim(item%levUnit) == 'hpa' .or. trim(item%levUnit)=='pa') then + item%havePressure = .true. + end if + if (item%havePressure) then + if (levFile(1)>levFile(size(levFile))) item%fileVDir="up" + else + positive => item%file_metadata%get_variable_attribute(levName,'positive',__RC__) + if (associated(positive)) then + if (MAPL_TrimString(positive)=='up') item%fileVDir="up" + end if + end if + + allocate(item%levs(item%lm),__STAT__) + item%levs=levFile + if (trim(item%fileVDir)/=trim(item%importVDir)) then + do i=1,size(levFile) + item%levs(i)=levFile(size(levFile)-i+1) + enddo + end if + if (trim(item%levunit)=='hpa') item%levs=item%levs*100.0 + if (item%isVector) then + item%units = item%file_metadata%get_variable_attribute(trim(item%fcomp1),"units",rc=status) + _VERIFY(status) + else + item%units = item%file_metadata%get_variable_attribute(trim(item%var),"units",rc=status) + _VERIFY(status) + end if + + else + item%LM=0 + end if + + _RETURN(ESMF_SUCCESS) + + end subroutine GetLevs + + subroutine CalcDerivedField(state,exportName,exportExpr,masking,rc) + type(ESMF_State), intent(inout) :: state + character(len=*), intent(in ) :: exportName + character(len=*), intent(in ) :: exportExpr + logical, intent(in ) :: masking + integer, optional, intent(out ) :: rc + + integer :: status + + type(ESMF_Field) :: field + + if (masking) then + call MAPL_ExtDataEvaluateMask(state,exportName,exportExpr,__RC__) + else + call ESMF_StateGet(state,exportName,field,__RC__) + call MAPL_StateEval(state,exportExpr,field,__RC__) + end if + _RETURN(ESMF_SUCCESS) + end subroutine CalcDerivedField + + subroutine MAPL_ExtDataInterpField(item,state,time,rc) + type(PrimaryExport), intent(inout) :: item + type(ESMF_State), intent(in) :: state + type(ESMF_Time), intent(in ) :: time + integer, optional, intent(out ) :: rc + + integer :: status + type(ESMF_Field) :: field + + call ESMF_StateGet(state,item%vcomp1,field,__RC__) + call item%modelGridFields%comp1%interpolate_to_time(field,time,__RC__) + if (item%vartype == MAPL_VectorField) then + call ESMF_StateGet(state,item%vcomp1,field,__RC__) + call item%modelGridFields%comp2%interpolate_to_time(field,time,__RC__) + end if + _RETURN(ESMF_SUCCESS) + end subroutine MAPL_ExtDataInterpField + + subroutine MAPL_ExtDataVerticalInterpolate(ExtState,item,filec,rc) + type(MAPL_ExtData_State), intent(inout) :: ExtState + type(PrimaryExport), intent(inout) :: item + integer, intent(in ) :: filec + integer, optional, intent(out ) :: rc + + integer :: status + integer :: id_ps + type(ESMF_Field) :: field, newfield,psF + + if (item%do_VertInterp) then + if (trim(item%importVDir)/=trim(item%fileVDir)) then + call MAPL_ExtDataFlipVertical(item,filec,rc=status) + _VERIFY(status) + end if + if (item%vartype == MAPL_fieldItem) then + call MAPL_ExtDataGetBracket(item,filec,newField,getRL=.true.,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataGetBracket(item,filec,Field,rc=status) + _VERIFY(STATUS) + id_ps = ExtState%primaryOrder(1) + call MAPL_ExtDataGetBracket(ExtState%primary%item(id_ps),filec,field=psF,rc=status) + _VERIFY(STATUS) + call vertInterpolation_pressKappa(field,newfield,psF,item%levs,MAPL_UNDEF,rc=status) + _VERIFY(STATUS) + + else if (item%vartype == MAPL_VectorField) then + + id_ps = ExtState%primaryOrder(1) + call MAPL_ExtDataGetBracket(ExtState%primary%item(id_ps),filec,field=psF,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataGetBracket(item,filec,newField,getRL=.true.,vcomp=1,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataGetBracket(item,filec,Field,vcomp=1,rc=status) + _VERIFY(STATUS) + call vertInterpolation_pressKappa(field,newfield,psF,item%levs,MAPL_UNDEF,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataGetBracket(item,filec,newField,getRL=.true.,vcomp=2,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataGetBracket(item,filec,Field,vcomp=2,rc=status) + _VERIFY(STATUS) + call vertInterpolation_pressKappa(field,newfield,psF,item%levs,MAPL_UNDEF,rc=status) + _VERIFY(STATUS) + + end if + + else if (item%do_Fill) then + if (item%vartype == MAPL_fieldItem) then + call MAPL_ExtDataGetBracket(item,filec,newField,getRL=.true.,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataGetBracket(item,filec,Field,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataFillField(item,field,newfield,rc=status) + _VERIFY(STATUS) + else if (item%vartype == MAPL_VectorField) then + call MAPL_ExtDataGetBracket(item,filec,newField,getRL=.true.,vcomp=1,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataGetBracket(item,filec,Field,vcomp=1,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataFillField(item,field,newfield,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataGetBracket(item,filec,newField,getRL=.true.,vcomp=2,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataGetBracket(item,filec,Field,vcomp=2,rc=status) + _VERIFY(STATUS) + call MAPL_ExtDataFillField(item,field,newfield,rc=status) + _VERIFY(STATUS) + end if + else + if (trim(item%importVDir)/=trim(item%fileVDir)) then + call MAPL_ExtDataFlipVertical(item,filec,rc=status) + _VERIFY(status) + end if + end if + + _RETURN(ESMF_SUCCESS) + end subroutine MAPL_ExtDataVerticalInterpolate + + subroutine GetMaskName(FuncStr,Var,Needed,rc) + character(len=*), intent(in) :: FuncStr + character(len=*), intent(in) :: Var(:) + logical, intent(inout) :: needed(:) + integer, optional, intent(out) :: rc + + integer :: status + integer :: i1,i2,i,ivar + logical :: found,twovar + character(len=ESMF_MAXSTR) :: tmpstring,tmpstring1,tmpstring2,functionname + + i1 = index(Funcstr,"(") + _ASSERT(i1 > 0,'Incorrect format for function expression: missing "("') + functionname = adjustl(Funcstr(:i1-1)) + functionname = ESMF_UtilStringLowerCase(functionname, __RC__) + if (trim(functionname) == "regionmask") twovar = .true. + if (trim(functionname) == "zonemask") twovar = .false. + if (trim(functionname) == "boxmask") twovar = .false. + tmpstring = adjustl(Funcstr(i1+1:)) + i1 = index(tmpstring,",") + _ASSERT(i1 > 0,'Incorrect format for function expression: missing ","') + i2 = index(tmpstring,";") + if (twovar) then + tmpstring1 = adjustl(tmpstring(1:i1-1)) + tmpstring2 = adjustl(tmpstring(i1+1:i2-1)) + else + tmpstring1 = adjustl(tmpstring(1:i1-1)) + end if + + found = .false. + do i=1,size(var) + if ( trim(tmpstring1) == trim(var(i)) ) then + ivar = i + found = .true. + exit + end if + end do + _ASSERT(found,'Var ' // trim(tmpstring1) // ' not found') + needed(ivar) = .true. + + if (twovar) then + found = .false. + do i=1,size(var) + if ( trim(tmpstring2) == trim(var(i)) ) then + ivar = i + found = .true. + exit + end if + end do + _ASSERT(found,'Secound Var ' // trim(tmpstring2) // ' not found') + needed(ivar) = .true. + end if + _RETURN(ESMF_SUCCESS) + end subroutine GetMaskName + + subroutine MAPL_ExtDataEvaluateMask(state,exportName,exportExpr,rc) + + type(ESMF_STATE), intent(inout) :: state + character(len=*), intent(in) :: exportName + character(len=*), intent(in) :: exportExpr + integer, optional, intent(out) :: rc + + integer :: status + + integer :: k,i + character(len=ESMF_MAXSTR) :: maskString,maskname,vartomask,functionname,clatS,clatN + character(len=ESMF_MAXSTR) :: strtmp + integer, allocatable :: regionNumbers(:), flag(:) + integer, allocatable :: mask(:,:) + real, pointer :: rmask(:,:) => null() + real, pointer :: rvar2d(:,:) => null() + real, pointer :: rvar3d(:,:,:) => null() + real, pointer :: var2d(:,:) => null() + real, pointer :: var3d(:,:,:) => null() + real(REAL64), pointer :: lats(:,:) => null() + real(REAL64), pointer :: lons(:,:) => null() + real(REAL64) :: limitS, limitN, limitE, limitW + real(REAL64) :: limitE1, limitW1 + real(REAL64) :: limitE2, limitW2 + type(ESMF_Field) :: field + type(ESMF_Grid) :: grid + integer :: rank,ib,ie,is,i1,nargs + integer :: counts(3) + logical :: isCube, twoBox + real, allocatable :: temp2d(:,:) + character(len=ESMF_MAXSTR) :: args(5) + + call ESMF_StateGet(state,exportName,field,__RC__) + call ESMF_FieldGet(field,rank=rank,grid=grid,__RC__) + i1 = index(exportExpr,"(") + _ASSERT(i1 > 0,'Expected "(" in expression: ' // trim(exportExpr)) + functionname = adjustl(exportExpr(:i1-1)) + functionname = ESMF_UtilStringLowerCase(functionname, __RC__) + + if (trim(functionname) == "regionmask") then + ! get mask string + ib = index(exportExpr,";") + ie = index(exportExpr,")") + maskString = trim(exportExpr(ib+1:ie-1)) + ! get mask name + ie = index(exportExpr,";") + is = index(exportExpr,"(") + ib = index(exportExpr,",") + vartomask = trim(exportExpr(is+1:ib-1)) + maskname = trim(exportExpr(ib+1:ie-1)) + call MAPL_GetPointer(state,rmask,maskName,__RC__) + if (rank == 2) then + call MAPL_GetPointer(state,rvar2d,vartomask,__RC__) + call MAPL_GetPointer(state,var2d,exportName,__RC__) + else if (rank == 3) then + call MAPL_GetPointer(state,rvar3d,vartomask,__RC__) + call MAPL_GetPointer(state,var3d,exportName,__RC__) + else + _ASSERT(.false.,'Rank must be 2 or 3') + end if + + k=32 + allocate(regionNumbers(k), flag(k), stat=status) + _VERIFY(STATUS) + regionNumbers = 0 + call MAPL_ExtDataExtractIntegers(maskString,k,regionNumbers,rc=status) + _VERIFY(STATUS) + flag(:) = 1 + WHERE(regionNumbers(:) == 0) flag(:) = 0 + k = SUM(flag) + deallocate(flag,stat=status) + _VERIFY(STATUS) + + ! Set local mask to 1 where gridMask matches each integer (within precision!) + ! --------------------------------------------------------------------------- + allocate(mask(size(rmask,1),size(rmask,2)),stat=status) + _VERIFY(STATUS) + mask = 0 + DO i=1,k + WHERE(regionNumbers(i)-0.01 <= rmask .AND. & + rmask <= regionNumbers(i)+0.01) mask = 1 + END DO + + if (rank == 2) then + var2d = rvar2d + where(mask == 0) var2d = 0.0 + else if (rank == 3) then + var3d = rvar3d + do i=1,size(var3d,3) + where(mask == 0) var3d(:,:,i) = 0.0 + enddo + end if + deallocate( mask) + elseif(trim(functionname) == "zonemask") then + + ib = index(exportExpr,"(") + ie = index(exportExpr,",") + vartomask = trim(exportExpr(ib+1:ie-1)) + ib = index(exportExpr,",") + is = index(exportExpr,",",back=.true.) + ie = index(exportExpr,")") + clatS = trim(exportExpr(ib+1:is-1)) + clatN = trim(exportExpr(is+1:ie-1)) + READ(clatS,*,IOSTAT=status) limitS + _VERIFY(status) + READ(clatN,*,IOSTAT=status) limitN + _VERIFY(status) + + call ESMF_GridGetCoord(grid, coordDim=2, localDE=0, & + staggerloc=ESMF_STAGGERLOC_CENTER, farrayPtr=lats, rc=status) + _VERIFY(status) + limitN=limitN*MAPL_PI_R8/180.0d0 + limitS=limitS*MAPL_PI_R8/180.0d0 + + if (rank == 2) then + call MAPL_GetPointer(state,rvar2d,vartomask,__RC__) + call MAPL_GetPointer(state,var2d,exportName,__RC__) + else if (rank == 3) then + call MAPL_GetPointer(state,rvar3d,vartomask,__RC__) + call MAPL_GetPointer(state,var3d,exportName,__RC__) + else + _ASSERT(.false.,'Rank must be 2 or 3') + end if + + if (rank == 2) then + var2d = 0.0 + where(limitS <= lats .and. lats <=limitN) var2d = rvar2d + else if (rank == 3) then + var3d = 0.0 + do i=1,size(var3d,3) + where(limitS <= lats .and. lats <=limitN) var3d(:,:,i) = rvar3d(:,:,i) + enddo + end if + + elseif(trim(functionname) == "boxmask") then + is=index(exportExpr,'(') + ie=index(exportExpr,')') + strtmp = exportExpr(is+1:ie-1) + do nargs=1,5 + is = index(strtmp,',') + if (is >0) then + args(nargs) = strtmp(:is-1) + else + args(nargs) = strtmp + end if + strtmp = strtmp(is+1:) + end do + + varToMask=args(1) + + READ(args(2),*,IOSTAT=status) limitS + _VERIFY(status) + READ(args(3),*,IOSTAT=status) limitN + _VERIFY(status) + READ(args(4),*,IOSTAT=status) limitW + _VERIFY(status) + READ(args(5),*,IOSTAT=status) limitE + _VERIFY(status) + _ASSERT(limitE > limitW,'LimitE must be greater than limitW') + _ASSERT(limitE /= limitW,'LimitE cannot equal limitW') + _ASSERT(limitN /= limitS,'LimitN cannot equal LimitS') + _ASSERT((limitE-limitW)<=360.0d0,'(LimitE - LimitW) must be less than or equal to 360') + + call ESMF_GridGetCoord(grid, coordDim=1, localDE=0, & + staggerloc=ESMF_STAGGERLOC_CENTER, farrayPtr=lons, rc=status) + _VERIFY(status) + call ESMF_GridGetCoord(grid, coordDim=2, localDE=0, & + staggerloc=ESMF_STAGGERLOC_CENTER, farrayPtr=lats, rc=status) + _VERIFY(status) + + ! do some tests if cube goes from 0 to 360, lat-lon -180 to 180 + call MAPL_GridGet(grid, globalCellCountPerDim=COUNTS,rc=status) + _VERIFY(STATUS) + if (counts(2)==6*counts(1)) then + isCube=.true. + else + isCube=.false. + end if + + twoBox = .false. + if (isCube) then + if (limitW < 0.0d0 .and. limitE >=0.0d0) then + ! need two boxes + twoBox=.true. + limitW1=0.0d0 + limitE1=limitE + limitW2=limitW+360.0d0 + limitE2=360.0d0 + + else if (limitW <0.0d0 .and. limitE <0.0d0) then + ! just shift + limitW1=limitW+360.d0 + limitE1=limitE+360.d0 + + else + ! normal case + limitW1=limitW + limitE1=limitE + end if + + else + + if (limitW <= 180.0d0 .and. limitE > 180.0d0) then + ! need two boxes + twoBox=.true. + limitW1=limitW + limitE1=180.0d0 + limitW2=-180.d0 + limitE2=limitE-360.0d0 + else if (limitW > 180.0d0 .and. limitE > 180.0d0) then + ! just shift + limitW1=limitW-360.d0 + limitE1=limitE-360.d0 + else + ! normal case + limitW1=limitW + limitE1=limitE + end if + + end if + + limitE1=limitE1*MAPL_PI_R8/180.0d0 + limitW1=limitW1*MAPL_PI_R8/180.0d0 + limitE2=limitE2*MAPL_PI_R8/180.0d0 + limitW2=limitW2*MAPL_PI_R8/180.0d0 + + limitN=limitN*MAPL_PI_R8/180.0d0 + limitS=limitS*MAPL_PI_R8/180.0d0 + if (rank == 2) then + call MAPL_GetPointer(state,rvar2d,vartomask,__RC__) + call MAPL_GetPointer(state,var2d,exportName,__RC__) + else if (rank == 3) then + call MAPL_GetPointer(state,rvar3d,vartomask,__RC__) + call MAPL_GetPointer(state,var3d,exportName,__RC__) + else + _ASSERT(.false.,'Rank must be 2 or 3') + end if + + if (rank == 2) then + var2d = 0.0 + where(limitS <= lats .and. lats <=limitN .and. limitW1 <= lons .and. lons <= limitE1 ) var2d = rvar2d + else if (rank == 3) then + var3d = 0.0 + do i=1,size(var3d,3) + where(limitS <= lats .and. lats <=limitN .and. limitW1 <= lons .and. lons <= limitE1 ) var3d(:,:,i) = rvar3d(:,:,i) + enddo + end if + + if (twoBox) then + allocate(temp2d(size(var2d,1),size(var2d,2)),stat=status) + _VERIFY(STATUS) + if (rank == 2) then + temp2d = 0.0 + where(limitS <= lats .and. lats <=limitN .and. limitW2 <= lons .and. lons <= limitE2 ) temp2d = rvar2d + var2d=var2d+temp2d + else if (rank == 3) then + do i=1,size(var3d,3) + temp2d = 0.0 + where(limitS <= lats .and. lats <=limitN .and. limitW2 <= lons .and. lons <= limitE2 ) temp2d = rvar3d(:,:,i) + var3d(:,:,i)=var3d(:,:,i)+temp2d + enddo + end if + deallocate(temp2d) + end if + + end if + + _RETURN(ESMF_SUCCESS) + + end subroutine MAPL_ExtDataEvaluateMask + + SUBROUTINE MAPL_ExtDataExtractIntegers(string,iSize,iValues,delimiter,verbose,rc) + +! !USES: + + IMPLICIT NONE + +! !INPUT/OUTPUT PARAMETERS: + + CHARACTER(LEN=*), INTENT(IN) :: string ! Character-delimited string of integers + INTEGER, INTENT(IN) :: iSize + INTEGER, INTENT(INOUT) :: iValues(iSize)! Space allocated for extracted integers + CHARACTER(LEN=*), OPTIONAL :: delimiter ! 1-character delimiter + LOGICAL, OPTIONAL, INTENT(IN) :: verbose ! Let me know iValues as they are found. + ! DEBUG directive turns on the message even + ! if verbose is not present or if + ! verbose = .FALSE. + INTEGER, OPTIONAL, INTENT(OUT) :: rc ! Return code +! !DESCRIPTION: +! +! Extract integers from a character-delimited string, for example, "-1,45,256,7,10". In the context +! of Chem_Util, this is provided for determining the numerically indexed regions over which an +! emission might be applied. +! +! In multiple passes, the string is parsed for the delimiter, and the characters up to, but not +! including the delimiter are taken as consecutive digits of an integer. A negative sign ("-") is +! allowed. After the first pass, each integer and its trailing delimiter are lopped of the head of +! the (local copy of the) string, and the process is started over. +! +! The default delimiter is a comma (","). +! +! "Unfilled" iValues are zero. +! +! Return codes: +! 1 Zero-length string. +! 2 iSize needs to be increased. +! +! Assumptions/bugs: +! +! A non-zero return code does not stop execution. +! Allowed numerals are: 0,1,2,3,4,5,6,7,8,9. +! A delimiter must be separated from another delimiter by at least one numeral. +! The delimiter cannot be a numeral or a negative sign. +! The character following a negative sign must be an allowed numeral. +! The first character must be an allowed numeral or a negative sign. +! The last character must be an allowed numeral. +! The blank character (" ") cannot serve as a delimiter. +! +! Examples of strings that will work: +! "1" +! "-1" +! "-1,2004,-3" +! "1+-2+3" +! "-1A100A5" +! Examples of strings that will not work: +! "1,--2,3" +! "1,,2,3" +! "1,A,3" +! "1,-,2" +! "1,2,3,4," +! "+1" +! "1 3 6" +! +! !REVISION HISTORY: +! +! Taken from chem utilities. +! +!EOP + CHARACTER(LEN=*), PARAMETER :: Iam = 'Chem_UtilExtractIntegers' + + INTEGER :: base,count,i,iDash,last,lenStr + INTEGER :: multiplier,pos,posDelim,sign + CHARACTER(LEN=255) :: str + CHARACTER(LEN=1) :: char,delimChar + LOGICAL :: Done + LOGICAL :: tellMe + +! Initializations +! --------------- + If (present(rc)) rc=0 + count = 1 + Done = .FALSE. + iValues(:) = 0 + base = ICHAR("0") + iDash = ICHAR("-") + +! Determine verbosity, letting the DEBUG +! directive override local specification +! -------------------------------------- + tellMe = .FALSE. + IF(PRESENT(verbose)) THEN + IF(verbose) tellMe = .TRUE. + END IF +#ifdef DEBUG + tellMe = .TRUE. +#endif +! Check for zero-length string +! ---------------------------- + lenStr = LEN_TRIM(string) + IF(lenStr == 0) THEN + If (present(rc)) rc=1 + PRINT *,trim(IAm),": ERROR - Found zero-length string." + RETURN + END IF + +! Default delimiter is a comma +! ---------------------------- + delimChar = "," + IF(PRESENT(delimiter)) delimChar(1:1) = delimiter(1:1) + +! Work on a local copy +! -------------------- + str = TRIM(string) + +! One pass for each delimited integer +! ----------------------------------- + Parse: DO + + lenStr = LEN_TRIM(str) + +! Parse the string for the delimiter +! ---------------------------------- + posDelim = INDEX(TRIM(str),TRIM(delimChar)) + IF(tellMe) PRINT *,trim(Iam),": Input string is >",TRIM(string),"<" + +! If the delimiter does not exist, +! one integer remains to be extracted. +! ------------------------------------ + IF(posDelim == 0) THEN + Done = .TRUE. + last = lenStr + ELSE + last = posDelim-1 + END IF + multiplier = 10**last + +! Examine the characters of this integer +! -------------------------------------- + Extract: DO pos=1,last + + char = str(pos:pos) + i = ICHAR(char) + +! Account for a leading "-" +! ------------------------- + IF(pos == 1) THEN + IF(i == iDash) THEN + sign = -1 + ELSE + sign = 1 + END IF + END IF + +! "Power" of 10 for this character +! -------------------------------- + multiplier = multiplier/10 + + IF(pos == 1 .AND. sign == -1) CYCLE Extract + +! Integer comes from remaining characters +! --------------------------------------- + i = (i-base)*multiplier + iValues(count) = iValues(count)+i + IF(pos == last) THEN + iValues(count) = iValues(count)*sign + IF(tellMe) PRINT *,trim(Iam),":Integer number ",count," is ",iValues(count) + END IF + + END DO Extract + + IF(Done) EXIT + +! Lop off the leading integer and try again +! ----------------------------------------- + str(1:lenStr-posDelim) = str(posDelim+1:lenStr) + str(lenStr-posDelim+1:255) = " " + count = count+1 + +! Check size +! ---------- + IF(count > iSize) THEN + If (present(rc)) rc=2 + PRINT *,trim(Iam),": ERROR - iValues does not have enough elements." + END IF + + END DO Parse + + _RETURN(ESMF_SUCCESS) + + END SUBROUTINE MAPL_ExtDataExtractIntegers + + function MAPL_ExtDataGridChangeLev(Grid,CF,lm,rc) result(NewGrid) + + type(ESMF_Grid), intent(inout) :: Grid + type(ESMF_Config), intent(inout) :: CF + integer, intent(in) :: lm + integer, optional, intent(out) :: rc + + integer :: status + character(len=ESMF_MAXSTR) :: Iam + + character(len=ESMF_MAXSTR) :: gname, comp_name + integer :: counts(3) + integer :: NX,NY + type(ESMF_Grid) :: newGrid + type(ESMF_Config) :: cflocal + character(len=*), parameter :: CF_COMPONENT_SEPARATOR = '.' + real :: temp_real + + IAM = "MAPL_ExtDataGridChangeLev" + + call MAPL_GridGet(grid,globalCellCountPerDim=counts,__RC__) + call ESMF_GridGet(grid,name=gName,__RC__) + call ESMF_ConfigGetAttribute(CF, value = NX, Label="NX:", __RC__) + call ESMF_ConfigGetAttribute(CF, value = NY, Label="NY:", __RC__) + + comp_name = "ExtData" + cflocal = MAPL_ConfigCreate(rc=status) + _VERIFY(status) + call MAPL_ConfigSetAttribute(cflocal,value=NX, label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"NX:",rc=status) + _VERIFY(status) + call MAPL_ConfigSetAttribute(cflocal,value=lm, label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"LM:",rc=status) + _VERIFY(status) + + if (counts(2) == 6*counts(1)) then + call MAPL_ConfigSetAttribute(cflocal,value="Cubed-Sphere", label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"GRID_TYPE:",rc=status) + _VERIFY(status) + call MAPL_ConfigSetAttribute(cflocal,value=6, label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"NF:",rc=status) + _VERIFY(status) + call MAPL_ConfigSetAttribute(cflocal,value=counts(1), label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"IM_WORLD:",rc=status) + _VERIFY(status) + call MAPL_ConfigSetAttribute(cflocal,value=ny/6, label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"NY:",rc=status) + _VERIFY(status) + call MAPL_ConfigSetAttribute(cflocal,value=trim(gname), label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"GRIDNAME:",rc=status) + _VERIFY(status) + call ESMF_AttributeGet(grid, name='STRETCH_FACTOR', value=temp_real, rc=status) + if (status == ESMF_SUCCESS) then + call MAPL_ConfigSetAttribute(cflocal,value=temp_real, label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"STRETCH_FACTOR:",rc=status) + _VERIFY(status) + endif + call ESMF_AttributeGet(grid, name='TARGET_LON', value=temp_real, rc=status) + if (status == ESMF_SUCCESS) then + call MAPL_ConfigSetAttribute(cflocal,value=temp_real*MAPL_RADIANS_TO_DEGREES, label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"TARGET_LON:",rc=status) + _VERIFY(status) + endif + call ESMF_AttributeGet(grid, name='TARGET_LAT', value=temp_real, rc=status) + if (status == ESMF_SUCCESS) then + call MAPL_ConfigSetAttribute(cflocal,value=temp_real*MAPL_RADIANS_TO_DEGREES, label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"TARGET_LAT:",rc=status) + _VERIFY(status) + endif + else + call MAPL_ConfigSetAttribute(cflocal,value=counts(1), label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"IM_WORLD:",rc=status) + _VERIFY(status) + call MAPL_ConfigSetAttribute(cflocal,value=counts(2), label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"JM_WORLD:",rc=status) + _VERIFY(status) + call MAPL_ConfigSetAttribute(cflocal,value=ny, label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"NY:",rc=status) + _VERIFY(status) + call MAPL_ConfigSetAttribute(cflocal,value=trim(gname), label=trim(COMP_Name)//CF_COMPONENT_SEPARATOR//"GRIDNAME:",rc=status) + _VERIFY(status) + end if + newgrid = grid_manager%make_grid(cflocal, prefix=trim(COMP_Name)//".", rc=status) + _VERIFY(status) + + _RETURN(ESMF_SUCCESS) + + end function MAPL_ExtDataGridChangeLev + + subroutine MAPL_ExtDataGetBracket(item,Bside,field,bundle,getRL,vcomp,rc) + + type(PrimaryExport), intent(inout) :: item + integer, intent(in ) :: bside + type(ESMF_Field), optional, intent(inout) :: field + type(ESMF_FieldBundle), optional, intent(inout) :: bundle + logical, optional, intent(in ) :: getRL + integer, optional, intent(in ) :: vcomp + integer, optional, intent(out ) :: rc + + character(len=ESMF_MAXSTR) :: Iam + integer :: status + + logical :: getRL_ + + Iam = "MAPL_ExtDataGetBracket" + + if (present(getRL)) then + getRL_=getRL + else + getRL_=.false. + end if + + if (present(vcomp)) then + + if (present(field)) then + + if (Bside == MAPL_ExtDataLeft .and. vcomp == 1) then + if (getRL_) then + call item%modelGridFields%auxiliary1%get_parameters('L',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + else + call item%modelGridFields%comp1%get_parameters('L',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + end if + else if (Bside == MAPL_ExtDataLeft .and. vcomp == 2) then + if (getRL_) then + call item%modelGridFields%auxiliary2%get_parameters('L',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + else + call item%modelGridFields%comp2%get_parameters('L',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + end if + else if (Bside == MAPL_ExtDataRight .and. vcomp == 1) then + if (getRL_) then + call item%modelGridFields%auxiliary1%get_parameters('R',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + else + call item%modelGridFields%comp1%get_parameters('R',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + end if + else if (Bside == MAPL_ExtDataRight .and. vcomp == 2) then + if (getRL_) then + call item%modelGridFields%auxiliary2%get_parameters('R',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + else + call item%modelGridFields%comp2%get_parameters('R',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + end if + end if + + else if (present(bundle)) then + _RETURN(ESMF_FAILURE) + end if + + else + + if (present(field)) then + if (Bside == MAPL_ExtDataLeft) then + if (getRL_) then + call item%modelGridFields%auxiliary1%get_parameters('L',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + else + call item%modelGridFields%comp1%get_parameters('L',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + end if + else if (Bside == MAPL_ExtDataRight) then + if (getRL_) then + call item%modelGridFields%auxiliary1%get_parameters('R',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + else + call item%modelGridFields%comp1%get_parameters('R',field=field,__RC__) + _RETURN(ESMF_SUCCESS) + end if + end if + else if (present(bundle)) then + !if (Bside == MAPL_ExtDataLeft) then + !bundle = item%binterp1 + !_RETURN(ESMF_SUCCESS) + !else if (Bside == MAPL_ExtDataRight) then + !bundle = item%binterp2 + !_RETURN(ESMF_SUCCESS) + !end if + + end if + + end if + _RETURN(ESMF_FAILURE) + + end subroutine MAPL_ExtDataGetBracket + + subroutine MAPL_ExtDataFillField(item,FieldF,FieldR,rc) + + type(PrimaryExport), intent(inout) :: item + type(ESMF_Field), intent(inout) :: FieldF + type(ESMF_Field), intent(inout) :: FieldR + integer, optional, intent(out) :: rc + + character(len=ESMF_MAXSTR) :: Iam + integer :: status + + real, pointer :: ptrF(:,:,:),ptrR(:,:,:) + integer :: lm_in,lm_out,i + + Iam = "MAPL_ExtDataFillField" + + call ESMF_FieldGet(FieldF,0,farrayPtr=ptrF,rc=status) + _VERIFY(STATUS) + call ESMF_FieldGet(FieldR,0,farrayPtr=ptrR,rc=status) + _VERIFY(STATUS) + ptrF = 0.0 + lm_in= size(ptrR,3) + lm_out = size(ptrF,3) + if (trim(item%importVDir)=="down") then + + if (trim(item%fileVDir)=="down") then + do i=1,lm_in + ptrF(:,:,lm_out-lm_in+i)=ptrR(:,:,i) + enddo + else if (trim(item%fileVDir)=="up") then + do i=1,lm_in + ptrF(:,:,lm_out-i+1)=ptrR(:,:,i) + enddo + end if + else if (trim(item%importVDir)=="up") then + if (trim(item%fileVDir)=="down") then + do i=1,lm_in + ptrF(:,:,lm_in-i+1)=ptrR(:,:,i) + enddo + else if (trim(item%fileVDir)=="up") then + do i=1,lm_in + ptrF(:,:,i)=ptrR(:,:,i) + enddo + end if + end if + + _RETURN(ESMF_SUCCESS) + + end subroutine MAPL_ExtDataFillField + + subroutine MAPL_ExtDataFlipVertical(item,filec,rc) + type(PrimaryExport), intent(inout) :: item + integer, intent(in) :: filec + integer, optional, intent(out) :: rc + + integer :: status + + type(ESMF_Field) :: Field,field1,field2 + real, pointer :: ptr(:,:,:) + real, allocatable :: ptemp(:,:,:) + integer :: ls, le + + if (item%isVector) then + + if (item%do_Fill .or. item%do_VertInterp) then + call MAPL_ExtDataGetBracket(item,filec,field=Field1,vcomp=1,getRL=.true.,__RC__) + call MAPL_ExtDataGetBracket(item,filec,field=Field2,vcomp=2,getRL=.true.,__RC__) + else + call MAPL_ExtDataGetBracket(item,filec,field=Field1,vcomp=1,__RC__) + call MAPL_ExtDataGetBracket(item,filec,field=Field2,vcomp=2,__RC__) + end if + + call ESMF_FieldGet(Field1,0,farrayPtr=ptr,rc=status) + _VERIFY(STATUS) + allocate(ptemp,source=ptr,stat=status) + _VERIFY(status) + ls = lbound(ptr,3) + le = ubound(ptr,3) + ptr(:,:,le:ls:-1) = ptemp(:,:,ls:le:+1) + + call ESMF_FieldGet(Field2,0,farrayPtr=ptr,rc=status) + _VERIFY(STATUS) + ptemp=ptr + ptr(:,:,le:ls:-1) = ptemp(:,:,ls:le:+1) + + deallocate(ptemp) + + else + + if (item%do_Fill .or. item%do_VertInterp) then + call MAPL_ExtDataGetBracket(item,filec,field=Field,getRL=.true.,__RC__) + else + call MAPL_ExtDataGetBracket(item,filec,field=Field,__RC__) + end if + + call ESMF_FieldGet(Field,0,farrayPtr=ptr,rc=status) + _VERIFY(STATUS) + allocate(ptemp,source=ptr,stat=status) + _VERIFY(status) + ls = lbound(ptr,3) + le = ubound(ptr,3) + ptr(:,:,le:ls:-1) = ptemp(:,:,ls:le:+1) + deallocate(ptemp) + end if + + _RETURN(ESMF_SUCCESS) + + end subroutine MAPL_ExtDataFlipVertical + subroutine MAPL_ExtDataPopulateBundle(item,filec,pbundle,rc) + type(PrimaryExport), intent(inout) :: item + integer, intent(in) :: filec + type(ESMF_FieldBundle), intent(inout) :: pbundle + integer, optional, intent(out) :: rc + + integer :: status + + type(ESMF_Field) :: Field,field1,field2 + type(ESMF_Grid) :: grid + + if (item%isVector) then + + if (item%do_Fill .or. item%do_VertInterp) then + call MAPL_ExtDataGetBracket(item,filec,field=Field1,vcomp=1,getRL=.true.,__RC__) + call MAPL_ExtDataGetBracket(item,filec,field=Field2,vcomp=2,getRL=.true.,__RC__) + else + call MAPL_ExtDataGetBracket(item,filec,field=Field1,vcomp=1,__RC__) + call MAPL_ExtDataGetBracket(item,filec,field=Field2,vcomp=2,__RC__) + end if + + call ESMF_FieldGet(Field1,grid=grid,rc=status) + _VERIFY(STATUS) + call ESMF_FieldBundleSet(pbundle,grid=grid,rc=status) + _VERIFY(STATUS) + call MAPL_FieldBundleAdd(pbundle,Field1,rc=status) + _VERIFY(STATUS) + call MAPL_FieldBundleAdd(pbundle,Field2,rc=status) + _VERIFY(STATUS) + + !block + !character(len=ESMF_MAXSTR) :: vectorlist(2) + !vectorlist(1) = item%fcomp1 + !vectorlist(2) = item%fcomp2 + !call ESMF_AttributeSet(pbundle,name="VectorList:", itemCount=2, & + !valuelist = vectorlist, rc=status) + !_VERIFY(STATUS) + !end block + + else + + if (item%do_Fill .or. item%do_VertInterp) then + call MAPL_ExtDataGetBracket(item,filec,field=Field,getRL=.true.,__RC__) + else + call MAPL_ExtDataGetBracket(item,filec,field=Field,__RC__) + end if + + call ESMF_FieldGet(Field,grid=grid,rc=status) + _VERIFY(STATUS) + call ESMF_FieldBundleSet(pbundle,grid=grid,rc=status) + _VERIFY(STATUS) + call MAPL_FieldBundleAdd(pbundle,Field,rc=status) + _VERIFY(STATUS) + + end if + + _RETURN(ESMF_SUCCESS) + + end subroutine MAPL_ExtDataPopulateBundle + + subroutine MAPL_ExtDataCreateCFIO(IOBundles, rc) + type(IOBundleNGVector), target, intent(inout) :: IOBundles + integer, optional, intent(out ) :: rc + + type (IOBundleNGVectorIterator) :: bundle_iter + type (ExtDataNG_IOBundle), pointer :: io_bundle + integer :: status + + bundle_iter = IOBundles%begin() + do while (bundle_iter /= IOBundles%end()) + io_bundle => bundle_iter%get() + call io_bundle%make_cfio(__RC__) + call bundle_iter%next() + enddo + + _RETURN(ESMF_SUCCESS) + + end subroutine MAPL_ExtDataCreateCFIO + + subroutine MAPL_ExtDataDestroyCFIO(IOBundles,rc) + type(IOBundleNGVector), target, intent(inout) :: IOBundles + integer, optional, intent(out ) :: rc + + type(IOBundleNGVectorIterator) :: bundle_iter + type (ExtDataNG_IOBundle), pointer :: io_bundle + integer :: status + + bundle_iter = IOBundles%begin() + do while (bundle_iter /= IOBundles%end()) + io_bundle => bundle_iter%get() + call io_bundle%clean(__RC__) + call bundle_iter%next + enddo + call IOBundles%clear() + + _RETURN(ESMF_SUCCESS) + + end subroutine MAPL_ExtDataDestroyCFIO + + subroutine MAPL_ExtDataPrefetch(IOBundles,rc) + type(IOBundleNGVector), target, intent(inout) :: IOBundles + integer, optional, intent(out ) :: rc + + integer :: n,nfiles + type(ExtDataNG_IOBundle), pointer :: io_bundle => null() + integer :: status + + nfiles = IOBundles%size() + + do n = 1, nfiles + io_bundle => IOBundles%at(n) + call io_bundle%cfio%request_data_from_file(io_bundle%file_name,io_bundle%time_index,rc=status) + _VERIFY(status) + enddo + + _RETURN(ESMF_SUCCESS) + + end subroutine MAPL_ExtDataPrefetch + + subroutine MAPL_ExtDataReadPrefetch(IOBundles,rc) + type(IOBundleNGVector), target, intent(inout) :: IOBundles + integer, optional, intent(out ) :: rc + + integer :: nfiles, n + type (ExtDataNG_IOBundle), pointer :: io_bundle + integer :: status + + + nfiles = IOBundles%size() + do n=1, nfiles + io_bundle => IOBundles%at(n) + call io_bundle%cfio%process_data_from_file(rc=status) + _VERIFY(status) + enddo + + _RETURN(ESMF_SUCCESS) + + end subroutine MAPL_ExtDataReadPrefetch + + subroutine createFileLevBracket(item,cf,rc) + type(PrimaryExport), intent(inout) :: item + type(ESMF_Config), intent(inout) :: cf + integer, optional, intent(out) :: rc + + integer :: status + type (ESMF_Grid) :: grid, newgrid + type(ESMF_Field) :: field,new_field + + call item%modelGridFields%comp1%get_parameters('L',field=field,__RC__) + newGrid = MAPL_ExtDataGridChangeLev(grid,cf,item%lm,__RC__) + new_field = MAPL_FieldCreate(field,newGrid,lm=item%lm,newName=trim(item%fcomp1),__RC__) + call item%modelGridFields%auxiliary1%set_parameters(left_field=new_field,__RC__) + new_field = MAPL_FieldCreate(field,newGrid,lm=item%lm,newName=trim(item%fcomp1),__RC__) + call item%modelGridFields%auxiliary1%set_parameters(right_field=new_field,__RC__) + if (item%vartype==MAPL_VectorField) then + new_field = MAPL_FieldCreate(field,newGrid,lm=item%lm,newName=trim(item%fcomp2),__RC__) + call item%modelGridFields%auxiliary2%set_parameters(left_field=new_field,__RC__) + new_field = MAPL_FieldCreate(field,newGrid,lm=item%lm,newName=trim(item%fcomp2),__RC__) + call item%modelGridFields%auxiliary2%set_parameters(right_field=new_field,__RC__) + end if + _RETURN(_SUCCESS) + + end subroutine createFileLevBracket + + + subroutine IOBundle_Add_Entry(IOBundles,item,entry_num,rc) + type(IOBundleNGVector), intent(inout) :: IOBundles + type(primaryExport), intent(inout) :: item + integer, intent(in) :: entry_num + integer, intent(out), optional :: rc + + integer :: status + + type (ExtDataNG_IOBundle) :: io_bundle + type (GriddedIOItemVector) :: items + logical :: update + character(len=ESMF_MAXPATHLEN) :: current_file + integer :: time_index + + call item%modelGridFields%comp1%get_parameters('L',update=update,file=current_file,time_index=time_index) + if (update) then + call items%push_back(item%fileVars) + io_bundle = ExtDataNG_IOBundle(MAPL_ExtDataLeft, entry_num, current_file, time_index, item%trans, item%fracval, item%file_template, & + item%pfioCollection_id,item%iclient_collection_id,items,rc=status) + _VERIFY(status) + call IOBundles%push_back(io_bundle) + call extdata_lgr%info('%a update L with with: %a %i2 ',item%name, current_file, time_index) + end if + call item%modelGridFields%comp1%get_parameters('R',update=update,file=current_file,time_index=time_index) + if (update) then + call items%push_back(item%fileVars) + io_bundle = ExtDataNG_IOBundle(MAPL_ExtDataRight, entry_num, current_file, time_index, item%trans, item%fracval, item%file_template, & + item%pfioCollection_id,item%iclient_collection_id,items,rc=status) + _VERIFY(status) + call IOBundles%push_back(io_bundle) + call extdata_lgr%info('%a update R with with: %a %i2 ',item%name,current_file, time_index) + end if + + _RETURN(ESMF_SUCCESS) + + end subroutine IOBundle_Add_Entry + + END MODULE MAPL_ExtDataGridComp2G diff --git a/gridcomps/ExtData2G/ExtDataLgr.F90 b/gridcomps/ExtData2G/ExtDataLgr.F90 new file mode 100644 index 000000000000..48654ffa982c --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataLgr.F90 @@ -0,0 +1,8 @@ +module MAPL_ExtDataLogger + use pFlogger + + public :: extdata_lgr + class(Logger), pointer :: extdata_lgr + +end module MAPL_ExtDataLogger + diff --git a/gridcomps/ExtData2G/ExtDataNode.F90 b/gridcomps/ExtData2G/ExtDataNode.F90 new file mode 100644 index 000000000000..3270f9868f9c --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataNode.F90 @@ -0,0 +1,73 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" +module MAPL_ExtDataNode + use ESMF + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_BaseMod, only: MAPL_UNDEF + implicit none + private + + type, public :: ExtDataNode + type(ESMF_Field) :: field + type(ESMF_Time) :: time + character(len=ESMF_MAXPATHLEN) :: file + integer :: time_index + logical :: was_set = .false. + contains + procedure :: set + procedure :: get + procedure :: equals + generic :: operator(==) => equals + end type + +contains + + subroutine set(this, unusable, field, time, file, time_index, was_set, rc) + class(ExtDataNode), intent(inout) :: this + class(KeywordEnforcer), optional, intent(in) :: unusable + type(ESMF_Time), optional, intent(in) :: time + type(ESMF_Field), optional, intent(in) :: field + character(len=*), optional, intent(in) :: file + integer, optional, intent(in) :: time_index + logical, optional, intent(in) :: was_set + integer, optional, intent(out) :: rc + + _UNUSED_DUMMY(unusable) + if (present(time)) this%time = time + if (present(field)) this%field = field + if (present(file)) this%file = trim(file) + if (present(time_index)) this%time_index = time_index + if (present(was_set)) this%was_set = was_set + _RETURN(_SUCCESS) + + end subroutine set + + subroutine get(this, unusable, field, time, file, time_index, was_set, rc) + class(ExtDataNode), intent(inout) :: this + class(KeywordEnforcer), optional, intent(in) :: unusable + type(ESMF_Time), optional, intent(out) :: time + type(ESMF_Field), optional, intent(out) :: field + character(len=*), optional, intent(out) :: file + integer, optional, intent(out) :: time_index + logical, optional, intent(out) :: was_set + integer, optional, intent(out) :: rc + + _UNUSED_DUMMY(unusable) + if (present(time)) time = this%time + if (present(field)) field = this%field + if (present(file)) file = trim(this%file) + if (present(time_index)) time_index = this%time_index + if (present(was_set)) was_set = this%was_set + _RETURN(_SUCCESS) + + end subroutine get + + logical function equals(a,b) + class(ExtDataNode), intent(in) :: a + class(ExtDataNode), intent(in) :: b + + equals = (trim(a%file)==trim(b%file)) .and. (a%time==b%time) .and. (a%time_index==b%time_index) + end function equals + +end module MAPL_ExtDataNode diff --git a/gridcomps/ExtData2G/ExtDataOldTypesCreator.F90 b/gridcomps/ExtData2G/ExtDataOldTypesCreator.F90 new file mode 100644 index 000000000000..bdca0eea4066 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataOldTypesCreator.F90 @@ -0,0 +1,204 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" + +module MAPL_ExtDataOldTypesCreator + use ESMF + use MAPL_BaseMod + use yafYaml + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_ExtDataTypeDef + use MAPL_ExtDataConfig + use MAPL_ExtDataFileStream + use MAPL_ExtDataFileStreamMap + use MAPL_ExtDataRule + use MAPL_ExtDataRuleMap + use MAPL_ExtDataDerived + use MAPL_ExtDataDerivedMap + use MAPL_RegridMethods + use MAPL_ExtDataAbstractFileHandler + use MAPL_ExtDataSimpleFileHandler + use MAPL_ExtDataClimFileHandler + use MAPL_ExtDataTimeSample + use MAPL_ExtDataTimeSampleMap + implicit none + public :: ExtDataOldTypesCreator + + type, extends(ExtDataConfig) :: ExtDataOldTypesCreator + private + contains + procedure :: fillin_primary + procedure :: fillin_derived + end type ExtDataOldTypesCreator + + interface ExtDataOldTypesCreator + module procedure :: new_ExtDataOldTypesCreator + end interface + + contains + + function new_ExtDataOldTypesCreator(config_file,current_time,unusable,rc ) result(ExtDataObj) + character(len=*), intent(in) :: config_file + type(ESMF_Time), intent(in) :: current_time + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(ExtDataOldTypesCreator) :: ExtDataObj + + integer :: status + + _UNUSED_DUMMY(unusable) + call ExtDataObj%ExtDataConfig%new_ExtDataConfig_from_yaml(config_file,current_time,rc=status) + _VERIFY(status) + + _RETURN(_SUCCESS) + end function new_ExtDataOldTypesCreator + + + subroutine fillin_primary(this,item_name,primary_item,time,clock,unusable,rc) + class(ExtDataOldTypesCreator), intent(inout) :: this + character(len=*), intent(in) :: item_name + type(PrimaryExport), intent(inout) :: primary_item + type(ESMF_Time), intent(inout) :: time + type(ESMF_Clock), intent(inout) :: clock + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(ExtDataRule), pointer :: rule + type(ExtDataFileStream), pointer :: dataset + type(ExtDataTimeSample), pointer :: time_sample + type(ExtDataTimeSample), target :: default_time_sample + type(ExtDataSimpleFileHandler) :: simple_handler + type(ExtDataClimFileHandler) :: clim_handler + integer :: status, semi_pos + logical :: disable_interpolation + + _UNUSED_DUMMY(unusable) + rule => this%rule_map%at(trim(item_name)) + time_sample => this%sample_map%at(rule%sample_key) + + if(.not.associated(time_sample)) then + call default_time_sample%set_defaults() + time_sample=>default_time_sample + end if + primary_item%isVector = allocated(rule%vector_partner) + ! name and file var + primary_item%name = trim(item_name) + if (primary_item%isVector) then + primary_item%vartype = MAPL_VectorField + primary_item%vcomp1 = trim(item_name) + primary_item%vcomp2 = trim(rule%vector_partner) + primary_item%var = rule%file_var + primary_item%fcomp1 = rule%file_var + primary_item%fcomp2 = rule%vector_file_partner + primary_item%fileVars%itemType = ItemTypeVector + primary_item%fileVars%xname = trim(rule%file_var) + primary_item%fileVars%yname = trim(rule%vector_file_partner) + else + primary_item%vartype = MAPL_FieldItem + primary_item%vcomp1 = trim(item_name) + primary_item%var = rule%file_var + primary_item%fcomp1 = rule%file_var + primary_item%fileVars%itemType = ItemTypeScalar + primary_item%fileVars%xname = trim(rule%file_var) + end if + + ! regrid method + if (trim(rule%regrid_method) == "BILINEAR") then + primary_item%trans = REGRID_METHOD_BILINEAR + else if (trim(rule%regrid_method) == "CONSERVE") then + primary_item%trans = REGRID_METHOD_CONSERVE + else if (trim(rule%regrid_method) == "VOTE") then + primary_item%trans = REGRID_METHOD_VOTE + else if (index(rule%regrid_method,"FRACTION;")>0) then + semi_pos = index(rule%regrid_method,";") + read(rule%regrid_method(semi_pos+1:),*) primary_item%fracVal + primary_item%trans = REGRID_METHOD_FRACTION + else + _ASSERT(.false.,"Invalid regridding method") + end if + + if (trim(time_sample%extrap_outside) =="clim") then + primary_item%cycling=.true. + else if (trim(time_sample%extrap_outside) == "persist_closest") then + primary_item%persist_closest=.true. + else if (trim(time_sample%extrap_outside) == "none") then + primary_item%cycling=.false. + primary_item%persist_closest=.false. + end if + + allocate(primary_item%source_time,source=time_sample%source_time) + ! new refresh + call primary_item%update_freq%create_from_parameters(time_sample%refresh_time, & + time_sample%refresh_frequency, time_sample%refresh_offset, time, clock, __RC__) + + disable_interpolation = .not.time_sample%time_interpolation + + call primary_item%modelGridFields%comp1%set_parameters(linear_trans=rule%linear_trans,disable_interpolation=disable_interpolation) + call primary_item%modelGridFields%comp2%set_parameters(linear_trans=rule%linear_trans,disable_interpolation=disable_interpolation) + call primary_item%modelGridFields%auxiliary1%set_parameters(linear_trans=rule%linear_trans, disable_interpolation=disable_interpolation) + call primary_item%modelGridFields%auxiliary2%set_parameters(linear_trans=rule%linear_trans, disable_interpolation=disable_interpolation) + + ! file_template + primary_item%isConst = .false. + if (index(rule%collection,"/dev/null")==0) then + dataset => this%file_stream_map%at(trim(rule%collection)) + primary_item%file_template = dataset%file_template + call dataset%detect_metadata(primary_item%file_metadata,time,get_range=(trim(time_sample%extrap_outside) /= "none"),__RC__) + else + primary_item%file_template = rule%collection + end if + + if (index(rule%collection,'/dev/null') /= 0) then + primary_item%isConst = .true. + primary_item%const=rule%linear_trans(1) + else + if (primary_item%cycling) then + call clim_handler%initialize(dataset,__RC__) + allocate(primary_item%filestream,source=clim_handler) + else + call simple_handler%initialize(dataset,persist_closest=primary_item%persist_closest,__RC__) + allocate(primary_item%filestream,source=simple_handler) + end if + end if + + _RETURN(_SUCCESS) + + end subroutine fillin_primary + + subroutine fillin_derived(this,item_name,derived_item,time,clock,unusable,rc) + class(ExtDataOldTypesCreator), intent(inout) :: this + character(len=*), intent(in) :: item_name + type(DerivedExport), intent(inout) :: derived_item + type(ESMF_Time), intent(inout) :: time + type(ESMF_Clock), intent(inout) :: clock + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(ExtDataDerived), pointer :: rule + integer :: status + type(ExtDataTimeSample), pointer :: time_sample + type(ExtDataTimeSample), target :: default_time_sample + + _UNUSED_DUMMY(unusable) + rule => this%derived_map%at(trim(item_name)) + derived_item%name = trim(item_name) + derived_item%expression = rule%expression + time_sample => this%sample_map%at(rule%sample_key) + + if(.not.associated(time_sample)) then + call default_time_sample%set_defaults() + time_sample=>default_time_sample + end if + call derived_item%update_freq%create_from_parameters(time_sample%refresh_time, & + time_sample%refresh_frequency, time_sample%refresh_offset, time, clock, __RC__) + derived_item%masking=.false. + if (index(derived_item%expression,"mask") /= 0 ) then + derived_item%masking=.true. + end if + + _RETURN(_SUCCESS) + + end subroutine fillin_derived + +end module MAPL_ExtDataOldTypesCreator diff --git a/gridcomps/ExtData2G/ExtDataRule.F90 b/gridcomps/ExtData2G/ExtDataRule.F90 new file mode 100644 index 000000000000..fa9ee35db272 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataRule.F90 @@ -0,0 +1,160 @@ +#include "MAPL_ErrLog.h" +module MAPL_ExtDataRule + use yaFyaml + use ESMF + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_TimeStringConversion + use MAPL_ExtDataTimeSample + use MAPL_ExtDataTimeSampleMap + implicit none + private + + type, public :: ExtDataRule + character(:), allocatable :: collection + character(:), allocatable :: file_var + character(:), allocatable :: sample_key + real, allocatable :: linear_trans(:) + character(:), allocatable :: regrid_method + character(:), allocatable :: vector_partner + character(:), allocatable :: vector_component + character(:), allocatable :: vector_file_partner + contains + procedure :: set_defaults + procedure :: split_vector + end type + + interface ExtDataRule + module procedure new_ExtDataRule + end interface + +contains + + function new_ExtDataRule(config,sample_map,key,unusable,rc) result(rule) + type(Configuration), intent(in) :: config + character(len=*), intent(in) :: key + type(ExtDataTimeSampleMap) :: sample_map + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(ExtDataRule) :: rule + logical :: is_present + integer :: status + type(Configuration) ::config1 + character(len=:), allocatable :: tempc + type(ExtDataTimeSample) :: ts + _UNUSED_DUMMY(unusable) + + if (allocated(tempc)) deallocate(tempc) + is_present = config%has("collection") + _ASSERT(is_present,"no collection present in ExtData export") + rule%collection = config%of("collection") + + if (allocated(tempc)) deallocate(tempc) + is_present = config%has("vname") + if (index(rule%collection,"/dev/null")==0) then + _ASSERT(is_present,"no vname present in ExtData export") + end if + if (is_present) then + tempc = config%of("vname") + rule%file_var=tempc + else + _ASSERT(.false.,"no variable name in rule") + end if + + if (config%has("sample")) then + config1=config%at("sample") + if (config1%is_mapping()) then + ts = ExtDataTimeSample(config1,_RC) + call sample_map%insert(trim(key)//"_sample",ts) + rule%sample_key=trim(key)//"_sample" + else if (config1%is_string()) then + rule%sample_key=config1 + else + _ASSERT(.false.,"sample entry unsupported") + end if + else + rule%sample_key = "" + end if + + if (allocated(rule%linear_trans)) deallocate(rule%linear_trans) + if (config%has("linear_transformation")) then + call config%get(rule%linear_trans,"linear_transformation") + else + allocate(rule%linear_trans,source=[0.0,0.0]) + end if + + if (allocated(tempc)) deallocate(tempc) + if (config%has("regrid")) then + tempc = config%of("regrid") + rule%regrid_method=tempc + else + rule%regrid_method="BILINEAR" + end if + + _RETURN(_SUCCESS) + end function new_ExtDataRule + + subroutine set_defaults(this,unusable,rc) + class(ExtDataRule), intent(inout), target :: this + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + _UNUSED_DUMMY(unusable) + this%collection='' + this%file_var='missing_variable' + this%regrid_method='BILINEAR' + _RETURN(_SUCCESS) + end subroutine set_defaults + + subroutine split_vector(this,original_key,ucomp,vcomp,unusable,rc) + class(ExtDataRule), intent(in) :: this + character(len=*), intent(in) :: original_key + type(ExtDataRule), intent(inout) :: ucomp,vcomp + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + integer :: semi_pos + character(len=:),allocatable :: uname,vname + + _UNUSED_DUMMY(unusable) + + semi_pos = index(this%file_var,";") + _ASSERT(semi_pos > 0,"vector rule does not have 2 variables in the file_var") + uname = this%file_var(1:semi_pos-1) + vname = this%file_var(semi_pos+1:len_trim(this%file_var)) + ucomp = this + vcomp = this + semi_pos = index(original_key,";") + ucomp%vector_partner = original_key(semi_pos+1:len_trim(original_key)) + vcomp%vector_partner = original_key(1:semi_pos-1) + ucomp%file_var = uname + vcomp%file_var = vname + ucomp%vector_file_partner = vname + vcomp%vector_file_partner = uname + ucomp%vector_component = "EW" + vcomp%vector_component = "NS" + _RETURN(_SUCCESS) + + end subroutine split_vector + +end module MAPL_ExtDataRule + +module MAPL_ExtDataRuleMap + use MAPL_ExtDataRule + +#include "types/key_deferredLengthString.inc" +#define _value type(ExtDataRule) +#define _alt + +#define _map ExtDataRuleMap +#define _iterator ExtDataRuleMapIterator + +#include "templates/map.inc" + +#undef _iterator +#undef _map + +#undef _alt +#undef _value + +end module MAPL_ExtDataRuleMap diff --git a/gridcomps/ExtData2G/ExtDataSample.F90 b/gridcomps/ExtData2G/ExtDataSample.F90 new file mode 100644 index 000000000000..ccf3d62c84dc --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataSample.F90 @@ -0,0 +1,114 @@ +#include "MAPL_ErrLog.h" +module MAPL_ExtDataTimeSample + use yaFyaml + use ESMF + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_TimeStringConversion + implicit none + private + + type, public :: ExtDataTimeSample + logical :: time_interpolation + type(ESMF_Time), allocatable :: source_time(:) + character(:), allocatable :: extrap_outside + character(:), allocatable :: refresh_time + character(:), allocatable :: refresh_frequency + character(:), allocatable :: refresh_offset + contains + procedure :: set_defaults + end type + + interface ExtDataTimeSample + module procedure new_ExtDataTimeSample + end interface + +contains + + function new_ExtDataTimeSample(config,unusable,rc) result(TimeSample) + type(Configuration), intent(in) :: config + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(ExtDataTimeSample) :: TimeSample + integer :: status + character(len=:), allocatable :: source_str + integer :: idx + _UNUSED_DUMMY(unusable) + + call TimeSample%set_defaults() + + if (config%has("extrapolation")) TimeSample%extrap_outside=config%of("extrapolation") + + if (config%has("time_interpolation")) then + TimeSample%time_interpolation = config%of("time_interpolation") + else + TimeSample%time_interpolation = .true. + end if + + if (config%has("update_reference_time")) TimeSample%refresh_time=config%of("update_reference_time") + + if (config%has("update_reference_time")) TimeSample%refresh_frequency=config%of("update_frequency") + + if (config%has("update_offset")) TimeSample%refresh_offset=config%of("update_offset") + + if (config%has("source_time")) then + call config%get(source_str,"source_time",rc=status) + _VERIFY(status) + if (allocated(TimeSample%source_time)) deallocate(TimeSample%source_time) + idx = index(source_str,',') + _ASSERT(idx/=0,'invalid specification of source_time') + allocate(TimeSample%source_time(2)) + TimeSample%source_time(1)=string_to_esmf_time(source_str(:idx-1)) + TimeSample%source_time(2)=string_to_esmf_time(source_str(idx+1:)) + else + if (.not.allocated(TimeSample%source_time)) allocate(TimeSample%source_time(0)) + end if + + _RETURN(_SUCCESS) + + end function new_ExtDataTimeSample + + + subroutine set_defaults(this,unusable,rc) + class(ExtDataTimeSample), intent(inout), target :: this + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + integer :: status + _UNUSED_DUMMY(unusable) + this%time_interpolation=.true. + this%extrap_outside='none' + this%refresh_time="00" + this%refresh_frequency="PT0S" + this%refresh_offset="PT0S" + if (allocated(this%source_time)) then + deallocate(this%source_time,stat=status) + _VERIFY(status) + end if + allocate(this%source_time(0),stat=status) + _VERIFY(status) + _RETURN(_SUCCESS) + end subroutine set_defaults + +end module MAPL_ExtDataTimeSample + +module MAPL_ExtDataTimeSampleMap + use MAPL_ExtDataTimeSample + +#include "types/key_deferredLengthString.inc" +#define _value type(ExtDataTimeSample) +#define _alt + +#define _map ExtDataTimeSampleMap +#define _iterator ExtDataTimeSampleMapIterator + +#include "templates/map.inc" + +#undef _iterator +#undef _map + +#undef _alt +#undef _value + +end module MAPL_ExtDataTimeSampleMap diff --git a/gridcomps/ExtData2G/ExtDataSimpleFileHandler.F90 b/gridcomps/ExtData2G/ExtDataSimpleFileHandler.F90 new file mode 100644 index 000000000000..7395aec3fb49 --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataSimpleFileHandler.F90 @@ -0,0 +1,159 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" +module MAPL_ExtdataSimpleFileHandler + use ESMF + use MAPL_ExtDataAbstractFileHandler + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_DataCollectionMod + use MAPL_CollectionVectorMod + use MAPL_DataCollectionManagerMod + use MAPL_FileMetadataUtilsMod + use MAPL_TimeStringConversion + use MAPL_StringTemplate + use MAPL_ExtDataBracket + use MAPL_ExtDataConstants + + implicit none + private + public ExtDataSimpleFileHandler + + type, extends(ExtDataAbstractFileHandler) :: ExtDataSimpleFileHandler + contains + procedure :: get_file_bracket + procedure :: get_file + end type + +contains + + subroutine get_file_bracket(this, input_time, source_time, bracket, rc) + class(ExtdataSimpleFileHandler), intent(inout) :: this + type(ESMF_Time), intent(in) :: input_time + type(ESMF_Time), intent(in) :: source_time(:) + type(ExtDataBracket), intent(inout) :: bracket + integer, optional, intent(out) :: rc + integer :: status + type(ESMF_TimeInterval) :: zero + + type(ESMF_Time) :: time + integer :: time_index + character(len=ESMF_MAXPATHLEN) :: current_file + logical :: get_left, get_right,in_range,was_set + type(ESMF_Time) :: target_time + + get_left=.true. + get_right=.true. + in_range=.true. + target_time=input_time + call bracket%set_parameters(intermittent_disable=.false.) + if (this%persist_closest) then + if (input_time <= this%valid_range(1)) then + target_time = this%valid_range(1) + get_right = .false. + in_range = .false. + call bracket%get_node('L',was_set=was_set) + if (was_set) get_left=.false. + call bracket%set_parameters(intermittent_disable=.true.) + else if (input_time >= this%valid_range(2)) then + target_time = this%valid_range(2) + get_right = .false. + in_range = .false. + call bracket%get_node('L',was_set=was_set) + if (was_set) get_left=.false. + call bracket%set_parameters(intermittent_disable=.true.) + end if + end if + if (bracket%time_in_bracket(target_time) .and. in_range) then + _RETURN(_SUCCESS) + end if + + call ESMF_TimeIntervalSet(zero,__RC__) + if (this%frequency == zero) then + current_file = this%file_template + if (get_left) then + call this%get_time_on_file(current_file,target_time,'L',time_index,time,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found in file") + call bracket%set_node('L',file=current_file,time_index=time_index,time=time,__RC__) + if (in_range .and. (bracket%left_node == bracket%right_node)) then + call bracket%swap_node_fields(rc=status) + _VERIFY(status) + else + bracket%new_file_left=.true. + call bracket%set_node('L',was_set=.true.) + end if + end if + if (get_right) then + call this%get_time_on_file(current_file,target_time,'R',time_index,time,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found in file") + call bracket%set_node('R',file=current_file,time_index=time_index,time=time,__RC__) + bracket%new_file_right=.true. + end if + else + if (get_left) then + call this%get_file(current_file,target_time,0,__RC__) + call this%get_time_on_file(current_file,target_time,'L',time_index,time,__RC__) + if (time_index == time_not_found) then + call this%get_file(current_file,target_time,-1,__RC__) + call this%get_time_on_file(current_file,target_time,'L',time_index,time,__RC__) + _ASSERT(time_index/=time_not_found,"Time not found in file") + end if + call bracket%set_node('L',file=current_file,time_index=time_index,time=time,__RC__) + if (in_range .and. (bracket%left_node == bracket%right_node)) then + call bracket%swap_node_fields(rc=status) + _VERIFY(status) + else + bracket%new_file_left=.true. + call bracket%set_node('L',was_set=.true.) + end if + end if + + if (get_right) then + call this%get_file(current_file,target_time,0,__RC__) + call this%get_time_on_file(current_file,target_time,'R',time_index,time,__RC__) + if (time_index == time_not_found) then + call this%get_file(current_file,target_time,1,__RC__) + call this%get_time_on_file(current_file,target_time,'R',time_index,time,__RC__) + _ASSERT(time_index /= time_not_found,"Time not found in file") + end if + call bracket%set_node('R',file=current_file,time_index=time_index,time=time,__RC__) + bracket%new_file_right=.true. + end if + + end if + _RETURN(_SUCCESS) + + end subroutine get_file_bracket + + subroutine get_file(this,filename,input_time,shift,rc) + class(ExtdataSimpleFileHandler), intent(inout) :: this + character(len=*), intent(out) :: filename + type(ESMF_Time) :: input_time + integer, intent(in) :: shift + integer, intent(out), optional :: rc + + type(ESMF_Time) :: ftime + integer :: n,status + logical :: file_found + integer(ESMF_KIND_I8) :: interval_seconds + + call ESMF_TimeIntervalGet(this%frequency,s_i8=interval_seconds) + if (interval_seconds==0) then + ! time is not representable as absolute time interval (month, year etc...) do this + ! brute force way. Not good but ESMF leaves no choice + ftime=this%reff_time + do while (ftime < input_time) + ftime = ftime + this%frequency + enddo + ftime=ftime -this%frequency + shift*this%frequency + else + n = (input_time-this%reff_time)/this%frequency + ftime = this%reff_time+(n+shift)*this%frequency + end if + call fill_grads_template(filename,this%file_template,time=ftime,__RC__) + inquire(file=trim(filename),exist=file_found) + _ASSERT(file_found,"get_file did not file a file using: "//trim(this%file_template)) + _RETURN(_SUCCESS) + + end subroutine get_file + +end module MAPL_ExtdataSimpleFileHandler diff --git a/gridcomps/ExtData2G/ExtDataTypeDef.F90 b/gridcomps/ExtData2G/ExtDataTypeDef.F90 new file mode 100644 index 000000000000..e1d2f953b5dd --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataTypeDef.F90 @@ -0,0 +1,80 @@ +module MAPL_ExtDataTypeDef + use ESMF + use MAPL_GriddedIOItemMod + use MAPL_ExtDataBracket + use MAPL_ExtDataPointerUpdate + use MAPL_ExtDataAbstractFileHandler + use MAPL_FileMetadataUtilsMod + implicit none + + public PrimaryExport + public DerivedExport + public BracketingFields + + integer, parameter :: MAPL_ExtDataNullFrac = -9999 + + type BracketingFields + ! fields to store endpoints for interpolation of a vector pair + type(ExtDataBracket) :: comp1 + type(ExtDataBracket) :: comp2 + ! if vertically interpolating vector fields + type(ExtDataBracket) :: auxiliary1 + type(ExtDataBracket) :: auxiliary2 + end type BracketingFields + + type PrimaryExport + character(len=ESMF_MAXSTR) :: name + character(len=ESMF_MAXSTR) :: units='' + integer :: Trans + character(len=ESMF_MAXSTR) :: var + character(len=ESMF_MAXPATHLEN) :: file_template ! remove + + logical :: isConst + real :: Const !remove + integer :: vartype ! MAPL_FieldItem or MAPL_BundleItem + + class(ExtDataAbstractFileHandler), allocatable :: filestream + + ! if primary export represents a pair of vector fields + logical :: isVector + type(BracketingFields) :: modelGridFields + + ! names of the two vector components in the gridded component where import is declared + character(len=ESMF_MAXSTR) :: vcomp1, vcomp2 + ! the corresponding names of the two vector components on file + character(len=ESMF_MAXSTR) :: fcomp1, fcomp2 + type(GriddedIOitem) :: fileVars + + integer :: pfioCollection_id + integer :: iclient_collection_id + + logical :: ExtDataAlloc + integer :: FracVal = MAPL_ExtDataNullFrac + ! do we have to do vertical interpolation + logical :: do_VertInterp = .false. + logical :: do_Fill = .false. + type(FileMetadataUtils) :: file_metadata + integer :: LM + real, allocatable :: levs(:) + character(len=4) :: importVDir = "down" + character(len=4) :: fileVDir = "down" + character(len=ESMF_MAXSTR) :: levUnit + logical :: havePressure = .false. + type(ExtDataPointerUpdate) :: update_freq + + ! new stuff + logical :: cycling + logical :: persist_closest + type(ESMF_Time), allocatable :: source_time(:) + end type PrimaryExport + + type DerivedExport + character(len=ESMF_MAXSTR) :: name + character(len=ESMF_MAXPATHLEN) :: expression + logical :: ExtDataAlloc + logical :: masking + type(ExtDataPointerUpdate) :: update_freq + end type DerivedExport + + +end module MAPL_ExtDataTypeDef diff --git a/gridcomps/ExtData2G/ExtDataUpdatePointer.F90 b/gridcomps/ExtData2G/ExtDataUpdatePointer.F90 new file mode 100644 index 000000000000..7b71faf2074e --- /dev/null +++ b/gridcomps/ExtData2G/ExtDataUpdatePointer.F90 @@ -0,0 +1,109 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" + +module MAPL_ExtDataPointerUpdate + use ESMF + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + use MAPL_TimeStringConversion + implicit none + private + + public :: ExtDataPointerUpdate + + type :: ExtDataPointerUpdate + private + logical :: disabled = .false. + type(ESMF_Alarm) :: update_alarm + type(ESMF_TimeInterval) :: offset + logical :: single_shot = .false. + contains + procedure :: create_from_parameters + procedure :: check_update + procedure :: is_disabled + procedure :: is_single_shot + procedure :: disable + end type + + contains + + subroutine create_from_parameters(this,update_time,update_freq,update_offset,time,clock,rc) + class(ExtDataPointerUpdate), intent(inout) :: this + character(len=*), intent(in) :: update_time + character(len=*), intent(in) :: update_freq + character(len=*), intent(in) :: update_offset + type(ESMF_Time), intent(inout) :: time + type(ESMF_Clock), intent(inout) :: clock + integer, optional, intent(out) :: rc + + type(ESMF_Time) :: reference_time + type(ESMF_TimeInterval) :: reference_freq + integer :: status,int_time,year,month,day,hour,minute,second + + if (update_freq == "-") then + this%single_shot = .true. + else if (update_freq /= "PT0S") then + int_time = string_to_integer_time(update_time) + hour=int_time/10000 + minute=mod(int_time/100,100) + second=mod(int_time,100) + call ESMF_TimeGet(time,yy=year,mm=month,dd=day,__RC__) + call ESMF_TimeSet(reference_time,yy=year,mm=month,dd=day,h=hour,m=minute,s=second,__RC__) + reference_freq = string_to_esmf_timeinterval(update_freq,__RC__) + this%update_alarm = ESMF_AlarmCreate(clock,ringTime=reference_time,ringInterval=reference_freq,sticky=.false.,__RC__) + end if + this%offset=string_to_esmf_timeinterval(update_offset,__RC__) + _RETURN(_SUCCESS) + + end subroutine create_from_parameters + + subroutine check_update(this,do_update,working_time,current_time,first_time,rc) + class(ExtDataPointerUpdate), intent(inout) :: this + logical, intent(out) :: do_update + type(ESMF_Time), intent(inout) :: working_time + type(ESMF_Time), intent(inout) :: current_time + logical, intent(in) :: first_time + integer, optional, intent(out) :: rc + type(ESMF_Time) :: previous_ring + + integer :: status + + if (this%disabled) then + do_update = .false. + _RETURN(_SUCCESS) + end if + if (ESMF_AlarmIsCreated(this%update_alarm)) then + if (first_time) then + call ESMF_AlarmGet(this%update_alarm,prevRingTime=previous_ring,__RC__) + working_time =previous_ring+this%offset + do_update = .true. + else + do_update = ESMF_AlarmIsRinging(this%update_alarm,__RC__) + working_time = current_time+this%offset + end if + else + do_update = .true. + if (this%single_shot) this%disabled = .true. + working_time = current_time+this%offset + end if + + end subroutine check_update + + function is_disabled(this) result(disabled) + class(ExtDataPointerUpdate), intent(in) :: this + logical :: disabled + disabled = this%disabled + end function is_disabled + + function is_single_shot(this) result(single_shot) + class(ExtDataPointerUpdate), intent(in) :: this + logical :: single_shot + single_shot = this%single_shot + end function is_single_shot + + subroutine disable(this) + class(ExtDataPointerUpdate), intent(inout) :: this + this%disabled = .true. + end subroutine + +end module MAPL_ExtDataPointerUpdate diff --git a/gridcomps/ExtData2G/ExtData_IOBundleMod.F90 b/gridcomps/ExtData2G/ExtData_IOBundleMod.F90 new file mode 100644 index 000000000000..1e116ee47a8d --- /dev/null +++ b/gridcomps/ExtData2G/ExtData_IOBundleMod.F90 @@ -0,0 +1,127 @@ +!#include "MAPL_Exceptions.h" +#include "MAPL_Generic.h" +#include "unused_dummy.H" + +!------------------------------------------------------------------------- +! NASA/GSFC, Global Modeling and Assimilation Office, Code 610.1 ! +!------------------------------------------------------------------------- + +module MAPL_ExtDataNG_IOBundleMod + use ESMF + use MAPL_BaseMod + use MAPL_GriddedIOMod + use MAPL_ExceptionHandling + use MAPL_GriddedIOItemMod + use MAPL_GriddedIOItemVectorMod + + public :: ExtDataNG_IOBundle + + type ExtDataNG_IOBundle + type (MAPL_GriddedIO) :: cfio + type (ESMF_FieldBundle) :: pbundle + character(:), allocatable :: template + integer :: regrid_method + + integer :: bracket_side + integer :: entry_index + character(:), allocatable :: file_name + integer :: time_index + integer :: fraction + integer :: metadata_coll_id + integer :: server_coll_id + type(GriddedIOItemVector) :: items + + contains + + procedure :: clean + procedure :: make_cfio + procedure :: assign + generic :: assignment(=) => assign + end type ExtDataNG_IOBundle + + + interface ExtDataNG_IOBundle + module procedure new_ExtDataNG_IOBundle + end interface ExtDataNG_IOBundle + +contains + + function new_ExtDataNG_IOBundle(bracket_side, entry_index, file_name, time_index, regrid_method, fraction, template, metadata_coll_id,server_coll_id,items,rc) result(io_bundle) + type (ExtDataNG_IOBundle) :: io_bundle + + integer, intent(in) :: bracket_side + integer, intent(in) :: entry_index + character(len=*), intent(in) :: file_name + integer, intent(in) :: time_index + integer, intent(in) :: regrid_method + integer, intent(in) :: fraction + character(len=*), intent(in) :: template + integer, intent(in) :: metadata_coll_id + integer, intent(in) :: server_coll_id + type(GriddedIOItemVector) :: items + integer, optional, intent(out) :: rc + + io_bundle%bracket_side = bracket_side + io_bundle%entry_index = entry_index + io_bundle%file_name = file_name + io_bundle%time_index = time_index + io_bundle%regrid_method = regrid_method + io_bundle%fraction = fraction + io_bundle%template = trim(template) + + io_bundle%metadata_coll_id=metadata_coll_id + io_bundle%server_coll_id=server_coll_id + io_bundle%items=items + + _RETURN(ESMF_SUCCESS) + end function new_ExtDataNG_IOBundle + + + subroutine clean(this, rc) + class (ExtDataNG_IOBundle), intent(inout) :: this + integer, optional, intent(out) :: rc + + integer :: status + call ESMF_FieldBundleDestroy(this%pbundle, noGarbage=.true.,rc=status) + _VERIFY(status) + + _RETURN(ESMF_SUCCESS) + + end subroutine clean + + + subroutine make_cfio(this, rc) + class (ExtDataNG_IOBundle), intent(inout) :: this + integer, optional, intent(out) :: rc + + this%cfio = MAPL_GriddedIO(output_bundle=this%pbundle,regrid_method=this%regrid_method, & + read_collection_id=this%server_coll_id, & + metadata_collection_id = this%metadata_coll_id, fraction = this%fraction, & + items=this%items) + + _RETURN(ESMF_SUCCESS) + + end subroutine make_cfio + + subroutine assign(to,from) + class(ExtDataNG_IOBundle), intent(out) :: to + type(ExtDataNG_IOBundle), intent(in) :: from + + to%bracket_side = from%bracket_side + to%entry_index = from%entry_index + to%file_name = from%file_name + to%time_index = from%time_index + to%regrid_method = from%regrid_method + to%fraction = from%fraction + to%template = from%template + + to%metadata_coll_id=from%metadata_coll_id + to%server_coll_id=from%server_coll_id + to%items=from%items + to%pbundle=from%pbundle + to%CFIO=from%CFIO + + end subroutine assign + +end module MAPL_ExtDataNG_IOBundleMod + diff --git a/gridcomps/ExtData2G/ExtData_IOBundleVectorMod.F90 b/gridcomps/ExtData2G/ExtData_IOBundleVectorMod.F90 new file mode 100644 index 000000000000..cdfc72c49b06 --- /dev/null +++ b/gridcomps/ExtData2G/ExtData_IOBundleVectorMod.F90 @@ -0,0 +1,10 @@ +module MAPL_ExtDataNG_IOBundleVectorMod + use MAPL_ExtDataNG_IOBundleMod + +#define _type type(ExtDataNG_IoBundle) +#define _vector IoBundleNGVector +#define _iterator IoBundleNGVectorIterator + +#include "templates/vector.inc" + +end module MAPL_ExtDataNG_IOBundleVectorMod diff --git a/gridcomps/ExtData2G/TimeStringConversion.F90 b/gridcomps/ExtData2G/TimeStringConversion.F90 new file mode 100644 index 000000000000..de5a527576de --- /dev/null +++ b/gridcomps/ExtData2G/TimeStringConversion.F90 @@ -0,0 +1,232 @@ +#include "MAPL_Exceptions.h" +#include "MAPL_ErrLog.h" +module MAPL_TimeStringConversion + use ESMF + use MAPL_KeywordEnforcerMod + use MAPL_ExceptionHandling + implicit none + private + + public :: string_to_integer_time + public :: string_to_integer_date + public :: string_to_esmf_time + public :: string_to_esmf_timeinterval + +contains + + function string_to_integer_time(time_string,unusable,rc) result(time) + character(len=*), intent(in) :: time_string + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + integer :: time + integer mpos(2), hpos(2), spos(2) + integer strlen + integer firstcolon, lastcolon + integer lastspace + integer hour,min,sec + + _UNUSED_DUMMY(unusable) + + strlen = LEN_TRIM (time_string) + + firstcolon = index(time_string, ':') + if (firstcolon .LE. 0) then + + ! If no colons, check for hour. + + ! Logic below assumes a null character or something else is after the hour + ! if we do not find a null character add one so that it correctly parses time + !if (time_string(strlen:strlen) /= char(0)) then + !time_string = trim(time_string)//char(0) + !strlen=len_trim(time_string) + !endif + lastspace = index(TRIM(time_string), ' ', BACK=.TRUE.) + if ((strlen-lastspace).eq.2 .or. (strlen-lastspace).eq.3) then + hpos(1) = lastspace+1 + hpos(2) = strlen-1 + read (time_string(hpos(1):hpos(2)), * ) hour + min = 0 + sec = 0 + else + hour = 0 + min = 0 + sec = 0 + endif + + else + hpos(1) = firstcolon - 2 + hpos(2) = firstcolon - 1 + lastcolon = index(time_string, ':', BACK=.TRUE.) + if ( lastcolon .EQ. firstcolon ) then + mpos(1) = firstcolon + 1 + mpos(2) = firstcolon + 2 + read (time_string(hpos(1):hpos(2)), * ) hour + read (time_string(mpos(1):mpos(2)), * ) min + sec = 0 + else + mpos(1) = firstcolon + 1 + mpos(2) = lastcolon - 1 + spos(1) = lastcolon + 1 + spos(2) = lastcolon + 2 + read (time_string(hpos(1):hpos(2)), * ) hour + read (time_string(mpos(1):mpos(2)), * ) min + read (time_string(spos(1):spos(2)), * ) sec + endif + endif + + time = hour*10000+min*100+sec + _RETURN(_SUCCESS) + + end function string_to_integer_time + + function string_to_integer_date(time_string,unusable,rc) result(date) + character(len=*), intent(in) :: time_string + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + integer :: date + integer ypos(2), mpos(2), dpos(2) + integer strlen + integer firstdash, lastdash + integer year,month,day + + _UNUSED_DUMMY(unusable) + + strlen = LEN_TRIM (time_string) + + firstdash = index(time_string, '-') + lastdash = index(time_string, '-', BACK=.TRUE.) + + if (firstdash .LE. 0 .OR. lastdash .LE. 0) then + _RETURN(_FAILURE) + endif + + ypos(2) = firstdash - 1 + mpos(1) = firstdash + 1 + ypos(1) = ypos(2) - 3 + + mpos(2) = lastdash - 1 + dpos(1) = lastdash + 1 + dpos(2) = dpos(1) + 1 + + read ( time_string(ypos(1):ypos(2)), * ) year + read ( time_string(mpos(1):mpos(2)), * ) month + read ( time_string(dpos(1):dpos(2)), * ) day + + date = year*10000+month*100+day + _RETURN(_SUCCESS) + + end function string_to_integer_date + + function string_to_esmf_time(input_string,unusable,rc) result(time) + character(len=*), intent(in) :: input_string + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(ESMF_Time) :: time + integer :: status + integer :: tpos + integer year,month,day,hour,min,sec + integer :: int_time, int_date + character(len=:), allocatable :: date_string,time_string + + _UNUSED_DUMMY(unusable) + + tpos = index(input_string,'T') + _ASSERT(tpos >0,"Invalid date/time format, missing date/time separator") + + date_string = input_string(:tpos-1) + time_string = input_string(tpos+1:) + int_time = string_to_integer_time(time_string,__RC__) + int_date = string_to_integer_date(date_string,__RC__) + + year=int_date/10000 + month=mod(int_date/100,100) + day=mod(int_date,100) + hour=int_time/10000 + min=mod(int_time/100,100) + sec=mod(int_time,100) + call ESMF_TimeSet(time,yy=year,mm=month,dd=day,h=hour,m=min,s=sec,__RC__) + _RETURN(_SUCCESS) + + end function string_to_esmf_time + + function string_to_esmf_timeinterval(time_interval_string,unusable,rc) result(time_interval) + character(len=*), intent(in) :: time_interval_string + class(KeywordEnforcer), optional, intent(in) :: unusable + integer, optional, intent(out) :: rc + + type(ESMF_TimeInterval) :: time_interval + integer :: status + + integer :: strlen,ppos,cpos,lpos,tpos + integer year,month,day,hour,min,sec + character(len=:), allocatable :: date_string,time_string + _UNUSED_DUMMY(unusable) + + year=0 + month=0 + day=0 + hour=0 + min=0 + sec=0 + strlen = len_trim(time_interval_string) + tpos = index(time_interval_string,'T') + ppos = index(time_interval_string,'P') + _ASSERT(time_interval_string(1:1) == 'P','Not valid time duration') + + if (tpos /= 0) then + if (tpos /= ppos+1) then + date_string = time_interval_string(ppos+1:tpos-1) + end if + time_string = time_interval_string(tpos+1:strlen) + else + date_string = time_interval_string(ppos+1:strlen) + end if + + if (allocated(date_string)) then + strlen = len_trim(date_string) + lpos = 0 + cpos = index(date_string,'Y') + if (cpos /= 0) then + read(date_string(lpos+1:cpos-1),*)year + lpos = cpos + end if + cpos = index(date_string,'M') + if (cpos /= 0) then + read(date_string(lpos+1:cpos-1),*)month + lpos = cpos + end if + cpos = index(date_string,'D') + if (cpos /= 0) then + read(date_string(lpos+1:cpos-1),*)day + lpos = cpos + end if + end if + if (allocated(time_string)) then + strlen = len_trim(time_string) + lpos = 0 + cpos = index(time_string,'H') + if (cpos /= 0) then + read(time_string(lpos+1:cpos-1),*)hour + lpos = cpos + end if + cpos = index(time_string,'M') + if (cpos /= 0) then + read(time_string(lpos+1:cpos-1),*)min + lpos = cpos + end if + cpos = index(time_string,'S') + if (cpos /= 0) then + read(time_string(lpos+1:cpos-1),*)sec + lpos = cpos + end if + end if + + call ESMF_TimeIntervalSet(time_interval,yy=year,mm=month,d=day,h=hour,m=min,s=sec,__RC__) + _RETURN(_SUCCESS) + + end function string_to_esmf_timeinterval + +end module MAPL_TimeStringConversion diff --git a/gridcomps/MAPL_GridComps.F90 b/gridcomps/MAPL_GridComps.F90 index fea600aa200a..daedebb7f624 100644 --- a/gridcomps/MAPL_GridComps.F90 +++ b/gridcomps/MAPL_GridComps.F90 @@ -1,6 +1,7 @@ module MAPL_GridCompsMod use mapl_CapMod use mapl_CapOptionsMod + use mapl_externalGCStorage #ifdef USE_FLAP use mapl_FlapCLIMod #endif diff --git a/include/MAPL_ErrLog.h b/include/MAPL_ErrLog.h index 36e8bb9a69fe..6c5dacb8a597 100644 --- a/include/MAPL_ErrLog.h +++ b/include/MAPL_ErrLog.h @@ -2,7 +2,7 @@ ! The error logging may eventually evolve into a module based ! on the ESMF logger. For now these macros provide simple -! traceback capability. +! traceback capability. #ifndef MAPL_ErrLog_DONE @@ -44,6 +44,12 @@ # ifdef _RC # undef _RC # endif +# ifdef _STAT +# undef _STAT +# endif +# ifdef _IOSTAT +# undef _IOSTAT +# endif # ifdef __return # undef __return # endif @@ -55,7 +61,7 @@ # ifdef I_AM_MAIN # define __return call MAPL_abort() -# define __rc(rc) +# define __rc(rc) # else # define __return return # define __rc(rc) ,rc @@ -92,6 +98,8 @@ # define _RC_(rc,status) rc=status);_VERIFY(status # define _RC _RC_(rc,status) +# define _STAT _RC_(stat,status) +# define _IOSTAT _RC_(iostat,status) # define _ASSERT_MSG_AND_LOC_AND_RC(A,msg,stat,file,line,rc) if(MAPL_Assert(A,msg,stat,file,line __rc(rc))) __return diff --git a/profiler/MeterNode.F90 b/profiler/MeterNode.F90 index b5fd5d0a360d..043b14340671 100644 --- a/profiler/MeterNode.F90 +++ b/profiler/MeterNode.F90 @@ -1,5 +1,5 @@ module MAPL_MeterNode - use, intrinsic :: iso_fortran_env, only: REAL64, REAL128 + use, intrinsic :: iso_fortran_env, only: REAL64 use MAPL_AbstractMeter use MAPL_AbstractMeterNode use MAPL_MeterNodeVector @@ -67,10 +67,10 @@ module MAPL_MeterNode interface MeterNodeIterator module procedure new_MeterNodeIterator end interface MeterNodeIterator - + integer, parameter :: NOT_FOUND = -1 - + contains @@ -79,7 +79,7 @@ function new_MeterNode(name, meter, depth) result(tree) character(*), intent(in) :: name class(AbstractMeter), intent(in) :: meter integer, optional, intent(in) :: depth - + tree%name = name tree%meter = meter @@ -99,14 +99,14 @@ function get_meter(this) result(meter) class (MeterNode), target, intent(in) :: this meter => this%meter end function get_meter - + function get_name(this) result(name) character(:), pointer :: name class (MeterNode), target, intent(in) :: this name => this%name end function get_name - + function get_inclusive(this) result(inclusive) real(kind=REAL64) :: inclusive @@ -121,11 +121,14 @@ function get_exclusive(this) result(exclusive) type (MeterNodevectorIterator) :: iter class (AbstractMeterNode), pointer :: child - real(kind=REAL128) :: tmp + real(kind=REAL64) :: tmp + + ! Subtract time of submeters from time of node meter. + ! Previously, this used 128-bit precision to avoid negative + ! exclusive times due to roundoff. But the GNU on M1 and NVHPC do + ! not allow REAL128. So tmp is now 64-bit and we use a max(tmp,0) + ! below to try and cap negatives - ! Subtract time of submeters from time of node meter. Note the - ! use of 128-bit precision to avoid negative exclusive times due - ! to roundoff. tmp = this%get_inclusive() iter = this%children%begin() @@ -135,7 +138,7 @@ function get_exclusive(this) result(exclusive) call iter%next() end do - exclusive = tmp + exclusive = max(tmp, 0.0_REAL64) end function get_exclusive @@ -178,7 +181,7 @@ function get_child(this, name) result(child) character(*), intent(in) :: name integer :: idx - + idx = this%find_child(name) if (idx /= NOT_FOUND) then child => this%children%at(idx) @@ -232,7 +235,7 @@ recursive integer function get_num_nodes(this) result(num_nodes) type (MeterNodeVectorIterator) :: iter class (AbstractMeterNode), pointer :: child - + num_nodes = 1 iter = this%children%begin() do while (iter /= this%children%end()) @@ -271,7 +274,7 @@ function begin(this) result(iterator) allocate(iterator, source=MeterNodeIterator(this)) end function begin - + function end(this) result(iterator) @@ -293,7 +296,7 @@ function end(this) result(iterator) end function end - + recursive subroutine next(this) class (MeterNodeIterator), intent(inout) :: this class (AbstractMeterNode), pointer :: current_child @@ -318,7 +321,7 @@ recursive subroutine next(this) deallocate(this%iterator_of_current_child) call this%iterator_over_children%next() if (this%iterator_over_children == this%reference%children%end()) then ! done - deallocate(this%iterator_over_children) + deallocate(this%iterator_over_children) else current_child => this%iterator_over_children%get() this%iterator_of_current_child = current_child%begin() ! always at least one node @@ -326,7 +329,7 @@ recursive subroutine next(this) end if end if end if - + end subroutine next @@ -360,7 +363,7 @@ logical function equals(a, b) type is (MeterNodeIterator) equals = associated(a%reference, b%reference) if (.not. equals) return - + equals = associated(a%current) .eqv. associated(b%current) if (.not. equals) return @@ -423,7 +426,7 @@ recursive subroutine accumulate(this, other) call t%reset() end if call t%accumulate(other%get_meter()) - + ! recurse over children of other iter = other%begin() call iter%next() ! skip top node (already handled) @@ -434,5 +437,5 @@ recursive subroutine accumulate(this, other) end subroutine accumulate - + end module MAPL_MeterNode