C $Header: /u/gcmpack/MITgcm/pkg/dic/dic_fields_load.F,v 1.36 2016/01/16 21:57:12 jmc Exp $ C $Name: $ #include "DIC_OPTIONS.h" CBOP C !ROUTINE: DIC_FIELDS_LOAD C !INTERFACE: ========================================================== SUBROUTINE DIC_FIELDS_LOAD ( I myIter, myTime, myThid ) C !DESCRIPTION: C Read in fields needed for CO2,O2 fluxterms, silica for pH calculation C !USES: =============================================================== IMPLICIT NONE #include "SIZE.h" #include "EEPARAMS.h" #include "PARAMS.h" #include "GRID.h" #include "DIC_VARS.h" #include "DIC_LOAD.h" C !INPUT PARAMETERS: =================================================== C myIter :: current timestep C myTime :: current time C myThid :: thread number INTEGER myIter _RL myTime INTEGER myThid #ifdef ALLOW_DIC c !LOCAL VARIABLES: =================================================== INTEGER bi, bj, i, j INTEGER intimeP, intime0, intime1 _RL aWght,bWght #ifdef READ_PAR CHARACTER*(MAX_LEN_MBUF) msgBuf #endif CEOP IF ( DIC_forcingCycle.GT.0. _d 0 ) THEN C-- Now calculate whether it is time to update the forcing arrays CALL GET_PERIODIC_INTERVAL( O intimeP, intime0, intime1, bWght, aWght, I DIC_forcingCycle, DIC_forcingPeriod, I deltaTClock, myTime, myThid ) bi = myBxLo(myThid) bj = myByLo(myThid) #ifdef ALLOW_DEBUG IF ( debugLevel.GE.debLevB ) THEN _BEGIN_MASTER(myThid) WRITE(standardMessageUnit,'(A,I10,A,4I5,A,2F14.10)') & ' DIC_FIELDS_LOAD (Light),', myIter, & ' : iP,iLd,i0,i1=', intimeP,DIC_ldRec(bi,bj), intime0,intime1, & ' ; Wght=', bWght, aWght _END_MASTER(myThid) ENDIF #endif /* ALLOW_DEBUG */ #ifdef ALLOW_AUTODIFF C- assuming that we call S/R DIC_FIELDS_LOAD at each time-step and C with increasing time, this will catch when we need to load new records; C But with Adjoint run, this is not always the case => might end-up using C the wrong time-records IF ( intime0.NE.intimeP .OR. myIter.EQ.nIter0 ) THEN #else /* ALLOW_AUTODIFF */ C- Make no assumption on sequence of calls to DIC_FIELDS_LOAD ; C This is the correct formulation (works in Adjoint run). C Unfortunatly, produces many recomputations <== not used until it is fixed IF ( intime1.NE.DIC_ldRec(bi,bj) ) THEN #endif /* ALLOW_AUTODIFF */ C-- If the above condition is met then we need to read in C data for the period ahead and the period behind myTime. IF ( debugLevel.GE.debLevZero ) THEN _BEGIN_MASTER(myThid) WRITE(standardMessageUnit,'(A,I10,A,2(2I5,A))') & ' DIC_FIELDS_LOAD (Light), it=', myIter, & ' : Reading new data, i0,i1=', intime0, intime1, & ' (prev=', intimeP, DIC_ldRec(bi,bj), ' )' _END_MASTER(myThid) ENDIF _BARRIER IF ( DIC_windFile .NE. ' ' ) THEN CALL READ_REC_XY_RS( DIC_windFile,dicwind0,intime0, & myIter,myThid ) CALL READ_REC_XY_RS( DIC_windFile,dicwind1,intime1, & myIter,myThid ) ENDIF C------ Shropshire--------- New Code----------- IF ( light_file .NE. ' ' ) THEN CALL READ_REC_XY_RS( light_file,light_field1,intime0, & myIter,myThid ) CALL READ_REC_XY_RS( light_file,light_field2,intime1, & myIter,myThid ) print*, '---------- Reading In Light Field ------------' ENDIF C-----Shropshire--------- End new code ------- IF ( DIC_atmospFile .NE. ' ' ) THEN CALL READ_REC_XY_RS( DIC_atmospFile,atmosp0,intime0, & myIter,myThid ) CALL READ_REC_XY_RS( DIC_atmospFile,atmosp1,intime1, & myIter,myThid ) ENDIF IF ( DIC_silicaFile .NE. ' ' ) THEN CALL READ_REC_XY_RS( DIC_silicaFile,silica0,intime0, & myIter,myThid ) CALL READ_REC_XY_RS( DIC_silicaFile,silica1,intime1, & myIter,myThid ) ENDIF IF ( DIC_iceFile .NE. ' ' ) THEN CALL READ_REC_XY_RS( DIC_iceFile,ice0,intime0, & myIter,myThid ) CALL READ_REC_XY_RS( DIC_iceFile,ice1,intime1, & myIter,myThid ) ENDIF #ifdef READ_PAR IF ( DIC_parFile .NE. ' ' ) THEN CALL READ_REC_XY_RS( DIC_parFile,par0,intime0, & myIter,myThid ) CALL READ_REC_XY_RS( DIC_parFile,par1,intime1, & myIter,myThid ) ENDIF #endif #ifdef LIGHT_CHL C-- Load chlorophyll climatology data, unit for chlorophyll : mg/m3 IF ( DIC_chlaFile .NE. ' ' ) THEN CALL READ_REC_XY_RS( DIC_chlaFile,chlinput,1, & myIter,myThid ) ENDIF #endif #ifdef ALLOW_FE IF ( DIC_ironFile .NE. ' ' ) THEN CALL READ_REC_XY_RS( DIC_ironFile,feinput0,intime0, & myIter,myThid ) CALL READ_REC_XY_RS( DIC_ironFile,feinput1,intime1, & myIter,myThid ) ENDIF #endif C-- fill-in overlap after loading temp arrays: _EXCH_XY_RS(dicwind0, myThid ) _EXCH_XY_RS(dicwind1, myThid ) C----- Shropshire------New Code -------- _EXCH_XY_RS(light_field1, myThid ) _EXCH_XY_RS(light_field2, myThid ) C----- Shropshire------End New Code--------- _EXCH_XY_RS(atmosp0, myThid ) _EXCH_XY_RS(atmosp1, myThid ) _EXCH_XY_RS(silica0, myThid ) _EXCH_XY_RS(silica1, myThid ) _EXCH_XY_RS(ice0, myThid ) _EXCH_XY_RS(ice1, myThid ) #ifdef READ_PAR _EXCH_XY_RS(par0, myThid ) _EXCH_XY_RS(par1, myThid ) #endif #ifdef LIGHT_CHL _EXCH_XY_RS(chlinput, myThid ) #endif #ifdef ALLOW_FE _EXCH_XY_RS(feinput0, myThid ) _EXCH_XY_RS(feinput1, myThid ) #endif C- save newly loaded time-record DO bj = myByLo(myThid), myByHi(myThid) DO bi = myBxLo(myThid), myBxHi(myThid) DIC_ldRec(bi,bj) = intime1 ENDDO ENDDO C- end if-bloc (time to load new fields) ENDIF DO bj = myByLo(myThid), myByHi(myThid) DO bi = myBxLo(myThid), myBxHi(myThid) IF ( DIC_windFile .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx WIND(i,j,bi,bj) = bWght*dicwind0(i,j,bi,bj) & + aWght*dicwind1(i,j,bi,bj) ENDDO ENDDO ENDIF C----- Shropshire------New Code -------- IF ( light_file .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx light(i,j,bi,bj) = bWght*light_field1(i,j,bi,bj) & + aWght*light_field2(i,j,bi,bj) ENDDO ENDDO ENDIF C Note: Rivers are not read in like the light fields are because C we only have river once a month and light fields we have daily. C----- Shropshire------End New Code--------- C calculate piston velocity C QQ: note - we should have wind speed variance in here C QQ also need to check units, and conversion factors c pisvel(i,j,bi,bj) =0.337*wind(i,j,bi,bj)**2/3.6d5 !QQQQ #ifndef USE_PLOAD IF ( DIC_atmospFile .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx AtmosP(i,j,bi,bj) = bWght*atmosp0(i,j,bi,bj) & + aWght*atmosp1(i,j,bi,bj) ENDDO ENDDO ENDIF #endif IF ( DIC_silicaFile .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx SILICA(i,j,bi,bj) = bWght*silica0(i,j,bi,bj) & + aWght*silica1(i,j,bi,bj) ENDDO ENDDO ENDIF IF ( DIC_iceFile .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx fIce(i,j,bi,bj) = bWght*ice0(i,j,bi,bj) & + aWght*ice1(i,j,bi,bj) ENDDO ENDDO ENDIF #ifdef READ_PAR IF ( DIC_parFile .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx PAR(i,j,bi,bj) = bWght*par0(i,j,bi,bj) & + aWght*par1(i,j,bi,bj) ENDDO ENDDO ELSE WRITE(msgBuf,'(2A)') & ' DIC_FIELDS_LOAD: You need to provide ', & ' a file if you want to use READ_PAR' CALL PRINT_ERROR( msgBuf, myThid ) STOP 'ABNORMAL END: S/R DIC_FIELDS_LOAD' ENDIF #endif #ifdef LIGHT_CHL IF ( DIC_chlaFile .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx CHL(i,j,bi,bj) = chlinput(i,j,bi,bj) ENDDO ENDDO ENDIF #endif #ifdef ALLOW_FE IF ( DIC_ironFile .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx InputFe(i,j,bi,bj) = bWght*feinput0(i,j,bi,bj) & + aWght*feinput1(i,j,bi,bj) ENDDO ENDDO ENDIF #endif ENDDO ENDDO C endif for DIC_forcingCycle ENDIF C -------Shropshire------New Code -------------: C Note: Blow just copies code above IF ( DIC_forcingCycle2.GT.0. _d 0 ) THEN C-- Now calculate whether it is time to update the forcing arrays CALL GET_PERIODIC_INTERVAL( O intimeP, intime0, intime1, bWght, aWght, I DIC_forcingCycle2, DIC_forcingPeriod2, I deltaTClock, myTime, myThid ) bi = myBxLo(myThid) bj = myByLo(myThid) #ifdef ALLOW_DEBUG IF ( debugLevel.GE.debLevB ) THEN _BEGIN_MASTER(myThid) WRITE(standardMessageUnit,'(A,I10,A,4I5,A,2F14.10)') & ' DIC_FIELDS_LOAD (River),', myIter, & ' : iP,iLd,i0,i1=', intimeP,DIC_ldRec2(bi,bj), intime0,intime1, & ' ; Wght=', bWght, aWght _END_MASTER(myThid) ENDIF #endif /* ALLOW_DEBUG */ #ifdef ALLOW_AUTODIFF C- assuming that we call S/R DIC_FIELDS_LOAD at each time-step and C with increasing time, this will catch when we need to load new C records; C But with Adjoint run, this is not always the case => might end-up C using C the wrong time-records IF ( intime0.NE.intimeP .OR. myIter.EQ.nIter0 ) THEN #else /* ALLOW_AUTODIFF */ C- Make no assumption on sequence of calls to DIC_FIELDS_LOAD ; C This is the correct formulation (works in Adjoint run). C Unfortunatly, produces many recomputations <== not used until it C is fixed IF ( intime1.NE.DIC_ldRec2(bi,bj) ) THEN #endif /* ALLOW_AUTODIFF */ C-- If the above condition is met then we need to read in C data for the period ahead and the period behind myTime. IF ( debugLevel.GE.debLevZero ) THEN _BEGIN_MASTER(myThid) WRITE(standardMessageUnit,'(A,I10,A,2(2I5,A))') & ' DIC_FIELDS_LOAD (River), it=', myIter, & ' : Reading new data, i0,i1=', intime0, intime1, & ' (prev=', intimeP, DIC_ldRec2(bi,bj), ' )' _END_MASTER(myThid) ENDIF _BARRIER IF ( riv_nutr_file .NE. ' ' ) THEN CALL READ_REC_XY_RS( riv_nutr_file,river_field1,intime0, & myIter,myThid ) CALL READ_REC_XY_RS( riv_nutr_file,river_field2,intime1, & myIter,myThid ) print*, '----- Reading River Nutrient Fields ---------' ENDIF IF ( riv_dis_file .NE. ' ' ) THEN CALL READ_REC_XY_RS( riv_dis_file,river_field3,intime0, & myIter,myThid ) CALL READ_REC_XY_RS( riv_dis_file,river_field4,intime1, & myIter,myThid ) print*, '----- Reading River Discharge Fields ---------' ENDIF C-- fill-in overlap after loading temp arrays: _EXCH_XY_RS(river_field1, myThid ) _EXCH_XY_RS(river_field2, myThid ) _EXCH_XY_RS(river_field3, myThid ) _EXCH_XY_RS(river_field4, myThid ) C- save newly loaded time-record DO bj = myByLo(myThid), myByHi(myThid) DO bi = myBxLo(myThid), myBxHi(myThid) DIC_ldRec2(bi,bj) = intime1 ENDDO ENDDO C- end if-bloc (time to load new fields) ENDIF IF ( riv_nutr_file .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx riv_nutr(i,j,bi,bj) = bWght*river_field1(i,j,bi,bj) & + aWght*river_field2(i,j,bi,bj) ENDDO ENDDO ENDIF IF ( riv_dis_file .NE. ' ' ) THEN DO j=1-OLy,sNy+OLy DO i=1-OLx,sNx+OLx riv_dis(i,j,bi,bj) = bWght*river_field3(i,j,bi,bj) & + aWght*river_field4(i,j,bi,bj) ENDDO ENDDO ENDIF C endif for DIC_forcingCycle ENDIF C------- Shropshire----- END New code --------- #endif /* ALLOW_DIC */ RETURN END