diff --git a/src/gsi/m_pm10Node.F90 b/src/gsi/m_pm10Node.F90 deleted file mode 100644 index 3536403fa9..0000000000 --- a/src/gsi/m_pm10Node.F90 +++ /dev/null @@ -1,252 +0,0 @@ -module m_pm10Node -!$$$ subprogram documentation block -! . . . . -! subprogram: module m_pm10Node -! prgmmr: j guo -! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 -! date: 2016-05-18 -! -! abstract: class-module of obs-type pm10Node (in-situ pm-10) -! -! program history log: -! 2016-05-18 j guo - added this document block for the initial polymorphic -! implementation. -! -! input argument list: see Fortran 90 style document below -! -! output argument list: see Fortran 90 style document below -! -! attributes: -! language: Fortran 90 and/or above -! machine: -! -!$$$ end subprogram documentation block - -! module interface: - use m_obsdiagNode, only: obs_diag - use m_obsdiagNode, only: obs_diags - use kinds , only: i_kind,r_kind - use mpeu_util, only: assert_,die,perr,warn,tell - use m_obsNode, only: obsNode - implicit none - private - - public:: pm10Node - - type,extends(obsNode):: pm10Node -! to avoid separate coding for pm10 profile e.g. from aircraft or -! soundings obs weights are coded as -! wij(8) even though for surface pm10 wij(4) would be sufficient. -! also surface pm10 may be treated differently than now for vertical -! interpolation - - !type(pm10_ob_type),pointer :: llpoint => NULL() - type(obs_diag), pointer :: diags => NULL() - real(r_kind) :: res ! pm10 residual - real(r_kind) :: err2 ! pm10 obs error squared - real(r_kind) :: raterr2 ! square of ratio of final obs error -! to original obs error - !real(r_kind) :: time ! observation time - real(r_kind) :: b ! variational quality control parameter - real(r_kind) :: pg ! variational quality control parameter - real(r_kind) :: wij(8) ! horizontal interpolation weights - integer(i_kind) :: ij(8) ! horizontal locations - !logical :: luse ! flag indicating if ob is used in pen. - !integer(i_kind) :: idv,iob ! device id and obs index for sorting - !real (r_kind) :: elat, elon ! earth lat-lon for redistribution - !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution - real (r_kind) :: dlev ! reference to the vertical grid - contains - procedure,nopass:: mytype - procedure:: setHop => obsNode_setHop_ - procedure:: xread => obsNode_xread_ - procedure:: xwrite => obsNode_xwrite_ - procedure:: isvalid => obsNode_isvalid_ - procedure:: gettlddp => gettlddp_ - - ! procedure, nopass:: headerRead => obsHeader_read_ - ! procedure, nopass:: headerWrite => obsHeader_write_ - ! procedure:: init => obsNode_init_ - ! procedure:: clean => obsNode_clean_ - end type pm10Node - - public:: pm10Node_typecast - public:: pm10Node_nextcast - interface pm10Node_typecast; module procedure typecast_ ; end interface - interface pm10Node_nextcast; module procedure nextcast_ ; end interface - - public:: pm10Node_appendto - interface pm10Node_appendto; module procedure appendto_ ; end interface - - character(len=*),parameter:: MYNAME="m_pm10Node" - -#include "myassert.H" -#include "mytrace.H" -contains -function typecast_(aNode) result(ptr_) -!-- cast a class(obsNode) to a type(pm10Node) - use m_obsNode, only: obsNode - implicit none - type(pm10Node),pointer:: ptr_ - class(obsNode),pointer,intent(in):: aNode - ptr_ => null() - if(.not.associated(aNode)) return - ! logically, typecast of a null-reference is a null pointer. - select type(aNode) - type is(pm10Node) - ptr_ => aNode - end select -return -end function typecast_ - -function nextcast_(aNode) result(ptr_) -!-- cast an obsNode_next(obsNode) to a type(pm10Node) - use m_obsNode, only: obsNode,obsNode_next - implicit none - type(pm10Node),pointer:: ptr_ - class(obsNode),target ,intent(in):: aNode - - class(obsNode),pointer:: inode_ - inode_ => obsNode_next(aNode) - ptr_ => typecast_(inode_) -return -end function nextcast_ - -subroutine appendto_(aNode,oll) -!-- append aNode to linked-list oLL - use m_obsNode , only: obsNode - use m_obsLList, only: obsLList,obsLList_appendNode - implicit none - type(pm10Node),pointer,intent(in):: aNode - type(obsLList),intent(inout):: oLL - - class(obsNode),pointer:: inode_ - inode_ => aNode - call obsLList_appendNode(oLL,inode_) - inode_ => null() -end subroutine appendto_ - -! obsNode implementations - -function mytype() - implicit none - character(len=:),allocatable:: mytype - mytype="[pm10Node]" -end function mytype - -subroutine obsNode_xread_(aNode,iunit,istat,diagLookup,skip) - use m_obsdiagNode, only: obsdiagLookup_locate - implicit none - class(pm10Node) , intent(inout):: aNode - integer(i_kind) , intent(in ):: iunit - integer(i_kind) , intent( out):: istat - type(obs_diags) , intent(in ):: diagLookup - logical,optional, intent(in ):: skip - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xread_' - logical:: skip_ - skip_=.false. - if(present(skip)) skip_=skip - - istat=0 - if(skip_) then - read(iunit,iostat=istat) - if(istat/=0) then - call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - endif - - else - read(iunit,iostat=istat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%dlev , & - aNode%wij , & - aNode%ij - if (istat/=0) then - call perr(myname_,'read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - end if - - aNode%diags => obsdiagLookup_locate(diagLookup,aNode%idv,aNode%iob,1_i_kind) - if(.not.associated(aNode%diags)) then - call perr(myname_,'obsdiagLookup_locate(), %idv =',aNode%idv) - call perr(myname_,' %iob =',aNode%iob) - call die(myname_) - endif - endif -_EXIT_(myname_) -return -end subroutine obsNode_xread_ - -subroutine obsNode_xwrite_(aNode,junit,jstat) - implicit none - class(pm10Node),intent(in):: aNode - integer(i_kind),intent(in ):: junit - integer(i_kind),intent( out):: jstat - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xwrite_' - - jstat=0 - write(junit,iostat=jstat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%dlev , & - aNode%wij , & - aNode%ij - if (jstat/=0) then - call perr(myname_,'write(%(res,err2,...)), iostat =',jstat) - _EXIT_(myname_) - return - end if -_EXIT_(myname_) -return -end subroutine obsNode_xwrite_ - -subroutine obsNode_setHop_(aNode) - use m_cvgridLookup, only: cvgridLookup_getiw - implicit none - class(pm10Node),intent(inout):: aNode - -! character(len=*),parameter:: myname_=MYNAME//'.obsNode_setHop_' -!_ENTRY_(myname_) - ! %dlev is defined, but not used. See setuppm10() for information. - call cvgridLookup_getiw(aNode%elat,aNode%elon,aNode%ij(1:4),aNode%wij(1:4)) - aNode% ij(5:8) = aNode%ij(1:4) - aNode%wij(5:8) = 0._r_kind -!_EXIT_(myname_) -return -end subroutine obsNode_setHop_ - -function obsNode_isvalid_(aNode) result(isvalid_) - implicit none - logical:: isvalid_ - class(pm10Node),intent(in):: aNode - -! character(len=*),parameter:: myname_=MYNAME//'::obsNode_isvalid_' -!_ENTRY_(myname_) - isvalid_=associated(aNode%diags) -!_EXIT_(myname_) -return -end function obsNode_isvalid_ - -pure subroutine gettlddp_(aNode,jiter,tlddp,nob) - use kinds, only: r_kind - implicit none - class(pm10Node), intent(in):: aNode - integer(kind=i_kind),intent(in):: jiter - real(kind=r_kind),intent(inout):: tlddp - integer(kind=i_kind),optional,intent(inout):: nob - - tlddp = tlddp + aNode%diags%tldepart(jiter)*aNode%diags%tldepart(jiter) - if(present(nob)) nob=nob+1 -return -end subroutine gettlddp_ - -end module m_pm10Node diff --git a/src/gsi/m_pm10node.F90 b/src/gsi/m_pm10node.F90 new file mode 100644 index 0000000000..5afd4ac5a2 --- /dev/null +++ b/src/gsi/m_pm10node.F90 @@ -0,0 +1,252 @@ +module m_pm10node +!$$$ subprogram documentation block +! . . . . +! subprogram: module m_pm10node +! prgmmr: j guo +! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 +! date: 2016-05-18 +! +! abstract: class-module of obs-type pm10node (in-situ pm-10) +! +! program history log: +! 2016-05-18 j guo - added this document block for the initial polymorphic +! implementation. +! +! input argument list: see Fortran 90 style document below +! +! output argument list: see Fortran 90 style document below +! +! attributes: +! language: Fortran 90 and/or above +! machine: +! +!$$$ end subprogram documentation block + +! module interface: + use m_obsdiagnode, only: obs_diag + use m_obsdiagnode, only: obs_diags + use kinds , only: i_kind,r_kind + use mpeu_util, only: assert_,die,perr,warn,tell + use m_obsnode, only: obsnode + implicit none + private + + public:: pm10node + + type,extends(obsnode):: pm10node +! to avoid separate coding for pm10 profile e.g. from aircraft or +! soundings obs weights are coded as +! wij(8) even though for surface pm10 wij(4) would be sufficient. +! also surface pm10 may be treated differently than now for vertical +! interpolation + + !type(pm10_ob_type),pointer :: llpoint => null() + type(obs_diag), pointer :: diags => null() + real(r_kind) :: res ! pm10 residual + real(r_kind) :: err2 ! pm10 obs error squared + real(r_kind) :: raterr2 ! square of ratio of final obs error +! to original obs error + !real(r_kind) :: time ! observation time + real(r_kind) :: b ! variational quality control parameter + real(r_kind) :: pg ! variational quality control parameter + real(r_kind) :: wij(8) ! horizontal interpolation weights + integer(i_kind) :: ij(8) ! horizontal locations + !logical :: luse ! flag indicating if ob is used in pen. + !integer(i_kind) :: idv,iob ! device id and obs index for sorting + !real (r_kind) :: elat, elon ! earth lat-lon for redistribution + !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution + real (r_kind) :: dlev ! reference to the vertical grid + contains + procedure,nopass:: mytype + procedure:: sethop => obsnode_sethop_ + procedure:: xread => obsnode_xread_ + procedure:: xwrite => obsnode_xwrite_ + procedure:: isvalid => obsnode_isvalid_ + procedure:: gettlddp => gettlddp_ + + ! procedure, nopass:: headerread => obsheader_read_ + ! procedure, nopass:: headerwrite => obsheader_write_ + ! procedure:: init => obsnode_init_ + ! procedure:: clean => obsnode_clean_ + end type pm10node + + public:: pm10node_typecast + public:: pm10node_nextcast + interface pm10node_typecast; module procedure typecast_ ; end interface + interface pm10node_nextcast; module procedure nextcast_ ; end interface + + public:: pm10node_appendto + interface pm10node_appendto; module procedure appendto_ ; end interface + + character(len=*),parameter:: myname="m_pm10node" + +#include "myassert.H" +#include "mytrace.H" +contains +function typecast_(anode) result(ptr_) +!-- cast a class(obsnode) to a type(pm10node) + use m_obsnode, only: obsnode + implicit none + type(pm10node),pointer:: ptr_ + class(obsnode),pointer,intent(in):: anode + ptr_ => null() + if(.not.associated(anode)) return + ! logically, typecast of a null-reference is a null pointer. + select type(anode) + type is(pm10node) + ptr_ => anode + end select + return +end function typecast_ + +function nextcast_(anode) result(ptr_) +!-- cast an obsnode_next(obsnode) to a type(pm10node) + use m_obsnode, only: obsnode,obsnode_next + implicit none + type(pm10node),pointer:: ptr_ + class(obsnode),target ,intent(in):: anode + + class(obsnode),pointer:: inode_ + inode_ => obsnode_next(anode) + ptr_ => typecast_(inode_) + return +end function nextcast_ + +subroutine appendto_(anode,oll) +!-- append anode to linked-list oll + use m_obsnode , only: obsnode + use m_obsllist, only: obsllist,obsllist_appendnode + implicit none + type(pm10node),pointer,intent(in):: anode + type(obsllist),intent(inout):: oll + + class(obsnode),pointer:: inode_ + inode_ => anode + call obsllist_appendnode(oll,inode_) + inode_ => null() +end subroutine appendto_ + +! obsnode implementations + +function mytype() + implicit none + character(len=:),allocatable:: mytype + mytype="[pm10node]" +end function mytype + +subroutine obsnode_xread_(anode,iunit,istat,diaglookup,skip) + use m_obsdiagnode, only: obsdiaglookup_locate + implicit none + class(pm10node) , intent(inout):: anode + integer(i_kind) , intent(in ):: iunit + integer(i_kind) , intent( out):: istat + type(obs_diags) , intent(in ):: diaglookup + logical,optional, intent(in ):: skip + + character(len=*),parameter:: myname_=myname//'.obsnode_xread_' + logical:: skip_ + skip_=.false. + if(present(skip)) skip_=skip + + istat=0 + if(skip_) then + read(iunit,iostat=istat) + if(istat/=0) then + call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + endif + + else + read(iunit,iostat=istat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%dlev , & + anode%wij , & + anode%ij + if (istat/=0) then + call perr(myname_,'read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + end if + + anode%diags => obsdiaglookup_locate(diaglookup,anode%idv,anode%iob,1) + if(.not.associated(anode%diags)) then + call perr(myname_,'obsdiaglookup_locate(), %idv =',anode%idv) + call perr(myname_,' %iob =',anode%iob) + call die(myname_) + endif + endif +_EXIT_(myname_) + return +end subroutine obsnode_xread_ + +subroutine obsnode_xwrite_(anode,junit,jstat) + implicit none + class(pm10node),intent(in):: anode + integer(i_kind),intent(in ):: junit + integer(i_kind),intent( out):: jstat + + character(len=*),parameter:: myname_=myname//'.obsnode_xwrite_' + + jstat=0 + write(junit,iostat=jstat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%dlev , & + anode%wij , & + anode%ij + if (jstat/=0) then + call perr(myname_,'write(%(res,err2,...)), iostat =',jstat) + _EXIT_(myname_) + return + end if +_EXIT_(myname_) + return +end subroutine obsnode_xwrite_ + +subroutine obsnode_sethop_(anode) + use m_cvgridlookup, only: cvgridlookup_getiw + implicit none + class(pm10node),intent(inout):: anode + +! character(len=*),parameter:: myname_=myname//'.obsnode_sethop_' +!_ENTRY_(myname_) + ! %dlev is defined, but not used. See setuppm10() for information. + call cvgridlookup_getiw(anode%elat,anode%elon,anode%ij(1:4),anode%wij(1:4)) + anode% ij(5:8) = anode%ij(1:4) + anode%wij(5:8) = 0._r_kind +!_EXIT_(myname_) + return +end subroutine obsnode_sethop_ + +function obsnode_isvalid_(anode) result(isvalid_) + implicit none + logical:: isvalid_ + class(pm10node),intent(in):: anode + +! character(len=*),parameter:: myname_=myname//'::obsnode_isvalid_' +!_ENTRY_(myname_) + isvalid_=associated(anode%diags) +!_EXIT_(myname_) + return +end function obsnode_isvalid_ + +pure subroutine gettlddp_(anode,jiter,tlddp,nob) + use kinds, only: r_kind + implicit none + class(pm10node), intent(in):: anode + integer(kind=i_kind),intent(in):: jiter + real(kind=r_kind),intent(inout):: tlddp + integer(kind=i_kind),optional,intent(inout):: nob + + tlddp = tlddp + anode%diags%tldepart(jiter)*anode%diags%tldepart(jiter) + if(present(nob)) nob=nob+1 + return +end subroutine gettlddp_ + +end module m_pm10node diff --git a/src/gsi/m_pm2_5Node.F90 b/src/gsi/m_pm2_5Node.F90 deleted file mode 100644 index 4a865e530d..0000000000 --- a/src/gsi/m_pm2_5Node.F90 +++ /dev/null @@ -1,254 +0,0 @@ -module m_pm2_5Node -!$$$ subprogram documentation block -! . . . . -! subprogram: module m_pm2_5Node -! prgmmr: j guo -! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 -! date: 2016-05-18 -! -! abstract: class-module of obs-type pm2_5Node (in-situ pm-2.5) -! -! program history log: -! 2016-05-18 j guo - added this document block for the initial polymorphic -! implementation. -! -! input argument list: see Fortran 90 style document below -! -! output argument list: see Fortran 90 style document below -! -! attributes: -! language: Fortran 90 and/or above -! machine: -! -!$$$ end subprogram documentation block - -! module interface: - use m_obsdiagNode, only: obs_diag - use m_obsdiagNode, only: obs_diags - use kinds , only: i_kind,r_kind - use mpeu_util, only: assert_,die,perr,warn,tell - use m_obsNode, only: obsNode - implicit none - private - - public:: pm2_5Node - - type,extends(obsNode):: pm2_5Node -! to avoid separate coding for pm2_5 profile e.g. from aircraft or -! soundings obs weights are coded as -! wij(8) even though for surface pm2_5 wij(4) would be sufficient. -! also surface pm2_5 may be treated differently than now for vertical -! interpolation - - !type(pm2_5_ob_type),pointer :: llpoint => NULL() - type(obs_diag), pointer :: diags => NULL() - real(r_kind) :: res ! pm2_5 residual - real(r_kind) :: err2 ! pm2_5 obs error squared - real(r_kind) :: raterr2 ! square of ratio of final obs error -! to original obs error - !real(r_kind) :: time ! observation time - real(r_kind) :: b ! variational quality control parameter - real(r_kind) :: pg ! variational quality control parameter - real(r_kind) :: wij(8) ! horizontal interpolation weights - integer(i_kind) :: ij(8) ! horizontal locations - !logical :: luse ! flag indicating if ob is used in pen. - !integer(i_kind) :: idv,iob ! device id and obs index for sorting - !real (r_kind) :: elat, elon ! earth lat-lon for redistribution - !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution - real (r_kind) :: dlev ! reference to the vertical grid - contains - procedure,nopass:: mytype - procedure:: setHop => obsNode_setHop_ - procedure:: xread => obsNode_xread_ - procedure:: xwrite => obsNode_xwrite_ - procedure:: isvalid => obsNode_isvalid_ - procedure:: gettlddp => gettlddp_ - - ! procedure, nopass:: headerRead => obsHeader_read_ - ! procedure, nopass:: headerWrite => obsHeader_write_ - ! procedure:: init => obsNode_init_ - ! procedure:: clean => obsNode_clean_ - end type pm2_5Node - - public:: pm2_5Node_typecast - public:: pm2_5Node_nextcast - interface pm2_5Node_typecast; module procedure typecast_ ; end interface - interface pm2_5Node_nextcast; module procedure nextcast_ ; end interface - - public:: pm2_5Node_appendto - interface pm2_5Node_appendto; module procedure appendto_ ; end interface - - character(len=*),parameter:: MYNAME="m_pm2_5Node" - -#include "myassert.H" -#include "mytrace.H" -contains -function typecast_(aNode) result(ptr_) -!-- cast a class(obsNode) to a type(pm2_5Node) - use m_obsNode, only: obsNode - implicit none - type(pm2_5Node),pointer:: ptr_ - class(obsNode ),pointer,intent(in):: aNode - ptr_ => null() - if(.not.associated(aNode)) return - ! logically, typecast of a null-reference is a null pointer. - select type(aNode) - type is(pm2_5Node) - ptr_ => aNode - end select -return -end function typecast_ - -function nextcast_(aNode) result(ptr_) -!-- cast an obsNode_next(obsNode) to a type(pm2_5Node) - use m_obsNode, only: obsNode,obsNode_next - implicit none - type(pm2_5Node),pointer:: ptr_ - class(obsNode ),target ,intent(in):: aNode - - class(obsNode),pointer:: inode_ - inode_ => obsNode_next(aNode) - ptr_ => typecast_(inode_) -return -end function nextcast_ - -subroutine appendto_(aNode,oll) -!-- append aNode to linked-list oLL - use m_obsNode , only: obsNode - use m_obsLList, only: obsLList,obsLList_appendNode - implicit none - type(pm2_5Node),pointer,intent(in):: aNode - type(obsLList),intent(inout):: oLL - - class(obsNode),pointer:: inode_ - inode_ => aNode - call obsLList_appendNode(oLL,inode_) - inode_ => null() -end subroutine appendto_ - -! obsNode implementations - -function mytype() - implicit none - character(len=:),allocatable:: mytype - mytype="[pm2_5Node]" -end function mytype - -subroutine obsNode_xread_(aNode,iunit,istat,diagLookup,skip) - use m_obsdiagNode, only: obsdiagLookup_locate - implicit none - class(pm2_5Node) , intent(inout):: aNode - integer(i_kind) , intent(in ):: iunit - integer(i_kind) , intent( out):: istat - type(obs_diags) , intent(in ):: diagLookup - logical,optional, intent(in ):: skip - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xread_' - logical:: skip_ -_ENTRY_(myname_) - skip_=.false. - if(present(skip)) skip_=skip - - istat=0 - if(skip_) then - read(iunit,iostat=istat) - if(istat/=0) then - call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - endif - - else - read(iunit,iostat=istat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%dlev , & - aNode%wij , & - aNode%ij - if (istat/=0) then - call perr(myname_,'read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - end if - - aNode%diags => obsdiagLookup_locate(diagLookup,aNode%idv,aNode%iob,1_i_kind) - if(.not.associated(aNode%diags)) then - call perr(myname_,'obsdiagLookup_locate(), %idv =',aNode%idv) - call perr(myname_,' %iob =',aNode%iob) - call die(myname_) - endif - endif -_EXIT_(myname_) -return -end subroutine obsNode_xread_ - -subroutine obsNode_xwrite_(aNode,junit,jstat) - implicit none - class(pm2_5Node),intent(in):: aNode - integer(i_kind),intent(in ):: junit - integer(i_kind),intent( out):: jstat - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xwrite_' -_ENTRY_(myname_) - - jstat=0 - write(junit,iostat=jstat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%dlev , & - aNode%wij , & - aNode%ij - if (jstat/=0) then - call perr(myname_,'write(%(res,err2,...)), iostat =',jstat) - _EXIT_(myname_) - return - end if -_EXIT_(myname_) -return -end subroutine obsNode_xwrite_ - -subroutine obsNode_setHop_(aNode) - use m_cvgridLookup, only: cvgridLookup_getiw - implicit none - class(pm2_5Node),intent(inout):: aNode - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_setHop_' -_ENTRY_(myname_) - ! %dlev is defined, but not used. See setuppm2_5() for information. - call cvgridLookup_getiw(aNode%elat,aNode%elon,aNode%ij(1:4),aNode%wij(1:4)) - aNode% ij(5:8) = aNode%ij(1:4) - aNode%wij(5:8) = 0. -_EXIT_(myname_) -return -end subroutine obsNode_setHop_ - -function obsNode_isvalid_(aNode) result(isvalid_) - implicit none - logical:: isvalid_ - class(pm2_5Node),intent(in):: aNode - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_isvalid_' -_ENTRY_(myname_) - isvalid_=associated(aNode%diags) -_EXIT_(myname_) -return -end function obsNode_isvalid_ - -pure subroutine gettlddp_(aNode,jiter,tlddp,nob) - use kinds, only: r_kind - implicit none - class(pm2_5Node), intent(in):: aNode - integer(kind=i_kind),intent(in):: jiter - real(kind=r_kind),intent(inout):: tlddp - integer(kind=i_kind),optional,intent(inout):: nob - - tlddp = tlddp + aNode%diags%tldepart(jiter)*aNode%diags%tldepart(jiter) - if(present(nob)) nob=nob+1 -return -end subroutine gettlddp_ - -end module m_pm2_5Node diff --git a/src/gsi/m_pm2_5node.F90 b/src/gsi/m_pm2_5node.F90 new file mode 100644 index 0000000000..02be31079d --- /dev/null +++ b/src/gsi/m_pm2_5node.F90 @@ -0,0 +1,254 @@ +module m_pm2_5node +!$$$ subprogram documentation block +! . . . . +! subprogram: module m_pm2_5node +! prgmmr: j guo +! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 +! date: 2016-05-18 +! +! abstract: class-module of obs-type pm2_5node (in-situ pm-2.5) +! +! program history log: +! 2016-05-18 j guo - added this document block for the initial polymorphic +! implementation. +! +! input argument list: see Fortran 90 style document below +! +! output argument list: see Fortran 90 style document below +! +! attributes: +! language: Fortran 90 and/or above +! machine: +! +!$$$ end subprogram documentation block + +! module interface: + use m_obsdiagnode, only: obs_diag + use m_obsdiagnode, only: obs_diags + use kinds , only: i_kind,r_kind + use mpeu_util, only: assert_,die,perr,warn,tell + use m_obsnode, only: obsnode + implicit none + private + + public:: pm2_5node + + type,extends(obsnode):: pm2_5node +! to avoid separate coding for pm2_5 profile e.g. from aircraft or +! soundings obs weights are coded as +! wij(8) even though for surface pm2_5 wij(4) would be sufficient. +! also surface pm2_5 may be treated differently than now for vertical +! interpolation + + !type(pm2_5_ob_type),pointer :: llpoint => null() + type(obs_diag), pointer :: diags => null() + real(r_kind) :: res ! pm2_5 residual + real(r_kind) :: err2 ! pm2_5 obs error squared + real(r_kind) :: raterr2 ! square of ratio of final obs error +! to original obs error + !real(r_kind) :: time ! observation time + real(r_kind) :: b ! variational quality control parameter + real(r_kind) :: pg ! variational quality control parameter + real(r_kind) :: wij(8) ! horizontal interpolation weights + integer(i_kind) :: ij(8) ! horizontal locations + !logical :: luse ! flag indicating if ob is used in pen. + !integer(i_kind) :: idv,iob ! device id and obs index for sorting + !real (r_kind) :: elat, elon ! earth lat-lon for redistribution + !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution + real (r_kind) :: dlev ! reference to the vertical grid + contains + procedure,nopass:: mytype + procedure:: sethop => obsnode_sethop_ + procedure:: xread => obsnode_xread_ + procedure:: xwrite => obsnode_xwrite_ + procedure:: isvalid => obsnode_isvalid_ + procedure:: gettlddp => gettlddp_ + + ! procedure, nopass:: headerread => obsheader_read_ + ! procedure, nopass:: headerwrite => obsheader_write_ + ! procedure:: init => obsnode_init_ + ! procedure:: clean => obsnode_clean_ + end type pm2_5node + + public:: pm2_5node_typecast + public:: pm2_5node_nextcast + interface pm2_5node_typecast; module procedure typecast_ ; end interface + interface pm2_5node_nextcast; module procedure nextcast_ ; end interface + + public:: pm2_5node_appendto + interface pm2_5node_appendto; module procedure appendto_ ; end interface + + character(len=*),parameter:: myname="m_pm2_5node" + +#include "myassert.H" +#include "mytrace.H" +contains +function typecast_(anode) result(ptr_) +!-- cast a class(obsnode) to a type(pm2_5node) + use m_obsnode, only: obsnode + implicit none + type(pm2_5node),pointer:: ptr_ + class(obsnode ),pointer,intent(in):: anode + ptr_ => null() + if(.not.associated(anode)) return + ! logically, typecast of a null-reference is a null pointer. + select type(anode) + type is(pm2_5node) + ptr_ => anode + end select + return +end function typecast_ + +function nextcast_(anode) result(ptr_) +!-- cast an obsnode_next(obsnode) to a type(pm2_5node) + use m_obsnode, only: obsnode,obsnode_next + implicit none + type(pm2_5node),pointer:: ptr_ + class(obsnode ),target ,intent(in):: anode + + class(obsnode),pointer:: inode_ + inode_ => obsnode_next(anode) + ptr_ => typecast_(inode_) + return +end function nextcast_ + +subroutine appendto_(anode,oll) +!-- append anode to linked-list oll + use m_obsnode , only: obsnode + use m_obsllist, only: obsllist,obsllist_appendnode + implicit none + type(pm2_5node),pointer,intent(in):: anode + type(obsllist),intent(inout):: oll + + class(obsnode),pointer:: inode_ + inode_ => anode + call obsllist_appendnode(oll,inode_) + inode_ => null() +end subroutine appendto_ + +! obsnode implementations + +function mytype() + implicit none + character(len=:),allocatable:: mytype + mytype="[pm2_5node]" +end function mytype + +subroutine obsnode_xread_(anode,iunit,istat,diaglookup,skip) + use m_obsdiagnode, only: obsdiaglookup_locate + implicit none + class(pm2_5node) , intent(inout):: anode + integer(i_kind) , intent(in ):: iunit + integer(i_kind) , intent( out):: istat + type(obs_diags) , intent(in ):: diaglookup + logical,optional, intent(in ):: skip + + character(len=*),parameter:: myname_=myname//'.obsnode_xread_' + logical:: skip_ +_ENTRY_(myname_) + skip_=.false. + if(present(skip)) skip_=skip + + istat=0 + if(skip_) then + read(iunit,iostat=istat) + if(istat/=0) then + call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + endif + + else + read(iunit,iostat=istat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%dlev , & + anode%wij , & + anode%ij + if (istat/=0) then + call perr(myname_,'read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + end if + + anode%diags => obsdiaglookup_locate(diaglookup,anode%idv,anode%iob,1) + if(.not.associated(anode%diags)) then + call perr(myname_,'obsdiaglookup_locate(), %idv =',anode%idv) + call perr(myname_,' %iob =',anode%iob) + call die(myname_) + endif + endif +_EXIT_(myname_) + return +end subroutine obsnode_xread_ + +subroutine obsnode_xwrite_(anode,junit,jstat) + implicit none + class(pm2_5node),intent(in):: anode + integer(i_kind),intent(in ):: junit + integer(i_kind),intent( out):: jstat + + character(len=*),parameter:: myname_=myname//'.obsnode_xwrite_' +_ENTRY_(myname_) + + jstat=0 + write(junit,iostat=jstat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%dlev , & + anode%wij , & + anode%ij + if (jstat/=0) then + call perr(myname_,'write(%(res,err2,...)), iostat =',jstat) + _EXIT_(myname_) + return + end if +_EXIT_(myname_) + return +end subroutine obsnode_xwrite_ + +subroutine obsnode_sethop_(anode) + use m_cvgridlookup, only: cvgridlookup_getiw + implicit none + class(pm2_5node),intent(inout):: anode + + character(len=*),parameter:: myname_=myname//'.obsnode_sethop_' +_ENTRY_(myname_) + ! %dlev is defined, but not used. See setuppm2_5() for information. + call cvgridlookup_getiw(anode%elat,anode%elon,anode%ij(1:4),anode%wij(1:4)) + anode% ij(5:8) = anode%ij(1:4) + anode%wij(5:8) = 0. +_EXIT_(myname_) + return +end subroutine obsnode_sethop_ + +function obsnode_isvalid_(anode) result(isvalid_) + implicit none + logical:: isvalid_ + class(pm2_5node),intent(in):: anode + + character(len=*),parameter:: myname_=myname//'::obsnode_isvalid_' +_ENTRY_(myname_) + isvalid_=associated(anode%diags) +_EXIT_(myname_) + return +end function obsnode_isvalid_ + +pure subroutine gettlddp_(anode,jiter,tlddp,nob) + use kinds, only: r_kind + implicit none + class(pm2_5node), intent(in):: anode + integer(kind=i_kind),intent(in):: jiter + real(kind=r_kind),intent(inout):: tlddp + integer(kind=i_kind),optional,intent(inout):: nob + + tlddp = tlddp + anode%diags%tldepart(jiter)*anode%diags%tldepart(jiter) + if(present(nob)) nob=nob+1 + return +end subroutine gettlddp_ + +end module m_pm2_5node diff --git a/src/gsi/m_pmslNode.F90 b/src/gsi/m_pmslNode.F90 deleted file mode 100644 index 2c41d0cf2a..0000000000 --- a/src/gsi/m_pmslNode.F90 +++ /dev/null @@ -1,244 +0,0 @@ -module m_pmslNode -!$$$ subprogram documentation block -! . . . . -! subprogram: module m_pmslNode -! prgmmr: j guo -! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 -! date: 2016-05-18 -! -! abstract: class-module of obs-type pmslNode (pmsl) -! -! program history log: -! 2016-05-18 j guo - added this document block for the initial polymorphic -! implementation. -! 2017-01-19 j guo - fixed a null reference bug in nextcast_(). -! -! input argument list: see Fortran 90 style document below -! -! output argument list: see Fortran 90 style document below -! -! attributes: -! language: Fortran 90 and/or above -! machine: -! -!$$$ end subprogram documentation block - -! module interface: - use m_obsdiagNode, only: obs_diag - use m_obsdiagNode, only: obs_diags - use kinds , only: i_kind,r_kind - use mpeu_util, only: assert_,die,perr,warn,tell - use m_obsNode, only: obsNode - implicit none - private - - public:: pmslNode - - type,extends(obsNode):: pmslNode - ! private - !type(pmsl_ob_type),pointer :: llpoint => NULL() - type(obs_diag), pointer :: diags => NULL() - real(r_kind) :: res ! pmsl residual - real(r_kind) :: err2 ! pmsl error squared - real(r_kind) :: raterr2 ! square of ratio of final obs error - ! to original obs error - !real(r_kind) :: time ! observation time in sec - real(r_kind) :: b ! variational quality control parameter - real(r_kind) :: pg ! variational quality control parameter - real(r_kind) :: wij(4) ! horizontal interpolation weights - integer(i_kind) :: ij(4) ! horizontal locations - !logical :: luse ! flag indicating if ob is used in pen. - - !integer(i_kind) :: idv,iob ! device id and obs index for sorting - !real (r_kind) :: elat, elon ! earth lat-lon for redistribution - !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution - contains - procedure,nopass:: mytype - procedure:: setHop => obsNode_setHop_ - procedure:: xread => obsNode_xread_ - procedure:: xwrite => obsNode_xwrite_ - procedure:: isvalid => obsNode_isvalid_ - procedure:: gettlddp => gettlddp_ - - ! procedure, nopass:: headerRead => obsHeader_read_ - ! procedure, nopass:: headerWrite => obsHeader_write_ - ! procedure:: init => obsNode_init_ - ! procedure:: clean => obsNode_clean_ - end type pmslNode - - public:: pmslNode_typecast - public:: pmslNode_nextcast - interface pmslNode_typecast; module procedure typecast_ ; end interface - interface pmslNode_nextcast; module procedure nextcast_ ; end interface - - public:: pmslNode_appendto - interface pmslNode_appendto; module procedure appendto_ ; end interface - - character(len=*),parameter:: MYNAME="m_pmslNode" - -#include "myassert.H" -#include "mytrace.H" -contains -function typecast_(aNode) result(ptr_) -!-- cast a class(obsNode) to a type(pmslNode) - use m_obsNode, only: obsNode - implicit none - type(pmslNode),pointer:: ptr_ - class(obsNode),pointer,intent(in):: aNode - ptr_ => null() - if(.not.associated(aNode)) return - ! logically, typecast of a null-reference is a null pointer. - select type(aNode) - type is(pmslNode) - ptr_ => aNode - end select -return -end function typecast_ - -function nextcast_(aNode) result(ptr_) -!-- cast an obsNode_next(obsNode) to a type(pmslNode) - use m_obsNode, only: obsNode,obsNode_next - implicit none - type(pmslNode),pointer:: ptr_ - class(obsNode),target ,intent(in):: aNode - - class(obsNode),pointer:: inode_ - inode_ => obsNode_next(aNode) - ptr_ => typecast_(inode_) -return -end function nextcast_ - -subroutine appendto_(aNode,oll) -!-- append aNode to linked-list oLL - use m_obsNode , only: obsNode - use m_obsLList, only: obsLList,obsLList_appendNode - implicit none - type(pmslNode),pointer,intent(in):: aNode - type(obsLList),intent(inout):: oLL - - class(obsNode),pointer:: inode_ - inode_ => aNode - call obsLList_appendNode(oLL,inode_) - inode_ => null() -end subroutine appendto_ - -! obsNode implementations - -function mytype() - implicit none - character(len=:),allocatable:: mytype - mytype="[pmslNode]" -end function mytype - -subroutine obsNode_xread_(aNode,iunit,istat,diagLookup,skip) - use m_obsdiagNode, only: obsdiagLookup_locate - implicit none - class(pmslNode), intent(inout):: aNode - integer(i_kind), intent(in ):: iunit - integer(i_kind), intent( out):: istat - type(obs_diags), intent(in ):: diagLookup - logical,optional,intent(in ):: skip - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_xread_' - logical:: skip_ -_ENTRY_(myname_) - - skip_=.false. - if(present(skip)) skip_=skip - - if(skip_) then - read(iunit,iostat=istat) - if(istat/=0) then - call perr(myname_,'skipping read(%(res,...)), iostat =',istat) - _EXIT_(myname_) - return - endif - - else - read(iunit,iostat=istat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%wij , & - aNode%ij - if (istat/=0) then - call perr(myname_,'read(%(res,...)), iostat =',istat) - _EXIT_(myname_) - return - end if - - aNode%diags => obsdiagLookup_locate(diagLookup,aNode%idv,aNode%iob,1_i_kind) - if(.not. associated(aNode%diags)) then - call perr(myname_,'obsdiagLookup_locate(), %idv =',aNode%idv) - call perr(myname_,' %iob =',aNode%iob) - call die(myname_) - endif - endif ! if(skip_); else -_EXIT_(myname_) -return -end subroutine obsNode_xread_ - -subroutine obsNode_xwrite_(aNode,junit,jstat) - implicit none - class(pmslNode),intent(in):: aNode - integer(i_kind),intent(in ):: junit - integer(i_kind),intent( out):: jstat - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_xwrite_' -_ENTRY_(myname_) - - write(junit,iostat=jstat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%wij , & - aNode%ij - if (jstat/=0) then - call perr(myname_,'write(%(res,...)), iostat =',jstat) - _EXIT_(myname_) - return - end if -_EXIT_(myname_) -return -end subroutine obsNode_xwrite_ - -subroutine obsNode_setHop_(aNode) - use m_cvgridLookup, only: cvgridLookup_getiw - implicit none - class(pmslNode),intent(inout):: aNode - -! character(len=*),parameter:: myname_=MYNAME//'::obsNode_setHop_' -!_ENTRY_(myname_) - call cvgridLookup_getiw(aNode%elat,aNode%elon,aNode%ij,aNode%wij) -!_EXIT_(myname_) -return -end subroutine obsNode_setHop_ - -function obsNode_isvalid_(aNode) result(isvalid_) - implicit none - logical:: isvalid_ - class(pmslNode),intent(in):: aNode - -! character(len=*),parameter:: myname_=MYNAME//'::obsNode_isvalid_' -!_ENTRY_(myname_) - isvalid_= associated(aNode%diags) -!_EXIT_(myname_) -return -end function obsNode_isvalid_ - -pure subroutine gettlddp_(aNode,jiter,tlddp,nob) - use kinds, only: r_kind - implicit none - class(pmslNode), intent(in):: aNode - integer(kind=i_kind),intent(in):: jiter - real(kind=r_kind),intent(inout):: tlddp - integer(kind=i_kind),optional,intent(inout):: nob - - tlddp = tlddp + aNode%diags%tldepart(jiter)*aNode%diags%tldepart(jiter) - if(present(nob)) nob=nob+1 -return -end subroutine gettlddp_ - -end module m_pmslNode diff --git a/src/gsi/m_pmslnode.F90 b/src/gsi/m_pmslnode.F90 new file mode 100644 index 0000000000..ca7f6e5db0 --- /dev/null +++ b/src/gsi/m_pmslnode.F90 @@ -0,0 +1,244 @@ +module m_pmslnode +!$$$ subprogram documentation block +! . . . . +! subprogram: module m_pmslnode +! prgmmr: j guo +! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 +! date: 2016-05-18 +! +! abstract: class-module of obs-type pmslnode (pmsl) +! +! program history log: +! 2016-05-18 j guo - added this document block for the initial polymorphic +! implementation. +! 2017-01-19 j guo - fixed a null reference bug in nextcast_(). +! +! input argument list: see Fortran 90 style document below +! +! output argument list: see Fortran 90 style document below +! +! attributes: +! language: Fortran 90 and/or above +! machine: +! +!$$$ end subprogram documentation block + +! module interface: + use m_obsdiagnode, only: obs_diag + use m_obsdiagnode, only: obs_diags + use kinds , only: i_kind,r_kind + use mpeu_util, only: assert_,die,perr,warn,tell + use m_obsnode, only: obsnode + implicit none + private + + public:: pmslnode + + type,extends(obsnode):: pmslnode + ! private + !type(pmsl_ob_type),pointer :: llpoint => null() + type(obs_diag), pointer :: diags => null() + real(r_kind) :: res ! pmsl residual + real(r_kind) :: err2 ! pmsl error squared + real(r_kind) :: raterr2 ! square of ratio of final obs error + ! to original obs error + !real(r_kind) :: time ! observation time in sec + real(r_kind) :: b ! variational quality control parameter + real(r_kind) :: pg ! variational quality control parameter + real(r_kind) :: wij(4) ! horizontal interpolation weights + integer(i_kind) :: ij(4) ! horizontal locations + !logical :: luse ! flag indicating if ob is used in pen. + + !integer(i_kind) :: idv,iob ! device id and obs index for sorting + !real (r_kind) :: elat, elon ! earth lat-lon for redistribution + !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution + contains + procedure,nopass:: mytype + procedure:: sethop => obsnode_sethop_ + procedure:: xread => obsnode_xread_ + procedure:: xwrite => obsnode_xwrite_ + procedure:: isvalid => obsnode_isvalid_ + procedure:: gettlddp => gettlddp_ + + ! procedure, nopass:: headerread => obsheader_read_ + ! procedure, nopass:: headerwrite => obsheader_write_ + ! procedure:: init => obsnode_init_ + ! procedure:: clean => obsnode_clean_ + end type pmslnode + + public:: pmslnode_typecast + public:: pmslnode_nextcast + interface pmslnode_typecast; module procedure typecast_ ; end interface + interface pmslnode_nextcast; module procedure nextcast_ ; end interface + + public:: pmslnode_appendto + interface pmslnode_appendto; module procedure appendto_ ; end interface + + character(len=*),parameter:: myname="m_pmslnode" + +#include "myassert.H" +#include "mytrace.H" +contains +function typecast_(anode) result(ptr_) +!-- cast a class(obsnode) to a type(pmslnode) + use m_obsnode, only: obsnode + implicit none + type(pmslnode),pointer:: ptr_ + class(obsnode),pointer,intent(in):: anode + ptr_ => null() + if(.not.associated(anode)) return + ! logically, typecast of a null-reference is a null pointer. + select type(anode) + type is(pmslnode) + ptr_ => anode + end select + return +end function typecast_ + +function nextcast_(anode) result(ptr_) +!-- cast an obsnode_next(obsnode) to a type(pmslnode) + use m_obsnode, only: obsnode,obsnode_next + implicit none + type(pmslnode),pointer:: ptr_ + class(obsnode),target ,intent(in):: anode + + class(obsnode),pointer:: inode_ + inode_ => obsnode_next(anode) + ptr_ => typecast_(inode_) + return +end function nextcast_ + +subroutine appendto_(anode,oll) +!-- append anode to linked-list oll + use m_obsnode , only: obsnode + use m_obsllist, only: obsllist,obsllist_appendnode + implicit none + type(pmslnode),pointer,intent(in):: anode + type(obsllist),intent(inout):: oll + + class(obsnode),pointer:: inode_ + inode_ => anode + call obsllist_appendnode(oll,inode_) + inode_ => null() +end subroutine appendto_ + +! obsnode implementations + +function mytype() + implicit none + character(len=:),allocatable:: mytype + mytype="[pmslnode]" +end function mytype + +subroutine obsnode_xread_(anode,iunit,istat,diaglookup,skip) + use m_obsdiagnode, only: obsdiaglookup_locate + implicit none + class(pmslnode), intent(inout):: anode + integer(i_kind), intent(in ):: iunit + integer(i_kind), intent( out):: istat + type(obs_diags), intent(in ):: diaglookup + logical,optional,intent(in ):: skip + + character(len=*),parameter:: myname_=myname//'::obsnode_xread_' + logical:: skip_ +_ENTRY_(myname_) + + skip_=.false. + if(present(skip)) skip_=skip + + if(skip_) then + read(iunit,iostat=istat) + if(istat/=0) then + call perr(myname_,'skipping read(%(res,...)), iostat =',istat) + _EXIT_(myname_) + return + endif + + else + read(iunit,iostat=istat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%wij , & + anode%ij + if (istat/=0) then + call perr(myname_,'read(%(res,...)), iostat =',istat) + _EXIT_(myname_) + return + end if + + anode%diags => obsdiaglookup_locate(diaglookup,anode%idv,anode%iob,1) + if(.not. associated(anode%diags)) then + call perr(myname_,'obsdiaglookup_locate(), %idv =',anode%idv) + call perr(myname_,' %iob =',anode%iob) + call die(myname_) + endif + endif ! if(skip_); else +_EXIT_(myname_) + return +end subroutine obsnode_xread_ + +subroutine obsnode_xwrite_(anode,junit,jstat) + implicit none + class(pmslnode),intent(in):: anode + integer(i_kind),intent(in ):: junit + integer(i_kind),intent( out):: jstat + + character(len=*),parameter:: myname_=myname//'::obsnode_xwrite_' +_ENTRY_(myname_) + + write(junit,iostat=jstat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%wij , & + anode%ij + if (jstat/=0) then + call perr(myname_,'write(%(res,...)), iostat =',jstat) + _EXIT_(myname_) + return + end if +_EXIT_(myname_) + return +end subroutine obsnode_xwrite_ + +subroutine obsnode_sethop_(anode) + use m_cvgridlookup, only: cvgridlookup_getiw + implicit none + class(pmslnode),intent(inout):: anode + +! character(len=*),parameter:: myname_=myname//'::obsnode_sethop_' +!_ENTRY_(myname_) + call cvgridlookup_getiw(anode%elat,anode%elon,anode%ij,anode%wij) +!_EXIT_(myname_) + return +end subroutine obsnode_sethop_ + +function obsnode_isvalid_(anode) result(isvalid_) + implicit none + logical:: isvalid_ + class(pmslnode),intent(in):: anode + +! character(len=*),parameter:: myname_=myname//'::obsnode_isvalid_' +!_ENTRY_(myname_) + isvalid_= associated(anode%diags) +!_EXIT_(myname_) + return +end function obsnode_isvalid_ + +pure subroutine gettlddp_(anode,jiter,tlddp,nob) + use kinds, only: r_kind + implicit none + class(pmslnode), intent(in):: anode + integer(kind=i_kind),intent(in):: jiter + real(kind=r_kind),intent(inout):: tlddp + integer(kind=i_kind),optional,intent(inout):: nob + + tlddp = tlddp + anode%diags%tldepart(jiter)*anode%diags%tldepart(jiter) + if(present(nob)) nob=nob+1 + return +end subroutine gettlddp_ + +end module m_pmslnode diff --git a/src/gsi/m_psNode.F90 b/src/gsi/m_psNode.F90 deleted file mode 100644 index 988edee011..0000000000 --- a/src/gsi/m_psNode.F90 +++ /dev/null @@ -1,260 +0,0 @@ -module m_psNode -!$$$ subprogram documentation block -! . . . . -! subprogram: module m_psNode -! prgmmr: j guo -! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 -! date: 2016-05-18 -! -! abstract: class-module of obs-type psNode (surface pressure) -! -! program history log: -! 2016-05-18 j guo - added this document block for the initial polymorphic -! implementation. -! 2019-09-20 X.Su - add new variational QC parameters -! -! input argument list: see Fortran 90 style document below -! -! output argument list: see Fortran 90 style document below -! -! attributes: -! language: Fortran 90 and/or above -! machine: -! -!$$$ end subprogram documentation block - -! module interface: - use kinds , only: i_kind,r_kind - use mpeu_util, only: assert_,die,perr,warn,tell - use m_obsdiagNode, only: obs_diag - use m_obsNode, only: obsNode - implicit none - private - - public:: psNode - - type,extends(obsNode):: psNode - !private - type(obs_diag), pointer :: diags => NULL() - real(r_kind) :: res =0._r_kind ! surface pressure residual - real(r_kind) :: err2 =0._r_kind ! surface pressure error squared - ! in reciprocal - real(r_kind) :: raterr2=0._r_kind ! square of ratio of final obs error - ! to original obs error - real(r_kind) :: b =0._r_kind ! variational quality control parameter - real(r_kind) :: pg =0._r_kind ! variational quality control parameter - real(r_kind) :: jb =0._r_kind ! variational quality control parameter - integer(i_kind) :: ib =0_i_kind ! new variational quality control parameter - integer(i_kind) :: ik =0_i_kind ! new variational quality control parameter - real(r_kind) :: wij(4) =0._r_kind ! horizontal interpolation weights - real(r_kind) :: ppertb =0._r_kind ! random number adding to the obs - integer(i_kind) :: ij(4) =0_i_kind ! horizontal locations - integer(i_kind) :: kx =0_i_kind ! ob type - - !type(ps_ob_type),pointer :: llpoint => NULL() - !real(r_kind) :: time ! observation time in sec - !logical :: luse ! flag indicating if ob is used in pen. - !integer(i_kind) :: idv,iob ! device id and obs index for sorting - !real (r_kind) :: elat, elon ! earth lat-lon for redistribution - !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution - contains - procedure,nopass:: mytype - procedure:: setHop => obsNode_setHop_ - procedure:: xread => obsNode_xread_ - procedure:: xwrite => obsNode_xwrite_ - procedure:: isvalid => obsNode_isvalid_ - procedure::gettlddp => gettlddp_ - - ! procedure, nopass:: headerRead => obsHeader_read_ - ! procedure, nopass:: headerWrite => obsHeader_write_ - ! procedure:: init => obsNode_init_ - ! procedure:: clean => obsNode_clean_ - end type psNode - - public:: psNode_typecast - public:: psNode_nextcast - interface psNode_typecast; module procedure typecast_ ; end interface - interface psNode_nextcast; module procedure nextcast_ ; end interface - - public:: psNode_appendto - interface psNode_appendto; module procedure appendto_ ; end interface - - character(len=*),parameter:: MYNAME="m_psNode" - -#include "myassert.H" -#include "mytrace.H" -contains -function typecast_(aNode) result(ptr_) -!-- cast a class(obsNode) to a type(psNode) - use m_obsNode, only: obsNode - implicit none - type (psNode ),pointer:: ptr_ - class(obsNode),pointer,intent(in):: aNode - ptr_ => null() - if(.not.associated(aNode)) return - ! logically, typecast of a null-reference is a null pointer. - select type(aNode) - type is(psNode) - ptr_ => aNode - end select -return -end function typecast_ - -function nextcast_(aNode) result(ptr_) -!-- cast an obsNode_next(obsNode) to a type(psNode) - use m_obsNode, only: obsNode,obsNode_next - implicit none - type (psNode ),pointer:: ptr_ - class(obsNode),target ,intent(in):: aNode - - class(obsNode),pointer:: inode_ - inode_ => obsNode_next(aNode) - ptr_ => typecast_(inode_) -return -end function nextcast_ - -subroutine appendto_(aNode,oll) -!-- append aNode to linked-list oLL - use m_obsNode , only: obsNode - use m_obsLList, only: obsLList,obsLList_appendNode - implicit none - type(psNode),pointer,intent(in):: aNode - type(obsLList),intent(inout):: oLL - - class(obsNode),pointer:: inode_ - inode_ => aNode - call obsLList_appendNode(oLL,inode_) - inode_ => null() -end subroutine appendto_ - -! obsNode implementations - -function mytype() - implicit none - character(len=:),allocatable:: mytype - mytype="[psNode]" -end function mytype - -subroutine obsNode_xread_(aNode,iunit,istat,diagLookup,skip) - use m_obsdiagNode, only: obsdiagLookup_locate - use m_obsdiagNode, only: obs_diags - implicit none - class(psNode), intent(inout):: aNode - integer(i_kind), intent(in ):: iunit - integer(i_kind), intent( out):: istat - type(obs_diags), intent(in ):: diagLookup - logical,optional,intent(in ):: skip - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_xread_' - logical:: skip_ -_ENTRY_(myname_) - - skip_=.false. - if(present(skip)) skip_=skip - - if(skip_) then - read(iunit,iostat=istat) - if(istat/=0) then - call perr(myname_,'skipping read(%(res,...)), iostat =',istat) - _EXIT_(myname_) - return - endif - - else - read(iunit,iostat=istat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%jb , & - aNode%ib , & - aNode%ik , & - aNode%ppertb , & - aNode%kx , & - aNode%wij , & - aNode%ij - if (istat/=0) then - call perr(myname_,'read(%(res,...)), iostat =',istat) - _EXIT_(myname_) - return - end if - - aNode%diags => obsdiagLookup_locate(diagLookup,aNode%idv,aNode%iob,1_i_kind) - if(.not. associated(aNode%diags)) then - call perr(myname_,'obsdiagLookup_locate(), %idv =',aNode%idv) - call perr(myname_,' %iob =',aNode%iob) - call die(myname_) - endif - endif ! if(skip_); else -_EXIT_(myname_) -return -end subroutine obsNode_xread_ - -subroutine obsNode_xwrite_(aNode,junit,jstat) - implicit none - class(psNode),intent(in):: aNode - integer(i_kind),intent(in ):: junit - integer(i_kind),intent( out):: jstat - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_xwrite_' -_ENTRY_(myname_) - - write(junit,iostat=jstat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%jb , & - aNode%ib , & - aNode%ik , & - aNode%ppertb , & - aNode%kx , & - aNode%wij , & - aNode%ij - if (jstat/=0) then - call perr(myname_,'write(%(res,...)), iostat =',jstat) - _EXIT_(myname_) - return - end if -_EXIT_(myname_) -return -end subroutine obsNode_xwrite_ - -subroutine obsNode_setHop_(aNode) - use m_cvgridLookup, only: cvgridLookup_getiw - implicit none - class(psNode),intent(inout):: aNode - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_setHop_' -_ENTRY_(myname_) - call cvgridLookup_getiw(aNode%elat,aNode%elon,aNode%ij,aNode%wij) -_EXIT_(myname_) -return -end subroutine obsNode_setHop_ - -function obsNode_isvalid_(aNode) result(isvalid_) - implicit none - logical:: isvalid_ - class(psNode),intent(in):: aNode - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_isvalid_' -_ENTRY_(myname_) - isvalid_= associated(aNode%diags) -_EXIT_(myname_) -return -end function obsNode_isvalid_ - -pure subroutine gettlddp_(aNode,jiter,tlddp,nob) - use kinds, only: r_kind - implicit none - class(psNode), intent(in):: aNode - integer(kind=i_kind),intent(in):: jiter - real(kind=r_kind),intent(inout):: tlddp - integer(kind=i_kind),optional,intent(inout):: nob - - tlddp = tlddp + aNode%diags%tldepart(jiter)*aNode%diags%tldepart(jiter) - if(present(nob)) nob=nob+1 -return -end subroutine gettlddp_ - -end module m_psNode diff --git a/src/gsi/m_psnode.F90 b/src/gsi/m_psnode.F90 new file mode 100644 index 0000000000..575e47af8a --- /dev/null +++ b/src/gsi/m_psnode.F90 @@ -0,0 +1,260 @@ +module m_psnode +!$$$ subprogram documentation block +! . . . . +! subprogram: module m_psnode +! prgmmr: j guo +! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 +! date: 2016-05-18 +! +! abstract: class-module of obs-type psnode (surface pressure) +! +! program history log: +! 2016-05-18 j guo - added this document block for the initial polymorphic +! implementation. +! 2019-09-20 X.Su - add new variational QC parameters +! +! input argument list: see Fortran 90 style document below +! +! output argument list: see Fortran 90 style document below +! +! attributes: +! language: Fortran 90 and/or above +! machine: +! +!$$$ end subprogram documentation block + +! module interface: + use kinds , only: i_kind,r_kind + use mpeu_util, only: assert_,die,perr,warn,tell + use m_obsdiagnode, only: obs_diag + use m_obsnode, only: obsnode + implicit none + private + + public:: psnode + + type,extends(obsnode):: psnode + !private + type(obs_diag), pointer :: diags => null() + real(r_kind) :: res =0._r_kind ! surface pressure residual + real(r_kind) :: err2 =0._r_kind ! surface pressure error squared + ! in reciprocal + real(r_kind) :: raterr2=0._r_kind ! square of ratio of final obs error + ! to original obs error + real(r_kind) :: b =0._r_kind ! variational quality control parameter + real(r_kind) :: pg =0._r_kind ! variational quality control parameter + real(r_kind) :: jb =0._r_kind ! variational quality control parameter + integer(i_kind) :: ib =0 ! new variational quality control parameter + integer(i_kind) :: ik =0 ! new variational quality control parameter + real(r_kind) :: wij(4) =0._r_kind ! horizontal interpolation weights + real(r_kind) :: ppertb =0._r_kind ! random number adding to the obs + integer(i_kind) :: ij(4) =0 ! horizontal locations + integer(i_kind) :: kx =0 ! ob type + + !type(ps_ob_type),pointer :: llpoint => null() + !real(r_kind) :: time ! observation time in sec + !logical :: luse ! flag indicating if ob is used in pen. + !integer(i_kind) :: idv,iob ! device id and obs index for sorting + !real (r_kind) :: elat, elon ! earth lat-lon for redistribution + !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution + contains + procedure,nopass:: mytype + procedure:: sethop => obsnode_sethop_ + procedure:: xread => obsnode_xread_ + procedure:: xwrite => obsnode_xwrite_ + procedure:: isvalid => obsnode_isvalid_ + procedure::gettlddp => gettlddp_ + + ! procedure, nopass:: headerread => obsheader_read_ + ! procedure, nopass:: headerwrite => obsheader_write_ + ! procedure:: init => obsnode_init_ + ! procedure:: clean => obsnode_clean_ + end type psnode + + public:: psnode_typecast + public:: psnode_nextcast + interface psnode_typecast; module procedure typecast_ ; end interface + interface psnode_nextcast; module procedure nextcast_ ; end interface + + public:: psnode_appendto + interface psnode_appendto; module procedure appendto_ ; end interface + + character(len=*),parameter:: myname="m_psnode" + +#include "myassert.H" +#include "mytrace.H" +contains +function typecast_(anode) result(ptr_) +!-- cast a class(obsnode) to a type(psnode) + use m_obsnode, only: obsnode + implicit none + type (psnode ),pointer:: ptr_ + class(obsnode),pointer,intent(in):: anode + ptr_ => null() + if(.not.associated(anode)) return + ! logically, typecast of a null-reference is a null pointer. + select type(anode) + type is(psnode) + ptr_ => anode + end select + return +end function typecast_ + +function nextcast_(anode) result(ptr_) +!-- cast an obsnode_next(obsnode) to a type(psnode) + use m_obsnode, only: obsnode,obsnode_next + implicit none + type (psnode ),pointer:: ptr_ + class(obsnode),target ,intent(in):: anode + + class(obsnode),pointer:: inode_ + inode_ => obsnode_next(anode) + ptr_ => typecast_(inode_) + return +end function nextcast_ + +subroutine appendto_(anode,oll) +!-- append anode to linked-list oll + use m_obsnode , only: obsnode + use m_obsllist, only: obsllist,obsllist_appendnode + implicit none + type(psnode),pointer,intent(in):: anode + type(obsllist),intent(inout):: oll + + class(obsnode),pointer:: inode_ + inode_ => anode + call obsllist_appendnode(oll,inode_) + inode_ => null() +end subroutine appendto_ + +! obsnode implementations + +function mytype() + implicit none + character(len=:),allocatable:: mytype + mytype="[psnode]" +end function mytype + +subroutine obsnode_xread_(anode,iunit,istat,diaglookup,skip) + use m_obsdiagnode, only: obsdiaglookup_locate + use m_obsdiagnode, only: obs_diags + implicit none + class(psnode), intent(inout):: anode + integer(i_kind), intent(in ):: iunit + integer(i_kind), intent( out):: istat + type(obs_diags), intent(in ):: diaglookup + logical,optional,intent(in ):: skip + + character(len=*),parameter:: myname_=myname//'::obsnode_xread_' + logical:: skip_ +_ENTRY_(myname_) + + skip_=.false. + if(present(skip)) skip_=skip + + if(skip_) then + read(iunit,iostat=istat) + if(istat/=0) then + call perr(myname_,'skipping read(%(res,...)), iostat =',istat) + _EXIT_(myname_) + return + endif + + else + read(iunit,iostat=istat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%jb , & + anode%ib , & + anode%ik , & + anode%ppertb , & + anode%kx , & + anode%wij , & + anode%ij + if (istat/=0) then + call perr(myname_,'read(%(res,...)), iostat =',istat) + _EXIT_(myname_) + return + end if + + anode%diags => obsdiaglookup_locate(diaglookup,anode%idv,anode%iob,1) + if(.not. associated(anode%diags)) then + call perr(myname_,'obsdiaglookup_locate(), %idv =',anode%idv) + call perr(myname_,' %iob =',anode%iob) + call die(myname_) + endif + endif ! if(skip_); else +_EXIT_(myname_) + return +end subroutine obsnode_xread_ + +subroutine obsnode_xwrite_(anode,junit,jstat) + implicit none + class(psnode),intent(in):: anode + integer(i_kind),intent(in ):: junit + integer(i_kind),intent( out):: jstat + + character(len=*),parameter:: myname_=myname//'::obsnode_xwrite_' +_ENTRY_(myname_) + + write(junit,iostat=jstat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%jb , & + anode%ib , & + anode%ik , & + anode%ppertb , & + anode%kx , & + anode%wij , & + anode%ij + if (jstat/=0) then + call perr(myname_,'write(%(res,...)), iostat =',jstat) + _EXIT_(myname_) + return + end if +_EXIT_(myname_) + return +end subroutine obsnode_xwrite_ + +subroutine obsnode_sethop_(anode) + use m_cvgridlookup, only: cvgridlookup_getiw + implicit none + class(psnode),intent(inout):: anode + + character(len=*),parameter:: myname_=myname//'::obsnode_sethop_' +_ENTRY_(myname_) + call cvgridlookup_getiw(anode%elat,anode%elon,anode%ij,anode%wij) +_EXIT_(myname_) + return +end subroutine obsnode_sethop_ + +function obsnode_isvalid_(anode) result(isvalid_) + implicit none + logical:: isvalid_ + class(psnode),intent(in):: anode + + character(len=*),parameter:: myname_=myname//'::obsnode_isvalid_' +_ENTRY_(myname_) + isvalid_= associated(anode%diags) +_EXIT_(myname_) + return +end function obsnode_isvalid_ + +pure subroutine gettlddp_(anode,jiter,tlddp,nob) + use kinds, only: r_kind + implicit none + class(psnode), intent(in):: anode + integer(kind=i_kind),intent(in):: jiter + real(kind=r_kind),intent(inout):: tlddp + integer(kind=i_kind),optional,intent(inout):: nob + + tlddp = tlddp + anode%diags%tldepart(jiter)*anode%diags%tldepart(jiter) + if(present(nob)) nob=nob+1 + return +end subroutine gettlddp_ + +end module m_psnode diff --git a/src/gsi/m_pwNode.F90 b/src/gsi/m_pwNode.F90 deleted file mode 100644 index 7ad1f390dc..0000000000 --- a/src/gsi/m_pwNode.F90 +++ /dev/null @@ -1,315 +0,0 @@ -module m_pwNode -!$$$ subprogram documentation block -! . . . . -! subprogram: module m_pwNode -! prgmmr: j guo -! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 -! date: 2016-05-18 -! -! abstract: class-module of obs-type pwNode (total column water) -! -! program history log: -! 2016-05-18 j guo - added this document block for the initial polymorphic -! implementation. -! -! input argument list: see Fortran 90 style document below -! -! output argument list: see Fortran 90 style document below -! -! attributes: -! language: Fortran 90 and/or above -! machine: -! -!$$$ end subprogram documentation block - -! module interface: - use m_obsdiagNode, only: obs_diag - use m_obsdiagNode, only: obs_diags - use kinds , only: i_kind,r_kind - use mpeu_util, only: assert_,die,perr,warn,tell - use m_obsNode, only: obsNode - implicit none - private - - public:: pwNode - - type,extends(obsNode):: pwNode - !type(pw_ob_type),pointer :: llpoint => NULL() - type(obs_diag), pointer :: diags => NULL() - real(r_kind) :: res ! precipitable water residual - real(r_kind) :: err2 ! precipitable water error squared - real(r_kind) :: raterr2 ! square of ratio of final obs error - ! to original obs error - !real(r_kind) :: time ! observation time in sec - real(r_kind) :: b ! variational quality control parameter - real(r_kind) :: pg ! variational quality control parameter - real(r_kind) :: wij(4) ! horizontal interpolation weights - real(r_kind),dimension(:),pointer :: dp => NULL() - ! delta pressure at mid layers at obs locations - integer(i_kind) :: ij(4) ! horizontal locations - !logical :: luse ! flag indicating if ob is used in pen. - - !integer(i_kind) :: idv,iob ! device id and obs index for sorting - !real (r_kind) :: elat, elon ! earth lat-lon for redistribution - !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution - contains - procedure,nopass:: mytype - procedure:: setHop => obsNode_setHop_ - procedure:: xread => obsNode_xread_ - procedure:: xwrite => obsNode_xwrite_ - procedure:: isvalid => obsNode_isvalid_ - procedure:: gettlddp => gettlddp_ - - procedure, nopass:: headerRead => obsHeader_read_ - procedure, nopass:: headerWrite => obsHeader_write_ - procedure:: init => obsNode_init_ - procedure:: clean => obsNode_clean_ - end type pwNode - - public:: pwNode_typecast - public:: pwNode_nextcast - interface pwNode_typecast; module procedure typecast_ ; end interface - interface pwNode_nextcast; module procedure nextcast_ ; end interface - - public:: pwNode_appendto - interface pwNode_appendto; module procedure appendto_ ; end interface - - character(len=*),parameter:: MYNAME="m_pwNode" - -#include "myassert.H" -#include "mytrace.H" -contains -function typecast_(aNode) result(ptr_) -!-- cast a class(obsNode) to a type(pwNode) - use m_obsNode, only: obsNode - implicit none - type(pwNode ),pointer:: ptr_ - class(obsNode),pointer,intent(in):: aNode - ptr_ => null() - if(.not.associated(aNode)) return - ! logically, typecast of a null-reference is a null pointer. - select type(aNode) - type is(pwNode) - ptr_ => aNode - end select -return -end function typecast_ - -function nextcast_(aNode) result(ptr_) -!-- cast an obsNode_next(obsNode) to a type(pwNode) - use m_obsNode, only: obsNode,obsNode_next - implicit none - type(pwNode),pointer:: ptr_ - class(obsNode),target,intent(in):: aNode - - class(obsNode),pointer:: inode_ - inode_ => obsNode_next(aNode) - ptr_ => typecast_(inode_) -return -end function nextcast_ - -subroutine appendto_(aNode,oll) -!-- append aNode to linked-list oLL - use m_obsNode , only: obsNode - use m_obsLList, only: obsLList,obsLList_appendNode - implicit none - type(pwNode),pointer,intent(in):: aNode - type(obsLList),intent(inout):: oLL - - class(obsNode),pointer:: inode_ - inode_ => aNode - call obsLList_appendNode(oLL,inode_) - inode_ => null() -end subroutine appendto_ - -! obsNode implementations - -function mytype() - implicit none - character(len=:),allocatable:: mytype - mytype="[pwNode]" -end function mytype - -subroutine obsHeader_read_(iunit,mobs,jread,istat) - use gridmod, only: nsig - implicit none - integer(i_kind),intent(in ):: iunit - integer(i_kind),intent(out):: mobs - integer(i_kind),intent(out):: jread - integer(i_kind),intent(out):: istat - - character(len=*),parameter:: myname_=myname//".obsHeader_read_" - integer(i_kind):: msig -_ENTRY_(myname_) - - read(iunit,iostat=istat) mobs,jread, msig - if(istat==0 .and. nsig/=msig) then - call perr(myname_,'unexpected dimension information, nsig =',nsig) - call perr(myname_,' but read msig =',msig) - call die(myname_) - endif -_EXIT_(myname_) -return -end subroutine obsHeader_read_ - -subroutine obsHeader_write_(junit,mobs,jwrite,jstat) - use gridmod, only: nsig - implicit none - integer(i_kind),intent(in ):: junit - integer(i_kind),intent(in ):: mobs - integer(i_kind),intent(in ):: jwrite - integer(i_kind),intent(out):: jstat - - character(len=*),parameter:: myname_=myname//".obsHeader_write_" -_ENTRY_(myname_) - write(junit,iostat=jstat) mobs,jwrite, nsig -_EXIT_(myname_) -return -end subroutine obsHeader_write_ - -subroutine obsNode_init_(aNode) - use gridmod, only: nsig - implicit none - class(pwNode),intent(out):: aNode - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_init_' -_ENTRY_(myname_) - aNode%llpoint => null() - aNode%luse = .false. - aNode%elat = 0._r_kind - aNode%elon = 0._r_kind - aNode%time = 0._r_kind - aNode%idv =-1 - aNode%iob =-1 - allocate(aNode%dp(nsig)) -_EXIT_(myname_) -return -end subroutine obsNode_init_ - -subroutine obsNode_clean_(aNode) - implicit none - class(pwNode),intent(inout):: aNode - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_clean_' -_ENTRY_(myname_) -!_TRACEV_(myname_,'%mytype() =',aNode%mytype()) - if(associated(aNode%dp)) deallocate(aNode%dp) -_EXIT_(myname_) -return -end subroutine obsNode_clean_ - -subroutine obsNode_xread_(aNode,iunit,istat,diagLookup,skip) - use m_obsdiagNode, only: obsdiagLookup_locate - implicit none - class(pwNode) , intent(inout):: aNode - integer(i_kind) , intent(in ):: iunit - integer(i_kind) , intent( out):: istat - type(obs_diags) , intent(in ):: diagLookup - logical,optional, intent(in ):: skip - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xread_' - logical:: skip_ -_ENTRY_(myname_) - skip_=.false. - if(present(skip)) skip_=skip - - istat=0 - if(skip_) then - read(iunit,iostat=istat) - if(istat/=0) then - call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - endif - - else - read(iunit,iostat=istat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%dp , & - aNode%wij , & - aNode%ij - if (istat/=0) then - call perr(myname_,'read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - end if - - aNode%diags => obsdiagLookup_locate(diagLookup,aNode%idv,aNode%iob,1_i_kind) - if(.not.associated(aNode%diags)) then - call perr(myname_,'obsdiagLookup_locate(), %idv =',aNode%idv) - call perr(myname_,' %iob =',aNode%iob) - call die(myname_) - endif - endif -_EXIT_(myname_) -return -end subroutine obsNode_xread_ - -subroutine obsNode_xwrite_(aNode,junit,jstat) - implicit none - class(pwNode),intent(in):: aNode - integer(i_kind),intent(in ):: junit - integer(i_kind),intent( out):: jstat - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xwrite_' -_ENTRY_(myname_) - - jstat=0 - write(junit,iostat=jstat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%dp , & - aNode%wij , & - aNode%ij - if (jstat/=0) then - call perr(myname_,'write(%(res,err2,...)), iostat =',jstat) - _EXIT_(myname_) - return - end if -_EXIT_(myname_) -return -end subroutine obsNode_xwrite_ - -subroutine obsNode_setHop_(aNode) - use m_cvgridLookup, only: cvgridLookup_getiw - implicit none - class(pwNode),intent(inout):: aNode - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_setHop_' -_ENTRY_(myname_) - call cvgridLookup_getiw(aNode%elat,aNode%elon,aNode%ij,aNode%wij) -_EXIT_(myname_) -return -end subroutine obsNode_setHop_ - -function obsNode_isvalid_(aNode) result(isvalid_) - implicit none - logical:: isvalid_ - class(pwNode),intent(in):: aNode - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_isvalid_' -_ENTRY_(myname_) - isvalid_=associated(aNode%diags) -_EXIT_(myname_) -return -end function obsNode_isvalid_ - -pure subroutine gettlddp_(aNode,jiter,tlddp,nob) - use kinds, only: r_kind - implicit none - class(pwNode), intent(in):: aNode - integer(kind=i_kind),intent(in):: jiter - real(kind=r_kind),intent(inout):: tlddp - integer(kind=i_kind),optional,intent(inout):: nob - - tlddp = tlddp + aNode%diags%tldepart(jiter)*aNode%diags%tldepart(jiter) - if(present(nob)) nob=nob+1 -return -end subroutine gettlddp_ - -end module m_pwNode diff --git a/src/gsi/m_pwnode.F90 b/src/gsi/m_pwnode.F90 new file mode 100644 index 0000000000..c73ff95069 --- /dev/null +++ b/src/gsi/m_pwnode.F90 @@ -0,0 +1,315 @@ +module m_pwnode +!$$$ subprogram documentation block +! . . . . +! subprogram: module m_pwnode +! prgmmr: j guo +! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 +! date: 2016-05-18 +! +! abstract: class-module of obs-type pwnode (total column water) +! +! program history log: +! 2016-05-18 j guo - added this document block for the initial polymorphic +! implementation. +! +! input argument list: see Fortran 90 style document below +! +! output argument list: see Fortran 90 style document below +! +! attributes: +! language: Fortran 90 and/or above +! machine: +! +!$$$ end subprogram documentation block + +! module interface: + use m_obsdiagnode, only: obs_diag + use m_obsdiagnode, only: obs_diags + use kinds , only: i_kind,r_kind + use mpeu_util, only: assert_,die,perr,warn,tell + use m_obsnode, only: obsnode + implicit none + private + + public:: pwnode + + type,extends(obsnode):: pwnode + !type(pw_ob_type),pointer :: llpoint => null() + type(obs_diag), pointer :: diags => null() + real(r_kind) :: res ! precipitable water residual + real(r_kind) :: err2 ! precipitable water error squared + real(r_kind) :: raterr2 ! square of ratio of final obs error + ! to original obs error + !real(r_kind) :: time ! observation time in sec + real(r_kind) :: b ! variational quality control parameter + real(r_kind) :: pg ! variational quality control parameter + real(r_kind) :: wij(4) ! horizontal interpolation weights + real(r_kind),dimension(:),pointer :: dp => null() + ! delta pressure at mid layers at obs locations + integer(i_kind) :: ij(4) ! horizontal locations + !logical :: luse ! flag indicating if ob is used in pen. + + !integer(i_kind) :: idv,iob ! device id and obs index for sorting + !real (r_kind) :: elat, elon ! earth lat-lon for redistribution + !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution + contains + procedure,nopass:: mytype + procedure:: sethop => obsnode_sethop_ + procedure:: xread => obsnode_xread_ + procedure:: xwrite => obsnode_xwrite_ + procedure:: isvalid => obsnode_isvalid_ + procedure:: gettlddp => gettlddp_ + + procedure, nopass:: headerread => obsheader_read_ + procedure, nopass:: headerwrite => obsheader_write_ + procedure:: init => obsnode_init_ + procedure:: clean => obsnode_clean_ + end type pwnode + + public:: pwnode_typecast + public:: pwnode_nextcast + interface pwnode_typecast; module procedure typecast_ ; end interface + interface pwnode_nextcast; module procedure nextcast_ ; end interface + + public:: pwnode_appendto + interface pwnode_appendto; module procedure appendto_ ; end interface + + character(len=*),parameter:: myname="m_pwnode" + +#include "myassert.H" +#include "mytrace.H" +contains +function typecast_(anode) result(ptr_) +!-- cast a class(obsnode) to a type(pwnode) + use m_obsnode, only: obsnode + implicit none + type(pwnode ),pointer:: ptr_ + class(obsnode),pointer,intent(in):: anode + ptr_ => null() + if(.not.associated(anode)) return + ! logically, typecast of a null-reference is a null pointer. + select type(anode) + type is(pwnode) + ptr_ => anode + end select + return +end function typecast_ + +function nextcast_(anode) result(ptr_) +!-- cast an obsnode_next(obsnode) to a type(pwnode) + use m_obsnode, only: obsnode,obsnode_next + implicit none + type(pwnode),pointer:: ptr_ + class(obsnode),target,intent(in):: anode + + class(obsnode),pointer:: inode_ + inode_ => obsnode_next(anode) + ptr_ => typecast_(inode_) + return +end function nextcast_ + +subroutine appendto_(anode,oll) +!-- append anode to linked-list oll + use m_obsnode , only: obsnode + use m_obsllist, only: obsllist,obsllist_appendnode + implicit none + type(pwnode),pointer,intent(in):: anode + type(obsllist),intent(inout):: oll + + class(obsnode),pointer:: inode_ + inode_ => anode + call obsllist_appendnode(oll,inode_) + inode_ => null() +end subroutine appendto_ + +! obsnode implementations + +function mytype() + implicit none + character(len=:),allocatable:: mytype + mytype="[pwnode]" +end function mytype + +subroutine obsheader_read_(iunit,mobs,jread,istat) + use gridmod, only: nsig + implicit none + integer(i_kind),intent(in ):: iunit + integer(i_kind),intent(out):: mobs + integer(i_kind),intent(out):: jread + integer(i_kind),intent(out):: istat + + character(len=*),parameter:: myname_=myname//".obsheader_read_" + integer(i_kind):: msig +_ENTRY_(myname_) + + read(iunit,iostat=istat) mobs,jread, msig + if(istat==0 .and. nsig/=msig) then + call perr(myname_,'unexpected dimension information, nsig =',nsig) + call perr(myname_,' but read msig =',msig) + call die(myname_) + endif +_EXIT_(myname_) + return +end subroutine obsheader_read_ + +subroutine obsheader_write_(junit,mobs,jwrite,jstat) + use gridmod, only: nsig + implicit none + integer(i_kind),intent(in ):: junit + integer(i_kind),intent(in ):: mobs + integer(i_kind),intent(in ):: jwrite + integer(i_kind),intent(out):: jstat + + character(len=*),parameter:: myname_=myname//".obsheader_write_" +_ENTRY_(myname_) + write(junit,iostat=jstat) mobs,jwrite, nsig +_EXIT_(myname_) + return +end subroutine obsheader_write_ + +subroutine obsnode_init_(anode) + use gridmod, only: nsig + implicit none + class(pwnode),intent(out):: anode + + character(len=*),parameter:: myname_=myname//'.obsnode_init_' +_ENTRY_(myname_) + anode%llpoint => null() + anode%luse = .false. + anode%elat = 0._r_kind + anode%elon = 0._r_kind + anode%time = 0._r_kind + anode%idv =-1 + anode%iob =-1 + allocate(anode%dp(nsig)) +_EXIT_(myname_) + return +end subroutine obsnode_init_ + +subroutine obsnode_clean_(anode) + implicit none + class(pwnode),intent(inout):: anode + + character(len=*),parameter:: myname_=myname//'.obsnode_clean_' +_ENTRY_(myname_) +!_TRACEV_(myname_,'%mytype() =',anode%mytype()) + if(associated(anode%dp)) deallocate(anode%dp) +_EXIT_(myname_) + return +end subroutine obsnode_clean_ + +subroutine obsnode_xread_(anode,iunit,istat,diaglookup,skip) + use m_obsdiagnode, only: obsdiaglookup_locate + implicit none + class(pwnode) , intent(inout):: anode + integer(i_kind) , intent(in ):: iunit + integer(i_kind) , intent( out):: istat + type(obs_diags) , intent(in ):: diaglookup + logical,optional, intent(in ):: skip + + character(len=*),parameter:: myname_=myname//'.obsnode_xread_' + logical:: skip_ +_ENTRY_(myname_) + skip_=.false. + if(present(skip)) skip_=skip + + istat=0 + if(skip_) then + read(iunit,iostat=istat) + if(istat/=0) then + call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + endif + + else + read(iunit,iostat=istat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%dp , & + anode%wij , & + anode%ij + if (istat/=0) then + call perr(myname_,'read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + end if + + anode%diags => obsdiaglookup_locate(diaglookup,anode%idv,anode%iob,1) + if(.not.associated(anode%diags)) then + call perr(myname_,'obsdiaglookup_locate(), %idv =',anode%idv) + call perr(myname_,' %iob =',anode%iob) + call die(myname_) + endif + endif +_EXIT_(myname_) + return +end subroutine obsnode_xread_ + +subroutine obsnode_xwrite_(anode,junit,jstat) + implicit none + class(pwnode),intent(in):: anode + integer(i_kind),intent(in ):: junit + integer(i_kind),intent( out):: jstat + + character(len=*),parameter:: myname_=myname//'.obsnode_xwrite_' +_ENTRY_(myname_) + + jstat=0 + write(junit,iostat=jstat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%dp , & + anode%wij , & + anode%ij + if (jstat/=0) then + call perr(myname_,'write(%(res,err2,...)), iostat =',jstat) + _EXIT_(myname_) + return + end if +_EXIT_(myname_) + return +end subroutine obsnode_xwrite_ + +subroutine obsnode_sethop_(anode) + use m_cvgridlookup, only: cvgridlookup_getiw + implicit none + class(pwnode),intent(inout):: anode + + character(len=*),parameter:: myname_=myname//'::obsnode_sethop_' +_ENTRY_(myname_) + call cvgridlookup_getiw(anode%elat,anode%elon,anode%ij,anode%wij) +_EXIT_(myname_) + return +end subroutine obsnode_sethop_ + +function obsnode_isvalid_(anode) result(isvalid_) + implicit none + logical:: isvalid_ + class(pwnode),intent(in):: anode + + character(len=*),parameter:: myname_=myname//'::obsnode_isvalid_' +_ENTRY_(myname_) + isvalid_=associated(anode%diags) +_EXIT_(myname_) + return +end function obsnode_isvalid_ + +pure subroutine gettlddp_(anode,jiter,tlddp,nob) + use kinds, only: r_kind + implicit none + class(pwnode), intent(in):: anode + integer(kind=i_kind),intent(in):: jiter + real(kind=r_kind),intent(inout):: tlddp + integer(kind=i_kind),optional,intent(inout):: nob + + tlddp = tlddp + anode%diags%tldepart(jiter)*anode%diags%tldepart(jiter) + if(present(nob)) nob=nob+1 + return +end subroutine gettlddp_ + +end module m_pwnode diff --git a/src/gsi/m_qNode.F90 b/src/gsi/m_qNode.F90 deleted file mode 100644 index b6b24d09f1..0000000000 --- a/src/gsi/m_qNode.F90 +++ /dev/null @@ -1,282 +0,0 @@ -module m_qNode -!$$$ subprogram documentation block -! . . . . -! subprogram: module m_qNode -! prgmmr: j guo -! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 -! date: 2016-05-18 -! -! abstract: class-module of obs-type qNode (moisture) -! -! program history log: -! 2016-05-18 j guo - added this document block for the initial polymorphic -! implementation. -! 2019-09-20 X.Su - add new variational QC parameters -! -! input argument list: see Fortran 90 style document below -! -! output argument list: see Fortran 90 style document below -! -! attributes: -! language: Fortran 90 and/or above -! machine: -! -!$$$ end subprogram documentation block - -! module interface: - use m_obsdiagNode, only: obs_diag - use m_obsdiagNode, only: obs_diags - use kinds , only: i_kind,r_kind - use mpeu_util, only: assert_,die,perr,warn,tell - use m_obsNode, only: obsNode - implicit none - private - - public:: qNode - - type,extends(obsNode):: qNode - !type(q_ob_type),pointer :: llpoint => NULL() - type(obs_diag), pointer :: diags => NULL() - real(r_kind) :: res ! moisture residual - real(r_kind) :: err2 ! moisture error squared - real(r_kind) :: raterr2 ! square of ratio of final obs error - ! to original obs error - !real(r_kind) :: time ! observation time in sec - real(r_kind) :: b ! variational quality control parameter - real(r_kind) :: pg ! variational quality control parameter - real(r_kind) :: jb ! variational quality control parameter - integer(i_kind) :: ib ! new variational quality control parameter - integer(i_kind) :: ik ! new variational quality control parameter - real(r_kind) :: wij(8) ! horizontal interpolation weights - real(r_kind) :: qpertb ! random number adding to the obs - integer(i_kind) :: ij(8) ! horizontal locations - integer(i_kind) :: k1 ! level of errtable 1-33 - integer(i_kind) :: kx ! ob type - !logical :: luse ! flag indicating if ob is used in pen. - - !integer(i_kind) :: idv,iob ! device id and obs index for sorting - !real (r_kind) :: elat, elon ! earth lat-lon for redistribution - !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution - real (r_kind) :: dlev ! reference to the vertical grid - - integer(i_kind) :: ich0=0 ! ich code to mark derived data. See - ! qNode_ich0 and qNode_ich0_PBL_Pseudo below - contains - procedure,nopass:: mytype - procedure:: setHop => obsNode_setHop_ - procedure:: xread => obsNode_xread_ - procedure:: xwrite => obsNode_xwrite_ - procedure:: isvalid => obsNode_isvalid_ - procedure:: gettlddp => gettlddp_ - - ! procedure, nopass:: headerRead => obsHeader_read_ - ! procedure, nopass:: headerWrite => obsHeader_write_ - ! procedure:: init => obsNode_init_ - ! procedure:: clean => obsNode_clean_ - end type qNode - - public:: qNode_typecast - public:: qNode_nextcast - interface qNode_typecast; module procedure typecast_ ; end interface - interface qNode_nextcast; module procedure nextcast_ ; end interface - - public:: qNode_appendto - interface qNode_appendto; module procedure appendto_ ; end interface - - ! Because there are two components in qNode for an ordinary wind obs, - ! ich values are set to (1,2). Therefore, ich values for PBL_pseudo_surfobsUV - ! are set to (3,4), and qNode_ich0_pbl_pseudo is set to 2. - - public:: qNode_ich0 - public:: qNode_ich0_PBL_pseudo - - integer(i_kind),parameter :: qNode_ich0 = 0 ! ich=0+1 - integer(i_kind),parameter :: qNode_ich0_PBL_pseudo = qNode_ich0+1 ! ich=1+1 - - character(len=*),parameter:: MYNAME="m_qNode" - -#include "myassert.H" -#include "mytrace.H" -contains -function typecast_(aNode) result(ptr_) -!-- cast a class(obsNode) to a type(qNode) - use m_obsNode, only: obsNode - implicit none - type(qNode ),pointer:: ptr_ - class(obsNode),pointer,intent(in):: aNode - ptr_ => null() - if(.not.associated(aNode)) return - ! logically, typecast of a null-reference is a null pointer. - select type(aNode) - type is(qNode) - ptr_ => aNode - end select -return -end function typecast_ - -function nextcast_(aNode) result(ptr_) -!-- cast an obsNode_next(obsNode) to a type(qNode) - use m_obsNode, only: obsNode,obsNode_next - implicit none - type(qNode ),pointer:: ptr_ - class(obsNode),target ,intent(in):: aNode - - class(obsNode),pointer:: inode_ - inode_ => obsNode_next(aNode) - ptr_ => typecast_(inode_) -return -end function nextcast_ - -subroutine appendto_(aNode,oll) -!-- append aNode to linked-list oLL - use m_obsNode , only: obsNode - use m_obsLList, only: obsLList,obsLList_appendNode - implicit none - type(qNode),pointer,intent(in):: aNode - type(obsLList),intent(inout):: oLL - - class(obsNode),pointer:: inode_ - inode_ => aNode - call obsLList_appendNode(oLL,inode_) - inode_ => null() -end subroutine appendto_ - -! obsNode implementations - -function mytype() - implicit none - character(len=:),allocatable:: mytype - mytype="[qNode]" -end function mytype - -subroutine obsNode_xread_(aNode,iunit,istat,diagLookup,skip) - use m_obsdiagNode, only: obsdiagLookup_locate - implicit none - class(qNode) , intent(inout):: aNode - integer(i_kind) , intent(in ):: iunit - integer(i_kind) , intent( out):: istat - type(obs_diags) , intent(in ):: diagLookup - logical,optional, intent(in ):: skip - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xread_' - logical:: skip_ -_ENTRY_(myname_) - - skip_=.false. - if(present(skip)) skip_=skip - - istat=0 - if(skip_) then - read(iunit,iostat=istat) - if(istat/=0) then - call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - endif - - else - read(iunit,iostat=istat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%jb , & - aNode%ib , & - aNode%ik , & - aNode%qpertb , & - aNode%k1 , & - aNode%kx , & - aNode%dlev , & - aNode%ich0 , & - aNode%wij , & - aNode%ij - if (istat/=0) then - call perr(myname_,'read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - end if - - aNode%diags => obsdiagLookup_locate(diagLookup,aNode%idv,aNode%iob,aNode%ich0+1_i_kind) - if(.not.associated(aNode%diags)) then - call perr(myname_,'obsdiagLookup_locate(), %idv =',aNode%idv) - call perr(myname_,' %iob =',aNode%iob) - call perr(myname_,' %ich0 =',aNode%ich0) - call die(myname_) - endif - endif -_EXIT_(myname_) -return -end subroutine obsNode_xread_ - -subroutine obsNode_xwrite_(aNode,junit,jstat) - implicit none - class(qNode),intent(in):: aNode - integer(i_kind),intent(in ):: junit - integer(i_kind),intent( out):: jstat - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xwrite_' -_ENTRY_(myname_) - - jstat=0 - write(junit,iostat=jstat) aNode%res , & - aNode%err2 , & - aNode%raterr2, & - aNode%b , & - aNode%pg , & - aNode%jb , & - aNode%ib , & - aNode%ik , & - aNode%qpertb , & - aNode%k1 , & - aNode%kx , & - aNode%dlev , & - aNode%ich0 , & - aNode%wij , & - aNode%ij - if (jstat/=0) then - call perr(myname_,'write(%(res,err2,...)), iostat =',jstat) - _EXIT_(myname_) - return - end if -_EXIT_(myname_) -return -end subroutine obsNode_xwrite_ - -subroutine obsNode_setHop_(aNode) - use m_cvgridLookup, only: cvgridLookup_getiw - implicit none - class(qNode),intent(inout):: aNode - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_setHop_' -_ENTRY_(myname_) - call cvgridLookup_getiw(aNode%elat,aNode%elon,aNode%dlev,aNode%ij,aNode%wij) -_EXIT_(myname_) -return -end subroutine obsNode_setHop_ - -function obsNode_isvalid_(aNode) result(isvalid_) - implicit none - logical:: isvalid_ - class(qNode),intent(in):: aNode - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_isvalid_' -_ENTRY_(myname_) - isvalid_=associated(aNode%diags) -_EXIT_(myname_) -return -end function obsNode_isvalid_ - -pure subroutine gettlddp_(aNode,jiter,tlddp,nob) - use kinds, only: r_kind - implicit none - class(qNode), intent(in):: aNode - integer(kind=i_kind),intent(in):: jiter - real(kind=r_kind),intent(inout):: tlddp - integer(kind=i_kind),optional,intent(inout):: nob - - tlddp = tlddp + aNode%diags%tldepart(jiter)*aNode%diags%tldepart(jiter) - if(present(nob)) nob=nob+1 -return -end subroutine gettlddp_ - -end module m_qNode diff --git a/src/gsi/m_qnode.F90 b/src/gsi/m_qnode.F90 new file mode 100644 index 0000000000..55ee986f54 --- /dev/null +++ b/src/gsi/m_qnode.F90 @@ -0,0 +1,282 @@ +module m_qnode +!$$$ subprogram documentation block +! . . . . +! subprogram: module m_qnode +! prgmmr: j guo +! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 +! date: 2016-05-18 +! +! abstract: class-module of obs-type qnode (moisture) +! +! program history log: +! 2016-05-18 j guo - added this document block for the initial polymorphic +! implementation. +! 2019-09-20 X.Su - add new variational QC parameters +! +! input argument list: see Fortran 90 style document below +! +! output argument list: see Fortran 90 style document below +! +! attributes: +! language: Fortran 90 and/or above +! machine: +! +!$$$ end subprogram documentation block + +! module interface: + use m_obsdiagnode, only: obs_diag + use m_obsdiagnode, only: obs_diags + use kinds , only: i_kind,r_kind + use mpeu_util, only: assert_,die,perr,warn,tell + use m_obsnode, only: obsnode + implicit none + private + + public:: qnode + + type,extends(obsnode):: qnode + !type(q_ob_type),pointer :: llpoint => null() + type(obs_diag), pointer :: diags => null() + real(r_kind) :: res ! moisture residual + real(r_kind) :: err2 ! moisture error squared + real(r_kind) :: raterr2 ! square of ratio of final obs error + ! to original obs error + !real(r_kind) :: time ! observation time in sec + real(r_kind) :: b ! variational quality control parameter + real(r_kind) :: pg ! variational quality control parameter + real(r_kind) :: jb ! variational quality control parameter + integer(i_kind) :: ib ! new variational quality control parameter + integer(i_kind) :: ik ! new variational quality control parameter + real(r_kind) :: wij(8) ! horizontal interpolation weights + real(r_kind) :: qpertb ! random number adding to the obs + integer(i_kind) :: ij(8) ! horizontal locations + integer(i_kind) :: k1 ! level of errtable 1-33 + integer(i_kind) :: kx ! ob type + !logical :: luse ! flag indicating if ob is used in pen. + + !integer(i_kind) :: idv,iob ! device id and obs index for sorting + !real (r_kind) :: elat, elon ! earth lat-lon for redistribution + !real (r_kind) :: dlat, dlon ! earth lat-lon for redistribution + real (r_kind) :: dlev ! reference to the vertical grid + + integer(i_kind) :: ich0=0 ! ich code to mark derived data. See + ! qnode_ich0 and qnode_ich0_pbl_Pseudo below + contains + procedure,nopass:: mytype + procedure:: sethop => obsnode_sethop_ + procedure:: xread => obsnode_xread_ + procedure:: xwrite => obsnode_xwrite_ + procedure:: isvalid => obsnode_isvalid_ + procedure:: gettlddp => gettlddp_ + + ! procedure, nopass:: headerread => obsheader_read_ + ! procedure, nopass:: headerwrite => obsheader_write_ + ! procedure:: init => obsnode_init_ + ! procedure:: clean => obsnode_clean_ + end type qnode + + public:: qnode_typecast + public:: qnode_nextcast + interface qnode_typecast; module procedure typecast_ ; end interface + interface qnode_nextcast; module procedure nextcast_ ; end interface + + public:: qnode_appendto + interface qnode_appendto; module procedure appendto_ ; end interface + + ! Because there are two components in qnode for an ordinary wind obs, + ! ich values are set to (1,2). Therefore, ich values for pbl_pseudo_surfobsuv + ! are set to (3,4), and qnode_ich0_pbl_pseudo is set to 2. + + public:: qnode_ich0 + public:: qnode_ich0_pbl_pseudo + + integer(i_kind),parameter :: qnode_ich0 = 0 ! ich=0+1 + integer(i_kind),parameter :: qnode_ich0_pbl_pseudo = qnode_ich0+1 ! ich=1+1 + + character(len=*),parameter:: myname="m_qnode" + +#include "myassert.H" +#include "mytrace.H" +contains +function typecast_(anode) result(ptr_) +!-- cast a class(obsnode) to a type(qnode) + use m_obsnode, only: obsnode + implicit none + type(qnode ),pointer:: ptr_ + class(obsnode),pointer,intent(in):: anode + ptr_ => null() + if(.not.associated(anode)) return + ! logically, typecast of a null-reference is a null pointer. + select type(anode) + type is(qnode) + ptr_ => anode + end select + return +end function typecast_ + +function nextcast_(anode) result(ptr_) +!-- cast an obsnode_next(obsnode) to a type(qnode) + use m_obsnode, only: obsnode,obsnode_next + implicit none + type(qnode ),pointer:: ptr_ + class(obsnode),target ,intent(in):: anode + + class(obsnode),pointer:: inode_ + inode_ => obsnode_next(anode) + ptr_ => typecast_(inode_) + return +end function nextcast_ + +subroutine appendto_(anode,oll) +!-- append anode to linked-list oll + use m_obsnode , only: obsnode + use m_obsllist, only: obsllist,obsllist_appendnode + implicit none + type(qnode),pointer,intent(in):: anode + type(obsllist),intent(inout):: oll + + class(obsnode),pointer:: inode_ + inode_ => anode + call obsllist_appendnode(oll,inode_) + inode_ => null() +end subroutine appendto_ + +! obsnode implementations + +function mytype() + implicit none + character(len=:),allocatable:: mytype + mytype="[qnode]" +end function mytype + +subroutine obsnode_xread_(anode,iunit,istat,diaglookup,skip) + use m_obsdiagnode, only: obsdiaglookup_locate + implicit none + class(qnode) , intent(inout):: anode + integer(i_kind) , intent(in ):: iunit + integer(i_kind) , intent( out):: istat + type(obs_diags) , intent(in ):: diaglookup + logical,optional, intent(in ):: skip + + character(len=*),parameter:: myname_=myname//'.obsnode_xread_' + logical:: skip_ +_ENTRY_(myname_) + + skip_=.false. + if(present(skip)) skip_=skip + + istat=0 + if(skip_) then + read(iunit,iostat=istat) + if(istat/=0) then + call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + endif + + else + read(iunit,iostat=istat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%jb , & + anode%ib , & + anode%ik , & + anode%qpertb , & + anode%k1 , & + anode%kx , & + anode%dlev , & + anode%ich0 , & + anode%wij , & + anode%ij + if (istat/=0) then + call perr(myname_,'read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + end if + + anode%diags => obsdiaglookup_locate(diaglookup,anode%idv,anode%iob,anode%ich0+1) + if(.not.associated(anode%diags)) then + call perr(myname_,'obsdiaglookup_locate(), %idv =',anode%idv) + call perr(myname_,' %iob =',anode%iob) + call perr(myname_,' %ich0 =',anode%ich0) + call die(myname_) + endif + endif +_EXIT_(myname_) + return +end subroutine obsnode_xread_ + +subroutine obsnode_xwrite_(anode,junit,jstat) + implicit none + class(qnode),intent(in):: anode + integer(i_kind),intent(in ):: junit + integer(i_kind),intent( out):: jstat + + character(len=*),parameter:: myname_=myname//'.obsnode_xwrite_' +_ENTRY_(myname_) + + jstat=0 + write(junit,iostat=jstat) anode%res , & + anode%err2 , & + anode%raterr2, & + anode%b , & + anode%pg , & + anode%jb , & + anode%ib , & + anode%ik , & + anode%qpertb , & + anode%k1 , & + anode%kx , & + anode%dlev , & + anode%ich0 , & + anode%wij , & + anode%ij + if (jstat/=0) then + call perr(myname_,'write(%(res,err2,...)), iostat =',jstat) + _EXIT_(myname_) + return + end if +_EXIT_(myname_) + return +end subroutine obsnode_xwrite_ + +subroutine obsnode_sethop_(anode) + use m_cvgridlookup, only: cvgridlookup_getiw + implicit none + class(qnode),intent(inout):: anode + + character(len=*),parameter:: myname_=myname//'::obsnode_sethop_' +_ENTRY_(myname_) + call cvgridlookup_getiw(anode%elat,anode%elon,anode%dlev,anode%ij,anode%wij) +_EXIT_(myname_) + return +end subroutine obsnode_sethop_ + +function obsnode_isvalid_(anode) result(isvalid_) + implicit none + logical:: isvalid_ + class(qnode),intent(in):: anode + + character(len=*),parameter:: myname_=myname//'::obsnode_isvalid_' +_ENTRY_(myname_) + isvalid_=associated(anode%diags) +_EXIT_(myname_) + return +end function obsnode_isvalid_ + +pure subroutine gettlddp_(anode,jiter,tlddp,nob) + use kinds, only: r_kind + implicit none + class(qnode), intent(in):: anode + integer(kind=i_kind),intent(in):: jiter + real(kind=r_kind),intent(inout):: tlddp + integer(kind=i_kind),optional,intent(inout):: nob + + tlddp = tlddp + anode%diags%tldepart(jiter)*anode%diags%tldepart(jiter) + if(present(nob)) nob=nob+1 + return +end subroutine gettlddp_ + +end module m_qnode diff --git a/src/gsi/m_radNode.F90 b/src/gsi/m_radNode.F90 deleted file mode 100644 index 33070e8382..0000000000 --- a/src/gsi/m_radNode.F90 +++ /dev/null @@ -1,453 +0,0 @@ -module m_radNode -!$$$ subprogram documentation block -! . . . . -! subprogram: module m_radNode -! prgmmr: j guo -! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 -! date: 2016-05-18 -! -! abstract: class-module of obs-type radNode (radiances) -! -! program history log: -! 2016-05-18 j guo - added this document block for the initial polymorphic -! implementation. -! 2016-07-19 kbathmann - add rsqrtinv and use_corr_obs to rad_ob_type -! 2019-04-22 kbathmann - change rsqrtinv to Rpred -! -! input argument list: see Fortran 90 style document below -! -! output argument list: see Fortran 90 style document below -! -! attributes: -! language: Fortran 90 and/or above -! machine: -! -!$$$ end subprogram documentation block - -! module interface: - use m_obsdiagNode, only: obs_diag,aofp_obs_diag => fptr_obsdiagNode - use m_obsdiagNode, only: obs_diags - use kinds , only: i_kind,r_kind - use mpeu_util, only: assert_,die,perr,warn,tell - use m_obsNode, only: obsNode - implicit none - private - - public:: radNode - - type,extends(obsNode):: radNode - type(aofp_obs_diag), dimension(:), pointer :: diags => NULL() - real(r_kind),dimension(:),pointer :: res => NULL() - ! obs-guess residual (nchan) - real(r_kind),dimension(:),pointer :: err2 => NULL() - ! error variances squared (nchan) - real(r_kind),dimension(:),pointer :: raterr2 => NULL() - ! ratio of error variances squared (nchan) - real(r_kind) :: wij(4) ! horizontal interpolation weights - real(r_kind),dimension(:,:),pointer :: pred => NULL() - ! predictors (npred,nchan) - real(r_kind),dimension(:,:),pointer :: dtb_dvar => NULL() - ! radiance jacobian (nsigradjac,nchan) - - real(r_kind),dimension(:,:),pointer :: Rpred => NULL() - ! square root of inverse of R, multiplied - ! by bias predictor jacobian - ! only used if using correlated obs - real(r_kind),dimension(: ),pointer :: rsqrtinv => NULL() - ! square root of inverse of R, only used - ! if using correlated obs - - integer(i_kind),dimension(:),pointer :: icx => NULL() - integer(i_kind),dimension(:),pointer :: ich => NULL() - integer(i_kind) :: nchan ! number of channels for this profile - integer(i_kind) :: ij(4) ! horizontal locations - logical :: use_corr_obs = .false. ! to indicate if correlated obs is implemented - integer(i_kind) :: iuse_PredOper_type = 0 ! indicate which type of correlated predictor operator is implemented - ! 0 uses none (diagonal) - ! 1 uses rsqrtinv - ! 2 uses Rpred - -!!! Is %isis or %isfctype ever being assigned somewhere in the code? -!!! They are used in intrad(). -!!! -!!! Now, they are not written to an obsdiags file, nor read from one. - - character(20) :: isis ! sensor/instrument/satellite id, e.g. amsua_n15 - !integer(i_kind) :: isfctype ! surf mask: ocean=0,land=1,ice=2,snow=3,mixed=4 - character(80) :: covtype ! surf mask: ocean=0,land=1,ice=2,snow=3,mixed=4 - contains - procedure,nopass:: mytype - procedure:: setHop => obsNode_setHop_ - procedure:: xread => obsNode_xread_ - procedure:: xwrite => obsNode_xwrite_ - procedure:: isvalid => obsNode_isvalid_ - procedure:: gettlddp => gettlddp_ - - procedure, nopass:: headerRead => obsHeader_read_ - procedure, nopass:: headerWrite => obsHeader_write_ - ! procedure:: init => obsNode_init_ - procedure:: clean => obsNode_clean_ - end type radNode - - public:: radNode_typecast - public:: radNode_nextcast - interface radNode_typecast; module procedure typecast_ ; end interface - interface radNode_nextcast; module procedure nextcast_ ; end interface - - public:: radNode_appendto - interface radNode_appendto; module procedure appendto_ ; end interface - - character(len=*),parameter:: MYNAME="m_radNode" - -#include "myassert.H" -#include "mytrace.H" -contains -function typecast_(aNode) result(ptr_) -!-- cast a class(obsNode) to a type(radNode) - use m_obsNode, only: obsNode - implicit none - type(radNode ),pointer:: ptr_ - class(obsNode),pointer,intent(in):: aNode - ptr_ => null() - if(.not.associated(aNode)) return - ! logically, typecast of a null-reference is a null pointer. - select type(aNode) - type is(radNode) - ptr_ => aNode - end select -return -end function typecast_ - -function nextcast_(aNode) result(ptr_) -!-- cast an obsNode_next(obsNode) to a type(radNode) - use m_obsNode, only: obsNode,obsNode_next - implicit none - type(radNode ),pointer:: ptr_ - class(obsNode),target ,intent(in):: aNode - - class(obsNode),pointer:: inode_ - inode_ => obsNode_next(aNode) - ptr_ => typecast_(inode_) -return -end function nextcast_ - -subroutine appendto_(aNode,oll) -!-- append aNode to linked-list oLL - use m_obsNode , only: obsNode - use m_obsLList, only: obsLList,obsLList_appendNode - implicit none - type(radNode),pointer,intent(in):: aNode - type(obsLList),intent(inout):: oLL - - class(obsNode),pointer:: inode_ - inode_ => aNode - call obsLList_appendNode(oLL,inode_) - inode_ => null() -end subroutine appendto_ - -! obsNode implementations - -function mytype() - implicit none - character(len=:),allocatable:: mytype - mytype="[radNode]" -end function mytype - -subroutine obsHeader_read_(iunit,mobs,jread,istat) - use radinfo, only: npred,nsigradjac - implicit none - integer(i_kind),intent(in ):: iunit - integer(i_kind),intent(out):: mobs - integer(i_kind),intent(out):: jread - integer(i_kind),intent(out):: istat - - character(len=*),parameter:: myname_=myname//'.obsHeader_read_' - integer(i_kind):: mpred,msigradjac -_ENTRY_(myname_) - - read(iunit,iostat=istat) mobs,jread, mpred,msigradjac - if(istat==0 .and. (npred/=mpred .or. nsigradjac/=msigradjac)) then - call perr(myname_,'unmatched dimension information, npred or nsigradjac') - if(npred/=mpred) then - call perr(myname_,' expecting npred =',npred) - call perr(myname_,' but read mpred =',mpred) - endif - if(nsigradjac/=msigradjac) then - call perr(myname_,'expecting nsigradjac =',nsigradjac) - call perr(myname_,' but read msigradjac =',msigradjac) - endif - call die(myname_) - endif -_EXIT_(myname_) -return -end subroutine obsHeader_read_ - -subroutine obsHeader_write_(junit,mobs,jwrite,jstat) - use radinfo, only: npred,nsigradjac - implicit none - integer(i_kind),intent(in ):: junit - integer(i_kind),intent(in ):: mobs - integer(i_kind),intent(in ):: jwrite - integer(i_kind),intent(out):: jstat - - character(len=*),parameter:: myname_=myname//'.obsHeader_write_' -_ENTRY_(myname_) - write(junit,iostat=jstat) mobs,jwrite, npred,nsigradjac -_EXIT_(myname_) -return -end subroutine obsHeader_write_ - -subroutine obsNode_clean_(aNode) - implicit none - class(radNode),intent(inout):: aNode - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_clean_' -_ENTRY_(myname_) -!_TRACEV_(myname_,'%mytype() =',aNode%mytype()) - if(associated(aNode%diags )) deallocate(aNode%diags ) - if(associated(aNode%ich )) deallocate(aNode%ich ) - if(associated(aNode%res )) deallocate(aNode%res ) - if(associated(aNode%err2 )) deallocate(aNode%err2 ) - if(associated(aNode%raterr2 )) deallocate(aNode%raterr2 ) - if(associated(aNode%pred )) deallocate(aNode%pred ) - if(associated(aNode%dtb_dvar)) deallocate(aNode%dtb_dvar) - if(associated(aNode%Rpred )) deallocate(aNode%Rpred ) - if(associated(aNode%rsqrtinv)) deallocate(aNode%rsqrtinv) - if(associated(aNode%icx )) deallocate(aNode%icx ) -_EXIT_(myname_) -return -end subroutine obsNode_clean_ - -subroutine obsNode_xread_(aNode,iunit,istat,diagLookup,skip) - use m_obsdiagNode, only: obsdiagLookup_locate - use radinfo, only: npred,nsigradjac - implicit none - class(radNode),intent(inout):: aNode - integer(i_kind),intent(in ):: iunit - integer(i_kind),intent( out):: istat - type(obs_diags),intent(in ):: diagLookup - logical,optional,intent(in ):: skip - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xread_' - integer(i_kind):: k,nchan - logical:: skip_ -_ENTRY_(myname_) - skip_=.false. - if(present(skip)) skip_=skip - - istat=0 - if(skip_) then - read(iunit,iostat=istat) - if (istat/=0) then - call perr(myname_,'skipping read(%(nchan,iuse_PredOper_type)), iostat =',istat) - _EXIT_(myname_) - return - end if - - read(iunit,iostat=istat) - if(istat/=0) then - call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - endif - - read(iunit,iostat=istat) - if(istat/=0) then - call perr(myname_,'skipping read(%(Rpred||rsqrtinv)), iostat =',istat) - _EXIT_(myname_) - return - endif - - else - read(iunit,iostat=istat) aNode%nchan,aNode%use_corr_obs,aNode%iuse_PredOper_type - if (istat/=0) then - call perr(myname_,'read(%(nchan,use_corr_obs,iuse_PredOper_type)), iostat =',istat) - _EXIT_(myname_) - return - end if - - if(associated(aNode%diags )) deallocate(aNode%diags ) - if(associated(aNode%ich )) deallocate(aNode%ich ) - if(associated(aNode%res )) deallocate(aNode%res ) - if(associated(aNode%err2 )) deallocate(aNode%err2 ) - if(associated(aNode%raterr2 )) deallocate(aNode%raterr2 ) - if(associated(aNode%pred )) deallocate(aNode%pred ) - if(associated(aNode%dtb_dvar)) deallocate(aNode%dtb_dvar) - if(associated(aNode%Rpred )) deallocate(aNode%Rpred) - if(associated(aNode%rsqrtinv)) deallocate(aNode%rsqrtinv) - if(associated(aNode%icx )) deallocate(aNode%icx ) - - nchan=aNode%nchan - allocate( aNode%diags(nchan), & - aNode%res (nchan), & - aNode%err2 (nchan), & - aNode%raterr2 (nchan), & - aNode%pred (npred,nchan), & - aNode%dtb_dvar(nsigradjac,nchan), & - aNode%ich (nchan), & - aNode%icx (nchan) ) - - read(iunit,iostat=istat) aNode%ich , & - aNode%res , & - aNode%err2 , & - aNode%raterr2 , & - aNode%pred , & - aNode%icx , & - aNode%dtb_dvar, & - aNode%wij , & - aNode%ij - if (istat/=0) then - call perr(myname_,'read(%(res,err2,...)), iostat =',istat) - _EXIT_(myname_) - return - end if - - if(.not.aNode%use_corr_obs) aNode%iuse_PredOper_type=0 - select case(aNode%iuse_PredOper_type) - case(1) - allocate(aNode%rsqrtinv(((nchan+1)*nchan)/2)) - read(iunit,iostat=istat) aNode%rsqrtinv - if (istat/=0) then - call perr(myname_,'read(%rsqrtinv), iostat =',istat) - _EXIT_(myname_) - return - end if - case(2) - allocate(aNode%Rpred(((nchan+1)*nchan)/2,npred)) - read(iunit,iostat=istat) aNode%Rpred - if (istat/=0) then - call perr(myname_,'read(%Rpred), iostat =',istat) - _EXIT_(myname_) - return - end if - - case default - read(iunit,iostat=istat) - end select - - do k=1,nchan - aNode%diags(k)%ptr => obsdiagLookup_locate(diagLookup,aNode%idv,aNode%iob,aNode%ich(k)) - if(.not.associated(aNode%diags(k)%ptr)) then - call perr(myname_,'obsdiagLookup_locate(k), k =',k) - call perr(myname_,' %idv =',aNode%idv) - call perr(myname_,' %iob =',aNode%iob) - call perr(myname_,' %ich(k) =',aNode%ich(k)) - call die(myname_) - endif - enddo - endif -_EXIT_(myname_) -return -end subroutine obsNode_xread_ - -subroutine obsNode_xwrite_(aNode,junit,jstat) - use radinfo, only: npred - implicit none - class(radNode),intent(in):: aNode - integer(i_kind),intent(in ):: junit - integer(i_kind),intent( out):: jstat - - character(len=*),parameter:: myname_=MYNAME//'.obsNode_xwrite_' - integer(i_kind):: k - integer(i_kind):: iuse_PredOper_type -_ENTRY_(myname_) - - jstat=0 - iuse_PredOper_type=0 - if(aNode%use_corr_obs) iuse_PredOper_type=aNode%iuse_PredOper_type - write(junit,iostat=jstat) aNode%nchan,aNode%use_corr_obs,iuse_PredOper_type - if (jstat/=0) then - call perr(myname_,'write(%(nchan,use_corr_obs, etc.)), iostat =',jstat) - _EXIT_(myname_) - return - end if - - write(junit,iostat=jstat) (/ (aNode%ich(k),k=1,aNode%nchan) /), & - aNode%res , & - aNode%err2 , & - aNode%raterr2 , & - aNode%pred , & - aNode%icx , & - aNode%dtb_dvar, & - aNode%wij , & - aNode%ij - if (jstat/=0) then - call perr(myname_,'write(%(ich,res,err2,...)), iostat =',jstat) - _EXIT_(myname_) - return - end if - - select case(iuse_PredOper_type) - case(1) - ASSERT(size(aNode%rsqrtinv)==((aNode%nchan+1)*aNode%nchan)/2) - write(junit,iostat=jstat) aNode%rsqrtinv - if (jstat/=0) then - call perr(myname_,'write(%rsqrtinv), iostat =',jstat) - _EXIT_(myname_) - return - end if - case(2) - ASSERT(size(aNode%Rpred,1)==((aNode%nchan+1)*aNode%nchan)/2) - ASSERT(size(aNode%Rpred,2)==npred) - write(junit,iostat=jstat) aNode%Rpred - if (jstat/=0) then - call perr(myname_,'write(%Rpred), iostat =',jstat) - _EXIT_(myname_) - return - end if - - case default - write(junit,iostat=jstat) - if (jstat/=0) then - call perr(myname_,'write as skip record, iostat =',jstat) - _EXIT_(myname_) - return - end if - end select - -_EXIT_(myname_) -return -end subroutine obsNode_xwrite_ - -subroutine obsNode_setHop_(aNode) - use m_cvgridLookup, only: cvgridLookup_getiw - implicit none - class(radNode),intent(inout):: aNode - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_setHop_' -_ENTRY_(myname_) - call cvgridLookup_getiw(aNode%elat,aNode%elon,aNode%ij,aNode%wij) -_EXIT_(myname_) -return -end subroutine obsNode_setHop_ - -function obsNode_isvalid_(aNode) result(isvalid_) - implicit none - logical:: isvalid_ - class(radNode),intent(in):: aNode - - character(len=*),parameter:: myname_=MYNAME//'::obsNode_isvalid_' - integer(i_kind):: k -_ENTRY_(myname_) - isvalid_=all( (/ (associated(aNode%diags(k)%ptr),k=1,aNode%nchan) /) ) -_EXIT_(myname_) -return -end function obsNode_isvalid_ - -pure subroutine gettlddp_(aNode,jiter,tlddp,nob) - use kinds, only: r_kind - implicit none - class(radNode), intent(in):: aNode - integer(kind=i_kind),intent(in):: jiter - real(kind=r_kind),intent(inout):: tlddp - integer(kind=i_kind),optional,intent(inout):: nob - - integer(kind=i_kind):: k - do k=1,aNode%nchan - tlddp = tlddp + aNode%diags(k)%ptr%tldepart(jiter)*aNode%diags(k)%ptr%tldepart(jiter) - enddo - if(present(nob)) nob=nob+aNode%nchan -return -end subroutine gettlddp_ - -end module m_radNode diff --git a/src/gsi/m_radnode.F90 b/src/gsi/m_radnode.F90 new file mode 100644 index 0000000000..bf473873ca --- /dev/null +++ b/src/gsi/m_radnode.F90 @@ -0,0 +1,453 @@ +module m_radnode +!$$$ subprogram documentation block +! . . . . +! subprogram: module m_radnode +! prgmmr: j guo +! org: NASA/GSFC, Global Modeling and Assimilation Office, 610.3 +! date: 2016-05-18 +! +! abstract: class-module of obs-type radnode (radiances) +! +! program history log: +! 2016-05-18 j guo - added this document block for the initial polymorphic +! implementation. +! 2016-07-19 kbathmann - add rsqrtinv and use_corr_obs to rad_ob_type +! 2019-04-22 kbathmann - change rsqrtinv to rpred +! +! input argument list: see Fortran 90 style document below +! +! output argument list: see Fortran 90 style document below +! +! attributes: +! language: Fortran 90 and/or above +! machine: +! +!$$$ end subprogram documentation block + +! module interface: + use m_obsdiagnode, only: obs_diag,aofp_obs_diag => fptr_obsdiagnode + use m_obsdiagnode, only: obs_diags + use kinds , only: i_kind,r_kind + use mpeu_util, only: assert_,die,perr,warn,tell + use m_obsnode, only: obsnode + implicit none + private + + public:: radnode + + type,extends(obsnode):: radnode + type(aofp_obs_diag), dimension(:), pointer :: diags => null() + real(r_kind),dimension(:),pointer :: res => null() + ! obs-guess residual (nchan) + real(r_kind),dimension(:),pointer :: err2 => null() + ! error variances squared (nchan) + real(r_kind),dimension(:),pointer :: raterr2 => null() + ! ratio of error variances squared (nchan) + real(r_kind) :: wij(4) ! horizontal interpolation weights + real(r_kind),dimension(:,:),pointer :: pred => null() + ! predictors (npred,nchan) + real(r_kind),dimension(:,:),pointer :: dtb_dvar => null() + ! radiance jacobian (nsigradjac,nchan) + + real(r_kind),dimension(:,:),pointer :: rpred => null() + ! square root of inverse of Rr multiplied + ! by bias predictor jacobian + ! only used if using correlated obs + real(r_kind),dimension(: ),pointer :: rsqrtinv => null() + ! square root of inverse of r, only used + ! if using correlated obs + + integer(i_kind),dimension(:),pointer :: icx => null() + integer(i_kind),dimension(:),pointer :: ich => null() + integer(i_kind) :: nchan ! number of channels for this profile + integer(i_kind) :: ij(4) ! horizontal locations + logical :: use_corr_obs = .false. ! to indicate if correlated obs is implemented + integer(i_kind) :: iuse_predoper_type = 0 ! indicate which type of correlated predictor operator is implemented + ! 0 uses none (diagonal) + ! 1 uses rsqrtinv + ! 2 uses rpred + +!!! Is %isis or %isfctype ever being assigned somewhere in the code? +!!! They are used in intrad(). +!!! +!!! Now, they are not written to an obsdiags file, nor read from one. + + character(20) :: isis ! sensor/instrument/satellite id, e.g. amsua_n15 + !integer(i_kind) :: isfctype ! surf mask: ocean=0,land=1,ice=2,snow=3,mixed=4 + character(80) :: covtype ! surf mask: ocean=0,land=1,ice=2,snow=3,mixed=4 + contains + procedure,nopass:: mytype + procedure:: sethop => obsnode_sethop_ + procedure:: xread => obsnode_xread_ + procedure:: xwrite => obsnode_xwrite_ + procedure:: isvalid => obsnode_isvalid_ + procedure:: gettlddp => gettlddp_ + + procedure, nopass:: headerread => obsheader_read_ + procedure, nopass:: headerwrite => obsheader_write_ + ! procedure:: init => obsnode_init_ + procedure:: clean => obsnode_clean_ + end type radnode + + public:: radnode_typecast + public:: radnode_nextcast + interface radnode_typecast; module procedure typecast_ ; end interface + interface radnode_nextcast; module procedure nextcast_ ; end interface + + public:: radnode_appendto + interface radnode_appendto; module procedure appendto_ ; end interface + + character(len=*),parameter:: myname="m_radnode" + +#include "myassert.H" +#include "mytrace.H" +contains +function typecast_(anode) result(ptr_) +!-- cast a class(obsnode) to a type(radnode) + use m_obsnode, only: obsnode + implicit none + type(radnode ),pointer:: ptr_ + class(obsnode),pointer,intent(in):: anode + ptr_ => null() + if(.not.associated(anode)) return + ! logically, typecast of a null-reference is a null pointer. + select type(anode) + type is(radnode) + ptr_ => anode + end select + return +end function typecast_ + +function nextcast_(anode) result(ptr_) +!-- cast an obsnode_next(obsnode) to a type(radnode) + use m_obsnode, only: obsnode,obsnode_next + implicit none + type(radnode ),pointer:: ptr_ + class(obsnode),target ,intent(in):: anode + + class(obsnode),pointer:: inode_ + inode_ => obsnode_next(anode) + ptr_ => typecast_(inode_) + return +end function nextcast_ + +subroutine appendto_(anode,oll) +!-- append anode to linked-list oll + use m_obsnode , only: obsnode + use m_obsllist, only: obsllist,obsllist_appendnode + implicit none + type(radnode),pointer,intent(in):: anode + type(obsllist),intent(inout):: oll + + class(obsnode),pointer:: inode_ + inode_ => anode + call obsllist_appendnode(oll,inode_) + inode_ => null() +end subroutine appendto_ + +! obsnode implementations + +function mytype() + implicit none + character(len=:),allocatable:: mytype + mytype="[radnode]" +end function mytype + +subroutine obsheader_read_(iunit,mobs,jread,istat) + use radinfo, only: npred,nsigradjac + implicit none + integer(i_kind),intent(in ):: iunit + integer(i_kind),intent(out):: mobs + integer(i_kind),intent(out):: jread + integer(i_kind),intent(out):: istat + + character(len=*),parameter:: myname_=myname//'.obsheader_read_' + integer(i_kind):: mpred,msigradjac +_ENTRY_(myname_) + + read(iunit,iostat=istat) mobs,jread, mpred,msigradjac + if(istat==0 .and. (npred/=mpred .or. nsigradjac/=msigradjac)) then + call perr(myname_,'unmatched dimension information, npred or nsigradjac') + if(npred/=mpred) then + call perr(myname_,' expecting npred =',npred) + call perr(myname_,' but read mpred =',mpred) + endif + if(nsigradjac/=msigradjac) then + call perr(myname_,'expecting nsigradjac =',nsigradjac) + call perr(myname_,' but read msigradjac =',msigradjac) + endif + call die(myname_) + endif +_EXIT_(myname_) + return +end subroutine obsheader_read_ + +subroutine obsheader_write_(junit,mobs,jwrite,jstat) + use radinfo, only: npred,nsigradjac + implicit none + integer(i_kind),intent(in ):: junit + integer(i_kind),intent(in ):: mobs + integer(i_kind),intent(in ):: jwrite + integer(i_kind),intent(out):: jstat + + character(len=*),parameter:: myname_=myname//'.obsheader_write_' +_ENTRY_(myname_) + write(junit,iostat=jstat) mobs,jwrite, npred,nsigradjac +_EXIT_(myname_) + return +end subroutine obsheader_write_ + +subroutine obsnode_clean_(anode) + implicit none + class(radnode),intent(inout):: anode + + character(len=*),parameter:: myname_=myname//'.obsnode_clean_' +_ENTRY_(myname_) +!_TRACEV_(myname_,'%mytype() =',anode%mytype()) + if(associated(anode%diags )) deallocate(anode%diags ) + if(associated(anode%ich )) deallocate(anode%ich ) + if(associated(anode%res )) deallocate(anode%res ) + if(associated(anode%err2 )) deallocate(anode%err2 ) + if(associated(anode%raterr2 )) deallocate(anode%raterr2 ) + if(associated(anode%pred )) deallocate(anode%pred ) + if(associated(anode%dtb_dvar)) deallocate(anode%dtb_dvar) + if(associated(anode%rpred )) deallocate(anode%rpred ) + if(associated(anode%rsqrtinv)) deallocate(anode%rsqrtinv) + if(associated(anode%icx )) deallocate(anode%icx ) +_EXIT_(myname_) + return +end subroutine obsnode_clean_ + +subroutine obsnode_xread_(anode,iunit,istat,diaglookup,skip) + use m_obsdiagnode, only: obsdiaglookup_locate + use radinfo, only: npred,nsigradjac + implicit none + class(radnode),intent(inout):: anode + integer(i_kind),intent(in ):: iunit + integer(i_kind),intent( out):: istat + type(obs_diags),intent(in ):: diaglookup + logical,optional,intent(in ):: skip + + character(len=*),parameter:: myname_=myname//'.obsnode_xread_' + integer(i_kind):: k,nchan + logical:: skip_ +_ENTRY_(myname_) + skip_=.false. + if(present(skip)) skip_=skip + + istat=0 + if(skip_) then + read(iunit,iostat=istat) + if (istat/=0) then + call perr(myname_,'skipping read(%(nchan,iuse_predoper_type)), iostat =',istat) + _EXIT_(myname_) + return + end if + + read(iunit,iostat=istat) + if(istat/=0) then + call perr(myname_,'skipping read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + endif + + read(iunit,iostat=istat) + if(istat/=0) then + call perr(myname_,'skipping read(%(rpred||rsqrtinv)), iostat =',istat) + _EXIT_(myname_) + return + endif + + else + read(iunit,iostat=istat) anode%nchan,anode%use_corr_obs,anode%iuse_predoper_type + if (istat/=0) then + call perr(myname_,'read(%(nchan,use_corr_obs,iuse_predoper_type)), iostat =',istat) + _EXIT_(myname_) + return + end if + + if(associated(anode%diags )) deallocate(anode%diags ) + if(associated(anode%ich )) deallocate(anode%ich ) + if(associated(anode%res )) deallocate(anode%res ) + if(associated(anode%err2 )) deallocate(anode%err2 ) + if(associated(anode%raterr2 )) deallocate(anode%raterr2 ) + if(associated(anode%pred )) deallocate(anode%pred ) + if(associated(anode%dtb_dvar)) deallocate(anode%dtb_dvar) + if(associated(anode%rpred )) deallocate(anode%rpred) + if(associated(anode%rsqrtinv)) deallocate(anode%rsqrtinv) + if(associated(anode%icx )) deallocate(anode%icx ) + + nchan=anode%nchan + allocate( anode%diags(nchan), & + anode%res (nchan), & + anode%err2 (nchan), & + anode%raterr2 (nchan), & + anode%pred (npred,nchan), & + anode%dtb_dvar(nsigradjac,nchan), & + anode%ich (nchan), & + anode%icx (nchan) ) + + read(iunit,iostat=istat) anode%ich , & + anode%res , & + anode%err2 , & + anode%raterr2 , & + anode%pred , & + anode%icx , & + anode%dtb_dvar, & + anode%wij , & + anode%ij + if (istat/=0) then + call perr(myname_,'read(%(res,err2,...)), iostat =',istat) + _EXIT_(myname_) + return + end if + + if(.not.anode%use_corr_obs) anode%iuse_predoper_type=0 + select case(anode%iuse_predoper_type) + case(1) + allocate(anode%rsqrtinv(((nchan+1)*nchan)/2)) + read(iunit,iostat=istat) anode%rsqrtinv + if (istat/=0) then + call perr(myname_,'read(%rsqrtinv), iostat =',istat) + _EXIT_(myname_) + return + end if + case(2) + allocate(anode%rpred(((nchan+1)*nchan)/2,npred)) + read(iunit,iostat=istat) anode%rpred + if (istat/=0) then + call perr(myname_,'read(%rpred), iostat =',istat) + _EXIT_(myname_) + return + end if + + case default + read(iunit,iostat=istat) + end select + + do k=1,nchan + anode%diags(k)%ptr => obsdiaglookup_locate(diaglookup,anode%idv,anode%iob,anode%ich(k)) + if(.not.associated(anode%diags(k)%ptr)) then + call perr(myname_,'obsdiaglookup_locate(k), k =',k) + call perr(myname_,' %idv =',anode%idv) + call perr(myname_,' %iob =',anode%iob) + call perr(myname_,' %ich(k) =',anode%ich(k)) + call die(myname_) + endif + enddo + endif +_EXIT_(myname_) + return +end subroutine obsnode_xread_ + +subroutine obsnode_xwrite_(anode,junit,jstat) + use radinfo, only: npred + implicit none + class(radnode),intent(in):: anode + integer(i_kind),intent(in ):: junit + integer(i_kind),intent( out):: jstat + + character(len=*),parameter:: myname_=myname//'.obsnode_xwrite_' + integer(i_kind):: k + integer(i_kind):: iuse_predoper_type +_ENTRY_(myname_) + + jstat=0 + iuse_predoper_type=0 + if(anode%use_corr_obs) iuse_predoper_type=anode%iuse_predoper_type + write(junit,iostat=jstat) anode%nchan,anode%use_corr_obs,iuse_predoper_type + if (jstat/=0) then + call perr(myname_,'write(%(nchan,use_corr_obs, etc.)), iostat =',jstat) + _EXIT_(myname_) + return + end if + + write(junit,iostat=jstat) (/ (anode%ich(k),k=1,anode%nchan) /), & + anode%res , & + anode%err2 , & + anode%raterr2 , & + anode%pred , & + anode%icx , & + anode%dtb_dvar, & + anode%wij , & + anode%ij + if (jstat/=0) then + call perr(myname_,'write(%(ich,res,err2,...)), iostat =',jstat) + _EXIT_(myname_) + return + end if + + select case(iuse_predoper_type) + case(1) + ASSERT(size(anode%rsqrtinv)==((anode%nchan+1)*anode%nchan)/2) + write(junit,iostat=jstat) anode%rsqrtinv + if (jstat/=0) then + call perr(myname_,'write(%rsqrtinv), iostat =',jstat) + _EXIT_(myname_) + return + end if + case(2) + ASSERT(size(anode%rpred,1)==((anode%nchan+1)*anode%nchan)/2) + ASSERT(size(anode%rpred,2)==npred) + write(junit,iostat=jstat) anode%rpred + if (jstat/=0) then + call perr(myname_,'write(%rpred), iostat =',jstat) + _EXIT_(myname_) + return + end if + + case default + write(junit,iostat=jstat) + if (jstat/=0) then + call perr(myname_,'write as skip record, iostat =',jstat) + _EXIT_(myname_) + return + end if + end select + +_EXIT_(myname_) + return +end subroutine obsnode_xwrite_ + +subroutine obsnode_sethop_(anode) + use m_cvgridlookup, only: cvgridlookup_getiw + implicit none + class(radnode),intent(inout):: anode + + character(len=*),parameter:: myname_=myname//'::obsnode_sethop_' +_ENTRY_(myname_) + call cvgridlookup_getiw(anode%elat,anode%elon,anode%ij,anode%wij) +_EXIT_(myname_) + return +end subroutine obsnode_sethop_ + +function obsnode_isvalid_(anode) result(isvalid_) + implicit none + logical:: isvalid_ + class(radnode),intent(in):: anode + + character(len=*),parameter:: myname_=myname//'::obsnode_isvalid_' + integer(i_kind):: k +_ENTRY_(myname_) + isvalid_=all( (/ (associated(anode%diags(k)%ptr),k=1,anode%nchan) /) ) +_EXIT_(myname_) + return +end function obsnode_isvalid_ + +pure subroutine gettlddp_(anode,jiter,tlddp,nob) + use kinds, only: r_kind + implicit none + class(radnode), intent(in):: anode + integer(kind=i_kind),intent(in):: jiter + real(kind=r_kind),intent(inout):: tlddp + integer(kind=i_kind),optional,intent(inout):: nob + + integer(kind=i_kind):: k + do k=1,anode%nchan + tlddp = tlddp + anode%diags(k)%ptr%tldepart(jiter)*anode%diags(k)%ptr%tldepart(jiter) + enddo + if(present(nob)) nob=nob+anode%nchan + return +end subroutine gettlddp_ + +end module m_radnode diff --git a/src/gsi/m_rerank.f90 b/src/gsi/m_rerank.f90 index fd905c570f..ee5573654d 100644 --- a/src/gsi/m_rerank.f90 +++ b/src/gsi/m_rerank.f90 @@ -1,5 +1,5 @@ module m_rerank -! From Jing Guo - use kinds for GSI consitency (Todling) +! From Jing Guo - use kinds for gsi consistency (todling) ! 01Aug2011 - Lueken - changed F90 to f90 (no machine logic ! 28Apr2011 - Todling - overload to handle single precision use kinds, only : r_double,r_single,i_kind @@ -24,7 +24,7 @@ module m_rerank end interface rerank character(len=*),parameter :: myname='m_rerank' -CONTAINS +contains function rerank_2in1r4_(i2) result(i1) implicit none @@ -267,15 +267,15 @@ end function rerank_hack_1in4r8_ end function rerank_1in4r8_ !>>> subroutine assert_eq_(lsize,lrank,who,what) - implicit none - integer(i_kind),intent(in)::lsize,lrank - character(len=*),intent(in)::who,what - if(lsize==lrank) return - write(*,*)' lsize= ',lsize - write(*,*)' lrank= ',lrank - write(*,*)' whois= ',who - write(*,*)' whats= ',what - call exit(2) + implicit none + integer(i_kind),intent(in)::lsize,lrank + character(len=*),intent(in)::who,what + if(lsize==lrank) return + write(*,*)' lsize= ',lsize + write(*,*)' lrank= ',lrank + write(*,*)' whois= ',who + write(*,*)' whats= ',what + call exit(2) end subroutine assert_eq_ end module m_rerank diff --git a/src/gsi/mpeu_mpif.F90 b/src/gsi/mpeu_mpif.F90 index 2888489858..305aac1596 100644 --- a/src/gsi/mpeu_mpif.F90 +++ b/src/gsi/mpeu_mpif.F90 @@ -6,7 +6,7 @@ module mpeu_mpif ! org: NASA/GSFC, Global Modeling and Assimilation Office, 900.3 ! date: 2010-03-22 ! -! abstract: a portable interface to include "mpif.h" for MPI. +! abstract: a portable interface to include "mpif.h" for mpi. ! ! program history log: ! 2010-03-22 j guo - added this document block @@ -26,88 +26,88 @@ module mpeu_mpif implicit none private ! except - public :: MPI_type - interface MPI_type; module procedure & - type_r0_of_integer1, & - type_r0_of_integer2, & - type_r0_of_integer4, & - type_r0_of_integer8, & - type_r0_of_real4 , & - type_r0_of_real8 , & - type_r0_of_real16 , & - type_r1_of_integer1, & - type_r1_of_integer2, & - type_r1_of_integer4, & - type_r1_of_integer8, & - type_r1_of_real4 , & - type_r1_of_real8 , & - type_r1_of_real16 , & - type_r2_of_integer1, & - type_r2_of_integer2, & - type_r2_of_integer4, & - type_r2_of_integer8, & - type_r2_of_real4 , & - type_r2_of_real8 , & - type_r2_of_real16 , & - type_r3_of_integer1, & - type_r3_of_integer2, & - type_r3_of_integer4, & - type_r3_of_integer8, & - type_r3_of_real4 , & - type_r3_of_real8 , & - type_r3_of_real16 + public :: mpi_type + interface mpi_type; module procedure & + type_r0_of_integer1, & + type_r0_of_integer2, & + type_r0_of_integer4, & + type_r0_of_integer8, & + type_r0_of_real4 , & + type_r0_of_real8 , & + type_r0_of_real16 , & + type_r1_of_integer1, & + type_r1_of_integer2, & + type_r1_of_integer4, & + type_r1_of_integer8, & + type_r1_of_real4 , & + type_r1_of_real8 , & + type_r1_of_real16 , & + type_r2_of_integer1, & + type_r2_of_integer2, & + type_r2_of_integer4, & + type_r2_of_integer8, & + type_r2_of_real4 , & + type_r2_of_real8 , & + type_r2_of_real16 , & + type_r3_of_integer1, & + type_r3_of_integer2, & + type_r3_of_integer4, & + type_r3_of_integer8, & + type_r3_of_real4 , & + type_r3_of_real8 , & + type_r3_of_real16 end interface - public :: MPI_IKIND - public :: MPI_ADDRESS_KIND + public :: mpi_ikind + public :: mpi_address_kind - public :: MPI_INTEGER1 - public :: MPI_INTEGER2 - public :: MPI_INTEGER4 - public :: MPI_INTEGER8 + public :: mpi_integer1 + public :: mpi_integer2 + public :: mpi_integer4 + public :: mpi_integer8 - public :: MPI_INTEGER - public :: MPI_REAL - public :: MPI_DOUBLE_PRECISION - public :: MPI_LOGICAL - public :: MPI_CHARACTER + public :: mpi_integer + public :: mpi_real + public :: mpi_double_precision + public :: mpi_logical + public :: mpi_character - public :: MPI_2INTEGER - public :: MPI_2REAL - public :: MPI_2DOUBLE_PRECISION + public :: mpi_2integer + public :: mpi_2real + public :: mpi_2double_precision - public :: MPI_REAL4 - public :: MPI_REAL8 - public :: MPI_REAL16 + public :: mpi_real4 + public :: mpi_real8 + public :: mpi_real16 - public :: MPI_COMM_WORLD - public :: MPI_COMM_NULL - public :: MPI_REQUEST_NULL + public :: mpi_comm_world + public :: mpi_comm_null + public :: mpi_request_null - public :: MPI_SUM - public :: MPI_PROD - public :: MPI_MIN - public :: MPI_MAX - public :: MPI_MINLOC - public :: MPI_MAXLOC + public :: mpi_sum + public :: mpi_prod + public :: mpi_min + public :: mpi_max + public :: mpi_minloc + public :: mpi_maxloc - public :: MPI_MAX_ERROR_STRING - public :: MPI_STATUS_SIZE - public :: MPI_ERROR + public :: mpi_max_error_string + public :: mpi_status_size + public :: mpi_error !#if !defined(sysLinux) - public :: MPI_OFFSET_KIND - public :: MPI_INFO_NULL - public :: MPI_MODE_RDONLY - public :: MPI_MODE_RDWR - public :: MPI_MODE_WRONLY - public :: MPI_MODE_CREATE - public :: MPI_SEEK_SET + public :: mpi_offset_kind + public :: mpi_info_null + public :: mpi_mode_rdonly + public :: mpi_mode_rdwr + public :: mpi_mode_wronly + public :: mpi_mode_create + public :: mpi_seek_set !#endif - public :: MPI_BYTE + public :: mpi_byte #ifdef MPICH_ - public :: MPIPRIV ! the common block name + public :: mpipriv ! the common block name #endif !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ @@ -115,14 +115,14 @@ module mpeu_mpif !----------------------------------------------------------------------- !BOP ! -! !MODULE: mpeu_mpif - a portable interface to the MPI "mpif.h" COMMONs. +! !MODULE: mpeu_mpif - a portable interface to the mpi "mpif.h" COMMONs. ! ! !DESCRIPTION: ! ! The purpose of \verb"mpeu_mpif" module is to provide a portable -! interface of \verb"mpif.h" with different MPI implementation. +! interface of \verb"mpif.h" with different mpi implementation. ! By combining module \verb"mpeu_mpif" and \verb"mpeu_mpif90", it may be -! possible to build a Fortran 90 MPI binding module graduately. +! possible to build a fortran 90 mpi binding module graduately. ! ! Although it is possible to use \verb'include "mpif.h"' directly ! in individual modules, it has several problems: @@ -130,7 +130,7 @@ module mpeu_mpif ! \item It may conflict with either the source code of a {\sl fixed} ! format or the code of a {\sl free} format; ! \item It does not provide the protection and the safety of using -! these variables as what a \verb"MODULE" would provide. +! these variables as what a \verb"module" would provide. ! \end{itemize} ! ! More information may be found in the module \verb"mpeu_mpif90". @@ -141,258 +141,286 @@ module mpeu_mpif ! !REVISION HISTORY: ! 01Apr98 - Jing Guo - initial prototype/prolog/code -! 16Feb05 - Todling - added a number of vars needed by GSI -! 18Mar05 - Todling - added a few more mpi vars used by GSI +! 16Feb05 - Todling - added a number of vars needed by gsi +! 18Mar05 - Todling - added a few more mpi vars used by gsi ! 20Mar05 - Todling - see no reason for sysLinux ifdef (commented out) ! 18Feb09 - Jing Guo -! - Copied from GMAO_mpeu/m_mpif.F here to avoid the -! dependency of GSI_GridComp/ on GMAO_mpeu. +! - Copied from gmao_mpeu/m_mpif.F here to avoid the +! dependency of gsi_gridcomp/ on gmao_mpeu. ! - Renamed to mpeu_mpif to avoid possible name conflict -! with GMAO_mpeu/. +! with gmao_mpeu/. ! 23Feb09 - Jing Guo ! - Renamed to *.F90, by assuming the system mpif.h is ! compatible with free-format Fortran. -! 22Feb17 - Todling - add dummy vars to satisfy EMC requirement +! 22Feb17 - Todling - add dummy vars to satisfy emc requirement ! of all vars being used - this is counter -! advanced FORTRAN interface practices. +! advanced fortran interface practices. !EOP !_______________________________________________________________________ - integer, parameter:: MPI_IKIND=kind(MPI_COMM_WORLD) + integer, parameter:: mpi_ikind=kind(mpi_comm_world) character(len=*),parameter :: myname='mpeu_mpif' contains function type_r0_of_integer1(mold) result(mtype) + use kinds, only: i_byte implicit none - integer(kind=MPI_IKIND):: mtype - integer*1,intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_byte),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER1 + mtype=mpi_integer1 end function type_r0_of_integer1 function type_r0_of_integer2(mold) result(mtype) + use kinds, only: i_short implicit none - integer(kind=MPI_IKIND):: mtype - integer*2,intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_short),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER2 + mtype=mpi_integer2 end function type_r0_of_integer2 function type_r0_of_integer4(mold) result(mtype) + use kinds, only: i_long implicit none - integer(kind=MPI_IKIND):: mtype - integer*4,intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_long),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER4 + mtype=mpi_integer4 end function type_r0_of_integer4 function type_r0_of_integer8(mold) result(mtype) + use kinds, only: i_llong implicit none - integer(kind=MPI_IKIND):: mtype - integer*8,intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_llong),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER8 + mtype=mpi_integer8 end function type_r0_of_integer8 function type_r0_of_real4(mold) result(mtype) + use kinds, only: r_single implicit none - integer(kind=MPI_IKIND):: mtype - real*4,intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_single),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real4 + mtype=mpi_real4 end function type_r0_of_real4 function type_r0_of_real8(mold) result(mtype) + use kinds, only: r_double implicit none - integer(kind=MPI_IKIND):: mtype - real*8,intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_double),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real8 + mtype=mpi_real8 end function type_r0_of_real8 function type_r0_of_real16(mold) result(mtype) + use kinds, only: r_quad implicit none - integer(kind=MPI_IKIND):: mtype - real*16,intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_quad),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real16 + mtype=mpi_real16 end function type_r0_of_real16 function type_r1_of_integer1(mold) result(mtype) + use kinds, only: i_byte implicit none - integer(kind=MPI_IKIND):: mtype - integer*1,dimension(:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_byte),dimension(:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER1 + mtype=mpi_integer1 end function type_r1_of_integer1 function type_r1_of_integer2(mold) result(mtype) + use kinds, only: i_short implicit none - integer(kind=MPI_IKIND):: mtype - integer*2,dimension(:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_short),dimension(:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER2 + mtype=mpi_integer2 end function type_r1_of_integer2 function type_r1_of_integer4(mold) result(mtype) + use kinds, only: i_long implicit none - integer(kind=MPI_IKIND):: mtype - integer*4,dimension(:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_long),dimension(:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER4 + mtype=mpi_integer4 end function type_r1_of_integer4 function type_r1_of_integer8(mold) result(mtype) + use kinds, only: i_llong implicit none - integer(kind=MPI_IKIND):: mtype - integer*8,dimension(:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_llong),dimension(:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER8 + mtype=mpi_integer8 end function type_r1_of_integer8 function type_r1_of_real4(mold) result(mtype) + use kinds, only: r_single implicit none - integer(kind=MPI_IKIND):: mtype - real*4,dimension(:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_single),dimension(:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real4 + mtype=mpi_real4 end function type_r1_of_real4 function type_r1_of_real8(mold) result(mtype) + use kinds, only: r_double implicit none - integer(kind=MPI_IKIND):: mtype - real*8,dimension(:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_double),dimension(:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real8 + mtype=mpi_real8 end function type_r1_of_real8 function type_r1_of_real16(mold) result(mtype) + use kinds, only: r_quad implicit none - integer(kind=MPI_IKIND):: mtype - real*16,dimension(:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_quad),dimension(:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real16 + mtype=mpi_real16 end function type_r1_of_real16 function type_r2_of_integer1(mold) result(mtype) + use kinds, only: i_byte implicit none - integer(kind=MPI_IKIND):: mtype - integer*1,dimension(:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_byte),dimension(:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER1 + mtype=mpi_integer1 end function type_r2_of_integer1 function type_r2_of_integer2(mold) result(mtype) + use kinds, only: i_short implicit none - integer(kind=MPI_IKIND):: mtype - integer*2,dimension(:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_short),dimension(:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER2 + mtype=mpi_integer2 end function type_r2_of_integer2 function type_r2_of_integer4(mold) result(mtype) + use kinds, only: i_long implicit none - integer(kind=MPI_IKIND):: mtype - integer*4,dimension(:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_long),dimension(:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER4 + mtype=mpi_integer4 end function type_r2_of_integer4 function type_r2_of_integer8(mold) result(mtype) + use kinds, only: i_llong implicit none - integer(kind=MPI_IKIND):: mtype - integer*8,dimension(:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_llong),dimension(:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER8 + mtype=mpi_integer8 end function type_r2_of_integer8 function type_r2_of_real4(mold) result(mtype) + use kinds, only: r_single implicit none - integer(kind=MPI_IKIND):: mtype - real*4,dimension(:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_single),dimension(:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real4 + mtype=mpi_real4 end function type_r2_of_real4 function type_r2_of_real8(mold) result(mtype) + use kinds, only: r_double implicit none - integer(kind=MPI_IKIND):: mtype - real*8,dimension(:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_double),dimension(:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real8 + mtype=mpi_real8 end function type_r2_of_real8 function type_r2_of_real16(mold) result(mtype) + use kinds, only: r_quad implicit none - integer(kind=MPI_IKIND):: mtype - real*16,dimension(:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_quad),dimension(:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real16 + mtype=mpi_real16 end function type_r2_of_real16 function type_r3_of_integer1(mold) result(mtype) + use kinds, only: i_byte implicit none - integer(kind=MPI_IKIND):: mtype - integer*1,dimension(:,:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_byte),dimension(:,:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER1 + mtype=mpi_integer1 end function type_r3_of_integer1 function type_r3_of_integer2(mold) result(mtype) + use kinds, only: i_short implicit none - integer(kind=MPI_IKIND):: mtype - integer*2,dimension(:,:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_short),dimension(:,:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER2 + mtype=mpi_integer2 end function type_r3_of_integer2 function type_r3_of_integer4(mold) result(mtype) + use kinds, only: i_long implicit none - integer(kind=MPI_IKIND):: mtype - integer*4,dimension(:,:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_long),dimension(:,:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER4 + mtype=mpi_integer4 end function type_r3_of_integer4 function type_r3_of_integer8(mold) result(mtype) + use kinds, only: i_llong implicit none - integer(kind=MPI_IKIND):: mtype - integer*8,dimension(:,:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + integer(i_llong),dimension(:,:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_INTEGER8 + mtype=mpi_integer8 end function type_r3_of_integer8 function type_r3_of_real4(mold) result(mtype) + use kinds, only: r_single implicit none - integer(kind=MPI_IKIND):: mtype - real*4,dimension(:,:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_single),dimension(:,:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real4 + mtype=mpi_real4 end function type_r3_of_real4 function type_r3_of_real8(mold) result(mtype) + use kinds, only: r_double implicit none - integer(kind=MPI_IKIND):: mtype - real*8,dimension(:,:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_double),dimension(:,:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real8 + mtype=mpi_real8 end function type_r3_of_real8 function type_r3_of_real16(mold) result(mtype) + use kinds, only: r_quad implicit none - integer(kind=MPI_IKIND):: mtype - real*16,dimension(:,:,:),intent(in):: mold + integer(kind=mpi_ikind):: mtype + real(r_quad),dimension(:,:,:),intent(in):: mold integer(mpi_ikind):: dummymold dummymold=kind(mold) - mtype=MPI_real16 + mtype=mpi_real16 end function type_r3_of_real16 end module mpeu_mpif diff --git a/src/gsi/mpimod.F90 b/src/gsi/mpimod.F90 index a249eadeac..baacfabfd1 100644 --- a/src/gsi/mpimod.F90 +++ b/src/gsi/mpimod.F90 @@ -3,7 +3,7 @@ !------------------------------------------------------------------------- !BOP ! -! !MODULE: mpimod --- GSI Module containing mpi related variables +! !MODULE: mpimod --- gsi module containing mpi related variables ! ! !INTERFACE: ! @@ -16,7 +16,7 @@ module mpimod #ifdef ibm_sp ! Include standard mpi includes file. - use mpi + use mpi, only: #else use mpeu_mpif, only : mpi_rtype4 => mpi_real4 #ifdef _REAL4_ @@ -64,13 +64,13 @@ module mpimod ! 2005-01-24 kleist - fix bug in array initialization ! 2005-02-15 todling - add use m_mpif, only for mpi_integer4, ! mpi_offset_kind, ... (only applies -! to non IBM SP machines) +! to non ibm sp machines) ! 2005-07-25 todling - add a couple more exports from m_mpif -! (only applies to non IBM SP machines) +! (only applies to non ibm sp machines) ! 2006-04-06 middlecoff - remove mpi_request_null since not used ! 2006-06-20 treadon - add mpi_itype -! 2006-06-28 da Silva - Added 2 integers represing a layout: nxPE and nyPE. -! 2009-02-19 jing guo - replaced m_mpif of GMAO_mpeu with gmaogsi_mpif. +! 2006-06-28 da Silva - Added 2 integers represing a layout: nxpe and nype. +! 2009-02-19 jing guo - replaced m_mpif of gmao_mpeu with gmaogsi_mpif. ! 2009-04-21 derber - add communications for strong balance constraint (bal) ! and unified uv (vec) transformation ! 2010-04-01 treadon - remove routines reorder, reorder2, strip_single, strip, @@ -78,7 +78,7 @@ module mpimod ! routines are now found in gridmod ! 2010-05-23 todling - nvarbal_id no longer wired to 1,2,3,4, rather linked ! to where fields are in control vector -! 2011-07-04 todling - allow proper setting of REAL*4 or REAL*8 +! 2011-07-04 todling - allow proper setting of real*4 or real*8 ! 2012-06-12 parrish - remove all communication variables, except for levs_id, nvar_id, and nvar_pe. ! Remove subroutines init_mpi_vars and destroy_mpi_vars. ! All communication variables that used to be here are now created in @@ -123,13 +123,13 @@ module mpimod #endif integer(i_kind) ierror - integer(i_kind) :: npe ! total num of MPI tasks - integer(i_kind) :: mype ! number of MPI task + integer(i_kind) :: npe ! total num of mpi tasks + integer(i_kind) :: mype ! number of mpi task -! Optional ESMF-like layout information: nxPE is the number of -! processors used to decompose the longitudinal dimensional, while nyPE +! Optional esmf-like layout information: nxpe is the number of +! processors used to decompose the longitudinal dimensional, while nype ! the number of processors used to decompose the latitudinal dimension. -! By construction, nPE = nxPE * nyPE. +! By construction, npe = nxpe * nype. ! integer(i_kind) :: nxpe=-1 ! optional layout information integer(i_kind) :: nype=-1 ! optional layout information diff --git a/src/gsi/mpl_allreduce.F90 b/src/gsi/mpl_allreduce.F90 index d01034cfe0..2d42130765 100644 --- a/src/gsi/mpl_allreduce.F90 +++ b/src/gsi/mpl_allreduce.F90 @@ -25,22 +25,22 @@ module mpl_allreducemod !$$$ end documentation block implicit none -PRIVATE -PUBLIC mpl_allreduce -PUBLIC mpl_allgather -PUBLIC mpl_reduce +private +public mpl_allreduce +public mpl_allgather +public mpl_reduce -INTERFACE mpl_allreduce -MODULE PROCEDURE rmpl_allreduce,qmpl_allreduce1d,qmpl_allreduce2d -END INTERFACE +interface mpl_allreduce +module procedure rmpl_allreduce,qmpl_allreduce1d,qmpl_allreduce2d +END interface -INTERFACE mpl_allgather -MODULE PROCEDURE mpl_allgatherq -END INTERFACE +interface mpl_allgather +module procedure mpl_allgatherq +END interface -INTERFACE mpl_reduce ! same as mpl_allreduce routines except reduces to 1 processor -MODULE PROCEDURE rmpl_reduce,qmpl_reduce1d,qmpl_reduce2d -END INTERFACE +interface mpl_reduce ! same as mpl_allreduce routines except reduces to 1 processor +module procedure rmpl_reduce,qmpl_reduce1d,qmpl_reduce2d +END interface contains @@ -87,7 +87,7 @@ subroutine rmpl_allreduce(klen,rpvals) call mpi_allgather(rpvals,klen,mpi_rtype, & zwork,klen,mpi_rtype, mpi_comm_world,ierror) -! Sum (note this is NOT reproducible across different numbers of processors) +! Sum (note this is not reproducible across different numbers of processors) do ii=1,klen rpvals(ii)=zwork(ii,1) end do @@ -100,7 +100,7 @@ subroutine rmpl_allreduce(klen,rpvals) endif ! ---------------------------------------------------------- -return + return end subroutine rmpl_allreduce ! ---------------------------------------------------------- subroutine qmpl_allreduce1d(klen,qpvals) @@ -176,7 +176,7 @@ subroutine qmpl_allreduce1d(klen,qpvals) endif ! ---------------------------------------------------------- -return + return end subroutine qmpl_allreduce1d ! ---------------------------------------------------------- subroutine qmpl_allreduce2d(ilen,klen,pvals,pvnew) @@ -243,9 +243,9 @@ subroutine qmpl_allreduce2d(ilen,klen,pvals,pvnew) if (present(pvnew)) then do kk=1,klen - do ii=1,ilen - pvnew(ii,kk)=pval2(ii,kk,1) - end do + do ii=1,ilen + pvnew(ii,kk)=pval2(ii,kk,1) + end do end do do nn=2,npe do kk=1,klen @@ -258,9 +258,9 @@ subroutine qmpl_allreduce2d(ilen,klen,pvals,pvnew) else do kk=1,klen - do ii=1,ilen - pvals(ii,kk)=pval2(ii,kk,1) - end do + do ii=1,ilen + pvals(ii,kk)=pval2(ii,kk,1) + end do end do do nn=2,npe do kk=1,klen @@ -273,7 +273,7 @@ subroutine qmpl_allreduce2d(ilen,klen,pvals,pvnew) endif ! ---------------------------------------------------------- -return + return end subroutine qmpl_allreduce2d subroutine mpl_allgatherq(idim,jdim,zloc,zall) @@ -381,23 +381,23 @@ subroutine rmpl_reduce(klen,iroot,rpvals) zwork,klen,mpi_rtype, iroot,mpi_comm_world,ierror) if(mype == iroot)then -! Sum (note this is NOT reproducible across different numbers of processors) - do ii=1,klen - rpvals(ii)=zwork(ii,1) - end do - do jj=2,npe - do ii=1,klen - rpvals(ii)=rpvals(ii)+zwork(ii,jj) - end do - end do +! Sum (note this is not reproducible across different numbers of processors) + do ii=1,klen + rpvals(ii)=zwork(ii,1) + end do + do jj=2,npe + do ii=1,klen + rpvals(ii)=rpvals(ii)+zwork(ii,jj) + end do + end do else - rpvals = 0._r_kind + rpvals = 0._r_kind end if endif ! ---------------------------------------------------------- -return + return end subroutine rmpl_reduce ! ---------------------------------------------------------- subroutine qmpl_reduce1d(klen,iroot,qpvals) @@ -472,13 +472,13 @@ subroutine qmpl_reduce1d(klen,iroot,qpvals) end do end do else - qpvals = 0._r_kind + qpvals = 0._r_kind end if endif ! ---------------------------------------------------------- -return + return end subroutine qmpl_reduce1d ! ---------------------------------------------------------- subroutine qmpl_reduce2d(ilen,klen,iroot,pvals,pvnew) @@ -490,7 +490,7 @@ subroutine qmpl_reduce2d(ilen,klen,iroot,pvals,pvnew) ! abstract: Reproducible (across different pe's) reduce ! ! program history log: -! 2008-12-09 todling - embed Derber's reproducible sum in subroutine +! 2008-12-09 todling - embed derber's reproducible sum in subroutine ! ! input argument list: ! ilen - first dimension of array pvals @@ -551,9 +551,9 @@ subroutine qmpl_reduce2d(ilen,klen,iroot,pvals,pvnew) if (present(pvnew)) then do kk=1,klen - do ii=1,ilen - pvnew(ii,kk)=pval2(ii,kk,1) - end do + do ii=1,ilen + pvnew(ii,kk)=pval2(ii,kk,1) + end do end do do nn=2,npe do kk=1,klen @@ -566,9 +566,9 @@ subroutine qmpl_reduce2d(ilen,klen,iroot,pvals,pvnew) else do kk=1,klen - do ii=1,ilen - pvals(ii,kk)=pval2(ii,kk,1) - end do + do ii=1,ilen + pvals(ii,kk)=pval2(ii,kk,1) + end do end do do nn=2,npe do kk=1,klen @@ -582,7 +582,7 @@ subroutine qmpl_reduce2d(ilen,klen,iroot,pvals,pvnew) end if ! ---------------------------------------------------------- -return + return end subroutine qmpl_reduce2d end module mpl_allreducemod