C
C Program GRDPS to Display GRiD values in Post Script C C Version: 5.60 C Date: 2001, September 4 C C Coded by Ludek Klimes C Department of Geophysics, Charles University Prague C Ke Karlovu 3, 121 16 Praha 2, Czech Republic C E-mail: klimes@seis.karlov.mff.cuni.cz C C....................................................................... C C Description of data files: C C Input data read from the standard input device (*): C The data are read by the list directed input (free format) and C consist of a single string 'SEP': C 'SEP'...String in apostrophes containing the name of the input C SEP parameter or history file with the input data. C No default, 'SEP' must be specified and cannot be blank. C C C Input data file 'SEP': C File 'SEP' has the form of the SEP C parameter file. The parameters, which do not differ from their C defaults, need not be specified in file 'SEP'. C Names of the input and output files: C GRD='string' ... Name of the input ASCII file with the grid values C to be plotted. C For general description of files with gridded data refer C to file forms.htm. C Default: GRD='grd.out' C PS='string' ... Name of the output PostScript file. C If non-blank PSHEX is given, output file PS contains just C the prolog part of the output PostScript file, see PSHEX. C If the value of parameter N4 is greater than 1, PS is C the output filename for the first spatial grid. C The filenames for the subsequent spatial grids are C generated in such a way that all digits contained within C the filename are considered a number which is increased C by 1 for each new time slice. C Default: PS='grd.ps' C PSHEX='string' ... Optional name of the output file containing C hexadecimal grid data and the PostScript trailer. C If PSHEX is blank (default), file PS contains the whole C PostScript file. C Otherwise, file PS contains just the prolog, setup, and C object definition parts, whereas the grid data are written C to file PSHEX. C The separation into two files is designed to facilitate C batch or manual editing of the PostScript definitions. C Resulting PostScript file is then obtained by merging C PS, possible new definitions and PSHEX. C If the value of parameter N4 is greater than 1, PSHEX is C the output filename for the first spatial grid. C The filenames for the subsequent spatial grids are C generated in such a way that all digits contained within C the filename are considered a number which is increased C by 1 for each new time slice. C Default: PSHEX=' ' C GRDPS2='string', GRDPS3='string', ... , GRDPS9='string' ... Names C of the optional input ASCII files with the grid values to C be plotted together with data of GRD. C For example, data GRD may control hue (colour) and data C GRDPS2 may control brightness or saturation. C For general description of files with gridded data refer C to file forms.htm. C Defaults: GRDSP2=' ', ... , GRDPS9=' ' C Data specifying grid dimensions: C N1=positive integer... Number of gridpoints along the X1 axis. C Default: N1=1 C N2=positive integer... Number of gridpoints along the X2 axis. C Default: N2=1 C N3=positive integer... Number of gridpoints along the X3 axis. C Default: N3=1 C N4=positive integer... Number of time slices (snapshots). C Individual time slices from the input file 'FILEIN' are C separated into different files. The name(s) of the first C file(s) is(are) given by input parameters 'PS' and C 'PSHEX' described above. Next filenames are generated in C such a way that all digits contained within the filenames C are considered as numbers which are increased by 1 for C each new time slice. C Default: N4=1 C D4=real... Time interval between two consecutive time slices. C Useful only if explicitly specifying time-dependent colour C scaling of snapshots (see input parameters VCIRC4, VSAT4 C and VBRI4 below). C Default: D4=1. C Additional data specific to this program: C (a) Dimensions and layout: C NH=positive integer... Layout of 2-D sections of a 3-D grid. C Individual 2-D sections of dimensions N1*N2 are merged C into rows NH sections long. Resulting (N3-1)/NH+1 rows C are merged together to form the final image. C Default: NH=INT(SQRT(N3)) C UNIT='string'... All lengths controlling the size and position of C the plot are assumed to be expressed in the units given C by the string. The units also influence the default C paper size, plot size and margins. Allowed values: C UNIT='cm': centimetres (default), C UNIT='in': inches (1in=2.54cm), C UNIT='pt': big points (1in=72pt). C Points (pt) are useful to generate plots for conversion to C bitmap forms, e.g., using GhostScript. C Default: UNIT='cm' C XSIGN=real... Determines the sign of the default value of HSIZE. C Default: XSIGN=1. C HSIZE=real... Size (in UNITs) of the image, corresponding to the X1 C axis (horizontal before a possible rotation). C If negative, the values will be displayed from the right C to the left. C Default: HSIZE=SIGN( 16.0,XSIGN) for UNIT='cm', C HSIZE=SIGN( 6.5,XSIGN) for UNIT='in', C HSIZE=SIGN(NH*N1,XSIGN) for UNIT='pt'. C YSIGN=real... Determines the sign of the default value of VSIZE. C Default: YSIGN=1. C VSIZE=real... Size (in UNITs) of the image, corresponding to the X2 C axis (vertical before a possible rotation). C If negative, the values will be displayed from the top to C the bottom. C Default (proportional display): C VSIZE=SIGN(HSIZE*((N3-1)/NH+1)*N2/(NH*N1),YSIGN), i.e., C VSIZE=SIGN(HSIZE*N2/N1,YSIGN) for N3=1. C HOFFSET=real... Distance (in UNITs) of the image from the leftmost C paper edge (before a possible rotation). Controls the C horizontal position of the figure. C Default: HOFFSET=2.5 for UNIT='cm', C HOFFSET=1.0 for UNIT='in', C HOFFSET=0.0 for UNIT='pt'. C VOFFSET=real... Distance (in UNITs) of the image from the bottom C paper edge (before a possible rotation). Controls the C vertical position of the figure. C Default: C if VSIZE.LE.HEIGHT-2*2.5: VOFFSET=HEIGHT-2.5-VSIZE C otherwise if VSIZE.LE.HEIGHT: VOFFSET=(HEIGHT-VSIZE)/2. C otherwise: VOFFSET=2.5 C HEIGHT=real... Height of the paper in a portrait position. C Default: HEIGHT=29.7 for UNIT='cm', C HEIGHT=11.0 for UNIT='in', C HEIGHT=((N3-1)/NH+1)*N2 for UNIT='pt'. C ROTATE=real... Enables to rotate the image by angle specified in C degrees (positive counterclockwise). The image is rotated C around the centre of the square paper of size HEIGHT. C If applied, the user will probably wish to specify the C value of ROTATE=90. C Parameters HSIZE,VSIZE,HOFFSET,VOFFSET apply to the image C before rotation. C Attention: BoundingBox is incorrect if ROTATE is not C multiple of 90 degrees. C Default: ROTATE=0. C (b) Range of displayed values: C VMIN=real, VMAX=real... Values less than or equal to VMIN, C or greater than or equal to VMAX will be deemed C undefined. For example, if VMIN=0 is set when plotting C velocities, zero velocities are rendered as undefined C values, that is usually desirable. C Defaults: VMIN=-999999, VMAX=999999000000. C R=real, G=real, B=real... Colour of the undefined C values. C Defaults: R=0.80, G=0.80, B=0.80 (light grey) C (c) Colour interpretation of values contained within file GRD: C VPLUS=real, VSIGN=real, VCIRC=real... VCIRC is the extent of C values corresponding to the whole colour circle RGB. C Negative VCIRC corresponds to the opposite direction C around the colour circle (BGR). C Undefined values are displayed in gray. C Defaults: VPLUS=0, VSIGN=1, C VCIRC=(GMAX1-GMIN1+VPLUS)*VSIGN, C where GMINi and GMAXi are the minimum and maximum grid C values defined in file 'GRDi'. C If the above value of VCIRC is zero (e.g., if C GMAX1=GMIN1), the default value of VCIRC is changed to C VCIRC=1. C Hint: Specify VPLUS=1 when plotting integer values. C Hint: To plot velocities with blue minimum velocity, C increasing through green, yellow and red again to blue, C specify VMIN=0 and VSIGN=-1 with default VREF and CREF. C VREF=real, CREF=real... Value VREF corresponds to colour CREF, C where Red=0, Yellow=1/6, Green=2/6, Cyan=3/6, Blue=4/6, C Magenta=5/6, Red=1, Yellow=7/6, Green=8/6, etc. C Default: VREF=GMIN1, CREF=4/6 (blue). C where GMIN1 is the minimum defined grid value. C Default VMIN, VMAX, VPLUS, VSIGN, VCIRC, VREF and CREF C render the central value (GMIN1+GMAX1)/2 in yellow. C Hint: To plot the values oscillating around 0, you may C wish to set VREF=0 and CREF=0.166667 to render 0 in C yellow. C VCIRC4=real... If VCIRC is specified, VCIRC4 enables to modify it C in dependence on progressing time levels I4=1,2,...,N4: C VCIRC(I4)=VCIRC*EXP(-VCIRC4*(I4-1)*D4) C Default: VCIRC4=0. C (d) Saturation interpretation of values contained within file GRDPS2: C VSAT=real... The extent of values corresponding to the whole C saturation range from 0 (white) to 1 (fully saturated C colours). C Default: VSAT=GMAX2-GMIN2 (or VSAT=1 for GMIN2=GMAX2). C VWHITE=real... Value corresponding to white (saturation=0). C Default: VWHITE=GMIN2 (strictly speaking, =GMAX2-VSAT). C VSAT4=real... If VSAT is specified, VSAT4 enables to modify it in C dependence on progressing time levels I4=1,2,...,N4: C VSAT(I4)=VSAT*EXP(-VSAT4*(I4-1)*D4) C Default: VSAT4=0. C (e) Brightness interpretation of values contained within file GRDPS3: C VBRI=real... The extent of values corresponding to the whole C brightness range from 0 (black) to 1 (bright colours). C Default: VBRI=GMAX3-GMIN3 (or VBRI=1 for GMIN3=GMAX3). C VBLACK=real... Value corresponding to black (brightness=0). C Default: VBLACK=GMIN3 (strictly speaking, =GMAX3-VBRI). C VBRI4=real... If VBRI is specified, VBRI4 enables to modify it in C dependence on progressing time levels I4=1,2,...,N4: C VBRI(I4)=VBRI*EXP(-VBRI4*(I4-1)*D4) C Default: VBRI4=0. C (f) Form of the output PostScript file: C SHOWPAGE=integer... PostScript command 'showpage' at the end of C file directs printer to print and delete the picture. C This is usually what we want. However, sometimes we may C wish to overlay the picture with another one. In this C case, we wish to remove the 'showpage'. C SHOWPAGE=0: Two last lines of the file containing strings C 'showpage' and '%%EOF' are not written. C SHOWPAGE=1: The lines are written. C Default: SHOWPAGE=1 C C======================================================================= C C Common block /RAMC/: INCLUDE 'ram.inc' C ram.inc C INTEGER MGRID PARAMETER (MGRID=MRAM) REAL GRID(MGRID) EQUIVALENCE (GRID,RAM) INTEGER IGRID(MGRID) EQUIVALENCE (GRID,IGRID) C C....................................................................... C CHARACTER*2 UNIT REAL UNITPT,HEIGHT,OFFSET,WIDTH C C UNIT... One of: 'cm', 'in', 'pt'. C UNITPT...Size of the length unit, in which input data controlling C the size and position of the plot are expressed, in big C points (pt). E.g., UNITPT=72./2.54 corresponds to plotting C in cm. C HEIGHT..Anticipated height of the paper sheet. C OFFSET..Left margin, and top or bottom margin for low or high C plots, respectively. C WIDTH...Default width of the plot. C INTEGER MFILE,LU PARAMETER (MFILE=9,LU=10) CHARACTER*80 FSEP,FPS,FPSHEX,FILE(MFILE) CHARACTER*1 HEX1(0:15) CHARACTER*2 HEX2(0:255) CHARACTER*6 FGRDPS INTEGER LUIN(MFILE) REAL UNDEF PARAMETER (UNDEF=-999999.) INTEGER NFILE,IFILE,N1,N2,N3,N31,N32,N4,J,J0,I,I0,I1,I2,I31,I32,I4 INTEGER ISHOW REAL XSIGN,YSIGN REAL VMIN,VMAX,VPLUS,VSIGN,VCIRC,VREF,CREF,ROTATE,RED,GREEN,BLUE REAL GMIN(MFILE),GMAX(MFILE),G,G0,GSTEP REAL BBOX1,BBOX2,BBOX3,BBOX4,BB1,BB2,BB3,BB4 REAL BB1P,BB2P,BB3P,BB4P,BB2DEF,BB4DEF,AUX,VAUX,C,S DATA LUIN /1,2,3,4,5,6,7,8,9/ DATA HEX1 * /'0','1','2','3','4','5','6','7','8','9','A','B','C','D','E','F'/ C C----------------------------------------------------------------------- C C Reading name of SEP file with input data: WRITE(*,'(A)') '+GRDPS: Enter input filename: ' FSEP=' ' READ(*,*) FSEP WRITE(*,'(A)') '+GRDPS: Working ... ' C C Reading all data from the SEP file into the memory: IF (FSEP.NE.' ') THEN CALL RSEP1(LU,FSEP) ELSE C GRDPS-04 CALL ERROR('GRDPS-04: SEP file not given') C Input file in the form of the SEP (Stanford Exploration Project) C parameter or history file must be specified. C There is no default filename. ENDIF C C Reading input parameters from the SEP file: CALL RSEP3T('GRD',FILE(1),'grd.out') CALL RSEP3T('PS',FPS,'grd.ps') CALL RSEP3T('PSHEX',FPSHEX,' ') FGRDPS='GRDPS0' DO 12 IFILE=1,MFILE FGRDPS(6:6)=CHAR(ICHAR('0')+IFILE) IF (IFILE.NE.1) CALL RSEP3T(FGRDPS,FILE(IFILE),' ') IF(FILE(IFILE).EQ.' ') THEN NFILE=IFILE-1 GO TO 13 END IF 12 CONTINUE 13 CONTINUE CALL RSEP3I('SHOWPAGE',ISHOW,1) C C Recalling the data specifying grid dimensions C (arguments: Name of value in input data, Variable, Default): CALL RSEP3I('N1',N1,1) CALL RSEP3I('N2',N2,1) CALL RSEP3I('N3',N3,1) CALL RSEP3I('N4',N4,1) C Transforming 2-D arrays to N1*N2*1 and 1-D arrays to N1*1*1 IF(N1.EQ.1) THEN N1=N2 N2=N3 N3=1 END IF IF(N1.EQ.1) THEN N1=N2 N2=1 END IF IF(N2.EQ.1) THEN N2=N3 N3=1 END IF C Determining the layout of 2-D slices of 3-D arrays CALL RSEP3I('NH',N31,INT(SQRT(FLOAT(N3))+.001)) N32=(N3-1)/N31+1 IF(N1*N2*N31*N32*NFILE.GT.MGRID) THEN C GRDPS-01 CALL ERROR('GRDPS-01: Too small array RAM(MRAM)') C Array RAM(MRAM) allocated in include file 'ram.inc' is too small C to contain the grid values to be plotted. You may wish to C increse the dimension MRAM in file 'ram.inc'. C ram.inc END IF C N3 slices will be arranged into N31*N32 slices. C C Recalling the plotting unit and setting default dimensions: CALL RSEP3T('UNIT',UNIT,'cm') CALL LOWER(UNIT) IF(UNIT.EQ.'cm') THEN UNITPT=72./2.54 HEIGHT=29.7 OFFSET=2.5 WIDTH=16.0 ELSE IF(UNIT.EQ.'in') THEN UNITPT=72. HEIGHT=11.0 OFFSET=1.0 WIDTH=6.5 ELSE IF(UNIT.EQ.'pt') THEN UNITPT=1. HEIGHT=FLOAT(N32*N2) OFFSET=0.0 WIDTH=FLOAT(N31*N1) ELSE C GRDPS-02 CALL ERROR('GRDPS-02: Unrecognized plotting units') C Allocated plotting units are UNIT='cm', UNIT='in' or UNIT='pt'. END IF C....................................................................... C C Loop over time slices: DO 90 I4=1,N4 IF(I4.EQ.1) THEN IF(N4.GT.1) THEN DO 15 IFILE=1,NFILE OPEN(LUIN(IFILE),FILE=FILE(IFILE), * FORM='FORMATTED',STATUS='OLD') FILE(IFILE)=' ' 15 CONTINUE END IF ELSE C Generating new output filenames: CALL NEWNAM(FPS) CALL NEWNAM(FPSHEX) END IF C C Reading input grid values: DO 16 IFILE=1,NFILE CALL RARRAY(LUIN(IFILE),FILE(IFILE),'FORMATTED',.TRUE.,UNDEF, * N1*N2*N3,GRID(N1*N2*N31*N32*(IFILE-1)+1)) 16 CONTINUE C C Rearranging 3-D input N1*N2*N3 grid into (N1*N31)*(N2*N32) grid: DO 29 IFILE=0,N1*N2*N31*N32*(NFILE-1),N1*N2*N31*N32 C Extending the grid to N1*N2*N31*N32 points DO 21 I=IFILE+N1*N2*N3+1,IFILE+N1*N2*N31*N32 GRID(I)=UNDEF 21 CONTINUE C Transposing the N1*N2*N31*N32 grid to N1*N31*N2*N32 grid DO 25 I32=0,N2*N31*(N32-1),N2*N31 DO 24 I0=0,N2*N31-1 J0=I0 22 CONTINUE I2=J0/N31 I31=J0-I2*N31 J0=I31*N2+I2 IF(J0.LT.I0) GO TO 22 I=IFILE+(I32+I0)*N1 J=IFILE+(I32+J0)*N1 DO 23 I1=1,N1 I=I+1 J=J+1 AUX=GRID(I) GRID(I)=GRID(J) GRID(J)=AUX 23 CONTINUE 24 CONTINUE 25 CONTINUE 29 CONTINUE C N1=N1*N31 N2=N2*N32 C C....................................................................... C C Recalling the data for the plotting area: CALL RSEP3R('XSIGN' ,XSIGN,1.) CALL RSEP3R('YSIGN' ,YSIGN,1.) AUX=HEIGHT CALL RSEP3R('HEIGHT' ,HEIGHT,AUX) CALL RSEP3R('HSIZE' ,BB3,SIGN(WIDTH,XSIGN)) CALL RSEP3R('HOFFSET',BB1,OFFSET) C Default height of the figure BB4DEF=ABS(BB3)*FLOAT(N2)/FLOAT(N1) CALL RSEP3R('VSIZE' ,BB4,SIGN(BB4DEF,YSIGN)) C Default vertical position of the figure IF(ABS(BB4).LE.HEIGHT-2.*OFFSET) THEN BB2DEF=HEIGHT-OFFSET-ABS(BB4) ELSE IF(ABS(BB4).LE.HEIGHT) THEN BB2DEF=(HEIGHT-ABS(BB4))/2. ELSE BB2DEF=OFFSET END IF CALL RSEP3R('VOFFSET',BB2,BB2DEF) IF(BB3.LT.0.) THEN BB1=BB1-BB3 END IF IF(BB4.LT.0.) THEN BB2=BB2-BB4 END IF CALL RSEP3R('ROTATE',ROTATE,0.) C C Transformation from plotting units (e.g. centimetres) to points: BB1P=BB1*UNITPT BB2P=BB2*UNITPT BB3P=BB3*UNITPT BB4P=BB4*UNITPT C Bounding box: BBOX1=AMIN1(BB1P,BB1P+BB3P) BBOX2=AMIN1(BB2P,BB2P+BB4P) BBOX3=AMAX1(BB1P,BB1P+BB3P) BBOX4=AMAX1(BB2P,BB2P+BB4P) IF(ROTATE.NE.0.) THEN C=COS(ROTATE*3.14159/180.) S=SIN(ROTATE*3.14159/180.) BBOX1=BBOX1-HEIGHT*UNITPT/2. BBOX2=BBOX2-HEIGHT*UNITPT/2. BBOX3=BBOX3-HEIGHT*UNITPT/2. BBOX4=BBOX4-HEIGHT*UNITPT/2. AUX =C*BBOX1-S*BBOX2 BBOX2=S*BBOX1+C*BBOX2 BBOX1=AUX AUX =C*BBOX3-S*BBOX4 BBOX4=S*BBOX3+C*BBOX4 BBOX3=AUX BBOX1=BBOX1+HEIGHT*UNITPT/2. BBOX2=BBOX2+HEIGHT*UNITPT/2. BBOX3=BBOX3+HEIGHT*UNITPT/2. BBOX4=BBOX4+HEIGHT*UNITPT/2. AUX =AMIN1(BBOX1,BBOX3) BBOX3=AMAX1(BBOX1,BBOX3) BBOX1=AUX AUX =AMIN1(BBOX2,BBOX4) BBOX4=AMAX1(BBOX2,BBOX4) BBOX2=AUX END IF C C....................................................................... C C Calculating minimum and maximum values GMIN(*) and GMAX(*): CALL RSEP3R('VMIN',VMIN,-999999.) CALL RSEP3R('VMAX',VMAX, 999999000000.) DO 33 IFILE=1,NFILE NNN=N1*N2*(IFILE-1) GMINI=VMAX GMAXI=VMIN DO 31 I=NNN+1,NNN+N1*N2 G=GRID(I) IF(G.NE.UNDEF) THEN IF(G.LE.VMIN.OR.VMAX.LE.G) THEN GRID(I)=UNDEF ELSE GMINI=AMIN1(GMINI,G) GMAXI=AMAX1(GMAXI,G) END IF END IF 31 CONTINUE C C Rescaling range GMIN--GMAX to 0--254 c (undefined values are 255): IF(GMINI.NE.GMAXI) THEN GSTEP=(GMAXI-GMINI)/254. ELSE GSTEP=1. END IF G0=GMINI-GSTEP/2. DO 32 I=NNN+1,NNN+N1*N2 IF(GRID(I).NE.UNDEF) THEN G=(GRID(I)-G0)/GSTEP IGRID(I)=INT(G) ELSE IGRID(I)=255 END IF 32 CONTINUE C GMIN(IFILE)=GMINI GMAX(IFILE)=GMAXI 33 CONTINUE C C Preparing hexadecimal coding table: I=0 DO 42 I2=0,15 DO 41 I1=0,15 HEX2(I)(1:1)=HEX1(I2) HEX2(I)(2:2)=HEX1(I1) I=I+1 41 CONTINUE 42 CONTINUE C C Parameters controlling colours: CALL RSEP3R('D4' ,D4 ,1.) TIME=FLOAT(I4-1)*D4 CALL RSEP3R('VPLUS' ,VPLUS ,0.) CALL RSEP3R('VSIGN' ,VSIGN ,1.) CALL RSEP3R('VCIRC4',VAUX ,0.) VAUX=EXP(VAUX*TIME) AUX=(GMAX(1)-GMIN(1)+VPLUS)*VSIGN*VAUX IF(AUX.EQ.0.) AUX=1. CALL RSEP3R('VCIRC' ,VCIRC ,AUX) VCIRC=VCIRC/VAUX CALL RSEP3R('VREF' ,VREF ,GMIN(1)) CALL RSEP3R('CREF' ,CREF ,2./3.) IF(NFILE.GE.2) THEN CALL RSEP3R('VSAT4' ,VAUX ,0.) VAUX=EXP(VAUX*TIME) AUX=(GMAX(2)-GMIN(2))*VAUX IF(AUX.EQ.0.) AUX=1. CALL RSEP3R('VSAT' ,VSAT ,AUX) VSAT=VSAT/VAUX CALL RSEP3R('VWHITE',VWHITE,GMAX(2)-AUX) END IF IF(NFILE.GE.3) THEN CALL RSEP3R('VBRI4' ,VAUX ,0.) VAUX=EXP(VAUX*TIME) AUX=(GMAX(3)-GMIN(3))*VAUX IF(AUX.EQ.0.) AUX=1. CALL RSEP3R('VBRI' ,VBRI ,AUX) VBRI=VBRI/VAUX CALL RSEP3R('VBLACK',VBLACK,GMAX(3)-AUX) END IF C C Writing PostScript prolog: WRITE(*,'(''+'',79('' ''))') WRITE(*,'(2A)') '+Writing: ',FPS(1:MIN0(LEN(FPS),70)) OPEN(LU,FILE=FPS) WRITE(LU,'(A/A,4I6,/(A))') * '%!PS-Adobe-3.0', * '%%BoundingBox:',INT(BBOX1+.5),INT(BBOX2+.5), * INT(BBOX3+.5),INT(BBOX4+.5), * '%%EndComments', * '%%BeginProlog', * '%%BeginProcSet: (grdps)', * '%%Creator: grdps', * '%-----------------------------------------------------------', * '% Default UNDEFINED procedure:' C520 WRITE(LU,'(A)') C520 * '/UNDEFINED {0 0 0.80 sethsbcolor} bind def' CALL RSEP3R('R',RED ,0.80) CALL RSEP3R('G',GREEN,0.80) CALL RSEP3R('B',BLUE ,0.80) WRITE(LU,'(A,3(F4.2,1X),A)') * '/UNDEFINED {',RED,GREEN,BLUE,'setrgbcolor} bind def' WRITE(LU,'(A)') * '% Default VALUEtoCOLOR procedure:' WRITE(LU,'(3(A,G13.6)/A,G13.6,A/A,G13.6,A/A)') * '/VCIRC ',VCIRC ,' def % Range of grid values:', * GMIN(1),' to',GMAX(1), * '/VREF ',VREF ,' def', * '/CREF ',CREF ,' def', * '/VRED VREF CREF VCIRC mul sub def' IF(NFILE.GE.2) THEN WRITE(LU,'(3(A,G13.6)/A,G13.6,A)') * '/VSAT ',VSAT ,' def % Range of grid values:', * GMIN(2),' to',GMAX(2), * '/VWHITE ',VWHITE,' def' END IF IF(NFILE.GE.3) THEN WRITE(LU,'(3(A,G13.6)/A,G13.6,A)') * '/VBRI ',VBRI ,' def % Range of grid values:', * GMIN(3),' to',GMAX(3), * '/VBLACK ',VBLACK,' def' END IF WRITE(LU,'(A)') * '/VALUEtoCOLOR{' IF(NFILE.LE.1) THEN WRITE(LU,'(A)') * ' VRED sub VCIRC div dup truncate sub dup 0 lt {1 add} if' ELSE WRITE(LU,'(A,I1,9A))') * ' ',NFILE,' -1 roll', * ' VRED sub VCIRC div dup truncate sub dup 0 lt {1 add} if' END IF IF(NFILE.LE.1) THEN WRITE(LU,'(A)') * ' 1' ELSE WRITE(LU,'(A,I1,9A))') * ' ',NFILE,' -1 roll', * ' VWHITE sub VSAT div' END IF IF(NFILE.LE.2) THEN WRITE(LU,'(A)') * ' 1' ELSE WRITE(LU,'(A,I1,9A))') * ' ',NFILE,' -1 roll', * ' VBLACK sub VBRI div' END IF DO 51 IFILE=4,NFILE WRITE(LU,'(A)') * ' pop' 51 CONTINUE WRITE(LU,'(A)') * ' sethsbcolor} bind def', * '%-----------------------------------------------------------', * '% User-defined VALUEtoCOLOR procedure may be inserted here:' IF(FPSHEX.NE.' ') THEN CLOSE(LU) WRITE(*,'(''+'',79('' ''))') WRITE(*,'(2A)') '+Writing: ',FPSHEX(1:MIN0(LEN(FPSHEX),70)) OPEN(LU,FILE=FPSHEX) END IF C C Writing the plotting procedure: WRITE(LU,'(A)') * '%-----------------------------------------------------------', * '%%EndProcSet', * '%%EndProlog', * '%-----------------------------------------------------------', * '%%BeginSetup', * '% Numerical values describing the image size and position:' WRITE(LU,'(A,I6,A)') '/N1',N1,' def' WRITE(LU,'(A,I6,A)') '/N2',N2,' def' WRITE(LU,'(A,F8.1,A,F8.3,A)') '/BB1',BB1P,' def %',BB1,'cm' WRITE(LU,'(A,F8.1,A,F8.3,A)') '/BB2',BB2P,' def %',BB2,'cm' WRITE(LU,'(A,F8.1,A,F8.3,A)') '/BB3',BB3P,' def %',BB3,'cm' WRITE(LU,'(A,F8.1,A,F8.3,A)') '/BB4',BB4P,' def %',BB4,'cm' WRITE(LU,'(A,F8.1,A)') '/PAPERSIZE',HEIGHT*UNITPT,' def' WRITE(LU,'(A,F8.1,A)') '/ROTATE',ROTATE,' def' WRITE(LU,'(A)') * '% Scaling of grid values:' DO 52 IFILE=1,NFILE WRITE(LU,'(A,I1,A,G13.6,A)') * '/GMIN',IFILE,' ',GMIN(IFILE),' def' WRITE(LU,'(A,I1,A,G13.6,A)') * '/GMAX',IFILE,' ',GMAX(IFILE),' def' 52 CONTINUE WRITE(LU,'(A)') * '%%EndSetup', * '%-----------------------------------------------------------', * '%%BeginObject: (grdps)', * 'PAPERSIZE 2 div dup translate ROTATE rotate', * 'PAPERSIZE -2 div dup translate', * '/SCALE1 N1 BB3 div def', * '/SCALE2 N2 BB4 div def' DO 53 IFILE=1,NFILE WRITE(LU,'(9(A,I1))') * '/VSTEP',IFILE,' GMAX',IFILE,' GMIN',IFILE,' sub 254 div def' 53 CONTINUE WRITE(LU,'(A)') * '/RGB 3 string def', * 'N1 N2 8 [ SCALE1 0 0 SCALE2 BB1 SCALE1 mul neg', * ' BB2 SCALE2 mul neg ]' WRITE(LU,'(A,I1,A)') * '{currentfile ',NFILE,' string readhexstring pop' DO 54 IFILE=1,NFILE WRITE(LU,'(9(A,I1))') * ' dup ',IFILE-1, * ' get VSTEP',IFILE,' mul GMIN',IFILE,' add exch' 54 CONTINUE WRITE(LU,'(11A)') * ' 0 get 255 eq {',('pop ',IFILE=1,NFILE), * 'UNDEFINED} {VALUEtoCOLOR} ifelse' WRITE(LU,'(A)') * ' currentrgbcolor', * ' 255 mul .5 add cvi RGB exch 2 exch put', * ' 255 mul .5 add cvi RGB exch 1 exch put', * ' 255 mul .5 add cvi RGB exch 0 exch put RGB}', * ' bind false 3 colorimage' C C Writing output hexadecimal values: N12=N1*N2 NNN=N1*N2*NFILE IF(NFILE.EQ.1) THEN WRITE(LU,'(40A2)') (HEX2(IGRID(I)),I=1,N12) ELSE DO 61 I=N12+1,NNN IF(IGRID(I).GE.255) THEN IGRID(MOD(I,N12))=255 END IF 61 CONTINUE WRITE(LU,'(40A2)') * ((HEX2(IGRID(IFILE)),IFILE=I,NNN,N12),I=1,N12) END IF C C Writing PostScript trailer: WRITE(LU,'(A)') * 'PAPERSIZE 2 div dup translate ROTATE neg rotate', * 'PAPERSIZE -2 div dup translate', * '%%EndObject' IF(ISHOW.NE.0) THEN WRITE(LU,'(A)') * 'showpage', * '%%EOF' END IF CLOSE(LU) C 90 CONTINUE C End of the loop over time slices. C IF(N4.GT.1) THEN DO 91 IFILE=1,NFILE CLOSE(LUIN(IFILE)) 91 CONTINUE END IF WRITE(*,'(A)') '+GRDPS: Done. ' STOP END C C======================================================================= C C C SUBROUTINE NEWNAM(NAME) CHARACTER*(*) NAME C C----------------------------------------------------------------------- C INTEGER N,I C IF(NAME.NE.' ') THEN N=LEN(NAME) 10 CONTINUE DO 11 I=N,1,-1 IF(LLE('0',NAME(I:I)).AND.LLE(NAME(I:I),'8')) THEN NAME(I:I)=CHAR(ICHAR(NAME(I:I))+1) GO TO 12 ELSE IF(NAME(I:I).EQ.'9') THEN NAME(I:I)='0' N=I-1 GO TO 10 END IF 11 CONTINUE C GRDPS-03 CALL ERROR('GRDPS-03: Too many output files') C The digits in the template name of the output files do not C allow for the generation of all output filenames. C The number of digits should be increased. 12 CONTINUE END IF RETURN END C C======================================================================= C INCLUDE 'error.for' C error.for INCLUDE 'sep.for' C sep.for INCLUDE 'forms.for' C forms.for INCLUDE 'length.for' C length.for C C======================================================================= C