diff --git a/src/offline/cable_driver_common.F90 b/src/offline/cable_driver_common.F90 index 2de00624f..272a5c8b2 100644 --- a/src/offline/cable_driver_common.F90 +++ b/src/offline/cable_driver_common.F90 @@ -104,6 +104,7 @@ MODULE cable_driver_common_mod PUBLIC :: cable_driver_init_default PUBLIC :: prepareFiles PUBLIC :: renameFiles + PUBLIC :: prepareFiles_princeton PUBLIC :: LUCdriver CONTAINS @@ -332,7 +333,8 @@ SUBROUTINE cable_driver_init_cru(dels, koffset, CRU) END SUBROUTINE cable_driver_init_cru SUBROUTINE prepareFiles(ncciy) - !* Select the correct files given the year for filenames following the gswp format + !* Select the correct files given the year for filenames following + ! the gswp format USE cable_IO_vars_module, ONLY: logn,gswpfile IMPLICIT NONE @@ -350,10 +352,11 @@ SUBROUTINE prepareFiles(ncciy) CALL renameFiles(logn,gswpfile%Tair,ncciy,'Tair') CALL renameFiles(logn,gswpfile%wind,ncciy,'wind') - END SUBROUTINE prepareFiles + END SUBROUTINE prepareFiles SUBROUTINE renameFiles(logn,inFile,ncciy,inName) - !! Replace the year in the filename with the value of ncciy. + !* Replace the year in the filename with the value of ncciy + ! for the gswp file format. IMPLICIT NONE INTEGER, INTENT(IN) :: logn !! Log file unit number @@ -370,6 +373,48 @@ SUBROUTINE renameFiles(logn,inFile,ncciy,inName) END SUBROUTINE renameFiles + SUBROUTINE prepareFiles_princeton(ncciy) + !* Select the correct files given the year for filenames following the + ! princeton format + USE cable_IO_vars_module, ONLY: logn,gswpfile + IMPLICIT NONE + INTEGER, INTENT(IN) :: ncciy + + WRITE(logn,*) 'CABLE offline global run using princeton forcing for ', ncciy + PRINT *, 'CABLE offline global run using princeton forcing for ', ncciy + + CALL renameFiles_princeton(logn,gswpfile%rainf,ncciy,'rainf') + CALL renameFiles_princeton(logn,gswpfile%LWdown,ncciy,'LWdown') + CALL renameFiles_princeton(logn,gswpfile%SWdown,ncciy,'SWdown') + CALL renameFiles_princeton(logn,gswpfile%PSurf,ncciy,'PSurf') + CALL renameFiles_princeton(logn,gswpfile%Qair,ncciy,'Qair') + CALL renameFiles_princeton(logn,gswpfile%Tair,ncciy,'Tair') + CALL renameFiles_princeton(logn,gswpfile%wind,ncciy,'wind') + + END SUBROUTINE prepareFiles_princeton + + SUBROUTINE renameFiles_princeton(logn,inFile,ncciy,inName) + !* Replace the year in the filename with the value of ncciy for + ! the princeton format + IMPLICIT NONE + INTEGER, INTENT(IN) :: logn,ncciy + INTEGER:: nn + CHARACTER(LEN=200), INTENT(INOUT) :: inFile + CHARACTER(LEN=*), INTENT(IN) :: inName + INTEGER :: idummy + + nn = INDEX(inFile,'19') + READ(inFile(nn:nn+3),'(i4)') idummy + WRITE(inFile(nn:nn+3),'(i4.4)') ncciy + nn = INDEX(inFile,'19', BACK=.TRUE.) + READ(inFile(nn:nn+3),'(i4)') idummy + WRITE(inFile(nn:nn+3),'(i4.4)') ncciy + READ(inFile(nn-5:nn-2),'(i4)') idummy + WRITE(inFile(nn-5:nn-2),'(i4.4)') ncciy + WRITE(logn,*) TRIM(inName), ' global data from ', TRIM(inFile) + + END SUBROUTINE renameFiles_princeton + !============================================================================== ! subroutine for SUBROUTINE LUCdriver( casabiome,casapool, & diff --git a/src/offline/cable_input.F90 b/src/offline/cable_input.F90 index 9c30208f7..debeb4ac7 100644 --- a/src/offline/cable_input.F90 +++ b/src/offline/cable_input.F90 @@ -400,7 +400,7 @@ SUBROUTINE open_met_file(dels,koffset,kend,spinup, TFRZ) PRINT*,'rainf' CALL handle_err( ok ) ENDIF - IF(.NOT. globalMetfile%l_gpcc) THEN + IF(.NOT. globalMetfile%l_gpcc .AND. cable_user%MetType .NE. "prin") THEN ok = NF90_OPEN(gswpfile%snowf,0,ncid_snow) IF (ok /= NF90_NOERR) THEN PRINT*,'snow' @@ -497,9 +497,12 @@ SUBROUTINE open_met_file(dels,koffset,kend,spinup, TFRZ) IF(ok/=NF90_NOERR) THEN ! if failed ! Try 'lon' instead of x ok = NF90_INQ_DIMID(ncid_met,'lon', xdimID) - IF(ok/=NF90_NOERR) CALL nc_abort & - (ok,'Error finding x dimension in '& - //TRIM(filename%met)//' (SUBROUTINE open_met_file)') + IF(ok/=NF90_NOERR) THEN ! MMY + ok = NF90_INQ_DIMID(ncid_met,'longitude', xdimID) ! MMY ! For princeton + IF(ok/=NF90_NOERR) CALL nc_abort & ! MMY + (ok,'Error finding x dimension in '& ! MMY + //TRIM(filename%met)//' (SUBROUTINE open_met_file)') ! MMY + END IF ! MMY END IF ok = NF90_INQUIRE_DIMENSION(ncid_met,xdimID,len=xdimsize) IF(ok/=NF90_NOERR) CALL nc_abort & @@ -510,9 +513,12 @@ SUBROUTINE open_met_file(dels,koffset,kend,spinup, TFRZ) IF(ok/=NF90_NOERR) THEN ! if failed ! Try 'lat' instead of y ok = NF90_INQ_DIMID(ncid_met,'lat', ydimID) - IF(ok/=NF90_NOERR) CALL nc_abort & - (ok,'Error finding y dimension in ' & - //TRIM(filename%met)//' (SUBROUTINE open_met_file)') + IF(ok/=NF90_NOERR) THEN ! MMY + ok = NF90_INQ_DIMID(ncid_met,'latitude', ydimID) ! MMY ! For princeton + IF(ok/=NF90_NOERR) CALL nc_abort & ! MMY + (ok,'Error finding y dimension in ' & ! MMY + //TRIM(filename%met)//' (SUBROUTINE open_met_file)') ! MMY + END IF ! MMY END IF ok = NF90_INQUIRE_DIMENSION(ncid_met,ydimID,len=ydimsize) IF(ok/=NF90_NOERR) CALL nc_abort & @@ -821,7 +827,10 @@ SUBROUTINE open_met_file(dels,koffset,kend,spinup, TFRZ) !===== done bug fixing for timevar in PALS met file =============== !===== gswp input file has bug in timeunits =========== - IF (ncciy > 0) WRITE(timeunits(26:27),'(i2.2)') 0 + IF (TRIM(cable_user%MetType) .NE. "prin") THEN ! MMY + IF (ncciy > 0) WRITE(timeunits(26:27),'(i2.2)') 0 + END IF ! MMY + !===== done bug fixing for timeunits in gwsp file ======== WRITE(logn,*) 'Time variable units: ', timeunits ! Get coordinate field: @@ -848,10 +857,17 @@ SUBROUTINE open_met_file(dels,koffset,kend,spinup, TFRZ) ! start time) from character to integer; calculate starting hour-of-day, ! day-of-year, year: IF (.NOT.cable_user%GSWP3) THEN - READ(timeunits(15:18),*) syear - READ(timeunits(20:21),*) smoy ! integer month - READ(timeunits(23:24),*) sdoytmp ! integer day of that month - READ(timeunits(26:27),*) shod ! starting hour of day + IF (cable_user%MetType .eq. "prin") THEN ! MMY + READ(timeunits(13:16),*) syear ! MMY + READ(timeunits(18:19),*) smoy ! integer month ! MMY + READ(timeunits(21:22),*) sdoytmp ! integer day of that month ! MMY + READ(timeunits(24:25),*) shod ! starting hour of day ! MMY + ELSE ! MMY + READ(timeunits(15:18),*) syear + READ(timeunits(20:21),*) smoy ! integer month + READ(timeunits(23:24),*) sdoytmp ! integer day of that month + READ(timeunits(26:27),*) shod ! starting hour of day + END IF ELSE syear=ncciy smoy=1 @@ -859,6 +875,7 @@ SUBROUTINE open_met_file(dels,koffset,kend,spinup, TFRZ) shod=0 END IF + ! if site data, shift start time to middle of timestep ! only do this if not already at middle of timestep ! vh_js ! @@ -1055,14 +1072,19 @@ SUBROUTINE open_met_file(dels,koffset,kend,spinup, TFRZ) IF(ok /= NF90_NOERR) CALL nc_abort & (ok,'Error finding Rainf units in met data file ' & //TRIM(filename%met)//' (SUBROUTINE open_met_file)') - IF(metunits%Rainf(1:8)=='kg/m^2/s'.OR.metunits%Rainf(1:6)=='kg/m2s'.OR.metunits%Rainf(1:10)== & - 'kgm^-2s^-1'.OR.metunits%Rainf(1:4)=='mm/s'.OR. & - metunits%Rainf(1:6)=='mms^-1'.OR. & - metunits%Rainf(1:7)=='kg/m^2s'.OR.metunits%Rainf(1:10)=='kg m-2 s-1'.OR.metunits%Wind(1:5)/='m s-1') THEN + IF(metunits%Rainf(1:8)=='kg/m^2/s' .OR. & + metunits%Rainf(1:7)=='kg/m2/s' .OR. & + metunits%Rainf(1:6)=='kg/m2s' .OR. & + metunits%Rainf(1:10)=='kgm^-2s^-1' .OR. & + metunits%Rainf(1:4)=='mm/s' .OR. & + metunits%Rainf(1:6)=='mms^-1' .OR. & + metunits%Rainf(1:7)=='kg/m^2s' .OR. & + metunits%Rainf(1:10)=='kg m-2 s-1' .OR. & + metunits%Wind(1:5)/='m s-1' ) THEN ! Change from mm/s to mm/time step: convert%Rainf = dels - ELSE IF(metunits%Rainf(1:4)=='mm/h'.OR.metunits%Rainf(1:6)== & - 'mmh^-1') THEN + ELSE IF(metunits%Rainf(1:4)=='mm/h' .OR. & + metunits%Rainf(1:6)== 'mmh^-1' ) THEN ! Change from mm/h to mm/time step: convert%Rainf = dels/3600.0 ELSE @@ -1817,19 +1839,39 @@ SUBROUTINE get_met_data(spinup,spinConv,met,soil,rad, & ! Get SWdown data for mask grid: IF (cable_user%GSWP3) ncid_met=ncid_sw ! since GSWP3 multiple met files - ok= NF90_GET_VAR(ncid_met,id%SWdown,tmpDat3, & - start=(/1,1,ktau/),count=(/xdimsize,ydimsize,1/)) - IF(ok /= NF90_NOERR) CALL nc_abort & - (ok,'Error reading SWdown in met data file ' & - //TRIM(filename%met)//' (SUBROUTINE get_met_data)') - ! Assign value to met data variable (no units change required): - !jhan:quick fix, use (1/2) dimension 1 here arbitrarily - DO i=1,mland ! over all land points/grid cells - met%fsd(landpt(i)%cstart:landpt(i)%cend,1) = & - 0.5 * REAL(tmpDat3(land_x(i),land_y(i),1)) - met%fsd(landpt(i)%cstart:landpt(i)%cend,2) = & - 0.5 * REAL(tmpDat3(land_x(i),land_y(i),1)) - ENDDO + + ! ______________________ MMY _________________________ + IF (TRIM(cable_user%MetType) .EQ. "prin") THEN + ok= NF90_GET_VAR(ncid_met,id%SWdown,tmpDat4, & + start=(/1,1,1,ktau/),count=(/xdimsize,ydimsize,1,1/)) + IF(ok /= NF90_NOERR) CALL nc_abort & + (ok,'Error reading SWdown in met data file ' & + //TRIM(filename%met)//' (SUBROUTINE get_met_data)') + DO i=1,mland ! over all land points/grid cells + met%fsd(landpt(i)%cstart:landpt(i)%cend,1) = & + 0.5 * REAL(tmpDat4(land_x(i),land_y(i),1,1)) + met%fsd(landpt(i)%cstart:landpt(i)%cend,2) = & + 0.5 * REAL(tmpDat4(land_x(i),land_y(i),1,1)) + ENDDO + ! PRINT *, "========== MMY ==========" + ! PRINT *, "met%fsd",met%fsd + ! ____________________________________________________ + ELSE ! MMY + + ok= NF90_GET_VAR(ncid_met,id%SWdown,tmpDat3, & + start=(/1,1,ktau/),count=(/xdimsize,ydimsize,1/)) + IF(ok /= NF90_NOERR) CALL nc_abort & + (ok,'Error reading SWdown in met data file ' & + //TRIM(filename%met)//' (SUBROUTINE get_met_data)') + ! Assign value to met data variable (no units change required): + !jhan:quick fix, use (1/2) dimension 1 here arbitrarily + DO i=1,mland ! over all land points/grid cells + met%fsd(landpt(i)%cstart:landpt(i)%cend,1) = & + 0.5 * REAL(tmpDat3(land_x(i),land_y(i),1)) + met%fsd(landpt(i)%cstart:landpt(i)%cend,2) = & + 0.5 * REAL(tmpDat3(land_x(i),land_y(i),1)) + ENDDO + END IF ! MMY ! inserted by rk4417 - phase2 ! Get Tair data for mask grid:- - - - - - - - - - - - - - - - - - IF(cable_user%GSWP3) ncid_met = ncid_ta ! since GSWP3 multiple met files @@ -2022,15 +2064,31 @@ SUBROUTINE get_met_data(spinup,spinConv,met,soil,rad, & ! Get Rainf and Snowf data for mask grid:- - - - - - - - - - - - - IF (cable_user%GSWP3) ncid_met = ncid_rain - ok= NF90_GET_VAR(ncid_met,id%Rainf,tmpDat3, & - start=(/1,1,ktau/),count=(/xdimsize,ydimsize,1/)) - IF(ok /= NF90_NOERR) CALL nc_abort & - (ok,'Error reading Rainf in met data file ' & - //TRIM(filename%met)//' (SUBROUTINE get_met_data)') - DO i=1,mland ! over all land points/grid cells - met%precip(landpt(i)%cstart:landpt(i)%cend) = & - REAL(tmpDat3(land_x(i),land_y(i),1)) ! store Rainf - ENDDO + +! ______________________ MMY _________________________ + IF (TRIM(cable_user%MetType) .eq. "prin") THEN ! MMY + ok= NF90_GET_VAR(ncid_met,id%Rainf,tmpDat4, & + start=(/1,1,1,ktau/),count=(/xdimsize,ydimsize,1,1/)) + IF(ok /= NF90_NOERR) CALL nc_abort & + (ok,'Error reading Rainf in met data file ' & + //TRIM(filename%met)//' (SUBROUTINE get_met_data)') + DO i=1,mland ! over all land points/grid cells + met%precip(landpt(i)%cstart:landpt(i)%cend) = & + REAL(tmpDat4(land_x(i),land_y(i),1,1)) ! store Rainf + ENDDO +! ____________________________________________________ + ELSE ! MMY + ok= NF90_GET_VAR(ncid_met,id%Rainf,tmpDat3, & + start=(/1,1,ktau/),count=(/xdimsize,ydimsize,1/)) + IF(ok /= NF90_NOERR) CALL nc_abort & + (ok,'Error reading Rainf in met data file ' & + //TRIM(filename%met)//' (SUBROUTINE get_met_data)') + DO i=1,mland ! over all land points/grid cells + met%precip(landpt(i)%cstart:landpt(i)%cend) = & + REAL(tmpDat3(land_x(i),land_y(i),1)) ! store Rainf + ENDDO + END IF + IF(exists%Snowf) THEN IF (cable_user%GSWP3) ncid_met = ncid_snow ok= NF90_GET_VAR(ncid_met,id%Snowf,tmpDat3, & @@ -2067,32 +2125,65 @@ SUBROUTINE get_met_data(spinup,spinConv,met,soil,rad, & ! Get LWdown data for mask grid: - - - - - - - - - - - - - - - - - IF (cable_user%GSWP3) ncid_met = ncid_lw - IF(exists%LWdown) THEN ! If LWdown exists in met file - ok= NF90_GET_VAR(ncid_met,id%LWdown,tmpDat3, & - start=(/1,1,ktau/),count=(/xdimsize,ydimsize,1/)) + + ! ______________________ MMY _________________________ + IF (TRIM(cable_user%MetType) .eq. "prin") THEN + ok= NF90_GET_VAR(ncid_met,id%LWdown,tmpDat4, & + start=(/1,1,1,ktau/),count=(/xdimsize,ydimsize,1,1/)) IF(ok /= NF90_NOERR) CALL nc_abort & - (ok,'Error reading LWdown in met data file ' & - //TRIM(filename%met)//' (SUBROUTINE get_met_data)') + (ok,'Error reading LWdown in met data file ' & + //TRIM(filename%met)//' (SUBROUTINE get_met_data)') DO i=1,mland ! over all land points/grid cells met%fld(landpt(i)%cstart:landpt(i)%cend) = & - REAL(tmpDat3(land_x(i),land_y(i),1)) + REAL(tmpDat4(land_x(i),land_y(i),1,1)) ENDDO - ELSE ! Synthesise LWdown based on temperature - ! Use Swinbank formula: - met%fld(:) = 0.0000094*0.0000000567*(met%tk(:)**6.0) + ! PRINT *, "========== MMY ==========" + ! PRINT *, "met%fld",met%fld + ! ____________________________________________________ + ELSE ! MMY + + IF(exists%LWdown) THEN ! If LWdown exists in met file + ok= NF90_GET_VAR(ncid_met,id%LWdown,tmpDat3, & + start=(/1,1,ktau/),count=(/xdimsize,ydimsize,1/)) + IF(ok /= NF90_NOERR) CALL nc_abort & + (ok,'Error reading LWdown in met data file ' & + //TRIM(filename%met)//' (SUBROUTINE get_met_data)') + DO i=1,mland ! over all land points/grid cells + met%fld(landpt(i)%cstart:landpt(i)%cend) = & + REAL(tmpDat3(land_x(i),land_y(i),1)) + ENDDO + ELSE ! Synthesise LWdown based on temperature + ! Use Swinbank formula: + met%fld(:) = 0.0000094*0.0000000567*(met%tk(:)**6.0) + END IF END IF ! Get CO2air data for mask grid:- - - - - - - - - - - - - - - - - - IF(exists%CO2air) THEN ! If CO2air exists in met file - ok= NF90_GET_VAR(ncid_met,id%CO2air,tmpDat4, & - start=(/1,1,1,ktau/),count=(/xdimsize,ydimsize,1,1/)) - IF(ok /= NF90_NOERR) CALL nc_abort & - (ok,'Error reading CO2air in met data file ' & - //TRIM(filename%met)//' (SUBROUTINE get_met_data)') - DO i=1,mland ! over all land points/grid cells - met%ca(landpt(i)%cstart:landpt(i)%cend) = & - REAL(tmpDat4(land_x(i),land_y(i),1,1))/1000000.0 - ENDDO +! ____________ MMY@23Apr2023 to read CO2air from PLUMBER2 __________ + ok = NF90_INQUIRE_VARIABLE(ncid_met,id%CO2air,ndims=ndims) + IF(ndims==3) THEN + ok= NF90_GET_VAR(ncid_met,id%CO2air,tmpDat3, & + start=(/1,1,ktau/),count=(/xdimsize,ydimsize,1/)) + IF(ok /= NF90_NOERR) CALL nc_abort & + (ok,'Error reading CO2air in met data file ' & + //TRIM(filename%met)//' (SUBROUTINE get_met_data)') + DO i=1,mland ! over all land points/grid cells + met%ca(landpt(i)%cstart:landpt(i)%cend) = & + REAL(tmpDat3(land_x(i),land_y(i),1))/1000000.0 + ENDDO + ELSE +! __________________________________________________________________ + ok= NF90_GET_VAR(ncid_met,id%CO2air,tmpDat4, & + start=(/1,1,1,ktau/),count=(/xdimsize,ydimsize,1,1/)) + IF(ok /= NF90_NOERR) CALL nc_abort & + (ok,'Error reading CO2air in met data file ' & + //TRIM(filename%met)//' (SUBROUTINE get_met_data)') + DO i=1,mland ! over all land points/grid cells + met%ca(landpt(i)%cstart:landpt(i)%cend) = & + REAL(tmpDat4(land_x(i),land_y(i),1,1))/1000000.0 + ENDDO + END IF ELSE ! Fix CO2 air concentration: met%ca(:) = fixedCO2 /1000000.0 diff --git a/src/offline/cable_mpimaster.F90 b/src/offline/cable_mpimaster.F90 index 973aa72d4..20b699edb 100644 --- a/src/offline/cable_mpimaster.F90 +++ b/src/offline/cable_mpimaster.F90 @@ -87,6 +87,7 @@ MODULE cable_mpimaster LALLOC, & prepareFiles, & renameFiles, & + prepareFiles_princeton, & LUCdriver USE cable_mpicommon USE cable_IO_vars_module, ONLY : NO_CHECK @@ -353,11 +354,17 @@ SUBROUTINE mpidrv_master (comm, trunk_sumbal, dels, koffset, kend, PLUME, CRU) ENDIF SELECT CASE (TRIM(cable_user%MetType)) - CASE ('gswp') + CASE ('gswp', 'prin') ncciy = CurYear WRITE(*,*) 'Looking for global offline run info.' - CALL prepareFiles(ncciy) + + IF ( TRIM(cable_user%MetType) == 'gswp' ) THEN + CALL prepareFiles(ncciy) + ELSE ! cable_user%MetType == 'prin' + CALL prepareFiles_princeton(ncciy) + END IF + CALL open_met_file( dels, koffset, kend, spinup, CTFRZ ) CASE ('plum') diff --git a/src/offline/cable_offline_driver.F90 b/src/offline/cable_offline_driver.F90 index d5e5d0f5f..c05b0bda4 100644 --- a/src/offline/cable_offline_driver.F90 +++ b/src/offline/cable_offline_driver.F90 @@ -44,6 +44,8 @@ PROGRAM cable_offline_driver CALL cable_driver_init_gswp(mpi_grp, GSWP_MID, NRRRR) CASE('gswp3') CALL cable_driver_init_gswp(mpi_grp) + CASE('prin') + CALL cable_driver_init_gswp(mpi_grp, GSWP_MID, NRRRR) CASE('plum') CALL cable_driver_init_plume(dels, koffset, PLUME) CASE('cru') diff --git a/src/offline/cable_serial.F90 b/src/offline/cable_serial.F90 index a919de3e8..88a920a56 100644 --- a/src/offline/cable_serial.F90 +++ b/src/offline/cable_serial.F90 @@ -76,6 +76,7 @@ MODULE cable_serial LALLOC, & prepareFiles, & renameFiles, & + prepareFiles_princeton, & LUCdriver USE cable_def_types_mod USE cable_IO_vars_module, ONLY: logn,gswpfile,ncciy,leaps, & @@ -297,10 +298,15 @@ SUBROUTINE serialdrv(trunk_sumbal, NRRRR, dels, koffset, kend, GSWP_MID, PLUME, ! Check for gswp run SELECT CASE (TRIM(cable_user%MetType)) - CASE ('gswp') + CASE ('gswp', 'prin') ncciy = CurYear - CALL prepareFiles(ncciy) + IF (cable_user%MetType == 'gswp') THEN + CALL prepareFiles(ncciy) + ELSE ! cable_user%MetType == 'prin' + CALL prepareFiles_princeton(ncciy) ! MMY + END IF + IF ( RRRR .EQ. 1 ) THEN CALL open_met_file( dels, koffset, kend, spinup, CTFRZ ) IF (leaps.AND.is_leapyear(YYYY).AND.kend.EQ.2920) THEN