00001
00002 c
00003 c Copyright (c) 1986,1987,1988,1989,1990,1991,1992,1993,
00004 c by Steve McMillan, Drexel University, Philadelphia, PA.
00005 c
00006 c All rights reserved.
00007 c
00008 c Redistribution and use in source and binary forms are permitted
00009 c provided that the above copyright notice and this paragraph are
00010 c duplicated in all such forms and that any documentation,
00011 c advertising materials, and other materials related to such
00012 c distribution and use acknowledge that the software was developed
00013 c by the author named above.
00014 c
00015 c THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR
00016 c IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED
00017 c WARRANTIES OF MERCHANTIBILITY AND FITNESS FOR A PARTICULAR PURPOSE.
00018 c
00019
00020 subroutine ez2dplot(a,na,xlen,ylen,iprompt,itype,*)
00021 save
00022 c
00023 c Plot a 2-D contour representation of the array a.
00024 c
00025 dimension a(1),xlevel(1000),temp(1000)
00026 character string*80,c*1,c2*1
00027 equivalence (string(1:1),c),(string(2:2),c2)
00028 character*50 list(20)
00029 c
00030 common /local offset/ xoff,yoff
00031 c
00032 save m,n,nl,xlevel,ipaint,invrt
00033 data m/0/n/0/nl/0/ipaint/0/invrt/0/
00034 c
00035 a2(i,j) = a(i+m*(j-1))
00036 c
00037 if (xlen.le.0..or.ylen.le.0.) return 1
00038 c
00039 call devoff
00040 if (iprompt.eq.1) write(6,'(''*** Two-dimensional plot ***'')')
00041 c
00042 call get_dimensions(m,n,na,iprompt,*99999)
00043 c
00044 call minmax(a,m*n,amin,amax)
00045 if (iprompt.eq.1) write(6,'(a,1p,e12.4,a,e12.4)')
00046 & 'Array minimum = ',amin,', maximum = ',amax
00047 if (amin.ge.amax) return
00048 c
00049 nl0 = nl
00050 ns = -1
00051 c
00052 xoff = 0.
00053 yoff = 0.
00054 c
00055 c Get the next sub-command (possibly prompting for input).
00056 c
00057 5 if (iprompt.eq.1) then
00058 call nextcmd(string,istart,ns,'Command(s) (h for help):')
00059 else
00060 call nextcmd(string,istart,ns,' ')
00061 end if
00062 c
00063 c Parse the command string.
00064 c
00065 if (c.eq.'q') then
00066
00067 if (nl.gt.0) then
00068 call sort(nl,xlevel)
00069 call devoff
00070 if (iprompt.eq.1) write(6,1000)(i,xlevel(i),i=1,nl)
00071 1000 format('Levels: ',1p,4(i4,': ',e12.4):/
00072 & (8x,4(i4,': ',e12.4)))
00073 end if
00074 return
00075
00076 else if (c.eq.'h') then
00077
00078 call devoff
00079 write(6,*)
00080 write(6,*)'b: redraw box'
00081 write(6,*)'co color: specify color'
00082 write(6,*)'e: erase'
00083 write(6,*)'d list: delete specified levels from the list'
00084 write(6,*)'g: draw contour through graphics cursor'
00085 write(6,*)'gd: delete contour through graphics cursor'
00086 write(6,*)'l: list current contour levels'
00087 write(6,*)'o: offset the contours'
00088 write(6,*)'p: paint the box'
00089 write(6,*)'p-: paint the box, inverted color map'
00090 write(6,*)'q: back to main mcdraw prompt'
00091 write(6,*)'r: repeat previously specified levels'
00092 write(6,*)'s list: select specified levels from the list'
00093 write(6,*)'t pattern: specify dashed lines'
00094 write(6,*)'w weight: set line weight'
00095 write(6,*)'x: delete contour list'
00096 write(6,*)'<number>: draw contour at this level'
00097 write(6,*)
00098
00099 else if (c.eq.'b') then
00100 c
00101 call box(0.,xlen,0.,ylen)
00102 c
00103 else if (c.eq.'c') then
00104 c
00105 c Set line color.
00106 c
00107 call readiq(string(istart:ns),1,
00108 & icolor,idum,idum,idum,*5)
00109 call color(icolor)
00110 c
00111 else if (c.eq.'d') then
00112 c
00113 c Delete specified levels, retain the rest.
00114 c
00115 if (nl.gt.0) then
00116 call gettokens(string(istart:ns),list,nlist)
00117 do 30 i=1,nlist
00118 read(list(i),*,err=30,end=30)ii
00119 call invert
00120 call xcontor(a,m,n,xlevel(ii),1,0.,xlen,0.,ylen,0,
00121 & amin,amax,ipaint,invrt)
00122 call invert
00123 30 continue
00124 c
00125 c Remove the levels from the list.
00126 c
00127 do 50 i=1,nlist
00128 read(list(i),*,err=50,end=50)ii
00129 do 40 j=ii+1,nl
00130 xlevel(j-1) = xlevel(j)
00131 40 continue
00132 if (nl.gt.0) nl = nl - 1
00133 50 continue
00134 end if
00135
00136 else if (c.eq.'e') then
00137 c
00138 c Erase the screen.
00139 c
00140 call clear
00141 call box(0.,xlen,0.,ylen)
00142 ipaint = 0
00143
00144 else if (c.eq.'g') then
00145 call devoff
00146 write(6,'(a)')'Entering graphics mode.'
00147 &
00148 call devon
00149 c
00150 c Graphic input.
00151 c
00152 60 call getgfx(r,s)
00153 if (r.le.0..or.s.le.0..or.r.ge.xlen.or.s.ge.ylen) go to 5
00154 if (nl.eq.0.and.string(2:2).eq.'d') go to 5
00155 c
00156 x = 1. + r*(m-1.)/xlen
00157 y = 1. + s*(n-1.)/ylen
00158 ix = x
00159 fx = x - ix
00160 jy = y
00161 fy = y - jy
00162 xx = (1.-fy)*((1.-fx)*a2(ix,jy) + fx*a2(ix+1,jy))
00163 & + fy*((1.-fx)*a2(ix,jy+1) + fx*a2(ix+1,jy+1))
00164 c
00165 ii = 0
00166 jtype = itype
00167 if (string(2:2).eq.'d') then
00168 c
00169 c Find the closest contour.
00170 c
00171 dxmin = 1.e30
00172 jmin = 0
00173 do 65 j=1,nl
00174 dx = abs(xx-xlevel(j))
00175 if (dx.lt.dxmin) then
00176 dxmin = dx
00177 jmin = j
00178 end if
00179 65 continue
00180 if (jmin.gt.0) then
00181 xx = xlevel(jmin)
00182 do 66 j=jmin+1,nl
00183 66 xlevel(j-1) = xlevel(j)
00184 if (nl.gt.0) nl = nl - 1
00185 c
00186 ii = -1
00187 jtype = 0
00188 call invert
00189 end if
00190 else
00191 ii = iaddlevel(nl,xlevel,xx)
00192 end if
00193 c
00194 c Always draw the contour (in case weight or color has changed).
00195 c
00196 call xcontor(a,m,n,xx,1,0.,xlen,0.,ylen,jtype,
00197 & amin,amax,ipaint,invrt)
00198 if (ii.lt.0) call invert
00199 go to 60
00200 c
00201 else if (c.eq.'l') then
00202 c
00203 c Sort and list current contour levels.
00204 c
00205 call devoff
00206 if (nl.gt.0) then
00207 call sort(nl,xlevel)
00208 write(6,1000)(i,xlevel(i),i=1,nl)
00209 else
00210 write(6,'(a)')'No points.'
00211 end if
00212
00213 else if (c.eq.'o') then
00214 c
00215 c Set contour offset.
00216 c
00217 call readrq(string(istart:ns),2,
00218 & xoff,yoff,dum,dum,*5)
00219 c
00220 else if (c.eq.'p') then
00221 c
00222 c Paint the box, according to the values in a.
00223 c
00224 if (c2.ne.'-') then
00225 invrt = 0
00226 else
00227 invrt = 1
00228 end if
00229 c
00230 call color2d(xlen,ylen,a,m,n,amin,amax,invrt)
00231 ipaint = 1
00232 c
00233 else if (c.eq.'r') then
00234 c
00235 c Replay previously saved levels.
00236 c
00237 if (nl0.gt.0) then
00238 call xcontor(a,m,n,xlevel,nl0,0.,xlen,0.,ylen,itype,
00239 & amin,amax,ipaint,invrt)
00240 nl = nl0
00241 end if
00242 else if (c.eq.'s') then
00243 c
00244 c Select specified levels, delete the rest.
00245 c
00246 if (nl.gt.0) then
00247 call gettokens(string(istart:ns),list,nlist)
00248 nl1 = 0
00249 do 70 i=1,nlist
00250 read(list(i),*,err=70,end=70)ii
00251 jj = iaddlevel(nl1,temp,xlevel(ii))
00252 70 continue
00253 c
00254 c Redraw:
00255 c
00256 call clear
00257 call box(0.,xlen,0.,ylen)
00258 nl = nl1
00259 do 80 i=1,nl
00260 80 xlevel(i) = temp(i)
00261 c
00262 call xcontor(a,m,n,xlevel,nl,0.,xlen,0.,ylen,itype,
00263 & amin,amax,ipaint,invrt)
00264 end if
00265 else if (c.eq.'t') then
00266 c
00267 c Specify line type.
00268 c
00269 do 90 j=istart,ns
00270 if (string(j:j).gt.' ') go to 100
00271 90 continue
00272 itype=0
00273 go to 5
00274 100 call readiq(string(istart:ns),4,
00275 & i1,i2,i3,i4,*5)
00276 itype=1
00277 call setpat(i1,i2,i3,i4)
00278 else if (c.eq.'w') then
00279 c
00280 c Set line weight.
00281 c
00282 call readiq(string(istart:ns),1,
00283 & iweight,idum,idum,idum,*5)
00284 call weight(iweight)
00285 c
00286 else if (c.eq.'x') then
00287 nl = 0
00288 else
00289 call gettokens(string,list,nlist)
00290 nl0 = nl
00291 do 110 i=1,nlist
00292 if (list(i)(1:1).eq.'#') then
00293 read(list(i)(2:len(list(i))),*,err=110,end=110)ix
00294 if (ix.le.0.or.ix.gt.nl) go to 110
00295 xx = xlevel(ix)
00296 else
00297 call readrtoken(list(i),xx,xx)
00298 ii = iaddlevel(nl,xlevel,xx)
00299 end if
00300 c
00301 c Always draw the contour (in case weight or color has changed).
00302 c
00303 call xcontor(a,m,n,xx,1,0.,xlen,0.,ylen,itype,
00304 & amin,amax,ipaint,invrt)
00305 110 continue
00306 end if
00307 c
00308 nl0 = nl
00309 go to 5
00310 c
00311 9999 return
00312 99999 return 1
00313 end
00314
00315
00316 subroutine get_dimensions(m,n,na,iprompt,*)
00317 c
00318 c Read in and check array dimensions.
00319 c
00320 character string*80
00321 character*50 list(20)
00322 c
00323 mm = 0
00324 nn = 0
00325 if (iprompt.eq.1) then
00326 call getstring('Array sub-dimensions:',1,21,string)
00327 else
00328 call getstring(' ',1,1,string)
00329 end if
00330 c
00331 call gettokens(string,list,nl)
00332 if (nl.eq.1) then
00333 call readitoken(list(1),mm,mm)
00334 if (mm.le.0) return 1
00335 nn = na/mm
00336 if (iprompt.eq.1) write(6,'(a,i6,a)')
00337 & 'Adopted ',nn,' as second array dimension.'
00338 else
00339 call readitoken(list(1),mm,mm)
00340 if (mm.le.0) return 1
00341 call readitoken(list(2),nn,nn)
00342 if (nn.le.0) return 1
00343 end if
00344 c
00345 if (mm.le.0.or.nn.le.0) then
00346 if (m.gt.0.and.n.gt.0) then
00347 if (iprompt.eq.1) write(6,'(a)')
00348 & 'Using previous values '
00349 & '-- may not replay correctly...'
00350 else
00351 if (iprompt.eq.1) write(6,'(a)')'No previous values.'
00352 return
00353 end if
00354 else
00355 m = mm
00356 n = nn
00357 end if
00358 c
00359 if (m*n.gt.na) then
00360 write(6,'(a)')'Too big!'
00361 return 1
00362 end if
00363 c
00364 end
00365
00366
00367 integer function iaddlevel(nl,xlevel,xx)
00368 save
00369 c
00370 c Add the specified level to the contour list, if it isn't on it.
00371 c
00372 dimension xlevel(1)
00373 c
00374 iaddlevel = 0
00375 do i=1,nl
00376 if (xx.eq.xlevel(i)) return
00377 end do
00378 c
00379 iaddlevel = 1
00380 nl = nl + 1
00381 xlevel(nl) = xx
00382 c
00383 end
00384
00385
00386 subroutine xcontor(a,m,n,v,nv,x1,x2,y1,y2,itype,
00387 & amin,amax,ipaint,invert)
00388 save
00389 c
00390 dimension v(1)
00391 common /local offset/ xoff,yoff
00392 c
00393 c Call contor or dcontor.
00394 c
00395 if (ipaint.eq.0) then
00396 if (itype.eq.0) then
00397 call contor(a,m,n,v,nv,x1+xoff,x2+xoff,y1+yoff,y2+yoff)
00398 else
00399 call dcontor(a,m,n,v,nv,x1+xoff,x2+xoff,y1+yoff,y2+yoff)
00400 end if
00401 else
00402 call getcolor(icsave)
00403 a2 = .5*(amin+amax)
00404 c
00405 do i=1,nv
00406 c
00407 if ((v(i).le.a2.and.invert.eq.0)
00408 & .or.(v(i).ge.a2.and.invert.eq.1))
00409 & call color(ncolors()-2)
00410 c
00411 if (itype.eq.0) then
00412 call contor(a,m,n,v(i),1,
00413 $ x1+xoff,x2+xoff,y1+yoff,y2+yoff)
00414 else
00415 call dcontor(a,m,n,v(i),1,
00416 $ x1+xoff,x2+xoff,y1+yoff,y2+yoff)
00417 end if
00418 call color(icsave)
00419 end do
00420 end if
00421 c
00422 end
00423
00424
00425 subroutine color2d(xlen,ylen,array,nx,ny,amin0,amax0,invert)
00426 save
00427 c
00428 c Color a rectangular region according to the values in array.
00429 c
00430 real array(nx,ny)
00431 common /color limits/ icmin,icmax
00432 c
00433 character*80 device
00434 real xp(4),yp(4)
00435 c
00436 icmin = 0
00437 icmax = 0
00438 c
00439 c Color Suns, X and PostScript only!
00440 c
00441 call devicename(device)
00442 if (device(1:1).ne.'s'.and.device(1:1).ne.'p'
00443 $ .and.device(1:1).ne.'x') return
00444 if (ncolors().le.2) return
00445 c
00446 c Set up the ranges of the colors to be used.
00447 c
00448 if (device(1:1).eq.'s') then
00449 c
00450 c Avoid strange colors at the ends of the SunCore color map!
00451 c
00452 icmin = 2
00453 icmax = ncolors() - 2
00454 else if (device(1:1).eq.'x'.and.ncolors().gt.128) then
00455 icmin = 1
00456 icmax = ncolors() - 2
00457 else
00458 icmin = 0
00459 icmax = ncolors() - 1
00460 end if
00461 nc1 = icmax - icmin
00462 c
00463 amin = amin0
00464 amax = amax0
00465 if (amin.ge.amax) call minmax(array,nx*ny,amin,amax)
00466 if (amin.ge.amax) return
00467 c
00468 c Assume uniform spacing along both axes.
00469 c
00470 dr = xlen/(nx-1.)
00471 ds = ylen/(ny-1.)
00472 c
00473 c Color the box.
00474 c
00475 do 100 i=1,nx-1
00476 do 100 j=1,ny-1
00477 rr = (i-1)*dr
00478 ss = (j-1)*ds
00479 xp(1) = rr
00480 yp(1) = ss
00481 xp(2) = rr
00482 yp(2) = ss + ds
00483 xp(3) = rr + dr
00484 yp(3) = ss + ds
00485 xp(4) = rr + dr
00486 yp(4) = ss
00487 c
00488 c Base color on average level in cell...
00489 c
00490 aa = .25*(array(i,j) + array(i+1,j)
00491 & + array(i,j+1) + array(i+1,j+1))
00492 c
00493 c Set up and clip the color, according to amin and amax:
00494 c
00495 ic = nint(icmin+nc1*(max(amin,aa)-amin)/(amax-amin))
00496 ic = max(icmin,min(icmax,ic))
00497 c
00498 c Inverse color map:
00499 c
00500 if (invert.ne.0) ic = icmax + icmin - ic
00501 c
00502 call setfill(ic)
00503 call polyfill(xp,yp,4)
00504 c
00505 100 continue
00506 c
00507 end
