C**** CPSIJ.FOR Color PostScript file of IxJ DataFiles 2006/05/31 C**** C**** Compile into Executable module: FCE CPSIJ.FOR C**** C**** Compilation notes: The directory /u/cmrun/ should contain the C**** files: CPIJ.I , COLORS , CPIJCBDT/* , and Z files. C**** /u/cmrun/ can be changed in this Fortran source listing. C**** The executable module a.out should be renamed to CPSIJ . C**** The PP command sends the PostScript file to a network color C**** printer. The printer's correct address should be part of the C**** argument to the CALL SYSTEM Fortran statment. C**** Parameter (MAXREC=1800, MAXIJR=72*46*MAXREC, IBACKC=0,IFOREC=1) C**** C**** MAXREC = dimension of TITLE array C**** MAXIJR = dimension of XYZD array C**** IBACKC = index of background color C**** IFOREC = index of foreground color C**** LOGICAL*4 QOPEN,QFOPEN,QPOPEN COMMON /PARMCB/ DATMIS,IM,JM,NRECM(4),QOPEN(3) DATA QFOPEN/.FALSE./, QPOPEN/.FALSE./ C**** C**** DATMIS = upper boundary for missing data C**** IM,JM = longitude and latitude horizontal resolution C**** NRECM = number of records read from data file so far C**** NRECM(4) = number of difference records calculated so far C**** QOPEN = whether data file is currently open C**** QFOPEN = whether CPSIJ.O and CPSIJ.IF have been opened C**** CHARACTER TITLE*80, JCBFN*8, JGRID*1,JLAND*1,JWHIT*1 COMMON /DATACB/ XYZD(MAXIJR,4) COMMON /TITLCB/ TITLE(MAXREC,4),CSCLJ(MAXREC,4),DSCLJ(MAXREC,4), * JCBDT(MAXREC,4),JROTA(MAXREC,4),JCBFN(MAXREC,4), * JGRID(MAXREC,4),JLAND(MAXREC,4),JWHIT(MAXREC,4) C**** C**** XYZD = longitude by latitude data read from files C**** TITLE = title describing each record C**** CSCLJ = center value of color scale C**** DSCLJ = scale difference between scale markers C**** JCBDT = index pointing to color bar distribution table C**** JROTA = color plot is rotated JROTA boxes to the left C**** JCBFN = file name of color bar distribution table C**** JGRID = color plot displays specified grid: A, B, U or V C**** JLAND = color plot boxes show: A = All, L = land, O = ocean C**** JWHIT = continental outline or unused areas: W = white, B = black C**** INTEGER*4 NREC(4), KMIN(3),KMAX(3) LOGICAL*4 QEXIST CHARACTER ARG*80, FILEIN*80, ZFILE*80, * CPIJCB*8, UNNN*4,UNN2*4,UNN3*4, CCSCL*8,CDSCL*8 COMMON /FIXDCB/ FOCEAN(720,360) EQUIVALENCE (NRX,NREC(1)),(NRY,NREC(2)),(NRZ,NREC(3)) DATA NREC/4*1/ C**** C**** FOCEAN = ocean fraction for each horizontal grid box C**** NREC = most recently displayed record of file C**** IARGS = IARGC() IF(IARGS.le.0) GO TO 800 C**** Read in CPIJ.I, COLORS and CPIJCBDT/ST12 CALL MATCH0 CALL COLOR CALL LOADCB ('SCT0 ',ICB,*750) C**** C**** Determine resolution C**** FOCEAN(1,1) = -1. IARG = 1 CALL GETARG (IARG,ARG) IF(ARG(1:1).lt.'0' .or. ARG(1:1).gt.'9') GO TO 50 LOCAX = SCAN(ARG,'Xx') IF(LOCAX.le.1) GO TO 50 READ (ARG(1:LOCAX-1) ,*,IOSTAT=IOSQ) IM IF(IOSQ.ne.0) GO TO 50 READ (ARG(LOCAX+1:80),*,IOSTAT=IOSQ) JM IF(IOSQ.ne.0) GO TO 50 IARG = 2 GO TO 200 C**** Resolution not specified on command line, derive it from file. 50 CALL RESOLU (ARG, IM,JM,ZFILE, *810) IF(ZFILE(1:4).eq.'None') GO TO 200 C**** C**** Read ocean fraction file C**** INQUIRE (FILE=ZFILE,EXIST=QEXIST) IF(.not.QEXIST) then CALL GETENV ('IFDIR',ARG) ZFILE = TRIM(ARG) // '/' // ZFILE ; endif OPEN (11,FILE=ZFILE,FORM='UNFORMATTED',STATUS='OLD',ERR=811) READ (11) ARG,(FOCEAN(I,1),I=1,IM*JM) CLOSE (11) WRITE (6,*) TRIM(ZFILE) // ' read in: ' // ARG(1:48) C**** C**** Open input files C**** C**** Open X file 200 CALL GETARG (IARG,ARG) OPEN (1,FILE=ARG,FORM='UNFORMATTED',STATUS='OLD',ERR=810) QOPEN(1) = .TRUE. NREC(1) = 1 NRX = 1 CALL READIN (1,NRX,*999) IMJMNR = 1 + IM*JM*(NRX-1) CALL MATCHI (IM*JM,XYZD(IMJMNR,1),TITLE(NRX,1),CSCLJ(NRX,1), * DSCLJ(NRX,1),JCBDT(NRX,1),JROTA(NRX,1),JCBFN(NRX,1), * JGRID(NRX,1),JLAND(NRX,1),JWHIT(NRX,1), DATMIS) IARG = IARG + 1 IF(IARG.gt.IARGS) GO TO 250 C**** Open Y file CALL GETARG (IARG,ARG) OPEN (2,FILE=ARG,FORM='UNFORMATTED',STATUS='OLD',ERR=810) QOPEN(2) = .TRUE. NREC(2) = 1 NRY = 1 CALL READIN (2,NRY,*999) IMJMNR = 1 + IM*JM*(NRY-1) CALL MATCHI (IM*JM,XYZD(IMJMNR,2),TITLE(NRY,2),CSCLJ(NRY,2), * DSCLJ(NRY,2),JCBDT(NRY,2),JROTA(NRY,2),JCBFN(NRY,2), * JGRID(NRY,2),JLAND(NRY,2),JWHIT(NRY,2), DATMIS) IARG = IARG + 1 IF(IARG.gt.IARGS) GO TO 250 C**** Open Z file CALL GETARG (IARG,ARG) OPEN (3,FILE=ARG,FORM='UNFORMATTED',STATUS='OLD',ERR=810) QOPEN(3) = .TRUE. NREC(3) = 1 NRZ = 1 CALL READIN (3,NRZ,*999) IMJMNR = 1 + IM*JM*(NRZ-1) CALL MATCHI (IM*JM,XYZD(IMJMNR,3),TITLE(NRZ,3),CSCLJ(NRZ,3), * DSCLJ(NRZ,3),JCBDT(NRZ,3),JROTA(NRZ,3),JCBFN(NRZ,3), * JGRID(NRZ,3),JLAND(NRZ,3),JWHIT(NRZ,3), DATMIS) C**** 250 NR = NRX IU = 1 WRITE (UNNN,923) 'X',NR C**** C**** Display color plot C**** 300 IMJMNR = 1 + IM*JM*(NR-1) CALL CPIJ (IM,JM, FOCEAN,XYZD(IMJMNR,IU),TITLE(NR,IU), * CSCLJ(NR,IU),DSCLJ(NR,IU),JCBDT(NR,IU),JROTA(NR,IU), * JGRID(NR,IU),JLAND(NR,IU),JWHIT(NR,IU), UNNN, DATMIS) C**** C**** Receive keyboard input C**** 400 READ (5,940,END=750) ARG IF(ARG(1:1).eq.'Q' .or. ARG(1:1).eq.'q') GO TO 750 IF(ARG(1:1).eq.'A' .or. ARG(1:1).eq.'a') GO TO 410 IF(ARG(1:1).eq.'L' .or. ARG(1:1).eq.'l') GO TO 420 IF(ARG(1:1).eq.'O' .or. ARG(1:1).eq.'o') GO TO 430 IF(ARG(1:1).eq.'W' .or. ARG(1:1).eq.'w') GO TO 440 IF(ARG(1:1).eq.'B' .or. ARG(1:1).eq.'b') GO TO 450 IF(ARG(1:1).eq.'G' .or. ARG(1:1).eq.'g') GO TO 460 IF(ARG(1:1).eq.'R' .or. ARG(1:1).eq.'r') GO TO 470 IF(ARG(1:1).eq.'C' .or. ARG(1:1).eq.'c') GO TO 480 IF(ARG(1:1).eq.'S' .or. ARG(1:1).eq.'s') GO TO 490 IF(ARG(1:1).eq.'M' .or. ARG(1:1).eq.'m') GO TO 500 IF(ARG(1:1).eq.'X' .or. ARG(1:1).eq.'x') GO TO 510 IF(ARG(1:1).eq.'Y' .or. ARG(1:1).eq.'y') GO TO 510 IF(ARG(1:1).eq.'Z' .or. ARG(1:1).eq.'z') GO TO 510 IF(ARG(1:1).eq.'D' .or. ARG(1:1).eq.'d') GO TO 550 IF(ARG(1:1).eq.'E' .or. ARG(1:1).eq.'e') GO TO 560 IF(ARG(1:1).eq.'T' .or. ARG(1:1).eq.'t') GO TO 590 IF(ARG(1:1).eq.'N' .or. ARG(1:1).eq.'n') GO TO 600 IF(ARG(1:1).eq.'F' .or. ARG(1:1).eq.'f') GO TO 640 IF(ARG(1:1).eq.'P' .or. ARG(1:1).eq.'p') GO TO 650 GO TO 400 C**** A: display All grid boxes 410 JLAND(NR,IU) = 'A' GO TO 300 C**** L: display only Land grid boxes 420 JLAND(NR,IU) = 'L' GO TO 300 C**** O: display only Ocean grid boxes 430 JLAND(NR,IU) = 'O' GO TO 300 C**** W: display continental outline or unused boxes in White 440 JWHIT(NR,IU) = 'W' GO TO 400 C**** B: display continental outline or unused boxes in Black 450 JWHIT(NR,IU) = 'B' GO TO 400 C**** G: display data on different Grid 460 IF(ARG(2:2).eq.'A' .or. ARG(2:2).eq.'a') JGRID(NR,IU) = 'A' IF(ARG(2:2).eq.'B' .or. ARG(2:2).eq.'b') JGRID(NR,IU) = 'B' IF(ARG(2:2).eq.'U' .or. ARG(2:2).eq.'u') JGRID(NR,IU) = 'U' IF(ARG(2:2).eq.'V' .or. ARG(2:2).eq.'v') JGRID(NR,IU) = 'V' GO TO 300 C**** R: Rotate color plot I grid boxes to the left 470 READ (ARG(2:10),*,IOSTAT=IOSQ) I IF(IOSQ.eq.0) JROTA(NR,IU) = MODULO (I+JROTA(NR,IU),IM) GO TO 300 C**** C**** C: use new Color bar distribution table C**** 480 CALL PARSE (ARG(2:80),79,1,KMIN,KMAX) IF(KMAX(1).le.0) GO TO 400 IF(KMAX(1).gt.KMIN(1)+7) KMAX(1) = KMIN(1)+7 CPIJCB = ARG(1+KMIN(1):1+KMAX(1)) CALL LOADCB (CPIJCB,ICBDT,*400) JCBFN(NR,IU) = CPIJCB JCBDT(NR,IU) = ICBDT GO TO 400 C**** C**** S: use new scaling parameters C**** 490 CALL PARSE (ARG(2:80),79,2,KMIN,KMAX) IF(KMIN(1).le.0) GO TO 491 READ (ARG(1+KMIN(1):1+KMAX(1)),*,IOSTAT=IOSQ) CSCLK IF(IOSQ.ne.0) GO TO 491 CSCLJ(NR,IU) = NINT(4.*CSCLK/DSCLJ(NR,IU))*DSCLJ(NR,IU)/4. IF(KMIN(2).le.0) GO TO 300 READ (ARG(1+KMIN(2):1+KMAX(2)),*,IOSTAT=IOSQ) DSCLK IF(IOSQ.ne.0 .or. DSCLK.eq.0.) GO TO 491 CSCLJ(NR,IU) = NINT(4.*CSCLK/DSCLK)*DSCLK/4. DSCLJ(NR,IU) = DSCLK GO TO 300 491 WRITE (6,949) CSCLJ(NR,IU),DSCLJ(NR,IU) GO TO 400 C**** C**** M: reset upper boundary of Missing data C**** 500 READ (ARG(2:80),*,IOSTAT=IOSQ) DATMIS IF(IOSQ.ne.0) GO TO 400 WRITE (6,950) DATMIS GO TO 300 C**** C**** Display another record C**** 510 CALL NEWREC (ARG,6, NREC, IU,NR,UNNN, *400) NREC(IU) = NR GO TO 300 C**** C**** D: display a difference record or calculate a different record C**** 550 IF(SCAN(ARG,'=').le.0) GO TO 510 C**** Calculate different record 551 ND = NRECM(4) + 1 IF(IM*JM*ND.gt.MAXIJR) GO TO 552 CALL DIFFER (ARG,IM*JM,NREC,ND,XYZD(1,4),DATMIS, *400) NR = ND IU = 4 NRECM(4) = NR NREC(4) = NR IMJMNR = 1 + IM*JM*(NR-1) CALL ESTSPP (IM*JM, XYZD(IMJMNR,4),CSCLJ(NR,4),DSCLJ(NR,4), * JCBDT(NR,4),JROTA(NR,4),JCBFN(NR,4),JGRID(NR,4), * JLAND(NR,4),JWHIT(NR,4), DATMIS) WRITE (UNNN,923) 'D',NR WRITE (6,930) UNNN // ARG(2:69) GO TO 300 C**** New D record exceeds internal dimension set by MAXIJR 552 WRITE (6,930) *'New D record exceeds internal dimension set by MAXIJR' GO TO 400 C**** C**** E: Estimate screen presentation parameters C**** 560 IMJMNR = 1 + IM*JM*(NR-1) + MAXIJR*(IU-1) CALL ESTSPP (IM*JM, XYZD(IMJMNR, 1),CSCLJ(NR,IU),DSCLJ(NR,IU), * JCBDT(NR,IU),JROTA(NR,IU),JCBFN(NR,IU),JGRID(NR,IU), * JLAND(NR,IU),JWHIT(NR,IU), DATMIS) GO TO 300 C**** C**** T: edit current Title C**** 590 READ (5,940,END=750) ARG DO 591 K=1,80 591 IF(ARG(K:K).ne.'"') TITLE(NR,IU)(K:K) = ARG(K:K) GO TO 300 C**** C**** N: load New data file C**** 600 KU = 0 IF(ARG(2:2).eq.'X' .or. ARG(2:2).eq.'x') KU = 1 IF(ARG(2:2).eq.'Y' .or. ARG(2:2).eq.'y') KU = 2 IF(ARG(2:2).eq.'Z' .or. ARG(2:2).eq.'z') KU = 3 IF(KU.le.0) GO TO 400 IF(QOPEN(KU)) CLOSE (KU) QOPEN(KU) = .FALSE. CALL PARSE (ARG(3:80),78,2,KMIN,KMAX) CALL FILEEX (ARG(2+KMIN(1):2+KMAX(1)),KMAX(1)-KMIN(1)+1,FILEIN) OPEN (KU,FILE=FILEIN,FORM='UNFORMATTED',STATUS='OLD',ERR=603) QOPEN(KU) = .TRUE. DO 601 N=1,NRECM(KU) 601 JCBDT(N,KU) = 0 NRECM(KU) = 0 KN = 1 IF(KMAX(2).gt.0) READ (ARG(2+KMIN(2):2+KMAX(2)),*,IOSTAT=IOSQ) KN 602 CALL READIN (KU,KN,*400) IMJMKN = 1 + IM*JM*(KN-1) CALL MATCHI (IM*JM,XYZD(IMJMKN,KU),TITLE(KN,KU),CSCLJ(KN,KU), * DSCLJ(KN,KU),JCBDT(KN,KU),JROTA(KN,KU),JCBFN(KN,KU), * JGRID(KN,KU),JLAND(KN,KU),JWHIT(KN,KU), DATMIS) IU = KU NR = KN NREC(IU) = NR WRITE (UNNN,923) CHAR(Z'57'+IU),MOD(NR,1000) IF(NR.gt.999) UNNN(2:2) = CHAR(Z'37'+NR/100) GO TO 300 603 WRITE (6,930) 'Unable to open file: ' // FILEIN(1:56) GO TO 400 C**** C**** F: write Files CPSIJ.O and CPSIJ.IF of title, data and parameters C**** 640 IF(QFOPEN) GO TO 641 OPEN (7,FILE='CPSIJ.O',FORM='UNFORMATTED') OPEN (8,FILE='CPSIJ.IF',FORM='FORMATTED') 641 WRITE (7) TITLE(NR,IU),(XYZD(I,IU),I=1+IM*JM*(NR-1),IM*JM*NR) WRITE (CCSCL,960) NINT(CSCLJ(NR,IU)) WRITE (CDSCL,960) NINT(DSCLJ(NR,IU)) IF(NINT(CSCLJ(NR,IU)*10.).ne.NINT(CSCLJ(NR,IU))*10 .and. * ABS(CSCLJ(NR,IU)).le.99999.9) then WRITE (CCSCL,961) CSCLJ(NR,IU) WRITE (CDSCL,961) DSCLJ(NR,IU) ; endif IF(NINT(CSCLJ(NR,IU)*100.).ne.NINT(CSCLJ(NR,IU)*10.)*10 .and. * ABS(CSCLJ(NR,IU)).le.9999.99) then WRITE (CCSCL,962) CSCLJ(NR,IU) WRITE (CDSCL,962) DSCLJ(NR,IU) ; endif IF(NINT(CSCLJ(NR,IU)*1000.).ne.NINT(CSCLJ(NR,IU)*100.)*10 .and. * ABS(CSCLJ(NR,IU)).le.999.999) then WRITE (CCSCL,963) CSCLJ(NR,IU) WRITE (CDSCL,963) DSCLJ(NR,IU) ; endif WRITE (8,964) TITLE(NR,IU)(1:32),CCSCL,CDSCL,JROTA(NR,IU), * JGRID(NR,IU),JLAND(NR,IU),JWHIT(NR,IU),JCBFN(NR,IU) WRITE (6,930) * 'Files CPSIJ.O and CPSIJ.IF written for record: ' // UNNN GO TO 400 C**** C**** P: Print 1 or 2 records on color printer C**** 650 CALL PARSE (ARG(3:80),78,2,KMIN,KMAX) CALL NEWREC (ARG(2+KMIN(1):80),6, NREC, KU,KN,UNN2, *400) IF(.not.QPOPEN) OPEN (9,FILE='CPSIJ.PS') QPOPEN = .TRUE. IF(KMIN(2).le.0) WRITE (9,965) IF(KMIN(2).gt.0) WRITE (9,966) IMJMKN = 1 + IM*JM*(KN-1) CALL WRITPS (IM,JM,FOCEAN,XYZD(IMJMKN,KU),TITLE(KN,KU), * CSCLJ(KN,KU),DSCLJ(KN,KU),JCBDT(KN,KU),JROTA(KN,KU), * JGRID(KN,KU),JLAND(KN,KU),JWHIT(KN,KU), DATMIS) UNN3 = ' ' IF(KMIN(2).le.0) GO TO 651 CALL NEWREC (ARG(2+KMIN(2):80),6, NREC, KU,KN,UNN3, *653) WRITE (9,967) IMJMKN = 1 + IM*JM*(KN-1) CALL WRITPS (IM,JM,FOCEAN,XYZD(IMJMKN,KU),TITLE(KN,KU), * CSCLJ(KN,KU),DSCLJ(KN,KU),JCBDT(KN,KU),JROTA(KN,KU), * JGRID(KN,KU),JLAND(KN,KU),JWHIT(KN,KU), DATMIS) 651 WRITE (9,968) IF(ARG(2:2).eq.'P' .or. ARG(2:2).eq.'p') GO TO 652 WRITE (6,930) 'PostScript file CPSIJ.PS written containing: ' * // UNN2 // ' ' // UNN3 GO TO 400 652 CLOSE (9) QPOPEN = .FALSE. CALL SYSTEM ('lp -d color6 CPSIJ.PS') WRITE (6,930) 'PostScript file CPSIJ.PS written and printed ' // * 'containing: ' // UNN2 // ' ' // UNN3 GO TO 400 653 WRITE (9,967) GO TO 400 C**** C**** Q: Quit CPSIJ program C**** 750 IF(.not.QFOPEN) GO TO 751 CLOSE (7) CLOSE (8) 751 GO TO 999 C**** 800 WRITE (6,*) 'Usage: CPSIJ [IMxJM] filex [filey [filez]] ', * 'keyboard IxJ data 2006/05/31' GO TO 999 810 WRITE (0,*) ' Unable to open file: ',ARG(1:56) STOP 810 811 WRITE (0,*) ' Ocean fraction file not found: ',ZFILE(1:20) STOP 811 C**** 923 FORMAT (A1,I3.3) 930 FORMAT (3X,A) 940 FORMAT (A) 949 FORMAT (' Current values:',F12.3,' = center value', * F12.3,' = marker difference') 950 FORMAT (' New upper boundary for missing data:',F12.2) 960 FORMAT (I8) 961 FORMAT (F8.1) 962 FORMAT (F8.2) 963 FORMAT (F8.3) 964 FORMAT (A32,2A8,I4,3(2X,A1),2X,A8) 965 FORMAT ('%!PS-Adobe-3.0 PostScript file CPSIJ.PS 2006/05/31'/ * '%%Orientation: Landscape' / * 'gsave 90 rotate 108 -306 translate' /) 966 FORMAT ('%!PS-Adobe-3.0 PostScript file CPSIJ.PS 2006/05/31'/ * '%%Orientation: Portrait' / * 'gsave .75 .75 scale 120 774 translate' /) 967 FORMAT (/'%%% Second record printed on bottom half of page' / * ' 0 -496 translate' /) 968 FORMAT (/ 'grestore' / 'showpage' /) 999 END SUBROUTINE READIN (IU,NR,*) C**** C**** READIN reads data from unit IU up to record NR C**** C**** NRECM = number of records read from file so far C**** Parameter (MAXREC=1800, MAXIJR=72*46*MAXREC) LOGICAL*4 QOPEN CHARACTER TITLE*80 COMMON /PARMCB/ DATMIS,IM,JM,NRECM(4),QOPEN(3) COMMON /DATACB/ XYZD(MAXIJR,4) COMMON /TITLCB/ TITLE(MAXREC,4) C**** IF(NR.le.NRECM(IU)) RETURN IF(.not.QOPEN(IU)) GO TO 800 IF(NR.gt.MAXREC) GO TO 810 IF(IM*JM*NR.gt.MAXIJR) GO TO 820 DO 10 N=NRECM(IU)+1,NR IMIN = IM*JM*(N-1) + 1 IMAX = IM*JM*N READ (IU,IOSTAT=IOSQ) TITLE(N,IU),(XYZD(I,IU),I=IMIN,IMAX) 10 IF(IOSQ.ne.0) GO TO 830 NRECM(IU) = NR RETURN C**** 800 WRITE (6,980) CHAR(Z'57'+IU),NRECM(IU) RETURN 1 810 WRITE (6,981) 'Requested record exceeds internal dimens' // * 'ion of the parameter MAXREC = 1800.' RETURN 1 820 WRITE (6,981) 'Requested record exceeds internal dimens' // * 'ion of the parameter MAXIJR.' RETURN 1 830 WRITE (6,983) N,CHAR(Z'57'+IU) CLOSE(IU) NRECM(IU) = N-1 QOPEN(IU) = .FALSE. RETURN 1 C**** 980 FORMAT (' File',A2,' was closed.',I6,' records were read in.') 981 FORMAT (3X,A) 983 FORMAT (' Error encountered reading record',I4,' of file',A2, * '. File is closed.') END SUBROUTINE MATCH0 C**** C**** MATCH0 loads data from the file CPIJ.I into CPIJoI common block C**** Parameter (MAXoI=1024) INTEGER*4 IJCBDT(10) LOGICAL*4 QEXIST REAL*4 XYZD(*), MEAN, CBDTK(9) CHARACTER TITLE*80, JCBFN*8, JGRID*1,JLAND*1,JWHIT*1, * TITLEI*32, ICBFN*8, IGRID*1,ILAND*1,IWHIT*1, FILEIN*80 COMMON /FIXDCB/ FOCEAN(720,360) COMMON /CPIJoI/ TITLEI(MAXoI),CSCLI(MAXoI),DSCLI(MAXoI), * IROTA(MAXoI),ICBFN(MAXoI),IGRID(MAXoI),ILAND(MAXoI), * IWHIT(MAXoI) C**** C**** TITLEI= matches first 32 characters of data record TITLE C**** CSCLI = center value of color scale C**** DSCLI = scale difference between scale markers C**** IROTA = color plot is rotated IROTA boxes to the left C**** ICBFN = file name of color bar distribution table C**** IGRID = color plot displays specified grid: A, B, U or V C**** ILAND = color plot boxes show: A = All, L = land, O = ocean C**** IWHIT = continental outline or unused areas: W = white, B = black C**** DATA IIMAX /0/ DATA CBDTK /-4.5,-4.,-3.,-2.,0.,2.,3.,4.,4.5/ C**** C**** IIMAX = number of lines read from the file CPIJ.I C**** C**** Read in CPIJ.I C**** C**** Locate file FILEIN = 'CPIJ.I' INQUIRE (FILE=FILEIN,EXIST=QEXIST) IF(QEXIST) GO TO 10 CALL GETENV ('IFDIR',FILEIN) FILEIN = TRIM(FILEIN) // '/CPIJ.I' INQUIRE (FILE=FILEIN,EXIST=QEXIST) IF(.not.QEXIST) GO TO 80 C**** Open and read file 10 OPEN (4,FILE=FILEIN,STATUS='OLD',ERR=90) WRITE (6,*) 'File CPIJ.I read in.' 20 DO 30 I=1,4 30 READ (4,904) DO 40 I=1,MAXoI 40 READ (4,904,END=50,ERR=50) TITLEI(I),CSCLI(I),DSCLI(I),IROTA(I), * IGRID(I),ILAND(I),IWHIT(I),ICBFN(I) 50 IIMAX = I-1 CLOSE (4) WRITE (6,*) TRIM(FILEIN),' contains',IIMAX,' lines of data.' IF(IIMAX.le.0) STOP 'MATCH0 50' RETURN C**** 80 WRITE (6,*) 'Neither CPIJ.I nor /u/cmrun/CPIJ.I was found.' RETURN 90 WRITE (6,*) 'Error opening: ',TRIM(FILEIN) RETURN C**** C**** ENTRY MATCHI (IMxJM, XYZD,TITLE,CSCLJ,DSCLJ,JCBDT,JROTA,JCBFN, * JGRID,JLAND,JWHIT, DATMIS) C**** C**** MATCHI matches the data record TITLE with TITLEIs in the file C**** CPIJ.I to get initial parameters for the color plot. C**** If the data record TITLE does not match any CPIJ.I TITLEI, C**** then the initial parameters are estimated. C**** IF(JCBDT.gt.0) RETURN DO 130 L=1,IIMAX DO 110 K=1,32 110 IF(TITLEI(L)(K:K).ne.TITLE(K:K) .and. TITLEI(L)(K:K).ne.'.') * GO TO 130 IF(ICBFN(L).eq.' ') GO TO 200 C**** TITLEI(L) matches TITLE CSCLJ = CSCLI(L) DSCLJ = DSCLI(L) JROTA = IROTA(L) JCBFN = ICBFN(L) JGRID = IGRID(L) JLAND = ILAND(L) JWHIT = IWHIT(L) CALL LOADCB (JCBFN,JCBDT,*120) IF(FOCEAN(1,1).eq.-1.) JLAND = 'A' RETURN C**** Desired color bar distribution table not found 120 JCBFN = 'SCT0 ' JCBDT = 1 RETURN 130 continue GO TO 200 C**** C**** C**** ENTRY ESTSPP (IMxJM, XYZD,CSCLJ,DSCLJ,JCBDT,JROTA,JCBFN,JGRID, * JLAND,JWHIT, DATMIS) C**** C**** ESTSPP estimates the screen presentation parmeters C**** 200 JGRID = 'A' JROTA = 0 JLAND = 'A' JWHIT = 'B' C**** Calculate mean and standard deviation IJXYZD = 0 SXYZD = 0. SXYZDQ = 0. DO 210 IJ=1,IMxJM IF(XYZD(IJ).le.DATMIS) GO TO 210 IJXYZD = IJXYZD + 1 SXYZD = SXYZD + XYZD(IJ) SXYZDQ = SXYZDQ + XYZD(IJ)**2 210 continue CSCLJ = SXYZD/IJXYZD STDV = SQRT(SXYZDQ/IJXYZD - CSCLJ**2) DSCLJ = ROUND (STDV/2.5) CSCLJ = NINT(4*CSCLJ/DSCLJ)*DSCLJ/4. C**** Determine Color Block Definition Table from difference between C**** median and mean DO 220 K=1,10 220 IJCBDT(K) = 0 DO 250 IJ=1,IMxJM IF(XYZD(IJ).le.DATMIS) GO TO 250 DO 230 K=1,9 230 IF(XYZD(IJ) .le. MEAN+DSCLJ*CBDTK(K)) GO TO 240 240 IJCBDT(K) = IJCBDT(K) + 1 250 continue SUMIJ = .5*IJCBDT(1) DO 260 K=1,8 SUMIJ = SUMIJ + .5*(IJCBDT(K)+IJCBDT(K+1)) 260 IF(SUMIJ.ge..5*IJXYZD) GO TO 270 C**** Load one of the Standard Color Block Definition Table 270 WRITE (JCBFN,927) K CALL LOADCB (JCBFN,JCBDT,*280) 280 RETURN C**** 904 FORMAT (A32,2F8.0,I4,3(2X,A1),2X,A8) 927 FORMAT ('SCT',I1) END SUBROUTINE COLOR C**** C**** COLOR reads the file COLORS containing the RGB color values C**** that are used for the PostScript output file. C**** Parameter (MAXC=63) CHARACTER FILEIN*80 COMMON /COLORS/ REDP(0:MAXC),GRNP(0:MAXC),BLUP(0:MAXC),NMCOLR C**** Open the file COLORS FILEIN = 'COLORS' OPEN (4,FILE=FILEIN,STATUS='OLD',ERR=10) GO TO 20 10 CALL GETENV ('IFDIR',FILEIN) FILEIN = TRIM(FILEIN) // '/' // 'COLORS' OPEN (4,FILE=FILEIN,STATUS='OLD',ERR=801) C**** Read in RGB color values 20 DO 30 N=1,5 30 READ (4,904) DO 40 N=0,MAXC READ (4,904,IOSTAT=IOSQ) REDS,GRNS,BLUS, REDP(N),GRNP(N),BLUP(N) 40 IF(IOSQ.ne.0) GO TO 50 50 NMCOLR = N-1 WRITE (6,906) NMCOLR,FILEIN RETURN C**** 801 WRITE (0,*) ' Unable to locate COLORS' STOP 801 C**** 904 FORMAT (3X,3F7.3,1X,3F7.3) 906 FORMAT (' Last line of RGB color values is:',I3, * '. Read from: ',A20) END SUBROUTINE LOADCB (CPIJCB,ICBDT,*) C**** C**** LOADCB determines whether the current color bar distribution table C**** has been loaded. If it has not, then LOADCB reads the color bar C**** distribution table and loads its contents into the next unused row C**** of KCBDT. In all cases, LOADCB returns the index ICBDT which C**** points to the current color bar distribution table. C**** C**** CPIJCB = requested color bar distribution table C**** ICBDT = index pointing to color bar distribution table C**** Parameter (MAXBLK=64, MAXCB=64) C**** C**** MAXBLK = maximum number of color blocks in any color bar table C**** MAXCB = maximum number of color bar distribution tables C**** CHARACTER CPIJCB*8, CBDTFN*8, FILEIN*80 COMMON /CBDTCB/ KCBDT(MAXBLK,MAXCB),NBLKS(MAXCB),CBDTFN(MAXCB) DATA MCBDT /0/ C**** C**** KCBDT = contents of color bar distribution tables C**** NBLKS = number of blocks in each color bar distribution table C**** CBDTFN = file names of color bar distribution tables C**** MCBDT = number of color bar distribution tables loaded so far C**** DO 10 J=1,MCBDT ICBDT = J 10 IF(CPIJCB.eq.CBDTFN(J)) RETURN C**** C**** Read in new color bar distribution table from disk C**** IF(MCBDT.ge.MAXCB) GO TO 800 MCBDT = MCBDT + 1 FILEIN = 'CPIJCBDT/' // CPIJCB OPEN (4,FILE=FILEIN,STATUS='OLD',ERR=110) GO TO 120 110 CALL GETENV ('IFDIR',FILEIN) FILEIN = TRIM(FILEIN) // '/CPIJCBDT/' // CPIJCB OPEN (4,FILE=FILEIN,STATUS='OLD',ERR=810) 120 DO 130 K=1,3 130 READ (4,*,ERR=820) C**** Contents of file is in reverse order DO 140 K=1,MAXBLK 140 KCBDT(K,MCBDT) = -1 READ (4,*,END=150,ERR=150) (KCBDT(K,MCBDT),K=MAXBLK,1,-1) 150 CLOSE (4) C**** Move loaded data into memory KCBDT(1,MCBDT) to KCBDT(NBLKS,MCBDT) KOFF = 0 DO 160 K=1,MAXBLK IF(KCBDT(K,MCBDT).ge.0) then KCBDT(K-KOFF,MCBDT) = KCBDT(K,MCBDT) ELSE KOFF = KOFF + 1 ; endif 160 continue NBLKS(MCBDT) = MAXBLK-KOFF CBDTFN(MCBDT) = CPIJCB ICBDT = MCBDT RETURN C**** 800 WRITE (6,980) 'Number of color bar distribution tables ' // * 'exceeds internal dimension MAXCB.' RETURN 1 810 WRITE (6,980) 'Neither CPIJCBDT/' // TRIM(CPIJCB) // ' nor ' * // TRIM(FILEIN) // ' was found.' RETURN 1 820 WRITE (6,980) 'Color bar distribution table is empty: ' // * FILEIN(1:38) RETURN 1 980 FORMAT (3X,A) END SUBROUTINE CPIJ (IM,JM,FOCEAN,XYZD,TITLE,CSCLJ,DSCLJ,JCBDT,JROTA, * JGRID,JLAND,JWHIT, UNNN, DATMIS) C**** C**** CPIJ write data minimum and maximum and scale values to unit 6 C**** C**** IM,JM = horizontal resolution C**** FOCEAN= ocean fraction array C**** XYZD = grid point data to be plotted C**** TITLE = title at top of plot C**** CSCLJ = center value of color scale C**** DSCLJ = scale difference between scale markers C**** JCBDT = index pointing to color bar distribution table C**** JGRID = type of grid, all data is assumed to be on A grid C**** JLAND = color plot boxes show: A = all, L = land, O = ocean C**** UNNN = character describing data file and record number C**** DATMIS= upper boundary for missing data C**** Parameter (MAXBLK=64, MAXCB=64) REAL*4 FOCEAN(IM,JM),XYZD(IM,JM) CHARACTER TITLE*80, JGRID*1,JLAND*1,JWHIT*1, UNNN*4 COMMON /CBDTCB/ KCBDT(MAXBLK,MAXCB),NBLKS(MAXCB) C**** Calculate data minimum and maximum XYZMIN = 1.d30 XYZMAX = -1.d30 DO 110 I=1,IM*JM IF((JLAND.eq.'O' .and. FOCEAN(I,1).le..5) .or. * (JLAND.eq.'L' .and. FOCEAN(I,1).gt..5)) GO TO 110 IF(XYZD(I,1).le.DATMIS) GO TO 110 IF(XYZD(I,1).lt.XYZMIN) XYZMIN = XYZD(I,1) IF(XYZD(I,1).gt.XYZMAX) XYZMAX = XYZD(I,1) 110 continue C**** Calculate scale minimum and maximum MBLK = NBLKS(JCBDT) BCEN = MBLK/2. N = NINT (4.*CSCLJ/DSCLJ) N4 = MODULO (N,4) B1 = BCEN - N4 - 4*INT((BCEN-N4-.5)/4.) BM = B1 + 4*INT((MBLK-.5-B1)/4.) VALMIN = CSCLJ + (B1-BCEN)*DSCLJ/4. VALMAX = CSCLJ + (BM-BCEN)*DSCLJ/4. C**** Display TITLE above the data WRITE (6,931) TITLE(1:79) WRITE (6,932) XYZMIN,XYZMAX,VALMIN,VALMAX,UNNN RETURN C**** 931 FORMAT (A) 932 FORMAT (' Data:',F8.2,' to',F8.2, * ' Scale:',F8.2,' to',F8.2,6X,A4) END SUBROUTINE WRITPS (IM,JM,FOCEAN,XYZD,TITLE,CSCLJ,DSCLJ,JCBDT, * JROTA,JGRID,JLAND,JWHIT, DATMIS) C**** C**** Write PostScript file for record KN from unit KU. C**** C**** IM,JM = horizontal resolution C**** FOCEAN= ocean fration array C**** XYZD = grid point data to be plotted C**** TITLE = title at top of plot C**** CSCLJ = center value of color scale C**** DSCLJ = scale difference between scale markers C**** JCBDT = index pointing to color bar distribution table C**** JROTA = color plot is rotated JROTA boxes to the left C**** JGRID = type of grid, all data is assumed to be on A grid C**** JLAND = color plot boxes show: A = all, L = land, O = ocean C**** JWHIT = continental outline or unused areas: W = white, B = black C**** UNNN = character describing data file and record number C**** DATMIS= upper boundary for missing data C**** Parameter (MAXC=63, MAXBLK=64, MAXCB=64, IBACKC=0, IFOREC=1) REAL*4 FOCEAN(IM,JM),XYZD(IM,JM) INTEGER*4 INDEXP(MAXBLK), KCOLOR(720),IEDGE(720) CHARACTER TITLE*80, JGRID*1,JLAND*1,JWHIT*1, ARG*80, CBAR*8 COMMON /COLORS/ REDP(0:MAXC),GRNP(0:MAXC),BLUP(0:MAXC),NMCOLR COMMON /CBDTCB/ KCBDT(MAXBLK,MAXCB),NBLKS(MAXCB) C**** C**** Write RGB color values to PostScript file C**** WRITE (9,910) WRITE (9,911) '/RED [',(REDP(N),N=0,NMCOLR) WRITE (9,*) ' ] def' WRITE (9,911) '/GREEN [',(GRNP(N),N=0,NMCOLR) WRITE (9,*) ' ] def' WRITE (9,911) '/BLUE [',(BLUP(N),N=0,NMCOLR) WRITE (9,*) ' ] def' C**** C**** Write title to PostScript file C**** LEN = LEN_TRIM(TITLE) NLP = 0 NRP = 0 DO 210 N=1,LEN IF(TITLE(N:N).eq.'(') NLP = NLP + 1 210 IF(TITLE(N:N).eq.')') NRP = NRP + 1 IF(NLP.le.NRP) GO TO 230 C**** Too many Left Parentheses: replace some with blanks DO 220 N=LEN,1,-1 IF(TITLE(N:N).ne.'(') GO TO 220 TITLE(N:N) = ' ' NLP = NLP - 1 IF(NLP.eq.NRP) GO TO 250 220 continue STOP 220 230 IF(NRP.le.NLP) GO TO 250 C**** Too many Right Parentheses: replace some with blanks DO 240 N=LEN,1,-1 IF(TITLE(N:N).ne.')') GO TO 240 TITLE(N:N) = ' ' NRP = NRP - 1 IF(NRP.eq.NLP) GO TO 250 240 continue STOP 240 250 WRITE (9,920) TITLE(1:LEN) C**** C**** Write longitude by latitude color plot to PostScript file C**** SCALEX = 72*9./IM SCALEY = 72*6./JM WRITE (9,930) IM,JM,SCALEX,SCALEY MBLK = NBLKS(JCBDT) BCEN = MBLK/2. DO 320 J=1,JM WRITE (9,931) J-JM/2.,J DO 310 IR=1,IM I = IR+JROTA IF(I.gt.IM) I = I-IM KCOLOR(IR) = 1 IF( JWHIT.eq.'B') KCOLOR(IR) = 0 IF((JLAND.eq.'O' .and. FOCEAN(I,J).le..5) .or. * (JLAND.eq.'L' .and. FOCEAN(I,J).gt..5)) GO TO 310 IF(XYZD(I,J).le.DATMIS) KCOLOR(IR) = 1-KCOLOR(IR) IF(XYZD(I,J).le.DATMIS) GO TO 310 K = CEILING (4.*(XYZD(I,J)-CSCLJ)/DSCLJ + BCEN) IF(K.le.0 ) K = 1 IF(K.gt.MBLK) K = MBLK KCOLOR(IR) = KCBDT(K,JCBDT) 310 continue WRITE (9,932) (KCOLOR(IR),IR=1,IM) 320 WRITE (9,933) C**** C**** Write continental outline to PostScript file C**** IF(FOCEAN(1,1).lt.0. .or. JLAND.ne.'A') GO TO 490 KPEDGE = 0 IF(JWHIT.eq.'W') KPEDGE = 1 WRITE (9,940) KPEDGE C**** Draw horizontal edges between north-south grid boxes DO 430 J=1,JM-1 WRITE (9,941) J-JM/2.,J NEDGES = 0 DO 420 IR=1,IM I = IR+JROTA IF(I.gt.IM) I = I-IM IF(FOCEAN(I,J).gt..5 .xor. FOCEAN(I,J+1).le..5) GO TO 420 NEDGES = NEDGES + 1 IEDGE(NEDGES) = IR 420 continue IF(NEDGES.le.0) GO TO 430 WRITE (9,942) (IEDGE(N),N=1,NEDGES) WRITE (9,943) NEDGES 430 continue C**** Draw vertical edges between east-west grid boxes DO 450 J=2,JM-1 WRITE (9,944) J-JM/2.,J NEDGES = 0 I = 1+JROTA DO 440 IR=1,IM-1 IP1 = IR+1+JROTA IF(IP1.gt.IM) IP1 = IP1-IM IF(FOCEAN(I,J).gt..5 .xor. FOCEAN(IP1,J).le..5) GO TO 440 NEDGES = NEDGES + 1 IEDGE(NEDGES) = IR 440 I = IP1 IF(NEDGES.le.0) GO TO 450 WRITE (9,942) (IEDGE(N),N=1,NEDGES) WRITE (9,945) NEDGES 450 continue 490 WRITE (9,*) ' grestore' C**** C**** Write colorbar to PostScript file C**** WRITE (9,950) 72*6./MBLK WRITE (9,932) (KCBDT(K,JCBDT),K=1,MBLK) WRITE (9,951) MBLK C**** Write colorbar labels WRITE (9,952) N = NINT (4.*CSCLJ/DSCLJ) N4 = MODULO (N,4) B1 = BCEN - N4 - 4.*INT((BCEN-N4-.5)/4.) DO 560 B=B1,MBLK-.5,4. YBAR = B*72*6./MBLK - 72*3. - 4 VALU = CSCLJ + (B-BCEN)*DSCLJ/4 WRITE (CBAR,953) NINT(VALU) IF(NINT(DSCLJ* 10).ne.NINT(DSCLJ )*10) WRITE (CBAR,954) VALU IF(NINT(DSCLJ* 100).ne.NINT(DSCLJ* 10)*10) WRITE (CBAR,955) VALU IF(NINT(DSCLJ*1000).ne.NINT(DSCLJ*100)*10) WRITE (CBAR,956) VALU 560 WRITE (9,957) YBAR,CBAR C**** C**** Write surrounding rectangle to PostScript file C**** WRITE (9,960) C**** C**** Write NASA/GISS logo to PostScript file C**** WRITE (9,970) RETURN C**** 910 FORMAT (/'%%% Load Color Table') 911 FORMAT (4X,A8,10F5.2 / (12X,10F5.2)) 920 FORMAT (/'%%% Produce TITLE' / * ' /Helvetica-Bold findfont 20 scalefont setfont' / * '(',A,')' / * ' gsave 0 226 moveto dup stringwidth pop' / * ' 648 exch div 1 scale show grestore') 930 FORMAT (/'/IM',I4,' def /JM',I4,' def %%% Resolution' / * '%%% Fill Each Grid Cell with Color' / * ' /COLORBOX { newpath Y moveto' / * ' -1 0 rlineto 0 -1 rlineto 1 0 rlineto closepath' / * ' dup RED exch get exch' / * ' dup GREEN exch get exch' / * ' BLUE exch get setrgbcolor fill } def' / * 'gsave ',2F8.4,' scale') 931 FORMAT (/'/Y',F6.1,' def %%% Data for J =',I4) 932 FORMAT (24I3) 933 FORMAT (' IM -1 1 {COLORBOX} for') 940 FORMAT (/'%%% Draw Continental Outline' / * ' /IEDGE { pop newpath Y moveto -1 0 rlineto stroke } def' / * ' /JEDGE { pop newpath Y moveto 0 -1 rlineto stroke } def' / * I5,' setgray 1 setlinecap .25 setlinewidth') 941 FORMAT (/'/Y',F6.1,' def %%% East-west edges for J =',I4) 942 FORMAT (18I4) 943 FORMAT (' 1 1',I4,' {IEDGE} for') 944 FORMAT (/'/Y',F6.1,' def %%% North-south edges for J =',I4) 945 FORMAT (' 1 1',I4,' {JEDGE} for') 950 FORMAT (/'%%% Draw Colorbar to left of data' / * '%%% COLORBAR procedure for each color table value' / * ' /COLORBAR { newpath 0 exch moveto' / * ' -1 0 rlineto 0 -1 rlineto 1 0 rlineto closepath' / * ' dup RED exch get exch' / * ' dup GREEN exch get exch' / * ' BLUE exch get setrgbcolor fill } def' / * 'gsave -10 -216 translate 8',F9.4,' scale') 951 FORMAT (I3,' -1 1 {COLORBAR} for grestore') 952 FORMAT ('%%% Label the Colorbar' / * ' /Helvetica-Bold findfont 16 scalefont setfont 0 setgray') 953 FORMAT (I7 ,'-') 954 FORMAT (F7.1,'-') 955 FORMAT (F7.2,'-') 956 FORMAT (F7.3,'-') 957 FORMAT (' -112',F8.2,' moveto (',A8,') dup stringwidth pop' / * ' 92 exch sub 0 rmoveto show') 960 FORMAT (/'%%% Draw surrounding rectangle' / * ' 2 setlinewidth newpath -1 -217 moveto' / * ' 650 0 rlineto 0 434 rlineto -650 0 rlineto', * ' closepath stroke') 970 FORMAT (/'%%% Draw NASA/GISS logo' / * ' /Helvetica-Bold findfont 12 scalefont setfont' / * ' newpath 582 -232 moveto (NASA/GISS) show') END SUBROUTINE PARSE (S,LENS,NTOK,KMIN,KMAX) C**** C**** PARSE parses a character string S of length LENS into its first C**** NTOK tokens. KMIN(n) and KMAX(n) are the first and last C**** characters of S for the n-th token. C**** INTEGER*4 KMIN(*),KMAX(*) CHARACTER S*(*) N = 0 K = 0 10 IF(N.ge.NTOK) RETURN N = N + 1 C**** Skip over blank characters between tokens 20 IF(K.ge.LENS) GO TO 50 K = K + 1 IF(S(K:K).eq.' ') GO TO 20 KMIN(N) = K C**** Contiguous nonblank characters are part of the token 30 IF(K.ge.LENS) GO TO 40 K = K + 1 IF(S(K:K).ne.' ') GO TO 30 KMAX(N) = K-1 GO TO 10 40 KMAX(N) = K N = N + 1 C**** Remaining tokens were not included in the string 50 DO 60 N0=N,NTOK KMIN(N0) = 0 60 KMAX(N0) = 0 RETURN END SUBROUTINE NEWREC (S,LENS, NREC, IU,NR,UNNN, *) C**** C**** NEWREC decodes the token S into a file and a record number, C**** and checks that the record is available C**** C**** Input: S = string containing X,Y,Z,D and/or record number C**** LENS = # of characters in string S C**** NREC = most recently dispaled records of files X, Y or Z C**** Output: IU = unit number 1,2,3,4 corresponding to X,Y,Z,D C**** NR = record numbr of unit IU determined from S and/or NREC C**** UNNN = 4 byte character string showing file and record numbr C**** Parameter (MAXREC=1800, MAXIJR=72*46*MAXREC) INTEGER*4 NREC(4) CHARACTER S*(LENS), UNNN*4 CHARACTER TITLE*80, JGRID*1,JLAND*1,JWHIT*1, JCBFN*8 COMMON /PARMCB/ DATMIS,IM,JM,NRECM(4),QOPEN(3) COMMON /DATACB/ XYZD(MAXIJR,4) COMMON /TITLCB/ TITLE(MAXREC,4),CSCLJ(MAXREC,4),DSCLJ(MAXREC,4), * JCBDT(MAXREC,4),JROTA(MAXREC,4),JCBFN(MAXREC,4), * JGRID(MAXREC,4),JLAND(MAXREC,4),JWHIT(MAXREC,4) C**** IF(S(1:1).eq.'D' .or. S(1:1).eq.'d') GO TO 10 C**** File is X, Y or Z KU = 0 IF(S(1:1).eq.'X' .or. S(1:1).eq.'x') KU = 1 IF(S(1:1).eq.'Y' .or. S(1:1).eq.'y') KU = 2 IF(S(1:1).eq.'Z' .or. S(1:1).eq.'z') KU = 3 IF(KU.eq.0) GO TO 810 IF(NRECM(KU).le.0) GO TO 820 KN = NREC(KU) IOSQ = 0 IF(LENS.gt.1 .and. S(2:2).ne.' ') * READ (S(2:LENS),*,IOSTAT=IOSQ) KN IF(IOSQ.ne.0) GO TO 830 CALL READIN (KU,KN,*800) IMJMKN = 1 + IM*JM*(KN-1) CALL MATCHI (IM*JM,XYZD(IMJMKN,KU),TITLE(KN,KU),CSCLJ(KN,KU), * DSCLJ(KN,KU),JCBDT(KN,KU),JROTA(KN,KU),JCBFN(KN,KU), * JGRID(KN,KU),JLAND(KN,KU),JWHIT(KN,KU), DATMIS) IU = KU NR = KN WRITE (UNNN,900) CHAR(Z'57'+IU),MOD(NR,1000) IF(NR.gt.999) UNNN(2:2) = CHAR(Z'37'+NR/100) RETURN C**** File is D 10 KN = NREC(4) IOSQ = 0 IF(LENS.gt.1 .and. S(2:2).ne.' ') * READ (S(2:LENS),*,IOSTAT=IOSQ) KN IF(IOSQ.ne.0) GO TO 830 IF(KN.gt.NRECM(4)) GO TO 840 IU = 4 NR = KN WRITE (UNNN,900) 'D',MOD(NR,1000) IF(NR.gt.999) UNNN(2:2) = CHAR(Z'37'+NR/100) RETURN C**** 800 RETURN 1 810 WRITE (6,981) 'File should be X, Y, Z or D. Not: ' // S RETURN 1 820 WRITE (6,981) S(1:1) // ' file was never opened.' RETURN 1 830 WRITE (6,981) 'Integer record number should follow file: ' // S RETURN 1 840 WRITE (6,981) 'Requested difference record is not defined: '// S RETURN 1 C**** 900 FORMAT (A1,I3.3) 981 FORMAT (3X,A) END SUBROUTINE FILEEX (FIN,LENF,FEX) C**** C**** FILEEX expands a filename FIN by environment variables defined C**** in the user's .profile file. C**** CHARACTER FIN*80, FEX*80 IF(FIN(1:1).eq.'$') GO TO 10 FEX = FIN(1:LENF) RETURN C**** Leading $ was found, determine / or . 10 DO 20 K=3,LENF 20 IF(FIN(K:K).eq.'/' .or. FIN(K:K).eq.'.') GO TO 30 FEX = FIN(1:LENF) RETURN C**** Internal / or . was found 30 KSLASH = K CALL GETENV (FIN(2:KSLASH-1),FEX) IF(FEX(1:1).gt.' ') GO TO 40 FEX = FIN(1:LENF) RETURN C**** Insert the rest of FIN into the output filename FEX 40 LENE = LEN_TRIM(FEX) FEX(LENE+1:80) = FIN(KSLASH:LENF) RETURN END FUNCTION ROUND (X) C**** C**** ROUND returns a rounded value of the positive number X C**** IF(X.lt.9000.) GO TO 10 ROUND = 80000. IF(X.ge.65000.) RETURN ROUND = NINT(X*.0001)*10000 RETURN C**** 10 IF(X.lt.900.) GO TO 20 ROUND = 8000. IF(X.ge.6500.) RETURN ROUND = NINT(X*.001)*1000 RETURN C**** 20 IF(X.lt.90.) GO TO 30 ROUND = 800. IF(X.ge.650.) RETURN ROUND = NINT(X*.01)*100 RETURN C**** 30 IF(X.lt.9.) GO TO 40 ROUND = 80. IF(X.ge.65.) RETURN ROUND = NINT(X*.1)*10 RETURN C**** 40 IF(X.lt..9) GO TO 50 ROUND = 8. IF(X.ge.6.5) RETURN ROUND = NINT(X) RETURN C**** 50 IF(X.lt..09) GO TO 60 ROUND = .8 IF(X.ge..65) RETURN ROUND = NINT(X*10.)*.1 RETURN C**** 60 IF(X.lt..015) GO TO 70 ROUND = .08 IF(X.ge..065) RETURN ROUND = NINT(X*100.)*.01 RETURN C**** 70 ROUND = .01 RETURN END SUBROUTINE DIFFER (ARG,IMxJM,NREC,ND,DDDD,DATMIS, *) C**** C**** D =: calculate a Different record C**** Parameter (MAXREC=1800, MAXIJR=72*46*MAXREC) INTEGER*4 NREC(4) CHARACTER ARG*80, UQQQ*4, * TITLE*80, JCBFN*8, JGRID,JLAND,JWHIT, TEXT*80 REAL*4 DDDD(IMxJM,*) COMMON /DATACB/ XYZD(MAXIJR,4) COMMON /TITLCB/ TITLE(MAXREC,4),CSCLJ(MAXREC,4),DSCLJ(MAXREC,4), * JCBDT(MAXREC,4),JROTA(MAXREC,4),JCBFN(MAXREC,4), * JGRID(MAXREC,4),JLAND(MAXREC,4),JWHIT(MAXREC,4) C**** DO 10 I=1,IMxJM 10 DDDD(I,ND) = 0. NQ = 0 C**** KEQU = SCAN(ARG,'=') KMAX = LEN_TRIM(ARG) 20 KMAX = LEN_TRIM(ARG(1:KMAX)) IF(KEQU.ge.KMAX) GO TO 80 KR = SCAN(ARG(KEQU+1:KMAX),'XYZDxyzd',.TRUE.) KS = SCAN(ARG(KEQU+1:KMAX),'+-',.TRUE.) KT = SCAN(ARG(KEQU+1:KMAX),'*',.TRUE.) IF(KT.gt.KS) GO TO 50 IF(KR.gt.KS) GO TO 40 C**** Current term is a signed constant READ (ARG(KEQU+KS+1:KMAX),*,ERR=90) FAC IF(ARG(KEQU+KS:KEQU+KS).eq.'-') FAC = -FAC DO 30 I=1,IMxJM 30 IF(DDDD(I,ND).gt.DATMIS) DDDD(I,ND) = DDDD(I,ND) + FAC KMAX = KEQU + KS - 1 GO TO 20 C**** Current term is a signed data record 40 FAC = 1. GO TO 60 C**** Current term is a signed number times a data record 50 READ (ARG(KEQU+KS+1:KEQU+KT-1),*,ERR=90) FAC 60 LENS = KMAX-KEQU-KR+1 CALL NEWREC (ARG(KEQU+KR:KMAX),LENS, NREC, IQ,NQ,UQQQ, *91) IF(ARG(KEQU+KS:KEQU+KS).eq.'-') FAC = -FAC DO 70 I=1,IMxJM IF(XYZD(I+IMxJM*(NQ-1),IQ).gt.DATMIS) then IF(DDDD(I,ND).gt.DATMIS) * DDDD(I,ND) = DDDD(I,ND) + FAC*XYZD(I+IMxJM*(NQ-1),IQ) else DDDD(I,ND) = DATMIS endif 70 continue KMAX = KEQU + KS - 1 GO TO 20 C**** 80 IF(NQ.gt.0.) TITLE(ND,4) = TITLE(NQ,IQ) TITLE(ND,4)(47:50) = 'DIFF' RETURN C**** Invalid format 90 WRITE (6,*) ' Invalid format for Different record.' 91 RETURN 1 END SUBROUTINE RESOLU (FILEIN, IM,JM,ZFILE, *) C**** C**** RESOLU reads in the first four byte integer of a datafile and C**** from this information it determines the horizontal resolution. C**** C**** Input: FILEIN = character file name of a DATAFILE C**** Output: IM,JM = horizontal resolution C**** ZFILE = character file name of the ocean fraction file C**** * = use RETURN 1 if resolution is not determined C**** CHARACTER FILEIN*(*), ZFILE*(*) C**** OPEN (1,FILE=FILEIN,ACCESS='DIRECT',RECL=4,STATUS='OLD',ERR=801) READ (1,REC=1) IJT CLOSE (1) IJRES = (IJT-80)/4 C**** IF(IJRES.eq.46*9) then IM = 46 JM = 9 ZFILE = 'None' RETURN endif C**** IF(IJRES.eq.46*12) then IM = 46 JM = 12 ZFILE = 'None' RETURN endif C**** If( IJRES == 72*46 .or. IJRES == 72*46*2 ) then IM = 72 JM = 46 ZFILE = 'Z5X4N' RETURN endif C**** If( IJRES == 72*46+46+1 ) then IM = 72 JM = 46 ZFILE = 'Z5X4N' RETURN endif C**** IF(IJRES.eq.90*60) then IM = 90 JM = 60 ZFILE = 'Z4X3N' RETURN endif C**** IF(IJRES.eq.144*72) then IM = 144 JM = 72 ZFILE = 'Z144X72' RETURN endif C**** IF(IJRES.eq.144*73) then IM = 144 JM = 73 ZFILE = 'None' RETURN endif C**** IF(IJRES.eq.120*90) then IM = 120 JM = 90 ZFILE = 'Z3X2N' RETURN endif C**** IF(IJRES.eq.144*90) then IM = 144 JM = 90 ZFILE = 'Z144X90N' RETURN endif C**** IF(IJRES.eq.360*180) then IM = 360 JM = 180 ZFILE = 'Z360X180' RETURN endif C**** IF(IJRES.eq.720*360) then IM = 720 JM = 360 ZFILE = 'Z720X360' RETURN endif C**** Unable to derive resolution from record length IM = IJRES JM = 1 ZFILE = 'None' RETURN C**** Unable to open datafile 801 RETURN 1 END BLOCK DATA C**** C**** Default values for PARMCB C**** LOGICAL*4 QOPEN COMMON /PARMCB/ DATMIS,IM,JM,NRECM(4),QOPEN(3) C**** C**** DATMIS = upper boundary for missing data C**** IM,JM = longitude and latitude horizontal resolution C**** NRECM = number of records read from data file so far C**** NRECM(4) = number of difference records calculated so far C**** QOPEN = whether data file is currently open C**** DATA DATMIS/-999999./, IM,JM/72,46/, NRECM/4*0/, QOPEN/3*.false./ END