PROGRAM amsr2lr IMPLICIT none !From Vera Gerald program of 18 Dec 2013 !F90 and other modifications R. Grumbine INTEGER iunt PARAMETER ( iunt = 12 ) CHARACTER(8) SUBSET CHARACTER(80) HS1, HS2, HS3, Htmp ! CHARACTER*(*) FLT CHARACTER(20) FLT PARAMETER ( FLT = 'bufr.tab' ) ! REAL HDR(9), XLOC(2,1), PRDLR(5,12), PRDHR(5,2) REAL(8) HDR(9), XLOC(2,1), PRDLR(5,12), PRDHR(5,2) ! REAL(8) HDR(9), XLOC(2,1), PRDLR(5,12), PRDHR(5,2) DATA HS1, HS2, HS3, Htmp & / 'SAID SIID YEAR MNTH DAYS HOUR MINU ORBN SLNM ', & 'CLATH CLONH ', & 'SCCF ALFR VIIRSQ ANPO TMBR ', & 'MTYP '/ ! INTEGER icount, iu, ntank, is REAL sccf, alfr, virsq, anpo, btmp DOUBLE PRECISION slon, slat INTEGER itk, lubfi, idate, iret, isad, issid, myr, mon, mda, mhr INTEGER minu, iorbn, iscanl, idatc INTEGER nchan, err INTEGER openout, outunit, creturns LOGICAL ok, hfok, lfok, weatherok !................................. OPEN ( UNIT = iunt, FILE = FLT ) CALL W3TAGB('DCODAMSR2',0095,0333,0077,'NP11 ') outunit = 52 OPEN (outunit, FORM="UNFORMATTED", STATUS="NEW") creturns = openout(outunit) icount = 0 iu = 11 ntank = 1 DO ITK = 1, NTANK !----------------------------------------------------------------------- ! DECODE BUFR FILE !----------------------------------------------------------------------- ! OPEN THE BUFR INPUT FILE ! ------------------------ LUBFI=IU CALL OPENBF(LUBFI,'IN',iunt) CALL READMG(LUBFI,SUBSET,IDATE,iret) IF(iret.NE.0) THEN CALL closbf(lubfi) IF (icount.lt.1) then WRITE(6, *) 'READMG ERROR W/ BUFR MSG:iret=', iret err = 23 STOP ! CALL exit(23) go to 99 ENDIF ENDIF ! PRINT*, 'READING DATA FOR ', IDATE ! ! READ A SUBSET - READ A MESSAGE WHEN NO SUBSETS - END WHEN NO MESSAGES ! --------------------------------------------------------------------- 150 CONTINUE CALL READSB(LUBFI, iret) IF(iret.NE.0) THEN CALL READMG(LUBFI, SUBSET, IDATE, iret) IF(iret.NE.0) GOTO 100 GOTO 150 ENDIF ! if(icount .LE. 5e8) go to 250 PRINT *, ' ***** NUMBER OF OBS EXCEEDED 5e8 *****' STOP 250 CONTINUE ! ----------------------------------------------------------------- ! CALL UFBINT, UFBREP TO GET THE DATA ! ! SAID SIID YEAR MNTH DAYS HOUR MINU ORBN SLNM ! CLATH CLONH ', ! SCCF ALFR VIIRSQ ANPO TMBR ! ----------------------------------------------- CALL UFBINT(LUBFI,HDR, 9, 1,iret, HS1) CALL UFBrep(LUBFI,XLOC,2, 1,iret, HS2) CALL UFBrep(LUBFI,PRDLR,5,12,iret, HS3) ! PRINT *,'iret hs3 = ',iret nchan = iret ! isad = NINT(hdr(1)) issid = NINT(hdr(2)) myr = NINT(hdr(3)) mon = NINT(hdr(4)) mda = NINT(hdr(5)) mhr = NINT(hdr(6)) minu = NINT(hdr(7)) iorbn = NINT(hdr(8)) iscanl = NINT(hdr(9)) ! SLAT = XLOC(1,1) SLON = XLOC(2,1) if(slon.lt.-180.0.or.slon.gt.180.0) go to 151 icount = icount + 1 ! Add some checking for quality and being in latitudes suitable for sea ! ice. 31 July 2017 Robert Grumbine ok = .TRUE. IF (nchan .EQ. 2) THEN ok = hfok(PRDLR) ELSE ok = lfok(PRDLR) ENDIF ok = ok .AND. (slat > 25.0 .OR. slat < -40.0) ! IF (ok) THEN CALL headerout(isad,myr,mon,mda,mhr,minu,slat,slon, nchan) DO is = 1,nchan sccf = PRDLR(1,is) alfr = PRDLR(2,is) virsq = PRDLR(3,is) anpo = PRDLR(4,is) btmp = PRDLR(5,is) CALL spotout(sccf, alfr, virsq, anpo, btmp) ENDDO ENDIF IF (MOD(icount,INT(1e5)) .EQ. 0) THEN WRITE(51,*)isad,myr,mon,mda,mhr,minu,slat,slon, nchan ENDIF 50 FORMAT(1x,i3,1x,i4,3i2,i3,1x,2f11.5,I3) ! WRITE(52+icount/5000000) sccf/1.e9,alfr,virsq,anpo,btmp ! IF (MOD(icount,INT(1e5)) .EQ. 0) THEN ! WRITE(51,54) sccf/1.e9,alfr,virsq,anpo,btmp ! ENDIF 54 FORMAT(1x,F6.2,1x,4F13.4) 151 continue ! get next report GO TO 150 100 CONTINUE ! CALL closbf(lubfi) IF (icount.lt.1) then WRITE(6,*) 'READMG ERROR W/ BUFR MSG:iret=',iret err = 23 ! CALL exit(23) STOP go to 99 endif ! WRITE(6,*) 'AFTER WOOLENS INTERFACE: icount= ',icount ! ENDDO ! !................. WRITE(6,692) icount 692 FORMAT(1x,'END OF JOB','count= ',i8) CALL W3TAGE('DCODNET3') 99 CONTINUE END !--------------------------------------------------------------- LOGICAL FUNCTION hfok(obs) IMPLICIT none REAL(8) obs(5,2) ! i = parameters, 5 = tmbr; j = channel LOGICAL tmp tmp = .TRUE. IF (obs(5,1) > 285 .OR. obs(5,2) > 285. ) tmp = .FALSE. IF (obs(5,1) > obs(5,2) ) tmp = .FALSE. hfok = tmp RETURN END LOGICAL FUNCTION lfok(obs) IMPLICIT none REAL(8) obs(5,12) ! i = parameters, 5 = tmbr; j = channel LOGICAL tmp INTEGER i tmp = .TRUE. DO i = 1, 6 IF (obs(5,2*i-1) > 285 .OR. obs(5,2*i) > 285 ) tmp = .FALSE. IF (obs(5,2*i-1) > obs(5,2*i) ) tmp = .FALSE. ENDDO lfok = tmp RETURN END LOGICAL FUNCTION weatherok(obs, lat) ! Weather filter. Should only apply after regressions to AMSRE have ! been performed. IMPLICIT none REAL AMSR_GR37LIM, AMSR_GR24LIM PARAMETER (AMSR_GR37LIM = 0.046) PARAMETER (AMSR_GR24LIM = 0.045) REAL(8) obs(5,12), amsre(5,12) REAL(8) lat LOGICAL tmp REAL gr3719, gr2419 amsre = obs CALL regress(amsre, lat) tmp = .FALSE. gr3719 = (amsre(5,12) - amsre(5,8) ) / (amsre(5,12)+amsre(5,8)) gr2419 = (amsre(5,10) - amsre(5,8) ) / (amsre(5,10)+amsre(5,8)) IF ((gr3719 < AMSR_GR37LIM) .AND. (gr2419 < AMSR_GR24LIM) ) THEN tmp = .TRUE. ENDIF weatherok = tmp RETURN END SUBROUTINE regress(amsre, lat) IMPLICIT none REAL(8) amsre(5,12) REAL(8) lat !Perform limited regression -- just those channels that are used in !weather filter IF (lat > 0) THEN amsre(5, 8) = amsre(5, 8)*1.031 - 9.710 amsre(5,10) = amsre(5,10)*0.999 - 1.706 amsre(5,12) = amsre(5,12)*0.997 - 2.610 ELSE amsre(5, 8) = amsre(5, 8)*1.032 -10.013 amsre(5,10) = amsre(5,10)*0.993 - 0.987 amsre(5,12) = amsre(5,12)*0.995 - 2.400 ENDIF RETURN END