subroutine statsconv(mype,& i_ps,i_uv,i_srw,i_t,i_q,i_pw,i_rw,i_dw,i_gps,i_sst,i_tcp,i_lag, & bwork,awork,ndata) !$$$ subprogram documentation block ! . . . . ! subprogram: statconv prints statistics for conventional data ! prgmmr: derber org: np23 date: 2003-05-22 ! ! abstract: The routine computes and prints statistics regarding the ! use of conventional observations. Printed information ! includes that about data counts, quality control decisions, ! statistics based on the innovations, and penalties - all ! as a observation type ! ! program history log: ! 2003-05-22 derber ! 2004-06-15 treadon - update documentation ! 2004-07-20 derber - add sst ! 2004-07-29 treadon - add only to module use, add intent in/out ! 2004-10-06 parrish - increase diagnostic array sizes and add ! output from nonlinear qc ! 2004-12-23 treadon - use module jfunc to pass jiter,first ! 2005-01-28 cucurull - modify summary output for refractivity to include ! the QC checks on incremental refractivity ! 2005-03-23 cucurull - cosmetic changes for refractivity print out ! 2005-04-20 treadon - correct error in wind ntot sum ! 2005-05-27 derber - level output changed ! 2005-07-27 derber - add print of monitoring and reject data ! 2005-12-02 cucurull - cosmetic changes for gps data ! 2006-02-03 derber - modify for new obs control and to clean up output ! 2006-02-24 derber - modify to take advantage of convinfo module ! 2006-04-02 derber - modify to eliminate dvast and move ob type printing to dtast ! 2008-04-11 safford - rm unused uses ! 2009-02-02 kleist - add synthetic tc-mslp ! 2009-03-05 meunier - add lagrangean data ! ! input argument list: ! mype - mpi task number ! i_ps - index in awork array holding surface pressure info ! i_uv - index in awork array holding wind info ! i_srw - index in awork array holding radar wind superobs info ! i_t - index in awork array holding temperature info ! i_q - index in awork array holding specific humidity info ! i_pw - index in awork array holding total precipitable water info ! i_rw - index in awork array holding radar radial winds info ! i_dw - index in awork array holding doppler lidar winds info ! i_gps - index in awork array holding gps info ! i_sst - index in awork array holding sst info ! i_tcp - index in awork array holding tcps info ! i_lag - index in awork array holding lag info ! bwork - array containing information for statistics ! awork - array containing information for data counts and gross checks ! ndata(*,1)- number of profiles retained for further processing ! ndata(*,2)- number of observations read ! ndata(*,3)- number of observations keep after read ! ! output argument list: ! ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ use kinds, only: r_kind,i_kind use constants, only: izero,ione,zero,three,five,izero use obsmod, only: iout_sst,iout_pw,iout_t,iout_rw,iout_dw,& iout_srw,iout_uv,iout_gps,iout_ps,iout_q,iout_tcp,iout_lag,& mype_dw,mype_rw,mype_srw,mype_sst,mype_gps,mype_uv,mype_ps,& mype_t,mype_pw,mype_q,mype_tcp,ndat,dtype,mype_lag use qcmod, only: npres_print,ptop,pbot,ptopq,pbotq use jfunc, only: first,jiter use gridmod, only: nsig use convinfo, only: nconvtype,ioctype implicit none ! Declare passed variables integer(i_kind) ,intent(in ) :: mype,i_ps,i_uv,& i_srw,i_t,i_q,i_pw,i_rw,i_dw,i_gps,i_sst,i_tcp,i_lag real(r_kind),dimension(7*nsig+100_i_kind,13) ,intent(in ) :: awork real(r_kind),dimension(npres_print,nconvtype,5,3),intent(in ) :: bwork integer(i_kind),dimension(ndat,3) ,intent(in ) :: ndata ! Declare local variables character(100) mesage integer(i_kind) numgrspw,numsst,nsuperp,nump,nhitopo,ntoodif integer(i_kind) numgrsq,numhgh integer(i_kind) ntot,numlow,k,numssm,i,j integer(i_kind) numgross,numfailqc,numfailqc_ssmi,nread,nkeep integer(i_kind) numfail1_gps,numfail2_gps,numfail3_gps,nreadspd,nkeepspd integer(i_kind),dimension(nsig)::num real(r_kind) grsmlt,tq,pw,rat,tgps,qmplty,tpw,tdw,rwmplty,trw real(r_kind) tmplty,tt,dwmplty,gpsmplty,umplty,tssm,qctssm,tu,tv,tuv real(r_kind) vmplty,uvqcplty,rat1,rat2,rat3 real(r_kind) dwqcplty,tqcplty,qctt,qctrw,rwqcplty,qctdw,qqcplty,qctgps real(r_kind) gpsqcplty,tpw3,pw3,qctq real(r_kind),dimension(1):: pbotall,ptopall logical,dimension(nconvtype):: pflag !********************************************************************************* ! Initialize constants and variables. ptopall(1)=zero; pbotall(1)=2000.0_r_kind ! Generate summary statistics as a function of observation type. ! Extensive comments are given for the winds. Similar comments ! apply for the remaining observation types in this routine pflag=.false. ! Summary report for winds if(mype==mype_uv) then ! Open output file so as to point to correct position in output file if(first)then open(iout_uv) else open(iout_uv,position='append') end if ! Compute and write counts, penalties, and ratio of penalty ! to data counts for each model level numssm=nint(awork(6,i_uv)); numgross=nint(awork(4,i_uv)) umplty=zero; vmplty=zero; uvqcplty=zero ; ntot=izero; tu=zero; tv=zero ; tuv=zero tssm=zero ; qctssm=zero nread=izero nkeep=izero nreadspd=izero nkeepspd=izero do i=1,ndat if(dtype(i)== 'uv')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) else if(dtype(i)== 'spd')then nreadspd=nreadspd+ndata(i,2) nkeepspd=nkeepspd+ndata(i,3) end if end do if(nkeep > izero .or. nkeepspd > izero)then ! Write header information mesage='current vfit of wind data, ranges in m/s$' ! Call routine to compute and write count, rms, and penalty information do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'spd' .or. trim(ioctype(j)) == 'uv' end do call dtast(bwork,npres_print,pbot,ptop,mesage,jiter,iout_uv,pflag) numlow = nint(awork(2,i_uv)) numhgh = nint(awork(3,i_uv)) write(iout_uv,900) 'wind',numhgh,numlow numfailqc=nint(awork(21,i_uv)) ! keep a seperate record of numfailqc for ssmi wind speeds numfailqc_ssmi=nint(awork(61,i_uv)) do k=1,nsig num(k)=nint(awork(6*nsig+k+100_i_kind,i_uv)) rat1=zero rat2=zero if(num(k) > izero)then rat1=awork(4*nsig+k+100_i_kind,i_uv)/float(num(k)) rat2=awork(5*nsig+k+100_i_kind,i_uv)/float(num(k)) end if umplty=umplty+awork(4*nsig+k+100_i_kind,i_uv) vmplty=vmplty+awork(5*nsig+k+100_i_kind,i_uv) ntot=ntot+num(k) write(iout_uv,241) 'w',num(k),k,awork(4*nsig+k+100_i_kind,i_uv),& awork(5*nsig+k+100_i_kind,i_uv),rat1,rat2 end do do k=1,nsig num(k)=nint(awork(6*nsig+k+100_i_kind,i_uv)) rat1=zero rat3=zero if(num(k) > izero)then rat1=(awork(4*nsig+k+100_i_kind,i_uv)+awork(5*nsig+k+100_i_kind,i_uv))/float(num(k)) rat3=awork(3*nsig+k+100_i_kind,i_uv)/float(num(k)) end if uvqcplty=uvqcplty+awork(3*nsig+k+100_i_kind,i_uv) write(iout_uv,240) 'w',num(k),k,awork(4*nsig+k+100_i_kind,i_uv)+awork(5*nsig+k+100_i_kind,i_uv), & awork(3*nsig+k+100_i_kind,i_uv),rat1,rat3 end do ! Write statistics gross checks write(iout_uv,*)' number ssm/i winds that fail nonlinear qc =',numfailqc_ssmi write(iout_uv,925) 'wind',numgross,numfailqc ! Write statistics regarding penalties if(ntot > izero)then tu=umplty/float(ntot) tv=vmplty/float(ntot) tuv=uvqcplty/float(ntot) end if if(numssm > izero)then tssm=awork(5,i_uv)/awork(6,i_uv) qctssm=awork(22,i_uv)/awork(6,i_uv) end if end if write(iout_uv,949) 'u',ntot,umplty,tu write(iout_uv,949) 'v',ntot,vmplty,tv write(iout_uv,950) 'uv',jiter,nread,nkeep,ntot*2 write(iout_uv,951) 'uv',umplty+vmplty,uvqcplty,tu+tv,tuv write(iout_uv,950) 'spd',jiter,nreadspd,nkeepspd,numssm write(iout_uv,951) 'spd',awork(5,i_uv),awork(22,i_uv),tssm,qctssm ! Close unit receiving summary output close(iout_uv) end if ! Summary report for radar wind superobs if(mype==mype_srw) then if(first)then open(iout_srw) else open(iout_srw,position='append') end if umplty=zero; vmplty=zero; uvqcplty=zero ; ntot=izero; tu=zero; tv=zero ; tuv=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'srw')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current fit of radar superob wind data, ranges in stderr$' do j=1,nconvtype pflag(j)= trim(ioctype(j)) == 'srw' end do call dtast(bwork,npres_print,pbot,ptop,mesage,jiter,iout_srw,pflag) do k=1,nsig num(k)=nint(awork(6*nsig+k+100_i_kind,i_srw)) rat1=zero rat2=zero rat3=zero if(num(k) > izero)then rat1=awork(4*nsig+k+100_i_kind,i_srw)/float(num(k)) rat2=awork(5*nsig+k+100_i_kind,i_srw)/float(num(k)) rat3=awork(3*nsig+k+100_i_kind,i_srw)/float(num(k)) end if umplty=umplty+awork(4*nsig+k+100_i_kind,i_srw) vmplty=vmplty+awork(5*nsig+k+100_i_kind,i_srw) uvqcplty=uvqcplty+awork(3*nsig+k+100_i_kind,i_srw) ntot=ntot+num(k) write(iout_srw,241) 's',num(k),k,awork(4*nsig+k+100_i_kind,i_srw),& awork(5*nsig+k+100_i_kind,i_srw),awork(3*nsig+k+100_i_kind,i_srw),rat1,rat2,rat3 end do numgross=nint(awork(4,i_srw)) numfailqc=nint(awork(21,i_srw)) write(iout_srw,925) 'srw',numgross,numfailqc tu=umplty/float(ntot) tv=vmplty/float(ntot) tuv=uvqcplty/float(ntot) numlow = nint(awork(2,i_srw)) numhgh = nint(awork(3,i_srw)) write(iout_srw,900) 'srw',numhgh,numlow end if write(iout_srw,950) 'srw1',jiter,nread,nkeep,ntot write(iout_srw,951) 'srw1',umplty,uvqcplty,tu,tuv write(iout_srw,950) 'srw2',jiter,nread,nkeep,ntot write(iout_srw,951) 'srw2',vmplty,uvqcplty,tv,tuv close(iout_srw) end if ! Summary report for gps if (mype==mype_gps)then if(first)then open(iout_gps) else open(iout_gps,position='append') end if gpsmplty=zero; gpsqcplty=zero ; ntot=izero tgps=zero ; qctgps=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'gps_ref' .or. dtype(i) == 'gps_bnd')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current fit of gps data in fractional difference$' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'gps_ref' .or. trim(ioctype(j)) == 'gps_bnd' end do call dtast(bwork,npres_print,pbot,ptop,mesage,jiter,iout_gps,pflag) do k=1,nsig num(k)=nint(awork(5*nsig+k+100_i_kind,i_gps)) rat=zero rat3=zero if(num(k)>0) then rat=awork(6*nsig+k+100_i_kind,i_gps)/float(num(k)) rat3=awork(3*nsig+k+100_i_kind,i_gps)/float(num(k)) end if ntot=ntot+num(k); gpsmplty=gpsmplty+awork(6*nsig+k+100_i_kind,i_gps) gpsqcplty=gpsqcplty+awork(3*nsig+k+100_i_kind,i_gps) write(iout_gps,240)'gps',num(k),k,awork(6*nsig+k+100_i_kind,i_gps), & awork(3*nsig+k+100_i_kind,i_gps),rat,rat3 end do numgross=nint(awork(4,i_gps)) numfailqc=nint(awork(21,i_gps)) numfail1_gps=nint(awork(22,i_gps)) numfail2_gps=nint(awork(23,i_gps)) numfail3_gps=nint(awork(24,i_gps)) write(iout_gps,925)'gps',numgross,numfailqc write(iout_gps,*)' number of gps obs failed stats qc in NH =',numfail1_gps write(iout_gps,*)' number of gps obs failed stats qc in SH =',numfail2_gps write(iout_gps,*)' number of gps obs failed stats qc in TR =',numfail3_gps numlow = nint(awork(2,i_gps)) numhgh = nint(awork(3,i_gps)) write(iout_gps,900) 'gps',numhgh,numlow tgps=gpsmplty/ntot qctgps=gpsqcplty/ntot end if write(iout_gps,950)'gps',jiter,nread,nkeep,ntot write(iout_gps,951)'gps',gpsmplty,gpsqcplty,tgps,qctgps close(iout_gps) endif ! Summary report for specific humidity if(mype==mype_q) then if(first)then open(iout_q) else open(iout_q,position='append') end if mesage='current fit of q data, units in per-cent of guess q-sat$' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'q' end do call dtast(bwork,npres_print,pbotq,ptopq,mesage,jiter,iout_q,pflag) qmplty=zero; qqcplty=zero ; ntot=izero tq=zero ; qctq=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'q')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then do k=1,nsig num(k)=nint(awork(k+6*nsig+100_i_kind,i_q)) rat=zero rat3=zero if(num(k) > izero)then rat=awork(5*nsig+k+100_i_kind,i_q)/float(num(k)) rat3=awork(3*nsig+k+100_i_kind,i_q)/float(num(k)) end if qmplty=qmplty+awork(5*nsig+k+100_i_kind,i_q) qqcplty=qqcplty+awork(3*nsig+k+100_i_kind,i_q) ntot=ntot+num(k) write(iout_q,240) 'q',num(k),k,awork(5*nsig+k+100_i_kind,i_q), & awork(3*nsig+k+100_i_kind,i_q),rat,rat3 end do grsmlt=five numgrsq=nint(awork(4,i_q)) numfailqc=nint(awork(21,i_q)) write(iout_q,*)' (scaled as precent of guess specific humidity)' write(iout_q,925) 'q',numgrsq,numfailqc write(iout_q,975) grsmlt,'q',awork(5,i_q) numlow = nint(awork(2,i_q)) numhgh = nint(awork(3,i_q)) write(iout_q,900) 'q',numhgh,numlow tq=qmplty/float(ntot) qctq=qqcplty/float(ntot) end if write(iout_q,950) 'q',jiter,nread,nkeep,ntot write(iout_q,951) 'q',qmplty,qqcplty,tq,qctq close(iout_q) end if ! Summary report for surface pressure if(mype==mype_ps) then if(first)then open(iout_ps) else open(iout_ps,position='append') end if nump=nint(awork(5,i_ps)) pw=zero ; pw3=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'ps')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current fit of surface pressure data, ranges in mb$' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'ps' end do call dtast(bwork,ione,pbotall,ptopall,mesage,jiter,iout_ps,pflag) numgross=nint(awork(6,i_ps)) numfailqc=nint(awork(21,i_ps)) write(iout_ps,925) 'psfc',numgross,numfailqc if(nump > izero)then pw=awork(4,i_ps)/float(nump) pw3=awork(22,i_ps)/float(nump) end if end if write(iout_ps,950) 'psfc',jiter,nread,nkeep,nump write(iout_ps,951) 'psfc',awork(4,i_ps),awork(22,i_ps),pw,pw3 close(iout_ps) end if ! Summary report for total precipitable water if(mype==mype_pw) then if(first)then open(iout_pw) else open(iout_pw,position='append') end if nsuperp=nint(awork(4,i_pw)) tpw=zero ; tpw3=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'pw')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current fit of precip. water data, ranges in mm$' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'pw' end do call dtast(bwork,ione,pbotall,ptopall,mesage,jiter,iout_pw,pflag) numgrspw=nint(awork(6,i_pw)) numfailqc=nint(awork(21,i_pw)) grsmlt=three tpw=zero tpw3=zero if(nsuperp > izero)then tpw=awork(5,i_pw)/nsuperp tpw3=awork(22,i_pw)/nsuperp end if write(iout_pw,925) 'p.w.',numgrspw,numfailqc write(iout_pw,975) grsmlt,'p.w.',awork(7,i_pw) end if write(iout_pw,950) 'pw',jiter,nread,nkeep,nsuperp write(iout_pw,951) 'pw',awork(5,i_pw),awork(22,i_pw),tpw,tpw3 close(iout_pw) end if ! Summary report for conventional sst if(mype==mype_sst) then if(first)then open(iout_sst) else open(iout_sst,position='append') end if numsst=nint(awork(5,i_sst)) pw=zero ; pw3=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'sst')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current fit of conventional sst data, ranges in C$' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'sst' end do call dtast(bwork,ione,pbotall,ptopall,mesage,jiter,iout_sst,pflag) numgross=nint(awork(6,i_sst)) numfailqc=nint(awork(21,i_sst)) if(numsst > izero)then pw=awork(4,i_sst)/numsst pw3=awork(22,i_sst)/numsst end if write(iout_sst,925) 'sst',numgross,numfailqc end if write(iout_sst,950) 'sst',jiter,nread,nkeep,numsst write(iout_sst,951) 'sst',awork(4,i_sst),awork(22,i_sst),pw,pw3 close(iout_sst) end if ! Summary report for temperature if (mype==mype_t)then if(first)then open(iout_t) else open(iout_t,position='append') end if tmplty=zero; tqcplty=zero ; ntot=izero tt=zero ; qctt=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 't')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current fit of temperature data, ranges in K $' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 't' end do call dtast(bwork,npres_print,pbot,ptop,mesage,jiter,iout_t,pflag) do k=1,nsig num(k)=nint(awork(5*nsig+k+100_i_kind,i_t)) rat=zero ; rat3=zero if(num(k) > izero) then rat=awork(6*nsig+k+100_i_kind,i_t)/float(num(k)) rat3=awork(3*nsig+k+100_i_kind,i_t)/float(num(k)) end if ntot=ntot+num(k); tmplty=tmplty+awork(6*nsig+k+100_i_kind,i_t) tqcplty=tqcplty+awork(3*nsig+k+100_i_kind,i_t) write(iout_t,240) 't',num(k),k,awork(6*nsig+k+100_i_kind,i_t), & awork(3*nsig+k+100_i_kind,i_t),rat,rat3 end do numgross=nint(awork(4,i_t)) numfailqc=nint(awork(21,i_t)) write(iout_t,925) 'temp',numgross,numfailqc numlow = nint(awork(2,i_t)) numhgh = nint(awork(3,i_t)) write(iout_t,900) 't',numhgh,numlow tt=tmplty/ntot qctt=tqcplty/ntot end if write(iout_t,950) 't',jiter,nread,nkeep,ntot write(iout_t,951) 't',tmplty,tqcplty,tt,qctt close(iout_t) endif ! Summary report for doppler lidar winds if(mype==mype_dw) then if(first)then open(iout_dw) else open(iout_dw,position='append') end if dwmplty=zero; dwqcplty=zero ; ntot=izero tdw=zero ; qctdw=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'dw')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current vfit of lidar wind data, ranges in m/s$' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'dw' end do call dtast(bwork,npres_print,pbot,ptop,mesage,jiter,iout_dw,pflag) do k=1,nsig num(k)=nint(awork(k+5*nsig+100_i_kind,i_dw)) rat=zero rat3=zero if(num(k) > izero) then rat=awork(6*nsig+k+100_i_kind,i_dw)/float(num(k)) rat3=awork(3*nsig+k+100_i_kind,i_dw)/float(num(k)) end if ntot=ntot+num(k) dwmplty=dwmplty+awork(6*nsig+k+100_i_kind,i_dw) dwqcplty=dwqcplty+awork(3*nsig+k+100_i_kind,i_dw) write(iout_dw,240) 'r',num(k),k,awork(6*nsig+k+100_i_kind,i_dw), & awork(3*nsig+k+100_i_kind,i_dw),rat,rat3 end do numgross=nint(awork(4,i_dw)) numfailqc=nint(awork(21,i_dw)) tdw=dwmplty/float(ntot) qctdw=dwqcplty/float(ntot) write(iout_dw,925) 'dw',numgross,numfailqc numlow = nint(awork(2,i_dw)) numhgh = nint(awork(3,i_dw)) write(iout_dw,900) 'dw',numhgh,numlow end if write(iout_dw,950) 'dw',jiter,nread,nkeep,ntot write(iout_dw,951) 'dw',dwmplty,dwqcplty,tdw,qctdw close(iout_dw) end if ! Summary report for radar radial winds if(mype==mype_rw) then if(first)then open(iout_rw) else open(iout_rw,position='append') end if rwmplty=zero; rwqcplty=zero ; ntot=izero trw=zero ; qctrw=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'rw')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current vfit of radar wind data, ranges in m/s$' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'rw' end do call dtast(bwork,npres_print,pbot,ptop,mesage,jiter,iout_rw,pflag) numgross=nint(awork(4,i_rw)) numfailqc=nint(awork(21,i_rw)) do k=1,nsig num(k)=nint(awork(k+5*nsig+100_i_kind,i_rw)) rat=zero rat3=zero if(num(k) > izero) then rat=awork(6*nsig+k+100_i_kind,i_rw)/float(num(k)) rat3=awork(3*nsig+k+100_i_kind,i_rw)/float(num(k)) end if ntot=ntot+num(k) rwmplty=rwmplty+awork(6*nsig+k+100_i_kind,i_rw) rwqcplty=rwqcplty+awork(3*nsig+k+100_i_kind,i_rw) write(iout_rw,240) 'r',num(k),k,awork(6*nsig+k+100_i_kind,i_rw), & awork(3*nsig+k+100_i_kind,i_rw),rat,rat3 end do trw=rwmplty/float(ntot) qctrw=rwqcplty/float(ntot) write(iout_rw,925) 'rw',numgross,numfailqc numlow = nint(awork(2,i_rw)) numhgh = nint(awork(3,i_rw)) nhitopo = nint(awork(5,i_rw)) ntoodif = nint(awork(6,i_rw)) write(iout_rw,900) 'rw',numhgh,numlow write(iout_rw,905) 'rw',nhitopo,ntoodif end if write(iout_rw,950) 'rw',jiter,nread,nkeep,ntot write(iout_rw,951) 'rw',rwmplty,rwqcplty,trw,qctrw close(iout_rw) end if if(mype==mype_tcp) then if(first)then open(iout_tcp) else open(iout_tcp,position='append') end if nump=nint(awork(5,i_tcp)) pw=zero ; pw3=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'tcp')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current fit of surface pressure data, ranges in mb$' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'tcp' end do call dtast(bwork,ione,pbotall,ptopall,mesage,jiter,iout_tcp,pflag) numgross=nint(awork(6,i_tcp)) numfailqc=nint(awork(21,i_tcp)) write(iout_tcp,925) 'psfc',numgross,numfailqc if(nump > izero)then pw=awork(4,i_tcp)/float(nump) pw3=awork(22,i_tcp)/float(nump) end if end if write(iout_tcp,950) 'psfc',jiter,nread,nkeep,nump write(iout_tcp,951) 'psfc',awork(4,i_tcp),awork(22,i_tcp),pw,pw3 close(iout_tcp) end if ! Summary report for lagrangian if (mype==mype_lag)then if(first)then open(iout_lag) else open(iout_lag,position='append') end if tmplty=zero; tqcplty=zero ; ntot=izero tt=zero ; qctt=zero nread=izero nkeep=izero do i=1,ndat if(dtype(i)== 'lag')then nread=nread+ndata(i,2) nkeep=nkeep+ndata(i,3) end if end do if(nkeep > izero)then mesage='current fit of lagangian data, ranges in m $' do j=1,nconvtype pflag(j)=trim(ioctype(j)) == 'lag' end do call dtast(bwork,npres_print,pbot,ptop,mesage,jiter,iout_lag,pflag) do k=1,nsig num(k)=nint(awork(6*nsig+k+100_i_kind,i_lag)) rat=zero ; rat3=zero if(num(k) > izero) then rat=awork(4*nsig+k+100_i_kind,i_lag)/float(num(k)) rat3=awork(3*nsig+k+100_i_kind,i_lag)/float(num(k)) end if ntot=ntot+num(k); tmplty=tmplty+awork(4*nsig+k+100_i_kind,i_lag) tqcplty=tqcplty+awork(3*nsig+k+100_i_kind,i_lag) write(iout_lag,240) 'lag',num(k),k,awork(4*nsig+k+100_i_kind,i_lag), & awork(3*nsig+k+100_i_kind,i_lag),rat,rat3 end do numgross=nint(awork(4,i_lag)) numfailqc=nint(awork(21,i_lag)) write(iout_lag,925) 'lag',numgross,numfailqc ! numlow = nint(awork(2,i_t)) ! numhgh = nint(awork(3,i_t)) ! write(iout_lag,900) 't',numhgh,numlow tt=tmplty/ntot qctt=tqcplty/ntot end if write(iout_lag,950) 'lag',jiter,nread,nkeep,ntot write(iout_lag,951) 'lag',tmplty,tqcplty,tt,qctt close(iout_lag) endif ! Format statements used above 111 format('obs lev num rms bias sumges sumobs cpen') 240 format(' num(',A1,') = ',i6,' at lev ',i4,' pen,qcpen,cpen,cqcpen = ',6(g11.5,1x)) 241 format(' num(',A1,') = ',i6,' at lev ',i4,' upen,vpen,cupen,cvpen = ',6(g11.5,1x)) 900 format(' number of ',a5,' obs extrapolated above',& ' top sigma layer=',i8,/,10x,' number extrapolated below',& ' bottom sigma layer=',i8) 905 format(' number of ',a5,' obs with station elevation > 2km = ',i8,/, & ' number with abs(guess topography-station elevation) > 200m = ',i8) 925 format(' number of ',a5,' obs that failed gross test = ',I5,' nonlin qc test = ',I5) 949 format(' number of ',a5,' obs = ',i6,' pen= ',e24.18,' cpen= ',g12.6) 950 format(' type ',a5,' jiter ',i3,' nread ',i7,' nkeep ',i7,' num ',i7) 951 format(' type ',a5,' pen= ',e24.18,' qcpen= ',e24.18,' r= ',g12.6,' qcr= ',g12.6) 952 format(t5,'it',t13,'sat',t21,'# read',t32,'# keep',t42,'# assim',& t52,'penalty',t67,'cpen') 975 format(' grsmlt=',f7.1,' number of bad ',a5,' obs=',f8.0) return end subroutine statsconv