!=============================================================================== ! SVN $Id: datm_physTN460.F90 4244 2007-05-02 21:04:16Z kauff $ ! SVN $URL: https://svn-ccsm-models.cgd.ucar.edu/datm7/trunk_tags/datm7_070531/datm_physTN460.F90 $ !=============================================================================== !BOP =========================================================================== ! ! !MODULE: datm_physTN460 - do physics necessary for TN460 mode ! ! !DESCRIPTION: ! Do physics necessary to advance the model state. Assumes a standard set ! of data as described in NCAR Tech Note #460 is provided as input. ! Assumes all reading of data stream data, mapping, filling, and ! time-interpolation has already been done. ! ! !REVISION HISTORY: ! 2007-Apr-18 - B. Kauffman - mods to support parallel datm7 ! 2005-Nov-02 - B. Kauffman - split code off from datm_phys.F90 ! 2005-Jul-15 - B. Kauffman - functional TN460 NFY mode ! ! !INTERFACE: ------------------------------------------------------------------ module datm_physTN460 ! !USES: use shr_sys_mod ! shared system calls use shr_mpi_mod ! mpi wrappers use shr_date_mod ! date data type and methods use shr_stream_mod ! stream data type and methods use shr_string_mod ! string methods use dshr_kind ! kinds for strong typing use dshr_const ! constants use dshr_domain ! domain data type and methods use dshr_bundle ! bundle data type and methods use dshr_map ! map data type and methods use datm_data ! holds global/public data use datm_control ! control flags/logic implicit none private ! !PUBLIC TYPES: ! none ! !PUBLIC MEMBER FUNCTIONS: public :: datm_physTN460_adv ! advance the model state public :: datm_physTN460_setDebug ! set internal debug level public :: datm_physTN460_getDebug ! get internal debug level ! !PUBLIC DATA MEMBERS: ! no public data members !EOP integer(IN) :: debug = 0 !=============================================================================== contains !=============================================================================== !=============================================================================== !BOP =========================================================================== ! ! !IROUTINE: datm_physTN460_adv -- do phsics as per NCAR Tech Note 460 ! ! !DESCRIPTION: ! Do physics necessary to advance the model state. Assumes a standard set ! of data as described in NCAR Tech Note #460 is provided as input. ! Assumes all reading of data stream data, mapping, filling, and ! time-interpolation has already been done. ! ! !REVISION HISTORY: ! 2005-Jul-11 - B. Kauffman - first version ! ! !INTERFACE: ------------------------------------------------------------------ subroutine datm_physTN460_adv(bunS2M,bunC2M,bunM2M,bunM2C,rc) use shr_file_mod ! file get/put methods use shr_ncread_mod ! netCDF reading methods use shr_map_mod ! mapping routines implicit none ! !INPUT/OUTPUT PARAMETERS: type(dshr_bundle_bundleType),intent(in) :: bunS2M ! ptr to s->m bundle type(dshr_bundle_bundleType),intent(in) :: bunC2M ! ptr to c->m bundle type(dshr_bundle_bundleType),intent(inout) :: bunM2M ! ptr to m->m bundle type(dshr_bundle_bundleType),intent(out) :: bunM2C ! ptr to m->c bundle integer(IN),optional ,intent(out) :: rc ! return code !EOP !--- local --- integer(IN) ,save :: i,j,n ! generic loop index integer(IN) ,save :: ni,nj ! dimensions of model grid integer(IN) ,save :: rCode ! return code logical ,save :: firstCall = .true. ! flags 1st envocation of this routine !--- pointers to bunM2C data: this data must be set --- real(R8),pointer,save :: z (:,:) ! bottom layer height real(R8),pointer,save :: u (:,:) ! u-velocity real(R8),pointer,save :: v (:,:) ! v-velocity real(R8),pointer,save :: tbot (:,:) ! bottom temperature real(R8),pointer,save :: shum (:,:) ! specific humidity real(R8),pointer,save :: ptem (:,:) ! potential temperature real(R8),pointer,save :: dens (:,:) ! density real(R8),pointer,save :: pslv (:,:) ! sea-level surface pressure real(R8),pointer,save :: pbot (:,:) ! bottom pressure real(R8),pointer,save :: rainc (:,:) ! net convective rain real(R8),pointer,save :: rainl (:,:) ! net large-scale rain real(R8),pointer,save :: snowc (:,:) ! net convective snow real(R8),pointer,save :: snowl (:,:) ! net large-scale snow real(R8),pointer,save :: lwdn (:,:) ! long-wave down real(R8),pointer,save :: swnet (:,:) ! short-wave, net real(R8),pointer,save :: swndr (:,:) ! short-wave down, near-inf, direct real(R8),pointer,save :: swndf (:,:) ! short-wave down, near-inf, diffuse real(R8),pointer,save :: swvdr (:,:) ! short-wave down, visable, direct real(R8),pointer,save :: swvdf (:,:) ! short-wave down, visable, diffuse real(R8),pointer,save :: co2prog(:,:) ! flux: co2, prognostic real(R8),pointer,save :: co2diag(:,:) ! flux: co2, diagnostic !--- pointers to bunS2M data not sent to cpl (used to compute bunM2C data) --- real(R8),pointer,save :: prec (:,:) ! net precip (snow+rain) real(R8),pointer,save :: swdn (:,:) ! short-wave down real(R8),pointer,save :: swup (:,:) ! short-wave up !--- pointers to bunC2M data (used to compute bunM2C data) --- real(R8),pointer,save :: anidr (:,:) ! albedo, near-inf, direct real(R8),pointer,save :: anidf (:,:) ! albedo, near-inf, diffuse real(R8),pointer,save :: avsdr (:,:) ! albedo, visable, direct real(R8),pointer,save :: avsdf (:,:) ! albedo, visable, diffuse !--- other data used read in and regrid data correction factors --- character(CL) ,save :: fileName = "unset" ! file name string real(R8) ,allocatable ,save :: windFactor(:,:) ! wind adjustment factors real(R8) ,allocatable ,save :: winddFactor(:,:) ! wind direction factor real(R8) ,allocatable ,save :: qsatFactor(:,:) ! rel humidty adjustment factors !--- data that describes the model domain --- real(R8) ,allocatable ,save :: yc (:,:) ! latitude of grid cell center !--- other data used to compute missing or incomplete bunM2C data --- real(R8) ,save :: cosFactor ! cosine factor real(R8) ,save :: factor ! generic/temporary correction factor real(R8) ,save :: avg_alb ! average albedo real(R8) ,save :: tMin ! minimum temperature real(R8),parameter :: tKFrz = dshr_const_tKFrz real(R8),parameter :: degtorad = SHR_CONST_PI/180.0_R8 real(R8),parameter :: avg_c0 = 61.846_R8 real(R8),parameter :: avg_c1 = 1.107_R8 real(R8),parameter :: amp_c0 = -21.841_R8 real(R8),parameter :: amp_c1 = -0.447_R8 real(R8),parameter :: phs_c0 = 0.298_R8 real(R8),parameter :: dLWarc = -5.000_R8 real(R8) ,save :: dTarc(12) data dTarc / 0.49_R8, 0.06_R8,-0.73_R8, -0.89_R8,-0.77_R8,-1.02_R8, & & -1.99_R8,-0.91_R8, 1.72_R8, 2.30_R8, 1.81_R8, 1.06_R8/ integer(IN) ,save :: eDay ! elapsed days real(R8) ,save :: rDay ! elapsed days, including fraction integer(IN) ,save :: year,month,day,sec !--- formats --- character(*),parameter :: subName = "(datm_physTN460_adv) " character(*),parameter :: F00 = "('(datm_physTN460_adv) ',4a) " character(*),parameter :: F01 = "('(datm_physTN460_adv) ',a,2i5)" character(*),parameter :: F02 = "('(datm_physTN460_adv) ',a,6e12.3)" character(*),parameter :: F03 = "('(datm_physTN460_adv) ',a,i8.8,i6,6es12.4)" !--- function statement --- real(R8) :: qSat,T qSat(T) = 640380.0_R8/ exp(5107.4_R8/T) !------------------------------------------------------------------------------- ! !------------------------------------------------------------------------------- if (present(rc)) rc=0 !---------------------------------------------------------------------------- ! one-time initializations !---------------------------------------------------------------------------- if (firstCall) then write(6,F00) "assuming input data stream is TN460 data" !------------------------------------------------------- ! get local dims and yc data !------------------------------------------------------- call dshr_domain_getdims(datm_data_domain,ni,nj) allocate(yc(ni,nj)) call dshr_domain_getData(datm_data_domain,yc,"lat" ) !------------------------------------------------------- ! read in correction factors for ncep wind & qsat !------------------------------------------------------- allocate(windFactor(ni,nj),winddFactor(ni,nj),qsatFactor(ni,nj)) fileName = datm_control_tn460CFfn call datm_physTN460getFactors(fileName,windFactor,winddFactor, & & qsatFactor) !------------------------------------------------------- ! assign pointers to bundle data !------------------------------------------------------- !--- assign pointers to bunM2C data --- call dshr_bundle_assignPtr(bunM2C,"z" ,z ,rCode) call dshr_bundle_assignPtr(bunM2C,"u" ,u ,rCode) call dshr_bundle_assignPtr(bunM2C,"v" ,v ,rCode) call dshr_bundle_assignPtr(bunM2C,"tbot" ,tbot ,rCode) call dshr_bundle_assignPtr(bunM2C,"ptem" ,ptem ,rCode) call dshr_bundle_assignPtr(bunM2C,"shum" ,shum ,rCode) call dshr_bundle_assignPtr(bunM2C,"dens" ,dens ,rCode) call dshr_bundle_assignPtr(bunM2C,"pslv" ,pslv ,rCode) call dshr_bundle_assignPtr(bunM2C,"pbot" ,pbot ,rCode) call dshr_bundle_assignPtr(bunM2C,"rainc" ,rainc ,rCode) call dshr_bundle_assignPtr(bunM2C,"rainl" ,rainl ,rCode) call dshr_bundle_assignPtr(bunM2C,"snowc" ,snowc ,rCode) call dshr_bundle_assignPtr(bunM2C,"snowl" ,snowl ,rCode) call dshr_bundle_assignPtr(bunM2C,"lwdn" ,lwdn ,rCode) call dshr_bundle_assignPtr(bunM2C,"swnet" ,swnet ,rCode) call dshr_bundle_assignPtr(bunM2C,"swndr" ,swndr ,rCode) call dshr_bundle_assignPtr(bunM2C,"swndf" ,swndf ,rCode) call dshr_bundle_assignPtr(bunM2C,"swvdr" ,swvdr ,rCode) call dshr_bundle_assignPtr(bunM2C,"swvdf" ,swvdf ,rCode) call dshr_bundle_assignPtr(bunM2C,"co2prog",co2prog ,rCode) call dshr_bundle_assignPtr(bunM2C,"co2diag",co2diag ,rCode) !--- assign pointers to bunS2M data not sent to cpl --- call dshr_bundle_assignPtr(bunS2M,"swdn" ,swdn ,rCode) call dshr_bundle_assignPtr(bunS2M,"swup" ,swup ,rCode) call dshr_bundle_assignPtr(bunS2M,"prec" ,prec ,rCode) !--- assign pointers to bunC2M data (data from cpl )--- call dshr_bundle_assignPtr(bunC2M,"anidr" ,anidr ,rCode) call dshr_bundle_assignPtr(bunC2M,"anidf" ,anidf ,rCode) call dshr_bundle_assignPtr(bunC2M,"avsdr" ,avsdr ,rCode) call dshr_bundle_assignPtr(bunC2M,"avsdf" ,avsdf ,rCode) end if !---------------------------------------------------------------------------- ! default bunM2C values (may get clobbered by bunS2M data) !---------------------------------------------------------------------------- co2prog(:,:) = 0.0_R8 co2diag(:,:) = 0.0_R8 !---------------------------------------------------------------------------- ! copy all matching fields from bunS2M into bunM2C !---------------------------------------------------------------------------- call dshr_bundle_copyFields(bunS2M,bunM2C,rc) !---------------------------------------------------------------------------- ! compute/correct/alter bunM2C fields as necessary !---------------------------------------------------------------------------- !--- phase-of-year factor used in air temperature correction --- call shr_date_getEDay(datm_data_date,eDay,sec) rDay = mod(eday,365) + real(sec)/dshr_const_cDay cosFactor = cos( (2.0_R8*dshr_const_pi*rDay)/365 - phs_c0) !--- get month for use wrt dTarc --- call shr_date_getYMD (datm_data_date,year,month,day,sec) do j=1,nj !dir$ concurrent !cdir nodep do i=1,ni !------------------------------------------------------------------------- ! STATE DATA !------------------------------------------------------------------------- !--- atm bottom layer height --- z(i,j) = 10.0_R8 !--- correction to NCEP winds based on QSCAT --- ! correct both magnitude and direction u(i,j) = u(i,j)*windFactor(i,j) v(i,j) = v(i,j)*windFactor(i,j) u(i,j) = u(i,j)*cos(winddFactor(i,j)*degtorad)- & & v(i,j)*sin(winddFactor(i,j)*degtorad) v(i,j) = u(i,j)*sin(winddFactor(i,j)*degtorad)+ & & v(i,j)*cos(winddFactor(i,j)*degtorad) !--- density, tbot, & pslv taken directly from input stream, set pbot --- ! dens(i,j) = ! tbot(i,j) = ! pslv(i,j) = pbot(i,j) = pslv(i,j) !--- correction to NCEP Arctic & Antarctic air T & potential T --- if ( yc(i,j) < -60.0_R8 ) then tMin = (avg_c0 + avg_c1*yc(i,j)) + (amp_c0 + amp_c1*yc(i,j))*cosFactor + tKFrz tbot(i,j) = max(tbot(i,j), tMin) else if ( yc(i,j) > 60.0_R8 ) then factor = MIN(1.0_R8, 0.1_R8*(yc(i,j)-60.0_R8) ) tbot(i,j) = tbot(i,j) + factor * dTarc(month) endif ptem(i,j) = tbot(i,j) !--- correction to NCEP relative humidity for heat budget balance --- ! from TN460: ! shum(i,j) = shum(i,j) - (qsat(tbot(i,j))*qsatFactor(i,j)/(100.0_R8*dens(i,j))) ! new Large&Yeager2007: shum(i,j) = shum(i,j) + qsatFactor(i,j) !------------------------------------------------------------------------- ! PRECIPITATION DATA !------------------------------------------------------------------------- prec(i,j) = prec(i,j)/86400.0_R8 ! convert mm/day to kg/m^2/s ! enhance GCGCS precipitation (only correct satellite products, not Serreze Arctic data) if ( yc(i,j) < 58. ) then prec(i,j) = prec(i,j)*1.14168_R8 endif if ( yc(i,j) >= 58. .and. yc(i,j) < 68. ) then factor = MAX(0.0_R8, 1.0_R8 - 0.1_R8*(yc(i,j)-58.0_R8) ) prec(i,j) = prec(i,j)*(factor*(1.14168_R8 - 1.0_R8) + 1.0_R8) endif if (tbot(i,j) < tKFrz ) then ! assign precip to rain/snow components rainc(i,j) = 0.0_R8 rainl(i,j) = 0.0_R8 snowc(i,j) = 0.0_R8 snowl(i,j) = prec(i,j) else rainc(i,j) = 0.0_R8 rainl(i,j) = prec(i,j) snowc(i,j) = 0.0_R8 snowl(i,j) = 0.0_R8 endif !------------------------------------------------------------------------- ! RADIATION DATA !------------------------------------------------------------------------- !--- fabricate required swdn components from net swdn --- swvdr(i,j) = swdn(i,j)*(0.28_R8) swndr(i,j) = swdn(i,j)*(0.31_R8) swvdf(i,j) = swdn(i,j)*(0.24_R8) swndf(i,j) = swdn(i,j)*(0.17_R8) !--- compute net short-wave based on albedo from cpl --- if (.not. firstCall) then !--- normal caluclation: net short-wave = swdn*(1-albedo) --- ! avg_alb = ( avsdr(i,j) + anidr(i,j) +avsdf(i,j) + anidf(i,j) ) * 0.25_R8 ! LY2007: avg_alb = ( 0.069 - 0.011*cos(2.0_R8*yc(i,j)*degtorad ) ) swnet(i,j) = swdn(i,j)*(1.0_R8 - avg_alb) else !--- alternate calculation: net short-wave = swdn - swup --- !--- albedo from cpl are not available on first call, but swup is --- swnet(i,j) = swdn(i,j) - swup(i,j) endif !--- corrections to GISS swdn for heat budget balancing --- ! TN460: factor = 1.0_R8 if ( -60.0_R8 < yc(i,j) .and. yc(i,j) < -50.0_R8 ) then factor = 1.0_R8 - (yc(i,j) + 60.0_R8)*(0.05_R8/10.0_R8) else if ( -50.0_R8 < yc(i,j) .and. yc(i,j) < 30.0_R8 ) then factor = 0.95_R8 else if ( 30.0_R8 < yc(i,j) .and. yc(i,j) < 40._R8 ) then factor = 1.0_R8 - (40.0_R8 - yc(i,j))*(0.05_R8/10.0_R8) endif ! LY2007: ! factor = 0.95_R8 swnet(i,j) = swnet(i,j)*factor swvdr(i,j) = swvdr(i,j)*factor swndr(i,j) = swndr(i,j)*factor swvdf(i,j) = swvdf(i,j)*factor swndf(i,j) = swndf(i,j)*factor !--- correction to GISS lwdn in Arctic --- if ( yc(i,j) > 60._R8 ) then factor = MIN(1.0_R8, 0.1_R8*(yc(i,j)-60.0_R8) ) lwdn(i,j) = lwdn(i,j) + factor * dLWarc endif end do end do firstCall = .false. end subroutine datm_physTN460_adv !=============================================================================== !=============================================================================== subroutine datm_physTN460getFactors(fileName,windFactor,winddFactor,qsatFactor) use shr_map_mod use shr_ncread_mod implicit none !--- arguments --- character(*),intent(in) :: fileName ! file name string real(R8) ,intent(inout) :: windFactor(:,:) ! wind adjustment factors real(R8) ,intent(inout) :: winddFactor(:,:) ! wind adjustment factors real(R8) ,intent(inout) :: qsatFactor(:,:) ! rel humidty adjustment factors !--- data that describes the local model domain --- integer(IN) :: ni0,nj0 ! dimensions of global bundle0 integer(IN) :: ni1,nj1,nf1 ! dimensions of global bundle1 integer(IN) :: i, j ! generic indicies type(dshr_domain_domainType) :: domain0 ! temp global domain (data domain) type(dshr_bundle_bundleType) :: bundle0 ! temp global bundle (data domain) type(dshr_bundle_bundleType) :: bundle1 ! temp global bundle (model domain) type(dshr_bundle_bundleType) :: bundle2 ! temp local bundle (model domain) logical :: doMap ! flags the need to map bundle0 -> bundle1 !--- temp arrays for data input --- real(R8) ,allocatable :: tempR4D(:,:,:,:) ! 4D data array real(R8) ,allocatable :: tempR2D(:,:) ! 2D data array integer(IN),allocatable :: tempI4D(:,:,:,:) ! 4D data array real(R8) ,pointer :: fldPtr (:,:) ! pointer to 2D data array real(R8) ,pointer :: x0 (:,:) ! pointer to 2D data array (longitude) real(R8) ,pointer :: x1 (:,:) ! pointer to 2D data array (longitude) !--- formats --- character(*),parameter :: subName = '(datm_physTN460getFactors) ' character(*),parameter :: F00 = "('(datm_physTN460getFactors) ',4a) " character(*),parameter :: F01 = "('(datm_physTN460getFactors) ',a,2i5)" character(*),parameter :: F02 = "('(datm_physTN460getFactors) ',a,6e12.3)" !------------------------------------------------------------------------------- ! Notes: ! bundle0: global bundle, data domain, on root PE only ! bundle1: global bundle, model domain, bcast this from root PE to all others ! bundle2: local bundle, model domain, extracted from bundle1 on all PEs !------------------------------------------------------------------------------- write(6,F00) "read in correction factors for ncep wind & qsat" !---------------------------------------------------------------------------- ! read in and map global correction factors !---------------------------------------------------------------------------- if (datm_data_commPID == 0) then !--- determine input file name --- write(6,F00) "reading correction factor data from file: ",trim(fileName) call shr_sys_flush(6) !--- verify necessary data is in input file --- if ( .not. shr_ncread_varExists(fileName,'lat' ) & .or. .not. shr_ncread_varExists(fileName,'lon' ) & .or. .not. shr_ncread_varExists(fileName,'mask' ) & .or. .not. shr_ncread_varExists(fileName,'windFactor') & .or. .not. shr_ncread_varExists(fileName,'winddFactor') & .or. .not. shr_ncread_varExists(fileName,'qsatFactor') ) then write(6,F00) "ERROR: invalid correction factor data file" call shr_sys_abort(subName//"invalid correction factor data file") end if !--- allocate temp arrays for input data --- call shr_ncread_varDimSizes(fileName,"windFactor",ni0,nj0) write(6,F01) "correction factor data size: ni,nj=",ni0,nj0 call dshr_domain_create(domain0,"TempDomain0",ni0,nj0) !--- read domain data: lon --- if (debug>1) write(6,F00) "read data: lon" allocate(tempR4D(ni0,1,1,1)) allocate(tempR2D(ni0,nj0)) call shr_ncread_field4dG(fileName,'lon' ,rfld=tempR4D) do j=1,nj0 tempR2D(:,j) = tempR4D(:,1,1,1) end do call dshr_domain_putData (domain0, tempR2D, "lon") deallocate(tempR4D) deallocate(tempR2D) !--- read domain data: lat --- if (debug>1) write(6,F00) "read data: lat" allocate(tempR4D(nj0,1,1,1)) allocate(tempR2D(ni0,nj0)) call shr_ncread_field4dG(fileName,'lat' ,rfld=tempR4D) do i=1,ni0 tempR2D(i,:) = tempR4D(:,1,1,1) end do call dshr_domain_putData (domain0, tempR2D, "lat") deallocate(tempR4D) deallocate(tempR2D) !--- read domain mask--- if (debug>1) write(6,F00) "read data: mask" allocate(tempI4D(ni0,nj0,1,1)) allocate(tempR2D(ni0,nj0)) call shr_ncread_field4dG(fileName,'mask',ifld=tempI4D) tempR2D(:,:) = real(tempI4D(:,:,1,1),R8) call dshr_domain_putData(domain0, tempR2D, "mask") deallocate(tempI4D) deallocate(tempR2D) !--- create temp bundle0 --- if (debug>0) write(6,F00) "create tempBundle0 (factors on global data domain)" call dshr_bundle_create (bundle0,domain0,"tempBundle0", & & "windFactor:winddFactor:qsatFactor") call dshr_bundle_putDate(bundle0,00010101,0) ! give temp bundle a bogus data allocate(tempR4D(ni0,nj0,1,1)) allocate(tempR2D(ni0,nj0)) !--- read bundle data: wind factor --- write(6,F00) "read data: windFactor" call shr_ncread_field4dG(fileName,'windFactor',rfld=tempR4D) tempR2D(:,:) = tempR4D(:,:,1,1) call dshr_bundle_putField (bundle0,tempR2D,"windFactor") !--- read bundle data: wind dir factor --- write(6,F00) "read data: winddFactor" call shr_ncread_field4dG(fileName,'winddFactor',rfld=tempR4D) tempR2D(:,:) = tempR4D(:,:,1,1) call dshr_bundle_putField (bundle0,tempR2D,"winddFactor") !--- read bundle data: qsat factor --- write(6,F00) "read data: qsatFactor" call shr_ncread_field4dG(fileName,'qsatFactor',rfld=tempR4D) tempR2D(:,:) = tempR4D(:,:,1,1) call dshr_bundle_putField (bundle0,tempR2D,"qsatFactor") deallocate(tempR4D) deallocate(tempR2D) !--- create temp bundle1 --- if (debug>0) write(6,F00) "create tempBundle1 (factors on global model domain)" call dshr_bundle_create(bundle1,datm_data_domainG,"tempBundle1", & & "windFactor:winddFactor:qsatFactor") !--- determine if mapping is needed --- doMap = .false. call dshr_bundle_getDims(bundle1,ni1,nj1,nf1) if ( ni0/=ni1 .or. nj0/=nj1 ) then doMap = .true. else call dshr_domain_assignPtr(domain0 ,"lon",x0) call dshr_domain_assignPtr(datm_data_domainG,"lon",x1) if ( maxval(abs(x1(:,:)-x0(:,:))) > 0.01_R8 ) doMap = .true. call dshr_domain_assignPtr(domain0 ,"lat",x0) call dshr_domain_assignPtr(datm_data_domainG,"lat",x1) if ( maxval(abs(x1(:,:)-x0(:,:))) > 0.01_R8 ) doMap = .true. end if !--- map onto global model domain --- if (.not. doMap) then write(6,F00) "map (ID map = copy) factors to global model grid" call dshr_bundle_copyFields(bundle0,bundle1) else write(6,F00) "map factors to global model grid" call dshr_map_mapData(bundle0,bundle1,type=shr_map_fs_remap , & & algo=shr_map_fs_bilinear, & & mask=shr_map_fs_srcmask ) end if !--- sanity check --- if (debug>0) then call dshr_bundle_assignPtr(bundle0,"windFactor",fldPtr) write(6,F02) "unmapped min/max windFactor=",minval(fldPtr),maxval(fldPtr) call dshr_bundle_assignPtr(bundle1,"windFactor",fldPtr) write(6,F02) " mapped min/max windFactor=",minval(fldPtr),maxval(fldPtr) call dshr_bundle_assignPtr(bundle0,"winddFactor",fldPtr) write(6,F02) "unmapped min/max winddFactor=",minval(fldPtr),maxval(fldPtr) call dshr_bundle_assignPtr(bundle1,"winddFactor",fldPtr) write(6,F02) " mapped min/max winddFactor=",minval(fldPtr),maxval(fldPtr) call dshr_bundle_assignPtr(bundle0,"qsatFactor",fldPtr) write(6,F02) "unmapped min/max qsatFactor=",minval(fldPtr),maxval(fldPtr) call dshr_bundle_assignPtr(bundle1,"qsatFactor",fldPtr) write(6,F02) " mapped min/max qsatFactor=",minval(fldPtr),maxval(fldPtr) end if call dshr_domain_clean(domain0) call dshr_bundle_clean(bundle0) end if !---------------------------------------------------------------------------- ! scatter global data to local domains !---------------------------------------------------------------------------- if (debug>0) write(6,F00) "create tempBundle2 (factors on local model domain)" call dshr_bundle_create(bundle2,datm_data_domain ,"tempBundle2", & & "windFactor:winddFactor:qsatFactor") call dshr_bundle_scatter (bundle1,bundle2,datm_data_comm,subName) call dshr_bundle_getField(bundle2,windFactor,"windFactor") call dshr_bundle_getField(bundle2,winddFactor,"winddFactor") call dshr_bundle_getField(bundle2,qsatFactor,"qsatFactor") if (debug>1) write(6,F00) "clean up temp bundles" call dshr_bundle_clean(bundle2) if (datm_data_commPID == 0) call dshr_bundle_clean(bundle1) end subroutine datm_physTN460getFactors !=============================================================================== !BOP =========================================================================== ! ! !IROUTINE: datm_physTN460_setDebug -- Set local debug level ! ! !DESCRIPTION: ! Set local debug level, 0 = production ! \newline ! General Usage: call datm\_physTN460\_setDebug(1) ! ! !REVISION HISTORY: ! 2005-May-10 - B. Kauffman - first version ! ! !INTERFACE: ------------------------------------------------------------------ subroutine datm_physTN460_setDebug(level) implicit none ! !INPUT/OUTPUT PARAMETERS: integer,intent(in) :: level !EOP !--- formats --- character(*),parameter :: subName = "(datm_physTN460_setDebug) " character(*),parameter :: F01 = "('(datm_physTN460_setDebug) ',a,i4) " !------------------------------------------------------------------------------- ! !------------------------------------------------------------------------------- debug = level write(6,F01) "debug level reset to ",level end subroutine datm_physTN460_setDebug !=============================================================================== !BOP =========================================================================== ! ! !IROUTINE: datm_physTN460_getDebug -- return local internal debug level ! ! !DESCRIPTION: ! Get local debug level, 0 = production ! \newline ! General Usage: call datm\_physTN460\_getDebug(level) ! ! !REVISION HISTORY: ! 2005-May-10 - B. Kauffman - first version ! ! !INTERFACE: ------------------------------------------------------------------ subroutine datm_physTN460_getDebug(level) implicit none ! !INPUT/OUTPUT PARAMETERS: integer,intent(out) :: level !EOP !--- formats --- character(*),parameter :: subName = "(datm_physTN460_getDebug) " character(*),parameter :: F00 = "('(datm_physTN460_getDebug) ',a) " !------------------------------------------------------------------------------- ! !------------------------------------------------------------------------------- level = debug end subroutine datm_physTN460_getDebug !=============================================================================== !=============================================================================== end module datm_physTN460