C $Header: /u/gcmpack/MITgcm/pkg/bulk_force/bulkf_fields_load.F,v 1.17 2011/06/07 20:45:30 jmc Exp $
C $Name: $
#include "BULK_FORCE_OPTIONS.h"
C !ROUTINE: BULKF_FIELDS_LOAD
C !INTERFACE:
SUBROUTINE BULKF_FIELDS_LOAD( myTime, myIter, myThid )
C !DESCRIPTION:
C *==========================================================*
C | SUBROUTINE BULKF_FIELDS_LOAD
C | o Control reading of fields from external source.
C *==========================================================*
C | Bulk formula External source field loading routine.
C | This routine is called every time we want to
C | load a a set of external fields. The routine decides
C | which fields to load and then reads them in.
C | This routine needs to be customised for particular
C | experiments.
C | Notes
C | =====
C | Two-dimensional and three-dimensional I/O are handled in
C | the following way under MITgcmUV. A master thread
C | performs I/O using system calls. This threads reads data
C | into a temporary buffer. At present the buffer is loaded
C | with the entire model domain. This is probably OK for now
C | Each thread then copies data from the buffer to the
C | region of the proper array it is responsible for.
C | =====
C | Conversion of flux fields are described in FFIELDS.h
C *==========================================================*
C !USES:
IMPLICIT NONE
C === Global variables ===
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "BULKF_PARAMS.h"
#include "BULKF.h"
C !INPUT/OUTPUT PARAMETERS:
C === Routine arguments ===
C myTime :: Simulation time
C myIter :: Simulation timestep number
C myThid :: Thread no. that called this routine.
_RL myTime
INTEGER myIter
INTEGER myThid
#ifdef ALLOW_BULK_FORCE
C !LOCAL VARIABLES:
C === Local arrays in common block ===
C bulkfRec :: time-record currently loaded (in temp arrays *[1])
C tair[01] :: Temp. for air temperature
C qair[01] :: Temp. for air specific humidity
C solar[01] :: Temp. for incoming solar radition
C flwdwn[01] :: Temp. for downward longwave flux
C cloud[01] :: Temp for cloud
C wspeed[01] :: Temp. for wind speed
C uwind[01] :: Temp. for zonal wind speed
C vwind[01] :: Temp. for meridional wind speed
C rain[01] :: Temp. for rain
C runoff[01] :: Temp. for runoff
C snow[01] :: Temp for snow
C thAir[01] :: Temp for Air potential temp. in the BL [K]
c qnetch[01] :: Temp for qnet (cheating)
C empch[01] :: Temp for empmr (cheating)
C [01] :: End points for interpolation
COMMON /BULKF_FIELDS_LOAD_I/ bulkfRec
COMMON /BULKF_FIELDS_LOAD_RS/
& tair0, tair1, qair0, qair1,
& solar0, solar1, flwdwn0, flwdwn1,
& cloud0, cloud1, wspeed0, wspeed1,
& uwind0, uwind1, vwind0, vwind1,
& rain0, rain1, runoff0, runoff1,
& snow0, snow1,
#ifdef ALLOW_FORMULA_AIM
& thAir0, thAir1,
#endif
& qnetch0, qnetch1, empch0, empch1
INTEGER bulkfRec(nSx,nSy)
_RS tair0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS tair1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS qair0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS qair1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS solar0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS solar1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS flwdwn0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS flwdwn1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS cloud0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS cloud1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS wspeed0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS wspeed1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS uwind0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS uwind1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS vwind0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS vwind1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS rain0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS rain1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS runoff0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS runoff1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS snow0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS snow1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
#ifdef ALLOW_FORMULA_AIM
_RS thAir0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS thAir1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
#endif
_RS qnetch0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS qnetch1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS empch0 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
_RS empch1 (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
C === Local variables ===
C aWght, bWght :: Interpolation weights
_RL aWght,bWght
INTEGER intimeP, intime0, intime1
INTEGER bi, bj, i, j
CEOP
IF ( periodicExternalForcing ) THEN
C First call requires that we initialize everything to zero for safety
IF ( myIter .EQ. nIter0 ) THEN
DO bj=myByLo(myThid),myByHi(myThid)
DO bi=myBxLo(myThid),myBxHi(myThid)
bulkfRec(bi,bj) = 0
DO j = 1-Oly,sNy+Oly
DO i = 1-Olx,sNx+Olx
tair0(i,j,bi,bj) = 0. _d 0
tair1(i,j,bi,bj) = 0. _d 0
qair0(i,j,bi,bj) = 0. _d 0
qair1(i,j,bi,bj) = 0. _d 0
rain0(i,j,bi,bj) = 0. _d 0
rain1(i,j,bi,bj) = 0. _d 0
solar0(i,j,bi,bj) = 0. _d 0
solar1(i,j,bi,bj) = 0. _d 0
flwdwn0(i,j,bi,bj) = 0. _d 0
flwdwn1(i,j,bi,bj) = 0. _d 0
cloud0(i,j,bi,bj) = 0. _d 0
cloud1(i,j,bi,bj) = 0. _d 0
wspeed0(i,j,bi,bj) = 0. _d 0
wspeed1(i,j,bi,bj) = 0. _d 0
uwind0(i,j,bi,bj) = 0. _d 0
uwind1(i,j,bi,bj) = 0. _d 0
vwind0(i,j,bi,bj) = 0. _d 0
vwind1(i,j,bi,bj) = 0. _d 0
runoff0(i,j,bi,bj) = 0. _d 0
runoff1(i,j,bi,bj) = 0. _d 0
snow0(i,j,bi,bj) = 0. _d 0
snow1(i,j,bi,bj) = 0. _d 0
#ifdef ALLOW_FORMULA_AIM
thAir0(i,j,bi,bj) = 0. _d 0
thAir1(i,j,bi,bj) = 0. _d 0
#endif
qnetch0(i,j,bi,bj) = 0. _d 0
qnetch1(i,j,bi,bj) = 0. _d 0
empch0(i,j,bi,bj) = 0. _d 0
empch1(i,j,bi,bj) = 0. _d 0
ENDDO
ENDDO
ENDDO
ENDDO
ENDIF
C-- Now calculate whether it is time to update the forcing arrays
CALL GET_PERIODIC_INTERVAL(
O intimeP, intime0, intime1, bWght, aWght,
I externForcingCycle, externForcingPeriod,
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)')
& ' BULKF_FIELDS_LOAD,', myIter,
& ' : iP,iLd,i0,i1=', intimeP,bulkfRec(bi,bj), intime0,intime1,
& ' ; Wght=', bWght, aWght
_END_MASTER(myThid)
ENDIF
#endif /* ALLOW_DEBUG */
C- Make no assumption on sequence of calls to BULKF_FIELDS_LOAD ;
C This is the correct formulation (works in Adjoint run).
IF ( intime1.NE.bulkfRec(bi,bj) ) THEN
_BARRIER
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))')
& ' BULKF_FIELDS_LOAD, it=', myIter,
& ' : Reading new data, i0,i1=', intime0, intime1,
& ' (prev=', intimeP, bulkfRec(bi,bj), ' )'
_END_MASTER(myThid)
ENDIF
IF ( AirTempFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( AirTempFile,tair0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( AirTempFile,tair1,intime1,
& myIter,myThid )
ENDIF
IF ( AirHumidityFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( AirhumidityFile,qair0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( AirhumidityFile,qair1,intime1,
& myIter,myThid )
ENDIF
IF ( SolarFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( SolarFile,solar0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( SolarFile,solar1,intime1,
& myIter,myThid )
ENDIF
IF ( LongwaveFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( LongwaveFile,flwdwn0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( LongwaveFile,flwdwn1,intime1,
& myIter,myThid )
ENDIF
IF ( CloudFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( CloudFile,cloud0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( CloudFile,cloud1,intime1,
& myIter,myThid )
ENDIF
IF ( WSpeedFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( WSpeedFile,wspeed0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( WSpeedFile,wspeed1,intime1,
& myIter,myThid )
ENDIF
IF ( UwindFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( UWindFile,uwind0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( UWindFile,uwind1,intime1,
& myIter,myThid )
ENDIF
IF ( VwindFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( VWindFile,vwind0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( VWindFile,vwind1,intime1,
& myIter,myThid )
ENDIF
IF ( RainFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( RainFile,rain0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( RainFile,rain1,intime1,
& myIter,myThid )
ENDIF
IF ( RunoffFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( RunoffFile,runoff0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( RunoffFile,runoff1,intime1,
& myIter,myThid )
ENDIF
IF ( SnowFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( SnowFile,snow0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( SnowFile,snow1,intime1,
& myIter,myThid )
ENDIF
#ifdef ALLOW_FORMULA_AIM
IF ( airPotTempFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( airPotTempFile, thAir0, intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( airPotTempFile, thAir1, intime1,
& myIter,myThid )
ENDIF
#endif
IF ( QnetFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( QnetFile,qnetch0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( QnetFile,qnetch1,intime1,
& myIter,myThid )
ENDIF
IF ( EmPFile .NE. ' ' ) THEN
CALL READ_REC_XY_RS( EmpFile,empch0,intime0,
& myIter,myThid )
CALL READ_REC_XY_RS( EmpFile,empch1,intime1,
& myIter,myThid )
c IF ( convertEmP2rUnit.EQ.mass2rUnit ) THEN
C- EmPmR is now (after c59h) expressed in kg/m2/s (fresh water mass flux)
_BARRIER
DO bj = myByLo(myThid), myByHi(myThid)
DO bi = myBxLo(myThid), myBxHi(myThid)
DO j=1-Oly,sNy+Oly
DO i=1-Olx,sNx+Olx
empch0(i,j,bi,bj) = empch0(i,j,bi,bj)*rhoConstFresh
empch1(i,j,bi,bj) = empch1(i,j,bi,bj)*rhoConstFresh
ENDDO
ENDDO
ENDDO
ENDDO
c ENDIF
ENDIF
C note: no need for extra barrier here since EXCH contains threads synchronization
c _BARRIER
_EXCH_XY_RS(tair0 , myThid )
_EXCH_XY_RS(tair1 , myThid )
_EXCH_XY_RS(qair0 , myThid )
_EXCH_XY_RS(qair1 , myThid )
_EXCH_XY_RS(solar0, myThid )
_EXCH_XY_RS(solar1, myThid )
_EXCH_XY_RS(flwdwn0, myThid )
_EXCH_XY_RS(flwdwn1, myThid )
_EXCH_XY_RS(cloud0, myThid )
_EXCH_XY_RS(cloud1, myThid )
_EXCH_XY_RS(wspeed0, myThid )
_EXCH_XY_RS(wspeed1, myThid )
c _EXCH_XY_RS(uwind0, myThid )
c _EXCH_XY_RS(uwind1, myThid )
c _EXCH_XY_RS(vwind0, myThid )
c _EXCH_XY_RS(vwind1, myThid )
CALL EXCH_UV_AGRID_3D_RS( uwind0, vwind0, .TRUE., 1, myThid )
CALL EXCH_UV_AGRID_3D_RS( uwind1, vwind1, .TRUE., 1, myThid )
_EXCH_XY_RS(rain0, myThid )
_EXCH_XY_RS(rain1, myThid )
_EXCH_XY_RS(runoff0, myThid )
_EXCH_XY_RS(runoff1, myThid )
_EXCH_XY_RS(snow0 , myThid )
_EXCH_XY_RS(snow1 , myThid )
#ifdef ALLOW_FORMULA_AIM
IF ( useFluxFormula_AIM ) THEN
_EXCH_XY_RS( thAir0, myThid )
_EXCH_XY_RS( thAir1, myThid )
ENDIF
#endif
_EXCH_XY_RS(qnetch0, myThid )
_EXCH_XY_RS(qnetch1, myThid )
_EXCH_XY_RS(empch0, myThid )
_EXCH_XY_RS(empch1, myThid )
C- save newly loaded time-record
DO bj = myByLo(myThid), myByHi(myThid)
DO bi = myBxLo(myThid), myBxHi(myThid)
bulkfRec(bi,bj) = intime1
ENDDO
ENDDO
C-- end if-block for loading new time-records
ENDIF
C-- Interpolate TAIR, QAIR, RAIN, SOLAR
DO bj = myByLo(myThid), myByHi(myThid)
DO bi = myBxLo(myThid), myBxHi(myThid)
DO j=1-Oly,sNy+Oly
DO i=1-Olx,sNx+Olx
Cswdblkf -- QQQQQ check if tair is K or C ------
C -- dasilva data in C, ncep data in K ---
TAIR(i,j,bi,bj) = bWght*tair0(i,j,bi,bj)
& +aWght*tair1(i,j,bi,bj)
c +273.15 _d 0
Cswdblkf -- QQQQQ set to kg.kg??? ---
C -- dasilva data in g, ncep in kg ---
QAIR(i,j,bi,bj) =(bWght*qair0(i,j,bi,bj)
& +aWght*qair1(i,j,bi,bj) )
SOLAR(i,j,bi,bj) = bWght*solar0(i,j,bi,bj)
& +aWght*solar1(i,j,bi,bj)
FLWDWN(i,j,bi,bj) = bWght*flwdwn0(i,j,bi,bj)
& +aWght*flwdwn1(i,j,bi,bj)
CLOUD(i,j,bi,bj) = bWght*cloud0(i,j,bi,bj)
& +aWght*cloud1(i,j,bi,bj)
WSPEED(i,j,bi,bj) = bWght*wspeed0(i,j,bi,bj)
& +aWght*wspeed1(i,j,bi,bj)
UWIND(i,j,bi,bj) = bWght*uwind0(i,j,bi,bj)
& +aWght*uwind1(i,j,bi,bj)
VWIND(i,j,bi,bj) = bWght*vwind0(i,j,bi,bj)
& +aWght*vwind1(i,j,bi,bj)
RAIN(i,j,bi,bj) = bWght*rain0(i,j,bi,bj)
& +aWght*rain1(i,j,bi,bj)
RUNOFF(i,j,bi,bj) = bWght*runoff0(i,j,bi,bj)
& +aWght*runoff1(i,j,bi,bj)
c#ifdef ALLOW_THSICE
c SNOW(i,j,bi,bj) = bWght*snow0(i,j,bi,bj)
c & +aWght*snow1(i,j,bi,bj)
c#endif
Qnetch(i,j,bi,bj) = bWght*qnetch0(i,j,bi,bj)
& +aWght*qnetch1(i,j,bi,bj)
EmPch(i,j,bi,bj) = bWght*empch0(i,j,bi,bj)
& +aWght*empch1(i,j,bi,bj)
ENDDO
ENDDO
#ifdef ALLOW_FORMULA_AIM
IF ( useFluxFormula_AIM ) THEN
DO j=1-Oly,sNy+Oly
DO i=1-Olx,sNx+Olx
thAir(i,j,bi,bj) = bWght*thAir0(i,j,bi,bj)
& + aWght*thAir1(i,j,bi,bj)
ENDDO
ENDDO
ENDIF
#endif
ENDDO
ENDDO
C- jmc: not needed since local fields have been exchanged.
c _EXCH_XY_RL(tair , myThid )
C-- Diagnostics
c IF (myThid.EQ.1 .AND. myTime.LT.62208000.) THEN
c write(*,'(a,1p5e12.4,2i6,2e12.4)')
c & 'time,TAIR,QAIR,RAIN,SOLAR,i0,i1,a,b = ',
c & myTime,
c & TAIR(1,sNy,1,1),QAIR(1,sNy,1,1),
c & RAIN(1,sNy,1,1),SOLAR(1,sNy,1,1),
c & intime0,intime1,aWght,bWght
c write(*,'(a,1p4e12.4,2e12.4)')
c & 'time,tair0,tair1,TAIR = ',
c & myTime,
c & tair0(1,sNy,1,1),tair1(1,sNy,1,1),TAIR(1,sNy,1,1),
c & aWght,bWght
c ENDIF
C endif for periodicForcing
ENDIF
#endif /*ALLOW_BULK_FORCE*/
RETURN
END