subroutine read_radarref_mosaic(nread,ndata,infile,obstype,lunout,twind,sis,nobs) !$$$ subprogram documentation block ! . . . . ! subprogram: read_radarref_mosaic Reading in reflectivity mosaic in RR grid ! ! PRGMMR: Ming Hu ORG: NP22 DATE: 2006-03-27 ! ! ABSTRACT: ! This routine read in reflectivity mosaic data. The data has already ! been interpolated into analysis grid and in form of BUFR. ! ! PROGRAM HISTORY LOG: ! 2008-12-20 Hu make it read in BUFR form reflectivity data ! 2010-04-09 Hu make changes based on current trunk style ! 2013-03-27 Hu add code to map obs from WRF mass H grid to analysis grid ! 2015-02-23 Rancic/Thomas - add l4densvar to time window logical ! 2015-03-23 Su -fix array size with maximum message and subset number from ! fixed number to dynamic allocated array ! ! input argument list: ! infile - unit from which to read mosaic information file ! obstype - observation type to process ! lunout - unit to which to write data for further processing ! twind - input group time window (hours) ! sis - observation variable name ! ! output argument list: ! nread - number of type "obstype" observations read ! ndata - number of type "obstype" observations retained for further processing ! nobs - array of observations on each subdomain for each processor ! ! USAGE: ! INPUT FILES: refInGSI ! ! OUTPUT FILES: ! ! REMARKS: ! ! ATTRIBUTES: ! LANGUAGE: FORTRAN 90 ! MACHINE: Linux cluster(Wjet) ! !$$$ ! !_____________________________________________________________________ ! use kinds, only: r_kind,r_double,i_kind use constants, only: zero,one use convinfo, only: nconvtype,ctwind,icuse,ioctype use gsi_4dvar, only: l4dvar,l4densvar,winlen use gridmod, only: nlon,nlat,nlon_regional,nlat_regional use mod_wrfmass_to_a, only: wrfmass_obs_to_a8 use mpimod, only: npe implicit none ! character(10), intent(in) :: infile,obstype integer(i_kind), intent(in) :: lunout integer(i_kind), intent(inout) :: nread,ndata integer(i_kind),dimension(npe) ,intent(inout) :: nobs real(r_kind), intent(in ) :: twind character(20), intent(in) :: sis ! ! For reflectiivty mosaic ! integer(i_kind) nreal,nchanl integer(i_kind) ifn,i integer(i_kind) :: ilon,ilat logical :: nsslrefobs ! ! for read in bufr ! real(r_kind) :: hdr(5),obs(1,35) character(80):: hdrstr='SID XOB YOB DHR TYP' character(80):: obsstr='HREF' character(8) subset integer(i_kind) :: lunin,idate integer(i_kind) :: ireadmg,ireadsb integer(i_kind) :: maxlvl integer(i_kind) :: numlvl,numref,numobsa,nmsgmax,maxobs integer(i_kind) :: k,iret integer(i_kind) :: nmsg,ntb real(r_kind),allocatable :: ref3d_column(:,:) ! 3D reflectivity in column integer(i_kind) :: ikx real(r_kind) :: timeo,t4dv !********************************************************************** ! ! END OF DECLARATIONS....start of program ! nsslrefobs = .false. ikx=0 do i=1,nconvtype if(trim(obstype) == trim(ioctype(i)) .and. abs(icuse(i))== 1) then nsslrefobs=.true. ikx=i endif end do nread=0 ndata=0 nchanl=0 ifn = 15 if(nsslrefobs) then !! get message and subset counts call getcount_bufr(infile,nmsgmax,maxobs) lunin = 10 maxlvl= 31 allocate(ref3d_column(maxlvl+2,maxobs)) open ( unit = lunin, file = trim(infile),form='unformatted',err=200) call openbf ( lunin, 'IN', lunin ) call datelen ( 10 ) nmsg=0 ntb = 0 msg_report: do while (ireadmg(lunin,subset,idate) == 0) nmsg=nmsg+1 if (nmsg>nmsgmax) then write(6,*)'read_radarref_mosaic: messages exceed maximum ',nmsgmax call stop2(50) endif loop_report: do while (ireadsb(lunin) == 0) ntb = ntb+1 if (ntb>maxobs) then write(6,*)'read_radarref_mosaic: reports exceed maximum ',maxobs call stop2(50) endif ! Extract type, date, and location information call ufbint(lunin,hdr,5,1,iret,hdrstr) ! check time window in subset if (l4dvar.or.l4densvar) then t4dv=hdr(4) if (t4dvwinlen) then write(6,*)'read_radarref_mosaic: time outside window ',& t4dv,' skip this report' cycle loop_report endif else timeo=hdr(4) if (abs(timeo)>ctwind(ikx) .or. abs(timeo) > twind) then write(6,*)'read_radarref_mosaic: time outside window ',& timeo,' skip this report' cycle loop_report endif endif ! read in observations call ufbint(lunin,obs,1,35,iret,obsstr) numlvl=min(iret,maxlvl) ref3d_column(1,ntb)=hdr(2)*10.0_r_kind ! observation location, grid index i ref3d_column(2,ntb)=hdr(3)*10.0_r_kind ! observation location, grid index j do k=1,numlvl ref3d_column(2+k,ntb)=obs(1,k) ! reflectivity (column 31 levels) enddo enddo loop_report enddo msg_report write(6,*)'read_radarref_mosaic: messages/reports = ',nmsg,'/',ntb numref=ntb ! ! covert BUFR value of missing (-64) and no echo (-63) to cloud analysis ! value of missing (-999.0) and no echo (-99.0) ! DO i=1,numref DO k=1,maxlvl if( abs(ref3d_column(k+2,i)+64.0_r_kind) <= 0.00001_r_kind) then ref3d_column(k+2,i)=-999.0_r_kind elseif( abs(ref3d_column(k+2,i)+63.0_r_kind) <= 0.00001_r_kind) then ref3d_column(k+2,i)=-99.0_r_kind endif enddo enddo ilon=1 ilat=2 nread=numref ndata=numref nreal=maxlvl+2 if(numref > 0 ) then if(nlon==nlon_regional .and. nlat==nlat_regional) then call count_obs(numref,maxlvl+2,ilat,ilon,ref3d_column,nobs) write(lunout) obstype,sis,nreal,nchanl,ilat,ilon write(lunout) ((ref3d_column(k,i),k=1,maxlvl+2),i=1,numref) else call wrfmass_obs_to_a8(ref3d_column,nreal,numref,ilat,ilon,numobsa) nread=numobsa ndata=numobsa call count_obs(numobsa,maxlvl+2,ilat,ilon,ref3d_column,nobs) write(lunout) obstype,sis,nreal,nchanl,ilat,ilon write(lunout) ((ref3d_column(k,i),k=1,maxlvl+2),i=1,numobsa) endif deallocate(ref3d_column) endif endif call closbf(lunin) close(lunin) return 200 continue write(6,*) 'read_radarref_mosaic, Warning : cannot find radar data file' end subroutine read_radarref_mosaic ! !