subroutine control2model(xhat,sval,bval) !$$$ subprogram documentation block ! . . . . ! subprogram: control2model ! prgmmr: tremolet ! ! abstract: Converts control variable to physical space ! ! program history log: ! 2007-04-13 tremolet - initial code ! 2007-04-27 tremolet - multiply by sqrt(B) (from ckerror D. Parrish) ! 2008-12-04 todling - update interface to ckgcov; add tsen/p3d ! 2008-12-29 todling - add call to strong balance contraint ! 2009-04-21 derber - modify call to getuv to getuv(*,0) ! 2009-08-14 lueken - update documentation ! 2010-03-15 zhu - make changes to ckbcov, add assign_cs2array ! 2010-04-28 todling - update to use gsi_bundle ! 2010-05-31 todling - better consistency checks; add co/co2 ! - ready to bypass analysis of (any) meteorological fields ! 2010-06-15 todling - generalized handling of chemistry ! 2011-05-15 auligne/todling - generalized cloud handling ! 2011-07-12 zhu - add cw_to_hydro and cwhydromod ! 2011-12-14 mkim - changed clouds4crtm to clouds in metguess ! 2012-06-25 parrish - modify wbundle by adding motley variables to control vector ! so will be in form expected by new version of ckgcov which uses general_grid2sub. ! 2012-10-09 wgu - update interface to tbalance (fut2ps) ! 2013-10-25 todling - nullify work pointers ! 2013-10-28 todling - rename p3d to prse ! 2013-05-23 zhu - add ntclen and predt for aircraft temperature bias correction ! ! input argument list: ! xhat - Control variable ! sval - State variable ! bval - Bias predictors ! ! output argument list: ! sval - State variable ! bval - Bias predictors ! ! attributes: ! language: f90 ! machine: ! !$$$ end documentation block use kinds, only: r_kind,i_kind use control_vectors, only: control_vector use control_vectors, only: cvars3d,cvars2d,cvarsmd use control_vectors, only: nc2d,mvars use bias_predictors, only: predictors use gsi_4dvar, only: nsubwin, l4dvar, lsqrtb use gridmod, only: lat2,lon2,nsig use jfunc, only: nsclen,npclen,ntclen use berror, only: varprd,fpsproj,fut2ps use balmod, only: balance use cwhydromod, only: cw2hydro_tl use gsi_bundlemod, only: gsi_bundlecreate use gsi_bundlemod, only: gsi_gridcreate use gsi_bundlemod, only: gsi_grid use gsi_bundlemod, only: gsi_bundle use gsi_bundlemod, only: gsi_bundlegetpointer use gsi_bundlemod, only: gsi_bundlegetvar use gsi_bundlemod, only: gsi_bundledestroy use gsi_chemguess_mod, only: gsi_chemguess_get use gsi_metguess_mod, only: gsi_metguess_get use mpeu_util, only: getindex use constants, only: max_varname_length implicit none ! Declare passed variables type(control_vector), intent(in ) :: xhat type(gsi_bundle) , intent(inout) :: sval(nsubwin) type(predictors) , intent(inout) :: bval ! Declare local variables character(len=*),parameter:: myname ='control2model' real(r_kind),dimension(lat2,lon2,nsig) :: workst,workvp,workrh type(gsi_bundle) :: wbundle type(gsi_grid) :: grid integer(i_kind) :: ii,jj,i,ic,id,ngases,nclouds,istatus character(len=10),allocatable,dimension(:) :: gases character(len=max_varname_length),allocatable,dimension(:) :: clouds ! Note: The following does not aim to get all variables in ! the state and control vectors, but rather the ones ! explicitly needed by this routine. ! Declare required local state variables logical :: ls_u,ls_v,ls_prse,ls_q,ls_tsen,ls_tv,ls_ps,ls_ql,ls_qi integer(i_kind), parameter :: nsvars = 9 integer(i_kind) :: isps(nsvars) character(len=4), parameter :: mysvars(nsvars) = (/ & ! vars from SV needed here 'u ', 'v ', 'prse', 'q ', 'tsen', 'tv ', 'ps ','ql ', 'qi ' /) character(len=max_varname_length),allocatable,dimension(:) :: cvars2dpm ! names of 2d fields including ! motley vars (if any) real(r_kind),pointer,dimension(:,:) :: sv_ps=>NULL() real(r_kind),pointer,dimension(:,:) :: sv_sst=>NULL() real(r_kind),pointer,dimension(:,:,:) :: sv_u=>NULL() real(r_kind),pointer,dimension(:,:,:) :: sv_v=>NULL() real(r_kind),pointer,dimension(:,:,:) :: sv_prse=>NULL() real(r_kind),pointer,dimension(:,:,:) :: sv_q=>NULL() real(r_kind),pointer,dimension(:,:,:) :: sv_tsen=>NULL() real(r_kind),pointer,dimension(:,:,:) :: sv_tv=>NULL() real(r_kind),pointer,dimension(:,:,:) :: sv_oz=>NULL() real(r_kind),pointer,dimension(:,:,:) :: sv_rank3=>NULL() real(r_kind),pointer,dimension(:,:) :: sv_rank2=>NULL() logical :: do_balance,do_getprs_tl,do_normal_rh_to_q,do_tv_to_tsen,do_getuv,cw_to_hydro !****************************************************************************** if (.not.lsqrtb) then write(6,*)trim(myname),': assumes lsqrtb' call stop2(104) end if if (nsubwin/=1 .and. .not.l4dvar) then write(6,*)trim(myname),': error 3dvar',nsubwin call stop2(105) end if ! create extended list of 2d variable names to include motley variables. ! NOTE: if mvars=0, there are no motley variables so cvars2dpm=cvars2d. allocate(cvars2dpm(nc2d+mvars)) do i=1,nc2d cvars2dpm(i)=cvars2d(i) end do do i=1,mvars cvars2dpm(nc2d+i)=cvarsmd(i) end do call gsi_gridcreate(grid,lat2,lon2,nsig) ! Inquire about chemistry call gsi_chemguess_get('dim',ngases,istatus) if (ngases>0) then allocate(gases(ngases)) call gsi_chemguess_get('gsinames',gases,istatus) endif ! Since each internal vector of xhat has the same structure, pointers are ! the same independent of the subwindow jj call gsi_bundlegetpointer (sval(1),mysvars,isps,istatus) ls_u =isps(1)>0; ls_v =isps(2)>0; ls_prse=isps(3)>0 ls_q =isps(4)>0; ls_tsen=isps(5)>0; ls_tv =isps(6)>0 ls_ps =isps(7)>0; ls_ql =isps(8)>0; ls_qi =isps(9)>0 ! Define what to do depending on what's in SV and ! what's explictly needed in this routine do_balance =ls_tv.and.ls_ps do_getprs_tl =ls_ps.and.ls_tv .and.ls_prse do_normal_rh_to_q=ls_tv.and.ls_prse.and.ls_q do_tv_to_tsen =ls_tv.and.ls_q .and.ls_tsen do_getuv =ls_u .and.ls_v ! Inquire about cloud-variables cw_to_hydro=.false. call gsi_metguess_get ('clouds::3d',nclouds,istatus) if (nclouds>0) then allocate(clouds(nclouds)) call gsi_metguess_get ('clouds::3d',clouds,istatus) if (getindex(cvars3d,'cw')>0 .and. ls_ql .and. ls_qi) cw_to_hydro=.true. end if ! Loop over control steps do jj=1,nsubwin ! create an internal structure w/ the same vars as those in the control vector, including motley vars call gsi_bundlecreate (wbundle,grid,'control2model work',istatus,names2d=cvars2dpm,names3d=cvars3d) if(istatus/=0) then write(6,*) trim(myname), ': trouble creating work bundle' call stop2(999) endif ! Multiply by sqrt of background error (ckerror) ! ----------------------------------------------------------------------------- ! Apply sqrt of variance, as well as vertical & horizontal parts of background ! error call ckgcov(xhat%step(jj)%values(:),wbundle,size(xhat%step(jj)%values(:))) ! Get pointers to required state variables call gsi_bundlegetpointer (sval(jj),'u' ,sv_u, istatus) call gsi_bundlegetpointer (sval(jj),'v' ,sv_v, istatus) call gsi_bundlegetpointer (sval(jj),'ps' ,sv_ps, istatus) call gsi_bundlegetpointer (sval(jj),'prse',sv_prse,istatus) call gsi_bundlegetpointer (sval(jj),'tv' ,sv_tv, istatus) call gsi_bundlegetpointer (sval(jj),'tsen',sv_tsen,istatus) call gsi_bundlegetpointer (sval(jj),'q' ,sv_q , istatus) call gsi_bundlegetpointer (sval(jj),'oz' ,sv_oz , istatus) call gsi_bundlegetpointer (sval(jj),'sst' ,sv_sst, istatus) ! Copy variables from CV to SV call gsi_bundlegetvar ( wbundle, 'sf' , workst, istatus ) call gsi_bundlegetvar ( wbundle, 'vp' , workvp, istatus ) call gsi_bundlegetvar ( wbundle, 'q' , workrh, istatus ) call gsi_bundlegetvar ( wbundle, 't' , sv_tv, istatus ) call gsi_bundlegetvar ( wbundle, 'oz' , sv_oz, istatus ) call gsi_bundlegetvar ( wbundle, 'ps' , sv_ps, istatus ) call gsi_bundlegetvar ( wbundle, 'sst', sv_sst, istatus ) ! Same one-to-one map for chemistry-vars; take care of them together do ic=1,ngases id=getindex(cvars3d,gases(ic)) if (id>0) then call gsi_bundlegetpointer (sval(jj),gases(ic),sv_rank3,istatus) call gsi_bundlegetvar (wbundle, gases(ic),sv_rank3,istatus) endif id=getindex(cvars2d,gases(ic)) if (id>0) then call gsi_bundlegetpointer (sval(jj),gases(ic),sv_rank2,istatus) call gsi_bundlegetvar (wbundle, gases(ic),sv_rank2,istatus) endif enddo ! Balance equation if(do_balance) call balance(sv_tv,sv_ps,workst,workvp,fpsproj,fut2ps) ! ----------------------------------------------------------------------------- ! Get 3d pressure if(do_getprs_tl) call getprs_tl(sv_ps,sv_tv,sv_prse) ! Convert input normalized RH to q if(do_normal_rh_to_q) call normal_rh_to_q(workrh,sv_tv,sv_prse,sv_q) ! Calculate sensible temperature if(do_tv_to_tsen) call tv_to_tsen(sv_tv,sv_q,sv_tsen) ! Convert streamfunction and velocity potential to u,v if(do_getuv) call getuv(sv_u,sv_v,workst,workvp,0) if (nclouds>0) then if (cw_to_hydro) then ! Case when cw is generated from hydrometeors call cw2hydro_tl(sval(jj),wbundle,clouds,nclouds) else ! Case when cloud-vars map one-to-one, take care of them together do ic=1,nclouds id=getindex(cvars3d,clouds(ic)) if (id>0) then call gsi_bundlegetpointer (sval(jj),clouds(ic),sv_rank3,istatus) call gsi_bundlegetvar (wbundle, clouds(ic),sv_rank3,istatus) endif enddo end if end if ! destroy temporary bundle call gsi_bundledestroy(wbundle,istatus) if(istatus/=0) then write(6,*) trim(myname), ': trouble destroying work bundle' call stop2(999) endif end do ! Bias correction terms do ii=1,nsclen bval%predr(ii)=xhat%predr(ii)*sqrt(varprd(ii)) enddo do ii=1,npclen bval%predp(ii)=xhat%predp(ii)*sqrt(varprd(nsclen+ii)) enddo if (ntclen>0) then do ii=1,ntclen bval%predt(ii)=xhat%predt(ii)*sqrt(varprd(nsclen+npclen+ii)) enddo end if ! Clean up if (ngases>0) then deallocate(gases) endif if (nclouds>0) deallocate(clouds) return end subroutine control2model